Comparative study of model updating methods using frequency response function data

2.1. The sensitivity method

Starting from the motion equations without considering damping in the frequency domain:

Assuming:

Eq. (1) can be transformed to:

where $\mathbf{x}\left(\omega \right)$ and $\mathbf{f}\left(\omega \right)$ are the displacement of the structure and the input force in frequency domain, respectively. $\mathbf{Z}\left(\omega \right)$ is the dynamic stiffness matrix in frequency domain and the inverse of which is the frequency response function $\mathbf{H}\left(\omega \right)$, the relationship between $\mathbf{Z}\left(\omega \right)$ and $\mathrm{}\mathbf{H}\left(\omega \right)$ gives:

in which:

The displacement vector can be expressed by:

When a unit force is applied on the $j$th degree-of-freedom (dof) of the structure, on the selected $m$ dofs, Eq. (6) can be simplified as:

where ${\mathbf{H}}_{mj}$ is the element at the $m$th row and the $j$th column in the frequency response matrix, the vector ${\mathbf{H}}_j$ denotes the $j$th column of the frequency response matrix.

The error is defined as the difference between the measured and computed response given by:

where $\mathbf{p}$ is the parameter vector. ${\mathbf{x}}^E$ and ${\mathbf{x}}^A$ is respectively the experimental and the analytical response. Substituting Eq. (6) and (7) into Eq. (8), the residual term gives:

The superscript $E$ and $A$ denote the FE model and the experimental model respectively. The sensitivity method is developed from a Taylor series expansion in which the high-order terms are truncated, after the truncation, the $j$th column of frequency response function matrix ${\mathbf{H}}_j^A\left(\omega ,\mathbf{p}\right)$ with respect to ${\mathbf{p}}^0$ gives:

where $\mathbf{S}$ is the sensitivity matrix that can be expressed as:

$\mathrm{\Delta }\mathbf{p}$ and $N_p$ are the vector of correction and the number of the parameters need to be updated:

Substituting Eq. (10) into Eq. (9), the residual Eq. (9) can be expressed in the following matrix form:

Assume:

where $\mathbf{b}={\left(b_1,b_2,\dots ,b_m\right)}^T$ is the error vector between the experimental and analytical response. After $n$th iterations, if the correlation between the computational and the experimental response is close to 1, the norm of the residual $‖\mathbf{R}\left(\mathbf{p}\right)‖$ will converged to zero, and the parameter vector ${\mathbf{p}}^0$ is updated to be ${\mathbf{p}}^n$. Then:

Model updating can be transformed to an iteration problem; the following equations should be solved exactly in every iteration step:

It is noted that the sensitivity matrix is varying during the iteration procedure, and the sensitivity matrix must be calculated at each step until the norm of $\mathrm{\Delta }\mathbf{p}$ converged to a small quantity.

The sensitivity method is based on the truncated Taylor expansion, and the solution of Eq. (16) can only identify the error direction of the parameters in a local constraint; however, due to the truncated errors, when the initial parameters contains large errors, the finite element model may not be updated with this method.

2.2. The system matrix error method

Derivation of the system matrix error method begins with the measured frequency response functions, which is briefly introduced in this section.

According to Eq. (1), when a unit force is applied as excitation, the Equations of motion of the experimental model and the analytical model of the system can respectively be written as:

The superscripts $A$ and $E$ denote FE model and experimental model, respectively.

Assuming that:

where $\mathrm{\Delta }\mathbf{K}$, $\mathrm{\Delta }\mathbf{M}$, $\mathrm{\Delta }\mathbf{H}$ are the stiffness, mass and frequency response function errors in the finite element model.

Substituting Eqs. (19) into (18), we have:

Rearranging Eq. (20), it gives:

Substituting Eq. (17) into (21), we have:

Eq. (23) can be employed to solve for the unknown modeling errors. However, before this can be achieved, some parameterization is needed. To parameterize the modeling errors, it is assumed that the error mass and stiffness matrices can be expressed as linear combinations of elemental mass and stiffness matrices, respectively:

where ${\left[{\mathbf{M}}^A\right]}_i^{ele}$ and ${\left[{\mathbf{K}}^A\right]}_i^{ele}$ are the $i$th element mass and stiffness matrices, respectively. $a_i$ and $b_i$ are the perturbation in design parameters corresponding to the $i\text{t}\text{h}$ element. Notation $\mathrm{\Sigma }$ represents the summation. Substituting Eq. (24) into Eq. (23), we can obtain the following equation:

Eq. (25) can be transformed into a set of linear algebraic equations in terms of unknown perturbations in design parameters which denoted as $a_i$ ($i=\text{1,}\text{2,…,}N$) and $b_i$ ($i=\text{1,}\text{2,…,}N$), respectively:

where ${\mathrm{s}}_i^a$, ${\mathrm{s}}_i^b$, $\mathbf{a}$ and $\mathbf{b}$ are:

Eq. (26) is formulated based on measured response function data with one measurement frequency selected. In practical vibration test, response function data are measured at many different measurement frequencies. When response function data at sufficient number $n$ of measurement frequencies are used, Eq. (26) can be transformed as a set of over-determined algebraic equations:

where:

The coefficient matrix $\mathbf{S}$ and the vector $\mathbf{q}$ can be calculated in advance which is formed using the analytical model and the measured response function data. Eq. (28) can be solved for $\mathbf{p}$ using linear least-square method and then the solution $\mathbf{p}$ is used to reconstruct the updated analytical model together with the original analytical model. However, when the test data is incomplete and contaminated with noise, Eq. (28) is needed to be transformed and constraints have to be added.

2.3. Frequency point selection

Frequency response function are very sensitivity to damping properties at resonance peaks, local modes influence is include and no modal analysis is required, so that it contains much more information than modal data, at each frequency point, the major system information are included. In the updating procedure, because of the huge number of frequency points, it is necessary to select frequency point for improving the convergence speed and the precision, the selected frequencies should ensure the rank of the coefficient matrix in Eqs. (16) or (28) is not less than the number of the parameters which need to be updated, and when the experimental FRF data is contaminated with measurement noises, frequency points selection seems more important.

There are two main principles of frequency are given in references [3, 12]:

1) Choose updating frequencies at the foot of experimental resonances peaks, and

2) Avoid frequencies between corresponding analytical and experimental resonances.

The first rule assures the coefficient matrix in Eqs. (16) or (28) to be well-posed. The second one avoids the problem of large residuals, as shown in Fig. 1, if frequencies between corresponding analytical and experimental resonances are selected, in the process of the analytical data approaching to the experimental frequency response function, the analytical resonances will definitely moving to the select frequency, and the amplitude at resonance is much bigger than the other positions, therefore, the right hand side (Residual 2) may become unbalance easily during the iterations, this will cause the difficulties in solving linear algebraic equations. At each iteration step, it is necessary to ensure the chosen frequencies comply with the rules.

Fig. 1Iteration problem if frequencies are selected between corresponding analytical and experimental resonances

2.4. Solution strategies

Finite element model updating can often be transformed into optimization problem, for accurately reflecting the practical situations and the physical significances of the updated model, and the solution of the optimum solution is often subject to a certain constraint. Due to the above reasons, based on the derivation of the updating theory, the linear algebraic Eq. (16) or (28) can be transformed to an optimization problem with constrains:

where $\mathbf{V}\mathbf{L}\mathbf{B}$ and $\mathbf{V}\mathbf{U}\mathbf{B}$ are respectively the lower limit and the upper limit of the parameter vector $\mathbf{x}$. Solving this problem by use of mathematical method to get the estimated vector $\mathrm{\Delta }\mathbf{p}$, the iteration procedure will continue until the norm of $\mathrm{\Delta }\mathbf{p}$ converged to a small quantity.

3.1. GARTEUR Truss

The simulation study using a GARTEUR Truss [14] is shown in Fig. 2, which is a cantilever truss. The numerical model of the Truss consists of 78 2-D beam elements, 74 nodes. Each node of the beam element has three dofs (two translations and one rotation) and hence, the total number of dofs in the FE model is 222. The material properties are used in the FE model: Young’s modulus is assumed to be $E=\text{0.75×}{\text{10}}^{\text{11}}\text{N/}{\text{m}}^{\text{2}}$ and density $\rho =\text{2800}\text{kg/}{\text{m}}^{\text{3}}$. The cross-sectional areas are $A_h=\text{0.004}{\text{m}}^{\text{2}}$ (horizontal), $A_v=0\text{.006}{\text{m}}^{\text{2}}$ (vertical), and $A_d=\text{0.003}{\text{m}}^{\text{2}}$ (diagonal), the bending moment is assumed to be the same for all the elements and is assumed to be $I=\text{0.0756}{\text{m}}^{\text{4}}$.

Fig. 2GARTEUR Truss model

Two simulated cases are carried out based on the GARTEUR Truss model. For Case-I, in order to generate the experimental data, stiffness modeling errors are introduced in the elements of the analytical model of the structure by changing the bending stiffness of elements as shown in Table 1, assuming that the stiffness of element shown in Table 1 are all overrated 20 % in the FE model. While in Case-II, the stiffness of elements in the FE model is overrated or underestimated have apparent discrepancies, details are shown in Table 2.

The two cases are employed to investigate the two methods introduced in Section 2, the updated results of coordinates incomplete and coordinates complete are discussed. Considering the small modeling error firstly, and then the large modeling errors.

3.2. Updated results evaluation criterions

To assess the progress of iterations and to compare different updated models certain model quality indices have been constructed. Percentage average error in natural frequencies (AENF) is calculated using following expressions:

where $f^A$ and $f^E$ are the natural frequencies, $m$ is the number of the considering modes. The updated model should have the abilities of reflect the test results in the considering frequency domain, outside of which the test results predicted.

Percentage average error in parameters (AEP) is calculated as error in the predicted parameters as a percentage of the known exact parameters:

where $np$ is the number of the parameters, $p_{exact}^i$ and $\mathrm{}{p}_{updated}^i$ are the value of the parameters.

3.3. Results and discussions

In order to demonstrate their practical application, the two model updating methods using frequency response function of GARTEUR Truss model which is introduced in Section 3.1. For comparative investigation on these methods, two cases including small modeling error and large modeling error are studied. In each case, the measured coordinates complete and incomplete are investigated, the complete case means the response of all coordinates of the structure can be tested, and the incomplete case meaning only translation degree of freedoms can be obtained.

3.3.1. Case-I small modeling error

For the sensitivity method, calculation of the frequency response function matrix ${\mathbf{H}}_j^A\left(\omega ,\mathbf{p}\right)$ with respect to the parameter vector $\mathbf{p}$ is very time-consuming, and a more important question is that the sensitivity matrix $\mathbf{S}$ will easily become ill-posed, therefore, the selection of parameters to be updated is also a very important aspect. In this case study, at first, we select the stiffness of all elements in the structure as parameters, the equation is an ill-posed problem, it is hard to get a credible solution. However, assume the position of the errors are known in advance, and select appropriate frequency points, iteration process is required during updating. The identification of element modeling errors is shown in Fig. 3. It can be seen that all the introduced modeling errors are well identified after convergence. The response function curves of the experimental, analytical and updated models are shown in Fig. 4.

Fig. 3Comparison of the exact and identified modeling errors (method 1, coordinates incomplete)

Fig. 4Comparison of the analytical, the experimental and the updated response function curves (method 1) (—o— updated, —— experimental, —*— analytical)

With regard to the second method, for the case of coordinates complete, it is need not to consider many limitations, the frequency point can be selected randomly, the identification of element modeling errors is shown in Fig. 5, and it is found that the agreement between the exact error and the identified error is excellent. But, in vibration test, it is not realistic that all coordinates which are specified in the analytical model have been measured. Therefore, the effect of coordinate incompleteness upon the updating procedure should be assessed. For the case of incomplete data, the model condensation and coordinates expansion technique should be used to solve this problem. In this work, those unmeasured elements of the response function vector are replaced by their analytical counterparts. The results for the identification of element modeling errors are shown in Fig. 6. The response function curves of the experimental, analytical and updated models are shown in Fig. 7. From this figure, there are a few differences between the regenerated and experimental response function data.

Fig. 5Comparison of the exact and identified modeling errors (method-2, coordinates complete)

Table 3 is the error indices comparison for the Case-I of small modeling error, from the table it is can be found that, under coordinates complete condition, method-2 can get an excellent solution.

Fig. 6Comparison of the exact and identified modeling errors (method-2, coordinates incomplete)

Fig. 7Comparison of the analytical, the experimental and the updated response function curves (method-2) (—o— updated, —— experimental, —*— analytical)

Table 3Error indices comparison for the Case-I of small modeling error (%)

Experimental coordinates	Percentage average error in natural frequencies (AENF)			Percentage average error in parameters (AEP)
	Before updating	Updated		Before updating	Updated
		Method-1	Method-2		Method-1	Method-2
Coordinates complete	2.82	–	0	20	–	0
Coordinates incomplete	2.82	0.01	0.18	20	0.08	0.59

3.3.2. Case-II large modeling error

For Case-II, the modeling error is much larger. It is hard to get a convergent solution by using of the sensitivity method; the reason is introduced in Section 2.1. But for the system matrix error method, in the case the coordinates is completely can be tested, it is also easy to get a reliable solution, the results can be seen in Fig. 8 and Fig. 9, it is shown that the solution is close to the exact error. However, for the case of incomplete coordinates, it is not easy get a convergent updating solution. But if we know the stiffness error position in advance, under this condition, the identification of stiffness error can be obtained. The comparison of the exact and identified modeling errors is shown in Fig. 10.

Fig. 8Comparison of the exact and identified modeling errors (method-2, coordinates complete)

Fig. 9Comparison of the analytical, the experimental and the updated response function curves (method-2) (—o— updated, —— experimental, —*— analytical)

Fig. 10Comparison of the exact and identified modeling errors (method-2, known error position in advance)