Frank-Kamenetskii theory

In combustion, Frank-Kamenetskii theory explains the thermal explosion of a homogeneous mixture of reactants, kept inside a closed vessel with constant temperature walls. It is named after a Russian scientist David A. Frank-Kamenetskii, who along with Nikolay Semenov developed the theory in the 1930s.

Problem description

Consider a vessel maintained at a constant temperature $T_{o}$, containing a homogeneous reacting mixture. Let the characteristic size of the vessel be $a$. Sine the mixture is homogeneous, the density $\rho$ is constant. During the initial period of ignition, the reactant concentration is negligible (see $t_{f}$ and $t_{e}$ below), thus the explosion is governed only by the energy equation. Assuming a one-step global reaction ${\text{fuel}}+{\text{oxidizer}}\rightarrow {\text{products}}+q$, where $q$ is the amount of heat released per unit mass of fuel consumed, and reaction rate governed by Arrhenius law, the energy equation becomes

$\rho c_{v}{\frac {\partial T}{\partial t}}=\lambda \nabla ^{2}T+q\rho BY_{Fo}e^{-E/(RT)}$

where

• $T$ is the temperature of the mixture
• $c_{v}$ is the specific heat at constant volume
• $\lambda$ is the thermal conductivity
• $B$ is the pre-exponential factor with dimension of one over time
• $Y_{Fo}$ is the initial fuel mass fraction
• $E$ is the activation energy
• $R$ is the universal gas constant

Non-dimensionalization

The non-dimensional activation energy $\beta$ and the heat-release parameter $\gamma$ are

$\beta ={\frac {E}{RT_{o}}},\quad \gamma ={\frac {qY_{Fo}}{c_{v}T_{o}}}$

The characteristic heat conduction time across the vessel is

$t_{c}={\frac {\rho c_{v}a^{2}}{\lambda }}$

the characteristic fuel consumption time is

$t_{f}=\left(Be^{-\beta }\right)^{-1}$

and the characteristic explosion/ignition time is

$t_{e}=\left(\beta \gamma Be^{-\beta }\right)^{-1}$

Note should be made that in combustion process, typically $\gamma \sim 6{-}8,\ \beta \sim 30{-}100$ so that $\beta \gamma \gg 1$. Therefore, $t_{f}=\beta \gamma t_{e}\gg 1$, i.e., the fuel is consumed at much longer times when compared with ignition time, the fuel consumption is essentially negligible to study ignition/explosion. That is the reason the fuel concentration is assumed to same as the initial fuel concentration $Y_{Fo}$.

The non-dimensional scales are

$\tau ={\frac {t}{t_{e}}},\quad \theta ={\frac {\beta (T-T_{o})}{T_{o}}},\quad \eta ^{j}={\frac {r^{j}}{a}},\quad \delta ={\frac {t_{c}}{t_{e}}}$

where $\delta$ is the Damköhler number and $r$ is the spatial coordinate with origin at the center, $j=0$ for planar slab, $j=1$ for cylindrical vessel and $j=2$ for spherical vessel. With this scale, the equation becomes

${\frac {\partial \theta }{\partial \tau }}={\frac {1}{\delta }}{\frac {1}{\eta ^{j}}}{\frac {\partial }{\partial \eta }}\left(\eta ^{j}{\frac {\partial \theta }{\partial \eta }}\right)+e^{\theta /(1+\theta /\beta )}$

Since $\beta \gg 1$, the exponential term can be linearized $e^{\theta /(1+\theta /\beta )}\approx e^{\theta }$, hence

${\frac {\partial \theta }{\partial \tau }}={\frac {1}{\delta }}{\frac {1}{\eta ^{j}}}{\frac {\partial }{\partial \eta }}\left(\eta ^{j}{\frac {\partial \theta }{\partial \eta }}\right)+e^{\theta }$

Semenov theory

Before Frank-Kamenetskii, his doctoral advisor Nikolay Semyonov(or Semenov) proposed a thermal explosion theory for a simple equation i.e., he assumed a linear function for heat conduction process instead of Laplacian operator. Semenove's equation reads as

${\frac {\partial \theta }{\partial \tau }}=e^{\theta }-{\frac {\theta }{\delta }},\quad \theta (0)=0$

For $e^{-1}<\delta <\infty$, the system explodes since the exponential term dominates. For $0<\delta <e^{-1}$, the system goes to a steady state, the system doesn't explode. In particular, Semenov found the critical Damköhler number, which is called as Frank-Kamenetskii parameter$\delta _{c}=e^{-1}$(where $\theta =1$) as a critical point where the system changes from steady state to explosive state. For $\delta \gg 1$, the solution is

${\frac {\partial \theta }{\partial \tau }}=e^{\theta },\quad \Rightarrow \quad \theta =\ln \left({\frac {1}{1-\tau }}\right)$

At time $\tau =1$, the system explodes.

Frank-Kamenetskii steady-state theory

The only parameter which characterizes the explosion is the Damköhler number $\delta$. When $\delta$ is very high, conduction time is longer than the chemical reaction time and the system explodes with high temperature since there is not enough time for conduction to remove the heat. On the other hand, when $\delta$ is very low, heat conduction time is much faster than the chemical reaction time, such that all the heat produced by the chemical reaction is immediately conducted to the wall, thus there is no explosion, it goes to an almost steady state, Amable Liñán coined this mode as slowly reacting mode. At a critical Damköhler number $\delta _{c}$ the system goes from slowly reacting mode to explosive mode. Therefore, $\delta <\delta _{c}$, the system is in steady state. Instead of solving the full problem to find this $\delta _{c}$, Frank-Kamenetskii solved the steady state problem for various Damköhler number until the critical value, beyond which no steady solution exists. So the problem to be solved is

${\frac {1}{\eta ^{j}}}{\frac {d}{d\eta }}\left(\eta ^{j}{\frac {d\theta }{d\eta }}\right)=-\delta e^{\theta }$

with boundary conditions

$\theta (1)=0,\quad {\frac {d\theta }{d\eta }}_{\eta =0}=0$

the second condition is due to the symmetry of the vessel. The above equation is called as Liouville–Bratu–Gelfand equation in mathematics.

Planar vessel

Frank-Kamenetskii explosion for planar vessel

For planar vessel, there is an exact solution and the equation was initially discussed by Joseph Liouville in 1853. Here $j=0$, then

${\frac {d^{2}\theta }{d\eta ^{2}}}=-\delta e^{\theta }$

If the transformations $\Theta =\theta _{m}-\theta$ and $\xi ^{2}=\delta e^{\theta _{m}}\eta ^{2}$, where $\theta _{m}$ is the maximum temperature which occurs at $\eta =0$ due to symmetry, are introduced

${\frac {d^{2}\Theta }{d\xi ^{2}}}=e^{-\Theta },\quad \Theta (0)=0,\quad {\frac {d\Theta }{d\xi }}_{\xi =0}=0$

Integrating once and using the second boundary condition, the equation becomes

${\frac {d\Theta }{d\xi }}={\sqrt {2(1-e^{-\Theta })}}$

and integrating again

$e^{(\theta _{m}-\theta )/2}=\cosh \left({\sqrt {\frac {\delta e^{\theta _{m}}}{2}}}\eta \right)$

The above equation is the exact solution, but $\theta _{m}$ maximum temperature is unknown, but we have not used the boundary condition of the wall yet. Thus using the wall boundary condition $\theta =0$ at $\eta =1$, the maximum temperature is obtained from an implicit expression,

$e^{\theta _{m}/2}=\cosh \left({\sqrt {\frac {\delta e^{\theta _{m}}}{2}}}\right)\quad {\text{or}}\quad \delta =2e^{-\theta _{m}}\left(\cosh ^{-1}e^{\theta _{m}/2}\right)^{2}$

Critical $\delta _{c}$ is obtained by finding the maximum point of the equation(look figure), i.e., $d\delta /d\theta _{m}=0$ at $\delta _{c}$.

${\frac {d\delta }{d\theta _{m}}}=0,\quad \Rightarrow \quad e^{\theta _{m,c}/2}-{\sqrt {e^{\theta _{m,c}}-1}}\cosh ^{-1}e^{\theta _{m,c}/2}=0$

$\theta _{m,c}=1.1868,\quad \Rightarrow \quad \delta _{c}=0.8785$

So the critical Frank-Kamentskii parameter is $\delta _{c}=0.8785$. The system has no steady state( or explodes) for $\delta >\delta _{c}=0.8785$ and for $\delta <\delta _{c}=0.8785$, the system goes to a steady state with very slow reaction.

Cylindrical vessel

Frank-Kamenetskii explosion for cylindrical vessel

For cylindrical vessel, there is an exact solution and the equation was initially discussed by G.Bratu in 1914. Though Frank-Kamentskii used numerical integration assuming there is no explicit solution, P.L. Chambre provided an exact soution in 1952. Here $j=1$, then

${\frac {1}{\eta }}{\frac {d}{d\eta }}\left(\eta {\frac {d\theta }{d\eta }}\right)=-\delta e^{\theta }$

If the transformations $\omega =\eta d\theta /d\eta$ and $\chi =e^{\theta }\eta ^{2}$ are introduced

${\frac {d\omega }{d\chi }}=-{\frac {\delta }{2+\omega }},\quad \omega (0)=0$

The general solution is

$\omega ^{2}+4\omega +C=-2\delta \chi$

But $C=0$ from the symmetry condition at the centre. Writing back in original variable, the equation reads,

$\eta ^{2}\left({\frac {d\theta }{d\eta }}\right)^{2}+4\eta {\frac {d\theta }{d\eta }}=-2\delta \eta ^{2}e^{\theta }$

But the the the original equation multiplied by $2\eta ^{2}$ is

$2\eta ^{2}{\frac {d^{2}\theta }{d\eta ^{2}}}+2\eta {\frac {d\theta }{d\eta }}=-2\delta \eta ^{2}e^{\theta }$

Now subtracting the last two equation from one another leads to

${\frac {d^{2}\theta }{d\eta ^{2}}}-{\frac {1}{\eta }}{\frac {d\theta }{d\eta }}-{\frac {1}{2}}\left({\frac {d\theta }{d\eta }}\right)^{2}=0$

This equation is easy to solve because it involves only the derivatives, so letting $g(\eta )=d\theta /d\eta$ transforms the equation

${\frac {dg}{d\eta }}-{\frac {g}{\eta }}-{\frac {g^{2}}{2}}=0$

This is a Bernoulli differential equation of order $2$, a type of Riccati equation. The solution is

$g(\eta )={\frac {d\theta }{d\eta }}=-{\frac {4\eta }{B'+\eta ^{2}}}$

Integrating once again, we have

$\theta =A-2\ln(B\eta ^{2}+1)\quad {\text{where}}\quad B={\frac {1}{B'}}$

We have used already one boundary condition, there is one more boundary condition left, but with two constants $A,\ B$. It turns out $A$ and $B$ are related to each other, which is obtained by substituting the above solution into the starting equation we arrive at

$A=\ln \left({\frac {8B}{\delta }}\right)$

Therefore the solution is

$\theta =\ln \left[{\frac {8B/\delta }{(B\eta ^{2}+1)^{2}}}\right]$

Now if we use the other boundary condition $\theta (1)=0$, we get an equation for $B$

$\delta (B+1)^{2}-8B=0$

The maximum value of $\delta$ for which solution is possible is when $B=1$, so the critical Frank-Kamentskii parameter is $\delta _{c}=2$. The system has no steady state( or explodes) for $\delta >\delta _{c}=2$ and for $\delta <\delta _{c}=2$, the system goes to a steady state with very slow reaction..The maximum temperature $\theta _{m}$ occurs at $\eta =0$

$\theta _{m}=\ln \left({\frac {8B}{\delta }}\right)\quad {\text{or}}\quad \delta =8Be^{-\theta _{m}}$

For each value of $\delta$, we have two values of $\theta _{m}$ since $B$ is multi-valued. The maximum critical temperature is $\theta _{m,c}=\ln 4$.

Spherical vessel

For spherical vessel, there is no known explicit solution, so Frank-Kamenetskii used numerical methods to find the critical value.Here $j=2$, then

${\frac {1}{\eta ^{2}}}{\frac {d}{d\eta }}\left(\eta ^{2}{\frac {d\theta }{d\eta }}\right)=-\delta e^{\theta }$

If the transformations $\Theta =\theta _{m}-\theta$ and $\xi ^{2}=\delta e^{\theta _{m}}\eta ^{2}$, where $\theta _{m}$ is the maximum temperature which occurs at $\eta =0$ due to symmetry, are introduced

${\frac {1}{\xi ^{2}}}{\frac {d}{d\xi }}\left(\xi ^{2}{\frac {d\Theta }{d\xi }}\right)=e^{-\Theta },\quad \Theta (0)=0,\quad {\frac {d\Theta }{d\xi }}_{\xi =0}=0$

The above equation is nothing but Chandrasekhar equation, which appears in astrophysics describing isothermal gas sphere.

From numerical solution, it is found that the critical Frank-Kamentskii parameter is $\delta _{c}=3.32$. The system has no steady state( or explodes) for $\delta >\delta _{c}=3.32$ and for $\delta <\delta _{c}=3.32$, the system goes to a steady state with very slow reaction..The maximum temperature $\theta _{m}$ occurs at $\eta =0$ and maximum critical temperature is $\theta _{m,c}=1.61$.

External links

• Planar solution in Chebfun solver http://www.chebfun.org/examples/ode-nonlin/BlowupFK.html