# Théorie de Frank-Kamenetskii

En combustion, la théorie de Frank-Kamenetskii est une explication de la première phase de l'explosion ou déflagration d'un mélange. Elle doit son nom à David Frank-Kamenetskii, qui, avec Nikolaï Semionov, a développé la théorie dans les années 1930$^,$$^,$$^,$.

## Description du problème

Considérons un récipient de dimension $a$ contenant un mélange réactif homogène de masse volumique $\rho$, maintenu à une température constante $T_{o}$. On met de côté tous les phénomènes liés à l'allumage du mélange : l'explosion est donc régie uniquement par l'équation d'énergie. En supposant une réaction globale en une étape ${\text{carburant}}+{\text{oxydant}}\rightarrow {\text{produits}}+q$, où $q$ est la quantité de chaleur libérée par unité de masse de carburant consommé, et une vitesse de réaction régie par la loi d'Arrhenius, l'équation d'énergie s'écrit$^,$$^,$$^,$$^,$ :
    $\rho c_{v}{\frac {\partial T}{\partial t}}=\lambda \nabla ^{2}T+q\rho BY_{Fo}e^{-E/(RT)}$
où :
| $T$ | est la température du mélange ; |
| $c_{v}$ | est la chaleur spécifique à volume constant ; |
| $\lambda$ | est la conductivité thermique ; |
| $B$ | est le facteur pré-exponentiel de la réaction ; |
| $Y_{Fo}$ | est la fraction massique initiale du combustible ; |
| $E$ | est l'énergie d'activation ; |
| $R$ | est la constante universelle des gaz parfaits. |

### Simplifications et adimensionalisation

Une augmentation de température de l'ordre de $T-T_{o}\sim RT_{o}^{2}/E$, où $RT_{o}^{2}/E$ est la « température de Frank-Kamenetskii » est suffisante pour augmenter la réaction chimique d'une quantité $e$ comme le montre le rapport :
    ${\frac {e^{-E/RT}}{e^{-E/RT_{o}}}}=\exp \left[{\frac {E}{RT_{o}}}\left(1-{\frac {T_{0}}{T}}\right)\right].$
Les échelles adimensionnelles de temps, de température, de longueur et de transfert de chaleur peuvent être définies comme suit :
    $\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}}}$
où
| $t_{c}=\rho c_{v}a^{2}/\lambda$ | est le temps caractéristique de conduction thermique à travers le récipient ; |
| $t_{f}=\left(Be^{-\beta }\right)^{\!-1}$ | est le temps caractéristique de consommation de carburant ; |
| $t_{e}=\left(\beta \gamma Be^{-\beta }\right)^{\!-1}$ | est le temps caractéristique d'allumage ; |
| $a$ | est la distance caractéristique, par exemple le rayon du récipient ; |
| $\beta =E/(R\,T_{o})$ | est l'énergie d'activation adimensionnelle ; |
| $\gamma =(qY_{Fo})/(c_{v}T_{o})$ | est le paramètre de dégagement de chaleur ; |
| $\delta$ | est le nombre de Damköhler ; |
| $r$ | est la coordonnée spatiale avec l'origine au centre ; |
| $j=0$ | pour une dalle plane ; |
| $j=1$ | pour un récipient cylindrique ; |
| $j=2$ | pour un récipient sphérique. |

#### Remarque

    Dans un processus de combustion typique, $\gamma \approx 6{-}8,\ \beta \approx 30{-}100$ de sorte que $\beta \gamma \gg 1$.
    Par conséquent, $t_{f}=\beta \gamma t_{e}\gg 1$. Autrement dit, le temps de consommation de carburant est beaucoup plus long que le temps d'allumage, de sorte que la consommation de carburant est essentiellement négligeable dans l'étude de l'allumage.
    C'est pourquoi la concentration de carburant est supposée rester voisine de la concentration initiale $Y_{Fo}$.
En remplaçant les variables non dimensionnelles dans l'équation d'énergie, on obtient :
    ${\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 )}$
Puisque $\beta \gg 1$, le terme exponentiel peut être linéarisé $e^{\theta /(1+\theta /\beta )}\approx e^{\theta }$, d'où la simplification :
    ${\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 }$
À $\tau =0$, nous avons $\theta (\eta ,0)=0$ et pour $\tau >0$, $\theta$ doit satisfaire $\theta (1,\tau )=0$ et $\partial \theta /\partial \eta |_{\eta =0}=0.$

## Théorie de Semionov

Avant David Frank-Kamenetskii, son directeur de thèse Nikolaï Semionov a proposé une théorie de l'explosion thermique avec un modèle plus simple dans lequel il suppose une fonction linéaire pour le processus de conduction thermique en lieu et place de l'opérateur Laplacien. L'équation de Semionov en géométrie plane s'écrit :
    ${\frac {d\theta }{d\tau }}=e^{\theta }-{\frac {\theta }{\delta }},\quad \theta (0)=0$
dans laquelle le terme exponentiel $e^{\theta }$ tend à augmenter la température réduite $\theta$ tandis que le terme linéaire $-\theta /\delta$ tend à la diminuer. L'importance relative entre les deux termes est déterminée par le nombre de Damköhler $\delta$. La solution numérique de l'équation ci-dessus pour différentes valeurs de $\delta$ est illustrée dans la figure.

### Régime stationnaire

Lorsque $0<\delta <e^{-1}$, le terme linéaire finit par dominer et le système atteint un état stationnaire lorsque $\tau \rightarrow \infty$. À l'état stationnaire ($d\theta /d\tau =0$), l'équilibre est donné par l'équation :
    $\delta e^{\theta }=\theta ,\quad \Rightarrow \quad \theta =-W(-\delta )$
où $W$ représente la fonction W de Lambert. À partir des propriétés de cette fonction il est facile de voir que la température à l'état stationnaire fournie par l'équation ci-dessus n'existe que lorsque $\delta \leq \delta _{c}=1/e$, où $\delta _{c}$ est appelé paramètre de Frank-Kamenetskii en tant que point critique où le système bifurque de l'existence d'un état stationnaire à un état explosif pour les grandes valeurs du temps.

### Régime explosif

Pour $\delta _{c}<\delta <\infty$, le système explose puisque le terme exponentiel domine lorsque la durée augmente. En raison du forçage exponentiel $\theta \rightarrow \infty$ pour une valeur finie de $\tau$. Ce temps est interprété comme le temps d'allumage ou le temps d'induction du système. Lorsque $\delta \gg \delta _{c}$ le terme de conduction thermique $-\theta /\delta$ peut être négligé, auquel cas le problème admet une solution explicite :
    ${\frac {d\theta }{d\tau }}=e^{\theta }\quad \Rightarrow \quad \theta =\ln \left({\frac {1}{1-\tau }}\right)$
À l'instant $\tau =1$, le système explose. Cet instant est également appelé « période d'induction adiabatique » puisque le terme de conduction thermique $-\theta /\delta$ est négligé.
Dans la condition quasi-critique, c'est-à-dire lorsque $\delta -\delta _{c}\rightarrow 0^{+}$, le système met très longtemps à exploser. L'analyse de cette limite a été réalisée pour la première fois par Frank-Kamenetskii bien que des asymptotiques appropriées n'aient été montrées que plus tard par D. R. Kassoy et Amable Liñán y compris la consommation de réactifs car celle-ci n'est pas négligeable lorsque $\tau \sim \beta \gamma$.
Une analyse simplifiée ignorant la consommation de réactifs est présentée ici. Définissons un petit paramètre $\epsilon \rightarrow 0^{+}$ tel que $\delta =\delta _{c}(1+\epsilon )$. Dans ce cas, l'évolution temporelle de $\theta$ est la suivante : elle croît d'abord jusqu'à la valeur de température à l'état stationnaire correspondant à $\delta =\delta _{c}$, qui est donnée par $\theta =-W(-\delta _{c})=1$ à des instants d'ordre $\tau \sim O(1)$, puis il reste très proche de cette valeur à l'état stationnaire pendant un temps long avant de finalement exploser à un temps long. La quantité qui nous intéresse est l'estimation à long terme de l'explosion. Pour trouver l'estimation on introduit les transformations $\zeta ={\sqrt {\epsilon }}\tau$ et $\theta =1+{\sqrt {\epsilon }}\psi (\zeta )-\epsilon$ utilisables dans la région où $\theta$ reste proche de $1$ dans l'équation de référence. En ne retenant que les termes d'ordre un on obtient :
    $\delta _{c}{\frac {d\psi }{d\zeta }}=1+{\frac {\psi ^{2}}{2}},\quad \psi (0)\rightarrow -(1-\theta )/{\sqrt {\epsilon }}=-\infty$
où la condition limite est dérivée en faisant correspondre avec la région initiale dans laquelle $\tau ,\theta \sim 1$. La solution de cette équation est donnée par :
    $\psi =-{\sqrt {2}}\cot \left({\frac {\zeta }{{\sqrt {2}}\delta _{c}}}\right)$
qui montre que $\psi \rightarrow \infty$ lorsque $\zeta ={\sqrt {2}}\pi \delta _{c}$. En écrivant cette condition en termes de $\tau$, le temps d'explosion dans la condition quasi-critique est :
    $\tau ={\sqrt {\frac {2\pi ^{2}\delta _{c}^{3}}{\delta -\delta _{c}}}}$
L'instant d'explosion $\tau \rightarrow \infty$ quand $\delta -\delta _{c}\rightarrow 0^{+}$.

## Théorie de l'état stationnaire de Frank-Kamenetskii

Le seul paramètre caractérisant l'explosion est le nombre de Damköhler $\delta$. Lorsque $\delta$ est très élevé, le temps de conduction est plus long que le temps de réaction chimique et le système explose quand il atteint une haute température. D'autre part, lorsque $\delta$ est très faible, le temps de conduction thermique est beaucoup plus rapide que le temps de réaction chimique, de sorte que toute la chaleur produite par la réaction chimique est immédiatement conduite vers la paroi, il n'y a donc pas d'explosion, on passe à un état presque stationnaire. Amable Liñán a qualifié ce phénomène de mode de réaction lente. À un nombre de Damköhler critique $\delta _{c}$, le système passe du mode de réaction lente au mode explosif. Au lieu de résoudre le problème complet pour trouver ce $\delta _{c}$, Frank-Kamenetskii a résolu le problème de l'état stationnaire pour différents nombres de Damköhler jusqu'à la valeur critique au-delà de laquelle aucune solution stationnaire n'existe. Le problème à résoudre est décrit par l'équation de Liouville – Bratu – Gelfand$^,$ :
    $\nabla ^{2}\theta +\delta e^{\theta }=0$
avec des conditions initiale $\theta (1)=0$ et la condition de symétrie ${\frac {d\theta }{d\eta }}_{\eta =0}=0$.

### Géométrie plane

Dans le cas plan, l'équation s'écrit :
    ${\frac {d^{2}\theta }{d\eta ^{2}}}=-\delta e^{\theta }$
On introduit les transformations $\Theta =\theta _{m}-\theta$ et $\xi ^{2}=\delta e^{\theta _{m}}\eta ^{2}$, où $\theta _{m}$ est la température maximale qui se produit à $\eta =0$ en raison de la symétrie, alors :
    ${\frac {d^{2}\Theta }{d\xi ^{2}}}=e^{-\Theta },\quad \Theta (0)=0,\quad {\frac {d\Theta }{d\xi }}_{\xi =0}=0$
En intégrant une fois et en utilisant la condition de symétrie l'équation devient :
    ${\frac {d\Theta }{d\xi }}={\sqrt {2(1-e^{-\Theta })}}$
En intégrant à nouveau :
    $e^{(\theta _{m}-\theta )/2}=\cosh \left(\eta {\sqrt {\frac {\delta e^{\theta _{m}}}{2}}}\right)$
L'équation ci-dessus est la solution exacte, mais la température maximale de $\theta _{m}$ est inconnue. En utilisant la condition limite au bord $\theta (1)=0$ la température maximale est obtenue à partir de l'expression implicite :
    $e^{\theta _{m}/2}=\cosh {\sqrt {\frac {\delta e^{\theta _{m}}}{2}}}\quad {\text{ou}}\quad \delta =2e^{-\theta _{m}}\left(\operatorname {arcosh} e^{\theta _{m}/2}\right)^{\!2}$
La valeur critique $\delta _{c}$ est obtenue en trouvant le maximum de l'équation $d\delta /d\theta _{m}=0$ à $\delta _{c}$ (voir figure) :
    ${\frac {d\delta }{d\theta _{m}}}=0\quad \Rightarrow \quad e^{\theta _{m,c}/2}-{\sqrt {e^{\theta _{m,c}}-1}}\operatorname {arcosh} e^{\theta _{m,c}/2}=0$
    $\theta _{m,c}=1,1868\quad \Rightarrow \quad \delta _{c}=0,8785$
Le système n'a pas d'état stationnaire (ou explose) pour $\delta >\delta _{c}=0.8785$.
Pour $\delta <\delta _{c}=0.8785$ le système est dans un état quasi-stationnaire avec une réaction très lente.

### Géométrie cylindrique

Paul L. Chambré a fourni une solution exacte du problème en 1952. H. Lemke a également fourni une solution sous une forme quelque peu différente en 1913.
L'équation s'écrit ici :
    ${\frac {1}{\eta }}{\frac {d}{d\eta }}\left(\eta {\frac {d\theta }{d\eta }}\right)=-\delta e^{\theta }$
Si les transformations $\omega =\eta d\theta /d\eta$ et $\chi =e^{\theta }\eta ^{2}$ sont introduites il vient :
    ${\frac {d\omega }{d\chi }}=-{\frac {\delta }{2+\omega }},\quad \omega (0)=0$
La solution générale est $\omega ^{2}+4\omega +C=-2\delta \chi$. $C=0$ d'après la condition de symétrie. En réécrivant dans la variable d'origine, l'équation est :
    $\eta ^{2}\left({\frac {d\theta }{d\eta }}\right)^{2}+4\eta {\frac {d\theta }{d\eta }}=-2\delta \eta ^{2}e^{\theta }$
Mais l'équation d'origine multipliée par $2\eta ^{2}$ est :
    $2\eta ^{2}{\frac {d^{2}\theta }{d\eta ^{2}}}+2\eta {\frac {d\theta }{d\eta }}=-2\delta \eta ^{2}e^{\theta }$
En soustrayant les deux dernières équations l'une de l'autre, on obtient :
    ${\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$
On introduit la fonction $g(\eta )=d\theta /d\eta$ qui vérifie donc :
    ${\frac {dg}{d\eta }}-{\frac {g}{\eta }}-{\frac {g^{2}}{2}}=0$
Il s'agit d'une équation différentielle de Bernoulli d'ordre $2$ dont la solution est :
    $g(\eta )={\frac {d\theta }{d\eta }}=-{\frac {4\eta }{B'+\eta ^{2}}}$
En intégrant celle-ci on obtient $\theta =A-2\ln(B\eta ^{2}+1)$ où $B=1/B'$. Il reste une condition pour les deux constantes $A,\ B$. Mais celles-ci sont liées comme l'on peut le constater en remplaçant la solution ci-dessus dans l'équation de départ. On trouve ainsi $A=\ln(8B/\delta )$. Par conséquent, la solution est :
    $\theta =\ln {\frac {8B/\delta }{(B\eta ^{2}+1)^{2}}}$
Si nous utilisons maintenant l'autre condition limite au bord du domaine $\theta (1)=0$, nous obtenons une équation pour $B$ telle que $\delta (B+1)^{2}-8B=0$. La valeur maximale de $\delta$ pour laquelle une solution est possible est lorsque $B=1$, donc le paramètre critique de Frank-Kamentskii est $\delta _{c}=2$. Le système n'a pas d'état stationnaire (ou explose) pour $\delta >\delta _{c}=2$ et pour $\delta <\delta _{c}=2$, le système passe à un état stationnaire avec une réaction très lente. La température maximale $\theta _{m}=\theta (\eta =0)$ vaut :
    $\theta _{m}=\ln {\frac {8B}{\delta }}\quad {\text{ou}}\quad \delta =8Be^{-\theta _{m}}$
Pour chaque valeur de $\delta$, nous avons deux valeurs de $\theta _{m}$ puisque $B$ est à valeurs multiples. La température critique maximale est $\theta _{m,c}=\ln 4$.

### Géométrie sphérique

Pour la géométrie sphérique, il n'existe pas de solution explicite connue pour l'équation :
    ${\frac {1}{\eta ^{2}}}{\frac {d}{d\eta }}\left(\eta ^{2}{\frac {d\theta }{d\eta }}\right)=-\delta e^{\theta }$
Si les transformations $\Theta =\theta _{m}-\theta$ et $\xi ^{2}=\delta e^{\theta _{m}}\eta ^{2}$, où $\theta _{m}$ est la température maximale qui se produit à $\eta =0$ en raison de la symétrie, sont introduites il vient :
    ${\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$
Cette équation est l'équation de Emden-Chandrasekhar. Contrairement au cas plan et cylindrique, le cas sphérique a une infinité de solutions pour $\delta <\delta _{c}$, oscillant autour du point $\delta =2$ au lieu de seulement deux solutions, ce qui a été montré par Israel Gelfand. La branche basse correspond au comportement explosif.
À partir de la résolution numérique, on constate que le paramètre critique de Frank-Kamenetskii est $\delta _{c}=3.3220$. Le système n'a pas d'état stationnaire (ou explose) pour $\delta >\delta _{c}=3.3220$. Pour $\delta <\delta _{c}=3.3220$ le système passe à un état stationnaire avec une réaction très lente. La température maximale $\theta _{m}$ se produit à $\eta =0$ et la température critique maximale est $\theta _{m,c}=1,6079$.