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)→ αH₂ (g) + βCO (g) + γCO₂ (g) + δCH₄ (g) + εC₂H₄ (g) + φH₂O (g) + λtar (g) + κchar (s)

Vapogasification :  commence à 600°C

• Char (s) + H₂O (g) → CO (g) + H₂ (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 H₂ 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.

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}$

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.