Réseaux d’écoulement

Traduction réalisée par : Eric Rosa (Canada) et Vincent Cloutier (Canada)

5.1 Construction graphique des réseaux d’écoulement

Nous avons vu au Chapitre 2 qu’un système d’écoulement d’eau souterraine peut être représenté par un ensemble de surfaces équipotentielles auxquelles correspondent des lignes d’écoulement leur étant orthogonales. Lorsqu’une coupe transversale significative peut être choisie au sein du système tridimensionnel, l’ensemble de lignes équipotentielles et de lignes d’écoulement y étant représenté constitue un réseau d’écoulement. La construction de réseaux d’écoulement représente l’un des outils analytiques les plus puissants pour l’analyse de l’écoulement de l’eau souterraine.

Dans la Section 2.11 et dans la Figure 2.25, nous avons vu qu’un réseau d’écoulement peut être considéré comme la solution à un problème de condition limite bidimensionnel en régime permanent. La solution requiert une connaissance de la région au sein de laquelle a lieu l’écoulement, des conditions limites aux pourtours de la région et de la distribution spatiale de la conductivité hydraulique au sein de la région. Dans l’Annexe III, une méthode de solution analytique mathématique est présentée. Dans cette section, nous allons apprendre que les réseaux d’écoulement peuvent aussi être construits graphiquement, sans appel aux mathématiques sophistiquées.

Les systèmes homogènes et isotropes

Considérons d’abord une région d’écoulement homogène, isotrope et pleinement saturée. Pour des conditions de régime permanent dans une telle région, trois types de limites peuvent exister : (1) des limites imperméables, (2) des limites à charge constante et (3) des limites de nappe d’eau souterraine. Considérons en premier lieu l’écoulement à proximité d’une limite imperméable [Figure 5.1 (a)]. Puisqu’il ne peut y avoir aucun écoulement à travers la limite, les lignes d’écoulement adjacentes à la limite doivent lui être parallèles, et les lignes équipotentielles doivent rencontrer la limite à angle droit. En invoquant la loi de Darcy et en imposant un débit spécifique nul à travers la limite, nous sommes conduits à l’énoncé mathématique de la condition limite. Pour des limites étant parallèles au plan xz:

$\frac{\partial h}{\partial x} = 0 \hspace{1cm} \text{ou} \hspace{1cm} \frac{\partial h}{\partial z} = 0$ (5.1)

**Figure 5.1** Écoulement souterrain à proximité (a) d’une limite imperméable, (b) d’une limite à charge constante et (3) d’une limite de nappe d’eau souterraine.

De fait, toute ligne d’écoulement au sein d’un réseau d’écoulement constitue une limite imperméable, en ce sens qu’il n’y a pas d’écoulement à travers une ligne d’écoulement. Lors de la construction des réseaux d’écoulement, il est souvent préférable de réduire la taille de la région d’écoulement en considérant uniquement les portions de la région se trouvant d’un côté ou de l’autre d’une ligne de symétrie. S’il est clair que la ligne de symétrie représente également une ligne d’écoulement, la condition limite devant être imposée à la limite de symétrie est celle définie par l’Éq. (5.1).

Une limite sur laquelle la charge hydraulique est constante [Figure 5.1 (b)] est une ligne équipotentielle. Les lignes d’écoulement doivent rencontrer cette limite à angle droit, et les lignes équipotentielles adjacentes doivent être parallèles à la limite. L’expression mathématique est :

$h = c$ (5.2)

Sur la nappe d’eau souterraine, la charge de pression, $\psi$ , est égale à zéro, et la simple relation de charge, $h = \psi +z$ , donne :

h = z (5.3)

pour la condition limite. Tel qu’illustré à la Figure 5.1 (c), dans le cas de la recharge, la nappe d’eau souterraine ne constitue ni une ligne d’écoulement, ni une ligne équipotentielle. Il s’agit simplement d’une ligne où est variable mais connu h.

Si nous connaissons la conductivité hydraulique K du matériel dans une région d’écoulement homogène et isotrope, il est possible de calculer le débit au sein du système à partir du réseau d’écoulement. La Figure 5.2 consiste en un réseau d’écoulement complété pour le cas simple illustré à la Figure 2.25 (a). La zone entre deux lignes d’écoulement adjacentes est désignée comme un tube d’écoulement. Si les lignes d’écoulement sont équidistantes, le débit s’écoulant au sein de chacun des tubes d’écoulement est le même. Considérons l’écoulement au sein de la région ABCD dans la Figure 5.2. Si les distances AB et BC correspondent à ds et dm, respectivement, et si la perte de charge hydraulique entre AD et BC correspond à dh, le débit s’écoulant dans cette région, à travers une section de profondeur unitaire dans l’axe perpendiculaire à la page, est:

$dQ = K\frac{dh}{ds}dm$ (5.4)

**Figure 5.2** Réseau d’écoulement quantitatif pour un système d’écoulement très simple.

En conditions de régime permanent, le débit s’écoulant à travers n’importe quel plan de profondeur unitaire (par exemple à AD, EH ou FG) au sein du tube d’écoulement doit aussi être égal à dQ. En d’autres termes, le débit s’écoulant au sein de n’importe quelle partie d’un tube d’écoulement peut être calculé en considérant l’écoulement s’observant dans une partie adjacente.

Si nous décidons arbitrairement de construire le réseau d’écoulement en carrés, avec $ds = dm$ , alors l’Éq. (5.4) devient

$dQ = K dh$ (5.5)

Pour un système comportan tubes d’écoulement m, le débit total est :

$Q = mK dh$ (5.6)

Si la perte de charge totale dans la région d’écoulement correspond à H et qu’il y a n subdivisions de charge dans le réseau d’écoulement (H = n dh), alors

$Q = \frac{mKH}{n}$ (5.7)

Pour la Figure 5.2, m = 3, n = 6, H = 60 m, et à partir de l’Éq. (5.7), Q = 30 K. Pour K = 10^-4 m/s, Q = 3 × 10^-3 m³/s (par mètre de section perpendiculaire au réseau d’écoulement).

L’équation (5.7) doit être employée avec précaution. Elle est applicable seulement à des systèmes d’écoulement simples ayant une limite de recharge et une limite de décharge. Pour des systèmes plus compliqués, il est préférable de simplement calculer dQ pour un des tubes d’écoulement et de multiplier par le nombre de tubes d’écoulements pour obtenir Q.

La Figure 5.3 illustre un réseau d’écoulement montrant la percolation sous un barrage au sein d’une formation reposant sur une limite imperméable. Elle peut être utilisée pour établir trois points additionnels par rapport à la construction des réseaux d’écoulement.

**Figure 5.3** Percolation sous un barrage au sein d’une formation homogène et isotrope.

Les « carrés » dans tous les réseaux d’écoulement, à l’exception des plus simples, sont en fait des carrés « curvilinéaires », c’est-à-dire qu’ils possèdent des dimensions centrales égales, ou, vu sous un autre angle, qu’ils comportent un cercle étant tangent à leurs arrêtes.
Il n’est pas nécessaire que les réseaux d’écoulement possèdent des limites finies de tous les côtés; les régions d’écoulement qui s’étendent à l’infini dans une direction ou plus, comme c’est le cas pour la couche horizontale infinie dans la Figure 5.3, sont traçables.
Un réseau d’écoulement peut être construit avec un tube d’écoulement « partiel » sur son côté.

Pour le réseau d’écoulement illustré à la Figure 5.3, $m = 3\frac{1}{2}$ . Si H = 100 m et K = 10^-4 m/s, alors, puisque n = 6, nous avons Q = 5.8 × 10^-3 m³/s (par mètre de section perpendiculaire au réseau d’écoulement).

Dans un milieu homogène et isotrope, la distribution des charges hydrauliques dépend seulement de la configuration des conditions limites. La nature qualitative du réseau d’écoulement est indépendante de la conductivité hydraulique du milieu. La conductivité hydraulique intervient uniquement lorsque des calculs quantitatifs de débit sont réalisés. Il incombe également de noter que les réseaux d’écoulements sont adimensionnels. Les réseaux d’écoulement des Figures 5.2 et 5.3 demeurent valides, nonobstant qu’il soit considéré que la région d’écoulement fasse quelques mètres carrés ou des milliers de mètres carrés.

Le traçage des réseaux d’écoulements est en quelque sorte un art. On peut habituellement réaliser la tâche sur la base d’une approche par essai-erreur. Certains hydrogéologues deviennent extrêmement talentueux pour arriver à construire des réseaux d’écoulement acceptables rapidement. Pour d’autres, il s’agit d’une source de frustration continuelle. Pour un réseau d’écoulement dans un milieu homogène et isotrope, les règles de construction graphique sont fort simples. Nous pouvons les résumer comme suit: (1) les lignes d’écoulement et les lignes équipotentielles doivent se croiser à angle droit au sein du système, (2) les lignes équipotentielles doivent rencontrer les limites imperméables à angle droit; (3) les lignes équipotentielles doivent être parallèles aux limites à charge constante; et (4) si le réseau d’écoulement est tracé de façon à ce que des carrés soient créés dans une portion du domaine, alors, à l’exception possible de tubes d’écoulement partiels sur les côtés, des carrés doivent exister sur l’ensemble du domaine.

Les systèmes hétérogènes et loi de la tangente

Lorsque les lignes d’écoulement traversent une limite géologique entre deux formations possédant des valeurs de conductivité hydraulique différentes, elles sont réfractées, un peu comme c’est le cas pour la lumière passant d’un milieu à un autre. Toutefois, en contradiction à la loi de Snell, qui est une loi sinusoïdale, la réfraction de l’eau souterraine obéit à une loi tangentielle.

Considérons le tube d’écoulement illustré à la Figure 5.4. L’écoulement s’effectue depuis un milieu possédant une conductivité hydraulique K₁ vers un milieu possédant une conductivité hydraulique K₂, avec K₂ > K₁.

**Figure 5.4** Réfraction des lignes d’écoulement à une limite géologique.

Le tube d’écoulement possède une profondeur unitaire dans l’axe perpendiculaire à la page, et les angles et distances sont tels que ceux illustrés dans la figure. Pour un écoulement permanent, le flux intrant Q₁ doit être égal au flux extrant Q₂, ce qui s’exprime selon la loi de Darcy,

$K_1a\frac{dh_1}{dl_1} = K_2c\frac{dh_2}{dl_2}$ (5.8)

où dh₁ représente la perte de charge hydraulique selon la distance dl₁, et dh₂ représente la perte de charge selon la distance dl₂. Sachant que dl₁ et dl₂ délimitent les mêmes lignes équipotentielles, il est clair que dh₁ = dh₂; et à partir de considérations géométriques, a = b cos θ₁ et c = b cos θ₂. Sachant que b/dl₁ = 1/sinθ₁, et que b/dl₂ = 1/sinθ₂, l’Éq. (5.8) devient

$K_1\frac{\cos \theta_1}{\sin \theta_1} = K_2\frac{\cos \theta_2}{\sin \theta_2}$ (5.9)

$\frac{K_1}{K_2} = \frac{\tan \theta_1}{\tan \theta_2}$ (5.10)

L’équation (5.10) représente la loi de la tangente pour la réfraction des lignes d’écoulement de l’eau souterraine à une limite géologique, pour un milieu hétérogène. Connaissant K₁, K₂, et θ₁, il est possible de résoudre l’Éq. (5.10) pour θ₂. La Figure 5.5 illustre la réfraction de lignes d’écoulement pour deux cas avec K₁/K₂ = 10. Les lignes d’écoulement, comme si elles possédaient un esprit leur étant propre, préfèrent utiliser les formations de forte perméabilité telles des conduits, alors qu’elles préfèrent traverser les formations de faible perméabilité par le chemin le plus court. Dans les systèmes aquifère-aquitard possédant des contrastes de perméabilité de deux ordres de grandeur ou plus, les lignes d’écoulement tendent à devenir quasi horizontales dans les aquifères et presque verticales dans les aquitards. Si l’on considère la vaste gamme de valeurs de conductivité hydraulique présentée à la Table 2.2, il est clair que des contrastes de deux ordres de grandeur ne sont absolument pas hors du commun.

**Figure 5.5** Réfraction des lignes d’écoulement dans un système en couches (selon Hubbert, 1940).

Si l’on tente de tracer les lignes équipotentielles pour compléter les réseaux d’écoulement des diagrammes illustrés à la Figure 5.5, il devient rapidement clair qu’il est impossible de construire des carrés dans toutes les formations. Dans les systèmes hétérogènes, les carrés dans une formation deviennent des rectangles dans une autre formation.

Nous pouvons résumer les règles de construction des réseaux d’écoulement pour des milieux hétérogènes et isotropes comme suit : (1) les lignes d’écoulement et les lignes équipotentielles doivent se croiser à angle droit au sein du système; (2) les lignes équipotentielles doivent rencontrer les imites imperméables à angle droit; (3) les lignes équipotentielles doivent être parallèles aux limites à charge constante; (4) la loi de la tangente doit être satisfaite aux limites géologiques; et (5) si le réseau d’écoulement est tracé de telle sorte que des carrés sont créés dans une portion d’une formation, des carrés doivent exister au sein de cette formation et au sein de toutes les formations ayant la même conductivité hydraulique. Des rectangles seront créés dans les formations ayant une conductivité hydraulique différente.

Ces deux dernières règles rendent le traçage de réseaux d’écoulements quantitatifs exacts extrêmement difficile dans les systèmes hétérogènes complexes. Néanmoins, les réseaux d’écoulement qualitatifs, au sein desquels l’orthogonalité est préservée mais où la création de carrés n’est pas tentée, peuvent grandement aider à la compréhension d’un système d’écoulement de l’eau souterraine. La Figure 5.6 illustre un réseau d’écoulement tracé de façon qualitative pour le problème de percolation sous un barrage introduit à la Figure 5.3, mais cette fois-ci pour une formation sous-jacente en deux couches.

**Figure 5.6** Percolation sous un barrage, au sein d’un remblai hétérogène et isotrope.

Les systèmes anisotropes et la section transformée

Dans un milieu homogène mais anisotrope, la construction du réseau d’écoulement est complexifiée par le fait que les lignes d’écoulement et les lignes équipotentielles ne sont pas orthogonales. Maasland (1957), Bear et Dagan (1965) et Liakopoulos (1965b) présentent des discussions relatives aux principes théoriques qui sous-tendent ce phénomène, alors que Bear (1972) présente une évaluation exhaustive du cadre théorique. Dans cette section, nous aborderons principalement la réponse pratique ayant été développée afin de contourner le problème de la non-orthogonalité. Elle implique la construction du réseau d’écoulement au sein de la section transformée.

Considérons l’écoulement dans une région bidimensionnelle au sein d’un milieu homogène et anisotrope dont les conductivités hydrauliques principales sont K_x et K_z. L’ellipse de conductivité hydraulique (Figure 5.7) aura des axes $\sqrt{K_x}$ et $\sqrt{K_z}$ .

**Figure 5.7** Ellipse de conductivité hydraulique pour un milieu anisotrope avec *K_x*/*K_z* = 5. Les cercles représentent deux transformations isotropes possibles.

Transformons l’échelle de la région d’écoulement de telle sorte que les coordonnées dans la région transformée avec les axes de coordonnées X et Z soient associées à celles du système original xy selon

$X = x$

$Z = \frac{z\sqrt{K_x}}{\sqrt{K_z}}$ (5.11)

Pour K_x > K_z, cette transformation engendrera une expansion de l’axe vertical de la région d’écoulement. Elle résultera aussi en une expansion de l’ellipse de conductivité hydraulique jusqu’à l’obtention d’un cercle de rayon $\sqrt{K_x}$ (le cercle externe dans la Figure 5.7); et la section d’écoulement étendue et fictive agira comme si elle était homogène et de conductivité K_x.

La validité de cette transformation peut être supportée sur la base de l’équation d’écoulement en régime permanent. Dans le système de coordonnées original xz, pour un milieu anisotrope, nous avons, selon l’Éq. (2.69),

$\frac{\partial}{\partial x}\left(K_x\frac{\partial h}{\partial x}\right) + \frac{\partial}{\partial z}\left(K_z\frac{\partial h}{\partial x}\right) = 0$ (5.12)

En divisant par K_x, on obtient

$\frac{\partial^2 h}{\partial x^2} + \frac{\partial}{\partial z} \left(\frac{K_z}{K_x}\frac{\partial h}{\partial x}\right) = 0$ (5.13)

Pour la section transformée, nous obtenons, à partir de la seconde expression de l’Éq. (5.11),

$\frac{\partial}{\partial z} + \frac{\sqrt{K_x}}{\sqrt{K_z}}\frac{\partial}{\partial Z}$ (5.14)

En considérant la première expression de l’Éq. (5.11), et en appliquant l’opération associée à l’Éq. (5.14) aux deux différenciations de l’Éq. (5.13), on obtient

$\frac{\partial^2 h}{\partial X^2} + \frac{\partial^2 h}{\partial Z^2} = 0$ (5.15)

ce qui correspond à l’équation d’écoulement pour un milieu homogène et isotrope au sein de la section transformée.

Une transformation tout aussi valide peut aussi être obtenue par une contraction de la région selon la x direction grâce aux relations

$X = \frac{x\sqrt{K_z}}{\sqrt{K_x}}$ (5.16)

$Z = z$

Dans ce cas, l’ellipse de conductivité hydraulique sera transformée et correspondra au petit cercle dans la Figure 5.7, et la section d’écoulement transformée et fictive agira comme si elle était homogène et de conductivité K_z.

Avec le concept de section transformée maintenant maîtrisé, les étapes de construction graphique d’un réseau d’écoulement dans un milieu homogène et anisotrope deviennent évidentes, il faut : (1) transformer les coordonnées en utilisant l’Éq. (5.11) ou l’Éq. (5.16); (2) construire un réseau d’écoulement dans la section transformée fictive, selon les règles prévalant pour un milieu homogène et isotrope; et (3) inverser le ratio d’échelle.

**Figure 5.8** (a) Problème d’écoulement au sein d’une région homogène et isotrope avec $\sqrt{K_x}/\sqrt{K_z} = 4$ . (b) Réseau d’écoulement au sein de la section transformée isotrope. (c) Réseau d’écoulement au sein de la section réelle anisotrope. T, transformation; I, inversion.

**Figure 5.8** (a) Problème d’écoulement au sein d’une région homogène et isotrope avec $\sqrt{K_x}/\sqrt{K_z} = 4$ . (b) Réseau d’écoulement au sein de la section transformée isotrope. (c) Réseau d’écoulement au sein de la section réelle anisotrope. T, transformation; I, inversion.

La Figure 5.8 présente un exemple de la technique. Le problème de valeur limite illustré à la Figure 5.8 (a) consiste en une section verticale qui représente l’écoulement depuis un étang à h = 100 vers un drain à h = 0. Le drain est considéré comme faisant partie d’un groupe de nombreux drains parallèles étant tous implantés à la même profondeur et orientés perpendiculairement à la page. Les limites imperméables verticales sont « imaginaires »; elles sont créées par la symétrie du système d’écoulement dans son ensemble. Le drain inférieur représente une limite réelle; il représente la base du sol superficiel, lequel est sus-jacent à un sol ou à une formation rocheuse possédant une conductivité hydraulique étant plus faible par plusieurs ordres de grandeur. Si l’axe vertical est arbitrairement défini afin que z = 0 au niveau du drain et que z = 100 à la surface, à partir de l’expression h = ψ + z, et des valeurs h de fournies, nous avons ψ = 0 aux deux limites. À la surface, cette condition implique que le sol est tout juste à saturation. L’étang est de profondeur nulle. Au drain, le fait que ψ = 0 implique des conditions d’écoulement libre. Le sol au sein du domaine d’écoulement possède une conductivité anisotrope de K_x/K_z = 16. La section transformée de la Figure 5.8 (b) présente donc une expansion verticale de $\sqrt{K_x}/\sqrt{K_z} = 4$ . La Figure 5.8 (c) illustre le résultat de la transformation inverse, au sein de laquelle le réseau d’écoulement homogène et isotrope issu de la section transformée est ramené dans la région d’écoulement d’échelle réelle. Suite à l’inversion, la charge hydraulique en tout point (X, Z) de la Figure 5.8 (b) devient la charge hydraulique au point (x, z) dans la Figure 5.8 (c).

La taille de la section transformée dépend évidemment de l’équation utilisée pour la transformation (Éq. (5.11) ou Éq. (5.16), mais la forme de la région et le réseau d’écoulement en résultant sont identiques dans les deux cas.

Si les débits relatifs aux vitesses d’écoulement sont requis, il est souvent plus simple de réaliser les calculs dans les sections transformées. Cela soulève la question visant à déterminer quelle valeur de conductivité hydraulique doit être utilisée pour de tels calculs. Il est clair qu’il serait incorrect d’utiliser K_x, pour une section ayant subi une expansion verticale et K_z, pour une section ayant subi une contraction horizontale, tel que cela puisse être déduit de la Figure 5.7, alors qu’il en résulterait l’établissement de deux ensembles de calculs quantitatifs différents pour deux représentations équivalentes du même problème. De fait, l’expression correcte devant être utilisée est

$K' = \sqrt{K_x \cdot K_z}$ (5.17)

La validité de l’Éq. (5.17) repose sur la condition voulant que l’écoulement dans chacune des deux représentations transformées de la région d’écoulement doit être égal. La preuve requiert l’application de la loi de Darcy à un tube d’écoulement unique dans chacune des deux transformations.

L’influence de l’anisotropie sur la nature des réseaux d’écoulement de l’eau souterraine est illustrée à la Figure 5.9 pour le même problème de valeur limite que celui impliqué à la Figure 5.8. L’aspect le plus important des réseaux d’écoulement anisotropes [Figure 5.9 (a) et 5.9 (c)] concerne leur manque d’orthogonalité. Il nous semble que les techniques de transformation introduites dans cette section nous fournissent une explication indirecte mais satisfaisante de ce phénomène.

**Figure 5.9** (a) Réseaux d’écoulement pour le problème illustré à la Figure 5.8 (a) pour $\sqrt{K_x}/\sqrt{K_z} =$ (a) ¼, (b) 1, (c) 4 (selon Maasland, 1957).

**Figure 5.9** (a) Réseaux d’écoulement pour le problème illustré à la Figure 5.8 (a) pour $\sqrt{K_x}/\sqrt{K_z} =$ (a) ¼, (b) 1, (c) 4 (selon Maasland, 1957).

Il existe plusieurs situations pour lesquelles il pourrait s’avérer souhaitable de construire un réseau d’écoulement à partir de données piézométriques de terrain. Si les formations géologiques sont connues comme étant anisotropes, une grande précaution doit être appliquée dans l’évaluation des directions d’écoulement à partir des données équipotentielles. Si le réseau d’écoulement complet est requis, une section transformée doit être produite, mais il existe une approche de construction graphique pouvant s’avérer utile lorsque seules les directions d’écoulement à des points spécifiques sont requises. Dans la Figure 5.10, la ligne pointillée représente la tendance directionnelle d’une ligne équipotentielle à un certain point d’intérêt au sein d’un domaine xy. Une ellipse de conductivité hydraulique inversée est alors construite autour du point. L’ellipse possède des axes principaux $1/\sqrt{K_x}$ et $1/\sqrt{K_z}$ (au lieu de $\sqrt{K_x}$ at $\sqrt{K_z}$ , comme dans la Figure 5.7). Une ligne est tracée selon la direction du gradient hydraulique, elle recoupe l’ellipse au point A. Si une tangente à l’ellipse est tracée à partir de A, la direction de l’écoulement est perpendiculaire à cette tangente. En guise d’exemple d’application de cette approche de construction, il est possible de comparer les résultats de la Figure 5.10 avec l’intersection des lignes d’écoulement et des lignes équipotentielles dans la portion de centre-droit de la Figure 5.9 (c).

**Figure 5.10** (a) Évaluation de la direction de l’écoulement pour une région anisotrope avec $\sqrt{K_x}/\sqrt{K_z} = 5$ .

**Figure 5.10** (a) Évaluation de la direction de l’écoulement pour une région anisotrope avec $\sqrt{K_x}/\sqrt{K_z} = 5$ .

5.2 Les réseaux d’écoulement par simulation analogue

Pour l’écoulement au sein d’un milieu homogène et isotrope dans un système de coordonnées xz, les lignes équipotentielles d’un réseau d’écoulement consistent en une illustration de la solution, h(x, z), du problème de valeur limite décrivant l’écoulement en régime permanent au sein de la région. La construction du réseau d’écoulement consiste en une solution indirecte de l’équation de Laplace :

$\frac{\partial^2 h}{\partial x^2} + \frac{\partial^2 h}{\partial z^2} = 0$ (5.18)

Cette équation est l’une des équations différentielles partielles les plus communes en physique mathématique. Parmi les différents phénomènes physiques qu’elle décrit, on compte le flux de chaleur au sein de solides et le flux de courant électrique au sein de milieux conducteurs. Dans ce dernier cas, l’équation de Laplace prend la forme

$\frac{\partial^2 V}{\partial x^2} + \frac{\partial^2 V}{\partial z^2} = 0$ (5.19)

où V représente le potentiel électrique ou le voltage.

La similarité entre les Éqs. (5.18) et (5.19) révèle une analogie mathématique et physique entre le flux de courant électrique et l’écoulement de l’eau souterraine. Les deux équations sont développées sur la base d’une équation linéaire de flux, soit d’une part l’équation de Darcy et d’autre part la loi d’Ohm; et sur la base d’une relation de continuité, soit d’une part la conservation de masse du fluide et d’autre part la conservation de la charge électrique. Une comparaison de la loi d’Ohm,

$l_x = -\sigma\frac{\partial V}{\partial x}$ (5.20)

et de la loi de Darcy,

$v_x = -K\frac{\partial h}{\partial x}$ (5.21)

clarifie spontanément cette analogie. Le débit spécifique, v_x (débit par unité de surface), est analogue à la densité de courant, I_x (courant électrique par unité de surface); la conductivité hydraulique, K, est analogue à la conductivité électrique spécifique σ; et la charge hydraulique, h, est analogue au potentiel électrique, V.

L’analogie entre le flux de courant électrique et l’écoulement de l’eau souterraine constitue la base de deux types de modèles analogues s’avérant utiles pour la génération de réseaux d’écoulement quantitatifs. Le premier type implique l’utilisation de papier conducteur et le second type utilise des réseaux de résistances électriques.

Papiers conducteurs analogues

Considérons à nouveau le problème hydraulique initialement illustré à la Figure 5.8 et maintenant reproduit à la Figure 5.11 (a). L’analogue électrique [Figure 5.11 (b)] consiste en une feuille de papier conducteur découpée selon la même forme géométrique que le domaine d’écoulement de l’eau souterraine. Une source de courant est utilisée afin de générer une différence de potentiel entre les limites, et un appareil de mesure branché au circuit grâce à un voltmètre est utilisé afin de mesurer la distribution du potentiel au sein de la feuille de papier conducteur. Des limites à charge constante, telle que celle correspondant à la limite V = 100 sur la Figure 5.11 (b), sont créées grâce à l’utilisation de peinture argentée très conductrice; alors que les limites imperméables sont simulées par des bords non reliés au sein du modèle en papier conducteur. Il est généralement possible de rechercher les lignes équipotentielles de façon relativement efficace, de telle sorte qu’un réseau d’équipotentielles peut être généré rapidement.

La méthode se limite à des systèmes bidimensionnels homogènes et isotropes, mais peut gérer des géométries et des conditions limites complexes. Des variations dans la conductivité du papier commercialement disponible peuvent engendrer des erreurs aléatoires limitant l’exactitude quantitative de la méthode. Deux des applications les plus détaillées de la méthode sont celle de Child (1963), associée à l’analyse théorique de systèmes d’écoulement superficiels au sein de terrains drainés et celle de Tóth (1968), associée à l’évaluation de réseaux d’écoulement régionaux pour un secteur de l’Alberta.

**Figure 5.11** Réseaux d’écoulement simulés par analogie électrique. (a) Problème hydrogéologique à valeur limite en régime permanent; (b) analogie avec le papier conducteur; (c) analogie avec un réseau de résistances.

Réseaux de résistance analogues

L’utilisation d’un réseau de résistances à titre d’analogue électrique s’appuie sur les mêmes principes que ceux employés pour le papier conducteur analogue. Dans cette approche [Figure 5.11 (c)], le domaine d’écoulement est remplacé par un réseau de résistances connectées les unes aux autres aux points nodaux d’une grille. Le flux de courant électrique à travers chacune des résistances est analogue à l’écoulement d’eau souterraine au sein d’un tube d’écoulement parallèle à la résistance et possédant une surface transversale égale à la distance entre les résistances telle que multipliée par une profondeur unitaire. Pour qu’il y ait un flux électrique au sein d’une résistance individuelle, le terme I de l’Éq. (5.20) doit maintenant être considéré comme le courant, alors que le terme σ est égal à 1/R, où R représente la résistance. Comme c’était le cas pour le papier conducteur analogue, une différence de potentiel est imposée entre les limites à charge constante du modèle. Un appareil de mesure est utilisé afin d’évaluer le voltage à chacun des points nodaux du réseau et ces valeurs, lorsqu’elles sont enregistrées et tracées sous forme de contours, créent un réseau équipotentiel.

En modifiant les résistances au sein du réseau, il est possible d’analyser des systèmes hétérogènes et anisotropes par la méthode des résistances analogues. L’approche possède une exactitude et une adaptabilité supérieure aux modèles de papiers conducteurs, mais elle n’est pas aussi flexible que les méthodes numériques décrites dans la prochaine section.

Karplus (1958) fournit un guide détaillé pour les simulations analogues. Les analogues par réseaux de résistances ont été utilisés par Luthin (1953) afin de générer des réseaux d’écoulement pour une application de drainage, alors que Bouwer et Little (1959) les ont utilisés pour un système d’écoulement en zones non saturée et saturée. Bouwer (1962) a employé cette approche afin d’analyser la configuration de dômes piézométriques se développant sous des étangs de recharge.

L’utilisation la plus répandu des méthodes électriques analogues en hydrogéologie est celle prenant la forme de réseaux de résistance-capacitance pour l’analyse de l’écoulement en régime transitoire au sein des aquifères. Cette application sera discutée à la Section 8.9.

5.3 Réseaux d’écoulement par simulation numérique

Le champ des charges hydrauliques, h(x, z), qui permet la construction d’un réseau d’écoulement peut être généré mathématiquement selon deux approches, à partir du problème de valeur limite pertinent. La première approche s’appuie sur l’utilisation de solutions analytiques telles que discutées à la Section 2.11 et à l’Annexe III; alors que la seconde approche s’appuie sur l’utilisation de solutions numériques. L’application des méthodes analytiques se limite aux problèmes d’écoulement pour lesquels la région d’écoulement, les conditions limites et la configuration géologique sont simples et régulières. Comme nous le verrons dans cette section, les méthodes numériques sont largement plus versatiles, bien que leur utilisation requiert inévitablement l’utilisation d’un ordinateur numérique.

Les approches numériques sont approximatives. Elles sont basées sur la discrétisation du continuum composant la région d’écoulement. Lors du processus de discrétisation, la région est divisée en un nombre fini de blocs, chacun possédant ses propres propriétés hydrogéologiques, et chacun ayant un point nodal central au niveau duquel une charge hydraulique est définie pour le bloc en entier. La Figure 5.12 (a) illustre une grille de 7 × 5 blocs à points nodaux centrés ( $i =1$ à i = 7 selon la direction x et j = 1 à j = 5 selon la direction z) pour une région d’écoulement rectangulaire.

Examinons maintenant le régime d’écoulement à proximité de l’un des points nodaux internes – par exemple, au sein du bloc nodal, i = 4, j = 3, et de ses quatre voisins périphériques. Afin de simplifier la notation, nous allons reclasser les points nodaux tel qu’indiqué à la Figure 5.12 (b). Si l’écoulement se produit depuis le point nodal 1 vers le point nodal 5, nous pouvons calculer le débit, Q₁₅, selon la loi de Darcy :

$Q_{15} = K_{15}\frac{h_1 - h_5}{\Delta z}\Delta x$ (5.22)

pour un écoulement au sein d’une coupe transversale de profondeur unitaire perpendiculaire à la page. En s’appuyant sur la prémisse voulant que l’écoulement soit dirigé vers le point nodal central dans chacun des cas, nous pouvons écrire des équations similaires pour Q₂₅, Q₃₅, et Q₄₅. Pour un écoulement en régime permanent, la prise en compte de la conservation de masse du fluide implique que la somme de ces quatre écoulements soit égale à zéro. Si le milieu est homogène et isotrope, K₁₅ = K₂₅ = K₃₅ = K₄₅, et si nous choisissons arbitrairement une grille nodale carré de telle sorte que Δx = Δz, la sommation des quatre termes donne

$h_5 = \frac{1}{4}(h_1 + h_2 + h_3 + h_4)$ (5.23)

Cette équation est connue à titre d’équation en différences finies. Si nous revenons à la notation ij de la Figure 5.12 (a), l’équation devient

$h_{i,j} = \frac{1}{4}(h_{i,j-1} + h_{i+1,j} + h_{i,j+1} + h_{i-1,j})$ (5.24)

**Figure 5.12** (a) Grille de blocs à points nodaux centraux pour la simulation numérique de réseau d’écoulement. (b) Point nodal interne et points nodaux avoisinants. (c) Schémas en différences finies pour un point nodal interne et pour des points nodaux se trouvant le long d’une limite imperméable basale et sur un coin imperméable.

Dans cette forme, l’Éq. (5.24) s’avère adéquate pour tous les points nodaux internes au sein de la grille nodale. Elle relate un constat d’une simplicité élégante : en régime permanent au sein d’un milieu homogène et isotrope, la charge hydraulique à tout point nodal interne correspond à la moyenne des quatre valeurs voisines.

Un exercice similaire révélera que l’équation en différences finies pour un point nodal se trouvant le long de la limite basale, en assumant cette limite comme étant imperméable, prendra la forme

$h_{i,j} = \frac{1}{4}(h_{i-1,j} + h_{i+1,j} + 2h_{i,j+1})$ (5.25)

et, à un point nodal de coin

$h_{i,j} = \frac{1}{4}(2h_{i-1,j} + 2h_{i,j+1})$ (5.26)

Les équations (5.24) à (5.26) sont représentées schématiquement, de manière explicite, par les trois schémas en différences finies de la Figure 5.12 (c).

En résumé, il est possible de développer une équation en différences finies pour chacun des points nodaux de la grille nodale. S’il y a points nodaux N, il y aura N équations en différences finies. Il y aura N aussi inconnues – soit les valeurs N de h aux points nodaux N. Nous avons donc produit N équations algébriques linéaires et N inconnues. Pour des valeurs très petites de N, nous pourrions directement résoudre les équations, en utilisant une technique telle la règle de Cramer, mais lorsque N est grand, comme c’est inévitablement le cas dans les simulations numériques de réseaux d’écoulement, nous devons introduire une approche plus efficace, connue sous le nom de relaxation.

Restons fidèles au problème d’écoulement rapporté à la Figure 5.11 (a) et assumons que nous souhaitons développer son réseau d’écoulement par un moyen numérique. Dans la grille nodale de la Figure 5.12 (a), les points nodaux où la charge hydraulique est connue sont encerclés: h = 0 à i = 1, j = 3, et h = 100 à tous les points de la rangée j = 5. La relaxation implique des balayages répétés au sein de la grille nodale, du bas vers le haut et de gauche à droite (ou de toute autre manière consistante), en appliquant l’équation de différences finies pertinente à chacun des points nodaux ou la charge hydraulique est inconnue. Il faut au préalable assumer une charge hydraulique initiale à chacun des points nodaux. Pour le problème à l’étude, pourrait être attribué à titre de valeur initiale h pour chacun des points nodaux non encerclés. Lors de l’application des équations en différences finies, la valeur de h calculée la plus récente est toujours utilisée à chacun des points nodaux. Chaque balayage au sein du système est désigné par le terme itération, et la valeur h calculée de approchera de plus en plus la réponse finale après chaque itération. La différence de h à tout point nodal entre deux itérations successives est désigné par le terme résiduel. Le résiduel maximal au sein du système décroitra à mesure que les itérations sont réalisées. Une solution est atteinte lorsque le résiduel maximal atteint une valeur inférieure à un seuil de tolérance préalablement établi.

Afin de tester sa compréhension du processus de relaxation, le lecteur pourrait réaliser un certain nombre d’itérations dans la portion haut-gauche de la grille. Si la valeur initiale attribuée au point nodal, i = 2, j = 4 est, par exemple, égale à 50, alors la valeur après la première itération devient 62.5 et après la seconde itération devient 64.0. La valeur finale, atteinte après plusieurs itérations, se trouve entre 65 et 66.

L’approche par simulation numérique est en mesure de gérer toute forme de région d’écoulement et n’importe quelle distribution de conditions limites. Il est simple de développer les équations en différences finies pour des grilles rectangulaires où Δx ≠ Δz, ainsi que pour des distributions hétérogènes et anisotropes de conductivité hydraulique où les valeurs de K_x et K_z et varient d’un point nodal à l’autre. Dans l’Éq. (5.22), la technique de moyennage classique définirait K₁₅ = (K₁ + K₅)/2, où K₁ et K₅ font ici référence aux conductivités hydrauliques verticales aux points nodaux 1 et 5, et ces dernières peuvent être distinctes l’une de l’autre et être différentes des conductivités hydrauliques horizontales pour ces même points nodaux. La simulation numérique permet la construction de réseaux d’écoulement pour des cas trop complexes pour la construction graphique ou l’utilisation de solutions analytiques. La simulation numérique est presque toujours programmée pour l’ordinateur numérique, et les codes informatiques sont généralement écrits selon une forme générale de telle sorte que seules de nouvelles sources de données stockées sont requises pour gérer des problèmes d’écoulement largement différents. Il s’agit d’un avantage distinctif par rapport aux réseaux de résistances analogues, qui requièrent un démontage complet du matériel assemblé pour la réalisation d’une nouvelle simulation.

Le développement des équations en différences finies présentées dans cette section a été plutôt informel. Il est possible de débuter avec l’équation de Laplace et de procéder avec une approche plus mathématique vers un résultat identique. Dans l’Annexe VI, nous présentons un bref développement en ce sens. Il convient probablement de noter ici que le développement informel s’appuie sur la loi de Darcy et sur la relation de continuité afin d’obtenir les expressions en différences finies. Il s’agit des mêmes deux étapes ayant mené à l’établissement de l’équation de Laplace dans la section 2.11.

La méthode que nous avons désignée par le terme relaxation (d’après Shaw et Southwell, 1941) présente de nombreux alias. Elle est reconnue de façon variée comme la méthode de Gauss-Seidel, la méthode de Liebmann et la méthode des déplacements successifs. Parmi les nombreuses méthodes permettant la résolution d’équations en différences finies, il s’agit de la méthode la plus simple, bien qu’elle soit loin d’être la plus efficace. À titre d’exemple, si les charges calculées sont corrigées lors du processus de relaxation selon

$h_{corr}^k = \omega h^k + (1 - \omega)h_{corr}^{k-1}$ (5.27)

où h^k représente la charge calculée à la itération kth, et $h_{corr}^{k-1}$ représente la charge corrigée selon l’itération précédente, alors la méthode est connue comme la sur-relaxation successive et le nombre d’itérations requises pour l’atteinte d’une solution par convergence est significativement réduit. Le paramètre ω est désigné paramètre de sur-relaxation, et il doit avoir une valeur comprise dans l’intervalle 1 ≤ ω ≤ 2.

Il existe plusieurs textes qui s’avéreront fort utiles au modélisateur. McCracken et Dord (1964) fournissent une introduction à la simulation par ordinateur dans leur manuel Fortran. Forsythe et Wasow (1960) présentent leur message à un niveau mathématique plus avancé. Remson et al. (1971) discutent un large éventail de techniques numériques avec des références spécifiques à l’écoulement de l’eau souterraine. Pinder et Gray (1977) traitent le sujet à un niveau plus avancé.

Les méthodes numériques ont été introduites à la littérature associée à l’hydrologie par Stallman (1956) dans le cadre d’une analyse des niveaux d’eau à l’échelle régionale. Fayers et Sheldon (1962) ont été parmi les premiers à promouvoir l’utilisation de simulations numériques en régime permanent pour la réalisation d’études hydrogéologiques régionales. Remson et al. (1965) ont utilisé cette méthode pour prédire l’effet d’un éventuel réservoir sur les niveaux d’eau souterraine au sein d’un aquifère de grès. Freeze et Witherspoon (1966) ont généré de nombreux réseaux d’écoulement numériques dans le cadre de leur étude théorique de l’écoulement de l’eau souterraine à l’échelle régionale. Cette méthode était déjà largement appliquée dans le domaine du drainage agricole (voir Luthin et Gaskell, 1950) et dans l’évaluation des patrons de percolation dans les barrages en remblais (Shaw et Southwell).

Au cours des années récentes, la méthode des différences finies a été égalée en termes de popularité par une autre méthode de solution numérique connue sous la désignation de méthode en éléments finis. Cette approche mène aussi à un ensemble de N équations avec N inconnus pouvant être résolues par relaxation, bien que les points nodaux de l’approche en différences finies correspondent ici aux sommets d’un maillage triangulaire ou quadrilatéral irrégulier pouvant être défini par le modélisateur pour chaque application spécifique, en remplacement à la grille rectangulaire de la méthode en différences finies. Dans plusieurs cas, une grille nodale plus petite suffit et il en résulte des économies au niveau de la tâche effectuée par l’ordinateur. La méthode en éléments finis peut aussi gérer une situation que la méthode en différences finies n’est pas en mesure de traiter. La méthode en différences finies requiert que les principales directions de l’anisotropie dans une formation anisotrope soient parallèles aux directions du système de coordonnées. S’il y a deux formations anisotropes au sein d’un domaine d’écoulement, chacune avec des directions principales différentes, la méthode en différences finies est cernée, alors que la méthode en éléments finis peut fournir une solution. Le développement des équations en éléments finis requiert une sophistication mathématique excédant le cadre de ce texte d’introduction. Le lecteur intéressé est référé à Pinder et Gray (1977).

Les méthodes numériques, tant en différences finies qu’en éléments finis, sont largement utilisées pour la simulation informatique numérique des écoulements d’eau souterraine en régime transitoire au sein des aquifères. Cette application est discutée à la Section 8.8.

5.4 Réseaux d’écoulement saturés – non saturés

Il existe un autre type de réseau d’écoulement qui s’avère extrêmement difficile à construire graphiquement. Pour des problèmes d’écoulement impliquant à la fois des écoulements en conditions saturée et non saturée, les réseaux d’écoulement en condition de régime permanent sont généralement obtenus par simulation numérique. Considérons le réseau d’écoulement illustré à la Figure 5.13. Ce dernier est similaire au problème que nous avons analysé de façon répétée dans les dernières sections, en ce sens qu’il implique de l’écoulement vers un drain dans un système possédant des limites imperméables sur trois de ses côtés. Il se distingue néanmoins par le fait que l’axe vertical a été ajusté de façon à ce que la charge hydraulique à la limite supérieure du modèle présente maintenant une charge de pression inférieure à la pression atmosphérique. Cela signifie que le sol n’est pas saturé en surface, bien que pour qu’il puisse y avoir un écoulement vers le drain, une zone saturée doit exister en profondeur. Le réseau d’écoulement qualitatif présenté à la Figure 5.13 a été développé pour un sol dont les fonctions caractéristiques non-saturées correspondent à celles illustrées dans les graphiques de la Figure 5.13.

**Figure 5.13** Réseau d’écoulement saturé – non saturé dans un sol homogène et isotrope. Les graphiques présentent les fonctions caractéristiques non-saturées pour ce sol.

Ces courbes de conductivité hydraulique, K, et de teneur en eau, θ, en fonction de $\psi$ , correspondent aux courbes caractéristiques issues de la Figure 2.13.

Comme c’était le cas pour l’exemple unidimensionnel ayant été illustré schématiquement à la Figure 2.12, il existe trois types d’extrants associés à un réseau d’écoulement bidimensionnel saturé-non saturé obtenu par simulation numérique, en régime permanent. En premier lieu, il y a le patron de distribution des charges hydrauliques, h(x, z), qui permet la construction du réseau équipotentiel (traits pointillés longs à la Figure 5.13). En second lieu, il y a le patron de distribution des charges de pression, ψ(x, z) (traits pointillés courts à la Figure 5.13), qui s’avère particulièrement utile afin de définir la position de la nappe d’eau souterraine (isobar ψ = 0). En troisième lieu, il y a le patron de distribution des teneurs en eau, θ(x, z), qui peut être évalué à partir du patron de distribution des charges de pression ψ(x, z) et de la courbe θ(ψ) du sol. À titre d’exemple, une teneur en eau, θ, de 27 % est évaluée le long de la ligne ψ = -50 à la Figure 5.13.

Les lignes d’écoulement et les lignes équipotentielles forment un réseau continu sur l’ensemble de la région saturée-non saturée. Elles se croisent à angle droit sur l’ensemble du système. Un réseau d’écoulement quantitatif pourrait être tracé à partir de carrés curvilinéaires au sein de la zone saturée homogène et isotrope, mais de tels tubes d’écoulement ne pourraient former des carrés lorsqu’ils traversent la zone non saturée, même pour un sol homogène et isotrope. Pour une décroissance de la charge de pression (et de la teneur en eau), une décroissance de la conductivité hydraulique s’observe, et des gradients hydrauliques accrus sont requis afin d’engendrer les mêmes débits au sein des tubes d’écoulement. Ce phénomène s’observe dans les tubes d’écoulement présentés dans la section haut-gauche du réseau d’écoulement de la Figure 5.13, où les gradients s’accentuent vers la surface.

Le concept d’un système d’écoulement saturé-non saturé intégré a été introduit dans la littérature hydrologique par Luthin et Day (1955). Ces auteurs ont utilisé la simulation numérique et une maquette expérimentale remplie de sable afin d’évaluer le patron de distribution h(x, z). Bouwer et Little (1959) ont utilisé un réseau de résistances électriques afin d’analyser des problèmes associés au drainage et à l’irrigation, selon une configuration similaire à celle illustrée à la Figure 5.13. Les réseaux d’écoulement saturés-non saturés sont requis afin d’expliquer les nappes perchées (Figure 2.15 et Figure 6.11), et pour comprendre le régime hydrologique en flanc de colline et le mécanisme de génération de ruissellement (Section 6.5). Reisenauer (1963) et Jeppson et Nelson (1970) ont utilisé des simulations numériques afin d’évaluer le régime d’écoulement saturé-non saturé sous des étangs et des canaux. Leurs solutions ont été appliquées à l’analyse de la recharge artificielle de l’eau souterraine (Section 8.11). Freeze (1971b) a évalué l’influence de la zone non saturée sur la percolation d’eau au sein de barrages en remblais.

5.5 La surface d’exfiltration et l’écoulement de Dupuit

Surface d’exfiltration, point d’exfiltration et surface libre

S’il existe un système d’écoulement saturé-non saturé à proximité d’une limite correspondant à un exutoire libre, telle la berge d’un cours d’eau ou la portion aval d’un barrage en remblai, une surface d’exfiltration se développera à la limite. Dans la Figure 5.14 (a), BC est une limite à charge constante et DC est une limite imperméable. S’il n’existe pas d’apport en eau à la surface, AB agira aussi comme une limite imperméable.

**Figure 5.14** Développement d’une surface d’exfiltration le long d’une limite d’exutoire libre. (a) Réseau d’écoulement saturé-non saturé; (b) réseau d’écoulement à la surface libre; (c) Réseau d’écoulement de Dupuit-Forchheimer.

La nappe d’eau souterraine EF recoupe la limite d’exutoire AD au point d’exfiltration E. En amont du point E, le long de la ligne AE, les charges de pression non saturées $\psi$ sont inférieures à la pression atmosphérique et les écoulements vers l’atmosphère sont impossibles. De fait, AE agit comme une limite imperméable. La condition le long du segment ED correspond à h = z, comme c’est le cas à la surface de la nappe d’eau souterraine. Le problème relatif à la préparation d’un réseau d’écoulement dans de tels cas réside dans le fait que la position du point d’exfiltration séparant les deux conditions limites à la limite d’exutoire n’est pas connu a priori. En simulation numérique, il est nécessaire de fournir une prédiction initiale de la position de ce point d’exfiltration. La position correcte du point d’exfiltration est ensuite évaluée grâce à une série de simulations en régime permanent conduite par une approche essai-erreur.

La construction de réseaux d’écoulement quantitatifs en condition saturée-non saturée requiert une connaissance de la conductivité hydraulique, K, en zone saturée et non saturée et de la courbe caractéristique K(ψ) pour le sol. Dans plusieurs applications en ingénierie, incluant l’analyse de la percolation sous les barrages, de telles données sont rarement disponibles. Dans de telles situations, l’hypothèse voulant que l’écoulement en zone non saturée soit négligeable est habituellement retenue. Vu sous un autre angle, cette hypothèse impliquerait que la conductivité hydraulique du sol pour des teneurs en eau inférieures à la saturation est négligeable par rapport à sa conductivité hydraulique à saturation. Dans un tel cas, la limite supérieure du réseau d’écoulement correspond à la position de la nappe d’eau souterraine, et cette dernière devient elle-même une ligne d’écoulement. Lorsque ces circonstances spéciales prévalent, cette limite supérieure est désignée comme une surface libre. Les réseaux d’écoulement dans des systèmes circonscrits par une surface libre peuvent être construits de la façon usuelle, mais il existe une complication. La position de la surface libre (et non pas seulement du point d’exfiltration) est inconnue à priori. Les conditions limites sur la surface libre doivent satisfaire à la fois celles de la nappe d’eau souterraine (h = z) et celles d’une ligne d’écoulement (les lignes équipotentielles doivent lui être perpendiculaire). Sa position est généralement déterminée graphiquement, par essai-erreur. Certains écrits relatifs à l’ingénierie, tels que Harr (1962) ou Cedergren (1967), fournissent des astuces pour la construction graphique et fournissent de nombreux exemples de réseaux d’écoulement avec surface libre, en régime permanent.

La Figure 5.14 (b) représente le réseau d’écoulement avec surface libre équivalent au réseau d’écoulement saturé-non saturé de la Figure 5.14 (a). Une comparaison des deux réseaux d’écoulement confirme que notre décision de définir la nappe d’eau souterraine à titre de ligne d’écoulement constitue une approximation plutôt valable pour ce système d’écoulement particulier. La limite d’exutoire ED demeure désignée comme une surface d’exfiltration. Au sens pratique, nous allons rencontrer des surfaces d’exfiltration lorsque nous travaillerons sur l’hydrologie en flanc de colline (Section 6.5) et lorsque nous évaluerons la percolation au sein d’un barrage en remblais.

La théorie de Dupuit-Forchheimer sur l’écoulement à la surface libre

Pour l’écoulement au sein de systèmes à nappe libre circonscrits par une surface libre, une approche initialement proposée par Dupuit (1863) et subséquemment développée par Forchheimer (1930) est souvent évoquée. Cette approche s’appuie sur deux hypothèses : (1) les lignes d’écoulement sont considérées horizontales et les lignes équipotentielles verticales et (2) le gradient hydraulique est considéré comme étant égal à la pente de la surface libre et indépendant de la profondeur (constant avec la profondeur). La Figure 5.14 (c) montre le réseau équipotentiel pour le même problème que celui illustré à la Figure 5.14 (a), mais en prenant en compte les hypothèses de Dupuit. La construction de lignes d’écoulement robustes est désormais impossible. Cette situation paradoxale présente la théorie de Dupuit-Forchheimer pour ce qu’elle est, une approximation empirique du domaine d’écoulement réel. De fait, la théorie néglige la composante verticale des écoulements. En pratique, son intérêt réside dans la possibilité de réduire à une dimension un système bidimensionnel pour des besoins analytiques. Les calculs appuyés sur les hypothèses de Dupuit s’avèrent comparables avec ceux appuyés sur des méthodes plus robustes lorsque la pente de la surface libre est faible et que la profondeur du domaine d’écoulement libre est de faible envergure. Le débit Q au sein d’une coupe transversale de profondeur unitaire par rapport à la page dans la Figure 5.14 (c) est fourni par

$Q = Kh(x)\frac{dh}{dx}$ (5.28)

où h(x) représente l’élévation de la surface libre par rapport à la base du système d’écoulement au point x, le gradient dh/dx étant donné par la pente de la surface libre Δh/Δx en x. Pour un écoulement en régime permanent Q, doit être constant sur l’ensemble du système, et cela peut s’avérer véridique uniquement si la surface libre définit une parabole.

L’équation d’écoulement pour la théorie de Dupuit-Forchheimer dans un milieu homogène et isotrope peut être développée sur la base de la relation de continuité dQ/dx = 0. À partir de l’Éq. (5.28), nous obtenons

$\frac{d^2(h^2)}{dx^2} = 0$ (5.29)

Si un domaine d’écoulement tridimensionnel à nappe libre est réduit à un domaine d’écoulement bidimensionnel xy horizontal par l’évocation de la théorie de Dupuit-Forchheimer, l’équation d’écoulement dans un milieu homogène et isotrope devient

$\frac{\partial^2(h^2)}{\partial x^2} + \frac{\partial^2(h^2)}{\partial y^2} = 0$ (5.30)

En d’autres termes, h² et non h pas doit satisfaire l’équation de Laplace. Il est possible de définir un problème de valeur limite sur la base de l’Éq. (5.30) et de le résoudre pour h(x, y) pour des domaines d’écoulement peu profonds et horizontaux par simulation analogue ou numérique. Il est aussi possible de développer une équation d’écoulement de type Dupuit-Forchheimer en régime transitoire pour des aquifères à nappe libre, avec qui h² remplace h du côté gauche de l’Éq. (2.77).

Harr (1962) discute les aspects pratiques de la théorie de Dupuit-Forchheimer avec un certain niveau de détail. Bear (1972) présente un développement théorique exhaustif. Kirkham (1967) examine les paradoxes de cette théorie et fournit des explications révélatrices. Cette approche est largement employée dans les applications d’ingénierie.

Lectures suggérées

CEDERGREN, H. R. 1967. Seepage, Drainage and Flow Nets. Chapter 4: Flow Net Construction, John Wiley & Sons, New York, pp. 148-169.

HARR, M. E. 1962. Groundwater and Seepage. Chapter 2: Application of the Dupuit Theory of Unconfined Flow, McGraw-Hill, New York, pp. 40-61.

KIRKHAM, D. 1967. Explanation of paradoxes in Dupuit-Forchheimer seepage theory. Water Resources Res., 3, pp. 609-622.

PRICKETT, T. A. 1975. Modeling techniques for groundwater evaluation. Adv. Hydrosci., 10, pp. 42-45, 66-75.

REMSON, I., G. M. HORNBERGER, and F. J. MOLZ. 1971. Numerical Methods in Subsurface Hydrology. Chapter 4: Finite-Difference Methods Applied to Steady-Flow Problems, Wiley Interscience, New York, pp. 123-156.

Problèmes

Considérez une coupe verticale rectangulaire ABCDA au sein d’un milieu homogène et isotrope, avec une limite supérieure AB, une limite inférieure DC, une limite AD du côté gauche et une limite BC du côté droit. La longueur du segment DC correspond au double du segment AD. Dessinez un réseau d’écoulement quantitatif adéquat pour chacun des cas suivants.
1. BC et DC sont imperméables. AB est une limite à charge constante avec . AD est divisé en deux segments égaux, la portion supérieure est imperméable et la portion inférieure est une limite à charge constante avec h = 40 m.
2. AD et BC sont imperméables, AB est une limite à charge constante avec . DC est divisé en trois segments égaux avec les extrémités de gauche et droite imperméables et le segment central est une limite à charge constante avec h = 40 m.

Prenez la coupe verticale ABCDA du problème 1 et transformez-la en trapèze en soulevant B verticalement de telle sorte que l’élévation des points C et D est de 0 m, A est à 100 m et B est à 130 m. les segments AD, DC et BC sont imperméables et AB représente la surface de la nappe d’eau souterraine, avec une pente constante (sur laquelle la charge hydraulique est égale à l’élévation).
1. Dessinez le réseau d’écoulement quantitatif adéquat pour cette situation. Identifiez les lignes équipotentielles avec les valeurs correspondantes de h.
2. Si la conductivité hydraulique de la région est de , calculez le débit total à travers le système en m^-5 m/s (pour une profondeur de 1 m perpendiculairement à la section).
3. Utilisez la loi de Darcy pour calculer la vitesse d’entrée ou de sortie de l’écoulement à chaque point où une ligne d’écoulement rencontre la limite supérieure.

1. Répétez les problèmes 2 (a) et 2 (b) pour un cas homogène et anisotrope où la conductivité hydraulique horizontale est de 10^-4 m/s, et où la conductivité hydraulique verticale est de 10^-5 m/s.
2. Dessinez l’ellipse de conductivité hydraulique pour la formation homogène et anisotrope rapportée en (a). Démontrez par des tracés sur l’ellipse que la relation entre la direction de l’écoulement et celle du gradient hydraulique indiquée par votre réseau d’écoulement est correcte pour deux points de votre réseau d’écoulement.

Répétez le problème 2 (a) pour le cas où un drain permettant un écoulement libre (i.e. à la pression atmosphérique) est situé au centre du segment BC. Le drain est orienté perpendiculairement au plan du réseau d’écoulement.

1. Répétez les problèmes 1 (a), 1 (b) et 2 (a) pour un système à deux couches où la moitié inférieure du domaine d’écoulement possède une conductivité hydraulique qui est cinq fois plus grande que celle de la moitié supérieure.
2. Répétez le problème 1 (b) pour un système à deux couches où la moitié supérieure du domaine d’écoulement possède une conductivité hydraulique qui est cinq fois plus grande que celle de la moitié inférieure.

Dessinez un piézomètre dont la base est installée près du centre des domaines d’écoulement de chacun des réseaux d’écoulement construits aux problèmes 2, 3, 4 et 5. Montrez les niveaux d’eau qui existeraient dans ces piézomètres selon vos réseaux d’écoulement.

1. Redessinez le réseau d’écoulement de la Figure 5.3 pour un barrage ayant une largeur de 150 m à sa base et étant sus-jacent à une couche de 120 m d’épaisseur. Imposez h₁, = 150 m y h₂ = 125 m.
2. Répétez le problème 7 (a) pour un système à deux couches au sein duquel la couche supérieure faisant 60 m d’épaisseur est 10 fois moins perméable que la couche inférieure, laquelle fait aussi 60 m d’épaisseur.

Deux piézomètres distants de 500 m se terminent à des profondeurs de 100 m et 120 m au sein d’un aquifère à nappe libre. L’élévation de l’eau est à 170 m au-dessus de la base imperméable horizontale dans le piézomètre le moins profond et à 150 m dans le piézomètre le plus profond. Utilisez l’approximation de Dupuit-Forchheimer pour calculer l’élévation de la nappe d’eau souterraine à mi-chemin entre les deux piézomètres et pour calculer le débit s’écoulant dans une tranche de 10 m au sein de laquelle K = 10^-3 m/s.

Dessinez des réseaux d’écoulement sur un plan horizontal au sein d’un aquifère à nappe captive horizontal :
1. Pour un écoulement en régime permanent vers un puits en pompage (puits au sein duquel le niveau d’eau demeure constant).
2. Pour deux puits pompés en régime permanent à des débits identiques (i.e. produisant des charges identiques au sein des puits).
3. Pour un puits situé à proximité d’une charge constante linéaire.

Freeze and Cherry Groundwater Book

Chapitre 5 : Réseaux d’écoulement