Performance evaluation method of solid hypersonic ramjet based on one-time combustion effect
Patent Information
- Application Number
- CN202611222033.X
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2026-08-12
- Publication Date
- 2026-10-13
AI Technical Summary
对于包含燃气发生器一次燃烧、超声速燃烧室二次燃烧两级燃烧过程的固体超燃冲压动力装置而言,二次燃烧的燃料入口状态由一次燃烧的产物参数直接决定,无法脱离一次燃烧过程独立设定;沿用现有评估方法直接预设燃料入口参数的处理方式,与该类动力装置的实际工作过程不匹配,会导致性能评估结果与实际工况存在偏差,无法准确反映动力装置的真实工作特性
[0016]本申请提供一种基于一次燃烧效应的固体超燃冲压动力装置性能评估方法,首先,通过压缩过程模型计算经进气道压缩后的超声速气流状态参数,准确得到第二反应区的主流来流边界;然后,基于热化学平衡原理求解第一反应区的平衡温度和产物组分,获得第一反应区燃烧后的真实产物状态,并以此作为第二反应区的燃料入口边界,避免直接预设燃料参数与实际工况不符的问题,客观体现第一反应区燃烧过程对后续燃烧过程的影响;接着,将上述两类参数共同作为第二反应区的入口边界条件,通过准一维控制方程组计算第二反应区内气流状态参数沿轴向的分布,能够还原两级燃烧作用下的燃烧室工作过程;再以第二反应区出口的气流状态参数为基础,通过膨胀过程模型求解尾喷管内气流状态的轴向分布,完成燃气膨胀加速过程的模拟;最后,根据尾喷管出口截面的气流状态参数计算动力装置的整体性能参数。由此,本申请通过将一次燃烧过程的求解结果作为二次燃烧的入口边界进行计算,能够适配固体超燃冲压动力装置的两级燃烧工作特性,使性能评估过程和动力装置实际工作过程的契合度更高,有效提高性能评估结果的准确性和可信度。
Smart Images

Figure CN122904794A_ABST
Abstract
Description
Technical Field
[0001] This application relates to the field of power unit testing technology, specifically to a performance evaluation method for solid scramjet power units based on the primary combustion effect. Background Technology
[0002] Scramjet propulsion systems are the core power systems of air-breathing hypersonic vehicles. Based on the form of fuel used, they can be divided into liquid-fuel scramjet propulsion systems and solid-fuel scramjet propulsion systems. Among them, solid scramjet propulsion systems employing a gas generator structure have broad engineering application prospects in the field of hypersonic propulsion technology due to their advantages such as simple structure, high volumetric impulse, long storage period, flexible gas flow adjustment, and good flame stability. In the preliminary design stage of the propulsion system, using quasi-one-dimensional analysis methods to conduct performance evaluation is a core technical means to verify the feasibility of the scheme and support rapid iteration of design parameters.
[0003] Currently, most existing quasi-one-dimensional performance evaluation methods are designed for liquid-fueled scramjet power plants. The combustion process in these plants occurs only within the supersonic combustion chamber. Therefore, the corresponding evaluation methods typically pre-define the initial state parameters and thermodynamic parameters of the fuel injected into the combustion chamber as known boundary conditions. Based on these, the distribution of airflow parameters along the combustion chamber is calculated, ultimately determining the overall performance of the power plant. However, for solid scramjet power plants, which involve a two-stage combustion process—primary combustion in the gas generator and secondary combustion in the supersonic combustion chamber—the fuel inlet state for secondary combustion is directly determined by the product parameters of primary combustion and cannot be set independently of the primary combustion process. The existing evaluation methods, which directly pre-define fuel inlet parameters, are incompatible with the actual operating process of these power plants, leading to deviations between performance evaluation results and actual operating conditions, and failing to accurately reflect the true operating characteristics of the power plant. In summary, existing quasi-one-dimensional performance evaluation methods cannot adapt to the two-stage combustion characteristics of solid scramjet power plants and are insufficient to meet the performance evaluation requirements during the design phase of such power plants.
[0004] The preceding description is intended to provide general background information and does not necessarily constitute prior art. Summary of the Invention
[0005] This application provides a performance evaluation method for solid scramjet power units based on the primary combustion effect, which can realize performance evaluation adapted to the two-stage combustion operating characteristics of solid scramjet power units and improve the accuracy of performance evaluation results.
[0006] Based on this, embodiments of this application provide a method for evaluating the performance of a solid scramjet power unit based on the primary combustion effect, including: Obtain flight operating parameters of the power unit and calculate supersonic airflow state parameters after compression through the inlet based on the compression process model; The initial composition of reactants and reaction pressure in the first reaction zone are obtained, and the equilibrium temperature and product composition of the first reaction zone are calculated based on the principle of thermochemical equilibrium. Using the supersonic airflow state parameters, the equilibrium temperature of the first reaction zone, and the product composition as the inlet boundary conditions of the second reaction zone, the distribution of the airflow state parameters along the axial direction in the second reaction zone is calculated based on the quasi-one-dimensional control equations. Using the airflow state parameters at the outlet of the second reaction zone as the inlet parameters of the expansion process, the distribution of the airflow state parameters along the axial direction in the tail nozzle is calculated based on the expansion process model, and the airflow state parameters at the outlet section of the tail nozzle are obtained. The overall performance parameters of the power unit are calculated based on the airflow state parameters of the tailpipe outlet section.
[0007] Furthermore, in some embodiments of this application, the step of obtaining the flight operating parameters of the power unit and calculating the supersonic airflow state parameters after compression through the inlet based on the compression process model includes: Based on the flight altitude of the power unit, standard atmospheric parameters are queried to calculate the static temperature, static pressure, atmospheric density, and air mass flow rate of the free flow. The intake outlet Mach number ratio and compression efficiency are preset, and the total temperature, static temperature, static pressure, Mach number and velocity of the airflow at the intake outlet are solved based on the adiabatic compression assumption and the total enthalpy conservation relationship. The loss of airflow state parameters during the compression process is calculated based on the preset compression efficiency, and the airflow state parameters of the inlet outlet section are output.
[0008] Furthermore, in some embodiments of this application, obtaining the initial components of the reactants and the reaction pressure of the first reaction zone, and calculating the equilibrium temperature and product components of the first reaction zone based on the principle of thermochemical equilibrium, includes: Obtain the initial components of the reactants in the first reaction zone and count the total atomic count of each chemical element in the reaction system; Gas-phase active species and candidate condensed phase species are selected to construct an active species set, and the thermodynamic property parameters of each species are calculated at the current temperature; Based on the conservation of total elements, chemical potential balance, conservation of total moles in the gas phase, and enthalpy balance, a set of governing equations was constructed, and the equilibrium temperature and the mole number of each species were solved iteratively using the descending Newton-Raphson method. The Gibbs stability of condensed phases not included in the set of active species is determined. If the inclusion condition is met, the condensed phase is included in the set of active species, and the iterative solution process is re-executed until all condensed phases outside the set of active species do not meet the inclusion condition. The equilibrium temperature, molar number of each species, and thermodynamic properties of each component obtained after iterative convergence are used as the output results of the first reaction zone.
[0009] Furthermore, in some embodiments of this application, the step of performing Gibbs stability determination on condensed phases not included in the set of active species, and if the inclusion condition is met, including the condensed phase in the set of active species, and re-executing the iterative solution process until all condensed phases outside the set of active species no longer meet the inclusion condition, includes: For each condensed phase species outside the set of active species, calculate the change in Gibbs free energy before and after the condensed phase species is added to the reaction system. If the change is less than zero, the condensed phase species is determined to meet the addition condition, the condensed phase species is added to the set of active species, and the iterative solution is performed again. If the number of moles of condensed phase species that have been included in the set of active species is negative, then the condensed phase species is removed from the set of active species and the iterative solution is performed again.
[0010] Furthermore, in some embodiments of this application, the step of using the supersonic airflow state parameters, the equilibrium temperature of the first reaction zone, and the product composition as the inlet boundary conditions of the second reaction zone, and calculating the axial distribution of the airflow state parameters in the second reaction zone based on the quasi-one-dimensional governing equations, includes: Obtain the geometric configuration parameters of the second reaction zone, the injection position and injection pattern of the reactants, and the equivalence ratio; set the integration start position and integration step size. At the current integration position, the mixing efficiency between the reactants and the mainstream air is calculated based on the mixing efficiency model to determine the mass flow rate of the reactants participating in the reaction and the heat released by the reaction. By combining the mass conservation equation, momentum conservation equation, energy conservation equation, and gas state equation, and combining the relationship between Mach number and speed of sound, a set of first-order ordinary differential equations concerning Mach number, pressure, and temperature is derived. The first-order ordinary differential equations were solved by axial integration using the fourth-order Runge-Kutta method to obtain the airflow Mach number, total pressure, static pressure, total temperature, static temperature and velocity at each integration position; Determine whether the current integration position has reached the outlet of the second reaction zone; if the current integration position has not reached the outlet of the second reaction zone, advance the integration position by one step, recalculate the mixing efficiency and heat release at the current integration position, and continue to perform the integration solution; if the current integration position has reached the outlet of the second reaction zone, output the airflow state parameters at the outlet of the second reaction zone; if the Mach number drops to 1 and supersonic flow cannot be maintained during the integration process, thermal blockage is determined to have occurred and the calculation is terminated.
[0011] Furthermore, in some embodiments of this application, the mass conservation equation includes an axial change term in mass flow rate caused by the addition of reactants; the momentum conservation equation also includes a cross-sectional area change term, a wall friction term, and a momentum change term caused by the addition of reactants; the energy conservation equation also includes an enthalpy change term caused by the addition of reactants, a heat change term caused by the heat release of the reaction, and a heat loss term caused by wall convection. The set of first-order ordinary differential equations includes the first-order differential equations of Mach number along the axial direction, pressure along the axial direction, and temperature along the axial direction. The first-order differential equations of Mach number along the axial direction, pressure along the axial direction, and temperature along the axial direction are solved by coupling the Mach number, pressure, temperature, rate of change of cross-sectional area, reactant addition rate, and reaction heat release rate at the current integration position.
[0012] Furthermore, in some embodiments of this application, the wall friction term in the momentum conservation equation is determined based on the wall friction coefficient, and the calculation method of the wall friction coefficient includes: The Reynolds number is temperature-corrected based on the dynamic viscosity at the mainstream airflow temperature and the dynamic viscosity at the wall temperature, and the corrected Reynolds number at the current integration position is calculated. The wall friction coefficient is calculated based on the modified Reynolds number and a semi-empirical formula.
[0013] Furthermore, in some embodiments of this application, the wall heat loss term in the energy conservation equation is determined based on the wall heat flux density, and the calculation method of the wall heat flux density includes: Based on the reference enthalpy method, the reference enthalpy is calculated according to the enthalpy of the main airflow, the enthalpy of the adiabatic wall, and the enthalpy at the wall temperature, and the reference temperature is determined according to the reference enthalpy; The convective heat flux density between the airflow and the wall is calculated based on the thermodynamic parameters and convective heat transfer coefficient at the reference temperature.
[0014] Furthermore, in some embodiments of this application, the step of using the airflow state parameters at the outlet of the second reaction zone as the inlet parameters of the expansion process, calculating the axial distribution of the airflow state parameters within the nozzle based on the expansion process model, and obtaining the airflow state parameters at the nozzle outlet cross-section includes: Ignoring the addition of reactants and the heat released by the reaction, and only considering the acceleration effect of cross-sectional area change and wall friction on the airflow, the mass conservation equation, momentum conservation equation, energy conservation equation and gas state equation are combined and derived into a set of first-order ordinary differential equations concerning Mach number, pressure and temperature. Using the airflow state parameters at the outlet of the second reaction zone as the initial values for integration, the fourth-order Runge-Kutta method is used to integrate along the axial direction to obtain the airflow Mach number, total pressure, static pressure, total temperature, static temperature and velocity at each integration position in the tail nozzle.
[0015] Furthermore, in some embodiments of this application, the method further includes: During multiple executions, the reactant composition or reaction pressure in the first reaction zone is changed; Based on the modified reactant composition and reaction pressure, the equilibrium temperature and product composition of the first reaction zone are recalculated based on the thermochemical equilibrium principle, and the product composition and equilibrium temperature of the first reaction zone and the overall performance parameters of the power unit are obtained for each execution. Based on the differences between the overall performance parameters obtained in each execution, the correspondence between the changes in reactant composition or reaction pressure in the first reaction zone and the overall performance of the power unit is determined.
[0016] This application provides a performance evaluation method for a solid scramjet power unit based on the primary combustion effect. First, the supersonic airflow state parameters after compression through the intake duct are calculated using a compression process model to accurately obtain the mainstream incoming flow boundary of the second reaction zone. Then, based on the thermochemical equilibrium principle, the equilibrium temperature and product composition of the first reaction zone are solved to obtain the actual product state after combustion in the first reaction zone, which is used as the fuel inlet boundary of the second reaction zone. This avoids the problem of directly pre-setting fuel parameters not matching actual operating conditions and objectively reflects the influence of the combustion process in the first reaction zone on the subsequent combustion process. Next, the above two types of parameters are used together as the inlet boundary conditions of the second reaction zone. The axial distribution of the airflow state parameters in the second reaction zone is calculated using a quasi-one-dimensional control equation set, which can recreate the working process of the combustion chamber under two-stage combustion. Then, based on the airflow state parameters at the outlet of the second reaction zone, the axial distribution of the airflow state in the tailpipe is solved using an expansion process model to complete the simulation of the gas expansion acceleration process. Finally, the overall performance parameters of the power unit are calculated based on the airflow state parameters at the outlet section of the tailpipe. Therefore, by using the solution of the primary combustion process as the inlet boundary for secondary combustion, this application can adapt to the two-stage combustion characteristics of solid scramjet power units, making the performance evaluation process more consistent with the actual working process of the power unit, and effectively improving the accuracy and reliability of the performance evaluation results. Attached Figure Description
[0017] To more clearly illustrate the technical solutions in the embodiments of this application, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments of this application. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.
[0018] Figure 1 This is a flowchart illustrating the performance evaluation method for a solid scramjet power unit based on the primary combustion effect provided in this application embodiment; Figure 2 This is a schematic diagram of the solid scramjet power unit provided in the embodiments of this application; Figure 3 The adiabatic compression efficiency η provided in the embodiments of this application c A schematic diagram showing the relationship between the incoming Mach number, the circulating static temperature ratio, and the number of oblique shock waves; Figure 4 This is a schematic diagram of the combustion chamber control unit of the solid scramjet power unit provided in the embodiments of this application; Figure 5 This is a schematic flowchart illustrating the solution of a single combustion process provided in an embodiment of this application; Figure 6 This is a schematic diagram of the quasi-one-dimensional analysis and calculation process provided in the embodiments of this application. Detailed Implementation
[0019] Exemplary embodiments will now be described in detail, examples of which are illustrated in the accompanying drawings. When the following description relates to the drawings, unless otherwise indicated, the same numerals in different drawings denote the same or similar elements. The embodiments described in the following exemplary embodiments do not represent all embodiments consistent with this application. Rather, they are merely examples of systems and methods consistent with those detailed in the appended claims or with some aspects of this application.
[0020] It should be noted that, in this embodiment, the terms "comprising," "including," or any other variations thereof are intended to cover descriptions such as non-exclusive inclusion, thereby ensuring that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Without further limitations, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes that element. Furthermore, components, features, and elements with the same names in different embodiments of this application may have the same meaning or different meanings, and their specific meanings need to be determined by their interpretation in that specific embodiment or by further considering the context of that specific embodiment.
[0021] It should be understood that the specific embodiments described herein are merely illustrative of this application and are not intended to limit this application.
[0022] In the following description, the use of suffixes such as "module," "part," or "unit" to denote elements is solely for the purpose of illustrative purposes and has no specific meaning in itself. Therefore, "module," "part," or "unit" may be used interchangeably.
[0023] To address the aforementioned technical problems and overcome the shortcomings of existing technologies, this application provides a performance evaluation method for solid scramjet power units based on the primary combustion effect. This method enables performance evaluation that adapts to the two-stage combustion characteristics of solid scramjet power units, thereby improving the accuracy of performance evaluation results.
[0024] Please see Figure 1 , Figure 1 This is a flowchart illustrating a method for evaluating the performance of a solid scramjet power unit based on the primary combustion effect, according to an embodiment of this application. Specifically, the method may include the following steps: S1. Obtain flight operating parameters of the power unit and calculate the supersonic airflow state parameters after compression through the air intake based on the compression process model; Specifically, for step S1, the flight operating parameters of the power unit are first obtained, including core operating conditions that determine the incoming flow state, such as flight altitude and Mach number, which are used as the calculation input for the inlet compression process. Based on the compression process model, the gas dynamics of the deceleration and pressurization process of the supersonic free flow passing through the inlet are solved, simulating the state changes of the airflow after compression in the inlet, quantifying the evolution of the thermodynamic and flow parameters of the airflow, and finally outputting the supersonic airflow state parameters after compression in the inlet, including key indicators such as airflow temperature, pressure, Mach number, and velocity, as the inlet boundary conditions for the mainstream incoming flow in the subsequent combustion chamber.
[0025] S2. Obtain the initial composition of reactants and reaction pressure in the first reaction zone, and calculate the equilibrium temperature and product composition of the first reaction zone based on the principle of thermochemical equilibrium. Specifically, for step S2, the initial composition of the reactants and the reaction pressure of the first reaction zone are obtained. The initial composition of the reactants refers to the chemical elemental composition and component ratio of the solid viscoelastic material participating in the first combustion, and the reaction pressure is the internal pressure condition of the first reaction zone during steady-state operation. Based on the principle of thermochemical equilibrium, the combustion process in the first reaction zone is regarded as an isobaric adiabatic combustion process. The thermodynamic equilibrium state of the reaction system is solved, the equilibrium temperature corresponding to the chemical equilibrium is calculated, and the component proportion and molar content of various combustion products under the equilibrium state are calculated. The complete thermodynamic state and component composition of the products after the first combustion are output.
[0026] S3. Using the supersonic gas flow state parameters, the equilibrium temperature of the first reaction zone, and the product composition as the inlet boundary conditions of the second reaction zone, the distribution of the gas flow state parameters along the axial direction in the second reaction zone is calculated based on the quasi-one-dimensional control equations. Specifically, in step S3, the obtained supersonic airflow state parameters after intake compression, the equilibrium temperature of the first reaction zone, and the product composition are used together as the inlet boundary conditions of the second reaction zone; where the supersonic airflow parameters correspond to the mainstream air inlet of the combustion chamber, and the product parameters of the first reaction zone correspond to the fuel injection inlet. Based on the quasi-one-dimensional governing equations, the supersonic flow and co-combustion process in the second reaction zone is solved. The calculation is performed segment by segment along the axial direction of the airflow to obtain the Mach number, temperature, pressure, and other state parameters of the airflow at each axial position, ultimately forming a complete axial distribution law of the airflow state parameters in the second reaction zone. Using the primary combustion products as the fuel input for secondary combustion can objectively reflect the influence of the state of the primary combustion products on the combustion and flow process in the second reaction zone, and restore the actual working process of the combustion chamber under the two-stage combustion series action.
[0027] S4. Using the airflow state parameters at the outlet of the second reaction zone as the inlet parameters of the expansion process, calculate the distribution of the airflow state parameters along the axial direction in the tail nozzle based on the expansion process model, and obtain the airflow state parameters at the outlet section of the tail nozzle. Specifically, for step S4, the airflow state parameters at the outlet section of the second reaction zone are used as the initial inlet conditions for the expansion process. Based on the expansion process model, the expansion and acceleration process of high-temperature and high-pressure gas in the expansion channel of the tail nozzle is simulated. The changes in thermodynamic and flow parameters of the airflow during the expansion and work process are solved. The airflow state parameters at each position are calculated segment by segment along the tail nozzle axis. Finally, the complete airflow state parameters of the tail nozzle outlet section are obtained, including core indicators such as outlet airflow velocity, static temperature, and static pressure, providing a data foundation for subsequent overall performance calculations.
[0028] S5. Calculate the overall performance parameters of the power unit based on the airflow state parameters of the tail nozzle exit section; Specifically, for step S5, based on the airflow state parameters of the tail nozzle exit section, combined with parameters such as the working flow rate and incoming flow conditions of the power unit, the core performance indicators such as the thrust and specific impulse of the power unit are quantitatively solved using the power performance calculation method, and finally the overall performance parameters of the power unit are obtained, thus completing a complete evaluation of the working performance of the solid scramjet power unit.
[0029] This embodiment calculates the inlet boundary of the secondary combustion by using the solution of the primary combustion process as the inlet boundary of the secondary combustion. This is suitable for the two-stage combustion working characteristics of the solid scramjet power unit, making the performance evaluation process more consistent with the actual working process of the power unit, thereby improving the accuracy and reliability of the performance evaluation results.
[0030] Furthermore, in some embodiments, step S1, "obtaining the flight operating parameters of the power unit and calculating the supersonic airflow state parameters after compression through the inlet based on the compression process model," may specifically include: S11. Based on the flight altitude of the power unit, query the standard atmospheric parameters and calculate the static temperature, static pressure, atmospheric density, and air mass flow rate of the free flow; Specifically, for step S11, the flight altitude of the power unit is obtained. Using the flight altitude as input, a standard atmospheric parameter database is queried to obtain basic atmospheric physical properties such as atmospheric static temperature, atmospheric static pressure, and atmospheric density at the corresponding altitude. Combined with the flight Mach number of the power unit, and based on ideal gas thermodynamic relations and aerodynamic function formulas, the flow and thermodynamic parameters such as airflow velocity, total temperature, and total pressure of the free flow are calculated. Then, combined with the airflow capture cross-sectional area of the inlet, the mass flow rate of the air entering the inlet is calculated, thus fully constructing the inlet boundary parameters of the inlet compression process.
[0031] S12. Preset the intake outlet Mach number ratio and compression efficiency, and solve the total temperature, static temperature, static pressure, Mach number and velocity of the airflow at the intake outlet based on the adiabatic compression assumption and the total enthalpy conservation relationship. Specifically, for step S12, the inlet outlet Mach number ratio and adiabatic compression efficiency are preset, where the Mach number ratio is the ratio of the inlet outlet Mach number to the incoming flow Mach number, used to constrain the Mach number level of the compressed airflow. Based on the adiabatic compression assumption, the airflow does not exchange heat with the outside environment during the inlet compression process, satisfying the total enthalpy conservation relationship; therefore, the total temperature of the compressed airflow remains consistent with the total temperature of the incoming flow. Combining the preset outlet Mach number, the outlet static temperature is solved through the aerodynamic function relationship between Mach number and total temperature and static temperature; then, based on the pressure evolution law of the adiabatic compression process and combined with the preset compression efficiency, the static pressure, Mach number, and airflow velocity at the inlet outlet are derived, completing the solution of the core parameters of the compression process.
[0032] S13. Solve for the loss of airflow state parameters during the compression process based on the preset compression efficiency, and output the airflow state parameters of the inlet outlet section; preferably, the compression efficiency ranges from 0.85 to 0.95.
[0033] Specifically, for step S13, based on the preset compression efficiency, the loss of airflow state parameters caused by factors such as shock wave action and boundary layer development during the compression process is quantitatively calculated. The calculation results of ideal adiabatic compression are corrected to obtain and output the complete airflow state parameters of the inlet outlet section under actual operating conditions. The lower the compression efficiency, the greater the aerodynamic loss in the compression process, and the greater the deviation of parameters such as the total pressure and static pressure of the outlet airflow from the ideal adiabatic compression results.
[0034] Based on standard atmospheric parameters and combined with the adiabatic compression assumption and compression efficiency loss correction, this embodiment can quickly and accurately solve for the supersonic airflow state parameters at the inlet outlet with low computational load, providing accurate and reliable mainstream inlet boundary conditions for subsequent processes.
[0035] Furthermore, in some embodiments, step S2, "obtaining the initial components of the reactants and the reaction pressure in the first reaction zone, and calculating the equilibrium temperature and product components of the first reaction zone based on the principle of thermochemical equilibrium," may specifically include: S21. Obtain the initial components of the reactants in the first reaction zone and count the total number of atoms of each chemical element in the reaction system; Specifically, for step S21, the initial components of the reactants in the first reaction zone are obtained, and the chemical composition and proportion of various reactants participating in the first combustion are clarified. Based on the molecular formula and proportion of each reactant, all types of chemical elements contained in the entire reaction system are counted, and the total molar amount of each element is calculated, which serves as the element conservation benchmark for the combustion reaction.
[0036] S22. Select gas-phase active species and candidate condensed phase species, construct an active species set, and calculate the thermodynamic property parameters of each species at the current temperature; Specifically, for step S22, the main gaseous active species present in the reaction process and the candidate condensed phase species that may be generated are selected to construct an initial set of active species; based on thermodynamic property data, the standard molar enthalpy, isobaric molar heat capacity, chemical potential and other thermodynamic property parameters of each species in the set at the current temperature are calculated to provide a basis for property calculation for subsequent equilibrium state solution.
[0037] S23. Based on the conservation of total elements, chemical potential balance, conservation of total moles in the gas phase, and enthalpy balance, a set of governing equations is constructed. The equilibrium temperature and the number of moles of each species are solved iteratively using the descending Newton-Raphson method. The set of governing equations includes the following constraint equations: total element conservation constraint, used to constrain the total amount of each chemical element before and after the reaction; chemical potential balance constraint, used to constrain the balance relationship between the chemical potential of the gas phase products and the Lagrange multipliers of the elements; total moles in the gas phase conservation constraint; and enthalpy balance constraint, used to constrain the total enthalpy conservation before and after the reaction. Specifically, for step S23, a complete set of governing equations is constructed based on four constraints: total element conservation, chemical potential balance, total mole number conservation in the gas phase, and enthalpy balance. The total element conservation constraint ensures the conservation of the total atomic quantity of each chemical element before and after the reaction; the chemical potential balance constraint describes the equilibrium relationship between the chemical potential of the gas-phase products and the corresponding Lagrange multipliers of the elements; the total mole number conservation constraint limits the total amount of substance in the gas-phase system; and the enthalpy balance constraint ensures the conservation of the total enthalpy of the reaction system during adiabatic combustion. The descent-type Newton-Raphson method is used to iteratively solve this set of governing equations, gradually correcting unknown quantities such as the equilibrium temperature and the mole number of each species until the calculation results meet the convergence accuracy requirements, thus obtaining the equilibrium temperature and the mole number of each species under the current set of active species.
[0038] S24. Perform Gibbs stability determination on condensed phases not included in the active species set. If the inclusion condition is met, include the condensed phase in the active species set and re-execute the iterative solution process until all condensed phases outside the active species set do not meet the inclusion condition. Specifically, for step S24, for all candidate condensed phase species not included in the current active species set, Gibbs stability is determined one by one, and the change in Gibbs free energy of the entire system before and after the addition of the condensed phase species to the reaction system is calculated. If the Gibbs free energy of the system decreases after the addition of the species, i.e., the addition condition is met, the condensed phase species is added to the active species set, and the iterative solution is re-executed based on the updated active species set. The above determination and iteration process is repeated until all candidate condensed phase species outside the active species set no longer meet the addition condition, confirming that the current set has covered all stable reaction products.
[0039] S25. The equilibrium temperature, molar number of each species, and thermodynamic properties of each component obtained after iterative convergence are used as the output results of the first reaction zone; Specifically, for step S25, after the iteration has fully converged and all candidate condensed phases have completed the stability determination, the final equilibrium temperature, the number of moles of each species, and the thermodynamic property parameters of each component at the corresponding equilibrium temperature are uniformly used as the output results of the first reaction zone to provide boundary conditions for subsequent combustion process calculations.
[0040] This embodiment constructs a multi-constraint control equation set based on the thermochemical equilibrium principle, and dynamically adjusts the product species set by combining the Gibbs stability judgment of condensed phase. It can accurately solve the temperature and product composition when the combustion reaches equilibrium in a single combustion. The calculation results closely match the actual composition of combustion products, providing accurate fuel inlet boundary parameters for subsequent combustion processes.
[0041] Furthermore, in some embodiments, when using the descent-type Newton-Raphson method for iterative solution, a descent control factor is used to limit the correction magnitude of each unknown quantity. The descent control factor is determined based on the temperature correction, the total number of moles in the gas phase correction, and the number of moles of each gas phase species correction, and the determination method of the descent control factor is updated synchronously with each correction. Specifically, this may include: The first descent control factor is calculated based on the temperature correction, the total gas phase mole count correction, and the mole count correction of the main gas phase species. The second decrease control factor is calculated based on the molar number correction of trace gaseous species; The minimum value among 1, the first descent control factor, and the second descent control factor is used as the actual descent control factor.
[0042] Specifically, in each iteration of the descent Newton-Raphson method for solving the thermochemical equilibrium problem, initial corrections for all unknowns are first obtained by solving the governing equations. From these initial corrections, temperature corrections, total gas phase mole count corrections, and corrections for the mole counts of major gas phase species with high molar percentages are extracted. Based on these three types of corrections, a first descent control factor is calculated. This control factor targets the core variables that significantly affect the thermodynamic state of the system, limiting the correction magnitude for each iteration and preventing excessively large step sizes in the core parameters from causing iterative oscillations and deviations from reasonable physical ranges.
[0043] The molar number correction values for all trace gaseous species are extracted. Based on the characteristics of trace species having small molar numbers and being prone to non-physical negative values, a second descent control factor is calculated separately. This control factor is specifically adapted to the correction process for low-content gaseous components, avoiding physically inconsistent results where the molar number of trace species is less than zero after iterative updates, thus ensuring the physical validity of the entire iteration process.
[0044] The three values—value 1, the first descent control factor, and the second descent control factor—are compared, and the minimum value is selected as the descent control factor actually used in this iteration. The correction amount for all unknowns is multiplied by this actual control factor before variable updates are performed. After obtaining a new correction amount in each iteration, both types of control factors are recalculated simultaneously, and their actual values are determined, achieving dynamic matching between the correction magnitude and the iteration state. Value 1 corresponds to the uncontrolled state where the initial correction amount is fully executed; the rule of taking the minimum value ensures that all types of unknowns satisfy their respective correction magnitude constraints.
[0045] This embodiment calculates the descent control factor by distinguishing between major and minor variables and using the most stringent minimum value as the actual control value. This allows for dynamic adaptation of the correction state in each round during the iterative solution process, effectively limiting the single correction amplitude of each unknown quantity. This avoids iterative divergence caused by large fluctuations in core parameters and prevents non-physical negative values of trace components, thereby improving the stability and convergence reliability of the thermochemical equilibrium iterative solution.
[0046] Furthermore, in some embodiments, step S24, "to determine the Gibbs stability of condensed phases not included in the active species set, and if the inclusion condition is met, to include the condensed phase in the active species set and to re-execute the iterative solution process until all condensed phases outside the active species set no longer meet the inclusion condition," may specifically include: S241. For each condensed phase species outside the active species set, calculate the change in Gibbs free energy before and after the condensed phase species is added to the reaction system. If the change is less than zero, the condensed phase species is determined to meet the addition condition, the condensed phase species is added to the active species set, and the iterative solution is performed again. Specifically, for step S241, for each candidate condensed phase species not included in the current active species set, the change in Gibbs free energy of the entire reaction system before and after the addition of that condensed phase species is calculated. Gibbs free energy is a core parameter characterizing the thermodynamic stability of a system. If the change in Gibbs free energy after adding a certain condensed phase species is less than zero, it indicates that the species can spontaneously and stably generate under the current reaction conditions, meaning that the condensed phase species meets the addition conditions. At this point, the condensed phase species is added to the active species set, and based on the updated active species set, the iterative solution process for thermochemical equilibrium is re-executed.
[0047] S242. For condensed phase species that have been included in the active species set, if the number of moles of condensed phase species that converges in the iteration is negative, then the condensed phase species shall be removed from the active species set and the iteration shall be performed again. Specifically, for step S242, for each condensed phase species included in the current active species set, the corresponding molar number calculation result is verified after iterative convergence. If the iterative convergence molar number of a certain condensed phase species is negative, the result has no physical meaning, indicating that the condensed phase species cannot stably exist in the reaction system under the current reaction conditions. At this time, the condensed phase species is removed from the active species set, and the iterative solution process for thermochemical equilibrium is re-executed based on the updated active species set.
[0048] This embodiment uses bidirectional verification of the inclusion of condensed phases outside the aggregate and the exclusion of condensed phases within the aggregate to dynamically adjust the range of product species participating in the equilibrium calculation. This ensures that all thermodynamically stable condensed phase products are included in the calculation, while excluding negative molar number results that have no physical meaning. As a result, the combustion product composition and equilibrium temperature obtained in the final solution are more consistent with the actual reaction state, thereby improving the accuracy and physical rationality of the thermochemical equilibrium calculation results.
[0049] Furthermore, in some embodiments, step S3, "using the supersonic gas flow state parameters, the equilibrium temperature of the first reaction zone, and the product composition as the inlet boundary conditions of the second reaction zone, and calculating the axial distribution of the gas flow state parameters in the second reaction zone based on the quasi-one-dimensional governing equations," may specifically include: S31. Obtain the geometric configuration parameters of the second reaction zone, the injection position and injection pattern of reactants, the equivalence ratio, and set the integration start position and integration step size; Specifically, for step S31, the basic input conditions such as the geometric configuration parameters of the second reaction zone, the injection position and injection pattern of reactants, and the equivalence ratio are obtained, while the starting position and integration step size of the axial integration are set. Among them, the geometric configuration parameters are used to determine the variation law of the cross-sectional area of the channel in the combustion chamber along the axial direction and the total length of the combustion chamber, forming the geometric boundary of the flow calculation; the injection position and injection pattern of reactants are used to clarify the axial position and mass flow distribution of fuel injected into the combustion chamber; the equivalence ratio is used to limit the overall ratio of fuel to incoming air, determining the total heat release level of the combustion process; the integration starting position usually corresponds to the inlet section of the combustion chamber, and the integration step size is the step distance of the axial calculation. The smaller the step size, the higher the calculation accuracy, and the corresponding increase in the amount of calculation.
[0050] S32. At the current integration position, calculate the mixing efficiency between the reactants and the mainstream air based on the mixing efficiency model, and determine the mass flow rate of the reactants participating in the reaction and the heat released by the reaction; the mixing efficiency model determines the functional relationship between the mixing efficiency and the axial position based on the reactant injection method; for the support plate injection method, the change of mixing efficiency with injection distance is described by a function curve including empirical coefficients, the empirical coefficients include: the value range of the first coefficient is 1.0~1.1, the value range of the second coefficient is 3.6~3.8, and the value range of the third coefficient is 0.8~0.9; Specifically, in step S32, at the current integration position, the mixing efficiency between the reactants and the mainstream air is calculated based on the mixing efficiency model, thereby determining the actual mass flow rate of the reactants participating in the combustion reaction and the corresponding heat released by the reaction. The mixing efficiency model is determined according to the injection method of the reactants and is used to describe the functional relationship between the mixing efficiency and the axial injection distance. For the support plate injection method, the change of the mixing efficiency with the injection distance is described by a function curve containing three sets of empirical coefficients. The three sets of coefficients correspond to different fitting characteristics, with value ranges of 1.0~1.1, 3.6~3.8, and 0.8~0.9, respectively. The appropriate coefficient value can be selected according to the actual injection conditions.
[0051] S33. By combining the mass conservation equation, momentum conservation equation, energy conservation equation and gas state equation, and combining the relationship between Mach number and speed of sound, a set of first-order ordinary differential equations about Mach number, pressure and temperature is derived. Specifically, in step S33, the mass conservation equation, momentum conservation equation, energy conservation equation, and gas state equation are combined. Simultaneously, considering the relationship between Mach number and sound speed, the governing equations are derived and transformed, ultimately resulting in a system of first-order ordinary differential equations with Mach number, pressure, and temperature as the variables to be solved and axial distance as the independent variable. This form of equation system can be directly solved stepwise using numerical integration methods, transforming the complex supersonic combustion flow problem into a differential problem that can be calculated piecewise.
[0052] S34. The first-order ordinary differential equations are solved by axial integration using the fourth-order Runge-Kutta method to obtain the Mach number, total pressure, static pressure, total temperature, static temperature and velocity of the airflow at each integration position. Specifically, for step S34, the fourth-order Runge-Kutta method is used to integrate and solve the constructed first-order ordinary differential equations along the axial direction. Using this high-precision numerical integration algorithm, complete state parameters of the airflow at the current integration position, such as Mach number, total pressure, static pressure, total temperature, static temperature, and flow velocity, are calculated. The fourth-order Runge-Kutta method improves integration accuracy through multi-order slope correction, ensuring the accuracy of the airflow parameter solution under a reasonable step size.
[0053] S35. Determine if the current integration position has reached the outlet of the second reaction zone; if the current integration position has not reached the outlet of the second reaction zone, advance the integration position by one step, recalculate the mixing efficiency and heat release at the current integration position, and continue to perform the integration solution; if the current integration position has reached the outlet of the second reaction zone, output the airflow state parameters at the outlet of the second reaction zone; if the Mach number drops to 1 and supersonic flow cannot be maintained during the integration process, thermal blockage is determined to have occurred and the calculation is terminated; Specifically, for step S35, it is determined whether the current integration position has reached the outlet section of the second reaction zone. If the outlet has not been reached, the integration position is moved forward along the axial direction, and the mixing efficiency and heat release at the new position are recalculated. The integration solution process continues, forming a cyclical calculation logic of "mixing calculation - heat release calculation - integration solution - position advancement". If the outlet section has been reached, the complete airflow state parameters of the outlet of the second reaction zone are output. If, during the integration advancement process, the airflow Mach number drops to 1 and supersonic flow cannot be maintained, thermal congestion is determined to have occurred, the calculation is terminated directly, and the condition is marked as a thermal congestion condition.
[0054] This embodiment correlates the combustion heat release and axial flow process through a blending efficiency model. It uses a quasi-one-dimensional ordinary differential equation system combined with a high-precision integral algorithm to solve the problem step by step. This can efficiently obtain the distribution of airflow state parameters along the path in the supersonic combustion chamber. At the same time, it has the ability to automatically determine the thermal congestion condition and accurately reconstruct the supersonic combustion flow process in the second reaction zone. This provides reliable combustion chamber outlet boundary data for subsequent expansion process calculations and overall performance evaluation.
[0055] Furthermore, in some embodiments, the mass conservation equation includes an axial change term in mass flow rate caused by the addition of reactants; the momentum conservation equation includes a change term in cross-sectional area, a wall friction term, and a momentum change term caused by the addition of reactants; the energy conservation equation includes an enthalpy change term caused by the addition of reactants, a heat change term caused by the heat release of the reaction, and a heat loss term caused by wall convection. The first-order ordinary differential equation system includes the first-order differential equations of Mach number along the axis, pressure along the axis, and temperature along the axis. The first-order differential equations of Mach number along the axis, pressure along the axis, and temperature along the axis are solved by coupling the Mach number, pressure, temperature, rate of change of cross-sectional area, reactant addition rate, and reaction heat release rate at the current integration position.
[0056] Specifically, a mass conservation equation suitable for supersonic combustion flow with fuel injection is constructed. Based on the conventional one-dimensional steady flow mass conservation relationship, an additional axial change term in mass flow rate caused by reactant addition is introduced. This term quantifies the increase in total mass flow rate caused by the continuous injection of fuel-rich gas into the mainstream during the axial propagation of the airflow. It accurately characterizes the actual law of gradual increase in total mass flow rate along the axial direction when fuel injection is present in the combustion chamber, distinguishing it from the isoflow flow model without external mass input.
[0057] When constructing the momentum conservation equation, three types of terms affecting the axial momentum change of the airflow are simultaneously included: the cross-sectional area change term, the wall friction term, and the momentum change term caused by reactant addition. Specifically, the cross-sectional area change term characterizes the momentum adjustment effect of the change in flow area on airflow pressure and velocity when the combustion chamber channel expands or contracts axially; the wall friction term characterizes the axial momentum loss caused by the viscous frictional resistance between the high-speed airflow and the combustion chamber wall; and the reactant addition momentum term characterizes the effect of the momentum carried by the injected fuel-rich gas itself on the change in the total momentum of the airflow after being injected into the mainstream.
[0058] When constructing the energy conservation equation, three types of terms affecting the total enthalpy level of the airflow are simultaneously included: the enthalpy change term caused by the addition of reactants, the heat change term caused by the heat release of the reaction, and the heat loss term caused by wall convection. Among them, the reactant enthalpy change term is used to characterize the direct contribution of the enthalpy carried by the injected fuel-rich gas itself to the total enthalpy of the airflow after entering the mainstream; the heat change term caused by the heat release of the reaction is used to characterize the effect of the chemical energy released and converted into heat energy during the combustion reaction of fuel and air, which increases the total enthalpy of the airflow; the wall convection heat loss term is used to characterize the total energy loss of the airflow caused by the transfer of heat from the high-temperature airflow to the combustion chamber wall through convection heat transfer.
[0059] By combining the aforementioned mass conservation equation, momentum conservation equation, energy conservation equation, and gas state equation, and using the thermodynamic correspondence between Mach number and speed of sound for mathematical derivation and variable substitution, we finally obtain the first-order differential equations for Mach number, pressure, and temperature along the axial direction. These three equations together constitute a complete set of first-order ordinary differential equations. These three differential equations are strongly coupled; the right-hand side of each equation contains multiple parameters at the current integration position, including Mach number, pressure, temperature, rate of change of cross-sectional area, reactant addition rate, and heat release rate. It is impossible to solve for a single variable individually; the axial rates of change of the three parameters must be obtained simultaneously through simultaneous calculations to fully characterize the axial evolution of the airflow state.
[0060] This embodiment incorporates various practical operating conditions such as mass addition, wall friction, combustion heat release, and wall heat dissipation into the conservation equations, and derives a set of first-order ordinary differential equations with strong coupling of multiple parameters. This allows for a more complete and accurate reproduction of the flow and combustion process within the supersonic combustion chamber, making the solution results of the airflow state parameters more consistent with the actual working characteristics of the combustion chamber, and effectively improving the physical reality and accuracy of the calculation of the combustion process within the combustion chamber.
[0061] Furthermore, in some embodiments, the wall friction term in the momentum conservation equation is determined based on the wall friction coefficient. The calculation method for the wall friction coefficient may specifically include: The Reynolds number is temperature-corrected based on the dynamic viscosity at the mainstream airflow temperature and the dynamic viscosity at the wall temperature, and the corrected Reynolds number at the current integration position is calculated. The wall friction coefficient is calculated based on the modified Reynolds number and a semi-empirical formula.
[0062] Specifically, when calculating the wall friction coefficient, a temperature correction is first applied to the Reynolds number to offset the influence of the temperature difference between the mainstream airflow and the wall on the gas viscosity parameters, thereby improving the accuracy of the boundary layer flow characteristics calculation. First, the mainstream airflow temperature and the working wall temperature at the current integration position are obtained. Based on the variation of gas dynamic viscosity with temperature, the dynamic viscosity corresponding to the mainstream airflow temperature and the dynamic viscosity corresponding to the wall temperature are calculated respectively. Using the initial Reynolds number calculated based on the physical properties at the mainstream temperature, the Reynolds number is corrected by combining the ratio of the dynamic viscosity at the two temperatures, ultimately obtaining a corrected Reynolds number that considers the wall temperature effect. The corrected Reynolds number can more realistically reflect the flow viscous state within the wall boundary layer, avoiding the distortion of flow characteristics caused by directly using the mainstream viscosity under large temperature difference conditions.
[0063] After obtaining the corrected Reynolds number, a semi-empirical friction coefficient formula adapted to supersonic turbulent flow is used. The corrected Reynolds number is substituted into the formula for calculation, and the wall friction coefficient at the current integration position is output. This semi-empirical formula is constructed based on a large amount of experimental data and flow law induction, and can establish a quantitative correspondence between the Reynolds number, wall shear stress, and friction coefficient. Using a temperature-corrected Reynolds number in the calculation allows the friction coefficient result to fully reflect the influence of wall temperature on gas viscosity, and is more in line with the actual wall resistance characteristics of high-temperature, high-speed airflow.
[0064] This embodiment introduces dynamic viscosity under dual operating conditions of mainstream temperature and wall temperature, performs temperature correction on the Reynolds number, and then solves for the wall friction coefficient. This effectively compensates for the viscosity calculation deviation caused by the large temperature difference between the airflow and the wall, making the calculated wall friction coefficient more consistent with the actual resistance characteristics of high-temperature supersonic flow, and improving the accuracy and physical reality of friction loss calculation during combustion chamber flow.
[0065] Furthermore, in some embodiments, the wall heat loss term in the energy conservation equation is determined based on the wall heat flux density. The calculation method for the wall heat flux density may specifically include: Based on the reference enthalpy method, the reference enthalpy is calculated according to the enthalpy of the main airflow, the enthalpy of the adiabatic wall, and the enthalpy at the wall temperature, and the reference temperature is determined based on the reference enthalpy. Calculate the convective heat flux density between the airflow and the wall based on the thermodynamic parameters and convective heat transfer coefficient at the reference temperature.
[0066] Specifically, the reference enthalpy method is used to correct for high-speed convective heat transfer when calculating the wall heat flux density. First, three basic parameters are obtained: the mainstream airflow enthalpy, the adiabatic wall enthalpy, and the enthalpy corresponding to the wall temperature at the current location. The mainstream airflow enthalpy corresponds to the thermodynamic enthalpy level of the airflow in the core region of the combustion chamber; the adiabatic wall enthalpy corresponds to the enthalpy level of the airflow after adiabatic stagnation at the wall; and the enthalpy at the wall temperature corresponds to the gas enthalpy at the wall's operating temperature. Based on the empirical calculation relationship of these three enthalpy values, a reference enthalpy is obtained for property correction. Then, based on the one-to-one correspondence between gas enthalpy and temperature, the corresponding reference temperature is determined by inversely using the reference enthalpy, and this temperature serves as the base temperature for subsequent convective heat transfer property calculations.
[0067] Based on a defined reference temperature, the thermodynamic and transport properties of the gas at that temperature are obtained. Combined with the airflow parameters and channel geometry at the current location, the convective heat transfer coefficient between the airflow and the wall is calculated. Finally, based on the fundamental principle of convective heat transfer, the convective heat flux density between the airflow and the wall is obtained by calculating the convective heat transfer coefficient and the enthalpy difference between the airflow and the wall. This heat flux density can quantitatively characterize the amount of heat transferred from the airflow to the wall per unit time and per unit area via convection.
[0068] This embodiment uses the reference enthalpy method to correct the physical properties and temperature of the convective heat transfer process under supersonic conditions. It can effectively consider the influence of the aerodynamic heating effect of high-speed airflow on wall heat transfer, making the calculated wall heat flux density more consistent with the actual heat transfer characteristics of the combustion chamber, accurately quantifying the wall heat dissipation level, and providing accurate data for combustion chamber heat loss assessment.
[0069] Furthermore, in some embodiments, step S4, "using the airflow state parameters at the outlet of the second reaction zone as the inlet parameters of the expansion process, calculating the axial distribution of the airflow state parameters within the nozzle based on the expansion process model, and obtaining the airflow state parameters at the nozzle outlet section," may specifically include: S41. Ignoring the addition of reactants and the heat released by the reaction, and only considering the acceleration effect of cross-sectional area change and wall friction on the airflow, the mass conservation equation, momentum conservation equation, energy conservation equation and gas state equation are combined to derive a set of first-order ordinary differential equations about Mach number, pressure and temperature. Specifically, for step S41, the calculation boundary for the nozzle expansion process is determined. This process involves no additional reactant injection and no further combustion heat release reaction; therefore, the effects of reactant addition and reaction heat release are ignored. Only the acceleration and loss effects of the airflow caused by changes in the channel cross-sectional area and wall friction are considered. Based on this, the mass conservation equation, momentum conservation equation, energy conservation equation, and gas state equation are combined. Mathematical transformations and variable substitutions are performed using the thermodynamic correspondence between Mach number and sound speed. Finally, a set of first-order ordinary differential equations is derived and simplified, with Mach number, pressure, and temperature as the variables to be solved and axial distance as the independent variable. The cross-sectional area change term corresponds to the driving effect of the change in the flow area of the nozzle expansion channel on the airflow expansion and acceleration, while the wall friction term corresponds to the momentum and total pressure losses caused by the viscous resistance between the high-speed airflow and the nozzle wall.
[0070] S42. Using the airflow state parameters at the outlet of the second reaction zone as the initial values for integration, the fourth-order Runge-Kutta method is used to integrate along the axial direction to obtain the airflow Mach number, total pressure, static pressure, total temperature, static temperature and velocity at each integration position in the tail nozzle. Specifically, for step S42, the airflow state parameters at the outlet section of the second reaction zone are used as the initial inlet values for the expansion process integral calculation, i.e., the initial airflow state at the nozzle inlet section, including core parameters such as Mach number, static pressure, static temperature, total pressure, total temperature, and velocity. The fourth-order Runge-Kutta method is used as the numerical integration algorithm, and the integration is performed step-by-step along the nozzle axis from the inlet to the outlet. For each integration step, based on the current channel cross-sectional area, airflow state, and wall friction characteristics, the airflow Mach number, total pressure, static pressure, total temperature, static temperature, and velocity at the next position are calculated. This step-by-step integration continues until the nozzle outlet section is reached, ultimately obtaining the distribution of airflow state parameters at each axial position within the entire nozzle flow channel.
[0071] The expansion process calculation method provided in this embodiment matches the working characteristics of the tail nozzle without combustion and without external mass input. It retains only the core influencing factors of cross-sectional area change and wall friction. By solving the quasi-one-dimensional ordinary differential equation system combined with a high-precision integral algorithm, it can efficiently and accurately simulate the expansion and acceleration process of the gas in the tail nozzle, and obtain the airflow state parameters of the entire flow channel and the outlet section of the tail nozzle, providing reliable outlet flow boundary data for the overall performance calculation of the power unit.
[0072] Furthermore, in some embodiments, the method may further include: S61. During multiple executions, the reactant composition or reaction pressure in the first reaction zone is changed; Specifically, for step S61, during the execution of multiple rounds of performance evaluation, the input boundary conditions of the first reaction zone are adjusted in a targeted manner. The object of the change is either the initial component ratio of the reactants in the first reaction zone or the working reaction pressure of the first reaction zone. Each parameter adjustment corresponds to an independent and complete performance evaluation calculation. The influence of a single parameter can be analyzed individually using a univariate adjustment method, or multiple parameters can be adjusted simultaneously for analysis. All evaluation rounds maintain consistent operating conditions except for the input parameters of the first reaction zone to ensure the effectiveness of the comparative analysis.
[0073] S62. Based on the changed reactant composition and reaction pressure, recalculate the equilibrium temperature and product composition of the first reaction zone according to the thermochemical equilibrium principle, and obtain the product composition and equilibrium temperature of the first reaction zone and the overall performance parameters of the power unit for each execution. Specifically, for step S62, for each operating condition after parameter changes, the adjusted reactant composition and reaction pressure are used as inputs. Based on the principle of thermochemical equilibrium, the combustion equilibrium state of the first reaction zone is re-solved, yielding the equilibrium temperature and product composition of the first reaction zone for each operating condition. Based on this, the subsequent combustion, expansion, and performance calculation processes are fully executed to obtain the overall performance parameters of the power unit corresponding to each operating condition. The calculation logic and solution standards for all rounds are kept consistent to ensure that the output results of each group are comparable.
[0074] S63. Based on the differences between the overall performance parameters obtained in each execution, determine the correspondence between the changes in reactant composition or reaction pressure in the first reaction zone and the overall performance of the power unit; Specifically, for step S63, the overall performance parameters obtained from multiple rounds of evaluation are compared, and the differences in performance indicators under different working conditions are analyzed. The differences in performance parameters are correlated with the component adjustment and pressure adjustment in the first reaction zone, and the direction and magnitude of the changes in reactant components and reaction pressure are clarified, as well as the quantitative correspondence between them and the overall performance trend and magnitude of the power unit.
[0075] This embodiment, through multiple rounds of adjusting the input parameters of the first reaction zone and conducting a complete performance evaluation, can quantify and clarify the influence law between reactant components, reaction pressure and the overall performance of the power unit. It can provide a quantitative basis for optimizing the viscoelastic material formulation and matching the working parameters of the gas generator, effectively supporting the iterative design and parameter optimization of the solid scramjet power unit.
[0076] To better understand the performance evaluation method for a solid scramjet power unit based on primary combustion effect provided in this embodiment, the following will describe each process of the method in detail. The actual working process of the solid scramjet power unit includes the pressure process in the inlet duct, the primary combustion process in the gas generator, the secondary combustion process in the supersonic combustion chamber, and the expansion process of the tail nozzle. The two combustion processes are the main power sources of the power unit and are crucial to its overall performance. The oxygen-deficient solid viscoelastic material first undergoes primary combustion in the gas generator to produce high-temperature, fuel-rich gas. Subsequently, the fuel-rich gas is injected into the supersonic combustion chamber, where it mixes with the high-speed incoming flow compressed through the inlet duct and undergoes secondary combustion. The high-temperature gas after combustion expands and accelerates through the tail nozzle, ultimately generating thrust. A schematic diagram of the solid scramjet power unit is shown in the attached diagram. Figure 2 The models for each working process are as follows.
[0077] Compression process: The operating environment of a scramjet propulsion system is constrained by atmospheric parameters, which are functions of flight altitude. Using the standard atmospheric conditions of the United States in 1976, parameters such as static temperature, static pressure, and atmospheric density at section 0 can be obtained at different flight altitudes. The relationship between air mass flow rate, incoming velocity, temperature, and pressure at section 0 can be given by the following formula: Therefore, when the flight altitude of the power unit is known H and Mach number Ma At time 0, the parameters at section 0 can all be obtained through calculation. The following calculation uses the assumption of adiabatic compression to solve for the airflow parameters at section 1. Since the total enthalpy remains constant during compression, the following relationship exists: According to thermodynamics, the specific heat capacity at constant pressure is... C p for: At this point, we can obtain the following formula for calculating the Mach number at section 1: Because the combustion process of a scramjet power unit has a maximum temperature limit at the inlet, typically... T Since 1 < 1600K, the Mach number at section 1 must also be limited accordingly, usually by a ratio greater than 0.4.
[0078] The atmospheric temperature variation is small within the operating altitude range of the scramjet propulsion system; the temperature change is less than 40K as the flight altitude changes from 20km to 40km. The primary objective of the compression process is to provide the required cycle static temperature ratio in a controllable and reliable manner within the operating range of the scramjet propulsion system. This minimizes aerodynamic losses. Cyclic static temperature ratio The definition of is: in T 0 represents the inlet airflow temperature. T 1 represents the outlet temperature during the air intake compression process. When the flight altitude is 30km, the cycle static temperature ratio is approximately 7.0.
[0079] The actual compression efficiency of a scramjet power unit is limited by a variety of factors, making it difficult to achieve the ideal state. (See attached image) Figure 3 The results of a systematic study of various compression system variables within the full expected range are presented. From these data, it can be concluded that the adiabatic compression efficiency is almost directly related to the cycle static temperature ratio and the incoming Mach number. Ma The efficiency is not related to 0, but is closely related to the number of oblique shock waves selected; and at least two, preferably three or four oblique shock waves must be used to achieve ideal compression performance. Based on this, the present invention assumes that the compression process consists of four oblique shock waves, at which point the adiabatic compression efficiency is... η c Always keep it above 0.9, take η c =0.9.
[0080] The four oblique shock waves have the same geometric rotation amplitude, and the flow rotation angle of each shock wave is... and the i Shock angle of the shock wave relative to its local upstream direction It can be solved by the following formula: The calculation results are taken from the weak shock root, and then the downstream Mach number is obtained using the normal shock wave relationship. Ma i Compared to temperature, through iteration Make: The mass conservation equation, energy conservation equation, and gas state equation involved in the compression process can be expressed as follows: In summary, when the power unit flies at an altitude H and Mach number Ma Given the Mach number ratio and compression efficiency, all parameters at section 1 can be solved by calculation.
[0081] One combustion process: Unlike liquid-fuel scramjet power units, solid-fuel scramjet power units involve two combustion processes during operation: primary combustion in the gas generator and secondary combustion in the supersonic combustion chamber. The fuel-rich gas composition and temperature required for secondary combustion are determined by the primary combustion process, making effective modeling of the primary combustion process essential. This method, based on the thermochemical equilibrium assumption, treats the primary combustion process during stable operation in the gas generator as isobaric adiabatic combustion. It employs the Gibbs minimum free energy method, using a descending Newton-Raphson method for iterative solution to obtain the equilibrium temperature and product composition of the fuel under given combustion conditions. The obtained equilibrium temperature and product composition represent the initial temperature and composition of the fuel-rich gas during secondary combustion. It is assumed that the primary combustion products contain NS components, with gaseous products consisting of components 1 to NG and condensed-phase products consisting of components NG+1 to NS. The governing equations involved in the calculation are as follows: Equations 12 to 15 represent the total element conservation equation, the condensed phase chemical potential equilibrium equation, the total mole number conservation equation in the gas phase, and the enthalpy equilibrium equation, respectively, for a single combustion process. For the first j The number of moles of each product species The molar number of the gas phase component. The dimensionless Lagrange multipliers corresponding to the elements. The number of chemical elements in the system. The first determined by the composition of the reactants k Total amount of all elements The first in the composition of the current iteration product k Total amount of all elements Components j Standard molar enthalpy, Components j Standard isobaric molar heat capacity Components j chemical potential, R This is the universal gas constant. T eq To balance the temperature.
[0082] The chemical potential of different components can be solved using the following formula: in: Dimensionless Lagrange multipliers It can be obtained from the following formula: in It corresponds to the first i Lagrange multipliers of elements.
[0083] When using the descent Newton-Raphson method for iterative solution, the descent control factor is calculated as follows: In the formula It is a Class I decrease control factor used to limit the correction magnitude for temperature, total gas phase moles, and major gas phase species; For gaseous species, the decrease control factor is the final control factor. Take 1, and The minimum value in.
[0084] Since a low-activity condensed phase may be generated during a single combustion process, the following judgment needs to be made when adding condensed phases with unknown activity during the iteration process: If this condition is met, it means that adding this condensate will reduce the Gibbs free energy, and it should be added to the active species set and iterated again.
[0085] The thermodynamic functions of the specific heat capacity, enthalpy, and entropy of each component can all be expressed as functions of temperature using least squares coefficients, where the coefficients are... a 1~ a 9 can be obtained from the NIST Chemistry WebBook, and the specific function form is as follows: The NIST Chemistry WebBook database contains thermodynamic data for commonly used components, which meets the requirements of this invention. When the temperature exceeds the range in the database, the temperature is limited to the range specified in the database by extrapolating from the highest temperature range.
[0086] In summary, when the initial composition of the viscoelastic material is known, the equilibrium temperature and product composition of the primary combustion process can be obtained by using the composition and initial reaction conditions as inputs, which can then be used as the initial conditions for the secondary combustion process. The specific solution process for the primary combustion process is shown in the appendix. Figure 5The solution process is based on the governing equations and employs a descending Newton-Raphson method. Unknowns include the molar number of gaseous components, the molar number of condensed phase components, temperature, and Lagrange multipliers. In the initial iterations, all gaseous components are evenly distributed; the initial liquid or solid phase content is assumed to be zero; and the initial temperature is set to 3800 K by default. When multiple condensates simultaneously meet the addition conditions in the Gibbs stability check, all stable phases meeting the conditions are added simultaneously to establish an extended component set, and then the Gibbs minimization is performed again.
[0087] Secondary combustion process The mass conservation equation at section 1 can be determined by the following equation: in The mass flow rate of the airflow. For airflow density, U For airflow velocity, A Let be the cross-sectional area of the channel through which the airflow passes. The corresponding differential form can be expressed as: Considering pressure changes, wall friction, and fuel addition, based on momentum conservation, the standard one-dimensional steady momentum equation for a scramjet power unit can be expressed as: in For airflow static pressure, is the coefficient of friction.
[0088] The coefficient of friction can be calculated using a semi-empirical formula: In the formula R e The Reynolds number is... The dynamic viscosity at the mainstream airflow temperature. is the dynamic viscosity at the wall temperature of the airflow.
[0089] The formula for calculating the Reynolds number is as follows: in D It is the hydraulic diameter of the channel through which the airflow passes.
[0090] The most commonly used empirical formula based on the kinetic theory of molecules is selected for solving the kinetic viscosity coefficient: in =1.7894×10 5 N˙s / m 2 The dynamic viscosity is at a reference temperature of 288.15 K. Take 100.4K.
[0091] The energy conservation equation takes into account the enthalpy change caused by fuel addition, the heat change caused by combustion heat release, and the heat loss caused by thermal convection between the main stream and the wall. Its specific differential form can be expressed as: in For the specific heat of airflow at constant pressure, The increase in enthalpy is due to fuel combustion. The change in enthalpy caused by fuel injection, The mass flow rate of the fuel. The wall heat flux density is given.
[0092] The wall heat flux density is calculated based on the Eckert reference enthalpy method, and the expression for the reference enthalpy is: According to the convection heat transfer formula, the wall heat flux formula can be obtained as follows: in h e The enthalpy of the main airflow. h w This represents the enthalpy of the airflow at the wall temperature. h aw This is the enthalpy of the adiabatic wall; the values marked with * are reference temperatures. T The variable * is used to determine the reference temperature, which can be obtained from the reference enthalpy. Pr* The Prandtl number at the reference temperature is taken as 0.71.
[0093] For an ideal gas, the molecules themselves do not occupy volume, the intermolecular forces are negligible, and its equation of state can be expressed as: The differential form can be expressed as: The relationship between Mach number and speed of sound is as follows: For an ideal gas, the formula for calculating the speed of sound is as follows: Solving the above equations of mass conservation, momentum conservation, energy conservation, and the ideal gas law simultaneously, and combining this with the conversion relationship between Mach number and speed of sound, we can obtain the following differential expressions for Mach number, pressure, and temperature: In the formula, From the differential expressions for Mach number, pressure, and temperature mentioned above, it can be seen that when the incoming flow parameters, combustion chamber geometry, mass addition rules, and combustion heat release rules are known, the differential equations can be solved using the fourth-order Runge-Kutta method to obtain the parameters within the combustion chamber as a function of axial position. x The changes.
[0094] One of the more tedious steps in the quasi-one-dimensional analysis of scramjet power units is calculating the heat distribution of the combustion process. For solid scramjet power units, the secondary combustion of fuel-rich gas in the supersonic combustion chamber is its main power source. From the governing equations and calculation results of the primary combustion process modeling, it is known that the fuel-rich gas produced by the primary combustion of a solid scramjet power unit has a complex composition, containing not only gaseous components but also a large number of solid particles and condensed phase components. Therefore, the actual secondary combustion process is extremely complex, involving the mixing and combustion of gas-solid two-phase flows. Since the current understanding of the mechanism of the secondary combustion process is still insufficient, this paper models the secondary combustion process based on the chemical equilibrium assumption. In this model, the chemical reaction rate is considered infinite, and the combustion efficiency mainly depends on the mixing efficiency, that is, the combustion efficiency of the fuel-rich gas and the incoming flow is approximately equal to the mixing efficiency. The final products include potentially generated condensed phases and unreacted condensed phases.
[0095] This scheme uses the following model to calculate the blending length Lmix and blending efficiency. The relevant calculation formulas are as follows: in The mass flow rate that can be used for combustion. Fuel mass flow rate, The nozzle diameter is... U f , and The speed, sound speed, and density of fuel injected into the combustion chamber. L in This refers to the distance from the fuel injection point to the combustion chamber inlet. U a , a a and for L in The velocity, sound speed, and density of the main current at the location, and constants. K The value of K is related to the fuel injection method and can be summarized through high-dimensional CFD simulation results. Cristian Birzer concluded in his research that for the case of rearward fuel injection using a thrust plate, K=390 can be taken; for other injection methods, KThe value needs further calibration through experiments. In this scheme, the injection process of the fuel-rich gas in the scramjet power unit adopts a support plate injection method, and the coefficient in the mixing efficiency formula is... a Take 1.06492, k Take 3.69639, d Take 0.80586.
[0096] Expansion process: As the main working component of the power unit, the exhaust nozzle converts the heat energy released from fuel combustion in the combustion chamber into the kinetic energy required for the power unit to operate. Compared with the secondary combustion process, the expansion process does not include the addition of rich fuel gas or the heat release from combustion. Therefore, in this scheme, the calculation of the exhaust nozzle adopts a quasi-one-dimensional treatment similar to the secondary combustion process, that is, considering the effects of exhaust nozzle area changes and wall friction. The differential expressions of Mach number, pressure, and temperature involved are as follows: Following the solution approach for the secondary combustion process, the airflow parameters along the nozzle can be obtained by using the fourth-order Runge-Kutta method to solve the differential equations.
[0097] Performance parameters: Compared to internally charged solid-fuel scramjet power units, solid-fuel scramjet power units can control the fuel flow rate for secondary combustion by adjusting the throat area of the combustion generator. Specific thrust I a Mass specific impulse I m and volumetric impulse I v These are important indicators for evaluating the overall performance of a scramjet power unit. Ignoring overflow losses between the intake manifold capture flow and the actual combustion chamber flow, their calculation formulas are as follows: In summary, the compression, primary combustion, secondary combustion, and expansion processes of the solid scramjet power unit were modeled, and the solution methods for the relevant performance parameters were also given.
[0098] Reference Figure 2 The solid scramjet propulsion system mainly consists of an intake, an isolation section, a gas generator, a supersonic combustion chamber, and a tail nozzle. Its basic working process is as follows: the oxygen-deficient solid viscoelastic material first undergoes a fuel-rich combustion process in the gas generator, producing high-temperature fuel-rich gas; subsequently, the fuel-rich gas is injected into the supersonic combustion chamber, where it mixes with the high-speed incoming flow compressed through the intake and undergoes secondary combustion; the high-temperature gas after combustion expands and accelerates through the tail nozzle, ultimately generating thrust.
[0099] Reference Figure 3 adiabatic compression efficiency η c The relationship between the incoming Mach number, the circulating static temperature ratio, and the number of oblique shock waves shows that the adiabatic compression efficiency is almost directly related to the circulating static temperature ratio and the incoming Mach number. Ma The number of oblique shock waves is irrelevant; rather, it is closely related to the number of oblique shock waves selected. At least two, and preferably three or four, oblique shock waves must be used to achieve the desired compression performance.
[0100] Reference Figure 4 The combustion chamber control unit of the solid scramjet power unit is shown in the figure, and the following assumptions are made about the flow field inside the combustion chamber: 1) The geometry of the combustion chamber and the airflow parameters of the quasi-one-dimensional flow are both functions of the axial distance x, and the physical property parameters are the same at all points on the same cross section; 2) The airflow in the combustion chamber is considered to be steady under given operating conditions, and the airflow parameters will not change with time; 3) The working fluid in the combustion chamber is regarded as a continuous medium, which approximately satisfies the ideal gas law. The molecules themselves do not occupy volume, and the intermolecular forces can be ignored.
[0101] Combination Figure 5 and 6 The specific work process of the implementation plan is explained as follows: Reference Figure 5 In the process of solving a combustion iteration, the pressure of the first combustion is first given. P reaction The total elemental composition of the viscoelastic material was calculated. bi 0 Next, active species from the gas phase and condensed phase are selected, and their current thermodynamic properties are calculated and output as input conditions for the governing equations. Unknown correction parameters, including solid-phase species, total gas phase molar number, temperature, and Lagrange multipliers, are obtained by solving the governing equations, and the gas phase species correction parameters and descent control factors are output. The variables are updated using the obtained correction parameters, and convergence is checked for each update result. If the convergence condition is not met, the data is re-input into the governing equations for solving, and the process is repeated until convergence. When the convergence condition is met, Gibbs stability checks are performed on condensed phases not included in the active set. If the iteration criterion is met, they are added to the candidate species set and iterated again. If the convergence molar number of an included condensed phase is negative, it is removed and iterated again. This process is repeated until the stability conditions of all condensed phases are met, and the final equilibrium temperature is output. T eq The product composition and thermodynamic properties serve as boundary conditions for the subsequent secondary combustion process.
[0102] Reference Figure 6In the quasi-one-dimensional analysis process, the flight altitude of the power unit is first given. H ,Mach number Ma 0. The Mach number ratio and compression efficiency can be calculated using the compression process formula to obtain all parameters at section 1, i.e., the inlet airflow conditions for secondary combustion. Next, given the oxygen-deficient viscoelastic material composition and reaction pressure, the reaction equilibrium temperature and product composition can be solved using the primary combustion solution process, i.e., the initial temperature and composition of the fuel-rich gas in secondary combustion. Then, the airflow conditions and fuel-rich gas conditions are input into the secondary combustion solution process, and the influence coefficient method is used to solve for the variation of airflow parameters along the path under given equivalence ratio and combustion chamber configuration conditions. Then, the combustion chamber outlet parameters are used as the inlet parameters for the expansion process, and the airflow parameters along the tailpipe are obtained through one-dimensional calculation of the expansion process. Finally, the overall performance parameters of the power unit are calculated using the performance solution formula.
[0103] In summary, the performance evaluation method for solid scramjet power units based on primary combustion effects provided in this embodiment proposes a quasi-one-dimensional performance analysis model for solid scramjet power units that considers the primary combustion process within the gas generator. Compared with existing one-dimensional models, this method considers the impact of primary combustion effects on the secondary combustion process in the combustion chamber and the overall performance of the power unit, improving the reliability of predictions in the initial design phase of the power unit. The method employs the Gibbs minimum free energy method to solve the primary combustion process, while the secondary combustion process is based on the chemical equilibrium assumption. The solution process is based on the influence coefficient method. Compared with two-dimensional and three-dimensional calculation models, this method involves less computation and is faster, better meeting the rapid evaluation needs in the initial design phase of solid scramjet power units. It can comprehensively consider the effects of factors including the addition of fuel-rich gas, combustion heat release, changes in the cross-sectional area of the combustion chamber, wall friction, and wall heat dissipation. Furthermore, it can analyze the impact of parameter changes during primary combustion on the power unit.
[0104] The technical features of the above embodiments can be combined in any way. For the sake of brevity, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, they should be considered to be within the scope of this specification.
[0105] The above provides a detailed description of a solid scramjet power device performance evaluation method based on primary combustion effect provided in the embodiments of this application. Specific examples are used in this embodiment to illustrate the principle and implementation of this application. The description of the above embodiments is only for the purpose of helping to understand the method and core idea of this application. At the same time, for those skilled in the art, there will be changes in the specific implementation and application scope based on the idea of this application. Therefore, the content of this specification should not be construed as a limitation of this application.
Claims
1. A method for evaluating the performance of a solid scramjet power unit based on the primary combustion effect, characterized in that, include: Obtain flight operating parameters of the power unit and calculate supersonic airflow state parameters after compression through the inlet based on the compression process model; The initial composition of reactants and reaction pressure in the first reaction zone are obtained, and the equilibrium temperature and product composition of the first reaction zone are calculated based on the principle of thermochemical equilibrium. Using the supersonic airflow state parameters, the equilibrium temperature of the first reaction zone, and the product composition as the inlet boundary conditions of the second reaction zone, the distribution of the airflow state parameters along the axial direction in the second reaction zone is calculated based on the quasi-one-dimensional control equations. Using the airflow state parameters at the outlet of the second reaction zone as the inlet parameters of the expansion process, the distribution of the airflow state parameters along the axial direction in the tail nozzle is calculated based on the expansion process model, and the airflow state parameters at the outlet section of the tail nozzle are obtained. The overall performance parameters of the power unit are calculated based on the airflow state parameters of the tailpipe outlet section.
2. The performance evaluation method for solid scramjet power units based on primary combustion effect according to claim 1, characterized in that, The acquisition of flight operating parameters of the power unit, based on the compression process model, calculates the supersonic airflow state parameters after compression through the inlet, including: Based on the flight altitude of the power unit, standard atmospheric parameters are queried to calculate the static temperature, static pressure, atmospheric density, and air mass flow rate of the free flow. The intake outlet Mach number ratio and compression efficiency are preset, and the total temperature, static temperature, static pressure, Mach number and velocity of the airflow at the intake outlet are solved based on the adiabatic compression assumption and the total enthalpy conservation relationship. The loss of airflow state parameters during the compression process is calculated based on the preset compression efficiency, and the airflow state parameters of the inlet outlet section are output.
3. The performance evaluation method for solid scramjet power units based on primary combustion effect according to claim 1, characterized in that, The process of obtaining the initial reactant composition and reaction pressure of the first reaction zone, and calculating the equilibrium temperature and product composition of the first reaction zone based on the principle of thermochemical equilibrium, includes: Obtain the initial components of the reactants in the first reaction zone and count the total atomic count of each chemical element in the reaction system; Gas-phase active species and candidate condensed phase species are selected to construct an active species set, and the thermodynamic property parameters of each species are calculated at the current temperature; Based on the conservation of total elements, chemical potential balance, conservation of total moles in the gas phase, and enthalpy balance, a set of governing equations was constructed, and the equilibrium temperature and the mole number of each species were solved iteratively using the descending Newton-Raphson method. The Gibbs stability of condensed phases not included in the set of active species is determined. If the inclusion condition is met, the condensed phase is included in the set of active species, and the iterative solution process is re-executed until all condensed phases outside the set of active species do not meet the inclusion condition. The equilibrium temperature, molar number of each species, and thermodynamic properties of each component obtained after iterative convergence are used as the output results of the first reaction zone.
4. The performance evaluation method for solid scramjet power units based on primary combustion effect according to claim 3, characterized in that, The step of performing Gibbs stability determination on condensed phases not included in the set of active species, and if the inclusion condition is met, then including the condensed phase in the set of active species, and re-executing the iterative solution process until all condensed phases outside the set of active species no longer meet the inclusion condition, includes: For each condensed phase species outside the set of active species, calculate the change in Gibbs free energy before and after the condensed phase species is added to the reaction system. If the change is less than zero, the condensed phase species is determined to meet the addition condition, the condensed phase species is added to the set of active species, and the iterative solution is performed again. If the number of moles of condensed phase species that have been included in the set of active species is negative, then the condensed phase species is removed from the set of active species and the iterative solution is performed again.
5. The performance evaluation method for a solid scramjet power unit based on primary combustion effect according to claim 1, characterized in that, The step of using the supersonic airflow state parameters, the equilibrium temperature of the first reaction zone, and the product composition as the inlet boundary conditions of the second reaction zone, and calculating the axial distribution of the airflow state parameters in the second reaction zone based on the quasi-one-dimensional governing equations, includes: Obtain the geometric configuration parameters of the second reaction zone, the injection position and injection pattern of the reactants, and the equivalence ratio; set the integration start position and integration step size. At the current integration position, the mixing efficiency between the reactants and the mainstream air is calculated based on the mixing efficiency model to determine the mass flow rate of the reactants participating in the reaction and the heat released by the reaction. By combining the mass conservation equation, momentum conservation equation, energy conservation equation, and gas state equation, and combining the relationship between Mach number and speed of sound, a set of first-order ordinary differential equations concerning Mach number, pressure, and temperature is derived. The first-order ordinary differential equations were solved by axial integration using the fourth-order Runge-Kutta method to obtain the airflow Mach number, total pressure, static pressure, total temperature, static temperature and velocity at each integration position; Determine whether the current integration position has reached the outlet of the second reaction zone; if the current integration position has not reached the outlet of the second reaction zone, advance the integration position by one step, recalculate the mixing efficiency and heat release at the current integration position, and continue to perform the integration solution; if the current integration position has reached the outlet of the second reaction zone, output the airflow state parameters at the outlet of the second reaction zone; if the Mach number drops to 1 and supersonic flow cannot be maintained during the integration process, thermal blockage is determined to have occurred and the calculation is terminated.
6. The performance evaluation method for a solid scramjet power unit based on primary combustion effect according to claim 5, characterized in that, The mass conservation equation includes an axial change term in mass flow rate caused by the addition of reactants; the momentum conservation equation also includes a change term in cross-sectional area, a wall friction term, and a momentum change term caused by the addition of reactants; the energy conservation equation also includes a change term in enthalpy caused by the addition of reactants, a change term in heat caused by the heat released from the reaction, and a term in heat loss caused by wall convection. The set of first-order ordinary differential equations includes the first-order differential equations of Mach number along the axial direction, pressure along the axial direction, and temperature along the axial direction. The first-order differential equations of Mach number along the axial direction, pressure along the axial direction, and temperature along the axial direction are solved by coupling the Mach number, pressure, temperature, rate of change of cross-sectional area, reactant addition rate, and reaction heat release rate at the current integration position.
7. The performance evaluation method for a solid scramjet power unit based on primary combustion effect according to claim 6, characterized in that, The wall friction term in the momentum conservation equation is determined based on the wall friction coefficient, which is calculated in the following ways: The Reynolds number is temperature-corrected based on the dynamic viscosity at the mainstream airflow temperature and the dynamic viscosity at the wall temperature, and the corrected Reynolds number at the current integration position is calculated. The wall friction coefficient is calculated based on the modified Reynolds number and a semi-empirical formula.
8. The performance evaluation method for a solid scramjet power unit based on primary combustion effect according to claim 6, characterized in that, The wall heat loss term in the energy conservation equation is determined based on the wall heat flux density, and the calculation method for the wall heat flux density includes: Based on the reference enthalpy method, the reference enthalpy is calculated according to the enthalpy of the main airflow, the enthalpy of the adiabatic wall, and the enthalpy at the wall temperature, and the reference temperature is determined according to the reference enthalpy; The convective heat flux density between the airflow and the wall is calculated based on the thermodynamic parameters and convective heat transfer coefficient at the reference temperature.
9. The performance evaluation method for a solid scramjet power unit based on primary combustion effect according to claim 1, characterized in that, The step of using the airflow state parameters at the outlet of the second reaction zone as the inlet parameters of the expansion process, calculating the axial distribution of the airflow state parameters within the nozzle based on the expansion process model, and obtaining the airflow state parameters at the nozzle outlet section includes: Ignoring the addition of reactants and the heat released by the reaction, and only considering the acceleration effect of cross-sectional area change and wall friction on the airflow, the mass conservation equation, momentum conservation equation, energy conservation equation and gas state equation are combined and derived into a set of first-order ordinary differential equations concerning Mach number, pressure and temperature. Using the airflow state parameters at the outlet of the second reaction zone as the initial values for integration, the fourth-order Runge-Kutta method is used to integrate along the axial direction to obtain the airflow Mach number, total pressure, static pressure, total temperature, static temperature and velocity at each integration position in the tail nozzle.
10. The performance evaluation method for a solid scramjet power unit based on primary combustion effect according to claim 1, characterized in that, The method further includes: During multiple executions, the reactant composition or reaction pressure in the first reaction zone is changed; Based on the modified reactant composition and reaction pressure, the equilibrium temperature and product composition of the first reaction zone are recalculated based on the thermochemical equilibrium principle, and the product composition and equilibrium temperature of the first reaction zone and the overall performance parameters of the power unit are obtained for each execution. Based on the differences between the overall performance parameters obtained in each execution, the correspondence between the changes in reactant composition or reaction pressure in the first reaction zone and the overall performance of the power unit is determined.