Next Article in Journal
Experimental Fitting of Efficiency Hill Chart for Kaplan Hydraulic Turbine
Previous Article in Journal
LSTM Networks for Home Energy Efficiency
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Dynamic Error Estimation in Higher-Order Finite Elements

Department of Mechanical, Energy and Management Engineering (DIMEG), University of Calabria, Cubo 45C, 87036 Rende, Italy
*
Author to whom correspondence should be addressed.
Designs 2024, 8(4), 79; https://doi.org/10.3390/designs8040079
Submission received: 20 June 2024 / Revised: 6 August 2024 / Accepted: 8 August 2024 / Published: 11 August 2024

Abstract

:
The Finite Element Method (FEM) has emerged as a powerful tool for predicting the behavior of industrial products, including those with complex geometries or uncommon materials. Finite Element Analysis (FEA) is widely used to study structural vibration-related aspects such as stress, displacement, and velocity. Modal analysis, a standard technique for characterizing the vibrational behavior of structures, is essential for identifying resonance frequencies, optimizing component design, and assessing structural integrity. Finite Elements (FE) modal analysis enables engineers to evaluate numerically the modal parameters, whereas model order reduction (MOR) schemes are exploited to achieve a balance between computational efficiency and accuracy, enabling a more efficient solution for computing transient dynamic analysis. Assessing the accuracy and reliability of FE solutions is a crucial aspect of the design cycle, and model-updating procedures are commonly employed to maximize the correlation between measured and predicted dynamic behavior. This study investigated the accuracy and computational efficiency of linear, quadratic, and cubic hexahedral FE formulations for modal analysis and transient dynamic solutions. More specifically, the documented results demonstrate the profitable use of the eigenenergy norm obtained in eigen solutions as a valid predictor of the accuracy reported using either the time response assurance criterion (TRAC) or the frequency response assurance criterion (FRAC), measured in transient dynamic cases. Moreover, our results also highlight the superior computational efficiency of higher-order formulations for both the eigen and transient dynamic solutions.

1. Introduction

In the broad spectrum of model-based design methodologies, the Finite Element Method (FEM) has surged as an invaluable tool, allowing the prediction of the behavior of an industrial manufact, even when dealing with very complex geometries or uncommon materials. Fueled by increasing computational power, the use of FEM in real-life 3D engineering projects allows the representation of complex physical phenomena to significantly reduce design costs [1].
Over the last few decades, Finite Element Analysis (FEA) has been affirming as one of the most widely used methodologies for studying structural vibration-related aspects such as stress, displacement, and velocity [2,3]. Finite Elements (FE) structural dynamics analysis is crucial for assessing the dynamic behavior of a product and verifying that its vibration patterns do not harm the expected functionality. Moreover, in the context of the vibration analysis of complex structures, the advent of FEA has complemented the use of several analytical or experimental techniques, such as modal analysis and wave propagation methods [4,5,6]. In highly competitive industrial fields, FEA shortens the time to market and supports designers in making informed choices, even during the concept modeling stages [7].
Modal analysis represents the de facto standard used to experimentally characterize the vibrational behavior of structures and machinery. The accurate estimation of modal parameters is essential for identifying resonance frequencies, optimizing component design, and assessing structural integrity [6,8]. Moreover, determining modal parameters, such as the natural frequencies, mode shapes, and related modal damping ratios, allows the formulation of a modal model [9], which is instrumental in more advanced applications such as damage detection and remaining useful life estimation [10]. FE modal analysis enables engineers to numerically evaluate the modal parameters [11,12,13,14]. For example, the use of modal model reduction schemes in the well-known Nastran Modal Transient Response (sol. 112) fueled the use of FEA for the direct solution of the transient dynamics of a vibrating structure. Although modern computing architectures allow the handling of complex transient dynamic problems, they still require large processing time and memory. Model order reduction (MOR) is regularly used to achieve a better balance between computational efficiency and accuracy [15]. For example, an advanced parametric-MOR scheme that exposes the material parameter dependencies out of the obtained reduction space was recently proposed: this enabled the exploitation of the reduced model for all cases where a canonical MOR would demand a recompilation of the reduction process, which could potentially lead to more advanced condition monitoring applications [16]. In this study, we streamlined a design procedure that allows engineers to increase the computational efficiency of a full-order FE model by selecting a proper order formulation, which in turn can have a beneficial impact on the existing MOR-based method.
In summary, being the use of FEA prominent in modern industrial practices, assessing the accuracy and reliability of FE solutions has become crucial, creating the need to adopt industrial procedures dedicated to classifying the errors and quantifying the level of uncertainty related to the FEA results. Interested readers may find more details regarding the FE model-updating process in previous studies [17,18]. The ultimate goal of model updating is to create a finite element model that yields reliable and accurate predictions of the dynamic behavior of a mechanical structure [6]. In practice, the model updating procedure requires two stages. First, dedicated FE model validation campaigns are generally performed [19], and then the selected FE model updating procedures are executed to maximize the correlation between the measured behavior of the structure under dynamic excitation and the corresponding FE predictions [20].
Several studies have been conducted to improve FE accuracy and reliability. The more commonly adopted strategy relies on adaptive mesh refinement [21,22,23], which consists of increasing the mesh density and/or the polynomial order of the interpolation functions. Moreover, attention has been dedicated to developing more effective FE formulations, proposing different interpolation schemes, and improving meshing procedures to guarantee tessellation regularity. However, a clear demonstration of the accuracy benefits of higher-order FE formulations has not yet been documented. Therefore, this study extends the preliminary work published in [24] by further investigating the accuracy and computational efficiency of linear, quadratic, and cubic hexahedral FE formulations for both the modal analysis solution and the transient dynamic one. For the modal accuracy, this study relies on a recently proposed error metric, namely the normalized eigen-energy (NEE) norm [25,26], which has the merit of considering both the eigenvalues relative errors and the eigenvector accuracy. For the transient dynamics, error analysis was performed using the time response assurance criterion (TRAC) and the frequency response assurance criterion (FRAC). Hence, this research aims to establish a relationship between the NEE metric, which can be evaluated by performing the FE eigenvalue solution, and the TRAC and FRAC metrics, which necessitate more intricate dynamic analysis to be evaluated.
The relationship between the modal and transient dynamic error analyses was assessed by considering a numerical case study for different hexahedral FE formulations ranging from linear to cubic order. After conducting a convergence study for modal analysis accuracy, a transient dynamic simulation was performed for three case studies opportunely selected for different levels of target accuracy (Cases 1 and 2), or by fixing the number of degrees of freedom (DOFs) and varying the order of the interpolant functions (Case 3). Based on the obtained results, this study empirically demonstrates that a strong correlation exists between the accuracy error measured using the NEE metric and the TRAC and FRAC values obtained in transient dynamics, irrespective of the order of the formulation used.
The remainder of this paper is structured as follows: Section 2 summarizes the background of the formulated design methodology, focusing on recalling the FE formulation of modal and transient dynamic problems. Section 3 describes the different metrics adopted for assessing the FE error accuracy in dynamic problems. Section 4 presents and discusses the results. Finally, concluding remarks are presented in Section 5.

2. Material and Methods

The FE design methodology proposed in this study involves higher-order FE formulations of structural dynamic problems. The necessary background is briefly recalled in this section while also discussing any impact on the overall FE model accuracy. Interested readers may find a complete and detailed derivation of the FE formulation in [24]. All the FE formulations presented in this study were implemented from scratch using the MATLAB coding environment.

2.1. Higher Order FE Formulations

The motion quantities (i.e., displacement, velocities, or accelerations) describing the vibration pattern of a structural component are modeled as continuous field quantities within the spatial domain of investigation, which can be limited to the space coincident with the volume of the considered object, assuming that the small-vibration hypothesis holds. In such circumstances, provided that the considered space can be divided into a countable finite set of cells, which are called finite elements, the continuous field motion quantities are mapped to a finite set of the corresponding nodal quantities by exploiting the interpolation functions given by the FE formulation. Thus, FEM allows the calculation of the corresponding field solutions over the entire domain [27,28]. For instance, within each element, the unknown displacement field u x , y , z , is approximated as a linear combination of a set of n interpolation functions n i ξ , η , ς , also known as shape functions [21]:
u x , y , z = N   q
where N is the matrix of the shape functions and q is the nodal displacement vector. It is worth highlighting that the shape function formulation plays an important role in how the numerical model captures rigid body motion and constant deformation conditions, while ensuring compatibility, completeness, and continuity across the inter-element boundaries; thus, it plays a crucial role in the accuracy of the solution [29].
The global mass and stiffness matrices were computed according to the isoparametric formulation, which in turn required the selection of a specific integration scheme. Gauss quadrature integration is widely used in FEA, Isogeometric analysis, and more generally within Galerkin discretization schemes [30] because it is known to improve accuracy and computational efficiency. Extended versions of the Gauss quadrature have recently been proposed, which allow the capture of discontinuities within elements [31], thus leading to more accurate integration in the case of discontinuous problems such as, for example, when investigating crack singularities [32]. Moreover, Gaussian integration is known to avoid the Runge phenomenon, which can increase the computational error and cause computational instabilities [33]. Based on the above considerations, the Gauss quadrature integration rules were selected for this study. The evaluation of the mass and stiffness matrices and of the force vector was obtained using the following weighted sums:
K e g w g   B g T D B g det J M e g w g   ρ   N g T N g det J F e g w g   N g T f b det J + g w g   N g T f s det J
where wg are the weight associated with each sampling point; Bg is the deformation tensor; D is the linear elasticity tensor; ρ is the material density; and J is the Jacobian matrix used to perform the transformation from the local coordinate system ξ , η , ς of each element to the global Cartesian system. Monitoring the Jacobian during the computations allows to check for excessive distortion of the element shape, thereby preventing for incurring into a negative volume error during the simulation, which would have led to incorrect results.
The locations and weights of the sampling points were selected to ensure exact integration of a polynomial of the maximum degree [34]. Motivated by the preliminary results discussed in [24], this study did not consider for the investigation the reduced integration schemes [35], which are commonly used in combination with the serendipity elements available in commercial codes, such as HEXA20 or HEXA32. In fact, by limiting the use to only full integration schemes, we prevented the surge of instability effects such as the “hourglass effects.” The computation of the global mass and stiffness matrices, M and K, as well as the global force vector, F, was performed by considering the contributions of every element (Equation (2)) and the inter-element connectivity.
K = e K e ; M = e K e ; F = e F e .
Consequently, the motion of a vibrating specimen can be investigated by considering the following discretized system of second-order ordinary differential equations (ODEs) [30]:
M   q ¨ + C   q ˙ + K   q = F
where M = m i j n × n , C = c i j n × n , K = k i j n × n are the mass, damping, and stiffness matrices of the system under investigation, respectively; n is the number of spatial DOFs; q, q ˙ , q ¨ are the vectors of nodal displacements, velocities, and accelerations, respectively; and F is the generalized nodal force vector. It can be demonstrated that for a physical system, these matrices are real-valued, symmetric, and positive definite [30]. After enforcing symmetry, to compensate for the truncation error introduced during the computation pipeline, the damping matrix C was obtained as a linear combination of the mass and stiffness matrices by adopting the Rayleigh Damping formulation (Appendix A.1).
The effects of structural constraints need to be considered before performing any form of FE solution. In this study, we limit ourselves to the Dirichlet boundary conditions, which generally consist of assigning a specific value to the constrained DOFs. In this context, the application of a set of kinematic constraints is formulated by splitting the DOFs, q, into the constrained DOFs, qc, and the independent ones, qi. After rearranging Equation (4), we get:
M ii M ic M ci M cc q ¨ i q ¨ c + C ii C ic C ci C cc q ˙ i q ˙ c + K ii K ic K ci K cc   q i q c = F i F c
Enforcing a given type of constraint on a node within a FE model is commonly obtained by imposing that the related DOFs are flagged as fixed, owing to imposing null displacements, q c = 0 , null velocities, q ˙ c = 0 , and null accelerations, q ¨ c = 0 . Rearranging the terms in Equation (5), and substituting the boundary conditions, we get:
M ii q ¨ i + C ii q ˙ i + K ii q i = F i
F c = M ci q ¨ i + C ci q ˙ i + K ci q i
It is worth noting that Equation (6) is structurally equivalent to the form of the unconstrained system expressed in Equation (4). The matrices to solve the dynamics of the constrained system can be obtained by removing the rows and columns corresponding to the fixed DOFs. To keep the notation simpler, in the remainder of the paper, pedices ii and i are dropped. After solving for the independent DOFS, Equation (7) can be used to quantify the reaction forces due to the constrains.
The Laplace transformation of Equation (6) [6] leads to the system of algebraic equations in Laplace variable p:
p 2 M + p C + K X p = Z p X p = F p
where Z(p) is the dynamic stiffness matrix.

2.2. FE Modal Analysis

Modal analysis is instrumental in assessing the main characteristics of the dynamic behavior of any vibrating structure. After rearranging the set of second-order ODEs (4) into an equivalent set of first-order ODEs, we obtain:
A   y ˙ + B   y = Q
where,
A = 0 M M C ; B = M 0 0 K ; y = q ˙ q ; Q = 0 F .
Modal analysis aims to solve the free vibration problem; thus, the external forces acting on the system are assumed to vanish, F = 0, leading to the formulation of the following generalized eigenvalue problem:
λ i   A   ϕ i = B   ϕ i ,   i   1 , , n  
solved by 2n complex conjugate eigenpairs (λi, ϕi) The complex conjugate roots ( λ n ,   λ n * ) constitute the system poles, which are often expressed in the Cartesian format as follows:
λ i = σ i + i ω i ,       λ i * = σ i i ω i ,
where σ i and ω i represent the damping factor and damped natural frequencies associated with the i-th eigenpair, respectively. For each generalized eigenpair, the modal vector ψi can be defined by satisfying the following relation:
ϕ i = λ i ψ i ψ i .
Substituting Equation (13) in (11),
λ i 2 M + λ i C + K ψ i = Z λ i ψ i = 0
It was demonstrated that each modal vector ψi, in combination with the associated system pole λi, leads to a general solution to the characteristic equation expressed in Equation (8).
Eigenvalue problem solvers are known to produce an arbitrarily scaled set of eigenpairs (λi, ϕi). However, the resulting modal vector ψi, can be arbitrarily changed by assuming a normalization strategy. For this work, we adopt modal mass normalization, which consists of normalizing the shape vectors with respect to the mass matrix, such that
ψ T M n × n ψ = I m × m ,
where ψ = ψ 1   ψ m n × m , represents a rectangular matrix constructed by considering a finite set of m modal vectors.

2.3. FE Frequency Response Functions

The solutions to the vibration Equation (6) can be expressed [6] as a function of the transfer function matrix, which is obtained as the inversion of the dynamic stiffness introduced in Equation (8):
H p = Z 1 p
Assuming that such an inversion can be performed, the knowledge of H(p) allows the description of the free evolution dynamics of the displacements, X(p), as a linear function of the input forces, F(p), according to
X p = H p   F p
By considering the substitution p = i ω , thus migrating from the complex frequency typical of the Laplace domain to the imaginary frequency of the Fourier domain, we obtain H(), which represents the well-known frequency response function (FRF). The FRF is widely used to characterize structural dynamic behavior, both numerically and experimentally. The FRF can be evaluated as the ratio between the output and input spectra.
H i ω = X i ω F i ω
Spectral analysis techniques can be used to reduce uncertainty and enhance the robustness of FRF estimations in the presence of noise. The most commonly used approaches are referred to as the H1 and H2 estimators, which are evaluated by computing auto- and cross-power spectra [36] according to:
H 1 i ω = X i ω · F * i ω F i ω · F * i ω H 2 i ω = X i ω · X * i ω F i ω · X * i ω
Because of the purely numerical nature of the benchmarks in this study, the use of either form is equivalent.

2.4. Transient Dynamic Using the Newmark Integration Scheme

For a damped system, predicting the dynamic response to an external excitation (Equation (4)) is a significant engineering problem that is posed as a multidimensional initial-value problem:
M   q ¨ k + 1 + C   q ˙ k + 1 + K   q k + 1 = F k + 1 q k = 0 = q 0 q ˙ k = 0 = q 0 ˙ ;
The solution to (20) can be obtained by performing discretization in the time domain. In this study, we adopted a Newmark implicit integration scheme with a fixed time step, Δt, which can be expressed in state-space form as follows:
M + γ Δ t   C γ Δ t   K β Δ t 2 C M + β Δ t 2   K   q ˙ q k + 1 = M 1 γ Δ t   C 1 γ Δ t   K Δ t   M 1 2 β Δ t 2   C M 1 2 β Δ t 2   K q ˙ q k + F γ   F β F γ = γ Δ t F k + 1 + 1 γ Δ t F k F β = β Δ t 2 F k + 1 + 1 2 β Δ t 2 F k
where γ and β are the integration factors that can be tuned to achieve unconditional stability.

3. FE Accuracy Estimation in a Dynamic Context

The literature extensively explores error estimation and mitigation techniques, focusing on error quantification and studying the convergence rate to enhance the mesh quality and achieve an optimal balance between accuracy and computational cost [37]. A preliminary investigation regarding the convergence of FE accuracy in relation to different mesh sizes and orders was recently discussed in [24]. In this study, we extend this research by considering a set of different error metrics to assess the accuracy of FE solutions. The proposed comprehensive FE error analysis comprises an evaluation of the accuracy of the FE modal analysis in relation to the transient dynamic performance.
Several error metrics can be considered for evaluating modal accuracy. Modal parameters can be used to express uncertainty in discrete approximations of linear boundaries and initial-value problems. Modal parameters are directly used within the perturbation method, which is a stochastic FE model-updating procedure that guarantees structural reliability throughout the lifecycle [38]. In this context, a method for estimating the modal parameters together with their covariance matrix was documented in [39]. The Modal Assurance Criterion (MAC) is a common approach for correlating mode shapes in modal analysis. MAC is a statistical indicator that provides a measure of the degree of resemblance between the two modal vectors [11,40]. More exactly, the MAC measures the correlation between two eigenvectors extracted from two different sets of eigenpairs. The MAC is computed as follows:
M A C ψ r ψ h = ψ r T ψ h 2 ψ r T ψ r ψ h T ψ h
where ψr and ψh are the reference and estimated modal vectors, respectively, and T denotes transpose operations. It is worth noting that MAC is not affected by the scaling of the modal vectors. A value of one indicates perfect correlation, whereas zero indicates no resemblance. MAC is a robust metric that is often applied to compare the modal vectors obtained through different experimental campaigns and to contrast the results between the FE predictions and experimental results.
In addition to the MAC, other metrics are also used, such as the L2-norm, which is obtained as the Root Mean Square Error and is typically adopted to minimize deviations in terms of displacements, elastic energy, or stress. As such, the L2-norm allows a comparison of the overall accuracy of the FE solutions. Another metric worth mentioning is the Energy Norm (H1-norm), which can be used to measure errors in gradient-dependent quantities such as stresses and strains. In summary, each metric is targeted to a specific need; different error norms are used depending on the target design goal [41].
After using the MAC to find a robust correspondence between two distinct sets of modal parameters, several error metrics were considered for evaluation: the relative errors measured on the eigenvalues, Eλ; eigenmodes error norm, EΦ, and normalized eigen energy norm, EΠ.
Eλ, is related to the residual errors due to the modal frequency shift and can be assessed according to
E λ = λ h λ r λ r .
EΦ, can be assessed according to the following L2-norm:
E Φ = | | ψ h ψ r | | 2 ,
where the differences between the eigen modes, are considered.
Reali and Hughes recently presented a comprehensive approach for estimating the error of FE approximations, demonstrating the existence of a relationship between the relative errors in the eigenvalues, the eigenmode error norm, and the eigen-energy norm [25,26]. In practice, they discussed the existence of the Pythagorean eigenvalue error theorem, demonstrating that an intimate relationship exists between Eλ, EΦ, and EΠ, as follows:
E Π = | | ψ h ψ r | | E 2 λ r = E λ + E Φ
Therefore, evaluating the normalized eigen energy norm, Eπ, has the interesting benefit of condensing the information of both the eigenvalues relative errors and the eigenvectors accuracy.
Two performance indicators were used in this study for assessing the transient dynamics error analysis. First, the frequency response assurance criterion (FRAC) allows the use of the FRFs as the criterion for model validation, thus increasing the reliability of the obtained results [42]. Specifically, FRAC is formulated as follows:
F R A C H r i ω , H h i ω = H r ω T H h ω 2 H r ω T H r ω H h ω T H h ω
where Hr and Hh denote the reference and approximated FRFs, respectively, and j ω denotes the complex angular frequency.
Second, we used the time response assurance criterion (TRAC) [43], which can be evaluated on two different time-history signals according to:
T R A C u r , u h = u r u h T 2 u r u r T u h u h T
where the u h and u r are the estimated and reference displacements, respectively. The TRAC is evaluated separately for each degree of freedom.

4. Results and Discussion

The numerical tests executed in this study aim to prove and assess the feasibility of using the accuracy measured on top of the FE modal results to predict or even design the accuracy of the same models in terms of transient time-domain dynamic performances. After introducing the numerical case study (Section 4.1), the modal analysis errors are discussed, along with the relevance of using the corresponding error metrics described in Section 3. Subsequently, a set of different case studies is proposed to highlight the relationship between the reported modal accuracy and the performance obtained using canonical transient dynamic metrics such as TRAC and FRAC (Section 4.3). It is worth noting that all the error metrics considered in this study provide non-dimensional results.

4.1. Numerical Case Study

This section describes a numerical case study designed to analyze the accuracy and performance of several C0 continuous hexahedral finite element formulations in the context of FE structural dynamic analysis.
As shown in Figure 1, a one-face-clamped cube was chosen as the numerical benchmark. The formulation of such boundary conditions is trivially realized by fixing all the DOFs related to the nodes belonging to the clamped face. The choice of a cube allowed us to adopt a very simple meshing strategy: the finite element meshes were easily derived by dividing the solid into equal cubic elements. The mechanical material properties of polytetrafluoroethylene (also known as Teflon) are listed in Table 1 and were utilized for the execution of the numerical benchmarks.

4.2. Accuracy of FE Modal Analysis

An h-refinement convergence study of the accuracy of the modal analysis solution was performed considering the linear (Hexa8), quadratic (Hexa27), and cubic (Hexa64) hexahedral element formulations. As discussed in Section 2, only the full Gauss integration schemes were applied; over-integration and under-integration schemes were not considered, based on the findings presented in [24]. This led to the adoption of the 2 × 2 × 2, 3 × 3 × 3, and 4 × 4 × 4 integration schemes for linear, quadratic, and cubic formulations, respectively.
As listed in Table 2, for each formulation, different testing conditions were considered by progressively increasing the number of elements in each of the three directions, thus maintaining a homogeneous discretization. The increasing trend of the three formulations was tuned to maintain the same number of DOFs.
In addition to the numerical tests planned in the numerical benchmark, a very fine cubic-order mesh featuring 18 × 18 × 18 equally spaced cubic hexahedral elements, corresponding to 499,125 DOFs, was adopted as the reference model for all validation scopes.
The testing protocol executed for performing the numerical evaluations for each case consisted of the following steps:
  • first, the mass and stiffness matrices, M and K, were formulated according to Equations (2) and (3);
  • the orthogonal undamped eigenvibration problem was solved using the Krylov-Schur method [44], resulting in a sparse subset of eigenpairs (λi, ψi);
  • MAC was used for disambiguation and matching of the computed modes with respect to those obtained for the reference case;
  • the different error metrics were computed and stored for later evaluation.
In particular, the first two steps were performed for the reference solution in order to identify the ten smallest eigenvalue modes, corresponding to a dynamic range of interest extending slightly above 1000 Hz. The modal frequencies obtained for the reference case are listed in Table 3. The Krylov-Schur method used within this study is part of the MATLAB package, available through the function ‘eigs’. For the scope of this study, every eigenproblem was solved, enforcing the symmetry and setting the convergence tolerance to
As discussed in Section 3, during the third step, the MAC was used for disambiguation purposes. For each testing case, the first 13 lower eigenvalue modes were extracted and compared with the 10 modes of the reference model. The obtained MAC matrix, as depicted in Figure 2a, was used to match the considered modes to the corresponding one in the reference set. Performing the disambiguation and correlation steps is crucial for a correct evaluation of the error metrics based on eigenvector differences, which could otherwise result in reporting excessively higher values in the case of mode switching. In fact, without disambiguation, the different modes could shift in frequency and change their positional order in the sequence, leading to a possibly erroneous comparison and excessively large errors in the modal vector-related metrics. Owing to the symmetries of the structure and double modes, such a disambiguation step could be instrumental in comparing, opportunely, the modes related to very close eigenvalues.
After disambiguation, only the ten matching modes were promoted to the following stages of the numerical experiments. Figure 2b depicts the MAC matrix resulting from the aforementioned disambiguation step.
For each testing case listed in Table 2, three error metrics were evaluated, according to the formulations discussed in Section 3: the relative errors measured on the eigenvalues, Eλ; the eigenmode error norm, EΦ, and the normalized eigen energy norm, Eπ. Figure 3 depicts the convergence of the maximum eigenvalue energy norm error, that is, the maximum value obtained using Equation (21), among the ten modes considered. It is worth noting that all the errors reported for different formulations are adimensional.
As expected, the obtained convergence rates indicate a favorable adoption of higher-order FE formulations. Using cubic Hexa64-Int64 and quadratic Hexa27-Int27 allowed us to obtain significantly smaller models while maintaining similar levels of accuracy.
Similarly, Figure 4 reports the h-convergence rates of the relative errors measured on the eigenvalues, Eλ, and of the eigenmode error norm, EΦ.
The findings presented in the results align with those obtained using the eigen-energy norm, suggesting that there is no discernible benefit in favoring one metric over the others. Except for a constant scaling factor, all metrics appeared to be equally effective in evaluating the relative importance of the various formulations.
Figure 5 details the relative prominence of each mode when the eigen-energy metric is used. The findings indicate the noteworthy stability of the proposed metric across the entire frequency range tested, demonstrating a consistent performance for all formulations. The representation of the relative contribution of each mode displays a characteristic profile on the error-mode plane, which remains unchanged even when the number of DOFs increases. This suggests that the relative significance of the modes does not depend on the resolution of the model or chosen formulation.
Table 4 details all the shape vectors and the related vector norm errors for the Hexa8-Int8 42 × 42 × 42 case. The scaling of all the modes was manually tuned for visualization purposes.

4.3. Accuracy of FE Transient Dynamic Analysis

As shown in Figure 6, three separate comparative case studies were conducted to assess the transient dynamic errors, which were measured in terms of the FRAC and TRAC. Cases 1 and 2 were chosen for two distinct levels of the error metric; thus, they were expected to lead to comparable levels of dynamic accuracy, regardless of the employed formulation. Instead, Case 3 was designed to compare the performances of different formulations while maintaining the same number of degrees of freedom.
Table 5 presents an overview of the six transient dynamic tests conducted to finalize the examination of the aforementioned case studies. In accordance with Equation (16), the initial value problem was formulated for each numerical test and then solved using the implicit beta Newmark integration scheme, as specified in Equation (17), and implemented by means of Matlab programming.
An impulse force is applied to one of the free vertices of the non-clamped side of the cube, as shown in Figure 7a. The displacement responses of the other vertices (labeled as A00, AB0, and A0C in Figure 7a) were recorded. The integration of all simulations was performed at a fixed rate of Δ t = 0.0004   s , for a total duration of 2 s. The peak value of the impulse was set to 4.6 N. The impulse was distributed over multiple time steps, resulting in a smoothly decreasing frequency of excitation, as shown in Figure 7c.
To perform the damped dynamic analysis, proportional damping was computed according to the Rayleigh damping model, C = α   K + μ   M , using coefficients α = 1.9 × 10 7 and μ = 6.33 . The damping tuning was performed according to the method proposed in [45] by setting a damping ratio, ζ , equal to 0.21 and 0.11 at the extremes of the frequency range (200–1250) Hz, which suffices to span over the dynamic content of the first ten eigenmodes (Table 3). Appendix A.1 provides details of the methodology adopted for the damping parameter determination.
Table 6 summarizes the results of the transient dynamic analysis, reported in terms of the FRAC and TRAC errors, in comparison with the eigen-energy error metric, Eπ, obtained for the modal analysis. Irrespective of the model order and number of DOFs, a clear and strong correlation was observed between the reported eigen-energy errors and the obtained levels of TRAC and FRAC. Details of each comparative case study are presented and discussed in the following subsections.

4.3.1. Case Study 1: EΠ ≈ 0.01

The displacement responses to the considered impulse loading scenario are shown in Figure 8 for the A00 vertex and for the three different FE models.
Figure 9 shows the displacement FRFs obtained for the A00 vertex, whereas the FRAC metric is shown in Figure 10.
The displacement responses and the corresponding FRFs were also measured for the other control points, labeled AB0 and A0C, and used to evaluate the corresponding TRAC and FRAC metrics, as reported in Table 6.

4.3.2. Case Study 2: EΠ ≈ 0.003

Similarly to case study 1, this case compares the performances of two FE models achieving a similar accuracy, as measured by means of the normalized eigen-energy norm, E Π , of approximately 0.003. As expected, the documented improvement in terms of modal errors was also reflected in the dynamic transient performance. Figure 11, Figure 12 and Figure 13 depict the time-displacement responses, FRFs, and FRACs, respectively.
It is important to recognize that the quadratic model demonstrated slightly greater accuracy than the cubic model when evaluated using the eigen energy metric, with a score of 0.0034 compared to 0.0037. Similarly, marginally better performance was also achieved in terms of FRACs and TRAcs. Both models achieved a good level of accuracy, which was high enough to produce a dynamic response qualitatively equivalent to the reference when looking at the time-domain displacement plots (Figure 11) or at the corresponding FRF results (Figure 12). The marginal difference can be assessed by looking at the FRAC(f) plot reported in Figure 13, where the axis scaling was tuned to maximize the differences, or by looking at the values reported in Table 6.

4.3.3. Case Study 3: Equally Dense Mesh

For case study 3, the goal was to compare different FE models formulated using increasing orders of the interpolants while maintaining an equal number of DOFs.
The different levels of accuracy provided by the different formulations are shown in Figure 14, where the FRAC levels are reported along the spectrum. It is worth noting how well the cubic formulation performed, even in the higher frequency range, as opposed to the linear elements, whose accuracy deteriorated irremediably.

4.4. Computational Efficiency of Higher-Order FE Formulations

The computational efficiency of the FE formulations is reported in this section, detailing the time required to assemble the matrices, solve the eigenvalue problem, and perform the time-domain integration. Figure 15 shows the computational performance in terms of the assembly time of the linear, quadratic, and cubic formulations.
For the assembly stage, the quadratic elements were found to be slightly more computationally efficient than the cubic formulation, offering the best Pareto-front solution. It is worth noting how the computational efficiency of the cubic formulation tends towards the slightly better results achieved by the quadratic one with an increasing number of DOFs. In particular, considering the performance corresponding to case study 1, the higher-order formulations Hexa27-Int27 and Hexa64-Int64 outperformed the linear Hexa8-Int8. While achieving a similar target error of about EΠ ≈ 0.01, the linear formulation required nearly 11.5 s, whereas the quadratic and cubic formulations required only 2.3 and 3.06 s, respectively. Regarding case study 2, the cubic element Hexa64-Int64 performed faster than the quadratic Hexa27-Int27, completing the task in 14 s instead of 16 s.
Considering case study 3, the linear and quadratic formulations reported similar performances, which is expected due to the equivalent number of DOfs to be considered; thus, size of the produced matrices is equal. The worse performance of the cubic formulation assembly stage seems to suggest some computational problems arising for the specific type, which may be due to several possible reasons. For example, a negative impact can be caused by the fact that the size of the element structures used to locally accumulate results in overflowing the local cache dimension, thus causing some memory overflow. Considering that the assembly stage is marginal compared to the time required for solving modal and transient dynamic problems, the solution to the documented issues was left for future investigations.
The computational efficiency of the compared formulations for solving eigenvalue problems is depicted in Figure 16. In general, the cubic formulation achieved more accurate results faster, dominating the optimal Pareto front. The results of case study 3 are emblematic: even though the formulated problems count the same number of DOFs, the linear FE required significantly more time than the other formulations. It took approximately 2000 s to achieve a mild accuracy level, while quadratic and cubic meshes took approximately 750 s, with the cubic FE model outperforming the quadratic in terms of accuracy. This could be explained by the fact that the increased accuracy of the higher-order formulations may positively affect the convergence of the iterative Krylov-Schur method [44] used to solve the sparse eigen problem.
As confirmed by the numerical results summarized in Table 7, higher-order finite elements, specifically Hexa27-Int27 and Hexa64-Int64, exhibit a greater aptitude for accurately capturing the modal behavior of structures than linear formulations. Within the examined frequency range, linear finite element models were found to significantly underperform in terms of accuracy and efficiency. The relatively lower TRAC values reported for the linear FE formulations correspond to relatively larger inaccuracies in the displacement response, suggesting that a target error of approximately 1% is insufficient for obtaining an accurate dynamic response. Moreover, the linear model exhibited a substantial decline in FRAC performance, especially at higher frequencies, whereas quadratic and cubic elements maintained superior accuracy, not influenced by the increasing frequencies.

5. Conclusions

This study thoroughly examines the efficacy of higher-order Finite Element (FE) methods in addressing dynamic structural problems, with a particular emphasis on both modal and transient dynamic analyses. The main objective was to evaluate the predictive accuracy and computational efficiency of various FE formulations. The results emphasize the value of the normalized eigen-energy (NEE) metric, which integrates errors resulting from eigenfrequency shifts and eigenshape aberrations, as a superior predictor of the FE model’s accuracy performance. The NEE metric demonstrated a strong relationship with the Time-Response Assurance Criterion (TRAC) and the Frequency-Response Assurance Criterion (FRAC), both of which are typically more computationally intensive. This relationship indicates that the NEE metric could be a reliable tool for fine-tuning FE models to achieve high-fidelity predictions, which are essential for contemporary industrial applications such as virtual sensing and digital twins.
In addition, the numerical investigation meticulously evaluated the computation time and accuracy of a 3D numerical model using three types of hexahedral meshes: cubic, quadratic, and linear. Our findings confirm that higher-order FE models can significantly decrease computational efforts while maintaining higher levels of accuracy. Quadratic and cubic elements not only achieve the desired error levels more rapidly but also demonstrate significantly reduced computation times compared to linear models. This efficiency makes higher-order formulations more favorable for applications that require precise dynamic response predictions.
Based on the reported findings, further investigation is motivated to improve the computational efficiency of the higher-order formulations, which may benefit larger speed-ups from properly using massive parallel accelerators, such as GPUs.

Author Contributions

Conceptualization, A.K. and F.C.; methodology, A.K., F.C. and D.M.; software, F.C.; validation, A.K., F.C. and D.M.; formal analysis, A.K., F.C. and D.M.; investigation, A.K., F.C. and D.M.; resources, D.M.; data curation, A.K. and F.C.; writing—original draft preparation, A.K.; writing—review and editing, F.C. and D.M.; visualization, A.K. and F.C.; supervision, F.C. and D.M.; project administration, F.C.; funding acquisition, D.M. All authors have read and agreed to the published version of the manuscript.

Funding

This research work was carried out in the context of the project “Centro Nazionale di Ricerca in High-Performance Computing, Big Data and Quantum Computing”, Project Number CN00000013, CUP H23C22000360005, and of the National Recovery. The views and opinions ex-pressed are only those of the authors and do not necessarily reflect those of the European Union or the European Commission. Neither the European Union nor the European Commission can be held responsible for them.

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 conflicts of interest.

Appendix A

Appendix A.1. Damped Systems

Real-life systems (Equation (4)) have energy dissipation mechanisms, that consequently, are not undamped. A linear system cannot be decoupled by modal analysis in the presence of damping unless it has all classical normal modes [46]. A common assumption made when applying modal analysis to dampened systems is the application of proportional viscous damping. This model expresses the damping matrix as a linear combination of the mass and stiffness matrices, as follows [45]:
C =   α K + μ M .
By transforming the damping into the modal-coordinate space, the damping ratio ζ can be calculated as
ζ = 1 2 α ω + μ ω ,
where ω is the circular frequency.
The coefficients α and μ can be determined once the damping ratios ζi and ζj are specified for the corresponding ω i and ω j frequencies [47], and by solving the following linear system of equations:
1 2 1 / ω i     1 / ω j     ω i       ω j μ α = ζ i ζ j

References

  1. Okereke, M.; Keates, S. Finite Element Mesh Generation. In Finite Element Applications; Springer International Publishing: Berlin/Heidelberg, Germany, 2018; pp. 165–186. [Google Scholar]
  2. Fish, J.; Belytschko, T. A First Course in Finite Elements; Wiley: Hoboken, NJ, USA, 2007. [Google Scholar]
  3. Logan, D. First Course in the Finite Element Method, 4th ed.; Thomson: Nelson, Toronto, ON, Canada, 2007; ISBN 0-534-55298-6. [Google Scholar]
  4. Cook, R.D.; Malkus, D.S.; Plhesha, M.E.; Witt, R.J. Concepts and Applications of Finite Element Analysis, 4th ed.; Wiley: Hoboken, NJ, USA, 2001; ISBN 978-0-471-35605-9. [Google Scholar]
  5. Yong, Y.; Lin, Y.K. Dynamic Response Analysis of Truss-Type Structural Networks: A Wave Propagation Approach. J. Sound Vib. 1992, 156, 27–45. [Google Scholar] [CrossRef]
  6. Heylen, W.; Lammens, S.; Sas, P. Modal Analysis Theory and Testing, 4th ed.; Katholieke Universiteit Leuven: Leuven, Belgium, 1998; ISBN 907380261X. [Google Scholar]
  7. Maressa, A.; Mundo, D.; Donders, S.; Desmet, W. A Wave-Based Substructuring Approach for Concept Modeling of Vehicle Joints. Comput. Struct. 2011, 89, 2369–2376. [Google Scholar] [CrossRef]
  8. Picard, C.; Frisson, C.; Faure, F.; Drettakis, G.; Kry, P.G. Advances in Modal Analysis Using a Robust and Multiscale Method. EURASIP J. Adv. Signal Process 2010, 2010, 392782. [Google Scholar] [CrossRef]
  9. Fu, Z.; He, J. Modal Analysis, 1st ed.; Butterworth-Heinemann: Oxford, UK, 2001; ISBN 0 7506 5079 6. [Google Scholar]
  10. Wang, L.; Zhong, H. A Time Finite Element Method for Structural Dynamics. Appl. Math. Model. 2017, 41, 445–461. [Google Scholar] [CrossRef]
  11. Mercer, J.F.; Aglietti, G.S.; Kiley, A.M. Modal and Frequency Domain Based Techniques for Finite Element Model Correlation. In Proceedings of the COMPDYN 2015 5th ECCOMAS Thematic Conference on Computational Methods in Structural Dynamics and Earthquake Engineering, Crete Island, Greece, 25–27 May 2015; pp. 191–208. [Google Scholar]
  12. El-Kafafy, M.; Peeters, B.; Geluk, T.; Guillaume, P. The MLMM Modal Parameter Estimation Method: A New Feature to Maximize Modal Model Robustness. Mech. Syst. Signal Process 2019, 120, 465–485. [Google Scholar] [CrossRef]
  13. Chandravanshi, M.L.; Mukhopadhyay, A.K. Modal Analysis of Structural Vibration. In Proceedings of the ASME International Mechanical Engineering Congress and Exposition, Volume 14: Vibration, Acoustics and Wave Propagation, San Diego, CA, USA, 15–21 November 2013. [Google Scholar]
  14. Franzoni, F.; Ferreira, C.A.E.; Weaver, P. Numerical and Experimental Dynamic Analyses of a Postbuckled Box Beam. AIAA J. 2016, 54, 1987–2003. [Google Scholar] [CrossRef]
  15. Yang, B.; Zhang, Y. A New Approach to Transient Vibration Analysis of Two-Dimensional Beam Structures at Medium and High Frequencies. J. Comput. Nonlinear Dyn. 2020, 15, 091006. [Google Scholar] [CrossRef]
  16. Capalbo, C.E.; De Gregoriis, D.; Tamarozzi, T.; Devriendt, H.; Naets, F.; Carbone, G.; Mundo, D. Parameter, Input and State Estimation for Linear Structural Dynamics Using Parametric Model Order Reduction and Augmented Kalman Filtering. Mech. Syst. Signal Process 2023, 185, 109799. [Google Scholar] [CrossRef]
  17. Ewins, D.J. Model Validation: Correlation for Updating. Sadhana—Acad. Proc. Eng. Sci. 2000, 25, 221–234. [Google Scholar] [CrossRef]
  18. Ereiz, S.; Duvnjak, I.; Fernando Jiménez-Alonso, J. Review of Finite Element Model Updating Methods for Structural Applications. Structures 2022, 41, 684–723. [Google Scholar] [CrossRef]
  19. Baghdasaryan, L.; Chen, W.; Buranathiti, T.; Cao, J. Model Validation via Uncertainty Propagation Using Response Surface Models. In Proceedings of the ASME 2002 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference. Volume 2: 28th Design Automation Conference, Montreal, QC, Canada, 29 September–2 October 2002; pp. 981–992. [Google Scholar]
  20. Zang, C.; Schwingshackl, C.W.; Ewins, D.J. Model Validation for Structural Dynamic Analysis: An Approach to the Sandia Structural Dynamics Challenge. Comput. Methods Appl. Mech. Eng. 2008, 197, 2645–2659. [Google Scholar] [CrossRef]
  21. Zienkiewicz, O.C.; Zhu, J.Z. A Simple Error Estimator and Adaptive Procedure for Practical Engineerng Analysis. Int. J. Numer. Methods Eng. 1987, 24, 337–357. [Google Scholar] [CrossRef]
  22. Zienkiewicz, O.C.; Taylor, R.L.; Zhu, J.Z. The Finite Element Method: Its Basis and Fundamentals; Butterworth-Heinemann: Oxford, UK, 2005. [Google Scholar]
  23. Babuška, I. Advances in the p and H-p Versions of the Finite Element Method. A Survey; Springer: Berlin/Heidelberg, Germany, 1988; pp. 31–46. [Google Scholar] [CrossRef]
  24. Karpik, A.; Cosco, F.; Mundo, D. Higher-Order Hexahedral Finite Elements for Structural Dynamics: A Comparative Review. Machines 2023, 11, 326. [Google Scholar] [CrossRef]
  25. Hughes, T.; Evans, J.; Reali, A. Finite Element and NURBS Approximations of Eigenvalue, Boundary-Value, and Initial-Value Problems. Comput. Methods Appl. Mech. Eng. 2014, 272, 290–320. [Google Scholar] [CrossRef]
  26. Bartoň, M.; Calo, V.; Deng, Q.; Puzyrev, V. Generalization of the Pythagorean Eigenvalue Error Theorem and Its Application to Isogeometric Analysis. In Numerical Methods for PDEs. SEMA SIMAI Springer Series, Vol 15; Springer: Cham, Switzerland, 2018; pp. 147–170. [Google Scholar]
  27. Anderson, R.; Andrej, J.; Barker, A.; Bramwell, J.; Camier, J.-S.; Cerveny, J.; Dobrev, V.; Dudouit, Y.; Fisher, A.; Kolev, T.; et al. A Modular Finite Element Methods Library. Comput. Math. Appl. 2021, 81, 42–74. [Google Scholar] [CrossRef]
  28. Arndt, D.; Bangerth, W.; Davydov, D.; Heister, T. The Deal. II Finite Element Library: Design, Features, and Insights. Comput. Math. Appl. 2021, 81, 407–422. [Google Scholar] [CrossRef]
  29. Neto, M.; Amaro, A.; Roseiro, L.; Cirne, J.; Leal, R. Engineering Computation of Structures: The Finite Element Method; Springer: Berlin/Heidelberg, Germany, 2015. [Google Scholar]
  30. Bathe, K. Finite Element Procedures, 2nd ed.; Bathe, K.J., Ed.; Massachusetts Institute of Technology: Watertown, MA, USA, 2014; ISBN 0979004950. [Google Scholar]
  31. Abedian, A.; Parvizian, J.; Düster, A.; Khademyzadeh, H.; Rank, E. Performance of Different Integration Schemes in Facing Discontinuities in the Finite Cell Method. Int. J. Comput. Methods 2013, 10, 1350002. [Google Scholar] [CrossRef]
  32. Mousavi, S.; Sukumar, N. Generalized Gaussian Quadrature Rules for Discontinuities and Crack Singularities in the Extended Finite Element Method. Comput. Methods Appl. Mech. Eng. 2010, 199, 3237–3249. [Google Scholar] [CrossRef]
  33. Venkateshan, S.; Swaminathan, P. Computational Methods in Engineering, 1st ed.; Elsevier: Amsterdam, The Netherlands, 2013; ISBN 0124167020. [Google Scholar]
  34. Wang, X.; Yuan, Z.; Jin, C. Weak Form Quadrature Element Method and Its Applications in Science and Engineering: A State-of-the-Art Review. Appl. Mech. Rev. 2017, 69, 030801. [Google Scholar] [CrossRef]
  35. Zienkiewicz, O.C.; Taylor, R.L.; Too, J.M. Reduced Integration Technique in General Analysis of Plates and Shells. Int. J. Numer. Methods Eng. 1971, 3, 275–290. [Google Scholar] [CrossRef]
  36. Mao, Z.; Todd, M. Statistical Modeling of Frequency Response Function Estimation for Uncertainty Quantification. Mech. Syst. Signal Process 2013, 38, 333–345. [Google Scholar] [CrossRef]
  37. Banz, L.; Hernández, O.; Stephan, E.P. A Priori and a Posteriori Error Estimates for Hp-FEM for a Bingham Type Variational Inequality of the Second Kind. Comput. Math. Appl. 2022, 126, 14–30. [Google Scholar] [CrossRef]
  38. Hua, X.; Wen, Q.; Ni, Y.; Chen, Z. Assessment of Stochastically Updated Finite Element Models Using Reliability Indicator. Mech. Syst. Signal Process 2017, 82, 217–229. [Google Scholar] [CrossRef]
  39. Pintelon, R.; Guillaume, P.; Schoukens, J. Uncertainty Calculation in (Operational) Modal Analysis. Mech. Syst. Signal Process 2007, 21, 2359–2373. [Google Scholar] [CrossRef]
  40. Brigante, D.; Rainieri, C.; Fabbrocino, G. The Role of the Modal Assurance Criterion in the Interpretation and Validation of Models for Seismic Analysis of Architectural Complexes. In Proceedings of the X International Conference on Structural Dynamics, EURODYN 2017, Rome, Italy, 10–13 September 2017; Volume 199, pp. 3404–3409. [Google Scholar]
  41. Talebi, H.; Zi, G.; Silani, M.; Samaniego, E.; Rabczuk, T. A Simple Circular Cell Method for Multilevel Finite Element Analysis. J. Appl. Math. 2012, 2012, 526846. [Google Scholar] [CrossRef]
  42. Manring, L.; Schultze, J.; Zimmerman, S. Improving Magnitude and Phase Comparison Metrics for Frequency Response Functions Using Cross-Correlation and Log-Frequency Shifting. J. Sound Vib. 2022, 539, 117255. [Google Scholar] [CrossRef]
  43. Iliopoulos, A.; Shirzadeh, R.; Weijtjens, R.; Guillaume, P.; Van Hemelrijck, D.; Devriendt, C. A Modal Decomposition and Expansion Approach for Prediction of Dynamic Responses on a Monopile Offshore Wind Turbine Using a Limited Number of Vibration. Mech. Syst. Signal Process 2016, 68–69, 84–104. [Google Scholar] [CrossRef]
  44. Kressner, D. Numerical Methods for General and Structured Eigenvalue; Springer: Berlin/Heidelberg, Germany, 2005. [Google Scholar]
  45. Adhikari, S.; Phani, A.S. Rayleigh’s Classical Damping Revisited. In Proceedings of the International Conference on Civil Engineering in the New Millennium: Opportunities and Challenges, Shibpur, India, 11–14 January 2007. [Google Scholar]
  46. Ma, F.; Imam, A.; Morzfeld, M. The Decoupling of Damped Linear Systems in Oscillatory Free Vibration. J. Sound Vib. 2009, 324, 408–428. [Google Scholar] [CrossRef]
  47. Chopra, A.K. Dynamics of Structures: Theory and Applications to Earthquake Engineering, 4th ed.; Pearson: London, UK, 2007. [Google Scholar]
Figure 1. 3D clamped cube model (a) and the corresponding FE reference mesh (b).
Figure 1. 3D clamped cube model (a) and the corresponding FE reference mesh (b).
Designs 08 00079 g001
Figure 2. Disambiguation step performed using MAC: original MAC matrix (a) and updated MAC matrix (b).
Figure 2. Disambiguation step performed using MAC: original MAC matrix (a) and updated MAC matrix (b).
Designs 08 00079 g002
Figure 3. Convergence rate of linear, quadratic, and cubic FE formulations as a function of the number of DOFs.
Figure 3. Convergence rate of linear, quadratic, and cubic FE formulations as a function of the number of DOFs.
Designs 08 00079 g003
Figure 4. Convergence rate of the eigenvalue Eλ and eigenvector EΦ errors for linear, quadratic, and cubic FE formulations as a function of the number of DOFs.
Figure 4. Convergence rate of the eigenvalue Eλ and eigenvector EΦ errors for linear, quadratic, and cubic FE formulations as a function of the number of DOFs.
Designs 08 00079 g004
Figure 5. Comprehensive convergence analysis as a function of the number of DOFs and modes.
Figure 5. Comprehensive convergence analysis as a function of the number of DOFs and modes.
Designs 08 00079 g005
Figure 6. Different case studies designed on top of the convergence rate obtained for the modal analysis results using linear, quadratic, and cubic FE formulations.
Figure 6. Different case studies designed on top of the convergence rate obtained for the modal analysis results using linear, quadratic, and cubic FE formulations.
Designs 08 00079 g006
Figure 7. 3D clamped cube model with applied force and control points (a) external impulse force represented in time domain (b) and frequency domain (c).
Figure 7. 3D clamped cube model with applied force and control points (a) external impulse force represented in time domain (b) and frequency domain (c).
Designs 08 00079 g007
Figure 8. Case study 1: displacement responses for point A00.
Figure 8. Case study 1: displacement responses for point A00.
Designs 08 00079 g008
Figure 9. Case study 1: displacement FRFs at the control point A00.
Figure 9. Case study 1: displacement FRFs at the control point A00.
Designs 08 00079 g009
Figure 10. Case study 1: FRAC(f) across frequency for the different FE models.
Figure 10. Case study 1: FRAC(f) across frequency for the different FE models.
Designs 08 00079 g010
Figure 11. Case study 2: Displacement responses for the point A00.
Figure 11. Case study 2: Displacement responses for the point A00.
Designs 08 00079 g011
Figure 12. Case study 2: Displacement FRFs at control point A00.
Figure 12. Case study 2: Displacement FRFs at control point A00.
Designs 08 00079 g012
Figure 13. Case study 2: FRAC(f) across frequency for the different FE models.
Figure 13. Case study 2: FRAC(f) across frequency for the different FE models.
Designs 08 00079 g013
Figure 14. FRAC(f) for all estimated frequency ranges across the entire model, Case 3.
Figure 14. FRAC(f) for all estimated frequency ranges across the entire model, Case 3.
Designs 08 00079 g014
Figure 15. Accuracy of linear, quadratic, and cubic FE formulations as a function of the assembly time.
Figure 15. Accuracy of linear, quadratic, and cubic FE formulations as a function of the assembly time.
Designs 08 00079 g015
Figure 16. Convergence analysis of linear, quadratic, and cubic FE formulations of the time for eigenvalue computations.
Figure 16. Convergence analysis of linear, quadratic, and cubic FE formulations of the time for eigenvalue computations.
Designs 08 00079 g016
Table 1. Mechanical properties of the material.
Table 1. Mechanical properties of the material.
Young’s Modulus E, KPaDensity ρ, kg/mm3Poisson’s Ratio ν
500,0002.13 × 10−60.46
Table 2. Design of the numerical experiment for the FE modal analysis accuracy.
Table 2. Design of the numerical experiment for the FE modal analysis accuracy.
Test NameShape TypeIntegration SchemeN. of Elements along the Edges
Hexa8-Int8Linear2 × 2 × 26121824303642
Hexa27-Int27Quadratic3 × 3 × 336912151821
Hexa64-Int64Cubic4 × 4 × 42468101214
N. of DOFs1029659120,57746,87589,373151,959238,521
Table 3. Natural frequencies obtained for the reference FE model.
Table 3. Natural frequencies obtained for the reference FE model.
Mode12345678910
ω [Hz]257.9258.0331.2630.6662.4662.4803.3943.6979.61020
Table 4. Mode shapes error of estimated model Hexa8-Int8 42 × 42 × 42 with respect to the reference model.
Table 4. Mode shapes error of estimated model Hexa8-Int8 42 × 42 × 42 with respect to the reference model.
Designs 08 00079 i001Designs 08 00079 i002Designs 08 00079 i003Designs 08 00079 i004Designs 08 00079 i005
Mode 1
257.94, Hz
Mode 2
257.99, Hz
Mode 3
331.19, Hz
Mode 4
630.60, Hz
Mode 5
662.35, Hz
Designs 08 00079 i006Designs 08 00079 i007Designs 08 00079 i008Designs 08 00079 i009Designs 08 00079 i010
Mode 6
662.36, Hz
Mode 7
803.32, Hz
Mode 8
943.59, Hz
Mode 9
979.63, Hz
Mode 10
1020, Hz
Table 5. Design of experiments for the transient dynamic case studies.
Table 5. Design of experiments for the transient dynamic case studies.
Test LabelShape
Order
Integration
Scheme
Elements
along the Edge
Free
DOFs
Case Study
123
Hexa64-Int64Cubic4 × 4 × 46 × 6 × 619,494X
Hexa27-Int27Quadratic3 × 3 × 312 × 12 × 1245,000X
Hexa64-Int64Cubic4 × 4 × 411 × 11 × 1186,490 X
Hexa8-Int8Linear2 × 2 × 242 × 42 × 42232,974X X
Hexa27-Int27Quadratic3 × 3 × 321 × 21 × 21232,974 XX
Hexa64-Int64Cubic4 × 4 × 414 × 14 × 14232,974 X
Table 6. Accuracy metrics for the transient dynamics case studies.
Table 6. Accuracy metrics for the transient dynamics case studies.
Test LabelCaseEπControl Points
AB0A0CA00
123TRACFRACTRACFRACTRACFRAC
Hexa64-Int64x 0.00970.510.850.510.850.520.86
Hexa27-Int27x 0.00970.500.850.500.850.510.87
Hexa8-Int8x x0.00990.510.850.510.860.520.87
Hexa64-Int64 x 0.00370.880.980.880.980.880.98
Hexa27-Int27 xx0.00340.900.980.900.980.900.98
Hexa64-Int64 x0.00150.981.000.981.000.981.00
Table 7. Computational efficiency for the transient dynamics.
Table 7. Computational efficiency for the transient dynamics.
Test LabelCaseEπAssembling
Time [s]
Eigenvalue
Time [s]
Dynamic
Time [s]
123
Hexa64-Int64x 0.00973.062.8596.42
Hexa27-Int27x 0.00972.3115.36257.20
Hexa8-Int8x x0.009911.451992.353187.15
Hexa64-Int64 x 0.003714.3740.23718.30
Hexa27-Int27 xx0.003416.18777.843015.62
Hexa64-Int64 x0.001546.15727.472576.63
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.

Share and Cite

MDPI and ACS Style

Karpik, A.; Cosco, F.; Mundo, D. Dynamic Error Estimation in Higher-Order Finite Elements. Designs 2024, 8, 79. https://doi.org/10.3390/designs8040079

AMA Style

Karpik A, Cosco F, Mundo D. Dynamic Error Estimation in Higher-Order Finite Elements. Designs. 2024; 8(4):79. https://doi.org/10.3390/designs8040079

Chicago/Turabian Style

Karpik, Anna, Francesco Cosco, and Domenico Mundo. 2024. "Dynamic Error Estimation in Higher-Order Finite Elements" Designs 8, no. 4: 79. https://doi.org/10.3390/designs8040079

APA Style

Karpik, A., Cosco, F., & Mundo, D. (2024). Dynamic Error Estimation in Higher-Order Finite Elements. Designs, 8(4), 79. https://doi.org/10.3390/designs8040079

Article Metrics

Back to TopTop