Augmented Tikhonov Regularization Method for Dynamic Load Identification

We introduce the augmented Tikhonov regularization method motivated by Bayesian principle to improve the load identification accuracy in seriously ill-posed problems. Firstly, the Green kernel function of a structural dynamic response is established; then, the unknown external loads are identified. In order to reduce the identification error, the augmented Tikhonov regularization method is combined with the Green kernel function. It should be also noted that we propose a novel algorithm to determine the initial values of the regularization parameters. The initial value is selected by finding a local minimum value of the slope of the residual norm. To verify the effectiveness and the accuracy of the proposed method, three experiments are performed, and then the proposed algorithm is used to reproduce the experimental results numerically. Numerical comparisons with the standard Tikhonov regularization method show the advantages of the proposed method. Furthermore, the presented results show clear advantages when dealing with ill-posedness of the problem.


Introduction
The dynamic loads acting on a structure cannot be directly measured in many cases, and it is necessary to use indirect methods. Among such methods, retrieving the actual dynamic loads by measuring the vibration response of a structure has attracted a wide attention. In this case, the dynamic load identification problem involves solving an ill-posed inverse problem. Hence, any small amount of measurement noise causes a huge disturbance of identification results. For example, it is common to have small singular values due to the measurement accuracy error. Although such values are often very small, the ill-posedness of the problem makes it difficult to get the correct solution. In general, the solution accuracy is often affected by an order of magnitude compared to the error in the measurement accuracy.
Many researchers investigated the ill-posedness problem, and several approaches were proposed [1][2][3][4][5][6][7][8]. In principle, these approaches are always based on regularization methods. For example, the Tikhonov regularization method can be used to identify the impact load, employing the frequency response function and the strain response of an experimental system [9]. Tikhonov regularization method, based on ordinary cross validation (OCV) [10], is also used to identify the excitation acting on a thin plate structure. The presented results show the robustness of Tikhonov regularization method when compared to the truncated singular value decomposition method. In another work, three regularization parameters, namely OCV, generalized cross validation (GCV), and L-curve method, are studied [11]. It is shown that the L-curve method combined with Tikhonov regularization method

Green Kernel Function
For a multi-degrees-of-freedom system, the general vibration differential equation is often written as M ..
Here, M, K,C are the mass, stiffness and damping matrix of the system, respectively, while x, . x, ..
x are the displacement, the velocity, and the acceleration vectors of the system response. The variable f represents external forces acting on the system. Assuming a zero initial state, i.e., x 0 and .
x 0 are both zero, the displacement response of the i-th degree of freedom is where x ij (t) is the displacement response of i-th degree of freedom caused by the external force acting on the j-th degree, and h ij is the unit impulse response function, while f j represents the external force acting on the j-th degree of freedom.
Assuming the sampling interval of the response is ∆t, x ij can be simplified using integral discretization Rewrite Equation (3) into matrix form as (1) x ij (2) . . .
Equation (4) can be abbreviated as Here, Y ij is a discrete vector form of x ij , while F j is a discrete vector form of f j . When only one load is acting on the system, i.e., on the j-th degree of freedom, the response of the i-th degree of freedom can be expressed as Equation (6) is the linear equations in the case of a single input and a single output. The unknown force F j can be obtained by solving Equation (6). Here, H ij is the kernel function matrix. When there are M loads acting on different degrees of freedom, the number of the responses to be collected N needs to meet the following formula M ≤ N. Thus, the relationship between M loads and N responses can be obtained: H ij F j (i = 1, 2, · · · , N) .
Rewrite Equation (7) in matrix form as Equation (8) is the linear equations in the case of multiple inputs and outputs, and the unknown force F 1 , F 2 , · · · , F M can be obtained by solving this equation. Clearly, Equations (6) and (8) both conform to the form of y = Hf. The inversion of equations in this form usually introduces a significant error into the solution due to the ill-posedness of the problem that is caused by the ill-posedness of the kernel function matrix H, especially for large systems, i.e., hundreds of degrees of freedom. Next, we discuss reducing the ill-posedness problem using Augmented Tikhonov regularization method.

Augmented Tikhonov Regularization Method
The augmented Tikhonov regularization is an application of the Bayesian principle in the field of inverse problems. The basic idea of the Bayesian theory is to obtain the posterior probability density function through the likelihood function and the prior probability density function, and then estimate the unknown quantity according to the posterior probability density.
In the model y = Hf mentioned above, the output signal y will always contain some noise µ.
According to the Bayesian principle, the posterior probability density function of the unknown quantity f can be obtained: where p(f y ) is the posterior probability density function, while p(f) is the prior probability density function. The variable p(y f ) is the likelihood function. Thus, the loads f can be estimated by the maximum posterior probability density function: The likelihood function represents the probability density function of y when the random variable f is given. The form of the likelihood function is generally determined by the noise µ. Usually, the white Gaussian noise with a zero mean and a standard deviation of σ can provide a proper representation for µ. The noise likelihood function can then be represent as where n is the data length of y.
The above probability density function needs to be selected according to the empirical method. In this study, Markov random field is selected as the prior model.
where γ is a magnification factor, while φ(γ( f i − f j )) is a non-negative function that determines the specific form of this Markov field, and W ij is a non-negative weight coefficient. The function φ is then set to the form φ(t) = 1/2t 2 , so that the Markov random field is where the variable λ denotes the scaling parameter, and m the dimension of vector f. W represents the symmetric matrix composed of W ij . The posterior probability density function can be obtained by substituting Equations (12) and (14) into Equation (10), and the following is obtained: Appl. Sci. 2020, 10, 6348

of 23
In the Bayesian hierarchical model, the two unknown parameters λ, σ in Equation (15) can be regarded as random variables; hence, they play the role of regularization parameters. For convenience, we set τ= 1/σ 2 . The prior distributions of λ, τ are assumed to be a Gamma distribution.
where α 1 , β 1 , α 2 , β 2 are the parameters of the Gamma distribution. Meanwhile, due to the introduction of two new random variables, the likelihood function and the prior probability density function are different from the previous ones.
Then, the joint posterior probability density function can be obtained Substitute Equations (16)- (19) into Equation (20), and we can obtain the specific expression The values of f, λ, τ can be estimated simultaneously by maximizing Equation (21). In order to facilitate the solution, we convert the Equation (21) into the minimum value of the following formula by using the negative form of the log-likelihood function: The alternating direction algorithm is used to solve the minimum problem in this study. The solution steps are as follows: • Calculate the partial derivatives of log-likelihood function with respect to the variables f, λ, τ respectively, and let them vanish.
• Three formulas can be obtained from the above formulas: (28) • The unknown force f can be calculated by setting appropriate initial values and repeating iterative Equations (26)-(28).
The flow chart are shown in the Figure 1, and specific iteration steps are:  (2) to repeat the iterative calculation until the termination condition is met. The termination condition in this study is

Initial Parameters Selection
As discussed above, the solution with augmented Tikhonov regularization method, six initial parameters need to be defined. Among these, the values of λ and τ changes after each iteration, while those of 1 1 2 2 , , , α β α β do not. In order to understand the effect of the initial values of , λ τ , we run hundreds of simulations, where we noticed that the initial values of these parameters have little effect on the load identification algorithm. We also noticed that setting the parameters 1 2 α α ， close to unity leads to stable results with accurate load identification. However, the choice of 1 2 , β β may
Calculate f k : f k = H T H + (λ k /τ k )I −1 H T y;

3.
Calculate new λ k+1 and τ k+1 : Hf k −y 2 2 +2β 2 ; 4. Set k = k + 1, and return to step (2) to repeat the iterative calculation until the termination condition is met. The termination condition in this study is

Initial Parameters Selection
As discussed above, the solution with augmented Tikhonov regularization method, six initial parameters need to be defined. Among these, the values of λ and τ changes after each iteration, while those of α 1 , β 1 , α 2 , β 2 do not. In order to understand the effect of the initial values of λ, τ, we run hundreds of simulations, where we noticed that the initial values of these parameters have little effect on the load identification algorithm. We also noticed that setting the parameters α 1 , α 2 close to unity leads to stable results with accurate load identification. However, the choice of β 1 , β 2 may lead to a significant change in the identification results. Therefore, in the following, we limit our discussion to the selection of these two parameters.
The correlation coefficient (r) and relative error (RE) are used to estimate the load identification error in this paper. These are calculated as follows: where f iden is the identified force, and f true is the true force, while f iden and f true are the means of f iden and f true , respectively. The values of β 1 , β 2 can have a significant impact on the load identification results f k 2 2 and the residual norm Hf k − y 2 2 . To analyze the impact of selecting the values of β 1 , β 2 on the load identification, we study, in this section, a single degree of freedom system. The considered parameters for this system are m = 1 kg , k = 1500 N/m and c = 15 N·s/m . The load applied on the system is periodic defined by F = 100 * sin(π * 50 * t). We also set the parameters α 1 , α 2 to 2, which is relatively a small value while the parameters λ, τ are set to 100 and 1 × 10 12 , respectively. To understand the effect of the parameters β 1 and β 2 on the identified load, their values will be varied, and the identified load will be recorded. Figure 2 shows the changes in the correlation coefficient and the relative error for varying the values of β 1 , β 2 when the signal noise ratio (SNR) of the response is 30 dB. It can be seen in the figure that the optimum range of values for β 1 is from 10 −5 to 10 5 , and, for β 2 , is a narrow band close to 10 4 . In this range, the correlation coefficient is close to 1, while the relative error is close to 0. This suggests that the load identification accuracy is higher at this specific range of parameters β 1 , β 2 . Clearly, the optimum range for β 1 is much wider compared to β 2 . Thus, choosing an optimum value for β 2 is our next concern.
The correlation coefficient (r) and relative error (RE) are used to estimate the load identification error in this paper. These are calculated as follows: where iden f is the identified force, and true f is the true force, while iden f and true f are the means of iden f and true f , respectively.
The values of 1 2 , β β can have a significant impact on the load identification results  , β β on the load identification, we study, in this section, a single degree of freedom system. The considered parameters for this system are 1 The load applied on the system is periodic defined by 100*sin( *50 * ) F t π = .
We also set the parameters 1 2 α α ， to 2, which is relatively a small value while the parameters , λ τ are set to 100 and 12 1 10 × , respectively. To understand the effect of the parameters 1 β and 2 β on the identified load, their values will be varied, and the identified load will be recorded. Figure 2 shows the changes in the correlation coefficient and the relative error for varying the values of 1 2 , β β when the signal noise ratio (SNR) of the response is 30 dB. It can be seen in the figure that the optimum range of values for 1 β is from 10 −5 to 10 5 , and, for 2 β , is a narrow band close to 10 4 . In this range, the correlation coefficient is close to 1, while the relative error is close to 0. This suggests that the load identification accuracy is higher at this specific range of parameters 1 2 , Clearly, the optimum range for 1 β is much wider compared to 2 β . Thus, choosing an optimum value for 2 β is our next concern.
(a) In real world applications, because the load acting on the structures is unknown, it is not possible to define the correlation coefficient and the relative error. Therefore, to determine the optimum values for the parameter 2 β , we use the relation between the parameter 2 β and the residual norm In real world applications, because the load acting on the structures is unknown, it is not possible to define the correlation coefficient and the relative error. Therefore, to determine the optimum values for the parameter β 2 , we use the relation between the parameter β 2 and the residual norm ||H f − y|| 2 2 of the identified load. Indeed, in the load identification process, we rebuild the dynamic responses according to the identified load for a given parameter β 2 and then calculate the residual norm ||H f − y|| 2 2 . When the noise level of the dynamic responses is known, we can find the optimum β 2 by comparing the real residual norm with the computation residual norm, which is shown in Figure 3a. When the real noise level is not available, we calculate the slope of the residual norm ||H f − y|| 2 2 with respect to β 2 and then take the local minimum value as the optimum point of the parameter β 2 , as shown in Figure 3b. In Figure 3, the blue curve represents the slope of the residual norm, while the red curve represents the relative error of the identification load. The black curve gives the real residual norm according to the response noise level. The pink vertical curve represents the selected β 2 value corresponding to the minimum slope point of the residual norm. The red curve shows that the load identification error is ideal with the selected β 2 value. It should be pointed out, in this case, that the response noise level is 30 dB (SNR). Based on the above parametric study, we choose the following values for α 1 , β 1 , α 2 , λ, τ, 2.1 × 10 −5 , 2100, 1 × 10 12 , respectively, while β 2 is determined according to the selection method shown in Figure 3.

Numerical Results
In order to study the efficiency and the accuracy of the proposed augmented Tikhonov regularization method on dynamic load identification, three numerical simulation examples are studied. These include two multi-degrees-of-freedom system examples and a cantilever beam example.

Single Force Identification Simulation
We first aim to identify the load applied on a four-degrees-of-freedom system. The system is shown in Figure 4, and the system parameters are listed underneath.
according to the response noise level. The pink vertical curve represents the selected 2 β value corresponding to the minimum slope point of the residual norm. The red curve shows that the load identification error is ideal with the selected 2 β value. It should be pointed out, in this case, that the response noise level is 30 dB (SNR). Based on the above parametric study, we choose the following values for 1 1 2 , , , , α β α λ τ

Numerical Results
In order to study the efficiency and the accuracy of the proposed augmented Tikhonov regularization method on dynamic load identification, three numerical simulation examples are studied. These include two multi-degrees-of-freedom system examples and a cantilever beam example.

Single Force Identification Simulation
We first aim to identify the load applied on a four-degrees-of-freedom system. The system is shown in Figure 4, and the system parameters are listed underneath.  The regularization parameters are obtained following the procedure described in Section 4. Three types of force are studied: simple harmonic excitation, triangular wave excitation, and impact load. The load is only applied to the first degree of freedom. We assume the system initially has a stationary state, i.e., zero initial conditions. The sampling time interval is 0.001 s. The considered time span is 0.5 s. The load identification results of the three loads with the SNR of 40 dB, 30 dB, and 20 dB are shown in Figures 5-13. On the left side of these figures, we show the results obtained with the Tikhonov regularization method, and, on the right, we show the result of the proposed augmented Tikhonov regularization method. The identification result errors are shown in Table 1.
Appl. Sci. 2020, 10, x FOR PEER REVIEW 10 of 24 The regularization parameters are obtained following the procedure described in Section 4. Three types of force are studied: simple harmonic excitation, triangular wave excitation, and impact load. The load is only applied to the first degree of freedom. We assume the system initially has a stationary state, i.e., zero initial conditions. The sampling time interval is 0.001 s. The considered time span is 0.5 s. The load identification results of the three loads with the SNR of 40 dB, 30 dB, and 20 dB are shown in Figures 5-13. On the left side of these figures, we show the results obtained with the Tikhonov regularization method, and, on the right, we show the result of the proposed augmented Tikhonov regularization method. The identification result errors are shown in Table 1.            Table 1 shows that the relative error and the correlation coefficient for the results of the augmented Tikhonov regularization method are significantly improved when compared to the Tikhonov method for both considered noise levels, i.e., SNR = 40 dB, 30 dB, and 20 dB. The figures show that the augmented Tikhonov regularization method leads to better identification accuracy for sinusoidal, impact and triangular wave load loads. They also show that the augmented method accurately identify the peak position of the triangular wave load, as well as the impact load. In general, the augmented Tikhonov regularization method shows much better accuracy for identifying sinusoidal loads under all considered noise conditions. For the impact load, the augmented Tikhonov regularization method tends to suppress the oscillation of the entire load, so the identification accuracy is very good when the noise is low, while the error in identifying the peak load level for higher noise, which is expected. This is also observed for other types of load where the augmented method shows a significant reduction in spurious oscillations compared to the standard method.

Two Force Identification Simulation
In this subsection, we study the proposed method when capturing loads applied at multiple degrees of freedom. To this end, two forces are applied on the same four-degrees-of-freedom system studied in Section 5.1. The loads are located at the first and the second degrees of freedom. The displacements are measured at the first three degrees of freedom in order to identify the unknown loads. Again, the identification results using traditional Tikhonov regularization method and the augmented Tikhonov regularization method are shown alongside in the Figures 14-21. The considered noise levels of the responses SNR are assumed to be 40 dB and 30 dB. We take two external loading cases sinusoidal and impact loads. The accuracy of the identified loads are shown in Tables 2 and 3. sinusoidal, impact and triangular wave load loads. They also show that the augmented method accurately identify the peak position of the triangular wave load, as well as the impact load. In general, the augmented Tikhonov regularization method shows much better accuracy for identifying sinusoidal loads under all considered noise conditions. For the impact load, the augmented Tikhonov regularization method tends to suppress the oscillation of the entire load, so the identification accuracy is very good when the noise is low, while the error in identifying the peak load level for higher noise, which is expected. This is also observed for other types of load where the augmented method shows a significant reduction in spurious oscillations compared to the standard method.

Two Force Identification Simulation
In this subsection, we study the proposed method when capturing loads applied at multiple degrees of freedom. To this end, two forces are applied on the same four-degrees-of-freedom system studied in Section 5.1. The loads are located at the first and the second degrees of freedom. The displacements are measured at the first three degrees of freedom in order to identify the unknown loads. Again, the identification results using traditional Tikhonov regularization method and the augmented Tikhonov regularization method are shown alongside in the Figures 14-21. The considered noise levels of the responses SNR are assumed to be 40 dB and 30 dB. We take two external loading cases sinusoidal and impact loads. The accuracy of the identified loads are shown in Tables  2 and 3.             Table 3. Identification effect of two impact loads.  The results again confirm the previous set of conclusions. The accuracy of the augmented Tikhonov regularization method is consistently better than the Tikhonov regularization method. The augmented Tikhonov regularization method can effectively capture multiple loads for the four-degrees-of-freedom system. It should be noted that the identification accuracy is relatively lower when compared to a single load case. This is caused by the interaction in between the two loads. However, the identification results still capture all the dynamics of the external loads.

Single Force Identification on the Cantilever Beam
Next, to investigate the proposed augmented method for continuous system, we study the load identification on cantilever beam. The cantilever model is shown in Figure 22, and the parameters of the cantilever beam are: length l = 1 m, cross-sectional area A = 0.0004 m 2 , moment of inertia I = 2.13 × 10 −9 m 4 , Modulus of Elasticity E = 209 Gpa, density ρ = 7800 kg/m 3 , and Poisson's ratio 0.3. The unknown load is applied at point 15, and the response is measured at point 17. Both points are displayed in Figure 22. Again, sinusoidal and impact loads are considered with the two previously considered noise levels, i.e., SNR = 30 and 40 dB. The system starts with a zero initial state. The time sampling interval is 1/1024 s, and the simulation time span is 1 s. For demonstration, the measured system response to the impulse load is shown in Figure 23. The identification results are listed in Figures 24-27, and the accuracy of the load identification results are shown in Table 4.           The results show that the proposed method can accurately identify the unknown loads on a continuous system, such as the considered cantilever beam. Consistent with the previous examples, the augmented Tikhonov method is more appropriate for solving the ill-posed problem of the load identification than the standard Tikhonov method. The sinusoidal load, under the two considered noise conditions, is identified with much better accuracy with the augmented method compared to the standard Tikhonov regularization method. The correlation coefficient and relative error of the results identified by tradition Tikhonov method are not as good as the proposed method, which suggests that the proposed method makes an improvement in the regularization parameters selected.

Experimental Examples
In this section, we aim to validate the proposed augmented method using experimental results. To this end, a load identification experiment for a simply supported bilateral beam is carried out. The The results show that the proposed method can accurately identify the unknown loads on a continuous system, such as the considered cantilever beam. Consistent with the previous examples, the augmented Tikhonov method is more appropriate for solving the ill-posed problem of the load identification than the standard Tikhonov method. The sinusoidal load, under the two considered noise conditions, is identified with much better accuracy with the augmented method compared to the standard Tikhonov regularization method. The correlation coefficient and relative error of the results identified by tradition Tikhonov method are not as good as the proposed method, which suggests that the proposed method makes an improvement in the regularization parameters selected.

Experimental Examples
In this section, we aim to validate the proposed augmented method using experimental results. To this end, a load identification experiment for a simply supported bilateral beam is carried out. The  The results show that the proposed method can accurately identify the unknown loads on a continuous system, such as the considered cantilever beam. Consistent with the previous examples, the augmented Tikhonov method is more appropriate for solving the ill-posed problem of the load identification than the standard Tikhonov method. The sinusoidal load, under the two considered noise conditions, is identified with much better accuracy with the augmented method compared to the standard Tikhonov regularization method. The correlation coefficient and relative error of the results identified by tradition Tikhonov method are not as good as the proposed method, which suggests that the proposed method makes an improvement in the regularization parameters selected.

Experimental Examples
In this section, we aim to validate the proposed augmented method using experimental results. To this end, a load identification experiment for a simply supported bilateral beam is carried out.
The model parameters are listed in Table 5. The natural frequency of the first four mode shapes are listed in Table 6. The frequencies are measured by modal experiment. To consider the model error and calculate the kernel function matrix, we build a finite element model and update it using the natural frequency obtained in the modal experiment. Table 6 shows the comparison of the natural frequency between the experimental results and the numerical results of the finite element model. The relative errors are evaluated based on the differences between the experimental and the numerical results. All the errors are smaller than 3%. Furthermore, for the first frequency, which is the most dominant, the error is less than 0.2%. These results verify that the finite element model is, in general, representative of the physical beam, and the numerical results are reliable. For demonstration, we show, in Figure 28, the displacements response obtained with the finite element model when an impulse load is applied. model parameters are listed in Table 5. The natural frequency of the first four mode shapes are listed in Table 6. The frequencies are measured by modal experiment. To consider the model error and calculate the kernel function matrix, we build a finite element model and update it using the natural frequency obtained in the modal experiment. Table 6 shows the comparison of the natural frequency between the experimental results and the numerical results of the finite element model. The relative errors are evaluated based on the differences between the experimental and the numerical results. All the errors are smaller than 3%. Furthermore, for the first frequency, which is the most dominant, the error is less than 0.2%. These results verify that the finite element model is, in general, representative of the physical beam, and the numerical results are reliable. For demonstration, we show, in Figure 28, the displacements response obtained with the finite element model when an impulse load is applied. Poisson's ratio μ 0.30   Figure 29 shows the experimental setting of the considered beam, while, in Table 7, we specify the used instruments. To apply a load on the beam, the exciter is installed at 0.42 m from one end. The force sensor is installed on the exciter connecting rod to measure the actual applied load. The  Figure 29 shows the experimental setting of the considered beam, while, in Table 7, we specify the used instruments. To apply a load on the beam, the exciter is installed at 0.42 m from one end. The force sensor is installed on the exciter connecting rod to measure the actual applied load. The laser displacement sensor is installed at 0.28 m from the beam end to measure the displacement response. A sampling rate of 4096 is used. Two different loads are applied. First, a sinusoidal and then a triangle wave signal is applied using the exciter. The load magnitude is adjusted through a power amplifier.
The load applied is then recovered using the proposed method and the measured displacement. The identified load is compared to the actual load measured as by the force sensor.
Appl. Sci. 2020, 10, x FOR PEER REVIEW 20 of 24 laser displacement sensor is installed at 0.28 m from the beam end to measure the displacement response. A sampling rate of 4096 is used. Two different loads are applied. First, a sinusoidal and then a triangle wave signal is applied using the exciter. The load magnitude is adjusted through a power amplifier. The load applied is then recovered using the proposed method and the measured displacement. The identified load is compared to the actual load measured as by the force sensor.   Figure 30 shows the measured displacement response of the sinusoidal, as well as the triangular, loads. Using the proposed augmented Tikhonov regularization method, 2 β selection results are shown in Figure 31. The load identification results are shown in Figures 32 and 33, respectively. The results accuracy is evaluated in Table 8.    Figure 30 shows the measured displacement response of the sinusoidal, as well as the triangular, loads. Using the proposed augmented Tikhonov regularization method, β 2 selection results are shown in Figure 31. The load identification results are shown in Figures 32 and 33, respectively. The results accuracy is evaluated in Table 8.
Appl. Sci. 2020, 10, x FOR PEER REVIEW 20 of 24 laser displacement sensor is installed at 0.28 m from the beam end to measure the displacement response. A sampling rate of 4096 is used. Two different loads are applied. First, a sinusoidal and then a triangle wave signal is applied using the exciter. The load magnitude is adjusted through a power amplifier. The load applied is then recovered using the proposed method and the measured displacement. The identified load is compared to the actual load measured as by the force sensor.   Figure 30 shows the measured displacement response of the sinusoidal, as well as the triangular, loads. Using the proposed augmented Tikhonov regularization method, 2 β selection results are shown in Figure 31. The load identification results are shown in Figures 32 and 33, respectively. The results accuracy is evaluated in Table 8.            The results in this section show experimentally that the proposed method can accurately identify the external loads applied on the considered beam. Figure 29 shows that the sinusoidal load is accurately recovered. Both the actual load period and amplitude match those obtained with the proposed approach with high accuracy. The zero initial stage and the initial change in the load amplitude for the first few cycles are all accurately reflected in the load identification. A similar observation can also be made for the triangular wave load identification. The load period, amplitude, and initial changes are all recovered with good accuracy, as can be seen in Figure 30. The correlation coefficients, for both types of considered loads, are above 99%. Similarly, the relative error is around 10% for both loads.

Conclusions
In this paper, we proposed the augmented Tikhonov regularization method in combination with the Green kernel function to identify the external dynamic loads acting on a structure. We also proposed a method for the initial values selection of the regularization parameters. An iterative algorithm for the load identification process was presented. The algorithm performance was evaluated using numerical, as well as experimental, studies. First, we used three numerical examples. The examples showed the efficiency and accuracy of the proposed method for different types of loading, including: an impact load, a sinusoidal wave, or a triangle wave load. The comparisons in these examples show that the proposed method outperforms the standard Tikhonov regularization method for load identification, even when considering noise in the input data. Moreover, the proposed method shows better stability than the standard method, which is attributed to the selection of the regularization parameters. After the numerical examples, we performed an experimental study to further validate the augmented method in a practical application. The method shows high accuracy in capturing different type of loading despite the inaccuracies in the measurement, as well as in the model, parameters. Dynamic load identification, in real conditions, is often affected by constant noise levels, which inevitably appear in the identification results. In the standard Tikhonov method, the noise can result in large fluctuations in the identified load. This is mitigated in the augmented Tikhonov method by a proposed approach to identify the initial parameters based on the noise characteristics. If the noise level does not change, the identification results are stable. Moreover, the paper also proposes a parameters selection approach where the noise level is unknown. Here, the initial parameter β 1 was chosen to be close to zero, while β 2 was chosen based on the slope of the residual norm ||H f − y|| 2 2 . Using this selection process, the values of β 1 and β 2 were chosen without knowing the noise level; hence, the process can be of great significance in practical engineering application. Moreover, the regularization parameters selected by the L-curve and GCV methods are not always optimal, which can lead to large errors. This is especially the case when the input data is noisy. However, the results show that the regularization parameters, selected by the proposed method, can still lead to accurate results even with relatively high level of noise. This paper presents a first attempt to utilize an augmented Tikhonov regularization method combined with the Green kernel function for load identification. Such a problem can be especially challenging due to the fact that they often are ill-posed. The proposed method to identify dynamic load can have a wide range of applications in several domains, such as design optimization, diagnostics, control, and monitoring of vibrating structures.