CN111241609A - Prediction method for blade tip clearance of rotor and stator assembly of aircraft engine - Google Patents

Prediction method for blade tip clearance of rotor and stator assembly of aircraft engine Download PDF

Info

Publication number
CN111241609A
CN111241609A CN202010065955.0A CN202010065955A CN111241609A CN 111241609 A CN111241609 A CN 111241609A CN 202010065955 A CN202010065955 A CN 202010065955A CN 111241609 A CN111241609 A CN 111241609A
Authority
CN
China
Prior art keywords
rotor
cas
blade tip
casing
clearance
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010065955.0A
Other languages
Chinese (zh)
Other versions
CN111241609B (en
Inventor
孙惠斌
王静
陈凯
夏华腾
袁剑
常智勇
王昱
邓王倩
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Northwestern Polytechnical University
AECC South Industry Co Ltd
Original Assignee
Northwestern Polytechnical University
AECC South Industry Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Northwestern Polytechnical University, AECC South Industry Co Ltd filed Critical Northwestern Polytechnical University
Priority to CN202010065955.0A priority Critical patent/CN111241609B/en
Publication of CN111241609A publication Critical patent/CN111241609A/en
Application granted granted Critical
Publication of CN111241609B publication Critical patent/CN111241609B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B21/00—Measuring arrangements or details thereof, where the measuring technique is not covered by the other groups of this subclass, unspecified or not relevant
    • G01B21/02—Measuring arrangements or details thereof, where the measuring technique is not covered by the other groups of this subclass, unspecified or not relevant for measuring length, width, or thickness
    • G01B21/04—Measuring arrangements or details thereof, where the measuring technique is not covered by the other groups of this subclass, unspecified or not relevant for measuring length, width, or thickness by measuring coordinates of points
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B21/00—Measuring arrangements or details thereof, where the measuring technique is not covered by the other groups of this subclass, unspecified or not relevant
    • G01B21/16—Measuring arrangements or details thereof, where the measuring technique is not covered by the other groups of this subclass, unspecified or not relevant for measuring distance of clearance between spaced objects

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Structures Of Non-Positive Displacement Pumps (AREA)

Abstract

The invention discloses a method for predicting blade tip clearance of rotor and stator assembly of an aircraft engine, which comprises the following steps of: s1, establishing a rotor blade tip calculation model; s2, establishing a support structure calculation model; s3, establishing a flow surface calculation model in the casing; s4, establishing a rotor and stator blade tip clearance calculation model, and inputting deviation data and size data of the rotor system, the supporting frame, the casing and the stator assembly into the model for calculation; the model outputs predicted values of clearances of rotor tips and stator tips of each stage under different rotor phases according to input data, gives corresponding clearance distribution curves, and obtains the maximum clearance, the minimum clearance and the average clearance of each stage of rotor and stator predicted by the model; and comparing the prediction result with the clearance requirement of the process standard to judge whether the blade tip clearance index of the rotor-stator assembly is qualified or not, or giving out an out-of-tolerance part according to the predicted value, namely completing the prediction of the blade tip clearance of the rotor-stator assembly of the aero-engine. The problem of current aeroengine rotor stator assembly apex clearance prediction difficulty is solved.

Description

Prediction method for blade tip clearance of rotor and stator assembly of aircraft engine
[ technical field ] A method for producing a semiconductor device
The invention belongs to the technical field of aero-engine assembly, and particularly relates to a prediction method for an aero-engine rotor and stator assembly blade tip clearance.
[ background of the invention ]
The rotor and stator blade tip clearance of the aero-engine is an important index influencing the safety and the performance of the aero-engine, and the problems of uneven circumferential distribution of the rotor and stator blade tip clearance and collision and abrasion of a rotor and an inner runner surface are a great problem which puzzles all major airlines and military enterprises in the world. A large number of errors exist in the manufacturing and assembling process of the aero-engine, blade tip gaps of aero-engine rotor and stator assembling cannot be accurately predicted due to accumulation and transmission of multiple errors, requirements can be met only through repeated assembly and adjustment, and accordingly assembling efficiency is low and quality consistency is poor. At present, the related disclosures of the blade tip clearance prediction method of the rotor and the stator of the aeroengine are less in China.
[ summary of the invention ]
The invention aims to provide a prediction method of blade tip clearance of an aircraft engine rotor and stator assembly, and aims to solve the problem that the blade tip clearance of the existing aircraft engine rotor and stator assembly is difficult to predict.
The invention adopts the following technical scheme: a method for predicting blade tip clearance of an aeroengine rotor-stator assembly, wherein the aeroengine comprises a rotor system, a supporting frame, a casing and a stator assembly, the blade tip clearance refers to the radial clearance of a rotor blade tip and a flow passage surface in the casing and the radial clearance of a stator blade tip and a rotor hub, and the method comprises the following steps:
s1, establishing a rotor blade tip calculation model;
s2, establishing a support structure calculation model;
s3, establishing a flow surface calculation model in the casing;
and S4, integrating the three calculation models obtained in the steps 1 to 3 into a rotor and stator blade tip clearance calculation model by using a homogeneous coordinate transformation matrix under an absolute coordinate system, and predicting the blade tip clearances under different rotor phase angles.
Further, in the step 1, a space coordinate system is established by taking the centroid of the axial reference plane of the primary rotor or the front journal as an origin; detecting a contour coordinate point of the rear end face of the first-stage rotor by using a part, establishing a space plane equation by using the coordinate point, and sequentially calculating rear end face plane equations of the rotors at all stages between the front support and the rear support; calculating unit normal vectors of rear end faces of the rotors at all levels by using a plane equation, namely completing geometric modeling of the end face run-out deviation; detecting the contour coordinate points of the outer side spoke plate surfaces of the front end surface and the rear end surface of the rotor by parts, establishing a circular equation in a cross section plane, estimating equation parameters by a least square method to obtain a centroid coordinate, namely completing radial eccentric deviation geometric modeling; calculating the direction vector of the inner centroid of the single-stage rotor according to the centroid coordinates of the front end face and the rear end face; and calculating the height and the initial phase angle of the blade tip of each blade so as to calculate the actual coordinates of the blade tip point, namely completing the rotor blade tip calculation model.
Further, the specific calculation process of step 1 is as follows:
establishing a coordinate system by taking a rotor datum plane as an assembly global datum, wherein an original point is a rotor datum plane center point, the axial stacking direction is the positive x direction, a cross-sectional plane is a yOz plane, and the vertical direction is the z-axis direction;
establishing a plane equation according to the profile data points of the rear end face of the ith-stage rotor acquired by a detection instrument, and fitting the parameters of the plane equation by a least square method to obtain an equation expression:
Apla,ix+Bpla,iy+Cpla,iz+Dpla,i＝0 (1)，
in the formula, Apla,Bpla,Cpla,DplaAll the coefficients are calibrated plane parameter coefficients, and a unit normal vector of the rear end face of the i-th-stage rotor is obtained by performing normalization processing on a plane normal vector obtained by the formula (1):
Figure BDA0002375971430000021
the method comprises the following steps of obtaining a section plane outer ring profile of a radial plate on the outer side of the front end face of the single-stage rotor at the axial x position by a rotary table detection instrument, and fitting according to a least square method to obtain a front end face outer side circular equation:
Acir,f,iy2+Bcir,f,iz2+Ccir,f,iy+Dcir,f,iz+Ecir,f,i＝0 (3)，
in the above formula, Acir,f,i,Bcir,f,i,Ccir,f,i,Dcir,f,i,Ecir,f,iCircle characteristic parameters of the front end face of the ith-stage rotor, which are calibrated by a least square method, are approximately estimated by a section plane circular equation, and the centroid position vector of the front end face of the ith-stage rotor is as follows:
Figure BDA0002375971430000031
in the formula (I), the compound is shown in the specification,
Figure BDA0002375971430000032
the approximate centroid position vector (x) of the rear end face can be obtained by the same principle of the formula (3) and the formula (4) for the accumulated coordinates of the axial dimension of the front i-1-stage rotorcir,l,ycir,l,zcir,l) Component nx,iThe unit direction vector of the ith stage rotor form axis is given by:
Figure BDA0002375971430000033
in the formula (I), the compound is shown in the specification,
Figure BDA0002375971430000034
is the die length of the axial vector of the ith-stage rotor inner form;
let the axial positioning dimension of the ith-stage rotor blade tip be liDesigned radial height dimension RiLet the total number of i-th-stage rotor blades be NiThe circumferential included angle of the adjacent blades is
Figure BDA0002375971430000035
Taking the positive direction of a z axis as an initial edge of a phase angle in a yOz plane, and numbering the circumferential blades from 1 to j in an increasing way along the positive direction of the phase angle, wherein j is more than 0 and is not more than NiIn which the phase angle of blade number 1 is phii,1The phase angle of the j blade is phii,j＝φi,1+(j-1)φi,Δ；
The unit direction vector from the centroid to the zero phase point of the blade tip in the section plane of the blade tip of the ith-stage rotor is set as
Figure BDA0002375971430000036
It satisfies the mutually perpendicular condition in space:
Figure BDA0002375971430000041
then the vector winding mandrel is formed from the cross section axis to the zero phase point of the corresponding unit direction vector of the j blade tip of the ith-stage rotor
Figure BDA0002375971430000042
Is rotated to obtain
Figure BDA0002375971430000043
Satisfies the formula:
Figure BDA0002375971430000044
the j-th blade tip of the ith-stage rotor has a designed height RiRun-out of deltaj(ii) a Thus, the actual height of the corresponding blade is Ri,jTherefore, the absolute coordinates of the j-th blade tip of the ith-stage rotor are as follows:
Figure BDA0002375971430000045
further, in the step 2, a coordinate system is established by taking the centroid of the left end face of the inner ring of the front bearing matched with the rotor as an origin, and the connecting line of the centroids of the adjacent end faces of the inner rings of the front bearing and the rear bearing is taken as the rotation axis of the rotor; calculating a unit vector of a rotor rotating axis, namely completing modeling of a rotor rotating shaft; and substituting the simulated deformation coordinate component of the elastic support and the bearing clearance value, and calculating a center offset vector of the reference surface of the casing under the bearing clearance deviation, namely completing the modeling of the supporting frame.
Further, the specific calculation process of step 2 is as follows: two bearings are arranged at the front and the back of the rotor structure, the bearings have radial play and axial play deviation, the left end face of the inner ring of the bearing is taken as a yOz plane, the centroid of the end face is taken as a coordinate origin O, the vertical direction is the z-axis direction, and the direction of the centroid towards the inner normal vector is the x direction; axial clearance of bearing is Deltaxcle,(xcle,min≤Δxcle≤xcle,max) Radial play of bearing of Deltarcle,(rcle,min≤Δrcle≤rcle,max)；
Respectively measuring the arc profiles of the rear end surfaces of the inner ring and the outer ring of the front bearing under a global coordinate system established by the rotor reference surface, and calculating to obtain the centroid coordinates of (x)f,in,cle,yf,in,cle,zf,in,cle) And (x)f,out,cle,yf,out,cle,zf,out,cle) (ii) a Similarly, the circular arc profiles of the front end surfaces of the inner ring and the outer ring of the rear bearing are respectively measured, and the centroid coordinates are calculated to be (x)l,in,cle,yl,in,cle,zl,in,cle) And (x)l,out,cle,yl,out,cle,zl,out,cle) And then the unit direction vector of the inner ring axis of the bearing is as follows:
Figure BDA0002375971430000051
the unit direction vector of the outer ring axis can be obtained by the same method
Figure BDA0002375971430000052
According to the structure simulation result of the elastic support, the deviation position vector of the reference plane of the casing caused by the deformation of the front elastic support and the rear elastic support is (delta x)sup,f,δysup,f,δzsup,f) And (δ x)sup,l,δysup,l,δzsup,l)。
Further, in step 3, establishing a coordinate system by taking the centroid of the left end face of the casing as an origin, measuring contour points at the matching part of the inner runner and the rotor blade tip, and recording the axial position and the coordinates of the contour points as model input; establishing a contour elliptic equation, calibrating parameters by using a least square method according to contour coordinates, and extracting the phases of the center, the long half shaft, the short half shaft and the long shaft of the elliptic equation; and calculating an interpolation centroid, substituting different phase flow channel surface contour jumping values in an elliptic equation, and establishing an inner channel surface contour equation to finish the casing calculation model.
Further, the specific establishment process of the casing calculation model is as follows:
establishing an inner runner surface elliptical equation of the single-section case, taking the centroid of the reference end surface of the case as the origin of coordinates, taking the axial direction as the x direction, setting the positive direction to be the same as the stacking direction of the rotors, setting a coordinate system for a section plane which is a yOz plane, and taking the vertical direction as the z-axis direction, and then, taking the inner runner surface elliptical equation matched with the i-th-stage rotor as the equation,
y2+Acas,iyz+Bcas,iz2+Ccas,iy+Dcas,iz+Ecas,i＝0 (10)，
in the formula, Acas,i,Bcas,i,Ccas,i,Dcas,i,Ecas,iFitting ellipse with least square method to obtain the centroid coordinate (x) of the flow channel cross sectioncas,i,ycas,i,zcas,i) The non-negative minimum phase angle of the major and minor axes of the elliptic equation is thetacas,iThe following conditions are satisfied:
Figure BDA0002375971430000061
when the casing assembly structure is horizontally placed, the relative sinking amount of the casing assembly structure in the global coordinate system of the rotor reference is
Figure BDA0002375971430000062
The centroid coordinate of the front datum plane of the casing is (x) in combination with the position deviation of the casing caused by the deformation of the supportcas,ben,ycas,ben,zcas,ben)；
According to the design shape of the flow channel surface in the casing, the change equation of the section profile radius of the flow channel surface in the casing along the axial direction is RcaR (x); and then the centroid curve is calculated.
Further, when the casing is a split casing, the centroid curve calculation method is as follows:
establishing cubic spline interpolation functions of coordinate components y and z relative to x for bisection of the centroid of the flow channel surface in the half casing respectively, wherein the interpolation node matrix formula is as follows:
Figure BDA0002375971430000063
in the formula, hiFor the step size of the node,
Figure BDA0002375971430000064
for interpolated node vectors, f [ x ]i-1,xi,xi+1]For third order difference, constraint m0,mnThe component y of the end face attitude vector of the left end face and the component z of the end face attitude vector of the right end face of the part can be calculated about the vector included angle of the axial direction x, and then the interpolation function is as follows:
Figure BDA0002375971430000071
wherein s (x) denotes ycas,c(x) Or zcas,c(x) Calculating to obtain (x)cas,c,ycas,c,zcas,c)。
Further, when the casing is an integral thin-wall annular casing, the centroid curve calculation mode is as follows:
for the integral thin-wall annular casing, the elliptical profile equation of the inner runner of the front end surface and the inner runner of the rear end surface of the ith section of the casing are given by the formula (10), and after undetermined parameters are evaluated by a least square method, the central coordinate of the inner ellipse of the front end surface and the deflection angle y of the long axis are calculated according to the formula (11)cas,f,i,zcas,f,i,θcas,f,iAnd the central coordinate of the ellipse and the deviation angle y of the long axis in the section plane of the rear end facecas,l,i,zcas,l,i,θcas,l,iAnd if the section i of the casing front end surface inner flow channel profile equation is obtained, the coordinates of any centroid points in the casing are as follows:
Figure BDA0002375971430000072
after the multiple sections of casings are stacked, the centroid coordinates of the front end face of the ith section of casing are calculated in the same way as in the formula (4):
Figure BDA0002375971430000073
in the formula (I), the compound is shown in the specification,
Figure BDA0002375971430000074
is a unit normal vector of the rear end face of the ith section of the casing, and is calculated by formula (1) and formula (2), Ecas,f,iCalibrating the section parameter of the front end surface ellipse of the ith section of the casing by a least square method;
then the centroid coordinate of any point in the ith section of the casing is the centroid coordinate (x) of the front end face and the rear end face of the single section of the casingcas,f,i,ycas,f,i,zcas,f,i) And (x)cas,l,i,ycas,l,i,zcas,l,i) Substituting formula (5) to obtain (x)cas,c,ycas,c,zcas,c)；
The total collection point number for measuring the radius value of the section of the flow channel matched with the blade tip of the ith-stage rotor is Ni,casAngle between collection points
Figure BDA0002375971430000075
The measured data can obtain the runout delta of the point to be measured of the runner surface relative to the design sizei,j,casEach data available axial position xi,casAnd phase angle phii,j,cas＝j·Δφi,casA unique representation; the coordinates of the inner flow channel surface corresponding to the axial arbitrary position and the phase angle of the casing are as follows:
Figure BDA0002375971430000081
in the formula (I), the compound is shown in the specification,
Figure BDA0002375971430000082
[g]is the rounding function, i.e. the largest integer not exceeding the value of g.
Further, the specific calculation process of step 4 is as follows:
under an absolute coordinate system, when the rotor does not generate phase deflection, the section of a flow channel matched with the j-th blade tip of the ith-stage rotor meets the relation:
Figure BDA0002375971430000083
substituting the condition satisfied by the formula (17) into the formula (16), and calculating the coordinate of the flow channel point matched with the formula; the assembly blade tip clearance cle of the i-th-stage rotor and the j-th blade in the non-rotating state with the blade tip rotating angle of 0 DEGi,jExpressed as:
Figure BDA0002375971430000084
when the rotor rotating shaft deflects by a phase angle thetarotIn the meantime, an initial blade tip position vector of the j blade of the i-th-stage rotor is set
Figure BDA0002375971430000085
Rotated tip position vector
Figure BDA0002375971430000086
And blade direction vector
Figure BDA0002375971430000087
Respectively satisfy the following formula：
Figure BDA0002375971430000088
Then the current phase angle thetarotUnder the condition, the corresponding point of a runner surface matched with the j-th stage rotor blade meets the relation:
Figure BDA0002375971430000091
calculating the coordinate of the corresponding point of the runner surface matched with the j-th blade of the i-th-stage rotor by using the conditional formula (20) back to the formula (16), and then, the blade tip clearance
Figure BDA0002375971430000092
Expressed as:
Figure BDA0002375971430000093
the maximum clearance cle of the mth stage rotorm,maxMinimum clearance clem,minRespectively as follows:
Figure BDA0002375971430000094
the number of times of measurement of the ith-stage rotor rotating for one circle at a certain rotation angle is set to be Nrot,iThen the phase angle of a single rotation is
Figure BDA0002375971430000095
Wherein, k is the corner number, the average clearance cle of the mth stage rotorm,meanComprises the following steps:
Figure BDA0002375971430000096
the invention has the beneficial effects that: according to the structural characteristics of the rotor and the stator and the actual assembly process, the blade tip clearance prediction model based on error transmission consists of three parts, namely a rotor blade tip, a supporting frame and a stator casing. In order to facilitate data acquisition and simplified calculation, local error transfer models are respectively established for the three parts according to an assembling procedure, and finally, a homogeneous coordinate transformation matrix is integrated into a complete blade tip clearance error transfer model in an absolute coordinate system, and blade tip clearances under different rotor phase angles are predicted. The method provided by the invention relates to the deviation data of each part, establishes a mathematical model through the rotor and stator assembly relation, provides a solution for the current situation that rotor and stator blade tip clearance is difficult to assemble and predict, helps an assembly engineer to analyze the accurate distribution of the rotor and stator blade tip clearance, and provides a basis for the optimization of the assembly process aiming at the rotor and stator blade tip clearance.
[ description of the drawings ]
FIG. 1 is a general frame diagram of a rotor stator blade tip clearance calculation model;
FIG. 2 is a schematic illustration of a rotor stator tip clearance forming mechanism;
FIG. 3 is a graph of the rotor blade tip and the corresponding flow surface in the casing;
FIG. 4 is a box diagram of the calculation results of different phase clearances of the tips of the rotors at different stages.
[ detailed description ] embodiments
The present invention will be described in detail below with reference to the accompanying drawings and specific embodiments.
The invention provides a prediction method of blade tip clearance in the assembly of an aero-engine rotor and stator, and an assembly structure influencing the blade tip clearance in the assembly of the aero-engine rotor and stator comprises three parts: rotor system, supporting frame and casing and stator assembly. The blade tip clearance refers to the radial clearance between the blade tip of the rotor and the inner flow channel surface of the casing and the radial clearance between the blade tip of the stator and the hub of the rotor. The method comprises the following steps:
s1, establishing a rotor blade tip calculation model:
establishing a space coordinate system by taking the centroid of the front end face of the primary rotor as an origin; detecting a contour coordinate point of the rear end face by using a part, establishing a space plane equation by using the coordinate point, and sequentially calculating rear end face plane equations of rotors at all stages between the front support and the rear support; calculating unit normal vectors of rear end faces of the rotors at all levels by using a plane equation, namely completing geometric modeling of the end face run-out deviation; detecting the outline coordinate points of the outer side spoke plate surfaces of the front end surface and the rear end surface of the rotor by parts, establishing an elliptic equation in a section plane, estimating equation parameters by a least square method to obtain a centroid coordinate, namely completing radial eccentric deviation geometric modeling; calculating the direction vector of the inner centroid of the single-stage rotor according to the centroid coordinates of the front end face and the rear end face; calculating the height and the initial phase angle of the blade tip of each blade so as to calculate the actual coordinate of each blade tip point, namely completing a rotor blade tip calculation model;
the specific calculation process is as follows:
establishing a coordinate system by taking a rotor datum plane as an assembly global datum, wherein an original point is a rotor datum plane center point, the axial stacking direction is the positive x direction, a cross-sectional plane is a yOz plane, and the vertical direction is the z-axis direction;
establishing a plane equation according to the profile data points of the rear end face of the ith-stage rotor acquired by a detection instrument, and fitting the parameters of the plane equation by a least square method to obtain an equation expression:
Apla,ix+Bpla,iy+Cpla,iz+Dpla,i＝0 (24)，
in the formula, Apla,Bpla,Cpla,DplaAll the coefficients are calibrated plane parameter coefficients, and a unit normal vector of the rear end face of the i-th-stage rotor is obtained by performing normalization processing on a plane normal vector obtained by the formula (1):
Figure BDA0002375971430000111
the method comprises the following steps of obtaining a section plane outer ring profile of a radial plate on the outer side of the front end face of the single-stage rotor at the axial x position by a rotary table detection instrument, and fitting according to a least square method to obtain a front end face outer side circular equation:
Acir,f,iy2+Bcir,f,iz2+Ccir,f,iy+Dcir,f,iz+Ecir,f,i＝0 (26)，
in the above formula, Acir,f,i,Bcir,f,i,Ccir,f,i,Dcir,f,i,Ecir,f,iCircle characteristic parameters of the front end face of the ith-stage rotor, which are calibrated by a least square method, are approximately estimated by a section plane circular equation, and the centroid position vector of the front end face of the ith-stage rotor is as follows:
Figure BDA0002375971430000112
in the formula (I), the compound is shown in the specification,
Figure BDA0002375971430000113
the approximate centroid position vector (x) of the rear end face can be obtained by the same principle of the formula (3) and the formula (4) for the accumulated coordinates of the axial dimension of the front i-1-stage rotorcir,l,ycir,l,zcir,l) Component nx,iThe unit direction vector of the ith stage rotor form axis is given by:
Figure BDA0002375971430000114
in the formula (I), the compound is shown in the specification,
Figure BDA0002375971430000115
is the die length of the axial vector of the ith-stage rotor inner form;
let the axial positioning dimension of the ith-stage rotor blade tip be liDesigned radial height dimension RiLet the total number of i-th-stage rotor blades be NiThe circumferential included angle of the adjacent blades is
Figure BDA0002375971430000116
Taking the positive direction of a z axis as an initial edge of a phase angle in a yOz plane, and numbering the circumferential blades from 1 to j in an increasing way along the positive direction of the phase angle, wherein j is more than 0 and is not more than NiIn which the phase angle of blade number 1 is phii,1The phase angle of the j blade is phii,j＝φi,1+(j-1)φi,Δ；
The unit direction vector from the centroid to the zero phase point of the blade tip in the section plane of the blade tip of the ith-stage rotor is set as
Figure BDA0002375971430000121
It satisfies the mutually perpendicular condition in space:
Figure BDA0002375971430000122
then the vector winding mandrel is formed from the cross section axis to the zero phase point of the corresponding unit direction vector of the j blade tip of the ith-stage rotor
Figure BDA0002375971430000123
Is rotated to obtain
Figure BDA0002375971430000124
Satisfies the formula:
Figure BDA0002375971430000125
the j-th blade tip of the ith-stage rotor has a designed height RiRun-out of deltaj(ii) a Thus, the actual height of the corresponding blade is Ri,jTherefore, the absolute coordinates of the j-th blade tip of the ith-stage rotor are as follows:
Figure BDA0002375971430000126
s2, establishing a support structure calculation model:
establishing a coordinate system by taking the centroid of the left end face of the front bearing inner ring matched with the rotor as an origin, and taking the connecting line of the centroids of the adjacent end faces of the front bearing inner ring and the rear bearing inner ring as the rotor rotation axis; calculating a unit vector of a rotor rotating axis, namely completing modeling of a rotor rotating shaft; and substituting the simulated deformation coordinate component of the elastic support and the bearing clearance value, and calculating a center offset vector of the reference surface of the casing under the bearing clearance deviation, namely completing the modeling of the supporting frame.
The specific calculation process is as follows: two bearings are arranged at the front and the back of the rotor structure, the bearings have radial play and axial play deviation, the left end face of the inner ring of the bearing is taken as a yOz plane, the centroid of the end face is taken as a coordinate origin O, the vertical direction is the z-axis direction, and the direction of the centroid towards the inner normal vector is the x direction; axial clearance of bearing is Deltaxcle,(xcle,min≤Δxcle≤xcle,max) Radial play of bearing of Deltarcle,(rcle,min≤Δrcle≤rcle,max)；
Respectively measuring the arc profiles of the rear end surfaces of the inner ring and the outer ring of the front bearing under a global coordinate system established by the rotor reference surface, and calculating to obtain the centroid coordinates of (x)f,in,cle,yf,in,cle,zf,in,cle) And (x)f,out,cle,yf,out,cle,zf,out,cle) (ii) a Similarly, the circular arc profiles of the front end surfaces of the inner ring and the outer ring of the rear bearing are respectively measured, and the centroid coordinates are calculated to be (x)l,in,cle,yl,in,cle,zl,in,cle) And (x)l,out,cle,yl,out,cle,zl,out,cle) And then the unit direction vector of the inner ring axis of the bearing is as follows:
Figure BDA0002375971430000131
the unit direction vector of the outer ring axis can be obtained by the same method
Figure BDA0002375971430000132
According to the structure simulation result of the elastic support, the deviation position vector of the reference plane of the casing caused by the deformation of the front elastic support and the rear elastic support is (delta x)sup,f,δysup,f,δzsup,f) And (δ x)sup,l,δysup,l,δzsup,l)。
S3, establishing a flow surface calculation model in the casing:
establishing a coordinate system by taking the centroid of the left end face of the casing as an origin, measuring a contour point at the matching part of the inner runner and the rotor blade tip, and recording an axial position and a contour point coordinate as model input; establishing a contour elliptic equation, calibrating parameters by using a least square method according to contour coordinates, and extracting the phases of the center, the long half shaft, the short half shaft and the long shaft of the elliptic equation; and calculating an interpolation centroid, substituting different phase flow channel surface contour jumping values in an elliptic equation, and establishing an inner channel surface contour equation to finish the casing calculation model.
The specific calculation process is as follows:
establishing an inner runner surface elliptical equation of the single-section case, taking the centroid of the reference end surface of the case as the origin of coordinates, taking the axial direction as the x direction, setting the positive direction to be the same as the stacking direction of the rotors, setting a coordinate system for a section plane which is a yOz plane, and taking the vertical direction as the z-axis direction, and then, taking the inner runner surface elliptical equation matched with the i-th-stage rotor as the equation,
y2+Acas,iyz+Bcas,iz2+Ccas,iy+Dcas,iz+Ecas,i＝0 (33)，
in the formula, Acas,i,Bcas,i,Ccas,i,Dcas,i,Ecas,iFitting ellipse with least square method to obtain the centroid coordinate (x) of the flow channel cross sectioncas,i,ycas,i,zcas,i) The non-negative minimum phase angle of the major and minor axes of the elliptic equation is thetacas,iThe following conditions are satisfied:
Figure BDA0002375971430000141
when the casing assembly structure is horizontally placed, the relative sinking amount of the casing assembly structure in the global coordinate system of the rotor reference is
Figure BDA0002375971430000142
The centroid coordinate of the front datum plane of the casing is (x) in combination with the position deviation of the casing caused by the deformation of the supportcas,ben,ycas,ben,zcas,ben)。
According to the design shape of the flow channel surface in the casing, the change equation of the section profile radius of the flow channel surface in the casing along the axial direction is RcaR (x). According to different structural characteristics of the casings, different centroid curve calculation modes are respectively adopted for the halved casing and the integral thin-wall annular casing.
(1) When the casing is a split casing, the centroid curve calculation method is as follows:
establishing cubic spline interpolation functions of coordinate components y and z relative to x for bisection of the centroid of the flow channel surface in the half casing respectively, wherein the interpolation node matrix formula is as follows:
Figure BDA0002375971430000143
in the formula, hiIs a node stepThe length of the utility model is long,
Figure BDA0002375971430000151
for interpolated node vectors, f [ x ]i-1,xi,xi+1]For third order difference, constraint m0,mnThe component y of the end face attitude vector of the left end face and the component z of the end face attitude vector of the right end face of the part can be calculated about the vector included angle of the axial direction x, and then the interpolation function is as follows:
Figure BDA0002375971430000152
wherein s (x) denotes ycas,c(x) Or zcas,c(x) Calculating to obtain (x)cas,c,ycas,c,zcas,c)。
(2) When the casing is an integral thin-wall annular casing, the centroid curve calculation mode is as follows:
for the integral thin-wall annular casing, the elliptical profile equation of the inner runner of the front end surface and the inner runner of the rear end surface of the ith section of the casing are given by the formula (10), and after undetermined parameters are evaluated by a least square method, the central coordinate of the inner ellipse of the front end surface and the deflection angle y of the long axis are calculated according to the formula (11)cas,f,i,zcas,f,i,θcas,f,iAnd the central coordinate of the ellipse and the deviation angle y of the long axis in the section plane of the rear end facecas,l,i,zcas,l,i,θcas,l,iAnd if the section i of the casing front end surface inner flow channel profile equation is obtained, the coordinates of any centroid points in the casing are as follows:
Figure BDA0002375971430000153
after the multiple sections of casings are stacked, the centroid coordinates of the front end face of the ith section of casing are calculated in the same way as in the formula (4):
Figure BDA0002375971430000154
in the formula (I), the compound is shown in the specification,
Figure BDA0002375971430000155
is the i-th section casing rear endA normal vector of surface unit calculated by the following formulae (1) and (2), Ecas,f,iCalibrating the section parameter of the front end surface ellipse of the ith section of the casing by a least square method;
then the centroid coordinate of any point in the ith section of the casing is the centroid coordinate (x) of the front end face and the rear end face of the single section of the casingcas,f,i,ycas,f,i,zcas,f,i) And (x)cas,l,i,ycas,l,i,zcas,l,i) Substituting formula (5) to obtain (x)cas,c,ycas,c,zcas,c)；
The total collection point number for measuring the radius value of the section of the flow channel matched with the blade tip of the ith-stage rotor is Ni,casAngle between collection points
Figure BDA0002375971430000161
The measured data can obtain the runout delta of the point to be measured of the runner surface relative to the design sizei,j,casEach data available axial position xi,casAnd phase angle phii,j,cas＝j·Δφi,casA unique representation; the coordinates of the inner flow channel surface corresponding to the axial arbitrary position and the phase angle of the casing are as follows:
Figure BDA0002375971430000162
in the formula (I), the compound is shown in the specification,
Figure BDA0002375971430000163
[g]is the rounding function, i.e. the largest integer not exceeding the value of g.
S4, establishing an absolute coordinate system by taking a primary rotor front end face assembly reference centroid as a global reference, taking the reference centroid as an original point, taking a reference normal as an x-axis direction, and taking a reference plane as a yOz plane; calculating the blade tip coordinates under an absolute coordinate system; calculating an internal runner surface equation under an absolute coordinate system, and substituting the blade tip coordinates into the equation to calculate the absolute coordinates of the runner fitting points; calculating the initial blade tip clearance of each stage of rotor; inputting a deflection phase angle, calculating blade tip clearances of different rotor phases, and recording numerical values to complete all modeling; after modeling is completed, inputting deviation data and size data required by the rotor, support and stator casing models into a model for calculation; the model outputs predicted values of clearances of rotor tips and stator tips of each stage under different rotor phases according to input data, gives corresponding clearance distribution curves, and obtains the maximum clearance, the minimum clearance and the average clearance of each stage of rotor and stator predicted by the model; and comparing the prediction result with the clearance requirement of the process standard to judge whether the blade tip clearance index of the rotor-stator assembly is qualified or not, or giving out an out-of-tolerance part according to the prediction value to provide a basis for reprocessing, namely completing the prediction of the blade tip clearance of the rotor-stator assembly of the aero-engine.
The method of the invention predicts the blade tip clearance after the assembly according to the measured data of the parts in the assembly process, because the blade tip clearance is invisible after the assembly. And fitting the measured data of different characteristics of the part according to a model method to obtain deviation characteristics, expressing the deviation characteristics by using a centroid and a unit direction vector of the model, then establishing a size chain transfer model of the deviation characteristics, and finally calculating the blade tip clearance. The method can obtain the prediction result of the blade tip clearance by introducing the measurement deviation data and the size data of the parts.
The specific calculation process is as follows:
the core part is a process that the relative position relation of 3 axis curves and 1 height vector shown in figure 2 changes along with the phase angle change of the rotor blade tip. Under an absolute coordinate system, when the rotor does not generate phase deflection, the section of a flow channel matched with the j-th blade tip of the ith-stage rotor meets the relation:
Figure BDA0002375971430000171
substituting the condition satisfied by the formula (17) into the formula (16), and calculating the coordinate of the flow channel point matched with the formula; the assembly blade tip clearance cle of the i-th-stage rotor and the j-th blade in the non-rotating state with the blade tip rotating angle of 0 DEGi,jExpressed as:
Figure BDA0002375971430000172
when the rotor rotating shaft deflects by a phase angle thetarotIn the meantime, an initial blade tip position vector of the j blade of the i-th-stage rotor is set
Figure BDA0002375971430000173
Rotated tip position vector
Figure BDA0002375971430000174
And blade direction vector
Figure BDA0002375971430000175
Respectively satisfy the following formula:
Figure BDA0002375971430000176
then the current phase angle thetarotUnder the condition, the corresponding point of a runner surface matched with the j-th stage rotor blade meets the relation:
Figure BDA0002375971430000177
calculating the coordinate of the corresponding point of the runner surface matched with the j-th blade of the i-th-stage rotor by using the conditional formula (20) back to the formula (16), and then, the blade tip clearance
Figure BDA0002375971430000181
Expressed as:
Figure BDA0002375971430000182
the maximum clearance cle of the mth stage rotorm,maxMinimum clearance clem,minRespectively as follows:
Figure BDA0002375971430000183
the number of times of measurement of the ith-stage rotor rotating for one circle at a certain rotation angle is set to be Nrot,iThen the phase angle of a single rotation is
Figure BDA0002375971430000184
Wherein, k is the corner number, the average clearance cle of the mth stage rotorm,meanComprises the following steps:
Figure BDA0002375971430000185
the invention relates to a prediction method of an aircraft engine rotor and stator assembly blade tip clearance, which comprises the steps of collecting the coordinates of contour characteristic points of front and rear end faces of each stage of rotor, the coordinates of circular contour points of a spoke plate face of an outer ring contour of the end face and blade tip jumping as input; inputting the axial size, the radial size, the tip height size, the number of blades and the circumferential phase angle of the No. 1 blade of the rotor according to an actual assembly structure; fitting the coordinate points of the deviation data according to characteristics, and expressing the matching characteristics with the deviation; the real end face is expressed by a plane equation, the real centroid and the spoke plate face profile are expressed by a space circle equation, and the axial curve of the rotor centroid is expressed in a vector form; when an enterprise assembles a blade disc on an aeroengine to adjust the circumferential installation angle, a certain position on the end face of the blade disc is marked so as to determine the circumferential phase angle when the blade disc is installed. To facilitate calculation of the actual height of each circumferential blade tip, it is defined that the zero phase positive first blade is denoted herein as blade number 1. After each step of assembly process is completed, for the rotor being assembled, the size data and the deviation data are substituted into the calculation model to complete one deviation transmission, and simultaneously, the initial phase tip coordinates of all the tip points of the rotor in the stage along the circumferential direction are generated in the model; the rotor assembly is a process of stacking parts, wherein the rotor of the next stage is assembled by taking the end face of the previous stage as a reference, the size deviation is transmitted among stages by influencing the centroid pose data of the matched surface, and finally the size deviation is mapped to the coordinates of the rotor blade tip point of each stage.
Detecting a bearing radial play part to obtain data as a rigidity deviation, simulating elastic support deformation by a structural finite element to obtain data as a flexibility deviation, and inputting the deviation data of the two parts as a calculation model of a support structure; the connecting line of the centroids of the front bearing and the rear bearing is the rotation axis of the rotor, and the rotor rotates around the rotation axis for a circle at a low speed, so that the blade tip gaps under different phase angles are calculated to present different distribution characteristics.
In the actual assembly and detection process of the casing structure, measuring the contour of a circle at an axially specified position of a flow channel surface in the casing in the circumferential direction, and acquiring contour coordinate point data of the cross section of the inner flow channel under an axial constant coordinate; establishing an inner runner surface elliptical equation with multiple sections, and determining the number and the axial position of the inner runner surface elliptical section equation according to the number of stages and the axial size of a rotor matched with a casing; the casing calculation model input deviation data comprises: elliptical deformation deviation, end face run-out deviation, inner runner face run-out deviation and stator blade tip run-out deviation of the casing; for the case, different case centroid estimation equations are established according to different structure types. The method comprises the following steps that a casing opening axial direction is a continuum, and a casing centroid at any axial position is calculated in a cubic spline interpolation mode; and for the integral thin-wall annular casing, the integral thin-wall annular casing is assembled in a multi-section stacking mode, and error transmission caused by the integral thin-wall annular casing is the same as the rotor stacking principle.
For two different casing structures, different centroid evaluation schemes are employed: for a half-split casing, establishing an interpolation function for calculating the centroid of any axial position by cubic spline interpolation; for the integral thin-wall annular casing, calculating a plane equation of the rear end face of the single-section casing to complete end face run-out deviation modeling; calculating the centroid by taking the contour detection data of the inner side flow channel of the front end surface and the rear end surface of the casing as input; and establishing a calculation function of the centroid coordinate of the axial arbitrary position in the single-section casing by using the centroid connecting line vector of the front end surface and the rear end surface.
Calculating in the three steps, and respectively simulating deviation transfer rules influencing the assembling processes of the rotor system, the stator system and the support structure in the rotor blade tip clearance; taking the front end face of the first-stage rotor as a global reference transmission structure of an assembly structure, respectively expanding the transmission of the size and the deviation chain in the calculation model along two paths, finally calculating the coordinates of the matching points of the blade tips and the corresponding flow channel faces, and calculating the coordinate difference to form a blade tip clearance value; because the centroid of the rotor is not concentric with the center of the rotating shaft, the blade tips of the rotor are different from the gaps formed by the flow channel surfaces under different phases; and setting a rotation deviation angle value in the calculation model, and rotating the deviation angle value around the rotation axis from the initial phase of all the characteristics of the rotor. And after one circle of circumferential direction is completed, the mean value, the maximum value and the minimum value of the blade tip clearances are obtained for all the blade tip clearances of the rotors at all stages.
Examples
The prediction method for the rotor and stator assembly blade tip clearance of the aero-engine forms a calculation program, a mathematical calculation model is constructed by taking the rotor end face runout, radial eccentricity, blade tip runout, the end face runout and the elliptical deformation of a casing, and the bearing clearance and the elastic support deformation in a supporting structure as input, the rotor blade tip coordinate, the cross-sectional equation of a flow channel surface in the casing, and the maximum clearance, the minimum clearance and the average clearance of each stage of rotor detected in the assembly process.
From input to output, the data stream undergoes four stages of processing, feature calculation, phantom building, global assembly and gap calculation, and inherits attributes in the form of classes in the first two stages. The method is described by taking the rotor tip clearance of a certain engine as a research object and taking the assembly key parameters of the aero-engine components related to a specific rotor and a specific casing as an example, and the structural key characteristic parameters are shown in the following table:
TABLE 1 rotor and casing structure key characteristic parameter table (Unit: mm)
Radial height of rotor Radial height of flow passage section of casing Blade tip clearance Number of blades Axial coordinate
One-stage rotor 110 110.8 0.8 20 25
Two-stage rotor 102.5 103.2 0.7 20 75
Three-stage rotor 95 95.6 0.6 20 125
Four-stage rotor 87.5 88 0.5 20 175
Five-stage rotor 80 80.4 0.4 20 225
The global error propagation path includes: the specific transmission paths of the three parts of the rotor, the supporting frame and the stator casing are expanded from a global reference as shown in figure 1, and are transmitted in two paths to finally form a blade tip gap.
Step 1, calculating the coordinate of the rotor blade tip:
in the example, there is a five-stage rotor blade disk, the structure data of which is referred to as a table, and the error characteristic data reference table is fitted in detail with the description of the characteristic data after fitting because the error detection data amount is too large.
TABLE 2 rotor deviation data after fitting to partial data
Rear end face runout Maximum run-out phase angle Maximum radial eccentricity Eccentric phase angle Maximum tip runout
One-stage rotor 0.019mm 97° 0.003mm 135° 0.031mm
Two-stage rotor 0.030mm 167° 0.005mm 247° 0.026mm
Three-stage rotor 0.025mm 295° 0.002mm 67° 0.021mm
Four-stage rotor 0.011mm 206° 0.003mm 94° 0.022mm
Five-stage rotor 0.018mm 341° 0.002mm 301° 0.018mm
And (3) calculating the rotor deviation data after the partial data fitting processing by respectively substituting equations (1) to (8) in the data in the table, and finally giving the absolute coordinates of each blade tip point of the rotor by using the equation (8).
Step 2, calculating the supporting structure:
according to the actual measurement parameters of the bearing structure, the maximum radial clearance of the front bearing structure is 0.005mm, the axial clearance is 0.008mm, the maximum radial clearance of the rear bearing structure is 0.006mm, and the axial clearance is 0.007mm, and the rotor rotating shaft vector can be obtained by the formula (9). According to design and deformation of elastic supportsSimulation analysis result, highest point deformation delta of front support in z-axis directionsupf(-0.001mm,0,0.0), and the highest point deformation delta of the rear support in the z-axis directionsup,l＝(0.001mm,0,0.002mm)。
Step 3, calculating the flow channel surface in the casing:
the calculation of the deviation of the inner runner surface is illustrated by taking a half casing as an example. And (3) processing the initial point cloud data of the internal flow channel surface deviation to obtain characteristic data shown as a table, wherein the axial coordinate in the table is the x-axis coordinate in the coordinate system of the case in the step 3.
TABLE 3 part data casing deviation data after fitting process
Figure BDA0002375971430000221
Respectively substituting the processed deviation characteristic data of the casing into equations (10) - (13), calculating cubic spline interpolation functions of an ellipse characteristic equation and a centroid of each section, finally giving out matching points of the inner flow channel surface of the casing matched with the rotor blade tip by an equation (14), and substituting the characteristic relation in the step (4) for calculation.
And 4, calculating the rotor blade tip clearance:
the core of a rotor and stator assembly blade tip clearance calculation model is a rotation deflection angle thetarotSuch that the rotor mandrel rotates in space around the axis of revolution, at different phase angles of the rotor as shown in figure 2, presenting a spatially different law of variation from that of the mandrel of the stator casing.
Substituting the deviation into a model for calculation to obtain a space coordinate diagram of the rotor blade tip coordinate in a single section and the corresponding flow channel section in the casing as shown in FIG. 3. And (3) calculating the coordinates of the corresponding blade tip points by substituting the coordinates of the rotor blade tip points into the formula (16), and according to the transformation relation of the formulas (17) - (20). The rotor blade tips rotate around the revolving shaft, the blade tip clearance value of each rotor blade tip point at any phase position is calculated by the formula (21), and a box line graph shown in the figure 4 is given by a program according to the calculation result. In the figure, one unit represents the fluctuation condition of 20 blade tips of the single-stage rotor at 20 phases in the circumferential direction respectively, and the fluctuation condition has 400 tip clearance values, the red line represents the median value of the tip clearance, the upper end and the lower end of the dotted line represent the maximum value and the minimum value respectively, and the red "+" represents the abnormal value with overlarge deviation. Meanwhile, the maximum clearance, the minimum clearance and the average clearance of the blade tips are calculated by the equations (22) and (23), respectively, and the calculation results are shown in table 4.
TABLE 4 rotor blade tip clearance calculation results (unit: mm)
Figure BDA0002375971430000231
From the results shown in table 4 and fig. 4, the method disclosed by the invention can be used for accurately evaluating the blade tip clearance of the rotor-stator assembly of the aero-engine. Through the total 400 blade tip clearance values of 20 blade tips of the single-stage rotor under 20 different phases, the maximum clearance, the minimum clearance and the average clearance of each-stage rotor can be estimated, the fluctuation range of the predicted result of each-stage rotor is within +/-0.15 mm, and the requirements of empirical expectation and tolerance are met.
If the result of assembling the tip clearance value exceeds the assembling process requirement of the related blade, an assembly engineer can predict the result according to the method, position the out-of-tolerance part and provide parameters according to the out-of-tolerance value, so that the engineer reprocesses the corresponding part, such as: grinding the blade tip, turning the runner surface or adjusting the process parameters.
In the conventional aircraft engine assembly process, only the parts are detected to meet the tolerance requirements of the parts, and the deviation data of the parts are not associated, and a size error transfer model of the rotor and stator assembly blade tip clearance is not established. For the blade tip clearance, the traditional process only detects the part with good openness, most of the blade tip clearances with poor openness cannot be detected, and other technical means are needed for evaluating and assembling the blade tip clearance. With the continuous improvement of the performance requirements of the aero-engine, the assembly process needs to be changed to be precise and digital, and the traditional technical means cannot meet the requirements. The invention provides a prediction method for rotor and stator assembly blade tip clearance of an aeroengine, which relates to deviation data of various parts, establishes a mathematical model through a rotor and stator assembly relation, can provide a solution for the current situation that rotor and stator blade tip clearance is difficult to assemble and predict, helps an assembly engineer analyze accurate distribution of rotor and stator blade tip clearance, and provides a basis for optimizing an assembly process aiming at rotor and stator blade tip clearance.

Claims (10)

1. A prediction method for an aeroengine rotor-stator assembly blade tip clearance is characterized in that the aeroengine comprises a rotor system, a supporting frame, a casing and a stator assembly, the blade tip clearance refers to the radial clearance between a rotor blade tip and a flow passage surface in the casing and the radial clearance between a stator blade tip and a rotor hub, and the prediction method comprises the following steps:
step 1, establishing a rotor blade tip calculation model;
step 2, establishing a support structure calculation model;
step 3, establishing a flow channel surface calculation model in the casing;
and 4, integrating the three calculation models obtained in the steps 1 to 3 into a rotor and stator blade tip clearance calculation model by using a homogeneous coordinate transformation matrix under an absolute coordinate system, and predicting the blade tip clearances under different rotor phase angles.
2. The method for predicting the blade tip clearance of the rotor-stator assembly of the aircraft engine as claimed in claim 1, wherein in the step 1, a space coordinate system is established by taking the centroid of the axial reference plane of the primary rotor or the front journal as an origin; detecting a contour coordinate point of the rear end face of the first-stage rotor by using a part, establishing a space plane equation by using the coordinate point, and sequentially calculating rear end face plane equations of the rotors at all stages between the front support and the rear support; calculating unit normal vectors of rear end faces of the rotors at all levels by using a plane equation, namely completing geometric modeling of the end face run-out deviation; detecting the contour coordinate points of the outer side spoke plate surfaces of the front end surface and the rear end surface of the rotor by parts, establishing a circular equation in a cross section plane, estimating equation parameters by a least square method to obtain a centroid coordinate, namely completing radial eccentric deviation geometric modeling; calculating the direction vector of the inner centroid of the single-stage rotor according to the centroid coordinates of the front end face and the rear end face; and calculating the height and the initial phase angle of the blade tip of each blade so as to calculate the actual coordinates of the blade tip point, namely completing the rotor blade tip calculation model.
3. The method for predicting the blade tip clearance of the rotor-stator assembly of the aircraft engine as claimed in claim 2, wherein the specific calculation process of the step 1 is as follows:
establishing a coordinate system by taking a rotor datum plane as an assembly global datum, wherein an original point is a rotor datum plane center point, the axial stacking direction is the positive x direction, a cross-sectional plane is a yOz plane, and the vertical direction is the z-axis direction;
establishing a plane equation according to the profile data points of the rear end face of the ith-stage rotor acquired by a detection instrument, and fitting the parameters of the plane equation by a least square method to obtain an equation expression:
Apla,ix+Bpla,iy+Cpla,iz+Dpla,i＝0 (1)，
in the formula, Apla,Bpla,Cpla,DplaAll the coefficients are calibrated plane parameter coefficients, and a unit normal vector of the rear end face of the i-th-stage rotor is obtained by performing normalization processing on a plane normal vector obtained by the formula (1):
Figure FDA0002375971420000021
the method comprises the following steps of obtaining a section plane outer ring profile of a radial plate on the outer side of the front end face of the single-stage rotor at the axial x position by a rotary table detection instrument, and fitting according to a least square method to obtain a front end face outer side circular equation:
Acir,f,iy2+Bcir,f,iz2+Ccir,f,iy+Dcir,f,iz+Ecir,f,i＝0 (3)，
in the above formula, Acir,f,i,Bcir,f,i,Ccir,f,i,Dcir,f,i,Ecir,f,iCircle characteristic parameters of the front end face of the ith-stage rotor, which are calibrated by a least square method, are approximately estimated by a section plane circular equation, and the centroid position vector of the front end face of the ith-stage rotor is as follows:
Figure FDA0002375971420000022
in the formula (I), the compound is shown in the specification,
Figure FDA0002375971420000023
the approximate centroid position vector (x) of the rear end face can be obtained by the same principle of the formula (3) and the formula (4) for the accumulated coordinates of the axial dimension of the front i-1-stage rotorcir,l,ycir,l,zcir,l) Component nx,iThe unit direction vector of the ith stage rotor form axis is given by:
Figure FDA0002375971420000024
in the formula (I), the compound is shown in the specification,
Figure FDA0002375971420000025
is the die length of the axial vector of the ith-stage rotor inner form;
let the axial positioning dimension of the ith-stage rotor blade tip be liDesigned radial height dimension RiLet the total number of i-th-stage rotor blades be NiThe circumferential included angle of the adjacent blades is
Figure FDA0002375971420000031
Taking the positive direction of a z axis as an initial edge of a phase angle in a yOz plane, and numbering the circumferential blades from 1 to j in an increasing way along the positive direction of the phase angle, wherein j is more than 0 and is not more than NiIn which the phase angle of blade number 1 is phii,1The phase angle of the j blade is phii,j＝φi,1+(j-1)φi,Δ；
The unit direction vector from the centroid to the zero phase point of the blade tip in the section plane of the blade tip of the ith-stage rotor is set as
Figure FDA0002375971420000032
It satisfies the mutually perpendicular condition in space:
Figure FDA0002375971420000033
then the vector winding mandrel is formed from the cross section axis to the zero phase point of the corresponding unit direction vector of the j blade tip of the ith-stage rotor
Figure FDA0002375971420000034
Is rotated to obtain
Figure FDA0002375971420000035
Satisfies the formula:
Figure FDA0002375971420000036
the j-th blade tip of the ith-stage rotor has a designed height RiRun-out of deltaj(ii) a Thus, the actual height of the corresponding blade is Ri,jTherefore, the absolute coordinates of the j-th blade tip of the ith-stage rotor are as follows:
Figure FDA0002375971420000037
4. the method for predicting the rotor-stator assembly tip clearance of the aircraft engine according to claim 1 or 2, wherein in the step 2, a coordinate system is established by taking a centroid of a left end face of the front bearing inner ring matched with the rotor as an origin, and a connecting line of centroids of adjacent end faces of the front bearing inner ring and the rear bearing inner ring is a rotor rotation axis; calculating a unit vector of a rotor rotating axis, namely completing modeling of a rotor rotating shaft; and substituting the simulated deformation coordinate component of the elastic support and the bearing clearance value, and calculating a center offset vector of the reference surface of the casing under the bearing clearance deviation, namely completing the modeling of the supporting frame.
5. The method for predicting the blade tip clearance of the rotor-stator assembly of the aircraft engine as claimed in claim 4, wherein the specific calculation process of the step 2 is as follows: two bearings are arranged at the front and the back of the rotor structure, the bearings have radial play and axial play deviation, the left end face of the inner ring of the bearing is taken as a yOz plane, the centroid of the end face is taken as a coordinate origin O, the vertical direction is the z-axis direction, and the direction of the centroid towards the inner normal vector is the x direction; axial clearance of bearing is Deltaxcle,(xcle,min≤Δxcle≤xcle,max) Radial play of bearing of Deltarcle,(rcle,min≤Δrcle≤rcle,max)；
Respectively measuring the arc profiles of the rear end surfaces of the inner ring and the outer ring of the front bearing under a global coordinate system established by the rotor reference surface, and calculating to obtain the centroid coordinates of (x)f，in,cle,yf,in,cle,zf,in,cle) And (x)f,out,cle,yf,out,cle,zf,out,cle) (ii) a Similarly, the circular arc profiles of the front end surfaces of the inner ring and the outer ring of the rear bearing are respectively measured, and the centroid coordinates are calculated to be (x)l,in,cle,yl,in,cle,zl,in,cle) And (x)l,out,cle,yl,out,cle,zl,out,cle) And then the unit direction vector of the inner ring axis of the bearing is as follows:
Figure FDA0002375971420000041
the unit direction vector of the outer ring axis can be obtained by the same method
Figure FDA0002375971420000042
According to the structure simulation result of the elastic support, the deviation position vector of the reference plane of the casing caused by the deformation of the front elastic support and the rear elastic support is (delta x)sup,f,δysup,f,δzsup,f) And (δ x)sup,l,δysup,l,δzsup,l)。
6. The method for predicting the blade tip clearance in the rotor-stator assembly of the aircraft engine as claimed in claim 1 or 2, wherein in the step 3, a coordinate system is established by taking the centroid of the left end face of the casing as an origin, the contour point of the matching part of the inner flow channel and the blade tip of the rotor is measured, and the axial position and the coordinates of the contour point are recorded as model input; establishing a contour elliptic equation, calibrating parameters by using a least square method according to contour coordinates, and extracting the phases of the center, the long half shaft, the short half shaft and the long shaft of the elliptic equation; and calculating an interpolation centroid, substituting different phase flow channel surface contour jumping values in an elliptic equation, and establishing an inner channel surface contour equation to finish the casing calculation model.
7. The method for predicting the blade tip clearance of the rotor-stator assembly of the aircraft engine as claimed in claim 6, wherein the specific establishment process of the casing calculation model is as follows:
establishing an inner runner surface elliptical equation of the single-section case, taking the centroid of the reference end surface of the case as the origin of coordinates, taking the axial direction as the x direction, setting the positive direction to be the same as the stacking direction of the rotors, setting a coordinate system for a section plane which is a yOz plane, and taking the vertical direction as the z-axis direction, and then, taking the inner runner surface elliptical equation matched with the i-th-stage rotor as the equation,
y2+Acas,iyz+Bcas,iz2+Ccas,iy+Dcas,iz+Ecas,i＝0 (10)，
in the formula, Acas,i,Bcas,i,Ccas,i,Dcas,i,Ecas,iFitting ellipse with least square method to obtain the centroid coordinate (x) of the flow channel cross sectioncas,i,ycas,i,zcas,i) The non-negative minimum phase angle of the major and minor axes of the elliptic equation is thetacas,iThe following conditions are satisfied:
Figure FDA0002375971420000051
when the casing assembly structure is horizontally placed, the relative sinking amount of the casing assembly structure in the global coordinate system of the rotor reference is
Figure FDA0002375971420000052
The centroid coordinate of the front datum plane of the casing is (x) in combination with the position deviation of the casing caused by the deformation of the supportcas,ben,ycas,ben,zcas,ben)；
According to the design shape of the flow channel surface in the casing, the change equation of the section profile radius of the flow channel surface in the casing along the axial direction is RcaR (x); and then the centroid curve is calculated.
8. The method for predicting the blade tip clearance in the rotor-stator assembly of the aircraft engine as claimed in claim 7, wherein when the casing is a split casing, the centroid curve calculation mode is as follows:
establishing cubic spline interpolation functions of coordinate components y and z relative to x for bisection of the centroid of the flow channel surface in the half casing respectively, wherein the interpolation node matrix formula is as follows:
Figure FDA0002375971420000061
in the formula, hiFor the step size of the node,
Figure FDA0002375971420000062
for interpolated node vectors, f [ x ]i-1,xi,xi+1]For third order difference, constraint m0,mnThe component y of the end face attitude vector of the left end face and the component z of the end face attitude vector of the right end face of the part can be calculated about the vector included angle of the axial direction x, and then the interpolation function is as follows:
Figure FDA0002375971420000063
wherein s (x) denotes ycas,c(x) Or zcas,c(x) Calculating to obtain (x)cas,c,ycas,c,zcas,c)。
9. The method for predicting the blade tip clearance in the rotor-stator assembly of the aircraft engine as claimed in claim 7, wherein when the casing is an integral thin-walled annular casing, the centroid curve calculation mode is as follows:
for the integral thin-wall annular casing, the elliptical profile equation of the inner runner of the front end surface and the inner runner of the rear end surface of the ith section of the casing are given by the formula (10), and after undetermined parameters are evaluated by a least square method, the central coordinate of the inner ellipse of the front end surface and the deflection angle y of the long axis are calculated according to the formula (11)cas,f,i,zcas,f,i,θcas,f,iAnd the central coordinate of the ellipse and the deviation angle y of the long axis in the section plane of the rear end facecas,l,i,zcas,l,i,θcas,l,iBefore the ith section of the casingThe end surface internal flow channel profile equation has the following coordinates of any centroid points inside the casing:
Figure FDA0002375971420000064
after the multiple sections of casings are stacked, the centroid coordinates of the front end face of the ith section of casing are calculated in the same way as in the formula (4):
Figure FDA0002375971420000071
in the formula (I), the compound is shown in the specification,
Figure FDA0002375971420000072
is a unit normal vector of the rear end face of the ith section of the casing, and is calculated by formula (1) and formula (2), Ecas,f,iCalibrating the section parameter of the front end surface ellipse of the ith section of the casing by a least square method;
then the centroid coordinate of any point in the ith section of the casing is the centroid coordinate (x) of the front end face and the rear end face of the single section of the casingcas,f,i,ycas,f,i,zcas,f,i) And (x)cas,l,i,ycas,l,i,zcas,l,i) Substituting formula (5) to obtain (x)cas,c,ycas,c,zcas,c)；
The total collection point number for measuring the radius value of the section of the flow channel matched with the blade tip of the ith-stage rotor is Ni,casAngle between collection points
Figure FDA0002375971420000073
The measured data can obtain the runout delta of the point to be measured of the runner surface relative to the design sizei,j,casEach data available axial position xi,casAnd phase angle phii,j,cas＝j·Δφi,casA unique representation; the coordinates of the inner flow channel surface corresponding to the axial arbitrary position and the phase angle of the casing are as follows:
Figure FDA0002375971420000074
in the formula (I), the compound is shown in the specification,
Figure FDA0002375971420000075
[g]is the rounding function, i.e. the largest integer not exceeding the value of g.
10. The method for predicting the blade tip clearance of the rotor-stator assembly of the aircraft engine as claimed in claim 1 or 2, wherein the specific calculation process of the step 4 is as follows:
under an absolute coordinate system, when the rotor does not generate phase deflection, the section of a flow channel matched with the j-th blade tip of the ith-stage rotor meets the relation:
Figure FDA0002375971420000076
substituting the condition satisfied by the formula (17) into the formula (16), and calculating the coordinate of the flow channel point matched with the formula; the assembly blade tip clearance cle of the i-th-stage rotor and the j-th blade in the non-rotating state with the blade tip rotating angle of 0 DEGi,jExpressed as:
Figure FDA0002375971420000081
when the rotor rotating shaft deflects by a phase angle thetarotIn the meantime, an initial blade tip position vector of the j blade of the i-th-stage rotor is set
Figure FDA0002375971420000082
Rotated tip position vector
Figure FDA0002375971420000083
And blade direction vector
Figure FDA0002375971420000084
Respectively satisfy the following formula:
Figure FDA0002375971420000085
then the current phase angle thetarotUnder the condition, the corresponding point of a runner surface matched with the j-th stage rotor blade meets the relation:
Figure FDA0002375971420000086
calculating the coordinate of the corresponding point of the runner surface matched with the j-th blade of the i-th-stage rotor by using the conditional formula (20) back to the formula (16), and then, the blade tip clearance
Figure FDA0002375971420000087
Expressed as:
Figure FDA0002375971420000088
the maximum clearance cle of the mth stage rotorm,maxMinimum clearance clem,minRespectively as follows:
Figure FDA0002375971420000089
the number of times of measurement of the ith-stage rotor rotating for one circle at a certain rotation angle is set to be Nrot,iThen the phase angle of a single rotation is
Figure FDA00023759714200000810
Wherein, k is the corner number, the average clearance cle of the mth stage rotorm,meanComprises the following steps:
Figure FDA0002375971420000091
CN202010065955.0A 2020-01-20 2020-01-20 Prediction method for blade tip clearance of rotor and stator assembly of aircraft engine Active CN111241609B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010065955.0A CN111241609B (en) 2020-01-20 2020-01-20 Prediction method for blade tip clearance of rotor and stator assembly of aircraft engine

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010065955.0A CN111241609B (en) 2020-01-20 2020-01-20 Prediction method for blade tip clearance of rotor and stator assembly of aircraft engine

Publications (2)

Publication Number Publication Date
CN111241609A true CN111241609A (en) 2020-06-05
CN111241609B CN111241609B (en) 2022-04-01

Family

ID=70864230

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010065955.0A Active CN111241609B (en) 2020-01-20 2020-01-20 Prediction method for blade tip clearance of rotor and stator assembly of aircraft engine

Country Status (1)

Country Link
CN (1) CN111241609B (en)

Cited By (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111609826A (en) * 2020-06-08 2020-09-01 清华大学 Method and device for measuring blade tip clearance of aircraft engine
CN112417593A (en) * 2020-11-18 2021-02-26 西北工业大学 Optimization method for connection assembly mounting angle of circular arc end teeth of aircraft engine
CN112880617A (en) * 2021-01-12 2021-06-01 西北工业大学 Free gap real-time monitoring method based on spectrum estimation
CN112989720A (en) * 2021-03-04 2021-06-18 上海交通大学 Method for realizing optimization of radial clearance of hub by using cantilever stator blades of axial-flow compressor
CN113158334A (en) * 2021-04-02 2021-07-23 北京航空航天大学 Digital false mounting quality prediction method for aero-engine case
CN113343357A (en) * 2021-06-23 2021-09-03 中国航发沈阳发动机研究所 Method for calculating radial clearance between high-pressure rotor and low-pressure rotor of aircraft engine
CN114510798A (en) * 2022-01-25 2022-05-17 中国航发沈阳发动机研究所 Pneumatic axial force error-proofing analysis method for aircraft engine
CN114627177A (en) * 2022-03-25 2022-06-14 南京航空航天大学 Aircraft skin gap and step difference measuring method based on image segmentation
CN114763990A (en) * 2021-01-15 2022-07-19 中国航发商用航空发动机有限责任公司 Method for measuring axial size of rotor and stator inner cavity
CN114812417A (en) * 2022-04-19 2022-07-29 天津大学 Rotor and stator gap error compensation method and device based on rotor position synchronization
CN115235404A (en) * 2022-09-22 2022-10-25 中国航发燃气轮机有限公司 Rotor and stator clearance correction method for heavy gas turbine

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104075677A (en) * 2014-06-30 2014-10-01 东北大学 Method and system for measuring concentricity and tip clearance of rotor and stator of aero-engine
US20150169811A1 (en) * 2013-12-18 2015-06-18 General Electric Company System and method for estimation and control of clearance in a turbo machine
CN106768302A (en) * 2017-02-08 2017-05-31 中国航发沈阳发动机研究所 A kind of critical vibration mode test method of aeroengine rotor

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20150169811A1 (en) * 2013-12-18 2015-06-18 General Electric Company System and method for estimation and control of clearance in a turbo machine
CN104075677A (en) * 2014-06-30 2014-10-01 东北大学 Method and system for measuring concentricity and tip clearance of rotor and stator of aero-engine
CN106768302A (en) * 2017-02-08 2017-05-31 中国航发沈阳发动机研究所 A kind of critical vibration mode test method of aeroengine rotor

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
JIA BINGHUI 等: "The development of aero-engine tip-clearance measurement technology: A simple review", 《2017 13TH IEEE INTERNATIONAL CONFERENCE ON ELECTRONIC MEASUREMENT & INSTRUMENTS (ICEMI)》 *
NKFMRD36: "航空发动机整机振动中的转静耦合机理分析及应用研究", 《豆丁》 *
丁健: "航空发动机叶尖间隙分析系统设计", 《中国优秀博硕士学位论文全文数据库(硕士)工程科技Ⅱ辑(月刊)》 *

Cited By (18)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111609826A (en) * 2020-06-08 2020-09-01 清华大学 Method and device for measuring blade tip clearance of aircraft engine
CN112417593A (en) * 2020-11-18 2021-02-26 西北工业大学 Optimization method for connection assembly mounting angle of circular arc end teeth of aircraft engine
CN112417593B (en) * 2020-11-18 2024-02-23 西北工业大学 Optimization method for connection assembly installation angle of arc end teeth of aero-engine
CN112880617A (en) * 2021-01-12 2021-06-01 西北工业大学 Free gap real-time monitoring method based on spectrum estimation
CN114763990A (en) * 2021-01-15 2022-07-19 中国航发商用航空发动机有限责任公司 Method for measuring axial size of rotor and stator inner cavity
CN114763990B (en) * 2021-01-15 2023-10-27 中国航发商用航空发动机有限责任公司 Method for measuring axial dimension of inner cavity of rotor and stator
CN112989720A (en) * 2021-03-04 2021-06-18 上海交通大学 Method for realizing optimization of radial clearance of hub by using cantilever stator blades of axial-flow compressor
CN112989720B (en) * 2021-03-04 2022-05-06 上海交通大学 Method for realizing optimization of radial clearance of hub by using cantilever stator blades of axial-flow compressor
CN113158334B (en) * 2021-04-02 2022-11-11 北京航空航天大学 Digital false mounting quality prediction method for aero-engine case
CN113158334A (en) * 2021-04-02 2021-07-23 北京航空航天大学 Digital false mounting quality prediction method for aero-engine case
CN113343357B (en) * 2021-06-23 2023-08-18 中国航发沈阳发动机研究所 Method for calculating radial clearance between high-low pressure rotors of aero-engine
CN113343357A (en) * 2021-06-23 2021-09-03 中国航发沈阳发动机研究所 Method for calculating radial clearance between high-pressure rotor and low-pressure rotor of aircraft engine
CN114510798A (en) * 2022-01-25 2022-05-17 中国航发沈阳发动机研究所 Pneumatic axial force error-proofing analysis method for aircraft engine
CN114510798B (en) * 2022-01-25 2023-07-21 中国航发沈阳发动机研究所 Aeroengine pneumatic axial force error-proofing analysis method
CN114627177A (en) * 2022-03-25 2022-06-14 南京航空航天大学 Aircraft skin gap and step difference measuring method based on image segmentation
CN114812417A (en) * 2022-04-19 2022-07-29 天津大学 Rotor and stator gap error compensation method and device based on rotor position synchronization
CN115235404A (en) * 2022-09-22 2022-10-25 中国航发燃气轮机有限公司 Rotor and stator clearance correction method for heavy gas turbine
CN115235404B (en) * 2022-09-22 2022-12-06 中国航发燃气轮机有限公司 Rotor and stator clearance correction method for heavy gas turbine

Also Published As

Publication number Publication date
CN111241609B (en) 2022-04-01

Similar Documents

Publication Publication Date Title
CN111241609B (en) Prediction method for blade tip clearance of rotor and stator assembly of aircraft engine
CN111460677B (en) Method for establishing rotor stacking precision prediction model based on geometric algebraic theory
CN109117460B (en) Method for calculating rotor assembly axis deflection based on end jump measurement
CN106446343B (en) Automatic extraction method for parameterized profile of straight-line blade of radial-flow impeller
CN106837426B (en) Method for optimizing eccentricity of mass center of rotor of engine core machine
CN111141518B (en) Model-based asymmetric rotor bearing system unbalance identification method
JP2011508153A (en) Rotor assembly system and method
CN104246635A (en) Method for machining the trailing edge of a turbine engine blade
JP7416674B2 (en) Turbine assembly method, turbine assembly support program, and turbine assembly support device
Huibin et al. A tip clearance prediction model for multistage rotors and stators in aero-engines
CN107908914B (en) Method for judging machinability of closed impeller of centrifugal compressor and calculating intermediate section
CN103995938B (en) A kind of compressor Middle casing end wall profile formative method
CN113343357A (en) Method for calculating radial clearance between high-pressure rotor and low-pressure rotor of aircraft engine
CN110177919B (en) Adaptive machining of cooled turbine airfoils
CN109960870B (en) Large-scale high-speed rotation equipment multi-stage part rigidity prediction method based on contact area maximization measurement and adjustment
CN110686640A (en) Method for measuring clearance between rotor part and stator part of centrifugal compressor
CN111076867B (en) Large-scale high-speed rotation equipment multistage part unbalance amount distribution method based on synchronous measurement and adjustment of mass center and inertia center
CN116595874A (en) Impeller mechanical performance prediction model parameter optimization method and device and storage medium
CN111241752A (en) Centrifugal impeller comprehensive optimization method based on digital twinning and reinforcement learning
CN113468788B (en) Deformation coordination design method for engine blisk
CN115164809A (en) Six-axis measurement and error compensation method for blisk
CN114510798A (en) Pneumatic axial force error-proofing analysis method for aircraft engine
Guo et al. Construction of digital twin for clamped near-net-shape blade in adaptive manufacturing
CN113761800B (en) Shafting dynamics parameter model scaling design method based on critical rotation speed correspondence
CN114763990B (en) Method for measuring axial dimension of inner cavity of rotor and stator

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant