A new regularization method for dynamic load identification

Abstract

Dynamic forces are very important boundary conditions in practical engineering applications, such as structural strength analysis, health monitoring and fault diagnosis, and vibration isolation. Moreover, there are many applications in which we have found it very difficult to directly obtain the expected dynamic load which acts on a structure. Some traditional indirect inverse analysis techniques are developed for load identification by measured responses. These inverse problems about load identification mentioned above are complex and inherently ill-posed, while regularization methods can deal with this kind of problem. However, most of regularization methods are only limited to solve the pure mathematical numerical examples without application to practical engineering problems, and they should be improved to exclude jamming of noises in engineering. In order to solve these problems, a new regularization method is presented in this article to investigate the minimum of this minimization problem, and applied to reconstructing multi-source dynamic loads on the frame structure of hydrogenerator by its steady-state responses. Numerical simulations of the inverse analysis show that the proposed method is more effective and accurate than the famous Tikhonov regularization method. The proposed regularization method in this article is powerful in solving the dyanmic load identification problems.

Introduction

It is well known that dynamic forces are very important boundary conditions in practical engineering applications, such as structural strength analysis, health monitoring and fault diagnosis, and vibration isolation. However, dynamic load identification theory in practical engineering is always very necessary and important in many cases. If dynamic load is exactly obtained, it is possible to exploit various advanced methods to ensure the safety and the reliability of engineering structures and to satisfy the requirement of modern industry. Moreover, there are many applications in which we have found it very difficult to directly obtain the expected dynamic load which acts on a structure. Therefore, it is very necessary and valuable to develop some indirect inverse analysis techniques for load identification by measured response. ¹

In the subject of structural dynamics, many experts exploited the frequency response functions of the structure to solve the force identification problems in the 1980s.^2–6 Chun et al. ⁷ identified the flight loads of an aircraft, exploiting the singular system theory to deal with the ill-posedness of the transfer function matrix. Thite and Thompson ⁸ presented the minimum condition numbers and sensor placement criterion for solving the force identification problem. Jiang et al. ⁹ eased the ill-posed problem of force identification using a weighted condition number method. Liu and Shepard ¹⁰ adopted the singular value decomposition (SVD) method to locate the small singular value affected by the noisy response. Liu and colleagues^11,12 presented an inverse procedure for identifying both concentrated and extended line load using Green’s function and Heaviside step function in time domain. Lourens et al. ¹³ developed an augmented Kalman filter method for solving inverse problem of load identification in practical engineering structure. Zhang and Ohsaki ¹⁴ exploited the optimal method to identify the member forces in a prestressed pin-jointed structure. Unfortunately, these inverse problems mentioned above are complex and inherently ill-posed. In addition, it is impossible to directly measure distributed dynamic loads acting on the engineering structure. However, we sometimes need to know the expected distributed dynamic loads which act on the engineering structure. Moreover, regularization methods can well deal with these difficulties in practical engineering structure. So far, many regularization methods have been developed to solve these ill-posed phenomena.^15,16 However, most of these methods are only limited to solve the pure mathematical numerical examples without application to practical engineering problems, and they should be improved to exclude jamming of noises in engineering. In order to solve these problems, a new regularization method is presented in this article to investigate the minimum of this minimization problem and applied to reconstructing multi-source dynamic loads on the frame structure of hydrogenerator by its steady-state responses.

This article is organized as follows. The determinate forward problem for a linear elastic structure is briefly introduced in section “Description of the forward problem.” In section “Reconstruction theory of distributed dynamic loads,” the regularization theory for load identification is detailedly analyzed; a new regularization method is established, and its validity and stability are strictly proved. The proposed method is applied to dynamic load identification of the frame structure of hydrogenerator in section “Application,” and we give a conclusion in section “Conclusion.”

Description of the forward problem

In order to assess the proposed method for use in reconstructing the expected unknown multi-source dynamic loads which act on the frame structure of hydrogenerator, it is very necessary for us to learn the following theory for a linear elastic structure.

For a linear and time-invariant dynamic system, we investigate the multi-source dynamic load identification problem. The response at an arbitrary receiving point in a structure can be expressed as a convolution integral of the forcing time-history and the corresponding Green’s kernel in time domain^{11,12,17–19}

$$y(t)=\underset{0}{\overset{t}{\int }}G(t-\tau )p(\tau )d\tau$$ (1)

where $y(t)$ is the measured response. It can be displacement, velocity, acceleration, strain, and so on. $G(t)$ is the corresponding Green’s function. It is the kernel of impulse response. $p(t)$ is the desired unknown dynamic load which needs to be identified.

Exploiting the discretization of this convolution integral, we separate the whole concerned time period into equally spaced intervals and transform equation (1) into a matrix form of

$$Y(t)=G(t)P(t)$$ (2)

or equivalently

$$\left(\begin{array}{c}{y}_1\hfill \\ \multicolumn{1}{c}{y_2}\\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{y_m}\end{array}\right)=\left(\begin{array}{cccc}{g}_1\hfill & \hfill & \hfill & \hfill \\ \multicolumn{1}{c}{g_2}& g_1\hfill & \hfill & \hfill \\ \multicolumn{1}{c}{⋮}& ⋮\hfill & \ddots \hfill & \hfill \\ \multicolumn{1}{c}{g_m}& g_{m-1}\hfill & \cdots \hfill & g_1\hfill \end{array}\right)\left(\begin{array}{c}{p}_1\hfill \\ \multicolumn{1}{c}{p_2}\\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{p_m}\end{array}\right)\mathrm{\Delta }t$$

in which $y_i,g_i,$ and $p_i$ are responses, Green’s function matrix, and input force at time $t=i\mathrm{\Delta }t,$ respectively. Moreover, $\mathrm{\Delta }t$ is the discrete time interval.

When the structure without applied force is static before force is applied, $y_0$ and $g_0$ are equal to zero. All the elements in the upper triangular parts of $G$ are zeros, so they are not shown. The special form of Green’s function matrix shows the special property of the convolution integral.

To reconstruct the time history $p(t)$ , the information of $y(t)$ and $G(t)$ are required. Actually, we can obtain the numerical Green’s function of a structure by finite-element method. The response at a receiving point can be measured. In addition, the problem of identifying the dynamic load $p(t)$ by $y(t)$ and $G(t)$ is usually ill-posed inverse problem. Next, we will establish a new regularization method and apply it to this kind of dynamic load identification problem.

Reconstruction theory of distributed dynamic loads

Usually, we exploit regularization methods to deal with ill-posed inverse problems, using the minimization of the discretization error, the data error and the round-off error. In addition, when we solve inverse problems about load identification in practical engineering structure, we will encounter bad results. Especially, the identified dynamic loads are very sensitive to the error of the output response which causes large error for the true solution in solving inverse problems without any regularization method. Therefore, we cannot deconvolute it by traditional direct inverse matrix method. In general, we adopt regularization methods to seek the expected solution which is numerically stable and efficient solution.

Considering the measured response is noisy, equation (2) can be transformed into the following equation

$$Y_{\mathrm{err}}=G{P}_{\mathrm{true}}+\mathrm{err}$$ (3)

where $\mathrm{err}$ denotes the noisy data of the response in the measurement. It is not difficult to obtain the SVD of $G,$ which is given as

$$G=\mathrm{US}{V}^T=\sum _{i=1}^m{u}_i{\mu }_i{v}_i^T$$ (4)

where $U=(u_1,u_2,\dots ,u_m)$ and $V=(v_1,v_2,\dots ,v_m)$ are, respectively, filled with orthonormal columns and satisfy $U^{T}U=V^{T}V=I_m.$ $S=\mathrm{diag}({\mu }_1,{\mu }_2,\dots ,{\mu }_m)$ has non-negative diagonal elements which appear in non-increasing order. In addition, it is easy to check that this decomposition generates the formulation of expected load

$$\tilde{P}=G^{-1}{Y}_{\mathrm{err}}=\mathrm{VDiag}({\mu }_i^{-1})U^T{Y}_{\mathrm{err}}=P_{\mathrm{true}}+\sum _{i=1}^m{\mu }_i^{-1}(u_i^T\mathrm{err})v_i$$ (5)

From equation (5), we can easily find that system equation (2) is ill-conditioned because the singular value ${\mu }_i$ is much smaller than maximal singular value ${\mu }_1$ during the process of inverse analysis. During the time that singular values decrease to zero, the noisy data of measured response will largely affect the effects of load identification. Therefore, it is necessary to filter small singular values and this can help to overcome this problem. The usual method for solving this difficulty is that ${\mu }_i^{-1}$ is multiplied by a regularization operator $f_{\alpha }({\mu }_i)$ , where $\alpha$ is the regularization parameter. It can make $f_{\alpha }({\mu }_i)/{\mu }_i$ approach to one when ${\mu }_i$ approaches to zero. Then we can have a good identified force $P_{\alpha }$

$$P_{\alpha }=\mathrm{VDiag}\left(\frac{f_{\alpha }({\mu }_i)}{{\mu }_i}\right)U^T{Y}_{\mathrm{err}}=\sum _{i=1}^m\frac{f_{\alpha }({\mu }_i)}{{\mu }_i}(u_i^T{Y}_{\mathrm{err}})v_i$$ (6)

When the regularization operator is given by

$$f_{\alpha }(\mu )=\frac{{\mu }^2}{\alpha +{\mu }^2}$$ (7)

the corresponding regularization method is just the famous Tikhonov regularization method. Its corresponding regularized solution is

$$P_{\alpha }=\sum _{i=1}^m\frac{{\mu }_i(u_i^T{Y}_{\mathrm{err}})}{{\mu }_i^2+\alpha }{v}_i=(G^{T}G+\alpha I{)}^{-1}{G}^T{Y}_{\mathrm{err}}$$ (8)

In fact, many experts have studied ill-posed problems which have attracted much attention in science and engineering.^20,21

Let $X$ and $Y$ be real Hilbert spaces and $G\in L(X,Y)$ that is., $G:X\to Y$ is a bounded linear operator. We consider the following equation

$$\mathrm{Gx}=y,(y\in R(G))$$ (9)

Throughout this article, we assume

• Hypothesis 1 (H1): $y^{\delta }\in X$ is the available noisy data with

$$\parallel y-y^{\delta }\parallel \le \delta$$

where $\delta$ represents the known noise level.

Then, we will need to solve

$$\mathrm{Gx}=y^{\delta }$$ (10)

If $R(G)$ is not closed, the problem equation (9) is ill-posed, then we have to use a regularization method to deal with it. In general terms, regularization is the approximation of an ill-posed problem by a family of neighboring well-posed problems. A regularization method consists of a regularization operator and a parameter choice rule which is convergent in the sense that if the regularization parameter is chosen according to that rule, then the regularized solutions will converge to the true solution in the norm as the noise level approaches zero. One of the most famous regularization methods is Tikhonov regularization method which exploits the following regularization operator

$$q(\alpha ,\mu )=\frac{{\mu }^2}{\alpha +{\mu }^2},\alpha >0,\mu >0$$ (11)

where the regularization parameter and the singular value will always be defined by $\alpha ,\mu ,$ respectively. Then the regularized solution $x_{\alpha }^{\delta }$ can be obtained by

$$x_{\alpha }^{\delta }=(G^{*}G+\alpha I{)}^{-1}{K}^{*}{y}^{\delta }$$ (12)

Theorem

Let $({\mu }_j,x_j,y_j{)}_{j\in N}$ be a singular system for the linear operator $G:X\to Y,$ and let $q:(0,+\infty )\times (0,\parallel G\parallel ]\to R.$ Then $q(\alpha ,\mu )$ is called a regularization operator, and the corresponding regularization method can be given as

$$x_{\alpha }:=R_{\alpha }y=\sum _{j=1}^{\infty }\frac{q(\alpha ,{\mu }_j)}{{\mu }_j}(y,y_j)x_j$$ (13)

$$\underset{\alpha \to 0}{lim}{R}_{\alpha }\mathrm{Gx}=x,x\in X$$ (14)

If the following conditions hold

• 1.  $|q(\alpha ,\mu )|\le 1$ for $\alpha >0$ and $\mu \in (0,\parallel G\parallel ].$

• 2. For any $\alpha >0,$ there exists $c(\alpha )>0$ such that

$$|q(\alpha ,\mu )|\le c(\alpha )\mu ,\mu \in (0,\parallel G\parallel ].$$

• 3.  $\underset{\alpha \to 0}{lim}q(\alpha ,\mu )=1$ for $\mu \in (0,\parallel G\parallel ].$

• Hypothesis 2 (H2): Let $G:X\to Y$ be a bi-univocal compact operator and $y\in R(G).$

We can easily obtain that equation (1) has unique solution $x$ under the condition of (H2). Exploiting the theory of singular system, we have

$$x=\sum _{j=1}^{\infty }\frac{1}{{\mu }_j}(y,y_j)x_j.$$ (15)

Since ${\mu }_j\to 0$ as $j\to +\infty ,$ as well as equations (13) and (15), we can obtain the convergent approximate solution by the attenuation of $q(\alpha ,\mu )$ to $1/\mu ,$ which can be performed by regularization operator. So we can obtain the corresponding regularization method if a proper regularization operator is established.

Next, we will propose a new regularization method based on new regularization operator, and prove its stability and regular property. Herein, we define $q(\alpha ,\mu ):R^{+}\times (0,\parallel G\parallel ]\to R^{+},$ given by

$$q(\alpha ,\mu )=1-e^{-\frac{{\mu }^2}{2{\alpha }^2}}$$ (16)

Actually, it is not difficult to check that: the function $q(\alpha ,\mu )$ defined by equation (16) is a regularization operator. Moreover

$$1-q(\alpha ,\mu )\le \frac{\sqrt{2}\alpha }{e\mu }$$ (17)

First, it follows from

$$0<e^{-\frac{x^2}{2}}\le 1$$

that $0\le q(\alpha ,\mu )<1$ for $\alpha >0$ and $\mu \in (0,\parallel G\parallel ].$

Second, it is also easy to obtain that

$$\underset{\alpha \to 0}{lim}q(\alpha ,\mu )=\underset{\alpha \to 0}{lim}(1-e^{-\frac{{\mu }^2}{2{\alpha }^2}})=1$$

In the following, we will show $1-e^{-\frac{{\mu }^2}{2{\alpha }^2}}\le \mu \alpha .$ In fact, we just need to prove that $1-e^{-\frac{x^2}{2}}<x$ for $x>0.$

First of all, we construct the following function

$$f(x)=e^{\frac{x^2}{2}}(1-x)-1.$$

Since

$$\begin{array}{c}\multicolumn{1}{c}{f\prime (x)=e^{\frac{x^2}{2}}x(1-x)-e^{\frac{x^2}{2}}}\\ =e^{\frac{x^2}{2}}(x-x^2-1)\hfill \\ =e^{\frac{x^2}{2}}[-{(x-\frac{1}{2})}^2-\frac{3}{4}]<0\hfill \end{array}$$

we have that the function $f(x)$ is a monotone decreasing function for $x>0.$ Then we can easily obtain $f(x)<f(0)=0.$ Then we can immediately obtain the assertion. At this time, we can assure that $q(\alpha ,\mu )$ defined by equation (16) is a regularization operator. Next, we will prove the second result.

Due to

$$e^{\frac{x^2}{2}}\ge \frac{\mathrm{ex}}{\sqrt{2}}(x>0)$$

we have

$$1-q(\alpha ,\mu )\le \frac{\sqrt{2}\alpha }{e\mu }$$

then we immediately complete the proof of the two results above.

According to theoretical analysis above, we can obtain the good approximate solution $P_{\alpha }$ of the identified force by the proposed method

$$P_{\alpha }=\sum _{i=1}^m\frac{(1-e^{-\frac{{\mu }^2}{2{\alpha }^2}})(u_i^T{Y}_{\mathrm{err}})}{{\mu }_i}{v}_i=R_{\alpha }{Y}_{\mathrm{err}}$$ (18)

where $R_{\alpha }:Y\to X$ is defined by

$$R_{\alpha }=[I-e^{-\frac{G^{*}G}{2{\alpha }^2}}](G^{*}G{)}^{-\frac{1}{2}},\alpha >0$$

Application

Dynamic load identification of hydrogenerator frame structure as shown in Figure 1 is investigated. Its corresponding practical engineering problem is to reconstruct the vertical loads which act on hydrogenerator frame structure. The material properties of the plate are as follows: $\rho =7.85\times {10}^3\mathrm{kg}/{\mathrm{m}}^3,E=210\mathrm{GPa},\nu =0.3$ . The vertical concentrated loads $F_1$ and $F_2$ act at nodes 167010 and 109870, respectively. The corresponding vertical displacement responses at nodes 166510 and 109710 can be obtained by finite-element method. Its bottom is fixed, and other parts of hydrogenerator frame structure are free. Its corresponding finite-element model is shown in Figure 1. The action point of dynamic load is denoted by the arrow in Figure 1. In this engineering example, the true loads are given as

$$\begin{array}{c}{F}_1(t)=\{\begin{array}{c}\begin{array}{cc}{q}_1\sin (\frac{2\pi t}{t_d}),\hfill & 0\le t\le 2t_d\hfill \end{array}\hfill \\ \multicolumn{1}{c}{0,t<0\mathrm{and}t>2t_d}\end{array}\hfill \\ F_2(t)=\{\begin{array}{c}\begin{array}{cc}4q_{2}t/t_d,\hfill & 0\le t\le t_d/4\hfill \end{array}\hfill \\ \multicolumn{1}{c}{\begin{array}{cc}2q_2-4q_{2}t/t_d,\hfill & t_d/4<t\le 3t_d/4\hfill \end{array}}\\ \multicolumn{1}{c}{\begin{array}{cc}4q_{2}t/t_d-4q_2,\hfill & 3t_d/4<t\le t_d\hfill \end{array}}\\ \multicolumn{1}{c}{\begin{array}{cc}0,\hfill & t>t_d\hfill \end{array}}\end{array}\hfill \end{array}$$

where $t_d$ is the period of sine force, and $q_i(i=1,2)$ is a constant amplitude of the force. When $t_d=0.004\mathrm{s},q_1=100,000\mathrm{N},$ and $q_2=80,000\mathrm{N},$ the sine force and triangle force are shown in Figures 2 and 3.

Figure 1.

The FEM model of the frame structure of hydrogenerator.

Figure 2.

The vertical concentrated sine load acting at node 167010.

Figure 3.

The vertical concentrated triangle load acting at node 109870.

We simulate the experimental data of measured response, using computational numerical solution. It can be easily obtained by the traditional finite-element method. The corresponding vertical displacement responses at nodes 166510 and 109710 are, respectively, shown in Figures 4 and 5. Moreover, in order to simulate the noise-contaminated measurement, we directly add a noise to the computer-generated response and define the noisy response as

$$Y_{\mathrm{err}}=Y_{\mathrm{cal}}+l_{\mathrm{noise}}·\mathrm{std}(Y_{\mathrm{cal}})·\mathrm{rand}(-1,1),$$

in which $Y_{\mathrm{cal}}$ is the computer-generated response; $\mathrm{std}(Y_{\mathrm{cal}})$ is the standard deviation of $Y_{\mathrm{cal}};$ $\mathrm{rand}(-1,1)$ represents the random number between $-1$ and +1; and $l_{\mathrm{noise}}$ is a parameter. It mainly controls the level of the noise contamination.

Figure 4.

The corresponding vertical displacement response at node 166510.

Figure 5.

The corresponding vertical displacement response at node 109710.

Next, the effect of measurement error on the accuracy of estimated values is studied by considering the case of noise level $5\%.$ The proposed regularization method is applied to multi-source dynamic load identification of hydrogenerator frame structure. Engl et al. ²⁰ showed that Tikhonov regularization method can stably and effectively solve ill-posed problems. So the regularized solution obtained by the proposed method will be compared with Tikhonov regularization method. The comparison will be made quantitatively by way of the relative estimation error

$$\tilde{F}=\left|\frac{F_{\mathrm{R}\mathrm{eal}}(i)-F_{\mathrm{I}\mathrm{dentified}}(i)}{max\{F_{\mathrm{j}}\}}\right|*100.$$ (19)

and the average error

$$F_{\mathrm{A}\mathrm{verage}}=1n\sum _{i=1}^n\left|\frac{F_{\mathrm{R}\mathrm{eal}}(i)-F_{\mathrm{I}\mathrm{dentified}}(i)}{max\{F_{\mathrm{j}}\}}\right|*100,$$ (20)

where $i=1,2,\dots ,n,j=1,2.$

To assess the stability and effectiveness of the proposed method and Tikhonov regularization method, we just select five time points. Moreover, we will compare the identified load and the corresponding true force for each point.

Numerical performances of the proposed method are given as follows: The proposed new regularization method (MRO) and the famous Tikhonov regularization method can both provide the stable and effective dynamic loads that act on the hydrogenerator frame structure, which can be shown in Figures 6 and 7. In addition, Table 1 shows more detailed results about the identified loads at five time points. From this table, we can find that most identified deviations by the famous Tikhonov regularization method are not greater than the proposed regularization method under the noise level $\pm 5\%,$ but the average deviations of the proposed method solving the identification of triangle force and sine force are smaller than the latter. This is mainly because of efficient identification of proposed new regularization method. We can also find that most identified deviations by the famous Tikhonov regularization method concentrate in the range of $13.9\%,$ while most identified deviations by the proposed method concentrate in the range of $10.6\%.$ Moreover, for the identification of sine force, the maximal deviation and average deviation by the Tikhonov regularization method are $13.89\%,5.02\%,$ respectively, while for the proposed method, its corresponding maximal and average deviation $10.54\%,4.04\%,$ respectively. In addition, the average deviation and maximal deviation of the identification of triangle force by the proposed new regularization method are $3.18\%,9.98\%,$ respectively, which are both smaller than those of the famous Tikhonov regularization method. From numerical results above, we can find that the proposed new regularization method is a stable, convenient and effective regularization method in solving multi-source dynamic load identification of practical engineering structure and gives satisfactory results.

Figure 6.

The identified sine force at noise level $5\%$ .

Figure 7.

The identified triangle force at noise level $5\%$ .

Table 1.

The identified force at five time points at noise level $5\%$ .

			Tikhonov		Present method
	Time point	Real force	Identified force	Error (%)	Identified force	Error (%)
Sine	0.001	100,000	91,797	8.20	91,110	8.89
Triangle	0.0006	48,000	47,068	1.16	46,108	2.37
Sine	0.003	−100,000	−92,499	7.50	−10,1870	1.87
Triangle	0.001	80,000	76,201	4.75	76,493	4.38
Sine	0.0045	70,711	69,454	1.26	64,579	6.13
Triangle	0.0016	32,000	32,729	0.91	36,629	5.79
Sine	0.0063	−45,399	−50,089	4.69	−49,569	4.17
Triangle	0.0033	−56,000	−53,575	3.03	−56,457	0.57
Sine	0.0073	−89,101	−90,738	1.64	−85,688	3.41
Triangle	0.0038	−16,000	−14,462	1.92	−18,880	3.6
Error (%)			Maximum	Average	Maximum	Average
Sine			13.89	5.02	10.54	4.04
Triangle			11.38	3.44	9.98	3.18

Conclusion

A new regularization method is proposed in this article, which is based on an improved regularization operator. It is applied to multi-source dynamic load identification of hydrogenerator frame structure. It has been found from the numerical simulations that the proposed method has better identified results than the famous Tikhonov regularization method in reconstructing the expected dynamic loads which act on the hydrogenerator frame structure. This shows that the proposed new regularization method can stably and effectively provides better solution in solving the load identification problems of practical engineering.