Numerical simulation of premixed hydrogen flame stabilized by heat transfer and differential diffusion blunt body

Abstract

The application discloses a numerical calculation method of bluff body stabilized premixed hydrogen flame considering heat transfer and differential diffusion, and belongs to the technical field of turbulent flame numerical calculation. The method is suitable for high-precision calculation of bluff body stabilized premixed hydrogen flame by coupling large eddy simulation, small flame model and artificial thickening flame model, and simultaneously processing multi-component differential diffusion effect and heat transfer effect. The method comprises the following steps: one-dimensional free propagation premixed small flame model calculation, three-dimensional small flame lookup table construction, large eddy simulation transport equation establishment containing differential diffusion correction, flame lookup table table lookup and interpolation, large eddy simulation solving and result output. Compared with the existing numerical method which does not consider or simplifies the differential diffusion and heat transfer effect, the application can more accurately predict the shape, stable position and temperature field distribution of the bluff body stabilized premixed hydrogen flame, and improve the design and optimization ability of the hydrogen combustion system under complex working conditions.

Description

Technical Field

[0001] This invention belongs to the field of numerical calculation technology of turbulent flames, specifically involving a numerical calculation method for bluff body stable combustion premixed hydrogen flames that takes into account heat transfer and differential diffusion. Background Technology

[0002] With the increasing application of hydrogen fuel in combustion devices such as gas turbines and industrial furnaces, bluff-body stabilized premixed hydrogen flames are widely used due to their excellent stability. However, hydrogen has extremely strong differential diffusion characteristics, which easily induce thermal diffusion instability under lean-burn conditions, manifesting as highly wrinkled flame surfaces and significant cellular structures, which have a significant impact on flame stability and combustion efficiency. Simultaneously, in actual bluff-body burners, the flame inevitably transfers heat to the bluff body wall and surrounding structures, leading to changes in the combustion state. Traditional small-flame models based on adiabatic assumptions cannot accurately reflect the impact of non-adiabatic effects on flame structure and stability.

[0003] Existing numerical methods for turbulent flames commonly employ the unified Lewis number assumption or neglect differential diffusion, using only a two-dimensional small flame lookup table with mixing fraction and progress variables as trajectory variables, and typically constructing the flame manifold under adiabatic conditions. While these methods offer high computational efficiency, they often struggle to accurately predict flame anchorage, temperature field distribution, and local reaction intensity for fuels with strong differential diffusion, such as hydrogen, and for bluff-body stable flames with significant wall heat transfer, thus impacting the design and safety assessment of hydrogen burners. Therefore, a high-precision numerical calculation method that can simultaneously consider differential diffusion and heat transfer effects is urgently needed. Summary of the Invention

[0004] The purpose of this invention is to provide a numerical calculation method for bluff body stable premixed hydrogen flame that takes into account heat transfer and differential diffusion. By introducing differential diffusion effect during the small flame manifold construction stage and adding total enthalpy coordinates to the manifold coordinates, the heat transfer process between the flame and the bluff body wall and the environment can be accurately described in large eddy simulation, thereby overcoming the problem of insufficient prediction accuracy caused by neglecting or simplifying differential diffusion and non-adiabatic effects in existing methods.

[0005] The technical solution of the present invention is as follows:

[0006] A numerical calculation method for bluff-body stabilized premixed hydrogen flames considering heat transfer and differential diffusion includes the following steps:

[0007] Step 1: Calculate thermochemical parameters using the one-dimensional premixed small flame control equation, map them to the mixing fraction Z, progress variable C, and total enthalpy He to form a three-dimensional manifold, thus forming a three-dimensional flame manifold;

[0008] Step 2: Introduce differential diffusion correction terms and artificially thickened flame models into the large eddy simulation framework to establish modified transport equations;

[0009] Step 3: Obtain the reaction source term and physical property parameters based on real-time trajectory variable values ​​through interpolation, iteratively update the flow field, and dynamically solve by looking up tables;

[0010] Step 4: Solve the governing equations using the pressure coupling algorithm;

[0011] Step 5: Output flame shape, temperature field and component distribution, and output combustion characteristics.

[0012] Compared with existing numerical methods based on the unified Lewis number assumption or adiabatic small flame manifolds, the present invention has the following advantages:

[0013] (1) By employing detailed chemical reaction mechanisms and average transport models of mixtures in the solution of one-dimensional free propagation small flames, this invention considers the differential diffusion effects of various species in the small flame manifold construction stage, which can more accurately describe the flame structure characteristics corresponding to the thermal diffusion instability of hydrogen.

[0014] (2) By changing the initial temperature of the premixed gas and introducing the total enthalpy trajectory variable to construct a small flame lookup table, the non-adiabatic effect caused by heat transfer on the bluff body wall and environmental heat exchange can be reflected in the large eddy simulation, thereby improving the prediction accuracy of the temperature field distribution of the actual bluff body stable combustion premixed hydrogen flame.

[0015] (3) A three-dimensional small flame lookup table with mixing fraction Z, progress variable C and total enthalpy He as independent variables is adopted. While maintaining high computational efficiency, it also takes into account computational accuracy and is suitable for the design, optimization and safety assessment of hydrogen burners on an engineering scale. Attached Figure Description

[0016] Figure 1 This is a flowchart of the numerical calculation method for a bluff body stable premixed hydrogen flame considering heat transfer and differential diffusion according to the present invention.

[0017] Figure 2 This is a schematic diagram of the computational domain geometric model and boundary conditions of the bluff body stable combustion premixed hydrogen burner in an embodiment of the present invention;

[0018] Figure 3 This is a cloud map comparing the time-averaged velocity field of numerical calculation results and experimental observation data in an embodiment of the present invention.

[0019] Figure 4 This is a cloud map comparing the morphology of the combustion reaction zone with the numerical calculation results and experimental observation data of the embodiments of the present invention;

[0020] Figure 5 A comparison diagram (I) of the radial distribution of time-averaged velocity at different flow directions based on numerical calculations and experimental measurements in an embodiment of the present invention.

[0021] Figure 6This is a comparison diagram (II) of the radial distribution of time-averaged velocity at different flow directions, based on numerical calculations and experimental measurements in an embodiment of the present invention.

[0022] Figure 7 This is a comparison diagram (III) of the radial distribution of time-averaged velocity at different flow directions, based on numerical calculations and experimental measurements in an embodiment of the present invention.

[0023] Figure 8 This is a comparison diagram (IV) of the radial distribution of time-averaged velocity at different flow directions, based on numerical calculations and experimental measurements in an embodiment of the present invention.

[0024] Figure 9 This is a comparison chart of the distribution of differential diffusion parameters calculated using the method of this invention and the distribution calculated using the differential diffusion neglect model;

[0025] Figure 10 This is a comparison chart of the component distributions calculated using the method of this invention and those calculated using the neglect of differential diffusion model;

[0026] Figure 11 This is a comparison diagram of the temperature distribution along the flow center axis calculated using the method of this invention and the temperature distribution calculated using an adiabatic model.

[0027] Figure 12 This is a comparison diagram of the near-wall component distribution calculated using the method of this invention and the distribution calculated using an adiabatic model. Detailed Implementation

[0028] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments. However, the following embodiments are only for explaining the present invention, and the scope of protection of the present invention should include all the contents of the claims. Moreover, through the description of the following embodiments, those skilled in the art can fully implement all the contents of the claims of the present invention.

[0029] Example 1:

[0030] This embodiment discloses a numerical calculation method for bluff body stabilized premixed hydrogen flames that takes into account heat transfer and differential diffusion, including the following steps:

[0031] Step 1: Calculation of a one-dimensional free-propagation premixed small flame model. Using a detailed chemical reaction mechanism and a mixed average transport model, the governing equations of the one-dimensional premixed small flame under different initial temperatures and equivalence ratios are solved to obtain thermochemical parameters such as the mass fraction, temperature, and density of each component along the normal coordinate of the flame surface. The small flame level includes multi-component differential diffusion and non-adiabatic effects.

[0032] Step 2: Construct a three-dimensional small flame lookup table. Based on the mixing fraction Z, progress variable C, and total enthalpy He calculated in Step 1, demap the one-dimensional small flame to a three-dimensional small flame lookup table with mixing fraction Z, progress variable C, and total enthalpy He as trajectory variables. Reconstruct parameters such as thermodynamic quantities, component mass fractions, reaction source terms, laminar flame velocity, and flame thickness in the space of mixing fraction Z, progress variable C, and total enthalpy He through sorting and interpolation to obtain a three-dimensional small flame lookup table for large eddy simulation lookup.

[0033] Step 3: Establish the large eddy simulation transport equations including differential diffusion corrections. Within the large eddy simulation framework of a bluff-body stabilized premixed flame, establish filtered continuity equations, momentum conservation equations, and pressure equations. For the trajectory variable mixing fraction Z, progress variable C, and total enthalpy He, corresponding transport equations are established respectively; the differential diffusion effect is introduced as a correction term into the transport equation for the mixing fraction Z. Simultaneously, an artificially thickened flame model is coupled, and the diffusion and reaction source terms in the transport equations are corrected by introducing a thickening factor F and an efficiency function E to achieve analytical analysis of the flame structure at the large eddy simulation grid scale; and the subgrid stress in the momentum equation is closed using a subgrid turbulence model.

[0034] Step 4: Flame lookup table lookup and interpolation. During the time progression of the large eddy simulation, for each time step and each control volume, using the current mixing fraction Z, progress variable C, and total enthalpy He value as lookup coordinates, local chemical reaction source terms, laminar flame parameters, thermophysical parameters, and other parameters are obtained by interpolation from the small flame lookup table constructed in Step 2. These parameters are then substituted into the control equations for the mixing fraction Z, progress variable C, and total enthalpy He to achieve real-time updates of the thermochemical state, turbulent combustion interaction, and differential diffusion effect.

[0035] Step 5: Large eddy simulation solution and result output. A pressure-based transient fluid solution algorithm, combined with an implicit time discretization scheme, is used to iteratively solve the governing equations coupled with detailed chemistry and differential diffusion effects. During the calculation, the time step is dynamically adjusted according to the Courant number to ensure computational stability. After the flow field passes through the initial transient and reaches a statistical quasi-steady state, the data statistics function is activated to output the instantaneous and statistically averaged velocity field, temperature field, component concentration field, and trajectory variables (mixing fraction Z, progress variable C, total enthalpy He) fields of the bluff body stable combustion premixed hydrogen flame, for subsequent combustion characteristic and differential diffusion effect analysis.

[0036] In this embodiment, the progress variable can monotonically characterize the degree of reaction advancement, and is preferably defined as the mass fraction of water, or as a linear combination of the mass fractions of multiple key components.

[0037] In this embodiment, the mixing fraction is used to characterize the local compositional state changes of the premixed system under differential diffusion, and satisfies the preset normalization definition and value range constraints.

[0038] In this embodiment, without changing the mixing fraction as a manifold coordinate used to characterize the composition state, the mixing fraction can also be obtained using other equivalent calculation methods.

[0039] In this embodiment, the total enthalpy is used to characterize the thermal state change caused by non-adiabatic effects, and is used as one of the coordinates in the three-dimensional small flame lookup table for lookup interpolation and state reconstruction.

[0040] In this embodiment, the mixed average transport model includes a modified velocity term to ensure that the diffusion flux of the multi-components meets the overall mass continuity, and the mixed average diffusion coefficient is calculated based on the binary diffusion coefficient.

[0041] In this embodiment, the detailed chemical reaction mechanism in step 1 refers to a mathematical model that provides a detailed description of combustion chemical kinetics based on a large number of elementary reactions, which is different from a single-step chemical reaction mechanism that only describes the overall changes of reactants and products.

[0042] In this embodiment, the equations for the hybrid average transport model in step 1 are as follows:

[0043] The diffusion flux ρYiVi of component i is calculated using the following formula:

[0044] ,

[0045] Among them, D i,m Let m be the average diffusion coefficient of component i, and m represent the mean. The formula is as follows:

[0046] ,

[0047] u cor To correct the speed and ensure continuity of quality, the calculation formula is as follows:

[0048] ,

[0049] In the formula, ρ is the density of the mixed gas; V is the average diffusion coefficient of component k, where m represents the mean; i Y is the diffusion coefficient of component i; i and X i These represent the mass fraction and mole fraction of component i, respectively; x is the spatial coordinate; D ji N is the binary diffusion coefficient between component i and component j; s This represents the total number of components.

[0050] Furthermore, the one-dimensional free-propagation premixed small flame equation in step 1 can be solved using open-source software such as Cantera and FlameMaster. The one-dimensional free-propagation premixed small flame equation is as follows:

[0051] mass conservation equation:

[0052] ,

[0053] Energy conservation equation:

[0054] ,

[0055] In the equation, x represents the spatial coordinates of a one-dimensional flame, ρ is the density, u is the flow velocity, and Y... i and ρY represents the mass fraction of component i and the chemical reaction rate, respectively. i V i That is the diffusion flux of component i, c p λ is the specific heat capacity of the mixture, T represents temperature, and λ is thermal conductivity. This represents the heat release rate, and in the last term of the energy equation, c p,k This refers to the specific heat capacity of component k.

[0056] In this embodiment, step 1, by solving the governing equations under different initial temperature conditions, effectively characterizes the impact of non-adiabatic effects on the flame structure. Meanwhile, given the significant differential diffusion characteristics of hydrogen, even under macroscopic premixing conditions, the local equivalence ratio in the microscopic reaction zone will still deviate from the initial set value due to differences in component diffusion rates. Therefore, it is necessary to solve the small flame equations under a wide range of equivalence ratios to fully cover the initial premixed equivalence ratio and the range of component changes caused by differential diffusion and stratification.

[0057] In this embodiment, the progress variable C in step 2 should be able to monotonically characterize the degree of reaction advancement, preferably defined as the mass fraction of water, or as a linear combination of the mass fractions of multiple key components; wherein the mass fraction of each component is obtained from the one-dimensional small flame calculation results in step 1. The mixing fraction Z in step 2 is used to characterize the degree of mixing between fuel and oxidizer, and is preferably calculated using the following formula:

[0058] ,

[0059] In the formula, Z is the mixing fraction; ν is the stoichiometric ratio; Y F and Y O These represent the mass fractions of fuel and oxidizer in the mixture, respectively; Y O,2 Y represents the mass fraction of the oxidant in the oxidant stream (e.g., air); F,1Z represents the mass fraction of fuel in the fuel stream. Optionally, provided that the normalized definition and range constraints of the mixing fraction are met, the mixing fraction Z can also be obtained using other equivalent calculation methods. The total enthalpy He mentioned in step 2 is used to characterize the change in thermal state caused by non-adiabatic effects, and can be calculated by weighted summation of the mass fraction and absolute enthalpy of each component, as expressed by:

[0060] ,

[0061] In the formula, He(T) j The initial temperature is T. j The total enthalpy of the mixture, Y j h represents the mass fraction of the j-th component. a,j It is the absolute enthalpy of the j-th component at the stated temperature, and n is the total number of components.

[0062] In this embodiment, step 2 is based on the one-dimensional small flame solution calculated in step 1 during the preprocessing procedure. The one-dimensional small flame solution is remapped from the original "physical coordinate-equivalence ratio-boundary temperature" space to an orthogonal "mixing fraction Z-progress variable C-total enthalpy He" phase space. During this process, linear interpolation and collapse algorithms are used to handle the alignment of the data grid.

[0063] In this embodiment, in step 2, for the mixing fraction region that exceeds the combustible limit of premixed combustion, the thermochemical state is supplemented by the inert mixing assumption or the linear extrapolation method based on chemical equilibrium, and the temperature, density and transport properties of the region are recalculated using a detailed mechanism library to ensure the continuity and integrity of the three-dimensional lookup table across the entire flow field.

[0064] In this embodiment, the continuity equation, momentum equation, and pressure equation established in step 3 within the framework of large eddy simulation are as follows:

[0065] ,

[0066] ,

[0067] ,

[0068] In the formula, t is time, x i and x j For Cartesian space coordinate components; ρ, p, and u i These represent the components of the density, pressure, and velocity vectors in the i-direction, respectively. This represents the molecular viscous stress tensor. The superscripts "¯" and "˜" denote spatial filtering and Favre filtering, respectively. The subgrid-scale (SGS) stress tensor reflects the influence of small-scale turbulent fluctuations filtered out by the grid on the large-scale flow field, and is defined as follows: ψ is the compressibility coefficient, and A p H represents the central coefficients after discretization of the momentum equation. j It includes neighbor contributions and source terms, excluding pressure gradients.

[0069] In this embodiment, the trajectory variable transport equations for the mixing fraction Z, progress variable C, and total enthalpy He established in step 3 within the framework of large eddy simulation combined with an artificially thickened flame model are as follows:

[0070] ,

[0071] ,

[0072] ,

[0073] In the formula, ρ is density, u j Ω represents the velocity component; E is the efficiency function of the artificially thickened flame model; F is the thickening factor; and Ω is the flame sensor. t The subgrid turbulent diffusion coefficient is... D is the source term for the schedule variable. C Let C be the mass diffusion coefficient, and D be the mass diffusion coefficient of the schedule variable. He D is the thermal diffusivity. Z The diffusion coefficient is the fractional diffusion coefficient. Let D be the cross-diffusion coefficient. Z and The definition is as follows:

[0074] ,

[0075] ,

[0076] In the formula, ν is the stoichiometric ratio; Y F,1 Y represents the mass fraction of fuel in the fuel stream. O,2 D represents the mass fraction of oxidant in the oxidant stream. F and D O These are the mixed average molecular diffusion coefficients of the fuel and oxidant, respectively, which were obtained through detailed chemical mechanism calculations.

[0077] In this embodiment, the subgrid turbulence model mentioned in step 3 can be any one of the Smagorinsky model, dynamic Smagorinsky model, wall adaptive local eddy viscosity model, etc., used to calculate the subgrid eddy viscosity and the subgrid stress term in the closed momentum conservation equation.

[0078] In this embodiment, step 4 is implemented based on a self-developed large eddy simulation solver. It constructs a dynamically coupled interface to call the lookup table and supports flexible setting of the lookup table update frequency to adapt to the calculation requirements under different time steps.

[0079] In this embodiment, the time-progression solution strategy described in step 5 can specifically employ pressure-velocity coupling algorithms such as the pressure implicit operator partitioning and pressure-coupled equation hybrid algorithm (PIMPLE) or the pressure implicit operator partitioning algorithm (PISO). The time discretization scheme can be selected from Euler implicit scheme (first-order accuracy) or Backward difference scheme (second-order accuracy), etc., to seek the best balance between numerical stability, computational efficiency and time analytical accuracy according to actual computational needs.

[0080] Example 2:

[0081] like Figure 1 and Figure 2 As shown, this embodiment takes a bluff burner as an example to illustrate the specific implementation process of the numerical calculation method of the present invention, which mainly includes the following steps:

[0082] Step 1: Determine the simulated operating conditions. In this embodiment, the operating conditions of the burner are set as follows: the equivalence ratio of hydrogen / air premixed gas is 0.4, the ambient pressure is 1 atm, and the initial temperature of the unburned premixed gas and the ambient temperature are 298.15K.

[0083] Step 2: Calculation of the one-dimensional free-propagating small flame model. In this embodiment, the open-source software FlameMaster is used to solve the one-dimensional premixed small flame control equation.

[0084] (1) Mechanism and model selection:

[0085] A detailed chemical reaction mechanism of hydrogen, involving 21 components and 109 elementary reactions, was employed, and a mixed-average transport model was used to calculate the component diffusion flux in order to accurately capture the differential diffusion characteristics of hydrogen.

[0086] (2) Parameter space construction:

[0087] Temperature dimension (considering non-adiabatic effects): In order to characterize the influence of wall heat transfer and recirculation zone heat transfer on flame structure, a series of different initial temperatures of unburned gas were selected for calculation, namely 298 K, 300 K, 400 K, 500 K, 600 K, 700 K, 800 K and 900 K.

[0088] Equivalent ratio dimension (considering differential diffusion effects): Given that the high diffusivity of hydrogen can cause local equivalent ratio deviations from the initial equivalent ratio, the equivalent ratio of the hydrogen / air mixture was calculated within the range of 0.3 to 0.8 under each of the initial temperature conditions described above. This range fully covers the local equivalent ratio fluctuations caused by differential hydrogen diffusion.

[0089] Step 3: Establish a small flame lookup table. Based on the series of one-dimensional small flame solutions obtained in Step 2, a three-dimensional small flame lookup table is constructed using a self-developed preprocessing program. This process first calculates the mixing fraction Z, the progress variable C, and the total enthalpy He. Then, using a manifold mapping algorithm, the one-dimensional solutions in the physical space are reconstructed into a three-dimensional phase space with the mixing fraction Z, the progress variable C (defined as the mass fraction of water), and the total enthalpy He as orthogonal coordinates. Subsequently, for each node in the lookup table, using detailed chemical mechanisms and a mixing-average transport model, the diffusion coefficients of each component and the chemical reaction source terms of the progress variable are explicitly calculated and stored, thereby generating a comprehensive database containing thermodynamic quantities and differential diffusion transport coefficients for subsequent large eddy simulation lookup.

[0090] Step 4: Large Eddy Simulation Solver Setup and Mesh Generation. Using the open-source fluid dynamics platform OpenFOAM, an unsteady solver coupled with the aforementioned lookup table and trajectory variable transport equations was developed, and a computational domain was established for numerical calculations. The computational domain was discretized using a structured mesh. To accurately capture the turbulent structure of the shear layer and flame front, the reaction core region was locally refined, resulting in approximately 7.7 million mesh cells. For boundary conditions, a velocity boundary condition was given at the burner inlet, with the average inlet velocity set to 8.56 m / s. The inlet component concentration and temperature were initialized according to the conditions set in Step 1. No-slip boundary conditions were used on the wall, with the temperature set to 298.15 K. For the turbulence model, the Large Eddy Simulation method was employed, and the turbulent stress term at the closed subgrid scale of the Smagorinsky model was selected to analyze the large-scale turbulent vortex structure and its interaction with the flame.

[0091] Step 5: Numerical Discretization Scheme and Solution Control. The finite volume method is used to discretize and solve the governing equations. For pressure-velocity coupling, the PIMPLE algorithm is employed, performing multiple internal iterations within each time step to ensure computational stability and convergence over large time steps. Second-order implicit backward Euler scheme is used for time discretization, and second-order central difference scheme (Gauss linear) is used for spatial discretization of the convection terms to minimize numerical dissipation. An adaptive strategy is used for time step control, dynamically adjusting the time step by limiting the maximum Coulomb number to no more than 0.9, while setting the maximum physical time step to 1e-05 s to meet numerical stability requirements. The solver is started for iterative calculations, and statistical averaging is activated after the flow field has passed the initial transient and reached a statistical quasi-steady state. The output includes instantaneous and averaged velocity, temperature, component, and trajectory variable (mixing fraction Z, progress variable C, total enthalpy He) distributions. Results analysis is then performed.

[0092] like Figure 3 As shown, the time-averaged velocity field cloud map calculated using the method of this invention is highly consistent with the average experimental results in both morphology and numerical value. Specifically, the numerically calculated height (H) of the inner circulation zone (IRZ) is... IRZ The actual value was 22.44 mm, while the experimentally measured value was 22.10 mm; the two are extremely close, with minimal error. Figure 4 The figure shows a comparison of the combustion reaction zone and velocity streamline distribution between numerical calculations and experimental results. The combustion reaction intensity in the experimental results is characterized by the normalized OH* chemiluminescence intensity, while the numerical simulation results are characterized by the normalized heat release rate (HRR). The results show that the combustion reaction zone is anchored above and surrounds the bluff body. Simultaneously, the velocity streamline distribution indicates the formation of a significant internal circulation zone above the bluff body, which is a key structure for maintaining stable hydrogen flame combustion. In this respect, the numerical calculations accurately reproduce the flow field characteristics observed experimentally. Figure 5 , Figure 6 , Figure 7 and Figure 8 The figures show the time-averaged radial velocity distribution curves at different flow directions. At all selected flow direction locations, the numerically calculated time-averaged velocity distribution trends and magnitudes agree well with the experimental data. In summary, the method of this invention can accurately predict complex flow field structures, verifying the effectiveness and accuracy of the model.

[0093] To further verify the ability of the method of this invention to capture the differential diffusion effect of hydrogen, a small flame lookup table that does not consider differential diffusion was established as Comparative Example 1. In the control group, the Lewis number of all components was set to 1, and the remaining solution steps were consistent with the method of this invention. In this embodiment, the differential diffusion parameter Z was calculated using the time-averaged component mass fraction.HO The calculation formula is as follows:

[0094] ,

[0095] Wherein, ξH and ξO are the mixing fractions calculated based on the time-averaged mass fractions of hydrogen and oxygen, respectively.

[0096] like Figure 9 As shown, in the results calculated using the method of this invention, it is evident that the calculation based on the time-averaged parameter... A significant deviation from 0 indicates a difference in the transport processes of hydrogen and oxygen, meaning the method of this invention successfully captured the differential diffusion effect of hydrogen. In contrast, in the control group where the Lewis number equals 1, The numerical value remains essentially zero, meaning that the diffusion rate of all components is forcibly assumed to be the same. Furthermore, as... Figure 10 As shown, considering differential diffusion, the mass fraction of component H is higher than when the Lewis number is equal to 1, which also proves that component H exhibits local aggregation, caused by differential diffusion. This comparison demonstrates the technical advantage of this invention in handling differential diffusion phenomena.

[0097] To further characterize the effect of non-adiabatic effects in the method of this invention, an adiabatic small flame lookup table was established as Comparative Example 2. In Comparative Example 2, the one-dimensional free-propagating small flame model was calculated only at the initial temperature of a single unburned gas (298.15 K in this example), and the remaining solution steps were consistent with those of this invention.

[0098] like Figure 11 The figure shows the time-averaged temperature distribution curve along the flow centerline, comparing the temperature changes along the flow direction under adiabatic and non-adiabatic conditions. It can be seen that under non-adiabatic conditions, the temperature along the flow direction is slightly higher than under adiabatic conditions. Although it is generally believed that the flame temperature should be lower under non-adiabatic conditions due to wall heat loss, the calculation results of this embodiment reveal a more complex flow field thermodynamic mechanism: during combustion, a high-temperature region is formed in the recirculation zone above the bluff body, which is crucial for maintaining stable hydrogen flame combustion. Due to the existence of this high-temperature recirculation zone, a thermal recirculation effect occurs, altering the initial thermal state of the hydrogen and unburned mixture involved in combustion (i.e., increasing the local total enthalpy). Specifically, the preheating effect of the high-temperature gas recirculation leads to an increase in the local initial temperature before hydrogen combustion. This heating effect numerically exceeds the cooling effect caused by wall heat loss, ultimately leading to an overall increase in the temperature along the flow direction. Under adiabatic conditions (Comparative Example 2), this heating phenomenon was not captured because this complex energy exchange and transfer mechanism was ignored. Figure 12The figure shows the mass fraction distribution of OH near the wall surface. The results indicate that the OH mass fraction is slightly higher under non-adiabatic conditions than under adiabatic conditions, suggesting a more vigorous combustion reaction under non-adiabatic conditions. This microscopic component distribution characteristic is consistent with the analysis results of the macroscopic temperature field, further confirming the above conclusions regarding the heat recirculation mechanism and demonstrating the accuracy of the method of this invention in handling non-adiabatic and complex thermal coupling problems.

[0099] The above description is merely a specific embodiment of this application, enabling those skilled in the art to understand or implement this application. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of this application. Therefore, this application is not to be limited to the embodiments shown herein, but is to be accorded the widest scope consistent with the principles and novel features claimed herein.

Claims

1. A numerical calculation method for bluff-body stable premixed hydrogen flame considering heat transfer and differential diffusion, characterized in that, Includes the following steps: Step 1: Calculate thermochemical parameters using the one-dimensional premixed small flame control equation, map them to the mixing fraction Z, progress variable C, and total enthalpy He to form a three-dimensional manifold, and construct a three-dimensional small flame lookup table; Step 2: Under the framework of large eddy simulation, establish the filtered continuity equation, momentum conservation equation, and pressure equation; For the trajectory variable mixing fraction Z, the progress variable C, and the total enthalpy He, corresponding transport equations are established respectively. Among them, the differential diffusion effect is introduced into the transport equation of mixing fraction Z in the form of a correction term, and coupled with the artificial thickened flame model, the diffusion term and reaction source term in the transport equation are corrected by introducing the thickening factor F and the efficiency function E. Step 3: During the time-progression process of the large eddy simulation, the current mixing fraction Z, progress variable C, and total enthalpy He value are used as lookup coordinates. Local chemical reaction source terms, laminar flame parameters, and thermophysical parameters are obtained by interpolation from the small flame lookup table constructed in Step 1. These parameters are substituted into the control equations of mixing fraction Z, progress variable C, and total enthalpy He to realize the real-time update of thermochemical state, turbulent combustion interaction, and differential diffusion effect. Step 4: The pressure coupling algorithm is used to iteratively solve the coupled control equation set consisting of the filtered continuity equation, momentum conservation equation, pressure equation, and control equations and transport equations for the mixing fraction Z, progress variable C, total enthalpy He. Step 5: Output flame morphology, temperature field, and component distribution; output combustion characteristics. In step 1, the solution of the one-dimensional premixed small flame control equation is remapped from the physical coordinate-equivalence ratio-initial temperature data space to the orthogonal mixing fraction Z-progress variable C-total enthalpy He phase space.

2. The numerical calculation method for bluff body stable premixed hydrogen flame considering heat transfer and differential diffusion according to claim 1, characterized in that, In step 1, the one-dimensional premixed small flame control equation is solved by covering the equivalence ratio range of 0.3 to 0.8 to cover the local equivalence ratio variation caused by differential diffusion.

3. The numerical calculation method for bluff body stable premixed hydrogen flame considering heat transfer and differential diffusion according to claim 1, characterized in that, Linear interpolation and collapse algorithms are used to align and reconstruct the data grid.

4. The numerical calculation method for bluff body stable premixed hydrogen flame considering heat transfer and differential diffusion according to claim 1, characterized in that, In step 1, for coordinate regions that exceed the premixed flammability limit or lookup table boundaries, the thermochemical state is supplemented using the inert mixing assumption or an extrapolation strategy based on chemical equilibrium. The temperature, density, and transport properties of the supplemented regions are then recalculated to ensure the continuity and integrity of the three-dimensional small flame lookup table.

5. The numerical calculation method for bluff body stable premixed hydrogen flame considering heat transfer and differential diffusion according to claim 1, characterized in that, In step 2, the differential diffusion correction term is set in the mixing fraction Z equation to characterize the local composition shift caused by the difference in diffusion rates of different components.

6. The numerical calculation method for bluff body stable premixed hydrogen flame considering heat transfer and differential diffusion according to claim 1, characterized in that, Step 2 employs the subgrid stress in the closed momentum equation of the subgrid turbulence model. The subgrid turbulence model is selected from any one of the Smagorinsky model, the dynamic Smagorinsky model, or the wall adaptive local eddy viscosity model.

7. The numerical calculation method for bluff body stable premixed hydrogen flame considering heat transfer and differential diffusion according to claim 1, characterized in that, In step 3, a dynamic coupling interface is constructed to call the three-dimensional small flame lookup table, and the table update frequency can be set to adapt to the calculation requirements under different time steps.

8. The numerical calculation method for bluff body stable premixed hydrogen flame considering heat transfer and differential diffusion according to claim 1, characterized in that, In step 4, a pressure-based transient fluid solution algorithm is used, combined with an implicit time discretization scheme, to iteratively solve the governing equations.

9. The numerical calculation method for bluff body stable premixed hydrogen flame considering heat transfer and differential diffusion according to claim 8, characterized in that, The implicit time discrete scheme can be either the implicit Euler scheme or the second-order Backward scheme.

Images