Numerical Analysis of the Main Wave Propagation Characteristics in a Steel-CFRP Laminate Including Model Order Reduction

: Guided ultrasonic waves are suitable for use in the context of structural health monitoring of thin-walled, plate-like structures. Hence, observing the wave propagation in the plates can provide an indication of whether damage has occurred in the structure. In this work, the wave propagation in ﬁber metal laminate consisting of thin steel foils and layers of carbon ﬁber-reinforced polymer is studied, focusing on the main propagation characteristics like dispersion diagrams and displacement ﬁelds. For this purpose, the dispersion diagrams derived from the analytical framework and numerical simulations are ﬁrst determined and compared to each other. Next, the displacement ﬁelds are computed using the global matrix method for two excitation frequencies. The results derived from the analytical framework is used to validate the numerically determined displacement ﬁelds based on a 2D and a 3D modeling approach. For both investigations the results of the analytical treatment and the numerical simulation show good agreement. Furthermore, the displacement ﬁeld reveals the typical and well-known characteristics of the propagation of guided waves in thin-walled structures. Since the use of full 3D models involves a very high computational cost, this work also successfully investigates the possibility for model order reduction to decrease the computational time and costs of the simulation without the loss of accuracy.


Introduction
It is only natural that the use of structures made of different materials is automatically associated with the fact that the defects in the material occur at different life cycles of these structures. Furthermore, not every defect is easily detectable by means of visual inspection. This is especially the case when dealing with composite materials and layered structures, which are increasingly being used in the aeronautical industry. One recent approach in the development of composite materials are fiber metal laminates (FML). They are of great interest because they combine the advantages of high ductility, often found in metals, with the high specific stiffness of fiber-reinforced polymers (FRP) [1,2]. However, this comes with the high risk of low-velocity impact induced internal damage like delaminations [3,4]. Therefore, it is of great importance to detect such hidden defects to ensure a safe operation of the component.
One technique here is the use of guided ultrasonic wave (GUW) based structural health monitoring (SHM) [5,6], which has been profoundly analyzed over the last decades for isotropic materials as well as thin-walled structures made from FRP. In both cases the stiffness and density are constant over the thickness. The complex multi-modal nature of Lamb waves, first described by Horace Lamb [7], derived by the solution of the analytical framework is discussed in Refs. [8][9][10][11] for isotropic materials. The adaption to single transversely-isotropic layers, orthotropic layers, and layered structures is provided in Refs. [12][13][14][15]. Since GUW propagation is not only of multi-modal nature but also strongly dispersive, the phase and group velocities are a function of the frequency [6]. A common way to depict the propagation characteristics are dispersion diagrams.
Besides their dispersive nature, the displacement fields are of great interest, because they are crucial for the interaction of the propagating wave with damage hidden in the structure. The displacement fields of different propagating modes are discussed in Refs. [5,11] for isotropic materials and [16,17] for wave guides made of FRP. It is shown that the constant material properties over the thickness lead to a displacement field covering the whole thickness of the component without a phase shift between the upper and lower surface. In contrast to these findings, the wave propagation in sandwich structures behaves significantly different [18]. It is shown that the impedance differences between the layers in sandwich structures lead to frequency-dependent wave propagation phenomena. Until a certain wavelength to thickness ratio is reached, the so called global Lamb waves propagate acting like the well-known Lamb wave propagation by covering the whole thickness without a phase shift [19,20]. This changes when the frequency is increased and hence, the wave length to thickness ratio passes a threshold value. In this case one speaks of true modes [21,22] propagating in the skin layers only and leaky Lamb waves [21]. In this case, large attenuation occurs due to waves propagating in the thickness direction, causing an energy dissipation into the core. It is shown in Ref. [23], that the leakage effect is strongly affected by geometry and acoustic impedance as well as wave mode and frequency. For more details about the leakage effects, the reader is kindly referred to [20,[22][23][24].
Due to the high impedance changes between the metal and FRP perpendicular to the fiber orientation, this might also be an issue for the wave propagation in FML. First investigations reveal that in structures formed of aluminum and glass fiber-reinforced polymer (GFRP) layers, the wave propagation meets the framework of GUW in thin-walled structures [25][26][27][28]. However, no further material combination is analyzed. In this work, the approach is extended to FML consisting of carbon fiber-reinforced polymers (CFRP) paired with thin steel foils. In contrast to aluminum-GFRP, which has a stiffness ratio of approximately 15, steel CFRP exhibits a stiffness ratio of up to 20 with respect to the Young's modulus. Based on the shear modulus, the ratios are 10 and 15, respectively. Therefore, the first objective of this work presented here is the analysis of the wave propagation characteristics in FML consisting of steel and CFRP layers by solving the analytical framework with respect to the well-known boundary conditions for the GUW propagation in thin-walled structures and numerical simulations. Subsequently, dispersion diagrams and the displacement field are obtained based on the analytical framework and numerical simulations before comparing the results.
The numerical models are used to calculate displacement fields in the wave propagation direction of the GUW's and in the thickness direction. Interlayer boundary conditions were considered by the use of the global matrix method [17,29]. The method used here for numerical determination of dispersion diagrams has been successfully applied in the evaluation of experimental data [30,31]. It is based on the evaluation of displacement fields along the wave propagation direction using a discrete 2D Fourier transform [32] and specially tuned excitation signals. Similar evaluations can also be found in Refs. [33][34][35].
To profoundly understand the propagation of GUW in the FML as well as to quantify the uncertainties of the material properties, it is required to solve the forward model several times with input parameter variations. As the multiple queries of the underlying high-fidelity system is associated with high computational cost, it is of great interest to investigate whether an order reduction of the numerical models can significantly reduce the computational effort. Model order reduction (MOR) seeks to swiftly capture the essential features of a higher order complex dynamical system by approximating it to a lower dimension without losing the accuracy of the system response. Depending on the nature of the system, reduced order modeling can be accomplished by any of the existing techniques that are classified into Krylov and singular value decomposition (SVD) methods [36]. The former methods are based on iterative approaches while the latter depends on the idea of expanding the highdimensional solution into a sequence of orthonormal basis functions, describing the most important features of spatial and temporal variation of the state of the system. One of the most renowned techniques in SVD methods is the proper orthogonal decomposition (POD) method [37]. The POD method has been successfully applied in several fields like geophysical fluid dynamics [38][39][40], meteorology, signal analysis, and pattern recognition [41,42]. Based on this, the second objective of the work is the application of POD based MOR on the wave propagation in undamaged FML to analyze the possibility of runtime reduction.
Following these two objectives presented above, Section 2 deals with the analytical framework of GUW in thin-walled structures. In Section 3, we define the material and models that are used in Sections 4 and 5 to analyze the main wave propagation characteristics of GUW in FML. Subsequently, Section 6 covers the application of MOR on the wave propagation in undamaged FML. Finally, Section 7 gives a summary and conclusion of the presented work.

Analytical Dispersion Relation in Fiber Metal Laminates
This section covers the analytical framework of GUW in layered structures. Hence, it first gives a brief description of the dispersion relation of transversely-isotropic layers and then a description of the layered structures. This is followed by the computation of displacement fields.

Dispersion Relation of a Single Layer
This work focuses on GUW in FML consisting of thin steel foils and CFRP layers. Since the analysis is limited to the wave propagation in fiber direction and perpendicular to it, it is sufficient to assume transversely-isotropic material symmetry. For this case, a detailed implementation can be found in Refs. [17,29,43]. The analytical consideration of the dynamic behavior of a solid is provided by the balance of momentum div σ + ρb = ρü. (1) here, σ is the Cauchy stress tensor, ρ is the density, b are the volume forces, u represents the displacement field, and the dots indicate the second derivative with respect to the time. Furthermore, to derive the equation of motion, linear strains are considered and the material behavior is described by Hooke's law In Equations (2) and (3), E denotes the linear Green Lagrange strain tensor and C is the elasticity tensor. Substituting the linear strain and Hooke's law into Equation (1) and neglecting the volume forces yields div(C : grad u) = ρü.
To solve Equation (4), a plane harmonic approach is used for the displacement field given by where k j is the wavenumber vector, the angular frequency is given by ω = kc p , c p is the phase velocity, and t stands for the time. Without loss of generality the wave propagation can be described in the x 1 -x 3 plane, with x 3 being perpendicular to the plate surface, see Figure 1. In this case, the wavenumber vector is k j = k x (1, 0, α) T = kn j . Introducing the harmonic approach in Equation (5) with the angular frequency and the provided wavenumber definition into the equation of motion in Equation (4) leads to with λ il = C ijkl n k n j , which is also known as the Christoffel equation [6]. Here, p l is the polarization vector, holding the information about the particle motion direction, and δ il is the Kronecker-delta.
In the case of isotropic and transversely-isotropic material behavior, the displacement fields of Lamb and shear horizontal waves are decoupled, which is indicated by the separation of the Christoffel equation into two independent systems of equations. Here, the analysis based on the existence of Lamb waves in FML is of interest. Therefore, we focus on the description in the x 1 -x 3 plane. To obtain the non-trivial solution the determinant of the matrix Λ il must be equal to zero This condition leads to a biquadratic equation for Lamb waves providing solutions for the parameter α Each individual solution of α corresponds to a partial wave. In order to describe the propagation of Lamb waves accurately, all these partial waves must be superimposed. Thus, the displacement field of Lamb waves is described by the following equations and where S j = By combining Equations (10) and (11) with (4), we get the following expressions describing the stress components result Here, H 1j and H 3j are defined as follows Finally, the dispersion relation is determined by assuming stress-free boundary conditions at the top and bottom surfaces of the plate. Summarizing the solution in matrix form yields with σ 33 = σ 33 ik , σ 13 = σ 13 ik . It is also often given in compact form Here the abbreviation UD stands for a unidirectional layer. The detailed consideration of the solution of the dispersion relation is described and discussed in Refs. [17,29].

Wave Propagation in Layered Structures
In this work, the dispersion relation for a layered plate is established utilizing the global matrix method [17,29]. In this method, a continuity condition for the stresses and displacements at the layer boundaries is introduced to describe the wave propagation in a layered structure. This condition is expressed by where "+" indicates an association with the top edge and "−" with the bottom edge of the layer, respectively. For materials with a symmetric stacking sequence, it is possible to formulate separate dispersion relations for the symmetric and antisymmetric wave mode, respectively. Therefore, only one half of the plate is described, which enables one to introduce new boundary conditions at the symmetry plane. For the symmetric wave modes the boundary conditions at the symmetry plane are given by For the antisymmetric wave modes, the boundary conditions read With these boundary conditions and the continuity condition, the following expression for the symmetric wave modes is obtained In the formulation of the dispersion relation of the antisymmetric modes, the indices in the symmetry plane are changed from 2 and 4 to 1 and 3 in Equation (21).

Displacement Fields
For each solution of the dispersion relation characterized by a frequency-phase velocity pair, the corresponding displacement field can be calculated by inserting the solution into the provided framework to determine the still missing amplitude. For this, the vector U 0 in Equation (5) is split into the amplitude A and the polarization vector p. In order to be able to determine the displacement field, first these two still unknown quantities are calculated. This is done by introducing the phase velocities, which is a solution of Equation (21), into Equation (6). The result is an eigenvalue problem. Hence, the eigenmodes are arbitrarily scalable. Therefore, to be able to determine the polarization vector, the component p 1 of the polarization vector is set to 1.
After the polarization vector is known, the last unknown quantity, the amplitude A of the displacement field is obtained by solving GHU 0 = 0. Again, one component is arbitrary. Following this procedure, the system of equations for the determination of the amplitude is given by The displacement field is finally obtained by introducing the polarization vector and the amplitude components into Equation (10). By doing so, the displacement field, which is of complex nature, can be derived at a certain position and a defined point in time. The real-valued solution is given by

Material and Model Definition
In this section, the numerical models are defined, which are used in this work for the analysis of the main wave propagation characteristics in FML. The section covers the specimen's layout and the numerical models used to compute dispersion diagrams and displacement fields, respectively. Furthermore, a model order reduction method is applied to the numerical model.

Material Layout
The FML plate analyzed here combines CFRP layers with thin steel foils. The stacking sequence is [steel/0 • 4 /steel/0 • 2 ] s , where the steel layers have a thickness of 0.12 mm each and the CFRP layers of 0.125 mm each. Thus, the total thickness of the plate is 1.98 mm. The material layout is depicted in Figure 1, including a coordinate system definition that is valid for all simulations. The x 1 -axis coincides with the fiber orientation, the x 2 -axis is perpendicular to it, and the x 3 -axis indicates the thickness direction of the specimen.
The corresponding material properties of the two constituents that are used throughout this work are provided in Table 1.

Model Definition
The numerical simulations provided in this work cover the determination of dispersion diagrams and displacement fields. Furthermore, a model order reduction algorithm is applied. All investigations are based on a two-dimensional model in two different directions, one coincides with the fiber direction, the second is perpendicular to it. Hence, the numerical model in fiber direction represents the x 1 -x 3 -cross section defined in Figure 1. Vice versa, the x 2 -x 3 plane is used to analyze the wave propagation perpendicular to the fiber orientation.
In addition to these simulations, the displacement field is also computed based on a full three-dimensional representation of the FML laminate to investigate any influences induced by the reduction from 3D to a 2D modeling approach.

Analysis of the Dispersion Diagram
This section covers the analysis of the dispersion diagram derived from the theoretical framework provided in Section 2 and computed by numerical simulations. The results are compared and discussed.

Analytical Treatment
There is no closed solution of the analytical framework provided in Section 2. Hence, to obtain the frequency-phase velocity pairs that fulfill the dispersion relation in Equation (21), numerical algorithms are used. In this case, the dispersion relation is solved by an iterative approach incorporating the bisection procedure up to a frequency of 250 kHz. Furthermore, only the fundamental S 0 and A 0 wave mode are of interest.

Numerical Simulations
In order to derive the dispersion diagrams from the numerical data, a method of the authors is used, which has already been successfully applied in the evaluation of experimental data [30,31]. The method is based on the evaluation of the displacement fields along the propagation direction of the waves, seen in Equation (10), using a discrete 2D-DFT [32]. Thereby, the property that the propagating part of the displacement field behaves sinusoidal in both time and spatial domain [10] is exploited. Similar evaluation methods can also be found in Refs. [33][34][35]. By transferring the data into the frequencywavenumber domain, a separation of the multifrequency signal components succeeds. Frequency-wavenumber pairs can be detected for each frequency, which describe the dispersion relations of the individual modes over the course of the frequency range.
The numerical simulations are performed using the finite element method within the commercial software COMSOL Multiphysics ® . The two-dimensional numerical model and the corresponding boundary conditions used are depicted in Figure 2. The total length of the model is 1 m and the discretization is realized by second-order Lagrangian elements under plane strain assumption. Subsequently, a structured mesh of rectangular elements with quadratic shape functions is generated. For the excitation of the wave field, a multifrequency signal is used, see Ref. [31]. It is tuned to be able to generate dispersion diagrams over a predefined frequency range and hence, to avoid multiple simulations. The same signal is used for both directions, covering a frequency range from 25 kHz to 245 kHz in the steps of 10 kHz. The signal is applied by a force boundary condition at a distance of 4 mm from the left edge at the upper surface of the plate. Subsequently, both fundamental A 0 -and S 0 -modes are excited. A mode-selective excitation is not necessary, because the evaluation is done at the center line of the plate. Due to the theoretical displacement field characteristics, the S 0 -mode shows only an in-plane displacement component at the center line, whereas the A 0 -mode is limited to an out-of-plane motion. Therefore, with the help of the displacement component u 1 the S 0 -mode and with the help of the component u 3 , the A 0 -mode can be evaluated separately. Details about the model size and the discretization are provided in Table 2.  To evaluate the required frequency range with the help of the 2D-DFT, the numerical model must fulfill requirements to the spatial and temporal resolution as well as simulation time and propagation distance. The choice of the model parameters ensures that the wave bodies can be clearly detected and the individual contained frequencies can be separated from each other during the evaluation. The specific model parameters are shown in Table 2.

Results and Comparison
In this subsection, the dispersion curves of GUW in FML are compared, which are obtained firstly from analyses based on the analytical framework and secondly, from finite element simulations. The corresponding results are collected in two diagrams, one holding the information of the phase velocity-frequency pairs for the wave propagation in fiber direction, see Figure 3, left, and the other one for the wave propagation perpendicular to the fiber direction (right).
In addition, Figure 4 gives the deviations between the results of the analytical and numerical approach. The comparison shows that the maximum deviation of the phase velocity is about 0.6%. Thus, the results prove that the numerical and analytical models provide almost identical results for the dispersion behavior over a larger frequency range. This agreement is valid for the fiber direction as well as perpendicular to it.
The dispersion diagrams reveal that the propagation of GUW in FML is still of multimodal and dispersive nature. At least two modes are generated even at low frequencies.
While the fundamental symmetric mode (S 0 ) shows almost a constant phase velocity over the frequency range of up to 250 kHz, the phase velocity of the fundamental antisymmetric mode (A 0 ) increases with increasing frequency. This meets the main characteristics of the wave propagation behavior in unidirectional composite [6]. In addition to this, the phase velocity of the S 0 -mode in the fiber direction is much higher than perpendicular to it. This is also observed for the wave propagation in FRP and can be assigned to the lower Young's modulus of the CFPR layer perpendicular to the fiber orientation.

Displacement Field Analysis
The analysis of the displacement field again comprises the solution of the analytical framework and two-dimensional numerical simulations. Furthermore, full 3D simulations are conducted to analyze the effect of the domain reduction on the numerical solution. The analysis focuses on a frequency of 120 kHz. In addition, although the dispersion plots in Figure 3 do not suggest significantly different findings, computations were performed at 80 kHz. The results are shown in the Appendix A.

Analytical Treatment
Based on the solution of the dispersion relation at a frequency of 80 kHz and 120 kHz the analytical results are determined by evaluating the displacement field given in Equation (15). Since the calculation of the analytical displacement fields involves solving an eigenvalue problem, the normalized displacements are computed.

Numerical Simulations
For the numerical simulations, 2D and a 3D modeling approaches are utilized. The 2D model with the corresponding boundary conditions is shown in Figure 5. In comparison to the dispersion diagram investigation, this model has a length of 0.3 m and the displacement field over the thickness is extracted only at one position. Furthermore, at this point, a mode selective excitation is realized to be able to analyze the displacement fields of the fundamental A 0 and S 0 -mode separately. Therefore, the excitation is implemented by a five cycle Hanning windowed sine burst at the upper and lower left corner of the wave guide. For the generation of the A 0 -mode, the displacements at these two points is oriented in the same direction. In contrast to this, they point in opposite directions for the excitation of the S 0 -mode. The signal is captured at a distance of 0.1 m from the excitation point. Again the model is discretized with second-order Lagrange elements by assuming a plane strain state.
For the full 3D numerical simulations, only a quarter of the steel CFRP plate is modeled to reduce the computational cost. A schematic of the model is provided in Figure 6. The model reduction requires the introduction of the symmetry boundary conditions at the sides marked in gray. Beside the adaption to the third dimension, the modeling approach is not changed. All parameters of the models and simulations for the two cases are summarized in Table 3.

Results and Comparison
Since the calculation of the analytical displacement fields involves solving an eigenvalue problem, the normalized displacements are evaluated and compared. First, the results for the wave propagation in fiber orientation for a frequency of 120 kHz are presented in Figure 7. Here, excellent agreement can be seen between the numerically and analytically determined results. In the first row the displacement fields for the symmetric wave mode are shown, where also the typical features of the displacement fields of the axial symmetry of the in-plane component (left) and the point symmetry for out-of-plane component (right) can be observed. The x-axis of the diagrams was strongly refined for the in-plane component in order to be able to represent the deviations. The second row shows the in-plane and out-of-plane components of the displacement field for the antisymmetric wave mode. Again, the already-known behavior of the displacement field can be seen. In contrast to the symmetric wave mode, the in-plane component is point symmetric and the out-of-plane component is axisymmetric.
Subsequently, Figure 8 shows the components of the displacement fields when considering the wave propagation perpendicular to the fiber orientation at 120 kHz. Here, the typical features such as point symmetry and axis symmetry can again be observed. Among other things, the analytical and numerical results also agree excellently.
In order to quantify the small observable deviations between the solution of the analytical framework and the numerical simulations, Tables 4 and 5 hold the mean error for each component and wave mode obtained by the 2D and 3D simulations at 120 kHz. Following these results, a 2D modeling approach for the wave propagation in FML is still sufficient to capture the main wave propagation characteristics adequately. The results are also confirmed for the wave propagation at 80 kHz. Here, Figure 9 holds the results for the displacement field of the S 0 -mode for the wave propagation perpendicular to the fiber orientation. All remaining results are summarized in Appendix A. Table 4. Calculated mean error for the 2D model at 120 kHz.

Out-of-Plane in-Plane Out-of-Plane
In fiber direction 0.072% 0.019% 0.042% 0.22% Perpendicular 0.57% 0.016% 0.06% 0.28% Table 5. Calculated mean error for the 3D model at 120 kHz.  Solution of: Analytical framework 2D simulation 3D simulation Solution of: Analytical framework 2D simulation 3D simulation Figure 9. In-plane (left) and out-of-plane (right) components of normalized displacement field for S 0 wave modes perpendicular to the fiber orientation at 80 kHz from analytical treatment and numerical 2D-and 3D-simulations.

Fundamentals
We performed several numerical simulations in Sections 4.2 and 5.2, by utilizing a finite element discretized undamped linear dynamic structural model as expressed in the following equation Equation (24) follows from the balance of momentum in Equation (1) by applying variational calculus and subsequent discretization of the equilibrium in its weak form, where M(θ) ∈ R N×N denotes the global mass matrix, K(θ) ∈ R N×N represents the global stiffness matrix, f ∈ R N×m the load matrix, and u(θ) ∈ R N×m the displacement matrix. The vector θ ∈ D ⊂ R N p contains N p parameters that define the system from parametric domain D and t ∈ [0, t max ] is the time variable. The total number of degrees of freedom in the system is given by N and the number of discretized time steps is indicated by m. The adapted model order reduction method is as described in Ref. [44] and summarized below.
The projection of basis functions, using a projection matrix Φ(θ) ∈ R N×n with n << N, from high-dimensional space R N allows the order reduction of the high-fidelity model to its lower dimensional space R n u ≈ Φu r ,ü ≈ Φü r .
Substituting Equation (25) in Equation (24) and projecting them on the reduced order space results in the following equations In Equation (26), M r , K r ∈ R n×n and f r ∈ R n×m denote the reduced system matrices. The projection matrix Φ can be obtained by the POD method. This involves the collection of observations from the HiFi model into a snapshot matrix Ψ(θ) Ψ = [u(θ, t 1 ), u(θ, t 2 ), . . . , u(θ, t m )], = [u 1 , u 2 , . . . , u m ], (27) and decomposing it by thin singular value decomposition method as follows: in Equation (28), the left and right-singular matrices are given by W(θ) = [w 1 , w 2 , . . . , w m ] ∈ R N×m and V, respectively. The right-singular matrix V serves no interest in this context, while the left-singular matrix W contains orthogonal basis functions also known as proper orthogonal models (POMs) of the system. The diagonal matrix Σ = diag(σ 1 , σ 2 , . . . , σ m ) ∈ R m×m contains singular values σ k m k=1 with σ 1 ≥ σ 1 ≥ . . . σ m > 0. The relative energy measure Ξ, which is used to determine the number of most influential modes and its approximation error ε are obtained as follows As the projection matrix Φ = [w 1 , w 2 , . . . , w n ] ∈ R N×n is computed, the reduced order Equation (26) can be deduced and solved for u r andü r . Eventually, the solution for the HiFi system can be evaluated using Equation (25). The choice of the observations extracted in the snapshot matrix largely influence the accuracy of the reduced-order model. As the reduced-order models produced through this technique are based on a certain parameter, their application at off-design parameters often fails to accurately predict the solution. In other words, whenever the system parameters changes, the reduced order basis should be rebuilt again. This is practically intractable in inverse problem analysis wherein the system should be solved for several thousand times. Therefore, a more sophisticated robust reduced-order modeling technique that preserves the parametric dependence of the system has to be developed. PMOR methods precisely satisfy this requirement. PMOR involves the extraction of global proper orthogonal modes (POMs) that accurately define the behavior of the underlying system for any given parameter set within the parametric domain.
The parameterized elliptic and parabolic PDE models are well approximated by PODbased global reduced-order models. However, the PMOR of a wide range of hyperbolic problems with non-linearity and discontinuity still remain a challenge. Therefore, the researchers within the MOR community are intensively seeking methodologies to reduce the Kolmogorov n-width of the solution manifold [45,46]. Boncoraglio et al. [47] efficiently solved multidisciplinary design optimization problem, which blends both linear and nonlinear constraints in aerodynamics using projection-based MOR along with an active manifold. The authors utilized a deep convolutional autoencoder to discover the relevant active manifold for dimensionality reduction of a high-dimensional design parameter space. An efficient adaptive algorithm was implemented by Bui-Thanh et al. [48] to achieve MOR for design optimization of a heat conduction fin to determine appropriate sample points over a large input parametric space. Based on Refs. [49,50] McBane et al. proposed a component-wise reduced-order model to optimize the topology of a lattice-type structure [51]. Moreover, they simplified the model to increase the speedup of the optimization process. A space-time MOR method built on least-squares Petrov-Galerkin projection was presented by Kim et al. to solve linear dynamical systems [52]. The approach was well demonstrated on 2D diffusion and 2D convection diffusion problems. Further contributions on PMOR span across the domains of contact in multibody nonlinear dynamics [53], nonlinear fluid-structure interaction problems [54], uncertainty quantification [55,56], and contact mechanics [57].

Application
In this research work, a POD-based adaptive parametric model order reduction technique along with a surrogate model for estimating the error indicator, based on the work of Paul-Dubois-Taine et al. [58], was employed to produce the global POMs. The adaptive sampling of the parameters was achieved in a greedy sense, wherein the locally best solution in a particular iteration was sampled, assuming it to be the global optimal solution. The training phase resulted in global reduced-order basis, which could approximate the HiFi solution with a very high accuracy. A detailed description of the procedure along with its algorithms can be found in Ref. [44], wherein the authors have applied the proposed PMOR method to approximate the GUW propagation in FML containing a damage. In this paper, the analysis is restricted to the GUW propagation in an intact plate. Since the PMOR method employed in this research is a model-intrusive approach, the Newmark's time integration technique is utilized to solve Equation (24). This required the excitation of the specimen to be applied in a force-controlled manner and the other boundary conditions are consistent with those from the study of the displacement fields on the 2D model. The discretization is realized by quadratic Lagrangian elements with a size of 0.5 mm in the length direction and one element over the thickness for each layer. The evaluation was performed at a distance of 57.5 mm from the excitation points in the center of the upper CFRP layer package.
The comparison plot of the out-of-plane displacement of FML over time obtained using the HiFi model and the reduced-order model, produced using the proposed method of adaptive PMOR approach with Galerkin projection, is shown in Figure 10. The displacement field was measured at 57.5 mm away from the left end of the laminate and 1.49 mm through the thickness from the bottom of the laminate. It can be seen that a precise approximation of the solution was achieved by the reduced-order model. The global reduced-order basis has gathered 250 modes to attain this accuracy of approximation of GUW propagation in an intact FML. This numerical experiment was conducted on a 4-core Intel(R) Core(TM) i7-10510U CPU @ 1.80 GHz processor with 16 GB RAM. Table 6 summarizes the computational costs of HiFi solution and its reduced-order solution. The application of a POD-based adaptive PMOR approach, which employed the surrogate model for estimating error indicators, resulted in a speedup factor of 45.72, corresponding to the setting mentioned in this research paper. This essential reduction of computation time affords an opportunity to conduct several parametric studies and analyze the GUW propagation in FML.

Conclusions
The investigations presented here are based on an FML plate with the stacking sequence [steel/0 4 /steel/0 2 ] s . As the first step, the dispersive nature of GUW in a steel-CFRP laminate was analyzed by computing the dispersion diagrams for the wave propagation in the fiber direction and perpendicular to it. This was done by utilizing the analytical framework and numerical 2D simulations. The analytically and numerically obtained dispersion diagrams for a frequency range up to 250 kHz coincide very well. Furthermore, it was found that the dispersive nature is similar to the propagation characteristics in FRP. In the following step, the phase velocity derived from the dispersion diagram for two excitation frequencies (80 kHz and 120 kHz) was used to compute the displacement field following the analytical framework. Furthermore, numerical 2D and 3D simulations were conducted to obtain the displacement field components for the wave propagation in the fiber direction and perpendicular to it. The resulting displacement fields show identical properties, as described in the literature for the GUW propagation in isotropic, transversely-isotropic, and orthotropic materials [16,17,29]. For the wave propagation in the fiber orientation as well as perpendicular to it, the in-plane component of the S 0 -wave mode showed an axisymmetric shape and the out-of-plane component showed a point-symmetric shape. For the A 0 -wave mode, the opposite features were observed-point symmetry for the in-plane component and axisymmetry for the out-of-plane component. The agreement between the qualitative analytical and numerical results both in the direction and perpendicular to the fiber orientation is excellent, and the deviations are difficult to detect, which was also shown by the overview of the percentage deviations of the displacement fields in Tables 4 and 5. The complex shape of the displacement fields indicate that finite shell elements are not suitable for modeling and 3D continuum elements should be used to accurately represent the displacement fields when a full 3D simulation is required.
Since the investigations presented here involve high computational effort, the possibility of model reduction is also investigated using a 2D model as an example. For this purpose, a PMOR method is applied to approximate the GUW propagation in an FML plate. In Figure 10 it can be seen that in the context of this model reduction, no disadvantages in the accuracy of the model occurred, and also a speed-up factor of the model of 45.72 was achieved.
In conclusion, the numerical investigation of the main wave propagation characteristics in a steel-CFRP laminate reveal the well-known wave propagation properties of GUW in thin-walled structures. As already shown for isotropic, transversely-isotropic, and orthotropic materials, the modeling approach can successfully be reduced to 2D for the wave propagation along a symmetry axis. Furthermore, the computation time can be significantly reduced by a PMOR method with a minimal loss of information regarding the wave propagation. In further investigations, the displacement fields on other FML plate configurations with different metal-UD layer ratios will be investigated.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Abbreviations
The following abbreviations are used in this manuscript:

Appendix A
Below, we summarize the results of the displacement field analysis at 80 kHz. Therefore, Figures A1 and A2 hold the results of the displacement field components of the A 0 -mode perpendicular to the fiber direction and the components of the displacement fields in fiber direction for both the A 0 -and S 0 -modes. The corresponding results of the mean error calculation are collected in Table A1 for the 2D simulation results and in Table A2 for the 3D simulation results. Solution of: Analytical framework 2D simulation 3D simulation Figure A1. In-plane (left) and out-of-plane (right) components of a normalized displacement field for A 0 wave mode perpendicular to the fiber orientation at 80 kHz derived from the analytical framework and numerical 2D-and 3D-simulations. Solution of: Analytical framework 2D simulation 3D simulation Figure A2. In-plane (left) and out-of-plane (right) components of a normalized displacement field for S 0 (top) and A 0 (bottom) wave modes in the fiber orientation at 80 kHz derived from the analytical framework and numerical 2D-and 3D-simulations.