PhreeqC est un logiciel de calcul géochimique de systèmes triphasique (liquide, solide, gaz) développé par l'USGS. C'est un logiciel libre. Les calculs chimiques peuvent être couplés avec des équations de transfert de matière par advection et/ou diffusion. La modélisation du transport s'effectue sous une dimension. Par ailleurs, le logiciel permet également de prendre en compte la cinétique des réactions. Les systèmes d'études sont fermés ou constitués d'un ensemble de systèmes fermés.
La simulation d'un système nécessite deux types de renseignement. Tout d'abord, une base de données référençant les données thermodynamiques des espèces chimiques. Plusieurs bases de données sont utilisables pour le logiciel. Le choix dépend, de la complexité de la simulation à effectuer mais aussi de l'objectif de la simulation. En exemple, la banque de données "iso.dat" renseigne des informations sur les isotopes les plus courants et sera utilisée en cas d'étude isotopique d'un système ou d'une réaction.
Le second renseignement consiste en les instructions qui permettent de réaliser les réactions recherchées. Le logiciel regroupent pas moins d'une vingtaine d'instructions permettant diverses études géochimiques.
A partir des différentes instructions, le logiciel permet de modéliser un grand nombre de réactions possibles:
De plus le système offre une certaine liberté pour définir un système:
Comme tout logiciel de modélisation PhreeqC possède certaines limites dans la modélisation des réactions chimiques. Le travail effectué par la suite est concerné principalement par les deux limites suivantes :
Notre étude se consacre à l'identification des dissolutions et précipitation ayant lieu au cours du mélange. L'étude des phénomènes de sorption et complexation en surface des minéraux nécessite une calibration expérimentale du logiciel. Par défaut de données, nous ne pouvons pas traiter ces réactions. De même la cinétique des réactions ne sera pas étudiée.
L'instruction "solution" permet de définir une solution aqueuse. Il est possible de choisir la température, le pH, le pe, et les éléments chimiques à introduire. Le logiciel chimique calcule la spéciation chimique des éléments selon le modèle d'équilibre thermodynamique, représenté à l'aide des lois d'action de masse (\ref{lois_masse}) :
\begin{equation}
\label{lois_masse}
K_{i} = a_{i}\prod_{m}^{M} a_{m}^{-C_{m,i}}
\end{equation}
Les activités des espèces en solution sont reliés aux concentrations par le modèle des solutions aqueuses (\ref{equa2}):
\begin{equation}
\label{equa2}
a_{i}=\gamma_{i}C_{i}
\end{equation}
Le coefficient d'activité de l'espèce i est déterminé suivant un modèle de solution aqueuse. Trois modèles sont ici disponible pour définir la solution aqueuse: le modèle de Davies (\ref{Davies}), utilisé par défaut car il ne nécessite pas d'entrée de données spécifiques à l'espèce; un modèle étendu de Debye Hückel ("Wataq model, \ref{Wateq}), ainsi que le modèle de Pitzer (\ref{Debye}) :
\begin{eqnarray}
log\gamma_{i} &= &-A\cdot z_{i}^2( \frac{\sqrt{I}}{1+\sqrt{I}}- 0.3 I) \qquad I < 0,5 \text{ mol.L$^{-1}$} \label{Davies} \\
log\gamma_{i} &= &-A\cdot z_{i}^2 \frac{\sqrt{I}}{1+rB\sqrt{I}}-\dot{B} I \qquad I < 1\text{ mol.L$^{-1}$} \label{Wateq} \\
log\gamma_{i} &= &-A\cdot z_{i}^2 \frac{\sqrt{I}}{1+\sqrt{I}} + \sum_{c}\sum_{a}B_{c,a}m_{c}m{a} \qquad I < 1\text{ mol.L$^{-1}$} \label{Debye}
\end{eqnarray}
Dans le cadre de notre étude, il sera important de vérifier la concordance du modèle utilisée avec la force ionique élevée de l'eau de mer étudiée. La figure 1, ci-dessous résume les domaines de validités des différents modèles de solutions aqueuses. Le modèle de Pitzer n'est pas représenté ci dessous. La modélisation établie par ce modèle est globalement supérieure à celles dérivées de la théorie de Debye Huckel, cependant il est plus difficile à appliquer à cause du nombre de paramètres à définir.
Figure 1 : Correspondance entre a force ionique de la solution et les coefficients d'activités des espèces (sources : d'après Garrels et Christ, 1965 et modilié par Merkel et Planer-Friedrich, 2008)
Outre la loi d'action de masse, le logiciel PhreeqC cherchera à satisfaire deux autres bilans chimiques : l'électroneutralité de la solution aqueuse, le bilan de matière pour chaque élément entré dans la solution aqueuse. Dans l'absolu, le système devrait être bouclé pour l'ensemble des éléments du système. Néanmoins, le logiciel offrent la possibilité d'utiliser certains éléments pour satisfaire le bilan en masse du système. Par exemple, le pH ou peut être utiliser pour satisfaire l'électroneutralité.
Pour prendre en compte les équilibres entre le gaz dissous dans les solutions aqueuses et le gaz présent à l'interface du liquide (par exemple équilibre O2 dissous et pression partiel de O2 dans l'air), la loi de Henry (\ref{Henry}) est utilisée. Elle relie la concentration de l'élément dissous en phase aqueuse à sa pression partielle dans la phase gaz :
\begin{equation}
\label{Henry}
P_{i} = C_{H} a_{i}
\end{equation}
PhreeqC permet de représenter deux types de réactions de dissolution / précipitation : les solides sont modélisés soit sous une forme de phases pure, soit sous la forme de solution solide. Les réactions sont alors modélisées par des lois d'action de masse telles que l'équation (\ref{lois_masse}). La modélisation diffère dans le calcul de l'activité du solide :
Lors de la réaction d'un solide, le logiciel dissout intégralement le minéral dans la solution et cherche ensuite à établir l'équilibre d'action de masse. A l'état final, il est possible de calculer l'indice de saturation (SI, \ref{SI}) à l'aide du produit d'activité ionique (PAI, \ref{PAI}) :
\begin{eqnarray}
SI &=& log(\frac{PAI}{K}) \label{SI} \\
PAI &=& a_{i}\prod_{m}^{M} a_{m}^{-C_{m,i}} \label{PAI}
\end{eqnarray}
Avec K, la constante d'équilibre de la réaction de dissolution. L'indice de saturation permet de déterminer si la solution est sous saturée ($ SI < 0 $); à l'équilibre ($SI = 0$) ou sursaturée ($SI > 0$). L'équilibre n'est pas obligatoirement atteint.
Le schéma suivant résume la modélisation de la solution aqueuse et des intéractions tri-phasiques :
Le modèle de transport de PhreeqC permet d'intégrer la diffusion et l'advection d'une solution aux calculs de réactions chimiques. Ces transferts sont modélisés à partir de l'équation générale une dimension d'advection diffusion avec terme source/puits :
\begin{equation}
\label{advec_diffu}
\frac{∂C_{i}}{∂t} = -\nu·\triangle C_{i} + D_{L}·\nabla^{2}C_{i} +\sum R_{i}
\end{equation}
Le modèle impose un coefficient de diffusion ainsi qu'une vitesse d'advection commune à toutes les espèces aqueuses. Concernant les bornes du milieu, trois types de conditions peuvent être imposés : Conditions aux limites de Cauchy, Dirichlet ou Neumann.
L'équation (\ref{advec_diffu}) est résolue par la méthode des différences finies avec un schéma centré dans l'espace et avancée temporellement. Les intéractions chimiques sont calculées séparément des termes de transport à chaque étape. Concrètement le logiciel calcule le transport advectif, puis l'équilibre chimique de la solution et enfin la diffusion. L'équilibre chimique est de nouveau calculé ensuite.
.