JPH085450A - Vibration characteristic analyzing device - Google Patents

Abstract

PURPOSE:To improve the convergence, and to prevent the quick increase of the computing time for computing reverse matrix by separately computing linear term and non-linear term with partial differential repetitive method, and providing a condition that a residual is not diverged even in the case where the linear term and the non-linear term are changed. CONSTITUTION:A computing device separately computes the linear term and the non-linear term with the partial differential repetitive method that the parameter is repetitively changed till a residual between the real measurement data and a theoretical value by a mode characteristic parameter becomes the minimum. The condition that the residual is not diverged even in the case where the linear term and the non-linear tent are changed is provided so as to identify the mode characteristic. As the condition that the residual is not diverged, partial differential by the linear term and the non-linear term in relation to the residual becomes zero, and at this stage, the secondary partial differential term can be omitted. Furthermore, as the condition that the residual is not diverged, the condition that the correcting vector of the non-linear term at the time of carrying the partial differential repetitive method is reduced to a value less than the threshold value is desirably added.

Description

Detailed Description of the Invention

[0001]

BACKGROUND OF THE INVENTION 1. Field of the Invention The present invention relates to a vibration characteristic analyzing apparatus, and more particularly to a vibration characteristic analysis for obtaining a response by vibrating an object under test such as an engine block and identifying a mode characteristic of the object under test from the response function. It relates to the device.

In order to analyze the dynamic characteristics of an actual mechanical structure, it is necessary to identify the dynamic characteristics of the structure from the vibration response characteristics obtained by subjecting the structure to a vibration test. In this case, the method for identifying the mode characteristics (modal parameters such as natural frequency, modal damping ratio, and eigenmode shape) that is a parameter that directly expresses the vibration phenomenon, is easy to understand the phenomenon directly, and has good convergence is currently used. It is the mainstream.

[0003]

2. Description of the Related Art FIG. 1 schematically shows the structure of a conventionally known vibration characteristic analyzing apparatus.
And an impulse hammer 2 for exciting the DUT 1 are prepared.

A force detecting load cell 3 is attached to the hammer 2, and the three-dimensional (x, y, z) direction acceleration of the vibration force by the hammer 2 is detected at an arbitrary position of the DUT 1. Attach the pickup 4 for.

Then, the output signal of the pickup 4 and the output signal of the load cell 3 are given to the arithmetic unit 5 and FFT processing is performed to perform a frequency response function (transfer function: compliance).
Is calculated and the mode characteristics are calculated and output.

In this case, the vibration test is usually carried out by fixing the vibration position of the DUT 1 by the hammer 2 and moving the pickup 4 sequentially to give a plurality of response functions to the arithmetic unit 5. The DUT 1 is fixed by a soft spring (not shown) or the like.

Many methods have been proposed as an algorithm of the arithmetic unit (experimental mode analysis unit) 5 of such a vibration characteristic analysis device. Among them, the "deviation method" in the frequency domain is accurate. It is currently most widely used due to its goodness.

In simple terms, the deviational iteration method is to identify the most probable mode characteristic by the least squares method by performing iterative calculation while testing the error between the experimental data and the theoretical value. Furthermore, weighting can be freely performed in the frequency domain, and the influence of mode characteristics other than the mode characteristic of interest can be considered.

<Divisional Iterative Method> The deviational iterative method will be described below. First, the response function of the mode characteristic given from the pickup 4 to the arithmetic unit 5 when the DUT 1 is vibrated by the hammer 2 ( Since the theoretical value) has an amplitude and a phase, it is well known that it is represented by a complex number composed of a real part G _R (ω) and an imaginary part G _I (ω) as in the following equation.

[0010]

[Equation 1]

In this case, the angular frequency is represented by ω, and the number of eigenmodes (degrees of freedom) within the angular frequency range of interest is n (corresponding to the number of peaks indicating resonance points in FIG. 7 described later). I am trying.
Note that U _r and V _r represent the real part and imaginary part of the r-th residue with respect to the measurement point.

Then, the sum of squares of the error between this complex number and the experimental data (data with the frequency after FFT calculation processing as the horizontal axis and the response function (compliance = displacement / force) as the vertical axis) obtained from the pickup 4 ( The parameter (constant) in Eq. (1) is determined so that the residual sum of squares is minimized, and the mode characteristic is determined from that parameter.

That is, the most probable parameter {γ} ^T = {..., σ _r , ω _dr , U _r , from the experimental data corresponding to the equation (1).
To obtain V _r , ..., C, D, E, F, ...} ({...} represents a vector), it is assumed that the error of the experimental data with respect to the theoretical value is normally distributed and the residual sum of squares S is If the parameter {γ} is sequentially corrected from the initial value (this is the value obtained by reading the frequency at the resonance point of the experimental data E of FIG. 7) so as to be the minimum, the experimental data E (solid line) is obtained as shown in FIG. The theoretical value T (dashed line) closest to (probably) the frequency response characteristic of
Will be obtained. This residual sum of squares S is expressed by the following equation.

[0014]

[Equation 2]

When simultaneously processing a plurality of response functions, the real part and the imaginary part thereof may be arranged in series, and the number of response functions (the number of measurement points) is m, and one response function ( If the number of data during one excitation is N, p = mN
Becomes W _i is a weight and is proportional to the reciprocal of the variance of each measured data.

In order to obtain the parameter {γ} that minimizes the residual sum of squares S, the denominator of the response function (theoretical expression) f in the equation (2) in the parameter {γ} is _expressed by σ _r , ω _dr, etc. Since the nonlinear term is included, a simplified modified iterative method (Gauss-Newton method) by linear approximation is applied as follows to find a correction vector {Δγ} for correcting the parameter {γ} and leave it. The calculation may be repeated until the sum of squared differences S converges to the minimum value.

[0017]

## EQU3 ## {Δγ} = ([A] ^T [W] [A]) ^-1 · [A] ^T [W] {y _i −f _i } Equation (3) where A _ij = ∂f _i / ∂γ _j , which is the theoretical formula
It shows how much the theoretical equation (1) changes when one of γ in (1) is changed. This shows that the denominator including the nonlinear term is forcibly converted to the numerator. Has achieved "approximation". [W] indicates a diagonal matrix having W _i as a diagonal component.

[0018]

However, the above-mentioned partial iteration method has the following two drawbacks. Especially when there are multiple response functions, the desired parameter {γ}
Also increases as the response function increases, so the calculation time of the inverse matrix of equation (3) increases sharply (see characteristic A in FIG. 3). Since the number of variables increases, iterative calculations are less likely to converge and more likely to diverge, and calculation may become impossible on the way (see FIG. 8).

Therefore, according to the present invention, a force detecting load cell attached to an impulse hammer for vibrating the DUT and a three-dimensional acceleration of the exciting force attached to an arbitrary position of the DUT. And a calculation device for receiving the output signal of the pickup and the output signal of the load cell to obtain a frequency response function and calculating a mode characteristic parameter to identify the mode characteristic of the DUT. It is an object of the present invention to realize an improvement related to the essence of an algorithm in which a vibration characteristic analyzing apparatus is provided, which is easy to diverge and takes a long calculation time.

[0020]

In order to achieve the above-mentioned object, the vibration characteristic analyzing apparatus according to the present invention is arranged such that the arithmetic unit operates until the residual between the measured data and the theoretical value by the mode characteristic parameter becomes minimum. The partial characteristic iterative method in which the parameters are repeatedly changed is divided into a linear term and a nonlinear term for calculation, and the mode characteristic is identified by providing a condition that the residual does not diverge even when the linear term and the nonlinear term change. It is characterized by doing.

Further, in the above computing device, the condition that the divergence does not occur is that the partial differential due to the linear term and the nonlinear term with respect to the residual becomes zero, and the second partial differential term at this time can be omitted. .

Further, it is preferable that the arithmetic unit adds a reduction condition within a threshold value to the correction vector of the nonlinear term when the partial iteration method is executed as the condition that the divergence does not occur.

[0023]

In the vibration characteristic analyzing apparatus according to the present invention, the mode characteristic parameter is divided into the linear term and the non-linear term.

That is, the number of non-linear term parameters is invariable because the position (frequency) of the peak shown in FIG. 7 does not change even if the number of response functions increases, and is not affected by the number of response functions. This is a feature, and if this non-linear term is found, the linear term can be easily found.

Therefore, the initial value (provisional value) of the nonlinear parameter is first set for only the nonlinear term, and the linear term parameter is calculated by this initial value.

The response functions G _R and G _I are calculated by the equation (1) using the linear term thus obtained and the nonlinear term of the initial value.

The response functions G _R , G calculated in this way
Calculate the residual sum of squares S between _I and the experimental value according to equation (2),
When it is determined that the residual sum of squares S is the minimum value (minimum value), it is assumed that the linear term and the nonlinear term obtained as described above are the final values of the mode characteristics to be identified.

When it is determined that the residual sum of squares S is not the minimum value, the residual sum of squares S finds a new nonlinear term in order to minimize the residual sum of squares S.

To obtain only this non-linear term, the correction vector {Δβ} is obtained in the same manner as in the above equation (3), and iterative calculation is performed until it converges. A constraint condition for preventing divergence is provided near the minimum point of the sum of squared differences S, and the modified vector {Δβ} of the nonlinear term is obtained by this constraint condition.

However, although the correction vector {Δβ} can be obtained,
If it is used as it is, it may diverge, so further side restrictions are set.

The correction vector {Δβ} thus obtained and the previous nonlinear term parameter {β} (initial value) is added to calculate the new nonlinear term parameter, and the new nonlinear term is calculated. By repeating the execution using the parameters, {α} and {β} showing the optimum mode characteristics can be determined.

[0032]

DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS The vibration characteristic analyzing apparatus according to the present invention may have the same structure as that of the conventional example shown in FIG. 1, but the difference from the conventional example is the arithmetic processing in the arithmetic unit 5. The calculation process will be described below with reference to the flowchart shown in FIG.

(1) Parameter division: When a plurality of response functions are considered at the same time, if the number of response functions is m and the number of eigenmodes to be identified is n as in the above, the number of parameters {γ} Is 2n + (2n + 4) m. Its parameters {γ}
In equation (1), the nonlinear term in the denominator is 2n, and the linear term in the numerator is (2n + 4) m. This parameter {γ} is divided into a linear term {α} and a non-linear term {β} as in the following equation.

[0034]

Equation 4] ^{{γ} T = {{α} } T, {β} T} {α} T = {‥, U r, V r, ‥, C, D, E, F, ‥} {β} T = {‥, σ _r , ω _dr , ...} Equation (4)

Here, σ _r in the nonlinear term {β} indicates each peak (n peaks shown in FIG. 7 that may exist), that is, the attenuation rate in each mode characteristic, and ω _dr indicates each peak, that is, each mode. The damping natural angular frequency in the characteristic is shown. This damping factor and natural frequency are response points (mounting position of pickup 4)
Does not change even when changes occur. In other words, the number 2n of the parameter {β} of the nonlinear term is invariant because the position (frequency) of the peak shown in FIG. 7 does not change even if the number m of the response functions increases. It turns out that it is not affected by the number of functions. Then, if the nonlinear term is obtained, only the linear term can be easily obtained.

Therefore, it is devised to first obtain only the nonlinear term. Therefore, first, the initial value of the nonlinear parameter {β} is set (step S1 in FIG. 2). This is done by reading the frequency of the peak of the frequency response as shown in FIG. 7, setting the value as ω _dr, and setting an appropriate value for σ _r .

If the non-linear term parameter {β} is initialized, the linear term parameter {α} can be obtained by one calculation by the following equation (step S2).

[0038]

Equation 5] {α} = ([D] T [W] [D]) -1 [D] T [W] {y i} ··· Equation (5) where, D _ij = ∂f _i / ∂ α _j (i = 1 to 2p, j = 1 to (2
n + 4) m)

Here, as shown in equation (4), {α} ^T
= {..., U _r , V _r , ..., C, D, E, F, ...} Therefore, D _ij
Is as follows.

First, when i = 1 to p, f _i = G
It is the _{R (ω i), D ij} = ∂G R (ω i) / ∂α j , and the following equation is obtained for each linear terms.

[0041]

(Equation 6)

Further, a _{r and} b _r in the above equation (6) are also _{expressed by}
As shown in (1), it can be expressed by nonlinear terms σ _r, ω _dr ,
Since the nonlinear terms σ _{r and} ω _dr are initialized as described above, D _ij can be calculated by inserting their initial values and ω _i .

Similarly, when i = p + 1 to 2p, f
_{Since i} = G _I (ω _i ), D _ij = ∂G _I (ω _i ) / ∂α _j , and the following equation is obtained for each linear term.

[0044]

(Equation 7)

Also in this equation, D _ij can be calculated by inserting the initial value of the nonlinear term and ω _i .

Using the linear term {α} and the non-linear term {β} thus obtained, the response functions G _R and G are obtained by the equation (1).
_I is calculated (step S3).

The response functions G _R and G calculated in this way
The residual sum of squares S of _I and the experimental value is calculated according to equation (2) (step S4).

Then, it is determined whether or not this residual sum of squares S is smaller than the previously obtained value (step S
5). In this case, the first residual sum of squares S may be set to an appropriate value.

As a result, when the residual sum of squares S is not smaller than the previous value, the residual sum of squares S has reached the minimum value. Therefore, the linear term {α} obtained as described above and The nonlinear term {β} is output as the final value of the mode characteristic to be identified (step S6), and this routine ends.

When the residual sum of squares S is smaller than the previous value, it means that the residual sum of squares S has not reached the minimum value yet. The non-linear term parameter {β} of is obtained as follows.

That is, in order to obtain only the nonlinear term, the correction vector {Δβ} in the above equation (3) is changed to {γ} in the above equation (3).
It is obtained in the same manner as in step (3), and iterative calculation is performed until it converges.

[0052]

[Formula 8] {Δβ} = ([C] ^T [W] [C]) ^-1 [C] ^T [W] · {y _i −f _i } ... Expression (8) where C _ij = ∂f _i / ∂β _j (i = 1 to 2p, j = 1 to 2
n)

However, in this state, the linear term {α}
Since the above is not taken into consideration, the following is devised in order to consider the linear term {α}.

The following equation holds at the minimum point of the residual sum of squares S.

[0055]

[Equation 9] ∂S / ∂α _j = 0 (j = 1 to (2n + 4) m) ・ ・ ・ Equation (9)

If the following equation is satisfied in this vicinity, divergence can be prevented.

[0057]

[Equation 10] ∂ ² S / ∂α _j ∂β _k = 0 (j = 1 to (2n + 4) m, k = 1 to 2n) Equation (10)

That is, the meaning of this equation (10) is ∂S / ∂α _j
Shows how much the residual sum of squares S changes when the linear term is changed, and shows that ∂S / ∂α _j does not change even if the nonlinear term is changed by partial differentiation of this with ∂β _k. (Because the right side is 0), it is a constraint condition for preventing divergence in consideration of the linear term {α}.

The number of the above equation (10) is 2n × (2n + 4) m, which is too many to eliminate the linear term {α}.
Since the objective is not to satisfy all of (0) but to reduce the variables, assuming that the differential coefficient of β _k of the response function f under the equation (10) is H _k , this differential coefficient H _k is The calculation represented by the formula is performed (step S7).

[0060]

[Equation 11]

To obtain this differential coefficient H _k , ∂α
_It is necessary to calculate _j / ∂β _k (j = 1 to (2n + 4) m, k = 1 to 2n).

Therefore, when the residual sum of squares S of the equation (2) is first substituted into the equation (10), the following equation is obtained.

[0063]

[Equation 12]

Substituting H _k of equation (11) for ∂f _i / ∂β _k of the first term on the left side of equation (12), and further rearranging it and displaying it in a matrix form the following equation.

[0065]

[Equation 13]

Therefore, by modifying the equation (13), {a} =
By setting {b} [B] ⁻¹ , ∂α _j / ∂β _k can be obtained.

[B] is a square matrix of (2n + 4) m × (2n + 4) m, but when B _st belongs to different response functions in B _st , [B] becomes As shown, since B _st = 0, it can be divided into m (2n + 4) × (2n + 4) matrices.

[0068]

[Equation 14]

Therefore, the size of the above-mentioned inverse matrix [B] ^-1 for obtaining {a} is irrelevant to the number m of response functions.

If the differential coefficient H _k of the equation (11) is obtained, the correction vector of the nonlinear term {Δ
β} is obtained (step S8).

[0071]

[Formula 15] {Δβ} = ([C] ^T [W] [C]) ^-1 [C] ^T [W] · {y _i −f _i } Equation (15) where C _ij = _{H j (f = f i)} (i = 1~2p, j = 1~2
n)

<Setting of Side Constraints> The correction vector {Δβ} of the non-linear term is obtained by the equation (15), but it may diverge if it is used as it is regardless of the constraint condition of the above equation (10). Further, set side constraints.

As a side constraint of divergence prevention, there is a method of applying a reduction factor to the correction vector to reduce the size in the same direction, but in this case, if the modes are different, the mutual effects of the nonlinear terms are reduced. Therefore, the constraint of the following equation is set for each nonlinear term (step S9).

[0074]

[Expression 16] | Δω _dr / ω _dr | ≦ K ₁ | Δσ _r / σ _r | ≦ K ₂ , σ _r > 0 ・ ・ ・ Equation (16)

This time, the coefficient K ₁ is 0.05 and the coefficient K ₂ is 0.8. If K ₂ <1 and the initial value of σ _r is positive, σ
_r remains positive after repeated calculations.

The correction vector {Δβ} thus obtained and the previous nonlinear term parameter {β} (initial value) are added to calculate the new nonlinear term parameter {β ^* } this time ( Step S10).

[0077]

## EQU17 ## {β ^* } = {β} + {Δβ} (17)

By using the new non-linear term parameter {β ^* } thus obtained, the above step S3 and thereafter are repeatedly executed, and {α} and {β} showing the optimum mode characteristics are output (step S6) It becomes.

In the present invention, the second expansion of the response function f, which is the second item of the equation (12), can be ignored and the same expansion can be performed.

[0080]

As described above, according to the vibration characteristic analyzing apparatus of the present invention, the deviation which repeatedly changes the parameter until the residual between the measured data and the theoretical value of the mode characteristic parameter is minimized. Since the iterative method is divided into a linear term and a non-linear term for calculation, and the mode characteristics are identified by providing a condition that the residual does not diverge even if the linear term and the non-linear term change, As shown in, the effect over the conventional technique can be obtained. (1) Calculation time: A comparison of the calculation times is shown in FIG. 3 and Table 1 below. The number of identified modes n is 5, and the number of data in one compliance is 573. For single point reference, the computer usage time is shorter than in the conventional example,
In the conventional example, when m of the frequency response function increases, the computer use time increases rapidly. This is because the size of the inverse matrix to be solved is (2n + 4) m. In the improved type, the size of the inverse matrix is constant at 2n, and since the number of calculations using this is only increased, it is increased almost linearly.

[0081]

[Table 1]

(2) Convergence: The compliance created by calculation is used to compare the convergence. In the model, the end of the 1st spring is fixed with 5 mass points and 5 springs. The specifications and undamped natural frequency are shown in Table 2 below. The number of modal characteristics to be identified is 5, and the data is 901 at 1 Hz intervals from 50 to 950 Hz. Mode characteristic damping ratio is 1 for all mode characteristics
%.

[0083]

[Table 2]

In order to compare the convergence of the present invention and the conventional example, the initial value of the nonlinear term parameter ω _dr was changed and repeated calculation was performed. The initial value of the modal damping ratio is 1 which is the correct answer.
/ 10 is 0.1%.

The results are shown in FIGS. 4 to 6.
The vertical axis represents the standardized residual Z, which is expressed by the following equation.

[0086]

(Equation 18)

Of these figures, FIG. 4 corresponds to the case where the second derivative term of the second term on the left side of the equation (12) is omitted, and FIG. 5 corresponds to the case where this second derivative term is not omitted. Then, FIG. 6 shows the characteristics when the second-order differential term is not omitted and the side surface constraint by the equation (16) is taken into consideration. The convergence is improved in this order, and the case of FIG. 6 is the conventional example. From this, it can be seen that it has been greatly improved.

[Brief description of drawings]

FIG. 1 is a block diagram showing a vibration characteristic analyzer common to the present invention and a conventional example.

FIG. 2 is a flowchart showing an arithmetic algorithm stored and executed in the arithmetic unit in the vibration characteristic analyzing apparatus according to the present invention.

FIG. 3 is a graph chart for comparing calculation times of the present invention and a conventional example.

FIG. 4 is a graph showing the convergence (No. 1) of the vibration characteristic analyzer according to the present invention.

FIG. 5 is a graph showing the convergence (No. 2) of the vibration characteristic analysis device according to the present invention.

FIG. 6 is a graph showing the convergence (part 3) of the vibration characteristic analysis device according to the present invention.

FIG. 7 is a graph diagram showing the measurement result of the compliance (response function) and the theoretical value when the DUT is vibrated with respect to the frequency.

FIG. 8 is a graph showing the convergence of the conventional example.

[Explanation of symbols]

 DESCRIPTION OF SYMBOLS 1 DUT 2 Impulse hammer 3 Load cell 4 Accelerometer 5 Calculation device In the figure, the same code | symbol shows the same or corresponding part.

Claims (3)

[Claims] 1. A force detecting load cell attached to an impulse hammer for vibrating a DUT, and a three-dimensional acceleration of the vibrating force attached to an arbitrary position of the DUT. And a calculation device for identifying the mode characteristic of the device under test by calculating the mode characteristic parameter by obtaining the frequency response function by receiving the actual measurement data from the pickup and the output signal of the load cell. In the vibration characteristic analyzer, the arithmetic unit divides the partial iteration method in which the parameter is repeatedly changed until the residual between the measured data and the theoretical value of the parameter is minimized into a linear term and a nonlinear term. The vibration characterized by identifying the modal characteristics by setting the condition that the residual does not diverge even if the linear term and the nonlinear term change Gender analysis device. 2. The arithmetic unit is characterized in that the condition that the divergence does not occur is that the partial derivative by the linear term and the nonlinear term with respect to the residual becomes zero, and the second partial derivative term at this time is omitted. The vibration characteristic analysis device according to claim 1. 3. The arithmetic unit adds a reduction condition within a threshold value to a correction vector of the non-linear term when the partial iteration method is executed as the non-divergence condition. Alternatively, the vibration characteristic analysis device according to item 2.