# Validation d'un formalisme cinétique par optimisation de paramètres : approche et méthode mathématiques

F. Alario, A. Bamberger, Z. Benjelloun-Dabaghi

► To cite this version:

F. Alario, A. Bamberger, Z. Benjelloun-Dabaghi. Validation d'un formalisme cinétique par optimisation de paramètres : approche et méthode mathématiques. Revue de l'Institut Français du Pétrole, 1998, 53 (2), pp.195-200. 10.2516/ogst:1998018. hal-02078961

HAL Id: hal-02078961

https://ifp.hal.science/hal-02078961v1

Submitted on 25 Mar 2019

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
VALIDATION D'UN FORMALISME
CINÉTIQUE PAR OPTIMISATION DE
PARAMÈTRES : APPROCHE ET
MÉTHODE MATHÉMATIQUES

F. ALARIO, A. BAMBERGER
et Z. BENJELLOUN-DABAGHI

Institut français du pétrole¹

Cet article est une synthèse du travail fait en commun entre la division Informatique et Mathématiques appliquées "DIMA" et la division Cinétique et Catalyse de l'Institut français du pétrole. On y expose une nouvelle approche mathématique pour la validation d'un formalisme cinétique : la déshydrocyclisation du n-heptane.

VALIDATION OF A KINETIC FORMALIZATION BY
MEANS OF PARAMETER OPTIMIZATION: APPROACH
AND MATHEMATICAL METHOD

This article is an overview of the work carried out together by
the "DIMA" and the "Cinétique et Catalyse" division of the Institut
français du pétrole. It is a new mathematical approach to validate a
kinetic formalization: the n-heptane dehydrocyclization.

VALIDACIÓN DE UN FORMALISMO CINÉTICO POR
OPTIMIZACIÓN DE PARÁMETROS: ENFOQUE Y
MÉTODO MATEMÁTICOS

Este estudio constituye una síntesis del trabajo llevado a cabo en común entre la "DIMA", y la división "Cinétique et Catalyse" del Institut français du pétrole. Se expone un nuevo enfoque matemático para la validación de un formalismo cinético: la deshidrociclización del n-heptano.

<footnote>(1) 1 et 4, avenue de Bois-Préau,
92852 Rueil-Malmaison Cedex - France</footnote>
# INTRODUCTION

Les schémas réactionnels de transformations chimiques catalysées intervenant dans les procédés de raffinage et pétrochimie font intervenir un certain nombre de composés tels que les réactifs présents dans la charge à traiter, les produits de la réaction chimique et les intermédiaires réactionnels. Il s'établit ainsi un véritable réseau constitué de réactions chimiques élémentaires pouvant être consécutives, parallèles, équilibrées ou non.

La modélisation de ces processus réactionnels est du plus grand intérêt puisqu'elle permettra de prédire la répartition des produits de réaction en fonction de la composition de la charge et des conditions opératoires pour un catalyseur, donc un procédé, donné.

Une condition nécessaire à l'obtention d'un modèle prédictif que l'on souhaite pouvoir faire évoluer au cours du temps est qu'il doit être basé sur des équations cinétiques possédant un sens physique, et non sur des formules mathématiques choisies pour représenter uniquement les allures des courbes expérimentales. Une telle approche consiste ainsi en la résolution cinétique du schéma réactionnel, c'est-à-dire en la détermination des vitesses des réactions de chaque étape élémentaire, maille du réseau.

Afin de mettre en équations ces vitesses, le cinéticien fait un certain nombre d'hypothèses sur la cinétique réactionnelle, qui vont conduire à l'écriture d'un mécanisme faisant intervenir notamment des intermédiaires de réaction dont l'existence est supposée. Ces hypothèses devront être confirmées ou infirmées en confrontant les résultats théoriques ainsi obtenus aux résultats expérimentaux. Les équations de vitesse, mises sous forme algébrique, vont contenir les variables expérimentales, qui sont dans la plupart des cas les pressions partielles en réactifs, et des paramètres cinétiques tels que la constante apparente de vitesse, les coefficients d'adsorption des composés en présence, les constantes d'équilibre le cas échéant. La difficulté pour le cinéticien résidera en la détermination numérique de ces paramètres, c'est-à-dire de répondre à la question « existe-t-il un jeu de paramètres permettant à l'équation de vitesse de représenter la réalité expérimentale et, si oui, ce jeu de paramètres est-il unique et a-t-il un sens physique ? ».

Cette étude expose la démarche qui a été retenue afin de répondre au problème posé. L'exemple choisi est celui de la déshydrocyclisation du *n*-heptane qui est une réaction importante du reformage catalytique.

# 1 POSITION DU PROBLÈME

Pour la réaction de déshydrocyclisation du *n*-heptane :

<figure><img src="image_2.png" /><figcaption>Figure 1
Déshydrocyclisation du *n*-heptane.</figcaption></figure>

un mécanisme réactionnel possible conduit à l'équation de vitesse [5] suivante :

$$V = \frac{K_{MFM} K'_{1} (b_{R} p_{R}) (b_{H} p_{H_{2}})^{\frac{5}{2}}}{\left[ \left( b_{H} p_{H_{2}} \right)^{2} + \left( b_{H} p_{H_{2}} \right)^{\frac{3}{2}} + \left( b_{R} p_{R} \right) \right] \left[ \left( b_{H} p_{H_{2}} \right) + K_{2} \left( b_{H} p_{H_{2}} \right)^{\frac{1}{2}} + K_{3} K_{2} \right]} \quad (1)$$

où les données expérimentales connues sont la vitesse *V* et les pressions partielles d'hydrogène et du *n*-heptane ($b_{H_2}$ et $p_R$ respectivement), et les constantes à déterminer sont au nombre de 5 : $b_H$, $b_R$, le produit $K_{MFM} K'_1$, $K_2$, $K_3$.

Les données expérimentales lissées se présentent sous la forme de courbes donnant *V* en fonction de $p_{H_2}$, et ce, à 5 valeurs de $p_R$ (fig. 2).

Le solveur d'Excel 4.0 a été utilisé dans un premier temps pour l'optimisation de ces 5 paramètres. Les différents essais réalisés ont montré que, d'une part, l'optimisation conduisait à différents jeux de paramètres en fonction des valeurs de départ fixées (solution non unique), et que d'autre part, les écarts entre les vitesses ainsi calculées et mesurées expérimentalement étaient trop importants pour être acceptés. Nous avions conclu, à ce stade, que l'outil d'optimisation fourni avec Excel 4.0 n'était pas adapté à ce problème, et qu'il fallait donc revoir la méthode d'optimisation, ainsi que le traitement mathématique de l'équation de vitesse.
<figure><img src="image_3.png" /><figcaption>Figure 2
Vitesse expérimentale à $p_{H_2}$ et $p_R$ fixé.</figcaption></figure>

## 2 APPROCHE RETENUE

Suite aux conclusions ci-dessus, un premier stage a été proposé conjointement par la DIMA et la division Cinétique et Catalyse pour l'étude mathématique et numérique du problème.

Ce premier stage [1] avait pour but de comparer différentes fonctions d'optimisation des logiciels MATLAB¹, MINPACK, IMSL et d'appliquer une méthode analytique sur laquelle on reviendra plus loin. Nous avons, grâce à ce stage, cerné différents problèmes:

- MINPACK et IMSL ont du mal à optimiser les 5 paramètres² $b_R$, $b_H$ le produit $K_{MFM} K_1$, $K_2$ et $K_3$ sans contraintes connues sur celles-ci ;

- avec la méthode analytique, des résultats satisfaisants ont été obtenus à $p_{H_2}$ fixé.

Un deuxième stage a été alors proposé [4], pour une analyse plus approfondie de la méthode analytique tout en s'appuyant cette fois-ci sur le « toolbox optimisation » de MATLAB [3]. Plus précisément, deux approches analytiques ont été étudiées grâce à deux écritures différentes de l'expression de vitesse : l'une non linéaire et l'autre linéaire qui permet d'appliquer des résultats mathématiques classiques sur le système matriciel qui en résulte.

Dans chacune des deux méthodes, l'indétermination du problème a été prouvée grâce à des arguments mathématiques simples et rigoureux.

### 2.1 Première expression analytique de la vitesse

Pour simplifier l'expression de vitesse définie plus haut, on pose :

$$ 
\left\{ 
\begin{array}{l} 
x_1 = b_H \\ 
x_2 = b_R \\ 
x_3 = K_{MFM} K_1' = K_1 \\ 
x_4 = K_2 \\ 
x_5 = K_3 
\end{array} 
\right. 
\quad 
\left\{ 
\begin{array}{l} 
p = p_{H_2} \\ 
q = p_R 
\end{array} 
\right. 
$$

La vitesse de réaction peut donc s'exprimer à $p$ et $q$ fixés sous la forme d'une fonction des 5 paramètres $x_i$ :

$$ V = f(x_1, x_2, x_3, x_4, x_5) = \frac{x_3 (x_2 q) (x_1 p)^{\frac{5}{2}}}{\left[ (x_1 p)^2 + (x_1 p)^{\frac{3}{2}} + x_2 q \left[ x_2 p + x_4 (x_1 p)^{\frac{1}{2}} + x_4 x_5 \right] \right]^2} \quad (2) $$

### 2.2 Seconde formulation de la vitesse sous forme linéaire

À partir de l'expression générale de la vitesse, on calcule $\frac{1}{\sqrt{V}}$ :

$$ 
\frac{1}{\sqrt{V}} = \frac{x_1^{\frac{3}{4}}}{(x_2 x_3)^{\frac{1}{2}}} \frac{p^{\frac{3}{4}}}{q^{\frac{1}{2}}} + \frac{x_1^{\frac{1}{4}}}{(x_2 x_3)^{\frac{1}{2}}} \frac{p^{\frac{1}{4}}}{q^{\frac{1}{2}}} + \frac{x_2^{\frac{1}{4}}}{(x_2 x_3)^{\frac{1}{2}}} \frac{q^{\frac{1}{4}}}{p^{\frac{1}{4}}} + \frac{x_2^{\frac{1}{2}}}{x_1^{\frac{1}{4}} x_3^{\frac{1}{2}}} \frac{q^{\frac{1}{2}}}{p^{\frac{1}{4}}} + \frac{x_2^{\frac{1}{2}}}{x_1^{\frac{1}{4}} x_3^{\frac{1}{2}}} \frac{q^{\frac{1}{2}}}{p^{\frac{1}{4}}} + \frac{x_2^{\frac{1}{2}}}{x_1^{\frac{1}{4}} x_3^{\frac{1}{2}}} \frac{x_4^{\frac{1}{2}}}{p^{\frac{1}{4}}} + \frac{x_2^{\frac{1}{2}}}{x_1^{\frac{1}{4}} x_3^{\frac{1}{2}}} \frac{x_4^{\frac{1}{2}}}{p^{\frac{1}{4}}} + \frac{x_2^{\frac{1}{2}}}{x_1^{\frac{1}{4}} x_3^{\frac{1}{2}}} \frac{x_5^{\frac{1}{2}}}{p^{\frac{1}{4}}} + \frac{x_2^{\frac{1}{2}}}{x_1^{\frac{1}{4}} x_3^{\frac{1}{2}}} \frac{q^{\frac{1}{2}}}{p^{\frac{1}{4}}} 
$$

En posant :

$$ z_1 = \frac{x_1^{\frac{3}{4}}}{(x_2 x_3)^{\frac{1}{2}}} \quad z_2 = \frac{x_1^{\frac{1}{4}}}{(x_2 x_3)^{\frac{1}{2}}} \quad z_3 = \frac{x_2^{\frac{1}{4}}}{x_1^{\frac{1}{4}} x_3^{\frac{1}{2}}} $$

$$ z_4 = \frac{x_2^2 x_4}{\frac{3}{x_1^4} \frac{1}{x_3^2}} \quad z_5 = \frac{x_2^2 x_4 x_5}{\frac{5}{x_1^4} \frac{1}{x_3^2}} $$

<footnote>(1) Le « toolbox optimisation » [3] n'étant pas à l'époque disponible, nous avons utilisé la fonction *fmin* [2], recherche directe du minimum sans information sur la dérivée.</footnote>

<footnote>(2) Nous ne pouvons pas mettre en cause les méthodes utilisées dans ces logiciels, il est bien évident que des contraintes sur les inconnues feraient mieux converger la solution qui, souvent pour des problèmes d'optimisation, dépend de l'initialisation du problème.</footnote>
alors :

$$
\frac{1}{\sqrt{V}}(p,q) = \left( \frac{\frac{3}{p^4}}{\frac{1}{q^2}} \right) z_1 + \left( \frac{\frac{1}{p^4}}{\frac{1}{q^2}} \right) z_2 + \left( \frac{\frac{1}{q^2}}{\frac{1}{p^4}} \right) z_3 \\
+ \left( \frac{\frac{1}{q^2}}{\frac{3}{p^4}} \right) z_4 + \left( \frac{\frac{1}{q^2}}{\frac{3}{p^4}} \right) z_5 \qquad (3)
$$

Nous remarquons que $\frac{1}{\sqrt{V}}$ a l'avantage d'être une
forme linéaire en de nouveaux paramètres $z_i^3$. On peut
ainsi traiter le problème par la résolution simple d'un
système linéaire, analyser les vecteurs propres du
système, et montrer que la projection de trois d'entre
eux était négligeable devant les deux autres, ce qui
explique bien l'indétermination du problème⁴.

3 RÉSULTATS

L'étude mathématique approfondie du problème a
permis de dégager deux résultats importants.

- Il y a un bon ajustement des données à q fixé (à 4 % en moyenne) sauf pour les valeurs de q extrêmes (q = 0,5 ; 5). Par contre, cet ajustement est moins bon à p fixé.

On représente sur la figure 3 les meilleures courbes
de vitesse obtenues à q fixé.

- L'indétermination des paramètres estimés a été démontrée par des arguments mathématiques : elle réside sur trois paramètres qui sont liés entre eux et qui nous donnent finalement un jeu de paramètres possible et améliorable (fig. 4).

Remarque :

L'étude générale de la vitesse ne donne évidem-
ment pas de résultats satisfaisants. Ceci résulte de
l'indétermination du problème que nous avons mise en
évidence.

<figure><img src="image_4.png" /><figcaption>Figure 3
Vitesse pour q = 0,5 ; 1 ; 2,2 ; 3,3 ; 5.</figcaption></figure>

<footnote>(3) En revanche, d'un point de vue pratique, les $z_i$ ne correspondent pas aux paramètres cinétiques, puisqu'il faudra déterminer tous les $z_i$ pour remonter aux valeurs des $x_i$.</footnote>

<footnote>(4) On pourra consulter [4] pour plus de précision.</footnote>
<figure><img src="image_5.png" /><figcaption>Figure 4

Comparaison de la vitesse théorique et de la vitesse perturbée de $\epsilon = 10$ sur les vecteurs $V_1, V_2, V_3, V_4, V_5$ à $q = 0,5$.</figcaption></figure>

## CONCLUSION

En résumé, nous avons une expression de vitesse qui est fonction de 5 paramètres à optimiser.

Une façon simple d'aborder le problème est de prendre le premier outil d'optimisation venu, et d'essayer de trouver la solution « miracle ». Cette approche ne s'est pas avérée satisfaisante.

Une étude analytique de l'expression à optimiser s'avère nécessaire, d'autant plus qu'aucune contrainte n'est connue sur un ou plusieurs paramètres.

La méthode proposée a eu le mérite de traiter le problème par une simple régression linéaire, et peut s'appliquer facilement à d'autres problèmes d'optimisation liés à la cinétique chimique. Deux solutions restent possibles à notre avis pour lever l'indétermination :

- donner d'autres relations entre paramètres cinétiques qui représentent des contraintes supplémentaires ;

- reprendre le formalisme chimique tel qu'il a été proposé au chapitre 1 pour déterminer une nouvelle expression de la vitesse.

Si l'on examine $z_i$, $i = 1 \text{ à } 5$, on remarque que leurs expressions dépendent toutes des paramètres $x_1 = b_H$,

$x_2 = b_R$, et $x_3 = K_{MFM} K_1'$, ce qui explique bien, comme on l'a vu plus haut, que grâce à un changement de variable, le problème initial d'optimisation à 5 paramètres, qui devient un problème d'optimisation à 2 paramètres, est indéterminé, et l'on pourra se demander si cette indétermination est liée au formalisme cinétique considéré. Il faut également considérer l'impact du « lissage » des courbes expérimentales et de l'incertitude des mesures sur les valeurs extrêmes de $p_R$.

## RÉFÉRENCES

1. Putot S. (1995) Optimisation des paramètres de l'équation de vitesse d'un catalyseur, partie I. *Rapport IFP 42461*.

2. Matlab user's guide. *The Math Works*, (décembre 1994).

3. Optimization toolbox user's guide, *The Math Works*, (mars 1994).

4. Belly C. et Kalfon P. (juin 1996) Optimisation des paramètres de l'équation de vitesse d'un catalyseur. Partie II. *Rapport de stage : Maîtrise d'ingénierie mathématique*.

5. Alvarez Herrera C.E. (1977) *Thèse de doctorat* de l'université de Poitiers.

Manuscrit définitif reçu en février 1998
# Annexe

## PREMIÈRE FORMULATION
DU PROBLÈME À OPTIMISER

On cherche $x_1, x_2, x_3, x_4, x_5$, tels que la somme des écarts quadratiques entre la vitesse expérimentale (obte- nue par mesures) et la vitesse théorique soit minimale.

On peut donc définir le problème de minimisation suivant :

Trouver $x = (x_1, x_2, x_3, x_4, x_5)$ tels que

$$ 
\begin{gathered}
J(x) \le J(y) \quad \forall y \in (R^+)^5 \\
\text{avec } J(x) = \sum_{p,q} (V_{\text{exp}} - V_{\text{théo}})^2
\end{gathered}
$$

on remarque que l'on peut écrire :

$$ V = \frac{\alpha q}{(1 + \gamma q)^2} \qquad (4) $$

$\alpha$ et $\gamma$ pouvant s'exprimer en fonction de $p$ et des paramètres recherchés :

$$ \alpha = \frac{x_2 x_3 x_1^2 p^2}{\left( (x_1 p)^2 + (x_1 p)^{\frac{3}{2}} \right)^2} \qquad (5) $$

$$ \gamma = \frac{x_2 x_1 p + x_2 x_4 x_1^2 p^2 + x_2 x_4 x_5}{(x_1 p)^2 + (x_1 p)^{\frac{3}{2}}} $$

Le critère $J$ à minimiser est alors :

$$ J(\alpha, \gamma) = \sum_{j=1}^{5} \left( \frac{\alpha q_j}{\left(1 + \gamma q_j\right)^2} - V_j \right)^2 \qquad (6) $$

où $V_{\text{exp}} = V_j$, $j$ variant de 1 à 5, pour chaque $p$ dans l'intervalle 6-22 bar.

En développant le carré de l'expression précédente et en dérivant par rapport à $\alpha$, on obtient la valeur optimale de $\alpha$ en résolvant l'équation :

$$ \frac{\partial J(\alpha, \gamma)}{\partial \alpha} = 0 $$

Cette valeur est unique car la courbe représentative de $J$ à $\gamma$ fixé est une parabole.

$$ \tilde{\alpha} = \frac{\sum_{j=1}^{5} \frac{q_j V_j}{(1 + \gamma q_j)^2}}{\sum_{j=1}^{5} \left( \frac{q_j}{1 + \gamma q_j} \right)^2} \qquad (7) $$

Il suffit de réinjecter dans le critère $J$ la valeur de $\tilde{\alpha}$ pour obtenir une fonction :

$$ F(\gamma) = J(\alpha, \gamma) \\ F(\gamma) = \frac{\left( \sum_{j=1}^{5} \frac{q_j V_j}{(1 + \gamma q_j)^2} \right)^2}{\sum_{j=1}^{5} \left( \frac{q_j}{1 + \gamma q_j} \right)^2} + \sum_{j=1}^{5} V_j^2 \qquad (8) $$

Le minimum $\tilde{\gamma}$ de $F$ est ensuite déterminé graphiquement à $p$ fixé.

## REFORMULATION SOUS FORME LINÉAIRE AU
SENS DES MOINDRES CARRÉS

Le critère à minimiser devient :

- déterminer $z_1, z_2, z_3, z_4$ et $z_5$ qui minimisent le critère $J$ suivant :

$$ J(z_1, z_2, z_3, z_4, z_5) = \sum_{p,q} \left( \frac{1}{\sqrt{V}} - \frac{1}{\sqrt{V_{\text{exp}}}} \right)^2 \qquad (9) $$

On peut écrire le critère de minimisation précédent sous une forme matricielle plus agréable⁵ :

$$ J(z_1, z_2, z_3, z_4, z_5) = Z' M Z - 2\Omega' Z + cte $$

et résoudre le système linéaire suivant pour minimiser le critère $j$ :

$$ MZ = \Omega $$

<footnote>(5) Les détails des calculs sont donnés dans [4].</footnote>