Modélisation de la diffusion anormale

Introduction 

Ce chapitre est consacré à la formulation mathématique du modèle que nous avons mis en oeuvre dans le cadre de cette thèse. D’une manière générale, les modèles de transport de masse dans les milieux poreux visent à comprendre et à reproduire la dynamique de la dispersion d’un soluté au sein de ces milieux. Dans le cadre des thématiques du laboratoire EMMAH concernant les transferts de solutés dans le sol, qui est en général un milieu poreux insaturé, un des objectifs finalisés de notre travail concerne la possibilité de prédiction de l’évolution d’un contaminant dans un tel milieu. Ceci nécessite le choix et la mise en oeuvre de modèles capables de fournir une interprétation cohérente des mesures en reproduisant les caractéristiques principales des observations. Comme il a été dit dans le premier chapitre de ce mémoire, aussi bien en laboratoire qu’en site naturel, une des propriétés principales du phénomène de transport de soluté dans les milieux poreux naturels, souvent mise en évidence expérimentalement, est le comportement non gaussien de la courbe de percée du traceur (qui représente la concentration du contaminant dans l’environnement). En particulier, un grand nombre d’observations rapportent un écart à la loi de Fick, matérialisé par une asymétrie marquée et par une traînée importante. Ces traînées peuvent traduire l’arrivée anticipée d’une fraction du panache de contaminant et l’arrivée tardive d’une autre fraction. Les interprétations avancées pour expliquer les longues traînées parfois observées s’appuient sur les propriétés intrinsèques de l’un ou l’autre des régimes supposés, et il y a peu de généralités à ce sujet. La difficulté d’avoir une connaissance complète des propriétés du milieu poreux conduit à le caractériser souvent de milieu hétérogène ou complexe. Cette hétérogénéité se trouvant à toutes les échelles et de manière variée, elle influence fortement le déplacement des par ticules et rend la migration sensible à ce qu’une moindre perturbation à petite échelle modifie les comportements à grande échelle. Dans la littérature, il est largement admis que la dynamique de transport dans un milieu possédant la propriété d’hétérogénéité peut être régie par une loi de diffusion anormale [75]. Dans ce contexte, Metzler et Klafter [75] ont introduit l’e.d.p appelée « Time Fractional Fokker Planck Equation » (TFFPE) qui est une extension de l’équation de Fokker-Planck classique. Des équations avec à la fois des dérivées fractionnaires en espace et en temps ont été proposées par Benson et al [76][77] comme par exemple le modèle « Space Fractional Fokker Planck Equation » (SFFPE). Une autre approche déja ancienne du transport dans les milieux poreux repose sur l’idée qu’à un instant donné, une fraction (inconnue) de la population de traceur est immobile et ne participe à la dynamique des déplacements qu’à un temps ultérieur. Cette idée très naturelle a donné lieu au modèle MIM [100] qui sépare la population de traceur en deux phases Mobile et Immobile. Ce modèle a permis de reproduire et de mieux comprendre les asymétries de nombreuses observations [41]. Néanmoins, le modèle MIM reste incapable de représenter certaines observations dans lesquelles les comportements des traînées semblent suivre des lois temporelles en inverse d’une puissance du temps t −γ (Haggerty et Gorelik [50] et Haggerty et al [49]). En particulier, c’est le cas d’observations en milieux poreux fortement insaturés. Le modèle que nous allons mettre en oeuvre constitue la version fractionnaire du modèle MIM, qui généralise ce dernier et possède des solutions se comportant aux grands temps comme des puissances. En effet, le modèle MIM fractionnaire vise à reproduire les comportements asymptotiques de la densité de probabilité de contaminants se comportant en loi de puissance du temps. L’equation de Fokker Planck fractionnaire appelée (FFPE) « Fractional Fokker Planck Equation » permet aussi de reproduire le comportement asymptotique de la densité de probabilité des contaminants se comportant comme une loi de puissance. De plus, elle est le sujet d’une littérature abondante. Cette équation de Fokker-Planck contient un paramètre de moins que le modèle MIM fractal, ce qui est intéressant, mais le degré de sophistication d’un modèle dépend aussi du nombre des paramètres. Le modèle MIM fractionnaire pour des temps courts est semblable au modèle d’advection dispersion (ADE) alors qu’aux grands temps il présente les mêmes comportements asymptotiques que l’équation de Fokker-Planck fractionnaire. Pour cette raison il correspond mieux que cette dernière à de nombreuses données, c’est pourquoi il a retenu notre attention. 

 Tout comme dans le cas du mouvement Brownien évoqué dans le chapitre précédent, le modèle MIM fractionnaire est la limite hydrodynamique d’un modèle de marche aléatoire dans lequel chaque particule effectue une marche au hasard marquée par des temps d’arrêt aléatoires qui modélisent la présence dans le milieu poreux de zones de stagnation. Un tel modèle de marche aléatoire est dit subordonné en temps. Plus précisément, pour chaque marcheur, le temps passé à bouger (ou encore le temps opérationnel) est différent du temps physique t. Ce phénomène fait partie de ce qu’on appelle la sous-diffusion ce qui est un terme général, puisqu’on l’utilise pour le transport chaque fois qu’un second moment n’évolue pas proportionnellement au temps. En effet, pour bien comprendre les mécanismes fondamentaux mis en jeu dans le phénomène de transport, il convient de se placer à trois niveaux d’échelle : – Echelle microscopique : à cette échelle, appelée aussi échelle de la particule, le phénomène de transport est modélisé par des marches aléatoires et dans la limite hydrodynamique, comme il y a une correspondance entre les phénomènes à l’échelle microscopique et les phénomènes à l’échelle macroscopique, la marche aléatoire tend vers la loi de transport macroscopique. Les modèles de diffusion normale correspondent à des marches aléatoires accumulant des sauts gaussiens régulièrement espacés dans le temps. – Echelle macroscopique ou les colonnes du laboratoire : à cette échelle le phénomène de transport est modélisé aussi bien par un processus stochastique que par des équations aux dérivées partielles, représentant l’évolution de sa densité. Ces e.d.p peuvent faire apparaître des opérateurs de dérivation et d’intégration d’ordre entier comme d’ordre fractionnaire, ainsi que des variables et des paramètres représentatifs des grandeurs moyennes au sein du milieu poreux. Ces modèles ont pour vocation d’être appliqués aussi à des échelles plus grandes, comme celles du terrain. Les modèles de diffusion normale donnent à cette échelle le mouvement Brownien et l’équation d’advection dispersion. – Echelle du terrain : cette échelle est plus grande que les deux autres. Lorsqu’on fait l’extension de l’approche macroscopique au milieu naturel, par exemple le sol, on se heurte à de sérieuses limitations. Le sol n’est jamais homogène et dans l’environnement naturel le caractère aléatoire des conditions aux limites renforce l’aspect stochastique des processus en jeu qu’il faut prendre compte pour adapter le modèle à cette échelle. On peut donc noter qu’un passage de l’échelle microscopique à l’échelle macroscopique permet d’établir que ces deux modèles sont liés de manière très étroite, puisque l’e.d.p régit l’évolution de la densité du processus stochastique auquel elle correspond. Cette étape de démonstration permettant de faire le passage à la limite hydrodynamique sera succinctement reprise, mais le lecteur intéressé aux détails peut consulter la référence .. Il est aussi question dans ce chapitre de la résolution des ces équations. La solution analytique et la solution numérique sont les deux méthodes de résolution de ces équations mais la première étant difficile voir impossible dans certains cas, la deuxième reste la meilleure possibilité pour avoir accès aux solutions des équations de manière générale. Après avoir présenté en détail ce modèle MIM, nous passerons aux modèles fractionnaires. 

Le modèle MIM classique

 Après avoir rappelé le principe du MIM, nous indiquerons une méthode de simulation numérique. 

Principe du MIM 

A l’origine, le modèle MIM décrit dans ce paragraphe a été développé par Coats et Smith [24] pour le génie pétrolier, puis adapté aux sols par Van Genuchten et Wierenga [109] et par bien d’autres parmi lesquels Gaudet et al.qui ont étudié les milieux poreux insaturés. L’idée de départ est qu’il existe dans un milieu poreux des pores en cul de sac, donc des zones d’eau stagnante. Comme on l’a dit plus haut, cette idée revient à considérer qu’une fraction de la population de contaminant est à l’état immobile cependant qu’une autre fraction est à l’état mobile. De plus ces deux fractions échangent continuement de la matière. Sous d’autres présentations, une variante du modèle MIM, le modèle à double milieu, consiste à rapporter ce concept de contaminants mobiles et immobiles aux propriétés du milieu poreux en considérant le milieu poreux comme la superposition de deux continua de porosité et de comportement différent : une de ces régions appelée macroporosité, permet l’écoulement du soluté, l’autre région, la microporosité est une zone de stagnation du fluide. La principale conséquence de ce concept de deux phases est la nécessité de considérer que la probabilité de trouver une particule donnée de contaminant au voisinage du point x, à l’instant t, se décompose en la probabilité de trouver cette particule mobile et la probabilité de trouver cette même particule immobile. En notant respectivement P, Pm et Pi les densités de probabilité correspondantes, on a alors P = Pm + Pi (3.1) Une autre implication du concept MIM est la nécessité de modéliser le passage de l’état immobile à l’état mobile. Dans le cadre du modèle MIM classique, ceci est fait en considérant une cinétique du premier ordre entre les populations mobiles et immobiles. Compte tenu de l’équation (3.1) et du fait que la fraction mobile de la population de contaminant est animé d’une dynamique d’advection dispersion, le modèle MIM classique peut se formuler de la manière suivante :    ∂tPm + ∂tPi = ∂x(D∂x(DPm − vPm)) + r ∂tPi = ω(KPm − Pi) (3.2) où D, v et r représentent respectivement le coefficient de dispersion, la vitesse d’advection et un terme source. Le paramètre ω est le coefficient d’échange dans la cinétique considérée. Ceci permet aussi d’écrire une relation entre les densités de probabilités Pm et Pi des populations mobiles et immobiles sous la forme Pm = Kω exp(−ωt) ∗ Pi P = (Id + exp(−ωt) ∗ Kω)Pm (3.3) où (∗) désigne un opérateur de convolution. Comme il s’applique à des fonctions causales, on parle de convolution « de Laplace » car dans ce cas cette convolution se traduit par la multiplication des transformées de Laplace. A partir des équations (3.1) et (3.3), on peut donner au modèle MIM une forme proche du modèle ADE (1.3) donné dans le premier chapitre. ∂tP = ∂x(D∂xHP − vHP) + r. (3.4) L’opérateur H est dans ce cas, H = (Id + Kω. exp(−ωt)∗) −1 (3.5) où Id est l’opérateur identité.