Improved complex mode theory and truncation and acceleration of complex mode superposition

Abstract

The traditional complex mode theory is substantially improved. By means of the conjugacy of the complex characteristic value and vector, the formulae to compute the accuracy solution to the time-domain dynamic response for a symmetrical system with generalized linear viscous damping are established. Based on the formulation constructed, the truncation of the complex mode superposition is investigated. By employing the complex mode expansion of the system flexibility matrix discovered in this work, an accelerating method for the complex mode superposition after truncation is proposed. According to such improved complex mode theory, the system dynamical parameters, for example, the natural frequency and the complex vibration shape, are all the functions with respect to the mass, stiffness, and damping. However, the influence of the damping on the natural frequency may always be relatively slight. The numerical simulating investigation shows that the precise dynamic response can be obtained simply and conveniently in line with the formulae established. Meanwhile, the truncating and accelerating methods are effective for both the narrow or wide band external excitation and can provide the approximate solution to the dynamic response with better precision even when the frequency of the excitation is higher. Under certain conditions, the accelerating approach can be used to boost the accuracy of the truncated response.

Introduction

Generally speaking, if the damping of a symmetrical multi-degree-of-freedom system is generalized linear viscous, the dynamic response of the system commonly needs to be worked out by means of the complex mode method in time domain. The reason essentially is that the damping matrix of such a system cannot be represented as the linear combination of its mass and stiffness matrices.

In the late 1950s, the study on the complex mode theory started from the work of Foss ¹ on the state space approach to solve the dynamic response of a damped linear system. Afterward, Schmitz ² developed the thought of Foss and preliminarily constructed the fundamental formulation for the complex mode theory. Nevertheless, he did not provide any numerical example. In the methods proposed by Foss and Schmitz, the n governing differential equations of motion of second order should be transformed to the $2n$ differential equations of first order. In order to avoid increasing the number of equation, some authors^3,4 straightforwardly established the orthogonal relationship to uncouple the differential equations of motion by virtue of solving a quadratic eigenvalue problem.

After coming into 1980s, Li ⁵ studied the relation between the complex and real modes. Meanwhile, the complex mode method was included in some monographs^6–8 and textbooks^9,10 for postgraduate students on mechanical vibration. In addition, the idea of the complex mode has been embedded into the component mode synthesis successfully. ¹¹ On the other hand, the complex mode theory has been applied in the experimental modal testing widely.^12,13

All these seem to show that the complex mode science has been mature. However, in the opinion of the author, it is not so at all. The deficiency mainly exists at least in three aspects. First, the quadratic eigenvalue problem generally cannot be solved using full-blown commercial computing software easily and directly. Second, the arithmetic of real number inextricably commingles with complex operation and this often results in an outcome with bad precision. Third, unlike the status of the real mode, the truncation and acceleration of the complex mode superposition have not been studied until recently. Unfortunately, in the past decade or so, there is almost no progress in such fields from the perspective of the literature in the mainstream journals.

In this work, the traditional complex mode theory is ameliorated considerably. The complex eigenvalue and eigenvector are proposed to be obtained by means of solving the generalized characteristic value problem in respect of the augmented matrices, which are the partitioned matrices made up to by the mass, stiffness, and damping matrices of the system. Hence, the trouble to solve the quadratic eigenvalue problem in the traditional complex mode theory can be kept away completely. In accordance with the classical complex mode philosophy, the dynamical response is generally attained by virtue of the integration transformation technique after knowing the complex characteristic value and vector. In the methodology suggested, the response is directly attained in time domain by making use of the artifice to separate the real and imaginary parts. Clearly, the new approach is more simple and convenient. On the basis of the orthogonal relation between the complex modes, the system flexibility matrix can be expanded as the superimposing summation of the complex modes. On the strength of this, the truncation and acceleration of the complex mode superposition are discussed in detail. An accelerating technique is put forward for the complex mode superposition after truncation. In addition, the utilization of the trick to separate the real and imaginary parts leads up to the avoidance of the mixture of the real and complex operations. This obviously is beneficial to improve the computing precision of the dynamical response.

Improved complex mode theory

Complex eigenvalue and complex mode of vibration

The governing differential equation of motion for a N-degree-of-freedom linear system can be represented as follows

$$\mathit{M}\stackrel{··}{\overrightarrow{\mathit{x}}}+\mathit{C}\stackrel{·}{\overrightarrow{\mathit{x}}}+\mathit{K}\overrightarrow{\mathit{x}}=\overrightarrow{\mathit{P}}\left(t\right)$$

(1)

whose initial conditions are

$$\overrightarrow{\mathit{x}}\left(0\right)={\overrightarrow{\mathit{x}}}_0,\stackrel{·}{\overrightarrow{\mathit{x}}}\left(0\right)=\stackrel{·}{\overrightarrow{\mathit{x}}}$$

(2)

where M , C , and K are the constant mass, damping, and stiffness matrices, respectively, and all of them are of order $N\times N$ . Concurrently, M , C , and K are all real symmetrical matrices, in which M is a positive definite matrix while K is a positive semi-definite matrix. What is more, $\overrightarrow{\mathit{P}}(t)$ is the excitation column matrix while $\overrightarrow{\mathit{x}}$ is the response column matrix.

Equation (1) may be transformed into the state equation for the system as below

$$\mathit{A}\stackrel{·}{\overrightarrow{\mathit{y}}}+\mathit{B}\overrightarrow{\mathit{y}}=\overrightarrow{\mathit{Q}}\left(t\right)$$

(3)

where $\overrightarrow{\mathit{y}}=\left[\begin{array}{c}\stackrel{·}{\overrightarrow{\mathit{x}}}\hfill \\ \multicolumn{1}{c}{\overrightarrow{\mathit{x}}}\end{array}\right]$ , $\overrightarrow{\mathit{Q}}\left(t\right)=\left[\begin{array}{c}0\hfill \\ \multicolumn{1}{c}{\overrightarrow{\mathit{P}}\left(t\right)}\end{array}\right]$ , $\mathit{A}=\left[\begin{array}{cc}0\hfill & \mathit{M}\hfill \\ \multicolumn{1}{c}{\mathit{M}}& \mathit{C}\hfill \end{array}\right]$ , and $\mathit{B}=\left[\begin{array}{cc}-\mathit{M}\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& \mathit{K}\hfill \end{array}\right]$ . It is easy to verify that ${\mathit{A}}^T=\mathit{A}$ and ${\mathit{B}}^T=\mathit{B}$ . Thereby both A and B are real symmetrical matrices.

Letting $\overrightarrow{\mathit{x}}=\overrightarrow{\mathit{\varphi }}{e}^{\lambda t}$ yields

$$\overrightarrow{\mathit{y}}=\overrightarrow{\mathit{\psi }}{e}^{\lambda t}$$

(4)

where $\overrightarrow{\mathit{\psi }}=\left[\begin{array}{c}\lambda \overrightarrow{\mathit{\varphi }}\hfill \\ \multicolumn{1}{c}{\overrightarrow{\mathit{\varphi }}}\end{array}\right]$ . Obviously, $\overrightarrow{\mathit{\varphi }}$ is a column matrix of order $N\times 1$ while $\overrightarrow{\mathit{\psi }}$ is a column matrix of order $2N\times 1$ .

Letting $\overrightarrow{\mathit{Q}}(t)=0$ and substituting equation (4) into equation (3) leads to

$$\mathit{B}\overrightarrow{\mathit{\psi }}=-\lambda \mathit{A}\overrightarrow{\mathit{\psi }}$$

(5)

Because both A and B are not positive definite commonly, their Cholesky decompositions cannot be implemented. ¹⁴ Therefore, in general, the generalized characteristic value problem (5) cannot be transformed into a characteristic value problem in regard to a real symmetrical matrix like the case of the real mode.

By taking advantage of the determinant computing approach for the partitioned matrix, ¹⁵ the determinant value of the matrix A can be obtained as follows

$$\left|\mathit{A}\right|={(-1)}^N\left|\begin{array}{cc}\mathit{M}\hfill & 0\hfill \\ \multicolumn{1}{c}{\mathit{C}}& \mathit{M}\hfill \end{array}\right|={(-1)}^N\left|\mathit{M}\right||\mathit{M}-\mathit{C}{\mathit{M}}^{-1}0|={(-1)}^N{\left|\mathit{M}\right|}^2\ne 0$$

(6)

Thereby the matrix A is reversible and the generalized proper value problem (5) can be transformed into a proper value problem in regard to a real square matrix as

$$\mathit{D}\overrightarrow{\mathit{\psi }}=\lambda \overrightarrow{\mathit{\psi }}$$

(7)

where $\mathit{D}=-{\mathit{A}}^{-1}\mathit{B}$ and generally the real square matrix D is not a symmetrical matrix.

Apparently the characteristic polynomial of the real square matrix D always has real coefficients. As a consequence, the $2N$ eigenvalues, ${\lambda }_1$ , ${\lambda }_2$ , …, and ${\lambda }_{2N}$ , belonging to the proper value problem (7), are usually pair-wise conjugate in complex number field. ¹⁶ This of course gives rise to that the corresponding $2N$ eigenvectors, ${\overrightarrow{\mathit{\psi }}}_1$ , ${\overrightarrow{\mathit{\psi }}}_2$ , …, and ${\overrightarrow{\mathit{\psi }}}_{2N}$ , are pair-wise conjugate as well. Parenthetically, the above eigenvalue ${\lambda }_j$ may be called the complex eigenvalue of the system while the corresponding eigenvector ${\overrightarrow{\mathit{\psi }}}_j$ may be called the state complex mode vector. The state complex mode vectors can compose the state complex vibration shape matrix, which can be written as

$$\mathit{\Psi }=\left[\begin{array}{ccc}{\overrightarrow{\mathit{\psi }}}_1\hfill & \cdots \hfill & {\overrightarrow{\mathit{\psi }}}_{2N}\hfill \end{array}\right]=\left[\begin{array}{ccc}{\lambda }_1{\overrightarrow{\mathit{\varphi }}}_1\hfill & \cdots \hfill & {\lambda }_{2N}{\overrightarrow{\mathit{\varphi }}}_{2N}\hfill \\ \multicolumn{1}{c}{{\overrightarrow{\mathit{\varphi }}}_1}& \cdots \hfill & {\overrightarrow{\mathit{\varphi }}}_{2N}\hfill \end{array}\right]=\left[\begin{array}{c}\mathit{\Phi }\mathit{\Lambda }\hfill \\ \multicolumn{1}{c}{\mathit{\Phi }}\end{array}\right]$$

(8)

where $\mathit{\Phi }$ is the complex vibration shape matrix and $\mathit{\Phi }=[\begin{array}{ccc}{\overrightarrow{\mathit{\varphi }}}_1\hfill & \cdots \hfill & {\overrightarrow{\mathit{\varphi }}}_{2N}\hfill \end{array}]$ , and $\mathit{\Lambda }=\mathrm{diag}(\begin{array}{c}{\lambda }_1,\cdots ,{\lambda }_{2N}\hfill \end{array})$ . Herein, ${\overrightarrow{\mathit{\varphi }}}_j$ $(j=1,2,\dots ,2N)$ may be called the complex mode vector. Obviously, the matrix $\mathit{\Psi }$ is a complex square matrix of order $2N\times 2N$ and the matrix $\mathit{\Phi }$ is a complex matrix of order $N\times 2N$ .

Provided the system in question does not have repeated complex proper value, it is easy to build the orthogonal relation between two different state complex mode vectors from equation (5), which can be expressed as

$${\overrightarrow{\mathit{\psi }}}_s^T\mathit{A}{\overrightarrow{\mathit{\psi }}}_k=0,{\overrightarrow{\mathit{\psi }}}_s^T\mathit{B}{\overrightarrow{\mathit{\psi }}}_k=0,s,k=1,2,\dots ,2N,s\ne k$$

(9)

Under the condition of no repeated complex proper value, ${\overrightarrow{\mathit{\psi }}}_1$ , ${\overrightarrow{\mathit{\psi }}}_2$ , …, and ${\overrightarrow{\mathit{\psi }}}_{2N}$ are the proper vector belonging to different proper value so that they are linearly independent. Consequently, the state complex vibration shape matrix, $\Psi$ , made up by them is invertible.

With the aid of equation (9), both the matrices A and B can be diagonalized and the results obtained can be expressed as follows

$${\mathit{\Psi }}^T\mathit{A}\mathit{\Psi }=\mathrm{diag}\left(\begin{array}{ccc}{a}_1,\hfill & \dots ,\hfill & a_{2N}\hfill \end{array}\right)={\mathit{A}}_p,{\mathit{\Psi }}^T\mathit{B}\mathit{\Psi }=\mathrm{diag}\left(\begin{array}{ccc}{b}_1,\hfill & \dots ,\hfill & b_{2N}\hfill \end{array}\right)={\mathit{B}}_p$$

(10)

where $a_j={\overrightarrow{\mathit{\psi }}}_j^T\mathit{A}{\overrightarrow{\mathit{\psi }}}_j=2{\lambda }_j{\overrightarrow{\mathit{\varphi }}}_j^T\mathit{M}{\overrightarrow{\mathit{\varphi }}}_j+{\overrightarrow{\mathit{\varphi }}}_j^T\mathit{C}{\overrightarrow{\mathit{\varphi }}}_j$ , $b_j={\overrightarrow{\mathit{\psi }}}_j^T\mathit{B}{\overrightarrow{\mathit{\psi }}}_j=-{\lambda }_j{\overrightarrow{\mathit{\varphi }}}_j^T\mathit{M}{\overrightarrow{\mathit{\varphi }}}_j+{\overrightarrow{\mathit{\varphi }}}_j^T\mathit{K}{\overrightarrow{\mathit{\varphi }}}_j$ , $j=1,2,\dots ,2N$ .

By right of equations (5) and (9), the relationship among $a_j$ , $b_j$ , and ${\lambda }_j$ can be attained as

$${\lambda }_j=-b_j/a_j={\mu }_j+i{\nu }_j,j=1,2,\dots ,2N$$

(11)

In order to ensure the stability of the system vibration, the real part of ${\lambda }_j$ should be negative, that is, ${\mu }_j<0$ .

Owing to the conjugacy of the $2N$ complex characteristic values, the natural circular frequency of the system can be defined as

$${\omega }_r=\left|{\lambda }_r\right|=\sqrt{{\mu }_r^2+{\nu }_r^2},r=1,2,\dots ,N$$

(12)

where ${\lambda }_1$ , ${\lambda }_2$ , …, and ${\lambda }_N$ denote the complex characteristic values extracted from the $2N$ complex characteristic values and they are not conjugate with each other.

Furthermore, corresponding to the natural frequency ${\omega }_r$ , the modal damping ratio in the complex mode theory can be defined as follows

$${\xi }_r=\frac{-{\mu }_r}{\left|{\lambda }_r\right|}=\frac{-{\mu }_r}{{\omega }_r}=\frac{-{\mu }_r}{\sqrt{{\mu }_r^2+{\nu }_r^2}}\in (0,1),r=1,2,\dots ,N$$

(13)

Equation (13) shows that, under the condition of the generalized viscous damping, the oscillation of the system can always pertain to the under-damped case in general.

On the basis of ${\xi }_r$ and ${\omega }_r$ , the complex characteristic value can as well be represented as

$${\lambda }_r=-{\xi }_r{\omega }_r+i{\nu }_r,r=1,2,\dots ,N$$

(14)

where ${\nu }_r={\omega }_r\sqrt{1-{\xi }_r^2}$ .

The complex mode vector ${\overrightarrow{\mathit{\varphi }}}_j$ $(j=1,2,\dots ,2N)$ can be written in the form of real and imaginary parts as

$${\overrightarrow{\mathit{\varphi }}}_j={\overrightarrow{\mathit{\varphi }}}_j^{\left(R\right)}+i{\overrightarrow{\mathit{\varphi }}}_j^{\left(I\right)}j=1,2,\dots ,2N$$

(15)

where $i=\sqrt{-1}$ . By employing equation (15), it is possible to have

$${\overrightarrow{\mathit{\varphi }}}_j^T\mathit{M}{\overrightarrow{\mathit{\varphi }}}_j=m_j^{\left(\mathrm{RR}\right)}-m_j^{\left(\mathrm{II}\right)}+2{\mathrm{im}}_j^{\left(\mathrm{RI}\right)},{\overrightarrow{\mathit{\varphi }}}_j^T\mathit{C}{\overrightarrow{\mathit{\varphi }}}_j=c_j^{\left(\mathrm{RR}\right)}-c_j^{\left(\mathrm{II}\right)}+2{\mathrm{ic}}_j^{\left(\mathrm{RI}\right)},j=1,2,\dots ,2N$$

(16)

where $m_j^{(\mathrm{RR})}=({\overrightarrow{\mathit{\varphi }}}_j^{(R)}{)}^T\mathit{M}{\overrightarrow{\mathit{\varphi }}}_j^{(R)}$ , $m_j^{(\mathrm{II})}=({\overrightarrow{\mathit{\varphi }}}_j^{(I)}{)}^T\mathit{M}{\overrightarrow{\mathit{\varphi }}}_j^{(I)}$ , $m_j^{(\mathrm{RI})}=m_j^{(\mathrm{IR})}=({\overrightarrow{\mathit{\varphi }}}_j^{(R)}{)}^T\mathit{M}{\overrightarrow{\mathit{\varphi }}}_j^{(I)}$ ; $c_j^{(\mathrm{RR})}=({\overrightarrow{\mathit{\varphi }}}_j^{(R)}{)}^T\mathit{C}{\overrightarrow{\mathit{\varphi }}}_j^{(R)}$ , $c_j^{(\mathrm{II})}=({\overrightarrow{\mathit{\varphi }}}_j^{(I)}{)}^T\mathit{C}{\overrightarrow{\mathit{\varphi }}}_j^{(I)}$ , and $c_j^{(\mathrm{RI})}=c_j^{(\mathrm{IR})}=({\overrightarrow{\mathit{\varphi }}}_j^{(R)}{)}^T\mathit{C}{\overrightarrow{\mathit{\varphi }}}_j^{(I)}$ .

After substituting equation (16) into the expression of $a_j$ in equation (10), it is easy to acquire the expression for $a_j$ in the form of real and imaginary parts as follows

$$a_j=a_j^{\left(R\right)}+{\mathrm{ia}}_j^{\left(I\right)},j=1,2,\dots ,2N$$

(17)

where $a_j^{(R)}=c_j^{(\mathrm{RR})}-c_j^{(\mathrm{II})}-2{\xi }_j{\omega }_j(m_j^{(\mathrm{RR})}-m_j^{(\mathrm{II})})-4{\nu }_j{m}_j^{(\mathrm{RI})}$ and $a_j^{(I)}=2c_j^{(\mathrm{RI})}-4{\xi }_j{\omega }_j{m}_j^{(\mathrm{RI})}+2{\nu }_j(m_j^{(\mathrm{RR})}-m_j^{(\mathrm{II})})$ .

On the other hand, from equation (15), it is easy to procure

$${\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T={\mathit{\Phi }}_j^{\left(1\right)}+i{\mathit{\Phi }}_j^{\left(2\right)},j=1,2,\dots ,2N$$

(18)

where ${\mathit{\Phi }}_j^{(1)}={\overrightarrow{\mathit{\varphi }}}_j^{(R)}({\overrightarrow{\mathit{\varphi }}}_j^{(R)}{)}^T-{\overrightarrow{\mathit{\varphi }}}_j^{(I)}({\overrightarrow{\mathit{\varphi }}}_j^{(I)}{)}^T$ and ${\mathit{\Phi }}_j^{(2)}={\overrightarrow{\mathit{\varphi }}}_j^{(R)}({\overrightarrow{\mathit{\varphi }}}_j^{(I)}{)}^T+{\overrightarrow{\mathit{\varphi }}}_j^{(I)}({\overrightarrow{\mathit{\varphi }}}_j^{(R)}{)}^T$ .

It is not difficult to understand that the intrinsic dynamical characteristics of the system, including the natural frequency and the complex mode, are all the functions regarding the mass, stiffness, and damping in accordance with the improved complex mode theory established here. Nonetheless, the situation in the real mode theory is far from this. There the intrinsic dynamical characteristics solely are functions regarding the mass and stiffness because the damping matrix always is dealt with to be the linear combination of the mass and stiffness matrices in terms of the Rayleigh assumption. From such a perspective, the improved complex mode theory proposed in this study is more impeccable and rigorous.

Treatment of initial condition

With the help of equation (10), letting $\overrightarrow{\mathit{y}}=\mathit{\Psi }\overrightarrow{\mathit{z}}$ and submitting it into equation (3) gives rise to

$${\stackrel{.}{z}}_j-{\lambda }_j{z}_j=\frac{1}{a_j}{\overrightarrow{\mathit{\varphi }}}_j^T\overrightarrow{\mathit{P}}\left(t\right),j=1,2,\dots ,2N$$

(19)

As already stated, the state complex vibration shape matrix, $\mathit{\Psi }$ , is reversible and its inverse matrix can be attained from equation (10) as ${\mathit{\Psi }}^{-1}={\mathit{A}}_p^{-1}{\mathit{\Psi }}^T\mathit{A}$ . Consequently, it is likely to have

$$\overrightarrow{\mathit{z}}={\mathit{\Psi }}^{-1}\overrightarrow{\mathit{y}}={\mathit{A}}_p^{-1}{\mathit{\Psi }}^T\mathit{A}\left[\begin{array}{c}\stackrel{·}{\overrightarrow{\mathit{x}}}\hfill \\ \multicolumn{1}{c}{\overrightarrow{\mathit{x}}}\end{array}\right]$$

(20)

In view of equations (2) and (20), the initial condition in respect of equation (19) can be established as below

$$z_j\left(0\right)=\frac{1}{{\left|a_j\right|}^2}{\overrightarrow{\mathit{\varphi }}}_j^T({\overrightarrow{\mathit{\eta }}}_j^{\left(R\right)}+i{\overrightarrow{\mathit{\eta }}}_j^{\left(I\right)}),j=1,2,\dots ,2N$$

(21)

where ${\overrightarrow{\mathit{\eta }}}_j^{(R)}=[({\nu }_j{a}_j^{(I)}-{\xi }_j{\omega }_j{a}_j^{(R)})\mathit{M}+a_j^{(R)}\mathit{C}]{\overrightarrow{\mathit{x}}}_0+a_j^{(R)}\mathit{M}{\stackrel{·}{\overrightarrow{\mathit{x}}}}_0$ and ${\overrightarrow{\mathit{\eta }}}_j^{(I)}=[({\xi }_j{\omega }_j{a}_j^{(I)}+{\nu }_j{a}_j^{(R)})\mathit{M}-a_j^{(I)}\mathit{C}]{\overrightarrow{\mathit{x}}}_0-a_j^{(I)}\mathit{M}{\stackrel{·}{\overrightarrow{\mathit{x}}}}_0$ .

Accurate solution to dynamic response

Equation (19) is a linear non-homogeneous ordinary differential equation of first order with regard to $z_j$ , whose solution is

$$z_j=z_j\left(0\right)e^{{\lambda }_{j}t}+\frac{1}{a_j}{\overrightarrow{\mathit{\varphi }}}_j^T{\displaystyle \underset{0}{\overset{t}{\int }}\overrightarrow{\mathit{P}}\left(\tau \right)e^{{\lambda }_j\left(t-\tau \right)}}d\tau ,j=1,2,\dots ,2N$$

(22)

After the coordinate transformation, the accurate solution of equation (1) can be acquired from equation (22) and the result is

$$\overrightarrow{\mathit{x}}=\sum _{j=1}^{2N}{z}_j{\overrightarrow{\mathit{\varphi }}}_j=\sum _{j=1}^{2N}{z}_j\left(0\right)e^{{\lambda }_{j}t}{\overrightarrow{\mathit{\varphi }}}_j+\sum _{j=1}^{2N}\frac{1}{a_j}\left({\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T\right)\underset{0}{\overset{t}{\int }}\overrightarrow{\mathit{P}}\left(\tau \right)e^{{\lambda }_j(t-\tau )}d\tau$$

(23)

From equation (23), making use of equation (18) and taking into account the conjugacy, the total response of the system in time domain can be represented as

$$\overrightarrow{\mathit{x}}\left(t\right)={\overrightarrow{\mathit{x}}}_{\mathrm{in}}\left(t\right)+2\sum _{r=1}^N\frac{1}{{\left|a_r\right|}^2}\underset{0}{\overset{t}{\int }}{e}^{-{\xi }_r{\omega }_r(t-\tau )}[{\overrightarrow{\mathit{Q}}}_r^{\left(A\right)}\left(\tau \right)\cos {\nu }_r(t-\tau )-{\overrightarrow{\mathit{Q}}}_r^{\left(B\right)}\left(\tau \right)\sin {\nu }_r(t-\tau )]d\tau$$

(24)

where ${\overrightarrow{\mathit{Q}}}_r^{(A)}(\tau )=(a_r^{(R)}{\mathit{\Phi }}_r^{(1)}+a_r^{(I)}{\mathit{\Phi }}_r^{(2)})\overrightarrow{\mathit{P}}(\tau )$ , ${\overrightarrow{\mathit{Q}}}_r^{(B)}(\tau )=(a_r^{(R)}{\mathit{\Phi }}_r^{(2)}-a_r^{(I)}{\mathit{\Phi }}_r^{(1)})\overrightarrow{\mathit{P}}(\tau )$ ; and ${\overrightarrow{\mathit{x}}}_{\mathrm{in}}\left(t\right)$ denotes the response gotten rise to by the initial condition. By virtue of equations (21) and (23), ${\overrightarrow{\mathit{x}}}_{\mathrm{in}}\left(t\right)$ can be achieved as

$${\overrightarrow{\mathit{x}}}_{\mathrm{in}}\left(t\right)=2\sum _{r=1}^N\frac{e^{-{\xi }_r{\omega }_{r}t}}{{\left|a_r\right|}^2}[{\overrightarrow{\mathit{\rho }}}_r^{\left(A\right)}\cos {\nu }_{r}t-{\overrightarrow{\mathit{\rho }}}_r^{\left(B\right)}\sin {\nu }_{r}t]$$

(25)

where ${\overrightarrow{\mathit{\rho }}}_r^{(A)}={\mathit{\Phi }}_r^{(1)}{\overrightarrow{\mathit{\eta }}}_r^{(R)}-{\mathit{\Phi }}_r^{(2)}{\overrightarrow{\mathit{\eta }}}_r^{(I)}$ and ${\overrightarrow{\mathit{\rho }}}_r^{(B)}={\mathit{\Phi }}_r^{(1)}{\overrightarrow{\mathit{\eta }}}_r^{(I)}+{\mathit{\Phi }}_r^{(2)}{\overrightarrow{\mathit{\eta }}}_r^{(R)}$ .

As shown by equation (24), the total dynamic response of the system consists of two parts. One is, ${\overrightarrow{\mathit{x}}}_{\mathrm{in}}\left(t\right)$ , which is the response caused by the initial condition and usually corresponds to the general solution. The other is, ${\overrightarrow{\mathit{x}}}_e\left(t\right)$ , the response caused by the external excitation and is generally analogous to the particular solution. Distinctly, the response caused by the initial conditions will decay in time domain following an exponential rule as indicated by equation (25). Therefore, the steady-state dynamic response is mainly determined by the response led to by the external excitation. This point is the same with the status in the real mode.

After transformation, the integral in equation (24) may have the roughly same form with the Duhamel integral, that is, the convolution integral, ¹⁷ in the real mode theory. Consequently, these integrals in equation (24) can be called the generalized Duhamel integral, which can be expressed as

$${\overrightarrow{\mathit{I}}}_r\left(t\right)=\underset{0}{\overset{t}{\int }}{\overrightarrow{\mathit{Q}}}_r^{\left(A\right)}\left(\tau \right)e^{-{\xi }_r{\omega }_r(t-\tau )}\cos {\nu }_r(t-\tau )d\tau ,{\overrightarrow{\mathit{J}}}_r\left(t\right)=\underset{0}{\overset{t}{\int }}{\overrightarrow{\mathit{Q}}}_r^{\left(B\right)}\left(\tau \right)e^{-{\xi }_r{\omega }_r(t-\tau )}\sin {\nu }_r(t-\tau )d\tau$$

(26)

In such integrals, ${\overrightarrow{\mathit{Q}}}_r^{(A)}(\tau )$ and ${\overrightarrow{\mathit{Q}}}_r^{(B)}(\tau )$ are associated with the excitation and they can be called the complex mode excitation.

Resorting to the preceding concept of the generalized Duhamel integral, the representation of the total dynamic response can further be reduced to

$$\overrightarrow{\mathit{x}}\left(t\right)={\overrightarrow{\mathit{x}}}_{\mathrm{in}}\left(t\right)+2\sum _{r=1}^N\frac{1}{{\left|a_r\right|}^2}[{\overrightarrow{\mathit{I}}}_r\left(t\right)-{\overrightarrow{\mathit{J}}}_r\left(t\right)]$$

(27)

Generally speaking, it is simpler to obtain the time-domain response for a symmetrical system directly from equation (27) because it need not be with the aid of the integral transformation. ¹⁸ Concurrently, the number of the equation essentially is not increased, nor is the quadratic eigenvalue problem required to be solved.

Truncation of complex mode superposition and related acceleration method

Truncation of complex modes superposition

At heart, the approach to obtain the response from equation (27) should pertain to the superposition of the complex modes in time domain. By employing such an approach, the accurate response can be worked out. However, if the number of the degree of freedom is relatively large, the preceding approach to search for accurate solution may be difficult or inconvenient to be implemented at least to some degree. To handle this challenge, one way is to truncate the summation of the complex modes in equation (27) and the outcome truncated can be represented as

$${\overrightarrow{\mathit{x}}}^{*}\left(t\right)={\overrightarrow{\mathit{x}}}_{\mathrm{in}}^{*}\left(t\right)+2\sum _{r=1}^S\frac{1}{{\left|a_r\right|}^2}[{\overrightarrow{\mathit{I}}}_r\left(t\right)-{\overrightarrow{\mathit{J}}}_r\left(t\right)],0<S<N$$

(28)

where ${\overrightarrow{\mathit{x}}}_{\mathrm{in}}^{*}\left(t\right)=2\sum _{r=1}^S\frac{e^{-{\xi }_r{\omega }_{r}t}}{{\left|a_r\right|}^2}[{\overrightarrow{\mathit{\rho }}}_r^{\left(A\right)}\cos {\nu }_{r}t-{\overrightarrow{\mathit{\rho }}}_r^{\left(B\right)}\sin {\nu }_{r}t]$ . Herein, S symbolizes the number of the terms kept in the superposition summation after the truncation.

In fact, equation (28) can be viewed as the truncation of the complex mode superposition, which can merely present an approximate solution to the dynamic response. However, the relevant computing workload is sure to be notably slighter. Given this, the investigation on how to truncate the complex mode superposition rationally is of great significance.

The truncation of the real mode superposition, covering the corresponding accelerating method, has been discussed in detail in the literature.^9,19 Nevertheless, the truncation of the complex mode superposition and the related accelerating method have never been reported until nowadays as far as the author has known.

Acceleration of complex mode superposition after truncation

So as to enhance the accuracy of the truncated approximate solution, the accelerating method of the complex mode superposition after truncation needs to be discovered. For this, the complex mode expansion of the flexibility matrix of the system ought to be established in advance.

From equation (10), the inverse matrix of the matrix B can be acquired and the result is

$$\begin{array}{c}{\mathit{B}}^{-1}=\left[\begin{array}{cc}-{\mathit{M}}^{-1}& 0\\ 0& {\mathit{K}}^{-1}\end{array}\right]=\mathit{\Psi }{\mathit{B}}_p^{-1}{\mathit{\Psi }}^T=\left[\begin{array}{ccc}{\lambda }_1{\overrightarrow{\mathit{\varphi }}}_1& \cdots & {\lambda }_{2N}{\overrightarrow{\mathit{\varphi }}}_{2N}\\ {\overrightarrow{\mathit{\varphi }}}_1& \cdots & {\overrightarrow{\mathit{\varphi }}}_{2N}\end{array}\right]\left[\begin{array}{ccc}{b}_1^{-1}& & 0\\ & \ddots & \\ 0& & b_{2N}^{-1}\end{array}\right]{\left[\begin{array}{ccc}{\lambda }_1{\overrightarrow{\mathit{\varphi }}}_1& \cdots & {\lambda }_{2N}{\overrightarrow{\mathit{\varphi }}}_{2N}\\ {\overrightarrow{\mathit{\varphi }}}_1& \cdots & {\overrightarrow{\mathit{\varphi }}}_{2N}\end{array}\right]}^T\\ =\left[\begin{array}{cc}{\displaystyle \sum _{j=1}^{2N}\frac{{\lambda }_j^2}{b_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T}& {\displaystyle \sum _{j=1}^{2N}\frac{{\lambda }_j}{b_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T}\\ {\displaystyle \sum _{j=1}^{2N}\frac{{\lambda }_j}{b_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T}& {\displaystyle \sum _{j=1}^{2N}\frac{1}{b_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T}\end{array}\right]\end{array}$$

(29)

By comparing with the counterparts of the partitioned matrices on both sides of equation (29), the following equivalence relations can be constructed as

$${\mathit{M}}^{-1}=\sum _{j=1}^{2N}\frac{{\lambda }_j}{a_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T,{\mathit{K}}^{-1}=\mathit{F}=-\sum _{j=1}^{2N}\frac{1}{a_j{\lambda }_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T,\sum _{j=1}^{2N}\frac{1}{a_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T=0$$

(30)

where F is the flexibility matrix of the system, which is the inverse matrix of the stiffness matrix K . The second equality in equation (30) is the so-called complex mode expansion of the flexibility matrix of the system.

It is assumed that the initial condition is zero and the external excitation is simple harmonic, that is, $\overrightarrow{\mathit{P}}(\tau )={\overrightarrow{\mathit{P}}}_0\sin \omega \tau$ , in which ${\overrightarrow{\mathit{P}}}_0$ is a constant column vector. Under these conditions, the steady-state dynamic response of the system can roughly be estimated from equation (22), whose consequence is

$$\begin{array}{c}\multicolumn{1}{c}{z_j=\frac{e^{{\lambda }_{j}t}}{a_j}{\overrightarrow{\mathit{\varphi }}}_j^T{\overrightarrow{\mathit{P}}}_0\underset{0}{\overset{t}{\int }}\sin \omega \tau e^{-{\lambda }_j\tau }d\tau =\frac{{\overrightarrow{\mathit{\varphi }}}_j^T{\overrightarrow{\mathit{P}}}_0}{a_j({\omega }^2+{\lambda }_j^2)}}\\ \multicolumn{1}{c}{(\omega e^{{\lambda }_{j}t}-\omega \cos \omega t-{\lambda }_j\sin \omega t)}\\ \multicolumn{1}{c}{\approx -\frac{{\overrightarrow{\mathit{\varphi }}}_j^T{\overrightarrow{\mathit{P}}}_0}{a_j({\omega }^2+{\lambda }_j^2)}({\lambda }_j\sin \omega t+\omega \cos \omega t)}\\ \multicolumn{1}{c}{\approx -\frac{{\overrightarrow{\mathit{\varphi }}}_j^T{\overrightarrow{\mathit{P}}}_0}{a_j{\lambda }_j^2}{\lambda }_j\sin \omega t=-\frac{1}{a_j{\lambda }_j}{\overrightarrow{\mathit{\varphi }}}_j^T\overrightarrow{\mathit{P}}\left(t\right),j=1,2,\dots ,2N}\end{array}$$

(31)

On the basis of equation (8), the relationship between the state vector, $\overrightarrow{\mathit{y}}$ , and the vector, $\overrightarrow{\mathit{z}}$ , can be built as

$$\overrightarrow{\mathit{y}}=\left[\begin{array}{c}\stackrel{·}{\overrightarrow{\mathit{x}}}\hfill \\ \multicolumn{1}{c}{\overrightarrow{\mathit{x}}}\end{array}\right]=\mathit{\Psi }\overrightarrow{\mathit{z}}=\left[\begin{array}{c}\mathit{\Phi }\mathit{\Lambda }\hfill \\ \multicolumn{1}{c}{\mathit{\Phi }}\end{array}\right]\overrightarrow{\mathit{z}}$$

(32)

such that the relationship between $\overrightarrow{\mathit{z}}$ and $\overrightarrow{\mathit{x}}$ is $\overrightarrow{\mathit{x}}=\mathit{\Phi }\overrightarrow{\mathit{z}}$ .

If both the complex vibration shape matrix $\mathit{\Phi }$ and the vector $\overrightarrow{\mathit{z}}$ are split into the higher and lower order parts, respectively, the above relation can be partitioned as

$$\overrightarrow{\mathit{x}}=\mathit{\Phi }\overrightarrow{\mathit{z}}=\left[\begin{array}{cc}{\mathbf{\Phi }}_{\mathbf{L}}\hfill & {\mathit{\Phi }}_H\hfill \end{array}\right]\left[\begin{array}{c}{\overrightarrow{\mathit{z}}}_L\hfill \\ \multicolumn{1}{c}{{\overrightarrow{\mathit{z}}}_H}\end{array}\right]={\mathit{\Phi }}_L{\overrightarrow{\mathit{z}}}_L+{\mathit{\Phi }}_H{\overrightarrow{\mathit{z}}}_H=\sum _{j=1}^{2S}{z}_j{\overrightarrow{\mathit{\varphi }}}_j+\sum _{j=2S+1}^{2N}{z}_j{\overrightarrow{\mathit{\varphi }}}_j$$

(33)

Substituting equation (31) into equation (33) yields an approximate solution to the response under the zero initial condition. By right of the complex mode expansion of the flexibility matrix in equation (30), the result acquired can be simplified to

$$\begin{array}{c}\multicolumn{1}{c}{\overrightarrow{\mathit{x}}\left(t\right)\approx \sum _{j=1}^{2S}{z}_j{\overrightarrow{\mathit{\varphi }}}_j-\left(\sum _{j=2S+1}^{2N}\frac{1}{a_j{\lambda }_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T\right)\overrightarrow{\mathit{P}}\left(t\right)}\\ \multicolumn{1}{c}{=\sum _{j=1}^{2S}{z}_j{\overrightarrow{\mathit{\varphi }}}_j-(\sum _{j=1}^{2N}\frac{1}{a_j{\lambda }_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T-\sum _{j=1}^{2S}\frac{1}{a_j{\lambda }_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T)\overrightarrow{\mathit{P}}\left(t\right)}\\ \multicolumn{1}{c}{=\sum _{j=1}^{2S}{z}_j{\overrightarrow{\mathit{\varphi }}}_j+(F+\sum _{j=1}^{2S}\frac{1}{a_j{\lambda }_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T)\overrightarrow{\mathit{P}}\left(t\right)}\end{array}$$

(34)

From equations (17) and (18), it is possible to have

$$\frac{1}{a_j{\lambda }_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T=\frac{-1}{{\omega }_j^2{\left|a_j\right|}^2}[h_j^{\left(1\right)}{\mathit{\Phi }}_j^{\left(1\right)}-h_j^{\left(2\right)}{\mathit{\Phi }}_j^{\left(2\right)}+i(h_j^{\left(1\right)}{\mathit{\Phi }}_j^{\left(2\right)}+h_j^{\left(2\right)}{\mathit{\Phi }}_j^{\left(1\right)})],j=1,2,\dots ,2N$$

(35)

where $h_j^{(1)}={\xi }_j{\omega }_j{a}_j^{(R)}+{\nu }_j{a}_j^{(I)}$ and $h_j^{(2)}={\nu }_j{a}_j^{(R)}-{\xi }_j{\omega }_j{a}_j^{(I)}$ .

In terms of the conjugacy, from equation (35), it is easy to obtain

$$\sum _{j=1}^{2S}\frac{1}{a_j{\lambda }_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T=2\mathrm{Re}\left[\sum _{j=1}^{2S}\frac{1}{a_j{\lambda }_j}{\overrightarrow{\mathit{\varphi }}}_j{\overrightarrow{\mathit{\varphi }}}_j^T\right]=-2\sum _{r=1}^S\frac{1}{{\omega }_r^2{\left|a_r\right|}^2}(h_r^{\left(1\right)}{\mathit{\Phi }}_r^{\left(1\right)}-h_r^{\left(2\right)}{\mathit{\Phi }}_r^{\left(2\right)})$$

(36)

Substituting equation (36) into equation (34) and making utilization of equation (28), an accelerated approximate solution to the response after the truncation during the complex mode superposition can be attained as

$${\overrightarrow{\mathit{x}}}_a\left(t\right)={\overrightarrow{\mathit{x}}}_{\mathrm{in}}^{*}\left(t\right)+2\sum _{r=1}^S\frac{1}{{\left|a_r\right|}^2}[{\overrightarrow{\mathit{I}}}_r\left(t\right)-{\overrightarrow{\mathit{J}}}_r\left(t\right)]+\Delta \mathit{F}\overrightarrow{\mathit{P}}\left(t\right)$$

(37)

where $\Delta F$ may be called the surplus flexibility matrix and

$$\Delta \mathit{F}=\mathit{F}-2\sum _{r=1}^S\frac{1}{{\omega }_r^2{\left|a_r\right|}^2}(h_r^{\left(1\right)}{\mathit{\Phi }}_r^{\left(1\right)}-h_r^{\left(2\right)}{\mathit{\Phi }}_r^{\left(2\right)})$$

(38)

Equation (37) is the computing formula for the approximate dynamic response by adopting the accelerating method of the complex mode superposition after truncation. In comparison to equation (28), $\Delta \mathit{F}\overrightarrow{\mathit{P}}\left(t\right)$ in equation (37) is the accelerating term, which is capable of compensating the truncation error at least in some measure. Hence, ${\overrightarrow{\mathit{x}}}_a\left(t\right)$ usually comes to be a better approximate solution to $\overrightarrow{\mathit{x}}\left(t\right)$ than ${\overrightarrow{\mathit{x}}}^{*}\left(t\right)$ , which is attained only by truncating the superimposed summation of the complex modes.

Numerical example investigation

The main goal of this section is to research the characteristics of the truncation and acceleration of the complex mode superposition via the numerical examples. To this end, a 5-degree-of-freedom system, as shown in Figure 1, is taken into consideration. The mass, stiffness, and damping matrices of the system are

$$\mathit{M}=\mathrm{diag}(3m,2m,m,2m,3m)$$

(39)

$$\mathit{K}=\left[\begin{array}{ccccc}3k\hfill & -2k\hfill & 0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{-2k}& 5k\hfill & -3k\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& -3k\hfill & 7k\hfill & -4k\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & -4k\hfill & 9k\hfill & -5k\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill & -5k\hfill & 11k\hfill \end{array}\right],\mathit{C}=\left[\begin{array}{ccccc}11c\hfill & -5c\hfill & 0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{-5c}& 9c\hfill & -4c\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& -4c\hfill & 7c\hfill & -3c\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & -3c\hfill & 5c\hfill & -2c\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill & -2c\hfill & 3c\hfill \end{array}\right]$$

(40)

in which $m=1\mathrm{kg}$ , $k=1000\mathrm{N}/\mathrm{m}$ , and $c=5\mathrm{kg}/\mathrm{s}$ .

Figure 1.

Physical structure of the system in question.

After solving eigenvalue problem equation (7), the natural frequencies of the system can be obtained, which are listed in Table 1. For the sake of contrast, the natural frequencies without damping are listed in the same table as well, which are attained by means of the real mode method. Table 1 shows that, strictly speaking, the damping effect makes the lower order natural frequencies rising while makes the higher order natural frequencies decreasing. However, the influence of the damping on the natural frequency is rather limited.

Table 1.

Natural frequencies of the system.

	${\omega }_1$	${\omega }_2$	${\omega }_3$	${\omega }_4$	${\omega }_5$
Un-damped natural frequency (rad/s)	18.4269	34.4389	49.8757	71.7138	97.5224
Natural frequency (rad/s)	18.6650	34.9984	49.1199	71.6522	96.2798

Because of the pair-wise conjugacy among the 10 complex modes, only 5 complex modes among them are given in Table 2. The subscripts reflect the corresponding relation between the complex mode and the related inherent frequency. In addition, the complex modal damping ratios are tabulated in Table 2.

Table 2.

Complex modes and complex modal damping ratios of the system.

${\overrightarrow{\mathit{\varphi }}}_1$	${\left[\begin{array}{ccccc}-0.0044-0.0309i,\hfill & -0.0001-0.0312i,\hfill & 0.0016-0.0239i,\hfill & 0.0020-0.0163i,\hfill & 0.0013-0.0081i\hfill \end{array}\right]}^T$
${\xi }_1$	0.1417
${\overrightarrow{\mathit{\varphi }}}_2$	${\left[\begin{array}{ccccc}-0.0018-0.0143i,\hfill & 0.0078+0.0017i,\hfill & 0.0058+0.0121i,\hfill & 0.0024+0.0157i,\hfill & 0.0000+0.0106i\hfill \end{array}\right]}^T$
${\xi }_2$	0.1487
${\overrightarrow{\mathit{\varphi }}}_3$	${\left[\begin{array}{ccccc}-0.0026+0.0070i,\hfill & -0.0029-0.0144i,\hfill & 0.0040-0.0055i\hfill & 0.0050+0.0047i,\hfill & 0.0018+0.0069i\hfill \end{array}\right]}^T$
${\xi }_3$	0.2003
${\overrightarrow{\mathit{\varphi }}}_4$	${\left[\begin{array}{ccccc}-0.0008+0.0004i,\hfill & 0.0014-0.0041i,\hfill & 0.0030+0.0071i,\hfill & -0.0018+0.0072i,\hfill & -0.0008-0.0077i\hfill \end{array}\right]}^T$
${\xi }_4$	0.1088
${\overrightarrow{\mathit{\varphi }}}_5$	${\left[\begin{array}{ccccc}0.0002-0.0002i,\hfill & -0.0003+0.0022i,\hfill & -0.0020-0.0090i,\hfill & 0.0021+0.0034i,\hfill & -0.0009-0.0006i\hfill \end{array}\right]}^T$
${\xi }_5$	0.2200

Hereinafter two cases of the external excitation are studied: the simple harmonic excitation and the slope step excitation.

Simple harmonic external excitation

In this case, the external excitation is

$$\overrightarrow{\mathit{P}}\left(t\right)={\overrightarrow{\mathit{P}}}_0\cos \omega t$$

(41)

where ${\overrightarrow{\mathit{P}}}_0=[\begin{array}{ccccc}10,000,\hfill & 0,\hfill & 0,\hfill & 0,\hfill & 0\hfill \end{array}{]}^T$ (N) without loss of generality.

Substituting equation (41) into equation (27), the response of the system caused by the simple harmonic external excitation can be worked out as

$${\overrightarrow{\mathit{x}}}_e\left(t\right)=2\sum _{r=1}^5\frac{1}{{\left|a_r\right|}^2}[f_A^{\left(r\right)}\left(t\right){\overrightarrow{\mathit{Q}}}_r^{\left(A\right)}-f_B^{\left(r\right)}\left(t\right){\overrightarrow{\mathit{Q}}}_r^{\left(B\right)}]$$

(42)

where in accordance with the formulae in equation (24), ${\overrightarrow{\mathit{Q}}}_r^{(A)}$ and ${\overrightarrow{\mathit{Q}}}_r^{(B)}$ can be figured out as

$$\left[\begin{array}{ccccc}{\overrightarrow{\mathit{Q}}}_1^{\left(A\right)},\hfill & {\overrightarrow{\mathit{Q}}}_2^{\left(A\right)},\hfill & {\overrightarrow{\mathit{Q}}}_3^{\left(A\right)},\hfill & {\overrightarrow{\mathit{Q}}}_4^{\left(A\right)},\hfill & {\overrightarrow{\mathit{Q}}}_5^{\left(A\right)}\hfill \end{array}\right]=\left[\begin{array}{ccccc}-0.4147\hfill & 0.0176\hfill & 0.0317\hfill & 0.0003\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{-0.1229}& 0.1061\hfill & -0.0322\hfill & -0.0019\hfill & -0.0001\hfill \\ \multicolumn{1}{c}{0.0229}& 0.0466\hfill & -0.0329\hfill & 0.0023\hfill & 0.0005\hfill \\ \multicolumn{1}{c}{0.0790}& -0.0134\hfill & -0.0060\hfill & 0.0032\hfill & -0.0002\hfill \\ \multicolumn{1}{c}{0.0611}& -0.0311\hfill & 0.0139\hfill & -0.0030\hfill & 0.0000\hfill \end{array}\right]$$

(43)

$$\left[\begin{array}{ccccc}{\overrightarrow{\mathit{Q}}}_1^{\left(B\right)},\hfill & {\overrightarrow{\mathit{Q}}}_2^{\left(B\right)},\hfill & {\overrightarrow{\mathit{Q}}}_3^{\left(B\right)},\hfill & {\overrightarrow{\mathit{Q}}}_4^{\left(B\right)},\hfill & {\overrightarrow{\mathit{Q}}}_5^{\left(B\right)}\hfill \end{array}\right]=\left[\begin{array}{ccccc}-2.0948\hfill & -0.2081\hfill & -0.0201\hfill & 0.0003\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{-2.1311}& 0.0475\hfill & 0.0668\hfill & 0.0001\hfill & -0.0000\hfill \\ \multicolumn{1}{c}{-1.6357}& 0.1884\hfill & 0.0100\hfill & -0.0024\hfill & -0.0002\hfill \\ \multicolumn{1}{c}{-1.1198}& 0.2294\hfill & -0.0342\hfill & -0.0005\hfill & 0.0001\hfill \\ \multicolumn{1}{c}{-0.5578}& 0.1502\hfill & -0.0337\hfill & 0.0017\hfill & -0.0001\hfill \end{array}\right]$$

(44)

On the other hand, the generalized Duhamel integrals can thus be reduced and the integral outcomes acquired in equation (42) are denoted by $f_A^{\left(r\right)}\left(t\right)$ and $f_B^{\left(r\right)}\left(t\right)$ . By employing equations (62) and (63) in Appendix 2, it is easy to have

$$\begin{array}{c}{f}_A^{\left(r\right)}\left(t\right)={\displaystyle \underset{0}{\overset{t}{\int }}{e}^{-{\xi }_r{\omega }_r\left(t-\tau \right)}\cos {\nu }_r\left(t-\tau \right)\cos \omega \tau d\tau }\\ =\frac{e^{-{\xi }_r{\omega }_{r}t}\sin \left({\nu }_{r}t-{\theta }_{A1}^{\left(r\right)}\right)}{{\omega }_r\sqrt{{\left(1-q_r^2\right)}^2+4{\xi }_r^2{q}_r^2}}-\frac{\sqrt{q_r^2{\left(1-2{\xi }_r^2-q_r^2\right)}^2+{\xi }_r^2{\left(1+q_r^2\right)}^2}}{{\omega }_r\left[{\left(1-q_r^2\right)}^2+4{\xi }_r^2{q}_r^2\right]}\sin \left(\omega t-{\theta }_{A2}^{\left(r\right)}\right)\end{array}$$

(45)

$$f_B^{\left(r\right)}\left(t\right)=\underset{0}{\overset{t}{\int }}{e}^{-{\xi }_r{\omega }_r(t-\tau )}\sin {\nu }_r(t-\tau )\cos \omega \tau d\tau =-\frac{{\omega }_r{e}^{-{\xi }_r{\omega }_{r}t}\cos ({\nu }_{r}t-{\theta }_{A1}^{\left(r\right)})-{\nu }_r\cos (\omega t-{\theta }_{B2}^{\left(r\right)})}{{\omega }_r^2\sqrt{{(1-q_r^2)}^2+4{\xi }_r^2{q}_r^2}}$$

(46)

where ${\theta }_{A1}^{(r)}=\arctan ({\xi }_r(1+q_r^2))/((1-q_r^2)\sqrt{1-{\xi }_r^2})$ , ${\theta }_{A2}^{(r)}=\arctan ({\xi }_r(1+q_r^2))/(q_r(1-2{\xi }_r^2-q_r^2))$ , ${\theta }_{B2}^{(r)}=\arctan (2{\xi }_r{q}_r)/(1-q_r^2)$ , $q_r=\omega /{\omega }_r$ . Hereinto ${\theta }_{A1}^{\left(r\right)},{\theta }_{A2}^{\left(r\right)},{\theta }_{B2}^{\left(r\right)}\in [0,\pi )$ .

The initial condition is presumed to be

$${\overrightarrow{\mathit{x}}}_0={\left[\begin{array}{ccccc}1,& 0,& 0,& 0,& 0\end{array}\right]}^T,{\stackrel{·}{\overrightarrow{\mathit{x}}}}_0={\left[\begin{array}{ccccc}0,& 1,& 1,& 1,& 1\end{array}\right]}^T$$

(47)

Taking advantage of the expressions of ${\overrightarrow{\mathit{\rho }}}_r^{(A)}$ and ${\overrightarrow{\mathit{\rho }}}_r^{(B)}$ in equation (25), it is possible to achieve

$$\left[\begin{array}{ccccc}{\overrightarrow{\mathit{\rho }}}_1^{\left(A\right)},\hfill & {\overrightarrow{\mathit{\rho }}}_2^{\left(A\right)},\hfill & {\overrightarrow{\mathit{\rho }}}_3^{\left(A\right)},\hfill & {\overrightarrow{\mathit{\rho }}}_4^{\left(A\right)},\hfill & {\overrightarrow{\mathit{\rho }}}_5^{\left(A\right)}\hfill \end{array}\right]=\left[\begin{array}{ccccc}0.0100\hfill & 0.0020\hfill & 0.0004\hfill & 0.0000\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{0.0110}& 0.0000\hfill & -0.0010\hfill & 0.0000\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{0.0087}& -0.0016\hfill & -0.0003\hfill & 0.0001\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{0.0061}& -0.0022\hfill & 0.0004\hfill & 0.0000\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{0.0031}& -0.0015\hfill & 0.0005\hfill & -0.0001\hfill & 0.0000\hfill \end{array}\right]$$

(48)

$$\left[\begin{array}{ccccc}{\overrightarrow{\mathit{\rho }}}_1^{\left(B\right)},\hfill & {\overrightarrow{\mathit{\rho }}}_2^{\left(B\right)},\hfill & {\overrightarrow{\mathit{\rho }}}_3^{\left(B\right)},\hfill & {\overrightarrow{\mathit{\rho }}}_4^{\left(B\right)},\hfill & {\overrightarrow{\mathit{\rho }}}_5^{\left(B\right)}\hfill \end{array}\right]=\left[\begin{array}{ccccc}-0.0078\hfill & -0.0006\hfill & 0.0002\hfill & 0.0000\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{-0.0064}& 0.0012\hfill & 0.0001\hfill & 0.0000\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{-0.0043}& 0.0012\hfill & -0.0003\hfill & 0.0000\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{-0.0026}& 0.0008\hfill & -0.0003\hfill & 0.0000\hfill & 0.0000\hfill \\ \multicolumn{1}{c}{-0.0012}& 0.0003\hfill & -0.0001\hfill & 0.0000\hfill & 0.0000\hfill \end{array}\right]$$

(49)

Then, the response caused by the initial condition, ${\overrightarrow{\mathit{x}}}_{\mathrm{in}}\left(t\right)$ , can be obtained by means of equation (25). Furthermore, the total response, $\overrightarrow{\mathit{x}}\left(t\right)$ , can be derived from equation (27) on this basis. The image of $\overrightarrow{\mathit{x}}\left(t\right)$ is shown in Figure 2. In order to validate the correctness of the theoretical solution of $\overrightarrow{\mathit{x}}\left(t\right)$ , its numerical solution is also acquired by virtue of the Newmark-β method $(\delta =1/2,\beta =1/4)$ , ²⁰ which is plotted as in Figure 2 by the intermittent curves formed by the mark symbols. Figure 2 indicates that the theoretical and numerical solutions tally with each other quite well.

Figure 2.

Accurate response curves under simple harmonic excitation $(\omega =({\omega }_1+{\omega }_2)/2)$ .

From equations (45) and (46), it is easy to make clear that the full response contains three parts: the damped vibration aroused by the initial condition, the damped vibration triggered by the external excitation, and the simple harmonic vibration under the same frequency with the simple harmonic excitation. Evidently the first two parts are both the attenuation vibration and the last one is the steady-state vibration caused by the external excitation. In consequence, only the third part vibration, namely, the steady-state simple harmonic vibration, can be reserved after an enough long time. Figure 2 reflects such characteristics definitely.

As an example for the lower frequency and narrow band external excitation, the frequency of the simple harmonic external excitation is selected to be $\omega =({\omega }_1+{\omega }_2)/2$ . At this time, let the truncation take place at the second-order complex mode to observe the effect. This means that only the complex modes of the first and second orders are kept in the complex mode superposition, that is, letting $S=2$ in equation (28). Thereupon the truncated response can be procured.

The truncated and accelerated responses, $x_1^{*}\left(t\right)$ and $x_{a1}\left(t\right)$ , are taken to be the examples to be investigated at this place since the external excitation locates at the mass 1 as shown by equation (41). Their images are plotted in Figure 3 by the black and red lines, respectively. With the purpose of contrast, the image of $x_1\left(t\right)$ , the precise response of the mass 1, is plotted together in Figure 3 by the blue line.

Figure 3.

Curves of $x_1\left(t\right)$ , $x_1^{*}\left(t\right)$ , and $x_{a1}\left(t\right)(\omega =({\omega }_1+{\omega }_2)/2)$ .

As manifested in Figure 3, if the external excitation mostly covers the low-frequency ingredient, truncating the complex vibration shape superposition at the lower order mode may give rise to a good approximate solution. Moreover, the precision of the truncated approximate solution can be further improved with the aid of the accelerating method given in equation (37). In fact, the remarkable difference between the precise response and the proximate response, including the truncated response and the accelerated response after truncation, may invariably take place at the peaks of the response curve as demonstrated in Figure 3.

Now investigate another instance with the higher frequency and narrow band external excitation, in which $\omega =5{\omega }_5$ . Under such condition, the displacement response curves of the mass 1, 3, and 5 are plotted in Figures 4–6, respectively. In these figures, only the accurate and related truncated response curves are considered and they are expressed by the blue and red lines in several. Because the frequency of the external excitation is higher, the truncation takes place at the fourth-order complex mode and the complex modes on the first four orders are reserved in the superposition.

Figure 4.

Curves of $x_1\left(t\right)$ and $x_1^{*}\left(t\right)$ $(\omega =5{\omega }_5)$ .

Figure 5.

Curves of $x_3\left(t\right)$ and $x_3^{*}\left(t\right)$ $(\omega =5{\omega }_5)$ .

Figure 6.

Curves of $x_5\left(t\right)$ and $x_5^{*}\left(t\right)$ $(\omega =5{\omega }_5)$ .

In Figure 4, the red and blue lines are almost coincident with each other. In Figure 6, the two lines are well coincident in most of the time except a short period of time at the beginning. In Figure 5, the two lines coincide slightly worse. Plainly, even under the circumstance of the high-frequency external excitation with narrow bandwidth, if the truncation is placed at the high-order mode, the precision of the truncated response may still be satisfactory. Meanwhile, if the truncation takes place at the higher order complex mode, the acceleration after the truncation may lose its own significance.

From equations (45) and (46), it is uncomplicated to know

$$\underset{\genfrac{}{}{0ex}{}{t\to \infty }{\omega \to \infty }}{lim}{f}_A^{\left(r\right)}\left(t\right)=\underset{\genfrac{}{}{0ex}{}{t\to \infty }{\omega \to \infty }}{lim}{f}_B^{\left(r\right)}\left(t\right)=0$$

such that $\underset{\genfrac{}{}{0ex}{}{t\to \infty }{\omega \to \infty }}{lim}\overrightarrow{\mathit{x}}\left(t\right)=0$ .

Figures 4–6 reflect such a tendency well.

Slope step external excitation

The expression of the slope step external excitation can be expressed as

$$\overrightarrow{\mathit{P}}(t)={\overrightarrow{\mathit{P}}}_{0}p\left(t\right),p\left(t\right)=\{\begin{array}{c}\frac{p_0}{t_0}t0\le t<t_0\\ p_{0}t\ge t_0\end{array}$$

(50)

where ${\overrightarrow{\mathit{P}}}_0=[\begin{array}{ccccc}0,\hfill & 1,\hfill & 0,\hfill & 0,\hfill & 0\hfill \end{array}{]}^T$ , $p_0=1000$ N, and $t_0=20\pi /{\omega }_1$ . The graph of the slope step function, $p\left(t\right)$ , is drawn in Figure 7.

Figure 7.

Graph of the slope step function $p\left(t\right)$ .

With the aid of equation (68) in Appendix 2, the amplitude spectrum of the slope step function, $p\left(t\right)$ , can be deduced as

$$\left|F_p\left(\omega \right)\right|=p_0\sqrt{\frac{2(1-\cos \omega t_0)}{{\omega }^4{t}_0^2}[1+\pi {\omega }^2{t}_0\delta \left(\omega \right)]+{\pi }^2\delta \left(\omega \right)}$$

(51)

where $\delta (•)$ is the so-called Dirac function. Equation (51) demonstrates the slope step excitation is a wideband excitation.

From equation (27), the response caused by the slope step external excitation over the time interval $[0,t_0]$ can be acquired as

$${\overrightarrow{\mathit{x}}}_e\left(t\right)=\frac{2p_0}{t_0}{\displaystyle \sum _{r=1}^5\frac{1}{{\left|a_r\right|}^2}\left[g_A^{\left(r\right)}\left(t\right){\overrightarrow{\mathit{Q}}}_r^{\left(A\right)}-g_B^{\left(r\right)}\left(t\right){\overrightarrow{\mathit{Q}}}_r^{\left(B\right)}\right]},0\le t<t_0$$

(52)

where by right of equations (76) and (79) in Appendix 2, it is undemanding to have

$$g_A^{\left(r\right)}\left(t\right)=\underset{0}{\overset{t}{\int }}\tau e^{-{\xi }_r{\omega }_r(t-\tau )}\cos {\nu }_r(t-\tau )d\tau =\frac{{\xi }_r{\omega }_{r}t+1-2{\xi }_r^2}{{\omega }_r^2}-\frac{e^{-{\xi }_r{\omega }_{r}t}}{{\omega }_r^2}\sin ({\nu }_{r}t+{\phi }_r)$$

(53)

$$g_B^{\left(r\right)}\left(t\right)=\underset{0}{\overset{t}{\int }}\tau e^{-{\xi }_r{\omega }_r(t-\tau )}\sin {\nu }_r(t-\tau )d\tau =\frac{\sqrt{1-{\xi }_r^2}}{{\omega }_r^2}({\omega }_{r}t-2{\xi }_r)+\frac{e^{-{\xi }_r{\omega }_{r}t}}{{\omega }_r^2}\cos ({\nu }_{r}t+{\phi }_r)$$

(54)

where ${\phi }_r=\arctan ((1-2{\xi }_r^2)/(2{\xi }_r\sqrt{1-{\xi }_r^2}))$ .

Similarly, the response caused by the slope step external excitation over the time interval $[t_0,+\infty ]$ can be acquired as

$${\overrightarrow{\mathit{x}}}_e\left(t\right)=\frac{2p_0}{t_0}\sum _{r=1}^5\frac{1}{{\left|a_r\right|}^2}[{\hat{g}}_A^{\left(r\right)}\left(t\right){\overrightarrow{\mathit{Q}}}_r^{\left(A\right)}-{\hat{g}}_B^{\left(r\right)}\left(t\right){\overrightarrow{\mathit{Q}}}_r^{\left(B\right)}],t_0\le t<+\infty$$

(55)

where with the aid of equations (71) and (75) in Appendix 2, it is possible to attain

$$\begin{array}{c}{\widehat{g}}_A^{\left(r\right)}\left(t\right)={\displaystyle \underset{0}{\overset{t_0}{\int }}\tau e^{-{\xi }_r{\omega }_r\left(t-\tau \right)}\cos {\nu }_r\left(t-\tau \right)d\tau }+t_0{\displaystyle \underset{t_0}{\overset{t}{\int }}{e}^{-{\xi }_r{\omega }_r\left(t-\tau \right)}\cos {\nu }_r\left(t-\tau \right)d\tau }\\ =t{g}_{A1}^{\left(r\right)}\left(t\right)-\left(t-t_0\right)g_{A1}^{\left(r\right)}\left(t-t_0\right)-g_{A2}^{\left(r\right)}\left(t\right)+g_{A2}^{\left(r\right)}\left(t-t_0\right)\\ =\frac{{\xi }_r{t}_0}{{\omega }_r}-\frac{e^{-{\xi }_r{\omega }_{r}t}}{{\omega }_r^2}\sin \left({\nu }_{r}t+{\phi }_r\right)+\frac{e^{-{\xi }_r{\omega }_r\left(t-t_0\right)}}{{\omega }_r^2}\sin \left[{\nu }_r\left(t-t_0\right)+{\phi }_r\right]\end{array}$$

(56)

$$\begin{array}{c}\multicolumn{1}{c}{{\hat{g}}_B^{\left(r\right)}\left(t\right)=\underset{0}{\overset{t_0}{\int }}\tau e^{-{\xi }_r{\omega }_r(t-\tau )}\sin {\nu }_r(t-\tau )d\tau }\\ \multicolumn{1}{c}{+t_0\underset{t_0}{\overset{t}{\int }}{e}^{-{\xi }_r{\omega }_r(t-\tau )}\sin {\nu }_r(t-\tau )d\tau }\\ \multicolumn{1}{c}{=\frac{t_0}{{\nu }_r}-\frac{{\xi }_r{\omega }_{r}t+1}{{\nu }_r}{g}_{A1}^{\left(r\right)}\left(t\right)+\frac{{\xi }_r{\omega }_r(t-t_0)t+1}{{\nu }_r}{g}_{A1}^{\left(r\right)}(t-t_0)}\\ \multicolumn{1}{c}{+\frac{{\xi }_r{\omega }_r}{{\nu }_r}[g_{A2}^{\left(r\right)}\left(t\right)-g_{A2}^{\left(r\right)}(t-t_0)]}\\ \multicolumn{1}{c}{=\frac{{\nu }_r{t}_0}{{\omega }_r^2}+\frac{e^{-{\xi }_r{\omega }_{r}t}}{{\omega }_r^2}\cos ({\nu }_{r}t+{\phi }_r)}\\ \multicolumn{1}{c}{-\frac{e^{-{\xi }_r{\omega }_r(t-t_0)}}{{\omega }_r^2}\cos [{\nu }_r(t-t_0)+{\phi }_r]}\end{array}$$

(57)

If presuming that ${\overrightarrow{\mathit{x}}}_0={\stackrel{·}{\overrightarrow{\mathit{x}}}}_0=0$ , the total response, $\overrightarrow{\mathit{x}}\left(t\right)$ , can be obtained with the help of equations (52) and (55), whose image is drawn in Figure 8. In order to verify the correctness of the theoretical solution of $\overrightarrow{\mathit{x}}\left(t\right)$ , its numerical solution is also attained by right of the Newmark-β method $(\delta =1/2,\beta =1/4)$ , which is together plotted as in Figure 8 by the intermittent curves formed by the mark symbols. Figure 8 illustrates that the theoretical and numerical solutions are in keeping with each other quite well.

Figure 8.

Response curves under slope step external excitation.

Equations (52) and (55) state clearly that the total dynamic response, $\overrightarrow{\mathit{x}}\left(t\right)$ , is a piecewise function. Equations (53) and (54) depict that the response over the time interval $[0,t_0]$ contains two parts: one is a linear motion and the other is a degenerative simple harmonic vibration. As depicted by equations (56) and (57), the response over the time interval $[t_0,+\infty ]$ is a damped simple harmonic vibration around a new equilibrium position. Moreover, the simple harmonic sections in the response suffer from attenuation with the passage of time following exponential law. Figure 8 reveals such properties very well.

Here, the truncating and accelerating characteristics of the responses of the masses 1 and 3 are selected to be the instances for investigation. The images of $x_1\left(t\right)$ , $x_1^{*}\left(t\right)$ , and $x_{a1}\left(t\right)$ are drawn together in Figure 9, among which $x_1^{*}\left(t\right)$ and $x_{a1}\left(t\right)$ are obtained under the condition of $S=2$ . Such a figure indicates that the curves of $x_1\left(t\right)$ and $x_{a1}\left(t\right)$ are nearly of coincidence. In consequence, for the response of the mass 1, the acceleration can observably boost the precision of the truncating solution and can supply a favorable approximate solution.

Figure 9.

Graphs of $x_1\left(t\right)$ , $x_1^{*}\left(t\right)$ , and $x_{a1}\left(t\right)$ under slope step external excitation.

In Figure 10, the graphs of $x_3\left(t\right)$ , $x_3^{*}\left(t\right)$ , and $x_{a3}\left(t\right)$ almost coincide mutually, among which $x_3^{*}\left(t\right)$ and $x_{a3}\left(t\right)$ are attained in the case of $S=2$ as well. It is thus evident that, under some certain conditions, merely the truncation without acceleration may lead up to a beneficial approximate solution to the response as well.

Figure 10.

Graphs of $x_3\left(t\right)$ , $x_3^{*}\left(t\right)$ , and $x_{a3}\left(t\right)$ under slope step external excitation.

Although the slope step excitation belongs to the broadband excitation, its leading energy chiefly distributes over the lower frequency ingredient as shown by the expression of its spectrum, that is, equation (51). Therefore, the truncation and acceleration can give rise to a nice proximate solution for the dynamical response like the case of the low-frequency excitation even if the truncation occurs at the lower order complex mode.

Conclusion

For a symmetrical system with general linear viscous damping, a method to solve its accurate dynamic response is proposed. In essence, the method is an improvement of the traditional complex mode technique. In the process of employing such an approach, the number of the differential equations to be solved need not be increased actually although the concept of state equation is quoted. Meanwhile, solving quadratic characteristic value problem is not required. In addition, by reason of taking advantage of the conjugacy of the complex characteristic value and vector, the complex operation is avoided completely. In comparison with the integration transformation method, the improved complex mode method suggested is more simple and convenient because it pertains to the time-domain method essentially and the integration transformation is not needed.

The truncation approach of the complex mode superposition is proposed and investigated using the improved complex mode theory established. The complex mode expansion of the flexibility matrix of the system and some other important and basic relations are discovered. Based on these, an accelerating method of the complex mode superposition after truncation is put forward.

For the low-frequency excitation with narrow band, only keeping the lower order modes during the superposition and accelerating the truncated solution can generally give a good approximate solution to the dynamic response. Under such circumstance, the accelerated solution after truncation is always more accurate than the pure truncated one.

With regard to the high-frequency excitation with narrow band, the higher order modes should usually be reserved in the superposition to guarantee the accuracy of the truncated solution. Under such a condition, the acceleration after truncation may not be necessary.

As to the broadband excitation, if the superimposition only comprises the complex modes of the lower orders, a nice approximate solution to the dynamical response can also be obtained in general. Furthermore, the acceleration after the truncation can be utilized to enhance the precision of the truncated solution.