Comparison of the CA/MCA and PROM methods.
Finite element analysis (FEA) is a well-established numerical simulation method for structural dynamics. It serves as the main computational tool for Noise, Vibration and Harshness (NVH) analysis in the low-frequency range. Because of developments in numerical methods and advances in computer software and hardware, FEA can now handle much more complex models far more efficiently than even a few years ago. However, the demand for computational capabilities increases in step with or even beyond the pace of these improvements. For example, automotive companies are constructing more detailed models with millions of degrees of freedom (DOFs) to study vibro-acoustic problems in higher frequency ranges. Although these tasks can be performed with FEA, the computational cost can be prohibitive even for high-end workstations with the most advanced software.
For large finite element (FE) models, a modal reduction is commonly used to obtain the system response. An eigenanalysis is performed using the system stiffness and mass matrices and a smaller in size modal model is formed which is solved more efficiently for the response. The computational cost is also reduced using substructuring (superelement analysis). Modal reduction is applied to each substructure to obtain the component modes and the system level response is obtained using Component Mode Synthesis (CMS).
When design changes are involved, the FEA analysis must be repeated many times in order to obtain the optimum design. Furthermore in probabilistic analysis where parameter uncertainties are present, the FEA analysis must be repeated for a large number of sample points. In such cases, the computational cost is even higher, if not prohibitive. Reanalysis methods are intended to analyze efficiently structures that are modified due to various changes. They estimate the structural response after such changes without solving the complete set of modified analysis equations. Several reviews have been published on reanalysis methods [1-3] which are usually based on local and global approximations. Local approximations are very efficient but they are effective only for small structural changes. Global approximations are preferable for large changes, but they are usually computationally expensive especially for cases with many design parameters. The well-known Rayleigh-Ritz reanalysis procedure [4, 5] belongs to the category of local approximation methods. The mode shapes of a nominal design are used to form a Ritz basis for predicting the response in a small parametric zone around the nominal design point. However, it is incapable of capturing relatively large design changes.
A parametric reduced-order modeling (PROM) method, developed by Balmes [6, 7], expands on the Rayleigh-Ritz method by using the mode shapes from a few selected design points to predict the response throughout the design space. The PROM method belongs to the category of local approximation methods and can handle relatively larger parameter changes because it uses multiple design points. An improved component-based PROM method has been introduced by Zhang et al. [8, 9] for design changes at the component level. The PROM method using a ‘parametric’ approach has been successfully applied to design optimization and probabilistic analysis of vehicle structures. However, the ‘parametric’ approach is only applicable to problems where the mass and stiffness matrices can be approximated by a polynomial function of the design parameters and their powers. A new ‘parametric’ approach using Kriging interpolation  has been recently proposed . It improves efficiency by interpolating the reduced system matrices without needing to assume a polynomial relationship of the system matrices with respect to the design parameters as in [6, 7].
The Combined Approximations (CA) method [12-14] combines the strengths of both local and global approximations and can be accurate even for large design changes. It uses a combination of binomial series (local) approximations, called Neumann expansion approximations, and reduced basis (global) approximations. The CA method is developed for linear static reanalysis and eigen-problem reanalysis [15-19]. Accurate results and significant computational savings have been reported. All reported studies on the CA method [12-19] use relatively simple frame or truss systems for static or dynamics analysis with a relatively small number of DOF and/or modes. For these problems, the computational efficiency was improved by a factor of 5 to 10. Such an improvement is beneficial for a single design change evaluation or even for gradient-based design optimization where only a limited number of reanalyses (e.g. less than 50) is performed. However, the computational efficiency of the CA method may not be adequate in simulation-based (e.g. Monte-Carlo) probabilistic dynamic analysis of large finite-element models where reanalysis must be performed hundreds or thousands of times in order to estimate accurately the reliability of a design.
A large number of modes must be calculated and used in a dynamic analysis of a large finite-element model with a high modal density, even if a reduced-order modeling approach (Section 2) is used. In such a case, the implicit assumption of the CA method that the cost of solving a linear system is dominated by the cost of matrix decomposition way no be longer valid (see Section 3.4) and the computational savings from using the CA method may not be substantial. For this reason, we developed a modified combined approximation (MCA) and integrated it with the PROM method to improve accuracy and computational efficiency. The computational savings can be substantial for problems with a large number of design parameters. Examples in this Chapter demonstrate the benefits of this reanalysis methodology.
The Chapter presents methodologies
for accurate and efficient vibration analysis methods of large-scale, finite-element models,
for efficient and yet accurate reanalysis methods for dynamic response and optimization, and
for efficient design optimization methods to optimize structures for best vibratory response.
The optimization is able to handle a large number of design variables and identify local and global optima. It is organized as follows. Section 2 presents an overview of reduced-order modeling and substructuring methods including modal reduction and component mode synthesis (CMS). Improvements to the CMS method are presented using interface modes and filtration of constraint modes. The section also overviews two Frequency Response Function (FRF) substructuring methods where two substructures, represented by FRFs or FE models, are assembled to form an efficient reduced-order model to calculate the dynamic response. Section 3 presents four reanalysis methods: the CDH/VAO method, the Parametric Reduced Order Modeling (PROM) method, the Combined Approximation (CA) method, and the Modified Combined Approximation (MCA) method. It also points out their strong and weak points in terms of efficiency and accuracy. Section 4 demonstrates how the reanalysis methods are used in vibration and optimization of large-scale structures. It also presents a new reanalysis method in Craig-Bampton substructuring with interface modes which is very useful for problems with many interface DOFs where the FRF substructuring methods cannot be used. Section 5 presents a vibration and optimization case study of a large-scale vehicle model demonstrating the value of reduced-order modeling and reanalysis methods in practice. Finally, Section 6 summarizes and concludes.
2. Reduced-Order Modeling for Dynamic Analysis
Computational efficiency is of paramount importance in vibration analysis of large-scale, finite-element models. Reduced-order modeling (or substructuring) is a common approach to reduce the computational effort. Substructuring methods can be classified in “mathematical” and “physical” methods. The “mathematical” substructuring methods include the Automatic Multi-level Substructuring (AMLS) and the Automatic Component Mode Synthesis (ACMS) in NASTRAN. The “physical” substructuring methods include the well known fixed-interface Craig-Bampton method. Both the AMLS and ACMS methods use graph theory to obtain an abstract (mathematical) substructuring using matrix partitioning of the entire finite-element model. The computational savings from the “mathematical” and “physical” methods can be comparable depending on the problem at hand.
2.1. Modal Reduction
For an undamped structure with stiffness and mass matrices and respectively, under the excitation force vector, the equations of motion (EOM) for frequency response are
where the displacement d is calculated at the forcing frequency. If the response is required at multiple frequencies, the repeated direct solution of Equation (1) is computationally very expensive and therefore, impractical for large scale finite-element models.
A reduced order model (ROM) is a subspace projection method. Instead of solving the original response equations, it is assumed that the solution can be approximated in a subspace spanned by the dominant mode shapes. A modal response approach can be used to calculate the response more efficiently. A set of eigen-frequencies and corresponding eigenvectors (mode shapes) are first obtained. Then, the displacement d is approximated in the reduced space formed by the first dominant modes as
The response d can be recovered by solving Equation (3) for the modal response U and projecting it back to the physical coordinates using Equation (2). If is the maximum excitation frequency, it is common practice to retain the mode shapes in the frequency range of. The system modes in the high frequency range can be safely truncated with minimal loss of accuracy.
Due to the modal truncation, the size of the ROM is reduced considerably, compared to the original model. However, the size increases with the maximum excitation frequency. An added benefit of the ROM is that Equations (3) are decoupled because of the orthogonality of the mode shapes and can be therefore, solved separately reducing further the overall computational effort.
Note that for a damped structure with a damping matrix, Equation (1) becomes and Equation (3) is modified as
2.2. Substructuring with Component Mode Synthesis (CMS)
To model the dynamics of a complex structure, a finite-element analysis of the entire structure can be very expensive, and sometimes infeasible, due to computer hardware and/or software constraints. This is especially true in the mid-frequency range, where a fine finite element mesh must be used in order to capture the shorter wavelengths of vibration. Component mode synthesis (CMS) was developed as a practical and efficient approach to modeling and analyzing the dynamics of a structure in such circumstances [20–23]. The structure is partitioned in component structures and the dynamics are described by the normal modes of the individual components and a set of modes that couple all component. Besides the significant computational savings, this component-based approach also facilitates distributed design. Components may be modified or redesigned individually without re-doing the entire analysis.
One of the most accurate and widely-used CMS methods is the Craig-Bampton method  where the component normal modes are calculated with the interface between connected component structures held fixed. Attachment at the interface is achieved by a set of “constraint modes.” A constraint mode shape is the static deflection induced in the structure by applying a unit displacement to one interface DOF while all other interface DOF are held fixed. The motion at the interface is thus completely described by the constraint modes. The Craig-Bampton reduced-order model (CBROM) results in great model size reduction by including only component normal modes within a certain frequency range. However, there is no size reduction for constraint modes because CBROM must have one DOF for each interface DOF. If the finite element mesh is sufficiently fine, the constraint-mode DOFs will dominate the CBROM mass and stiffness matrices, and increase the computational cost.
We address this problem by using interface modes (also called characteristic constraint –CC- modes) in order to reduce the number of retained interface DOFs of the Craig-Bampton approach. For that, a secondary eigenvalue analysis is performed using the constraint-mode partitions of the CMS mass and stiffness matrices. The CC modes are the resultant eigenvectors. The basic formulation is described in Sections 2.2.1 and 2.2.2. The interface modes represent more “natural” physical motion at the interface. Because they capture the characteristic motion of the interface, they may be truncated as if they were traditional modes of vibration, leading to a highly compact CC-mode-based reduced order model (CCROM). In addition, the CC modes provide a significant insight into the physical mechanisms of vibration transmission between the component structures. This information could be used, for example, to determine the design parameters that have a critical impact on power flow. Figure 1 compares a conventional constraint mode used in Craig-Bampton analysis with an interface mode, for a simple cantilever plate which is subdivided in two substructures.
It should be noted that the calculation of the CC modes is essentially a secondary modal analysis. Therefore, the benefits are the same as those of a traditional modal analysis. For instance, refining the finite-element mesh increases the accuracy of a CCROM without introducing any additional degrees of freedom. The ability of the CC mode approach to produce CCROM whose size does not depend on the original level of discretization makes it especially well suited for finite-element based analysis of mid-frequency vibration.
2.2.1. Craig-Bampton Fixed Interface CMS
This section provides the basics of Craig-Bampton method using the fixed-interface assumption. The method is commonly used in CMS algorithms. The finite-element model of the entire structure is partitioned into a group of substructures. The DOFs in each substructure are divided into interface () DOF and interior () DOF. The equations of motion for the i th substructure are then expressed as
The fixed-interface Craig-Bampton CMS method utilizes two sets of modes to represent the substructure motion; substructure normal (N) modes, and constraint (C) modes, where i denotes the i th substructure. The size reduction of the Craig-Bampton method comes from the truncation of the normal modes which are calculated by fixing all interface DOFs and solving the following eigenvalue problem
The static constraint modesare calculated by enforcing a set of static unit constraints to the interface DOFs as
The original physical DOFs
Using the above Craig-Bampton transformation, the original EOM of Equation (7), can be expressed as
where the superscripts C and N are used to indicate partition related to static constraint mode DOFs and fixed-interface normal mode DOFs, respectively. The matrix partitions of Equation (8) are
The matrices of each substructure are then assembled by applying displacement continuity and force balance along the interface to obtain the EOM of the reduced system. A secondary modal analysis is finally carried out using the mass and stiffness matrices of the reduced system to obtain the eigenvalues and eigenvectors.
Note that constraint mode matrix
solving for the normal modes,
solving for the constraint modes, and
the transformation calculation in Equation (9).
2.2.2. Craig-Bampton CMS with Interface Modes
In Craig-Bampton CMS, the matrices from all substructures are assembled into a global CBROM with substructures coupled at interfaces by enforcing displacement compatibility. This synthesis yields the modal displacement vector
where is the number of substructures in the global structure. The corresponding synthesized CMS mass and stiffness matrices are as follows
where the component modal matrices
The eigenvectors are transformed into the finite-element DOFs for the i th component structure using the following transformation
is a selected set of interface eigenvectors which are few compared to the number of the constraint-mode DOFs. The matrix
Finally, the CMS matrices are transformed using the CC modes and the reduced-order CMS matrices are obtained similarly to Equations (18). Now, the unknown displacement vector dROM is partitioned as
where superscript CC indicates the partition associated with the CC modes. The equations of motion of the reduced order CMS model (ROM) are expressed by
The mass matrix MROM, the stiffness matrix KROM, and the applied force vector fROM, are explicitly written as
2.2.3. Filtration of Constraint Modes
Figure 2 shows a typical constraint mode for a plate substructure. The non-zero displacement field (indicated by red color) is usually limited to a very small region close to the perturbed interface DOF.
If the small-displacement part of the constraint mode shape is explicitly replaced by zero, the density of the resulting “filtered” constraint mode will be significantly reduced. Consequently, the computational cost in Equation (9) will be considerably reduced. To filter the constraint modes, the following criterion is used
where is the p th element of the q th constraint mode
2.3. Frequency Response Function (FRF) Substructuring and Assembling
If the number of interface nodes (or DOFs) between connected substructures is small, a reduced-order model of small size can be developed using an FRF representation of each substructure. This is known as FRF substructuring. The FRF representation can be easily obtained from a finite element (FE) model or even experimentally. If the FE model of one substructure is very small (e.g. a vehicle suspension model), it can be easily coupled directly to another substructure which is represented using FRFs. This section provides the fundamentals of FRF substructuring for both FRF-FE and FRF-FRF coupling.
2.3.1. Algorithm for FRF/FE Coupling
The numerical algorithm is explained using the two-substructure example of Figure 3. Sub1 is an FRF type substructure, meaning that its dynamic behavior is described using FRF matrices which are denoted by H (see Equation 30 for notation). The elements of H are frequency dependent and complex if damping is present. A bold letter indicates a matrix or vector. According to Equation (30), HAC for example, represents the displacement XA of DOF A due to a unit force FC on DOF C. Sub2 is a finite element (FE) type substructure. Its dynamic behavior is described using the stiffness K, mass M and damping B matrices.
The FRF matrix of Sub1 can be calculated by either a direct frequency response method or a modal response method. In the former case, the original FE equations are used in the frequency domain while in the latter a modal model is first developed and then used to calculate the FRF matrix. The size of the FRF matrix is small and depends on the number of DOFs of the excitation, response and interface DOFs. Usually FRFs are calculated between excitation and response DOFs. However in order to assemble two substructures, FRFs are also calculated between interface DOFs and excitation/response DOFs. The Sub1 FRF matrix H in Equation (30) is thus partitioned into interior (A) DOFs and interface/coupling (C) DOFs. The interior DOFs include all excitation and all response DOFs (Figure 3).
The second substructure Sub2 is expressed in FE format. The system FE matrices K, M, and B form the frequency dependent dynamic matrix which is then partitioned according to the interior (B) and interface (C) DOFs. The interface DOFs for Sub1 and Sub2 have the same node IDs so that they can be assembled to obtain the system FRF matrix. The procedure to assemble the H matrix of Sub1 with the Z matrix of Sub2 and calculate (solve for) the system matrix H is described below.
The equations of motion for Sub1 are expressed as
where subscript A indicates the interior (excitation plus response) DOFs of Sub1, and subscript C indicates the connection/common/coupling DOFs between Sub1 and Sub2. The equations of motion for Sub2 are expressed as
where subscript B indicates the interior DOFs of Sub2, and subscript C indicates the connection/common/coupling DOFs between Sub2 and Sub1. Because of displacement compatibility at the interface, appears on the left-hand side of Equation (30) for Sub1 and on the right-hand side of Equation (31) for Sub2. Superscripts 1 and 2 are used to differentiate the interface forces FC at Sub1 and Sub2.
To couple Sub1 and Sub2, compatibility of forces at the interface is applied as
From Equation (31),
Solving Equation (30) for
resulting in the following FRF of the assembled system
where S = [Φ Θ I] and
2.4.2. Algorithm for FRF/FRF Coupling
Figure 4 shows the coupling of two FRF type substructures. The equations of motion for Sub1 and Sub2 and are expressed as
To couple Sub1 and Sub2, we enforce displacement compatibility at the interface and apply the interface force relationship
The FRF/FRF coupling is a special case of the FRF/FE coupling.
3. Reanalysis Methods for Dynamic Analysis
3.1. CDH/VAO Method
The CDH/VAO method, developed by CDH AG for vibro-acoustic analysis, is a Rayleigh-Ritz type of approximation. If the stiffness and mass matrices of the baseline design structure are K0 and M0, the exact mode shapes in Φ0 are obtained by solving the eigen-problem
where Λ0 is the diagonal matrix of the baseline eigenvalues. A new design (subscript p) has the following stiffness and mass matrices
For a modest design change where ΔK and ΔM are small, it is assumed that the change in mode shapes is small and the new response can be therefore, captured in the sub-space spanned by the mode shapes Φ0 of the initial design. The new stiffness and mass matrices are then condensed asand and the following reduced eigen-value problem is solved to calculate the eigen-vector
The approximate eigenvalues of the new design are given by Λp and the approximate eigenvectors Φp are
where. Thus, the modal response of the modified structure can be easily obtained and the actual response can be recovered using the eigenvectors.
3.2. Parametric Reduced-Order Modeling (PROM) Method
The PROM method approximates the mode shapes of a new design in the subspace spanned by the dominant mode shapes of some representative designs, which are selected such that the formed basis captures the dynamic characteristics in each dimension of the parameter space. Balmes et al. [6, 7] suggested that these representative designs should correspond to the middle points on the faces of a box in the parameter space representing the range of design parameters. For a structure with design variables, Zhang  suggested that the representative designs include a baseline design for which all parameters are at their lower limits plus m designs obtained by perturbing the design variables from their lower limits to their upper limits, one at a time. The points representing these designs in the space of the design variables are called corner points (see Figure 5). This selection of representative designs results in a more accurate PROM algorithm.
The mode shapes of a new design are approximated in the space of the mode shapes of the corner points as
where the modal matrix P includes the basis vectors as in Equation (49) and represents the participation factors of these vectors. The columns of P are the dominant mode shapes of the above designs,
whereis the modal matrix composed of the dominant mode shapes of the baseline design, and is the modal matrix of the corner point. The mode shapes of the new design satisfy the following eigenvalue problem,
where Λ is a diagonal matrix of the first s eigenvalues.
A reduced eigenvalue problem is obtained by pre-multiplying both sides of Equation (50) by as
where the reduced stiffness and mass matrices are
Thus, the matrix in Equation (48) consists of the eigenvectors of the reduced stiffness and mass matrices and.
For design variables, eigenvalue problems must be solved in order to form the basis of Equation (49). Therefore, both the cost of obtaining the modal matrices and the size of matrix increase linearly with m. The PROM approach uses the following algorithm to compute the mode shapes of a new design:
Calculate the mode shapes of the baseline design and the designs corresponding to the corner points in the design space, and form subspace basis.
Calculate the reduced stiffness and mass matrices and from Equation (52).
Solve eigenproblem (51) for matrix.
Reconstruct the approximated eigenvectors in using Equation (48).
Step 1 is performed only once. A reanalysis requires only steps 2 to 4. For a small number of mode shapes and a small number of design variables, the cost of steps 2 to 4 is much smaller than the cost of a full analysis. The computational cost of PROM consists of
the cost of performing full eigen-analyses to form subspace basis in Equation (49), and
the cost of reanalysis of each new design in steps 2 to 4.
The former is the fixed cost of PROM because it does not depend on the numbers of reanalyses and the latter is the variable cost of PROM because it is proportional to the number of reanalyses. The fixed cost is not attributed to the calculation of the response for a particular design. It is simply required to obtain the information needed to apply PROM. The variable cost (cost of reanalysis of a new design in part b) is small compared to the fixed cost. The fixed cost of PROM is proportional to the number of design variables m because it consists of the dominant eigenvectors of the baseline design and the dominant eigenvectors of the m corner design points (see Equation 49). When the size of basis increases so does the fixed cost because more eigenvalue problems and mode shapes must be calculated to form basis. The PROM method results in significant cost savings when applied to problems that involve few design variables (less than 10) and require many analyses (e.g. Monte Carlo simulation or gradient-free optimization using genetic algorithms).
3.3. Combined Approximations (CA) Method
The PROM method requires an eigenvalue analysis for multiple designs (corner points) to form a basis for approximating the eigenvectors at other designs. It is efficient only when the number of design parameters is relatively small. On the contrary, the CA method of this section does not have such a restriction because the reanalysis cost is not proportional to the number of design parameters. The CA method is thus more suitable than the PROM method, when the number of reanalyses is less than or comparable to the number of design parameters, such as in gradient-based design optimization.
The fundamentals of the combined approximations (CA) method [15-19] are given below. A subspace basis is formed through a recursive process for calculating the natural frequencies and mode shapes of a system. If and are the stiffness and mass matrices of the original (baseline) design, the exact mode shapes are obtained by solving the eigen-problem. We want to approximate the mode shapes of a modified design (subscript p) with stiffness and mass matrices and where and represent large perturbations. The CA method estimates the new eigenvalues and eigenvectors without performing an exact eigenvalue analysis.
The eigen-problem for the modified design can be expressed as
leading to the following recursive equation
which produces a sequence of approximations of the mode shapes. The CA method uses the changes to form a subspace basis to approximate the modes of the new design. In order to simplify the calculations, in Equation (53) is replaced with and Equation (54) becomes showing that the basis vectors satisfy the following recursive equation
where the first basis vector is assumed to be.
The CA method forms a subspace basis
Condense the stiffness and mass matrices as
Solve the reduced eigen-problem (using matrices KR and MR) to calculate the eigenvector matrix.
Reconstruct the approximate eigenvectors of the new designas
The eigenvalues of the new design are approximated by the eigenvalues of the reduced eigen-problem.
The CA method has three main advantages:
it only requires a single matrix decomposition of stiffness matrix in Equation (55) to calculate the subspace basis,
it is accurate because the basis is updated for every new design, and
the eigenvectors of a new design are efficiently approximated by Equation (58) where the eigenvectors correspond to a much smaller reduced eigen-problem.
However for a large number of reanalyses, the computational cost can increase substantially because a new basis and the condensed mass and stiffness matrices in Equation (57) must be calculated for every reanalysis. Examples where many analyses are needed are optimization problems in which a Genetic Algorithm is employed to search for the optimum, and probabilistic analysis problems using Monte-Carlo simulation. The PROM method can be more suitable for these problems because the subspace basis does not change for every new design point. Note that steps 1 and 3 (Equations 57and 58) are similar to steps 2 and 4 of PROM. CA uses basis and PROM uses basis.
The CA method is more efficient than PROM for design problems where few reanalyses are required for two reasons. First, it does not require calculation of the eigenvectors, of the corner design points, and second the cost of matrix condensation of Equation (57) is much lower than that of Equation (52), because the size (number of columns) of basis is not proportional to the number of parameters as in basis P . For problems with a large number of design parameters, the PROM approach is efficient only when a ‘parametric’ relationship is established  because a large overhead cost, proportional to the number of design parameters, is required. In contrast, the CA method does not require such an overhead cost because the reanalysis cost is not proportional to the number of design parameters. The CA method is therefore, more suitable than PROM, if the number of reanalyses is less than or comparable to the number of design parameters.
3.4. Modified Combined Approximations (MCA) Method
In the literature, the accuracy and efficiency of the CA method has been mostly tested on problems involving structures with up to few thousands of DOFs, such as frames or trusses [12-19]. We have tested the CA method using, among others, the structural dynamics response of a medium size (65,000 DOFs) finite-element model (Figure 7 of Section 4.3). Due to its high modal density, there were more than two hundred dominant modes in the frequency range of zero to 50 Hz. It was observed that the computational savings of the CA method, using the recursive process of Equation (55), were not substantial. For this reason, we developed a modified combined approximations method (MCA) by modifying the recursive process of Equation (55) which is much more efficient than the original CA method for large size models.
The cost of calculating the subspace basis in Equation (55) consists of one matrix decomposition (DCMP) and one forward-backward substitution (FBS). The DCMP cost is only related to the size and density of the symmetric stiffness matrix, while the FBS cost depends on both the size and density of the stiffness matrix and the number of columns of. As the frequency range of interest increases, more modes are needed to predict the structural response accurately. In such a case, although a single DCMP is needed in Equation (55), the number of columns in may increase considerably, thereby increasing the cost of the repeated FBS. When the number of dominant modes becomes very large, the cost of performing the calculations in Equation (55) can be dominated by the FBS cost. For example, the vehicle model of Section 4.3 (Figure 7) has 65,000 degrees of freedom and 1050 modes in the frequency range of 0-300 Hz. The cost of one DCMP is 1.1 seconds (using a SUN ULTRA workstation and NASTRAN v2001) and the cost of one FBS is less than 0.1 seconds if has only one mode. In this case, the total cost is dominated by the DCMP, and the CA method reduces the cost of one reanalysis considerably. However, if has 1050 modes, the cost of FBS increases to 29 seconds dominating the cost of the DCMP. The CA method can therefore, improve the efficiency only when the number of retained modes is small. Otherwise, the computational savings do not compensate for the loss of accuracy from using (stiffness matrix of baseline design) instead of (stiffness matrix of new design). The modified combined approximations (MCA) method of this section addresses this issue.
The MCA method uses a subspace basis T whose columns are constructed using the recursive process
instead of that of Equation (55). The selection of the appropriate number of basis vectors is discussed later in this section. The only difference between Equations (55) and (59) is that matrix is inverted in the former while matrix is inverted in the latter. The DCMP of must be repeated for every new design. However, the cost of the repeated DCMP does not significantly increase the overall cost in Equation (59) because the latter is dominated by the FBS cost. The iterative process of Equation (59) provides a continuous mode shape updating of the new design. If the process converges in s iterations, the mode shapes will be the exact mode shapes. Equation (55) does not have the same property. The vectors provide therefore, a more accurate approximation of the exact mode shapes than the vectors of the original CA method. This is an important advantage of MCA.
Because the mode shapes in Equation (59) can quickly converge to the exact mode shapes, for many practical problems only one iteration (i.e. s = 1) may be needed, resulting in
If the convergence is slow, multiple sets of updated mode shapes can be used so that
For better accuracy, the above basis also includes the mode shapes of the baseline design. Because the approximate modes are more accurate than the CA vectors in approximating the exact mode shapes, MCA can achieve similar accuracy to the CA method using fewer modes. The example of Section 4.3 demonstrates that the MCA method achieves good accuracy with only 1 basis vector whereas the CA method requires 3 to 6 basis vectors [13-17].
In summary, the proposed MCA method involves four steps in calculating the approximate eigenvectors as follows
Calculate the condensed stiffness and mass matrices KR and MR as
Solve the following reduced eigen-problem to calculate the eigenvalues and the projections of the modes in the reduced space spanned by T
Reconstruct the approximate eigenvectors as
The slightly increased cost of using Equation (61) instead of Equation (60) is usually smaller than the realized savings in steps 2 through 4 of Equations (62) through (64) due to the smaller size of the reduced basis. The bases of Equations (60) and (61) are smaller in size than the CA basis of Equation (56) for comparable accuracy. The MCA method requires therefore, less computational effort for steps 2 through 4. The computational savings compensate for the increased cost of DCMP for dynamic reanalysis of large finite-element models with a large number of dominant modes.
All mode shapes in Equation (63) must be calculated simultaneously in order to ensure that the approximate mode shapes are orthogonal with respect to the mass and stiffness matrices. However, the cost of estimating the mode shapes using Equations (62) to (64) may increase quickly (quadratically) with the number of modes, and as a result, the MCA method may become more expensive than a direct eigen-solution when the number of dominant modes exceeds a certain limit. One way to circumvent this problem is to divide the frequency response into smaller frequency bands and calculate the frequency response in each band separately instead of solving for the frequency response in one step. The modal basis in Equation (61) is divided into k groups as
Each group contains roughly original modes from, and their corresponding improved modes. The eigenvectors of the new design are calculated using instead of T in Equations (62) to (64). The process is repeated times using a modal basis that is of the size of the original modal basis. All groups of eigenvectors are then collected together to calculate the frequency response of the new design. As demonstrated in Section 4.3.2, this reduces the cost considerably with minimal loss of accuracy.
3.5. Comparison of CA/MCA and PROM Methods
As we have discussed, a large overhead cost which is proportional to the number of design parameters is required before the PROM reanalysis is carried out. However, the CA/MCA method does not require this overhead cost because it does not need the basis of Equation (49) (see Section 3.2). Therefore, the CA/MCA method is more suitable, when the number of reanalyses is comparable to the number of design parameters. This is usually true in gradient-based design optimization. The CA and MCA methods can become expensive however, when many reanalyses are needed because, for each reanalysis, they require a new basis or (see Equations 56and 61, respectively) and new condensed mass and stiffness matrices in Equations (57) and (62). This is the case in gradient-free optimization problems employing a Genetic Algorithm for example, and in simulation-based probabilistic analysis problems employing the Monte-Carlo method. For these problems, the PROM method is more suitable because the subspace basis does not change for every new design point. Table 1 summarizes the main characteristics, advantages and application areas of the CA/MCA and PROM methods.
|CA/MCA Method||PROM Method|
|Overhead Cost||None||Required cost to construct.|
Cost proportional to the number of design parameters.
|Basis Vector||Variable basis/.|
Size proportional toand.
Size proportional toand.
Relatively small size of/.
Must recalculate/at every new design.
|High if no parametric relationship exists due to the condensation of large size and dense.|
Low if a parametric relationship exists.
|Best Application||Small number of reanalyses compared to the number of design parameters.|
Evaluation of few design alternatives and gradient-based optimization.
|Very large number of reanalyses.|
Gradient-free optimization (e.g. genetic algorithms) and probabilistic analysis.
4. Reanalysis Methods in Dynamic Analysis and Optimization
The reanalysis methods of Section 3 can be used in different dynamic analyses such as modal or direct frequency response and free or forced vibration in time domain. Depending on the problem and the type of analysis, a particular reanalysis method may be preferred considering how many times it will be performed and how many design parameters will be allowed to change. This section demonstrates the computational efficiency and accuracy of reanalysis methods in dynamic analysis and optimization. It also introduces a new reanalysis method in Craig-Bampton substructuring with interface modes which is very useful for problems with many interface DOFs where FRF substructuring is not practically applicable.
4.1. Integration of MCA Method in Optimization
We have mentioned that the MCA method provides a good balance between accuracy and efficiency for problems that require a moderate number of reanalyses, as in gradient-based optimization. For problems where a large number of reanalyses is necessary, as in probabilistic analysis and gradient-free (e.g. genetic algorithms) optimization, a combination of the MCA and PROM methods is more suitable.
Figure 6 shows a flowchart of the optimization process for modal frequency response problems. The DMAP (Direct Matrix Abstraction Program) capabilities in NASTRAN have been used to integrate the MCA method and the NASTRAN modal dynamic response and optimization (SOL 200). The highlighted boxes indicate modifications to the NASTRAN optimizer. Starting from the original design, the code first calculates the design parameter sensitivities in order to establish a local search direction and determine an improved design along the local direction. At the improved design, an eigen-solution is obtained to calculate a modal model and the corresponding modal response. The dynamic response at certain physical DOFs is then recovered from the modal response. At this point, a convergence check is performed to decide if the optimal design is obtained. If not, further iterations are needed and the above procedure is repeated. Many iterations are usually needed for practical problems to obtain the final optimal design. Section 4.3.2 demonstrates how this process was used to optimize the vibro-acoustic behavior of a 65,000 DOF, finite-element model of a truck. Using the MCA method, the computational cost of the entire optimization process was reduced in half compared with the existing NASTRAN approach.
As for a stand alone modal frequency response, the eigen-solution accounts for a large part of the overall optimization cost for vibratory problems where a modal model is used. A reanalysis method can be inserted into the procedure as shown in Figure 6 to provide an approximate eigen-solution saving therefore, substantial computational cost. Other reanalysis methods such as the CDH/VAO, CA or PROM can also be used depending on the number of design variables and the number of expected iterations.
4.2. Integration of MCA and PROM Methods
The PROM method requires exact calculation of the mode shapes for all designs corresponding to the corner points of the parameter space in order to calculate the subspace basis of Equation (49). The required computational effort can be prohibitive for a large number of parameters (optimization design variables). This effort can be reduced substantially if the modes of each corner point are approximated by the MCA method. In this case, an exact eigen-solution is required only for the baseline design. The following steps describe an algorithm to integrate the MCA and PROM methods:
Perform an exact eigen-analysis at the baseline design point p0 all parameters are at their lower limit, to obtain the baseline mode shapes.
Use the MCA method at design point pi with all parameters at their low limit except the parameter which is set at its upper limit. Obtain approximate mode shapes for the corner point using the following recursive process
Form the subspace basis T as
where is the total number of parameters.
The modal basis can be subsequently used in a modal dynamic response solution. Only step 4 is repeated in reanalysis. The computational savings can be substantial especially for problems where many reanalyses are needed.
4.3. Combined MCA and PROM Methods: Vibro-Acoustic Response of a Vehicle
The pickup truck vehicle model with 65,000 DOFs of Figure 7 is used in this section to demonstrate the advantages of the combined MCA and PROM method in optimizing the vibro-acoustic response of a vehicle. The model has 78 components and roughly 11,000 nodes and elements. The example is performed on a SUN ULTRA workstation using NASTRAN v2001. The MCA and PROM methods have been implemented in NASTRAN DMAP.
The sound pressure level at the driver’s ear location is calculated using a vibro-acoustic analysis. The structural forced vibration response due to unit harmonic forces in x, y, and z directions at the engine mount locations, is coupled with an interior acoustic analysis. The first and second eigen-frequencies of the acoustic volume inside the cabin are 95.9 Hz and 128.3 Hz. The sound pressure level is calculated in the 80 to 140 Hz frequency range. The structure and fluid domains are coupled through boundary conditions ensuring continuity of vibratory displacement and acoustic pressure. A finite-element formulation of the coupled undamped problem yields the following system equations of motion .
where the vibratory displacement and the acoustic pressure are the primary variables. The finite-element representation of the two domains consists of stiffness and mass matrix pairs and, respectively. The air density and wave speed are and. The right hand side of Equation (69) denotes the external forces.
The spatial coupling matrix indicates coupling between the two domains which is usually referred to as “two-way coupling.” Due to this coupling, the combined structural-acoustic system of equations is not symmetric. If the acoustic effect on the structural response is small, the coupling term can be omitted, resulting in the so-called “one-way coupling,” where the structural response is first calculated and then used as input (in Equation 69) to solve for the acoustic response. The coupled structure-acoustic system can be solved either by a direct method, or more efficiently by a modal response method which can be applied to both the structural and acoustic domains.
4.3.1. Combined MCA and PROM Methods
To demonstrate the computational effectiveness and accuracy of integrating MCA in PROM, a reanalysis was performed for a modified design where five plate thickness parameters vary; chassis and its cross links, cabin, truck bed, left door, and right door. All parameters were increased by 100% from their baseline values. The sound pressure at the driver’s ear was calculated using “two-way” coupling. A structural modal frequency response was used. The acoustic response was calculated using a direct method because the size of the acoustic model is relatively small. For the structural analysis, 1050 modes were retained in the 0 to 300 Hz frequency range. The combined MCA and PROM approach was compared against the NASTRAN direct solution for a modified design where all five parameters were at their upper limits. Only one iteration was used in Equation (59) in order to get the set of once-updated mode shapes for each corner design point. The subspace basis, which includes information for all five design parameters, is therefore, represented by
The maximum error in natural frequencies as predicted by the combined MCA and PROM method and NASTRAN, is less than 0.45% in the entire frequency range. Figure 8 indicates that the sound pressures calculated by both methods are almost identical. The computational effort for the MCA method to obtain approximate mode shapes at each corner design point is about 30 seconds. In contrast, it takes about 180 seconds for an exact eigen-solution using NASTRAN. The computational cost to construct the reduced basis (in PROM and in PROM+MCA) is compared in Table 2. The total cost was reduced from 1080 seconds to 330 seconds. The computational saving is more significant if the number of design parameter increases.
|Method||Solving for mode shape Φ 0 at baseline design||Solving for mode shapes at 5 corner design points||Total Cost|
|PROM||180 sec||180*5=900 sec||1080 sec|
|PROM+MCA||180 sec||30*5=150 sec||330 sec|
4.3.2. Optimization using MCA Method
The goal here is to minimize the sound pressure at the driver’s ear. A total of 41 design parameters are used representing the thickness of all vehicle components modeled with plate elements. All thicknesses are allowed to change by 100% from their baseline values. Table 3 describes all design parameters. At the initial point of the optimization process, all parameters are at their low bound.
|Prm.#||Description(thickness of)||Prm.#||Description(thickness of)||Prm.#||Description(thickness of)|
|1||Bumper||15||Radiator mtg.||29||Tire, front right|
|2||Rails||16||Radiator mtg., mid.||30||Tire, rear left|
|3||A-arm, low left||17||Fan cover, low||31||Tire, rear right|
|4||A-arm, low right||18||Fan cover, up||32||Engine outer|
|5||A-arm, up left||19||Cabin||33||A-arm conn., up left|
|6||A-arm, up right||20||Cabin mtg. reinf.||34||A-arm conn., up right|
|7||Tire rim||21||Door, left||35||A-arm conn., low left|
|8||Engine Oil-box||22||Door, right||36||A-arm conn., low right|
|10||Hood||24||Brake, front left||38||Glass, right|
|11||Fender, left||25||Brake, front right||39||Glass, rear|
|12||Fender, right||26||Rail conn., rear||40||Glass, front|
|13||Wheel house, left||27||Rail mount||41||Rail conn., front|
|14||Wheel house, right||28||Tire, front left|
Because of the large number of design parameters, the combined MCA and PROM approach Section 4.3.1 is not computationally efficient because the size of the PROM basis is very large (see Equation 70). For this reason, we use the MCA reanalysis method and demonstrate its capability to handle a large number of parameters. It approximates the mode shapes at intermediate design points using only in Equation (59). The subspace basis at each optimization step is thus. Because 1050 modes exist in the frequency range of 0 to 300 Hz of the initial design, the size of the MCA modal basis is 2*1050 = 2100.
|Eq. (59)||31 sec||31 sec|
|Eq. (62)||258 sec||50 sec|
|Eq. (63)||48 sec||6 sec|
|Eq. (64)||67 sec||10 sec|
|Total Cost||404 sec||97 sec|
The cost of solving for 1050 modes directly from NASTRAN is 180 seconds (see Table 2). In the MCA method, the cost of solving the linear system of equations in Equation (59) is 31 seconds, and the additional combined cost of Equations (62) to (64) is 373 seconds, resulting in a total cost of 404 seconds (see Table 4). To reduce this cost, the 1050 modes are divided into 21 groups and the modes in each group are obtained separately as explained in the last paragraph of Section 3.4. This reduces the cost of Equations (62) to (64) to 66 seconds for a total cost of 97 seconds, which is about half the cost of the direct NASTRAN method.
The gradient-based optimizer in NASTRAN (SOL 200) using the optimization process of Figure 6 needed three iterations to calculate the optimal design. Figure 9 compares the sound pressure at the driver’s ear between the optimal and initial designs. Figure 10 shows the percentage increase of optimal values relative to the initial values for all 41 design parameters. In the frequency range of 80-140Hz, the maximum sound pressure is slightly reduced from 7.9E-7 to 7.2E-7 Pascal. Most parameters are minimally changed. The largest increase is 20% for the rail mount thickness (parameter #27).
The Sequential Quadratic Programming (SQP) algorithm of NASTRAN can only find a local optimum. To obtain a more significant design improvement, two additional studies were performed using a Genetic Algorithm with the MCA method. The first study used 20 initial populations and 4 generations, and the second study used 100 initial populations and 6 generations. Figures 11and 12 show that the number of initial populations and the number of generations, affect the optimization results. While a higher number of initial populations and generations results in a slightly better result, both studies produced a much better optimum than the SQP algorithm. In the case of 100 initial populations and 6 generations, the sound pressure is reduced from 7.9E-7 Pascal to 2.0E-7 Pascal, which is equivalent to about 15 dB in sound pressure level (SPL).
To verify the accuracy of the MCA approximation, the sound pressure response at the optimal design from MCA+GA with 100 initial populations and 6 generations was evaluated by both direct NASTRAN and MCA. Figure 13 shows that the MCA method is very accurate. For a similar to MCA accuracy, the original CA method needed three sets of mode shapes to form the subspace basis, requiring 90 seconds to solve the linear equations. The much larger mode basis R in CA increases the computational cost to calculate the triple matrix products of Equation (57). Therefore for large scale, finite-element models with a high modal density, the proposed MCA method can be more efficient compared to either a complete NASTRAN analysis or the original CA method.
4.4. Reanalysis in Craig-Bampton Substructuring with Interface Modes
The FRF substructuring of Section 2.3 couples two structures using FRF information between the coupling (interface) DOFs, and the excitation and/or response DOFs. Although this approach is very efficient, it is practical only if we have a few coupling DOFs; e.g. connection of a vehicle suspension to chassis or connection of the exhaust system to body through a few hangers. If the physical substructures have interfaces with many DOFs, a different reduced-order modeling (ROM) approach must be used such as the Craig-Bampton ROM of Section 2.2.1. The Craig-Bampton ROM can be large however, if the number of retained interface DOFs is large. We address this problem by performing a secondary eigenvalue analysis which yields the so-called interface modes (see Section 2.2.2). The following section describes a reanalysis methodology for physical substructuring with Craig-Bampton ROMs using interface modes. We show that its accuracy is very good and the computational savings are substantial.
4.4.1. Craig-Bampton with Interface Modes and Reanalysis
In the Craig-Bampton CMS method (Craig-Bampton reduced-order model or CBROM), the mass and stiffness matrices of each substructure are partitioned into interface sub-matrices, interior (omitted DOF) sub-matrices, and their coupling sub-matrices. The dynamics of a structure are then described by the normal modes of its individual components, plus a set of modes called constraint modes that couple the components. In CBROM, there is no size reduction for constraint modes since all of them are kept in the reduced equations. If the finite element mesh is sufficiently fine, the constraint-mode DOFs will dominate the size of CBROM mass and stiffness matrices and result in a large computational cost. This issue is addressed by using interface modes (also called characteristic constraint –CC- modes). For that, a secondary eigenvalue analysis is performed using the constraint-mode partitions of the CMS mass and stiffness matrices. The CC modes are the resultant eigenvectors. Details are provided in Sections 2.2.1 and 2.2.2.
The number of constraint modes equals to the number of interface DOF. For many FE models of large structures, the number of interface DOF can be rather large. The calculation of constraint modes in Equation (6) involves a decomposition step and a FBS step. The cost of FBS is proportional to. For any matrix multiplication that involves, the cost is proportional to. For any triple-product that involves the cost is proportional to.
The matrices from all substructures are assembled into a global CBROM with substructures coupled at interfaces by enforcing displacement compatibility. If and are the component (substructure) matrices, the global matrices and are assembled as
and a secondary eigenvalue analysis is performed to calculate the interface modes as
where the matrices miCC , miCCN and kiCC are of much smaller size than matrices miC , miCN and kiC .
The interface modes reduce the interface size producing a smaller reduced order model (ROM) compared with the traditional Craig-Bampton ROM (CBROM). However, they are calculated from the assembled interface K and M matrices. Thus, the calculation of constraint modes and all matrix multiplications related to constraint modes are still necessary. The interface mode method reduces the size of ROM but it does not reduce the computational cost related to the constraint modes.
The following observations can be made:
In Equation (73), the number of columns of matrix is equal to the number of interface modes which is usually smaller than. Therefore, the FBS cost of solving for is proportional toand it is much smaller than the FBS cost of solving for.
Because both and are of size the cost of matrix multiplication and triple-product in Equations (74) to (76) are now proportional to and. Therefore, the cost is much smaller than the corresponding cost in Equations (9), (10) and (12).
In this CCROM method which is based on CBROM, the interface modes are obtained using the assembled interface partitions of the CBROM formulation. Thus, it is impossible to know before hand for a new design. For this reason, Equations (73) to (76) can not be theoretically implemented to improve efficiency. For this reason, we propose a reanalysis approach where the calculated interface modes Φ C C for original (baseline) design can be used as an approximation of the new interface modes at any modified design. In this case, Equations (73) to (76) are applied to improve the computational efficiency.
4.4.2. A Car Door Example
The car door model of Figures 14 and 15 is used to demonstrate the proposed reanalysis method for substructuring with Craig-Bampton method using interface modes. It has 25,800 nodes and 25,300 elements and is divided into two substructures. The first substructure includes the outer door shell and a bar attached to it. The second substructure includes the rest of the door. There are 293 nodes (1758 DOFs) on the interface. Therefore, the CBROM or CCROM method must calculate 1758 constraint modes according to Equation (6) for both substructures. The 1758 constraint modes are involved in matrix multiplication or triple-products in Equations (9), (10) and (12). Figure 16 shows the interface nodes.
For the initial design using the CCROM method, 52 interface modes are calculated below 600 Hz. A modified design is created where the shell thicknesses for the outer door (substructure 1) and inner door (substructure 2) are doubled. To provide baseline numbers, the CCROM method is used on the new design to solve for the system natural frequencies. The new reanalysis approach is used on the new design to calculate approximate natural frequencies which are then compared with the baseline numbers. The interface modes calculated at the original design are used as an initial guess for the interface modes of the new design.
Figure 17 compares the natural frequencies of the new (modified) design between the original CCROM (Craig-Bampton with Interface modes) method and the new approach where reanalysis is used in CCROM to approximate the interface modes. We observe that the natural frequencies of the modified design are very different from those of the original design. Also, the accuracy of the proposed reanalysis method is excellent. The frequencies for the modified design calculated by the original CCROM and the proposed new approach are almost identical. The percentage error of the new approach versus the original CCROM approach is less than 1% on average. The computation cost is summarized in Table 5.
|CPU Time||Normal Modes||Constraint Modes||Multiplication||Other Cost||Total Cost|
|CCROM||8 sec||61 sec||65 sec||3 sec||137 sec|
|New Approach||7 sec||2 sec||0.3 sec||2 sec||11 sec|
|CPU Time||Normal Modes||Constraint Modes||Multiplication||Other Cost||Total Cost|
|CCROM||108 sec||282 sec||927 sec||10 sec||1327 sec|
|New Approach||110 sec||16 sec||3 sec||10 sec||139 sec|
In the new approach to reduce the cost related to constraint modes, the total remaining cost is dominated by the cost of calculating the normal modes for each substructure. For example, the calculation of the normal modes for Substructure 2 took 110 seconds out of a total of 139 seconds (see Table 5). It should be noted that the normal modes cost can be further reduced by applying another reanalysis method such as CDH/VAO, CA or MCA to approximate the normal modes. Therefore, the overall cost of substructuring based on Craig-Bampton with interface modes, can be drastically reduced by using the proposed reanalysis to approximate the constraint modes and a CDH/VAO or MCA reanalysis to approximate the normal modes at a new design.
5. Optimization of a Vehicle Model
A detailed optimization study is presented using a large-scale FE model of a vehicle. For simplicity, we call it “BETA” car model. It is composed of approximately 7.1 million DOFs and 1.1 million elements. Figure 18 shows all modeling details.
We form an optimization problem in terms of the maximum vibratory displacement at any location of the outer shell of the two front doors by minimizing the maximum displacement among ten locations of the two front doors (Figure 19) due to a hypothetical engine excitation in the vertical (up-down) direction. The engine is represented by a lumped mass connected rigidly to the engine mounts (Figure 20). The powertrain-exhaust model has about 1.3 million DOFs and is composed of 29 PSHELL components and 12 PSOLID components. There are also some RBE2 and PBUSH elements which are used as connectors. The maximum displacement at each of the ten door locations is observed in the y direction (lateral direction – perpendicular to the door plane).
Fifteen design variables are chosen; five structural elements of each door (thickness of door shell, front frame, rear frame, top panel, middle pipe), vertical stiffness of each of the four engine mounts, and vertical stiffness of each of the six exhaust system supports. All design variables are schematically indicated in Figures 20 and 21.
The optimization problem is stated as follows:
The optimal value of each of the fifteen design variables is calculated in order to minimize the maximum response among the ten locations on the doors while the mass of the vehicle remains less or equal to the mass of the initial (nominal) vehicle. The response is calculated in the 100 Hz to 200 Hz frequency range and a 3% structural damping is used. The optimization problem is numerically very challenging because of
the many local optima and
the computational cost of each dynamic analysis.
The former was handled by using a hybrid optimization algorithm which first explores the entire design space using a Niching Genetic Algorithm (GA)  and then switches to a gradient-based optimizer (fmincon in MATLAB) using the best estimate of the optimal point from the GA as initial point. This ensures a rapid convergence to the final optimum because although all GA optimizers can move quickly to the vicinity of the final optimum, they have a very slow convergence rate in pinpointing the final optimum.
FRF substructuring is used to assemble all components of the vehicle (body, doors, and engine-exhaust) into a small reduced-order model. This keeps the computational cost of each dynamic analysis low (4 minutes per analysis). A modal model is created only once for the body subsystem and then used to generate an FRF representation. This modal model does not change during the optimization because the chosen design variables are not associated with the body. However, the modal models of the doors change during the optimization. The final model for the entire vehicle is created by assembling the FRF models of each component. The FRF assembly operation is repeated during optimization because the FRF models of the two doors keep changing.
The Niching GA optimizer maximizes a fitness function by modifying all design variables. A proper fitness function which minimizes the maximum response among the ten door locations while satisfying the vehicle mass constraint is chosen as follows
The ratio of the nominal maximum response over the actual maximum response is used so that the fitness value increases when the actual response is reduced. This ratio is multiplied by where p = 10 is a penalty value and. Thus, c is positive if is less than satisfying the constraint and the value of is equal to one.
Otherwise, c becomes negative if is greater thanand the term assumes a large negative value which reduces the fitness value considerably. As a result, the GA optimizer always satisfies the mass constraint while maximizing the value of the fitness function.
Figure 22 summarizes the optimization results by comparing the maximum door response between the optimal and initial designs. The optimizer determined that the maximum response occurs at location 9 (center of left front door of Figure 19) at approximately 105 Hz. Figure 23 shows that this represents a vehicle local mode involving motion of the doors only. At the optimal design the maximum response was reduced from the initial 10-3 m to 0.47*10-3 m (Table 6).
Table 6 compares the value of each design variable between the initial (nominal) and final optimal designs. It also indicates that all designed variables were allowed to vary within a lower and upper bound. The values of the five door design variables changed considerably between the initial and optimal designs. This is expected because the optimizer tried to suppress the local door mode. The stiffness of the four engine mounts and the six exhaust supports also changed. Although we intuitively expect the stiffness of the engine mounts to change but not the stiffness of the exhaust supports, this is not the case in this example. Table 6 also indicates that at the optimum we not only reduced the maximum response from 10-3 m to 0.47*10-3 m but the vehicle mass was also reduced from the initial 55.12 units to the final 51.92 units.
Figure 24 shows the actual function evaluations (design points where the vehicle dynamic response was calculated) in the X1-X2-X3 space and indicates the vicinity of the optimal design point. The GA optimizer needed only 359 function evaluations and used a population size of 5*(15+1) = 80 and a maximum of 10 generations. The population size and the number of allowed generations were kept at a minimum in order to locate the vicinity of the optimum quickly without “wasting” valuable computational effort.
Considering that the computational cost for each function evaluation was 4 minutes, the total computational time was (359 function evaluations) * (4 minutes per evaluation) = 1436 minutes or 23.9 hours. This is acceptable considering the size and type of performed analysis. The computational cost, in terms of number of function evaluations, was kept low by coupling the Niching GA with a Lazy Learning metamodeling technique [26, 27]. The latter estimates the value of the fitness function from existing values at close by designs without calculating the actual response. It uses an error measure to figure out if the estimation is accurate. The error is small if enough previous designs, for which the fitness value was evaluated, are close to the new design. In this case, the metamodel estimates the current fitness value without running an actual dynamic response. If the error is large, an actual response is calculated and the fitness value of this new design is added to the “pool” of previous designs the Lazy Learning metamodeling technique will use downstream.
6. Conclusions and Future Work
Reduced-order models and reanalysis methodologies were presented for accurate and efficient vibration analysis of large-scale, finite element models, and for efficient design optimization of structures for best vibratory response. The optimization is able to handle a large number of design variables and identify local and global optima.
For large FE models, it is common to solve for the system response through modal reduction in order to improve computational efficiency. An eigenanalysis is performed using the system stiffness and mass matrices and a modal model is formed which is then solved for the response. The computational cost can be also reduced using substructuring (or reduced-order modeling) methods. A modal reduction is applied to each substructure to obtain the component modes and the system level response is then obtained using component mode synthesis. In optimization of dynamic systems involving design changes (e.g. thicknesses, material properties, etc) the FEA analysis must be repeated many times in order to obtain the optimum design. Also in probabilistic analysis where parameter uncertainties are present, the FEA analysis must be repeated for a large number of sample points. In such cases, the computational cost is very high, if not prohibitive.
To drastically reduce the computational cost without compromising accuracy beyond an acceptable level, we developed and used various reanalysis methods in conjunction with reduced-order modeling, in optimization of vibratory systems. Reanalysis methods are intended to efficiently calculate the structural response of a modified structure without solving the complete set of modified analysis equations. We presented a variety of reanalysis methods including the CDH/VAO method, the Combined Approximations (CA) and Modified Combined Approximations (MCA) method, and the Parametric Reduced-Order Modeling (PROM) method. Their advantages and limitations were fully described and demonstrated with practical examples.
Future work will concentrate on developing reanalysis methodologies for shape and topology optimization of vibratory systems and extend the presented work in optimization under uncertainty where efficient deterministic reanalysis methods will be combined with efficient probabilistic reanalysis methods.