# Evaluating the Effect of Modelling Errors in Load Identification Using Classical Identification Methods

Authors: Michael Vigsø; Rune Brincker; Christos Georgakis
Journal: Shock and Vibration
 (2019)

---

## Abstract

Load identification, or input identification as the more general term, is a field of study that requires a wide set of disciplines, which suffers from uncertainties caused by the challenges within each discipline. When making load identification, several different approaches exist. For all (or at least most) methods, however, some sort of system model is required. This model may be simple or complex, depending on the system at hand. Typically, if the identification process is vibration fed, the system model will be created from modal parameters. These parameters, however, are often subject to uncertainty and thus may be considered as stochastic variables. In this paper, the root causes of uncertainty for load identification are demonstrated using classical identification techniques. From a numerical perspective, uncertainty is quantified through Monte Carlo simulations. Two results are outlined: one where the identification process is completely blindfolded in its most naive form, and one where the spatial distribution of the load is predefined. In general, it is found that fixing the spatial distribution of the load can compensate for truncation errors in the modal parameters.


---

## Body

## 1. Introduction

Modal analysis is a convenient and efficient method for retrieving the dynamic properties of a civil engineering structure in a condensed form. Traditionally, two methods exist for modal analysis—experimental modal analysis (EMA) and operational modal analysis (OMA). The return from a modal analysis consists of the following modal parameters: mode shapes, natural frequencies, and modal damping [1]. Modal analysis is a commonly used approach, yet the confidence in the extracted modal parameters will depend on the quality of excitation and response measurement. Although new methods are continuously surfacing to overcome some of the challenges, a remaining uncertainty on the estimates persists. In the application for civil engineering structures, another range of challenges may arise, for instance, when the stationary assumption is violated through variations in environmental conditions such as temperature or mass loading.The exertion of indirect load identification may be a successive step to the modal analysis. If the retrieved modal parameters accurately and in whole can describe the dynamics of the system, they can be used in the inverse calculation needed when doing indirect load measurements. The principle is that if a dynamic model is available (from EMA or OMA) and the response of the system is measured, then an estimate of the load can be obtained through deconvolution. Load identification has been in focus in the recent decade, and many attempts have been made in order to obtain a stable and accurate method. In recent years, identification techniques in real time using Kalman filters have proven successful in certain cases [2–7]. Besides the recursive model description, another benefit from these methods is the capability to merge different sensor type information. For these filters to perform well, however, several noise models must be determined, which may be difficult. Aucejo et al. [8] thoroughly examined how the different fitting parameters may influence the performance of identification though they omitted model errors from their analysis. Using identification techniques in the time domain, Wang et al. [9–11] studied the interval envelope on the load estimation given uncertainty in the response measurements and integration errors.Another approach is to deconvolute the response in the frequency domain [12–15]. This was the main approach from the early days in the automotive and aviation industry [16, 17]. The method is based on establishing the relation between the system response and force at different frequencies. The technique requires a recording time, which is why the method does not perform well in real-time applications. Comparative studies between the two approaches are also seen in [18]. The response measurements are often recorded in terms of displacements, velocities, accelerations, or strains. These can be obtained through attached sensors (invasive) or from noninvasive techniques such as lasers, digital image correlation, or even by acoustics [19, 20].For every method referred to above, a system model is needed. This model may be derived from a finite element (FE) model or from modal parameters directly obtained from the EMA or OMA. Now, an obvious question arises: how accurate must the modal parameters be in order for the load identification to be satisfying? Although some authors already have addressed some of these questions [21–24], this paper will revisit some of the fundamental challenges from a visual point of view. By studying a practical example, the sensitivity of estimation errors will be demonstrated through a Monte Carlo simulation. In order not to disappoint the reader, we note that this paper deals merely with the consequence and not mitigation of the estimation errors in modal parameters.


## 2. Simulation Setup

A plane cantilever beam will host the basis for this study. The beam is made from Euler–Bernoulli beam elements with two degrees of freedom (DOFs) at each node—a translational and a rotational as indicated in Figure1. We assume that the beam is proportionally damped through the Rayleigh coefficients α and β. The mass is assumed to be distributed along the elements. The beam is discretized by three elements and fixed at the bottom node which yields, in total, six free DOFs and hence six modes of vibration. We note that the FE model has not yet fully converged for all modes at this nodal resolution; however, this is not of importance for the given study.Figure 1
Simulated system. Left shows the input: 10 s load history applied at each node. Right shows the corresponding response. Dashed lines are the moment and corresponding nodal rotation. Solid lines are the force and the corresponding nodal displacement.The beam model is then subjected to a loading history from which the dynamic response is simulated at a sampling frequency of 1 kHz. The loading time history is selected so that both dynamic response and quasistatic response will be visible within the output. The load is composed of a square impulse followed by a low-frequency wave, again followed by a superimposed high-frequency wave train (15 Hz). Figure1 shows the system in focus for this current study. The four-noded cantilever beam constitutes the system model. The loading history (input) is shown in the left-hand side of the figure. Each line corresponds to a load time series for the corresponding node, i.e., a flat line resembles zero loading. The solid lines are the horizontal force, whereas the dashed lines are the moment. The loading is set to act at one node only—the second node from the top. The right-hand side shows the corresponding response (output) for each node. Solid lines are the horizontal displacements, and the dashed lines show the rotation of the nodes. Amplitude scaling is irrelevant as the system is linear.The process of load identification will hence be the reverse order of the process shown in Figure1. First, the response of the structure is measured. Then, if the dynamic properties of the system are known, we can make an inverse calculation and thus estimate the load causing the measured response. Before proceeding to the load identification, a small amount of noise is added to the system response and hereafter considered as the measured signal (100 dB Gaussian white noise is added).


## 3. Load Identification Process

In order to assess the sensitivity of the load estimate given the variations in the estimated modal parameters, a method must be chosen. For this paper, two methods will be outlined—both operating in the frequency domain. We have selected two methods that easily can be applied using the results from an EMA or OMA (Typically, if modal parameters are obtained from an OMA, problems persist as the modal mass and hence the mode shapescaling is unknown. For now, we will ignore this challenge and assume that the mode shapes obtained are correctly scaled). The methods originate from the automotive and aviator industry [22, 23] and are often referred to as the “transfer path analysis matrix inversion” or just “pseudoinverse technique.”


### 3.1. Method 1

The first method, referred to as Method 1, is the most simple and straight forward. This method might be a tempting first approach, but it shall be demonstrated later on why it could lead to unfortunate results. Consider the response of a linear, time-invariant system which, in the frequency domain, is given as follows:(1)Y¯ω=H¯¯ωF¯ω,where Y¯ω is the system response, H¯¯ω is the frequency response function (FRF), and F¯ω is the load, all given at the discrete frequency ω. This process resembles very well what has been illustrated in Figure 1. The FRF matrix can be computed as a sum of N modal contributions as follows:(2)H¯¯ω=∑r=1NQrϕr¯ϕr¯Tiω−λr+Qr∗ϕr¯∗ϕr¯Hiω−λr∗,where ϕr¯ is the mode shape vector, Qr is the mode scaling function, and λr is the complex pole, all for mode r. T is the transpose, ∗ is the complex conjugate, and H is the Hermitian transpose. See [25] for the derivation of this formula.Now, we presume that the outputY¯ω is available through measurements—both the translation and the rotation of all DOFs. All that is left to do is to calculate the inverse of the FRF and premultiply on both sides in order to obtain an estimate on load F^¯ω.(3)F^¯ω=H¯¯†ωY¯ω.Since errors or modal truncation can cause the inversion of the FRF to become singular, the inversion is performed using the Moore–Penrose pseudoinverse [26]. Method 1 has been used in [15, 27, 28] with the minor change that the inversion of the FRF matrix is performed through a singular value decomposition technique. Regularization techniques have also shown to reduce the sensitivity in the matrix inversion [29].


### 3.2. Method 2

Another method, Method 2, will be added for comparison. This method resembles Method 1 to a high degree with only one minor, yet important, difference. In Method 1, there were no restrictions to the solution, and the load estimateF^¯ω could be distributed to any DOF in the system. If the spatial distribution of the load is known (and unchanging), this can be incorporated into the identification process through a separation of variables as follows:(4)F¯ω=f0¯gω,where f0¯ constitutes the spatial distribution of the load and gω is its scaling function. If the given loading scenario allows this separation and if f0¯ is known, equation (1) and (3) can be rewritten as(5)Y¯ω=H¯¯ωF¯ω,=H¯¯ωf0¯gω,=c¯ωgω,which then can be solved for the scaling constant gω:(6)g^ω=c¯†ωY¯ω,and the final load estimate then reads(7)F^¯ω=g^ωf0¯.Methods to derive the distributionf0¯ in a similar system have been shown in [13, 30–32], but for this study let us assume that the distribution is a known quantity; hence, f0¯=001000T.The equations, (1)–(7), describing Method 1 and Method 2 directly yield a frequency domain estimate of the load. All the following figures in this paper will show the load estimate as a time history in line with the input (left-hand side) in Figure 1. The equivalent time-domain representation is obtained by means of the inverse Fourier transformation.


### 3.3. Condition of the Frequency Response Function

The condition number may be a measure of how sensitive the matrix is to inversion, where a large condition number indicates an ill-conditioned matrix and vice versa [33]. Singular value decomposition may separate the FRF matrix into singular vectors and singular values. Letting the condition number be defined as the ratio between the largest and smallest singular value, we thus obtain a quantitative measure on the sensitivity. Since the FRF matrix, H¯¯ω is defined for a range of frequencies, each frequency will be inherent to a condition number. The condition number(s) for the full rank FRF matrix is shown in Figure 2 as a function of frequency. Here, it is seen how the condition number increases at the resonance frequencies.Figure 2
Condition number for the full rank frequency response function at different frequencies.There is no sharp limit to assess whether a condition number will yield a stable inversion, but a condition number larger than 100 should call for immediate concern which is the case for most frequencies in our case [33]. Since the FRF matrix is constructed from a set of modal contributions, the rank of the FRF matrix will depend on the number of independent modes. If modal truncation is present, this causes the FRF matrix to be rank deficient, which means that the smallest singular value becomes close to zero and consequently yields a condition number going to infinity.(8)condH¯¯ω⟶∞rkH¯¯ω<minm,n,where m and n are the number of rows and columns in H¯¯ω while cond⋅ and rk⋅ indicates the condition number and rank, respectively.If the FRF matrix is constrained by a spatial distribution vectorf0¯ as it is done in Method 2, the product reduces to the vector c¯ω. Singular values of an m×1 nonnull vector is a single nonzero element arranged in an m×1 vector [34]. This consequently provides a condition number of unity, which means that the “inversion” will be stable for all frequencies.(9)condH¯¯ωf0¯≡1.


## 3.1. Method 1

The first method, referred to as Method 1, is the most simple and straight forward. This method might be a tempting first approach, but it shall be demonstrated later on why it could lead to unfortunate results. Consider the response of a linear, time-invariant system which, in the frequency domain, is given as follows:(1)Y¯ω=H¯¯ωF¯ω,where Y¯ω is the system response, H¯¯ω is the frequency response function (FRF), and F¯ω is the load, all given at the discrete frequency ω. This process resembles very well what has been illustrated in Figure 1. The FRF matrix can be computed as a sum of N modal contributions as follows:(2)H¯¯ω=∑r=1NQrϕr¯ϕr¯Tiω−λr+Qr∗ϕr¯∗ϕr¯Hiω−λr∗,where ϕr¯ is the mode shape vector, Qr is the mode scaling function, and λr is the complex pole, all for mode r. T is the transpose, ∗ is the complex conjugate, and H is the Hermitian transpose. See [25] for the derivation of this formula.Now, we presume that the outputY¯ω is available through measurements—both the translation and the rotation of all DOFs. All that is left to do is to calculate the inverse of the FRF and premultiply on both sides in order to obtain an estimate on load F^¯ω.(3)F^¯ω=H¯¯†ωY¯ω.Since errors or modal truncation can cause the inversion of the FRF to become singular, the inversion is performed using the Moore–Penrose pseudoinverse [26]. Method 1 has been used in [15, 27, 28] with the minor change that the inversion of the FRF matrix is performed through a singular value decomposition technique. Regularization techniques have also shown to reduce the sensitivity in the matrix inversion [29].


## 3.2. Method 2

Another method, Method 2, will be added for comparison. This method resembles Method 1 to a high degree with only one minor, yet important, difference. In Method 1, there were no restrictions to the solution, and the load estimateF^¯ω could be distributed to any DOF in the system. If the spatial distribution of the load is known (and unchanging), this can be incorporated into the identification process through a separation of variables as follows:(4)F¯ω=f0¯gω,where f0¯ constitutes the spatial distribution of the load and gω is its scaling function. If the given loading scenario allows this separation and if f0¯ is known, equation (1) and (3) can be rewritten as(5)Y¯ω=H¯¯ωF¯ω,=H¯¯ωf0¯gω,=c¯ωgω,which then can be solved for the scaling constant gω:(6)g^ω=c¯†ωY¯ω,and the final load estimate then reads(7)F^¯ω=g^ωf0¯.Methods to derive the distributionf0¯ in a similar system have been shown in [13, 30–32], but for this study let us assume that the distribution is a known quantity; hence, f0¯=001000T.The equations, (1)–(7), describing Method 1 and Method 2 directly yield a frequency domain estimate of the load. All the following figures in this paper will show the load estimate as a time history in line with the input (left-hand side) in Figure 1. The equivalent time-domain representation is obtained by means of the inverse Fourier transformation.


## 3.3. Condition of the Frequency Response Function

The condition number may be a measure of how sensitive the matrix is to inversion, where a large condition number indicates an ill-conditioned matrix and vice versa [33]. Singular value decomposition may separate the FRF matrix into singular vectors and singular values. Letting the condition number be defined as the ratio between the largest and smallest singular value, we thus obtain a quantitative measure on the sensitivity. Since the FRF matrix, H¯¯ω is defined for a range of frequencies, each frequency will be inherent to a condition number. The condition number(s) for the full rank FRF matrix is shown in Figure 2 as a function of frequency. Here, it is seen how the condition number increases at the resonance frequencies.Figure 2
Condition number for the full rank frequency response function at different frequencies.There is no sharp limit to assess whether a condition number will yield a stable inversion, but a condition number larger than 100 should call for immediate concern which is the case for most frequencies in our case [33]. Since the FRF matrix is constructed from a set of modal contributions, the rank of the FRF matrix will depend on the number of independent modes. If modal truncation is present, this causes the FRF matrix to be rank deficient, which means that the smallest singular value becomes close to zero and consequently yields a condition number going to infinity.(8)condH¯¯ω⟶∞rkH¯¯ω<minm,n,where m and n are the number of rows and columns in H¯¯ω while cond⋅ and rk⋅ indicates the condition number and rank, respectively.If the FRF matrix is constrained by a spatial distribution vectorf0¯ as it is done in Method 2, the product reduces to the vector c¯ω. Singular values of an m×1 nonnull vector is a single nonzero element arranged in an m×1 vector [34]. This consequently provides a condition number of unity, which means that the “inversion” will be stable for all frequencies.(9)condH¯¯ωf0¯≡1.


## 4. Modal Truncation

The number of DOFs in a real structure will, of course, be approaching infinity, but the modes which can be identified are depending on the sensor distribution and capability. Consequently, truncation of the modal space will always be present in real-life situations, which is why this section is included [35, 36].From equation (2), it is seen that the response of a system is a linear contribution from all modes in the system. If modes are omitted from this sum, this will truncate the response function. Figures 3 and 4 show how this modal truncation affects the estimate on the load estimate from the inverse calculations. One additional mode is removed for every case. All modal parameters are treated as deterministic and at their true values, i.e., the only error in the FRF is the higher modes being omitted.Figure 3
Load estimate using Method 1. The figure shows the time-domain load estimate based on the number of modes included in the FRF matrix. Each time series shown has a 10 sec duration. Solid lines show the force (N), and the dashed lines show the moment (Nm).Figure 4
Load estimate using Method 2. The figure shows the time-domain load estimate based on the number of modes included in the FRF matrix. Each time series shown has a 10 sec duration. Solid lines show the force (N), and the dashed lines show the moment (Nm).As seen in Figure3, the load estimate is acceptable when all six modes are used despite the small content of noise on the response measurements. However, when omitting modes from the FRF matrix, the load estimate from Method 1 quickly becomes erroneous, both in terms of distribution and magnitude. Even though the resonance frequency of the modes omitted may be far from the frequency of the load, this truncation still makes a great impact. Meanwhile, in Figure 4, it is seen that Method 2 retrieves more consistent estimates on the load. Any DOF, besides DOF 3, obtains a zero-load estimate following the definition of the spatial distribution given in equation (7). For both methods, it is seen that the estimate for the square impulse overshoots at the discontinuities caused by the Gibbs phenomena for the truncated modal space [37].If we consider the resynthesized responseY^¯ω from the load estimates F^¯ω from either Method 1 or Method 2 and the truncated response function H^¯¯ω, we observe(10)Y^¯ω=H^¯¯ωF^¯ω,εY¯ω,Y^1¯ω<ϵY¯ω,Y^2¯ω.Here,ε is an arbitrary error function between the measured response Y¯ω and the resynthesized response Y^¯ω. The solution from Method 1 will yield a better response approximation than Method 2. This means that Method 1 is a more mathematically accurate solution. Yet, when comparing Figures 3 and 4, it is clear that in a physical sense Method 2 is more consistent. This is one of the major challenges when making load identification solely based on the response using a least-square approach. Consequently, the load estimate may return as an equivalent loading rather than the actual, as seen in Figure 3.For the given simulation, the load is conveniently acting at a point where the response is being “measured.” In other cases one might not be this lucky and expansion is needed in order to estimate the response at locations that were not recorded originally. Several expansion techniques (/virtual sensing techniques) exist for this [38, 39], but effectively, similar truncation errors will be introduced as no new modes are added during the expansion.


## 5. Stochastic Modelling

The modal parameters needed to establish the response function may be obtained through an experimental identification process—either EMA or OMA. The parameters obtained in this process will always be subject to variations. The variations can originate from physical causes such as influence from temperature, operational mass loading, or other nonanticipated nonlinearities violated by the model description. Also, nonphysical uncertainties related to data postprocessing and pole extraction are known to exist [40–43]. These magnify if the excitation of the structure is unfortunate or if the sensor resolution and location is poorly chosen. Noise and other limitations on the sensor may also be a cause of uncertainty.The variations in the modal parameters can be implemented from many different uncertainty models. The most common is probabilistic, where a probability density function is used to describe the statistical variations. Other formulations such as the fuzzy-set model or a nonprobabilistic interval method could also have been utilized. For this study, however, we will rely on the probabilistic approach and assume that all these aforementioned causes of uncertainty effectively can be modelled as an uncorrelated stochastic variation in each of the estimated modal parameters. The uncertainty in natural frequencies is assumed to follow a normal distribution given a coefficient of variation; hence, estimates for the lower frequencies are the most certain. The coefficient of variance is chosen ascv=2% for all frequency estimates. The stochastic model for the natural frequency of mode r hence reads(11)Ωr∼Nμ,σ,where the expectation μ is assumed to be in line with the true natural frequency for mode r, i.e., μ=ωr. The standard deviation is derived from the coefficient of variation as σ=cvωr.Damping estimates are known to be the most uncertain parameters to quantify. Since the damping typically involves fairly low values—say a few percents—and the estimates are subjected to a vast amount of uncertainty, it must have some skewness in order to avoid negative estimates. The study [44] also noted a positive skewness in damping estimates, yet the distribution fit was not studied. For this simulation, we assume that the distribution can be modelled as a standard gamma distribution with shape parameter α. The stochastic damping parameter model for mode r hence reads(12)Zr∼Γα,where the shape parameter is assumed to be equal to the true damping ratio for mode r (in percent). That is, α=ζr which yields a mean EZr≈ζr.We assume that mode shapes always will be real following the proportionally damped system. The variations in mode shapes are implemented as a stochastic process described by a spatially uncorrelated variation in individual DOFs and a variation in scaling of the mode shape. For moder, this means(13)Φr¯=Σ¯¯ϕr¯κ,where ϕr¯ is the true mode shape for mode r, Σ¯¯ is a diagonal matrix whose entries contain the uncorrelated noise on individual DOFs, and κ is the scaling of the mode shape. Both noise models are assumed to be normally distributed with a mean μ and standard deviation σ as(14)Σ¯¯ii∼N1,5×10−3,κ∼N1,1×10−2.The mode shapes are assumed to be mass normalized. Thus, the alteration in mode shapes will affect the modal mass as well. The spatially uncorrelated alterations may violate the mode shape orthogonality, but only to a small degree as the error standard deviation for the mode shapes’ DOFs is set at 0.5%. The consequences of omitting modes have already been shown; hence, only small variations are now included.Using the principles of Monte Carlo simulation, a realization of the stochastic modelling is shown in Figure5 for the first three and the last modes. The uncertainty band for the mode shapes is shown using cubic interpolation in order to show the rotational DOFs. The modal assurance criterion (MAC) [45] between the true mode shape and the stochastically altered is shown for each mode. We note that the variations in MAC values are small compared to what may be experienced from experimental work.Figure 5
Realization of the different modal parameters from 100 000 Monte Carlo simulations. The mode shape uncertainty band is shown for the 98% quantile using cubic interpolation.


## 6. Results

Before turning to the identification process, let us have a look at the FRF matrix. The stochastic variables shown in Figure5 are fed into equation (2) through a Monte Carlo simulation, and a sample of the FRF matrix is shown in Figure 6. For the 100 000 Monte Carlo simulations, the upper and lower 98% quantiles are indicated in the figure as the hatched green area. The black solid line is the mean value. The consequences of missing modes in the FRF matrix were studied by Maes et al. [46] and shall not be repeated here.Figure 6
Sample frequency response function (dynamic flexibility) for DOF 3, given an input in DOF 3. The figure shows the upper and lower 98% uncertainty band, given the modal parameters presented in Figure5. Solid black line indicates the mean value.In order to evaluate the sensitivity in load estimates, the frequency response function is now considered as a stochastic process given the natural variation in modal parameters. Keep in mind that the “measured” system response is kept the same for all simulations and at its true value with a minimum of noise added. In the following three sections, the parameter variations are introduced a little at a time so that it will be more clear what happens to the estimates. Since the mean value of the modal parameters is in line with the true values, the averaged value of the estimates converges towards the true value as well. All six modes will be used in the making of the FRF matrix.


### 6.1. Natural Frequencies

First, only the natural frequencies are treated as stochastic variables, and the remaining modal parameters are kept as deterministic and exact. For the 100 000 Monte Carlo simulations, the FRF is synthesized, and a load estimate for each simulation is obtained through both methods. Figure7 shows the upper and lower 98% quantiles for the time-domain load estimates. It is seen that both methods yield a reasonably stable result despite variation in natural frequencies. Three points in time (a, b, and c) are extracted to highlight the error distribution on the estimate. The relative error for the two methods is given by the error standard deviation and summarized in Table 1.Figure 7
Natural frequencies are introduced as stochastic variables. The figure consequently shows the time history of the identified load using Method 1 and Method 2. Each time series shown has a 10 sec duration. The solid lines are the mean values of the estimated force (N), and the dashed line is the mean value of the estimated moment (Nm). The hatched green area indicates uncertainty on the force, while the hatched orange area is the uncertainty of the moment. Both are shown using the upper and lower 98% quantiles. Three selected points in time show the error probability density (in percent) for the estimated load in DOF 3.Table 1
Error standard deviation for different levels of uncertainty and for three selected points in time. All are based on the force estimate in DOF 3.
Uncertainties introduced
Method 1
Method 2
Error std (%)
Error std (%)
a
b
c
a
b
c
Natural frequencies
2.73
2.29
2.13
6.08
4.39
3.64
Natural frequencies and damping ratios
3.28
2.31
2.12
7.67
4.45
3.64
Natural frequencies, damping ratios, and mode shapes
53.1
46.3
36.9
7.96
4.94
4.10


### 6.2. Natural Frequencies and Damping Ratios

Next, in addition to variations in natural frequency, the damping ratios are also treated as stochastic variables. That means only the mode shapes are left deterministic. The variation in the estimated load is shown for the two methods in Figure8. As expected, the uncertainty bound for the estimates increases as the damping ratios are also treated as stochastic parameters. Method 1 seems to outperform Method 2, given this level of uncertainty in the modal parameters and noise in the output signal. A ringing effect in the estimate is observed at the discontinuities near the square impulse. In Method 1, however, it is seen how the load estimates are spreading to different nodes to compensate for the errors in the FRF matrix.Figure 8
Natural frequencies and damping ratios are introduced as stochastic variables. The figure consequently shows the time history of the identified load using Method 1 and Method 2. Each time series shown has a 10 sec duration. The solid lines are the mean value of the estimated force (N), and the dashed line is the mean value of the estimated moment (Nm). The hatched green area indicates uncertainty on the force, while the hatched orange area is the uncertainty of the moment. Both are shown using the upper and lower 98% quantiles. Three selected points in time show the error probability density (in percent) for the estimated load in DOF 3.


### 6.3. Natural Frequencies, Damping Ratios, and Mode Shapes

Finally, all modal parameters are now considered as stochastic processes, and the consequent result for the estimated load is given in Figure9. Although the mode shapes are altered by values down to one or two percent for each node, the uncertainty in the load estimates for Method 1 shows an exponential growth. When mode shape errors occur, Method 1 fails to predict the load distribution and consequently the scaling of the load. For reference, the sample FRF shown in Figure 6 hosts the basis for the load estimate given in Figure 9. Note that in Figure 9, the error probability density is unequally scaled for Method 1 and Method 2.Figure 9
Natural frequencies, damping ratios, and mode shapes are introduced as stochastic variables. The figure consequently shows the time history of the identified load using Method 1 and Method 2. Each time series shown has a 10 sec duration. The solid lines are the mean value of the estimated force (N), and the dashed line is the mean value of the estimated moment (Nm). The hatched green area indicates uncertainty on the force, while the hatched orange area is the uncertainty of the moment. Both are shown using the upper and lower 98% quantiles. Three selected points in time show the error probability density (in percent) for the estimated load in DOF 3. Note that the probability density is scaled differently between Method 1 and Method 2.Each of the load estimates shown by Figures7 to 9 has three highlighted points in time, where an error probability density is shown. The corresponding standard deviation (std) is given in Table 1. We note that Method 1 seems to perform better than Method 2 as long as the mode shapes are intact. However, when introducing minor changes to the mode shapes, Method 1 cannot fully predict to which node the load is applied, which results in a poor estimate on DOF 3.


## 6.1. Natural Frequencies

First, only the natural frequencies are treated as stochastic variables, and the remaining modal parameters are kept as deterministic and exact. For the 100 000 Monte Carlo simulations, the FRF is synthesized, and a load estimate for each simulation is obtained through both methods. Figure7 shows the upper and lower 98% quantiles for the time-domain load estimates. It is seen that both methods yield a reasonably stable result despite variation in natural frequencies. Three points in time (a, b, and c) are extracted to highlight the error distribution on the estimate. The relative error for the two methods is given by the error standard deviation and summarized in Table 1.Figure 7
Natural frequencies are introduced as stochastic variables. The figure consequently shows the time history of the identified load using Method 1 and Method 2. Each time series shown has a 10 sec duration. The solid lines are the mean values of the estimated force (N), and the dashed line is the mean value of the estimated moment (Nm). The hatched green area indicates uncertainty on the force, while the hatched orange area is the uncertainty of the moment. Both are shown using the upper and lower 98% quantiles. Three selected points in time show the error probability density (in percent) for the estimated load in DOF 3.Table 1
Error standard deviation for different levels of uncertainty and for three selected points in time. All are based on the force estimate in DOF 3.
Uncertainties introduced
Method 1
Method 2
Error std (%)
Error std (%)
a
b
c
a
b
c
Natural frequencies
2.73
2.29
2.13
6.08
4.39
3.64
Natural frequencies and damping ratios
3.28
2.31
2.12
7.67
4.45
3.64
Natural frequencies, damping ratios, and mode shapes
53.1
46.3
36.9
7.96
4.94
4.10


## 6.2. Natural Frequencies and Damping Ratios

Next, in addition to variations in natural frequency, the damping ratios are also treated as stochastic variables. That means only the mode shapes are left deterministic. The variation in the estimated load is shown for the two methods in Figure8. As expected, the uncertainty bound for the estimates increases as the damping ratios are also treated as stochastic parameters. Method 1 seems to outperform Method 2, given this level of uncertainty in the modal parameters and noise in the output signal. A ringing effect in the estimate is observed at the discontinuities near the square impulse. In Method 1, however, it is seen how the load estimates are spreading to different nodes to compensate for the errors in the FRF matrix.Figure 8
Natural frequencies and damping ratios are introduced as stochastic variables. The figure consequently shows the time history of the identified load using Method 1 and Method 2. Each time series shown has a 10 sec duration. The solid lines are the mean value of the estimated force (N), and the dashed line is the mean value of the estimated moment (Nm). The hatched green area indicates uncertainty on the force, while the hatched orange area is the uncertainty of the moment. Both are shown using the upper and lower 98% quantiles. Three selected points in time show the error probability density (in percent) for the estimated load in DOF 3.


## 6.3. Natural Frequencies, Damping Ratios, and Mode Shapes

Finally, all modal parameters are now considered as stochastic processes, and the consequent result for the estimated load is given in Figure9. Although the mode shapes are altered by values down to one or two percent for each node, the uncertainty in the load estimates for Method 1 shows an exponential growth. When mode shape errors occur, Method 1 fails to predict the load distribution and consequently the scaling of the load. For reference, the sample FRF shown in Figure 6 hosts the basis for the load estimate given in Figure 9. Note that in Figure 9, the error probability density is unequally scaled for Method 1 and Method 2.Figure 9
Natural frequencies, damping ratios, and mode shapes are introduced as stochastic variables. The figure consequently shows the time history of the identified load using Method 1 and Method 2. Each time series shown has a 10 sec duration. The solid lines are the mean value of the estimated force (N), and the dashed line is the mean value of the estimated moment (Nm). The hatched green area indicates uncertainty on the force, while the hatched orange area is the uncertainty of the moment. Both are shown using the upper and lower 98% quantiles. Three selected points in time show the error probability density (in percent) for the estimated load in DOF 3. Note that the probability density is scaled differently between Method 1 and Method 2.Each of the load estimates shown by Figures7 to 9 has three highlighted points in time, where an error probability density is shown. The corresponding standard deviation (std) is given in Table 1. We note that Method 1 seems to perform better than Method 2 as long as the mode shapes are intact. However, when introducing minor changes to the mode shapes, Method 1 cannot fully predict to which node the load is applied, which results in a poor estimate on DOF 3.


## 7. Discussion

We have studied the performance of the two methods given the natural variation in the modal parameters used for the model description. The load estimates when using Method 1, shown in Figures3 and 9, are not as diverse as they might appear at first glance. If the resulting force is considered instead of every single entry in the estimated load vector, this leads to a more appealing result. First, reexamining the estimates shown in Figure 9 and by summarizing the contributions from each DOF, an equivalent baseline load is obtained. This resulting load estimate is compared with the true baseline load and shown in Figure 10. Now, it is seen how the uncertainty band is narrowed down. For this static system, an equivalent global loading can be estimated using Method 1 despite the stochastic variations in modal parameters as long as all modes are represented.Figure 10
Estimated resulting baseline forces given stochastic modelling used for Figure9 for Method 1. The black line is the true baseline force, and the hatched area shows the upper and lower 98% uncertainty bounds for the estimated moment and force, respectively.
(a)
(b)Returning to the modal truncation example in Section4, we once again compute the resulting load estimate from Method 1 as previously shown in Figure 3. Using the same analogy as explained above, a resulting baseline load is estimated for different number of modes included in the FRF. As the estimates from the truncation study are deterministic, we choose to plot each result on top of each other in Figure 11 together with the true baseline load. We see that the overturning moment is fairly accurate regardless of the number of modes included, while the shear force is drifting when modes are being omitted. When modes are omitted, Method 1 fails to estimate the point of attack and thus the scaling of the load, as it is seen in Figure 11. It should be noted that the degree of overshooting at the impulse is increased as modes are being omitted. Whether the resulting moment estimate remains stable is likely dependent on the static system.Figure 11
Resulting baseline forces given the different levels of modal truncation for Method 1. Every coloured line corresponds to an estimate given a number of modes, while the black lines are the true input.
(a)
(b)There are some items which have not been covered by this study and which should be mentioned as they limit the conclusions.(i)
We have chosen a simple static system with plane deformation. The system is blessed by not having any closely spaced modes and is behaving perfectly linear. Now, if the system is more complex, which may be the case for most civil engineering structures, the generalization of the number of modes needed and the precision on the modal parameters may differ.(ii)
We have not touched upon multiple load sources or a moving load.(iii)
For this paper, the response function is driven by displacement measurements. It has not been assessed whether acceleration or strain measurements would have yielded different results.(iv)
Only one noise level on the output signal has been studied.(v)
From Figure8, it may seem like the variation in damping ratios and frequency yields nearly no effect on the load estimate. If the system is subjected to a harmonic load at a frequency near a natural frequency, these modal values may be vital for the load estimate.


## 8. Conclusion

The sensitivity in load identification following the uncertainty in modal parameters has been studied for two different methods. The uncertainties have been introduced stepwise into the system to demonstrate the influence of different parameters. It was found that fixing the spatial distribution of the load can compensate for the modal truncation in the response function. If one mode is dominating in the response, parameter estimation following this mode is crucial for the following load identification. However, modes with a resonant frequency above the frequency of the load still contribute to the quasistatic response, and omitting these modes will cause an error. A possible idea is to compensate through static deflection shapes/Ritz vector in the FRF matrix, but this has not been included. The spatial distribution of the load was fixed in Method 2 through the distribution vectorf0¯. If the assumption regarding the load distribution f0¯ is not correct, this will introduce some systematic errors in the load estimate. Demonstrating this has been omitted from this paper.Whether or not the stochastic modelling outlined in Section5 is appropriate for a real-life structure will be left for the reader to decide. However, it has been demonstrated how variations in these parameters may alter the estimates on the loading, and shown how a sensitivity analysis may reveal flaws in the algorithm.As a general observation, when doing vibration-based load identification, the algorithm for estimating the input must be supported by additional information besides the system response. This additional information may be in terms of the spatial distribution of the load—which was shown with Method 2—or it may be other load models which are driven by wind speed measurements, wave gauge readings, local pressure measurements, or possibly information about the frequency content of the loading. Without aiding the inverse calculations with one of these, the uncertainty of the estimate will exceed what is acceptable. For most of the successful methods available in the literature, the spatial distribution of the load is also defined as a fixed measure, which leaves room for further research in this field.



---

*Source: 9490760-2019-10-21.xml*