Formula of Laminar Flame Speed Coupled with Differential Equation for
Temperature in Low-Mach-Number Model

Keigo Wada111Institute of Philosophy in Interdisciplinary Sciences, Kanazawa
University, Kakuma, Kanazawa, Ishikawa 920-1192, Japan 222Additional Post
Member of OCAMI, Osaka Metropolitan University, Sugimoto, Sumiyoshi-ku, Osaka
558-8585, Japan

<EMAIL_ADDRESS>
###### Abstract

The non-monotonic profile of temperature is to be considered in the context of
combustion inside tubes or thermonuclear flames, which may accelerates to
become detonation waves. This transition is known as deflagration-to-
detonation transition. Focusing on this point, we extend the study of
adiabatic flames, which assumes zero temperature gradient at the burned-side
edge of a flame front, to the case of non-zero temperature gradient by
introducing compression zone behind the front. Such a treatment has been
performed in the previous work based on the ignition temperature
approximation. In that case, the effect of gas compression on flames was
described by asymptotic solutions with respect to small squared Mach numbers:
for example, in the compression zone, inner solutions are expressed by use of
Lambert $W$ function. However, due to the aim of obtaining analytical
solutions inside a flame front, the temperature dependence of reaction term
was ignored and, as a result, the cold boundary difficulty was not treated
seriously. This issue is resolved in this study by use of large activation
energy asymptotics with the reaction term expressed by the Arrhenius type. The
temperature distribution is computed numerically in the reaction zone and
calculated analytically in the compression zone. The burning-rate eigenvalue
problem is also considered to derive a formula of laminar flame speed, which
is coupled with differential equation for temperature in the reaction zone.
Then, the dependence of laminar flame speed on thermal-diffusive parameters,
accompanied with hydrodynamic properties brought by the compression effect,
becomes clear. The non-zero gradient condition of temperature on the burned
side of a flame front is obtained for non-zero values of Mach number. Opposed
to the monotonic profile of temperature in the zero-Mach-number limit, the
condition of negative temperature gradient brings the non-monotonic behavior
of temperature. It is found that the smallest value of laminar flame speed is
achieved for some negative small value of temperature gradient.

Keywords: gas compression, temperature gradient, laminar flame speed, singular
perturbations, large activation energy

##  1 INTRODUCTION

The adiabatic profile of a premixed flame front has been studied usually based
on the assumption of zero temperature gradient at the burned-side edge of the
front [52]. In the zero-Mach-number limit, this assumption is consistent with
the results obtained in experiments [38]. In this case, the large activation
energy asymptotics is a powerful method to understand the flow properties
inside and outside the reaction zone [65, 66, 51, 53, 45]. The boundary
conditions across the reaction zone have been derived based on the method of
matched asymptotic expansions [21]. These conditions are in turn used for the
study of preheat zone inside a flame front [12, 47, 44, 43], which gives the
correction of Darrieus-Landau instability (DLI) [15, 36] due to the effects of
curvature of a flame front [41] or tangential velocity [19, 20].

In the study of reaction zone, the exothermic reaction term is usually assumed
to be the Arrhenius type, which has an explicit dependence on temperature [54,
55]. In this case, it is possible to derive the parameter dependence of
laminar flame speed by considering the burning-rate eigenvalue problem [22,
58, 6]. Due to the low-Mach-number limit, the thermal-diffusive properties are
decoupled completely from the hydrodynamic flow inside the reaction zone [29].
Therefore, only the heat-conduction and mass-fraction equations are calculated
to give the asymptotic solution of burning-rate eigenvalue under the
assumption of zero temperature gradient on the burned side. On the other hand,
the non-zero gradient of temperature has also been studied [39, 37]. In this
study, we consider the relation between the temperature gradient and laminar
flame speed. This relation may be important in some situations related to the
deflagration-to-detonation transition (DDT), where the non-zero values of
temperature gradient are observed for a fast propagating deflagration or a
reaction wave behind the shock of detonation [11].

The DDT is a serious problem in some scientific fields: for example, the
terrestrial combustion in narrow channels [34] and mines [16], or the
thermonuclear explosions which are related to the understanding of local
expansion rate of the universe [25]. The mechanism of DDT remains unclear and
the numerical analysis of flow properties has been performed in channels [30,
35]. In confined circumstances, the acceleration of flame propagation speed is
induced by the effects of thermal expansion and wall frictions [4, 17]. The
increase of area of flame surface is essential for the rise of burning rate
[32, 5, 42]. This is mainly due to the corrugation of front invoked by the DLI
[10, 46]. The theoretical study about the effect of gas compression on the DLI
in tubes has been performed by small Mach number expansions [8, 57].

The DDT is also important when we consider the Type Ia supernova (SN Ia) in
the process of thermonuclear explosions [62, 49]. The SN Ia has attracted many
researchers’ interest to its theoretical analysis. This is because the SN Ia
may be modelled by the classical thermal theory of flames [40]: for example,
the burning of premixed carbon-oxygen fuel [28, 2]. At the early stage of
explosions, the propagation of flames is observed with the Mach number in the
range of $0.01$-$0.1$ [64, 1, 48]. In this case, the spatial distribution of
temperature has been studied numerically based on the deflagration model [24],
detonation model [33, 50] and reconstruction model [56]. In every work, the
non-monotonic profile of temperature is reported with its maximum value is
achieved inside the reaction zone. As for the acceleration of thermonuclear
flames, the Rayleigh-Taylor instability due to buoyancy effect is an essential
factor inducing the corrugated fronts [9]. Therefore, the multidimensional
structure, which lead to the formation of turbulence, must be considered in
the process of DDT. Although the curvature effect is neglected in this study,
our results may be adopted to the local estimation of properties of flames
propagating with small Mach numbers at the early stage of DDT. Indeed, our
analysis is performed in the microscopic scale expressed by the order of
thickness of reaction zone.

In the previous work [61], the compression zone, whose thickness is
$O({Ma}^{2})$, was introduced behind a flame front under ${Ma}^{2}\ll 1$ with
${Ma}$ denoting the Mach number, which is defined by the ratio of laminar
flame speed to local adiabatic sound speed: similar approaches are found, for
example, in [31, 7]. This approach enables us to capture the non-adiabatic
temperature profile with the temperature gradient at the burned-side edge
being no longer zero. Especially, such a non-adiabatic profile is described by
inner solutions inside the compression zone, which are expressed by use of
Lambert $W$ function [13]. However, the strong ignition temperature
approximation was employed in order to obtain analytical solutions inside a
flame front, which consists of the preheat and reaction zones [14, 26, 27, 3].
Therefore, the temperature dependence of reaction term was ignored and the
cold boundary difficulty was unresolved. Instead, in this work, the reaction
term is expressed by the Arrhenius type and its temperature dependence is
accounted for [23, 60]. By employing the large activation energy asymptotics,
the cold boundary difficulty is resolved [63].

In this study, the order of Mach number is assumed to be squared root of the
inverse of non-dimensional activation energy. Thanks to this specific order,
under the large values of activation energy, the effect of gas compression is
able to be captured in the reaction zone as well as in the compression zone.
As a result, the condition of temperature gradient at the burned-side edge of
a flame front is simply described by the difference between the flame
temperature and its adiabatic value. Besides, the burning-rate eigenvalue
contains the influence of not only thermal-diffusive properties but also
hydrodynamic flow. In this case, the equations of heat conduction and burning-
rate eigenvalue need to be computed simultaneously to obtain the solutions of
temperature and laminar flame speed.

The non-zero temperature gradient for small Mach numbers brings two specific
features of spatial temperature distribution. For the positive temperature
gradient, the monotonic increase of temperature is observed, which is similar
with that in the adiabatic case. In this case, the augmentation of viscous
effect, which is denoted by the Prandtl number $Pr$, induces the drop of
temperature. On the other hand, in the case of negative temperature gradient,
the non-monotonic temperature profile is captured with its maximum value taken
inside the reaction zone. Then, the increase of $Pr$ causes that of
temperature.

##  2 PRELIMINARIES

###  2.1 Governing Equations

We consider a planar premixed flame front freely propagating in the negative
z-direction. The coordinate $z$ is non-dimensionalised by the thermal
thickness of a flame front, or the length scale of preheat zone [52]. The
density $R$, velocity of gases $V$, pressure $P$, temperature $T$ and mass
fraction of deficient species $Y$ are made dimensionless by use of their
values at a position far from a flame front, $z\to-\infty$, denoted by
$\tilde{R}_{-\infty}$, $\tilde{S}_{L}$, $\tilde{P}_{-\infty}$,
$\tilde{T}_{-\infty}$ and $\tilde{Y}_{-\infty}$. Especially, $\tilde{S}_{L}$
is called the laminar flame speed, or the burning velocity, which is the
propagation velocity of a planar flame front into a quiescent fuel mixture
[63]. Hereafter, the variables with tilde symbol stand for dimensional ones.

$\displaystyle R$ $\displaystyle=\frac{\tilde{R}}{\tilde{R}_{-\infty}},$
$\displaystyle V$ $\displaystyle=\frac{\tilde{V}}{\tilde{S}_{L}},$
$\displaystyle P$ $\displaystyle=\frac{\tilde{P}}{\tilde{P}_{-\infty}},$
$\displaystyle T$ $\displaystyle=\frac{\tilde{T}}{\tilde{T}_{-\infty}},$
$\displaystyle Y$ $\displaystyle=\frac{\tilde{Y}}{\tilde{Y}_{-\infty}},$
$\displaystyle\tilde{l}_{d}$
$\displaystyle=\frac{\tilde{D}_{th}}{\tilde{S}_{L}},\quad$ $\displaystyle Le$
$\displaystyle=\frac{\tilde{D}_{th}}{\tilde{D}},\quad$ $\displaystyle Pr$
$\displaystyle=\frac{\tilde{\mu}\tilde{c}_{p}}{\tilde{\lambda}},\quad$
$\displaystyle\tilde{D}_{th}$
$\displaystyle=\frac{\tilde{\lambda}}{\tilde{R}_{-\infty}\tilde{c}_{p}},\quad$
$\displaystyle Ma$ $\displaystyle=\frac{\tilde{S}_{L}}{\tilde{c}_{s}},\ \ $
(1)

where $\tilde{l}_{d}$, $\tilde{D}_{th}$, $\tilde{D}$, $\tilde{\mu}$,
$\tilde{c}_{p}$ and $\tilde{\lambda}$ are the length scale of flame front, the
thermal diffusivity, the diffusion coefficient of a deficient species, the
coefficient of viscosity, the specific heat at constant pressure of the
mixture and the thermal conductivity. We assume that $\tilde{D}$,
$\tilde{\mu}$, $\tilde{c}_{p}$ and $\tilde{\lambda}$ are all constant. The
adiabatic sound speed is defined by
$\tilde{c}_{s}=(\gamma\tilde{P}_{-\infty}/\tilde{R}_{-\infty})^{1/2}$. The
specific heat ratio, Mach number, Prandtl number and Lewis number are denoted
by $\gamma$, $Ma$, $Pr$ and $Le$, respectively.

Because of (1), we have the following far-field boundary conditions.

$\displaystyle R=V=P=T=Y=1\quad\mbox{as}\quad z\to-\infty.$ (2)

For a steady one-dimensional flow, the conservation law of mass tells us that
the mass flux $M$ is constant in the entire region including the inside of
flame front. Then, the boundary condition (2) leads to

$\displaystyle M=RV=1.$ (3)

Accounting for (3), we have the following non-dimensional governing equations.

$\displaystyle\frac{dV}{dz}$
$\displaystyle=-\frac{1}{\gamma{Ma}^{2}}\frac{dP}{dz}+\frac{4}{3}Pr\frac{d^{2}V}{dz^{2}},$
(4) $\displaystyle\frac{dT}{dz}$
$\displaystyle=\frac{d^{2}T}{dz^{2}}+\frac{\gamma-1}{\gamma}V\frac{dP}{dz}+\frac{4}{3}(\gamma-1)Pr{Ma}^{2}\left(\frac{dV}{dz}\right)^{2}+qQ,$
(5) $\displaystyle\frac{dY}{dz}$
$\displaystyle=\frac{1}{Le}\frac{d^{2}Y}{dz^{2}}-Q,$ (6) $\displaystyle P$
$\displaystyle=RT,$ (7)

where $q=\tilde{q}\tilde{Y}_{-\infty}/(\tilde{c}_{p}\tilde{T}_{-\infty})$ and
$Q$ are the non-dimensional heat release and reaction rate with $\tilde{q}$
being the dimensional heat release. We assume that the deficient species is
depleted inside a flame front and does not exist on the burned side.

$\displaystyle Y\equiv 0\quad(z\geq 0).$ (8)

We note that the assumption (8) contains the equality in the range of $z$.
This means that the burned-side edge of a flame front is located at $z=0$.

The reaction rate $Q$ is assumed to be the Arrhenius type:

$Q=RY\frac{\tilde{A}\tilde{D}_{th}}{\tilde{S}_{L}^{2}}\exp\left(-\frac{N}{T}\right),$
(9)

where $\tilde{A}$ is the frequency factor and assumed to be constant. The non-
dimensional activation energy is expressed by
$N=\tilde{E}/\tilde{R}_{g}\tilde{T}_{-\infty}$ with the dimensional activation
energy $\tilde{E}$ and the universal gas constant $\tilde{R}_{g}$. In order to
avoid the cold boundary difficulty in the limit of large activation energy
$N\gg 1$, $Q$ is transformed into the following form [44, 63].

$\displaystyle Q$ $\displaystyle=\frac{\Lambda
RY}{\delta}\left(\frac{N}{T_{b}^{2}}\right)^{2}\exp\left(N\frac{T-T_{b}}{T_{b}T}\right),$
(10)

where

$\displaystyle\Lambda=\frac{\tilde{A}\tilde{D}_{th}}{\tilde{S}_{L}^{2}}\left(\frac{T_{b}^{2}}{N}\right)^{2}\exp\left(-\frac{N}{T_{b}}\right).$
(11)

The coefficient $\Lambda$ plays a role of the eigenvalue which determines the
laminar flame speed $\tilde{S}_{L}$ and is called the burning-rate eigenvalue
[65, 66, 22, 58, 6, 63]. Under the large-activation-energy asymptotics, $N\gg
1$, we introduce a small parameter $\epsilon$, which characterizes the
thickness of reaction zone.

$\displaystyle\epsilon=T_{b}^{2}/N=(T_{b}-1)/\beta,$ (12)

where the Zel’dovich number $\beta$ [23] is defined by

$\displaystyle\beta=\frac{\tilde{E}}{\tilde{R}_{g}\tilde{T}_{b}}\frac{\tilde{T}_{b}-\tilde{T}_{-\infty}}{\tilde{T}_{b}}.$
(13)

In comparing our result with that of previous works, it may be helpful to
rewrite $\epsilon$ into the form expressed by the Zel’dovich number. The flame
temperature $T_{b}$ is evaluated at the burned-side edge of a flame front
under $\epsilon\ll 1$, and is different from the adiabatic flame temperature
$T_{ad}$.

######  Remark 1.

The essential assumption, which determines the adiabatic profile of a flame
front, is that the temperature gradient is zero at $z=0$. Then, the adiabatic
flame temperature is given by $T_{ad}=1+q$ [53, 44].

In this study, the above assumption of zero-gradient of temperature at $z=0$
is removed. As a result, the non-monotonic behavior of temperature is captured
for its negative gradients with the maximum value of temperature taken inside
the reaction zone. Such a temperature profile may be observed in the context
of thermonuclear explosions [24, 33, 50, 56].

By use of (12), the burning-rate eigenvalue (11) is rewritten as

$\displaystyle\Lambda=\frac{\tilde{A}\tilde{D}_{th}}{\tilde{S}_{L}^{2}}\epsilon^{2}\exp\left(-\frac{T_{b}}{\epsilon}\right).$
(14)

We seek for the burning-rate eigenvalue in the form of power series with
respect to $\epsilon\;(\ll 1)$.

$\Lambda=\Lambda_{0}+O(\epsilon).$ (15)

All quantities are also supposed to be expanded asymptotically with respect to
$\epsilon$. In the outer regions ($z<0,z>0$),

$\displaystyle R$ $\displaystyle={R}_{0}+\epsilon{R}_{1},\quad$ $\displaystyle
V$ $\displaystyle={V}_{0}+\epsilon{V}_{1},\quad$ $\displaystyle P$
$\displaystyle={P}_{0}+\epsilon{P}_{1},$ (16a) $\displaystyle T$
$\displaystyle={T}_{0}+\epsilon{T}_{1},$ $\displaystyle Y$
$\displaystyle={Y}_{0}+\epsilon{Y}_{1},$ (16b)

in the reaction zone ($-\infty<\eta<0$),

$\displaystyle R$ $\displaystyle=\hat{\rho}_{0}+\epsilon\hat{\rho}_{1},$
$\displaystyle V$ $\displaystyle=\hat{v}_{0}+\epsilon\hat{v}_{1},\;$
$\displaystyle P$ $\displaystyle=\hat{p}_{0}+\epsilon\hat{p}_{1},$ (17a)
$\displaystyle T$
$\displaystyle=\hat{\theta}_{0}+\epsilon\hat{\theta}_{1}+\epsilon^{2}\hat{\theta}_{2},\;$
$\displaystyle Y$
$\displaystyle=\hat{c}_{0}+\epsilon\hat{c}_{1}+\epsilon^{2}\hat{c}_{2},$ (17b)

and in the compression zone ($0<\xi<+\infty$),

$\displaystyle R$ $\displaystyle=\check{\rho}_{0}+\epsilon\check{\rho}_{1},$
$\displaystyle V$ $\displaystyle=\check{v}_{0}+\epsilon\check{v}_{1},\;$
$\displaystyle P$ $\displaystyle=\check{p}_{0}+\epsilon\check{p}_{1},$ (18a)
$\displaystyle T$
$\displaystyle=\check{\theta}_{0}+\epsilon\check{\theta}_{1}+\epsilon^{2}\check{\theta}_{2},\;$
$\displaystyle Y$
$\displaystyle=\check{c}_{0}+\epsilon\check{c}_{1}+\epsilon^{2}\check{c}_{2}.$
(18b)

###  2.2 The Order of Mach Number

In this study, we consider the following order of Mach number with positive
constant $b=O(1)$.

$\displaystyle{Ma}=b\sqrt{\epsilon}.$ (19)

The assumption (19) enables us to comprehend analytically the effect of
compressibility, or finite Mach number, on the temperature profile of a flame
front. Indeed, the pressure and dissipation terms in the heat-conduction
equation emerge as a consequence of non-zero values of Mach number, or $b>0$.
Then, the adiabatic profile of a flame front is greatly changed: for example,
the non-monotonic behavior of temperature is observed. Especially, the
variation of temperature in the compression zone is described by the inner
asymptotic solution expressed by Lambert $W$ function. The similar result is
obtained in [61] based on the strong ignition temperature approximation, which
ignored the temperature dependence of reaction rate. Contrary to the previous
work, we systematically deal with the reaction rate represented by the
Arrhenius type [54, 55] by use of the large activation energy asymptotics. As
a result, the burning-rate eigenvalue is found to be affected by the
compressibility. In other words, the laminar flame speed reflects the non-
adiabatic feature of a flame front.

###  2.3 Notation in Reaction and Compression Zones

We introduce the reaction zone inside a flame front and the compression zone
behind it based on the large activation energy asymptotics under $\epsilon\ll
1$. The inner variables, or the stretched variables, are introduced as

$\displaystyle z=\epsilon\eta\quad(-\infty<\eta<0),\quad
z={Ma}^{2}\xi=b^{2}\epsilon\xi\quad(0<\xi<+\infty),$ (20)

where $\eta$ and $\xi$ represent the coordinates in the reaction and
compression zones respectively, and $z$ is the variable in outer regions,
which consist of the unburned ($z<0$) and burned ($z>0$) sides. The study of
inner regions is made to calculate inner solutions, which are given by (17)
and (18), and the boundary conditions across them, which are in turn used to
seek for outer solutions.

####  2.3.1 Matching Conditions Between Outer and Inner Solutions

For convenience in the later calculation, we recall the matching conditions
between the inner and outer solutions in the overlapping regions, where
$\eta\to-\infty$ and $z\to 0-$ or $\xi\to+\infty$ and $z\to 0+$ as
$\epsilon\to 0$. Let $\Phi$ be any function and its asymptotic expansions with
respect to $\epsilon$ are assumed to be

$\displaystyle\Phi$
$\displaystyle={\Phi}_{0}(z)+\epsilon{\Phi}_{1}(z)+\epsilon^{2}{\Phi}_{2}(z)+O(\epsilon^{3})\quad$
$\displaystyle(z<0,z>0),$ $\displaystyle\Phi$
$\displaystyle=\hat{\phi}_{0}(\eta)+\epsilon\hat{\phi}_{1}(\eta)+\epsilon^{2}\hat{\phi}_{2}(\eta)+O(\epsilon^{3})\quad$
$\displaystyle(-\infty<\eta<0),$ $\displaystyle\Phi$
$\displaystyle=\check{\phi}_{0}(\xi)+\epsilon\check{\phi}_{1}(\xi)+\epsilon^{2}\check{\phi}_{2}(\xi)+O(\epsilon^{3})\quad$
$\displaystyle(0<\xi<+\infty).$

We assume that the inner and outer solutions and their derivatives are
continuous in the overlapping regions. Then, according to [21], the matching
conditions are given by, between the reaction and unburned zones,

$\displaystyle\hat{\phi}_{0}|_{\eta\to-\infty}$
$\displaystyle={\Phi}_{0}|_{z\to
0-},\quad\hat{\phi}_{1}|_{\eta\to-\infty}={\Phi}_{1}|_{z\to
0-}+\eta\frac{d{\Phi}_{0}}{dz}\Big{|}_{z\to 0-},$ (21a)
$\displaystyle\frac{d\hat{\phi}_{0}}{d\eta}\Big{|}_{\eta\to-\infty}$
$\displaystyle=0,\quad\quad\quad\;\frac{d\hat{\phi}_{1}}{d\eta}\Big{|}_{\eta\to-\infty}=\frac{d{\Phi}_{0}}{dz}\Big{|}_{z\to
0-},$ (21b)
$\displaystyle\frac{d\hat{\phi}_{2}}{d\eta}\Big{|}_{\eta\to-\infty}$
$\displaystyle=\frac{d{\Phi}_{1}}{dz}\Big{|}_{z\to
0-}+\eta\frac{d^{2}{\Phi}_{0}}{dz^{2}}\Big{|}_{z\to 0-},$ (21c)
$\displaystyle\frac{d^{2}\hat{\phi}_{1}}{d\eta^{2}}\Big{|}_{\eta\to-\infty}$
$\displaystyle=0,\quad\quad\quad\;\frac{d^{2}\hat{\phi}_{2}}{d\eta^{2}}\Big{|}_{\eta\to-\infty}=\frac{d^{2}{\Phi}_{0}}{dz^{2}}\Big{|}_{z\to
0-}$ (21d)

and between the compression and burned zones,

$\displaystyle\check{\phi}_{0}|_{\xi\to+\infty}$
$\displaystyle={\Phi}_{0}|_{z\to
0+},\quad\check{\phi}_{1}|_{\xi\to+\infty}={\Phi}_{1}|_{z\to
0+}+b^{2}\xi\frac{d{\Phi}_{0}}{dz}\Big{|}_{z\to 0+},$ (22a)
$\displaystyle\frac{d\check{\phi}_{0}}{d\xi}\Big{|}_{\xi\to+\infty}$
$\displaystyle=0,\quad\quad\quad\;\frac{d\check{\phi}_{1}}{d\xi}\Big{|}_{\xi\to+\infty}=b^{2}\frac{d{\Phi}_{0}}{dz}\Big{|}_{z\to
0+},$ (22b)
$\displaystyle\frac{d\check{\phi}_{2}}{d\xi}\Big{|}_{\xi\to+\infty}$
$\displaystyle=b^{2}\frac{d{\Phi}_{1}}{dz}\Big{|}_{z\to
0+}+b^{4}\xi\frac{d^{2}{\Phi}_{0}}{dz^{2}}\Big{|}_{z\to 0+},$ (22c)
$\displaystyle\frac{d^{2}\check{\phi}_{2}}{d\xi^{2}}\Big{|}_{\xi\to+\infty}$
$\displaystyle=b^{4}\frac{d^{2}{\Phi}_{0}}{dz^{2}}\Big{|}_{z\to 0+}.$ (22d)

In addition, we impose the following assumptions.

######  Remark 2.

We assume that all quantities are continuous at the burned-side edge of a
flame front, $z=0$.

$\displaystyle\hat{\phi}_{i}|_{\eta\to 0-}=\check{\phi}_{i}|_{\xi\to
0+}\quad(i=0,1).$ (23)

######  Remark 3.

We further assume that the derivatives of inner solutions at the order higher
than $\epsilon$ are continuous at $z=0$. This condition is not imposed on
leading-order solutions, because we assume inner solutions of temperature and
mass fraction is constant inside of reaction zone, and in this case, only
constant solutions are admissible in reaction zone.

$\displaystyle\frac{d\hat{\phi}_{i}}{d\eta}\Big{|}_{\eta\to
0-}=\frac{1}{b^{2}}\frac{d\check{\phi}_{i}}{d\xi}\Big{|}_{\xi\to
0+},\quad(i=1,2).$ (24)

For simplicity, we will use the following notation.

$\displaystyle{\Phi}_{0}|_{z\to
0\pm}={\Phi}_{0\pm},\quad\frac{d{\Phi}_{0}}{dz}\Big{|}_{z\to
0\pm}=\frac{d{\Phi}_{0}}{dz}\Big{|}_{\pm}.$

The same things are true for the higher order terms.

##  3 INNER SOLUTIONS

At first, we study the inner structure of compression zone by employing the
assumption of Mach number (19) and the matching conditions (22). After that,
the inner solutions of reaction zone are calculated with the help of matching
conditions (21). Then, the solutions in reaction and compression zones are
connected by use of continuity conditions (23) and (24) to derive the boundary
conditions for outer solutions.

According to previous works [58, 6, 45], we further impose the following
assumption.

######  Remark 4.

We assume the following solutions inside reaction zone ($-\infty<\eta<0$).

$\displaystyle\hat{\theta}_{0}=T_{b},\quad\hat{c}_{0}=0.$ (25)

This means that the temperature is sufficiently high in the entire region of
reaction zone and nearly equal to its value at the burned-side edge of a flame
front. Correspondingly, at the leading order of $\epsilon$, the mass fraction
also has its value on the burned side.

###  3.1 Compression Zone $(0<\xi<+\infty)$

####  3.1.1 $O(\epsilon^{0})$ Solutions in Compression Zone

From the assumption (8), we have

$\displaystyle\check{c}_{i}=0\quad(i=0,1,2\ldots).$ (26)

Collecting $O(\epsilon^{0})$ terms, governing equations (3)-(5) and (7) become

$\displaystyle\check{\rho}_{0}\check{v}_{0}$ $\displaystyle=1,$ (27)
$\displaystyle 0$
$\displaystyle=-\frac{1}{\gamma}\frac{d\check{p}_{0}}{d\xi}+\frac{4}{3}Pr\frac{d^{2}\check{v}_{0}}{d\xi^{2}},$
(28) $\displaystyle 0$
$\displaystyle=\frac{d^{2}\check{\theta}_{0}}{d\xi^{2}},$ (29)
$\displaystyle\check{p}_{0}$
$\displaystyle=\check{\rho}_{0}\check{\theta}_{0}.$ (30)

From (29), accounting for (21b), (23) and (25), we readily find that

$\displaystyle\check{\theta}_{0}=T_{b}.$ (31)

Integrating (28), with the help of (22), we get

$\displaystyle\check{p}_{0}=\nu\frac{d\check{v}_{0}}{d\xi}+P_{0+},$ (32)

where we put

$\displaystyle\nu=\frac{4}{3}\gamma Pr.$ (33)

We find that the pressure is constant for inviscid flow, $Pr=0$. Then, the
velocity also becomes constant and only the adiabatic profile of a premixed
flame front is observed.

From (22a), (27), (30) and (31), in the limit of $\xi\to+\infty$, we find

$\displaystyle\check{\theta}_{0}=T_{b}={P}_{0+}{V}_{0+}.$ (34)

Substituting (30) into (32), with the help of (34), we obtain the following
differential equation for $\check{v}_{0}$ ($>0$).

$\displaystyle\frac{\check{v}_{0}}{\check{v}_{0}-{V}_{0+}}\frac{d\check{v}_{0}}{d\xi}=-\frac{{P}_{0+}}{\nu}.$
(35)

By use of the upper branch of Lambert $W$ function [14], which is expressed by
$y=W_{0}\left(x\right)$ as the inverse relation of $y=x\mathrm{e}^{x}$ for
real numbers $x$ and $y$ ($x\geq-1/e,y\geq-1$), we get the solution of (35) as
follows [61].

$\displaystyle\check{v}_{0}={V}_{0+}\left(1+W\right),$ (36)

where

$\displaystyle W$
$\displaystyle=W_{0}\left(\alpha\exp\left(-\frac{1}{\nu}\frac{{P}_{0+}}{{V}_{0+}}\xi+\alpha\right)\right),$
(37) $\displaystyle\alpha$
$\displaystyle=W|_{\xi=0+}=\frac{\check{v}_{0+}}{{V}_{0+}}-1.$ (38)

Substituting (36) into (32), we obtain

$\displaystyle\check{p}_{0}={P}_{0+}(1+W)^{-1}.$ (39)

####  3.1.2 $O(\epsilon)$ Solutions in Compression Zone

At $O(\epsilon)$, the heat-conduction equation (5) reads

$\displaystyle 0$
$\displaystyle=\frac{d^{2}\check{\theta}_{1}}{d\xi^{2}}+b^{2}\frac{\gamma-1}{\gamma}\left(\check{v}_{0}\frac{d\check{p}_{0}}{d\xi}+\nu\left(\frac{d\check{v}_{0}}{d\xi}\right)^{2}\right).$
(40)

Accounting for (22b), (36) and (39), the integration of (40) leads to

$\displaystyle\frac{d\check{\theta}_{1}}{d\xi}=b^{2}\frac{\gamma-1}{\gamma}\check{\theta}_{0}W+b^{2}\frac{d{T}_{0}}{dz}\Big{|}_{+}.$
(41)

Integrating (41), by use of (22a), we get

$\displaystyle\check{\theta}_{1}={T}_{1+}+\xi
b^{2}\frac{d{T}_{0}}{dz}\Big{|}_{+}-b^{2}\frac{\gamma-1}{\gamma}\nu{V}_{0+}^{2}W\left(1+\frac{W}{2}\right).$
(42)

Next, we calculate $O(\epsilon)$ solutions of velocity and pressure in the
compression zone. Collecting $O(\epsilon)$ terms in governing equations (3),
(4) and (7), we have

$\displaystyle 0$
$\displaystyle=\check{\rho}_{1}\check{v}_{0}+\check{\rho}_{0}\check{v}_{1},$
(43) $\displaystyle\gamma b^{2}\frac{d\check{v}_{0}}{d\xi}$
$\displaystyle=-\frac{d\check{p}_{1}}{d\xi}+\nu\frac{d^{2}\check{v}_{1}}{d\xi^{2}},$
(44) $\displaystyle\check{p}_{1}$
$\displaystyle=\check{\rho}_{1}\check{\theta}_{0}+\check{\rho}_{0}\check{\theta}_{1}.$
(45)

From (27), (30), (43) and (45), $\check{p}_{1}$ is calculated by

$\displaystyle\check{p}_{1}=-\frac{\check{p}_{0}}{\check{v}_{0}}\check{v}_{1}+\frac{\check{\theta}_{1}}{\check{v}_{0}}.$
(46)

Substituting (46) into the integration of (44), by use of the matching
conditions (22a) and (22b) with the matching relation $dP_{0}/dz|_{+}=0$ which
is found by taking $\xi\to+\infty$ in (44), we obtain the following
differential equation for $\check{v}_{1}$.

$\displaystyle\nu\frac{d\check{v}_{1}}{d\xi}+\frac{\check{p}_{0}}{\check{v}_{0}}\check{v}_{1}=\frac{\check{\theta}_{1}}{\check{v}_{0}}+\gamma
b^{2}\check{v}_{0}+b^{2}\nu\frac{dV_{0}}{dz}\Big{|}_{+}-P_{1+}-\gamma
b^{2}V_{0+},$ (47)

Solving (47), by use of solutions (36), (39) and (42), we get

$\displaystyle\check{v}_{1}$
$\displaystyle=\frac{1}{{P}_{0+}}\frac{W}{1+W}\left\\{\check{v}_{1c}+({T}_{1+}-{P}_{1+}{V}_{0+})\left(\frac{1}{W}-\ln{|W|}\right)\right.$
$\displaystyle\left.\quad-b^{2}{V}_{0+}^{2}\frac{1+\gamma-B}{4}W(2+W)+\frac{{P}_{0+}}{\nu}(b^{2}{V}_{0+}(1-B)-{P}_{1+})\xi\right.$
$\displaystyle\left.\quad+b^{2}\frac{dT_{0}}{dz}\Big{|}_{+}\left(\left(\frac{1}{W}-1-\ln{|W|}\right)\xi-\nu\frac{V_{0+}}{P_{0+}}\left(W+\frac{\ln{|W|}}{2}\right)\ln{|W|}\right)\right.$
$\displaystyle\left.\quad-b^{2}V_{0+}\nu\frac{dV_{0}}{dz}\Big{|}_{+}(W+2\ln{|W|})\right\\},$
(48)

where

$\displaystyle B$ $\displaystyle=\frac{\gamma-1}{\gamma}(\nu-\gamma),$ (49)
$\displaystyle\check{v}_{1c}$
$\displaystyle=\frac{1+\alpha}{\alpha}{P}_{0+}\check{v}_{1+}-({T}_{1+}-{P}_{1+}{V}_{0+})\left(\frac{1}{\alpha}-\ln{|\alpha|}\right)$
$\displaystyle\quad+b^{2}{V}_{0+}^{2}\frac{1+\gamma-B}{4}\alpha(2+\alpha)+b^{2}\frac{V_{0+}}{P_{0+}}\nu\frac{dT_{0}}{dz}\Big{|}_{+}\left(\alpha+\frac{\ln{|\alpha|}}{2}\right)\ln{|\alpha|}$
$\displaystyle\quad+b^{2}V_{0+}\nu\frac{dV_{0}}{dz}\Big{|}_{+}(\alpha+2\ln{|\alpha|}).$
(50)

Solutions (46) and (48) are used to calculate the heat-conduction equation at
$O(\epsilon^{2})$ in Section 3.1.3. Besides, they are matched with inner
solutions of reaction zone to determine the values of $T_{1+}$ and $P_{1+}$.

####  3.1.3 Heat-Conduction Equation at $O(\epsilon^{2})$ in Compression Zone

We consider the heat-conduction equation (5) at $O(\epsilon^{2})$ as follows.

$\displaystyle
b^{2}\frac{d\check{\theta}_{1}}{d\xi}=\frac{d^{2}\check{\theta}_{2}}{d\xi^{2}}+b^{2}\frac{\gamma-1}{\gamma}\left(\check{v}_{1}\frac{d\check{p}_{0}}{d\xi}+\check{v}_{0}\frac{d\check{p}_{1}}{d\xi}+2\nu\frac{d\check{v}_{0}}{d\xi}\frac{d\check{v}_{1}}{d\xi}\right).$
(51)

The integration of (51), by use of solutions (36), (39), (46) and (48), leads
to

$\displaystyle\frac{\gamma}{\gamma-1}\check{\theta}_{1}-\frac{1}{b^{2}}\frac{\gamma}{\gamma-1}\frac{d\check{\theta}_{2}}{d\xi}+\check{\theta}_{1c}$
$\displaystyle=\frac{\check{v}_{1c}}{1+W}+T_{1+}\frac{W}{1+W}\left(1+\ln{|W|}\right)-P_{1+}V_{0+}\left\\{\frac{W\ln{|W|}-1}{1+W}+2(1+W)\right.$
$\displaystyle\left.\quad+\ln{|W|}+\frac{P_{0+}}{\nu
V_{0+}}\frac{\xi}{1+W}\right\\}+b^{2}V_{0+}^{2}(1-B)\left\\{\frac{1}{4(1+W)}+\frac{3}{4}W(3+W)\right.$
$\displaystyle\left.\quad+\ln{|W|}+\frac{P_{0+}}{\nu
V_{0+}}\frac{\xi}{1+W}\right\\}+b^{2}V_{0+}^{2}\frac{\gamma}{4}\left(\frac{1}{1+W}+(W+1)(W-4)\right)$
$\displaystyle\quad+b^{2}\frac{V_{0+}}{P_{0+}}\nu\frac{dT_{0}}{dz}\Big{|}_{+}\left\\{\frac{1}{1+W}-W+\frac{W\ln{|W|}}{1+W}\left(W+\frac{\ln{|W|}}{2}\right)-\ln{|W|}\right.$
$\displaystyle\left.\quad+\frac{P_{0+}}{\nu
V_{0+}}\frac{W\ln{|W|}+W-1}{1+W}\xi\right\\}+b^{2}V_{0+}\nu\frac{dV_{0}}{dz}\Big{|}_{+}\left(3+2W+\frac{2W\ln{|W|}}{1+W}\right).$
(52)

The integral constant $\check{\theta}_{1c}$ is determined by taking
$\xi\to+\infty$ in (52), with the help of matching condition (22) while being
careful about $dT_{0}/dz|_{+}=d^{2}T_{0}/dz^{2}|_{+}$ which is obtained by
taking $\xi\to+\infty$ in (51), as

$\displaystyle\check{\theta}_{1c}$
$\displaystyle=\frac{1+\alpha}{\alpha}P_{0+}\check{v}_{1+}-{T}_{1+}\left(\frac{\gamma}{\gamma-1}+\frac{1}{\alpha}-\ln{|\alpha|}\right)+\frac{\gamma}{\gamma-1}\frac{dT_{1}}{dz}\Big{|}_{+}$
$\displaystyle\quad-
P_{1+}V_{0+}\left(1+\alpha-\frac{1}{\alpha}+2\ln{|\alpha|}\right)+b^{2}{V}_{0+}^{2}\frac{1-B}{4}\left(1+6\alpha+\alpha^{2}+4\ln{|\alpha|}\right)$
$\displaystyle\quad+b^{2}{V}_{0+}^{2}\frac{\gamma}{4}(\alpha-1)(\alpha+3)+b^{2}\frac{V_{0+}}{P_{0+}}\nu\frac{dT_{0}}{dz}\Big{|}_{+}\left((1-\alpha)(1-\ln{|\alpha|})+\frac{(\ln{|\alpha|})^{2}}{2}\right)$
$\displaystyle\quad+b^{2}V_{0+}\nu\frac{dV_{0}}{dz}\Big{|}_{+}(3+\alpha+2\ln{|\alpha|}),$
(53)

where we note that $\lim_{\xi\to+\infty}\xi W=0$ and $\lim_{\xi\to+\infty}\xi
W\ln{|W|}=0$, which are readily confirmed by l’Hôpital’s theorem.

######  Remark 5.

The terms $\ln|W|+\xi{P}_{0+}/(\nu V_{0+})$ in (52) converge to
$\alpha+\ln{|\alpha|}$ as $\xi\to+\infty$. In fact, the following relation
holds.

$\displaystyle\ln{|W|}+\frac{{P}_{0+}}{\nu V_{0+}}\xi=-W+\ln|\alpha|+\alpha.$
(54)

The relation (54) is obtained by taking logarithm of both sides of the
following equation.

$\displaystyle
W\mathrm{e}^{W}=\alpha\exp\left(-\frac{1}{\nu}\frac{{P}_{0+}}{{V}_{0+}}\xi+\alpha\right),$
(55)

which is the definition of Lambert $W$ function for the solution (36) [61].

###  3.2 Reaction Zone $(-\infty<\eta<0)$

####  3.2.1 $O(\epsilon^{0})$ Solutions in Reaction Zone

Collecting $O(\epsilon^{0})$ terms, governing equations (3), (4) and (7)
become

$\displaystyle\hat{\rho}_{0}\hat{v}_{0}$ $\displaystyle=1,$ (56)
$\displaystyle 0$
$\displaystyle=-\frac{d\hat{p}_{0}}{d\eta}+b^{2}\nu\frac{d^{2}\hat{v}_{0}}{d\eta^{2}},$
(57) $\displaystyle\hat{p}_{0}$
$\displaystyle=\hat{\rho}_{0}\hat{\theta}_{0}.$ (58)

In the same procedure as Section 3.1.1, we have the following system of
differential equations in the reaction zone.

$\displaystyle\frac{\hat{v}_{0}}{\hat{v}_{0}-{V}_{0-}}\frac{d\hat{v}_{0}}{d\eta}$
$\displaystyle=-\frac{{P}_{0-}}{b^{2}\nu},$ (59) $\displaystyle\hat{p}_{0}$
$\displaystyle=b^{2}\nu\frac{d\hat{v}_{0}}{d\eta}+{P}_{0-}.$ (60)

Because $\hat{v}_{0}$ should be bounded in space, we find that only constant
solutions are admissible for (59) and (60) as follows.

$\displaystyle\hat{v}_{0}={V}_{0-},\quad\hat{p}_{0}={P}_{0-}.$ (61)

Due to (25), (58) and (61), the density at $O(\epsilon^{0})$ is also constant.

$\displaystyle\hat{\rho}_{0}={P}_{0-}/T_{b}.$ (62)

####  3.2.2 $O(\epsilon)$ Solutions in Reaction Zone

At $O(\epsilon)$, the governing equations (3)-(7) are written down, while
being careful about the constancy of leading-order solutions (61), as

$\displaystyle 0$
$\displaystyle=\hat{\rho}_{1}\hat{v}_{0}+\hat{\rho}_{0}\hat{v}_{1},$ (63)
$\displaystyle 0$
$\displaystyle=-\frac{d\hat{p}_{1}}{d\eta}+b^{2}\nu\frac{d^{2}\hat{v}_{1}}{d\eta^{2}},$
(64) $\displaystyle 0$
$\displaystyle=\frac{d^{2}\hat{\theta}_{1}}{d\eta^{2}}+q\hat{\Lambda}_{0}\hat{\rho}_{0}\hat{c}_{1}\mathrm{e}^{\hat{\theta}_{1}},$
(65) $\displaystyle 0$
$\displaystyle=\frac{1}{Le}\frac{d^{2}\hat{c}_{1}}{d\eta^{2}}-\hat{\Lambda}_{0}\hat{\rho}_{0}\hat{c}_{1}\mathrm{e}^{\hat{\theta}_{1}},$
(66) $\displaystyle\hat{p}_{1}$
$\displaystyle=\hat{\rho}_{1}\hat{\theta}_{0}+\hat{\rho}_{0}\hat{\theta}_{1}.$
(67)

At first, we calculate the temperature gradient at $O(\epsilon)$,
$d\hat{\theta}_{1}/d\eta$. To that end, we combine (65) and (66) to eliminate
the reaction terms.

$\displaystyle
0=\frac{d^{2}}{d\eta^{2}}\left(\hat{\theta}_{1}+\frac{q}{Le}\hat{c}_{1}\right).$
(68)

Integrating (68), with the help of (21b), we get

$\displaystyle\frac{d}{d\eta}\left(\hat{\theta}_{1}+\frac{q}{Le}\hat{c}_{1}\right)=\frac{d}{dz}\left({T}_{0}+\frac{q}{Le}{Y}_{0}\right)\Big{|}_{-}.$
(69)

Accounting for (21a), further integration of (69) leads to

$\displaystyle\hat{\theta}_{1}+\frac{q}{Le}\hat{c}_{1}={T}_{1-}+\frac{q}{Le}{Y}_{1-}+\eta\frac{d}{dz}\left({T}_{0}+\frac{q}{Le}{Y}_{0}\right)\Big{|}_{-}.$
(70)

By use of (70), the heat-conduction equation (65) is rewritten as

$\displaystyle
0=\frac{d^{2}\hat{\theta}_{1}}{d\eta^{2}}+Le\hat{\Lambda}_{0}\hat{\rho}_{0}\left({T}_{1-}+\frac{q}{Le}{Y}_{1-}-\hat{\theta}_{1}+\eta\frac{d}{dz}\left({T}_{0}+\frac{q}{Le}{Y}_{0}\right)\Big{|}_{-}\right)\mathrm{e}^{\hat{\theta}_{1}}.$
(71)

Multiplying $2d\hat{\theta}_{1}/d\eta$ to (71), we obtain

$\displaystyle
0=\frac{d}{d\eta}\left(\frac{d\hat{\theta}_{1}}{d\eta}\right)^{2}-2Le\hat{\Lambda}_{0}\hat{\rho}_{0}\frac{d}{dz}\left({T}_{0}+\frac{q}{Le}{Y}_{0}\right)\Big{|}_{-}\mathrm{e}^{\hat{\theta}_{1}}$
$\displaystyle+2Le\hat{\Lambda}_{0}\hat{\rho}_{0}\frac{d}{d\eta}\left\\{\left(1+{T}_{1-}+\frac{q}{Le}{Y}_{1-}-\hat{\theta}_{1}+\eta\frac{d}{dz}\left({T}_{0}+\frac{q}{Le}{Y}_{0}\right)\Big{|}_{-}\right)\mathrm{e}^{\hat{\theta}_{1}}\right\\}.$
(72)

Integrating (72), accounting for matching conditions (21a) and (21b), we
obtain

$\displaystyle\left(\frac{d\hat{\theta}_{1}}{d\eta}\right)^{2}=\left(\frac{d{T}_{0}}{dz}\Big{|}_{-}\right)^{2}+2Le\hat{\Lambda}_{0}\hat{\rho}_{0}\frac{d}{dz}\left({T}_{0}+\frac{q}{Le}{Y}_{0}\right)\Big{|}_{-}\int_{-\infty}^{\eta}\mathrm{e}^{\hat{\theta}_{1}}d\eta$
$\displaystyle-2Le\hat{\Lambda}_{0}\hat{\rho}_{0}\left(1+{T}_{1-}+\frac{q}{Le}{Y}_{1-}-\hat{\theta}_{1}+\eta\frac{d}{dz}\left({T}_{0}+\frac{q}{Le}{Y}_{0}\right)\Big{|}_{-}\right)\mathrm{e}^{\hat{\theta}_{1}},$
(73)

where we note $\hat{\theta}_{1}|_{-\infty}=T_{1-}+\eta
dT_{0}/dz|_{-}\to-\infty$ as $\eta\to-\infty$, because it is natural to
consider that there exists temperature gradient on the unburned side of
reaction front, $dT_{0}/dz|_{-}\neq 0$. Taking $\eta=0-$ in (73), with
$\hat{\theta}_{1-}={T}_{1-}+q{Y}_{1-}/Le$ due to (23), (26) and (70), we get

$\displaystyle 0$
$\displaystyle=\left(\frac{\gamma-1}{\gamma}T_{b}\alpha+\frac{d{T}_{0}}{dz}\Big{|}_{+}\right)^{2}-\left(\frac{d{T}_{0}}{dz}\Big{|}_{-}\right)^{2}$
$\displaystyle\quad+2Le\hat{\Lambda}_{0}\hat{\rho}_{0}\left(\mathrm{e}^{\hat{\theta}_{1-}}-\frac{d}{dz}\left({T}_{0}+\frac{q}{Le}{Y}_{0}\right)\Big{|}_{-}I\right),$
(74)

with

$\displaystyle I=\int_{-\infty}^{0-}\mathrm{e}^{\hat{\theta}_{1}}d\eta.$ (75)

The first term in (74) is due to (24), (31), (41) and (38). By substituting
outer solutions calculated in Section 4 into (74), the leading order term of
burning-rate eigenvalue $\hat{\Lambda}_{0}$ is calculated in Section 5.

Next, we calculate $\hat{v}_{1}$ and $\hat{p}_{1}$. Once $\hat{v}_{1}$ is
obtained, $\hat{\rho}_{1}$ is calculated from (63) as

$\displaystyle\hat{\rho}_{1}=-\frac{\hat{\rho}_{0}}{\hat{v}_{0}}\hat{v}_{1}.$
(76)

Accounting for $dP_{0}/dz|_{-}=0$, which is found by taking $\eta\to-\infty$
in (64) with the help of (21d), the integration of (64) tells us $\hat{p}_{1}$
is calculated by

$\displaystyle\hat{p}_{1}=b^{2}\nu\frac{d\hat{v}_{1}}{d\eta}+P_{1-}-b^{2}\nu\frac{dV_{0}}{dz}\Big{|}_{-},$
(77)

where (21a) and (21b) are employed. Differentiating (67) with respect to
$\eta$, by use of (58) and (76), we have

$\displaystyle\frac{d\hat{p}_{1}}{d\eta}=-\frac{\hat{p}_{0}}{\hat{v}_{0}}\frac{d\hat{v}_{1}}{d\eta}+\hat{\rho}_{0}\frac{d\hat{\theta}_{1}}{d\eta}.$
(78)

Substituting (78) into (64), we obtain the following differential equation for
$\hat{v}_{1}$.

$\displaystyle\frac{d^{2}\hat{v}_{1}}{d\eta^{2}}+\frac{1}{b^{2}\nu}\frac{\hat{p}_{0}}{\hat{v}_{0}}\frac{d\hat{v}_{1}}{d\eta}=\frac{\hat{\rho}_{0}}{b^{2}\nu}\frac{d\hat{\theta}_{1}}{d\eta}.$
(79)

Putting $\lambda=\hat{p}_{0}/(b^{2}\nu\hat{v}_{0})$, we get

$\displaystyle\hat{v}_{1}=V_{1-}-\frac{T_{1-}}{P_{0-}}+\frac{1}{\lambda
P_{0-}}\frac{dT_{0}}{dz}\Big{|}_{-}+\frac{\hat{\theta}_{1}}{P_{0-}}-\frac{\mathrm{e}^{-\lambda\eta}}{P_{0-}}\int_{-\infty}^{\eta}\mathrm{e}^{\lambda\eta}\frac{d\hat{\theta}_{1}}{d\eta}d\eta,$
(80)

where we required $\hat{v}_{1}$ is bounded in space. In determining integral
constant, we note matching relation
$dV_{0}/dz|_{-}=(\hat{\rho}_{0}\hat{v}_{0}/\hat{p}_{0})dT_{0}/dz|_{-}$, which
is obtained by taking $\eta\to-\infty$ in (78), and the limit
$\lim_{\eta\to-\infty}\mathrm{e}^{-\lambda\eta}\int_{-\infty}^{\eta}\mathrm{e}^{\lambda\eta}(d\hat{\theta}_{1}/d\eta)d\eta=\lambda^{-1}dT_{0}/dz|_{-}$,
which is confirmed by l’Hôpital’s theorem. Then, from (77) and (80),
$\hat{p}_{1}$ is calculated as

$\displaystyle\hat{p}_{1}=P_{1-}-b^{2}\nu\frac{dV_{0}}{dz}\Big{|}_{-}+\frac{\mathrm{e}^{-\lambda\eta}}{V_{0-}}\int_{-\infty}^{\eta}\mathrm{e}^{\lambda\eta}\frac{d\hat{\theta}_{1}}{d\eta}d\eta.$
(81)

####  3.2.3 Heat-Conduction and Mass-Fraction Equations at $O(\epsilon^{2})$
in Reaction Zone

At $O(\epsilon^{2})$, with the help of (61), we consider the combination of
heat-conduction equation (5) and mass-fraction equation (6), which is written
as

$\displaystyle\frac{d}{d\eta}(\hat{\theta}_{1}+q\hat{c}_{1})=\frac{d^{2}}{d\eta^{2}}\left(\hat{\theta}_{2}+\frac{q}{Le}\hat{c}_{2}\right)+\frac{\gamma-1}{\gamma}\hat{v}_{0}\frac{d\hat{p}_{1}}{d\eta}.$
(82)

The integration of (82) leads to

$\displaystyle\hat{\theta}_{1}+q\hat{c}_{1}$
$\displaystyle=\frac{d}{d\eta}\left(\hat{\theta}_{2}+\frac{q}{Le}\hat{c}_{2}\right)+\frac{\gamma-1}{\gamma}\hat{v}_{0}\hat{p}_{1}+{T}_{1-}+q{Y}_{1-}$
$\displaystyle\quad-\frac{d}{dz}\left({T}_{1}+\frac{q}{Le}{Y}_{1}\right)\Big{|}_{-}-\frac{\gamma-1}{\gamma}{V}_{0-}{P}_{1-},$
(83)

where the matching conditions (21a) and (21c) are used to determine the
integral constant in the limit of $\eta\to-\infty$, while being careful about
the matching relation
$d(T_{0}+qY_{0})/dz|_{-}=d^{2}(T_{0}+qY_{0}/Le)/dz^{2}|_{-}$, which is
obtained by taking $\eta\to-\infty$ in (82) by use of (21d) and
$dP_{0}/dz|_{-}=0$.

###  3.3 Boundary Conditions Across Reaction and Compression Zones

We summarize the boundary conditions across the reaction and compression
zones, which are obtained subject to assumptions (23) and (24). From (25),
(26), (31), (41), (42), (69), (47), (70) and (77), by taking $\eta=0-$ or
$\xi=0+$, we find

$\displaystyle[{T}_{0}]=[{Y}_{0}]$ $\displaystyle=0,$ (84)
$\displaystyle\left[\frac{d}{dz}\left({T}_{0}+\frac{q}{Le}{Y}_{0}\right)\right]$
$\displaystyle=-\frac{\gamma-1}{\gamma}T_{b}\alpha,$ (85)
$\displaystyle\left[{T}_{1}+\frac{q}{Le}{Y}_{1}\right]$
$\displaystyle=b^{2}\frac{\gamma-1}{\gamma}\nu{V}_{0+}^{2}\alpha\left(1+\frac{\alpha}{2}\right),$
(86) $\displaystyle[P_{1}]$
$\displaystyle=b^{2}\left(\nu\left[\frac{dV_{0}}{dz}\right]+\gamma\alpha\right),$
(87)

where the brackets $[\;]$ means the difference of solutions evaluated at the
burned-side edge of compression zone and the unburned-side edge of reaction
zone, $[\Phi(z)]=\Phi_{+}-\Phi_{-}$. The constant $\alpha$ is defined by (38)
and now found to be $\alpha=V_{0-}/V_{0+}-1$.

From (23), (26), (70) and (80), we find

$\displaystyle\check{v}_{1+}=\hat{v}_{1-}=V_{1-}+\frac{q}{Le}\frac{Y_{1-}}{P_{0-}}+\frac{1}{\lambda
P_{0-}}\frac{dT_{0}}{dz}\Big{|}_{-}-\frac{1}{P_{0-}}\int_{-\infty}^{0-}\mathrm{e}^{\lambda\eta}\frac{d\hat{\theta}_{1}}{d\eta}d\eta,$
(88)

Taking $\xi=0+$ in (52) and $\eta=0-$ in (83), by use of (23) and (24), we
obtain

$\displaystyle 0$
$\displaystyle=\frac{\gamma}{\gamma-1}\left\\{\left[\frac{d}{dz}\left(T_{1}+\frac{q}{Le}Y_{1}\right)\right]+q\left(1-\frac{1}{Le}\right)Y_{1-}\right\\}$
$\displaystyle\quad+(1+\alpha)[V_{0}P_{1}]+\alpha
T_{1+}+\frac{\alpha}{\lambda}\frac{dT_{0}}{dz}\Big{|}_{-}-\alpha\int_{-\infty}^{0-}\mathrm{e}^{\lambda\eta}\frac{d\hat{\theta}_{1}}{d\eta}d\eta$
$\displaystyle\quad-b^{2}V_{0+}^{2}\alpha\left\\{\frac{\gamma}{2}\alpha+\nu\left(1+\frac{\gamma-1}{\gamma}\alpha\right)\left(1+\frac{\alpha}{2}\right)+\frac{\nu}{V_{0+}}\frac{dV_{0}}{dz}\Big{|}_{+}\right\\}.$
(89)

##  4 OUTER REGIONS

In the previous section, we have obtained the inner solutions in the reaction
and compression zones and the boundary conditions across them. In this
section, we seek for the solutions in the outer regions, which consist of the
unburned and burned zones. Now, the governing equations (4)-(6) in the outer
regions ($z<0,z>0$) are written, by use of (19), as

$\displaystyle b^{2}\gamma\frac{dV}{dz}$
$\displaystyle=-\frac{1}{\epsilon}\frac{dP}{dz}+b^{2}\nu\frac{d^{2}V}{dz^{2}},$
(90) $\displaystyle\frac{dT}{dz}$
$\displaystyle=\frac{d^{2}T}{dz^{2}}+\frac{\gamma-1}{\gamma}\left(V\frac{dP}{dz}+\epsilon
b^{2}\nu\left(\frac{dV}{dz}\right)^{2}\right),$ (91)
$\displaystyle\frac{dY}{dz}$
$\displaystyle=\frac{1}{Le}\frac{d^{2}Y}{dz^{2}},$ (92)

where the reaction term $Q$ is omitted from (91) and (92) due to the large-
activation-energy limit and $T<T_{b}$ for $z<0$ and due to $Y=0$ for $z>0$.

###  4.1 Outer Solutions

We seek for solutions of governing equations (3), (7), (90)-(92) in the form
of asymptotic expansions with respect to $\epsilon$ as introduced by (16).
From (90), with the help of (2), we find that the pressure at
$O(\epsilon^{0})$ is constant.

$\displaystyle{P}_{0}=\left\\{\begin{array}[]{lll}1&(z<0)\\\
{P}_{0+}&(z>0)\end{array}\right..$ (95)

Taking account of (95), the solutions of other quantities at $O(\epsilon^{0})$
and pressure at $O(\epsilon)$ are obtained by solving the following equations.

$\displaystyle{R}_{0}{V}_{0}$ $\displaystyle=1,$ (96) $\displaystyle
b^{2}\gamma\frac{d{V}_{0}}{dz}$
$\displaystyle=-\frac{d{P}_{1}}{dz}+b^{2}\nu\frac{d^{2}{V}_{0}}{dz^{2}},$ (97)
$\displaystyle\frac{d{T}_{0}}{dz}$
$\displaystyle=\frac{d^{2}{T}_{0}}{dz^{2}},$ (98)
$\displaystyle\frac{d{Y}_{0}}{dz}$
$\displaystyle=\frac{1}{Le}\frac{d^{2}{Y}_{0}}{dz^{2}},$ (99)
$\displaystyle{P}_{0}$ $\displaystyle={R}_{0}{T}_{0}.$ (100)

Boundary conditions (2), (8), (84) and (87) are imposed on (96)-(100). Then,
the solutions in the outer regions are given as follows.

$\displaystyle{T}_{0}$
$\displaystyle=\left\\{\begin{array}[]{lll}1+(T_{b}-1)\mathrm{e}^{z}&(z<0)\\\
T_{b}&(z>0)\end{array}\right.,\quad{Y}_{0}=\left\\{\begin{array}[]{lll}1-\mathrm{e}^{Lez}&(z<0)\\\
0&(z>0)\end{array}\right.,$ (105) $\displaystyle{V}_{0}$
$\displaystyle=\left\\{\begin{array}[]{lll}{T}_{0}&(z<0)\\\
V_{0+}&(z>0)\end{array}\right.,\quad\quad{R}_{0}={V}_{0}^{-1},$ (108)
$\displaystyle{P}_{0}$ $\displaystyle=\left\\{\begin{array}[]{lll}1&(z<0)\\\
T_{b}V_{0+}^{-1}&(z>0)\end{array}\right.,\quad{P}_{1}=\left\\{\begin{array}[]{lll}b^{2}(\nu-\gamma)(T_{b}-1)\mathrm{e}^{z}&(z<0)\\\
b^{2}\gamma(1+\alpha-T_{b})&(z>0)\end{array}\right.,$ (113)

We find that the following relation holds due to (85) and (105).

$\displaystyle\frac{\gamma-1}{\gamma}T_{b}\alpha=T_{b}-(1+q)=T_{b}-T_{ad}.$
(114)

Then, by use of (38), we get

$\displaystyle
V_{0+}=\frac{(\gamma-1)T_{b}^{2}}{(2\gamma-1)T_{b}-\gamma(1+q)}.$ (115)

We note that if $T_{b}$ takes its adiabatic value $T_{ad}=1+q$, the adiabatic
profile of a premixed flame front is recovered. However, in the present case,
its value is not restricted to adiabatic one, and the non-zero temperature
gradient is described by the deviation of $T_{b}$ from $T_{ad}$ as will be
shown in (128).

At $O(\epsilon)$, governing equations (91) and (92) are written down as

$\displaystyle\frac{d{T}_{1}}{dz}$
$\displaystyle=\frac{d^{2}{T}_{1}}{dz^{2}}+\frac{\gamma-1}{\gamma}\left({V}_{0}\frac{d{P}_{1}}{dz}+b^{2}\nu\left(\frac{dV_{0}}{dz}\right)^{2}\right),$
(116) $\displaystyle\frac{d{Y}_{1}}{dz}$
$\displaystyle=\frac{1}{Le}\frac{d^{2}{Y}_{1}}{dz^{2}}.$ (117)

Solving (116) and (117), by use of (2) and (8), we obtain

$\displaystyle{T}_{1}$
$\displaystyle=\left\\{\begin{array}[]{lll}C_{T}\mathrm{e}^{z}-b^{2}\frac{\gamma-1}{\gamma}(T_{b}-1)\big{\\{}(\nu-\gamma)z+\frac{T_{b}-1}{2}(2\nu-\gamma)\mathrm{e}^{z}\big{\\}}\mathrm{e}^{z}&(z<0)\\\
{T}_{1+}&(z>0)\end{array}\right.,$ (120) $\displaystyle{Y}_{1}$
$\displaystyle=\left\\{\begin{array}[]{lll}C_{Y}\mathrm{e}^{Lez}&(z<0)\\\
0&(z>0)\end{array}\right.,$ (123)

where $C_{T}$ and $C_{Y}$ are integral constants.

The value of $T_{1+}$ is calculated from (89), with (108), (113), (120) and
(123) substituted, as follows.

$\displaystyle T_{1+}$
$\displaystyle=\frac{\gamma-1}{\gamma}\frac{T_{b}}{T_{ad}}\left\\{\frac{b^{2}}{2}\left(\gamma+T_{b}^{2}\frac{\gamma(\alpha^{3}-1)-(B-1)\alpha^{2}(\alpha+2)}{(1+\alpha)^{2}}\right)\right.$
$\displaystyle\left.\quad-\alpha\int_{-\infty}^{0-}\mathrm{e}^{\lambda\eta}\frac{d\hat{\theta}_{1}}{d\eta}d\eta\right\\},$
(124)

##  5 MAIN RESULTS

###  5.1 Formula of Laminar Flame Speed

The formula of laminar flame speed is obtained by calculating the burning-rate
eigenvalue [22, 58, 6]. The leading-order term of burning-rate eigenvalue
$\hat{\Lambda}_{0}$ is calculated by substituting (105) into (74), with the
help of (114), as

$\displaystyle\hat{\Lambda}_{0}=\frac{q}{Le\hat{\rho}_{0}}\left(T_{b}-1-\frac{q}{2}\right)\left(\mathrm{e}^{\hat{\theta}_{1-}}-(T_{b}-1-q)I\right)^{-1}.$
(125)

Recalling (12), (14) and (62), the parameter dependence of laminar flame speed
is obtained from (125) as

$\displaystyle\tilde{S}_{L}=\left(2\tilde{A}\tilde{D}_{th}Le\frac{T_{b}^{3}}{N^{2}}\mathrm{e}^{-N/T_{b}}\frac{\mathrm{e}^{\hat{\theta}_{1-}}-(T_{b}-1-q)I}{2q(T_{b}-1)-q^{2}}\right)^{1/2}.$
(126)

###  5.2 Differential Equation for Temperature

In the reaction zone, $\eta<0$, the temperature at $O(\epsilon^{0})$ has a
constant value $T_{b}$ as mentioned in (25). At the next order, the value of
$\hat{\theta}_{1}$ is computed from (71). Using (23), (26) and (70) as
$\eta=0-$ and (105), (125), we transform (71) into the following form.

$\displaystyle
0=\frac{d^{2}\hat{\theta}_{1}}{d\eta^{2}}+q\left(T_{b}-1-\frac{q}{2}\right)\frac{\hat{\theta}_{1-}-\hat{\theta}_{1}+\eta(T_{b}-1-q)}{\mathrm{e}^{\hat{\theta}_{1-}}-(T_{b}-1-q)I}\mathrm{e}^{\hat{\theta}_{1}},$
(127)

subject to the boundary conditions

$\displaystyle\frac{d\hat{\theta}_{1}}{d\eta}\Big{|}_{-}$
$\displaystyle=T_{b}-T_{ad},$ (128) $\displaystyle\hat{\theta}_{1-}$
$\displaystyle={T}_{1+}-\frac{1}{2}{Ma}^{2}\frac{N}{T_{b}^{2}}\nu{V}_{0+}^{2}\left(1-\frac{1+q}{T_{b}}\right)\left(1+\frac{T_{b}}{V_{0+}}\right).$
(129)

The first condition (128) is obtained by taking $\xi=0+$ in (41), by use of
(24), (25), (105) and (114). From (128), we find that the temperature gradient
at the burned-side edge of a flame front is simply described by the deviation
of flame temperature from its adiabatic value. The second condition (129) is
obtained by taking $\xi=0+$ in (42), with the help of (12), (19) and (114).

In the compression zone, $\xi>0$, the asymptotic solution of temperature is
obtained from (12), (19), (31), (42) and (105) as

$\displaystyle
T=T_{b}+\frac{T_{b}^{2}}{N}{T}_{1+}-{Ma}^{2}\frac{\gamma-1}{\gamma}\nu{V}_{0+}^{2}W\left(1+\frac{W}{2}\right)+O(N^{-2}).$
(130)

In the subsequent subsections, by use of (125), the flow properties are
investigated for adiabatic, inviscid and non-adiabatic cases.

###  5.3 Adiabatic case

In the adiabatic case, $T_{b}=T_{ad}=1+q$, we find $V_{0-}/V_{0+}-1=\alpha=0$
due to (114). Then, from (124), we find that the temperature on the burned
side is decreased by the compressibility effect brought by $O({Ma}^{2})$ term,
because $T_{ad}>1$ for $q>0$.

$\displaystyle T_{1+}=b^{2}\frac{\gamma-1}{2}(1-T_{ad}^{2}).$ (131)

Now, outer solutions (105)-(113), (120) and (123) read

$\displaystyle{T}_{0}$
$\displaystyle=\left\\{\begin{array}[]{lll}1+q\mathrm{e}^{z}&(z<0)\\\
1+q&(z>0)\end{array}\right.,\quad{Y}_{0}=\left\\{\begin{array}[]{lll}1-\mathrm{e}^{Lez}&(z<0)\\\
0&(z>0)\end{array}\right.,$ (136) $\displaystyle{V}_{0}$
$\displaystyle={T}_{0},\quad{R}_{0}={V}_{0}^{-1},\quad{P}_{0}=1\quad(z\gtrless
0),$ (137) $\displaystyle{P}_{1}$
$\displaystyle=\left\\{\begin{array}[]{lll}b^{2}(\nu-\gamma)q\mathrm{e}^{z}&(z<0)\\\
-b^{2}\gamma
q&(z>0)\end{array}\right.,\quad{Y}_{1}=\left\\{\begin{array}[]{lll}C_{Y}\mathrm{e}^{Lez}&(z<0)\\\
0&(z>0)\end{array}\right.,$ (142) $\displaystyle{T}_{1}$
$\displaystyle=\left\\{\begin{array}[]{lll}C_{T}\mathrm{e}^{z}-b^{2}\frac{\gamma-1}{\gamma}q\left((\nu-\gamma)z+\frac{q}{2}(2\nu-\gamma)\mathrm{e}^{z}\right)\mathrm{e}^{z}&(z<0)\\\
b^{2}\frac{\gamma-1}{2}(1-T_{ad}^{2})&(z>0)\end{array}\right.,$ (145)

and inner solutions in compression zone (26), (27), (31), (36) and (39) are

$\displaystyle\check{c}_{i}=0\;(i=0,1\ldots),\quad\check{\theta}_{0}=\check{v}_{0}=\check{\rho}_{0}^{-1}=T_{ad},\quad\check{p}_{0}=1,$
(146)
$\displaystyle\check{\theta}_{1}=b^{2}\frac{\gamma-1}{2}(1-T_{ad}^{2}),\quad\check{v}_{1}=b^{2}q\left(1+\frac{q}{2}(1+\gamma)\right),\quad\check{p}_{1}=-b^{2}\gamma
q,$ (147)

and those in reaction zone (25), (61), (62), (70), (80) and (81) are

$\displaystyle\hat{c}_{0}=0,\quad\hat{\theta}_{0}=\hat{v}_{0}=\hat{\rho}_{0}^{-1}=T_{ad},\quad\hat{p}_{0}=1,$
(148)
$\displaystyle\hat{c}_{1}=C_{Y}+\frac{Le}{q}\left(C_{T}-b^{2}\frac{\gamma-1}{\gamma}\frac{q^{2}}{2}(2\nu-\gamma)-\hat{\theta}_{1}\right),$
(149) $\displaystyle\hat{v}_{1}=b^{2}\gamma
qT_{ad}+\hat{\theta}_{1}-\mathrm{e}^{-\lambda\eta}\int_{-\infty}^{\eta}\mathrm{e}^{\lambda\eta}\frac{d\hat{\theta}_{1}}{d\eta}d\eta,$
(150) $\displaystyle\hat{p}_{1}=-b^{2}\gamma
q+\frac{\mathrm{e}^{-\lambda\eta}}{T_{ad}}\int_{-\infty}^{\eta}\mathrm{e}^{\lambda\eta}\frac{d\hat{\theta}_{1}}{d\eta}d\eta,$
(151)

The value of $\hat{\theta}_{1}$ is computed from the differential equation
(127) for the adiabatic case:

$\displaystyle
0=\frac{d^{2}\hat{\theta}_{1}}{d\eta^{2}}-\frac{q^{2}}{2}\left(\hat{\theta}_{1}-T_{1+}\right)\mathrm{e}^{\hat{\theta}_{1}-T_{1+}}.$
(152)

Accounting for (114), (124) and (129), the value of laminar flame speed for
the adiabatic case, $\tilde{S}_{Lad}$, is expressed as follows.

$\displaystyle\tilde{S}_{Lad}=\left(2\tilde{A}\tilde{D}_{th}Le\frac{T_{ad}^{3}}{q^{2}N^{2}}\mathrm{e}^{T_{1+}-N/T_{ad}}\right)^{1/2}.$
(153)

Taking $b\to 0$ in (153), which means $T_{1+}\to 0$ due to (124), we recover
the formula of laminar flame speed in the zero-Mach-number model [29, 63]:

$\displaystyle\tilde{S}_{Lzero}=\left(2\tilde{A}\tilde{D}_{th}Le\frac{T_{ad}^{3}}{q^{2}N^{2}}\mathrm{e}^{-N/T_{ad}}\right)^{1/2}.$
(154)

Because it is convenient to treat non-dimensional quantities, we introduce
${Ma}_{ad}=\tilde{S}_{Lad}/\tilde{c}_{s}$ and
${Ma}_{zero}=\tilde{S}_{Lzero}/\tilde{c}_{s}$, which are the values in the
adiabatic low-Mach-number and zero-Mach-number cases, respectively. Then, from
(153) and (154), we have

$\displaystyle{Ma}_{ad}={Ma}_{zero}\exp\left(\frac{T_{1+}}{2}\right).$ (155)

From (19), (131) and (155), the adiabatic value of $T_{1+}$ is expressed by
use of the upper branch of Lambert $W$ function as

$\displaystyle
T_{1+}=-W_{0}\left(\frac{{Ma}_{zero}^{2}}{\epsilon}\frac{\gamma-1}{2}(T_{ad}^{2}-1)\right).$
(156)

###  5.4 Inviscid Limit $Pr\to 0$

For inviscid limit $Pr\to 0$, or $\nu\to 0$, inner solutions are approximated
by constant values because $W\to 0$ for $\nu\to 0$ in the range of
$0<\xi<+\infty$ as confirmed from (37): for example, $\check{v}_{0}={V}_{0+}$
which implies ${V}_{0-}={V}_{0+}$, or $\alpha=0$, due to (23) and (61).
Therefore, from (114), the temperature on the burned side of a flame front is
equal to its adiabatic value.

$\displaystyle T_{b}=T_{ad}=1+q.$ (157)

As a consequence, if the flow is inviscid, then it is adiabatic and its
solutions are given by those in Section 5.3 with $\nu=0$.

###  5.5 Non-adiabatic case

In this section, we reveal the relation between the laminar flame speed and
the temperature profile for a viscous compressible flow. Because it is
convenient to use ${Ma}$ rather than $\tilde{S}_{L}$, we transform (126) into
the following non-dimensional form, by use of (1) and (153):

$\displaystyle{Ma}={Ma}_{zero}\left(q\left(\frac{T_{b}}{T_{ad}}\right)^{3}\mathrm{e}^{N(1/T_{ad}-1/T_{b})}\frac{\mathrm{e}^{\hat{\theta}_{1-}}-(T_{b}-1-q)I}{2(T_{b}-1)-q}\right)^{1/2}.$
(158)

Requiring
$\frac{\mathrm{e}^{\hat{\theta}_{1-}}-(T_{b}-1-q)I}{2(T_{b}-1)-q}>0$, we find
the value of $T_{b}$ is limited to

$\displaystyle
1+\frac{q}{2}<T_{b}<1+q+\frac{\mathrm{e}^{\hat{\theta}_{1-}}}{I}.$ (159)

The effect of non-zero temperature gradient denoted by the deviation of flame
temperature from its adiabatic value, on the laminar flame speed, or the Mach
number, is investigated by solving (127) and (158) subject to boundary
conditions (128) and (129). The computation of $\hat{\theta}_{1}$ and ${Ma}$
from these equations is made by iterative method: at first, we calculate the
adiabatic temperature solution by solving (152). This solution is, in turn,
used as an initial value in studying the non-adiabatic case. For example, by
fixing some value of ${Ma}$ with given parameters, we change the value of
$d\hat{\theta}_{1}/d\eta|_{-}\neq 0$ so as to $I$ defined by (75) converges
numerically, where $\hat{\theta}_{1}$ is computed from (127). Next, we check
the numerical convergence of ${Ma}$ from (158), and if it does not converge,
then we repeat the previous step on $d\hat{\theta}_{1}/d\eta|_{-}\neq 0$. In
this study, we set parameters as $q=2$, $\gamma=1.4$, $N=90$ and
${Ma}_{zero}=0.1$. We investigate the effect of non-zero temperature gradient
at the burned-side edge of a flame front, $d\hat{\theta}_{1}/d\eta|_{-}\neq 0$
in (128), due to the gas compression. Besides, the influence of viscosity,
expressed by Prandtl number $Pr$, is also revealed.

For the fixed value of Prandtl number $Pr=3/4$, the temperature profiles are
shown in Fig. 2 for $d\hat{\theta}_{1}/d\eta|_{-}>0$ and in Fig. 2 for
$d\hat{\theta}_{1}/d\eta|_{-}<0$. Figure 2 shows that the positive value of
temperature gradient at the burned-side edge of a flame front brings the
monotonic increase of temperature. This behavior is similar to that in
adiabatic case except for non-zero values of temperature gradient. On the
other hand, the negative gradient of temperature induces the non-monotonic
profile with the maximum value of temperature achieved inside the reaction
zone, as plotted in Fig. 2. Such a profile is to be observed, for example, in
the numerical study of thermonuclear flames based on the deflagration model
[24], detonation model [33, 50] or reconstruction model [56]. Our results are
consistent with those obtained in previous works qualitatively. In the present
study, the asymptotic solution of temperature on the burned side, $z>0$, is
described by Lamber $W$ function as obtained in (130).

How the variation of $Pr$ influences the temperature profile is investigated
in Figs. 4 and 4. The opposed effects are observed in positive or negative
values of temperature gradient, respectively. In the case of
$d\hat{\theta}_{1}/d\eta|_{-}>0$, the increase of $Pr$ leads to the decrease
of temperature as shown in Fig. 4. On the other hand, for
$d\hat{\theta}_{1}/d\eta|_{-}<0$, Fig. 4 indicates that the temperature is
enhanced by the viscous effect.

The effect of non-adiabaticity on the Mach number ${Ma}$ is found in Table 1,
accompanied with the value of $I$. Clearly, the rise of temperature gradient
induces the decrease of ${Ma}$. Therefore, we may consider that fast
deflagration waves possess the non-adiabatic temperature profile as
classically indicated by ZND Theory for one dimensional detonation wave [67,
59, 18].

Figure 1: Temperature distribution with $q=2$, $N=90$, $\gamma=1.4$,
${Ma}_{zero}=0.1$ and $Pr=3/4$ for positive values of temperature gradient
given in (128): $d=d\hat{\theta}_{1}/d\eta|_{-}$.

Figure 2: The same as Fig. 2 but for negative values of temperature gradient.

Figure 3: The same as Fig. 2 but for several values of Prandtl number with
${Ma}=0.08$.

Figure 4: The same as Fig. 2 but for several values of Prandtl number with
${Ma}=0.11$.

Table 1: The values of $I$ and $d\hat{\theta}_{1}/d\eta|_{-}$, which is the
temperature gradient at the burned-side edge of a flame front given in (128)
are computed for several values of ${Ma}$. Parameters are set as $q=2$,
$N=90$, $\gamma=1.4$, ${Ma}_{zero}=0.1$ and $Pr=3/4$.

${Ma}$ | $I$ | $d\hat{\theta}_{1}/d\eta|_{-}$
---|---|---
$0.12$ | $41.540$ | $-0.0257$
$0.11$ | $34.046$ | $-0.0184$
$0.1$ | $28.041$ | $-0.00835$
$0.095$ | $23.829$ | $-0.0024$
$0.09$ | $21.364$ | $0.005$
$0.085$ | $18.385$ | $0.0142$
$0.08$ | $15.399$ | $0.0264$
$0.075$ | $12.796$ | $0.042$

##  6 CONCLUSIONS

We have employed the method of matched asymptotic expansions with respect to
$\epsilon$ ($=T_{b}^{2}/N$) under the assumption of large activation energy
$N\gg 1$. By assuming that the order of Mach number $Ma$ is squared root of
$\epsilon$, the compression zone, whose thickness is $O({Ma}^{2})$, has been
introduced behind a premixed flame front. Inside of the compression zone, the
effect of compressibility originates from the pressure variation due to the
existence of viscosity in the Navier-Stokes equation. The Navier-Stokes
equation is solved analytically inside the compression zone to reveal the
inner solutions of velocity and pressure are expressed by Lambert $W$
function. Then, the consideration of pressure term in the heat-conduction
equation leads to the non-adiabatic variation of temperature with its
asymptotic solution obtained. In the reaction zone, the analytical solution of
temperature is not obtained due to the existence of reaction term. Instead,
the temperature profile is captured by solving the heat-conduction equation
numerically subject to non-zero gradient condition of temperature at the
burned-side edge of a flame front.

Unlike the previous works, the derivative of combination of temperature and
mass fraction, or enthalpy, is no longer zero at the unburned-side edge of
reaction zone. This is because the temperature gradient at the burned-side
edge of reaction zone is not restricted to zero and there is a possibility of
positive or negative gradient. The non-adiabatic feature of a flame front is
observed by examining non-zero values of temperature gradient
$d\hat{\theta}_{1}/d\eta|_{-}\neq 0$. If $d\hat{\theta}_{1}/d\eta|_{-}>0$, the
monotonic increase of temperature is observed, which is similar with that in
the adiabatic case, except for non-zero temperature gradient on the burned
side of a flame front. On the other hand, in the case of
$d\hat{\theta}_{1}/d\eta|_{-}<0$, the temperature shows its non-monotonic
behavior and has a single peak structure inside the reaction zone. This
temperature profile is consistent with that obtained in previous works on
thermonuclear flames for small Mach numbers. Besides, the viscosity brings the
decrease of temperature for $d\hat{\theta}_{1}/d\eta|_{-}>0$ and the rise of
it for $d\hat{\theta}_{1}/d\eta|_{-}<0$, though it has a slight impact.

The main assumption of our study is ${Ma}=O(1/\sqrt{N})$. Due to this specific
order of Mach number, our results strongly depend on the value of non-
dimensional activation energy $N$. The removal of such a specific dependence
and the study of general values of Mach number need to be treated in the
future work.

## Acknowledgements

This work was supported by the Sasakawa Scientific Research Grant from The
Japan Science Society (Grant Number 2022-2002). The author is grateful to
Professor Makoto Hirota for fruitful discussions and kind support.

## Disclosure statement

No potential conflict of interest was reported by the author.

## References

  * [1] A. S. Almgren, J. B. Bell, C. A. Rendleman, and M. Zingale.  Low Mach number modeling of type Ia supernovae. i. hydrodynamics.  Astrophys. J., 637(2):922–936, feb 2006.
  * [2] John R. Boisseau, J. Craig Wheeler, Elaine S. Oran, and Alexei M. Khokhlov.  The multidimensional structure of detonations in type Ia supernovae.  Astrophys. J., 471(2):L99–L102, nov 1996.
  * [3] J. R. Bowen.  Approximations to the structure of a laminar detonation wave by the method of matched asymptotic expansions.  Phys. Fluids, 10(2):290–296, 1967.
  * [4] I. Brailovsky and G. I. Sivashinsky.  Hydraulic resistance as a mechanism for deflagration–to–detonation transition.  Combust. Flame, 122(4):492–499, 2000.
  * [5] J. Buckmaster.  The quenching of two–dimensional premixed flames.  Acta Astronaut., 6:741–769, 1979.
  * [6] W. B. Bush and F. E. Fendell.  Asymptotic analysis of laminar flame propagation for general Lewis numbers.  Combust. Sci. Tech., 1:421–428, 1970.
  * [7] W. B. Bush and L. Krishnamurthy.  Asymptotic analysis of the structure of a steady planar detonation: review and extension.  Math. Probl. Eng., 5:223–254, 1999.
  * [8] V. Bychkov, V. Akkerman, D. Valiev, and C. K. Law.  Role of compressibility in moderating flame acceleration in tubes.  Phys. Rev. E, 81:026309, Feb 2010.
  * [9] V. Bychkov, M. Modestov, V. Akkerman, and L.-E. Eriksson.  The Rayleigh–Taylor instability in inertial fusion, astrophysical plasma and flames.  Plasma Phys. Control. Fusion, 49(12B):B513–B520, Nov 2007.
  * [10] V. Bychkov, M. Modestov, and M. Marklund.  The Darrieus–Landau instability in fast deflagration and laser ablation.  Phys. Plasmas, 15:032702, 2008.
  * [11] P. Clavin and G. Searby.  Combustion Waves and Fronts in Flows.  Cambridge University Press, Cambridge, 2016.
  * [12] P. Clavin and F. A. Williams.  Effects of molecular diffusion and of thermal expansion on the structure and dynamics of premixed flames in turbulent flows of large scale and low intensity.  J. Fluid Mech., 116:251–282, 1982.
  * [13] R. M. Corless, G. H. Gonnet, D. E. G. Hare, D. J. Jeffrey, and D. E. Knuth.  On the Lambert W function.  Adv. Comput. Math., 5:329–359, 1996.
  * [14] C. F. Curtiss, J. O. Hirschfelder, and M. P. Barnett.  Theory of detonations. iii. ignition temperature approximation.  J. Chem. Phys., 30(2):470–492, 1959.
  * [15] G. Darrieus.  Propagation d’un front de flamme.  Unpublished work presented at La Technique Moderne (Paris), 1938\.
  * [16] Sinan Demir, Anish Raman Calavay, and V’yacheslav Akkerman.  Influence of gas compressibility on a burning accident in a mining passage.  Combust. Theory Model., 22(2):338–358, 2018.
  * [17] B. Demirgok, O. Ugarte, D. Valiev, and V. Akkerman.  Effect of thermal expansion on flame propagation in channels with nonslip walls.  Proc. Combust. Inst, 35:929–936, 2015.
  * [18] W. Döring.  Über detonationsvorgang in gasen.  Ann. Physik, 43:421–436, 1943.
  * [19] W. Eckhaus.  On the stability of laminar flame–fronts.  M.I.T Fluid Dynamics Research Group Report, (59-4), 1959.
  * [20] W. Eckhaus.  Theory of flame–front stability.  J. Fluid Mech., 10:80–100, 1961.
  * [21] P. C. Fife.  Dynamics of internal layers and diffusive interfaces.  SIAM, Philadelphia, 1988.
  * [22] R. Friedman and E. Burke.  A theoretical model of a gaseous combustion wave governed by a first–order reaction.  J. Chem. Phys., 21:710–714, 1953.
  * [23] I. Glassman and R. A. Yetter.  Combustion.  Academic Press, Burlington, 4th edition, 2008.
  * [24] S. I. Glazyrin, S. I. Blinnikov, and A. Dolgov.  Flame fronts in type Ia supernovae and their pulsational stability.  Mon. Notices Royal Astron. Soc., 433(4):2840–2849, 06 2013.
  * [25] Wolfgang Hillebrandt and Jens C. Niemeyer.  Type Ia supernova explosion models.  Annu. Rev. Astron. Astrophys., 38(1):191–230, Sep 2000.
  * [26] J. O. Hirschfelder.  Effect of chemical reaction order on flame propagation.  Phys. Fluids, 4:253–259, 1961.
  * [27] J. O. Hirschfelder, C. F. Curtiss, and M. P. Barnett.  Ignition temperature approximation for bimolecular detonations.  Phys. Fluids, 4:260–261, 1961.
  * [28] F. Hoyle and William A. Fowler.  Nucleosynthesis in supernovae.  Astrophys. J., 132:565–590, nov 1960.
  * [29] G. Joulin and P. Vidal.  Hydrodynamics and Nonlinear Instabilities, chapter 5.  Cambridge University Press, Cambridge, 1998.
  * [30] L. Kagan and G. I. Sivashinsky.  Deflagration–to–detonation transition in narrow channels: Hydraulic resistance vs. flame folding.  25th ICDERS Leeds, UK, 2015.
  * [31] A. K. Kapila, B. J. Matkowsky, and A. van Harten.  An asymptotic theory of deflagrations and detonations i. the steady solutions.  SIAM J. Appl. Math., 43(3):491–519, 1983.
  * [32] B. Karlovitz, D.W. Denniston, D.H. Knapschaefer, and F.E. Wells.  Studies on turbulent flames: A. flame propagation across velocity gradients B. turbulence measurement in flames.  Fourth Symposium (International) on Combustion, pages 613–620, 1953\.
  * [33] A. M. Khokhlov.  The structure of detonation waves in supernovae.  Mon. Notices Royal Astron. Soc., 239(3):785–808, 08 1989.
  * [34] V. N. Kurdyumov and M. Matalon.  Effects of gas compressibility on the dynamics of premixed flames in long narrow adiabatic channels.  Combust. Theory Model., 20:6:1046–1067, 2016.
  * [35] V. N. Kurdyumov and M. Matalon.  Critical conditions for flame acceleration in long adiabatic channels closed at their ignition end.  Proc. Combust. Inst, 36:1549–1557, 2017.
  * [36] L. D. Landau.  On the theory of slow combustion.  Acta Phys. U.R.S.S., 19(1):77–85, 1944.
  * [37] C. K. Law.  Combustion Physics.  Cambridge University Press, Cambridge, 2006.
  * [38] B. Lewis and G. von Elbe.  Combustion flames and explosions of gases.  Academic Press, San Diego, 3rd edition, 1987.
  * [39] A. Liñán.  The asymptotic structure of counterflow diffusion flames for large activation energies.  Acta Astronaut., 1:1007–1039, 1974.
  * [40] E. Mallard and H. L. Le-Chatelier.  Thermal model for flame propagation.  Ann. Mines, 4:379, 1883.
  * [41] G. H. Markstein.  Experimental and theoretical studies of flame–front stability.  J. Aero. Sci., 18:199–209, 1951.
  * [42] M. Matalon.  On flame stretch.  Combust. Sci. Tech., 31:169–181, 1983.
  * [43] M. Matalon, C. Cui, and J. K. Bechtold.  Hydrodynamic theory of premixed flames: effects of stoichiometry, variable transport coefficients and arbitrary reaction orders.  J. Fluid Mech., 487:179–210, 2003.
  * [44] M. Matalon and B. J. Matkowsky.  Flames as gasdynamic discontinuities.  J. Fluid Mech., 124:239–259, 1982.
  * [45] B. J. Matkowsky and G. I. Sivashinsky.  An asymptotic derivation of two models in flame theory associated with the constant density approximation.  SIAM J. Appl. Math., 37:686–699, 1979.
  * [46] M. Modestov, V. Bychkov, D. Valiev, and M. Marklund.  Growth rate and the cutoff wavelength of the Darrieus–Landau instability in laser ablation.  Phys. Rev. E, 80:046403, 2009.
  * [47] P. Pelce and P. Clavin.  Influence of hydrodynamics and diffusion upon the stability limits of laminar premixed flames.  J. Fluid Mech., 124:219–237, 1982.
  * [48] J. M. M. Pfannes, J. C. Niemeyer, W. Schmidt, and C. Klingenberg.  Thermonuclear explosions of rapidly rotating white dwarfs – i. deflagrations.  Astron. Astrophys., 509:A74, Jan 2010.
  * [49] Friedrich K. Röpke and Stuart A. Sim.  Models for type Ia supernovae and related astrophysical transients.  Space Sci. Rev., 214:72, Apr 2018.
  * [50] Gary J. Sharpe.  The structure of steady detonation waves in type Ia supernovae: pathological detonations in C–O cores.  Mon. Notices Royal Astron. Soc., 310(4):1039–1052, 12 1999.
  * [51] G. I. Sivashinsky.  Structure of bunsen flames.  J. Chem. Phys., 62:638–643, 1975.
  * [52] G. I. Sivashinsky.  On a distorted flame front as a hydrodynamic discontinuity.  Acta Astronaut., 3:889–918, 1976.
  * [53] G. I. Sivashinsky.  Nonlinear analysis of hydrodynamic instability in laminar flames – i. derivation of basic equations.  Acta Astronaut., 4:1177–1206, 1977.
  * [54] D. B. Spalding.  I. Predicting the laminar flame speed in gases with temperature–explicit reaction rates.  Combust. Flame, 1:287–295, 1957.
  * [55] D. B. Spalding.  II. One–dimensional laminar flame theory for temperature–explicit reaction rates.  Combust. Flame, 1:296–307, 1957.
  * [56] Dean M. Townsley, Broxton J. Miles, F. X. Timmes, Alan C. Calder, and Edward F. Brown.  A tracer method for computing type Ia supernova yields: Burning model calibration, reconstruction of thickened flames, and verification for planar detonations.  Astrophys. J. Suppl. Ser., 225(1):3, jul 2016.
  * [57] Damir M. Valiev, V’yacheslav Akkerman, Mikhail Kuznetsov, Lars-Erik Eriksson, Chung K. Law, and Vitaly Bychkov.  Influence of gas compression on flame acceleration in the early stage of burning in tubes.  Combust. Flame, 160(1):97–111, 2013.
  * [58] T. von Kármán.  The present status of the theory of laminar flame propagation.  Sixth Symposium on Combustion, Reinhold, New York, pages 1–11, 1957\.
  * [59] J. von Neumann.  Theory of detonation waves.  O.S.R.D. Rept. No. 549, 1942.
  * [60] Keigo Wada.  Effect of small heat release and viscosity on thermal-diffusive instability.  Scientific Reports, 11(20225):1–22, 2021.
  * [61] Keigo Wada.  Lambert W function for effect of gas compression on premixed flames with ignition temperature approximation.  Combustion Theory and Modelling, 25(3):542–569, 2021.
  * [62] J. C. Wheeler and R. P. Harkness.  Type I supernovae.  Rep. Prog. Phys., 53(12):1467–1557, dec 1990.
  * [63] F. A. Williams.  Combustion Theory: The Fundamental Theory of Chemically Reacting Flow Systems.  Addison-Wesley, Menlo Park, 2nd edition, 1985.
  * [64] S. E. Woosley.  Models for type Ia supernovae.  Nucl. Phys. A, 688(1):9–16, 2001.
  * [65] Y. B. Zel’dovich and D. A. Frank-Kamenetskii.  Theory of thermal flame propagation.  Zhu. Fiz. Khim., 12:100–105, 1938.
  * [66] Y. B. Zel’dovich and D. A. Frank-Kamenetskii.  A theory of thermal propagation of flame.  Acta Physicochim. U.R.S.S., 9(2):341–350, 1938.
  * [67] Y. B. Zel’dovich and A. S. Kompaneets.  Theory of Detonation.  Academic Press, New York, 1960.