Calculation method for demand quantity of turbine cooling air
By constructing a heat transfer geometric model and using an explicit differential calculation method, the problems of high time consumption and low reliability in calculating the cooling gas demand of supercritical carbon dioxide radial turbines were solved, achieving efficient and accurate calculation of cooling gas demand and ensuring the safety and lifespan of the turbine.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- INST OF ENGINEERING THERMOPHYSICS - CHINESE ACAD OF SCI
- Filing Date
- 2026-01-22
- Publication Date
- 2026-05-08
AI Technical Summary
Existing technologies for calculating the cooling gas demand of supercritical carbon dioxide radial turbines suffer from high computational cost, long processing time, and poor numerical convergence in three-dimensional simulation methods, resulting in low reliability of the results and failing to meet the requirements of actual operating conditions.
The method for calculating turbine cooling gas demand using electronic devices constructs a heat transfer geometric model, calculates steady-state and transient temperature fields based on initial conditions and multiple algorithms, determines the cooling gas demand, and uses an explicit differential calculation method for iterative calculation, simplifying it into a simplified model of solid and fluid regions, avoiding direct handling of supercritical region property changes.
It improves the accuracy and consistency of calculations, reduces computing resources and time, and ensures the reliability and safety of cooling gas requirements after turbine shutdown.
Smart Images

Figure CN121997585A_ABST
Abstract
Description
Technical Field
[0001] This disclosure relates to the field of turbine cooling technology, specifically to a method for calculating the turbine cooling gas demand. Background Technology
[0002] In supercritical carbon dioxide radial turbines, to meet the cooling requirements of the turbine shaft end, dry sealing gas needs to be introduced as cooling gas into the rotor-stator gap of the turbine to cool the rotor. For calculating the cooling gas requirement after the turbine stops, a three-dimensional simulation method can be used to transiently simulate the temperature field changes of the rotor and casing over time to determine whether cooling gas injection is necessary.
[0003] Three-dimensional numerical simulation methods have high computational requirements, long computation time, and poor numerical convergence. Especially when dealing with the physical properties of supercritical carbon dioxide working fluid, local non-physical solutions or solution failures are prone to occur, resulting in low reliability of turbine cooling gas demand results, excessively long calculation cycles, and inability to meet actual operating conditions. Summary of the Invention
[0004] In view of the above problems, this disclosure provides a method for calculating the turbine cooling air demand performed by an electronic device.
[0005] According to a first aspect of this disclosure, a method for calculating the cooling air demand of a turbine, executed by an electronic device, is provided, comprising: constructing a heat transfer geometric model based on the structural parameters of the turbine, wherein the heat transfer geometric model includes a solid region for simulating temperature changes of solid components of the turbine, and a fluid region for simulating temperature changes of cooling airflow through a cooling channel; constructing a first initial condition based on the basic rotational speed of the turbine rotor, the initial temperature of the solid components of the turbine, and the basic parameters of the cooling air, and calculating the steady-state temperature fields of the solid and fluid regions under the cooling air supply condition during turbine operation using a first algorithm based on the first initial condition; constructing a second initial condition based on the adjusted rotational speed of the turbine rotor and the temperature of the corresponding solid components of the turbine in the steady-state temperature field, and calculating the basic transient temperature field of the solid and fluid regions under the cooling air supply condition after the turbine stops using a second algorithm based on the second initial condition, wherein the value of the adjusted rotational speed is zero; and outputting a calculation result indicating that cooling air supply is not required after the turbine stops, provided that the temperature of the target region in the basic transient temperature field does not exceed a predetermined threshold.
[0006] According to embodiments of this disclosure, the method for calculating the turbine cooling gas demand executed by an electronic device further includes: when the temperature of the target region in the basic transient temperature field exceeds a predetermined threshold, constructing multiple sets of third initial conditions based on the turbine rotor's adjusted rotational speed, the temperature of the corresponding turbine solid component in the steady-state temperature field, and multiple sets of cooling gas adjustment parameters; and calculating multiple sets of reference transient temperature fields of the solid region and fluid region under the condition of continued cooling gas supply after turbine shutdown using a third algorithm based on the multiple sets of third initial conditions, wherein the multiple sets of reference transient temperature fields correspond to multiple sets of cooling gas adjustment parameters; determining the cooling gas adjustment parameter with the minimum flow rate that minimizes the temperature of the target region from exceeding the predetermined threshold based on the calculation results of the multiple sets of reference transient temperature fields, and outputting the calculation result of the cooling gas demand after turbine shutdown based on the minimum cooling gas adjustment parameter.
[0007] According to embodiments of this disclosure, constructing a heat transfer geometric model based on the structural parameters of the turbine includes: constructing a heat transfer geometric model based on the structural parameters of the turbine's rotor, casing, and wheel; the heat transfer geometric model includes multiple solid regions used to simulate the temperature changes of the rotor, casing, and wheel wall, respectively, and multiple fluid regions used to simulate the temperature changes of the cooling airflow through the rotor gap and wheel back gap, respectively, wherein the rotor gap is the gap between the bottom of the casing and the rotor, and the wheel back gap is the gap between the side of the casing and the wheel.
[0008] According to embodiments of this disclosure, constructing the first initial condition includes: constructing the first initial condition based on the basic rotational speed of the turbine rotor, the basic disk wall temperature, the initial temperature of the rotor and casing, the first cooling gas inlet temperature, the first cooling gas inlet pressure, and the first cooling gas mass flow rate, wherein the basic disk wall temperature is equal to the temperature of the working fluid inside the turbine; constructing the second initial condition includes: constructing the second initial condition based on the adjusted rotational speed of the turbine rotor and the temperatures corresponding to the disk wall, rotor, and casing in the steady-state temperature field; constructing multiple sets of third initial conditions includes: constructing multiple sets of third initial conditions corresponding to the multiple second cooling gas mass flow rates based on the adjusted rotational speed of the turbine rotor, the temperatures corresponding to the disk wall, rotor, and casing in the steady-state temperature field, the second cooling gas inlet temperature, the second cooling gas inlet pressure, and multiple increasing second cooling gas mass flow rates.
[0009] According to the embodiments of this disclosure, the calculation of the steady-state temperature field, the basic transient temperature field, and the reference transient temperature field of the solid region and the fluid region all include multiple iterative calculations. Specifically, for any of the first algorithm, the second algorithm, and the third algorithm, the calculation at any k+1 time includes: calculating the temperature calculation results of the solid region and the fluid region at the k+1 time based on the temperature calculation results at the k-th time of the solid region and the fluid region, where k is a positive integer.
[0010] According to embodiments of this disclosure, the method further includes meshing the solid region and the fluid region; based on the temperature calculation result at time k, calculating the temperature calculation result at time k+1 of the solid region and the fluid region includes: at time k+1, calculating the temperature of any target mesh with two-dimensional coordinates i and j in the solid region and the fluid region. The following methods are used:
[0011] ;
[0012] Where A, B, C, D, E, and F are calculation coefficients. , , , , These refer to the temperature of the target grid with two-dimensional coordinates i and j at time k, and the temperature of multiple adjacent grids surrounding the target grid.
[0013] According to embodiments of this disclosure, the fluid region includes an axial flow region of cooling gas corresponding to the rotor clearance and a radial flow region of cooling gas corresponding to the wheel back clearance; the solid region includes a rotor solid region and a casing solid region; wherein, the rotor solid region includes a rotor center region and multiple rotor boundary surfaces, and the casing solid region includes a casing center region and multiple casing boundary surfaces; in any of the first algorithm, the second algorithm, and the third algorithm, the calculation methods for the coefficients A, B, C, D, E, and F for the axial flow region of cooling gas, the radial flow region of cooling gas, the rotor center region, the multiple rotor boundary surfaces, the casing center region, and the multiple casing boundary surfaces are different.
[0014] According to embodiments of this disclosure, in the first algorithm:
[0015] For the axial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are obtained based on the following parameters: the convective heat transfer coefficient of the first rotor surface, the axial heat transfer coefficient of the first casing, the mass flow rate of the cooling gas, the isobaric heat capacity of the cooling gas, the density of the cooling gas, the rotor radius, and the rotor clearance dimension; wherein, the convective heat transfer coefficient of the first rotor surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the rotor surface during turbine operation, and the axial heat transfer coefficient of the first casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the bottom of the casing during turbine operation;
[0016] For the radial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are obtained based on the following parameters: radial heat transfer coefficient of the first casing, convective heat transfer coefficient of the first wheel surface, mass flow rate of the cooling gas, isobaric heat capacity of the cooling gas, density of the cooling gas, and wheel back clearance dimension; the radial heat transfer coefficient of the first casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the casing side during turbine operation, and the convective heat transfer coefficient of the first wheel surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the wheel surface during turbine operation;
[0017] For the target rotor boundary surface among multiple rotor boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the first rotor surface, the rotor thermal diffusivity, and the rotor thermal conductivity; where the target rotor boundary surface is the area of the rotor surface in direct contact with the cooling gas.
[0018] For the rotor center region and the other rotor boundary surfaces except the target rotor boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the rotor thermal diffusivity.
[0019] For two first target casing boundary surfaces among multiple casing boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: first casing radial heat transfer coefficient, first casing axial heat transfer coefficient, casing thermal diffusivity, and casing thermal conductivity; wherein, the two target casing boundary surfaces are the areas of the casing side and the casing bottom that are in direct contact with the cooling gas, respectively.
[0020] For the second target casing boundary surface, which is in direct contact with the working fluid inside the turbine among multiple casing boundary surfaces, A=D=1, B=C=E=F=0; the temperature of the target casing boundary surface is equal to the temperature of the working fluid inside the turbine.
[0021] For the central region of the casing and the remaining casing boundary surfaces other than the first target casing boundary surface and the second target casing boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the casing thermal diffusivity. The temperature at the boundary surface where the rotor and the wheel contact, and the temperature at the boundary surface where the cooling gas and the wheel contact, are equal to the temperature of the working fluid inside the turbine.
[0022] According to embodiments of this disclosure, the convective heat transfer coefficient of the first rotor surface and the axial heat transfer coefficient of the first casing are calculated based on the first Nusselt number, the thermal conductivity of the cooling gas, and the rotor length; the first Nusselt number is calculated based on the rotating Reynolds number at the rotor gap, the axial Reynolds number at the rotor gap, the rotor length, the rotor radius, the rotor gap size, the dynamic viscosity of the cooling gas, the isobaric heat capacity of the cooling gas, and the thermal conductivity of the cooling gas; the radial heat transfer coefficient of the first casing and the convective heat transfer coefficient of the first disc surface are calculated based on the second Nusselt number, the thermal conductivity of the cooling gas, and the disc radius; the second Nusselt number is calculated based on the rotating Reynolds number at the turbine disc, the radial Reynolds number at the turbine disc, the back clearance size, and the disc radius.
[0023] According to embodiments of this disclosure, in the third algorithm:
[0024] For the axial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are obtained based on the following parameters: the convective heat transfer coefficient of the second rotor surface, the axial heat transfer coefficient of the second casing, the mass flow rate of the cooling gas, the isobaric heat capacity of the cooling gas, the density of the cooling gas, the rotor radius, and the rotor clearance dimension; wherein, the convective heat transfer coefficient of the second rotor surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the rotor surface after the turbine is shut down, and the axial heat transfer coefficient of the second casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the bottom of the casing after the turbine is shut down;
[0025] For the radial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are obtained based on the following parameters: radial heat transfer coefficient of the second casing, convective heat transfer coefficient of the second wheel surface, mass flow rate of the cooling gas, isobaric heat capacity of the cooling gas, density of the cooling gas, and wheel back clearance dimension; the radial heat transfer coefficient of the second casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the side of the casing after the turbine stops, and the convective heat transfer coefficient of the second wheel surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the wheel surface after the turbine stops;
[0026] For the target rotor boundary surface among multiple rotor boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the second rotor surface, the rotor thermal diffusivity, and the rotor thermal conductivity; where the target rotor boundary surface is the area of the rotor surface in direct contact with the cooling gas.
[0027] For the rotor center region and the other rotor boundary surfaces except the target rotor boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the rotor thermal diffusivity.
[0028] For two first target casing boundary surfaces among multiple casing boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: second casing radial heat transfer coefficient, second casing axial heat transfer coefficient, casing thermal diffusivity, and casing thermal conductivity; wherein, the two first target casing boundary surfaces are the areas of the casing side and the casing bottom that are in direct contact with the cooling gas, respectively.
[0029] For the second target casing boundary surface, which is in direct contact with the working fluid inside the turbine among multiple casing boundary surfaces, it is set as an adiabatic boundary. The calculation coefficients A, B, C, D, E, and F are obtained based on the casing thermal diffusivity.
[0030] For the central region of the casing and the remaining casing boundary surfaces other than the first target casing boundary surface and the second target casing boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the casing thermal diffusivity.
[0031] The boundary surfaces where the rotor and the disk contact each other, as well as the boundary surfaces where the cooling gas and the disk contact each other, are set as adiabatic boundaries. The calculation coefficients A, B, C, D, E, and F are obtained based on the rotor thermal diffusivity.
[0032] According to embodiments of this disclosure, the convective heat transfer coefficient of the second rotor surface and the axial heat transfer coefficient of the second casing are calculated based on the third Nusselt number, the thermal conductivity of the cooling gas, and the rotor length; the third Nusselt number is calculated based on the axial Reynolds number at the rotor gap, the dynamic viscosity of the cooling gas, the isobaric heat capacity of the cooling gas, and the thermal conductivity of the cooling gas; the radial heat transfer coefficient of the second casing and the convective heat transfer coefficient of the second wheel surface are calculated based on the fourth Nusselt number, the thermal conductivity of the cooling gas, and the wheel back gap size; the fourth Nusselt number is calculated based on the radial Reynolds number at the turbine wheel, the dynamic viscosity of the cooling gas, the isobaric heat capacity of the cooling gas, and the thermal conductivity of the cooling gas.
[0033] According to an embodiment of this disclosure, in the second algorithm: for the axial flow region and the radial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: thermal conductivity of the cooling gas, constant pressure heat capacity of the cooling gas, and density of the cooling gas.
[0034] For the rotor center region, the target rotor boundary surface among multiple rotor boundary surfaces, and the other rotor boundary surfaces other than the target rotor boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the rotor thermal diffusivity.
[0035] For the central region of the casing, two first target casing boundary surfaces among multiple casing boundary surfaces, the second target casing boundary surface on the side of the casing that is in direct contact with the working fluid inside the turbine, and the remaining casing boundary surfaces other than the first and second target casing boundary surfaces, the calculation coefficients A, B, C, D, E, and F are obtained based on the casing thermal diffusivity.
[0036] The boundary surfaces where the rotor and the disk contact each other, as well as the boundary surfaces where the cooling gas and the disk contact each other, are set as adiabatic boundaries. The calculation coefficients A, B, C, D, E, and F are calculated based on the rotor thermal diffusivity.
[0037] According to embodiments of this disclosure, the target area is the solid wall surface at the inlet of the cooling channel. Attached Figure Description
[0038] The foregoing contents, as well as other objects, features, and advantages of this disclosure, will become clearer from the following description of embodiments with reference to the accompanying drawings, in which:
[0039] Figure 1 A schematic diagram of a turbine shaft end cooling structure according to an embodiment of the present disclosure is shown.
[0040] Figure 2 A flowchart illustrating a method for calculating turbine cooling gas demand according to an embodiment of the present disclosure is shown schematically.
[0041] Figure 3 A schematic diagram of a heat transfer geometry model constructed from the various components of a turbine according to an embodiment of the present disclosure is shown.
[0042] Figure 4 A schematic diagram of a heat transfer geometry model mesh according to an embodiment of the present disclosure is shown. Detailed Implementation
[0043] The embodiments of the present disclosure will now be described with reference to the accompanying drawings. However, it should be understood that these descriptions are exemplary only and are not intended to limit the scope of the disclosure. In the following detailed description, numerous specific details are set forth to provide a thorough understanding of the embodiments of the present disclosure for ease of explanation. However, it will be apparent that one or more embodiments may be practiced without these specific details. Furthermore, descriptions of well-known structures and techniques are omitted in the following description to avoid unnecessarily obscuring the concepts of the present disclosure.
[0044] The terminology used herein is for the purpose of describing particular embodiments only and is not intended to limit this disclosure. The terms “comprising,” “including,” etc., as used herein indicate the presence of the stated features, steps, operations, and / or components, but do not exclude the presence or addition of one or more other features, steps, operations, or components.
[0045] All terms used herein (including technical and scientific terms) have the meanings commonly understood by those skilled in the art, unless otherwise defined. It should be noted that the terms used herein are to be interpreted in a manner consistent with the context of this specification, and not in an idealized or overly rigid way.
[0046] When using expressions such as "at least one of A, B and C", they should generally be interpreted in accordance with the meaning that is commonly understood by those skilled in the art (e.g., "a system having at least one of A, B and C" should include, but is not limited to, a system having A alone, a system having B alone, a system having C alone, a system having A and B, a system having A and C, a system having B and C, and / or a system having A, B and C, etc.).
[0047] During the shutdown of a supercritical carbon dioxide radial turbine, the rotor, due to its large heat capacity, dissipates heat slowly, and residual heat continues to transfer to the shaft end area. Without effective cooling control, components such as shaft end seals may remain in a high-temperature environment for extended periods, affecting their functional integrity and service life. Therefore, an accurate assessment of the rotor temperature field decay process is necessary to determine the appropriate cooling timing and gas flow rate.
[0048] Three-dimensional transient coupled simulation is typically used to perform detailed modeling and solution of the entire turbine domain, including the fluid-structure interface. This method needs to handle the drastic changes in physical properties such as density and specific heat capacity of supercritical carbon dioxide near the critical point, which places extremely high demands on the stability of the numerical solution. To capture property gradients and heat transfer details, extremely fine meshes and extremely small time steps must be used, leading to a sharp increase in computational load and excessively long simulation times. At the same time, abrupt changes in physical properties can easily cause difficulties in solving the equations, often resulting in computational divergence or non-physical solutions, making it difficult to guarantee the reliability of the simulation results.
[0049] In view of this, embodiments of the present disclosure provide a method for calculating the cooling air demand of a turbine executed by an electronic device, comprising: constructing a heat transfer geometric model based on the structural parameters of the turbine, wherein the heat transfer geometric model includes a solid region for simulating temperature changes of the turbine solid components and a fluid region for simulating temperature changes of the cooling airflow through the cooling channel; constructing a first initial condition based on the basic rotational speed of the turbine rotor, the initial temperature of the turbine solid components, and the basic parameters of the cooling air, and calculating the steady-state temperature fields of the solid and fluid regions under the cooling air supply condition during turbine operation using a first algorithm based on the first initial condition; constructing a second initial condition based on the adjusted rotational speed of the turbine rotor and the temperature of the corresponding turbine solid components in the steady-state temperature field, and calculating the basic transient temperature fields of the solid and fluid regions under the cooling air supply condition after the turbine stops using a second algorithm based on the second initial condition, wherein the value of the adjusted rotational speed is zero; and outputting the calculation result that cooling air supply is not required after the turbine stops, provided that the temperature of the target region in the basic transient temperature field does not exceed a predetermined threshold.
[0050] Figure 1 A schematic diagram of a turbine shaft end cooling structure according to an embodiment of the present disclosure is shown.
[0051] like Figure 1 As shown, a two-dimensional cross-sectional view illustrates the main structure and cooling airflow path at the shaft end of a supercritical carbon dioxide radial turbine.
[0052] In such Figure 1 In the turbine shown, the rotor is the rotating component used to convert the energy of the high-temperature, high-pressure carbon dioxide working fluid into mechanical energy and output it. The impeller disk is mounted on the rotor shaft, and its blades directly capture the kinetic energy of the high-speed working fluid, thereby driving the entire rotor to rotate. When the turbine stops after high-speed operation, the impeller stores a large amount of residual heat due to direct contact with the high-temperature working fluid. If the impeller disk and rotor are not cooled, the high temperature continuously accumulating at the rotor shaft end will be transferred along the rotor axis to the other shaft end with a lower temperature, for example... Figure 1 The dry gas seal device shown is shown in the image. When the temperature of the transferred heat accumulates to a level exceeding the temperature threshold of the dry gas seal device, it may affect the safety and service life of the dry gas seal device.
[0053] Therefore, cooling is required for the impeller disk and the rotor shaft end near the impeller disk. Forced convection cooling is achieved by introducing a low-temperature dry sealing gas (such as carbon dioxide) into the shaft end gap as the cooling medium to reduce the temperature of the impeller disk and rotor.
[0054] like Figure 1As shown, the cooling channel is formed by a narrow annular gap between the turbine disk wall and the rotor shaft surface, and another narrow annular gap between the bottom wall of the casing and the axial surface of the rotor. Cooling gas is injected from the channel inlet point 0 and flows sequentially through the axial section. and radial segment Finally, it is discharged from outlet point 3. During this flow process, the cooling gas undergoes convective heat transfer with the high-temperature shaft surface, the bottom wall of the casing, and the wall of the impeller disk near the rotor, continuously carrying away heat and achieving turbine cooling, thereby ensuring the turbine's safety during shutdown.
[0055] Figure 2 A flowchart illustrating a method for calculating turbine cooling gas demand according to an embodiment of the present disclosure is shown. Figure 3 A schematic diagram of a heat transfer geometry model constructed from the various components of a turbine according to an embodiment of the present disclosure is shown.
[0056] like Figure 2 As shown, the turbine cooling gas demand calculation method 100 of this embodiment includes operations S110 to S140, and this transaction processing method can be executed by an electronic device.
[0057] In operation S210, a heat transfer geometry model is constructed based on the turbine's structural parameters. The heat transfer geometry model includes a solid region for simulating temperature changes in the turbine's solid components and a fluid region for simulating temperature changes in the cooling airflow as it passes through the cooling channels.
[0058] In operation S220, a first initial condition is constructed based on the basic speed of the turbine rotor, the initial temperature of the turbine solid components, and the basic parameters of the cooling gas. Based on the first initial condition, the steady-state temperature field of the solid region and the fluid region under the condition of cooling gas flow during turbine operation is calculated using a first algorithm.
[0059] In operation S230, a second initial condition is constructed based on the adjusted speed of the turbine rotor and the temperature of the corresponding solid components of the turbine in the steady-state temperature field. Based on the second initial condition, the basic transient temperature field of the solid region and the fluid region is calculated using the second algorithm when the cooling air supply stops after the turbine stops. The value of the adjusted speed is zero.
[0060] When operating S240, if the temperature of the target region in the basic transient temperature field does not exceed a predetermined threshold, the calculation results will be output, and cooling gas will not need to be continued after the turbine is shut down.
[0061] According to embodiments of this disclosure, constructing a heat transfer geometric model based on the turbine's structural parameters includes: constructing a heat transfer geometric model based on the structural parameters of the turbine's rotor, casing, and disk.
[0062] according to Figure 1The positions and structural parameters of the turbine components shown are used to construct the following... Figure 3 The heat transfer geometry model shown. Figure 3 As shown, the heat transfer geometry model includes multiple solid regions used to simulate the temperature changes of the rotor, casing, and disk walls, such as... Figure 3 The dark areas in the image. And, several fluid regions used to simulate temperature changes as cooling airflow passes through rotor clearance and wheel back clearance, such as... Figure 3 The light-colored area in the image. The rotor clearance is the gap between the bottom of the casing and the rotor, corresponding to... Figure 1 In Axial section; wheel back clearance is the clearance between the side of the casing and the wheel disc, corresponding to Figure 1 In Radial segment.
[0063] According to embodiments of this disclosure, the structural parameters of the turbine include at least: the structural parameters of the rotor and the structural parameters of the disk. The structural parameters of the rotor may include: rotor radius. Rotor clearance dimension, i.e., axial annular clearance Rotor length ; Cross-sectional area of the annular gap The structural parameters of the wheel include: the wheel back clearance dimension. ; Roulette radius Receiver length The turbine's structural parameters may also include: the casing extending beyond the impeller height. Total height of the casing See Table 1 for reference.
[0064] Table 1
[0065] ;
[0066] Figure 4 A schematic diagram of a heat transfer geometry model mesh according to an embodiment of the present disclosure is shown.
[0067] Each component is divided into an independent computing area.
[0068] like Figure 4 As shown, the computational mesh is formed by spatially discretizing the solid region, the fluid region, and the coupling interface between the two, with the axial direction being x and the radial direction being r. Each mesh node represents a location, and variables such as temperature and pressure are solved at these nodes using this mesh. The rotor solid region is used to simulate the temperature field of the rotor body; the casing solid region is used to simulate the temperature field of the casing wall; and the impeller disk solid region is used to simulate the temperature field of the impeller disk wall.
[0069] Based on the different gaps through which the cooling gas flows, it is divided into independent fluid channel regions. For example... Figure 3, Figure 4 As shown, the first fluid region is the rotor clearance flow channel, referring to the gap between the bottom of the casing and the outer surface of the rotor, corresponding to... Figure 1 In In the first fluid region, cooling gas flows axially and exchanges heat, directly cooling the rotor shaft. The second fluid region is the impeller back clearance channel, referring to the gap between the inner wall of the casing near the impeller and the back of the impeller disc. Figure 1 In In the second fluid region, the cooling gas flows radially, directly cooling the impeller disk with the highest temperature and largest heat capacity.
[0070] Specifically, the rotor region is (1:n, 1:m_1), where the axial nodes are i=1~n and the radial nodes are j=1~m_1. The casing region is (1:n-2, m_1+2:m_3), where the axial nodes are i=1~(n-2) and the radial nodes are j=(m_1+2)~m_3. The cooling gas radial flow region is (n-1, m_1+1:m_2), where the axial nodes are i=(n-1)~n and the radial nodes are j=(m_1+1)~m_2. The cooling gas axial flow region is (1:n-1, m_1+1), where the axial nodes are i=1~n-1 and the radial nodes are j=(m_1+1)~(m_1+2). The impeller wall is (n, 1:m_2), where the axial nodes are i=n and the radial nodes are j=1~m_2.
[0071] The axial mesh length of each node in this two-dimensional computational domain is The radial grid length is .
[0072] According to embodiments of this disclosure, constructing the first initial conditions includes: constructing the first initial conditions based on the turbine rotor's base rotational speed, base disc wall temperature, initial temperatures of the rotor and casing, first cooling gas inlet temperature, first cooling gas inlet pressure, and first cooling gas mass flow rate, wherein the base disc wall temperature is equal to the temperature of the working fluid inside the turbine. See Table 2.
[0073] Table 2
[0074] ;
[0075] Specifically, the initial temperatures of the rotor and casing are referenced to formula (1).
[0076] i=1~n; j=1~m_3 Formula (1);
[0077] The first cooling gas inlet temperature reference formula (2);
[0078] i=1; j=m_1+1 Formula (2);
[0079] Reference formula for basic disk wall temperature (3a);
[0080] i=n; j=m_1+1~m_2 Formula (3a);
[0081] The initial temperature of the contact surface between the rotor and the disk is referenced by formula (3b).
[0082] i=n; j=1~m_1 Formula (3b).
[0083] In addition, the boundary surface on the side of the casing that is in direct contact with the working fluid inside the turbine corresponds to the position inside the turbine volute that is in contact with the high-temperature mainstream working fluid, and can be set as the wall temperature condition, as shown in formula (3c).
[0084] i=(n-2), j=(m_2+1)~m_3 Formula (3c).
[0085] According to embodiments of this disclosure, constructing the second initial conditions includes: constructing the second initial conditions based on the adjusted rotational speed of the turbine rotor and the temperatures at the corresponding disk wall, rotor, and casing in the steady-state temperature field.
[0086] According to embodiments of this disclosure, constructing multiple sets of third initial conditions includes: based on the adjusted rotational speed of the turbine rotor, the temperatures at the corresponding disk wall, rotor, and casing in the steady-state temperature field, the second cooling gas inlet temperature, the second cooling gas inlet pressure, and multiple increasing second cooling gas mass flow rates, constructing multiple sets of third initial conditions corresponding to the multiple second cooling gas mass flow rates.
[0087] According to embodiments of this disclosure, the target area is the solid wall surface at the inlet of the cooling channel, which can be a measuring point area in the rotor or casing, such as... Figure 4 Point A (1, m_1) and / or Figure 4 Point B (1, m_1+2) in the equation.
[0088] According to embodiments of this disclosure, the temperature change curve of the target area over time is monitored within the calculated basic transient temperature field. A predetermined threshold is set, representing the highest safe temperature allowed by components such as dry gas seals. If the temperature of the target area does not exceed this safe threshold throughout the entire simulation's resting time, it indicates that the component will not overheat even without cooling air supply. The conclusion is output that cooling air supply is not required after shutdown. Conversely, if the temperature of the target area exceeds the threshold, it indicates that active cooling is required, and the minimum required cooling air volume is calculated.
[0089] According to embodiments of this disclosure, by abstracting the complex turbine cooling process into a simplified heat transfer geometric model consisting of solid and fluid regions, the steady-state temperature fields of the solid and fluid regions are first calculated based on operating parameters. Then, using this steady-state temperature field as accurate initial conditions, transient heat transfer simulations are performed under no-cooling conditions after shutdown. Based on the physical principles of solid heat conduction and convection heat transfer, the temperature change of the target region under natural cooling conditions can be calculated. Furthermore, by comparing the transient calculation results with a preset safe temperature threshold, a judgment result on whether cooling gas injection is necessary can be directly output. The step-by-step calculation strategy, using the steady-state results as the transient initial conditions, ensures the continuity of the physical process and avoids errors caused by assuming a non-uniform temperature field at the beginning of the transient calculation, thus improving the accuracy and consistency of the analysis results. Simultaneously, the NIST database is used to calculate the working fluid properties, avoiding direct solution to the strongly nonlinear property changes of supercritical carbon dioxide working fluid during the transient process, thereby improving computational stability. Compared to traditional simulation methods that require a complete solution of the three-dimensional transient fluid-structure interaction field and direct handling of drastic property changes in the supercritical region, this method reduces computational resources and simulation time.
[0090] According to embodiments of this disclosure, when the temperature of the target region in the basic transient temperature field exceeds a predetermined threshold, multiple sets of third initial conditions are constructed based on the turbine rotor's adjusted rotational speed, the temperature of the corresponding turbine solid component in the steady-state temperature field, and multiple sets of cooling gas adjustment parameters. Based on these multiple sets of third initial conditions, a third algorithm is used to calculate multiple sets of reference transient temperature fields for the solid and fluid regions under the condition of continued cooling gas supply after turbine shutdown. These multiple sets of reference transient temperature fields correspond to multiple sets of cooling gas adjustment parameters. Based on the calculation results of the multiple sets of reference transient temperature fields, the cooling gas adjustment parameter that minimizes the flow rate to ensure the target region temperature does not exceed the predetermined threshold is determined. Based on this cooling gas adjustment parameter, the calculation result of the cooling gas demand after turbine shutdown is output.
[0091] According to embodiments of this disclosure, the natural cooling process without cooling air supply during shutdown revealed that the temperature in the target area exceeded a safe threshold. Further determination is needed to determine the minimum amount of cooling air required to ensure the target area does not exceed the temperature limit. Specifically, multiple sets of cooling air adjustment parameters are introduced, including cooling gas flow rate, inlet pressure, and temperature. For example, the cooling gas flow rate after shutdown can be set in a stepped manner, such as setting the flow rate to 10%, 20%, 30%, 40%, 50%, 60%, 70%, etc., of the reference flow rate for multiple trial calculations. The reference flow rate can be the flow rate of cooling gas supplied to the turbine during its stable operation phase. Based on the constructed third initial conditions, a third algorithm is used to simulate the temperature change over time in the solid and fluid regions under the condition of shutdown but with a specific flow rate of cooling air. Each calculation yields a temperature-time change curve for the target area, i.e., a reference transient temperature field. The initial temperature of the solid components is the steady-state temperature field calculated in operation S120, i.e., the thermal state at the moment of shutdown, to ensure that all calculations start from the same point.
[0092] Analyze the above simulation results. Check whether the temperature curve of the target area remains below the safe threshold throughout the entire static period under each set of cooling air flow rates. Select the set with the smallest flow rate value from all sets of cooling air parameters that meet the temperature control requirements (i.e., do not exceed the threshold). For example, if the simulation shows that a flow rate of 60% of the rated value is just enough to avoid overheating, while 50% would cause overheating, then 60% is the minimum cooling air adjustment parameter. Finally, output the calculated cooling air demand after turbine shutdown.
[0093] According to embodiments of this disclosure, a parametric scanning method is employed to construct multiple sets of boundary conditions based on different cooling gas adjustment parameters and perform parallel transient simulations. This allows for the determination of the influence of key operational variables, such as cooling gas flow rate, on the temperature drop rate of the target region. Secondly, by comparing and analyzing the temperature curves of the target region in multiple sets of reference transient temperature fields, the critical cooling conditions that satisfy predetermined temperature threshold constraints can be accurately identified. Finally, the minimum cooling gas adjustment parameters determined based on these critical conditions provide a quantitative control target for optimal energy consumption while meeting thermal safety requirements in engineering operations.
[0094] The first, second, and third algorithms all use the explicit difference calculation method as an example to describe the calculation process, that is, the temperature field at time k+1 is calculated based on the temperature field at time k.
[0095] According to the embodiments of this disclosure, the calculation of the steady-state temperature field, the basic transient temperature field, and the reference transient temperature field of the solid region and the fluid region all involve multiple iterative calculations. Specifically, for any of the first algorithm, the second algorithm, and the third algorithm, the calculation at any time k+1 includes: based on the temperature calculation results of the solid region and the fluid region at time k, calculating the temperature calculation results of the solid region and the fluid region at time k+1, where k is a positive integer.
[0096] According to embodiments of this disclosure, the solid region and the fluid region are meshed; based on the temperature calculation result at time k, the temperature calculation result at time k+1 of the solid region and the fluid region includes: at time k+1, calculating the temperature of any target mesh with two-dimensional coordinates i and j in the solid region and the fluid region. It is calculated using formula (4).
[0097] Formula (4);
[0098] Where A, B, C, D, E, and F are calculation coefficients. , , , , These refer to the temperature of the target grid with two-dimensional coordinates i and j at time k, and the temperature of multiple adjacent grids surrounding the target grid.
[0099] Boundary conditions characterize the external constraints and initial internal state of the entire heat transfer geometry model at the initial time (i.e., k=1). Specifically, these include the mass flow rate of the cooling gas flowing through the cooling channels. The turbine rotor's rotational speed N and the initial temperature of the disk wall surface. Initial temperature distribution of the rotor (Including its central region and multiple boundary surfaces), initial temperature distribution of the casing (Including its central region and multiple boundary surfaces).
[0100] According to embodiments of this disclosure, the fluid region includes an axial flow region of cooling gas corresponding to the rotor clearance and a radial flow region of cooling gas corresponding to the wheel back clearance; the solid region includes a rotor solid region and a casing solid region; wherein, the rotor solid region includes a rotor center region and a plurality of rotor boundary surfaces, and the casing solid region includes a casing center region and a plurality of casing boundary surfaces.
[0101] In any of the first, second, and third algorithms, the calculation methods for coefficients A, B, C, D, E, and F differ for the axial flow region of cooling gas, the radial flow region of cooling gas, the rotor center region, multiple rotor boundary surfaces, the casing center region, and multiple casing boundary surfaces.
[0102] According to embodiments of this disclosure, due to differences in the physical conditions, dominant heat transfer mechanisms, and physical properties of the solid and fluid regions, different calculation coefficients A, B, C, D, E, and F are required for different regions. Specifically, the governing equations for the solid region are primarily based on heat conduction, and their calculation coefficients are mainly determined by the thermal conductivity, density, and specific heat capacity of the solid material. In contrast, the governing equations for the fluid region require coupling flow and convective heat transfer, and their calculation coefficients need to incorporate fluid dynamic parameters such as flow velocity and convective heat transfer coefficient. Furthermore, under the same operating conditions, the thermal boundary conditions of the central and boundary regions of the same component (such as a rotor or casing) also differ. Specifically, the central region involves only pure conductive heat exchange with adjacent solid meshes, and its calculation coefficients are relatively simple and symmetrical. However, the boundary region is in direct contact with the fluid or other solids, and its calculation coefficients must be coupled with complex boundary condition parameters such as external convective heat transfer coefficients or other thermal resistances, resulting in a different calculation method compared to the central region. Therefore, during the calculation, it is also necessary to distinguish the calculation coefficients of the rotor center region, multiple rotor boundary surfaces, casing center region, and multiple casing boundary surfaces, which are all solid regions.
[0103] The physical properties of the solid materials used in the calculation process are shown in Table 3.
[0104] Table 3
[0105] ;
[0106] According to embodiments of this disclosure, the method of calculation using the first algorithm will be described in detail below.
[0107] (1) For the axial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the first rotor surface. Axial heat transfer coefficient of the first casing Cooling gas mass flow rate Cooling gas constant pressure heat capacity Cooling gas density Rotor radius Rotor clearance dimensions Among them, the convective heat transfer coefficient of the first rotor surface The axial heat transfer coefficient of the first casing represents the convective heat transfer between the cooling gas and the rotor surface during turbine operation. It represents the heat transfer coefficient that represents the convective heat transfer between the cooling gas and the bottom of the casing during turbine operation.
[0108] Specifically, such as Figure 4 As shown, for nodes with axial flow regions of cooling gas i=1~(n-1); j=m_1+1, the calculation method of coefficients A, B, C, D, E, and F is referenced in formula (5).
[0109] (5);
[0110] (2) For the radial flow zone of cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: radial heat transfer coefficient of the first casing. The convective heat transfer coefficient of the first disc surface Cooling gas mass flow rate Cooling gas constant pressure heat capacity Cooling gas density Wheel back clearance dimensions The radial heat transfer coefficient of the first casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the side of the casing during turbine operation, and the convective heat transfer coefficient of the first wheel surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the wheel surface during turbine operation.
[0111] Specifically, such as Figure 4 As shown, for nodes with radial flow regions of cooling gas i=n-1; j=(m_1+1)~m_2, the calculation method of coefficients A, B, C, D, E, and F is referred to formula (6).
[0112] (6);
[0113] Where, r j The radius at the current grid node; The radial dimension of the grid; The radial heat transfer coefficient of the first casing; This refers to the wheel back clearance dimension.
[0114] (3) For the rotor solid region, it is specifically divided into the rotor center region and multiple rotor boundary surfaces. The multiple rotor boundary surfaces are further divided into the target rotor boundary surface and the remaining rotor boundary surfaces. The target rotor boundary surface is the region where the rotor surface is in direct contact with the cooling gas.
[0115] Specifically, such as Figure 4 As shown, the rotor solid region is defined as the node range of 1≤i≤n and 1≤j≤m_1.
[0116] (3.1) For the target rotor boundary surface among multiple rotor boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the first rotor surface. Rotor thermal diffusivity Rotor thermal conductivity The target rotor boundary surface is the area where the rotor surface is in direct contact with the cooling air.
[0117] For the target rotor boundary surface, that is, when the node position is i=1~n and j=m_1, corresponding to the axial convection heat transfer boundary of the cold air on the rotor surface, the calculation method of the coefficients A, B, C, D, E and F is referred to formula (7).
[0118] (7).
[0119] (3.2) For the rotor center region and the other rotor boundary surfaces except the target rotor boundary surface, the calculated coefficients A, B, C, D, E, and F are based on the rotor thermal diffusivity. Calculated.
[0120] For the rotor center region, when the node position is 2≤i≤(n-1) and 2≤j≤(m_1-1), the calculation method for coefficients A, B, C, D, E, and F is based on formula (8):
[0121] (8).
[0122] For the two rotor boundary surfaces other than the target rotor boundary surface, the calculation methods for the coefficients A, B, C, D, E, and F are as follows (9) and (10).
[0123] Among them, when the node position is i=1, j=1~m_1, the calculation method of coefficients A, B, C, D, E, F is referred to formula (9).
[0124] (9);
[0125] When the node position is i=1~n and j=1, the calculation method of coefficients A, B, C, D, E and F is based on formula (10).
[0126] (10).
[0127] (4) For the solid area of the casing, it is specifically divided into the central area of the casing and multiple casing boundary surfaces. The multiple casing boundary surfaces are further divided into two first target casing boundary surfaces (the areas on the side and bottom of the casing that are in direct contact with the cooling gas), the second target casing boundary surface on the side of the casing that is in direct contact with the working fluid inside the turbine, and the remaining casing boundary surfaces other than the first target casing boundary surface and the second target casing boundary surface.
[0128] Specifically, such as Figure 4 As shown, when the node range is 1≤i≤(n-2) and (m_1+2)≤j≤m_3, it is the solid region of the casing.
[0129] (4.1) For two first target casing boundary surfaces among multiple casing boundary surfaces, the calculated coefficients A, B, C, D, E, and F are obtained based on the following parameters: first casing radial heat transfer coefficient Axial heat transfer coefficient of the first casing Casing thermal diffusivity Thermal conductivity of the casing Among them, the two target casing boundary surfaces are the areas on the side and bottom of the casing that are in direct contact with the cooling air.
[0130] Among them, in the two first target casing boundary surfaces mentioned above, the area where the casing side directly contacts the cooling gas, that is, when the node position is i=(n-2) and j=(m_1+2)~m_2, corresponds to the radial convection boundary of the cooling gas. The calculation method of the coefficients A, B, C, D, E, and F is referred to formula (11):
[0131] (11);
[0132] Among them, in the two first target casing boundary surfaces mentioned above, the area where the bottom of the casing is in direct contact with the cooling gas, when the node position is i=1~(n-2) and j=(m_1+2), corresponds to the axial convection boundary of the cooling gas. The calculation method of the coefficients A, B, C, D, E, and F is referred to formula (12):
[0133] (12).
[0134] (4.2) For the second target casing boundary surface in which the casing side directly contacts the working fluid inside the turbine among multiple casing boundary surfaces, that is, when the node position is i=(n-2) and j=(m_2+1)~m_3, corresponding to the position of the main working fluid in contact with the high temperature inside the turbine casing, it can be set as the wall temperature condition. In this case, the calculation method of the calculation coefficients A, B, C, D, E, F refers to formula (13):
[0135] A= D=1, B=C=E=F=0 (13).
[0136] Furthermore, the temperature at the boundary of the second target casing is equal to the temperature of the working fluid inside the turbine, that is... It should be noted that the temperature of the second target casing boundary surface is defined here as a boundary condition for the temperature field calculation.
[0137] (4.3) For the central region of the casing and the remaining casing boundary surfaces other than the first target casing boundary surface and the second target casing boundary surface, the calculated coefficients A, B, C, D, E, and F are based on the casing thermal diffusivity. Calculated.
[0138] For the central region of the casing, i.e., when the node position is 2≤i≤(n-3) and (m_1+3)≤j≤(m_3-1), the calculation method for coefficients A, B, C, D, E, and F is based on formula (14):
[0139] (14).
[0140] For the remaining casing boundary surfaces, the calculation methods for coefficients A, B, C, D, E, and F are referenced in formulas (15) and (16).
[0141] Specifically, when the node positions are i=1 and j=(m_1+2)~m_3:
[0142] (15);
[0143] Specifically, when the node positions are i=1~(n-2) and j=m_3:
[0144] (16).
[0145] (5) Furthermore, in the first algorithm, the temperature at the interface between the rotor and the disk, and the temperature at the interface between the cooling gas and the disk, are equal to the temperature of the working fluid inside the turbine.
[0146] The following sections introduce the methods for calculating the surface convective heat transfer coefficient at the fluid-solid interface in the axial flow section and at the fluid-solid interface in the radial flow section, respectively.
[0147] (1) For the fluid-solid interface of the axial flow section, the convective heat transfer coefficient of the first rotor surface can be assumed. Axial heat transfer coefficient of the first casing They are equal. It should be noted that the convective heat transfer coefficient on the surface of the first rotor is... Axial heat transfer coefficient of the first casing They can also be unequal.
[0148] The convective heat transfer coefficient of the first rotor surface Axial heat transfer coefficient of the first casing Based on the first Nusselt number Thermal conductivity of cooling gas Rotor length Calculated.
[0149] First Nusel number Based on the first rotating Reynolds number at the rotor clearance Axial Reynolds number at rotor clearance Rotor length Rotor radius Rotor clearance dimensions (length-to-width ratio) gap ratio Cooling aerodynamic viscosity Cooling gas constant pressure heat capacity Thermal conductivity of cooling gas The calculated values are: cooling gas dynamic viscosity, cooling gas isobaric heat capacity, cooling gas thermal conductivity, and Prandtl number. related.
[0150] Specifically, the convective heat transfer coefficient of the first rotor surface The convective heat transfer coefficient between the rotor and the cooling fluid (in W / (m²·K)) can be obtained through the first Nusselt number. The calculation yields the result, referring to formula (17):
[0151] (17);
[0152] First Nusel number
[0153] The calculation can be performed using formula (18):
[0154] (18);
[0155] in, The first effective Reynolds number is represented by formula (19), which is calculated using the following method:
[0156] Formula (19).
[0157] in, The first axial Reynolds number at the rotor clearance. This is the first rotating Reynolds number at the rotor clearance.
[0158] The formula is applicable to Taylor numbers. ( );
[0159] First axial Reynolds number at rotor clearance ;
[0160] Aspect Ratio gap ratio .
[0161] If the scope of application is outside this range, other calculation methods may be used.
[0162] The convective heat transfer coefficient of the first rotor surface It can be assumed that the axial heat transfer coefficient with the first casing is... The values are the same, refer to formula (20):
[0163] (20).
[0164] For the parameters and calculation methods included in formulas (17) to (20), please refer to Table 4.
[0165] Table 4
[0166]
[0167] (2) For the flow-solid interface of the radial flow section, the radial heat transfer coefficient of the first casing can be assumed. The convective heat transfer coefficient of the first disc surface They are equal. It should be noted that the radial heat transfer coefficient of the first casing is equal. The convective heat transfer coefficient of the first disc surface They can also be unequal.
[0168] According to embodiments of this disclosure, the radial heat transfer coefficient of the first casing is... The convective heat transfer coefficient of the first disc surface Based on the second Nusselt number Thermal conductivity of cooling gas Wheel radius The calculated second Nusselt number Based on the first rotational Reynolds number at the turbine disk The first radial Reynolds number at the turbine disk Wheel back clearance dimensions roulette radius Calculated.
[0169] Specifically, the radial heat transfer coefficient of the first casing The calculation method for the convective heat transfer coefficient (in W / (m²·K)) between the impeller disk and the cooling fluid is given by formula (21):
[0170] (twenty one).
[0171] Second Nusel number The calculation can be performed using formula (22):
[0172] (twenty two).
[0173] in, The first radial Reynolds number at the turbine disk, is the Reynolds number of the first rotation at the roulette wheel.
[0174] For the other parameters and calculation methods in formulas (21) to (22), please refer to Table 5.
[0175] Table 5
[0176]
[0177] The embodiments disclosed herein are merely exemplary recommendations for using equations (17) to (22) to calculate the convective heat transfer coefficient during the turbine's rotating operation phase, but are not limited to these calculation methods. Any applicable method can be selected according to the actual operating conditions.
[0178] Furthermore, calculating the convective heat transfer coefficient requires knowledge of the thermal properties of the cooling gas fluid. Since the properties of supercritical carbon dioxide change drastically with temperature, the fluid properties are updated in each iteration to ensure simulation accuracy.
[0179] Specifically, in k= Within the time frame, The time step begins to be calculated iteratively hour by hour.
[0180] In this process, the physical property parameters are updated once after each iteration. For example, after calculating the temperature field of the fluid region at time k, the physical property parameters at time k+1 are determined based on the temperature of the fluid region at time k.
[0181] The physical properties can be updated by calling REFPROP, see formula (23).
[0182] (twenty three)
[0183] in, The density updated at time k+1 (kg / m³) 3 ), isobaric heat capacity (J / (kg·K)), thermal conductivity (W / (m·K)) and dynamic viscosity (Pa·s); Let K be the temperature of each node in the fluid region at time k, in K.
[0184] The specific method for calculating using the first algorithm includes the following steps:
[0185] (1) In Within the time frame, During the hourly iteration calculation of the time step, the initial condition is Equation (1), and the boundary conditions are Equations (2), (3a), (3b), and (3c).
[0186] (2) Fluid region temperature update calculation:
[0187] (2.1) Calculation of radial flow cooling gas fluid: Based on the fluid temperature in the radial gap at time k, the thermal properties of the cooling gas are retrieved by equation (23); the local convective heat transfer coefficient at each radial position is calculated, as shown in equations (21) and (22); then the fluid temperature in the radial gap at time k+1 is updated according to equations (4) and (6).
[0188] (2.2) Calculation of axial flow fluid: Based on the fluid temperature in the axial gap at time k, the heat properties of the cooling gas are retrieved by equation (23); the local convective heat transfer coefficient at each axial position is calculated, as shown in equations (17) to (20); then the fluid temperature in the axial gap at time k+1 is updated according to equations (4) and (5).
[0189] (3) Calculation of heat conduction process in the rotor solid region:
[0190] (3.1) Rotor boundary update: According to equations (4), (7), (9), (10), update the temperature of each node on the rotor boundary at time k+1.
[0191] (3.2) Heat conduction inside the rotor: According to equations (4) and (8), update the temperature of each node inside the rotor at time k+1.
[0192] (4) Calculation of heat conduction process in the solid region of the casing:
[0193] (4.1) Casing boundary conditions: According to equations (4), (11), (12), (13), (15), (16), update the temperature of each node on the casing boundary at time k+1.
[0194] (4.2) Heat conduction inside the casing: According to equations (4) and (14), update the temperature of each node inside the casing at time k+1.
[0195] (5) Repeat steps (2) to (4) above until convergence is achieved.
[0196] If the temperature field difference between two adjacent time points satisfies equation (24), the calculation is considered convergent, and the iteration can be exited.
[0197] (twenty four).
[0198] (6) Output the calculation results. After the calculation converges, the two-dimensional temperature field of the entire computational domain can be output, including the temperature distribution of the fluid domain, rotor and casing solid domain, to provide initial values for subsequent calculations.
[0199] After the calculation converges, the two-dimensional temperature field T of the entire computational domain can be output, including the temperature distribution of the fluid domain, rotor and casing solid domain, providing initial values for subsequent calculations.
[0200] According to an embodiment of this disclosure, if cooling continues after the turbine is shut down, the temperature field is calculated according to the third algorithm.
[0201] The third algorithm and the aforementioned first algorithm both involve convective heat transfer in the fluid-structure interaction region and conduction heat transfer within the solid regions of the rotor and casing. Therefore, the third algorithm and the first algorithm share the same algorithm principle for calculating coefficients A, B, C, D, E, and F, and still use equations (4) to (16) to calculate the transient temperature field of the two-dimensional fluid-structure interaction computational domain. However, the calculation method for the convective heat transfer coefficient needs to be changed to an empirical formula under static conditions. In addition, the temperature conditions on some boundary surfaces need to be updated from isothermal conditions to adiabatic conditions. These will be explained separately below.
[0202] Several boundary conditions need to be updated.
[0203] First, in the third algorithm, after the turbine stops, there is still cold air exchanging heat between the rotor and stator. The temperature boundary at the cooling gas inlet node needs to be updated from formula (2) to the following formula (25-1).
[0204] i=1; j=m_1+1 (25-1);
[0205] Second, for the second target casing boundary surface in which the casing side directly contacts the working fluid inside the turbine among multiple casing boundary surfaces, i.e., when the node positions are i=(n-2) and j=(m_2+1)~m_3, the condition changes from isothermal conditions to adiabatic boundaries. The calculation coefficients A, B, C, D, E, and F are obtained based on the casing thermal diffusivity. The specific calculation method is shown in the following formula (25-2):
[0206] (25-2).
[0207] Third, for the boundary surfaces of the rotor and the disk, and the boundary surfaces of the cooling gas and the disk, i.e. when i=n, j=1~m_2, the condition changes from isothermal conditions to adiabatic boundaries. The calculation coefficients A, B, C, D, E, and F are obtained based on the rotor thermal diffusivity.
[0208] For specific calculation methods, please refer to the following formula (25-3):
[0209] (25-3).
[0210] Compared to the first algorithm, except for the updated algorithm for the boundary conditions, the third algorithm calculates the coefficients A, B, C, D, E, and F in other regions as follows:
[0211] Specifically, for the axial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the second rotor surface, the axial heat transfer coefficient of the second casing, the mass flow rate of the cooling gas, the isobaric heat capacity of the cooling gas, the density of the cooling gas, the rotor radius, and the rotor clearance dimension; wherein, the convective heat transfer coefficient of the second rotor surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the rotor surface after the turbine is shut down, and the axial heat transfer coefficient of the second casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the bottom of the casing after the turbine is shut down.
[0212] For the radial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are obtained based on the following parameters: radial heat transfer coefficient of the second casing, convective heat transfer coefficient of the second wheel surface, mass flow rate of the cooling gas, isobaric heat capacity of the cooling gas, density of the cooling gas, and wheel back clearance dimension. The radial heat transfer coefficient of the second casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the side of the casing after the turbine is shut down, and the convective heat transfer coefficient of the second wheel surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the wheel surface after the turbine is shut down.
[0213] For the target rotor boundary surface among multiple rotor boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the second rotor surface, the rotor thermal diffusivity, and the rotor thermal conductivity; where the target rotor boundary surface is the area where the rotor surface is in direct contact with the cooling gas.
[0214] For the rotor center region and the other rotor boundary surfaces other than the target rotor boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the rotor thermal diffusivity.
[0215] For two first target casing boundary surfaces among multiple casing boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: second casing radial heat transfer coefficient, second casing axial heat transfer coefficient, casing thermal diffusivity, and casing thermal conductivity; wherein, the two first target casing boundary surfaces are the areas of the casing side and the casing bottom that are in direct contact with the cooling gas, respectively.
[0216] For the central region of the casing and the remaining casing boundary surfaces other than the first target casing boundary surface and the second target casing boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the casing thermal diffusivity.
[0217] In the above algorithms for calculating coefficients A, B, C, D, E, and F, the third algorithm and the first algorithm have the same algorithm principle for calculating coefficients A, B, C, D, E, and F. They still use equations (4) to (16) to calculate the transient temperature field of the two-dimensional fluid-structure interaction computational domain, which will not be repeated here. It should be noted that equation (13) needs to be replaced with equation (25-2). However, the calculation method of the convective heat transfer coefficient needs to be replaced with an empirical formula under static conditions. Therefore, the calculation method of the surface convective heat transfer coefficient related to the heat transfer between the rotor and the fluid, and between the casing and the fluid, is updated from equations (17) to (22) to the following equations (26) to (28).
[0218] The following sections introduce the methods for calculating the surface convective heat transfer coefficient at the fluid-solid interface in the axial flow section and at the fluid-solid interface in the radial flow section, respectively.
[0219] (1) For the fluid-solid interface of the axial flow section, the convective heat transfer coefficient of the second rotor surface can be assumed. axial heat transfer coefficient of the second casing equal.
[0220] Specifically, the convective heat transfer coefficient of the second rotor surface axial heat transfer coefficient of the second casing Based on the third Nusselt number Thermal conductivity of cooling gas Rotor length
[0221] The third Nusselt number is calculated based on the axial Reynolds number at the rotor clearance, the dynamic viscosity of the cooling gas, the constant pressure heat capacity of the cooling gas, and the thermal conductivity of the cooling gas.
[0222] Among them, the first Nusselt number Updated to the third Nusselt number Specifically, formula (18) is updated to formula (26).
[0223] Formula (26);
[0224] in, The second axial Reynolds number at the updated rotor clearance.
[0225] Depending on the actual working conditions, other convective heat transfer calculation methods applicable to static states and annular gaps can be selected.
[0226] The convective heat transfer coefficient of the second rotor surface axial heat transfer coefficient of the second casing Calculate using formulas (17) and (20).
[0227] (2) For the flow-solid interface of the radial flow section, the radial heat transfer coefficient of the second casing can be considered as... 2. Convection heat transfer coefficient of the second disc surface equal.
[0228] Specifically, the radial heat transfer coefficient of the second casing and the convective heat transfer coefficient of the second disc surface are based on the fourth Nusselt number. Thermal conductivity of cooling gas Wheel back clearance dimensions The fourth Nusselt number is calculated based on the radial Reynolds number at the turbine disk and the viscosity of the cooling aerodynamics. Cooling gas constant pressure heat capacity Thermal conductivity of cooling gas Calculated.
[0229] Among them, the second Nusselt number Updated to the fourth Nusselt number Specifically, formula (22) is updated to formula (27).
[0230] Formula (27);
[0231] in, The second radial Reynolds number at the updated turbine disk.
[0232] Depending on the actual working conditions, other calculation methods applicable to convective heat transfer in static states and disc-shaped gaps can be selected.
[0233] Second casing radial heat transfer coefficient 2. Convection heat transfer coefficient of the second disc surface It can be calculated according to formula (28).
[0234] Formula (28);
[0235] According to an embodiment of this disclosure, if the cooling air supply is stopped after the turbine is shut down, the temperature field is calculated according to the second algorithm.
[0236] In the second algorithm, if the cooling air supply stops after the turbine shuts down, the heat exchange mode between the cooling air and the rotor and casing will change from convection to conduction. This requires updating the calculation boundary conditions and the algorithms for calculating coefficients A, B, C, D, E, and F. These will be explained below.
[0237] (1) For the axial flow region and the radial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: thermal conductivity of cooling gas, dynamic viscosity of cooling gas, constant pressure heat capacity of cooling gas, and density of cooling gas.
[0238] Specifically, such asFigure 4 As shown, for nodes with cooling gas axial flow regions of i=1~(n-1); j=m_1+1, compared to the first algorithm, the calculation method of coefficients A, B, C, D, E, and F is updated to formula (29).
[0239] (29).
[0240] Specifically, such as Figure 4 As shown, for nodes in the radial flow region of cooling gas i=n-1, j=(m_1+1)~m_2, the calculation method of coefficients A, B, C, D, E, and F is the same as that of formula (29).
[0241] (2) For the rotor solid region, it is specifically divided into the rotor center region and multiple rotor boundary surfaces. The multiple rotor boundary surfaces are further divided into the target rotor boundary surface and the remaining rotor boundary surfaces. The target rotor boundary surface is the region where the rotor surface is in direct contact with the cooling gas.
[0242] Specifically, for the rotor center region, the target rotor boundary surface among multiple rotor boundary surfaces, and the remaining rotor boundary surfaces other than the target rotor boundary surface, the calculated coefficients A, B, C, D, E, and F are based on the rotor thermal diffusivity. Calculated.
[0243] (2.1) Among them, the target rotor boundary surface between the outer periphery of the rotor and the cold air: i=1~n, j=m_1, is transformed from a convective heat transfer boundary to a conductive heat transfer boundary. The calculation methods of the coefficients A, B, C, D, E, and F are also referred to formula (30).
[0244] (30).
[0245] (2.2) For the rotor center region and the other rotor boundary surfaces except the target rotor boundary surface, the calculated coefficients A, B, C, D, E, and F are based on the rotor thermal diffusivity. The calculation is based on the algorithm of the aforementioned formula (8).
[0246] (2.3) The boundary surfaces of the rotor and the disk, as well as the boundary surfaces of the cooling gas and the disk, are set as adiabatic boundaries. That is, when i=n and j=1~m_2, compared with the first algorithm, the condition of constant wall temperature is changed to adiabatic boundary. The calculation coefficients A, B, C, D, E, and F are based on the rotor thermal diffusivity. The calculation is obtained by referring to formula (31). The algorithm for formula (31) is the same as that for the aforementioned formula (25-3).
[0247] (31);
[0248] (3) For the solid area of the casing, it is specifically divided into the central area of the casing and multiple casing boundary surfaces. The multiple casing boundary surfaces are further divided into two first target casing boundary surfaces (the areas of the casing side and bottom that are in direct contact with the cooling gas), a second target casing boundary surface that is in direct contact with the working fluid inside the turbine, and the remaining casing boundary surfaces other than the first and second target casing boundary surfaces. These are described in detail below:
[0249] (3.1) For the boundary where the side of the casing contacts the high-temperature working fluid, that is, the second target casing boundary surface in which the side of the casing directly contacts the working fluid inside the turbine among multiple casing boundary surfaces, that is, when the node positions are i=(n-2) and j=(m_2+1)~m_3, it does not contact the high-temperature working fluid after shutdown. The calculated coefficients A, B, C, D, E, and F are based on the casing thermal diffusivity. The calculation is obtained by referring to formula (32); the algorithm of formula (32) is the same as that of the aforementioned formula (25-2).
[0250] (32).
[0251] (3.2) For the two first target casing boundary surfaces (the areas on the casing side and the bottom of the casing that are in direct contact with the cooling air), the calculated coefficients A, B, C, D, E, and F are based on the casing thermal diffusivity. Calculated.
[0252] Specifically, for the interface between the left side of the casing and the cold air, i.e., when i=(n-2) and j=(m_1+2)~m_2, corresponding to the radial convection boundary of the cold air, after shutdown, the heat transfer needs to be changed from convection to internal conduction; and for the interface between the bottom of the casing and the cold air, i.e., when i=1~(n-2) and j=(m_1+2), corresponding to the axial convection boundary of the cold air, after shutdown, the heat transfer needs to be changed from convection to internal conduction, the algorithm for calculating coefficients A, B, C, D, E, and F is referenced in formula (33);
[0253] (33).
[0254] (3.3) For the central region of the casing, i.e. when the node position is 2≤i≤(n-3) and (m_1+3)≤j≤(m_3-1), the coefficients A, B, C, D, E and F are calculated based on the thermal diffusivity of the casing. The calculation method is the same as the formula (14) mentioned above.
[0255] (3.4) For the other casing boundary surfaces besides the first target casing boundary surface and the second target casing boundary surface, i.e. when the node position is i=1, j=(m_1+2)~m_3, the calculation coefficients A, B, C, D, E, and F are calculated based on the casing thermal diffusivity. The calculation method is the same as the above formulas (15) and (16).
[0256] The general process of calculation using the first, second, and third algorithms described above is explained below.
[0257] (1) First, the steady-state temperature fields of the solid region and the fluid region under the cooling gas condition during turbine operation are calculated using the first algorithm, using formulas (4) to (22). Whether the second rotational speed N is 0 is used to determine whether the shutdown state has been reached. After entering the shutdown state, the initial temperature field of the calculation domain is the result calculated using the first algorithm.
[0258] (2) Then, calculate the transient temperature field when the cooling air supply is stopped. If the cooling air mass flow rate... If no cooling gas is supplied after shutdown, the second algorithm is used to calculate the transient temperature field of the two-dimensional fluid-structure interaction computational domain, and the updated formulas are (29) to (33).
[0259] Furthermore, it outputs the temperature of a specific point of interest (target region) in the solid domain. Data or curves changing over time, to determine whether it is in The temperature exceeded the limit within a certain time period, such as exceeding the temperature that the sealing ring can withstand, which is 200℃.
[0260] If there is no overheating, there is no need to run the cooling system again after shutdown.
[0261] (3) If overheating occurs, subsequent steps can be performed to estimate the cooling demand. The transient temperature field of the two-dimensional fluid-structure interaction computational domain is then calculated according to the third algorithm.
[0262] When using the third algorithm, the boundary conditions need to be updated according to formulas (25-1), (25-2), and (25-3), and the cold air mass flow rate needs to be set. Import pressure and temperature The temperature field is calculated using formulas (4) to (16) and (26) to (28) (where formula (13) is replaced by formula (25-2)). Formulas (26) to (28) indicate that the local convective heat transfer coefficient of the fluid domain should be calculated using a method applicable under shutdown conditions.
[0263] If, in the calculated temperature field, the temperature at a certain point of interest (target area) in the solid domain still exceeds the set cooling air volume, then the parameters such as the cooling air flow rate are updated. This process is repeated until the cooling requirements are met and the overheating phenomenon is eliminated.
[0264] (4) After the calculation is completed, the output of the turbine shutdown status can be generated. The two-dimensional temperature field T at any given time includes the temperature distribution of the fluid domain, rotor, and casing solid domain. It can also output data or curves showing the temperature change of the target region in the solid domain over time after the turbine stops operating without cooling, and under different cooling flow rates. Specifically, this can be: the steady-state temperature field T of the fluid, rotor, and casing after sufficient heat exchange during turbine operation; and the temperature field T after turbine shutdown. Within a given timeframe, without cooling gas supply, the transient temperature field in the two-dimensional fluid-solid computational domain, and the temperature variation data or curves of the points of interest in the solid domain over time; after turbine shutdown. Within a given time interval, under conditions where cooling gas flows at different rates, the transient temperature field of the fluid-solid two-dimensional computational domain, and the temperature change data or curves of the points of interest in the solid domain over time.
[0265] The embodiments of this disclosure have been described above. However, these embodiments are for illustrative purposes only and are not intended to limit the scope of this disclosure. Although various embodiments have been described above, this does not mean that the measures in the various embodiments cannot be used advantageously in combination. Various substitutions and modifications can be made by those skilled in the art without departing from the scope of this disclosure, and all such substitutions and modifications should fall within the scope of this disclosure.
Claims
1. A method for calculating turbine cooling air demand executed by electronic equipment, comprising: A heat transfer geometric model is constructed based on the structural parameters of the turbine, wherein the heat transfer geometric model includes a solid region for simulating the temperature change of the solid components of the turbine, and a fluid region for simulating the temperature change of the cooling airflow as it passes through the cooling channel. The first initial conditions are constructed based on the basic speed of the turbine rotor, the initial temperature of the turbine solid components, and the basic parameters of the cooling gas. Based on the first initial conditions, the steady-state temperature fields of the solid and fluid regions under the condition of cooling gas flow during turbine operation are calculated using the first algorithm. Based on the adjusted rotational speed of the turbine rotor and the temperature of the corresponding solid turbine component in the steady-state temperature field, a second initial condition is constructed. Based on the second initial condition, the basic transient temperature field of the solid region and the fluid region under the condition of stopping cooling air after the turbine stops is calculated using a second algorithm. The value of the adjusted rotational speed is zero. If the temperature of the target region in the basic transient temperature field does not exceed a predetermined threshold, the calculation results will be output so that cooling gas does not need to be continued after the turbine is shut down.
2. The method according to claim 1, wherein, The method further includes: When the temperature of the target region in the basic transient temperature field exceeds a predetermined threshold, multiple sets of third initial conditions are constructed based on the turbine rotor's adjusted speed, the temperature of the corresponding turbine solid component in the steady-state temperature field, and multiple sets of cooling gas adjustment parameters. Based on these multiple sets of third initial conditions, a third algorithm is used to calculate multiple sets of reference transient temperature fields for the solid and fluid regions under the condition of continued cooling gas supply after turbine shutdown. These multiple sets of reference transient temperature fields correspond to multiple sets of cooling gas adjustment parameters. Based on the calculation results of multiple sets of reference transient temperature fields, the cooling gas adjustment parameters that minimize the flow rate to ensure that the temperature in the target area does not exceed a predetermined threshold are determined, and the cooling gas demand calculation results after turbine shutdown are output based on the cooling gas adjustment parameters.
3. The method according to claim 1 or 2, wherein, The heat transfer geometric model based on the turbine's structural parameters includes: constructing a heat transfer geometric model based on the structural parameters of the turbine's rotor, casing, and disk; The heat transfer geometry model includes multiple solid regions for simulating temperature changes on the rotor, casing, and wheel disc walls, and multiple fluid regions for simulating temperature changes as cooling airflow passes through the rotor gap and wheel back gap. The rotor gap is the gap between the bottom of the casing and the rotor, and the wheel back gap is the gap between the side of the casing and the wheel disc.
4. The method according to claim 3, wherein, The first initial conditions are constructed based on the turbine rotor's base speed, base wheel disk wall temperature, initial temperatures of the rotor and casing, first cooling gas inlet temperature, first cooling gas inlet pressure, and first cooling gas mass flow rate. The base wheel disk wall temperature is equal to the temperature of the working fluid inside the turbine. The second initial conditions are constructed based on the adjusted rotational speed of the turbine rotor and the temperatures at the corresponding disk wall, rotor, and casing in the steady-state temperature field. The construction of multiple sets of third initial conditions includes: based on the adjusted rotational speed of the turbine rotor, the temperatures at the corresponding disk wall, rotor and casing in the steady-state temperature field, the second cooling gas inlet temperature, the second cooling gas inlet pressure, and multiple increasing second cooling gas mass flow rates, constructing multiple sets of third initial conditions corresponding to the multiple second cooling gas mass flow rates.
5. The method according to claim 3, wherein, The calculation of the steady-state temperature field, the basic transient temperature field, and the reference transient temperature field of the solid region and the fluid region all involve multiple iterative calculations. Specifically, for any of the first, second, and third algorithms, the calculation at any time k+1 includes: Based on the temperature calculation results of the solid region and the fluid region at time k, the temperature calculation results of the solid region and the fluid region at time k+1 are calculated, where k is a positive integer.
6. The method according to claim 5, wherein, The method further includes meshing the solid region and the fluid region; Based on the temperature calculation result at time k, the temperature calculation results for the solid region and the fluid region at time k+1 include: at time k+1, calculating the temperature of any target grid with two-dimensional coordinates i and j in the solid region and the fluid region. The following methods are used: ; Where A, B, C, D, E, and F are calculation coefficients. , , , , These refer to the temperature of the target grid with two-dimensional coordinates i and j at time k, and the temperature of multiple adjacent grids surrounding the target grid.
7. The method according to claim 6, wherein, The fluid region includes an axial flow region of cooling gas corresponding to the rotor clearance and a radial flow region of cooling gas corresponding to the wheel back clearance. The solid region includes a rotor solid region and a casing solid region; wherein, the rotor solid region includes a rotor center region and multiple rotor boundary surfaces, and the casing solid region includes a casing center region and multiple casing boundary surfaces; In any of the first, second, and third algorithms, the calculation methods for the calculation coefficients A, B, C, D, E, and F are different for the axial flow region of cooling gas, the radial flow region of cooling gas, the rotor center region, multiple rotor boundary surfaces, the casing center region, and multiple casing boundary surfaces.
8. The method according to claim 7, wherein, In the first algorithm: For the axial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the first rotor surface, the axial heat transfer coefficient of the first casing, the mass flow rate of the cooling gas, the isobaric heat capacity of the cooling gas, the density of the cooling gas, the rotor radius, and the rotor clearance dimension; wherein, the convective heat transfer coefficient of the first rotor surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the rotor surface during turbine operation, and the axial heat transfer coefficient of the first casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the bottom of the casing during turbine operation; For the radial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: radial heat transfer coefficient of the first casing, convective heat transfer coefficient of the first wheel surface, mass flow rate of the cooling gas, isobaric heat capacity of the cooling gas, density of the cooling gas, and wheel back clearance dimension; the radial heat transfer coefficient of the first casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the side of the casing during turbine operation, and the convective heat transfer coefficient of the first wheel surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the wheel surface during turbine operation; For a target rotor boundary surface among multiple rotor boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the first rotor surface, the rotor thermal diffusivity, and the rotor thermal conductivity; wherein, the target rotor boundary surface is the region of the rotor surface that is in direct contact with the cooling gas. For the rotor center region and the other rotor boundary surfaces other than the target rotor boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the rotor thermal diffusivity. For two first target casing boundary surfaces among multiple casing boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: first casing radial heat transfer coefficient, first casing axial heat transfer coefficient, casing thermal diffusivity, and casing thermal conductivity; wherein, the two target casing boundary surfaces are the areas of the casing side and the casing bottom that are in direct contact with the cooling gas, respectively. For the second target casing boundary surface, which is in direct contact with the working fluid inside the turbine among multiple casing boundary surfaces, A=D=1, B=C=E=F=0; the temperature of the second target casing boundary surface is equal to the temperature of the working fluid inside the turbine. For the central region of the casing and the remaining casing boundary surfaces other than the first target casing boundary surface and the second target casing boundary surface, the calculation coefficients A, B, C, D, E, and F are calculated based on the casing thermal diffusivity. The temperature at the interface between the rotor and the disk, and the temperature at the interface between the cooling gas and the disk, are equal to the temperature of the working fluid inside the turbine.
9. The method according to claim 8, wherein: The first rotor surface convective heat transfer coefficient and the first casing axial heat transfer coefficient are calculated based on the first Nusselt number, the cooling gas thermal conductivity, and the rotor length; the first Nusselt number is calculated based on the rotating Reynolds number at the rotor gap, the axial Reynolds number at the rotor gap, the rotor length, the rotor radius, the rotor gap size, the cooling gas dynamic viscosity, the cooling gas isobaric heat capacity, and the cooling gas thermal conductivity. The radial heat transfer coefficient of the first casing and the convective heat transfer coefficient of the first wheel surface are calculated based on the second Nusselt number, the thermal conductivity of the cooling gas, and the wheel radius; the second Nusselt number is calculated based on the rotating Reynolds number at the turbine wheel, the radial Reynolds number at the turbine wheel, the wheel back clearance dimension, and the wheel radius.
10. The method according to claim 7, wherein, In the third algorithm: For the axial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the second rotor surface, the axial heat transfer coefficient of the second casing, the mass flow rate of the cooling gas, the isobaric heat capacity of the cooling gas, the density of the cooling gas, the rotor radius, and the rotor clearance dimension; wherein, the convective heat transfer coefficient of the second rotor surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the rotor surface after the turbine is shut down, and the axial heat transfer coefficient of the second casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the bottom of the casing after the turbine is shut down; For the radial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: radial heat transfer coefficient of the second casing, convective heat transfer coefficient of the second wheel surface, mass flow rate of the cooling gas, isobaric heat capacity of the cooling gas, density of the cooling gas, and wheel back clearance dimension; the radial heat transfer coefficient of the second casing represents the heat transfer coefficient of convective heat transfer between the cooling gas and the side of the casing after the turbine stops, and the convective heat transfer coefficient of the second wheel surface represents the heat transfer coefficient of convective heat transfer between the cooling gas and the wheel surface after the turbine stops; For a target rotor boundary surface among multiple rotor boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: the convective heat transfer coefficient of the second rotor surface, the rotor thermal diffusivity, and the rotor thermal conductivity; wherein, the target rotor boundary surface is the region of the rotor surface that is in direct contact with the cooling gas; For the rotor center region and the other rotor boundary surfaces other than the target rotor boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the rotor thermal diffusivity. For two first target casing boundary surfaces among multiple casing boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: second casing radial heat transfer coefficient, second casing axial heat transfer coefficient, casing thermal diffusivity, and casing thermal conductivity; wherein, the two first target casing boundary surfaces are the areas of the casing side and the casing bottom that are in direct contact with the cooling gas, respectively. For the second target casing boundary surface, which is in direct contact with the working fluid inside the turbine among multiple casing boundary surfaces, it is set as an adiabatic boundary. The calculation coefficients A, B, C, D, E, and F are calculated based on the casing thermal diffusivity. For the central region of the casing and the remaining casing boundary surfaces other than the first target casing boundary surface and the second target casing boundary surface, the calculation coefficients A, B, C, D, E, and F are calculated based on the casing thermal diffusivity. The boundary surfaces where the rotor and the disk contact each other, as well as the boundary surfaces where the cooling gas and the disk contact each other, are set as adiabatic boundaries. The calculation coefficients A, B, C, D, E, and F are calculated based on the rotor thermal diffusivity.
11. The method of claim 10, wherein: The second rotor surface convective heat transfer coefficient and the second casing axial heat transfer coefficient are calculated based on the third Nusselt number, the cooling gas thermal conductivity, and the rotor length; the third Nusselt number is calculated based on the axial Reynolds number at the rotor gap, the cooling gas dynamic viscosity, the cooling gas isobaric heat capacity, and the cooling gas thermal conductivity. The radial heat transfer coefficient of the second casing and the convective heat transfer coefficient of the second wheel surface are calculated based on the fourth Nusselt number, the thermal conductivity of the cooling gas, and the wheel back clearance size; the fourth Nusselt number is calculated based on the radial Reynolds number at the turbine wheel, the dynamic viscosity of the cooling gas, the isobaric heat capacity of the cooling gas, and the thermal conductivity of the cooling gas.
12. The method according to claim 7, wherein, In the second algorithm: For the axial flow region and the radial flow region of the cooling gas, the calculation coefficients A, B, C, D, E, and F are calculated based on the following parameters: thermal conductivity of the cooling gas, constant pressure heat capacity of the cooling gas, and density of the cooling gas. For the rotor center region, the target rotor boundary surface among multiple rotor boundary surfaces, and the remaining rotor boundary surfaces other than the target rotor boundary surface, the calculation coefficients A, B, C, D, E, and F are obtained based on the rotor thermal diffusivity. For the central region of the casing, two first target casing boundary surfaces among multiple casing boundary surfaces, a second target casing boundary surface that is in direct contact with the working fluid inside the turbine on the side of the casing, and the remaining casing boundary surfaces other than the first and second target casing boundary surfaces, the calculation coefficients A, B, C, D, E, and F are calculated based on the thermal diffusivity of the casing. The boundary surfaces where the rotor and the disk contact each other, as well as the boundary surfaces where the cooling gas and the disk contact each other, are set as adiabatic boundaries. The calculation coefficients A, B, C, D, E, and F are calculated based on the rotor thermal diffusivity.
13. The method according to claim 1, wherein: The target area is the solid wall surface at the inlet of the cooling channel.