A Copula-Based and Monte Carlo Sampling Approach for Structural Dynamics Model Updating with Interval Uncertainty

As the uncertainty is widely existent in the engineering structure, it is necessary to study the finite element (FE) modeling and updating in consideration of the uncertainty. A FE model updating approach in structural dynamics with interval uncertain parameters is proposed in this work. Firstly, the mathematical relationship between the updating parameters and the output interesting qualities is created based on the copula approach and the vast samples of inputs and outputs are obtained by the Monte Carlo (MC) sampling technology according to the copula model. Secondly, the samples of updating parameters are rechosen by combining the copula model and the experiment intervals of the interesting qualities. Next, 95% confidence intervals of updating parameters are calculated by the nonparameter kernel density estimation (KDE) approach, which is regarded as the intervals of updating parameters. Lastly, the proposed approach is validated in a two degree-of-freedommass-spring system, simple plates, and the transport mirror system. The updating results evidently demonstrate the feasibility and reliability of this approach.


Introduction
Finite element (FE) models that numerically solve various engineering problems can aid virtual prototyping, reduce product development cycle, and cut down the cost of performing the physical tests.However, the reliability of the simulation results by finite element modeling is not always guaranteed since FE models are the approximations of real world phenomena based on various assumptions.These assumptions may detract from the quality and accuracy of simulation results.In order to improve the accuracy of FE simulating to serve the structural design better, the FE model updating techniques are needed to develop.In the past few decades, various kinds of FE model updating approaches have been widely investigated based on the actually observed behaviors of the system.Additionally, experimental modal and vibration data are often used in FE model updating in the field of structural dynamics [1][2][3][4][5][6].
In most model updating approaches, the simulations are usually deterministic where each of the updating parameters is considered to have one "true" value and the purpose of the updating procedure is to provide a deterministic estimation.
In reality, there are always uncertainties in nominally identical structures, such as the structural parameter uncertainty (physical material properties, geometric parameters), the assembly joints uncertainty, and the experiment uncertainty (measurement noise, modal identification techniques, etc.).As a result, the FE model updating approaches with uncertainty have received great attentions recently.Studies have shown that the simulation results are more reliable when the uncertainties are taken into account [7], suggesting that it is necessary to consider the uncertainties during modeling and simulating [8].
In FE model updating approaches with uncertainty, the updated parameters are no longer deterministic and are described as random variables.Usually, the FE model updating approaches with uncertainty can be classified into two major categories: probabilistic and nonprobabilistic approaches.In the earlier works, a probabilistic approach proposed incorporated the measurement noise into model updating [9].Subsequently, Bayesian statistical frameworks were adopted to estimate the posterior probabilities of uncertain parameters [10][11][12].However, high computational costs due to a large amount of samples required for a satisfactory estimation greatly restrain the applications of Bayesian updating approaches.As a result, surrogate models such as the Gaussian process model with the perturbation approaches and sensitivity analysis approaches have been employed in stochastic model updating to improve the efficiency [13][14][15][16].Though, the surrogate model approaches own the superiority of computational efficiency over Monte Carlo (MC) based methods.Nevertheless, the prerequisite of small uncertainties, together with the Gaussian distribution assumption, also limits the applications to complex problems.Moreover, perturbation based predictions are sensitive to the initial estimates of parameters.Recently, an approach with the response surface models and MC simulation has been developed, which decomposed a stochastic updating process into a series of deterministic ones [17].On the other hand, the accuracy of the probabilistic approaches depends on the estimation of the probability distribution characteristics of the structural parameters and the responses.The establishment of an accurate probability distribution function (PDF) needs lots of experiment data in the probabilistic approaches, which greatly limits its application in engineering.
In nonprobabilistic approaches, the interval approach has been intensively investigated.By comparison, the experiment samples are not strictly needed in the FE model updating with interval analysis as was proposed.In the field of interval model updating (IMU), the inclusion theorem was employed to establish an interval inverse problem.And the convergence was achieved when measured responses fall into numerically predicted intervals [18][19][20][21][22]. Considering the easy implementation, IMU problems are usually solved within a deterministic framework where the upper and lower bounds of parameters are sought separately.For example, an IMU problem was decomposed into two deterministic constrained optimization processes where the midpoints and interval radii of parameters were separately estimated [19].Alternatively, the vertex solution theorem is effective and cost-efficient for IMU due to its easy implementation [20], particularly in the solution of Eigen value problems [21].But the vertex solution was valid only for particular parameterization of an FE model without the involvement of eigenvectors, which highly limits its further applications.Due to this drawback, global optimization algorithms were taken into account for more general solutions.Surrogate models such as the Kriging predictor and interval response surface were used to improve the efficiency of gradient computation and facilitate the convergence [23,24].So far most of IMU problems are solved within a deterministic framework since direct interval arithmetic operations are difficult to implement during inverse solutions.Therefore the upper and lower bounds of parameters should be sought separately through a deterministic inverse procedure.Additionally, global optimization of interval variables is difficult to realize due to the fact that the interval arithmetic is quite different with the traditional mathematical arithmetic.
Though several probabilistic and interval model updating approaches have been developed in the past years, most of them are still complicated for implementation.Additionally, these approaches with uncertainty suffer from the challenges such as ill-condition, nonuniqueness and local optimal solution, etc.To overcome such inconvenience, an IMU approach is developed in this work based on the copula model and MC sampling.In the proposed approach, the copula model between the updating parameters and the interesting qualities is constructed firstly.Then a large amount of samples is obtained according to the copula model, and the samples are rechosen based on experiment intervals of interesting qualities.Lastly, the updating intervals of parameters are obtained by estimating on the rechosen samples with kernel density estimation (KDE).The remainder of the paper is organized as follows.In Section 2, the copula-based FE model updating approach and procedure with interval uncertainty are presented.In Section 3, three examples are provided to validate the accuracy and reliability of the proposed approach.Conclusions are presented in Section 4.

Identification of Interval Parameters
The FE model updating problems are classic inverse problems in structural mechanics where the standard "forward" relationship between input and output variables of a model is inverted.The key in solving a FE model updating problem is to construct the mathematical relationship between the updating parameters and the output interesting qualities.The copula function is one of the most effective mathematical tools to determine this relationship, which expediently characterizes the correlation between the marginal functions of multivariables and the joint distribution function.

Brief Introduction of the Copula Function.
A copula function is a general function in statistics to formulate a multivariate distribution with various statistical dependence patterns, which was presented by Sklar in 1959 [25].Formally, a copula is a joint distribution function of standard uniform random variables.According to the Sklar's theorem, there exists a two-dimensional copula C such that variables  1 and  2 in a real random space. where Equation ( 1) can be spread for m-dimensional variables easily, that is, Consequently, the m-dimensional PDF is as follows: where ( 1 ( 1 ),  2 ( 2 ), . . .,   (  )) is the copula PDF, ( 1 ,  2 , . . .,   ) is the united PDF for m-dimensional random variables, and   (  ) is the PDF of the th random variable.
At present, the general copula function types include the Gaussian copula function, t-copula function, and Archimedean copula function [26,27].Among them, the Gaussian copula is widely utilized because most of the parameters in the engineering satisfy the normal distribution.In the study, the Gaussian copula is adopted for the FE model.Specifically, the Gaussian copula function is constructed by multidimensional Gaussian distribution and the linear correlation parameters, and its distribution function is as follows: where Φ  is the distribution function of the standard normal function for d-dimensional with the correlation matrix , Φ −1 is the inverse function of the distribution function of the standard normal function, and   =   (  ),  = 1, 2, . . ., .

Copula-Based Approach for Model Updating with Interval
Uncertainty.Firstly, the original design spaces of updating parameters   ,  = 1, 2, ⋅ ⋅ ⋅ , , are assumed, and a few samples are obtained by design of experiment (DOE) approach and subsequent deterministic FE analysis on samples according to DOE.Then, the samples of the output interesting qualities   ,  = 1, 2, ⋅ ⋅ ⋅ , , are obtained from the FE analysis results.Secondly, the copula model is constructed according to the samples of updating parameters and output interesting qualities, and resampling is performed to get large samples with number N for updating parameters and interesting response qualities based on the copula model.The samples falling into the experiment data space are considered to characterize the input-output relationship of the physical structure believably, and unuseful samples are needed to remove.Next, the samples of updating parameters   are rechosen according to the experiment intervals of interesting qualities, as follows: where  is the real space and  +  and  −  are the upper and lower bounds of the th output interesting quality which can be obtained from the experiment results.
In practical model updating, the measured data are only a few samples in general.Reasonable interval estimation on experiment data is the precondition to obtain the reliable updated FE model.However, the KDE allows for the capture of the observed distributional structure for the random variables, without having to assume a particular parametric distribution form.
Following [28], the kernel density estimator for variable x has the form where  is the number of observations used to construct the estimate, (⋅) is a kernel function,   is the th observation, and ℎ is the window width, or bandwidth.A typical choice for the kernel (⋅) is the standard normal density and is implemented here.The choice of the window width ℎ is usually based on the optimization of some scoring function.A least-square cross-validation score function is adopted for this work [28].
The empirical CDF and 95% confidence interval (CI) [ −  95  +  95 ]of the random variable  can be obtained by KDE in Matlab that is regarded as the interval of the random variable in this work.
Considering the fact that the estimation on the original intervals of the updating parameters may be inaccurate, the reliable intervals are not identified through one copula-based FE model updating procedure.In order to overcome this problem, the idea of the adaptive response surface technique is adopted for this work [29].That is, in order to get the final updating results, multiloop on the copula-based model updating procedures may be performed.
In the FE model updating procedure, the convergent criterion is that the difference of the updating parameter intervals between the ( + 1)th iteration step and the th iteration step is less than a small value, or the difference of the output interesting quality intervals between the th iteration step results and the experiment results is less than a critical value.
In order to improve the efficiency and the validity of model updating, the renewal strategy of updating parameters is as follows in each iteration step.The current intervals of updating parameters are updated according to the results of the previous iteration step, and the interval medians of the previous step are regarded as the current interval medians, and about 80% of the interval width of the previous step is regarded as the current interval width.Also, the Latin Hypercube Sample (LHS) method is suggested in the DOE, and the number of samples is not less than ten.
The copula-based model updating procedure is repeated until the convergent criterion is satisfied.The flow chart for the copula-based FE model updating is outlined in Figure 1.

Example 1: A Two-Degree-of-Freedom Mass-Spring System.
A two degree-of-freedom mass-spring system is shown in Figure 2. The deterministic parameters in the system are  For simplicity, it is assumed that the uncertain parameters are uniformly distributed.To create such kind of uncertainty, the LHS method is used to generate twenty experiment samples.Afterwards, the experiment results of the first two natural frequencies are obtained according   Assume that the original intervals of  2 and  3 are the same as [1.2 2.2] N/m, and the first two natural frequencies are regarded as the output interesting qualities.The intervals of  2 and  3 are identified according to the copula-based model updating flow with interval uncertainty in Figure 1.The updating results are convergent after three iteration steps.
Figure 3 is the scatter map between the updating parameters and the output interesting qualities when the updating results are convergent, which shows that there is strong correlation between  3 and  1 and between  2 and  2 and weak correlation between  2 and  1 and between  3 and  2 .
Then, the updated intervals of  2 and  3 are obtained by KDE for the resample in Figure 3, which are [0.7751.230] N/m and [0.882 1.133] N/m, respectively.The comparison between the original uncertain interval and the updated interval of updating parameters is shown in Figure 4, which shows that the updated interval matches the real interval better.Because the effect of small samples on experiment data is considered in the model updating procedure, the updated interval is bigger than the real interval and the result is reasonable.
In order to validate the updating results on the interval uncertain parameters, the copula models are reconstructed according to the original and the updated intervals of  2 and  3 .As a result, 5000 samples of  2 ,  3 , and the first two frequencies of the system are resampled by the copula models constructed just now.The scatter map between the simulating and the experimental results is shown in Figure 5, which indicates that the frequencies of the updated model are agreement with the experimental results better.

Example 2: Interval Model Updating in Simple Plates.
Impact hammer modal testing with free-free boundary conditions was conducted on thirty-three nominally identical steel plates in [30].The nominal geometric dimensions of the plates are 564 mm(length) × 110 mm(width) × 1.45 mm(thickness).And the nominal material properties are Young's modulus of 210 GPa, the shear modulus of 83 GPa, and the mass density of 7860 kg/m 3 .The statistical properties of the first five measured natural frequencies of the plates are given in Table 1.
The FE model of the plate is created by the SHELL181 element in ANSYS, which has 300 shell elements.The boundary is free-free in the model.The material properties of the original FE model are the nominal values.The FE model and the first five mode shapes of simulation are shown in Figure 6.
After investigation, the uncertainty of rectangular plates can be characterized by Young's modulus , the shear modulus , and the thick  in the FE model.The model updating  is to identify the interval of these three uncertain parameters , , and  by the experiment results.Firstly, the intervals of the first five modal frequencies are estimated by KDE according to the experiment data in Ref. [30].Secondly, it is assumed that the original intervals of E, G, and T are [195 220] GPa, [78 87] GPa, and [1.30 1.60] mm, respectively, and intervals of E, G, and T are identified according to the copula-based model updating flow with interval uncertainty in Figure 1.There are 20 experimental design data in each iteration step, i.e., 20 determined FE simulations in each iteration step.The model updating of the plate is convergent after four iteration steps with the convergent indices  1 = 0.002 and  2 = 0.002.
After updating, the estimated intervals of the three parameters were  = [203.21205.98]GPa,  = [82.7585.20] GPa, and  = [1.44 1.46] mm, respectively.The comparison between the original uncertain interval and the updated interval of the three parameters is shown in Figure 7, which shows that the updated intervals of these three parameters are much less than their original uncertain interval.
In order to validate the updating results on the interval uncertain parameters, the copula models are reconstructed on the original and the updated intervals of , , and .The 5000 samples of , , , and the first five natural frequencies of plates are obtained by MC sampling.The scatter plots for the simulation and the experimental results the first five natural frequencies of plates are shown in Figure 8.The comparisons on the natural frequencies between the simulation and the experimental results are listed in Table 2.It can be seen from Table 2 and Figure 8 that the frequencies of the updated model are in better agreement with the experimental results and the mean errors of frequencies decrease from the initial [12.51 12.62]% to [0.72 0.24]%.Thus, this experimental validation has also confirmed the feasibility of the proposed IMU approach in estimating interval parameters of steel plates.

Example 3: Interval Model Updating in Transport Mirror
System.The ShenGuangIII (SGIII) facility is designed for inertial confinement fusion (ICF) high energy experiments with 48 laser beams exactly transported and oriented to target.Figure 9 is the view of the beam transport system in SGIII facility target area.There are 276 transport mirror systems in the facility, and the dynamic response under ambient vibration is a key factor to affect the stability of the SGIII facility [31].A classic transport mirror system is shown in Figure 10, which consists of mirror component and support frame.
The material of the support frame is steel with nominal Young's modulus 200 Gpa, nominal Poisson's ratio 0.3, and the nominal density 7850 kg/m 3 .The material of the mirror is K9 class with nominal Young's modulus 80 Gpa, nominal Poisson's ratio 0.21, and the nominal density 2510 kg/m 3 .The transport mirror system is about 1.32 m in height.
The uncertainty exists in the transport mirror system because of the installation fluctuation and the welding technology fluctuation about the support frame, which induces the uncertainty of the vibration characteristics.In order to study how these uncertain factors affect the natural frequencies of the transport mirror system, impact hammer modal testing with fixed boundary conditions was conducted on ten nominally identical transport mirror systems.The experiment mount of the transport mirror system is shown in Figure 11.Three natural frequencies in the interesting frequency range were found to significantly influence the response of the transport mirror system under the work condition.Then, the first three natural frequencies of the transport mirror system should be updated before calculating the response.The corresponding natural frequencies were obtained by the modal experiments, as listed in Table 3.It is observed that the frequency variations become more obvious with the increase of the mode order.Generally, the bolts are ignored in the structural dynamic analysis model generally, and the FE model of the transport mirror system is established by the SOLID185 and SHELL181 element in ANSYS, which is shown in Figure 12.A fully fixed bottom of the support frame is used as the boundary condition.The material properties of the original FE model are the nominal values.The first three simulation mode shapes are shown in Figure 13.Theoretical analysis finds that the installation uncertainty of the transport mirror system can be characterized by Young's modulus  1 of the bottom part of the support frame, and the weld uncertainty of the support frame can be characterized by Young's modulus  2 of the upper part of the support frame.The goal of FE model updating is to identify the interval of  1 and  2 .
Firstly, the intervals of the first three modal frequencies are estimated by KDE according to the experiment data in Secondly, assuming that the original intervals of E 1 and E 2 are [100 150] GPa and [160 200] GPa, respectively, the intervals of E 1 and E 2 are identified according to the copula-based   In order to validate the updating results on the interval uncertain parameters, the copula models are reconstructed on the original and the updated intervals of E 1 and E 2 , Next, the 5000 samples of E 1 , E 2 , and the first three natural frequencies of the transport mirror system are obtained by MC sampling.The scatter plots for the simulation and the experimental results of the first three natural frequencies of  the transport mirror system are shown in Figure 15.The comparisons on the natural frequencies between the simulation and the experimental results are listed in Table 4.It can be seen from Table 4 and Figure 15 that the frequencies of the updated model agree with the experimental results better and the mean errors of frequencies decrease from the initial [2.4910.83]% to [0.47 0.40]%.As a result, this experimental validation has also proved the feasibility of the proposed

Figure 1 :
Figure 1: Flow chart of the copula-based FE model updating.

Figure 2 :
Figure 2: A two degree-of-freedom mass-spring system.

Figure 3 :
Figure 3: Scatter map between the updating parameters and the output quantities of mass-spring system.

Figure 4 :
Figure 4: Patch plots for the updating parameters of mass-spring system.

Figure 5 :
Figure 5: Scatter plots for the first two natural frequencies of massspring system.

Figure 7 :
Figure 7: Patch plots for the updating parameters of steel plates.
f 1 vs f 2 vs f 3 f 3 vs f 4 vs f 5

Figure 8 :
Figure 8: Scatter plots for the first five natural frequencies of steel plates.

Figure 9 :Figure 10 :
Figure 9: A view of the beam transport system in SGIII facility target area, which is reproduced from Chen X J et al. (2014).

Figure 11 :
Figure 11: Experiment setup of the transport mirror system.

Figure 12 :
Figure 12: FE model of the transport mirror system.

Figure 13 :Figure 14 :
Figure 13: Mode shapes of the transport mirror system.

Figure 15 :
Figure 15: Scatter plots for the first three natural frequencies of transport mirror system.

Table 1 :
Statistical properties of the measured frequencies of the plates.

Table 3 :
Experiment results of the first three natural frequencies of the transport mirror system.