Adaptive Reduced Basis Methods Applied to Structural Dynamic Analysis

The reduced basis methods (RBM) have been demonstrated as a promising numerical technique for statics problems and are extended to structural dynamic problems in this paper. Direct stepby-step integration and mode superposition are the most widely used methods in the field of the finite element analysis of structural dynamic response and solid mechanics. Herein these two methods are both transformed into reduced forms according to the proposed reduced basis methods. To generate a reduced surrogate model with small size, a greedy algorithm is suggested to construct sample set and reduced basis space adaptively in a prescribed training parameter space. For mode superposition method, the reduced basis space comprises the truncated eigenvectors from generalized eigenvalue problem associated with selected sample parameters. The reduced generalized eigenvalue problem is obtained by the projection of original generalized eigenvalue problem onto the reduced basis space. In the situation of direct integration, the solutions of the original increment formulation corresponding to the sample set are extracted to construct the reduced basis space. The reduced increment formulation is formed by the same method as mode superposition method. Numerical example is given in Section 5 to validate the efficiency of the presented reduced basis methods for structural dynamic problems.


Introduction
Nowadays structural dynamic problems are usually solved by the finite element technique.Solution of dis-placement responses of all the nodes requires great effort.The scale and complexity of dynamics problems of practical engineering structure are ever increasing such that it requests more memory and computing time than before.Despite of the continuing advances in computer speeds and hardware capabilities, the dimension for numerical simulation is too large to provide real-time response in the design, optimization, control and characterization of engineering components or systems.Thus there are many motivations to develop methods that can not only reduce significantly the problem size and computational cost but also retain the accuracy of the solution and the physics of the structures.
Model order reduction techniques [1]- [11] have been proposed to reduce the size of a large-sized model before a detailed analysis performed.They are widely used in global-local analysis, reanalysis and structural dynamic optimization, eigenvalue problem, structural vibration and buckling, sensitivity studies and control parameter design, model update, and damage detection.A detailed review on model reduction techniques can be found in Noor [12].These reduction methods usually include two steps.The first step is the classic finite element discretization; the second is the computation of some basis vectors in order to perform a Rayleigh-Ritz analysis.Clearly, the success of the method depends chiefly on the proper selection of the basis vectors.
However, order reduction has long been focused on control problems [4] [13] [14]; most of the reduction methods in that field are designed for small or moderate-size systems and cannot be directly applied in the largescale case.Nevertheless, the reduced model cannot retain all features of the full model due to the truncated errors.Even for features within an interested frequency range, they may not be exactly kept in the reduced model resulting from most of the model reduction techniques.In recent years, the requirement of reduction techniques for large-scale systems has triggered a revival of research activities in model order reduction [13] [15] [16].Many powerful reduction techniques have been devised, in particular for linear time-invariant systems.Despite this progress, there are still many open problems.
Different from the traditional reduction methods, the reduced basis method (RBM) [17]- [20] is a very promising method which requires a projection onto the parameter-induced reduced basis space, as makes it very suitable for the analysis of large-scale system.The RBM has first been introduced for single-parameter problems in nonlinear structural analysis in the late 1970s and subsequently developed for multi-parameter problems.However, RBM rarely has been extended to perform model reduction in the structural dynamic problems yet.
In this paper we adopt the reduced basis method to perform the dynamic analysis of structures based on mode superposition method and direct integration method, respectively.A greedy algorithm is suggested to perform the adaptively selection of reduced basis vectors.Numerical example of a simplified one-dimensional seismic model is presented to demonstrate the feasible application of reduced basis method in structural dynamic problems.The error of the reduced system is evaluated numerically.

Theoretical Background
In structural dynamic analysis, the equations of motion are generally written as a set of linear second-order differential equations.The matrix form of these equations may be expressed by: where d , d  and d  are the acceleration, velocity, and displacement response vectors of the nodes, respec- tively, in the total Cartesian coordinate system.The upper dot means derivative with respect to time; F is the equivalent force vector acting on the nodes; the total mass matrix , where N is the shape function, B is the strain matrix, and D is the elasticity matrix.In the following analysis, the structure is subjected to initial conditions given by 0 0 ,

The Reduced Basis Method Applied to Dynamic Problems
In the following analysis, the stiffness and mass matrices are assumed as parameter-decomposition forms ( ) In Equation (3) and Equation (4), , P Q are the numbers of stiffness matrix and mass matrix that can be decomposed, respectively.They are determined by the problem itself.
The damping matrix is considered to be proportional.

Brief Introduction of Mode Superposition Technique
The mode superposition method can be used to perform a time history analysis to obtain the response of structure due to a transient loading as a function of time.It requires the solution of Equation ( 6) for the frequencies and mode shapes.
where mode shapes Φ can be shown to be orthogonal to the mass and stiffness matrices, as permit the equa- tions of motion to be uncoupled.
where ( ) i r t is the loading of the ith order mode, i ξ is the damping ratio for the ith mode and i ω is the fre- quency of the ith order mode as following., 2 Equation ( 11) can be solved by a procedure for solving single-degree-of-freedom dynamic problems.
It should be mentioned that the higher mode shapes of the system are unimportant for a practical engineering structure or component.Neglecting the higher frequencies and mode shapes of the system generally does not introduce significant errors.Thus modal truncation is often considered to reduce the computational effort when the number of DOF is large.

Reduced Basis Method to Generalized Eigenvalue Problem
Before the application of reduced basis method, a sample set of parameter domain is selected in a training space, which comprised of parameters spanning the parameter domain roughly.

{ }
1 , , The truncated eigenvectors corresponding to the parameters in the sample set are extracted to construct the reduced basis space where m is the number of mode is retained in terms of the required accuracy.
It should be noted that the basis vectors are the solutions of the system equations at different parameters.They are perhaps nearly oriented in the same direction.Consequently, the resulted reduced system is very ill-posed especially for large N , i.e. the condition number of the reduced stiffness matrix grows exponentially with N .To guarantee the basis vectors' linearly independence and make the reduced system well-posed, QR decomposition is used to generate a new basis which is orthogonal and able to retain all approximation properties of the original basis , , , The corresponding transform matrix is ( ) Then, the eigenvectors corresponding to a new parameter can be expressed as a linear combination of the basis vectors The above equation also can be rewritten in a matrix form ( ) To get the reduced system, the parameter-independent matrices are projected onto the reduced basis in terms of a Galerkin form.

( )
T 1, 2, , From this parameter-decomposition expression, the reduced system can be easily obtained and explored in the whole parameter domain.
( ) Obviously, the reduced eigenvalue problem can be solved more efficiently for each new parameter in test parameter-space.
The truncated eigenvectors can be regenerated approximately by The approximation of eigenvalues can be demonstrated in terms of Rayleigh's quotient.

The Adaptively Selection of Basis Vectors for Reduced Generalized Eigenvalue Problem
The basis vectors selection is critical for the efficiency and accuracy of the reduced basis method.Too many or too few vectors selected should be avoided.The former results in computational inefficiency, while the latter in unacceptable error.To obtain an appropriate basis space, a greedy algorithm is suggested to select the vectors adaptively.
At first, the error in approximated eigenvalues is presented.
( ) The maximum error in the training space is definite as The performing procedure of greedy algorithm is summarized as follows.
Step 1.One parameter in the training space is selected as the start point; the associated truncated eigenvectors are extracted as the vectors of the reduced basis space.
Step 2. QR decomposition is applied to perform orthogonalization of basis vectors.
Step 3. The reduced generalized eigenvalue problem is solved in the training space to yield the approximated modes ˆi λ , ˆi Φ .
Step 4. The maximum error max e is determined.
Step 5.The truncated eigenvectors corresponding to the maximum error will be selected as the next basis vectors and added to the reduced basis space.Then steps 2 to 4 are repeated.The greedy algorithm will terminate when the maximum error is less than a prescribed tolerance ε .

Brief Introduction of Direct Integration Technique
Direct integration provides a step-by-step numerical procedure to solve the equations of motion in Equation ( 1) directly without prior transformation of the equations to a different form.It can compute an approximate solution at discrete time intervals 0, ,2 ,3 , , , , , where T is duration of the input motion or loading and t ∆ is the time step.The widely used explicit methods are only conditionally stable such that some re- strictions over the size of the selected time step.On the other hand implicit methods may be unconditionally stable, but the computational work and storage requirement per time step can be much greater than explicit methods because a solution of a coupled algebraic system is always involved.Generally, explicit methods are most preferred for wave propagation problems, while implicit methods are widely employed and advocated for structural dynamic problems.
Newmark method is considered as the example.It is a widely employed linear one-step implicit method with two basic assumptions ( ) ( ) ( ) ( ) ( ) and is unconditionally stable under the following parameter limitation.
From the initial condition given by Equation ( 2), the initial acceleration given by ( )

Reduced Basis Method Based on Direct Integration Technique
Just as the same in mode superposition, a sample set will be introduced from a training space comprised of a span of parameters and all time steps.The reduced basis space is defined as the span of N finite element dis- placement response and doesn't change with time for a new arbitrary parameter.
, , , , For the same reason mentioned in foregoing section, QR decomposition is applied to generate an orthogonal reduced basis space.

{ }
1 , , The transform matrix for projection can be written as: ( ) The displacement response corresponding to new parameter and new time step can be approximated as the linear combination of the vectors in the reduced basis space.
It also can be expressed as a matrix form.

( ) ( ) ( )
The approximated velocity and acceleration can be obtained by first order and second order derivatives of the approximated displacement response with respect to time, respectively.

( ) ( ) ( )
The reduced Newmark formulation can be obtained by Galerkin projection of original space onto the reduced basis space.
The reduced stiffness, mass and damping matrices are respectively given by ( ) where the reduced parameter-independent matrices are ( ) T 1, 2, , The reduced load vector is The initial condition corresponding to the reduced system is ( ) ( )

The Adaptively Selection of Basis Vectors for Reduced Newmark Formulation
Similarly, the greedy algorithm is adopted to select the vectors adaptively and subsequently obtain an appropriate reduced basis space.At first, the projection error is defined in the training space as , , , , The maximum norm of the projection error is defined as The perform procedure of greedy algorithm is summarized as follows.
Step 1.To span the training space, the displacement of the last time step is selected as the first basis vector, corresponding to one source in the training space.
Step 2. QR decomposition is applied to perform orthogonalization of basis vectors.
Step 3. The reduced Newmark's formula is solved in the training space to yield the reduced basis displacements N d .
Step 4. The maximum norm of the projection error Π e is determined.Step 5.The displacement corresponding to the maximum norm of the projection error will be selected as the next basis vector and steps 2 to 4 are repeated.The greedy algorithm will terminate when the maximum norm of projection error is less than a prescribed tolerance ε .

Numerical Model
A simplified one-dimensional seismic model [21] is presented to numerically validate the application of reduced basis method to structural dynamic problem.The pressure variable P during an earthquake is governed by dynamic equilibrium equation: The earthquake source S and occurring time T are considered as system parameters, and denoted by = .The pressure P changes with the occurring time T and doesn't change in spatial dis- tribution such that it is fixed as 0.5.The pressure is zero in the earth's crust and the pressure gradient is zero on the earth's surface.The initial condition is given by To obtain parameter-decomposition forms of stiffness and mass matrices, namely, to express the stiffness and  mass matrices as the combination form of product of system parameter function and the matrix independent of system parameters.The original x-domain is decomposed into the left zone ( ) , right zone ( ) 3 S Ω and output zone ( )

S Ω
. A standard y-domain is introduced as reference and decomposed into where M and K are 200 200 × matrices, i M and i K are both independent of the system parameters, the parameter-dependent coefficients are: ( )

The Numerical Results
As the reduced structural dynamic analysis performed by using mode superposition, the 12 th truncation of mode is considered.It can be found from Figure 4 that the maximum error of approximated eigenvalue decreases rapidly with the increasing of the basis vectors.For Newmark integration case, the numerical parameters are selected as 200, 0.01 and the same convergence phenomenon as mode superposition can be found in Figure 5.However, it is obvious that the former converges sooner than the later.
The resulted reduced eigenvalue problem is 60 in dimensional, while the reduced Newmark formulation is 85 in dimensional for a prescribed error tolerance .Figure 6 and Figure 7 show that both reduced mode superposition and reduced Newmark integration approximate the original algorithms very well for engineering analysis.The CPU time for the reduced system and the original system are given in Table 1.The original Newmark costs expensively CPU time, while the reduced Newmark gives the lowest cost for test parameter       space.Despite the original mode superposition more effectively executed than the original Newmark method, the reduced form of the former costs more CPU time than the reduced form of the later.It can be concluded that the dynamic analysis have been performed much more effectively by either reduced mode superposition or reduced Newmark method.
It should be point out that the dimensional of the reduced system is determined by the reduced basis space and independent of the original system.For larger dynamic system, the efficiency of the reduced basis methods can be further enhanced.

Conclusion
Two kinds of reduced basis methods for dynamic problems are proposed in this paper.In the numerical example, the direct integration for the dynamic analysis is not numerically efficient as compared with the mode superposition method using eigenvectors due to the linear property of the seismic problem.However, it proves that the reduced basis method is available for structural dynamic analysis based on either mode superposition or direct integration.Though the undamped case studied, the reduced basis method can be applied to damped structures without any more effort.Furthermore, although the reduced Newmark method is only considered here, the reduced basis method can be very easily extended to other direct integration techniques.
displacement, velocity, and acceleration at time m t , the algorithmic values at time 1 m t + are the solution of the linear algebraic equations.

(
distribution and the temporal characteristics of earthquake source respectively, are showed in Figure1and Figure2, respectively.The spatial domain is divided into linear elements and normalized to unit length,

Figure 1 .
Figure 1.The time history of earthquake source.

Figure 2 .
Figure 2. The spatial distribution of earthquake source.

3 Ω
mapping from the standard y-domain to the original x-domain is given in Figure3: ; and the identity mapping from4  Ω to4  Ω .The resultant parameter-decomposition matrices are

Figure 3 .
Figure 3.The affine mapping from y-domain to x-domain.

Figure 4 .
Figure 4. Maximum norm of projection error changing with the increasing of number of basis vectors.

Figure 5 .
Figure 5. Maximum norm of projection error changing with the increasing of number of basis vectors.

Figure 6 .
Figure 6.Maximum norm of projection error changing with the increasing of number of basis vectors.

Figure 7 .
Figure 7.The pressure in the output domain corresponding to S = 0.5.

Table 1 .
CPU time comparison of mode superposition and implicit Newmark method.