A method for predicting tip clearance with cross-scale multi-physics coupling

Abstract

The present invention discloses a method for predicting tip clearance with cross-scale multi-physical field coupling. During the tip clearance prediction process, the heat conduction in the solid domain and the parameters of the main flow path and the secondary air system are all unsteady. The present invention combines the one-dimensional fluid network method of the secondary air system, the two-dimensional axisymmetric finite element method, the engineering algorithms for blade temperature and displacement, and the main flow path volume method. By establishing a cross-scale multi-physical field coupling calculation model and adopting a real-time information interaction method of two-way coupling of gas, heat, and solid across time scales and weak coupling of the primary and secondary flows at the same time scale, accurate prediction of the tip clearance considering the changes in the main flow path parameters, the secondary air system parameters, and solid deformation is achieved. The present invention combines the advantages of one-dimensional clearance prediction and three-dimensional clearance simulation, achieving a balance between calculation accuracy and efficiency.

Description

Technical Field

[0001] The present invention relates to a method for predicting tip clearance with cross-scale multi-physical field coupling, belonging to the technical field of overall structure design of aeroengines. Background Art

[0002] The tip clearance is an important content in the design of an aero-gas turbine engine, which directly affects the operation safety and performance of the engine. A reasonable hot-state clearance can improve the performance of the engine on the premise of ensuring no rubbing. Therefore, a method for predicting tip clearance needs to be established during the engine design.

[0003] Research shows that accurate prediction and analysis of the tip clearance must consider the fluid-thermal-solid coupling effect between the gas-thermal environment and the solid domain. Considering the modeling and calculation costs comprehensively, the analysis considering the fluid-thermal-solid coupling effect in engineering often adopts a calculation method that couples a one-dimensional air system fluid network model and a two-dimensional axisymmetric finite element model. In the exploration of steady-state problems, Muller et al. carried out fluid-thermal-solid coupling calculations for a low-pressure turbine rotor. During the process, the boundary conditions of the main flow path and the air system boundary were given, and the interaction between the air system, heat transfer, and structural deformation was considered. Comparing the decoupled calculation results, the difference in the air system flow rate obtained by considering the fluid-thermal-solid coupling reached more than 90%, and the difference in the solid temperature field reached more than 10K. Accurate prediction of the transient clearance of an aeroengine is a more critical issue in refined design. Ganine et al. carried out transient fluid-thermal-solid coupling calculations for a gas turbine using a strong coupling method. During the calculation, the flow boundaries of the secondary air system and the main flow path were simplified as functions of time for input. Comparing the decoupled calculation results, it was found that there were significant differences in the variation laws of the fluid and solid temperatures in the transition state. Guo Naxian used the "cross-time-scale fluid-thermal-solid coupling" calculation method to explore the transition state change law of the tip clearance of a certain type of core engine. During the calculation, the linear interpolation boundary of the steady-state parameters in different states of the main flow path was adopted. The calculation results better captured the tip clearance change characteristics during the whole process.

[0004] In the existing simulation studies related to tip clearance or solid domain deformation mentioned above, the boundaries of the air system and the main flow path are both given by considering them as functions of time or interpolating based on the steady-state parameters under different engine operating conditions, without considering the coupled effects caused by the dynamic interaction among the clearance change, main flow path parameters, and air system parameters. During the engine operation, the main flow path parameters are affected by parameters such as tip clearance and the air extraction volume of the secondary air system, and do not show a simple linear relationship with the rotational speed. Moreover, during the engine R & D stage, it is also impossible to obtain the functions of the secondary air system boundary and the relevant cross-sections of the main flow path with respect to time through experimental measurements. If the interpolated steady-state parameters obtained by calculation are used as boundary conditions, it will cause certain deviations in the calculation and analysis of tip clearance, especially the analysis of thermal deformation. Therefore, some scholars have conducted relevant research on the dynamic interaction among tip clearance, secondary air system, and main flow path parameters.

[0005] Regarding the research on the interactive effects of tip clearance and main flow path parameters at the level of the whole engine or core engine, Li et al., Chapman et al., Wang Jiamei et al., and Sheng et al. all adopted the method of coupling calculation by bringing a one-dimensional simplified model into the engine performance model to achieve the dynamic interaction between tip clearance and main flow path parameters. The above methods can be summarized as follows: after each step of performance calculation, the parameters at the inlet section of the high-pressure turbine rotor in the main flow path are used as the boundary of the main flow path, and the temperature of the secondary flow cooling gas is predicted according to the empirical formula and the heat transfer coefficient at the corresponding position is calculated, and then the tip clearance is calculated. After the calculation is completed, the newly obtained clearance is used to modify the component characteristics at the next moment, and then the main flow path parameters are calculated, and so on. However, although this method can ensure the calculation efficiency, the tip clearance model it uses is too simplified. Although it can reflect the change trend of tip clearance, the calculation accuracy is relatively limited.

[0006] Regarding the coupling of primary and secondary flows, for the dynamic interaction between the secondary air system and the main flow path parameters, sufficient research has been carried out both at home and abroad. Liu Chuankai et al. constructed an integrated primary-secondary flow solution model using a modular approach and solved it using the N + 1 residual method to analyze the dynamic interaction mechanism and its effects between the air system and the gas path of the main flow path during the transient process. Nikolaidis et al. constructed an integrated primary-secondary flow solution model considering the rotational effect and cavity effect, and carried out calculations for normal operating conditions and strong transient operating conditions respectively, believing that the quasi-steady method can be used for air system analysis under normal operating conditions. However, the influence of the deformation of the engine solid structure, including changes in tip clearance, labyrinth clearance, etc., on the fluid domain has not been considered in this type of whole-engine flow simulation, which brings certain limitations to the final calculation results.

[0007] In summary, during the current analysis of the transitional tip clearance, the dynamic coupling interactions among tip clearance changes, the mainstream, and the secondary air system have not been considered simultaneously, which limits the prediction accuracy of the engine tip clearance. Therefore, it is of great significance to establish a tip clearance prediction method that couples multi-physical fields across scales with a two-dimensional model as the main body. Summary of the Invention

[0008] The technical problem to be solved by the present invention is to provide a tip clearance prediction method that couples multi-physical fields across scales. The fluid-structure weak coupling method is adopted, a multi-dimensional calculation model is established during the calculation process, and the bidirectional transfer of information between physical fields is considered, so as to improve the calculation accuracy of the tip clearance while reducing the calculation time.

[0009] The present invention adopts the following technical solutions to solve the above technical problems:

[0010] A tip clearance prediction method that couples multi-physical fields across scales, and the method is specifically as follows:

[0011] Establish a one-dimensional network model of the engine air system, a two-dimensional axisymmetric model of the disk and the casing, a parametric model of the blade, and a quasi-one-dimensional model of the main flow path;

[0012] Use the volume method to solve the quasi-one-dimensional model of the main flow path, and at the same time use the one-dimensional transient fluid network method to solve the one-dimensional network model of the air system to obtain the heat transfer boundary between the engine disk and the casing;

[0013] Considering the thermal load and centrifugal load in the unsteady heat conduction process of the engine solid domain, use the two-dimensional axisymmetric finite element method to solve the two-dimensional axisymmetric model of the disk and the casing to obtain the radial displacements of the disk and the casing;

[0014] Use the engineering algorithms of temperature and displacement to solve the parametric model of the blade to obtain the temperature and radial displacement of the blade;

[0015] Integrate the radial displacements of the disk, the casing and the blade and the initial tip clearance value to obtain the tip clearance value during the transient process;

[0016] In the process of solving the radial displacements of the disk and the casing above, a weak fluid-structure coupling method is adopted for two-way information interaction between the fluid domain and the solid domain, that is: the heat transfer boundary conditions are transferred from the main flow channel and the air system flow field to the solid domain temperature field, the thermal load and the qualitative temperature of the engine material properties are transferred from the solid domain temperature field to the solid domain displacement field, and the geometric parameters of the throttling elements in the air system flow field are updated by using the solid domain displacement field; conversely, the change of the grid node coordinates of the engine solid domain is provided from the solid domain displacement field to the solid domain temperature field, the change of the solid domain wall temperature and the labyrinth clearance is provided from the solid domain temperature field to the air system flow field, and the change of the tip clearance is provided from the solid domain displacement field to the main flow channel, which is used as the basis for correcting the main flow channel parameters.

[0017] As a preferred embodiment of the present invention, for the radial displacements of the disk and the casing, a weak fluid-structure coupling method is adopted for two-way information interaction between the fluid domain and the solid domain. During a solid domain time step Δt s multiple fluid domain calculations are carried out, during which it is considered that the solid wall temperature remains unchanged. Δt s = n·Δt f where Δt f is the fluid domain time step, specifically:

[0018] 1) During the current fluid domain time step of the current solid domain time step, the tip clearance calculated in the previous solid domain time step is fed back to the main flow channel, and the total temperature, total pressure and flow rate of the fluid at each section of the engine main flow channel, the engine speed, and the pressures and temperatures at the inlet and outlet of the air system are calculated according to the volume method;

[0019] 2) During the same fluid domain time step as 1), the labyrinth clearance calculated in the previous solid domain time step is fed back to the air system, and the flow rate, pressure and temperature of the fluid in the internal components of the air system are calculated by using the one-dimensional transient fluid network method according to the pressures and temperatures at the inlet and outlet of the air system calculated in 1);

[0020] 3) During the next fluid domain time step of the current solid domain time step, 1) and 2) are repeated until the iterative calculations for all fluid domain time steps within the current solid domain time step are completed, and then enter 4);

[0021] 4) Based on the results calculated in the last fluid domain time step of the current solid domain time step, the third type of heat transfer boundary conditions are converted, and the transient thermal analysis of the engine solid domain is carried out to obtain the solid domain temperature field;

[0022] 5) The centrifugal load of the solid domain is determined according to the rotational speed calculated in the last fluid domain time step of the current solid domain time step, and the solid domain temperature field calculated in 4) is used as the temperature load to solve the deformation amount of the solid domain, and the deformation amounts of the disk and the casing and the deformation amount of the labyrinth are output.

[0023] 6) Obtain the radial displacements of the disk and the casing based on the deformation amounts of the disk and the casing, and obtain the tip clearance according to the radial displacements of the disk and the casing, the radial displacement of the blade, and the initial tip clearance value; obtain the labyrinth clearance according to the deformation amount of the labyrinth.

[0024] 7) Update the two-dimensional axisymmetric model using the deformation amounts of the disk and the casing, and return to 1) to enter the calculation of the next solid domain time step.

[0025] As a preferred solution of the present invention, the specific solution process of the radial displacement of the blade is as follows:

[0026] (1) Divide the blade body radially into several segments according to the blade height, extract the geometric information of each cross-section, and establish a parametric model of the blade.

[0027] (2) Establish the energy conservation equation for each segment as follows:

[0028] Q = hA1ΔT1 + λA2ΔT2 - λA3ΔT3

[0029] QΔt s = c p mΔT

[0030] where Q is the net heat flux of the blade segment, h is the heat transfer coefficient between the gas flow and the wall surface, A1 is the gas-solid convection heat transfer area, ΔT1 is the heat transfer temperature difference between the gas and the solid, λ is the thermal conductivity of the blade segment, A2 is the end face area for introducing heat into the blade segment, ΔT2 is the average temperature difference between the upper and lower segments of the end face for introducing heat, A3 is the end face area for exporting heat from the blade segment, ΔT3 is the average temperature difference between the upper and lower segments of the end face for exporting heat, Δt s is the solid domain time step, c p is the specific heat capacity of the blade segment, m is the mass of the blade segment, and ΔT is the change in the average temperature of the blade segment within a single solid domain time step;

[0031] (3) Calculate the thermal load of the blade segment from (2), and superimpose the thermal displacements of each segment as follows:

[0032]

[0033] where ΔR is the overall thermal displacement of the blade, α is the linear expansion coefficient of the blade segment, ΔT j is the temperature change of the j-th segment, Δl j is the length of the j-th segment, and n is the number of segments;

[0034] Superimpose the displacements caused by the centrifugal loads of each blade segment as follows:

[0035]

[0036] Among them, ΔL is the overall centrifugal displacement of the blade, and σ j is the centrifugal stress borne by the j-th segment, and E is the elastic modulus of the j-th segment;

[0037] The total radial displacement of the blade is obtained by superimposing and summing the overall thermal displacement and the overall centrifugal displacement of the blade.

[0038] As a preferred solution of the present invention, in the above (2), during the solution process of the one-dimensional transient fluid network method, the volume effect of the cavity is considered, and the relationship between the inlet and outlet flow rates and pressures is as follows:

[0039]

[0040] Among them, m in is the inlet flow rate of the cavity, m out is the outlet flow rate of the cavity, A in is the inlet cross-sectional area of the cavity, A out is the outlet cross-sectional area of the cavity, is the total pressure of the inlet air flow, R is the ideal gas constant, is the total temperature of the inlet air flow, k is the adiabatic index, P0 is the static pressure in the middle of the cavity, C D is the flow rate correction coefficient, P r is the static pressure at the outer edge of the cavity, is the average static temperature of the air flow in the cavity, P out is the static pressure of the air flow at the outlet of the cavity.

[0041] As a preferred solution of the present invention, in the above (4), during the transient thermal analysis of the solid domain of the engine, the two-dimensional axisymmetric finite element method is used for unsteady heat conduction calculation. There is no internal heat source in the calculation model, and the thermal conductivity of the solid domain is isotropic. Therefore, the unsteady heat conduction differential equation in the cylindrical coordinate system is as follows:

[0042]

[0043] Among them, ρ is the density of the solid domain, c is the specific heat capacity at constant pressure of the solid domain, T is the temperature of the solid domain, τ is the time, λ is the thermal conductivity of the solid domain, r is the radial coordinate, and z is the axial coordinate.

[0044] As a preferred solution of the present invention, in the above (1), when using the volume method for calculation, the change rates of the inlet and outlet parameters of the combustion chamber in the main flow channel are:

[0045]

[0046] Among them, T3 is the outlet temperature of the combustion chamber, p3 is the outlet pressure of the combustion chamber, R is the gas constant, V3 is the volume of the main combustion chamber, C v is the constant volume heat capacity, gto is the fuel supply amount, H uis the fuel calorific value, η3 is the main combustion efficiency, i c is the enthalpy value carried by the fuel when it enters the combustion chamber, i3 is the unit enthalpy at the combustion chamber outlet, k3 is the heat capacity ratio, g2 is the flow rate at the combustion chamber inlet, i2 is the unit enthalpy at the combustion chamber inlet, and g3 is the flow rate at the combustion chamber outlet

[0047] A computer device includes a memory, a processor, and a computer program stored in the memory and capable of running on the processor. When the processor executes the computer program, the steps of the cross-scale multi-physical field coupling tip clearance prediction method are implemented.

[0048] A computer-readable storage medium stores a computer program. When the computer program is executed by a processor, the steps of the cross-scale multi-physical field coupling tip clearance prediction method are implemented.

[0049] Compared with the prior art, the present invention adopts the above technical solutions and has the following technical effects:

[0050] The present invention combines the one-dimensional fluid network method, the two-dimensional axisymmetric finite element method, the engineering algorithms for temperature and displacement, and the main flow path volume method, and uses a weak coupling method to achieve two-way interaction of parameters between the fluid domain and the solid domain within each coupling iteration step. The displacements of the disk and the casing are obtained through coupled calculations, the displacements of the main flow blades are obtained through engineering algorithms, and the tip clearance value is obtained by synthesizing the three displacements and the initial clearance. Further improving the calculation accuracy of the tip clearance, compared with the traditional unidirectional coupling tip clearance calculation method, the calculation error of the tip clearance under the test run condition is less than 10%, and the calculation accuracy is improved by more than 30%, which has important engineering reference value for the prediction of the tip clearance of aeroengines. Description of the Drawings

[0051] Figure 1 is a schematic diagram of the main components and cross-sections of an aeroengine;

[0052] Figure 2 is a schematic diagram of the two-dimensional solid domain calculation model of the aeroengine core;

[0053] Figure 3 is a schematic diagram of the main flow blade calculation model of the core;

[0054] Figure 4 is the one-dimensional fluid network calculation model of the core air system;

[0055] Figure 5 is a schematic diagram of the finite element mesh of the core solid domain;

[0056] Figure 6 is the throttle and speed change curve;

[0057] Figure 7It is a schematic diagram of data interaction and time-step coupling between physical fields;

[0058] Figure 8 It is a flow chart of the cross-scale multi-physical field coupling algorithm;

[0059] Figure 9 It is a comparison chart of the calculated value and the experimental value of the tip clearance. Specific implementation manner

[0060] The following details the implementation manner of the present invention. The examples of the implementation manner are shown in the accompanying drawings. The implementation manner described below with reference to the accompanying drawings is exemplary and is only used to explain the present invention, and cannot be construed as a limitation to the present invention.

[0061] The present invention proposes a method for predicting the tip clearance of cross-scale multi-physical field coupling, which is aimed at the problem of predicting the transient tip clearance of aero-engines. First, a one-dimensional network model of the engine air system, a two-dimensional axisymmetric model of the disk and the casing, a parametric model of the blade, and a quasi-one-dimensional model of the main flow path are established; then, after calculating the fluid pressure, temperature, and flow rate by using the main flow path cavity method and the one-dimensional transient fluid network method of the air system, the heat transfer boundary between the engine disk and the casing is converted, considering the thermal load and centrifugal load in its unsteady heat conduction process, the radial displacement of the disk and the casing is calculated by using the two-dimensional axisymmetric finite element method, and the temperature and radial displacement of the blade are calculated by using engineering algorithms; finally, the radial displacements of the disk, the blade, and the casing, as well as the initial tip clearance value, are comprehensively considered to obtain the tip clearance value of each stage in the transient process.

[0062] During the calculation process, the heat transfer boundary condition is transferred from the main flow path and the air system flow field to the solid domain temperature field, the qualitative temperature of the thermal load and material properties is transferred from the solid domain temperature field to the solid domain displacement field, and the geometric parameters of the throttling elements in the air system flow field are updated by the solid domain displacement field; conversely, the change of the grid node coordinates is provided from the solid domain displacement field to the solid domain temperature field, the solid domain wall temperature is provided from the solid domain temperature field to the air system flow field, the change of the labyrinth clearance is provided from the solid domain displacement field to the air system flow field, and the change of the tip clearance is provided from the solid domain displacement field superimposed with the blade deformation to the main flow path, which is used as the basis for correcting the main flow path parameters. Thus, a two-way transfer of information between different physical fields is formed.

[0063] Figure 1 It is a schematic diagram of the main components and cross-section division of the aero-engine. The following takes Figure 2 the core engine of the shown aero-engine as an example to detail the specific implementation process of the present invention for calculating the tip clearance.

[0064] First, establish the solid domain model of the core engine. For the casing and the disk of the core engine, considering the symmetry of its structure, a two-dimensional axisymmetric model of the casing and the disk is established, as shown in Figure 2As shown. For the mainstream blades, due to their axial discrete arrangement, a parametric model is adopted. For example, Figure 3 as shown, it is divided into several segments along the radial direction according to the blade height, and parameters such as the perimeter, area, and chord length on each interface at the corresponding height are extracted. Secondly, a calculation model of the core engine air system is established. For example, Figure 4 as shown, the air system is simplified into a fluid network composed of components and nodes, and is divided into four independent flow paths according to the structural characteristics of each network.

[0065] Since the temperature and displacement of the blades in the solid domain adopt engineering algorithms, and the air system in the fluid domain is linearly solved by the one-dimensional fluid network method, and the main flow path adopts the quasi-one-dimensional cavity method, the above models do not need to be meshed. However, the casing and disk are calculated for temperature and displacement by the finite element method, and it is necessary to mesh them with finite elements. The casing and disk are meshed by ANSYS MECHANICAL APDL, and the mesh is refined at the labyrinth teeth. The final number of mesh elements obtained is 33988, and the overall meshing effect is as shown in Figure 5 as shown. In addition, during the calculation process, the material types of each component are considered, and the mesh elements are divided into several material zones, and the overall number of material types exceeds 14.

[0066] To better illustrate the implementation process of the present invention, the tip clearance of the core engine in the engine process is calculated. Figure 6 is the curve of the throttle lever angle (PLA) change during the test run. During the entire process, the engine maintains the inlet and ambient temperature T1 = 300K and the total pressure P1 = 101325Pa. The entire test run process includes the 0 - t1 start process; the throttle angle is maintained at 12° for ground idle speed from t1 to t2, and the high-pressure target speed is 72% of the relative physical speed; the warm-up process from t2 to t3, the throttle angle is maintained at 35°, and the high-pressure target speed is 90% of the relative physical speed; the throttle is closed to 12° and maintained at ground idle speed from t3 to t4, and the throttle angle is maintained at 45° from t4 to t5, and the high-pressure target speed is 95% of the relative physical speed. Among them, the high-pressure relative physical speed N 2R is defined as the ratio of the current speed to the maximum physical speed.

[0067] During the clearance calculation process, the mainstream blade calculation is first carried out. Select a segment as shown in Figure 3 as the research object, and the following control equation is established according to the law of conservation of energy for temperature calculation:

[0068] Q = hA1ΔT1 + λA2ΔT2 - λA3ΔT3

[0069] QΔt s = c p mΔT

[0070] The first formula represents three parts: the heat transferred from the air flow to the blade through convective heat transfer, the heat input from one end of the segment through heat conduction, and the heat output from the other end of the segment. The average convective heat transfer coefficient of the blade adopts the empirical relationship:

[0071] Nu = 0.032Re 0.8

[0072] The heat transfer coefficient of the blade takes the chord length as the characteristic dimension and is calculated according to the flat plate with non-separated flow around the air flow. Considering the influence of rotation, the following correction coefficient is established:

[0073]

[0074] where u is the linear velocity of the rotation of the segment; l is the length of the segment; v is the relative velocity of the air flow inlet; D is the average diameter of the segment.

[0075] When calculating the blade displacement, the total thermal displacement ΔR of the blade is the sum of the displacements of each segment, as shown in the following formula:

[0076]

[0077] When calculating the centrifugal displacement, the displacement Δl caused by the centrifugal load of each segment j is as follows:

[0078] ΔP i = ρω 2 A i r i Δr i

[0079]

[0080] Δl j = l·(σ j / E)

[0081] where ΔP i is the centrifugal force received by this segment; σ j is the stress of the low-radius end face of this segment; E is the elastic modulus of the blade. Then the total centrifugal displacement ΔL of the blade is as follows:

[0082]

[0083] The total radial displacement of the blade is obtained by superimposing and summing the total thermal displacement and the total centrifugal displacement of the blade.

[0084] After obtaining the blade displacement, the displacement of the casing and disk is calculated. At this time, the problem of cross-time-scale information transfer between different physical fields needs to be considered, so the weak fluid-structure coupling method as Figure 7 shown is adopted, that is, in a solid domain time step Δt sWithin, multiple fluid domain calculations are performed. During this period, it is assumed that the temperature of the solid wall remains unchanged. After the number of fluid domain calculation steps reaches the set value, the solid domain parameter calculation is carried out, thus completing the calculation of a coupled iteration step. For the main flow and the secondary air system, both belong to the fluid domain and use the same time step Δt f Perform the calculation.

[0085] Through the above weak coupling method, two-way information interaction between the fluid and the solid is carried out to calculate the displacement of the disk and the casing. Refer to Figure 8 The specific process is as follows:

[0086] 1) The main flow path parameter module reads in relevant preset information such as the engine component characteristics parameters and corresponding control laws. The secondary air system module reads in preset information such as the size parameters and fluid parameters of relevant components. The solid domain module uses the Ansys APDL finite element analysis module to process the geometric information and divide the mesh, and outputs the relevant dimension coordinate information to the finite element calculation module.

[0087] 2) Start the calculation from the main flow path module, and calculate the parameters of the main flow path at time Δt f The total temperature T, total pressure p, and flow rate g of the fluid at each key cross-section, and output relevant information such as the high-pressure physical rotation speed n2. Δt f Is the time step of the fluid domain.

[0088] 3) The one-dimensional transient fluid network calculation model of the secondary air system reads the aerodynamic boundary at the moment and starts the calculation, and outputs the flow rate g, pressure p, and fluid temperature T of the internal components of the secondary air system at the moment.

[0089] 4) Since the response time of the fluid domain is relatively fast, in the typical transition process, the iterative calculation between the main flow and the secondary flow is carried out first, and the loop steps 2 - step 3 are repeated, and the number of iterations is n.

[0090] 5) When the number of iterations n satisfies Δt s = n·Δt f After that (where Δt s Is the time step of the solid domain), suspend the fluid calculation. Based on the calculation results of the main and secondary fluids at the current moment, convert the third-kind heat transfer boundary condition, and perform the transient thermal analysis of the engine components, and output the component temperature field.

[0091] 6) Determine the centrifugal load of the components based on the rotation speed n2 calculated by the main flow path parameter module, and use the component temperature field calculated in step 5) as the temperature load to solve the deformation amount of the components, and output the newly calculated deformation amount and new coordinates.

[0092] 7) Update the node parameter information of the geometric model according to the deformation amount of the key parts, recalculate the aerodynamic parameters, and perform the next round of fluid calculation with the number of times n.

[0093] 8) By repeating steps 1) to 7), the rotor-stator displacements at each moment can be obtained, and then the tip clearance can be calculated. In the above process, the transient temperature calculation is carried out by establishing a system of linear equations in the following matrix form:

[0094] [K] τ [t] τ =[F] τ

[0095] where [K] τ is the overall temperature stiffness matrix, [t] τ is the nodal temperature matrix, and [F] τ is the load matrix.

[0096] The two-dimensional axisymmetric finite element method is used to solve the disk-casing displacement. The thermal loads and centrifugal loads of each element are summed separately and superimposed to obtain the overall load matrix, and combined with the overall stiffness matrix to form the following system of linear equations for displacement solution:

[0097] [K][δ]=[L]

[0098] where [K] is the overall displacement stiffness matrix, [δ] is the displacement matrix of each grid node, and [L] is the overall load matrix of the grid nodes.

[0099] After calculating the displacements of the casing, disk, and blades, the variation law of the tip clearance at each stage in the transient condition can be obtained by combining the initial tip clearance value.

[0100] Taking Figure 2 the core engine model as an example, the tip clearance calculation and the engine test results are compared and analyzed.

[0101] As can be seen from Figure 9 , the method of the present invention can better capture the variation trend of the high-pressure turbine clearance with the test run time, and the numerical value is also close to the test results. The quantitative comparison results of the tip clearance variation and the clearance prediction error at each moment are shown in Table 1. It can be seen that the errors of the tip clearance variation at the end moment and during the stable rotation speed of each rotation speed calculated by the method of the present invention are all within 10%.

[0102] Table 1 Comparison of test and calculation results

[0103]

[0104] In terms of details, during the test process, the experimental results of the high-pressure turbine clearance were affected by vibrations, showing a fluctuating pattern. At the same time, since manual throttle control cannot achieve transient steps like a computational program, there are slight deviations between the acceleration and deceleration processes in actual engine tests and those calculated in this invention, resulting in a difference of several seconds in the clearance response time between the two. In addition, before t1 is the startup process, and the internal flow and heat transfer state of the engine are quite different from those at idle speed and above, so the error in model simulation is relatively large, which further leads to a relatively large error in clearance prediction in the first few minutes after idle speed starts. However, compared with the stable operating stage at idle speed and above, the startup process has a lower temperature, a smaller heat transfer coefficient, and a shorter duration. Therefore, after several minutes of continuous ground idle speed, the impact caused by the simulation deviation in the startup state gradually weakens to an ignorable level.

[0105] Based on the same inventive concept, an embodiment of the present application provides a computer device, including a memory, a processor, and a computer program stored in the memory and executable on the processor. When the processor executes the computer program, it implements the steps of the aforementioned cross-scale multi-physical field coupling tip clearance prediction method.

[0106] Based on the same inventive concept, an embodiment of the present application provides a computer-readable storage medium storing a computer program, which, when executed by a processor, implements the steps of the aforementioned cross-scale multi-physical field coupling tip clearance prediction method.

[0107] Those skilled in the art should understand that the embodiments of the present invention can be provided as a method, a system, or a computer program product. Therefore, the present invention can take the form of a completely hardware embodiment, a completely software embodiment, or an embodiment combining software and hardware aspects. Moreover, the present invention can take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk memories, CD-ROMs, optical memories, etc.) containing computer-usable program code.

[0108] The present invention is described with reference to the flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to embodiments of the present invention. It should be understood that each flow and / or block in the flowcharts and / or block diagrams, as well as the combination of flows and / or blocks in the flowcharts and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to the processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing devices to generate a machine, such that the instructions executed by the processor of the computer or other programmable data processing devices generate means for implementing the specified functions in Figure 1 one process or multiple processes and / or blocks Figure 1 one block or multiple blocks.

[0109] These computer program instructions can also be stored in a computer-readable memory that can direct a computer or other programmable data processing apparatus to work in a particular manner, such that the instructions stored in the computer-readable memory produce a manufacture including an instruction device that implements the functions specified in one or more processes and / or blocks Figure 1 in the flow Figure 1 or boxes.

[0110] These computer program instructions can also be loaded onto a computer or other programmable data processing apparatus, such that a series of operational steps are performed on the computer or other programmable apparatus to produce a computer-implemented process, whereby the instructions executed on the computer or other programmable apparatus provide steps for implementing the functions specified in one or more processes and / or blocks Figure 1 in the flow Figure 1 or boxes.

[0111] The above embodiments are only for illustrating the technical idea of the present invention, and the protection scope of the present invention cannot be limited thereby. Any modification made on the basis of the technical solution according to the technical idea proposed by the present invention shall fall within the protection scope of the present invention.

Claims

1. A method for predicting tip clearance with cross-scale multi-physics field coupling, characterized in that The method is specifically as follows: Establish a one-dimensional network model of the engine air system, a two-dimensional axisymmetric model of the impeller and casing, a parametric model of the blades, and a quasi-one-dimensional model of the main flow channel; The quasi-one-dimensional model of the main channel is solved using a volumetric method, and the one-dimensional network model of the air system is solved using a one-dimensional transient fluid network method to obtain the heat exchange boundary between the engine wheel and the casing; Considering the thermal load and centrifugal load during the unsteady heat conduction process in the solid domain of the engine, the two-dimensional axisymmetric finite element method is used to solve the two-dimensional axisymmetric model of the wheel and casing to obtain the radial displacement of the wheel and casing; Solving the parameterized model of the blade using an engineering algorithm of temperature and displacement to obtain the temperature and radial displacement of the blade; The tip clearance value in the transient process is obtained by combining the radial displacement of the wheel, casing and blades and the initial tip clearance value; The radial displacement of the wheel and the casing is solved by using the fluid-solid weak coupling method to perform two-way information interaction between the fluid domain and the solid domain. During the calculation, multiple fluid domain calculations are performed, during which the solid wall temperature is considered to be constant. , is the fluid domain time step, is the number of iterations, specifically: 1) During the current fluid-domain time step of the current solid-domain time step, the tip clearance calculated in the previous solid-domain time step is fed back to the main flow channel. The total temperature, total pressure, and flow rate of the fluid at each section of the engine main flow channel, the engine speed, and the pressure and temperature of the air system inlet and outlet are calculated using the volumetric method. 2) Within the same fluid domain time step as in 1), the grate gap calculated in the previous solid domain time step is fed back to the air system. Based on the pressure and temperature of the air system inlet and outlet calculated in 1), the flow rate, pressure, and temperature of the fluid in the air system components are calculated using a one-dimensional transient fluid network method. 3) Repeat 1) and 2) in the next fluid domain time step of the current solid domain time step until all fluid domain time steps within the current solid domain time step are iteratively calculated and proceed to 4); 4) Based on the results of the last fluid domain time step of the current solid domain time step, the third type of heat transfer boundary conditions are converted to perform transient thermal analysis of the engine solid domain to obtain the solid domain temperature field; 5) Determine the centrifugal load in the solid domain based on the rotational speed calculated in the last fluid domain time step of the current solid domain time step. Use the solid domain temperature field calculated in 4) as the temperature load to solve for the solid domain deformation, including the deformation of the output wheel and casing, as well as the grate teeth. 6) The radial displacement of the wheel disc and casing is obtained based on the deformation of the wheel disc and casing. The tip clearance is obtained based on the radial displacement of the wheel disc and casing, the radial displacement of the blades, and the initial tip clearance value; the grate clearance is obtained based on the grate deformation; 7) Update the 2D axisymmetric model using the deformations of the wheel and casing, and return to step 1) to enter the calculation of the next solid domain time step.

2. The tip clearance estimation method based on cross-scale multi-physics field coupling according to claim 1, characterized in that: The specific solution process for the radial displacement of the blade is as follows: (1) Divide the blade into several sections along the radial direction according to the blade height, extract the geometric information of each section, and establish a parametric model of the blade; (2) The energy conservation equation for each segment is established as follows: ， ， Among them, is the net heat flux of the blade segment, is the heat transfer coefficient between the gas flow and the wall surface, is the gas-solid convective heat transfer area, is the temperature difference of heat transfer between the gas and the solid, is the thermal conductivity of the blade segment, is the end face area for the heat introduced into the blade segment, is the average temperature difference between the two segments above and below the end face where the heat is introduced, is the end face area for the heat exported from the blade segment, is the average temperature difference between the two segments above and below the end face where the heat is exported, is the time step of the solid domain, is the specific heat capacity of the blade segment, is the mass of the blade segment, is the change in the average temperature of the blade segment within a single time step of the solid domain; (3) The blade segment thermal load is calculated from (2), and the thermal displacement of each segment is superimposed as follows: ， in, is the overall thermal displacement of the blade, is the blade segment linear expansion coefficient, For the Segment temperature variation, For the Segment length, is the number of segments; The displacements caused by the centrifugal loads of each blade segment are superimposed as follows: ， Among them, is the overall centrifugal displacement of the blade, is the centrifugal stress borne by the th segment, and is the elastic modulus of the th segment; The total radial displacement of the blade is obtained by superimposing the overall thermal displacement of the blade and the overall centrifugal displacement of the blade.

3. The tip clearance prediction method for cross-scale multi-physical field coupling according to claim 1, characterized in that In the above 2), the volume effect of the cavity is considered in the solution process of the one-dimensional transient fluid network method, and the relationship between the inlet and outlet flow and pressure is: ， ， in, is the cavity inlet flow rate, is the cavity outlet flow rate, is the cross-sectional area of the cavity inlet, is the cross-sectional area of the cavity outlet, is the total pressure of the inlet airflow, is the ideal gas constant, is the total temperature of the inlet airflow, is the adiabatic index, is the static pressure in the middle of the cavity, is the flow correction coefficient, is the static pressure outside the cavity, is the average static temperature of the airflow in the cavity, is the static pressure of the airflow at the cavity outlet.

4. The tip clearance prediction method for cross-scale multi-physical field coupling according to claim 1, wherein In the above 4), when analyzing the transient thermal performance of the engine solid domain, the two-dimensional axisymmetric finite element method is used to calculate the unsteady heat conduction. There is no internal heat source in the calculation model, and the thermal conductivity of the solid domain is isotropic. Therefore, the unsteady heat conduction differential equation in the cylindrical coordinate system is as follows: ， in, is the solid domain density, is the specific heat capacity of the solid domain at constant pressure, is the solid domain temperature, For time, is the thermal conductivity of the solid domain, is the radial coordinate, is the axial coordinate.

5. The tip clearance prediction method for cross-scale multi-physical field coupling according to claim 1, wherein In the above 1), when the volume method is used for calculation, the change rate of the inlet and outlet parameters of the combustion chamber in the main channel is: ， ， in, is the combustion chamber outlet temperature, is the combustion chamber outlet pressure, is the ideal gas constant, is the volume of the main combustion chamber, is the constant volume heat capacity, is the oil supply amount, is the calorific value of fuel oil, The main combustion efficiency, is the enthalpy of the fuel when it enters the combustion chamber. is the unit enthalpy at the combustion chamber outlet, is the heat capacity ratio, is the combustion chamber inlet flow rate, is the unit enthalpy at the combustion chamber inlet, is the combustion chamber outlet flow rate.

6. A computer device, comprising a memory, a processor, and a computer program stored in the memory and capable of running on the processor, characterized in that, When the processor executes the computer program, the steps of the cross-scale multi-physics field coupled tip clearance estimation method according to any one of claims 1 to 5 are implemented.

7. A computer-readable storage medium storing a computer program, characterized in that: When the computer program is executed by a processor, the steps of the cross-scale multi-physics field coupled blade tip clearance estimation method according to any one of claims 1 to 5 are implemented.

Images