Simulation numérique de l'étape de pyrolyse lors de la gazéification de la biomasse en lit fluidisé gaz-particule

Photo de l'installation au LGC Toulouse:

Réalisé par :

Vincent GROSJEAN (élève ingénieur à l'ENSIACET, Option Fluide, Energétique et Procédé)

Zhiya DUAN (élève ingénieur à l'ENSIACET, Option Fluide, Energétique et Procédé)

Encadré par :

Renaud ANSART (Chercheur au Laboratoire de Génie Chimique de Toulouse)

Hervé NEAU (Chercheur à l'Institut de Mécanique des Fluides de Toulouse)


Introduction

Contexte industriel

 

Projet GAYA

Le projet Gaya est un projet de recherche et développement qui vise à valider les choix technologiques innovant et les applications du biométhane. Financé par GDF SUEZ et l'ADEME, il a été lancé en 2010 à l'initiative de 11 partenaires issus du monde industriel, institutionnel et académique, en France et en Europe. 

LGC

Parmis ces partenaires, le laboratoire de génie chimique de l’ENSIACET à Toulouse (LGC) apporte au projet GAYA ses connaissances des lits fluidisés, ses compétances dans l'étude et la modélisation des phénomènes de pyrolyse et de gazéification au niveau de la particule. 

BEI

Dans notre cas, ce BEI vise à étudier l'effet de la pyrolyse dans un lit fluidisé qui joue le rôle d'un gazéifieur de biomasse. Cette étude est basée sur une unité pilote du procédé située au LGC. Nous utilisons le logiciel NEPTUNE_CFD pour réaliser la modélisation et la simulation. 

 

 

 

Réactions principales

Pyrolyse : commence à 400°C
  • Biomasse (g)→ αH2 (g) + βCO (g) + γCO2 (g) + δCH4 (g) + εC2H4 (g) + φH2O (g) + λtar (g) + κchar (s)
Vapogasification :  commence à 600°C
  • Char (s) + H2O (g) → CO (g) + H2 (g)        

La pyrolyse consiste en la dévolatilisation partielle de la biomasse. Par traitement thermique uniquement (T>400°C) jusqu'à 80% de la masse de biomasse est convertie en un mélange de gaz valorisable par différents procédés (méthanation, cogénération ...). La composition de ce gaz et la quantité produite sont variables selon la technologie du réacteur et la température de pyrolyse [2]. Cette réaction s'opère en un temps très court et induit un dégagement de gaz assez violent au sein du système réactionnel.

En outre le résidu solide non convertit en gaz consiste en un matériaux très riche en carbone et très poreux appelé le char. Ce char peut réagir avec de la vapeur d'eau pour se convertir de façon complète en CO et H2 supplémentaire lors d'une réaction plus lente appelée vapogasification.

Ce système de réactions est globalement endothermique [2] : il faut apporter de l'énergie au système pour entretenir les réactions. Le problème est que le procédé doit rester autonome en énergie pour avoir une raison d'être.

Procédé

Solution à l'endothermicité :

La solution choisie lors au sein du projet GAYA est de brûler une partie du char produit lors de la pyrolyse dans un autre réacteur afin de produire l'énergie nécessaire à la gasification. Pour transporter la chaleur du combusteur vers le gazéifieur on lie ces deux réacteurs par un lit fluidisé circulant dans lequel se déplace un média caloporteur (sable ou olivine suivant les besoins catalytiques). Le média "froid" (850°C) sort du gazéifieur par une déverse en entraînant du char non convertit vers le combusteur. Le char est alors brûlé et réchauffe le média environnant (950°C) qui retourne vers le gazéifeur en sortant par le haut du combusteur. L'unité pilote du LGC fait circuler un débit de média de 378 kg/h pour 20 kg/h de biomasse.

Schéma du procédé :

Autre rôle du média :

La biomasse utilisée dans l'unité pilote est sous forme de granulés de bois aggloméré d'environs s 2,5 cm de longs.

 

Longueur caractéristique

Masse volumique Vitesse minimale de fluidisation

Média

(sable ou olivine)

300µm 3040 kg/m3 4 cm/s

Biomasse

(granulés de bois)

2.5mm 600 kg/m3 5 m/s

Dans ce tableau on a calculé les vitesses minimales de fluidisation du média et de la biomasse grâce à la relation de Thonglimp [1].

  • $Re_{mf}=(31,6^2+0,0425Ga)^{\frac{1}{2}}-31,6$

On remarque que la vitesse minimale de fluidisation de la biomasse est considérable ce qui donne un autre rôle au média : mettre en mouvement la biomasse grâce aux multiples collisions des grains de média sur les granulés.

Mise en place du problème

Pyrolyse

Cinétique :

Le modèle cinétique choisit est un modèle de cœur rétrécissant [2]. Il est basé sur l'hypothèse que les particules de biomasse sont sphériques et ne changent pas de forme lors de la pyrolyse. L'avancement de la réaction est suivit par le déplacement d'un front de réaction à l'intérieur de la particule : le cœur composé de biomasse non convertie rétrécit pour laisser place au char. La cinétique est donc entièrement déterminée par la vitesse de déplacement de ce front :

$U_l= |\frac{dr_c}{dt}|=0,00144e^{\frac{-17180}{RT}}$

Ce modèle devrait toutefois être amélioré pour tout travaux ultérieurs car si il est bien adapté pour les particules sphériques de diamètre inférieur à 5mm il est beaucoup moins précis pour les granulés utilisés dans le procédé.

Stoechiométrie:

La stoechiométrie de la réaction a été simplifiée pour limiter le nombre de scalaire à traiter. Au final la pyrolyse ne produit que du char, du $H_2$, du $CO$, du $CH_4$ et du $H_2O$. Pour des besoins de modélisation on considèrera que chaque mole de gaz produite a une composition donnée par l'équation de la réaction simplifiée.

Nous avons retenu :

  • $x_{H_2}=0,158$
  • $x_{CO}=0,455$
  • $x_{CH_4}=0,137$
  • $x_{H_2O}=0,250$

Ce qui donne une masse molaire de gaz de pyrolyse à $W_{pyro}=0,026kg/mol$, valeur utilisée dans les bilans suivants.

 

 

Equations bilan

Cette partie se base en grande partie sur différents cas de calculs similaires au notre correspondant aux références [3], [4] et [6].

Transfert de masse :

Pour chacune des phases la conservation de la masse s'écrit de la manière suivante :

$\frac{\partial \alpha_k\rho_k}{\partial t}+\frac{\partial \alpha_k\rho_kU_{k,i}}{\partial x_i}=\Gamma_k$ avec $\sum \Gamma_k​=0$

On désignera les différente phase par leurs indice : 1, 2 et 3 pour la phase gaz, le média et la biomasse + le char respectivement. Le transfert de masse ne pouvant avoir lieu, vu la réaction en jeu, que de la biomasse vers le gaz on a :

$\Gamma_1=-\Gamma_3$ et $\Gamma_2=0$

$\Gamma_1$, le taux de transfert de la biomasse vers le gaz en $kg/m³/s$, s'écrit de la façon suivante :

$\Gamma_1=-n_3\frac{dm_3}{dt}$ avec $n_3$ le nombre de particules de biomasse par unité de volume en $m^{-3}$.

En écrivant $n_3=\frac{\rho_3\alpha_3}{m_3}$ on obtient $\Gamma_1=\frac{\alpha_3\rho_3(\rho_{biom}-\rho_{char})r_c^2U_l}{r_3^3\rho_{char}+r_c^3(\rho_{biom}-\rho_{char})}$

Il apparaît au final que cette expression est peut être plus adaptée aux réactions à masse volumique de particule constante et non à volume constant car l'expression de $\Gamma_1$ fait intervenir $\rho_3$. Pour tout travaux qui suivraient il serait peut être plus intéressant d'écrire $n_3=\frac{\alpha_3}{V_3}$ et ainsi obtenir $\Gamma_1=\frac{\alpha_3(\rho_{biom}-\rho_{char})r_c^2U_l}{r_3^3}$

En effet dans cette dernière formule $\Gamma_1$ ne dépend plus de $\rho_3$ qui est susceptible de varier. La formule utilisée dans ce BEI est la première des deux versions. $\Gamma_1$ est codé directement dans le fichier user ustrmv.F.

 

Transfert de quantité de mouvement :

Pour chaque phase k la forme non conservative de la conservation de la quantité de mouvement projetée sur la direction i s'écrit de la manière suivante :

$\alpha_{k}\rho_{k}\frac{\partial​ U_{k,i}}{\partial t}+\alpha_k \rho_k U_{k,i}\frac{\partial U_{k,i}}{\partial x_j} $

$=$

$-\alpha_k \frac{\partial​ P_g}{\partial x_i}+\alpha_k \rho_k g_i+\frac{\partial}{\partial x_i}[-\alpha_k \rho_k <u'_{k,i}u'_{k,j}>+\Theta_{k,ij}]+I'_{k,i}+[U_{\sigma,i}-U_{k,i}]\Gamma_k+\delta_{kp}\sum_{n \neq k}S_{kn,i}$

Tout ces termes à l'exeption de $[U_{\sigma,i}-U_{k,i}]\Gamma_k$ sont déjà pris en charge par NEPTUNE et ne nécessitent pas de traitement particulier. Le terme $[U_{\sigma,i}-U_{k,i}]\Gamma_k$ représente le gain de quantité de mouvement induit par le transfert de matière.  $U_{\sigma,i}$ est ici supposée égale à la vitesse de particules de la phase 3 au point considéré. C'est une hypothèse assez réductrice mais nous ne disposons pas de modèle plus élaboré, ce point est une piste d'amélioration de ce travail. Sous cette hypothèse on remarque que ce terme est non nul pour la phase gaz uniquement. En effet quand on considère la phase média $\Gamma_k$ est nul et quand on s'intéresse à la phase biomasse c'est $(U_{\sigma,i}-U_{k,i})$ qui est nul.

Sous NEPTUNE dans l'équation de conservation $\alpha_k \rho_k$ passe du coté des termes sources et le terme résultant, $\frac{1}{\alpha_{k}\rho_{k}}[U_{\sigma,i}-U_{k,i}]\Gamma_k$ se code dans le fichier ustsns.F sous la forme $\sum_{l=1}^3TSA(l)U_{l,i}$

Il en découle pour la phase gaz :

  • $TSA(3)=\frac{\Gamma_1}{\alpha_1 \rho_1}$
  • $TSA(2)=0$
  • $TSA(1)=-TSA(3)$
  • $TSB=0$

 

Transport des scalaires :

Pour bien décrire la composition des gaz produits il est nécessaire d'introduire dans la simulation trois scalaires représentant la fraction massique de chaque composés : $Y_{H_2}$, $Y_{CO}$, et $Y_{CH_4}$. On a pas besoin d'introduire $Y_{H_2O}$ car il peut être calculé à partir des trois autres.

Le transport d'un scalaire S quelconque porté par une phase k s'écrit de la façon suivante sous forme conservative :

$\frac{\partial}{\partial x_i}\alpha_k \rho_{k}S+\frac{\partial}{\partial x_i} \alpha_{k} \rho_{k} U_{k,i} S = \frac {\partial}{\partial x_i}(\alpha_{k}\rho_{k}D_{k,S}\frac{\partial S}{\partial x_i}) +\Psi_{k,S}$

Avec $\Psi_{k,S}$ terme source de S dans k.

Pour la fraction massique de l'espèce m on a : $\Psi_{k,S}=x_m\frac{W_m}{W_{pyro}}\Gamma_1$

En développant les termes du premier membre, en identifiant les termes relatifs à la conservation de la masse et en passant $\alpha_k$ on obtient l'équation sous forme non conservative telle qu'elle est traitée dans NEPTUNE :

$\rho_k[\frac{\partial S}{\partial t}+U_{k,i} \frac{\partial}{\partial x_i}S]=\frac{1}{\alpha_k}[\frac{\partial}{\partial x_i} (\alpha{k} \rho_{k} D_{k,S} \frac{\partial}{\partial x_i} S)+\Psi_{k,S}-\Gamma_k]$

$\frac{\partial}{\partial x_i} (\alpha{k} \rho_{k} D_{k,S} \frac{\partial}{\partial x_i} S)$ correspond à la diffusion du scalaire dans la phase porteuse. Pour les fractions massiques, le coefficient $D_{k,S}$ a été pris égal à $10^{-5}$ pour chaque espèces. On pourrait améliorer cette valeur en la calculant grâce à la théorie de Chapman-Enskog.

$\frac{1}{\alpha_k }[\Psi_{k,S}-\Gamma_k]$ est le seul terme nécessitant un traitement dans le fichier ustssp.F. Pour les fractions massiques, vu la forme de $\Psi_{k,S}$ se terme se simplifie en : $\frac{\Gamma_1}{\alpha_1}(x_m\frac{W_m}{W_{pyro}}-Y_m)$

Dans le fichier user ce terme s'identifie à TSA et on doit aussi coder $TSB=\frac{\partial TSA}{\partial S}$

Au final pour les fractions massiques on obtient :

  • $TSA_{m}=\frac{\Gamma_1}{\alpha_1}(x_m\frac{W_m}{W_{pyro}}-Y_m)$
  • $TSB_m=-\frac{\Gamma_1}{\alpha_1}$

Pour résoudre le problème on a besoin d'un dernier scalaire $\chi_d$ qui régit la conservation du nombre de particule par unité de masse. Il est définit comme $\chi_d=\frac{m_3(t=0)}{m_3(t)}$ ce qui dans notre cas se simplifie en $\chi_d=\frac{\rho_3(t=0)}{\rho_3(t)}$. Son transport se traite exactement de la même manière si ce n'est qu'il est porté par la phase 3, qu'il n'est pas sujet à la diffusion et que $\Psi_{3,\chi_d}=0$ vu qu'il n'y a ni attrition ni agglomération. On a au final :

  • $TSA_{\chi_d}=\Gamma_1\frac{\chi_d}{\alpha_3}$
  • $TSB_{\chi_d}=\frac{\Gamma_1}{\alpha_3}$

Propriétés physiques

Toute les propriétés physiques suivantes sont calculées dans le fichier usphyv.F

Masse volumique des particules de biomasse :

Cette grandeur découle directement de la définition de $\chi_d$ : $\rho_3=\frac{rho_3(t=0)}{\chi_d}$. On peut toutefois noter que cette masse masse volumique n'est qu'apparente, elle résulte de la structure des particule : un cœur de biomasse entouré d'une croute de char.

 

Rayon de cœur :

Ce rayon délimitant la frontière entre biomasse et char dans les particules de la phase 3 est utilisé pour calculer $\Gamma_1$. Son expression découle de l'écriture de la masse d'une particule de biomasse.

$r_c=[(\frac{\rho_{biom}}{\chi_d}-\rho_{char})\frac{r_0^3}{\rho_{biom}-\rho_{char}}]^\frac{1}{3}$

 

Masse volumique du gaz :

La masse volumique du gaz se calcule grâce à sa composition donnée par les fractions massiques.
$\rho_1=\frac{P}{RT}(\frac{Y_{H_2}}{W_{H_2}}+\frac{Y_{CO}}{W_{CO}}+\frac{Y_{CH_4}}{W_{CH_4}}+\frac{Y_{H_2O}}{W_{H_2O}})^{-1}$

Avec $Y_{H_2O}=1-Y_{H_2}-Y_{CO}-Y_{CH_4}$

Présentation de la simulation et des résultats

Géométrie et maillage

Géométrie 

La géométrie utilisée est une simplification du pilote du LGC. C'est un cylindre de 2,5m de haut et de 0,214m de diamètre avec un convergent en tête qui assure l'évacuation des gaz et l'absence de backflow dans la simulation. La géométrie est réalisée sur le logiciel Simail en 2D. Nous avons décidé dans un premier temps de traiter un cas de procédé batch pour simplifier les simulations et la géométrie ne présente donc pas d'entrée ou de sortie de solide.

Schéma de la géométrie 

Maillage 

Lorsque l'objectif de ce BEI est plutôt de comprendre la physique derrière les équations  plutôt que de donner un résultat parfait, nous décidons donc de réaliser la simulation en 2D afin d'économiser du temps de calcul. Pour les géométrie 2D, le logiciel NEPTUNE_CFD travaille en pseudo-2D, c'est-à-dire en 3D sur une épaisseur d'une maille. Nous utilisons au total 5880 mailles quadrangles pour la simulation. La taille caractéristique des mailles est de 1cm.

Conditions aux limites et initialisation

Conditions aux limites

Schéma des conditions aux limites de la géométrie 

              Tableau des types de conditions limites en pression
Phase Entrée Symétrie Paroi Sortie
CL en pression extr.P<i,w> dP=0 <i,s,w> extr.P <i,w> ddp=0 <o>

                                                

               Tableau des conditions limites des phases 
Position Phase Type de limite Vitesse Taux de présence
Bas gaz inlet (0;0;0,298) 1
média wall (0;0;0) 0
biomasse wall (0;0;0) 0
Haut gaz outlet (0;0;0) 1
média outlet (0;0;0) 10-12
biomasse outlet (0;0;0) 10-12

La vitesse d'injection de la vapeur est de 0,298 m/s qui correspond à 7 fois la vitesse minimale de fluidisation du média. Ce choix de vitesse assure une bonne fluidisation du lit.

      Tableau des conditions limites des scalaires
Scalaire Condition en entrée Condition en sortie Valeur
 YH2 dirichlet flux 10-12
YCO dirichlet flux 10-12
YCO2 dirichlet flux 10-12
Xd flux flux 10-12

 

 

 

Définition des paramètres

Paramètres physiques

Les propriétés physiques de l'eau ont été récupérées sur la base de données NIST [5].

On initialise les différentes phases particulaires à $10^{-12}$ éviter des erreurs lors du calcul. Le réacteur est ensuite chargé en biomasse et média dans le fichier usiniv.F avec 1,07kg de média et 0,057 kg de biomasse, soit un rapport média / biomasse de 20 ce qui correspond au point de fonctionnement de l'unité pilote du LGC (378 kg/h de média pour 20 kg/h de biomasse).

Tableau des propriétés physiques  
Phase Vapeur d'eau Média Biomasse
Etat  gas particle particle
Masse volumique (kg/m3) 0,19295 3040 600
Tréf (°C) 1123,5 1123,5 1123,5
Viscosité dynamique (Pa.s) 4,22.10-5 0 0
Diamètre (m) 0 0,0003 0,025
Cp (J/(kg.K)) 2378,1 1000 1400
Alpha init ---- 10-12 10-12

Paramètres des modèles

Tableau des paramètres des modèles 
Phase Vapeur d'eau Média Biomasse
Modèle de turbulence k-eps q2-q12 q2-q12
Couplage inverse turbulent non Small.inc. Small.inc.
Modèle de friction non fluxes fluxes
Modèle granulaire non fluxes fluxes
Modèle cinétique non uncorr.coll. uncorr.coll.
Modèle polydispersé Oui avec 0,64 en compacité maximale
CL paroi friction part.no slip part.no slip
Loi de trainée non Wen&Yu-Ergun Wen&Yu-Ergun

Le modèle de friction sur les phases particulaires impose une force de répulsion lorsque le système atteint un seuil de densité. Cela facilite la convergence du calcul car le système présentera moins souvent des points de trop forte compacité.

La CL en paroi sur ces mêmes phases est adaptée aux lits denses comme dans notre cas.

Vidéo

Cette vidéo correspond à 60 secondes de réaction sur un lit préfluidisé sans réaction pendant 5 secondes. La simulation a tourné pendant 26h sur un cœur et a réalisé 85193 itérations. La préfluidisation du lit permet d'éviter d'avoir un phénomène de "bouchon" lorsque la réaction commence.

Sur la partie correspondant au taux de présence du média on peut voir que le dégagement gazeux dû à la pyrolyse parvient à entrainer une partie du média vers la sortie du réacteur, la biomasse par contre reste dans le réacteur. On voit aussi qu'au cour de la réaction l'écoulement devient de moins en moins agité. Ceci s'explique par le fait que le dégagement gazeux ($\Gamma_1$) décroit avec le temps. En effet $r_c$ diminue linéairement lors de la réaction et donc les particules offrent un front de réaction de plus en plus réduit.

Analyse des résultats

Perte de média

Avec l'évaluation de la masse de média dans le réacteur, nous constatons qu'il y a une chute de la masse tout au début. Cette perte de média est à relier avec la vidéo où l'on peut voir du média s'échapper par le haut du réacteur. Cet envol d'une partie de la masse de média est dû à la grosse production de gaz pendant la pyrolyse. Presque 13 % de la masse de média fuit du lit.   

Évaluation de la masse gazeuse dans le réacteur

En traçant la variation de la masse gaz dans le réacteur, on voit que l'augmentation se présente essentiellement  dans les premières 5 - 6 secondes. Ceci s'explique par le fait que le gaz produit par la pyrolyse a une plus grande masse volumique que la vapeur d'eau. Après un certain temps la réaction n'est plus assez rapide pour compenser l'effet de vidange par la vapeur de fluidisation. A la fin de la simulation on remarque qu'il y a une masse résiduelle dans le réacteur :  la simulation n'a pas duré assez longtemps pour permettre à la vapeur de vidanger entièrement le réacteur.

Évaluation de la production de gaz

Nous avons aussi tracée l'évaluation de la production de gaz. Cette quantité s'obtient en faisant la différence entre le gaz sortant du réacteur et la vapeur injectée pour fluidiser.






 

H2

CO

CH4

H2O

Fraction massique

0,016

0,645

0,111

0,228

Production totale (kg)

4,02E-4

1,62E-2

2,79E-3

5,72E-3

 On peut calculer la fraction massique de chaque espèce dans le gaz produit grâce à la stoechiométrie de la réaction afin de récupérer la production de chaque espèce séparément.

Bibliographie

Références:

[1] Mehrdji HEMATI, Contacteurs fluide-solide (Caractérisaton, lits fixes et fluidisés), support de cours à l'ENSIACET 

[2] Mehrdji HEMATI, Gazéification de la biomasse, support de cours à l'ENSIACET   

[3] Jérôme LAVIEVILLE et al. (2007), NEPTUNE_CFD V1.0 - USER GUIDE, mode d'emploi de NEPTUNE_CFD EDF-CEA

[4] Hervé NEAU; Renaud ANSART et Olivier SIMONIN, Mise en œuvre d'un cas de calcul NEPTUNE_CFD V1.07@Tlse à l'ENSEEIHT RELEASE 2014-2015, mode d'emploi de NEPTUNE_CFD à l'ENSEEIHT     

[5] WebBook de Chimie NIST: Base de données standard de référence NIST numéro 69 Disponible sur: http://webbook.nist.gov/chemistry/

[6] Renaud ANSART; Olivier SIMONIN, réaction de carbochloration du zircon, AREVA, IMFT, aout 2011

Conclusion et perspectives

Nous avons pu voir dans ce BEI que la réaction de pyrolyse produisait un fort dégagement gazeux susceptible de faire s'envoler le média. Dans le cas du procédé réel ce phénomène devrait être moins visible car le procédé étant continu la biomasse est introduite progressivement et ne réagit pas d'un seul coup comme c'est la cas dans un procédé batch tel que celui simulé ici. Nous avons essayé de simuler un cas plus complexe se rapprochant du procédé réel mais nous ne sommes pas parvenu à ajuster les conditions limites et les propriétés du siphon (voir vidéo suivante) pour éviter que le réacteur se vide.

Hormis le fait que le dégagement gazeux serait moins violent la simulation d'un procédé continu présenterait aussi l'avantage de pouvoir accéder à différentes moyennes (taux de présence par exemple). On pourrait ainsi étudier l'influence de la pyrolyse sur la ségrégation de la biomasse, facteur déterminant pour son temps de séjour dans le réacteur et donc pour la qualité des gaz produits. En effet plus le temps de séjour est long plus la vapogasification à le temps d'agir sur le char. L'image suivante est le résultat d'une simulation de 30 secondes de fluidisation sur un réacteur non réactif. On peut voir la moyenne des taux de présence du média et de la biomasse. On vois bien que la biomasse se situe plutôt en surface. Il serait intéressant dans un prochain travail de simuler le procédé continu pour comparer les deux moyennes.

Un autre axe d'amélioration serait aussi de trouver un modèle cinétique de la pyrolyse plus adapté aux particules utilisées.