From CFD-Wiki

Jump to: navigation, search

The power of Fire, or Flame, for instance, which we designate by some trivial chemical name, thereby hiding from ourselves the essential character of wonder that dwells in it as in all things, is with these old Northmen, Loke, a most swift subtle Demon of the brood of the J\"otuns... From us too no Chemistry, if it had not Stupidity to help it, would hide that Flame is a wonder. What is Flame?

Carlyle on Heroes Odin and Scandinavian Mythology., [1].

What is combustion -- physics versus modelling

Combustion phenomena consist of many physical and chemical processes which exhibit a broad range of time and length scales. A mathematical description of combustion is not always trivial, although some analytical solutions exist for simple situations of laminar flame. Such analytical models are usually restricted to problems in zero or one-dimensional space.

Fundamental Aspects

Main Specificities of Combustion Chemistry

Combustion can be split into two processes interacting with each other: thermal, and chemical.

The chemistry is highly exothermal (this is the reason of its use) but also highly temperature dependent, thus highly self-accelerating. In a simplified form, combustion can be represented by a single irreversible reaction involving 'a' fuel and 'an' oxidizer:

$\frac{\nu_F}{\bar M_F}Y_F + \frac{\nu_O}{\bar M_O} Y_O \rightarrow Product + Heat$

Althgough very simplified compared to real chemistry involving hundreds of species (and their individual transport properties) and elemental reactions, this rudimentary chemistry has been the cornerstone of combustion analysis and modelling.

The most widely used form for the rate of the above reaction is the Arrhénius law:

$\dot\omega = \rho^2 A \left ( \frac{Y_F}{\bar M_F} \right )^{n_F} \left (\frac{Y_O}{\bar M_O}\right )^{n_O} \exp^{-T_a/T}$

$T_a$ is the activation temperature, high in combustion, consistently with the temperature dependence. This is where the high non-linearity in temperature is modelled. A is the pre-exponential constant. One of the interpretation of the Arrhénius law comes from gas kinetic theory: the number of molecules whose kinetic energy is larger than the minimum value allowing a collision to be energetic enough to trig a reaction is proportional to the exponential term introduced above divided by the square root of the temperature. This interpretation allows one to think that the temperature dependence of A is very weak compared to the exponential term. A is eventually considered as constant. The reaction rate is also naturally proportional to the molecular density of each of the reactant. Nonetheless, the orders of reaction $n_i$ are different from the stoichiometric coefficients as the single-step reaction is global, not governed by collision for it represents hundreds of elementary reactions. If one goes into the details, combustion chemistry is based on chain reactions, decomposed into three main steps: (i) generation (where radicals are created from the fresh mixture), (ii) branching (where products and new radicals appear from interaction of radicals with reactants), and (iii) termination (where radicals collide and turn into products). The branching step tends to accelerate the production of active radicals (autocatalytic). The impact is nevertheless small compared to the high non-linearity in temperature. This explains why single-step chemistry has been sufficient for most of the combustion modelling work up to now.

The fact that a flame is a very thin reaction zone separating, and making the transition between, a frozen mixture and an equilibrium is explained by the high temperature dependence of the reaction term, modelled by a large activation temperature, and a large heat release (the ratio of the burned and fresh gas temperatures is about 7 for typical hydrocarbon flames) leading to a sharp self-acceleration in a very narrow area. To evaluate the order of magnitude of the quantities, the terms in the exponential argument are normalized:

$\beta=\alpha\frac{T_a}{T_s} \qquad \alpha=\frac{T_s-T_f}{T_s}$

$\beta$ is named the Zeldovitch number and $\alpha$ the heat release factor. Here, $T_s$ has been used instead of $T_b$ , the conventional notation for burned gas temperature (at final equilibrium). $T_s$ is actually $T_b$ for a mixture at stoichiometry and when the system is adiabatic, i.e. this is the reference highest temperature that can be obtained in the system. $T_f$ is the ambient temperature of the fresh gases. That said, typical value for $\beta$ and $\alpha$ are 10 and 0.9, giving a good taste of the level of non-linearity of the combustion process with respect to temperature. Actually, the reaction rate is rewritten as:

$\dot\omega = \rho^2 A \left ( \frac{Y_F}{\bar M_F} \right )^{n_F} \left (\frac{Y_O}{\bar M_O}\right )^{n_O} \exp^{-\frac{\beta}{\alpha}}\exp^{-\beta\frac{1-\theta}{1-\alpha(1-\theta)}}$

where the non-dimensionalized temperature is:

$\theta=\frac{T-T_f}{T_s-T_f}$

The non-linearity of the reaction rate is seen from the exponential term:

${\mathcal O}(\exp^{-\beta})$ for $\theta$ far from unity (in the fresh gas)
${\mathcal O}(1)$ for $\theta$ close to unity (in the reaction zone close to the burned gas whose temperature must be close to the adiabatic one $T_s$ ), more exactly $1-\theta \sim {\mathcal O}(\beta^{-1})$

Temperature Non-Linearity of the Source Term: the Temperature-Dependent Factor of the Reaction Term for Some Values of the Zeldovitch and Heat Release Parameters

Note that for an infinitely high activation energy, the reaction rate is piloted by a $\delta(\theta)$ function. The figure, beside, illustrates how common values of $\beta$ around 10 tend to make the reaction rate singular around $\theta$ of unity. Two set of values are presented: $\beta = 10$ and $\beta = 8$ . The first magnitude is the representative value while the second one is a smoother one usually used to ease numerical simulations. In the same way, two values for the heat release $\alpha$ , 0.9 and 0.75, are explored. The heat release is seen to have a minor impact on the temperature non-linearity.

Transport Equations

Additionally to the Navier-Stokes equations, at least with variable density, the transport equations for a reacting flow are the energy and species transport equations. In usual notations, the specie i transport equation is written as:

$\frac{D \rho Y_i}{D t} = \nabla\cdot \rho D_i\vec\nabla Y_i - \nu_i\bar M_i\dot\omega$

and the temperature transport equation:

$\frac{D\rho C_p T}{Dt} = \nabla\cdot \lambda\vec\nabla T + Q\nu_F \bar M_F \dot\omega$

The diffusion is modelled thanks to Fick's law that is a (usually good) approximation to the rigorous diffusion velocity calculation. Regarding the temperature transport equation, it is derived from the energy transport equation under the assumption of a low-Mach number flow (compressibility and viscous heating neglected). The low-Mach number approximation is suitable for the deflagration regime (as it will be demonstrated below), which is the main focus of combustion modelling. Hence, the transport equation for temperature, as a simplified version of the energy transport equation, is usually retained for the study of combustion and its modelling.

Note: The species and temperature equations are not closed as the fields of velocity and density also need to be computed. Through intense heat release in a very small area (the jump in temperature in typical hydrocarbon flames is about seven and so is the drop in density in this isobaric process, and the thickness of a flame is of the order of the millimetre), combustion influences the flow field. Nevertheless, the vast majority of combustion modelling has been developed based on the species and temperature equations, assuming simple flow fields, sometimes including hydrodynamics perturbations.

Low-Mach Number Equations

In compressible flows, when the motion of the fluid is not negligible compared to the speed of sound (which is the speed at which the molecules can reorganize themselves), the heap of molecules results in a local increase of pressure and temperature moving as an acoustic wave. It means that, in such a system, a proper reference velocity is the speed of sound and a proper pressure reference is the kinetic pressure. A contrario, in low-Mach number flows, the reference speed is the natural representative speed of the flow and the reference pressure is the thermodynamic pressure. Hence, the set of reference quantities to characterize a low-Mach number flow is given in the table below:

Density $\rho_o$ A reference density (upstream, average, etc.)

Velocity $U_o$ A reference velocity (inlet average, etc.)

Temperature $T_o$ A reference temperature (upstream, average, etc.)

Pressure (static) $P_o=\rho_o \bar r T_o$ From Boyle-Mariotte

Length $L_o$ A reference length (representative of the domain)

Time $L_o/U_o$

Energy $C_p T_o$ Internal energy at constant reference pressure. $C_p$ must also be chosen in a reference thermodynamical state.

The equations for fluid mechanics properly adimensionalized can be written:

Mass conservation:

$\frac{D\rho}{Dt} =0$

Momentum:

$\frac{D\rho\vec U}{Dt}=-\frac{1}{\gamma M}\vec\nabla P+\nabla\cdot\frac{1}{Re}\bar\bar\Sigma$

Total energy:

$\frac{D\rho e_T}{Dt}-\frac{1}{RePr}\nabla\cdot\lambda\vec\nabla T=-\frac{\gamma-1}{\gamma}\nabla\cdot P\vec U +\frac{1}{Re}M^2(\gamma-1)\nabla\cdot \vec U\bar\bar\Sigma + \frac{\rho}{C_pT_o(U_o/L_o)}Q\nu_F \bar M_F\dot\omega$

Specie:

$\frac{D\rho Y}{Dt}-\frac{1}{ScRe}\nabla\cdot\rho D\vec\nabla Y=-\frac{\rho}{U_o/L_o}\nu\bar M\dot\omega$

State law:

$P=\rho T$

The low-Mach number equations are obtained considering that $M^2$ is small. 0.1 is usually taken as the limit, which recovers the value of a Mach number of 0.3 to characterize the incompressible regime.

Considering the energy equation, in addition to the terms with $M^2$ in factor in the equation, the total energy reduces to internal energy as: $e_T = T/\gamma+M^2(\gamma-1)\vec U^2/2$ . Moreover, the work of pressure is considered as negligible because the gradient of pressure is negligible (low-Mach number approximation is indeed also named isobaric approximation) and the flow is assumed close to a divergence-free state. For the same reason, volumic energy and enthalpy variations are assumed equal as they only differ through the addition of pressure. Hence, redimensionalized, the low-Mach number energy equation leads to the temperature equation as used in combustion analysis:

$\frac{D\rho C_p T}{Dt} = \nabla\cdot \lambda\vec\nabla T + Q\nu_F \bar M_F \dot\omega$

The Damköhler Number

A flame is a reaction zone. From this simple point of view, two aspects have to be considered: (i) the rate at which it is fed by reactants, let call $\tau_d$ the characteristic time, and the strength of the chemistry to consume them, let call the characteristic chemical time $\tau_c$ . In combustion, the Damköhler number, Da, compares these both time scales and, for that reason, it is one of the most integral non-dimensional groups:

$Da=\frac{\tau_d}{\tau_c}$ .

If Da is large, it means that the chemistry has always the time to fully consume the fresh mixture and turn it into equilibrium. Real flames are usually close to this state. The characteristic reaction time, $(Ae^{-T_a/T_s})^{-1}$ , is estimated of the order of the tenth of a ms. When Da is low, the fresh mixture cannot be converted by a too weak chemistry. The flow remains frozen. This situation happens with ignition or misfire, for instance.

The picture of a deflagration lends itself to a description based on the Damköhler number. A reacting wave progresses towards the fresh mixture through preheating of the upstream closest layer. The elevation of the temperature strengthens the chemistry and reduces its characteristic time such that the mixture changes from a low-Da region (far upstream, frozen) to a high-Da region in the flame (intense reaction to equilibrium).

Conservation Laws

The processus of combustion transforms the chemical enthalpy into sensible enthalpy (i.e. rise the temperature of the gases thanks to the heat released). Simple relations can be drawn between species and temperature by studying the source terms appearing in the above equations:

$\frac{Y_F}{\nu_F\bar M_F} - \frac{Y_O}{\nu_O\bar M_O} = \frac{Y_{F,u}}{\nu_F\bar M_F} - \frac{Y_{O,u}}{\nu_O\bar M_O}$

$Y_F + \frac{Cp T}{Q} = Y_{F,u} + \frac{CpT_u}{Q}$

Those coupling functions are named Schwalb-Zeldovich variables.

Hence $T_b = T_u + \frac{Q Y_{F,u}}{Cp}$ , $Y_{O,b} = Y_{O,u} - \frac{\nu_O \bar M_O}{\nu_F \bar M_F} Y_{F,u} = Y_{O,u} -sY_{F,u}$ and $Y_{F,b} = 0$ . Here, the example has been taken for a lean case.

As mentioned in Sec. Main Specificities, the stoichiometric state is used to non-dimensionalize the conservation equations:

$Y_i^* = Y_{i,u}^* - \theta\frac{Y_{F,s}}{Y_{i,s}}s^{\delta_{i,O}} \qquad ; \qquad Y_i^*=Y_i/Y_{i,s}$ .

A comprehensive form of the reaction rate can be reconstituted to understand the difficulty of numerically resolving the reaction zone:

$\dot\omega = B \prod_{i=O,F}(Y_{i,u}^*-\theta\frac{Y_{F,s}}{Y_{i,s}}s^{\delta_{i,O}} )^{n_i} \exp{\left ( -\beta\frac{1-\theta}{1-\alpha(1-\theta)}\right)}$

where $B$ stands for all the constant terms present in this reaction rate, plus density.

Source Term versus Temperature

For the stoichiometric case and a global order of two, the reaction rate is graphed versus the reduced temperature for different values of the heat release and Zeldovitch parameter. A high value of $\beta$ makes the reaction rate very sharp, versus temperature. It means that reaction is significant beyond a temperature level (sometimes called ignition temperature) that is close to one (the exponential term above is non-negligible for $1-\theta \sim \beta^{-1}$ ). The heat release has qualitatively the same impact but not so strong. Transposed to the case of a flame sheet, it effectively shows that the reaction exists only in a fraction of the thermal thickness of the flame (the region close to the flame that the latter preheats, hence, where the reduced temperature rises from 0 to 1 here) where the temperature deviates few from the maximal one (density can be assumed as constant and equal to its burned-gas value). Numerically capturing such a sharp reaction zone can be costly and the lower values of $\beta$ and $\alpha$ as presented here are usually preferred whenever possible.

Most problems in combustion involve turbulent flows, gas and liquid fuels, and pollution transport issues (products of combustion as well as for example noise pollution). These problems require not only extensive experimental work, but also numerical modelling. All combustion models must be validated against the experiments as each one has its own drawbacks and limits. In this article, we will address the modeling fundamentals only.

The combustion models are often classified on their capability to deal with the different combustion regimes.

Three Combustion Regimes

Depending on how fuel and oxidizer are brought into contact in the combustion system, different combustion modes or regimes are identified. Traditionally, two regimes have been recognized: the premixed regime and the non-premixed regime. Over the last two decades, a third regime, sometime considered as a hybrid of the two former ones to a certain extend, has risen. It has been named partially-premixed regime.

The Non-Premixed Regime

Sketch of a diffusion flame

This regime is certainly the easiest to understand. Everybody has already seen a lighter, candle or gas-powered stove. Basically, the fuel issues from a nozzle or a simple duct into the atmosphere. The combustion reaction is the oxidization of the fuel. Because fuel and oxidizer are in contact only in a limited region but are separated elsewhere (especially in the feeding system) this configuration is the safest. The non-premixed flame is also named diffusion flame as diffusion of oxidant and fuel has to occur simultaneously to reaction to sustain combustion, the flame being a surface of separation of fuel and oxidant streams. The non-premixed flame has some other advantages. By controlling the flows of both reactants, it is (theoretically) possible to locate the stoichiometric interface, and thus, the location of the flame sheet. Moreover, the strength of the flame can also be controlled through the same process. Depending on the width of the transition region from the oxidizer to the fuel side, the species (fuel and oxidizer) feed the flame at different rates. This is because the diffusion of the species is directly dependent on the unbalance (gradient) of their distribution. A sharp transition from fuel to oxidizer creates intense diffusion of those species towards the flame, increasing its burning rate. As said above, this burning rate control through the diffusion process is certainly one of the reasons of the alternate name of such a flame and combustion mode: diffusion flame and diffusion regime.

Because a diffusion flame is fully determined by the inter-penetration of the fuel and oxidizer streams, it is convenient to introduce a tracer of the state of the mixture. This is the role of the mixture fraction, usually called Z or f. Z is usually taken as unity in the fuel stream and is null in the oxidizer stream. It varies linearly between this two bounds such that at any point of a frozen flow the fuel mass fraction is given by $Y_F=ZY_{F,o}$ and the oxidizer mass fraction by $Y_O = (1-Z)Y_{O,o}$ . $Y_{F,o}$ and $Y_{O,o}$ are the fuel and oxidizer mass fractions in the fuel and oxidizer streams, respectively. The mixture fraction posses a transport equation that is expected to not have any source term as a tracer of a mixture must be a conserved scalar. First, the fuel and oxidizer mass fraction transport equations are written in usual notations:

$\frac{D \rho Y_F}{D t} = \nabla\cdot \rho D\vec\nabla Y_F - \nu_F\bar M_F\dot\omega$

$\frac{D \rho Y_O}{D t} = \nabla\cdot \rho D\vec\nabla Y_O - \nu_O\bar M_O\dot\omega$

The two above equations are linearly combined in a single one in a manner that the source term disappears:

$\frac{D \rho (\nu_O\bar M_O Y_F-\nu_F\bar M_F Y_O)}{D t} = \nabla\cdot \rho D\vec\nabla (\nu_O\bar M_O Y_F-\nu_F\bar M_F Y_O)$

The quantity $(\nu_O\bar M_O Y_F-\nu_F\bar M_F Y_O)$ is thus a conserved scalar, one recognizes one of the Schwalb-Zeldovich quantities introduced in Sec. Conservation Laws. One also remarks that this transport equations combination is made in the equidiffusional approximation, i.e. all the scalars, including temperature, diffuse at the same rate. The last step is to normalize it such that it equals unity in the pure fuel stream ( $Y_F=Y_{F,o}$ and $Y_O=0$ ) and is null in the pure oxidizer stream ( $Y_F=0$ and $Y_O=Y_{O,o}$ ). The resulting normalized passive scalar is the mixture fraction:

$Z=\frac{\Phi(Y_F/Y_{F,o})-(Y_O/Y_{O,o})+1}{\Phi+1}\qquad \Phi=\frac{\nu_O\bar M_O Y_{F,o}}{\nu_F\bar M_F Y_{O,o}}$

governed by the transport equation:

$\frac{D \rho Z}{D t} = \nabla\cdot \rho D\vec\nabla Z$

The stoichiometric interface location (and thus the approximate location of the flame if the flow is reacting) is where $\Phi(Y_F/Y_{F,o})-(Y_O/Y_{O,o})$ vanishes (or $Y_F$ and $Y_O$ are both null in the reacting case). This leads to a stoichiometry definition:

$Z_s=\frac{1}{1+\Phi}$

As the mixture fraction qualifies the degree of inter-penetration of fuel and oxidizer, the elements originally present in these molecules are conserved and can be directly traced back to the mixture fraction. This has led to an alternate defintion of the mixture fraction, based on element conservation. First, the elemental mass fraction $X_{j}$ of element j is linked to the species mass fraction $Y_i$ :

$X_j = \sum_{i=1}^{n} \frac{a_{i,j} \bar M_j}{\bar M_i} Y_i$

where $a_{i,j}$ is a matrix counting the number of element j atoms in specie molecule named i and n is the number of species in the mixture. The group pictured by the summation above is a linear combination of $Y_i$ . Because the transport equations of species mass fraction, few lines earlier, are also linear, a transport equation for the elemental mass fraction can be written:

$\frac{D \rho X_j}{D t} = \nabla\cdot \sum_{i=1}^n\frac{a_{i,j}\bar M_j}{\bar M_i}\rho D_i\vec\nabla Y_i$

For mass is conserved, the linear combination of the source terms vanishes. Furthermore, by taking the same diffusion coefficient $D_i$ for all the species, the elemental mass fraction transport equation has exactly the same form as the specie transport equation (except the source term). Notice that the assumption of equal diffusion coefficient was also made in the previous definition of the mixture fraction and is justified in turbulent combustion modelling by the turbulence diffusivity flattening the diffusion process in high Reynolds number flows. Hence, the elemental mass fraction transport equation has the same structure as the mixture fraction transport equation seen above. Properly renormalized to reach unity in the fuel stream and zero in the oxidizer stream, the elemental mass fraction is a convenient way of determining the mixture fraction field in a flow. Indeed, it is widely used in practice for this purpose.

Simplified Diffusion Flame Solution

Because Schwalb-Zeldovitch variables rid off any source, they allow a first approximation of the flame description independent of the kinetics. From the definition of $Z$ , the Schwalb-Zeldovich variables have the following expression:

$\left \{ \begin{array}{lll} Y_O-sY_F & = & Y_{O,o}-Z(sY_{F,o}+Y_{O,o}) \\ \frac{Q Y_F}{C_p} + T & = & T_{O,o} + Z(T_{F,o}-T_{O,o}+\frac{Q Y_{F,o}}{C_p}) \end{array} \right.$

Two simplified limit cases for the diffusion flame equation can be built from those relations:

In the frozen flow, the chemical source term is null everywhere. Hence, $Y_F$ and $Z$ follow

strictly the same linear transport equation. Then, $Y_F = Y_{F,o} Z$ as anticipated above. Follow: $Y_O=Y_{O,o}(1-Z)$ , $T=T_{O,o}+Z(T_{F,o}-T_{O,o})$ .

The equilibrium solution, also with zero source term, the 'historical' diffusion flame description laid out by Burke and Schumann in 1928:
- Rich side: $Y_O = 0$ , $Y_F = -Y_{O,o}/s + Z(Y_{F,o}+Y_{O,o}/s)$ , $T=T_{F,o}+(1-Z)(T_{O,o}-T_{F,o}+Y_{O,o}\frac{Q}{sC_p})$
- Lean side: $Y_F = 0$ , $Y_O = Y_{O,o} - Z(sY_{F,o}+Y_{O,o})$ , $T=T_{O,o}+Z(T_{F,o}-T_{O,o}+Y_{F,o}\frac{Q}{C_p})$

From the value of $Z$ at stoichiometry above, the adiabatic flame temperature is obtained:

$T_s=T_{O,o}+Z_s(T_{F,o}-T_{O,o}+Y_{F,o} \frac{Q}{C_p})$

Remark 1: The Burke-Schumann solution for a diffusion flame, despite its age, is still one of the most used model, under some refreshments such as Mixed-Is-Burned, or Equilibrium combustion model where the burned value of the quantities is obtained by Gibbs function extrema instead of the full adiabatic consumption.

Remark 2: In the Burke-Schuman solution, as the chemistry is infinitely fast, the flame is reduced to an interface separating the rich and lean sides. Therefore, quantities must be continuous at the flame and their gradients observe the jump conditions, from Sec. Conservation Laws:

$\left \{ \begin{array}{lll} s\rho D_F||\vec\nabla Y_F|| & = & \rho D_O||\vec\nabla Y_O|| \\ \frac{Q}{C_p}\rho D_F||\vec\nabla Y_F|| & = & \lambda ||\vec\nabla T|| \end{array} \right .$

Remark 3 From Remark 2, it is seen that the reactants fluxes must reach the flame in stoichimetric proportion (the fluxes and not the values themselves). Therefore, the relations developed in Remark 2 are true even without the equidiffusional approximation, whereas the relation based on $Z$ , no. On the other hand, the heat released by chemistry is conducted away from the flame thanks to the temperature gradients on both sides.

Dissipation Rate and Non-Equilibrium Effect

A very important quantity, derived from the mixture fraction concept, is the scalar dissipation rate, usually noted: $\chi$ . In the above introduction to non-premixed combustion, it has been said that a diffusion flame is fully controlled through: (i) the position of the stoichiometric line, dictating where the flame sheet lies; (ii) the gradients of fuel on one side and oxidizer on the other side, dictating the feeding rate of the reaction zone through diffusion and thus the strength of combustion. According the the mixture fraction definition, the location of the stoichiometric line is naturally tracked through the $Z_s$ iso-line and it is seen here how the mixture fraction is a convenient tracer to locate the flame. In the same manner, the mixture fraction field should also be able to give information on the strength of the chemistry as the gradients of reactants are directly linked to the mixture fraction distribution. The feeding rate of the reaction zone is characterized through the inverse of a time. Because it is done through diffusion, it must be obtained through a combination of mixture fraction gradient and diffusion coefficient (dimensional analysis):

$\chi_s = \frac{(\rho D)_s}{\rho_s} ||\vec \nabla Z||^2_s \qquad (\rho D)_s = \left ( \frac{\lambda}{C_p} \right )_s$

where the subscript s refers to quantities taken effectively where the reacting sheet is supposed to be, close to stoichiometry. This deduction of the scalar dissipation rate, scaling the feeding rate of the flame, is obtained here through physical arguments and is to be derived from equations below in a more mathematical manner. Note that the transport coefficient for the mixture fraction is identified to the one for temperature. This notation is usually used in the literature to emphasize that the rate of temperature diffusion (that is commensurable to the rate of species diffusion) is the retained parameter (as introduced in the following approach highlighting the role of the scalar dissipation rate). Depending if the diffusion process is free (mixing layer, boundary layer) or controlled (stagnation), the dissipation rate and then the non-premixed combustion may be of unsteady or steady nature.

Ignition / Extinction Curve

Let consider a mixing layer between fuel and oxidizer whose strain (and then intensity of reactant inter-diffusion) is carefully controlled. Let start from a flame already existing. Let also use the Damk\"ohler number as defined above with the inverse of the scalar dissipation rate being the characteristic mechanical time of flow of reactants feeding the reaction zone. For high Da (low dissipation rate), the reactants diffuse slowly. The reaction is not very intense but the chemistry is fully achieved such that the maximum temperature is reached. Now, let decrease slowly (to avoid any unsteady effect on the chemistry) the Da. On one hand, the rate of feeding of flame through diffusion increases, and so does the reaction rate. On the other hand, the chemistry may not have the time to `eat' every reactant molecules that begin to leak through the flame: the fully sensible value is not realized and a lower temperature is achieved. With still increasing the dissipation rate, this mechanism leads to lower temperature in the flame zone down to a level that cannot be sustained by combustion (that strongly depends on temperature). The diffusion flame extinguishes suddenly. It is said to be quenched. This is experienced in real life when blowing off a small wood fire. Slightly blowing increases the transfer between reactants and strengthens the reaction but blowing too much extinguishes the fire.

S-Curve Diffusion Flame.

If now the start state is a frozen mixing layer between fuel and oxidant, at low Da (high dissipation rate), the flow remains frozen. When increasing Da, there will be a point where the chemistry will self-accelerate and the flame will light up through a sudden increase in temperature. Because the starting temperature is low, the ignition Da number is higher that the extinction Da above, exhibiting an hysteresis phenomenon.

When looking at the trend of the maximum temperature only versus Da, a curve with a shape in `S' appears, named as S-Curve for diffusion flames (see beside figure).

Flamelet Equation

In turbulent combustion modelling, the flamelet regime is identified as a regime that saves the integrity of the flame structure (Sec. The Wrinkled Regime). Turbulent eddies do not enter the structure and only contord the flame at a large scale. Therefore, it is usual to call the flamelet equation, the equation that describes the diffusion flame structure when not perturbed by turbulence. The model problem is usually retained as the counterflow diffusion flame. The strain imposed to the flame mimicks the one due to inhomogeneities in a real turbulent flow.

Sketch of the Stretched Mixing Layer Created by the Counter-Flow Configuration.

A counter-flow diffusion flame is basically made of two ducts, front-to-front, one issuing fuel at a mass fraction $Y_{F,o}$ (may be diluted by an inert) and temperature $T_{F,o}$ , and the other issuing oxidizer at a mass fraction $Y_{O,o}$ (the inert is usually di-nitrogen in air such that $Y_{O,o} =.23310$ is a value commonly encountered) and temperature $T_{O,o}$ . A stretched mixging layer develops at the location the flows oppose. This stretch can be modulated by the flow rate from the feeding ducts. The diffusion flame sits around the iso-surface corresponding to stoichiometry. The figure beside sketches this widely used configuration, with the coordinate frame origin taken at the stagnation location in the middle of the mixing layer.

The Howarth-Dorodnitzyn transform and the Chapman approximation are applied to the above mixture fraction transport equation. In the Chapman approximation, the thermal dependence of $(\rho D)$ is approximated as $\rho^{-1}$ . The Howarth-Dorodnitzyn transform introduces $\rho$ in the space coordinate system: $\vec\nabla \rightarrow \rho\vec\nabla \quad ; \quad \nabla\cdot\rightarrow \rho\nabla\cdot$ . The effect of these both mathematical operations is to `digest' the thermal variation of quantities such as density or transport coefficient. Hence, the mixture fraction equation comes in a simpler mathematical shape:

$\rho\vec U \cdot \vec\nabla Z = \rho_s ( \rho D )_s \Delta Z$

Here the references are taken in the stoichiometric zone (s subscript).

For the sake of simplicity, we shall use the potential flow result for this counter-flow problem with two identical jets although the presence of the flame (the associated density change) jeopardizes this theory, strictly speaking. The velocity solution is thus:

$\left \{ \begin{array}{lll} u & = & -a\int\frac{dx}{\rho} \\ v & = & \frac{a}{2}\int\frac{dr}{\rho} \end{array} \right.$

$a$ is named the stretch of the mixing layer and is the inverse of a time.

In such a flow, the transport equation for the mixture fraction may be written as a purely 1-D problem (assuming constant properties):

$-ax\frac{dZ}{dx} = \rho_s (\rho D)_s \frac{d^2Z}{dx^2}$

from which a `length' scale emerges for the mixing layer: $\sqrt{\rho_s (\rho D)_s /a}$ . Non-dimensionalized, with $\eta = x \sqrt{a/\rho_s (\rho D)_s}$ :

$\eta \frac{dZ}{d\eta} + \frac{d^2 Z}{d \eta^2} = 0$

The boundary conditions may be approximated as $\eta \rightarrow +\infty\quad ; \quad Z\rightarrow 1$ and $\eta \rightarrow -\infty\quad ; \quad Z\rightarrow 0$ . Note that this approximation requires a thin layer compared to the counter-flow distance, and thus a relatively high strain rate. Then, the solution for $Z$ :

$Z=\frac{1}{2}(1+erf{\frac{\eta}{\sqrt{2}}})$

The scalar dissipation rate is deduced from this equation as:

$\chi = \frac{1}{2\pi}\exp\left( - [erf^{-1} (2Z)]^2 \right) \frac{a\rho^2}{\rho_s^2}$

Because combustion is highly temperature-dependent, T is certainly the scalar to which attention must be paid for. The temperature equation in the Low-Mach Number regime (Sec. Low-Mach Number Equations) is written below in steady-state:

$\rho\vec U \cdot \vec\nabla T = \nabla\cdot \frac{\lambda}{C_p}\vec\nabla T + \frac{Q}{C_p}\nu_F \bar M_F \dot\omega$

After the Chapman approximation and the Howarth-Dorodnitzyn transform:

$\frac{\rho}{\rho_s}\vec U \cdot \vec\nabla T = \left ( \frac{\lambda}{C_p}\right )_s \nabla\cdot \vec\nabla T + \frac{Q}{\rho\rho_s C_p}\nu_F \bar M_F \dot\omega$

Here the references are taken in the flame, i.e. close to the stoichiometric line (s subscript).

In a non-premixed system, strictly speaking, $\vec U$ , is not really relevant as the flame is fully controlled by the diffusion process. Notwithstanding, in practice, non-premixed flames must be stabilized by creating a strain in the direction of diffusion. This is the reason why the velocity is left in the equation. Because a diffusion flame is fully described by the mixture fraction field, a change in coordinate can be applied:

$\left ( \begin{array}{c} r \\ x \end{array}\right ) \longrightarrow \left ( \begin{array}{c} r \\ Z \end{array}\right )$

The Jacobian of the transform is given as:

$\left [ \begin{array}{cc} \frac{\partial r}{\partial r} & \frac{\partial Z}{\partial r} \\ \frac{\partial r}{\partial x} & \frac{\partial Z}{\partial x} \end{array}\right ] = \left [ \begin{array}{cc} 1 & 0 \\ 0 & (\rho l_d)^{-1} \end{array}\right ]$

Note that the diffusive layer of thickness $l_d$ is defined as the region of transition between fuel and oxidizer and is thus given by the gradient of Z along the x direction. This transform is applied to the vectorial operators:

$\nabla\cdot = \nabla_r\cdot+\nabla_x\cdot=\nabla_r\cdot+\nabla_x\cdot Z\nabla_Z\cdot$

$\vec\nabla = \vec\nabla_r+\vec\nabla_x=\vec\nabla_r+\vec\nabla_x Z\nabla_Z\cdot$

With this transform, the above temperature equation looks like:

$\frac{\rho}{\rho_s}\vec U \cdot (\vec\nabla_r T + \vec\nabla Z \nabla_Z\cdot T) = \left (\frac{\lambda}{C_p}\right )_s \left ( \nabla_r\cdot \vec\nabla_r T + \nabla_x\cdot Z \nabla_Z\cdot \vec\nabla_x Z \nabla_Z\cdot T + \nabla_x\cdot Z \nabla_Z\cdot \vec\nabla_r T + \nabla_r\cdot \vec\nabla_x Z \nabla_Z\cdot T \right ) + \frac{Q}{\rho\rho_s C_p}\nu_F \bar M_F \dot\omega$

As mentioned above, the velocity and the variation along the tangential direction to the main flame structure r are not supposed to play a major role. By emphasizing the role of the gradient of Z along x as a key parameter defining the configuration the following equation is obtained:

$0 = \left (\frac{\lambda}{C_p}\right )_s ||\vec\nabla_x Z||^2 \Delta_Z T + \frac{Q}{\rho\rho_s C_p}\nu_F \bar M_F \dot\omega$

This equation (sometimes named the flamelet equation) serves as the basic framework to study the structure of diffusion flames. It highlights the role of the dissipation rate with respect to the strength of the source term and shows that the dissipation rate calibrates the combustion intensity. From Sec. Simplified Solution of Diffusion Flame, the temperature jump at stoichiometry is identified as $Q Y_{F,o} Z_s/Cp$ .

Although chemistry in a diffusion flame is usually very fast, one expects a very residual amount of reactants in the reaction zone laying along the stoichiometric line. It thus means that equilibrium is not reached in practice and that the temperature is somewhat lower than the equilibrium temperature. Thermodynamics is moderated by kinetics effects.

To describe the structure of the diffusion flame, the reduced mixture fraction is set:

$\xi = \frac{Z-Z_s}{Z_s(1-Z_s)\varepsilon}$

The utility of the reduced mixture fraction is to focus on the reaction zone. This reaction zone is supposed located on the stoichiometric line (this is why the reduced mixture fraction is centred on $Z_s$ ) and to be very thin (reason of the introduction of the magnifying factor $\varepsilon$ ).

Ignition

The flamelet equation may be used to explore the behaviour of the mixing-layer for not too high strain rate. As mentioned earlier, the Damk&oumlaut;hler may be high enough to allow a self-ignition. In this regime, the following assumptions may be made:

the exponential argument in the reaction term, $T_a/T$ , is linearized around the frozen temperature value (indeed, as the activation temperature is large, a small change in the temperature from its frozen value is magnified in the reaction rate magnitude): $T_a/T = T_a/T_o-T_a/T_o^2(T-T_o) = T_a/T_o - \phi$ . It can be seen that a small change of temperature leading to $T_a/T_o^2(T-T_o)\quad {\mathcal O} (1)$ impacts the exponential term of the reaction rate.
the reactants may be approximated by their frozen flow values: $Y_F=Y_{F,o}Z-\frac{Cp T_o^2}{QT_a}\phi\approx Y_{F,o}Z \quad ; \quad Y_O=Y_{O,o}(1-Z)-s\frac{Cp T_o^2}{QT_a} \phi \approx Y_{O,o}(1-Z)$ , under the hypothesis of large activation energy.

The flamelet equation is recast as: $\Delta_Z\phi = - Z^{n_F}(1-Z)^{n_O} e^{\phi}\overbrace{\frac{Q\nu_F}{\rho_s\lambda_s}\frac{A}{\nabla^2_x Z \bar M_F^{n_F-1}\bar M_O^{n_O}} \frac{T_a}{T_o^2} e^{-\frac{T_a}{T_o}}Y_{F,o}^{n_F}Y_{O,o}^{n_O}}^{\mathcal D}$ In this expression $\mathcal D$ is a Damk\"ohler number.

This expression must be solved numerically but, with the help of simplifications, an analytical solution may be found. It must be noticed that this analytical solution has no quantitative prediction capabilities. However, it is sufficient to illustrate the critical Damk\"ohler number for ignition for pedagogical purpose.

First, the equation is integrated from $Z_s$ to $Z$ 'somewhere in the rich side'.

$\int_{Z_s}^{Z} dZ \Delta_Z\phi = \left [\frac{d\phi}{dZ} \right]_{Z_s}^{Z} = - \int_{Z_s}^{Z} dZ Z^{n_F}(1-Z)^{n_O} e^{\phi}{\mathcal D}$

$\phi$ may be approximated as, on the asymptote $\hat{\phi}(Z-1)/(Z_s-1)$ with a maximum $\hat{\phi}$ on the stoichiometric line, where ignition is supposed the strongest. Then, on the stoichiometric line, the derivative of $\phi$ vanishes and, sufficiently far from the stoichiometric line it tends to $\hat{\phi}/(Z_s-1)$ . To keep the solution simple, we choose the integration limit to be on the border of the region where $\phi$ remains within 10% of its maximum such that it may be approximated as constant on the right hand side. The value of $Z$ that fits this trade-off of being sufficiently far from $Z_s$ to allow the derivative of $\phi$ to reach its asymptote and sufficiently close to prevent $\phi$ from decreasing is estimated to 40% of $Z_s-Z$ . This estimation, quite disputable, is based on temperature profile observed in diffusion flames at low Da when the diffusive effect strongly smoothes it.

Ignition Da curve in a cold diffusion layer.

The end equation has the form:

$\hat{\phi}e^{-\hat{\phi}}=-{\mathcal D} \overbrace{(Z_s-1)\left (\frac{Z^2-Z_s^2}{2}-\frac{Z^3-Z_s^3}{3}\right )}^{a}$

This makes $\mathcal D$ a function of $\hat{\phi}$ , maximum temperature in the ignition zone at stoichiometry, that has a maximum for $\hat{\phi}=1$ . This gives the bottom tipping point of the S-curve presented above. For the academical case $Z_s = .5$ and $n_F=n_O=1$ , this gives $\mathcal D = 16.$ which is in a satisfactory agreement with the full numerical solution. It is interesting to note that the numerical solution follows the physics process of ignition of the diffusion layer by slowly increasing the Damk\"ohler number. Therefore, it cannot catch properly the turning point and the doubled-value relation. Instead, once the maximum Damk\"ohler number corresponding to ignition is reached, the reduced temperature soars sharply. Further increasing D in the forbidden range, no solution can be found. Physically, for larger values of the Damköhler number, the configuration jumps from this nearly frozen situation to an equilibrium situation with a well developed diffusion flame, see the S-Curve above. The simplified equation for ignition is thus unadapted and a specific solution for a developed diffusion flame is now to be presented.

Diffusion Burning Process

To a certain extend, the flamelet equation Sec. Flamelet Equation may be used to present a solution for the non-premixed flame when the flame is already ignited, and the chemistry strength controlled by diffusion (we have stressed the author that an interesting aspect of diffusion flames is the ability to control the reaction rate through the level of mixing in the neighbourhood of the flame). As already emphasized, due to a strong non-linearity in temperature, the structure of a burning flame can never be far from the Burke-Schumann simplified solution.

Non-dimensionalized, the temperature solution in the flamelet may be written as:

$\theta=1+(Z-Z_s)\Delta T_o + \frac{Z_s - Z}{1-Z_s} -\varepsilon\Gamma$

in the rich side of the mixing layer, and:

$\theta=1+(Z-Z_s)\Delta T_o + \frac{Z-Z_s}{Z_s} -\varepsilon\Gamma$

in the lean side. On recognizes the structure of the Burke-Schumann solution (Sec. Simplified Diffusion Flame Solution). $\Delta T_o$ is the non-dimensionalized temperature difference between the reactants feed streams:

$\Delta T_o = \frac{T_{F,o}-T_{O,o}}{Q Y_{F,o} Z_s/C_p}$

Those solutions are amended by the term $\varepsilon\Gamma$ with $\Gamma\sim{\mathcal O}(1)$ in the reaction zone and vanishing on the sides. This amendment corresponds to the departure from the equilibrium solution due to kinetics effect. The quantity $\varepsilon$ is taken small for the physical reasons explained above.

The departure in temperature from equilibrium in the reaction zone is reflected into the species concentrations, as it corresponds to still unburned reactants. With the help of the above Conservation Laws one finds:

Lean side: $Y_F=\varepsilon\Gamma Y_{F,o}Z_s$ , $Y_O=Y_{O,o}(1-Z)-s(Y_{F,o}Z-\varepsilon\Gamma Y_{F,o}Z_s)$
Rich side: $Y_F=Y_{F,o}Z+\frac{1}{s}(s Y_{F,o}Z_s\varepsilon\Gamma-Y_{O,o}(1-Z))$ , $Y_O=sY_{F,o}Z_s \varepsilon\Gamma$

The Premixed Regime

Sketch of a premixed flame

In contrast to the non-premixed regime above, the reactants are here well mixed before entering the combustion chamber. Chemical reaction can occur everywhere and this flame can propagate upstream into the feeding system as a subsonic (deflagration regime) chemical wave. This presents lots of safety issues. Some situations prevent them: (i) the mixture is made too rich (lot of fuel compared to oxidizer) or too lean (too much oxidizer) such that the flame is close to its flammability limits (it cannot easily propagate); (ii) the feeding system and regions where the flame is not wanted are designed such that they impose strong heat loss to the flame in order to quench it. For a given thermodynamical state of the mixture (composition, temperature, pressure), the flame has its own dynamics (speed, heat release, etc) on which there is few control: the wave exchanges mass and energy through diffusion process in the fresh gases. On the other hand, those well defined quantities are convenient to describe the flame characteristics. The mechanism of spontaneous propagation towards fresh gas through the thermal transfer from the combustion zone to the immediate slice of fresh gas such that the ignition temperature is eventually reached for this latter was highlighted as early as by the end of the 19th century by Mallard and LeChatelier. The reason the chemical wave is contained in a narrow region of reaction propagating upstream is the consequence of the discussion on the non-linearity of the combustion with temperature in the Sec. Fundamental Aspects. It is of interest to compare the orders of magnitude of the temperature dependent term $\exp{(-\beta(1-\theta)/(1-\alpha(1-\theta)))}$ of the reaction source upstream in the fresh gas ( $\theta\rightarrow 0$ ) and in the reaction zone close to equilibrium temperature ( $\theta\rightarrow 1$ ) for the set of representative values: $\beta = 10$ and $\alpha=0.9$ . It is found that the reaction is about $10^{43}$ times slower in the fresh gas than close to the burned gas. It is known that the chemical time scale is about 0.1 ms in the reaction zone of a typical flame, then the typical reaction time in the fresh gas in normal conditions is about $10^{39} s$ . To be compared with the order of magnitude of the estimated Universe age: $1 0^{17} s$ . Non-negligible chemistry is only confined in a thin reaction zone stuck to the hot burned gas at equilibrium temperature. In this zone, the #Damköhler number is high, in contrast to in the fresh mixture. It is natural and convenient to consider that the reaction rate is strictly zero everywhere except in this small reaction zone (one recovers the Dirac-like shape of the reaction profile, provided that one can see the upstream flow as a region of increasing temperature towards the combustion zone and the downstream flow as in fully equilibrium).

As the premixed flame is a reaction wave propagating from burned to fresh gases, the basic parameter is known to be the progress variable. In the fresh gas, the progress variable is conventionally put to zero. In the burned gas, it equals unity. Across the flame, the intermediate values describe the progress of the reaction to turn into burned gas the fresh gas penetrating the flame sheet. A progress variable can be set with the help of any quantity like temperature, reactant mass fraction, provided it is bounded by a single value in the burned gas and another one in the fresh gas. The progress variable is usually named c, in usual notations:

$c=\frac{T-T_f}{T_b-T_f}$

It is seen that c is a normalization of a scalar quantity. As mentioned above, the scalar transport equations are assumed linear such that the transport equation for c can be obtained directly. Actually, the transport equation for T (Sec. Transport Equations) is linear if constant heat capacity is further assumed (combustion of hydrocarbon in air implies a large excess of nitrogen whose heat capacity is only slightly varying) and the progress variable equation is directly obtained (here for a default of fuel - lean combustion):

$\frac{D\rho c}{Dt}=\nabla\cdot\rho D\vec\nabla c + \frac{\nu_F\bar M_F}{Y_{F,u}}\dot\omega \qquad \rho D = \frac{\lambda}{C_p}$

The fact that the default or excess of fuel has been discussed above leads to the introduction of another quantity: the equivalence ratio. The equivalence ratio, usually noted $\Phi$ , is the ratio of two ratios. The first one is the ratio of the mass of fuel with the mass of oxidizer in the mixture. The second one is the same ratio for a mixture at stoichiometry. Hence, when the equivalence ratio equals unity, the mixture is at stoichiometry. If it is greater than unity, the mixture is named rich as there is an excess of fuel. In contrast, when it is smaller than unity the mixture is named lean. The equivalence ratio presented here for premixed flames has little connection with the equivalence ratio introduced earlier regarding the non-premixed regime. Basically, the equivalence ratio as defined for non-premixed flames gives the equivalence ratio of a premixed mixture with the same mass of fuel and oxidizer. Moreover, the equivalence ratio as defined for a premixed mixture can be obtained based on the mixture fraction (it is thus the local equivalence ratio at a point in the non-homogeneous mixture described by the mixture fraction). From the definitions given above:

$\Phi=\frac{Z}{1-Z}\frac{1-Z_s}{Z_s}$

Premixed Flame Péclet Number

Earlier in this section, it has been said that a premixed flame posses its own dynamics, as a free propagating surface, and has thus characteristic quantities. For this reason, a Péclet number may be defined, based on these quantities. The Péclet number has the same structure as the Reynolds number but the dynamical viscosity is replaced by the ratio of the thermal conductivity and the heat capacity of the mixture. The thickness $\delta_L$ of a premixed flame is essentially thermal. It means that it corresponds to the distance of the temperature rise between fresh and burned gases. This thickness is below the millimetre for conventional flames. The width of the reaction zone inside this flame is even smaller, by about one order of magnitude. This reaction zone is stuck to the hot side of the flame due to the high thermal dependency of the combustion reactions, as seen above. Hence, the flame region is essentially governed by a convection-diffusion process, the source term being negligible in most of it.

It is convenient to write the progress variable transport equation in a steady-state framework. The quantities at flame temperature ( $_f$ ) are used to non-dimensionalize the equation:

$\overbrace{(\rho ||\vec S_L||)}^{||\vec M||}\frac{\vec S_L}{S_L}\vec\nabla c = (\rho D)_f \nabla\cdot (\rho D)^* \vec\nabla c$

Note that the source term is neglected, consistently with what has been said above. This convection-diffusion equation makes appear a first approximation of a flame Péclet number:

$Pe_f = \frac{||\vec M|| \delta_L}{(\rho D)_f} \approx 1$

From the Péclet number, it is possible to obtain an expression for the flame velocity (remembering that $\delta_L/S_{L,f} \approx \tau_c$ , vid. inf. Sec. Three Turbulent-Flame Interaction Regimes):

$S_{L,f}^2\approx \frac{(\rho D)_f}{\rho_f \tau_c}$

For typical hydrocarbon flames, the speed is some tens of centimetres per second and the diffusivity is some $10^{-5}$ square metres per second. The chemical time in the reaction zone of about one tenth of a millisecond is recovered.

Details of the Premixed Unstrained Planar Flame

A plane combustion wave propagating in a homogeneous fresh mixture is the reference case to describe the premixed regime. At constant speed, it is convenient to see the flame at rest with a flowing upstream mixture. This is actually the way propagating flames are usually stabilized. In the frame of description, the physics is 1-D, steady with a uniform (the flame is said unstrained) mass flowing across the system. Two types of equations are thus sufficient to describe the problem, the temperature transport equation and the species transport equations, as in Sec. Transport Equations. The transport coefficients will be chosen as equal: $\rho D_i = \lambda / C_p$ (unity Lewis numbers). Suppose the 1-D domain is described thanks to a conventional (Ox) axis with a flame propagating towards negative x (this is the conventional usage), the boundary conditions are:

in the frozen mixture:
- $Y_i \rightarrow Y_{i,u} \qquad ; \qquad x\rightarrow-\infty$
- $T \rightarrow T_u \qquad ; \qquad x\rightarrow-\infty$
in the burned gas region supposed at equilibrium:
- $Y_i \rightarrow Y_{i,b} \qquad ; \qquad x\rightarrow+\infty$
- $T \rightarrow T_b \qquad ; \qquad x\rightarrow+\infty$

$Y_{i,b}$ and $T_b$ are obtained from Sec. Conservation Laws.

The quantities that have been mentioned just above (scalar and temperature profiles, mass flow rate through the system) are the solution to be sought.

According to the discussions above, the temperature transport equation in its full normalized form may be written as (lean /stoichiometric case):

$||\vec M|| \frac{\partial \theta}{\partial x} = \frac{\partial \ }{\partial x}\frac{\lambda}{Cp} \frac{\partial \theta}{\partial x} + \frac{\nu_F \bar M_F B}{Y_{F,s}} \prod_{i=O,F}(Y_{i,u}^*-\theta\frac{Y_{F,s}}{Y_{i,s}}s^{\delta_{i,O}} )^{n_i} \exp{-\beta\frac{1-\theta}{1-\alpha(1-\theta)}}$

This equation is further simplified by the variable change $d\xi=||\vec M||/(\lambda/Cp)dx$ :

$\frac{\partial \theta}{\partial \xi}=\frac{\partial^2 \theta}{\partial \xi^2} + \overbrace{\frac{\lambda/Cp}{||\vec M||^2}\frac{\nu_F \bar M_F B}{Y_{F,s}}}^{\Lambda} \prod_{i=O,F}(Y_{i,u}^*-\theta)^{n_i} \exp{-\beta\frac{1-\theta}{1-\alpha(1-\theta)}}$

Although somewhat out of scope, the existence and unicity of the solution of this type of equation are usually demonstrated with the help of the Schauder Theorem and Maximum Principle. From the point of view of physicists and engineers, the solution that is found analytically is de facto considered as the unique solution of the equation.

Scenarii of Combustion Process in the Phase Portrait

In the frame moving with the flame, both phase variables are the reduced temperature and its gradient. To ease the reading with usual notations, it is written: $X_1 = \theta \quad ; \quad X_2 = \partial \theta / \partial \xi = \dot \theta$ . The system arises:

$\left\{ \begin{array}{lll} \dot X_1 & = & X_2 \\ \dot X_2 & = & X_2 - \varpi(X_1) \end{array} \right.$

with $\varpi$ being the full source term in the above equation.

In the frame moving with the flame, two singular nodes are found in the frozen flow $(X_1,X_2) = (0,0)$ and the equilibrium region $(\theta_b,0)$ , i.e when $\varpi(X_1)$ vanishes. $x_1,x_2$ are defined as small departures from the singular nodes such that the linearized system in their neighbourhood is:

$\left\{ \begin{array}{lll} \dot x_1 & = & l_{1,1} x_1 + l_{1,2} x_2 \\ \dot x_2 & = & l_{2,1} x_1 + l_{2,2} x_2 \end{array} \right.$

provided: $l_{1,1} = 0 \quad l_{1,2} = 1 \quad l_{2,1} = -\varpi'_{X_1 = 0,\theta_b} \quad l_{2,2} = 1$ . The characteristic polynom is, in usual notations: $s^2 - s + \varpi'$ such that the eigenvalues are:

$s^{\pm}=\frac{1\pm\sqrt{1-4\varpi'}}{2}$

A priori, those eigenvalues may be (i) real distinct, (ii) real identical, or (iii) conjugated complex. In the first case, the orbits in the phase diagram are organized, in the immediate neighbourhood of the singular node, with respect to the eigenvectors directions associated to the eigenvalues. The following task is to identify the nature of those eigenvalues and of the corresponding nodes. Because $0 < X_1 < \theta_b$ is bounded, complex eigenvalues are excluded as they would lead to a spiral node. This remark is important for the node on the cold side as it imposes a bound:

$\varpi'_{X_1=0} \le \frac{1}{4}$

As the mass flow rate through the flame is included into $\varpi$ , it imposes a minimum value on the flame speed to tackle with the cold boundary difficulty (rise of the chemical rate in the frozen flow). In this condition, it is an unstable node (improper in case of equality). On the other hand, because $\varpi'|_{X_1=\theta_b}$ is not positive, the node on the hot side is found as a saddle point. The overall scenario of combustion within the flame is thus an orbit leaving the cold node to join the hot node by branching on a trajectory compatible with the negative eigenvalue of the saddle.

Sketch of orbits for a combustion process across a premixed front. Dashed lines represent forbidden orbits from the physics. The red line describes the orbit expected in an idealized combustion process.

It must be noted that the associated eigenvectors are of the form:

$\left ( \begin{array}{l} t \\ t s^{\pm} \end{array} \right ); \quad t\in \Re^*$

that is, on the cold node, a positive departure on $X_1$ following any of the two eigendirections, leads to a consistent creation of positive temperature gradient, while on the hot node, only the stable direction will allow a consistent creation of a positive temperature gradient for any departure of the temperature towards region where it is inferior to $\theta_b$ . Another remark is the structure of the eigendirections. The leaving directions on the cold node have a slope larger than the one of $\varpi$ while the stable direction of the hot node has a slope smaller than the one of $\varpi$ . It means that there is some point where the orbit must cross the profile of the chemical term versus temperature. For that temperature, the gradient equals the reaction rate through construction of the phase space. When looking back to the equation of the premixed flame, it happens in a region of inflexion for the temperature (the second order derivative must vanish). Furthermore, at this intersection, the orbit is horizontal (if the frame of reference for $(X_1, X_2)$ is Cartesian) due to the shape of the premixed flame equation above that can be recast into $X_{2,X_1}' = (X_2-\varpi)/X_2$ . Close to the cold node, the orbits have a shape of parabola whose axis is the direction with the largest eigenvalue magnitude. Close to the hot node, the orbits have an hyperboloid shape with asymptots as the eigendirections. Now the ingredients are here to draw a sketchy scenario of the combustion in a premixed flame. It will be superimposed on the reaction rate graph studied in Sec. Fundamentals. Some typical orbits from the above analysis are drawn in the figure on the right. The basic geometrical arguments developed are reproduced. In particular, the dashed lines represent forbidden orbits by the physics (boundedness of $X_1$ , irreversibility). Orbits must be travelled from left to right, corresponding to increasing free parameter $\xi$ . It is of integral importance to remember that, in combustion in conventional conditions, the source term is highly non-linear. Therefore, it is localized in a very thin sheet and, upstream, $\varpi$ and $\varpi' \rightarrow 0$ . Only the most unstable eigendirection of the cold node is compatible as an orbit. The other trajectories, being paraboloidal, are tangent to the other direction that is flat at the limit. It physically means an elevation of temperature in bulk, that is contradictory to what is expected from a highly non-linear combustion term. Hence, the additional orbit in red is the one expected in idealized combustion conventional conditions.

Sketch of orbits for a combustion process across a premixed front, DNS simulations for the usual combinations of and . Note that the case with does not mean an athermal reaction but the mathematical simplification of the denominator of the exponential argument of the source term.

Sketch of orbits for a combustion process across a premixed front, DNS simulations for the usual combinations of $\beta$ and $\alpha$ . Note that the case with $\alpha = 0$ does not mean an athermal reaction but the mathematical simplification of the denominator of the exponential argument of the source term.

The following picture is the result of actual computations of the above 1-D flame equation (stoichiometry and $n_i = 1$ ) with the help of a high-order (6) code. The already presented combinations of $\beta$ and $\alpha$ are used and the orbits are retrieved. For most of the cases, as predicted above, the system selects a solution leaving the cold node with the most unstable direction (identifiable with its slope close to unity for vanishing $\varpi'_{X_1=0}$ ). There is also an additional curve, for $\beta = 10$ with the denominator of the exponential argument suppressed. This curve is remarkable as the solution selected by the system leaves the cold node in a manner fully controlled by the $\varpi'_{X_1=0}$ . This is a very singular solution, not expected in combustion in conventional conditions, as explained above. The purpose of this remark is to question the well-posedness of considering simplifying the exponential argument for $\beta$ "sufficiently high", as it is usually proposed for this type of modelling. As observed, the dynamical system analysis demonstrates a switch in the nature of the solution selected. At the physical level, when the orbit follows the most unstable direction with a slope close to unity, it means that $X_2$ "follows" $X_1$ , which is a signature of a diffusion process. In other words, the preheating mechanism of a premixed flame propagation as proposed for more than one century is in work. On the other hand, when $X_2$ is dependent on the evolution of $\varpi'$ , it shows that "cold" chemistry drives the solution in the frozen flow and not the acknowledged mechanism of deflagration.

Flame Solution

As already mentioned, the flame system may be split into three zones. Upstream, the conventional mechanism of deflagration is supported by diffusion of heat. Downstream, the mixture is at equilibrium after combustion. In between, there exists the reaction layer. For large $\beta$ the reaction layer is very thin such that it can be seen as a discontinuity between the fresh and burned gases. This is this difference in scales that introduces the use of the asymptotic method to resolve some of the flame characteristics, such as speed, time, heat region thickness, or reaction zone thickness.

The domain is partitioned, according to this zoning defined by the scales driving the physics with an outer domain, driven by large scales and an inner domain refining the description within the discontinuity. If the discontinuity (flame reaction zone) is at $\xi=0$ , everywhere but 0, the equation is simplified as:

$\frac{\partial^2 \theta}{\partial \xi^2} = \frac{\partial\theta}{\partial \xi}$

For $\xi>0$ , it is expected that the mixture has reached equilibrium chemistry, such that: $\forall \xi > 0, \quad \theta=\theta_b \quad ; \quad \partial \theta / \partial \xi =0$ .
For $\xi < 0$ , this is the preheat zone and the solution is $\theta = \theta_b e^{\xi}$ with $\theta$ reaching $\theta_b$ at the disconstinuity and vanishing, together with its gradient, far upstream in the frozen mixture.

The solution for the species and the value of $\theta_b$ are obtained from Conservation Equations above. The 'big picture' is thus an exponential variation in the thermal thickness matched with a plateau in the downstream region, the line of matching being the discontinuity (flame) that has no thickness at this scale of description.

To refine the analysis in the discontinuity region, a magnifying factor $\varepsilon$ is used to stretch the coordinates: $\xi = \varepsilon \Xi$ . The inner solution is thus a slowly-varying function of $\Xi$ . Hence, in this inner region, the equation for the premixed flame becomes:

$\frac{1}{\varepsilon}\frac{\partial \theta}{\partial \Xi}=\frac{1}{\varepsilon^2}\frac{\partial^2 \theta}{\partial \Xi^2} +\varpi$

In order to stretch and 'look inside' a discontinuity, $\varepsilon$ is very small. It yields two remarks:

convection is negligible compared to diffusion. The heat losses from the reaction zone are essentially diffusion driven.
The reaction zone is governed by a diffusion-reaction budget and the reaction term $\varpi$

Combustion