CN105843780A - Sparse deconvolution method for impact load identification of mechanical structure - Google Patents

Abstract

The invention relates to a sparse deconvolution method for impact load identification of a mechanical structure. The method is used for solving the ill-posed nature of the impact load identification inverse problem and comprises steps as follows: 1) a frequency response function between an impact load acting point and a response measurement point of the mechanical structure is measured with a hammering method, a unit impulse response function is obtained through inverse fast fourier transform, discretization is further performed, and a transfer matrix is obtained; 2) an acceleration signal generated by impact load of the mechanical structure is measured with an acceleration sensor; 3) an L1-norm-based sparse deconvolution convex optimization model for impact load identification is established; 4) the sparse deconvolution optimization model is resolved with a primal-dual interior point method, and a sparse deconvolution solution, namely, a to-be-identified impact load, is obtained. The sparse deconvolution method is suitable for identifying the impact load acting on the mechanical structure. Compared with conventional Tikhonov regularization methods based on an L2 norm, the sparse deconvolution method has the advantages of high identification accuracy, high computation efficiency and high stability.

Description

A kind of sparse solution convolution method of frame for movement shock loading identification

Technical field

The invention belongs to system state machine monitoring field, be specifically related to the sparse of a kind of frame for movement shock loading identification Deconvolution method.

Background technology

Shock loading is a very important dynamic loading of class, particularly in the monitoring structural health conditions of composite.Such as Fan blade for wind-power electricity generation is running and in maintenance process, is inevitably being subjected to dust storm, flying bird, hail, maintenance The impact of the exotics such as instrument, and occurrence frequency is higher, impact injury accumulate can to the integrity of fan blade and Bearing capacity causes potential safety hazard；Act on sudden change air-flow, flying bird, the maintenance tool shock to wing of aircraft wing.In real time Monitor this kind of impact signal to be very important, and it is extremely difficult for directly measuring them.

Shock loading recognition methods is all inside L2 norm space such as classical Tikhonov, truncated singular value decomposition etc. Solving unknown load, relate to matrix inversion or decomposition operation, (transfer matrix dimension is more than not to be suitable for extensive indirect problem 10⁴) computing；In recent years, according to the pattern of shock loading, functional approaching was (such as Db small echo, Chebshev multinomial, three B Batten etc.) it is applied to by numerous scholars in the single-impact load identification of data volume on a small scale.This method is by controlling basic function Number reach regularization and approach the purpose of dynamic loading, maximum predicament of its application is its regularization parameter basic function number Mesh is difficult to determine, too much or very few basic function number all can cause the invalid of result；And around the god of shock loading identification Needing substantial amounts of sample to go to learn all kinds of impact event through network method, this is very limited in practical engineering application.Pass The regularization method of system is difficult to tackle repetitive shocks identification under big data scale, and now transfer matrix conditional number is very big, Numerical stability is excessively poor；Additionally when response noises level is higher, traditional regularization method hydraulic performance decline, even incapability are Power.

Summary of the invention

Based on this, the invention discloses the sparse solution convolution method of a kind of frame for movement shock loading identification, described method Comprise the following steps:

S100, the frequency response function measured between machinery structural impact load application point and frame for movement response measuring point also calculate biography Pass matrix；

S200, frame for movement is applied shock loading measure shock response；

S300, construct the convex optimization of sparse solution convolution of shock loading identification based on L1 norm based on step S100 and S200 Model；

S400, solve the convex Optimized model of sparse solution convolution, it is thus achieved that the sparse solution convolution solution of shock loading.

The present invention compared with prior art has the advantage that

1. it is different from traditional truncated singular value decomposition based on L2 norm, Tikhonov regularization method, based on L1 model The sparse solution convolution method of the shock loading identification of number makes full use of the time domain sparse features of shock loading, greatly inhibits sound Answer noise amplification in the shock loading identified；

2., compared with traditional Tikhonov regularization algorithm, sparse solution convolution iterative algorithm accuracy of identification is high, calculate effect Rate is high and stability is strong；

3. the prim al-dual interior point m ethod used has and is not related to transfer matrix inversion operation and need not clear and definite regularization parameter Advantage, its iterative process is i.e. the process of regularization；

4. the sparse solution convolution model that the present invention is given and corresponding prim al-dual interior point m ethod, in high precision and solve efficiently (transfer matrix dimension is more than 10 to big data scale⁴A shock loading identification difficult problem under).

Accompanying drawing explanation

Fig. 1 is the sparse solution convolution method flow process of a kind of frame for movement shock loading identification in one embodiment of the invention Figure；

Fig. 2 is 300W composite wind turbine blade shock loading identification device schematic diagram in an embodiment；

Fig. 3 (a), 3 (b), 3 (c) are force transducer and the 300W blower fan leaves of acceleration transducer actual measurement in an embodiment Sheet data, wherein, Fig. 3 (a) shock loading signal, the acceleration responsive of Fig. 3 (b) measuring point R1；The acceleration of Fig. 3 (c) measuring point R2 Response；

When Fig. 4 (a), 4 (b) are that in an embodiment, PDIPM identifies fan blade shock loading, the change of duality gap becomes Gesture, wherein, Fig. 4 (a) response signal is R1；Fig. 4 (b) response signal is R2；

Fig. 5 (a), 5 (b) are the sparse solution convolution results of blower fan blade percussion load in an embodiment, wherein, Fig. 5 (a) Response signal is R1；Fig. 5 (b) response signal is R2；

Fig. 6 (a), 6 (b) they are the Tikhonov regularization results of blower fan blade percussion load in an embodiment, wherein, and figure 6 (a) response signal is R1；Fig. 6 (b) response signal is R2.

Detailed description of the invention

The invention will be further described for 1-6 and a specific embodiment below in conjunction with the accompanying drawings, it should be emphasised that, following Illustrate to be merely exemplary, and the application of the present invention does not limit to following example.

In one embodiment, the invention discloses the sparse solution convolution method of a kind of frame for movement shock loading identification, Said method comprising the steps of:

S100, the frequency response function measured between machinery structural impact load application point and frame for movement response measuring point also calculate biography Pass matrix；

S200, frame for movement is applied shock loading measure shock response；

S300, construct the convex optimization of sparse solution convolution of shock loading identification based on L1 norm based on step S100 and S200 Model；

S400, solve the convex Optimized model of sparse solution convolution, it is thus achieved that the sparse solution convolution solution of shock loading.

In the present embodiment, described step S200 can not apply shock loading, and not in any position of frame for movement The measurement point being limited in step S100.

Method described in the present embodiment is used for solving that conventional impact load recognition method precision is low, efficiency is low and poor stability Predicament, use the advanced prim al-dual interior point m ethod in convex optimization field to impact under big data scale with high accuracy and having solved efficiently A load identification difficult problem.The sparse theory based on L1 norm of contemporary scientific circle and engineering circles extensive concern is applied to by described method Load identification field, the interior point method optimized algorithm using optimization field advanced solves sparse solution convolution model, greatly inhibits Response noises amplification in the shock loading identified.

In one embodiment, described step S100 specifically includes following steps:

S1001, the frequency response function measured between machinery structural impact load application point and frame for movement response measuring point；

S1002, described frequency response function is carried out inverse fast Fourier transform obtain unit impulse response function, so discrete Change and obtain transfer matrix.

In the present embodiment, the measuring method of frequency response function mainly includes hammering method and vibrator advocate approach, wherein hammering method The most convenient, in the present embodiment, prioritizing selection hammering method measures frequency response function.

In one embodiment, described step S200 utilize sensor measurement to put on the shock loading institute of frame for movement The shock response produced.

In the present embodiment, the corresponding letter that the dynamic loading using acceleration transducer measurement to act on frame for movement produces Number, it is possible to use speed, displacement or strain transducer to measure vibratory response.

In one embodiment, the convex Optimized model of sparse solution convolution in described step S300 is:

$\begin{array}{cc}\underset{f}{\min imize}& \left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1\end{array};$

Wherein, | | g | |₂Represent the L2 norm of vector, | | g | |₁Representing the L1 norm of vector, f represents the sparse of shock loading Deconvolution solution, λ represents regularization parameter, H systems communicate matrix, and y is shock loading response vector.

In one embodiment, employing prim al-dual interior point m ethod (primal-dual interior point method, PDIPM) the convex Optimized model of sparse solution convolution in solution procedure S400, specifically includes following steps:

S401, combine sparse solution convolution convex Optimized model structure duality gap function and center path object function, and just Beginningization；

S402, the direction of search of calculating center path object function；

S403, determine the iteration step length of center path function；

S404, the barrier parameter updated in center path function；

S405, update current solution based on step S401-404；

S406, setting PDIPM Stopping criteria, if current solution meets Stopping criteria, then terminate iterative process, institute State sparse solution convolution solution f of the current i.e. shock loading of solution；Otherwise, iterative process returns step S401 and continues iterative computation, until Meet PDIPM Stopping criteria.

In the present embodiment, described duality gap function is used as to judge the condition of prim al-dual interior point m ethod iteration ends, center Path object function is the object function that prim al-dual interior point m ethod is to be optimized.

In one embodiment, the duality gap function in described step S401 is:

Wherein, antithesis feasible variable v=2 (Hf-y), wherein, | | g | |₂Represent the L2 norm of vector, | | g | |₁Represent vector L1 norm, f is sparse solution convolution solution, and λ represents regularization parameter, and H is systems communicate matrix, and y is load response vector, G (v) Represent Lagrange dual problem, $\begin{array}{cc}\underset{v}{\max imize}G\left(v\right)=-\frac{1}{4}{v}^{T}v-v^{T}y;subjectto& {\left(A^{T}v\right)}_i\≤\mathrm{\λ}\end{array};$

Center path function is:

$\begin{array}{cc}\underset{(f,u)}{\min imize}& {\mathrm{\φ}}_t(f,u)=\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}\underset{i=1}{\overset{n}{\Σ}}{u}_i+\frac{1}{t}\mathrm{\Φ}(f,u)\end{array}$

Wherein, t ∈ (0, ∞) is barrier parameter, u ∈ RⁿFor Obstacles Constraints variable, u_iFor Obstacles Constraints variable u ∈ RⁿIn I-th element, (f u) is logarithm barrier function to Φ；

Initialization in described step S401 specifically includes: initialization non-negative initial solution f=[0 ...., 0]^T∈Rⁿ, obstacle Bound variable u=[1 ...., 1]^T∈Rⁿ, terminate threshold epsilon ＞ 0, regularization parameter λ=0.02 | | H^T y||_∞, initial obstacle ginseng Number t=1/ λ.

In one embodiment, described step S402 use following formula calculate the direction of search [the Δ f of center path function^T, Δu^T]^T:

${\▿}^2{\mathrm{\φ}}_t(f,u)\left[\begin{array}{c}\mathrm{\Δ}f\\ \mathrm{\Δ}u\end{array}\right]=\▿{\mathrm{\φ}}_t(f,u)$

Wherein,It it is center path function phi_t(f, Hessian matrix u),It it is Center Road Footpath function phi_t(f, gradient u)；

Set iteration step length in described step S403 as s, make s=β^j, wherein j represents the smallest positive integral meeting following formula:

${\mathrm{\φ}}_t(f+{\mathrm{\β}}^j,u+{\mathrm{\β}}^j)\≤{\mathrm{\φ}}_t(f,u)+{\mathrm{\α\β}}^j\▿{\mathrm{\φ}}_t{(f,u)}^T{\[{\mathrm{\Δf}}^T,{\mathrm{\Δu}}^T\]}^T$

Wherein, α ∈ (0,1/2) and β ∈ (0,1) is linear search parameter.

In the present embodiment, more excellent parameter is set to: α=0.01 and β=0.5.

In one embodiment, described step S404 specially utilizes following formula regeneration barrier parameter t:

$t=\left\{\begin{array}{cc}max\left\{min\right\{2\mathrm{\μ}n/\mathrm{\η},\mathrm{\μ}t\},t\},& s\≥s_{min}\\ t,& otherwise\end{array}\right.$

Wherein, n is the dimension of transfer matrix H, constant μ=2, s_min=0.5.

In one embodiment, the current solution in described S405 includes present percussion load sparse solution convolution solution f_newAnd barrier Hinder bound variable currently solves u_new, the most more new formula is:

f_new=f_old+sΔf

u_new=u_old+sΔu

Wherein, f_oldAnd u_oldFor last PDIPM iteration result.

In one embodiment, in described step S406, PDIPM Stopping criteria is:

Wherein, terminate threshold epsilon and represent acceptable error.

In the present embodiment, the value terminating threshold epsilon constantly to be attempted, although the least precision of value is somewhat improved, But the iterative steps more and more calculating time can be caused；On the contrary, ε value is too big, can cause premature end iteration；Typical case takes Value is such as ε=0.01.

In one embodiment, the invention discloses the sparse solution convolution method of a kind of frame for movement shock loading identification, Described method specifically includes following steps:

1) measure frequency response function and calculate transfer matrix.Hammering method is used to measure machinery structural impact load application point and machine Frequency response function H (ω) between tool structural response measuring point, by inverse fast Fourier transform (Inverse Fast Fourier Transform, IFFT) obtain unit impulse response function h (t), and then discretization obtains transfer matrix H, wherein, ω represents round Frequency variable, t express time variable；

2) apply shock loading and measure acceleration shock response, using acceleration transducer to measure by acting on machinery knot The response signal y that the shock loading of structure produces.Note: step 1) and step 2) acceleration transducer can change speed, position into Move and strain-ga(u)ge transducer；

3) consider that shock loading signal is that (within the testing time, impulsive force is only in punching for a sparse sequence in time domain Hitting load phase has bigger numerical value, non-load district to can be considered zero), construct sparse solution convolution based on L1 norm convex optimization mould Type:

$\begin{array}{cc}\underset{f}{\min imize}& \left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1\end{array};---\left(1\right)$

Wherein, | | g | |₂Represent the L2 norm of vector, | | g | |₁Representing the L2 norm of vector, f represents that shock loading, λ represent Regularization parameter；Minimum criteria energetic from regularization model based on L2 norm is different, regularization model based on L1 norm The energy that suppression is understood, in the diffusion of all time points, remains the openness feature treating reconstruction signal.

4) there is explicit solution f=(H with traditional Tikhonov method based on L2 regularization^TH+λI)^-1H^TY is different, and L1 is just Then changing model and do not have explicit solution, utilize prim al-dual interior point m ethod to solve sparse solution convolution Optimized model, it comprises steps that:

Initialize 4.1): the sparse solution convolution model in convolution (1), structure duality gap function:

$\mathrm{\η}=\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1-G\left(v\right)---\left(2\right)$

Wherein, antithesis feasible variable v=2 (Af-y) and Lagrange dual problem G (v):

$\begin{array}{cc}\underset{v}{\max imize}G\left(v\right)=-\frac{1}{4}{v}^{T}v-v^{T}y;subjectto& {\left(A^{T}v\right)}_i\≤\mathrm{\λ}\end{array}---\left(3\right)$

Initialize 4.2): structure center path function:

$\begin{array}{cc}\underset{(f,u)}{\min imize}& {\mathrm{\φ}}_t(f,u)=\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}\underset{i=1}{\overset{n}{\Σ}}{u}_i+\frac{1}{t}\mathrm{\Φ}(f,u)\end{array}---\left(4\right)$

Wherein, t ∈ (0, ∞) is barrier parameter, u ∈ RⁿFor Obstacles Constraints variable and logarithm barrier function Φ (f, u):

$\mathrm{\Φ}(f,u)=-\underset{i=1}{\overset{n}{\Σ}}ln(u_i+f_i)-\underset{i=1}{\overset{n}{\Σ}}ln(u_i-f_i)---\left(5\right)$

Initialize 4.3): non-negative initial solution f=[0 ...., 0]^T∈Rⁿ, Obstacles Constraints variable u=[1 ...., 1]^T∈ Rⁿ, terminate terminate threshold epsilon ＞ 0, regularization parameter λ=0.02 | | H^T y||_∞, initial obstacle parameter t=1/ λ；Wherein, | | | |_∞ Represent infinity norm；

Step 4.1): calculate the direction of search [the Δ f of center path function^T, Δ u^T]^T:

${\▿}^2{\mathrm{\φ}}_t(f,u)\left[\begin{array}{c}\mathrm{\Δ}f\\ \mathrm{\Δ}u\end{array}\right]=\▿{\mathrm{\φ}}_t(f,u)---\left(6\right)$

Wherein,It is the Hessian matrix of center path function,It it is the ladder of center path function Degree.

Step 4.2): calculate iteration step length s, make s=β^j, wherein j is the smallest positive integral meeting following formula:

${\mathrm{\φ}}_t(f+{\mathrm{\β}}^j,u+{\mathrm{\β}}^j)\≤{\mathrm{\φ}}_t(f,u)+{\mathrm{\α\β}}^j\▿{\mathrm{\φ}}_t{(f,u)}^T{\[{\mathrm{\Δf}}^T,{\mathrm{\Δu}}^T\]}^T---\left(7\right)$

Wherein, α ∈ (0,1/2) and β ∈ (0,1) is linear search parameter, and typical value is α=0.01 and β=0.5.

Step 4.3): regeneration barrier parameter:

$t=\left\{\begin{array}{cc}max\left\{min\right\{2\mathrm{\μ}n/\mathrm{\η},\mathrm{\μ}t\},t\},& s\≥s_{min}\\ t,& otherwise\end{array}\right.---\left(8\right)$

Wherein, n is the dimension of transfer matrix H, constant μ=2 and s_min=0.5.

Step 4.4): update and currently solve:

$\begin{array}{c}{f}_{new}=f_{old}+s\mathrm{\Δ}f\\ u_{new}=u_{old}+s\mathrm{\Δ}u\end{array}---\left(9\right)$

Step 4.5): duality gap function and dual objective function be used for PDIPM Stopping criteria:

$\mathrm{\ϵ}\≤\frac{\mathrm{\η}}{G\left(v\right)}---\left(10\right)$

Wherein, terminate threshold epsilon and represent that one can accept error, typical value ε=0.01.

If currently solving f_newMeet above formula Stopping criteria, then terminate iterative process, it is thus achieved that the sparse uncoiling of shock loading Long-pending solution f；Otherwise, iterative process returns step 4.1) continue iterative computation, until meeting above formula.

Fig. 1 is the flow chart of the sparse solution convolution method of a kind of frame for movement shock loading identification that the present invention completes, should The time domain that method makes full use of shock loading is openness, according to the frequency response function measured and acceleration signal, builds shock loading Sparse solution convolution model, the prim al-dual interior point m ethod utilizing optimization field advanced solves, it is achieved the purpose of shock loading identification, specifically Step is as follows:

1) measure frequency response function and calculate transfer matrix.Use hammering method measure frame for movement dynamic loading apply location point with Frequency response function H (ω) between frame for movement response measuring point, by inverse fast Fourier transform (Inverse Fast Fourier Transform, IFFT) obtain unit impulse response function h (t), and then discretization obtains transfer matrix H, wherein, ω represents round Frequency variable, t express time variable；Wherein, described hammering method is Modal Test method of testing commonly used in the art；

11) test schematic diagram is as in figure 2 it is shown, a length of 702mm of composite wind turbine blade, the Breadth Maximum selected are 112mm and thickness are 11mm.In order to simulate the real mounting condition of fan blade, with clamped section of clamping blower fan of cantilever beam bearing Root of blade 0～60mm position, two pieces of models are that PCB 333B50 acceleration transducer is arranged in away from root and end 120 ～160mm position, it is respectively labeled as measuring point R1 and R2.The position of simulation impact force action is at distance root of blade 70%.Adopt Hammer (tup top is embedded with force transducer) into shape with the impulsive force of model PCB 086C02, repeat to tap application point five times, simultaneously by LMS SCADASIII data collection system synchronizing record impulsive force and acceleration signal, five Secondary Shocks load application points to acceleration are surveyed Frequency response function between point is H₁(ω)、H₂(ω)、H₃(ω)、H₄(ω) and H₅(ω), LMS IMPACT module it is calculated it Meansigma methods is H (ω)；

12) sample frequency when measuring system frequency response function is 2048Hz, and the sampling time is 1s, and data length is 2050. Excitation point F is respectively 4.92E+17 and 2.66E+18 (conditional number with the conditional number of response point R1 and the high-dimensional transfer matrix of R2 It is the index weighing matrix Degree of Ill Condition).Understanding, this fan blade shock loading identification problem is Very Ill-conditioned.

2) apply shock loading and measure acceleration shock response, using acceleration transducer to measure by acting on machinery knot The response signal y that the shock loading of structure produces；

21) use impulsive force hammer to tap continuously fan blade to be easily hit position five times, and simultaneously by LMS SCADASIII Data acquisition with the sample frequency of 2048Hz, synchronous recording acceleration signal and shock loading signal (such as Fig. 3 (a), 3 (b) 3 (c) Shown in).Wherein, actual measurement force signal is as the comparison other of shock loading sparse convolution method recognition result.Notice that this step is executed Add application point and the step 1 of impulsive force) measure frequency response function application point consistent, remain constant with brief acceleration position.

22), during impact test, impulsive force hammer continuously taps fan blade point five times, be respectively labeled as F1, F2, F3, F4 and F5, as shown in Fig. 3 (a).In order to verify the high efficiency of sparse solution convolution method of the present invention, four times will simultaneously identified in Fig. 3 (a) Repetitive shocks.The bump persistent period is 14s, and sample frequency is 2050Hz, response data a length of 26638.

23) reconstructed due to bump time history simultaneously, cause the dimension of Solve problems to become the hugest, original The dimension of transmission function is much smaller than the response data length of actual measurement, and load identification governing equation now is defined as:

$\left[\begin{array}{c}y\left(\mathrm{\Δ}t\right)\\ y\left(2\mathrm{\Δ}t\right)\\ M\\ y\left(n\mathrm{\Δ}t\right)\\ y\left(\right(n+1\left)\mathrm{\Δ}t\right)\\ M\\ y\left(m\mathrm{\Δ}t\right)\end{array}\right]=\mathrm{\Δ}t\left[\begin{array}{ccccc}h\left(\mathrm{\Δ}t\right)& 0& L& 0& 0\\ h\left(2\mathrm{\Δ}t\right)& h\left(\mathrm{\Δ}t\right)& L& 0& 0\\ M& M& L& M& M\\ h\left(n\mathrm{\Δ}t\right)& h\left(\right(n-1)\mathrm{\Δ}t)& L& 0& 0\\ 0& h\left(n\mathrm{\Δ}t\right)& L& 0& 0\\ M& M& L& h\left(\mathrm{\Δ}t\right)& 0\\ 0& 0& L& h\left(2\mathrm{\Δ}t\right)& h\left(\mathrm{\Δ}t\right)\end{array}\right]\left[\begin{array}{c}f\left(\mathrm{\Δ}t\right)\\ f\left(2\mathrm{\Δ}t\right)\\ M\\ f\left(n\mathrm{\Δ}t\right)\\ f\left(\right(n+1\left)\mathrm{\Δ}t\right)\\ M\\ f\left(m\mathrm{\Δ}t\right)\end{array}\right]---\left(1\right)$

Wherein: m is the data length of actual measurement response, and n is unit impulse response function, and m ＞ n.After discrete, formula (1) Deconvolution model can again be write as the compact form of matrix-vector:

Y=Hf； (2)

Wherein: new transfer matrix H ∈ R^m×mIt it is a banding Toeplitz matrix blocked.Under sparse framework, transmission Matrix is referred to as observing matrix H.In present case, for the response data dimension m=26638, then observation matrix H of inverting dynamic loading Dimension is 26638 × 26638.

3) consider that shock loading signal is that (within the testing time, impulsive force is only in punching for a sparse sequence in time domain Hitting load phase has bigger numerical value, non-load district to can be considered zero), construct sparse solution convolution based on L1 norm convex optimization mould Type:

$\begin{array}{cc}\underset{f}{\min imize}& \left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1\end{array};---\left(4\right)$

Wherein, | | g | |₂Represent the L2 norm of vector, | | g | |₁Representing the L1 norm of vector, f represents the sparse of shock loading Deconvolution solution, λ represents regularization parameter；

4) there is explicit solution f=(H with traditional Tikhonov method based on L2 regularization^TH+λI)^-1H^TY is different, and L1 is just Then changing model and do not have explicit solution, utilize prim al-dual interior point m ethod to solve sparse solution convolution Optimized model, it comprises steps that:

Initialize 4.1): combine the sparse solution convolution model in claim 1, structure duality gap function:

$\mathrm{\η}=\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1-G\left(v\right)---\left(5\right)$

Wherein, antithesis feasible variable v=2 (Af-y) and Lagrange dual problem G (v):

$\begin{array}{cc}\underset{v}{\max imize}G\left(v\right)=-\frac{1}{4}{v}^{T}v-v^{T}y;subjectto& {\left(A^{T}v\right)}_i\≤\mathrm{\λ}\end{array}---\left(6\right)$

Initialize 4.2): structure center path function:

$\begin{array}{cc}\underset{(f,u)}{\min imize}& {\mathrm{\φ}}_t(f,u)=\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}\underset{i=1}{\overset{n}{\Σ}}{u}_i+\frac{1}{t}\mathrm{\Φ}(f,u)\end{array}---\left(7\right)$

Wherein, t ∈ (0, ∞) becomes barrier parameter, u ∈ RⁿFor Obstacles Constraints variable and logarithm barrier function Φ (f, u):

$\mathrm{\Φ}(f,u)=-\underset{i=1}{\overset{n}{\Σ}}\ln (u_i+f_i)-\underset{i=1}{\overset{n}{\Σ}}\ln (u_i-f_i)---\left(8\right)$

Initialize 4.3): non-negative initial solution f=[0 ...., 0]^T∈Rⁿ, Obstacles Constraints variable u=[1 ...., 1]^T∈ Rⁿ, terminate threshold epsilon ＞ 0, regularization parameter λ=0.02 | | H^Ty||_∞, initial obstacle parameter t=1/ λ；Wherein, | | | |_∞Represent Infinity norm；.

Present case calculates the exploitativeness of high-dimensional indirect problem in view of Tikhonov method, and computing environment is a job Standing, it is configured to 2 Intel Xeon E5580 CPU and 1 48GB internal memory；

Step 4.1): calculate the direction of search [the Δ f of center path function^T, Δ u^T]^T:

${\▿}^2{\mathrm{\φ}}_t(f,u)\left[\begin{array}{c}\mathrm{\Δ}f\\ \mathrm{\Δ}u\end{array}\right]=\▿{\mathrm{\φ}}_t(f,u)---\left(9\right)$

Wherein,It is the Hessian matrix of center path function,It it is the ladder of center path function Degree.

Step 4.2): calculate iteration step length s, make s=β^j, wherein j is the smallest positive integral meeting following formula:

${\mathrm{\φ}}_t(f+{\mathrm{\β}}^j,u+{\mathrm{\β}}^j)\≤{\mathrm{\φ}}_t(f,u)+{\mathrm{\α\β}}^j\▿{\mathrm{\φ}}_t{(f,u)}^T{\[{\mathrm{\Δf}}^T,{\mathrm{\Δu}}^T\]}^T---\left(10\right)$

Wherein, α ∈ (0,1/2) and β ∈ (0,1) is linear search parameter, and typical value is α=0.01 and β=0.5.

Step 4.3): regeneration barrier parameter:

$t=\left\{\begin{array}{cc}max\left\{min\right\{2\mathrm{\μ}n/\mathrm{\η},\mathrm{\μ}t\},t\},& s\≥s_{min}\\ t,& otherwise\end{array}\right.---\left(11\right)$

Wherein, n is the dimension of transfer matrix H, constant μ=2 and s_min=0.5.

Step 4.4): what the current solution of renewal included shock loading currently solves f_newCurrent solution with Obstacles Constraints variable u_new:

$\begin{array}{c}{f}_{new}=f_{old}+s\mathrm{\Δ}f\\ u_{new}=u_{old}+s\mathrm{\Δ}u\end{array}---\left(12\right)$

What current solution included shock loading currently solves f_newWith Obstacles Constraints variable currently solve u_new。

Step 4.5): duality gap function and dual objective function be used for PDIPM Stopping criteria:

$\mathrm{\ϵ}\≤\frac{\mathrm{\η}}{G\left(v\right)}---\left(13\right)$

Wherein, terminate threshold epsilon and represent that one can accept error, typical value ε=0.01.

If currently solving f_newMeet above formula Stopping criteria, then terminate iterative process, it is thus achieved that the sparse uncoiling of shock loading Long-pending solution f；Otherwise, iterative process returns step 4.1) continue iterative computation, until meeting above formula.

5) for quantitative assessment sparse solution convolution algorithm PDIPM and the essence of the identified load of Tikhonov regularization method Degree, definition time domain overall situation relative error and shock loading peak value relative error are respectively:

$\frac{\left|\right|f_{measured}-f_{identified}|{|}_2}{\left|\right|f_{measured}|{|}_2}\×100\%---\left(14\right)$

$\frac{\left|\right|\max \left(f_{measured}\right)-\max \left(f_{identified}\right)|{|}_2}{\left|\right|\max \left(f_{measured}\right)|{|}_2}\×100\%---\left(15\right)$

Wherein, f_measuredAnd f_identifiedIt is shock loading and the application regularization method reconstruct of force transducer actual measurement respectively Shock loading.

51) Fig. 4 (a), 4 (b) are that sparse solution convolution algorithm PDIPM is respectively with responding signal R1 and R2 recognition reaction in shell During the repetitive shocks of structure, duality gap is along with the increase of iterative steps, and rapid decrease the most gradually tends to mild.Terminate Iterative steps is respectively 23 and 26.

52) Fig. 5 (a), 5 (b) are that the bump that sparse solution convolution algorithm PDIPM response signal R1 and R2 reconstructs carries Lotus.Understanding, PDIPM all coincide with actual measurement load height by five Secondary Shocks load of response signal R1 and R2 reconstruct respectively, and Shock duration is interval, and the shock loading reconstructed is the most sparse.Table 1 gives PDIPM response signal R1 and R2 reconstruct Repetitive shocks relative error respectively be only 13.02% and 12.93%.The peak force of five Secondary Shocks load of actual measurement divides Wei 11.19N, 36.24N, 23.88N, 38.90N and 26.00N.Peak by the load (see Fig. 5 (a)) of response signal R1 inverting Value power is respectively 11.43N, 36.81N, 24.61N, 39.22N and 26.20N, relative error is respectively 2.14%, 1.51%, 3.06%, 0.82% and 0.77% (being shown in Table 1).

53) Fig. 6 (a), figure (b) are that the bump that Tikhonov regularization method response signal R1 and R2 reconstructs carries Lotus.Understanding, Tikhonov solution is exaggerated in the non-load district of shock loading, is used to inverting dynamic loading when response point R2 especially Time (see Fig. 6 (b)).Table 1 gives the relative of repetitive shocks of Tikhonov method response signal R1 with R2 reconstruct and misses Difference Gao Da 58.09% and 290.24%.On the contrary, by the prior information that shock loading time domain is openness, PDIPM reconstruct Repetitive shocks is impacting non-load district almost nil (see Fig. 5 (a), 5 (b)).It should be noted that: although present case The relative error of Tikhonov method is relatively big, but its peak value relative error still can accept (less than 10%) even than PDIPM peak value relative error is little.Time needed for reconstructing as repetitive shocks, Tikhonov is time-consumingly far more than PDIPM. Therefore, the sparse solution convolution method of the present invention can high accuracy, high efficiency, stably solve high-dimensional, the load of Very Ill-conditioned Identify indirect problem.

Table 1 compares for fan blade repetitive shocks, PDIPM and Tikhonov recognition result

The foregoing is only presently preferred embodiments of the present invention, not in order to limit the present invention, all essences in the present invention Any amendment, equivalent and the improvement etc. made within god and principle, should be included within the scope of the present invention.

Claims (10)

1. the sparse solution convolution method of a frame for movement shock loading identification, it is characterised in that described method includes following step Rapid:

S100, the frequency response function measured between machinery structural impact load application point and frame for movement response measuring point also calculate transmission square Battle array；

S200, frame for movement is applied shock loading measure shock response；

S300, construct the sparse solution convolution convex optimization mould of shock loading identification based on L1 norm based on step S100 and S200 Type；

S400, solve the convex Optimized model of sparse solution convolution, it is thus achieved that the sparse solution convolution solution of shock loading to be identified.

Method the most according to claim 1, it is characterised in that preferably, described step S100 specifically includes following steps:

S1001, the frequency response function measured between machinery structural impact load application point and frame for movement response measuring point；

S1002, described frequency response function is carried out inverse fast Fourier transform obtain unit impulse response function, and then discretization obtains Obtain transfer matrix.

Method the most according to claim 1, it is characterised in that: described step S200 utilize sensor measurement to put on machine Shock response produced by the shock loading of tool structure.

Method the most according to claim 1, it is characterised in that the convex Optimized model of sparse solution convolution in described step S300 For:

$\begin{array}{cc}\underset{f}{\min imize}& \left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1\end{array};$

Wherein, | | g | |₂Represent the L2 norm of vector, | | g | |₁Representing the L1 norm of vector, f represents shock loading sparse solution convolution Solving, λ represents regularization parameter, and H is transfer matrix, and y is shock loading response vector.

Method the most according to claim 4, it is characterised in that described step S400 utilizes prim al-dual interior point m ethod (PDIPM) Solve the sparse solution convolution convex Optimized model sparse solution convolution solution with acquisition shock loading, and specifically include following steps:

S401, combine sparse solution convolution convex Optimized model structure duality gap function and center path object function, and initialize；

S402, the direction of search of calculating center path object function；

S403, the iteration step length of calculating center path function；

S404, the barrier parameter updated in center path function；

S405, update current solution based on step S401-404 iteration；

S406, set PDIPM Stopping criteria, if current solution meets Stopping criteria, then terminate iterative process, described work as Front solution is the sparse solution convolution solution of shock loading f；Otherwise, iterative process returns step S401 and continues iterative computation, until full Foot PDIPM Stopping criteria.

Method the most according to claim 5, it is characterised in that

Duality gap function in described step S401 is:

Wherein, antithesis feasible variable v=2 (Hf-y), f is shock loading sparse solution convolution solution, and λ represents regularization parameter, and H is for being System transfer matrix, y is shock loading response vector, and G (v) represents Lagrange dual function,

$\begin{array}{ccc}\underset{v}{\max imize}& G\left(v\right)=-\frac{1}{4}{v}^{T}v-v^{T}y;subjectto& {\left(A^{T}v\right)}_i\≤\mathrm{\λ}\end{array};$

Center path object function is:

$\begin{array}{cc}\underset{(f,u)}{\min imize}& {\mathrm{\φ}}_t(f,u)=\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}\underset{i=1}{\overset{n}{\Σ}}{u}_i+\frac{1}{t}\mathrm{\Φ}(f,u)\end{array}$

Wherein, t ∈ (0, ∞) is barrier parameter, u ∈ RⁿFor Obstacles Constraints variable, u_iFor Obstacles Constraints variable u ∈ RⁿIn i-th Individual element, (f u) is logarithm barrier function to Ф；

Initialization in described step S401 specifically includes: initialization shock loading sparse solution convolution solution f=[0 ...., 0]^T∈ Rⁿ, Obstacles Constraints variable u=[1 ...., 1]^T∈Rⁿ, terminate threshold epsilon ＞ 0, regularization parameter λ=0.02 | | H^Ty||_∞, initial Barrier parameter t=1/ λ.

Method the most according to claim 6, it is characterised in that

Described step S402 use following formula calculate the direction of search [the Δ f of center path function^T, Δ u^T]^T:

${\▿}^2{\mathrm{\φ}}_t(f,u)\left[\begin{array}{c}\mathrm{\Δ}f\\ \mathrm{\Δ}u\end{array}\right]=\▿{\mathrm{\φ}}_t(f,u)$

Wherein,It it is center path function phi_t(f, Hessian matrix u),It it is center path function φ_t(f, gradient u)；

Iteration step length in described step S403 is s, makes s=β^j, wherein j represents the smallest positive integral meeting following formula:

${\mathrm{\φ}}_t(f+{\mathrm{\β}}^j,u+{\mathrm{\β}}^j)\≤{\mathrm{\φ}}_t(f,u)+{\mathrm{\α\β}}^j\▿{\mathrm{\φ}}_t{(f,u)}^T{\[{\mathrm{\Δf}}^T,{\mathrm{\Δu}}^T\]}^T$

Wherein, α ∈ (0,1/2) and β ∈ (0,1) is linear search parameter.

Method the most according to claim 7, it is characterised in that

Described step S404 utilizes following formula regeneration barrier parameter t:

$t=\left\{\begin{array}{cc}max\left\{min\right\{2\mathrm{\μ}n/\mathrm{\η},\mathrm{\μ}t\},t\},& s\≥s_{min}\\ t,& otherwise\end{array}\right.$

Wherein, n is the dimension of transfer matrix H, constant μ=2, S_min=0.5.

Method the most according to claim 8, it is characterised in that

Current solution in described S405 includes present percussion load sparse convolution solution f_newWith Obstacles Constraints variable currently solve u_new, Concrete more new formula is:

f_new=f_old+sΔf

u_new=u_old+sΔu

Wherein, f_oldAnd u_oldFor shock loading sparse solution convolution solution in last PDIPM iteration result and Obstacles Constraints variable Result.

Method the most according to claim 9, it is characterised in that

In described step S406, PDIPM Stopping criteria is: