A numerical method for phase-change transpiration cooling coupled with chemical reactions at supersonic conditions

Abstract

The application discloses a numerical calculation method for phase transition sweating cooling of coupling cracking reaction under supersonic speed conditions, comprising the following steps: a geometric model containing a high-speed outer flow field and a porous field is established, and the interface of the two is defined; the high-speed outer flow field and the porous field are solved by alternately iterating until the global calculation converges; wherein, when the porous field is solved, a sequential solving mechanism of 'phase transition first and then cracking' is adopted, that is, in each iteration step, for a single grid cell, the phase transition condition is judged based on the fluid temperature and pressure, and after the phase transition process is completed and the liquid phase volume fraction is reduced to close to zero, the cracking reaction model is triggered to calculate the chemical source term; in the process of alternately iterating, the point-to-point transmission of physical quantities such as temperature, mass flow and component mass fraction between the high-speed outer flow field and the porous field is realized in the form of a boundary data file through the interface. The numerical divergence problem is solved, the cross-speed domain oscillation risk is eliminated, and high-fidelity global coupling simulation is realized.

Description

Technical Field

[0001] This invention relates to the field of thermal protection technology for aerospace vehicles, and in particular to a numerical calculation method for phase change sweating cooling of coupled pyrolysis reactions under supersonic conditions. Background Technology

[0002] During high-speed flight, critical components of an aircraft, such as the nose cone, air intake, wing leading edge, and combustion chamber walls, experience extremely high aerodynamic heat flux densities. Phase change evaporation cooling technology utilizes the large specific surface area of ​​porous structures and the characteristic of the cooling medium to absorb the latent heat of phase change, achieving efficient thermal protection with relatively low medium consumption. Liquid water is often used as a cooling medium due to its easy availability and large latent heat of phase change. However, pure water has limitations in engineering applications. For example, it is prone to freezing and blockage in low-temperature and low-pressure environments, and it can easily induce gas-liquid two-phase flow instability and temperature oscillations under high-speed flow conditions. Furthermore, when used in combustion chambers, liquid water may adversely affect combustion efficiency and combustion stability.

[0003] To address the aforementioned issues, existing research proposes using a mixed working fluid, specifically by incorporating crackable components such as propylene glycol or hydrocarbon fuels into water, to improve the cooling fluid's adaptability to complex operating conditions. Under high-temperature environments, the cooling fluid incorporating crackable components not only undergoes a liquid-gas phase transition but may also undergo further cracking reactions, absorbing a significant amount of heat during the reaction and thus providing additional "chemical heat sinking." However, existing numerical calculation methods typically only consider the liquid-gas phase transition process, primarily assessing the "physical heat sinking" effect generated by the latent heat of phase change, making it difficult to simulate the "chemical heat sinking" effect brought about by cracking reactions and its impact on overall cooling performance. Furthermore, after the coolant seeps out through the porous structure, it enters the external supersonic flow field and may continue to undergo cracking reactions, thereby affecting the wall heat flux distribution and overall cooling effect.

[0004] Therefore, accurate modeling of the sweating cooling process using a mixed working fluid requires establishing a strongly coupled multiphysics mathematical system. This system must not only achieve strong coupling between gas-liquid phase change and decomposition reaction within the porous medium, but also further couple it globally with the external supersonic airflow accompanying the coolant decomposition reaction. Currently, such strongly coupled calculations involving multiphase flow, chemical reactions, and high-speed compressible flow face significant numerical convergence difficulties within the existing computational fluid dynamics (CFD) framework, easily leading to computational divergence or solution distortion, and failing to meet the computational stability and accuracy requirements of engineering predictions. Summary of the Invention

[0005] This invention aims to propose a numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions, which can solve the problems of calculation divergence and accuracy distortion that easily occur when the gas-liquid phase change and pyrolysis reaction in porous structures are strongly coupled with the external high-speed compressible flow in the prior art.

[0006] This invention proposes a numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions, comprising the following steps: establishing a geometric model including a high-speed external flow domain and a porous domain, and defining the interface between the two; performing alternating iterative solutions for the high-speed external flow domain and the porous domain until the global calculation converges; wherein, when solving the porous domain, a sequential solution mechanism of "phase change first, then pyrolysis" is adopted, that is, in each iteration step, for a single grid cell, the phase change conditions are first determined based on the fluid temperature and pressure, and after the phase change process is completed and the liquid phase volume fraction drops to near zero, the pyrolysis reaction model is triggered to calculate the chemical source term; during the alternating iteration process, through the interface, point-to-point transfer of physical quantities such as temperature, mass flow rate and component mass fraction between the high-speed external flow domain and the porous domain is realized in the form of boundary data files.

[0007] Furthermore, the sequential solution mechanism, including the "phase transformation first, then pyrolysis" method, specifically involves extracting the true fluid temperature of the current mesh cell. With local pressure P; when Reaching the saturated boiling point of the working fluid under the current pressure P When the phase transition begins, the phase transition source term is activated, and the liquid volume fraction within the unit is monitored. When the liquid volume fraction is detected to drop to near 0, the local phase transition process is determined to be complete. Then, the pyrolysis reaction model is called to calculate and update the pyrolysis chemical source term.

[0008] Furthermore, the phase transition calculation within the porous domain is based on the modified Local Thermal Nonequilibrium Two-Phase Mixing Model (LTNE-TPMM).

[0009] Furthermore, when solving the high-speed external flow domain, a density-based solver is used to simultaneously couple the supersonic compressible flow, convective heat transfer, and the accompanying gas-phase pyrolysis reaction.

[0010] Furthermore, the boundary data file is a Profile file, which contains the coordinates of each grid node at the interface and the corresponding physical quantity information, used to achieve precise point-to-point transmission of data between different computational domains.

[0011] Furthermore, the convergence criterion for the alternating iterative solution is: in two adjacent iterations, the relative errors of coolant mass flow rate, temperature, and component mass fraction at each grid node of the interface are all less than 0.1‰.

[0012] Furthermore, the mixed working fluid used is an aqueous solution containing a pyrolyzable component, preferably an aqueous solution of propylene glycol.

[0013] The numerical calculation method for phase transition sweating cooling of coupled pyrolysis reaction under supersonic conditions of the present invention has the following advantages: The numerical divergence problem has been solved: by establishing a sequential solution mechanism of "phase change first, then pyrolysis", the computational divergence problem that is easily caused by the simultaneous strong coupling of multiple nonlinear source terms of gas-liquid phase change and pyrolysis reaction in porous media is effectively avoided, and the stability of numerical calculation is significantly improved.

[0014] Eliminates the risk of cross-velocity oscillation: By adopting a point-to-point data transfer strategy based on profile files, the physical quantities of the internal and external flow fields at the interface are accurately matched, eliminating the risk of numerical oscillation caused by the mismatch in data transfer across the velocity domain (low-speed flow inside the pore and supersonic flow outside).

[0015] Achieving high-fidelity global coupled simulation: Breaking through the stability limitations of traditional CFD methods when dealing with multiple nonlinear source terms, it is the first to realize global steady-state coupled calculation covering phase change in porous structures, chemical reaction kinetics, and external high-speed compressible flow. It can more accurately evaluate the combined effect of "physical heat sink" and "chemical heat sink" in sweating cooling systems, and provides a high-fidelity numerical simulation tool for efficient heat protection design. Attached Figure Description

[0016] The accompanying drawings are provided to further illustrate the invention and form part of the specification. They are used in conjunction with embodiments of the invention to explain the invention and do not constitute a limitation thereof. In the drawings: Figure 1 The diagram shows a flowchart illustrating the numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions in this invention.

[0017] Figure 2 The diagram shows a physical model used for calculations in this invention, which includes a high-speed external flow domain and a porous domain.

[0018] Figure 3 The diagram shows the results calculated for a specific working condition in this invention. Detailed Implementation

[0019] The technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments of the present invention. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments.

[0020] like Figures 1 to 3 As shown, this invention provides a numerical calculation method for phase change sweating cooling of coupled pyrolysis reactions under supersonic conditions, comprising the following steps: Step 1: Preprocessing and Initialization A geometric model comprising a high-speed external flow domain and a porous domain is established, and meshing and global flow field initialization are performed separately. The common boundary between the high-speed external flow domain and the porous domain is defined as the data interface.

[0021] Step 2: Multi-region iterative coupled solution An alternating iterative solution process involving high-speed external flow domain calculations and porous domain calculations is initiated. During the iteration process, the two computational domains exchange physical quantity information through the data interface until the global multiphysics coupling reaches a preset convergence criterion.

[0022] Step 3: Solving the porous domain Within the porous domain's internal cells, a sequential solution mechanism of "phase transformation first, then fragmentation" is implemented. Specifically, within each iteration step, the solution is first determined based on the actual fluid temperature within the current mesh cell. Determine whether the phase transition condition has been met by comparing the local pressure P; when Reaching the saturated boiling point of the working fluid under the current pressure When the phase transition model is activated to calculate the phase transition source term, the liquid phase volume fraction within the unit is monitored. When the liquid phase volume fraction drops to near zero, it is determined that the local phase transition process within the unit is complete, and then the pyrolysis reaction model is called to calculate the pyrolysis chemical source term. In a preferred embodiment, in step 3, the phase transition calculation within the porous domain is based on the modified Local Thermal Non-equilibrium Two-Phase Mixing Model (LTNE-TPMM).

[0023] Step 4: Solving the high-speed external watershed In the high-speed external flow domain, a density-based solver is used to simultaneously couple the solutions for supersonic compressible flow, convective heat transfer, and the accompanying gas-phase pyrolysis reaction.

[0024] Step 5: Interface data transfer and convergence determination After each single-step or multi-step iteration of a computational domain, an interface boundary data file for transfer is generated or updated. This data file contains the interface grid node coordinates, temperature, mass flow rate, and multi-component mass fraction information. After the high-speed outer flow domain calculation is completed, the wall pressure and local heat flux density distribution are transferred to the porous domain as boundary conditions. After the porous domain calculation is completed, the temperature, mass flow rate, and component mass fraction at the outlet boundary are transferred to the high-speed outer flow domain as boundary conditions. By comparing the interface physical quantity data of two adjacent iterations, when the error of the corresponding data for each node is less than a preset threshold, the calculation is considered converged, and the calculation results of the global component field, flow field, and temperature field are output. In a preferred embodiment, the data transfer in step 5 is specifically achieved by compiling a boundary data file (Profile) to realize the precise point-to-point transfer of mixture data between interface grid nodes.

[0025] The novel phase change-pyrolysis coupled sweating cooling mathematical model constructed in this invention can uniformly characterize the coupled control equations of the mixed working fluid phase change process and the pyrolysis reaction process, as shown in equation [equation missing]. As shown:

[0026]

[0027]

[0028] . In equation (1), ,in Density, unit: u represents velocity, unit: m / s; In equation (2), K is the permeability, in units of: ; The viscosity is the mixed kinematic viscosity, a function of the saturation s, in units of: ; Pressure, unit: Pa; In equation (3), The specific heat capacity of the coolant liquid phase at constant pressure, in J / (kg·K); Corrected temperature for fluid, unit: K; h is the fluid diffusion coefficient, in W / (m·K); h is the enthalpy, in J / kg; the subscript i indicates the type of component; J represents the mass diffusion flux, in J / kg. Subscript Indicates the types of components; This represents the heat transfer rate between the solid skeleton and the fluid per unit volume, in units of: ; Represents the source term of a chemical reaction, unit: ; In equation (4), The effective thermal conductivity of the porous framework is expressed in W / (m·K). Temperature of porous framework solids, unit: K; This represents the heat transfer rate between the solid skeleton and the fluid per unit volume, in units of: .

[0029] In the numerical solution of the above equations, this invention employs a "phase transformation first, then fragmentation" solution mechanism. Within each iteration step, the solver first extracts the true fluid temperature of the current mesh element. With local pressure P. When Reaching the saturated boiling point of the working fluid under the current pressure When the phase transition process begins, the phase transition source term is activated, and the liquid volume fraction (or liquid saturation s) within the unit is monitored. When s drops to near 0 (i.e., the fluid has completely transformed into the gas phase), the local phase transition process is considered complete, and the pyrolysis reaction model is then called to update the pyrolysis chemical source term. This mechanism avoids strong coupling of multiple nonlinear source terms, ensuring computational convergence. A density-based solver is used in the high-speed external flow domain to simultaneously and coupledly solve for supersonic compressible flow, convective heat transfer, and the accompanying gas-phase pyrolysis reaction, accurately capturing the complex interaction between the supersonic boundary layer and the coolant outflow. In the single-step iterative solution of the high-speed external flow, the computational domain interface directly reads the profile file transferred from the porous domain (containing temperature, mass flow rate, and component mass fraction from the porous domain outlet boundary). These boundary conditions serve as inputs to the component transport and energy equations of the high-speed external flow domain, ultimately outputting the pressure and local heat flux density distribution at the wall for use in the next iteration step of the porous domain.

[0030] The following is an illustration using a specific example: See Figure 1 The calculation flowchart shown and Figure 2 The physical model shown in this embodiment provides a numerical calculation method for phase change sweating cooling under supersonic conditions with coupled pyrolysis reactions. The physical model consists of two parts: a high-speed outer flow domain and a porous domain. The coolant flows into the porous domain, undergoes phase change and pyrolysis reactions sequentially during flow, and then seeps out of the porous structure into the high-speed outer flow domain, forming a gas film layer near the wall to prevent direct erosion of the thermal protection structure by high-temperature gases.

[0031] In this embodiment, the porous structure has an outer diameter of 5 mm, a thickness of 3 mm, and an opening angle of 8°. The computational domain is divided using a structured grid with a total of approximately 560,000 grid cells. The adjacent boundary between the high-speed external flow domain and the porous domain is defined as the data interface.

[0032] Taking propylene glycol aqueous solution as a mixed cooling medium as an example, the main cracking reaction pathway of propylene glycol was first calculated and established using the KiSThelP program based on transition state theory, and the corresponding Arrhenius kinetic parameters were obtained by fitting. The main cracking reaction pathway of propylene glycol is shown in the reaction equation. As shown, the specific parameters are shown in Table 1 of the instruction manual:

[0033]

[0034] .

[0035] Table 1 Pyrolysis reaction parameters

[0036] The specific calculation steps in this embodiment are as follows: (1) Initialization and start-up calculation: Before the solution begins, the temperature field, velocity field and component field in the high-speed external flow domain and the porous domain are initialized respectively. First, the high-speed external flow domain example is started, and the initial wall temperature, coolant flow rate and component mass fraction are preset at the interface.

[0037] (2) Initial data transfer from the outer basin to the porous basin: After the high-speed outer basin calculation has initially converged, the heat flow and pressure distribution at the interface are obtained and exported as a boundary data file in Profile format according to the node coordinates.

[0038] (3) Solving the porous domain: Start the porous domain example, import the Profile file containing pressure and heat flux distribution as boundary conditions, and set the porosity, particle size, pore size, propylene glycol aqueous solution concentration, and related thermophysical parameters of the porous domain. Within each iteration step of the porous domain, execute the core solution mechanism of this invention, namely the coupled method of "phase transformation first, then pyrolysis": Phase transition determination and calculation: The solver first extracts the true fluid temperature of the current mesh element. With local pressure P. When Reaching the saturated boiling point of the working fluid under the current pressure When the phase change process is determined to have started, the phase change source term based on the Modified Local Thermal Nonequilibrium Two-Phase Mixing Model (LTNE-TPMM) is activated, and the liquid volume fraction within the unit is continuously monitored.

[0039] Cracking reaction triggering: When the liquid volume fraction (or liquid saturation) in the unit is monitored to drop to near 0, indicating that the fluid has been completely converted into the gas phase, the local phase transition process is determined to be complete. At this time, the propylene glycol cracking reaction model is immediately invoked, the current gas phase state parameters are substituted, and the chemical composition source term and energy source term brought about by the cracking reaction are calculated and updated.

[0040] (4) Data transfer from the porous domain to the external flow domain: After the porous domain calculation is converged, the coordinates of each grid node at the interface and the corresponding coolant mass flow rate, temperature and mass fraction of the components after pyrolysis are extracted using post-processing code, and then compiled into a Profile boundary data file.

[0041] (5) Iterative convergence determination: Compare the interface data output by the porous domain in this iteration with the data of the corresponding nodes in the previous iteration. When the relative errors of the coolant mass flow rate, temperature and component mass fraction of all nodes are less than 0.1‰, the global calculation is determined to be converged and the process jumps to step (6); otherwise, the updated Profile file is passed to the high-speed external flow domain as the inlet boundary condition for its next iteration, and the process returns to step (2) to continue the alternating iterative calculation.

[0042] (6) Output results: After the calculation converges, the final global component field, flow field and temperature field distribution results are output.

[0043] To verify the effectiveness of this method, calculations were performed using a case study of an altitude of 30 km, an incoming flow Mach number of 6, and a propylene glycol mass flux of 0.6 kg / m² / s. The calculated temperature contour map, component distribution, and other results are shown below. Figure 3 As shown, the calculation process is stable, without divergence or non-physical oscillations, demonstrating the superiority and reliability of the method of this invention in handling such multi-nonlinear coupling problems.

[0044] Therefore, the numerical calculation method for phase change sweating cooling under supersonic conditions coupled with pyrolysis provided by this invention solves the problem of numerical divergence: by establishing a sequential solution mechanism of "phase change first, then pyrolysis", it effectively avoids the calculation divergence problem that is easily caused by the simultaneous strong coupling solution of multiple nonlinear source terms of gas-liquid phase change and pyrolysis reaction in porous media, and significantly improves the stability of numerical calculation. It eliminates the risk of cross-velocity domain oscillation: by adopting a point-to-point data transfer strategy based on profile files, it achieves accurate matching of physical quantities at the interface between the internal and external flow fields, eliminating the risk of numerical oscillation caused by the mismatch of data transfer across velocity domains (low-speed flow inside the porous structure and external supersonic flow). It achieves high-fidelity full-domain coupled simulation: breaking through the stability limitations of traditional CFD methods in handling multiple nonlinear source terms, it is the first to achieve full-domain steady-state coupled calculation covering phase change inside the porous structure, chemical reaction kinetics, and external high-speed compressible flow. This allows for a more accurate evaluation of the combined effects of "physical heat sink" and "chemical heat sink" in the sweating cooling system, providing a high-fidelity numerical simulation tool for efficient heat protection design.

[0045] The above description is only a preferred embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any equivalent substitutions or modifications made by those skilled in the art within the scope of the technology disclosed in the present invention, based on the technical solution and inventive concept of the present invention, should be covered within the scope of protection of the present invention.

Claims

1. A numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions, characterized in that, Including the following steps: Establish a geometric model that includes a high-speed external flow domain and a porous domain, and define the interface between the two; Perform alternating iterative solutions for the high-speed external flow domain and the porous domain until the global computation converges; In solving the porous domain, a sequential solution mechanism of "phase transformation first, then pyrolysis" is adopted. That is, in each iteration step, for a single grid cell, the phase transformation conditions are first determined based on the fluid temperature and pressure. After the phase transformation process is completed and the liquid volume fraction drops to near zero, the pyrolysis reaction model is triggered to calculate the chemical source term. During the alternating iteration process, the interface enables point-to-point transfer of physical quantities such as temperature, mass flow rate, and component mass fraction between the high-speed external flow domain and the porous domain in the form of boundary data files.

2. The numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions according to claim 1, characterized in that, The sequential solution mechanism, including the "phase transition first, then pyrolysis" method, is as follows: Extract the true fluid temperature of the current mesh cell With local pressure P; when Reaching the saturated boiling point of the working fluid under the current pressure P When the phase transition begins, the phase transition source term is activated, and the liquid volume fraction within the unit is monitored. When the liquid phase volume fraction is detected to drop to near 0, the local phase transition process is determined to be complete. Then, the pyrolysis reaction model is called to calculate and update the pyrolysis chemical source term.

3. The numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions according to claim 1, characterized in that, The phase transition calculations within the porous domain are based on a modified Local Thermal Nonequilibrium Two-Phase Mixing Model (LTNE-TPMM).

4. The numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions according to claim 1, characterized in that, When solving the high-speed external flow domain, a density-based solver is used to simultaneously couple the supersonic compressible flow, convective heat transfer, and the accompanying gas-phase pyrolysis reaction.

5. The numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions according to claim 1, characterized in that, The boundary data file is a Profile file, which contains the coordinates of each grid node at the interface and the corresponding physical quantity information, and is used to realize the precise point-to-point transfer of data between different computational domains.

6. The numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions according to claim 1, characterized in that, The convergence criterion for the alternating iterative solution is: in two adjacent iterations, the relative errors of coolant mass flow rate, temperature, and component mass fraction at each grid node of the interface are all less than 0.1‰.

7. The numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions according to claim 1, characterized in that, The mixed working fluid used is an aqueous solution containing pyrolyzable components.

8. The numerical calculation method for phase change sweating cooling of coupled pyrolysis reaction under supersonic conditions according to claim 7, characterized in that, The aqueous solution is a propylene glycol aqueous solution.

Images