The Applications of Order Reduction Methods in Nonlinear Dynamic Systems

Two different order reduction methods of the deterministic and stochastic systems are discussed in this paper. First, the transient proper orthogonal decomposition (T-POD) method is introduced based on the high-dimensional nonlinear dynamic system. The optimal order reduction conditions of the T-POD method are provided by analyzing the rotor-bearing system with pedestal looseness fault at both ends. The efficiency of the T-POD method is verified via comparing with the results of the original system. Second, the polynomial dimensional decomposition (PDD) method is applied to the 2 DOFs spring system considering the uncertain stiffness to study the amplitude-frequency response. The numerical results obtained by the PDD method agree well with the Monte Carlo simulation (MCS) method. The results of the PDD method can approximate to MCS better with the increasing of the polynomial order. Meanwhile, the Uniform-Legendre polynomials can eliminate perturbation of the PDD method to a certain extent via comparing it with the Gaussian-Hermite polynomials.


Introduction
Most of the actual engineering systems are high-dimensional and complex, the qualitative analysis for these systems is difficult and the numerical calculations are expensive. So, the original systems should be reduced to low-DOF ones, then the reduced models are used to represent the original ones. The common order reduction methods in the deterministic systems contain center manifold method [1,2], inertial manifold method [3,4], Lyapunov-Schmidt (L-S) method [5,6], Galerkin method [7,8], POD method [9][10][11][12][13] and other order reduction methods [14,15], Rega et al. [16,17] reviewed the order reduction methods in nonlinear dynamics systems.
Lu [20] applied the T-POD method to reduce the 7-DOF rotor model supported by a pair of ball bearings to a 2-DOF one, the reduced model reserves the dynamical characteristics of the original one, meanwhile, the stability of the reduced model is studied. Lu [21] applied the T-POD method to the rotor system with pedestal looseness at both ends, the order reduction condition and the optimal reduced model are provided, and the dynamical characteristics are also discussed in comparison to the model with looseness at one end. The T-POD method was also applied to the rotor system supported by cubically nonlinear stiffness and the bifurcation behaviors of the universal unfolding are analyzed [22]. However, the engineering systems usually involve uncertainties, such as design uncertainty, essential uncertainty, computational uncertainty, etc. The POD method can only be applied in the deterministic system, and it will lose efficiency to study high-dimensional problems with uncertainties. The MCS is an exact method to study the uncertain systems, but the computational time is very expensive, so order reduction in the stochastic systems should be discussed.
The numerical methods for order reduction in stochastic system include perturbation method [23][24][25], polynomial chaos expansion (PCE) method [26,27], PDD method [28], etc. In the recent twenty years, the PCE method is widely applied in a series of research areas [29,30]. The calculation quantity of the PCE method is more expensive than the PDD method and the error of the PCE method is larger [31]. The formulas of the first and second moment of the PDD method are provided by Rahman [32]. At present, the PDD method is widely applied in reliability-based design optimization, stochastic fracture mechanics [33,34], etc. Nevertheless, the PDD method is never used to study the dynamical problem, and it is a main research topic to apply the PDD method in the dynamical models with uncertainties for order reduction in stochastic systems.
The motivation of this paper is to summary the application of two order reduction methods in different types of systems. In Section 2, the T-POD method is introduced briefly, and the T-POD method is applied in the 24-DOF rotor model with pedestal looseness at both ends. In Section 3, the PDD method is applied in the 2-DOF spring system for the first time, the efficiency of the PDD method is verified via comparing it with the MCS method, Legendre polynomial is used to reduce the perturbation around the resonant frequency based on the discussion of amplitude-frequency response. At last, the conclusions and future work are drawn in Section 4.

Remark 1
One of the novelty of this paper is to generalize the T-POD method to the more complex rotor system model (looseness fault at both ends of the bearings), this 24-DOFs rotor system model is new and more complex than the 16-DOFs rotor system model in [21], the new model can further verify the accuracy and efficiency of the POD method. The second is to apply the PDD method of different polynomials to the mass spring system model and efficiency of the PDD method is also verified.

The T-POD Method in Deterministic System
In this section, the T-POD method in the multi-DOF dynamical system is introduced first, and a 24-DOF rotor model is established by the Newton's second law. The T-POD method is applied to the rotor system and provides the optimal order reduction condition. Finally, the bifurcation behavior is used to verify the efficiency of the T-POD method.

The T-POD Method
In general, the dynamical equation of a multi-DOF system can be written as: whereZ represents the vector of vibration displacement, C is the equivalent damping matrix, K is the equivalent stiffness matrix, Fðt;ZÞ is the nonlinear function about time t andZ.
The basic processes of the T-POD method are expressed as follows: (1) Provide the initial conditions, extract displacement and velocity information from transient process of each DOF, which is denoted by z 1 ðtÞ; z 2 ðtÞ; …z m t ð Þ. m is the DOF number of the system, transient time interval displacement sequence of all DOFs is The interval is equal to each other and the number of time interval is N , which can be denoted as t 1 ; t 2 …; t N . The time series form the matrix X ¼ z 1 t ð Þ; z 2 t ð Þ; …; z m t ð Þ ½ , the order of matrix X is N Â m, the matrix X can be expressed as formula (2): . . .
The self-correlation matrix T ¼ X T X can be obtained and the order is m Â m. The eigenvectors and the corresponding eigenvalues of the matrix can be denoted as (2)Q is the matrix formed by the first n order eigenvectors of the matrix T, and it can also be indicated thatQ contains the first n largest eigenvalues of T ¼ X T X. The order of matrixQ is m Â n, so the order of Q TQ is n Â n. Taking the coordinates transformation on system coordinatesZ, the new coordinates P are obtained,Z ¼QP, substitutingZ into Eq. (1), we can get Eq. (3): The column vector ofQ is nonzero and orthogonal, the matrixQ TQ is n order diagonal and full rank, so the inverse matrixQ TQ exists. Taking ðQ TQ Þ À1Q T left multiplication on both sides of Eq. (3), and we obtain Eq. (4): Setting C n ¼ ÀðQ TQ Þ À1Q T CQ; K n ¼ ÀðQ TQ Þ À1Q T KQ; F n ¼ ðQ TQ Þ À1Q T F, we can get Eq. (5): As mentioned above, the T-POD method is used to reduce the m DOF system to an n DOF one, where n < m. The reduced dynamical equation is expressed as Eq. (5).

Dynamical Equation of the Rotor Model
The 24-DOF rotor model with looseness at both ends is shown in Fig. 1. The pedestal looseness fault model is discussed in the author's previous work [19][20][21][22] and the work of other researchers [35,36]. The dynamical equation of the 24-DOF system is established by the Newton's second law, the equation form and the corresponding parameters are similar as the 23-DOF rotor model with looseness at one end [19], Figure 1: 24-DOF rotor model with looseness at both ends only one more loose mass at right end is considered, so we don't provide the detailed equation here. All the parameters of the rotor system are authentic and available in engineering.
The optimal order reduction conditions are confirmed by proper orthogonal mode (POM) energy method. We consider that order reduction efficiency is best when the energy of the first few order POMs is maximum and approximate to 100%. The energy of the POMs varies with the rotating speed, the POM energy formula is denoted in formula (6). The number of l is the optimal modal number when the energy is approximate to 100%. This method can confirm the DOF number of the reduced model and provide the rotating speed corresponding to the largest POM energy, see details in [21]. Fig. 2 reflects the rotating speed of y 1 DOF varies with the energy. When the rotating speed is x ¼ 1120 rad=s, the corresponding five order POMs hold the highest energy. The optimal order reduction condition is obtained and the reduced model is a 5-DOF one.

The Results
Provide the initial conditions of the system: the integral is p=256, the initial displacement and velocity are Þ , the rotating speed is x ¼ 1120 rad=s. The time history curve of y 1 DOF is shown in Fig. 3, when s 2 0; 120p ½ , the system is in transient process, and the system is in the periodic motion state after 120p.
The T-POD method is applied to reduce the original rotor model to a 5-DOF one, and the bifurcation behavior is discussed. In Fig. 4, the bifurcation point of the order reduction system is consistent with the original system, and the topological structure of the reduced system is similar to that of the original system. So the reduced system reserves the bifurcation behavior of the original one, which verifies the efficiency of the T-POD method in complex rotor system.   The application of the T-POD method in the high-dimensional rotor system is discussed in this section. This method extracts the transient information of dynamical system first, and then the POM energy method is used to confirm the POM number of the reduced model, finally, the original model is projected to the reduced one. The 24-DOF rotor model is only an example, the T-POD method can be applied to the rotor systems with different faults.

The PDD Method in Stochastic System
This section introduces the general theories of the PDD method first, and this method is applied to the 2-DOF spring system. The numerical results indicate that the PDD method can provide accurate, convergent and computationally efficient estimation.

The PDD Method
Let ; F; P ð Þ be a complete probability space, is a sample space, F is a rÀfield on , and P : F ! 0; 1 ½ is a probability measure. B N represents the Borel rÀfield on R N , define the independent random vector X ¼ Þ that describes input of the complex random system. The joint probability density function (PDF) of X is defined as Let y X ð Þ be a real, square-integrabel and measurable transformation on ; F ð Þ, define a relevant performance function of a stochastic system. In general the multi-variate function: y : R N ! R is implicit, the analytical form can't be obtained, and can only be considered as a high-dimensional input-output mapping.

Dimensional Decomposition
The common dimension decomposition methods include referential dimensional decomposition (RDD), analysis of variance dimensional decomposition (ADD), polynomial dimensional decomposition (PDD), etc. They are introduced as follows: Þat a reference point c 2 R N as the kernel function, which leads to: X v ; c Àv ð Þis an N dimensional vector, when i 2 v, ith component is X i , if i = 2 v, it is c i . Both the recursive and explicit forms of Eq. (7) exist: recursive form presents as cut-HDMR [37], and explicit form presents as dimensional decomposition [38], the independent two forms were proved to be equivalent [39]. Even so, the (2) ADD ADD regards the PDF f X x ð Þ of X as the kernel function, which results in [40]: where y 0 is the expansion coefficient. There exists the recursive form of Eq. (8). The ADD has some synonyms, notably, Sobol decomposition, which was generalized by Sudret [41], etc. The ADD can provide desirable orthogonal characteristics, the component function is difficult to obtain, which requires the calculations of high-dimensional integrals.
is a set of orthogonal polynomial basis function in the Hilbert space L 2 i ; F i ; P i ð Þ , and it is consistent with the probability measure P i of X i , the ADD can be expanded as PDD as follows: is an additional expansion coefficient, the high-dimensional integral is required, 1 u j j N . PDD entails the orthogonal component function and considers the smoothness of y, for efficiently calculating the probabilistic characteristics.

Orthogonal Polynomial Basis
Let w ij x i ð Þ; j ¼ 0; 1; Á Á Á n o be univariate set of orthogonal polynomial basis function in the Hilbert space. L 2 i ; F i ; P i ð Þ , and is consistent with the probability measure P i . Let X i follow the PDF f i x i ð Þ with support where dx i is the Lebesgue measure. Define the space of square integrable functions with measure P i as: The inner product space associated with the random variable X i can be defined as: where y i x i ð Þ and z i x i ð Þ are two univariate functions of x i , E P i y i 2 X i ð Þ ½ < 1, E P i represents the expectation operator with respect to the probability measure P i . L 2 i ; F i ; P i ð Þcan also be interpreted as the space of functions of random variables with finite second moments.
Consider w ij x i ; j ¼ 0; 1; . . . ð Þas a set of orthonormal basis in the Hilbert space L 2 i ; F i ; P i ð Þ , which satisfies: Eqs. (12) and (13) represent that the basis functions have a zero mean, unit variance, and are mutually orthogonal when j; j 1 ; j 2 ! 1. Many traditional polynomials satisfy these characteristics, such as Hermite, Legendre polynomials, etc.

Parameter Calculation
On the basis of previous work of Xu et al., the low-dimensional integrals can be applied to estimate the high-dimensional ones, and to calculate the corresponding parameters [38]. Consider y i 1 ;ÁÁÁ;i s x i 1 ; Á Á Á ; x i s ð Þis SÀvariate component function, the Fourier expansion exists: the Fourier approximation can be obtained: Eq. (15) approximates to y i 1 ;ÁÁÁ;i s x i 1 ; Á Á Á ; x i s ð Þwhen m ! 1. If the component function y i 1 ;ÁÁÁ;i s x i 1 ; Á Á Á ; x i s ð Þ is smooth, the approximate convergence speed will be fast. The convergent rate is complex in normal conditions, such as the multi-DOF dynamical system, the rotor system, etc.
Setting S ¼ 1, 2, and 3, the univariate, bivariate, and trivariate component functions, are respectively shown as follows: The dimension reduction method is applied to calculate the corresponding parameters, the formulas and parametric expressions can refer to the work of Xu et al. [28,38]. We only provide the parametric calculation formula of univariate component function, the bivariate and trivariate expressions will not be discussed in details here.

Statistic Moments Analysis
The PDD of y X ð Þ can be expressed as follows [42,43]: which can be viewed as a finite, hierarchical expansion of an output function in terms of its input variables with increasing dimensions. y 0 , C i 1 ;ÁÁÁ;i s j 1 ;ÁÁÁ;j s ; s ¼ 1; 2; Á Á Á ; N are the expansion coefficients. Eq. (21) can be approximated as the sum of SÀvariate component functions, which contain m order orthogonal polynomials, where 1 m 1, 1 S N are all integers, so SÀvariate approximation is: y S X ð Þ will converge to y X ð Þ when S ¼ N ; m ! 1. The formulas of first and second order moments are provided as follows: Formulas of high order moments are discussed in details in [32].

Uncertainties in Dynamical System
In this section, the method to calculate the amplitude-frequency response of the general multi-DOF system will be provided first, and uncertainties in the dynamical system will be discussed. A general dynamical equation of system can be expressed as: where M, C, K, F t ð Þ are the matrixes of mass, damping, stiffness and external force, y t ð Þ is the steady response of the system. We consider the external force is harmonic, it can be denoted as F t ð Þ ¼ F 0 e ixt , so the response can be expressed as y t ð Þ ¼ Ye ixt , Y is the solution of Eq. (26): The amplitude of response Y is Y 1 þ iY 2 j j , Y 1 and Y 2 are the real and imaginary part of Y.
Consider the stiffness K is uncertain, which can be expressed as: where K is the deterministic matrix, which represents the mean stiffness matrix; K is controlled parameter d K , the deterministic system corresponds to d K ¼ 0. This paper studies the case that n K meets the Gaussian distribution (mean is 1, standard deviation is 0.05), we should guarantee that d K is very small, the uncertain stiffness is positive, the actually physical significance is guaranteed. This uncertainty is very simple, the highlight is not to find a new method to derive the moments, but to apply the PDD method to discuss the dynamical problem. Actually, the mass, damping, and external force can all be uncertain, we only apply the PDD method to study the statistical characteristics of the amplitude-frequency curves in the case of uncertain stiffness. In this paper, the solutions of the PDD method are compared with the accurate ones obtained by the MCS method to verify the efficiency of the PDD method.

Two-DOF Spring Model
2-DOF spring model with stochastic stiffness coefficients is shown in Fig. 5, where m ¼ 1 kg, In Fig. 5, the uncertain stiffness is k ¼ k 1 þ d k n k ð Þ . The mass, damping, stiffness and external force are as follows The frequency response equation is expressed in formula (26), the univariate method is used to solve this problem, for there is only one uncertain quantity in the system. We choose the Hermite polynomial to approximate the response, for the uncertain quantity meets the Gaussian distribution. The first and second order moments can be calculated by the PDD method, see the formulas and parameters in details in Sections 3.1.3 and 3.1.4.

The Results
We study the first and second order moments of amplitude-frequency response of 2-DOF spring system based on the formulas in Section 3.1.4. The MCS method is used for 10000 independent samples of random variable, the first two order moments are shown in Fig. 6, where Fig. 6a also contains the case of deterministic response.
In Figs. 7 and 8, the amplitude-frequency responses of order 2 and 9 of polynomials are calculated, and compared with the MCS method. The PDD results agree very well with the MCS results except for the  The Hermite polynomial based on the Gaussian distribution is used, the perturbations occur near the resonant frequencies. Then the Legendre polynomial based on Uniform distribution will be applied to discuss the statistical characteristics of the amplitude-frequency responses.
Figs. 9 and 10 show the results of PDD method when the order of Legendre polynomial is 2 and 9, respectively. In Figs. 9a and 10a, the first moment of amplitude-frequency response meets well with the MCS results. The perturbations disappear when the order increases to 9, the PDD results reserve the behaviors of the first moment. The second moment has smaller perturbations in comparison to Hermite polynomial, but the approximation is less efficient. As a whole, the perturbations around the resonant frequencies are solved preliminarily by the PDD method based on the Legendre polynomials.

Remark 3
The PDD results can agree well with the MCS results, the results are better as the polynomial order increases. Legendre polynomials are used to compare with the Hermite polynomials based on the amplitude-frequency characteristics. In a certain extent, Legendre polynomials eliminate the perturbations around the resonant frequencies. It is not clear that why the Legendre polynomial can reduce the perturbations, we will finally solve the problem in the future.

Conclusions and Outlooks
Two order reduction methods have been studied in this paper, one is the T-POD method in the deterministic system and the other is the PDD method in the stochastic system. The T-POD method has been applied to the 24-DOF rotor system with looseness at both ends, POM energy method has provided the optimal order reduction condition of the rotor system. The reduced model has reserved the bifurcation behaviors of the original one, which verifies the efficiency of the T-POD method in the deterministic system. The PDD method has been applied to the 2-DOF spring system with uncertain stiffness, the accuracy of the PDD method has been verified by comparing with the MCS method. The PDD results approximate to the accurate solutions better as the polynomial order increases. The perturbations around the resonant frequencies can be reduced by using the Legendre polynomials. Further studies on this research topic will be carried out by the present authors in two aspects: one is to use other methods to eliminate the perturbations around the resonant frequencies of the amplitude frequency response completely, and the other is to apply the PDD method in the complex rotor systems.