Infiltration depuis la surface plane du sol

Description de l’expérience

L’expérience consiste en une infiltration depuis la surface plane du sol, sur une aire circulaire (8 cm de diamètre) avec une pression d’eau négative (Figure 2). Les points d’injection (électrodes) sont disposés tous les 5 centimètres. Les mesures électriques sont obtenues en injectant un courant électrique via 2 électrodes et en mesurant le champ électrique induit via 2 autres électrodes. Les lieux d’injection et de mesure du potentiel (électrodes) sont, dans notre modèle, positionnés sur des nœuds du maillage. Le protocole d’acquisition comprend un mélange de mesures gradient et dipôle-dipôle, incluant des mesures directes et des mesures réciproques (Tableau 1). Les mesures dipôle-dipôle permettent d’obtenir une profondeur d’investigation verticale importante et une bonne sensibilité sous les électrodes, alors que dispositif gradient possède une très bonne sensibilité sous les électrodes de mesure. 

Modélisation hydraulique

Les flux d’eau dans un milieu poreux variablement saturé sont décrits par l’équation de Richards (1931) (équation 3). L’équation est résolue numériquement par le code SWMS_3D (Simunek et al., 1995) pour les conditions initiales et aux limites suivantes : ( , , ) ( , , 0) 1 2 3 i 1 2 3 θ x x x t, = θ x x x , (45) h(x , x , x t, )= h ( ,0 x , x ,t) 1 2 3 0 2 3 2 2 3 2 2 x + x ≤ r (46) ( ) 1 , , , 1 1 2 3 − ∂ ∂ = x h x x x t x1 = 0 2 2 3 2 2 x + x ≥ r (47) où θi [L3 L -3] est la teneur en eau volumique initiale, h0 [L] est le potentiel auquel l’infiltration est réalisée, et r [L] est le rayon du disque d’infiltration. Les équations (46) et (47) spécifient respectivement une condition de potentiel imposé pour les nœuds placés en surface du sol sous l’infiltromètre, et une condition de flux nul pour les nœuds de surface qui ne sont pas sous l’infiltromètre (pas d’évaporation envisagée). Des conditions de flux nul, similaires à l’équation (47), ont été appliquées aux faces latérales et dans le fond du domaine de simulation, en faisant l’hypothèse que l’évaporation, le drainage ou les transports latéraux hors du domaine sont négligeables pendant l’infiltration. Des dimensions du domaine de simulation très larges ont été choisies afin de s’affranchir d’effets de bords, et une discrétisation spatiale très fine (85731 nœuds) a été choisie pour obtenir une bonne précision sur les calculs (Tableau 2).

 Modélisation électrique

Le champ de potentiel électrique dans le profil de sol résultant de l’injection du courant électrique par un dipôle d’électrodes a été calculé en résolvant l’équation de Poisson suivante : ( ) ( ) ( ) ( ) n c i n i j b i j i = Iδ x x Iδ x x x x t, σ x t, x − − −         ∂ ∂ ∂ ∂ ϕ (48) où ϕn est le champ de potentiel électrique [V] du à l’injection du courant à l’électrode n au temps tj [T], δ est la fonction Dirac, xn et xc [L] sont respectivement les positions de l’électrode n et de l’électrode puits virtuelle (placée au centre du domaine), I [A] est l’intensité du courant, et σb [A V-1 L -1] est la conductivité électrique volumique du sol. L’équation (48) a été résolue en utilisant des conditions aux limites de flux nuls. Afin d’éviter des effets de bords, les limites du domaine de simulation électrique ont été placées à une distance largement supérieure à celle du domaine utilisé pour l’hydraulique (Tableau 2). Le domaine électrique comprend ainsi 8,5×105 nœuds arrangés dans un maillage écossais (accroissement des distances entre les nœuds au fur et à mesure de l’éloignement du centre). En raison de l’analogie qui existe entre le potentiel hydraulique et le potentiel électrique, l’équation (48) a été résolue, elle aussi, avec le code SWMS_3D (Annexe). Pour évaluer la précision de la solution numérique, une solution analytique a été calculée pour un champ uniforme de σb grâce à l’équation suivante qui donne la différence de potentiel électrique ∆V entre deux électrodes de potentiel P1 et P2 pour une injection réalisée entre les électrodes C1 et C2, (respectivement l’électrode source et l’électrode puits) :       − − − ) 1 1 ) ( 1 1 ( 2 1 2 3 4 π r r r r I ∆V = σ b (49) où r1 [L] est la distance entre C1 et P1, r2 [L] est la distance entre C2 et P1, r3 [L] est la distance entre C1 et P2, et r4 [L] est la distance entre C2 et P2. La différence entre les solutions numériques et analytiques a été trouvée toujours inférieure à 2 %.

 Unicité et stabilité de la solution

Comme l’a observé Russo et al. (1991), l’optimisation d’un nombre important de paramètres de Mualem-van Genuchten peut mener vers une non-unicité de la solution et vers une instabilité de la solution inverse. Ces auteurs reportent que, lors d’une infiltration 1D d’eau à saturation, la courbe d’infiltration cumulée au cours du temps n’apporte pas assez d’information pour obtenir un jeu unique de paramètres. Des résultats similaires ont été obtenus pour des expériences d’infiltration à une dimension lorsque les données d’infiltration étaient utilisées seules pour l’inversion des paramètres (van Dam et al., 1992, 1994 ; Eching et Hopmans 1993 ; Toorman et al., 1992). Dans le cas d’une infiltration à disque, Simunek et van Genuchten (1996) ont montré que la courbe de volume d’eau infiltré en fonction du temps ne contient pas assez d’information pour déduire une solution unique dans l’espace à 3 dimensions (α, n, Ks). Il est généralement recommandé, pour tester l’unicité de la solution, de résoudre le problème direct en partant de différentes estimations des paramètres initiaux. Pour vraiment s’assurer de l’unicité de la solution, il est fort utile de calculer la topographie de la fonction objective. Cependant, pour notre problème, la fonction objective est à 5 dimensions (autant que de paramètres) et, pour une raison de temps de calcul, il est impossible de calculer la valeur de la fonction objective sur l’ensemble du domaine des paramètres. On peut opter dans ce cas pour une solution intermédiaire, qui est de calculer les valeurs de la fonction objective dans des plans du domaine des paramètres, c’est-à-dire en ne faisant varier que 2 paramètres, les autres étant fixés à la valeur exacte (les plans contiennent la solution). De telles analyses de surface de la fonction objective donnent des indications sur l’unicité de la solution, mais ne peuvent permettre de s’en assurer complètement. Pour tester la stabilité de la solution, il est courant de réaliser des scénarios avec une contamination des données par un bruit aléatoire, en faisant l’hypothèse que l’erreur peut être représentée par une loi normale de type N(µ,σ). 

Paramétrisation

Pour des raisons de non-unicité, il est préférable de ne pas optimiser des paramètres de Mualem-van Genuchten qui ne jouent pas un rôle majeur dans la description des courbes de teneur en eau et de conductivité hydraulique en fonction du potentiel matriciel de l’eau, ou qui peuvent être mesurés facilement avec une bonne précision. De même que Simunek et van Genuchten (1996) ou Lambot et al. (2004), nous avons décidé de ne pas optimiser θs, lors des expérimentations numériques, pour les raisons suivantes : – On peut mesurer ce paramètre directement. – La saturation du milieu n’est pas atteinte par définition lors de l’infiltration avec le tensiomètre à disque (le paramètre correspond à un état physique non atteint lors de l’infiltration). Le paramètre θr peut être décrit comme le contenu en eau du sol pour une succion infinie. Mais, comme il n’y a pas de raison qu’une eau adsorbée soit présente dans le sol pour une succion infinie, ce paramètre n’a probablement aucun sens physique. C’est un paramètre d’ajustement, et son existence semble plutôt relever d’une fonctionnalité pratique qui permet pour certaines expériences de décrire avec une meilleure précision les mesures. De plus ce paramètre est très sensible au type d’expérience utilisé pour son estimation. Par exemple, certains auteurs ont montré que l’estimation de θr à partir de courbe de rétention en eau et de Chapitre 2 : Expérimentation numérique 55 fonction de pédo-transfert pouvaient mener vers des valeurs négatives, qui n’ont pas de sens physique (Haverkamp et al., 2005), ou vers des valeurs très élevées (autour de 0,2) pour lesquelles, comme le souligne Nimmo (1991), aucune preuve de la prétendue immobilité de l’eau n’a été apportée. Cependant il apparaît très probable que l’eau fortement liée aux constituants soit immobile à l’échelle de temps considérée lors de l’infiltration contrôlée. Ce paramètre ne sera donc pas ajusté lors de la procédure inverse, du moins dans le cadre des expérimentations numériques. Dans le cas des expérimentations réelles, on pourra l’estimer par inversion ou le laisser fixé à 0. En ce qui concerne l, il a été choisi, afin de simplifier le processus d’inversion et d’éviter la non-unicité de la solution, de ne pas l’optimiser. Sa valeur a été fixée à 0,5. Finalement, notre procédure d’estimation concerne les 3 paramètres Ks, α et n. Des valeurs typiques pour différents types de sols ont été choisies (Tableau 3) en se référant à la base de données de Carsel et Parrish (1988). Des valeurs des paramètres pétrophysiques de la relation de Rhoades issus de la littérature (Rhoades et al., 1976 ; Amente et al., 2000) ont été rapportées également (Tableau 4) et ont permis d’envisager l’ensemble du domaine possible pour les valeurs de ces paramètres. Les paramètres initiaux ont été choisis au milieu de l’espace des paramètres (Tableau 5). On a choisi de travailler sur 2 types de sol : un sol sablo-limoneux et un sol limoneux. Les paramètres hydrauliques retenus pour nos 2 sols sont ceux de Carsel et Parrish (1988).