A hybrid expansion method for frequency response functions of non-proportionally damped systems

Mechanical Systems and Signal Processing 42 (2014) 31–41

Contents lists available at ScienceDirect

Mechanical Systems and Signal Processing journal homepage: www.elsevier.com/locate/ymssp

A hybrid expansion method for frequency response functions of non-proportionally damped systems Li Li, Yujin Hu n, Xuelin Wang, Lei Lü School of Mechanical Science and Engineering, Huazhong University of Science and Technology, Wuhan 430074, People's Republic of China

a r t i c l e i n f o

abstract

Article history: Received 18 April 2013 Received in revised form 28 July 2013 Accepted 31 July 2013 Available online 17 September 2013

This study is aimed at eliminating the influence of the higher-order modes on the frequency response functions (FRFs) of non-proportionally viscously damped systems. Based on the Neumann expansion theorem, two power-series expansions in terms of eigenpairs and system matrices are derived to obtain the FRF matrix. The relationships satisfied by eigensolutions and system matrices are established by combining the two power-series expansions. By using the relationships, an explicit expression on the contribution of the higher-order modes to FRF matrix can be obtained by expressing it as a sum of the lower-order modes and system matrices. A hybrid expansion method (HEM) is then presented by expressing FRFs as the explicit expression of the contribution of the higher-order modes and the modal superposition of the lower-order modes. The HEM maintains original-space without having to use the state-space equation of motion such that it is efficient in computational effort and storage capacity. Finally, a two-stage floating raft isolation system is used to illustrate the effectiveness of the derived results. & 2013 Elsevier Ltd. All rights reserved.

Keywords: Modal truncation error Viscous damping FRF Neumann expansion Modal analysis Mode superposition method

1. Introduction The frequency response functions (FRFs) of mechanical and structural systems are of interest in dynamic problems subjected to harmonically varying loading that may be caused by reciprocating or rotating machine parts including motors, fans, compressors, and forging hammers [1,2]. FRFs are of fundamental importance and play a very important role in many areas such as model updating [3,4], structural damage detection [5,6], vibration and noise control [7,8], system identification [9–11], dynamic optimization [12,13] and many other applications. Two kinds of methods, i.e., direct frequency response method (DFRM) and modal superposition method, are usually used to calculate the FRFs. The DFRM is based on the direct frequency results in an exact calculation by solving the inversion of the dynamic stiffness matrix directly. The modal superposition method calculates the FRFs by expressing them as the summation of the contributions of all the modes. The modal superposition method has been extensively used in many dynamic analyses, and also programmed in some commercial software, e.g., NASTRAN, ANSYS or ABAQUS. However, the modal superposition method requires that all the frequencies and mode shapes should be available. As we know, often it is difficult, or even unnecessary, to obtain all the eigenpairs of a large-scaled model, which means that the modal truncation scheme is generally used and the modal truncation error is therefore introduced. As a result, the quality of the calculated FRFs may be adversely affected. The corrections to the modal truncation scheme have been investigated by several authors. Model acceleration methods [14–18] approximates the contribution of higher-order (unavailable) modes to dynamic response in terms of a pseudo-static

n

Corresponding author. Tel.: +86 27 87543972. E-mail addresses: [email protected], [email protected] (L. Li), [email protected] (Y. Hu).

0888-3270/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ymssp.2013.07.020

32

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

term, which is a particular solution of the dynamic equation of motion (the particular solution can be considered as the exciting frequency equals zero). Since the mode acceleration method neglects the contribution of both velocity and acceleration terms to the dynamic response, it can be considered as a state approximation method. Dynamic correction methods [19–21], which were developed to improve the accuracy of mode acceleration method, include the contribution of higher-order modes by a sum of the particular solutions of both the equation of motion and the reduced differential equation of motion. Force derivative methods [22,23] reduce the modal truncation error by considering the higher-order derivatives of the forcing function, which means the forcing function should be described by analytical laws. The force derivative methods take advantage of the fact that successive integrations by parts of the convolution integral produce terms, which can be expressed as the combination form of system matrices, the forcing function and its derivatives. In addition, the correction methods for stochastic systems have been studied by Refs. [24–26] and can give a corrected response for both the deterministic and random forcing function. Bilello et al. [27,28] studied the corrected methods for continuous systems. For the past two decades, high accurate modal superposition methods [29–33] were developed to the problem on the correction to the modal truncation scheme by combining the mode superposition of the lower (available) modes and a power-series expansion of dynamic response in terms of system matrices. These methods have been applied widely in the sensitivity of mode shapes [34–36] and the sensitivity of responses in the frequency domain [37–39]. Recently, Qu [40] presented an adaptive mode superposition and acceleration technique to solve the problem how many items of the convergent power-series expansion of dynamic response should be considered to satisfy the necessary accuracy. However, the high accurate modal superposition methods [36,39] approximate the FRFs using 2N-space (state-space) formulation, where N is the system dimension. Although these 2N-space correction methods are exact in nature, the 2N-space correction method usually needs heavy computational cost for real-life multiple degrees-of-freedom (DOF) systems since the size of system matrices of state-space equations is double. More recently, Li et al. [41] developed N-space correction methods to calculate the FRFs of nonviscously (viscoelastically) damped systems. The N-space correction methods attempt to approximate the influence of higher-order modes in terms of the lower-order modes and system matrices by using the first one or two terms of Neumann expansion of the contribution of higher-order modes. However, these procedures cannot be extended to further high-order terms since all of them will be affected by the nonviscous damping matrix which is frequency-dependent. In this study, the correction problem of the modal truncation scheme of non-proportionally viscously damped systems is studied. The aim of this paper is to propose an N-space power-series expansion method to the problem on the correction of FRFs. Based on the Neumann expansion theorem, two power-series expansions in terms of eigenpairs and system matrices are derived to obtain the FRF matrix. By using the two power-series expansions, an explicit expression on the contribution of the higher-order modes can be expressed as a sum of the lower-order modes and system matrices. Then, a hybrid expansion method is presented by expressing the FRFs as the explicit expression of the contribution of the higher-order modes and the modal superposition of the lower-order modes. The second section of this paper simply reviews the dynamic of non-proportionally viscously damped systems. The third section presents two power-series expansions to compute the FRF matrix. The fourth section gives some relationships between eigensolutions and system matrices, and presents an N-space power-series expansion method to the problem on the correction of FRF matrix and the displacement vectors. And the fifth section presents illustrate the engineering application, accuracy and efficiency of the presented method by a two-stage floating raft isolation system.

2. Dynamic of viscously damped systems The dynamic equation of motion for a viscously damped system with N DOF in Laplace domain can be expressed as ðs2 M þ sC þ KÞXðsÞ ¼ FðsÞ or DðsÞXðsÞ ¼ FðsÞ

ð1Þ

where M, C and K A ℝ are, respectively, the mass, damping and stiffness matrices (only consider symmetric system matrices in this study), FðsÞ is the forcing vector and XðsÞ is the displacement vector. The matrix DðsÞ ¼ s2 M þ sC þ K is pffiffiffiffiffiffiffi so-called the dynamic stiffness matrix. In the context of structural dynamics, s¼iω, where i¼ 1 and ω denotes the exciting frequency in rad/s. The viscously damped system cannot be simultaneously decoupled by modal analysis unless it also possesses a full set of classical normal modes. The condition of viscously damped systems to possess classical normal modes (known as the proportionally damped system), originally introduced by Rayleigh [42] in 1877, is still extensively used. It shows that a viscous damping is proportionally damping if the damping matrix is a linear combination of inertia and stiffness matrices. This damping is routinely assumed in engineering applications. Later, Caughey and O’Kelly [43] and Adhikari [44] gave some more restrictive conditions which make damped systems possess normal modes as well. Generally speaking, proportional damping means that energy dissipation is almost uniformly distributed throughout the mechanical system [45]. However, there is no reason why these mathematical conditions must be satisfied. In practical, systems with two or more parts with significantly different levels of energy dissipation are encountered frequently in engineering designs. To this end, the non-proportionally damped system is considered in this study, i.e., the concern of this study is when these mathematical conditions are not met, the most general case in engineering applications. NN

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

33

From the complex mode superposition theorem, the complex FRF (receptance) matrix HðsÞ ¼ DðsÞ1 can be expressed in terms of the eigenvalues and eigenvectors of viscously damped systems [46,47] φj φTj

2N

HðsÞ ¼ ∑

ð2Þ

j ¼ 1 θ j ðsλj Þ

where θj ¼ φTj ð2λj M þ CÞφj

ð3Þ

where φi denotes the eigenvector corresponding to the ith eigenvalueλi . It should be noted that Eq. (2) is not working in the case of systems with repeated eigenvalues (in that case, the Jordan form is needed in place of a set of independent eigenvectors, see Newland [48] or Lancaste [49] for further reading). It means that in this study the eigenvalues should be distinct. In addition, we assume the number of eigenvalues is 2N (i.e., the system is underdamped and has N complex conjugate pairs of eigenvalues), the most likely case for many practical applications. By using HðsÞ ¼ DðsÞ1 , one can obtain the displacement vector from Eq. (1) as the following form: 2N

XðsÞ ¼ ∑

j¼1

φTj FðsÞφj

ð4Þ

ðsλj Þθj

Increasing the number of DOF used in finite element analysis, exact eigensolutions of large-scaled structures are very computationally time-consuming. Practically, often only a few lower mode shapes are considered in the dynamic analysis of large-scale systems. In that case, the lower-order (available) modes can be only used to approximate FRF matrix and displacement vector as the following forms: L

HMDM ðsÞ ¼ ∑

φj φTj

j ¼ 1 θ j ðsλj Þ

L

and XMDM ðsÞ ¼ ∑

φTj FðsÞφj

j¼1

ð5Þ

θj ðsλj Þ

Eq. (5) is called modal displacement method (MDM). Since some higher-order modes are ignored, the modal truncation error may be introduced. When L«2N, the results obtained by MDM may be not only inaccurate, but also misleading. 3. Two power-series expansions for FRF matrix Neumann expansion is of fundamental importance and is often used to approximate the inverse matrix [50]: ðIN þAÞ1 ¼ IN þ AA2 þ A3 A4 þ ⋯

ð6Þ

Here A is a (N  N) matrix and Ie denotes the identity matrix of size e. If all the eigenvalues of A have absolute values less than one, the previous series-expansion converges to the exact results. In this section, two power-series expansion methods will be presented to compute FRF matrix based on the Neumann expansion theorem. These two power-series expansions allow us to establish some relationships between eigensolutions and system matrices, which can be used to correct the mode truncation error of non-proportionally viscously damped systems. 3.1. Power-series expansion in terms of eigensolutions By casting Eq. (2) into the matrix form, one obtains HðsÞ ¼ UΘ1 Λ1 ðI2N sΛ1 Þ1 UT

ð7Þ

where Λ ¼ diag½ λ1 ;

λ2 ;

⋯;

λ2N ; U ¼ ½ φ1 ;

φ2 ;

⋯;

 φ2N  and Θ ¼ diag θ1 ;

θ2 ;

⋯;

θ2N



A Neumann series-expansion can be used to expand the inverse term ðI2N sΛ1 Þ1 1

ðI2N sΛ1 Þ1 ¼ I2N þ sΛ1 þ s2 Λ2 þ s3 Λ3 þ ⋯ ¼ ∑ sr1 Λrþ1 r¼1

ð8Þ

It should be emphasized that the previous series-expansion converges only at the region |s| o|λ1|. That is, the seriesexpansion has the convergence results if and only if the frequency range is chosen at 0–|λ1| rad/s. Substituting Eq. (8) into Eq. (7), one obtains r1 UΘ1 Λr UT HðsÞ ¼ ∑1 r ¼ 1s

ð9Þ

As can be seen, the FRF matrix derived here is in terms of complex eigensolutions. As mentioned before, when the frequency range is chosen at 0–|λ1| rad/s, the series-expansion in Eq. (9) can give convergence results. For this reason, the series-expansion is not interested in engineering, but it can allow us to establish some relationships between eigensolutions and system matrices.

34

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

3.2. Power-series expansion in terms of system matrices The expression of the dynamic stiffness matrix D(s) can be reformed as DðsÞ ¼ K ½IN þ sK1 ðC þ sMÞ

ð10Þ

Taking the inverse of the previous equation yields DðsÞ1 ¼ ½IN þ sK1 ðC þ sMÞ1 K1

ð11Þ

By using Neumann expansion, we expand the inverse of matrix in the bracket of the previous equation ½IN þ sK1 ðC þ sMÞ1 ¼ IN þ ðsÞK1 ðC þ sMÞ þ s2 K1 ðC þ sMÞK1 ðC þ sMÞ þ ⋯ 2 3 þðsÞk 4K1 ðC þ sMÞ⋯K1 ðC þ sMÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

5þ⋯

ð12Þ

k terms

which can be simplified as 1

½IN þ sK1 ðC þ sMÞ1 ¼ IN þ ðsÞK1 C þ s2 ðK1 CK1 CK1 MÞ þ ⋯ ¼ ∑ ðsÞk Ξk k¼0

ð13Þ

where Ξ0 ¼ IN ; Ξ1 ¼ K1 C and Ξk ¼ K1 CΞk1 K1 MΞk2

ð14Þ

Note that matrices Ξk only need to be obtained once for different exciting frequencies since Ξk are independent to frequency. More importantly, Ξk can be determined by a simple iteration process (the proof of the iteration process can be seen in Appendix A). Using relationship H(s) ¼D  1(s), and substituting Eqs. (13) and (14) into Eq. (11), one obtains 1

HðsÞ ¼ ∑ ðsÞk Ξk K1

ð15Þ

k¼0

From the previous equation, the FRF matrix can be simplified as 1

HðsÞ ¼ ∑ sk Γk

ð16Þ

Γ0 ¼ K1 ; Γ1 ¼ K1 CK1 and Γk ¼ K1 CΓk1 K1 MΓk2

ð17Þ

k¼0

with

As can be seen from Eq. (16), the FRF matrix derived here only involves system matrices. 4. Hybrid expansion method for frequency response functions 4.1. Relationships between eigensolutions and system matrices The purpose of the section is to develop some relationships satisfied by eigensolutions and system matrices. Combining Eqs. (9) and (16), one has 1

1

r¼1

k¼0

 ∑ sr1 UΘ1 Λr UT ¼ ∑ sk Γk

ð18Þ

Equating the coefficients of s0, …, s3, several relationships involving eigensolutions and system matrices can be obtained: UΘ1 Λ1 UT ¼ Γ0 ¼ K1

ð19Þ

UΘ1 Λ2 UT ¼ Γ1 ¼ K1 CK1

ð20Þ

UΘ1 Λ3 UT ¼ Γ2 ¼ K1 CK1 CK1 K1 MK1

ð21Þ

UΘ1 Λ4 UT ¼ Γ3 ¼ K1 CðK1 CK1 CK1 K1 MK1 Þ þ K1 MK1 CK1

ð22Þ

It should be noted that Fawzy and Bishop [51] and Adhikari [52,53] obtained the same relationships expressed by Eqs. (19), (20) and (21). Here we give the more general relationships satisfied by eigensolutions and system matrices by equating the coefficients of sk for k¼ 0, 1, 2,…, 1 UΘ1 Λk1 UT ¼ Γk

ð23Þ

So far, the relationships satisfied by eigensolutions and system matrices are established. Next, we use these relationships to correct the modal truncation problem of FRF matrix and the displacement vector.

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

35

4.2. Corrected methods to calculate frequency response functions Separating the right side of Eq. (7) in terms of the lower-order (available) modes and the higher-order (unavailable) modes, the FRF matrix can be expressed as 1 1 1 T 1 1 1 1 T HðsÞ ¼ UL Θ1 L ΛL ðIL sΛL Þ UL UH ΘH ΛH ðIH sΛH Þ UH

ð24Þ

where ΛL ¼ diag½ λ1 ;

λ2 ;

⋯;

λL ; UL ¼ ½ φ1 ;

φ2 ;

⋯;

ΛH ¼ diag½ λLþ1 ;

λLþ2 ;

⋯;

λ2N ; UH ¼ ½ φLþ1 ;

ΘH ¼ diag½ θLþ1 ;

θLþ2 ;

⋯;

θ2N 

φL ; ΘL ¼ diag½ θ1 ; φLþ2 ;

⋯;

θ2 ;

⋯;

θL ;

φ2N ; and

Based on Neumann expansion (applying a similar procedure in Section 3.1), the second term of the right side of Eq. (24), i.e., the contribution of the higher-order modes to FRF matrix can be expressed as 1

1 1 1 T r1 r T UH Θ1 UH Θ1 H ΛH ðIH sΛH Þ UH ¼  ∑ s H ΛH UH r¼1

ð25Þ

It should be noted that the previous series-expansion has the convergence results when the frequency range is chosen at 0–|λL+1| rad/s and enlarge the convergence region given in Section 3.1. Using the relationships given by Eq. (23), one has k1 T k1 T UH Θ1 UH ¼ Γk þ UL Θ1 UL for k ¼ 0; 1; 2; …; 1 H ΛH L ΛL

ð26Þ

Then, Eq. (25) can be expressed in terms of the lower-order modes and system matrices 1

1 1 1 T r1 r T ½Γr1 þ UL Θ1 UH Θ1 H ΛH ðIH sΛH Þ UH ¼ ∑ s L Λ L UL  r¼1

ð27Þ

Substituting the previous equation into Eq. (24), a novel mode superposition method to calculate FRF matrix can be obtained " # L 1 L φj φT φj φTj j þ ∑ sr1 Γr1 þ ∑ HHEM ðsÞ ¼ ∑ ð28Þ r r¼1 j ¼ 1 θ j ðsλj Þ j ¼ 1 θ j λj Eq. (28) is proposed based on the two power-series expansions and the mode superposition of the lower-order modes, therefore Eq. (28) is called the hybrid expansion method (HEM). Since the effect of the higher-order modes is also considered, the HEM can improve the accuracy of the approximate FRF matrix obtained by MDM. For engineering application, a few terms in power-series expressed by the second term of the right-hand of Eq. (28) are considered for suitable accuracy requirements. Suppose the first h terms is retained, the FRF matrix computed by HEM is " # L h L φj φT φj φTj j þ ∑ ðiωÞr1 Γr1 þ ∑ HHEM ðiωÞ ¼ ∑ ð29Þ r r¼1 j ¼ 1 θ j ðiωλj Þ j ¼ 1 θj λj By using X(iω)¼H(iω)F(iω), the displacement vectors approximated by HEM can be easily obtained as " # L φT FðiωÞφj h L φT FðiωÞφj j j þ ∑ ðiωÞr1 Γr1 FðiωÞ þ ∑ XHEM ðiωÞ ¼ ∑ θj λrj r¼1 j ¼ 1 θ j ðiωλj Þ j¼1

ð30Þ

As can be noted, the accuracy of displacement vectors can be also improved compared to that calculated by Eq. (5) (MDM) since the FRF matrix is computed by HEM. As can be observed from the previous equation, the first term of HEM is identical with the MDM, and therefore the second term (the power-series expansion) can be considered as some corrections to the MDM. As it can be known from the theorem of Neumann expansion, when the number h is increased, the errors of FRFs can be decreased, but the calculations of HEM will suffer from more computational effort. As mentioned before, the HEM has the convergence region at 0–|λL+1| rad/s. It is implied that, although not all of the modes are required in the HEM, it needs to retain, at a minimum, all the modes whose resonant frequencies lie within the range of exciting frequencies. It should be noted that matrices Γr1 can be determined by an iteration process and only need to be obtained once for different exciting frequencies since Γr1 are independent to frequency. In addition, the HEM maintains original-space without having to use the state-space equation of motion. In this sense, the HEM can be suitable in dynamic analysis of engineering structures. 5. Example study To illustrate the accuracy, efficiency and engineering application, of the proposed method, a two-stage floating raft isolation system [54], shown in Fig. 1 is considered. The machine vibration is attempted to be isolated by mounting two machines on a single intermediate raft structure. The two machines have the values: m1 ¼200 kg, m2 ¼250 kg. The length– width–thickness of the raft and foundation plates is 1200–800–20 mm and 2000–800–40 mm, respectively. The two short sides of the foundation plate are clamped and two long sides of it are free, while the four sides of the raft plate are all free. The raft and foundation plates have Young's modulus E¼ 2.0  1011 N/m2, material density ρ¼7800 kg/m3 and Poisson's

36

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

Fig. 1. The schematic of a two-stage floating raft isolation system.

Fig. 2. The finite element model of the raft isolation system.

ratio μ¼0.3. Two springs and dashpots are used to mount two machines on a single intermediate raft structure, and eight springs and dashpots are fixed between the raft and foundation plates. These spring stiffness and damping coefficients have the following values: k1 ¼ 1:0  105 N=m; k2 ¼ 5:0  105 N=m; c1 ¼ 1:0  102 Ns=m; c2 ¼ 5:0  102 Ns=m The finite element model, shown in Fig. 2, is modeled with 1258 DOF using the finite element software-NASTRAN. System matrices of the finite element model are assembled by NASTRAN and exported by means of DAMP language. By using these system matrices, the proposed HEM is programmed in MATLAB. The distribution of the nonzero terms of these sparse system matrices obtained using the MATLAB function spy(  ) is shown in Fig. 3. It has been checked that this floating raft isolation system is a non-proportionally damped system. The lower 20 complex eigenvalues are listed in Table 1. It should be noted that the imaginary part of each complex eigenvalue is the damped natural circle frequency in rad/s. The pffiffiffiffiffiffiffiffi undamped natural frequency ωin is the norm of the complex eigenvalue λi , i.e., ωin ¼ λi λni where the superscript n denote the complex conjugate. The real part of each complex eigenvalue is a measure of the decay rate of a damped system. The system is stable since all the real parts of eigenvalues are negative. FRFs are by definition the amplitude of the response of the system under study to the harmonic excitation with a unit force amplitude at each frequency point. Resonance and antiresonance frequencies can be determined by phase changes and the minima of the FRFs. To show the capacity of the proposed method to calculate both FRF matrix and displacement vectors, we assume that a unit harmonic force is located at the machine which has the mass m2. Suppose that the first 9 pair complex modes are available to calculate the FRFs. This means that the modes of the 19th and the higher-order modes are truncated when the complex modal truncation is adopted. In order to show the accuracy of the presented HEM, the FRF matrix has been also calculated by MDM (i.e., Eq. (5)) and DFRM which calculates the FRFs by solving the inversion of the dynamic stiffness matrix directly. As we know, the DFRM is an exact method. Assume the frequency points of concern are considered at the frequency ranges, 500–700 rad/s. The approximation results of one FRF are plotted in Fig. 4. As can be observed, the errors of MDM may be significant if only a few lower modes are retained especially for the FRF around the antiresonance frequencies, and HEM can make the errors reduce. The errors of the amplitude of the FRF computed in terms of MDM and HEM are shown in Fig. 5. As can be seen from Fig. 5, the errors of the amplitude of the FRF, which computed using HEM when the number h equals 1 and 2, are almost the same. The same phenomenon can be observed when the number h equals 3 and 4. To show the errors of the FRF can be decreased when the number h considered in HEM is increased, the phase of the approximated results in terms of HEM with h¼1 and h¼2 is shown in Fig. 6. It is shown that the phase of the FRF computed by HEM with h ¼2 are more accurate than that calculated by HEM with h¼ 1. In addition, as can be seen from Fig. 7, the approximated phase of the FRF computed by HEM with h¼4 are more accurate than that calculated by HEM with h ¼3. As can be observed from Figs. 4 and 7, the approximated results calculated in terms of HEM with h¼4 can give almost the same accuracy with the exact results computed by DFRM. For the proposed HEM, 2001 times steps, i.e., Δω¼0.1 rad/s, are applied to calculate the FRFs. The computed time of HEM with h¼4 is 34.85 s, while the DFRM takes 170.09 s for this problem. As can be seen, the latter is much more computationally expensive than the former. The reason is that when DFRM is used for this problem, 2001 times of decomposition, forward and back-substitutions of the dynamic

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

0

0

200

200

400

DOF

400

DOF

37

600

600

800

800

1000

1000

1200

1200 0

200

400

600

800

1000

1200

0

200

400

DOF

600

800

1000

1200

DOF

0 200

DOF

400 600 800 1000 1200 0

200

400

600

800

1000

1200

DOF Fig. 3. The distribution of the nonzero terms of system matrices. (a) Mass matrix; (b) Damping matrix; (c) Stiffness matrix. Table 1 The lower 20 eigenvalues and frequencies of the floating raft isolation system. Mode 1, 2 3, 4 5, 6 7, 8 9, 10 11, 12 13, 14 15, 16 17, 18 19, 20

Eigenvalue

Damped natural frequency (Hz) 4

 1.8816  10 7 19.543i  2.3768  10  4 7 21.923i  9.3877  10  3 7 152.17i  1.7662  10  2 7 198.96i  2.9400  10  2 7 259.33i  6.2425  10  3 7 355.05i  4.7305  10  2 7 521.52i  2.5394  10  2 7 559.31i  1.9205  10  2 7640.00i  4.6052  10  3 7 941.99i

3.1103 3.4891 24.219 31.666 41.273 56.507 83.002 89.018 101.86 149.92

stiffness matrix D(iω) should be required to solve algebraic Eq. (1). If the number of frequency points of concern or the DOF becomes larger, the computational time of the DFRM will increase rapidly. 6. Concluding remarks FRFs have been extensively applied to the areas of model updating, structural damage detection, vibration and noise control, system identification, dynamic optimization and many other applications. This study is aimed at eliminating the influence of the higher-order modes on the frequency response functions (FRFs) of non-proportionally viscously damped systems. But here assume that the eigenvalues of the system are distinct and the system is underdamped such that it has N complex conjugate pairs of eigenvalues. Based on the Neumann expansion theorem, two power-series expansions in terms

38

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

−160

DFRM MDM h=1 h=2 h=3 h=4

Amplitude (dB)

−180

−200

−220

−240

−260

−280 500

520

540

560

580

600

620

640

660

680

700

Frequency (rad/s) Fig. 4. Amplitude of one FRF of the floating raft isolation system.

102 MDM h=1 h=2 h=3 h=4

Error of Amplitude (%)

101

100

10−1

10−2

10−3

10−4 500

520

540

560

580

600

620

640

660

680

700

Frequency (rad/s) Fig. 5. Errors of the amplitude of the FRF computed in terms of MDM and HEM.

4 3

DFRM MDM h=1 h=2

Phase (rad)

2 1 0 −1 −2 −3 −4 500

520

540

560

580

600

620

640

660

680

Frequency (rad/s) Fig. 6. Phase of the FRF in terms of HEM with h ¼ 1 and h ¼ 2.

700

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

39

4 DFRM MDM h=3 h=4

3

Phase (rad)

2 1 0 −1 −2 −3 −4 500

550

600

650

700

Frequency (rad/s) Fig. 7. Phase of the FRF in terms of HEM with h ¼3 and h ¼4.

of eigenpairs and system matrices are derived to obtain the FRF matrix. By using the two power-series expansions, the relationships satisfied by eigensolutions and system matrices are established, and an explicit expression on the contribution of the higher-order modes to FRF matrix can be expressed as a sum of the lower-order modes and system matrices. Then, a HEM is presented by expressing FRF matrix as the explicit expression of the contribution of the higher-order modes and the modal superposition of the lower-order modes. The HEM maintains original-space without having to use the state-space equation of motion such that it is efficient in computational effort and storage capacity. The HEM has the convergence region at 0–|λL+1| rad/s. It is implied that, although not all of the modes are required in the HEM, it needs to retain, at a minimum, all the modes whose resonant frequencies lie within the range of exciting frequencies. The HEM also requires that the system matrices are available. These pose some serious limitations in the experimental and computational calculation of FRFs. Finally, a two-stage floating raft isolation system is used to illustrate the effectiveness of the derived results. It is shown the errors of MDM may be significant if only a few lower modes are retained and the MDM may be not a good approach to perform the frequency response analysis, especially for the FRF around the antiresonance frequencies. The presented HEM can make the errors reduce. When the number h considered in HEM is increased, the errors of FRFs can be decreased, but the calculations of HEM will suffer from more computational effort. As can be seen from the two-stage floating raft isolation system, the approximated results calculated in terms of HEM with h ¼4 can give almost the same accuracy with the exact results computed by DFRM, but save much computational cost.

Acknowledgments This research was supported by the National Science and Technology Major Project of China (Grant no. 2012ZX04003021), the National Natural Science Foundation of China (Grant nos. 30870605 and 50675077), and the Fundamental Research Funds for Universities supported by central authority (Grant no. 2011QN126). Appendix A. Proof of the iterative procedure The Appendix will prove the following iteration problem using the principle of mathematical induction. 2 3 6 7 1 6 7 IN þ ðsÞK1 ðC þ sMÞ þ ⋯ þ ðsÞk 6K1 ðC þ sMÞ⋯K1 ðC þ sMÞ 7 þ ⋯ ¼ ∑ ðsÞk Ξk 4|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}5 k¼0 k terms 1

1

1

where Ξ0 ¼ IN , Ξ1 ¼ K C and Ξk ¼ K CΞk1 K MΞk2 . By equating the coefficients of s2 and s3, we can check the first two iterative terms Ξ2 ¼ K1 CK1 CK1 M ¼ K1 CΞ1 K1 MΞ0 Ξ3 ¼ K1 CðK1 CK1 CK1 MÞK1 MK1 C ¼ K1 CΞ2 K1 MΞ1 The iterative problem for k43 will be proved in terms of the following two distinct cases.

ðA1Þ

40

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

Case 1. k is the even number The series-expansion terms involving the coefficient matrix of sk  2 and sk  1 can be found as 2 3 2 3 ðsÞ

k2 2

6 1 1 4K ðC þ sMÞ⋯K ðC þ sMÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} 2

ðsÞ

k1 4

ðk2Þ=2 terms

k 6 1 7 1 5 þ ðsÞ2 4K ðC þ sMÞ⋯K ðC þ sMÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

k=2 terms

3

1

5 ¼ ðsÞ

1

K ðC þ sMÞ⋯K ðC þ sMÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

7 5þ⋯

k2

k1

Ξk2 þ ðsÞ

ðA2Þ

Ξk1 þ ⋯

k1 terms

The series-expansion terms involving the coefficient matrix of sk are 2 3 2 k 6 1 7 1 1 1 k ðsÞ2 4K ðC þ sMÞ⋯K ðC þ sMÞ 5 þ ⋯ þ ðsÞ 4K ðC þ sMÞ⋯K ðC þ sMÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} k=2 terms

3 5 k terms

k

¼ ðsÞ Ξk þ ⋯

ðA3Þ

Comparing Eq. (A2) with Eq. (A3), one obtains ðsÞK1 ðC þ sMÞ½ðsÞk2 Ξk2 þ ðsÞk1 Ξk1 þ ⋯ ¼ ðsÞk Ξk þ ⋯

ðA4Þ

which can be further simplified as ðsÞk1 K1 CΞk2 ðsÞk K1 MΞk2 þ ðsÞk K1 CΞk1 ðsÞkþ1 K1 MΞk1 þ ⋯ ¼ ðsÞk Ξk þ ⋯

ðA5Þ k

By equating the coefficients of s , one has Ξk ¼ K

1

CΞk1 K

1

MΞk2 .

Case 2. k is the odd number We find the series-expansion terms involving the coefficient matrix of sk  2 and sk  1 2 3 2 ðsÞ

k1 2

6 1 1 4K ðC þ sMÞ⋯K ðC þ sMÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} 2

k1 4

7 5 þ ðsÞ ðk1Þ=2terms

6 1 1 4K ðC þ sMÞ⋯K ðC þ sMÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

3

1

5 ¼ ðsÞ

1

K ðC þ sMÞ⋯K ðC þ sMÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

þðsÞ

kþ1 2

k2

Ξk2 þ ðsÞ

k1

3 7 5þ⋯ ðkþ1Þ=2terms

ðA6Þ

Ξk1 þ ⋯

k1terms

In this case, the series-expansion terms involving the coefficient matrix of sk, i.e., Ξk as 2 3 2 kþ1 6 1 7 1 1 1 k ðsÞ 2 4K ðC þ sMÞ⋯K ðC þ sMÞ 5 þ ⋯ þ ðsÞ 4K ðC þ sMÞ⋯K ðC þ sMÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} ðkþ1Þ=2terms

3 5 kterms

ðA7Þ

k

¼ ðsÞ Ξk þ ⋯ By comparing Eq. (A6) with Eq. (A7), one obtains the iterative procedure Ξk ¼ K1 CΞk1 K1 MΞk2 when k is the odd number. So far, the proof is complete. References [1] K.K. Choi, N.H. Kim, Structural sensitivity analysis and optimization 1: linear systems, in: F.L. Frederic (Ed.), Mechanical Engineering Series, SpringerVerlag, New York, 2005. [2] M.D. Nastran, Dynamic Analysis User's Guide, MSC Software Corporation, 2010. [3] J.E. Mottershead, L. Michael, M.I. Friswell, The sensitivity method in finite element model updating: a tutorial, Mechanical Systems and Signal Processing 25 (2011) 2275–2296. [4] J.E. Mottershead, M.I. Friswell, Model updating in structural dynamics: a survey, Journal of Sound and Vibration 167 (1993) 347–375. [5] Z. Wang, R. Lin, M. Lim, Structural damage detection using measured FRF data, Computer Methods in Applied Mechanics and Engineering 147 (1997) 187–197. [6] Q. Huang, Y. Xu, J. Li, Z. Su, H. Liu, Structural damage detection of controlled building structures using frequency response functions, Journal of Sound and Vibration 331 (2012) 3476–3492. [7] H. Ouyang, D. Richiedei, A. Trevisani, G. Zanardo, Eigenstructure assignment in undamped vibrating systems: A convex-constrained modification method based on receptances, Mechanical Systems and Signal Processing 27 (2012) 397–409. [8] K. Sanliturk, O. Cakar, Noise elimination from measured frequency response functions, Mechanical Systems and Signal Processing 19 (2005) 615–631. [9] S. Adhikari, J. Woodhouse, Identification of damping: Part 1, viscous damping, Journal of Sound and Vibration 243 (2001) 43–61. [10] R. Lin, Identification of modal parameters of unmeasured modes using multiple FRF modal analysis method, Mechanical Systems and Signal Processing 25 (2011) 151–162. [11] S. Pradhan, S. Modak, A method for damping matrix identification using frequency response data, Mechanical Systems and Signal Processing 33 (2012) 69–82.

L. Li et al. / Mechanical Systems and Signal Processing 42 (2014) 31–41

41

[12] D.-H. Lee, W.-S. Hwang, C.-M. Kim, Design sensitivity analysis and optimization of an engine mount system using an FRF-based substructuring method, Journal of Sound and Vibration 255 (2002) 383–397. [13] Y.-H. Park, Y.-S. Park, Structure optimization to enhance its natural frequencies based on measured frequency response functions, Journal of Sound and Vibration 229 (2000) 1235–1255. [14] J.M. Dickens, J.M. Nakagawa, M.J. Wittbrodt, A critique of mode acceleration and modal augmentation methods for modal response analysis, Computers and Structures 62 (1997) 985–998. [15] P. Leger, E.L. Wilson, Modal summation methods for structural dynamic computations, Earthquake Engineering and Structural Dynamics 16 (1988) 23–27. [16] S.H. Kulkarni, S.F. Ng, Inclusion of higher modes in the analysis of non-classically damped systems, Earthquake Engineering and Structural Dynamics 21 (1992) 543–549. [17] N.R. Maddox, On the number of modes necessary for the accurate response and resulting forces in dynamic analyses, Journal of Applied MechanicsASME 42 (1975) 516–517. [18] O.E. Hansteen, K. Bell, On the accuracy of mode superposition analysis in structural dynamics, Earthquake Engineering and Structural Dynamics 7 (1979) 405–411. [19] G. Borino, G. Muscolino, Mode-superposition methods in dynamic analysis of classically and non-classically damped linear systems, Earthquake Engineering and Structural Dynamics 14 (1986) 705–717. [20] A. D’Aveni, G. Muscolino, Improved dynamic correction method in seismic analysis of both classically and non-classically damped structures, Earthquake Engineering and Structural Dynamics 30 (2001) 501–517. [21] M. Di Paola, G. Failla, A correction method for dynamic analysis of linear systems, Computers and Structures 82 (2004) 1217–1226. [22] C.J. Camarda, R.T. Haftka, M.F. Riley, An evaluation of higher-order modal methods for calculating transient structural response, Computers and Structures 27 (1987) 89–101. [23] M.A. Akgun, A new family of mode-superposition methods for response calculations, Journal of Sound and Vibration 167 (1993) 289–302. [24] S. Benfratello, G. Muscolino, Mode-superposition correction method for deterministic and stochastic analysis of structural systems, Computers and Structures 79 (2001) 2471–2480. [25] P. Cacciola, N. Maugeri, G. Muscolino, A modal correction method for non-stationary random vibrations of linear systems, Probabilistic Engineering Mechanics 22 (2007) 170–180. [26] A. Palmeri, M. Lombardo, A new modal correction method for linear structures subjected to deterministic and random loadings, Computers and Structures 89 (2011) 844–854. [27] C. Bilello, M. Di Paola, S. Salamone, A correction method for dynamic analysis of linear continuous systems, Computers and Structures 83 (2005) 662–670. [28] C. Bilello, M. Di Paola, S. Salamone, A correction method for the analysis of continuous linear one-dimensional systems under moving loads, Journal of Sound and Vibration 315 (2008) 226–238. [29] Z.S. Liu, S.H. Chen, W.T. Liu, Y.Q. Zhao, An accurate modal method for computing responses to harmonic excitation, International Journal of Analytical and Experimental Modal Analysis 9 (1994) 1–13. [30] Z.S. Liu, K. Wang, C. Huang, S.H. Chen, An extended hybrid method for contribution due to truncated lower- and higher-frequency modes in modal summation, Engineering Structures 18 (1996) 558–563. [31] C. Huang, Z.S. Liu, S.H. Chen, An accurate modal method for computing response to periodic excitation, Computers and Structures 63 (1997) 625–631. [32] Z.-Q. Qu, Z.-F. Fu, A highly power series expansion method for calculating frequency response function, Chinese Journal of Computational Mechanics 17 (1998) 91–95. [33] W. Ning, J. Wang, J.-H. Zhang, An improved series expansion method of frequency response function under medium and high frequency excitations, International Journal of Mechanics and Materials in Design 6 (2010) 11–15. [34] S.H. Chen, Matrix Perturbation Theory in Structural Dynamic Design, Science Press, Beijing, 2007. [35] Z.S. Liu, S.H. Chen, Contrubution of the truncated modes to eigenvector derivatives, AIAA Journal 32 (1994) 1551–1553. [36] Q.H. Zeng, Highly accurate modal method for calculating eigenvector derivatives in viscous damping systems, AIAA Journal 33 (1995) 746–751. [37] Z.-Q. Qu, Accurate methods for frequency responses and their sensitivities of proportionally damped systems, Computers and Structures 79 (2001) 87–96. [38] Z.-Q. Qu, Hybrid expansion method for frequency responses and their sensitivities, Part I: undamped systems, Journal of Sound and Vibration 231 (2000) 175–193. [39] Z.-Q. Qu, R.P. Selvam, Hybrid expansion method for frequency responses and their sensitivities, Part II: viscously damped systems, Journal of Sound and Vibration 238 (2000) 369–388. [40] Z.-Q. Qu, Adaptive mode superposition and acceleration technique with application to frequency response function and its sensitivity, Mechanical System and Signal Processing 21 (2007) 40–57. [41] L. Li, Y.J. Hu, X.L. Wang, Improved approximate methods for calculating frequency response function matrix and response of MDOF systems with viscoelastic hereditary terms, Journal of Sound and Vibration 332 (2013) 3945–3956. [42] L. Rayleigh, Theory of Sound (two volumes), Dover Publications, New York, 1877. [43] T.K. Caughey, M.E.J. O’Kelly, Classical normal modes in damped linear dynamic systems, Journal of Applied Mechanics-ASME 32 (1965) 583–588. [44] S. Adhikari, Damping modelling using generalized proportional damping, Journal of Sound and Vibration 293 (2006) 156–170. [45] F. Ma, M. Morzfeld, A. Imam, The decoupling of damped linear systems in free or forced vibration, Journal of Sound and Vibration 329 (2010) 3182–3202. [46] F. Tisseur, K. Meerbergen, The quadratic eigenvalue problem, SIAM Review 43 (2001) 235–286. [47] S. Adhikari, Modal analysis of linear asymmetric nonconservative systems, Journal of Engineering Mechanics 125 (1999) 1372–1379. [48] D.E. Newland, Mechanical vibration analysis and computation, Courier Dover Publications, New York, 2006. [49] P. Lancaster, Lambda-matrices and vibrating systems, Dover Publications, New York, 2002. [50] J.H. Wilkinson, The algebraic eigenvalue problem, Oxford Univ Press, New York, 1965. [51] I. Fawzy, R.E.D. Bishop, On the dynamics of linear non-conservative systems, Proceedings of the Royal Society of London, Series-A 352 (1976) 25–40. [52] S. Adhikari, Eigenrelations for nonviscously damped systems, AIAA Journal 39 (2001) 1624–1630. [53] S. Adhikari, Damping Models for Structural Vibration, University of Cambridge, Cambridge, 2000. [54] L. Li, Y.J. Hu, X.L. Wang, L. Ling, Eigensensitivity analysis of damped systems with distinct and repeated eigenvalues, Finite Elements in Analysis and Design 72 (2013) 21–34.