CN100510697C - Quantitative diagnosis method for rotor transverse crack by B-spline wavelet on the interval - Google Patents

Abstract

This invention relates to a quantitative diagnosis method for rotor transverse crack by B-spline wavelet on the interval. <0}Wherein, the wavelet unit comprises the effect to rotor natural frequency from structural moment of inertia to improve calculation precision; further, it creates transverse crack rotor interval B-spline wavelet finite element crack quantitative diagnosis model database to quantitative diagnose the relative position and depth of the crack by contour-line method. The corresponding calculation results and experimental research show that this invention needs a few wavelet units, and has high precision and well robustness for wide application.

Description

The interval B-spline wavelet unit is used for the method for rotor transverse crack quantitative Diagnosis

Technical field

The invention belongs to structural crack quantitative Diagnosis field, be specifically related to a kind of interval B-spline wavelet B-spline wavelet on the interval that is used for the rotor transverse crack quantitative Diagnosis, the BSWI unit.

Background technology

Rotor crack is one of the most outstanding hidden danger of rotor operation safety, if untimely detection, diagnosis, eliminating, can upset normal production run, especially modern electromechanical equipment develops towards maximization, high speed and intelligentized direction day by day, and but physical construction develops towards light-duty, exquisite direction, make the rotor crack accident that causes owing to crack propagation in recent years constantly take place, cause great economic loss even casualties.As: in April, 2002, first homemade 600MW shaft system of unit abnormal vibration of No. 3 units in the 3rd generating plant, Harbin is through the visual examination no problem.Yet find that in being sent to producer maintenance there are the crackle of dark 180mm in generator amature body and footstalk knuckle place along 165 ° of scopes of periphery of rotor, and the crackle area reaches 1/3 of axle journal cross section, cause the badly damaged accident of generator amature; In October, 2003 Hebei province's southern power grid 2# of factory turbine rotor generation crackle accident, these accidents have been brought enormous economic loss to commercial production.Therefore, adopt simple technology, determine the crack defect of inside configuration, become the problem that engineering circle extremely is concerned about and is constantly sought.Correctly the position of quantitative Diagnosis crackle and the degree of depth are the important decision foundations of estimating its life-span.Though have a large amount of flaw detection new technologies and new instrument equipment to occur in recent years, with these local defect-detecting equipments large scale structure is carried out complete detection, then to expend great amount of manpower, time and funds.

Because any power system can be regarded the mechanical system of being made up of quality, damping and stiffness matrix as, in case the crackle damage occurs, structural parameters just change thereupon, thereby cause the change of system mode parameter (natural frequency, damping, the vibration shape), so the change of modal parameters can be considered the sign that the structure earlier damage takes place.In the method for many detection crackles, utilize occur on the structure crackle and the damage after, will reduce the local stiffness of structure, thereby change this principle of natural frequency of structure, by the test natural frequency, especially in recent years adopt and measure first three rank natural frequency of structure easily, set up structure tradition finite element model, draw out crackle parameter (two parameters of relative position and relative depth change) in advance to first three rank natural frequency influence curve of structure, utilize contouring method, relative position and relative depth that quantitative Diagnosis goes out the structural crack existence are the methods of using always.Yet these methods are based on traditional finite element diagnosis crackle, and problem such as have that robustness is strong, counting yield and precision are not high must adopt a large amount of higher-dimension unit for obtaining higher diagnostic accuracy.The wavelet finite element method is a kind of new numerical analysis method that developed recently gets up, substitute traditional polynomial expression as approximating function with scaling function or wavelet function, utilize the characteristic of wavelet multiresolution, can obtain to be used for the multiple basis function of structure analysis, at the accuracy requirement of finding the solution problem, adopt different basis functions.Spline wavelets is to be defined on the whole real number axis R or the square integrable real number space L of one-period function ²(R) the complete base on when finding the solution boundary value problem with it, numerical oscillation can occur on the border.Interval B-spline wavelet has good local character in spatial domain, when finding the solution boundary value problem, can overcome this defective of numerical oscillation on the border.

At present, during both at home and abroad to static rotor structure dynamic analysis, generally with rotor equivalent for not considering traditional free beam of cross section rotary inertia influence, often can not get the result who conforms to actual measured value.

Summary of the invention

The object of the present invention is to provide a kind of interval B-spline wavelet unit to be used for the method for rotor transverse crack quantitative Diagnosis, this method can solve the horizontal fatigue crack quantitative Diagnosis of the rotor problem that extensively exists in the engineering efficiently, reliably, can go out the position and the degree of depth that the horizontal fatigue crack of rotor-support-foundation system exists by simple exciting experiment quantitative Diagnosis.

Technical scheme of the present invention is to solve like this:

1) node of configuration interval B spline wavelet, structure comprises the interval B-spline wavelet unit of structure moment of inertia;

Rotor transverse crack interval B-spline wavelet finite element model is set up in the unit of being constructed 2) employing 1), obtains crackle quantitative Diagnosis database, utilizes the contouring method quantitative Diagnosis to go out relative position and relative depth that the horizontal fatigue crack of rotor exists.

The node of described configuration interval B spline wavelet, structure comprises the interval B-spline wavelet unit of structure moment of inertia, may further comprise the steps:

The employing exponent number is m, and yardstick is the interval B-spline wavelet scaling function of j, is designated as BSWIm _jScaling function, as interpolation basis function tectonic element, near configuration node singular point, unit are divided into the n=2 of unequal interval ^j+ m-4 section, cell node number are 2 ^j+ m-3, total number of degrees of freedom, is 2 ^j+ m-1, boundary node are 1, n+1; Internal node is 2,3, Λ, and n can obtain

w(ξ)＝N ^ew ^e

Shape function N ^eFor

$N^e=\mathrm{\&Phi;}{T}_b^e$

In the formula For

$T_b^e={\left\{{\left[{\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_1\right)\frac{1}{l_e}\frac{\mathrm{d\&Phi;}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_1}{\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_2\right)K{\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_n\right){\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_{n+1}\right)\frac{1}{l_e}\frac{\mathrm{d\&Phi;}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_{n+1}}\right]}^{-1}\right\}}^T$

$\mathrm{\&Phi;}=\left\{{\mathrm{\&phi;}}_{m,-m+1}^j\left(\mathrm{\&xi;}\right){\mathrm{\&phi;}}_{m,-m+2}^j\left(\mathrm{\&xi;}\right)K{\mathrm{\&phi;}}_{m,2^j-1}^j\left(\mathrm{\&xi;}\right)\right\}$ Expression BSWIm _jThe capable vector of scaling function;

w _eFor physics degree of freedom column vector, be expressed as

w ^e＝{w ₁?θ ₁?w ₂Λw _nw _n+1θ _n+1} ^T

In the formula, ${\mathrm{\&theta;}}_1=\frac{1}{l_e}\frac{\mathrm{dw}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_1}$ With ${\mathrm{\&theta;}}_{n+1}=\frac{1}{l_e}\frac{\mathrm{dw}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_{n+1}}$ Corner on the expression elementary boundary node, l _eBe element length;

The unit Vibration Frequency Equations is

$|{\mathrm{\&omega;}}^2(M_r^e+M_b^e)-K_b^e|=0$

In the formula, ω is circular frequency (rad/s), makes vibration frequency f=ω/2 π (Hz);

Unit rotational inertia consistent Mass Matrix

For

$M_r^e=\frac{\mathrm{\&rho;I}}{l_e}{\left(T_b^e\right)}^T{\mathrm{\&Gamma;}}^{\mathrm{1,1}}{T}_b^e$

The unit bending consistent Mass Matrix

For

$M_b^e=\mathrm{\&rho;}{\mathrm{Al}}_e{\left(T_b^e\right)}^T{\mathrm{\&Gamma;}}^{\mathrm{0,0}}{T}_b^e$

The unit bending stiffness matrix

For

$K_b^e=\frac{\mathrm{EI}}{l_e^3}{\left(T_b^e\right)}^T{\mathrm{\&Gamma;}}^{\mathrm{2,2}}{T}_b^e$

More than various in, E is an elastic modulus, I is a cross sectional moment of inertia, ρ is a density of material, A is an area of section, l _eBe element length,

For be transformed into the cell translation matrix the finite element space, Γ from wavelet space ^2,2, Γ ^0,0And Γ ^1,1Be respectively interval B-spline wavelet scaling function respective derivative integration matrix on the same yardstick.

The unit of being constructed described employing 1), set up rotor transverse crack interval B-spline wavelet finite element model, obtain crackle quantitative Diagnosis database, utilize the contouring method quantitative Diagnosis to go out relative position and relative depth that the horizontal fatigue crack of rotor exists, may further comprise the steps:

Utilize fracturing mechanics knowledge, can obtain the relevant torsion line rigidity k of crackle relative depth α _tAnd corresponding crackle stiffness matrix K _s, promptly

$k_t=\frac{\mathrm{\&pi;E}{r}^8}{32(1-\mathrm{\&mu;})}\&times;\frac{1}{{\&Integral;}_{-r\sqrt{1-{(1-2\mathrm{\&alpha;})}^2}}^{r\sqrt{1-{(1-2\mathrm{\&alpha;})}^2}}(r^2-{\mathrm{\&xi;}}^2)\left[{\&Integral;}_0^{a\left(\mathrm{\&xi;}\right)}\mathrm{\&eta;}{F}^2(\mathrm{\&eta;}/H)\mathrm{d\&eta;}\right]\mathrm{d\&xi;}}$

And

$H_s=\left[\begin{array}{cc}{k}_t& -k_t\\ -k_t& k_t\end{array}\right]$

With crackle stiffness matrix K _sAdd in the overall finite element matrix, obtain the overall vibration frequency equation of rotor transverse crack system, i.e. rotor transverse crack interval B-spline wavelet finite element model

|ω ²(M _r+M _b)-K _b|＝0

In the formula, K _bBe stack of interval B-spline wavelet unit and adding crackle stiffness matrix K _sThe global stiffness matrix that the back forms, M _r+ M _bBe interval B-spline wavelet unit rotational inertia consistent Mass Matrix

With the unit bending consistent Mass Matrix

The overall consistent Mass Matrix that the stack back forms;

Find the solution rotor transverse crack interval B-spline wavelet finite element model by the eigenwert method for solving, obtain crackle quantitative Diagnosis database, utilize the contouring method quantitative Diagnosis in the quantitative identification of crackle to go out relative position and the relative depth that the horizontal fatigue crack of rotor exists.

The present invention influences owing to the BSWI unit of being constructed comprises the structure moment of inertia, and is applied in the rotor transverse crack quantitative Diagnosis, has the following significant advantage that is different from traditional finite element method diagnosis rotor transverse crack:

1) comprise that the BSWI unit of rotor moment of inertia influence is suitable finds the solution singularity problem and structure is had excellent adaptability, quantitatively diagnosing with degree of depth high precision for rotor horizontal fatigue crack position provides a kind of new unit;

2) the present invention can be good with calculation cost acquisition seldom very high computational accuracy, robustness, improved the reliability and the adaptability of rotor transverse crack quantitative Diagnosis;

3) the present invention can be according to the rotor-support-foundation system of different structure, the crackle quantitative Diagnosis that is consistent with actual test result database is provided in advance, the horizontal natural frequency in first three rank of rotor structure that only needs input measurement to obtain in the practical application, but just quantitative Diagnosis goes out relative position and the relative depth that crackle exists.Therefore, can in engineering practice, extensively promote the use of.

Description of drawings

Fig. 1 is BSWIm of the present invention _jWavelet unit is separated territory Ω _eLast node and degree of freedom Pareto diagram;

Fig. 2 is a transverse cracking rotor sketch of the present invention;

Fig. 2 (a) is a rotor-support-foundation system model of the present invention;

Fig. 2 (b) is a rotor transverse crack section of the present invention;

Fig. 3 is rotor transverse crack quantitative Diagnosis BSWI finite element model solving result figure of the present invention;

The f that Fig. 3 (a) expression surface fitting obtains ₁=F ₁(α, β) figure;

The f that Fig. 3 (b) expression surface fitting obtains ₂=F ₂(α, β) figure;

The f that Fig. 3 (c) expression surface fitting obtains ₃=F ₃(α, β) figure;

Fig. 4 tests the horizontal fatigue crack quantitative Diagnosis of rotor figure as a result for the present invention;

Fig. 4 (a) expression crackle operating mode is β=0.78889, α=0.3 o'clock crackle quantitative Diagnosis result, and among the figure, the horizontal ordinate that the A point is corresponding is crackle relative depth α, ordinate is crackle relative position β;

Fig. 4 (b) expression crackle operating mode is β=0.78889, α=0.4 o'clock crackle quantitative Diagnosis result.

Embodiment

Accompanying drawing is specific embodiments of the invention;

Below in conjunction with accompanying drawing content of the present invention is described in further detail:

1) node of configuration interval B spline wavelet, structure comprises the interval B-spline wavelet unit of structure moment of inertia;

Rotor transverse crack interval B-spline wavelet finite element model is set up in the unit of being constructed 2) employing 1), obtains crackle quantitative Diagnosis database, utilizes the contouring method quantitative Diagnosis to go out relative position and relative depth that the horizontal fatigue crack of rotor exists.

With reference to shown in Figure 1, the employing exponent number is m, and yardstick is the interval B-spline wavelet scaling function of j, notes by abridging to be BSWIm _jScaling function, as interpolation basis function tectonic element, because the interval B-spline wavelet scaling function is complete base in the interval, therefore, on the reversible basis of the cell translation matrix that guarantees to be transformed into physical space from wavelet space, more node can be near singular point, disposed, thereby singularity problems such as gradient is big, displacement is discontinuous can be described.BSWIm _jWavelet unit is separated territory Ω _eLast node and degree of freedom are arranged, and the unit is divided into the n=2 of unequal interval ^j+ m-4 section, cell node number are 2 ^j+ m-3, total number of degrees of freedom, is 2 ^j+ m-1.Boundary node is 1, n+1; Internal node is 2,3, Λ, n.

Unknown field function w (ξ) can be expressed as

$w\left(\mathrm{\&xi;}\right)=\underset{k=-\mathrm{<figure-callout\; id="1"\; label="m\; +"\; filenames="C200610042652E00171.png,C200610042652E00181.png"\; state="\{\{state\}\}">m</figure-callout>}+1}{\overset{2^j-1}{\mathrm{\&Sigma;}}}{a}_{m,k}^j{\mathrm{\&phi;}}_{m,k}^j\left(\mathrm{\&xi;}\right)=\mathrm{\&Phi;}{a}^e---\left(1\right)$

In the formula

$a^e={\left\{\begin{array}{ccc}{a}_{m,-m+1}^j& a_{m,-m+2}^{j}K& a_{m,2^j-1}^j\end{array}\right\}}^T$ Expression small echo interpolation coefficient column vector;

$\mathrm{\&Phi;}={\left\{\begin{array}{ccc}{\mathrm{\&phi;}}_{m,-m+1}^j\left(\mathrm{\&xi;}\right)& {\mathrm{\&phi;}}_{m,-m+2}^j\left(\mathrm{\&xi;}\right)K& {\mathrm{\&phi;}}_{m,2^j-1}^j\left(\mathrm{\&xi;}\right)\end{array}\right\}}^{}$ Expression BSWIm _jThe capable vector of scaling function.

Definition physics degree of freedom column vector is

w _e＝{w ₁?θ ₁?w ₂Λw _n?w _n+1?θ _n+1} ^T (2)

In the formula, ${\mathrm{\&theta;}}_1=\frac{1}{l_e}\frac{\mathrm{dw}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_1}$ With ${\mathrm{\&theta;}}_{n+1}=\frac{1}{l_e}\frac{\mathrm{dw}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_{n+1}}$ Represent the corner on the elementary boundary node respectively, l _e=x _N+1-x ₁Be element length, ${\mathrm{\&xi;}}_i=\frac{x_i-x_1}{l_e},$ i＝1，2，Λ，n+1。

W (ξ with different nodes in the formula (1) _i) difference substitution formula (2),

$a^e=T_b^e{w}^e---\left(3\right)$

Matrix in the formula

For

$T_b^e={\left\{{\left[{\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_1\right)\frac{1}{l_e}\frac{\mathrm{d\&Phi;}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_1}{\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_2\right)K{\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_n\right){\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_{n+1}\right)\frac{1}{l_e}\frac{\mathrm{d\&Phi;}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_{n+1}}\right]}^{-1}\right\}}^T---\left(4\right)$

With formula (3) substitution formula (1),

$w\left(\mathrm{\&xi;}\right)=\mathrm{\&Phi;}{T}_b^e{w}^e=N^e{w}^e---\left(5\right)$

In the formula, shape function N ^eFor

$N^e=\mathrm{\&Phi;}{T}_b^e---\left(6\right)$

It is pointed out that and guaranteeing matrix

Under nonsingular prerequisite, the unit internal node can dispose arbitrarily.Therefore, can near singularity point (as crackle), dispose more internal node, thereby can obtain from integral body to the part all solving result preferably with less calculation cost.

When considering that moment of inertia influences, beam deflection problem potential energy functional is

${\mathrm{\&Pi;}}_p={\&Integral;}_a^b\frac{\mathrm{EI}}{2}{\left(\frac{d^{2}w}{{\mathrm{dx}}^2}\right)}^2-{\&Integral;}_a^b\frac{\mathrm{\&rho;A}{\mathrm{\&omega;}}^2}{2}{\mathrm{<figure-callout\; id="2"\; label="w"\; filenames="C200610042652E00171.png"\; state="\{\{state\}\}">w</figure-callout>}}^2\mathrm{dx}-{\&Integral;}_a^b\frac{\mathrm{\&rho;I}{\mathrm{\&omega;}}^2}{2}{\left(\frac{\mathrm{dw}}{\mathrm{dx}}\right)}^2\mathrm{dx}---\left(7\right)$

In the formula, ρ is a density of material, and ω is the vibration circular frequency, and A is the beam cross-sectional area, and E is an elastic modulus, and I is a cross sectional moment of inertia, l _e=b-a is an element length.

Make Unknown Displacement field function w BSWIm _jThe scaling function interpolation representation, then

$w=\mathrm{\&Phi;}{T}_b^e{w}^e---\left(8\right)$

To potential energy functional formula (7), at first the unit is found the solution territory Ω _eBe mapped to the unit standard and find the solution territory Ω _s,, and, make δ Π by variational principle then with unknown field functional expression (8) the substitution formula (7) of beam _p=0, can obtain unit mode of oscillation equation

In the formula, ω is circular frequency (rad/s), makes vibration frequency f=ω/2 π (Hz).

The unit Vibration Frequency Equations that comprises the moment of inertia influence that formula (9) is corresponding is

$|{\mathrm{\&omega;}}^2(M_r^e+M_b^e)-K_b^e|=0---\left(10\right)$

In the formula, unit rotational inertia consistent Mass Matrix

For

$M_r^e=\frac{\mathrm{\&rho;I}}{l_e}{\left(T_b^e\right)}^T{\mathrm{\&Gamma;}}^{\mathrm{1,1}}{T}_b^e---\left(11\right)$

The unit bending consistent Mass Matrix

For

$M_b^e=\mathrm{\&rho;}{\mathrm{Al}}_e{\left(T_b^e\right)}^T{\mathrm{\&Gamma;}}^{\mathrm{0,0}}{T}_b^e---\left(12\right)$

The unit bending stiffness matrix

For

$K_b^e=\frac{\mathrm{EI}}{l_e^3}{\left(T_b^e\right)}^T{\mathrm{\&Gamma;}}^{\mathrm{2,2}}{T}_b^e---\left(12\right)$

More than various in, E is an elastic modulus, I is a cross sectional moment of inertia, ρ is a density of material, A is an area of section, l _eBe element length,

For be transformed into the cell matrix the finite element space from wavelet space.Γ ^2,2, Γ ^0,0And Γ ^1,1Be respectively interval B-spline wavelet scaling function respective derivative integration matrix on the same yardstick, be expressed as respectively

${\mathrm{\&Gamma;}}^{\mathrm{2,2}}={\&Integral;}_0^1\frac{{\mathrm{<figure-callout\; id="2"\; label="d"\; filenames="C200610042652E00171.png"\; state="\{\{state\}\}">d</figure-callout>}}^2{\mathrm{\&Phi;}}^T}{d{\mathrm{\&xi;}}^2}\frac{{\mathrm{<figure-callout\; id="2"\; label="d"\; filenames="C200610042652E00171.png"\; state="\{\{state\}\}">d</figure-callout>}}^2\mathrm{\&Phi;}}{d{\mathrm{\&xi;}}^2}\mathrm{d\&xi;}---\left(14\right)$

${\mathrm{\&Gamma;}}^{1,1}={\&Integral;}_0^1\frac{{\mathrm{d\&Phi;}}^T}{d\mathrm{\&xi;}}\frac{d\mathrm{\&Phi;}}{d\mathrm{\&xi;}}\mathrm{d\&xi;}---\left(15\right)$

${\mathrm{\&Gamma;}}^{\mathrm{0,0}}={\&Integral;}_0^1{\mathrm{\&Phi;}}^T\mathrm{\&Phi;d\&xi;}---\left(16\right)$

With reference to shown in Figure 2, be typical transverse cracking rotor sketch.For determining relational expression f _j=F _j(α, β), (j=1,2,3), wherein, crackle relative depth α=δ/d ₁, crackle relative position β=e/L ₂At first determine the torsion line rigidity k relevant with crackle relative depth α _tAnd corresponding crackle stiffness matrix K _sFor

$H_s=\left[\begin{array}{cc}{k}_t& -k_t\\ -k_t& k_t\end{array}\right]---\left(17\right)$

In the formula, k _tCan calculate by following formula

$k_t=\frac{\mathrm{\&pi;E}{r}^8}{32(1-\mathrm{\&mu;})}\&times;\frac{1}{{\&Integral;}_{-r\sqrt{1-{(1-2\mathrm{\&alpha;})}^2}}^{r\sqrt{1-{(1-2\mathrm{\&alpha;})}^2}}(r^2-{\mathrm{\&xi;}}^2)\left[{\&Integral;}_0^{a\left(\mathrm{\&xi;}\right)}\mathrm{\&eta;}{F}^2(\mathrm{\&eta;}/H)\mathrm{d\&eta;}\right]\mathrm{d\&xi;}}---\left(18\right)$

In the formula

$a\left(\mathrm{\&xi;}\right)=2\mathrm{r\&alpha;}-(r-\sqrt{r^2-{\mathrm{\&xi;}}^2})---\left(19\right)$

F(η/H)＝1.122-1.40(η/H)+7.33(η/H) ²-13.08(η/H) ³+14.0(η/H) ⁴?(20)

$H=2\sqrt{r^2-{\mathrm{\&xi;}}^2}---\left(21\right)$

More than various in, r=d ₁/ 2 is the cracked rotor radius, and E is an elastic modulus, and μ is a Poisson's Ratio.

Then with crackle stiffness matrix K _sAdd in the overall finite element matrix, it adds the implantation site by crackle relative position β decision, obtains the overall vibration frequency equation of implicit crackle relative position β and relative depth α

|ω ²(M _r+M _b)-K _b|＝0 (22)

In the formula, K _bBe stack of interval B-spline wavelet unit and adding crackle stiffness matrix K _sThe global stiffness matrix that the back forms, M _r+ M _bBe interval B-spline wavelet unit rotational inertia consistent Mass Matrix

With the unit bending consistent Mass Matrix

The overall consistent Mass Matrix that the stack back forms.

Under the prerequisite of given different crackle relative position β and relative depth α, find the solution the overall vibration frequency equation formula (22) relevant with different beta and α, can obtain the corresponding relation formula of crackle relative position β and relative depth α and first three horizontal natural frequency in rank

f _j＝F _j(α，β)，(j＝1，2，3) (23)

Because funtcional relationship F _jTherefore the unknown, can be obtained by surface fitting technology by the discrete value that calculates, and is rotor crack quantitative Diagnosis database.

In order from formula (23), to pass through known f _jSolve α and β, promptly determine relational expression $(\mathrm{\&alpha;},\mathrm{\&beta;})=F_j^{-1}\left(f_j\right),(j=\mathrm{1,2,3}),$ At first employing power hammer knocks exciting, picks up impulse response signal with the Polytec laser vibration measurer, by to the response signal spectrum analysis, obtains the horizontal natural frequency f in first three rank of rotor ₁, f ₂And f ₃

Then, the horizontal natural frequency in first three rank of rotor that test obtains is distinguished substitution rotor crack quantitative Diagnosis database, and make corresponding level line, be plotted on the same plane, corresponding three isocontour intersection points of different frequency can unique definite crackle relative position β and relative depth α, promptly determines relational expression

$(\mathrm{\&alpha;},\mathrm{\&beta;})=F_j^{-1}\left(f_j\right),(j=\mathrm{1,2,3})---\left(24\right)$

Embodiment 1: present embodiment verifies that mainly the BSWI unit finds the solution the horizontal calculation on Natural Frequency precision of rotor.For flawless single disk rotor system shown in Figure 2, get each shaft part length of rotor-support-foundation system and be respectively: L=300mm, L ₁=8mm, L ₂=133mm, L ₃=18mm, the armature spindle diameter d ₁=9.5mm, disk diameter d ₂=76mm, elastic modulus E=2.06 * 10 ¹¹N/m ², density p=7860kg/m ³, Poisson ratio μ=0.3.

Consider following four kinds of situations:

(1), five BSWI4 that do not comprise the moment of inertia influence ₃Beam element;

(2), 400 traditional beam elements that do not comprise the moment of inertia influence;

(3), five BSWI4 that comprise the moment of inertia influence ₃Beam element;

(4), 400 traditional beam elements that comprise the moment of inertia influence are found the solution.

Table 1 is found the solution horizontal calculation on Natural Frequency result of rotor and experiment test result for the variety classes unit.

The horizontal natural frequency result of rotor is found the solution and surveyed in table 1 variety classes unit

As shown in Table 1, adopt five BSWI43 beam element and 400 traditional beam element frequency solving results that do not comprise the moment of inertia influence that do not comprise the moment of inertia influence very approaching.But the second order frequency calculated value and experiment test value error are very big, if adopt these two kinds of unit to carry out rotor transverse crack quantitative Diagnosis based on model, will certainly draw false result.

Five BSWI43 beam elements that comprise the moment of inertia influence are consistent with 400 traditional beam element solving results that comprise the moment of inertia influence, and coincide with the experiment test frequency resultant, can be used for the rotor transverse crack quantitative Diagnosis.Simultaneously, can adopt five BSWI43 unit that comprise the moment of inertia influence to reach and 400 identical solving precision of traditional beam element that comprise the moment of inertia influence, calculated amount obviously reduces, the counting yield height.

Embodiment 2:

Present embodiment is mainly by rotor crack high precision quantitative diagnostic method in actual applications the validity of single span fatigue crack rotor checking based on the BSWI unit.For single disk rotor system shown in Figure 2, crackle appears at L ₂Shaft part is got each shaft part length of rotor-support-foundation system and is respectively: L=300mm, L ₁=8mm, L ₂=188mm, L ₃=18mm, the armature spindle diameter d ₁=9.5mm, disk diameter d ₂=76mm, elastic modulus E=2.06 * 10 ¹¹N/m ², density p=7860kg/m ³, Poisson ratio μ=0.3, crackle relative depth α=δ/d ₁, relative position β=e/L ₂

Horizontal fatigue crack operating mode of table 2 rotor and diagnostic result

With reference to shown in Figure 3, the frequency response curved surface that expression is obtained by the match of crackle diagnostic data base, promptly the rotor crack system frequency is with respect to the funtcional relationship of β and α.Suppose that crackle relative position and degree of depth span are α, β ∈ [0.1,0.9], the f that Fig. 3 (a) expression surface fitting obtains ₁=F ₁(α, β), the f that Fig. 3 (b) expression surface fitting obtains ₂=F ₂(α, β), the f that Fig. 3 (c) expression surface fitting obtains ₃=F ₃(α, β).

Employing power hammer knocks exciting, picks up impulse response signal with the Polytec laser vibration measurer, by to the response signal analysis, obtains the horizontal natural frequency in first three rank of rotor, and as the input of rotor crack quantitative Diagnosis, 14 BSWI4 are adopted in substitution respectively ₃This rotor crack high-precision diagnosis database that the unit is set up, and make corresponding level line, be plotted on the same plane.

With reference to shown in Figure 4, be the horizontal fatigue crack high-precision diagnosis of rotor experimental result.As shown in Table 2, the crackle operating mode is β=0.78889, α=0.3 and β=0.78889, α=0.4.Diagnostic result is β=0.78, α=0.29 and β=0.795, α=0.41.The relative error of crack position diagnosis is no more than 1.13%, and the relative error of crack depth diagnosis is no more than 3.33%.This shows: under different crackle operating modes, adopt 14 BSWI4 ₃The horizontal fatigue crack quantitative Diagnosis of the rotor database that the unit is set up has very high diagnostic accuracy in the horizontal fatigue crack diagnosis of actual rotor, can carry out the quantitative Diagnosis of the horizontal fatigue crack of rotor reliably.

Claims (1)

1. the interval B-spline wavelet unit is used for the method for rotor transverse crack quantitative Diagnosis, it is characterized in that,

1) node of configuration interval B spline wavelet, structure comprises the interval B-spline wavelet unit of structure moment of inertia; May further comprise the steps:

The employing exponent number is m, and yardstick is the interval B-spline wavelet scaling function of j, is designated as BSWIm _jScaling function, as interpolation basis function tectonic element, near configuration node singular point, unit are divided into the n=2 of unequal interval ^j+ m-4 section, cell node number are 2 ^j+ m-3, total number of degrees of freedom, is 2 ^j+ m-1, boundary node are 1, n+1; Internal node is 2,3 ..., n can obtain

w(ξ)＝N ^ew ^e

Shape function N ^eFor

$N^e=\mathrm{\&Phi;}{T}_b^e$

In the formula

For

$T_b^e={\left\{{[{\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_1\right)\frac{1}{l_e}\frac{\mathrm{d\&Phi;}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_1}{\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_2\right)...{\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_n\right){\mathrm{\&Phi;}}^T\left({\mathrm{\&xi;}}_{n+1}\right)\frac{1}{l_e}\frac{\mathrm{d\&Phi;}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_{n+1}}]}^{-1}\right\}}^T$

$\mathrm{\&Phi;}=\{{\mathrm{\&Phi;}}_{m,-m+1}^j\left(\mathrm{\&xi;}\right){\mathrm{\&phi;}}_{m,-m+2}^j\left(\mathrm{\&xi;}\right)\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&phi;}}_{m,2^j-1}^j\left(\mathrm{\&xi;}\right)\}$ Expression BSWIm _jThe capable vector of scaling function;

w _eFor physics degree of freedom column vector, be expressed as

w ^e＝{w ₁?θ ₁?w ₂…w _n?w _n+1?θ _n+1} ^T

In the formula, ${\mathrm{\&theta;}}_1=\frac{1}{l_e}\frac{\mathrm{dw}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_1}$ With ${\mathrm{\&theta;}}_{n+1}=\frac{1}{l_e}\frac{\mathrm{dw}\left(\mathrm{\&xi;}\right)}{\mathrm{d\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_{n+1}}$ Corner on the expression elementary boundary node, l _eBe element length;

The unit Vibration Frequency Equations is

$|{\mathrm{\&omega;}}^2(M_r^e+M_b^e)-K_b^e|=0$

In the formula, ω is circular frequency rad/s, makes vibration frequency f=ω/2 π Hz;

Unit rotational inertia consistent Mass Matrix

For

$M_r^e=\frac{\mathrm{\&rho;I}}{l_e}{\left(T_b^e\right)}^T{\mathrm{\&Gamma;}}^{1.1}{T}_b^e$

The unit bending consistent Mass Matrix

For

$M_b^e=\mathrm{\&rho;}{\mathrm{Al}}_e{\left(T_b^e\right)}^T{\mathrm{\&Gamma;}}^{\mathrm{0,0}}{T}_b^e$

The unit bending stiffness matrix For

$K_b^e={\frac{\mathrm{EI}}{l_e^3}\left(T_b^e\right)}^T{\mathrm{\&Gamma;}}^{2,2}{T}_b^e$

More than various in, E is an elastic modulus, I is a cross sectional moment of inertia, ρ is a density of material, A is an area of section, l _eBe element length,

For be transformed into the cell translation matrix the finite element space, Γ from wavelet space ^2,2, Γ ^0,0And Γ ^1,1Be respectively interval B-spline wavelet scaling function respective derivative integration matrix on the same yardstick.

The unit of being constructed 2) employing 1), set up rotor transverse crack interval B-spline wavelet finite element model, obtain crackle quantitative Diagnosis database, utilize the contouring method quantitative Diagnosis to go out relative position and relative depth that the horizontal fatigue crack of rotor exists, may further comprise the steps:

Utilize fracturing mechanics knowledge, can obtain the relevant torsion line rigidity k of crackle relative depth α _tAnd corresponding crackle stiffness matrix K _s, promptly

$k_t=\frac{\mathrm{\&pi;}{\mathrm{Er}}^8}{32(1-\mathrm{\&mu;})}\&times;\frac{1}{{\&Integral;}_{-r\sqrt{1-{(1-2a)}^2}}^{r\sqrt{1-{(1-2a)}^2}}(r^2-{\mathrm{\&xi;}}^2)\left[{\&Integral;}_0^{a\left(\mathrm{\&xi;}\right)}\mathrm{\&eta;}{F}^2(\mathrm{\&eta;}/H)\mathrm{d\&eta;}\right]\mathrm{d\&xi;}}$

And

$K_s=\left[\begin{array}{cc}{k}_t& -k_t\\ -k_t& k_t\end{array}\right]$

With crackle stiffness matrix K _sAdd in the overall finite element matrix, obtain the overall vibration frequency equation of rotor transverse crack system, i.e. rotor transverse crack interval B-spline wavelet finite element model

|ω ²(M _r+M _b)-K _b|＝0

In the formula, K _bBe stack of interval B-spline wavelet unit and adding crackle stiffness matrix K _sThe global stiffness matrix that the back forms, M _r+ M _bBe interval B-spline wavelet unit rotational inertia consistent Mass Matrix

With the unit bending consistent Mass Matrix

The overall consistent Mass Matrix that the stack back forms;

Find the solution rotor transverse crack interval B-spline wavelet finite element model by the eigenwert method for solving, obtain crackle quantitative Diagnosis database, utilize the contouring method quantitative Diagnosis in the quantitative identification of crackle to go out relative position and the relative depth that the horizontal fatigue crack of rotor exists.

Images