Improved FEM natural frequency calculation for structural frames by local correction procedure

Citation: Urruzola, J.; Garmendia, I. Improved FEM Natural Frequency Calculation for Structural Frames by Local Correction Procedure. Buildings 2024,14, 1195. https://doi.org/ 10.3390/buildings14051195 Academic Editors: Shaohong Cheng and Haijun Zhou Received: 18 March 2024 Revised: 12 April 2024 Accepted: 14 April 2024 Published: 23 April 2024 Copyright: © 2024 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/). buildings Article Improved FEM Natural Frequency Calculation for Structural Frames by Local Correction Procedure Javier Urruzola and Iñaki Garmendia * Mechanical Engineering Department, Engineering School of Gipuzkoa, University of the Basque Country UPV/EHU, Plaza de Europa, 1, E-20018 Donostia-San Sebastián, Spain; javier[email protected] *Correspondence: [email protected] Abstract: The accurate calculation of natural frequencies is important for vibration and earthquake analyses of structural frames. For this purpose, it is necessary to discretize each beam or column of the frame into one or more smaller elements. The required number of elements per member increases when the frame’s modal shapes have wavelengths similar to the beam lengths. This paper presents a method that reduces the number of elements needed for a precise calculation. This is achieved by implementing a straightforward local correction to the kinetic and elastic energy of certain elements, resulting in a substantial decrease in error. The validity of this method is demonstrated through a range of examples, from simple canonical cases to more realistic ones. Additionally, the paper discusses the unique features of this method and examines its relationship with other approaches. Keywords: structure; frame; mechanical; vibration; natural frequency; finite element; beam; column 1. Introduction Modal analysis and natural frequency calculation by the FEM are very valuable tools to study the dynamic behavior of building structures [ 1 ]. For example, the Spanish Structural Code [ 2 ], the Eurocode [ 3 ] and the American Code ASCE [ 4 , 5 ] accept their validity for seismic analysis and wind load induced vibration analysis. Therefore, the most popular structural analysis programs such as ETABS, Robot or Staad implement these numerical techniques. Frames are usually modelled using one element per member (beam or column), which is accurate enough for linear structural analysis, but it falls short for vibration eigenvalue problems [ 6 , 7 ] because of the inadequacy of polynomials to represent localized modal shapes. Therefore, the need arises to develop methods to perform modal analysis in a more accurate way with the least numerical cost and implementation effort. Consequently, several approaches have been proposed in the scientific literature to estimate the incurred error and possibly reduce it, including correction formulas, the superconvergent patch recovery technique (SPR), the hierarchical FEM (HFEM), the smoothed FEM (SFEM), the massredistributed FEM (MRFEM) and the use of various higher-order beam finite elements. Correction formulas were applied by Xie and Steven [ 8 ] to improve the accuracy of the FEM calculation of natural frequencies in beam/column elements. Their approach stems from a previous study by Mackie [ 9 ] on the topic of numerical dispersion error reduction. A similar technique [ 10 ] can be applied to linear structural buckling critical load calculations. Their method can also be applied to structural frames by means of a weighted average of single beam/column correction terms. FEM error estimates [ 11 , 12 ] deal with the problem of numerical inaccuracy induced by the discretization of the continuum of differential equations. They are more detailed than convergence graphs and can be used to refine the mesh where necessary to achieve a certain level of precision. Residual-based estimators measure the error on the exact differential equations [ 11 ] while recovery-based estimators build a better approximation of Buildings 2024,14, 1195. https://doi.org/10.3390/buildings14051195 https://www.mdpi.com/journal/buildings Buildings 2024,14, 1195 2 of 19 the displacement or stress field [ 13 – 15 ] that can be used to obtain a more precise value of the natural frequencies. The superpatch recovery technique (SPR) [ 16 ] uses a patch of neighboring elements to adjust a higher order polynomial to approximate the stress in a finite element using the element values as well as those in the conveniently weighted patch. It originates [ 14 , 17 ] from the idea of fitting an improved stress distribution field to a set of so-called superconvergent points, when they exist, where stresses are calculated with a higher accuracy. Wiberg et al. [ 18 ] fitted the polynomial to displacements instead of stresses (SPRD) in order to improve the calculated value of natural frequencies, which depend not only on the displacement derivatives, but also on the displacements themselves. The smoothed finite element method (SFEM) [ 19 ] uses a gradient smoothing technique to reduce the overstiffening of the FEM. This method is based on the G space theory [ 20 ] that makes it possible to use discontinuous shape functions in the element formulation while maintaining stability and convergence to the exact solution. The node-based smoothed finite element method (NS-FEM) [ 21 ] improves accuracy and gives an upper bound of the elastic energy, whereas the edge-based smoothed finite element method (ES-FEM) [ 22 ] provides a lower bound. The hierarchical FEM [ 23 ] employs nested polynomial shape functions of different orders to increase the accuracy of the elements when necessary. Therefore, it can be used for error estimation and adaptative mesh refinement. Early application of the method to dynamic analysis focused on Bernoulli–Euler beams [ 24 ]. More recently, the method has been applied to various types of beams such as Timoshenko beams [ 25 ], three-dimensional sandwich beams [26], etc. Modifying the element mass matrix is another strategy to improve natural frequency calculations. Fried and Chavez [ 27 ] used a weighted average of the consistent matrix and the lumped matrix to model strings and membranes. A more economical alternative was developed by Fried and Leong [ 28 ] using the consistent matrix for the modal shape calculation and a weighted mass matrix for a Rayleigh quotient correction. Li and He [ 29 ] changed the location of the Gaussian points used to integrate the element mass matrix. Their approach stems from previous work in acoustics by Gudatti [30]. Higher-order Euler–Bernoulli elements [ 31 ] and the Timoshenko element [ 32 – 34 ] provide improved accuracy due to the better representation of displacements. This leads to a reduction in the number of elements required to calculate natural frequencies. However, their implementation is more complex, and the resulting equations have more unknown variables and will be worse conditioned [ 35 ]. Thin-walled beams [ 36 ] also require enriched sets of modelling variables because of their complex geometrical and deformation patterns. The present paper shows a new method for improving the accuracy of the calculation of the natural frequencies of structural frames when beams and columns are modelled with a small number of elements (possibly one or two). Sway frames can often be analyzed accurately with one element per member but non sway ones usually require a finer mesh or a higher precision technique like ours. The algorithm proceeds in two stages. In the first stage, a coarse solution is calculated whereas in the second one, local corrections are added at a finer level. If necessary, some elements are subdivided if the local correction excessively distorts the modal shape. Concerning the novelty of this work, the authors recently wrote a closely related paper [ 37 ] about calculating the critical buckling loads of structural frames using one element per member. This latest work presents some fundamentally novel developments. First, preventing structural buckling requires knowing just the lowest critical load, but in order to model structural dynamics accurately, multiple natural frequencies are needed and the algorithm has to be modified accordingly. Second, the interplay between multiple frequencies coupled with the limited accuracy of individual elements leads to using two subelements instead of four. Third, because of the same reasons, some members will have to be modelled with more than one element per member according to a novel specific element Buildings 2024,14, 1195 3 of 19 distortion criterion that we will later introduce. Lastly, the derivation of the equations and algorithms has been optimized for clarity, ease of implementation and performance. The subsequent sections of this paper are outlined next. First, natural frequency calculations are carried out for five fundamental cases of one-bar (beam/column) structural elements with the aim of assessing the error associated with coarse meshes and the problems arising from multiple modal interplay. Second, the vibration modes of these bars are modified by a local correction procedure and the elements are subdivided in two according to a distortion criterion. Third, the method is extended to structural frames made up of more than one bar. Next, the devised algorithms are validated using 2D and 3D cases representative of realistic building structures taken from [ 10 , 37 ]. Finally, the results, discussion and conclusions are presented. 2. Natural Frequency Analysis of Some Fundamental Cases Using One Element Per Bar We have selected a set of fundamental cases [ 38 ] (see Figure 1) to test our method against the standard FEM. Various support conditions such as clamped (C), pinned (P) and free (F) are considered. We have added the case of the second mode of the pinned–pinned beam (PP2) because it will help us better explain the algorithm. Figure 1. Five fundamental beam vibration cases. Natural frequencies of a structural frame can be calculated by the FEM as the solution of the eigenvalue problem K−ω2Mϕ=0 (1) where K and M are the stiffness and mass matrices, ω is any natural frequency and ϕ is its corresponding modal shape. Table 1shows the accuracy of the FEM calculation with Nel cubic elements. The relative errors using a single element are excessive in all cases except CF. Errors larger than 1% are considered excessive from a structural engineering point of view [37]. Table 1. Relative error 1in natural frequency computation for some fundamental cases. Nel 2CC CP PP CF PP2 1 - 32.92% 10.99% 0.48% 27.14% 2 1.62% 0.93% 0.39% 0.05% 10.98% 3 0.41% 0.20% 0.08% 0.01% 1.17% 4 0.13% 0.06% 0.03% 0.00% 0.38% 1 Relative error refers to “nearly exact” values calculated with Abaqus and Nel = 10. 2 Nel: number of elements in the discretization. Looking at Figure 1and by analogy with the column buckling problem, we can interpret that a single element is not accurate enough to model more than a quarter of Buildings 2024,14, 1195 4 of 19 a sinusoidal deformation wavelength. We will see how to reduce these errors in the next section. 3. Corrected Calculation of Natural Frequencies in Some Fundamental Cases Using One Element Per Bar In this section, we will improve the quality of the displacements inside the structural element in two ways: (1) we will use an auxiliary discretization of the bar elements (see Figure 2) with two subelements and three nodes (1–3) to obtain a local correction of the coarse mesh solution and (2) we will split some elements in half when necessary (adaptive mesh refinement). This approximation results in acceptable errors near those obtained in Table 1with four elements. Figure 2. Bar element “global” displacement ug(from 1–3 to 1′–3′). In our previous work on buckling [ 37 ] we used an auxiliary discretization of four subelements instead of two, but we will see that this approach is not possible when calculating multiple eigenvalues (or multiple frequencies) because of the override problem that is later explained. We express the nodal displacements ur (r: 1–3) as the sum of a “global” term derived from the coarse solution and a “local” correction term. The global displacement ug (see Figure 2) results from a static analysis in which we fix the external displacements ug 1 and ug 2 and obtain the value of the internal nodal displacement ug 3by condensation [7]. The nodal displacements of the element nodes in the local reference frame of the bar, ul 1and ul 2, can be expressed as ul 1=ϕ1ηul 2=ϕ2η(2) where ϕ1 and ϕ2 are the coarse modal shapes at nodes 1 and 2 expressed in the local reference frame and ηis a modal amplitude variable. In order to find the inner nodal displacement ul 3 we will need the element stiffness matrix. Assuming a uniform beam, each subelement (1-3 and 3-2) will have the same stiffness and mass matrices, KS and MS , which we can express in terms of their nodal submatrices: KS=KAA KAB KBA KBB MS=MAA MAB MBA MBB (3) Buildings 2024,14, 1195 5 of 19 Assembling these subelement matrices, we obtain the stiffness and mass matrices of the refined element, Krand Mr: Kr=  KAA 0KAB 0KBB KBA KBA KAB KAA +KBB Mr=  MAA 0MAB 0MBB MBA MBA MAB MAA +MBB (4) Therefore, the sought internal displacement results in ul 3=(KAA +KBB)−1(KBAul 1+KABul 2=ϕ3η(5) where we define ϕ3as ϕ3=(KAA +KBB)−1(KBAϕ1+KABϕ2(6) The local displacement term ∆ul (see Figure 3), increases the internal node displacement ( ∆ul 3 ) without modifying the external nodal displacements (for economy of notation, we group nodal rotations and displacements in one term). ∆ul=   0 0 ∆ul 3    (7) Figure 3. Element incremental local displacements ∆ul(from 1–3 to 1′–3′). Figure 4shows the total nodal displacement of the refined element ur resulting from both the global and the local term. Figure 4. Refined element total displacements ur (global from 1–3 to 1 ′ 3 ′ and local from 1’–3’ to 1”–3”). Buildings 2024,14, 1195 6 of 19 Now, we define the refined element nodal displacements ur as a function of the modal coordinate ηand the internal incremental displacement ∆ul 3using a projection matrix P ur=  ϕ10 ϕ20 ϕ3I η ∆ul 3=Pη ∆ul 3(8) where I is the identity matrix. As a result, we can obtain a corrected natural frequency ωp by solving the projected eigenvalue problem PTKrPϕp=ω2 pPTMrPϕp(9) We summarize the procedure to calculate the corrected natural frequency in Algorithm 1. Algorithm 1. Correction of the natural frequency of a one-element bar Evaluate Krand Mras 2-subelement refinements of Keand Meusing Equation (3) Evaluate the projection matrix P with Equation (8) Calculate ωpas the lowest natural frequency in Equation (9) After the correction process the errors in natural frequencies change as shown in Table 2. Table 2. Relative error in corrected natural frequency calculation (Nel = 1). Nel CC CP PP CF PP2 1 1.61% 0.93% 0.39% 0.05% 42.42% Looking at these results, we can see that correcting the one element per member model in the CC and PP2 cases does not reduce the error up to a level that is acceptable in engineering. We can understand what is happening if we study what we will designate as the distortion factor: the maximum change in V or T after applying the correction (see Table 3), where V and T are the stiffness and mass quadratic forms, respectively. Table 3. Distortion factor (maximum percentage change in V or T) after correction (Nel = 1). Nel CC CP PP CF PP2 1 -% 211.33% 49.66% 1.73% 1.4 ×1029% What we can see here is that the local correction has largely distorted V and/or T in both cases (most notably for PP2). First, we will study the cause for the most troubling case, PP2. We can see a graphical depiction of its distortion with altered scales in Figure 5. The softer CC mode (K = 22.4) has almost completely overridden the stiffer PP2 mode (K = 39.5) thereby suppressing an existing mode and replacing it with a rough duplicate of a previously calculated one (CC). Second, we turn our attention to the origin of the slightly unacceptable error of the CC case. We can attribute it to the fact that the correction process can never surpass the accuracy of doubling the element at the coarse level. In order to solve both problems, we propose splitting elements in half when the distortion factor (from now on called γ ) surpasses the 100% threshold. The PP2 spurious modal override problem will be solved because the half element corrections cannot roughly represent the CC mode in isolation. In turn, the CC slightly unacceptable error will be reduced because of the superior quality of the refined coarse mesh. It Is clear now that using a four subelement discretization for local element correction is not acceptable because it will be plagued by the override problem in the same way as the two-subelement one. Buildings 2024,14, 1195 7 of 19 Figure 5. Locally corrected PP2 modal shape using one element and two subelements. Therefore, we will allow one or two elements per member at the coarse level and two subelements at the local level, which allows for a total of four subelements per member, roughly equivalent to what we did for buckling in [37]. In the next section, we extend this correction/refinement process to general structural frames made up of multiple bars and one or two elements per bar. 4. Correction of Natural Frequencies for Multiple-Element Structures Similarly to what we did in [ 37 ], we are going to generalize the procedure for single elements based on four main ideas: 1. The local element corrections can be combined additively into an overall modal correction. 2. When calculating local corrections for an element, the rest of the structure can be sufficiently represented by the frame modal shape ϕand an amplitude variable η. 3. The corrected natural frequency for the whole frame can be calculated using Rayleigh’s quotient with the corrected modal shape. 4. Local corrections for different natural frequencies can be calculated in isolation from each other once the distortion factor has been introduced to solve the override problem. Let us examine how the whole procedure would work for our most problematic case, PP2. Figure 6shows the modal shape of the PP2 case beam discretized at the coarse level with two elements per member. In order to improve the quality of the modal shape, we will fix the end nodal displacements of each coarse element and subsequently correct its inner displacements, and, as a result, we will obtain the corrected vibration shape in the same figure. Figure 6. PP2 modal shape using two elements after local correction. Buildings 2024,14, 1195 8 of 19 The local corrections for each element will be calculated separately (as shown in Figure 7). For this purpose, we maintain all the elements in the frame mesh except the one to be corrected, which is replaced with two subelements (see Figure 7). Figure 7. Structural discretization used to correct the upper part of two-element PP2 modal shape. The local modal shape correction is the solution of a projected eigenvalue problem that is derived below. The structure stiffness quadratic form calculated with the coarse mesh can be expressed as: V=uTKu (10) where uand K are the nodal displacement vector and stiffness matrix of the frame. We can modify Vby replacing element e contribution with its refined counterpart V=uTKu −uT eKeue+uT rKrur(11) where ue and Ke are the nodal displacement vector and stiffness matrix of element e, while ur and Kr are their refined versions (using two subelements and three nodes) calculated in the local reference frame of the element. Next, the nodal displacements can be expressed as a function of the modal coordinate η and the incremental inner nodal displacement of the element being corrected ∆ul 3 , similarly to what we did for a single-element bar in Equation (8). u=ϕη (12) ue=ϕeη(13) ur=   ϕ1η ϕ2η ϕ3η+∆ul 3    (14) and the same operations can be performed on the mass quadratic form: T=uTMu −uT eMeue+uT rMrur(15) As a result, we obtain a projected eigenvalue problem whose solution contains the local modal shape correction: Kpϕp=ω2 pMpϕp(16) where Kp=V−Ve+Vr0 0KAA +KBB(17) Buildings 2024,14, 1195 9 of 19 Mp=T−Te+TrϕT 1MAB +ϕT 2MBA +ϕT 3(MAA +MBB symmetric MAA +MBB (18) ϕr=   ϕ1 ϕ2 ϕ3   (19) Vr=ϕT 1KAAϕ1+ϕT 2KBBϕ2+ϕT 1KABϕ3+ϕT 3KABϕ2(20) Tr=ϕT 1MAAϕ1+ϕT 2MBBϕ2+ϕT 3(MAA +MBB)ϕ3+2ϕT 1MABϕ3+2ϕT 3MABϕ2(21) We can partition the projected mode in terms of its modal amplitude part ϕp0 and its internal correction ϕp3 ϕpe =ϕp0 ϕp3(22) and dividing the right-hand side by ϕp0, we generate a projected mode ϕ∗ p=(1 ϕp3 ϕp0)(23) which represents the sum of the overall modal shape ϕ plus the local internal correction term ∆ϕc3 ∆ϕc3=ϕp3 ϕp0 (24) and this modal shape can be used to calculate the corrected mass and stiffness quadratic forms of the element. Vce =Vr+∆ϕT c3(KAA +KBB∆ϕc3(25) Tce =Tr+2∆ϕT c3(MBAϕ1+MABϕ2+∆ϕT c3(MAA +MBB(2ϕ3+∆ϕc3)(26) We show in Algorithm 2 the complete procedure for computing Vce and Tce. Algorithm 2. Computation of the corrected quadratic forms of an element Vce and Tce Overall frame inputs: ϕTKϕ,ϕTMϕ Frame element inputs: ϕT eKeϕe,ϕT eMeϕe,ϕe Subelement inputs: KAA,KAB,KBB,MAA,MAB,MBB Evaluate Krand Mras 2-subelement refinements of Keand Meusing Equation (3) Convert ϕeto local element coordinates by the following operations: ϕ1=RT eϕe1ϕ2=RT eϕe2(Re: element rotation matrix) Calculate ϕ3,ϕrin Equations (6) and (19) Calculate Kp,Mpin Equations (17) and (18) Obtain ϕpas the first eigenvector of Equation (16) Evaluate Vce,Tce by applying Equations (24)–(26) The element-corrected quadratic forms of the whole structure can be collected in Rayleigh’s quotient to obtain an improved value of the structure’s natural frequency ωc ω2 c=∑eVce ∑eTce (27) The complete procedure to calculate N natural frequencies is given in Algorithm 3. It should be pointed out that the denominators of the distortion factor γe have been modified to cope with the possibility of elements with very small Ve or Te , which would lead to near division by zero. Therefore, one hundredth of the frame V or T is distributed equally among all elements when measuring relative change, while the other 99% comes from the element itself. Buildings 2024,14, 1195 16 of 19 Table 12. Corrected calculation statistics (3D braced building structure with 2 elems./member). Mode # Exact ω(rad/s) 2-Elem. ω(rad/s) Relative Error (%) Corrected ω(rad/s) Relative Error (%) Distortion Factor γ(%) Distorted Elements # 1 124.49 124.54 0.03 124.50 0.00 2.35 0 2 126.75 126.79 0.03 126.75 0.00 4.65 0 3 146.27 146.33 0.04 146.28 0.00 3.94 0 4 151.65 151.71 0.04 151.66 0.00 4.33 0 5 184.03 184.13 0.05 184.04 0.00 6.87 0 6 212.50 212.66 0.08 212.52 0.01 4.66 0 7 245.42 245.67 0.10 245.44 0.01 7.20 0 8 276.69 277.13 0.16 276.73 0.01 17.03 0 9 366.95 368.89 0.53 367.08 0.04 9.62 0 10 384.71 387.27 0.67 384.89 0.05 16.12 0 11 392.67 395.59 0.74 392.87 0.05 3.82 0 12 403.95 407.28 0.82 404.18 0.06 4.83 0 After examining the results, we can conclude that our models achieve the same level of accuracy as Abaqus standard FEM models with twice as many bars (on the condition that distortion factors lie below 100%), therefore they can work with stiffness and mass matrices twice smaller. In addition, we have attained calculation times 15% smaller measuring the main components of required processing power, i.e., the whole frame eigenvalue problem and the local beam eigenvalue problems. 6. Discussion The method by Xie and Steven [ 10 ] provides local natural frequency updates for single bars rather than local modal shape corrections, which leads to using a weighted criterion with a weaker physical foundation than Rayleigh’s quotient. In addition, structural members are discretized with four or five elements while in our case one or two (in a few bars) elements are needed. The SPRD technique by Wiberg et al. [ 18 ] bears some similarities with our approach. It relies on a polynomial fitting to the existing mode at some superconvergent points while our procedure completely updates the mode at inner points without being constrained by the coarse calculation on the inside. Plus, the whole structure rather than a single element and its neighboring patch participates in the adjustment by means of the coarse modal amplitude η . It can also be noted that SPR techniques use an external patch of elements while our method relies on an inner set of refined elements. Gradient smoothing methods such as [ 19 ] also rely on neighboring elements to improve the quality of the solution but they do it before solving the system of equations of the whole structure, thereby increasing connectivity and the enhanced element matrix computation time. In addition, there is no straightforward way of applying the smoothing concept to beam elements of different sections and orientations sharing a node. Modified stiffness and mass matrices [ 29 ], while improving the accuracy of dynamic analysis, are limited by the fact that they do not depend on the natural frequency being studied (like in the case of dynamic stiffness methods) or on the actual modal shape, as happens in our method. Higher-order finite elements [ 31 – 34 ] provide better accuracy but are more complex to implement, have to solve larger systems of equations and lead to worse conditioned matrices. In contrast, our method works well with the standard finite element method and could work with higher-order finite elements as well to improve their accuracy. As for thin-walled beams [ 36 ], they require higher-order models in order to represent complex deformation patterns, but they are fully compatible with our correction algorithm. The hierarchical FEM [ 23 – 25 ] relies on error estimators to refine the structural mesh and improve accuracy. In contrast, our correction does not need a full reanalysis to increase precision but an array of concurrent element-centered corrections. If necessary, Buildings 2024,14, 1195 17 of 19 the distortion factor indicates what bars require two elements instead of one. Like the higher-order elements discussed above, hierarchical elements of higher order can be used for the corrections instead of our two-subelement set. As mentioned in the introduction, the authors recently wrote a closely related paper [37] about calculating the critical buckling loads of structural frames using one element per member. For this purpose, a local correction procedure was applied using four subelements. This latest work presents some fundamental differences. First, preventing structural buckling requires knowing just the lowest critical load, but in order to model structural dynamics accurately, multiple natural frequencies are needed. Second, the interplay between multiple frequencies coupled with the limited accuracy of individual elements makes it necessary to use two subelements instead of four. Third, because of the smaller number of subelements involved in the local corrections, some members have to be split into half beams according to a novel element distortion criterion which measures the relative change in kinetic and elastic energy caused by the correction of modal shapes. As far as the efficiency of our method is concerned, most of what was stated in [ 37 ] remains applicable: the algorithm’s greatest source of efficiency comes from its fully parallelizable nature and the local eigenvalue problems can be solved with minimal computational resources by the power method. In addition, after shrinking the submodel to two elements, all the matrices that appear in the local eigenvalue problem can be easily programmed with scalar operators and functions, thereby reducing the cost to solving the eigenvalue problem. Our technique can be applied to enhance the standard FEM analysis of any structure made up of beam/column elements. The structure could also contain shell, plate or lumped elements but the gain in accuracy would only occur for the bar elements. Our approach can be used to calculate natural frequencies with errors acceptable in engineering (below 1%) using one or two elements per member or to increase the accuracy of a calculation with any number of elements per member. Therefore, our approach offers the potential for a reduction in memory requirements and calculation speed when compared with the standard FEM at the cost of some additional coding. In order to fully exploit the advantages of the algorithm, it is advisable to distribute the correction calculations to the GPU. Concerning future research developments based on the present work, we have selected a few areas of interest. First, there is always a significant discrepancy between numerical vibrational properties and experimental measurements [ 39 , 40 ] because of approximate modelling, nonlinearities, temperature effects, etc., which could be addressed advantageously with the proposed numerical technique or an enhanced version of it. Likewise, in the field of health structural monitoring, there is also the need to deal with discrepancies caused by structural failure or deterioration, and to diagnose their nature and location [ 41 , 42 ]. However, the present study is only concerned with numerical efficiency and has no direct application in these areas in its present form. 7. Conclusions A new method for improving the accuracy of the standard FEM natural frequency calculation of structural frames made up of beam/column elements has been presented. The algorithms are based on previous work by the authors on structural frame buckling, but significant novel modifications have been made to improve efficiency and ease implementation, and to account for the challenges of calculating several eigenvalues instead of the lowest one. The fully parallel nature of the method makes it very convenient to take advantage of the current trend towards GPU-based architectures. For this purpose, the main calculation cost driver is a small individual nodal centered eigenvalue problem solvable with a few power iterations. Structural members are modelled with one or two elements following a novel subdivision criterion based on the degree of distortion caused by the correction of the original modal shape. As a result, enough accuracy for engineering applications is achieved with a modest increase in computation time and storage requirements. The approach is very Buildings 2024,14, 1195 18 of 19 flexible and can accommodate different beam types, higher-order models and finer meshes and target precision levels, even though it has been demonstrated with simple cubic elements. Algorithm inputs are readily available data on FEM codes such as frame and element stiffness and mass matrices, rotations and modal shapes that can be processed with scalar functions and operators before the local eigensolver correction step. Author Contributions: Conceptualization, J.U.: methodology, J.U.; software, J.U.; validation, J.U. and I.G.; writing, J.U. and I.G. All authors have read and agreed to the published version of the manuscript. Funding: This research received no external funding. Data Availability Statement: The original contributions presented in the study are included in the article. Further inquiries can be directed to the corresponding author. Conflicts of Interest: The authors declare no conflict of interest. References 1. Chopra, A.K. Dynamics of Structures: Theory and Applications to Earthquake Engineering, 6th ed.; Pearson: London, UK, 2023. 2. Spanish Structural Code (Código Estructural); Spanish Ministry of Transport, Mobility and the Urban Agenda (Ministerio de Transportes, Movilidad y Agenda Urbana): Madrid, Spain, 2021. 3. EN1993-1-1; Eurocode 3: Design of Steel Structures. European Committee for Standardization: Brussels, Belgium, 2005. 4. Staff, American Society of Civil Engineers (ASCE). Minimum Design Loads for Buildings and Other Structures, 3rd ed.; American Society of Civil Engineers: Reston, VA, USA, 2013. 5. Biswas, P.; Peronto, J. Design and Performance of Tall Buildings for Wind, 1st ed.; American Society of Civil Engineers (ASCE): Reston, VA, USA, 2020. 6. Petyt, M. Introduction to Finite Element Vibration Analysis, 2nd ed.; Cambridge University Press: Cambridge, UK, 2015. 7. Zienkiewicz, O.C.; Taylor, R.L.; Zhu, J.Z. The Finite Element Method: Its Basis and Fundamentals; Butterworth-Heinemann: Cambridge, UK, 2013; Volume 3. 8. Xie, Y.M.; Steven, G.P. Improving finite element predictions of buckling loads of beams and frames. Comput. Struct. 1994, 52, 381–385. [CrossRef] 9. Mackie, R.I. Improving finite element predictions of modes of vibration. Int. J. Numer. Methods Eng. 1992,33, 333–344. [CrossRef] 10. Xie, Y.M.; Steven, G.P. Explicit formulas for correcting finite-element predictions of natural frequencies. Commun. Numer. Methods Eng. 1993,9, 671–680. [CrossRef] 11. Babuška, I.; Rheinboldt, W.C. A-posteriori error estimates for the finite element method. Int. J. Numer. Methods Eng. 1978, 12, 1597–1615. [CrossRef] 12. Babuška, I. Accuracy Estimates and Adaptive Refinements in Finite Element Computations; Wiley: Chichester, UK, 1986. 13. Zienkiewicz, O.C.; Zhu, J.Z. The superconvergent patch recovery (SPR) and adaptive finite element refinement. Comput. Methods Appl. Mech. Eng. 1992,101, 207–224. [CrossRef] 14. Boroomand, B.; Zienkiewicz, O.C. Recovery by equilibrium in patches (REP). Int. J. Numer. Methods Eng. 1997,40, 137–164. [CrossRef] 15. Sun, H.; Yuan, S. An improved local error estimate in adaptive finite element analysis based on element energy projection technique. Eng. Comput. 2023,40, 246–264. [CrossRef] 16. Wiberg, N.-E.; Abdulwahab, F.; Ziukas, S. Improved element stresses for node and element patches using superconvergent patch recovery. Commun. Numer. Methods Eng. 1995,11, 619–627. [CrossRef] 17. Wiberg, N.; Abdulwahab, F. Patch recovery based on superconvergent derivatives and equilibrium. Int. J. Numer. Methods Eng. 1993,36, 2703–2724. [CrossRef] 18. Wiberg, N.; Bausys, R.; Hager, P. Improved eigenfrequencies and eigenmodes in free vibration analysis. Comput. Struct. 1999, 73, 79–89. [CrossRef] 19. Liu, G.R. A generalized gradient smoothing technique and the smoothed bilinear form for Galerkin formulation of a wide class of computational methods. Int. J. Comput. Methods 2008,5, 199–236. [CrossRef] 20. Liu, G.R. On G space theory. Int. J. Comput. Methods 2009,6, 257–289. [CrossRef] 21. Liu, G.R.; Nguyen-Thoi, T.; Nguyen-Xuan, H.; Lam, K.Y. A node-based smoothed finite element method (NS-FEM) for upper bound solutions to solid mechanics problems. Comput. Struct. 2009,87, 14–26. [CrossRef] 22. Liu, G.R.; Nguyen-Thoi, T.; Lam, K.Y. An edge-based smoothed finite element method (ES-FEM) for static, free and forced vibration analyses of solids. J. Sound Vib. 2009,320, 1100–1130. [CrossRef] 23. Zienkiewicz, O.C.; De, S.R.; Gago, J.P.; Kelly, D.W. The hierarchical concept in finite element analysis. Comput. Struct. 1983, 16, 53–65. [CrossRef] 24. Ganesan, N.; Engels, R.C. Hierarchical Bernoulli-Euler beam finite elements. Comput. Struct. 1992,43, 297–304. [CrossRef] 25. Tai, C.-Y.; Chan, Y.J. A hierarchic high-order Timoshenko beam finite element. Comput. Struct. 2016,165, 48–58. [CrossRef] Buildings 2024,14, 1195 19 of 19 26. Hui, Y.; Giunta, G.; Belouettar, S.; Huang, Q.; Hu, H.; Carrera, E. A free vibration analysis of three-dimensional sandwich beams using hierarchical one-dimensional finite elements. Compos. Part B Eng. 2017,110, 7–19. [CrossRef] 27. Fried, I.; Chavez, M. Superaccurate finite element eigenvalue computation. J. Sound Vib. 2004,275, 415–422. [CrossRef] 28. Fried, I.; Leong, K. Superaccurate finite element eigenvalues via a Rayleigh quotient correction. J. Sound Vib. 2005,288, 375–386. [CrossRef] 29. Li, E.; He, Z.C. Development of a perfect match system in the improvement of eigenfrequencies of free vibration. Appl. Math. Model. 2017,44, 614–639. [CrossRef] 30. Guddati, M.N.; Yue, B. Modified integration rules for reducing dispersion error in finite element methods. Comput. Methods Appl. Mech. Eng. 2004,193, 275–287. [CrossRef] 31. Shang, H.Y.; Machado, R.D.; Abdalla Filho, J.E. Dynamic analysis of Euler–Bernoulli beam problems using the Generalized Finite Element Method. Comput. Struct. 2016,173, 109–122. [CrossRef] 32. Hsu, Y.S. Enriched finite element methods for Timoshenko beam free vibration analysis. Appl. Math. Model. 2016,40, 7012–7033. 33. Cornaggia, R.; Darrigrand, E.; Le Marrec, L.; Mahé, F. Enriched finite elements and local rescaling for vibrations of axially inhomogeneous Timoshenko beams. J. Sound Vib. 2020,474, 115228. [CrossRef] 34. Necib, B.; Sun, C.T. Analysis of truss beams using a high order Timoshenko beam finite element. J. Sound Vib. 1989,130, 149–159. [CrossRef] 35. Wang, T.; Mikkola, A.; Matikainen, M.K. An Overview of Higher-Order Beam Elements Based on the Absolute Nodal Coordinate Formulation. J. Comput. Nonlinear Dynam. 2022,17, 091001. [CrossRef] 36. Nguyen, T.; Nguyen, N.; Lee, J.; Nguyen, Q. Vibration analysis of thin-walled functionally graded sandwich beams with non-uniform polygonal cross-sections. Compos. Struct. 2021,278, 114723. [CrossRef] 37. Urruzola, J.; Garmendia, I. Calculation of Linear Buckling Load for Frames Modeled with One-Finite-Element Beams and Columns. Computation 2023,11, 109. [CrossRef] 38. Meirovitch, L. Fundamentals of Vibrations; McGraw-Hill: New York, NY, USA, 2001. 39. Luo, J.; Huang, M.; Lei, Y. Temperature Effect on Vibration Properties and Vibration-Based Damage Identification of Bridge Structures: A Literature Review. Buildings 2022,12, 1209. [CrossRef] 40. Huang, M.; Zhang, J.; Hu, J.; Ye, Z.; Deng, Z.; Wan, N. Nonlinear modeling of temperature-induced bearing displacement of long-span single-pier rigid frame bridge based on DCNN-LSTM. Case Stud. Therm. Eng. 2024,53, 103897. [CrossRef] 41. Deng, Z.; Huang, M.; Wan, N.; Zhang, J. The Current Development of Structural Health Monitoring for Bridges: A Review. Buildings 2023,13, 1360. [CrossRef] 42. Huang, M.; Ling, Z.; Sun, C.; Lei, Y.; Xiang, C.; Wan, Z.; Gu, J. Two-stage damage identification for bridge bearings based on sailfish optimization and element relative modal strain energy. Struct. Eng. Mech. Int. J. 2023,86, 715–730. Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.