A Modal Perturbation Method for Eigenvalue Problem of Non-Proportionally Damped System

: The non-proportionally damped system is very common in practical engineering structures. The dynamic equations for these systems, in which the damping matrices are coupled, are very time consuming to solve. In this paper, a modal perturbation method is proposed, which only requires the first few lower real mode shapes of a corresponding undamped system to obtain the complex mode shapes of non-proportionally damped system. In this method, an equivalent proportionally damped system is constructed by taking the real mode shapes of a corresponding undamped system and then transforming the characteristic equation of state space into a set of nonlinear algebraic equations by using the vibration modes of an equivalent proportionally damped system. Two numerical examples are used to illustrate the validity and accuracy of the proposed modal perturbation method. The numerical results show that: (1) with the increase of vibration modes of the corresponding undamped system, the eigenvalues and eigenvectors monotonically converge to exact solutions; (2) the accuracy of the proposed method is significantly higher than the first-order perturbation method and proportional damping method. The calculation time of the proposed method is shorter than the state space method; (3) the method is particularly suitable for finding a few individual orders of frequency and mode of a system with highly non-proportional damping.


Introduction
Damping is one of critical factors affecting a structure's responses under dynamic excitations. In practical engineering structures, damped systems are commonly characterized non-proportional damping, such as the soil-structure system [1][2][3], the steel-concrete structures [4][5][6], and structures with supplemental dampers [7,8]. Even though the dynamic response of a non-proportionally damped system can be solved by the direct integration method, the mode superposition method is frequently used. One reason is that the natural vibration properties can help to realize the dominant frequency region. Another, is that few modes usually provide sufficiently accurate results, making the mode superposition method highly efficient. The important key contribution of mode superposition method is decoupling the coupled governing equations of motion in physical coordinates into a series of generalized single degree-of-freedom (DOF) equations of motion by mode shape coordinates. It is well known that traditional mode shapes based on the modal analysis of undamped systems can decouple undamped or proportional damped systems [9]. For the nonproportionally damped system, the traditional mode shapes cannot decouple the equations of motion, and an additional technique must be adopted [10].
Foss [11] proposed a state-space method to transform the equation of a damped linear dynamic system to an uncoupled set of complex modal system, but the increased computational cost hindered the wider application of this damped modal analysis. The equation of motion of the state-space method is double in dimensions compared to the traditional method, and it is necessary to solve the eigenvalue problem with complex-valued operation. Consequently, calculating the complex-value mode shapes in state space requires approximately eight times the numerical effort as the calculation of undamped real mode shapes [12]. Cronin [13] forced diagonalization of the transformed damping matrix by ignoring the off-diagonal terms. However, these solutions would cause unpredictable errors [14,15]. Furthermore, the equivalent uniform damping ratio method is a simple and approximate procedure for irregularly damped cases [16,17] such as the modal strain energy (MSE) method [18,19], but that method could significantly overestimate damping ratios for certain systems [20]. Ibrahimbegovic and Wilson [12], Udwadia and Esfandiari [21], Adhikari [22] proposed iterative solution methods which account for the effect of off-diagonal terms, and acceptable accuracy and efficiency have been observed.
The non-proportionally damped system can be viewed as having been modified from a damped system, so the dynamic characteristics and responses can be determined by the perturbation method [23][24][25][26]. Cha [27] used the modes of an undamped system as a basis to find the first-order perturbation solution for an arbitrary but weakly non-proportionally damped system. Lou and Chen [28,29] proposed a direct modal perturbation method (MPM) that treated the new system as a minor change based on the original system. Hračov [30], using the first-order perturbation method to evaluate the complex eigensolution of a linearly proportional system supplemented with a viscous damper. Tang and Wang [31] proposed a perturbation method which can deal with the undamped system with repeated eigenvalues. Those perturbation methods are usually suitable to analyze the weakly non-proportionally damped system [27,[31][32][33], and need the complete eigensolution set of the unperturbed system [27,30].
In this paper, a modal perturbation method for non-proportional damping problems applied to the vibrations of a highly non-proportionally damped system was proposed. In the method, an equivalent proportionally damped system, which is taken as the unperturbed system, is constructed from the real mode shapes of a corresponding undamped system and a non-proportional damping matrix. The equivalent proportionally damped system is close to the non-proportionally damped system to reduce the effect of the deviation damping matrix so that the method can deal with the highly non-proportionally damped system. Then, the characteristic equation of state space can be transformed into a set of nonlinear algebraic equations based on the MPM. The method does not require the complete eigenvector set of the undamped system, which is convenient and efficient for solving the large system. Several numerical examples are used to illustrate the validity and accuracy of the proposed modal perturbation method.

Theoretical Background
For a viscously damped linear system with N degrees of freedom (DOFs), the equation of motion for the forced vibration can be expressed as: It is well-known that the modal superposition method is a powerful method for evaluating the dynamic response of a linear system. For non-proportionally damped structures, the undamped mode shapes cannot uncouple Equation (1). To uncouple the equation of motion with nonproportional damping, the N second order differential Equation (1) has to be reformulated as 2N first order differential equations in state space [11]. By adding the identify the following state-space equation is derived from Equation (1):

Corresponding Proportionally Damped System
For a corresponding undamped system, the real eigenvalue problem can be expressed as: in which ωj and { } φ j are the jth undamped circular natural frequency and its associated real mode shape. Assume the undamped circular natural frequency is distinct, then the real modes satisfy the following orthogonality relation: in which δ is Kronecker's Delta function, Ml is the modal mass.
The non-proportional damping matrix can be split into the sum of two matrices: in which ζ l is the modal damping ratio. In fact,     p c is the superposition of modal damping matrix, which can be solved but need not to be explicitly determined in the proposed method. Therefore, the system composed of [m], [cp], and [k] can be called as equivalent proportionally damped system, which is close to the non-proportionally damped system to reduce the effect of the deviation damping matrix. By defining: in which ≤ n N is the first n modes determined from Equation (5). "  " denotes complex conjugation, On introducing the symbol Clk, defined by and making use of the orthogonality conditions of real mode { } φ , the complex mode { } η will satisfy the following relationships: for ≤ > , for > ≤ , for > > ,

Modal Perturbation Method
In the modal perturbation method (MPM), the jth eigenvalue and its associated mode of vibration with non-proportionally damped system is related to those of the corresponding proportionally damped system by: For the jth mode ( ≤ j n ), substituting Equations (22) and (23) into Equation (4) gives Premultiplying Equation (24)  gives: After has been introduced and l is taken from 1 to 2n, Equations (25)- (27) can be written in the matrix form: x are 2n × 1 vectors. The matrices and vectors are given below: denotes that the matrix     * is a diagonal matrix, the jth element of which is zero. For example, [2 ] l l diag s M denotes that the diagonal elements are 2 l l s M , . Equation (28) can be solved by the modified-Newton-Raphson iteration method. By assuming: in which The iteration solution of {x} can be obtained by: The initial solution can be set as: Iteration can be terminated when: where ε is a predetermined tolerance and usually is set to 1 × 10 −6 .
Once { } x is determined from Equation (28), the jth eigenvalue and associated complex mode of vibration of non-proportionally damped system can be obtained by Equations (22) and (23). Then, by defining: in which Re stands for the real part of the quantity that follows. The jth eigenvalue may be expressed as Equation (34) is the same as that eigenvalue of a viscously damped single DOF system with an undamped natural frequency ω Equation (23). So, the perturbation solution will converge to the exact solution. Usually, the effect of the higher order mode is lighter than the first lower order modes, and the mode superposition method only needs the first lower modes. Therefore, the n of Equation (23) is far less than N, which cause the dimension of Equation (28) to be small and reduce the time of calculation time significantly.

Two Degrees-Of-Freedom System
The first example is a 2 degrees-of-freedom system as shown in Figure 1, which had been extensively studied by many researchers [10]. The mechanical parameters of the system are chosen as follows:   (5) and the damping ratios ζ j from Equation (9) based on the equivalent proportional damping system. Meirovitch and Ryland [23] used the first-order perturbation method (FPM) to determine the complex eigenvalues. Cha [27] constructed a proportional damping matrix by least squares approach, then used the first-order perturbation method to calculate the complex eigenvalues. Those results are also showed in Table 1 and Table 2.
For non-proportionally damped systems, the matrix [C] is not diagonal. The relative maximum magnitude of the off-diagonal elements of [C] to the diagonal element can be expressed as [34]:  The coupling index α indicates how much the damping is non-proportional is shown in Table 3. For case 1 and case 2, coupling index α is equal to 0.237 and 0.708, respectively. Therefore, the damping of case 2 is highly non-proportional. The relative error e of the jth modal eigenvalues can be expressed into:

Frame with Concentrated Damper
Although the modal perturbation solution would converge to the exact solution based on the complete eigenvector set of the undamped system, usually only the first lower order eigenvectors can be obtained, especially for large-scale system with hundreds and thousands of DOF. Relatively few modes can provide sufficiently accurate modal superposition results in engineering applications, therefore the number n of undamped eigenvalues in Equation (23) is only a little larger than the truncated mode number. Similar to the subspace iteration method, to calculate the jth eigenvalue of the non-proportionally damped system, the number of modes n must be: As the additional number of modes, Δn, increases, the result gradually converges to the exact solution. The second example, a simple frame structure [12,35], is used to investigate the effect of Δn on the accuracy of eigenvalues. As schematically shown in Figure 2, the frame structure is idealized as 9 beam elements with 24 DOFs. The Young's modulus and density of material are 500 N/m 2 , and 1 kg/m 3 . The cross-section inertia and area are 1m 4 and 1m 2 . The length of each beam element is 1 m. Using the lumped-mass, the frame has 16 possible real mode shapes. The non-proportional damping arises from a damper attached to node 3 and the other attached to node 8. The damping coefficients of the two dampers are The coupling index α of the frame structure is equal to 1, which indicates that system is highly non-proportionally damped. The relative errors of the first three natural frequencies and damping ratios by MPM for various number of modes Δn are presented in Figure 3. It can be seen that the MPM leads to monotonically convergent solutions and gradually converges to the exact solution with the increase of Δn even if it is a highly non-proportionally damped system. It represents that the equivalent proportionally damped system reduces the effect of the deviation damping matrix so that the method can deal with the highly non-proportionally damped system. The curves of error are almost merging into one when Δn is greater than 8.  Then, the first three natural frequencies and damping ratios of the frame system by the proposed perturbation method with n = 11 and n = 16 are shown in Tables 4-6. As a comparison, those results obtained by exact method, PDM, FPM [23], and FPM [27] are also listed.
Obviously, the first three natural frequencies and damping ratios of perturbation method are same as the exact solution when n = 16. This is because all of the 16 modes of the frame system are the complete eigenvector set. The relative errors of the first 3 natural frequencies and damping ratios obtained by MPM are within 0.146% when n = 11, which are obviously less than the PDM and FPM. Therefore, the additional number of modes Δn = 8 is enough to obtain sufficiently accurate results. The errors of damping ratios of the other three methods may exceed 25%. Therefore, the proposed method can deal with the highly non-proportionally damped system for the equivalent proportionally damped system reduces the effect of the deviation damping matrix.   Table 6. The errors of natural frequencies and damping ratios of the frame system.  The main calculation time for the proposed method contains two parts, one is for solving the characteristic equation of undamped system and the other one is for solving the nonlinear algebraic equations. For the frame with concentrated damper, the first part consumed 0.057 ms, and the second part consumed 0.309 ms, the total time is 0.366 ms. The PDM only needed to solve the characteristic equation of undamped system which consumed 0.057 ms. To solve the characteristic equation of state space, the exact method consumed 0.470 ms. The FPM [23] needs to solve the characteristic equation of undamped system and to calculate the perturbation of undamped system which totally consumed 0.132 ms. The FPM [27] needs to optimize proportional damping coefficient and to solve the characteristic equation of undamped system. It totally consumed 0.123 ms. It can be seen that the time consumed by proposed method is shorter than state space method and longer than the other three existing methods. The advantages of the proposed method will become more prominent as the degrees of the system increase because it only needs to solve the characteristic equation of undamped system. Although the other three methods consumed less time, the accuracy of the proposed method is much higher, especially for high non-proportional damping systems. The proposed method outperforms the others when considering both the efficiency and the accuracy.

Conclusions
Non-proportional damping structures widely exist in practical engineering. In this paper, a modal perturbation method is proposed for evaluation of natural frequencies and damping ratios to significantly simplifies the solution process. Based on extensive analyses and numerical results, the following conclusions can be drawn: (1) The modal perturbation method transforms the characteristic equation of state space into a set of nonlinear algebraic equations by using the real mode shapes of the corresponding undamped system. The method only requires the first few lower eigenvectors of the undamped system, the calculation time is shorter than the state space method, which will simplify the solution and be convenient to solve the large system. (2) With the increase of vibration modes of the corresponding undamped system, the eigenvalues and eigenvectors monotonically converge to exact solutions. Usually j + 8 is enough to obtain sufficiently accurate natural frequencies and damping ratios. When all of the mode shapes of the undamped system are used, the modal perturbation method will obtain the exact solution. (3) The equivalent proportionally damped system is close to the non-proportionally damped system to reduce the effect of the deviation damping matrix so that the method can deal with the highly non-proportionally damped system. The method is particularly suitable for finding a few individual orders of frequency and mode for system with highly non-proportional damping.