Load Localization and Reconstruction Using a Variable Separation Method

Load identification is very important in engineering practice. In this paper, a novel method for load reconstruction and localization is proposed. In the traditional load localization method, location information is coupled to the impulse response matrix. +e inversion of the impulse response matrix leads the process of load localization to be time-consuming. So we propose a variable separation method to separate the load location information from the impulse response matrix. An error optimization function of load histories in different modes is employed to determine the true load location. After locating the external load, the load time history can be easily reconstructed by the measurement responses and determinate impulse response matrix. +is method is verified by simulations of a simply supported beam acted by a sine load and an impact separately. An experiment is also carried out to validate the feasibility and accuracy of the proposed method.


Introduction
Knowledge of external loads and load locations is crucial in various fields, such as structural dynamic design, noise reduction, and fault diagnosis.However, it is usually difficult to directly measure the structure external loads and the locations due to some physical or economical limitations.
is is why indirect methods must be developed to identify structural loads by using measured structural responses, such as acceleration, velocity, displacement, and strain.
In recent years, with the deepening of the research on structure dynamics, the techniques of load identification have developed rapidly.Currently, there are two mature methods of load identification: frequency-domain methods and time-domain methods.e frequency-domain method converts the kinetic equation into linear equation in the frequency domain.Bartlett and Flannelly [1] firstly used the frequency-domain methods to identify the hub forces in a helicopter model.Starkey and Merrill [2] used direct inversion of the frequency response function (FRF) to determine the load.In their research, they found the frequency response function (FRF) was ill-posed near the resonance zone, and with the number of the loads increasing, the accuracy of the identification result was reduced.Doyle [3] established a waveguide model to describe the dynamic response and employed the spectral analysis to reconstruct the impacting force acting on a two-dimensional bimaterial beam system.Liu and Shepard [4] used enhanced least squares schemes to reduce the random errors of structural response signals.In their researches, a total least squares scheme was also employed to solve the errors associated with the FRF matrix.e applications of the frequency-domain methods are often restricted since sufficient long data are required to apply the Fourier transformation or other harmonic transformation which is often used in frequencydomain methods.us, the time-domain method of load identification was developed.e research of time-domain methods started relatively late and still has many problems need to be solved.Desanghere and Snoeys [5] established the time-domain method of load identification and firstly used the modal coordinate transformation to identify the external load.Chang and Sun [6] employed a deconvolution method to reconstruct the time history of structure loads.e load identification is an inverse problem, where the structural properties and responses are known while the external loads need to be determined.However, this inverse problem is mathematically ill-posed in most case, which means that the uniqueness and the stability of the solution are lost.A little noise in the measured responses of the structure, which is inevitable, would result in the great change of the estimated load.To solve this ill-posed problem, some regularization techniques are introduced to stabilize the identified load.Jacquelin et al. [7] compared the efficiency of two widely used regularization techniques, the Tikhonov method and the truncated singular value decomposition (TSVD), by identifying an impact acting on a circular Kirchhoff plate.In their research, they also compared the effect of twoparameter choice criterion, L-curve method, and generalized cross validation (GCV).Besides these two methods, there are also some commonly used regularization methods, such as the modified TSVD [8], the damped singular value decomposition [9], and the iterative regularization methods [10].Liu et al. [11,12] introduced several new regularization filter function and confirmed the effectiveness and the accuracy of proposed methods for solving load identification problems.In many researches including those mentioned above, the main purpose is to reconstruct the time history of the excitation and the load location is assumed as known information, which is often not the case and needs to be identified before the reconstruction of load time history.In fact, the process of load identification consists of two separate parts: localization of the external load and reconstruction of the load time history.In order to identify the load location, wavelets were introduced to locate the external loads.Gaul and Hurlebaus [13] employed the wavelet transform to determine the arrival time of impact wave to different sensors and used an optimization method to identify the impact location.As for impact load, there are several localization methods by analyzing wave propagation, such as time difference of arrival (TDOA) [14,15], direction of arrival (DOA) [16,17], and wave's energy loss with traveled distance [18,19].Ginsberg and Fritzen [20] created a sample-force-dictionary as the prior knowledge to transform the impact identification into a sparse recovery task.Li and Lu [21] adopted a complex method to determine the location of the impact and then identify the impact history by using a constrained optimization scheme.
ey [22] used a two-step iterative approach to both localize and reconstruct a single point force acting on a structure.In their research, they proposed the stabilization diagram of identified locations to determine the appropriate regularization range and the true force location.However, these two studies have the same problem: the location information is coupled to the impulse response matrix, which generally needs to be inverted in the process of load identification.us, to identify the load position, a considerable amount of matrix inversion is needed and a lot of time will be consumed.
In this paper, a new localization method is proposed, which separates the location variable from the impulse response matrix.In this method, a modal decomposition method is adopted and only once matrix inversion needs to be done to obtain several-order modal force.e load position is identified by optimizing an error function of the modal force.Numerical simulation and identification test on a simply supported beam structure are implemented to demonstrate the accuracy and effectiveness of the proposed method.is paper is organized as follows: in Section 2, the inverse model for load localization and reconstruction is discussed, with introducing regularization techniques.e applied method for fast localization is also proposed.In Section 3, both numerical simulation and experiment are carried out to validate the proposed method.Finally, concluding remarks are presented in Section 4.

Load Localization and Reconstruction Scheme
2.1.Inverse Model.First, we construct the relationship between a concentrated force and responses of a structure.Without loss of generality, a multiple degree-of-freedom system is considered as the structure model, and the degree of freedom of the system is assumed as N. Generally, the dynamic equilibrium equation of the MDOF system is expressed as where M, C, and K represent the mass, damping, and stiffness matrixes, respectively.u(t) is the response vector of the system, and P(t) is the external load vector.Under the zero initial condition, the relation between the responses and the load can be written as the form of a convolution: where h(t) is the unit impulse response matrix, and it can be expressed through the following modal superposition: where m r , ξ r , and ω r are the rth-order modal mass, damping ratio, and natural circular frequency, respectively.φ r is the rth-order mode of vibration.For a concentrated force, s f is assumed as the load location, and f(t) is the load time history.e response of the ith degree of freedom is obtained: where h r (t) is the rth-order impulse response function: e summation part in equation ( 4) can be represented by 2 Shock and Vibration e subscript of h i,s f (t) denotes the ith DOF unit impulse response with load acting on the location s f .e convolution in equation ( 4) can be taken into a discrete form: where Δt denotes the sampling interval, Q denotes the number of sampling points, and Equation ( 7) can be simplified as For a certain system, the only unknown variables in equation ( 8) are the load location s f and load time history f.To determine the location s f , two or more measurement points are chosen to obtain corresponding load time history f.A minimum optimization problem is established to find the suitable value of location s f , which makes the error function value of several load vectors f minimum.In this process, there are lots of matrix inversion because the location variable s f is included in the matrix H. Generally, the size of the matrix H is large, due to the number of sampling points Q.
us, the process of load localization will cost plenty of time which may cause unknown problems during researches and practical engineering situations.

New Load Localization Method.
To reduce the time of load localization, a new method is considered which separates the load location s f from the matrix H.In this method, the number of matrix inversions is reduced significantly so that the load location can be identified faster.
From equation ( 6), h i,s f (t) is expressed as the form of the modal superposition.By substituting equations ( 6) and ( 8), the following equation can be derived: where Combining the component φ r (s f ) and load time history f, equation ( 9) is turned into where S r � φ r (s f ) • f(r � 0, 1, 2, . . ., N). en, summation in equation (11) can be transformed to a matrix form: where From equation ( 12), we can see that the load location s f is included in the vector S, and the matrix W i only contains the information of measurement point and properties of the system.us, the matrix W i is defined when the measurement point is settled.However, the vector S cannot be calculated by directly inversing equation ( 12) because equation ( 12) is an undetermined system of equation.us, to obtain the vector S uniquely, the modal truncation method is employed and the m modes which contribute the most to the response are chosen.Equation ( 12) is transformed to where In addition, the responses of other measurement points are necessary.Assuming that responses of m points are known, then we have where Equation ( 16) can be inversed directly because W trun is a mQ × mQ matrix.
en S i (i � 1, 2, . . ., m) are calculated, and m load vectors can be obtained: For the true load location s f , the m load vectors us, an optimization function is introduced: From equation (19), the optimization function η(s f ) gets the minimum value when the variable s f is the true load location.
us, the load localization is transformed to a minimum optimization problem.
In this method of load localization, only once matrix inversion is computed, which reduces the operation time greatly and increases the efficiency of load localization.For a certain system, the matrix W trun is determined when the measurement points, the sampling interval, and the number Shock and Vibration of sampling points are all defined.erefore, the load location can be quickly identified, as the responses are measured by the sensors.In addition, the selection of modes is a very important part in this method.e identified location may not be accurate, if a specific mode is not kept.us, the selected modes should reflect the response as far as possible.As the number of selected modes is limited by the number of measurement points, the identified result of a load with a wide range of frequencies would not be effective when the number of measuring points is not large enough.Furthermore, the normalized modal shape values must be measured, which is difficult to realize in the real work.is problem needs to be solved with the help of finite element techniques.
Due to the modal truncation, the load time history which is obtained in this method may not be accurate enough.
us, after identifying the load location, the load time history can be reconstructed by equation (8).

Regularization Method.
Since the matrix W trun is generally ill-posed, the direct inversion of equation ( 16) would not give a stable result.e identified result is sensitive to the errors including signal sampling errors, modal truncation errors, and rounding errors.A little error will cause the direct inverse solution change a great deal.
Considering the errors, the measurement response can be expressed as where δ is the error data of responses and Y δ is the measurement response including errors data.e singular value decomposition of W trun is where U � [u 1 , u 2 , . . ., u mQ ] and V � [v 1 , v 2 , . . ., v mQ ] are matrixes with orthonormal columns and Ω � Diag(σ i )(i � 1, 2, . . ., mQ) has nonnegative diagonal singular values appearing in the nonincreasing order.us, equation (20) gives a formulation of the identified result As shown in equation ( 22), the ill-posed factor of the load identification model is controlled by using the smaller singular value σ i rather than the maximal singular value σ 1 .
With the singular values going down to zero, the error of measurement responses would strongly influence the stability of identified results.To reduce this effect, the regularization method is introduced to filter the small singular values.σ −1 i in equation ( 22) is coupled with a regularization operator p(α, σ i ), and the parameter α can make p(α, σ i ) decay to zero when σ i approaches zero.us, the stable identified result can be obtained: e regularization operator is selected as follows: is regularization method has turned into the famous regularization method called the Tikhonov regularization method.
e regularization parameter α is generally determined with a certain criterion, such as L-curve criterion and GCV criterion.

Numerical Simulation.
In this section, two numerical simulations of a simply supported beam model are implemented to verify the proposed method.In this beam model, the length is 1 m long, and the cross section is 6 cm × 1 cm.e material of the beam has the density 7800 kg/m 3 , Young's module 210 GPa and Poisson's ratio 0.3.e beam is evenly divided into 1000 elements, and there are 1001 nodes in total.e first two-order natural frequencies of the beam are f 1 � 23.5 Hz and f 2 � 94.1 Hz.Two displacement sensors are installed at x 1 � 0.23 m and x 2 � 0.45 m from one simply supported end, and the beam model is excited by a sine load and an impact separately, as shown in Figure 1.

Identification of a Sine Load.
e beam is excited by a load F 1 , which is shown in Figure 2. F 1 is a sine load with frequency f � 60 Hz acting at the point a 1 � 0.35 m from one simply supported end.e measurement responses of the two points x 1 and x 2 are computed through the analytical response expression of simply supported beams, and 5% (5% of the maximum of the measured signals) Gaussian white noise has been added to the responses.e responses of the two points are shown in Figures 3 and 4.
To localize the load, the first two-order modes are adopted.ere is only one matrix inversion needing to be done in the identification procedure.Furthermore, the Tikhonov regularization method is used to overcome the illcondition of the matrix inversion.
e regularization parameter α is determined by the L-curve method, which is shown in Figure 5. e red dot indicates the regularization parameter α � 0.0028.e two vectors S 1 and S 2 are calculated through equation (23).To identify the true load location, the 999 nodes are assumed as the load location, except the two ends.At each assumed load location, the two load vectors f 1 and f 2 can be obtained by using equation (18) and corresponding optimization function value η is computed.Figure 6 shows the curve of η and the load location with minimum value of η is 0.356 m from one simply supported end, which is very close to the true load location a 1 � 0.35 m. f 1 and f 2 in the identified load location are shown in Figure 2. At the same time, the relative 4 Shock and Vibration error and the correlation coefficient of f 1 and f 2 are listed in Table 1.
In Table 1, the identified f 1 and f 2 can be used to describe the external load, but they still have some noticeable oscillatory components.To reconstruct the load history more accurately, the responses of measure points x 1 and x 2 are, respectively, used to identify the load history by regularizing equation ( 8) as the load location has already been determined.e reconstruction of the load history is shown in Figure 7. e relative error and the correlation coefficient of the reconstruction are listed in Table 1.It can be found that the reconstruction by measurement response of point 2 is closest to the real load.is is because the response value and the SNR (signal-to-noise ratio) of point x 2 are larger to these of point x 1 .e final result is considered satisfactory.Shock and Vibration

Identification of an Impact.
In this section, a simulation of the beam excited by using an impact load F 2 is carried out to verify the proposed method.e impact load F 2 is shown in Figure 8, acting at the point a 2 � 0.26 m from one simply supported end.e measurement responses of the two points x 1 and x 2 are computed through the analytical response expression of simply supported beams, and 5% (5% of the maximum of the measured signals) Gaussian white noise has been added to responses.e responses of the two points are shown in Figures 9 and 10.
In the same way, the first two-order modes are used to identify the impact location, and the Tikhonov regularization method is employed to reduce the influence of the errors on the identified result.e regularization parameter α is determined by the L-curve method, which is shown in Figure 11.e red dot indicates the regularization parameter α � 0.0054.
To localize the impact load, the optimization function value η at each assumed load location is computed, and the curve of η is shown in Figure 12. e identified load location

6
Shock and Vibration is 0.259 m from one simply supported end. is identified result is very close to the actual impact location a 2 � 0.26 m. e corresponding f 1 and f 2 in the identified load location are compared in Figure 8.We also use the response of measure point x 2 to reconstruct the impact history, which is shown in Figure 13.e relative error and the correlation coefficient of the reconstruction are listed in Table 2. e reconstruction can describe the external load well, so the identified result is accepted.

Experiment.
For experimental validation of the proposed localization method, an experimental test of a steel simply supported beam was set up.As shown in Figure 14, the steel rectangular beam is simply supported at the two ends and the geometric dimension of this beam is measured.e length, width, and thickness are l � 0.695 m, w � 0.04 m, and h � 0.007 m, respectively.e sensors used in this experiment are two PCB piezoelectric accelerometers and one force transducer of Model Number 208C02.An impact hammer with a rubber head is used to excite the beam, and the sampling frequency is 4096 Hz.
e modal parameters are obtained by the system identification techniques.e experiment values of the first three-order natural frequencies are recorded in Table 3. e experiment values of the first three-order modal damping ratio are ξ 1 � 0.0026, ξ 2 � 0.0013, and ξ 3 � 0.0009.According to the natural frequencies obtained by the experiment, the finite element model of the beam is established and the modal shape values φ can be calculated from this finite element model.e beam is evenly divided into 695 elements, and there are 696 nodes in total.e simulation values of the first three-order natural frequencies are also recorded in Table 3.
Eleven mark points are distributed on the beam, and two accelerometers are installed at points 3 and 6, shown in Figure 14.
e beam is impacted at points A, B, and C separately. e real location of these impacts is shown in Table 4. e measurement of the accelerometers is integrated twice to obtain the displacement responses of the measure points.As shown in Figure 15, an FRF analysis shows that   Shock and Vibration 7 the displacement responses are mainly involved in the first two modes, which are used to identify the impact location.After computing the optimization function value η at each node, the impact locations are identified and shown in Figures 16-18 separately.e details of the identified locations are shown in Table 4.It can be seen that the identified locations coincide with the true impact locations well for all three impact points A, B, and C.After determining the impact location, the force histories can be reconstructed.By using the responses of sensor 2 installed at point 6, the comparison of the measured and identified time histories of the three impacts are shown in Figures 19-21, respectively.e relative error and the correlation coefficient of the reconstructions are listed in Table 5.For all three impact test, the identified results have good consistency.us, the experiment result is satisfactory and verifies the credibility of the proposed method.

Conclusion
In this paper, the problem of load localization and reconstruction from structural responses addressed.In order to improve the problem that a lot of matrix inversions in the identification process consume much operation time, a variable separation method is proposed.is method separates the load location variable from the impulse response matrix which needs to be inversed in the localization process.We take several order modes in which the response is mainly included and identify the stable time history vectors of the corresponding modes.In this step, there is     Shock and Vibration only once matrix inversion and it can be replaced by singular value decomposition method.e regularization method is employed to overcome the ill-posed problem and obtain the stable time history vector.en we compute the force history through dividing the time history vector by corresponding mode shape value.e load location is determined with error function using the minimum optimization method.After identifying the location, this problem is transformed into the classic reconstruction of force history.By using the regularization method and L-curve criterion, the force history is inversely computed.e proposed method is fully demonstrated and verified with simulations of a simply supported beam separately acted by a sine load and an impact.An experiment is also implemented to prove the validity of the method.

Supplementary Materials
e codes and data used in this manuscript are included in the "supplementary information.zip"file.e specific description of the files is listed as follows: Simulation 1: the M file "simulation _sin load.m" is the Matlab code to simulate the identification of a sine load.Simulation 2: the txt files "u1err_sin2d.txt"and "u2err_sin2d.txt"are the responses data of the sine load simulation.Simulation 3: the M file "simulation_impact.m" is the Matlab code to simulate the identification of an impact load.Simulation 4: the txt files "u1err_impact2d.txt"and "u2err_impact2d.txt"are the responses data of the impact load simulation.Experiment 5: the M file "experiment_impactA.m" is the Matlab code to verify the identification of impactA.Experiment 6: the txt files "experiment_impactA_response1.txt" and "exper-iment_impactA_response2.txt" are the measured responses of the experiment of impactA.Experiment 7: the txt file "experiment_impactA_Force.txt" is the measured data of impactA.Experiment 8: the M file "experiment_impactB.m"  is the Matlab code to verify the identification of impactB.Experiment 9: the txt files "experiment_impactB_res-ponse1.txt" and "experiment_impactB_response2.txt" are the measured responses of the experiment of impactB.Experiment 10: the txt file "experiment_impactB_Force.txt" is the measured data of impactB.Experiment 11: the M file "experiment_impactC.m" is the Matlab code to verify the identification of impactC.Experiment 12: the txt files "experiment_impactC_response1.txt" and "exper-iment_impactC_response2.txt" are the measured responses of the experiment of impactC.13. e txt file "exper-iment_impactC_Force.txt" is the measured data of impactC.(Supplementary Materials)

Figure 6 :
Figure 6: e optimization function value of F 1 .

Figure 11 :Figure 12 : e optimization function value of F 2 Figure 13 :
Figure 11: e L-curve plot of F 2 .

Figure 14 :
Figure 14: Experimental model of the simply supported beam.

FrequencyFigure 15 :Figure 17 :
Figure 15: FRF of the displacement responses by impact A.

Figure 16 :Figure 18 :
Figure 16: e optimization function value of impact A.

Figure 21 :
Figure 21: e reconstruction of the impact C.

Table 5 :Figure 19 :Figure 20 :
Figure 19: e reconstruction of the impact A.

Table 1 :
e relative error and the correlation coefficient of the identified F 1 .

Table 2 :
e relative error and the correlation coefficient of the identified F 2 .

Table 4 :
Identified locations of the experiment.

Table 3 :
Natural frequency of the beam (unit: Hz).