L’environnement de modélisation

Cast3M (site internet www-cast3m.cea.fr) est un code de calcul généraliste développé depuis le début des années 80 principalement au Commissariat à l’Energie Atomique (CEA). Ce code éléments finis permet la résolution d’équations aux dérivées partielles. D’abord développé pour résoudre des problèmes de mécanique, il a ensuite été adapté pour traiter des problèmes de mécanique des fluides, de thermique, d’acoustique et bien sûr d’écoulements en milieux poreux saturés et non saturés.Ce code est une boite à outils constituée d’un grand nombre d’opérateurs permettant de discrétiser et de résoudre des problèmes aux dérivées partielles [Verpeaux et al, 1989]. Ce code est constitué de deux niveaux de programmation. Deux langages de programma- tion sont donc utilisés. Les opérateurs sont codés en Esope, langage dérivé du fortran. Les jeux de données sont programmés en Gibiane, langage spécifique à Cast3M. La syntaxe de base d’opération dans Cast3M est la suivante :où Opérateur est un opérateur Cast3M programmé en Esope, Objet1 l’argument de l’opé- rateur et Objet2 le résultat de l’opérateur. Il existe dans la bibliothèque Cast3M des opé- rateurs de maillage, des opérateurs permettant de définir les propriétés des milieux ou des opérateurs de résolution, soit au total 500 opérateurs environ. Les jeux de données Cast3M sont donc une succession de commandes Gibiane permettant de définir le pro- blème, de le résoudre et de réaliser les post-traitements adéquats.

Une des particularités de Cast3M est que chaque utilisateur a la possibilité de dé- velopper ses outils pour résoudre son problème. Ce développement peut se situer aux deux niveaux évoqués précédemment et conduire à des développements d’opérateurs ou de procédures. Une fois acceptés par un comité de contrôle qui vérifie que les nouveaux outils sont cohérents avec les développements passés et la philosophie générale, ils peuvent être intégrés directement dans le code commercial. Cast3M est donc un code de calcul approprié pour la recherche puisque l’utilisateur contrôle l’intégralité de ses opérationset peut librement réaliser des développements si nécessaire. Dans le cadre de ce travail, aucun opérateur Esope n’a été développé. Les procédures nécessaires à la résolution de problèmes d’écoulements saturé et non saturé existaient déjà. Le travail a donc consisté à développer une procédure Gibiane pour modéliser les processus de surface et à implé- menter le couplage entre ces processus et les processus de subsurface.Plusieurs formulations numériques sont disponibles pour la résolution. D’abord déve- loppées en éléments finis classiques, des formulations éléments finis mixtes hybrides ou vo- lumes finis peuvent aussi être utilisées [Bernard Michel et al, 2004] [Le Potier et al, 1998] [Le Potier, 2005a] [Le Potier, 2005b]. La partie suivante présente la formulation EFMH utilisée pour la résolution des équations d’écoulements dans notre modèle. En revanche, l’équation de transport implémentée dans notre modèle est résolue en volumes finis. La partie la plus importante de ce travail ayant été réalisée sur les processus d’écoulements, la formulation volumes finis ne sera pas détaillée dans ce manuscrit.

Les éléments finis mixtes hybrides

La résolution des équations d’écoulements s’appuie sur la méthode des éléments fi-nis mixtes hybrides (EFMH) [Chavent et Jaffré, 1986] [Mosé et al, 1994] [Dabbene, 1998] [Le Potier et al, 1998]. Dans une formulation éléments finis classiques, le champ de charge et de vitesse sont calculés successivement. L’équation de Darcy est introduite dans l’équa- tion de conservation de la masse. On obtient alors une équation ne contenant plus que la charge piézométrique. Après résolution de cette équation, le champ de vitesse est dé- terminé à partir des charges calculées et d’une estimation de la valeur de conductivité hydraulique. Au cours de cette opération, le champ de conductivité hydraulique est sou-La formulation EFMH permet la résolution simultanée de l’équation de Darcy et de l’équation de conservation de la masse. On ajoute au système classique composé des équations de Darcy et de conservation de la masse des contraintes de continuité de la composante normale de la vitesse et de la charge aux interfaces entre éléments. Dans le cas d’une résolution en régime permanent, le système d’équations résolues est le suivant :Un multiplicateur de Lagrange, appelé trace de charge et égal à l’intégrale de la charge selon une face, est introduit pour résoudre ce problème sous contrainte. On reformule alors ce système sous forme matricielle en fonction de la seule inconnue trace de charge. Une fois celle-ci calculée, les champs de charge et de vitesse contraints sont déduits. La continuité des flux et des charges aux interfaces entre éléments est donc imposée naturellement dans cette formulation. Les points de calcul en EFMH sont différents des points de calcul en EF classiques. En EFMH, les traces de charge et les flux sont calculées aux centres des faces et la charge moyenne est calculé aux centres des éléments (voir Fig 2.1). Cette méthode est particulièrement intéressante dans le traitement de milieux très hétérogènes. Elle assure la conservation de la masse à l’échelle d’une maille et donc dans l’ensemble du domaine de calcul. Pour plus de détails sur les développements numériques et mathématiques, nous renvoyons à [Dabbene, 1993]. Les intérêts d’une telle méthode par rapport à une approche différences finis pourront être trouvés dans [Durlovsky, 1994].