# Model Updating for Systems with General Proportional Damping Using Complex FRFs

Authors: Xinhai Wu; Huan He; Yang Liu; Guoping Chen
Journal: Shock and Vibration
 (2020)

---

## Abstract

In this paper, we propose a model updating method for systems with nonviscous proportional damping. In comparison to the traditional viscous damping model, the introduction of nonviscous damping will not only reduce the vibration of the system but also change the resonance frequencies. Therefore, most of the existing updating methods cannot be directly applied to systems with nonviscous damping. In many works, for simplicity, the Rayleigh damping model has been applied in the model updating procedure. However, the assumption of Rayleigh damping may result in large errors of damping for higher modes. To capture the variation of modal damping ratio with frequency in a more general way, the diagonal elements of the modal damping matrix and relaxation parameter are updated to characterize the damping energy dissipation of the structure by the proposed method. Spatial and modal incompleteness are both discussed for the updating procedure. Numerical simulations and experimental examples are adopted to validate the effectiveness of the proposed method. The results show that the systems with general proportional damping can be predicted more accurately by the proposed updating method.


---

## Body

## 1. Introduction

Damping is a significant factor which describes the energy dissipation from vibration. To predict the dynamic behavior accurately, the effects of damping cannot be overlooked, especially if the external forced frequency is in the vicinity of the natural frequency of the main structure. In order to describe different energy dissipation mechanisms, many damping models have been proposed by scholars. One of the representative damping models is the viscous damping, and this model has been widely used. However, as engineering systems become more and more sophisticated and the use of composite materials is much greater than before, the application of viscous damping may result in significant errors in dynamic analysis of light space structures. Thus, several nonviscous damping models have been proposed to illustrate the dissipative forces in a more general way, such as the Biot model [1, 2], the Golla–Hughes–McTavish (GHM) model [3, 4], the Anelastic Displacement Field (ADF) model [5], and the fractional derivative model [6, 7]. The nonviscous damping assumes that the damping force depends on any variable other than the instantaneous velocity. In this paper, the exponential damping model, which describes the damping force via convolution integrals over the exponential function and velocity, will be chosen as the nonviscous damping model to be discussed. Assuming that the kernel function is composed of several exponential functions, the equations of motion of N-DOF linear system can be expressed as follows:(1)Mx¨t+∫0t∑i=1Nμμie−μit−τCx˙τdτ+Kxt=ft,where M,K∈ℝN×N are the mass and stiffness matrices, respectively; the constant μ is the relaxation parameter; Nμ is the number of exponential kernel functions; and C∈ℝN×N is the exponential damping coefficient matrix.By transforming exponential damping into the frequency domain, the equations can be written as follows:(2)−ω2M+∑i=1Nμjωμijω+μiC+KXω=Fω,where j=−1 denotes the imaginary unit. It is noted that if the relaxation parameters go to infinity, the exponential damping model will degenerate into the well-known viscous damping model. In other words, the viscous damping model is a special case of the exponential damping model.The dynamic stiffness matrix (DSM)Z∈ℂN×N and frequency response function matrix H∈ℂN×N of such a system can be expressed as follows:(3)Z=−ω2M+∑i=1Nμjωμijω+μiC+K,(4)H=−ω2M+∑i=1Nμjωμijω+μiC+K−1.In order to accurately describe the vibration characteristic and help reduce the vibrations of engineering structures, it is essential to obtain reliable mathematical models. Through the unremitting efforts of many researchers, the finite element (FE) theory has been greatly developed, and now the FE method has become the most widely used method to obtain mathematical models. If the material and geometry parameters of the structure are obtained, the mass and stiffness matrices can be calculated rather accurately. However, until now, the understanding of damping is still very vague. The underlying damping mechanisms of structures are usually too complicated, and the damping matrix in the FE model cannot be obtained as accurate as the mass and stiffness matrices. With the assumption of the viscous or structural damping model, a lot of identification methods have been proposed to estimate the corresponding damping coefficient matrix. Usually, according to the data used in the identification process, the damping identification methods can be divided into two categories, namely, the matrix method [8–11] and the modal method [12–16]. For some methods, the analytical mass and stiffness matrices need to be obtained before the identification procedure. Many scholars have compared the performances of these identification methods from different aspects, for example, Phani and Woodhouse [17], Prandina et al. [18], and Pradhan and Modak [19]. Although the finite element method can provide a general estimation of the mass and stiffness matrices, usually, the simulation results computed by the FE model are not in good agreement with the corresponding experimental results. On the contrary, it is rather difficult to build accurate mathematical simulation models of boundary conditions, joints, and damping for large complex structures. Besides, the modelling errors of material and geometry properties can also lead to the differences between the analytical and experimental models. To minimize the differences, many scholars make great efforts to the update the FE models. The purpose of model updating is to obtain a more reliable analytical model of a structure and the corresponding model updating method can be regarded as the combination of the FE method and experiment method. Usually, model-updating methods can be classified into two categories, including the direct method and the iterative method. Based on the experimental modal date, the direct updating method [20] can minimize the error between analytical and experimental modal results by modifying mass and stiffness matrices. Recently, Bagha et al. [21] use the direct updating method to update the model of structures with composite materials. For the iterative method, both the modal data and frequency response functions can be utilized to update the analytical model. By the use of frequency response functions (FRFs), Lin and Ewins [22] firstly proposed the iterative method to improve the analytical model, and this method was referred as the response function method (RFM). Imregun et al. [23, 24] discussed several solution strategies for the RFM. Lu and Tu [25] proposed the two-level neural network method to update structural parameters and damping ratios. The mass and stiffness matrices are updated during the first step and the updating of the damping ratios follows after that, which means that this is the two-step method. Lin and Zhu [26] updated FE models by incorporating damping matrices. Based on the frequency response function, Arora et al. modified the model parameters in the complex field [27, 28]. Pradhan and Modak [29] proposed the model-updating method based on the normal frequency response functions (nFRFs). The nFRFs are deduced from experimental complex FRFs and can be utilized to update the damping matrix separately from mass and stiffness matrices. Gang et al. [30] considered the effect of model reduction and frequency shift and proposed a new iterative model updating method. The convergence of the updating method can be improved by using the so called “pseudomaster DOFs” as additional reference. Lin [31] proposed the improved model-updating method based on the function-weighted sensitivity method. Mondal and Chakraborty [32] identified nonproportional viscous damping matrix by matching the imaginary parts of complex mode shapes. In that method, the mass and stiffness have to be accurately obtained firstly. Many scholars have conducted extensive studies of different model updating methods, e.g., Arora et al. [33, 34], Jiang et al. [35], and Modak et al. [36].For most real-life structures, complex modes can always be obtained from experimental results. Thus, nonproportional damping model is more plausible. However, for complex structures, the spatial distribution of damping and the underlying damping mechanism are unclear. While updating the damping of structure, the location and underlying damping mechanism usually need to be assumed in priority, which may differ from the real damping. Currently, the proportional damping model is still preferred in most complex structures, and many researchers implement proportional damping to perform dynamic analysis from different aspects. With the assumption of the proportional damping model, the multidegree-of-freedom system can be completely decoupled; therefore, the dynamic analysis can be simplified significantly. Liu et al. [37] evaluated the damping nonproportionality by the identified modal information. Adhikari [38, 39] proposed a generalized proportional damping model represented by smooth continuous functions, in which variables are represented by mass and stiffness matrices of special arrangements. For the simplicity of dynamic analysis, Bilbao et al. [40] proposed two methods to build the analytical model with proportional damping for structures with additional viscoelastic dampers. Lázaro et al. [41] modelled the systems with viscoelastic damping by using the equivalent proportional viscous damping model. Li and Law [42] identified the damping ratios based on the sensitivity of acceleration response of structural with respect to the damping parameters, and the viscous damping matrix is obtained by inverse modal transformation. For the model updating of systems with proportional viscous damping, Lin and Zhu [26] proposed the model updating method for the Rayleigh and Cauchy damping model. However, it is noted that the mass and stiffness matrices should be firstly corrected before the identification of proportional damping matrix. Besides, for the Cauchy damping model, the exponential power of frequency exists in the calculation of sensitivity, and this may lead to ill-conditioned coefficient matrix. With the transformation of exponential damping into frequency domain, it can be found that both the damping and stiffness matrix depend on frequency. Accordingly, the identification and model updating of such models are more complex. As most existing identification and updating methods are based on the assumption of viscous damping model, some corresponding improved methods have to be proposed for these kinds of damping models. Adhikari and Woodhouse [43] identified the exponential damping by modal parameters based on the first-order perturbation method. Arora et al. [44] updated the nonproportional and nonviscous damping model by using the complex FRFs through modifying material damping parameters of different components. Pan and Wang [45] updated the exponentially damped systems with the assumption of the Rayleigh damping model; however, the damping ratios of higher modes are obviously high in that model, which may not be valid for practical engineering structures. Consequently, it is necessary to propose a model updating method for systems with more general proportional damping.As the derivation of modal data may lead to extra errors, especially for the imaginary parts of the mode-shapes, to overcome the drawback, this paper will focus on the model updating based on complex FRFs. The main contributions of this paper is to derive the sensitivity of the dynamic stiffness matrix with respect to the physical and damping parameters for systems with nonviscous proportional damping, which can been seen in Section2.2 and 2.3. To characterize the damping energy dissipation of the structure, in addition to the material and geometry parameters, the diagonal elements of the modal damping matrix and the relaxation parameters are chosen as the updating parameters. In this way, without the assumption of Rayleigh or Cauchy damping models, the arbitrary variation of modal damping ratios with respect to frequency can be captured within the interested frequency range. The proposed method is validated by numerical and experimental examples. The results show that the dynamic response can be predicted accurately by the updated analytical model.


## 2. Methodology

For the proportionally damped system, the damping matrix can be diagonalized as(5)C′=ΦTCΦ=diagc1′,c2′,…,cN′,where C′∈ℝN×N is the diagonal modal damping matrix, ck′ is the kth diagonal element, and Φ=φ1,φ2,…,φN is the normalized undamped mode-shape matrix of the analytical model. With the mass normalization condition, the inverse of Φ and ΦT can be calculated as follows:(6)Φ−T=MΦ,Φ−1=ΦTM.Once the modal damping matrix is updated using experimental results, the damping matrix in the physical coordinate can be obtained accurately by the inverse transformation:(7)C=Φ−TC′Φ−1=MΦC′ΦTM.For most real-life structures, the modes of high frequency are hard to be excited; hence, the modal information of the experimental data is often incomplete. With the assumption that only the firstNc elements of the modal damping matrix are considered in the updating procedure, then the damping matrix in the physical coordinate can be expressed as follows:(8)C^=Φ^−TC^′Φ^−1=MΦ^C^′Φ^TM,where(9)C^′=diagc1′,c2′,…,cNc′,Φ^=φ1,φ2,…,φNc.In this paper, the updating parameters are divided into two categories: the physical parameterpk (geometry, material, boundary condition, etc.) and the damping parameters (diagonal elements of modal damping matrix ck′ and the relaxation parameter μk). The detailed updating formulas are given in the following sections.


### 2.1. Mode-Shape Sensitivity with Respect to the Updating Parameter

Firstly, the mode-shape sensitivity with respect to the physical parameter is derived [46], which will be utilized later in Section 2.2. The characteristic equations of the undamped dynamic system are represented as follows:(10)−λiM+Kφi=0,i=1,2,…,N,where λi and φi are the corresponding ith eigenvalue and eigenvector. Taking the derivative of equation (10) with respect to the kth physical updating parameter pk, we obtain(11)−λiM+K∂φi∂pk+−λi∂M∂pk−∂λi∂pkM+∂K∂pkφi=0,k=1,2,…,Np,where Np is the number of physical parameters that are chosen to be updated.Generally, the derivative of the mode shape can be expressed by the linear combination of the mode-shape matrix:(12)∂φi∂pk=∑l=1Nβlikφl,where βlik is the coefficient to be solved. The superscript k indicates that the coefficient is related to the kth physical updating parameter. Substituting equation (12) into equation (11) yields:(13)−λiM+K∑l=1Nβlikφl+−λi∂M∂pk−∂λi∂pkM+∂K∂pkφi=0.Premultiplying equation (13) by φr and considering the orthogonality of the eigenvector with respect to the mass and stiffness matrices, one has(14)∂Φ∂pk=Φβk,∂ΦT∂pk=ΦβkT,where βk can be calculated by(15)βlik=φlT−λi∂M/∂pk+∂K/∂pkφiλi−λl,l≠i,−φlT∂M/∂pkφi2,l=i.Combining equation (6) and equation (14), the derivative of the inverse of mode-shape matrix with respect to the physical parameters can be calculated by(16)∂Φ−T∂pk=∂MΦ∂pk=∂M∂pkΦ+M∂Φ∂pk,(17)∂Φ−1∂pk=∂ΦTM∂pk=∂ΦT∂pkM+ΦT∂M∂pk.Based on the assumption of proportional damping, the mode shapes are the same in either damped or undamped systems. Regarding the mode-shape sensitivity with respect to the damping parameters, one can obtain the following relationship:(18)∂φi∂ck′=0,i,k=1,2,…,N;(19)∂φi∂μk=0,i=1,2,…,N;k=1,2,…,Nμ.


### 2.2. Dynamic Stiffness Sensitivity with Respect to the Physical Parameter

From equation (3), the sensitivity of dynamic stiffness with respect to the kth physical parameter pk can be expressed as follows:(20)∂Z∂pk=−ω2∂M∂pk+∑i=1Nμjωμijω+μi∂C∂pk+∂K∂pk.The mass and stiffness sensitivity in equation (20) can be calculated easily while constructing the global mass and stiffness matrix by the FE method. As the updating of physical parameters will lead to the changes of mode shapes, from equation (7), it can be found that the damping matrix in the physical coordinate is also dependent on the physical parameters. The damping matrix sensitivity with respect to physical parameters can be calculated as follows:(21)∂C∂pk=∂Φ−TC′Φ−1∂pk=∂Φ−T∂pkCΦ−1+Φ−T∂C′∂pkΦ−1+Φ−TC′∂Φ−1∂pk.Since the modal damping parameters are independent updating variables, regarding the sensitivity of the modal damping matrix with respect to the physical parameter, one can obtain the following relationship:(22)∂C′∂pk=0.Substituting equation (6), (16), (17), and (22) into equation (21), we can obtain(23)∂C∂pk=∂M∂pkΦ+M∂Φ∂pkC′ΦTM+MΦC′∂ΦT∂pkM+ΦT∂M∂pk.For modal incomplete case, using equation (8), equation (23) can be expressed as follows:(24)∂C^∂pk=∂M∂pkΦ^+M∂Φ^∂pkC^′Φ^TM+MΦ^C^′∂Φ^T∂pkM+Φ^T∂M∂pk.Substituting equation (24) into equation (20), we can obtain the dynamic stiffness sensitivity with respect to the physical parameter.


### 2.3. Dynamic Stiffness Sensitivity with Respect to the Damping Parameter

For the exponential damping model, the damping parameters that need to be updated include the diagonal elementsck′ of the modal damping matrix and the relaxation parameter μk.Since the mass and stiffness sensitivity with respect to the damping parameter are equal to zero, the dynamic stiffness sensitivity with respect to thekth element of the modal damping matrix can be calculated by(25)∂Z∂ck′=∑i=1Nμjωμijω+μi∂C∂ck′=∑i=1Nμjωμijω+μi∂Φ−TC′Φ−1∂ck′.For the proportional damping model, substituting equation (18) into equation (25), one has(26)∂Z∂ck′=∑i=1Nμjωμijω+μiΦ−T∂C′∂ck′Φ−1.Since the modal damping matrix is the diagonal matrix constructed byck′, one can obtain(27)∂C′∂ck′=Ik=diag0,…,1,…,0,where Ik∈ℝN×N is a zero matrix except for the unit element in the kth diagonal. Substituting equation (27) into equation (26), one has(28)∂Z∂ck′=∑i=1Nμjωμijω+μi∂C∂ck′=∑i=1Nμjωμijω+μiΦ−TIkΦ−1=∑i=1Nμjωμijω+μiMΦIkΦTM.For thekth relaxation parameter μk, the dynamic stiffness sensitivity can be computed easily by(29)∂Z∂μk=jωjω+μk−jωμijω+μk2C=jωjω+μk−jωμkjω+μk2MΦC′ΦTM.


### 2.4. FRF-Based Model-Updating Method for Proportionally Damped Systems

For the analytical and experimental model, the following relationship always holds:(30)ZAHA=I,(31)ZEHE=I,where the superscripts A and E denote the analytical and experimental model, respectively.Express the experimental dynamic stiffness matrix as(32)ZE=ZA+ΔZ.Substituting equations (32) and (30) into equation (31), one has(33)ΔZHE=ZAHA−HE.Premultiply equation (33) with the analytical FRF HA,(34)HAΔZHE=HA−HE.If only one column of the experimental FRF matrix is measured, then the above equation can be reduced as(35)HAΔZHjE=HjA−HjE,where j denotes the jth column of the experimental FRF.LinearizingΔZ with respect to the updating parameters αi, we have(36)ΔZ=∑i=1Nt∂Z∂αiΔαi,Nt=Np+Nc+Nμ,where α=p1,p2,…,pNp,c1′,c2′,…,cNc′,μ1,μ2,…,μNμT∈ℝN×1 and Nt is the total number of updating parameters.The dynamic stiffness sensitivities with respect to the updating parameters are all given in Sections2.2 and 2.3. Substituting equation (36) into equation (35), one obtains(37)∑i=1NtHA∂ZA∂αiHjEΔαi=HjA−HjE.If FRFs at several frequency points (Nf) are considered, it is easy to construct the following algebraic equations:(38)SΔα=ΔH,where(39)S=HAω1∂ZA∂α1HjEω1⋯HAω1∂ZA∂αNtHjEω1⋮⋱⋮HAωNf∂ZA∂α1HjEωNf⋯HAωNf∂ZA∂αNtHjEωNf,(40)ΔH=HjAω1−HjEω1⋮HjAωNf−HjEωNf.Since the updating parameters are all real variables, separating the real and imaginary parts of equation (38), we can obtain(41)ReSImSΔα=ReΔHImΔH.Once enough proper points are chosen within the interested frequency range, equation (41) can be solved by the SVD technique. As equation (36) is just a first-order approximation, the iterative method of the above process can be applied to update the parameters. Assuming that at the kth iteration, the updating variable is αk, then Δαk can be calculated with the updated sensitivity matrix Sk and residual vector ΔHk by using equation (39)–(41). The updated parameter can be obtained as(42)αk+1=αk+Δαk.The iteration procedure will not stop until the updating parameter satisfies either one of the following conditions:(43)Δαkαk≤ε,k≥Nstep,where ⋅ denotes the norm of a vector and ε and Nstep are the given convergence threshold and the maximum number of iterations, respectively.Usually, the incompleteness of experimental data is inevitable and some approximation has to be introduced during the updating procedure. For higher modes outside the interested frequency range, the corresponding modal damping is neglected when constructing the damping matrix in this paper. It should be pointed out that the number of the diagonal elements of the modal damping matrix which need to be updated should be larger than the number of modes considered in the interested frequency range, as the damping of each mode is described by the corresponding diagonal element. The proposed model-updating method can be regarded as an application of the mode superposition. To solve the problem of spatial incompleteness, the unmeasured FRFs ofHjE in equation (37) are replaced by the analytical counterparts of HjA [46].Since that the damping effect is more pronounced around the resonance peaks, to have a better estimation of the damping parameter, the FRFs around the resonance peaks are used in the updating procedure. Besides, to improve the robustness of the updating procedure, some weighting techniques proposed by Lin [47] are adopted to balance the FRFs at different frequencies and locations.


## 2.1. Mode-Shape Sensitivity with Respect to the Updating Parameter

Firstly, the mode-shape sensitivity with respect to the physical parameter is derived [46], which will be utilized later in Section 2.2. The characteristic equations of the undamped dynamic system are represented as follows:(10)−λiM+Kφi=0,i=1,2,…,N,where λi and φi are the corresponding ith eigenvalue and eigenvector. Taking the derivative of equation (10) with respect to the kth physical updating parameter pk, we obtain(11)−λiM+K∂φi∂pk+−λi∂M∂pk−∂λi∂pkM+∂K∂pkφi=0,k=1,2,…,Np,where Np is the number of physical parameters that are chosen to be updated.Generally, the derivative of the mode shape can be expressed by the linear combination of the mode-shape matrix:(12)∂φi∂pk=∑l=1Nβlikφl,where βlik is the coefficient to be solved. The superscript k indicates that the coefficient is related to the kth physical updating parameter. Substituting equation (12) into equation (11) yields:(13)−λiM+K∑l=1Nβlikφl+−λi∂M∂pk−∂λi∂pkM+∂K∂pkφi=0.Premultiplying equation (13) by φr and considering the orthogonality of the eigenvector with respect to the mass and stiffness matrices, one has(14)∂Φ∂pk=Φβk,∂ΦT∂pk=ΦβkT,where βk can be calculated by(15)βlik=φlT−λi∂M/∂pk+∂K/∂pkφiλi−λl,l≠i,−φlT∂M/∂pkφi2,l=i.Combining equation (6) and equation (14), the derivative of the inverse of mode-shape matrix with respect to the physical parameters can be calculated by(16)∂Φ−T∂pk=∂MΦ∂pk=∂M∂pkΦ+M∂Φ∂pk,(17)∂Φ−1∂pk=∂ΦTM∂pk=∂ΦT∂pkM+ΦT∂M∂pk.Based on the assumption of proportional damping, the mode shapes are the same in either damped or undamped systems. Regarding the mode-shape sensitivity with respect to the damping parameters, one can obtain the following relationship:(18)∂φi∂ck′=0,i,k=1,2,…,N;(19)∂φi∂μk=0,i=1,2,…,N;k=1,2,…,Nμ.


## 2.2. Dynamic Stiffness Sensitivity with Respect to the Physical Parameter

From equation (3), the sensitivity of dynamic stiffness with respect to the kth physical parameter pk can be expressed as follows:(20)∂Z∂pk=−ω2∂M∂pk+∑i=1Nμjωμijω+μi∂C∂pk+∂K∂pk.The mass and stiffness sensitivity in equation (20) can be calculated easily while constructing the global mass and stiffness matrix by the FE method. As the updating of physical parameters will lead to the changes of mode shapes, from equation (7), it can be found that the damping matrix in the physical coordinate is also dependent on the physical parameters. The damping matrix sensitivity with respect to physical parameters can be calculated as follows:(21)∂C∂pk=∂Φ−TC′Φ−1∂pk=∂Φ−T∂pkCΦ−1+Φ−T∂C′∂pkΦ−1+Φ−TC′∂Φ−1∂pk.Since the modal damping parameters are independent updating variables, regarding the sensitivity of the modal damping matrix with respect to the physical parameter, one can obtain the following relationship:(22)∂C′∂pk=0.Substituting equation (6), (16), (17), and (22) into equation (21), we can obtain(23)∂C∂pk=∂M∂pkΦ+M∂Φ∂pkC′ΦTM+MΦC′∂ΦT∂pkM+ΦT∂M∂pk.For modal incomplete case, using equation (8), equation (23) can be expressed as follows:(24)∂C^∂pk=∂M∂pkΦ^+M∂Φ^∂pkC^′Φ^TM+MΦ^C^′∂Φ^T∂pkM+Φ^T∂M∂pk.Substituting equation (24) into equation (20), we can obtain the dynamic stiffness sensitivity with respect to the physical parameter.


## 2.3. Dynamic Stiffness Sensitivity with Respect to the Damping Parameter

For the exponential damping model, the damping parameters that need to be updated include the diagonal elementsck′ of the modal damping matrix and the relaxation parameter μk.Since the mass and stiffness sensitivity with respect to the damping parameter are equal to zero, the dynamic stiffness sensitivity with respect to thekth element of the modal damping matrix can be calculated by(25)∂Z∂ck′=∑i=1Nμjωμijω+μi∂C∂ck′=∑i=1Nμjωμijω+μi∂Φ−TC′Φ−1∂ck′.For the proportional damping model, substituting equation (18) into equation (25), one has(26)∂Z∂ck′=∑i=1Nμjωμijω+μiΦ−T∂C′∂ck′Φ−1.Since the modal damping matrix is the diagonal matrix constructed byck′, one can obtain(27)∂C′∂ck′=Ik=diag0,…,1,…,0,where Ik∈ℝN×N is a zero matrix except for the unit element in the kth diagonal. Substituting equation (27) into equation (26), one has(28)∂Z∂ck′=∑i=1Nμjωμijω+μi∂C∂ck′=∑i=1Nμjωμijω+μiΦ−TIkΦ−1=∑i=1Nμjωμijω+μiMΦIkΦTM.For thekth relaxation parameter μk, the dynamic stiffness sensitivity can be computed easily by(29)∂Z∂μk=jωjω+μk−jωμijω+μk2C=jωjω+μk−jωμkjω+μk2MΦC′ΦTM.


## 2.4. FRF-Based Model-Updating Method for Proportionally Damped Systems

For the analytical and experimental model, the following relationship always holds:(30)ZAHA=I,(31)ZEHE=I,where the superscripts A and E denote the analytical and experimental model, respectively.Express the experimental dynamic stiffness matrix as(32)ZE=ZA+ΔZ.Substituting equations (32) and (30) into equation (31), one has(33)ΔZHE=ZAHA−HE.Premultiply equation (33) with the analytical FRF HA,(34)HAΔZHE=HA−HE.If only one column of the experimental FRF matrix is measured, then the above equation can be reduced as(35)HAΔZHjE=HjA−HjE,where j denotes the jth column of the experimental FRF.LinearizingΔZ with respect to the updating parameters αi, we have(36)ΔZ=∑i=1Nt∂Z∂αiΔαi,Nt=Np+Nc+Nμ,where α=p1,p2,…,pNp,c1′,c2′,…,cNc′,μ1,μ2,…,μNμT∈ℝN×1 and Nt is the total number of updating parameters.The dynamic stiffness sensitivities with respect to the updating parameters are all given in Sections2.2 and 2.3. Substituting equation (36) into equation (35), one obtains(37)∑i=1NtHA∂ZA∂αiHjEΔαi=HjA−HjE.If FRFs at several frequency points (Nf) are considered, it is easy to construct the following algebraic equations:(38)SΔα=ΔH,where(39)S=HAω1∂ZA∂α1HjEω1⋯HAω1∂ZA∂αNtHjEω1⋮⋱⋮HAωNf∂ZA∂α1HjEωNf⋯HAωNf∂ZA∂αNtHjEωNf,(40)ΔH=HjAω1−HjEω1⋮HjAωNf−HjEωNf.Since the updating parameters are all real variables, separating the real and imaginary parts of equation (38), we can obtain(41)ReSImSΔα=ReΔHImΔH.Once enough proper points are chosen within the interested frequency range, equation (41) can be solved by the SVD technique. As equation (36) is just a first-order approximation, the iterative method of the above process can be applied to update the parameters. Assuming that at the kth iteration, the updating variable is αk, then Δαk can be calculated with the updated sensitivity matrix Sk and residual vector ΔHk by using equation (39)–(41). The updated parameter can be obtained as(42)αk+1=αk+Δαk.The iteration procedure will not stop until the updating parameter satisfies either one of the following conditions:(43)Δαkαk≤ε,k≥Nstep,where ⋅ denotes the norm of a vector and ε and Nstep are the given convergence threshold and the maximum number of iterations, respectively.Usually, the incompleteness of experimental data is inevitable and some approximation has to be introduced during the updating procedure. For higher modes outside the interested frequency range, the corresponding modal damping is neglected when constructing the damping matrix in this paper. It should be pointed out that the number of the diagonal elements of the modal damping matrix which need to be updated should be larger than the number of modes considered in the interested frequency range, as the damping of each mode is described by the corresponding diagonal element. The proposed model-updating method can be regarded as an application of the mode superposition. To solve the problem of spatial incompleteness, the unmeasured FRFs ofHjE in equation (37) are replaced by the analytical counterparts of HjA [46].Since that the damping effect is more pronounced around the resonance peaks, to have a better estimation of the damping parameter, the FRFs around the resonance peaks are used in the updating procedure. Besides, to improve the robustness of the updating procedure, some weighting techniques proposed by Lin [47] are adopted to balance the FRFs at different frequencies and locations.


## 3. Numerical Examples

In this section, as shown in Figure1, the cantilever beam is considered to evaluate the effectiveness of the proposed method. The length of the beam is 1000mm and the cross-sectional area is 50×20mm. The elastic modulus and density are E=7.2×1010N/m2 and ρ=2780kg/m3, respectively. The FE model of the Euler–Bernoulli beam consists of 10 elements, and each node has two degrees of freedom (translational and rotational). Consequently, the size of the constructed mass and stiffness matrices is 20×20.Figure 1
Geometry and FE model of the cantilever beam.Proportional exponential damping model is utilized to model the damping of the beam. For simplicity, only one exponential kernel function is considered here. The degree of nonviscosity can be described by comparing the relaxation parameter with the maximum mode frequency within the interested frequency range. While the relaxation parameter is much larger than the interested frequency, the nonviscosity is so weak that the exponential damping model can be regarded as the traditional viscous damping model. In this example, the damping coefficient matrix is assumed to have the following form as(44)C=2ξMM−1K,where ξ is the constant variable, which denotes the damping level of the system. It is obvious that the given damping is proportional in the system. For viscous damping model, the damping matrix obtained from equation (44) ensures that the damping ratios of all the modes are the same. Considering different levels of nonviscosity, two cases listed in Table 1 are discussed in this section. To simulate experimental FRFs, the frequency response functions are firstly transformed into time domain by inverse Fourier Transform (IFFT). Then, random noise is added to the calculated time functions to simulate the experimental impulse response functions. Finally, the impulse responses are transformed back to frequency domain to simulate the FRFs contaminated with noise. The noise level are determined by the signal-noise ratio (SNR) based on the RMS value. In this example, noise is added to make SNR = 50 dB. Considering the spatial incompleteness in practical situations, only the FRFs of translational degrees of freedom are assumed to be predetermined in this example. That is to say, the percentage of spatial incompleteness is 50%.Table 1
Values of damping parameters for different levels of nonviscous effect.
Caseμ (s−1)ξ15600.2525.6 × 1040.025To simulate the errors between the real structure and the corresponding FE model, as shown in Table2, some discrepancies are added in the thicknesses of each beam element. Before the updating procedure, the initial values of the damping should be provided first. In this example, the initial relaxation parameter is set as the first natural frequency. That is to say, the nonviscous effect is assumed to be strong at first. With the modal damping ratios estimated from the half-power bandwidth, the initial diagonal elements of the modal damping matrix can be calculated easily by using the mass and stiffness matrices. In this section, for the comparison with the proposed method, we also update the damping coefficient matrix which is expressed by the Rayleigh damping model αMM+βKK. Usually, for the Rayleigh damping model, the damping coefficients can be estimated using just 2 modes. In this paper, “1 and 2” denotes that the coefficients for the Rayleigh model are updated by the FRFs around the 1st and 2nd modes. “Full modes” denote that the coefficients are obtained by the FRFs around all the modes within the interested frequency range. Figure 2 presents the comparison of FRFs H1,1 between the simulated and initial analytical model for Case 1. It can be observed that the results do not match with each other at all.Table 2
Discrepancies between the finite element and the simulated experimental model.
Element no.1, 7, 9, 102, 63, 548Deviation in thickness (%)0+20+30+50−10Figure 2
Comparison of simulated and initial analytical FRFs H1,1 for case 1.The FRFs within the frequency range from 0 to 1600 Hz are assumed to be measured to simulate experimental data. Here, 10 thicknesses of the beam elements are chosen as the physical parameters to update the FE model. For the proposed updating method, as only 6 modes are contained in the interested frequency range, then the first 6 diagonal elements of the modal damping matrix and the relaxation parameter are selected as the damping parameters that need to be updated. For the Rayleigh damping model, the damping parameters that need to be updated areαM,βK,andμ.First, for Case 1, the relaxation parameter is smaller than the maximum natural frequency, which means the nonviscous effect is strong. Figure3(a) presents the convergence of the thicknesses of each beam element by the proposed method. Figure 3(b) presents the convergence of the thicknesses of each beam element, while the Rayleigh damping model is updated by the FRFs around the 1st and 3rd modes. It can be observed that the thicknesses can be updated accurately by the both two methods. Tables 3–5 present the comparisons of different parameters before and after updating (letter “N” denotes no convergence is reached). As demonstrated in Table 3, it can be found that both the proposed and Rayleigh damping model can give good estimation of the relaxation parameter. The max error of the updated relaxation parameter is less than 0.5%. As shown in Table 5, for the Rayleigh damping model, it can be found that different coefficients will be obtained if different FRFs are chosen during the updating procedure. However, if full modes within the interested frequency range are considered to update the Rayleigh damping model, no convergence is reached for the updating parameter. It is concluded that no results will be obtained if the wrong damping model is chosen for the updating procedure. This is also one of the motivations for the proposed method.Figure 3
Convergence of the thicknesses of each element for case 1. (a) Proposed updating method; (b) updating with the Rayleigh damping model using the FRFs around the 1st and 3rd modes.
(a)(b)Table 3
Comparison of the relaxation parameterμ before and after updating for case 1.
Initial (s−1)Accurate (s−1)Updated (s−1)Error (%)Proposed16560562.690.48Rayleigh(1 and 2)560.680.12Rayleigh(1 and 3)558.81−0.21Full modesN—Table 4
Comparison of the diagonal elements of the modal damping matrix before and after updating for case 1.
ModeDamping ratioInitial (Ns/m)Updated (Ns/m)Accurate (Ns/m)Error (%)10.2960.5951.6951.650.0820.084109.20322.87323.72−0.2631.79 × 10–265.05897.10906.63−1.0545.45 × 10–338.761728.961777.87−2.7552.04 × 10–323.982863.022943.55−2.7469.32 × 10–416.434218.504409.77−4.34Table 5
Comparison of the parameters of the Rayleigh damping model before and after updating for case 1.
VariablesInitialRayleigh (1 and 2)Rayleigh (1 and 3)Full modesαM59.7244.5648.86NβK1.19 × 10–46.64 × 10–42.62 × 10–4NFigure4 illustrates the comparison of the simulated experimental and updated analytical FRFs H1,1 for Case 1. Figure 5 presents the relative error between the amplitudes of the reconstructed FRF and the accurate one. It can be observed that the FRFs predicted by the proposed updating method are in good agreement with the simulated ones. The results outside the frequency range are not discussed here. However, it should be pointed out that, for those modes, as no damping has been considered, the FRFs would be larger than the exact ones. For the Rayleigh damping model, it can be found that only the modes considered in the updated procedure are in good accordance with the experimental ones, and the damping effects of the higher modes are larger than the real values. As presented in Figure 4(b), the resonance peaks for the first 4–6 modes are smaller than the simulated ones, while the resonance peak for the 2nd mode is larger than the simulated one. This is consistent with the variation of damping with respect to frequency for the Rayleigh damping model. Besides, the resonance peaks of higher modes shift toward right, and this is quite different from the viscous damping model. This may also be the reason why no convergence is reached while updating the Rayleigh damping by all the interested modes.Figure 4
Comparison of simulated and updated analytical FRFs H1,1 for case 1. (a) Rayleigh (1 and 2). (b) Rayleigh (1 and 3).
(a)(b)Figure 5
Relative error between the amplitudes of the reconstructed FRF and the accurate one for case 1.For Case 2, the relaxation parameter is much larger than the maximum natural frequency, which means that the nonviscous effect is very weak. Figure6 presents the updated thickness of each element, and it can be found that all the thicknesses converge to the real values. Tables 6 and 7 illustrate the comparisons of the initial and updated damping parameters. It can be observed that the diagonal elements of the modal damping matrix can be updated accurately by the proposed method. The maximum error of the diagonal elements for the first 6 modes is less than 5%. Besides, as shown in Table 7, the initial diagonal elements of the modal damping matrix which are estimated by the half-power bandwidth are very close to the real ones. Thus, it can be concluded that, for viscously damped system, the damping can be estimated by the half-power bandwidth even though small discrepancies exist in the mass and stiffness matrix. For the relaxation parameter, the initial value is also set as the first natural frequency. However, for Case 2, either updating the diagonal elements of the modal damping matrix or the coefficients of Rayleigh damping, the magnitude of the updated relaxation parameter always tends to be larger than 1010, which is much larger than the real value. As mentioned above, the nonviscous effect is too weak that the exponential damping model can be regarded as the viscous damping model. Figure 7 gives the comparison of the decoupled FRFs for the first six modes. “Equivalent viscous” means that the FRFs are obtained by the updated mass, stiffness, and damping matrices, and “Accurate” means that the FRFs are obtained by the accurate mass, stiffness, and damping model. It can be found that all the modal FRFs match well with the accurate ones within the interested frequency range. That is to say, for Case 2, even though the relaxation parameter is different from the exact value, the updated results can still capture the damping characteristics within the specific frequency range. As shown in Figure 8(b), the FRFs predicted by the proposed method are in good accordance with the simulated experimental ones. With the assumption of the Rayleigh damping model, only the peaks of the first two modes match with the simulated ones, and the resonance peaks of the higher modes are smaller than the experimental ones as expected.Figure 6
Updated thicknesses of each element for case 2.Table 6
Comparison of the parameters of the Rayleigh damping model before and after updating for case 2.
VariablesInitialRayleigh (1 and 2)αM0.834.45βK4.02 × 10–46.66 × 10–5Table 7
Comparison of the diagonal elements of the modal damping matrix before and after updating for case 2.
ModeDamping ratioInitial (Ns/m)Updated (Ns/m)Accurate (Ns/m)Error (%)10.02865.925.165.17−0.1920.024732.0332.3232.377−0.1830.024990.2990.9290.6670.2840.0251178.66178.26177.790.2650.0253297.77289.20294.36−1.7560.0249438.47424.97440.98−3.63Figure 7
Comparison of the modal FRFs for the first six modes.Figure 8
Comparison of simulated and analytical FRFs H1,1 for case 2. (a) Before updating and (b) after updating.
(a)(b)


## 4. Experimental Verification

In this section, for the validation of the proposed updating method, as shown in Figure9, the model of the steel cantilever beam is considered. The length of the beam is 1760 mm, the cross-sectional area is 45×8mm, and the density is ρ = 7784 kg/m3. The FE model of the Euler–Bernoulli beam is divided into 10 elements. Each node has two degrees of freedom (translational and rotational). To make a difference from the pure beam structure, a concentrated mass is fixed at the tip of the beam, and the mass is 0.224 kg. The accelerometers are attached at the 10 specific nodes to measure all the acceleration responses simultaneously. The accelerometers are too light that their effects on the dynamic response can be neglected. The impulsive excitation is repeatedly applied to the free end of the beam by the hammer for 4 times, and the FRFs are averaged for the further analysis.Figure 9
FE model (a) and experimental set up (b) for the cantilever beam.
(a)(b)One of the strategies to obtain an accurate model is to choose the updating parameters on the basis of engineering judgment. In the experimental apparatus, the end of the beam is clamped tightly enough that the constraint can be regarded as the ideal state. The values of the geometry parameters of the beam are accurately measured; thus, there is no need to update these parameters. The elastic modulus of the beam is estimated to beE=1.92×1011N/m2 by other experimental analysis. For the updating procedure, the initial values of the elastic modulus and the added mass are set as EA=0.8E and m=0.1kg, respectively. The elastic modulus and the added mass are chosen as the physical parameters that need to be updated. Before the updating procedure, the analytical undamped model can be established for the beam system by all the above given parameters. For this experimental beam model, the sources of the damping include junction, material damping, and viscous damping. However, the primary dissipative mechanism can be regarded as viscous on the whole. In this paper, the analytical model is updated based on the assumption of the proportional viscous damping model. Figure 10 presents the comparison of experimental and analytical unupdated FRFs H1,1 in the frequency range 0–120 Hz. It can be found that the analytical model do not match well with the real structure.Figure 10
Comparison of experimental and initial analytical FRFs H1,1.For the damping parameters, the initial values are calculated as discussed in Section3. Tables 8–11 present the comparisons of different parameters before and after updating, and it can be observed that the values of both the added mass and the elastic modulus can be updated rather accurately. The good agreement proves that the updated results are meaningful. As shown in Table 11, with the updated coefficients of the Rayleigh damping model, the diagonal elements of the modal damping matrix are also calculated. Figure 11 shows the overlay of experimental and updated analytical FRFs. Since all the physical parameters have been updated accurately, both the proposed damping model and the Rayleigh damping model can provide good estimation of the structure. However, for the responses around the 3rd and 4th resonance frequencies, the proposed method can give a better prediction as the modal damping has been captured in a more general way.Table 8
Values of the added mass before and after updating.
Initial (kg)Accurate (kg)Updated (kg)Error (%)Proposed0.10.2240.209−6.25Rayleigh0.219−11.16Table 9
Values of the elastic modulus (N/m2) before and after updating.
InitialExperimentalUpdatedError (%)Proposed1.536 × 10111.92 × 10111.887 × 1011−1.72Rayleigh1.879 × 1011−2.14Table 10
Comparison of the parameters of the Rayleigh damping model before and after updating.
VariablesInitialUpdatedαM0.0870.25βK2.394 × 10–58.74 × 10–6Table 11
Values of the diagonal elements of modal damping matrix for the first 5 modes before and after updating.
ModeProposedDamping ratioInitial (Ns/m)Updated (Ns/m)Rayleigh (Ns/m)10.00380.08510.11650.253220.00150.20830.38810.301830.00361.42421.52070.648240.00120.90661.34761.793150.00344.32074.54304.5164Figure 11
Comparison of experimental and updated analytical FRFs.
(a)(b)(c)(d)


## 5. Conclusions

In this paper, the FRF- based model- updating method is extended to update the system model with non-viscous proportional damping. The difficulties of updating non-viscously damped system lie in that the introduction of non-viscous damping will not only reduce the vibration, but also shift the resonance frequencies. Thus, the traditional updating methods which update the analytical model in two steps are inappropriate for such systems. Usually, the damping is non-proportional in the complex structure and the spatial distribution is unclear. For the sake of simplicity, the proportional damping is assumed for the updating procedure. Different from the updating of the coefficients of the Rayleigh damping model, the diagonal elements of the modal damping matrix and the relaxation parameter are updated to capture the damping characteristics. The number of the diagonal elements of the modal damping matrix contained in the updating procedure depends on the considered modes within the frequency range. To construct the updating formula, the dynamic stiffness matrix sensitivities with respect to the updating parameters are firstly derived. The identification of the damping matrix can be regarded as the application of the mode-superposition method. The accuracy of the proposed method is verified by the simulated and experimental examples. Spatial and modal incompleteness are also discussed for the updating procedure. The results show that, the proposed updating method can provide a better prediction for systems with general proportional damping. If the degree of non-viscosity is very small, the equivalent viscously damped system can be obtained.



---

*Source: 8886082-2020-11-20.xml*