La notion d' UNITE HYDROGEOLOGIQUE essai de définition THESE présentée ò la Faculté des Sciences de l'Université de Neuchâtel pour obtenir le grade de Docteur ès-Scïences par LASZLO KIRALY Jury: A. BURGER, J.-B. GRIZE, E. RECORDON, J.-P. SCHAER Extrait du Bulletin du Centre d'Hydrogéologîe No 2 Neuchâtel 1978 IMPRIMATUR POUR LA THÈSE ..La...notion...d.\uni^ Jura (essai de définition )___ de M .P.n.s.ieur....LasZlo^^Kiralx UNIVERSITÉ DE NEUCHATEL FACULTÉ DES SCIENCES La Faculté des sciences de l'Université de Neuchâtel, sur le rapport des membres du jury, J.-Bl. Grize et E. Recordon autorise l'impression de la présente thèse sans exprimer d'opinion sur les propositions qui y sont contenues. Neuchâtel, Je...........3...mars.,.1.9I.8.......................................................... Le doyen: ir^,^ 'j J.-P. Schaer 83 BULLETIN DU CENTRE D'HYDROGEOLOGIE N° 2, 1977 LA NOTION D'UNITE HYDROGEOLOGIQUE ESSAI.DE DEFINITION par Laszlo KIRALY CHAPITRE PREMIER INTRODUCTION ET DEFINITION DU PROBLEME 1.1. Le but initial Le Centre d'Hydrogeologie de l'Université de NeuchStel étudie surtout l'nydrogéologie dea roches fissurées et karati- fiées du Jura. Cn même temps que la mise en route de travaux hydrogéologiques détaillés concernant des régions relativement restreintes, il fallait définir un cadre plus général facili- tant la corrélation des travaux particuliers et permettant une certaine planification des études à venir. L'organisation des régions déjà étudiées ou devant être étudiées en unités hydro- géologiques exigeait la définition explicite de cette notion. 84 Au début de notre travail la notion "unité hydrogéologique" (UHG) semblait intuitivement claire: il s'agissait d'un "piège aquifère" (AURDUZ, 1966), c'est-à-dire d'une série géologique "perméable" et "poreuse", délimitée par une ou plusieurs séries "imperméables", l'ensemble ayant une structure permettant l'ali- mentation et la formation, au moins temporaire, d'une nappe d'eau souterraine dans la série "perméable" (par exemple: une série calcaire reposant sur une série marneuse dans une structure syn- clinale, avec possibilité d'alimentation et de drainage des cal- caires) . Dans ce cas les relations entre qéoloqie, morphologie et "pièges aquifères" sont relativement simples (par exemple: calcaire = perméable; marne = imperméable; talweg = exutoire; affleurement des calcaires = réqian alimentaire) et l'on trouve qu'à une "unité géologique" (lithologique et structurale) cor- respond, généralement, une "unité hydroqéoloqique" (bassin ver- sant souterrain ou "piège aquifère"). La localisation et la dé- limitation de ces "pièges aquifères" peut se faire aisément à l'aide de cartes topographiques, de cartes géologiques et de cartes représentant la surface structurale des "imperméables", donc par des méthodes" surtout géologiques et, pour ainsi dire, en l'absence de renseignements détaillés sur les conditions hy- drauliques de la nappe contenue dans ces pièges aquifères. Une telle définition (ou essai de définition] correspondait, d'ailleurs, à notre connaissance de 1'hydrogéologie du Jura. En effet, si les études géologiques sont nombreuses dans le Jura, la connaissance des nappes d'eaux souterraines dans les calcaires est extrêmement fragmentaire. Lea forages et puits profonds sont peu nombreux dans les calcaires et très souvent les études hydro- ?éoloqiques de détail se réduisent a des essais de coloration, à quelques essais de pompaqe et à l'analyse chimique des sources. Il était donc compréhensible que, dans la mesure du possible, nous voulions utiliser l'information géologique abondante pour notre travail. Notre intention initiale était, par conséquent, 1'établis- sement d'une sorte de répertoire des différents types de pièges aquifères rencontrés dans le Jura, ainsi que leur délimitation 85 en utilisant surtout les caractéristiques géologiques et morpho- logiques des terrains. 1.2. Changement et redéfinition du but En coura de route nous nous sommes, toutefois, vite aperçu que la notion "unité hydraqéaloqique" ne pouvait pas être réduite à la notion "piège aquifère". En effet, en définissant formellement le terme général "unité hydrogéologique" par l'expression "classe d'équivalence dans le champ d'une variable caractérisant la qualité, la quan- tité ou le mouvement de l'eau souterraine", il est devenu possi- ble de définir plusieurs sortes d'unités hydrogéologiques; pra- tiquement autant que le nombre des variables envisagées {compo- sition chimique, température, vitesse d'écoulement, perméabilité, etc) et le nombre des relations d'équivalence définies dans le champ de ces variables (voir chapitre II). Cette définition est d'autant plus importante qu'elle prescrit explicitement la détermination du champ des variables hydrogéolo- giques, c'est-à-dire l'attribution d'une valeur bien définie (numérique ou nonil de cette variable à chaque "point" d'un vo- lume de terrain. Par conséquent, le premier problème, éminemment pratique et antérieur à la subdivision d'une région en unités hydroqéoloqiques, doit être la détermination concrète du champ de certaines variables, autrement dit la détermination de la dis- tribution spatiale de leurs valeurs dans la région étudiée. D'après ces quelques remarques on comprendra aisément que la détermination, directe ou indirecte, des champs hydroqéolo- qiques représentera, par la force des choses et apparemment en contradiction avec le titre, le problème central de notre travail. Est-ce que cela signifie qu'avec la nouvelle définition générale nous alIons nous désintéresser des facteurs géologiques ? Certainement pas, mais il faut préciser que si nous les prenons en considération, c'est seulement "faute de mieux", comme le mon- tre l'exemple suivant: 66 Les "pièges aquifères" déjà mentionnés sont, en fait, dé- terminés dans le champ des perméabilités et des porosités effi- caces (voir paragraphe 1.1.). Toutefois, si nous ne connaissons pas assez bien la distribution réelle de leurs valeurs dans une région alorB, "faute de mieux", nous utiliserons les facteurs géologiques grâce à une "salto mortale" méthodologique: il s'agit de "corré1er" les classes de perméabilité avec les classes litho- logiques et d'utiliser cette relation pour transformer le "champ des lithologies" (supposé connu) en un "champ des perméabilités".' Cette approximation du champ des perméabilités est très grossière et qualitative seulement, mais elle montre assez bien le rôle réel des relations entre facteurs géologiques et variables hydro- géologiques: que ces relations soient "théoriques-causales" ou "empiriques-statistiques", numériques ou non-numériques, on les utilise dans la mesure seulement où elles permettent la détermi- nation indirecte des champs hydragéologiques è partir des "champs géologiques" supposés connus. Elles jouent, en quelque sorte, le rôle des "fonctions d'interpolation" et leur analyse occupera une partie importante de notre travail. Nous arrivons maintenant à un choix fondamental: parmi toutes les "unités" possibles, quelle est l'unité hydrogéologi- que que nous jugeons la plus importante et quelles sont les va- riables hydrologiques dont il faut déterminer le champ en prio- rité ? L'unité hydrogéologique "idéale" serait celle qui représente le mieux le comportement d'une nappe (distribution des vitesses d'écoulement, distribution des potentiels, distribution des dé- bits positifs ou négatifs) dans des conditions naturelles et sous l'effet perturbateur de l'intervention humaine. Qr, TOTH (1962, 1963) définit le "hydrodynemicai flow system" comme une classe d'équivalence dans l'ensemble des lignes de courant et, è notre avis, ce "système d'écoulement hydrodynamique" représente l'unité hydrogéologique la plus importante, aussi bien pour le praticien que pour le chercheur (voir figure t ï. Il ne fait, en effet, aucun doute que seule la détermination du champ des potentiels hydrauliques et des lignes de courant 87 FIgM Systèmes d ' écoulements (d'après TOTH 1963 < GC < CS D 86 permet de comprendre véritablement ce qui Be passe dans une ¦nappe d'eau souterraine et seule la partition d'une nappe en systèmes d'écoulement juxtaposés ou superposés fournit un,ca- dre rationnel pour l'étude du bilan hydrique, du temps de sé- jour des eaux aboutissant aux exutoires, du champ thermique, de la composition chimique des eaux, de la propagation des pol- lutions et des influences mutuelles des ouvrages de captage. THOMAS et LEOPOLD (1964) vont jusqu'à affirmer que les proJBtB de recherche en hydrogéologie devraient être considé- rés comme "supportinq activities of the overall objective of defining numerically the regional flow patterns and super- imposed chemical systems". MEYBDOM (1966), plus modéré, admet pourtant la bonne correspondance entre composition chimique des eaux souterraines et leur position dans les systèmes d'é- coulement: "It is safe to say that the size of the flow sys- tem determines the prevailing hydrochemical facies as much, or even more sd, as does the minerelogical composition of the host rock". Etant donné que les systèmes d'écoulement permettent de "structurer", en quelque sorte, les nappes d'eau souterraines et forment un cadre idéal pour la corrélation des travaux hydro- géologiques particuliers ou spéciaux, nous fixons comme nouveau but du présent travail la détermination des systèmes d'écoule- rnen^ dans IeB roches carbonetées. fissurées et karstifiées du Jura. S'il est relativement facile de définir le but, l'exécu- tion du orojet pour y arriver présente de très grandes diffi- cultés et, dans la plupart des cas, la solution n'est obtenue qu'eu prix de "terribles simplifications" qui sont, toutefois, des simplifications "conscientes". Parmi les principales diffi- cultés nous devons mentionner: a) La détermination du champ des potentiels hydrauliques

. Soit une relation binaire R définie dans E. 5i R est re- flexive, symétrique et transitive, elle est une relation d'é- quivalence. Dans un ensemble E toute relation d'équivalence déter- mine une partition de E en sous-ensembles disjoints £., Ep...E appelés classes d'équivalence. L'ensemble des classes d'équivalence par rapport à R est un "ensemble-quotient" et il est désigné: E/R - {Er e2....En} Deux éléments X , X appartenant à la même classe d'é- quivalence sont équivalents par rapport à la relation R et la paire (X , X) est dite élément de la relation R : ( X1, X2 ) C R Considérons maintenant 1'ensemble-quotient E* . E/R » (E1. E2...En} Dans E nous pouvons définir une autre relation d* équi- ¦ i valence R qui partitionne E en classes d'équivalence E /R = JE1, E2. ..E1n } Chaque classe (ou "unité") E. appartiendra è une seule classe E . . Il est évident que les classes E . auront un ordj^e de grandeur supérieur è celui des classes E.: E/R est une classification plus "fine" que E /R De cette façon nous pouvons réaliser une classification "emboîtée" de E. c'est-à-dire nous pouvons définir des "unités" {classes d'équivalence) de plusieurs ordres de grandeur. I1 s'aqit là d'un processus d'abstraction très important qui consiste à considérer les classes d'équivalence {ou "unités") comme de nouveaux "individus" ou éléments et à définir des relations d'équivalence dans l'ensemble de ces nouveaux élé- ments . Si nous avons réalisé une classification emboîtée de E (par E/R, E'/R', etc), alors l'ensemble des "unités" E/R, E /R',... est partiellement ordonné par la relation d'inclu- sion (2 ¦ Ainsi, la définition des "unités" dans un ensemble E revient à la définition d'une relation d'équivalence R dans 94 2.2. Les relations d 'équivalence engendrées par des applica- tions 5oit un ensemble E (par exemple un volume dans 1'es- pace ponctuel tridimensionnel). Soit une application f qui fait correspondre à chaque élément d'un ensemble de E (ensem- ble de départ) un et un seul élément d'un ensemble d'arrivée A (ensemble numérique ou non). Dans ce cas, nous pouvons définir une relation d'équivalence R dans E de la façon suivante: deux éléments X , X £ E sont équivalents si leurs images parfsontégales. (X1, X2)CR ? Hx1) = f(x2);i HX1), f{x2) C a Si 1'ensemble d'arrivée A est partitionné en classes d'équivalence {par exemple, en intervalles disjoints dans un ensemble numérique) par une relation d'équivalence R , alors deux éléments X., X £ E sont équivalent si leurs images par f sont dans la même classe d* équivalence de A : (X1, X2) C R ? f(x1), f(x2) C R* Ainsi, à chaque "unité" de E correspond une classe d'équivalence de A. Cette définition de la relation d'équivalence assure une très grande souplesse dans son utilisation pratique. En variant l'application f ou l'ensemble d'arrivée A , nous pouvons partitionner ou "classifier" I'ensemble de départ de plusieurs "pointe de vue". Soit E un espace ponctuel et R l'ensemble des nom- bres réels. En définissant n applications f., i ~ 1..-n, de E dens R , ou une application F de E dans R, nous pouvons attacher un vecteur n-dimensionnel à chaque élément de E, c'est-à-dire nous pouvons transformer E 95 en un champ vectoriel (si n = 1, il s'agit d'un champ scalaire). Si R est une relation d'équivalence dans R , alors la relation d'équivalence P dans E est définie par (*r X2 ) C R ? T(X1), F(X2) C R* X1, X2 C E Les classes d'équivalence ("unités") de E sont des "volumes" dans 1'espace considéré et représentent une classi- fication de n "points de vue" à la fois- Les "volumes" de E appartenant à la même classe d'équivalence de E/R peuvent être toutefois fort éloignés les uns des autres dans 1'es- pace considéréI Pour certains problèmes {géologiques et hydrogéologiques) nous considérons les "volumes" disjoints comme des "unités" séparées, c'est-à-dire nous exigeons qu' une "unité" soit spatialement connexe ( "continue"). Dans ce cas nous devons combiner la relation d'équivalence R avec une relation d'équivalence topoloqique R- et la notion d'"uni- té" est définie comme suit : Deux éléments X , X £ E appartiennent à la même "unité" Si (X1, X2)O ? MX1), F(X2) £ R* .et S'il existe un chemin de X à X situé entièrement à l'intérieur de la même classe d'équivalence de E/R . 2.3. Résumé du schéma abstrait définissant la notion d'"unité" Une unité est une classe d'équivalence. Pour définir une classe d'équivalence nous devons avoir: - un ensemble référentiel E bien défini - une relation d'équivalence R dans E - une opération pratique permettant de décider si deux élé- ments X , X de E forment, oui ou non, un élément de R E/R = ^E ...E > est l'ensemble des unités de E par rapport à R . 96 En employant une application f pour engendrer la re- lation d'équivalence R , nous devons avoir: - un ensemble référentiel E bien défini ("départ") - un ensemble A • bien défini ("arrivée") - une application f de E dans (ou sur} A - une relation d'équivalence R* dans A - une opération pratique permettant de décider si un élément de A est, oui ou non, 1'image d'un élément de E par f. R est définie par: (X1, X?) £fl ^ fUj), f(X2) C R' et E/R a E ... .E_ est l'ensemble des unités de E par rapport à R . Pour certains problèmes nous exigeons l'existence d'un chemin de X à X situé entièrement dans la même classe d'é- quivalence de E/R, c'est-à-dire nous exigeons que l'unité soit connexe dans E : (X1, X2) C P et (X1, X2) C RT CHAPITRE III 97 SCHEMA GENERAL DE5 UNITES HYDROGE0LOGIQUE5 3.1. Introduction Dans le chapitre précédent nous avons vu que la défini- tion des "unités" dans un ensemble E revient à la défini- tion d'une relation d'équivalence R dans E , c'est-à-dire à une partition de E par R . La définition des "unités hydrogéologiques" (UGH) dans un volume de terrain E revient à la définition d'une rela- tion d*équivalence hydrogéologique R. dans E , c'est-à-dire à une partition de E par R. ng L'étude des unités Hydroaéalpgigués se réduit donc h une analyse empirique et, théorique des relations d'équiva- lence R. dans l'ensemble des volumes de terrains du point de vue hydrogéologique. Avant de continuer plus loin il faut insister sur deux faits importants: 1. Nous créerons les relations d'équivalence Ru selon nos hg besoins, pour la résolution de problèmes pratiques et théoriques. Cela veut dire que 1'interprétation concrète et 1'extension des unités hydrogéologiques changera avec nos problèmes théoriques et pratiques à résoudre. Seul le schéma de construction des unités hydrogéologiques reste invariant : il s'agit d'organiser un volume de terrain donné en classes d'équivalence à l'aide d'une relation d'équivalence hydrogéologique Rl . 2. L'affirmation "un volume de terrain E. (^ E est une unité1 hydrogéologique", n'a de sens que par rapport a une relation d ' équivalence R. bien définie et par rapport a un volume référentiel E bien délimité. 98 3.2. Quelques conceptions récentes des UHG G. CASTANY (1968, p. 76 et suite) avait inclu dans son ouvrage un paragraphe intitulé "Unités hydragéologiques. Bassins des eaux souterraines". Il écrit: "Le tracé des li- gnes de partage des eaux souterraines (crêtes hydrologiques) permet de délimiter les grandes unités hydrologiques, les bas- sins des eaux souterraines (bassins hydrogéologiques) et les grandes zones hydrogéologiques. Il repose sur l'étude des nap- pes elliptiques et du niveau de base souterrain ou superfi- ciel (émergences) de 1'écoulement souterrain. Les limites sont portées sur les cartes en courbes isopièzes. Ces délimitations sont la base des études hydrogéologiques ,régionales et des bassins versants représentatifs ou expérimentaux". L'UHG est définie à 1'aide de la surface piézométriqué d'une nappe d* eau : les crêtes de la surface piézométrique séparent les UHG et les eaux souterraines d'une UHG s'écoulent vers une zone collectri- ce commune. L'idée d'une nappe d'eau constituée est impliquée. L'extension verticale d'une UHG n'est pas précisée. Mademoiselle 0. DELARDZIERE (1968) a déterminé des UHG dans le Bassin du Doubs (Jura français et suisse). Elle ne définit pas explicitement l'UHG, mais on comprend que pour elle l'UHG est un bassin versant (BV) superficiel (et souter- rain?) pour lequel le bilan hydrique est équilibré: alimenta- tion = évapotranspiration + décharge à l'exutoire du BV super- ficiel, Autrement dit: dans une UHG il n'y a pas de pertes souterraines dans la direction d'un autre BV mj il y a autant de pertes souterraines que de gains. L'UHG est donc définie d'après le bilan hydrique et non pas d'après l'allure du ni- veau piézomètri que. A. GEIS5LER (1957, p.138) définit les unités oéohvdrolo- oiques; "Ein bezüglich seiner oberirdischen und unter- irdischen Wasserscheiden erforderliches und im Hinblick auf Bodenfiltrationsfä'higkeit bzw. Ver- 99 sickerung und Grundwasseraufnähme sowie Grund- wasserleitfähigkeit bekanntes Gebiet einheit- lichen geohydrologischen Charakters ist als geo- hydrologische Einheit anzusprechen". Il ajoute que la difficulté principale dans la délimita- tion d'une UHG provient du fait que le BV superficiel et le BV souterrain ne se recouvrent pas. Le nombre de facteurs dé- finissant l'UGH est plus grand que dans les définitions pré- cédentes: il faut connaître non seulement le BV superficiel et le BV souterrain, mais aussi les propriétés hydrauliques des terrains (porosité, perméabilité, transmissivité, etc} à l'intérieur des bassins. G. B. MAXEY (196-1, p. 124 et suite) définit l'unité géo- hydrologique (UGH) de la façon suivante : "An aquifer, a confining unit, dt a combination of aquifers and confining units that compose a framework for a reasonable distinct hydraulic system is a geohydrologic unit". 11 souligne que la lithologie, l'épaisseur, 1'extension et la structure de la série géologique ne suffisent pas pour délimiter les UGH: "movement, storage, and occurence of water, as well as the mode of locating, drilling, and completing wells (well logs and water level and pumpaqe records from the wells are virtually our only source of information except the sur- ficial deposits ) are also elements in the definition of the units". La définition de MAXEY est très importante car, outre les facteurs géologiques et .hydrologiques habituels, il con- sidère aussi des facteurs pragmatiques (méthodes de recher- ches et d'exploitation nécessaires) comme étant décisifs pour la définition des UGH. Les conditions géologiques sont importantes dans la mesure où elles forment un cadre pour IQD un ayatèrne d'écoulement "raisonnablement distinct", c 'est-è- dire separable des autres systèmes d'écoulement. 11 ne dit, toutefois, pas ce qu'il entend par système hydraulique "rai- sonnablement distinct" . J. TDTH (1962, 1963, 1966) va plus loin que MAXEY en ce sens qu'il analyse et définit explicitement la structure des systèmes hydrodynamiques ("flow systems"), c'est-à-dire, qu'il explicite ce qu'on doit entendre par un système hydraulique "raisonnablement distinct". Il analyse, en outre, l'influence de certains paramètres du cadre qéologique et géomorphologi- que sur la géométrie des systèmes hydrauliques, c'est-à-dire, qu'il explicite certaines relations entre le cadre géologique et les systèmes hydrauliques. Il définit le système d'écoulement souterrain par une propriété topologique des liqnes d'écoulement: "a flow system is a set of flow lines in which any two flow lines adjacent at one point of the flow region remain adjacent through the whole region; they can be intersected anywhere by an uninter- rupted surface across which flow takes place in one direction only". La région alimentaire et le région*d'exutoire d'un même système d'écoulement {"flow system") sont reliées par des lignes d'écoulement qui sont, topologiquement, des chemins homotopes. La figure 1 représente la con- ception de TOTH en deux dimensions, dans un milieu homo- gène et isotrope. Dans des communications personnelles (correspondance 1968) J, TOTH a bien voulu nous faire part de sa conception des UHG, Selon lui on peut parler de deux sortes d'UHG: l) Le premier type d'UHb est un volume de terrain où le régime des eaux souterraines ("ground water regime") est homogène. TOTH (196B) définit le ré- gime des eaux souterraines à l'aide de six facteurs: 101 - la quantité d'eau contenue dans un volume de terrain, - la direction de 1'écoulement, - le débit de 1'écoulement, - le chimisme de l'eau, - la température de 1'eau, - la variation dans le temps des facteurs précédents. Les six facteurs peuvent être reliés, do moins au ni- veau conceptuel, au cadre géologique, morphologique et cli- matique par six équations. 2 ) Le deuxième type d ' UHG est, pour TDTH, le bassin hydroqéoloqique ("qround-water basin") caractérisé par la géométrie du système d'écoulement souter- rain ("flow system"). Il s'aqit, donc, d'une unité hydraulique, d'un système hydraulique "raisonnable- ment distinct". 3.3. Commentaires et remarques Les quelques conceptions sur les UHG citées dans le paragraphe précédent montrent que le nombre de facteurs déterminant l'UHG varie d'un auteur à 1'autre, 1'importance des facteurs varie de même. D'après les conceptions récentes et plus anciennes ("classifications des nappes d'eau souter- raines" ) nous avons pu dégager trois tendances principales dans les définitions des UHG: 1. Pour le premier type d'unités on exige que le "régime" des eaux souterraines soit homogène dans 1'UHG, quels que soient, d'ailleurs, les facteurs déterminant le "régime". On reconnaît facilement que ces unités sont des classes d'équivalence con- formes au schéma abstrait du paragraphe 2.3. 102 Quand la nature des facteurs est explicitée, ces fac- teurs concernent les propriétés de L1 eau souterraine (tem- pérature, chimisme, mouvement, quantité de l'eau, etc.) et non pas celles des roches contenant l'eau. Les classes d'é- quivalence sont donc définies dans le champ des facteurs hy- drologiques et théoriquement il n'est pas nécessaire de con- naître les facteurs géologiques. Il n'est pas nécessaire que ce type d'UHG s'étende sur toute la région d'écoulement, depuis les régions alimentaires jusqu'aux régions d'exutoires. 2. La deuxième conception exige que l'UHG.forme "un tout" fermé et bien délimité de façon "naturelle". Dans ce ces 1' UHG est un système hydrodynamique souterrain, avec région alimentaire, région d'écoulement et région d'exutoire (H. SCHÛLLER, 1962, p. 158) quel que soit le nom qu'on lui donne; "bassin hydrologique" de CA5TANY; "appareil hydrau- lique souterrain" de SCHOLLERi "système hydraulique raison- nablement distinct" de HAXEY; "flow system" de TOTH; "unterir- disches Einzugsgebiet" de GIES5LER, etc. Ces termes signi- fient, tous, des classes d'équivalence dans 1'ensemble des liqnes de courants (les lignes de courant sont des abstractions: ce sont des courbes fictives relient un "point" de la région alimentaire à un "point" de la région d'exutoire). Le calcul des liqnes de courant implique: - une théorie sur 1'écoulement de 1'eau dans les milieux poreux; - la connaissance de la conductivité hydraulique du milieu; - des hypothèses sur les conditions eux limites; - la connaissance de méthodes mathématiques ou physiques pour déterminer la direction des lignes de courant en tous les points du volume considéré (modèles). Ce type d'UHG n'aura donc de sens qu'à l'intérieur d'une théorie et d'un système d'hypothèses explicitées: il se si- tue à un niveau d'abstraction plus élevé que le premier. 103 3. La troisième conception des UHG consiste à définir une "unité" du cadre géologique, morphologique et climatique et à poser que le volume de terrain ainsi délimité forme uneUHG. Cette troisième conception, apparemment très dif- férente des autres, peut, toutefois, être ramenée aux deux premières. En effet, en définissant une unité du cadre géo- logique, morphologique et climatique on admet, explicitement ou implicitement, un des deux cas suivants: - le cadre-type défini contient un système hydraulique sou- terrain "raisonnablement distinct" - a l'intérieur du cadre-type le régime des eaux souter- raines est "homoqène". La troisième conception consiste donc aussi à définir des classes d'équivalence dans le champ des facteurs hydro- logiques et hydrauliques, mais en admettant Is validité des hypothèses suivantes : - il existe des transformations transformant le champ des facteurs géologiques, morphologiques et climatiques en le champ des facteurs hydrologiques - ces transformations sont isomorphes ou homomorphes par rapport aux relations d'équivalence dans les deux champs. Par conséquent le champ des facteurs géologiques, mor- pholooioues et climatiques n'est utilisable pour la déter- mination des UHG gue si l'on peut le transformer en un champ des facteurs hydroIogigues. En définitive, l'unité hvdroaéologique est, sur le plan conceptuel, toujours une classe d'équivalence dans le champ de facteurs hydrologioues et hydrauligues. Que 1'on obtienne le valeur de ces facteurs è partir d'hypothèses, d'estima- tions statistiques, de mesures directes ou par la transforma- tion des facteurs géologiques, morphologiques et climatiques {ou mime géophysiques), cela ne change pas le fait fondamental. 104 Nous proposons, par conséquent, la définition suivante : L'unité hydrogéologique {UHG) est une classe d'équiva- lence spatialement connexe dans le champ des variables carac- térisant la qualité, la quantité ou le mouvement de l1eau souterraine. Cette définition nous parait particulièrement importante car elle exige explicitement la détermination du champ des va- riables hydrologiques, c'est-à-dire lpattribution d'une valeur bien définie (numérique ou non I ) de ces variables à chaque "point" de la région étudiée. Les méthodes employées pour la détermination des champs hydrologiques doivent être explicitées car l'adéquation à la réalité, c'est-à-dire l'utilité de ces champs "re-construits" dépendra des techniques utilisées pour obtenir les données de départ (mesures ou observations isolées et ponctuelles) et des transformations que l'on fait subir à ces données de départ (interpolation entre les points de mesure, extrapolation des valeurs pour toute la région étudiée, corréla- tion avec d'autres variables, etc. ) . Ainsi on satisfait aux exigences "pragmatiques" de HAXEY (voir paragraphe 3.2.) et l'on pourrait "classifier" les UHG non seulement selon la na- ture du champ dans lequel elles sont définies et selon la na- ¦ ture de la relation d'équivalence qui partitionne ce champ, mais aussi suivent les méthodes qui ont permis la "re-construction" de ce champ. Or, les principales différences entre ces méthodes concer- nent, le plus sDuvent, 1'utilisation des "fonctions d'interpola- tion" dans la détermination indirecte deB champs. On comprend, en effet, aisément que dans la pratique il est impassible de mesurer ou observer directement la valeur d'une variable en tous les points d'un volume de terrain et que le champ de la variable doit être déterminé, en grande partie, indirectement: 105 - par la mesure directe de la valeur du champ en des endroits isolés et généralement éloignés les uns des autres ; - par interpolation entre les point9 de mesure à l'eide d'une fonction d'interpolation. On utilise, principalement, trois types de fonctions d'interpolation: a) des théories "génétiques", "causales" ou "déterministes" qui décrivent, parfois à l'aide d'équations différentielles, le comportement d'un système en interaction avec son entou- rage. La construction de modèles simulant le comportement du système permet, qénéralement, d'obtenir le champ de certaines variables b) des fonctions d* interpolation purement "statistiques"(krigeage ponctuel, d'autres types de moyennes pondérées,"trends", etc) permettant l'estimation de la valeur du champ entre IeB points de mesure c) enfin, soit une variable H dont il faut déterminer le champ et une variable G dont on connaît déjà le champ. On peut utiliser comme fonction d'interpolation toutes les applica- tions, relations statistiques ou simples correspondances qualitatives qui permettent d'estimer la valeur (numérique ou non) de H à partir de la valeur de G (par exemple: utili- sation de la lithologie pour 1'estimation approximative des perméabilités). Par conséquent, dans la détermination des unités hydrogéo- logiques le véritable problème consiste à trouver de bonnes fonc- tions d'interpolation pour la détermination indirecte des champs hydrologiques choisis. Nous arrivons ainsi à la question fonda- mentale: quelle est l'unité hydrogéologique que nous jugeons la plus importante et quelles sont lea variables dont nous aimerions déterminer le champ en priorité ? 106 3.4, L1UHE fondamentale: le système d'écoulement A notre avis les systèmes d1écoulement définis par TQTH (1962, 1963) représentent le mieux le comportement d'une nappe dans des conditions naturelles et Sous 1'effet perturbateur de 1'intervention humaine. Ils offrent un cadre idéal [jour la syn- thèse de tous leg renseignements que l'on possède sur une nappe: distribution des potentiels, des vitesses d'écoulement, des ré- gions alimentaires et des régimes d'exutoires; distribution des températures, variation de la composition chimique de l'eau, influence mutuelle des ouvrages de captaqe, etc. Espérant que dans le cadre dea systèmes d'écoulement nous pourrions aborder les problèmes les plus importants de l'hydro- géologie des roches carbonatées, nous définissons l'unité hydro- géologique fondamentale, désignée par UHG-I: UHG-I = df. "système d'écoulement hydrodynamique" Les systèmes d'écoulement sont des classes d* équivalence dans l'ensemble des lignes de courant, donc la connaissance du champ des potentiels hydrauliques *P et du champ des vecteurs- vitesses de filtration q sera nécessaire pour leur détermination. Etant donné que la détermination du champ q est impossible sans la connaissance des perméabilités K dans l'aquifère, nous nous fixons comme but la détermination des champs K. 1P et g . Dans les chapitres suivants nous examinons la possibi- lité de déterminer ces champs dans les aquifères karstiques et, en particulier, nous examinons s'il est possible de déterminer la structure caractéristique du champ (très hétérogène) des perméabilités à l'aide des facteurs géologiques. CHAPITRE IV 107 RELATIONS ENTRE SYSTEMES D'ECOULEMENT, CARACTERES-PHYSIQUES ET FACTEUR5 GE0L0GIQUE5 DANS LE5 RDCHES KARSTIQUES Dans ce chapitre nous présentons un schéma conceptuel des relations entre facteurs géologiques, caractères physiques des roches karstiques et systèmes d'écoulement dans le karst, la connaissance de ces derniers étant le but final de notre étude. Par la suite nous examinerons, entre autres, dans quelle mesure on peut utiliser les relations ainsi structurées pour la déter- mination indirecte des champs hydrologiques et de certains champs géologiques (par exemple: champ de la fissuration). 4.1. Autorégulation partielle entre systèmes d'écoulement et caractères physiques des aquifères dans le karst Dans l'état actuel de nos connaissances et moyennant quel- ques hypothèses simplificatrices, nous pouvons admettre que le comportement de la zone saturée ("nappe d'eau souterraine") est décrit, en tous ses points et avec une approximation suffisante, par les équations générales suivantes (voir, par exemple: BEAR - ZASLAVSKY - IRKAY, 1968): - dans les nappes captives et à 1'intérieur des nappes libres: h

Ss,S Equations hydrodyn. MODELES ^J Partition par R. Systèmes d'écoulement ? q = -K grad Fig:2 L'emploi des modèles.de simulation permet une autre consta- tation intéressante: le comportement des nappes ne dépend que du champ des caractères physiques de l'aquifère (principalement des perméabilités K) et des "conditions aux limites" imposées naturellement ou artificiellement (en particulier de 1'altitude des exutoires, de l'alimentation et du débit des prélèvements) . Autrement dit, les facteurs géologiques, morphp^logicjues et climatiques exercent leur influence sur les écoulements souter- rains uniquement par, l'intermédiaire des champs de la perméabi- lité K, de la porosité m , du coefficient d'emmagasinement spéci- fiaue S et des "conditions aux limites". Si nous connaissions, —3----s---------------------------------- par exemple, la distribution de la perméabilité at de la porosité dans le sous-sol, nous n'aurions pas besoin de la géologie pour la recherche d'eau et pour la détermination des systèmes d'écou- lement . Ainsi, si nous voulons donner une signification hydrogéolo- gique ou hydraulique aux facteurs géologiques, morphologiques et climatiques, si nous voulons examiner leur influence sur l'évo- lution des systèmes d'écoulement, nous devons les "traduire" en 110 termes de "conditions aux limites" et en termes de "caractères physiques des squifères" (perméabilité, porosité). 11 semble évident que les facteurs géologiques influencent ou déterminent la perméabilité et la porosité par la distribu- tion des vides (pores interstitiels, microfissures, fractures, etc.): plus la fréquence. 1'ouverture et la connexité (exten- sion) des vides sont grandes, plus la valeur de la perméabilité et de la porosité est élevée. Si les microfissures et les frac- tures montrent une orientation préférentielle, la perméabilité devient anisotrope, c'est-à-dire la roche "conduit" l'eau sou- terraine mieux dans une direction que dans d'autres. Examinons si les facteurs géologiques sont les seuls à influencer le ré- partition des vides. Dans les roches meubles la distribution des vides dépend des fractions granulométriques, de la forme des grains, de la structure sedimentai re, etc. Dans les roches sedimentaireb con- solidées la distribution des vides est modifiée par les proces- sus diâgénétiques et, surtout, par l'apparition des fissures dues aux déformations tectoniques, mais toutes ces modifications dans la distribution des vides sont encore déterminées par des facteurs géoloqiques. La situation se complique dans le cas des roches solubles dans l'eau, comme les roches carbonatées, car la circulation des eaux souterraines peut modifier la répartition des vides (et, par conséquent, la perméabilité) en modifiant l'ouverture des fissures par la dissolution des parois (karstification) ou par la formation de dépôts. L'élargissement des fissures par dissolution dépend, bien sûr, de la composition chimique des roches carbonatées et de l'eau, mais la karstification relative de certains systèmes de fissures dépendra surtout de la grandeur et de la direction du vecteur vitesse de filtration ~q" (BEDINGER, 1966). Etant donné que selon la loi de DARCY-Q*= -K ¦ grad 0 (voir équation 1), nous arrivons dans une "feed back loop" ou "boucle de retour" très caractéristique de 1'hydrogéoloqie des roches karstiques: 111 - q dépend de la perméabilité K et du gradient hydraulique grad O i - la perméabilité K est fortement influencée par 1'ouverture des fissures karstifiées; - 1'ouverture des fissures karstifiées est fortement influen- cée par la direction et la grandeur du vecteur vitesse q pendant les états antérieurs de 1 'écoulement souterrain. Autrement dit, dans le karat, le champ de l'ouverture des vides et, par conséquent, le champ des perméabilités sont le résultat non seulement de l'histoire géologique des roches, mais de toute 1'histoire, de toute 1'évolution des systèmes d'écoule- ment souterrains : 1 ' état actuel des systèmes d'écoulement et du champ des perméabilités est le résultat d'auto-réglages succes- sifs (à court et a long terme) entre le champ des vecteurs vites- ses g. le champ des caractères physiques des aquifères (K, m , 5 ) et les"conditions aux limites" (alimentation de la nappe, alti- tude des exutoires et géométrie de la région d'écoulement) . Etant donné qu'il y a interaction entre les systèmes d'écou- lement et le champ des perméabilités par l'intermédiaire de la dissolution et de la formation de dépôts, on comprend aisément que 1'estimation, qualitative ou quantitative de la structure du champ K ne peut se faire uniquement à l'aide des facteurs géo- logiques, mais encore il est nécessaire d'estimer 1'orientation générale du gradient pour les états antérieurs de la nappe, ainsi que 1'importance relative des débits antérieurs parallèlement aux principaux groupes de fissures. En d'autres termes: il fau- drait prendre en considération les conditions "paléo-hydrauliques" comme l'ont proposé MANDEL (1966) et KIRALY - MATHEY - TRIPET, (1971). A cause du caractère transitif des influences dans la chaîne (conditions aux limites) ^ (champs q) + (perméabilités) nous pou- vons raisonnablement supposer que les corrélations entre facteurs géologiques et caractères physiques établies pour un "bassin karstique" bien défini ne seront pas forcément valables pour un 112 OtICNIfTtOM BlQMU - 8ilm - Hfititiimmit mi tiatoirtt - Mttrrts, tt ultra t -Stitiitiints, intimas ItICHIfHOl IBCJILE -fi* ci»vittilt — Dìrtttia* ,finititi d'itotttmttt — Timpènltrit, cbtmiimt -Eptiisiir it li mpn -Pirmiioiliti, untiti ITtTEMlI 0 -[COUlEMtMT IDDTfRIIAINI CIIHM 7 il * EOUAT I OHt HTDROOVNAMIflllfS DEI ECOULEMENT! [HOD(IEt DE tIMULATIONI 1 t CONDITIOMt IUK LIMIIES, * L IME NTATION. f IE L E VE Mf NTS1 POtIIION DEt EXUTOIREt CHAMFS DEt CARACTENEt PHTtIQUEt K . !,.«,. DISTIIIIUTION DEt 'VIDEI" -Orientatiti!, (/étante. ot/vtrtvn itltBtìon T----------------- » (B Q. i CT fi) O JT « a o -tu MOCEtIUt CMIMIOUEt -Oittoliitìtn. iipttt EACT EURI Gf OLO GICUEt — Slmc tort, dé fermât iòni -lithùlaait, fteitt - Gitehimit, pttigfiaçripnît FACTEURt GEOMORF1NOtOGIOUES - KtUtI Itlritviti. piatti. phittmèfitt ktrttitattl -Oritm'tÊlie» tu riniti it itimtgt stiptrficitl tt taittrrtin h---« FACTEURS UIRUTIQUEt ET RIOLDGIOUEt -Ttmpiritiw. pritipititiiitt -Sili.tipftitie», iltirttitt. ivtßttrittpi/itit» Fig 3— Schèma des relations entre facteurs hydrologiques, propriétés physiques de l'a- quifère et caractères géologiques dans le karst. 113 autre "bassin karstique" géographiquement distinct et appartenant à un autre réseau de drainage superficiel. En effet, les posi- tions géographiques des exutoires naturels et des régions ali- mentaires représentent des "conditions aux limites" et leur évo- lution dans le temps (paléogéographie, géomorphologie) pourrait influencer les écoulements souterrains, la karstifieation et la distribution des perméabilités tout aussi intensément que les facteurs géologiques {faciès, structure, etc). Toutes ces relations conceptuelles entre systèmes d'écoule- ment, caractères physiques, distribution des vides et facteurs géologiques des roches karstiques sont représentées schémetique- ment par le système partiellement autoréqulateur de la figure 3 La "feed back loop" de ce schéma fait toute la différence entre 1'hydrogéoloqie des roches karstiques et 1'hydrogéologie des roches non karstiques dû 1'influence des écoulements sur la per- méabilité est relativement faible, Il faut souliqner que l'action du champ q sur le champ des perméabilités est une action relativement "lente" et produira un effet à lonq terme seulement, introduisant ainsi la "durée" è 1'échelle géologique dans le schéma général. En simplifiant à 1'extrême les relations entre les "conditions aux limites", le champ des perméabilités et le champ q nous obtenons le schéma suivant : conditions aux limites champ g" champ K -I I I r-l = influence à court terme = influence à long terme Qn imagine facilement qu'un changement dans les conditions aux limites (alimentation de la nappe, altitude des exutoires) produira un changement très rapide dans la distribution des vec- teurs vitesse q, mais le champ des perméabilités ne sera pas im- médiatement affecté: pour que le champ K soit "réajusté" aux nou- velles conditions aux limites, il faut que ces dernières persis- tent pendant un temps "asse? lonq". IU Puisque la perméabilité occupe une position clé dans le schéma général, puisqu'elle est reliée, plus ou moins directe- ment, à toutes les autres variables, nous pouvons considérer le système partiellement autorégulateur de la figure 3 comme la représentation schématique d'une théorie sur le développement de la perméabilité et de la porosité dans les roches carbonetées. à .2. La structure générale du champ des perméabilités dans les roches karstiques Le schéma général de la figure 3 est basé sur les tra- vaux de RHODES et SINACORI {19dl), 5WINNERT0N (1919), LEGRAND et STRINGFIELD (1966), MANDEL (1966) et BEDINGER (1966) concer- nant le développement de la perméabilité et de la porosité dans le karst. Dans tous ces travaux on admet, au départ, l'existence de systèmes d'écoulement dans un aquifère de roches carbonatees, chaque système possédant une réqion alimentaire et une région d'exutoire. Les "vides", représentés par des pores et fissures interconnectés, ne sont pas encore élargis par la dissolution. A ce stade les hétérogénéités de la perméabilité sont relative- ment faibles, elles sont influencées surtout par les lithofaciès et les déformations tectoniques. Lors de l'écoulement des eaux souterraines vers les exu- toires (probablement encore "diffus"), la dissolution affectera surtout les fissures subparallèles au gradient hydraulique local. Les vitesses de filtration q et l'agressivité des eaux étant plus grandes près de la surface de la nappe qu'en profondeur, la dis- solution et la karstification seront aussi plus importantes dans la zone de battement de la nappe, ainsi qu'au voisinage des exu- toires (HEDINGER, 1966). L *hétéroqénéité et l'anisotropie de la perméabilité augmentent et l'on est en droit de supposer que les zones devenues plus perméables par la dissolution sont con- nexes entre elles et aboutissent aux exutoires de la nappe. 115 Il faut souligner que les zones plus karstifiées fonction- nent comme des exutoires locaux par rapport aux parties restées moins perméables, pouvant ainsi changer, localement, la direction d'écoulement dans les blocs à perméabilités plus faibles. La "concurrence" entre zones fortement karstifiées aboutira, for- cément , à des "captures": une zone évoluant très rapidement drai- nera les autres zones perméables, contribuant ainsi à l'unifica- tion du réseau karstique et â la concentration des exutoires de l'aquifère. Dans ce karst "évolué" on est en droit de schématiser le champ des perméabilités comme étant composé d'un réseau de drainage très perméable, mais de faible volume drainant deB blocs peu perméables, maiB de volumes très importants. L1hétérogénéité spatiale de 1'alimentation de la nappe aug- mente aussi, soit par la suite d'un drainage déjà dans la zone non saturée, soit à cause d'un drainage superficiel par IeB do- lines et les entonnoirs dans la région alimentaire. L'hétérogénéité organisée de la perméabilité et l'alimentation hétérogène de la zone saturée rendent le régime des exutoires (sources) de plus en plus "karstique", avec des crues violentes et des décrues non- exponentielles . Grâce aux travaux de BURGER (1957), SCHOELLER (1967) BERKALDFF (1967), FDRKASIEWICZ et PALOC (1967) et DROGUE (1967 et 1973) on peut considérer les courbes de décrue des sources kars- tiques comme les indications indirectes les plus importantes sur la structure du champ des perméabilités d'un aqui fere karstique. La perméabilité du réseau karstique étant de plusieurs ordres de qrandeur supérieure aux perméabilités rencontrées dans les blocs, la structure générale du champ des perméabilités d'un aquifère karstique pourrait itre caractérisée par: - la distribution statistique des valeurs de perméabilités dans les "blocs" - la densité at le degré d'organisation du réseau karstique - 1* effet d'échelle sur la perméabilité, c'est-è-dire la comparaison deB perméabilités moyennes obtenues a partir des échantillons d'ordres de grandeur différents (essais de labora- toire, essais dans les forages, valeur moyenne globale pour l'aquifère) qui nous renseigne sur le "degré de karstification". 116 d. 3. Les UHG nécessaires pour la simulation des systèmes d'écoulement La détermination dea systèmes d'écoulement (UHG-I) dans un volume de terrain donné exige la connaissance du champ des propriétés physiques et la connaissance du champ des conditions imposées (débits ou potentiels) à l'intérieur et sur les limites de la région d'écoulement. Les classes d'équivalence définies dans ces champs sont des unités hydrogéologiques importantes parce que nécessaires à le détermination des systèmes d'écoule- ment. Par la suite nous utiliserons les définitions suivantes : UHG-2 = classes d'équivalence dans le champ des perméabi- lités UHG-3 = classes d'équivalence dans le champ des porosités efficaces (nappes libres) ou des coefficients d* emmsgasinement (nappe captive) UHG-4 = classe d'équivalence dans le champ des conditions imposées. Il s'agit de la répartition spatiale des débits ou des potentiels imposés à 1'intérieur ou sur les limites de la réqion d*écoulement, c1est- à-dire des limites imperméables, des régions ali- mentaires et des régions d'exutoire. Pratiquement il est très difficile de déterminer les UHG-4 sur toutes les limites d'un volume arbitrairement découpé dans un aquifère et il nous paraît intéressant de définir un cinquième type UHG plus naturelle, équivalent du "bassin hydrogéologique": UHG-5 = volume de terrain délimité, latéralement et vers le bas, par des frontières imperméables et, vers le haut, par la surface de la nappe. Dans un tel volume toutes les régions alimentaires et toutes les régions d'exutoire se situent à la surface de la nappe, les conditions imposées sont relativement faciles à déterminer (infiltrations à la zone saturée, altitude des 117 exutoires, pompages dans les puits) et la séparation des systèmes d'écoulement se fait sans trop de difficultés. Les limites "im- perméables" latérales sont, le plus souvent, des limites hydrau- liques (lignes d'écoulement verticales), et non des limites géo- logiques (roche "imperméable"). La représentation cartographique des unités hydrogéologi- ques définies, à commencer par l'UHG-5, devrait être un des buts importants des cartes hydrogéologiques. 118 CHAPITRE V LA SIMULATION DES ECOULEMENTS DANS LE5 AQUIFERES KARSTIQUES 5.1. Introduction La technique des modèles de simulation est employée depuis fort longtemps en hydrogéologie, mais généralement pour des equi- feres i perméabilités d'interstices où ces dernières varient de façon relativement "continue" ou restent constantes. Lors de la simulation des aquifères karstiques on utilise des perméabilités ou des transmissivitéB."équivalentes" (TRIPET, 1971) qui tiennent compte de 1'effet du réseau karstique à l'é- chelle régionale, mais ce réseau, de volume très faible et de perméabilité très grande, n'est jamais introduit explicitement dans le modèle. Ces modèles (analogiques électriques ou mathéma- tiques à différences finies) permettent de simuler 1'allure gé- nérale de la surface d'une nappe karstique dans la mesure où 1'on admet que cette surface varie de façon continue, Sans brus- ques changements du gradient hydraulique. Souvent la surface de la nappe est observée (piézomètres, forages) dans les parties relativement peu perméables de l'aquifere et l'on "étalonne" le modèle (on choisit les transmissivités équivalentes) de façon à simuler ces observations pour une alimentation donnée. En bref: on emploie un champ de perméabilités équivalentes continu et "lissé" pour obtenir un champ de potentiels continu et "lissé". L'accent est mis, dans ce cas, sur la simulation des potentiels (surface approximative de la nappe, par exemple), gé- néralement dans le but de prévoir les effets à long terme de l'exploitation d'un aquifère. Ce genre de modèle convient très bien pour des régions sans exutaires ponctuels importants, c'est- à-dire pour des régions où la simulation de 1'hydrogramme des sources karstiques (vauclusiennes) n* est pas une obligation. 119 En effet, les premières contradictions dans 1'emploi des modèles apparaissent avec la simulation de 1'hydrogramme des sources karstiques. Si, par exemple, Dn utilise un champ de perméabilités équivalentes continu et lissé, il sera pour ainsi dire impossible de simuler 1'hydrogramme d'une source karstique "nerveuse", à réactions rapides, malgré le fait que ce même champ de perméabilité permet de simuler, en écoulement permanent, par exemple, la surface de la nappe d'une façon acceptable. Si, par contre, on veut simuler correctement une crue com- plète de la source karstique (et pas seulement la courbe de ta- rissement exponentielle), alorB il faut introduire explicitement un réseau karstique et des "blocs" peu perméables dans le modèle comme l'ont suggéré KIRALY et MOREL (1976a et 1976b). Mais, dans ce cas, il n'y a pas moyen de cacher ou "adoucir" la "terrible simplification" que 1'on introduit dans le modèle, c'est-è-dire il faut expliciter l'hypothèse que l'on fait sur la densité et 1'organisation du réseau karstique très perméable dans 1'aquifère. Or, ce dernier n'est jamais connu et tout ce que nous pou- vons faire, c'est de choisir un réseau "équivalent" théorique, une "structure équivalente" du champ des perméabilités qui au- rait le même effet sur 1'hydrogramme de la source que le ré- seau réel inconnu. La forme du réseau réel n 'aura, probablement, rien de commun avec la forme du réseau dans le modèle et la sur- face de là nappe simulée n'aura, probablement pas "bonne façon" (pas comme on la dessine sur les cartes hydrogéologiques): il y aura dee "bosses" dans les blocs peu perméables et des "creux" sur le réseau très perméable alors que 1'emplacement même de ces zones dans le modèle est plus ou moins illusoire. Malgré les inconvénients cités, 1'étude des "structures équivalentes" du champ des perméabilités nous paraît particuliè- rement intéressante et, laissant de côté les modèles analogiques et les modèles mathématiques a différences finieB, nous utilise- rons la technique qui permet la simulation des réseaux karsti- ques, c'est-ô-dire les modèles mathématiques à éléments finis. 120 Fig 4 1 5.2. Le principe ties modèles à éléments finis Les modèles à éléments finis sont décrits en détail par ZIENKIEWICZ (1971), JAVANDEL et WITHER5PDDN (1966), DESAI et ABEL {19 72), etc. Nous ne mentionnons ici que les principes généraux en insistant un peu plus sur quelques difficultés dans l'emploi du modèle pour le karst. Pour simplifier les problèmes de représentation, prenons comme exemple l'équation qénérale 3 décrivant les écoulements laminaires dans une nappe "bidimensionnelle" (voir pagelOS): 5 ¦¦ ¦-- 4 div (-T * g"rad h) + q = 0 ^3J S'il s'aqit d'un écoulement bidimensionnel dans le plan vertical ("profil", "coupe"), 1'équation 3 est modifiée en conséquence (les limites de la région d'écoulement restant fixes) 5—^---+ div t-K - qrad 0 > + q = 0 [4J où : h= hauteur de la nappe = potentiel hydraulique K= perméabilité Tu e * K = transmissivité 5= coefficient d'emmagasinement e= épaisseur de 1'aquifere qe densité de source (alimentation ou prélèvement Lb rénion ri'écoulement bien délimitée est subdivisée en éléments bidimensionnels, chaque élément étant caractérisé par sa forme géométrique (voir figure 4 ) et par ses propriétés physiques (transmissivité T, perméabilité K ou coefficient d ' emmagasinement S). 5ur chaque élément la hauteur h est approximée par une fonction h* que l'on choisit à l'avance comme une fonction de premier, de deuxième ou de troisième degré. Sur un élément, la 122 valeur de la fonction approximative h peut être déterminée en chaque point, de coordonnées "globales" (x.v) ou de coordonnées "IqcsIbb" (s. t) ¦ par l'interpolation des valeurs nodales K. (définies aux "noeuds" de 1'élément) à l'aide de fonctions d'in- terpolation M. qui ne dépendent que des coordonnées spatiales, généralement des coordonnées locales (s , t) : 't h (S, t) = { N; (S, t) N. Jt h. > ou simplement h = Si la correspondance entre coordonnées globales (x, y) et coordonnées locales (s, t) est définie à l'aide des mimes fonctions d'interpolation N. que l'on utilise pour la descrip- tion de la fonction h , nous parlons d'éléments isoparamétriques. Dans le cas des éléments isoparamétriques la fonction h et la forme de l'élément sont parfaitement déterminées si nous con- naissons les fonctions d'interpolation N., les valeurs nodales h. et les coordonnées qlobales x. et y. aux noeuds de 1'élément : = N. h. î et 1 yi Pour les fonctions d'interpolation de quelques éléments isoparamétriques linéaires et quadratiques , voir: ZIENKIEWICZ (1971), ERGATDUDIS et alii (196B). Pour trouver une solution approximative des équations 3 ou 4 sur chaque, élément, nous employons la méthode des résidus pondérés, plus particulièrement la méthode de 6ALERKINE (voir ZIENKIEWICZ. 1971): représentons l'équation 3 sous la forme A (h) = 0 où A est l'opérateur différentiel. Une première approximation consiste à remplacer la fonction de type inconnu h par une fonc- tion approximative h linéaire ou quadratique, sur chaque élément It M h. 1 Avec la méthode de GALERKINE Xe problème revient à trouver une fonction h qui vérifie la relation N. [K Y -M]" *-° © Après avoir choisi les fonctions d'interpolation et la géométrie des éléments noua pouvons intégrer l'équation 5 sur la réqîon d'écoulement pour obtenir un système d'équations liné- aires de la forme générale : © où h, = hauteur'de Xe nappe aux noeuds du modèle q, = débits (alimentation ou prélèvements) aux noeuds La "matricP des porosités" Dj et la "matrice des di- vergences" f A I sont définies, pour chaque élément, par : avec KbJ = —-----— = "matrice des gradients" ò* ò\ et H = tenseur des transmissivités Après assemblage des matrices élémentaires on obtient les matrices globales pour le modèle : [°]=£[°h -N-EH I3 faut souligner que le vecteur des débits q. est la somme de trois vecteurs et nous avons, pour chaque élément ; q, - QD SN.I + QL ' JLN. } + { QN, 124 où DD = débit "distribué" par unité de surface de l'élément QL = débit (entrant ou sortant) par unité de longueur de limite ?N. = débit "concentré" impose au noeud i (par exemple: pompage! |SN.\ = // jfJ. > dxdy = "surface d'influence" du noeud i {lN.I = /{n.| dl u "longueur d'influence" du noeud i situé sur la limite de la région d'écoulement Le système d'équationslinéaires 6 contient autant d'équa- tions qu'il y a de noeuds dans le modèle et pour chaque noeud nous pouvons calculer soit la hauteur h. (si l'on impose le débit q.)i soit le débit q. (si l'on impose la hauteur h.) à chaque moment de 1'écoulement. Pour ce faire, 1'équation 6 doit être intégrée encore "dans le temps", en choisissant des "paB de temps" At convenables . L ' équation 7 ci-après assure la conver- gence de la solution (JAVANDEL et WITHERSPOÛN, 19fiB): At D+-----A --< h.(t + AtJ + h 2 ' * .U)) = [d] .jhi(tij*îL-i 'q.( t+ At) 4 q. (t) 0 Cela revient à calculer d'abord une solution "moyenne" au temps (t +----) et, en partant de cette solution "moyenne", a estimer les hauteurs (ou lea débits) au temps (t + At). L*utilisation d'éléments tridimensionnels ne change pas la forme des équations 6 et 7, seuls le calcul des "matrices élémen- taires" et 1'imposition des "conditions aux limites" deviennent plus compliqués. Etant donné que la signification concrète du système d'é- quBtions 6 et des matrices I BJ , | DJ et [ AJ devrait Être claire même pour les hydrogéologues qui détestent les équations diffé- rentielles, il n'est, peut-être, pas inutile d'ajouter les re- marques suivantes: La"matrice des gradients" I Bj qui dépend uniquement de la forme de 1*élément, des coordonnées des noeuds et des fonctions 1 25 d'interpolation N. , permet de calculer, en tous les points d'un élément, les composantes du qradient hydraulique par simple mul- tiplication des valeurs nodales h. : - l°T: |m g rad h Le vecteur vitesse de filtration se calcule simplement par: "q*" -T ¦ qrTd h = -T • f BJ l • jb. L'anisotropie de la transmissivité ne pose pas de problème: on multiplie le vecteur qrad h par le tenseur - T I Il faut souligner, ensuite, qu'une "région d'influence" est attribuée à chaque noeud d'un élément et cela en fonction de la nature et de la forme de l'élément isoparamétrique. Les équa- tions linéaires du système 6 représentent, au fond, des "équa- tions de continuité" qui permettent de calculer le "bilan" pour la région d'influence de chaque noeud. Etant donné que la somme des différents termes du bilan doit être îéro, il sera toujours possible de calculer un terme inconnu par noeud si l'on connaît la valeur des autres termes. Dans l'équation 6 le "bilan" se compose de trois termes: Le premier terme se calcule à l'aide de la "matrice des porosités" I D I qui a la dimension d'une surface (m ) et _p,ui dépend uniquement du champ 5 et de la géométrie des éléments. Si l'on multiplie les hauteurs nodales h. d'un élément par la matrice ^l on obtient les volumes d'eau attribués à la ré- gion d'influence de chacun des noeuds pour cet état de la nappe. Le premier terme de l'équation fi, le produit M St donne, par conséquent, le chanqement de volume d'eau par unité de temps dans la réqion d'influence de chaque noeud, ce change- ment étant déterminé par la variation du "niveau de la nappe", Il s'agit, au fond, du débit entrant dans la réqion d'influence d'un noeud ou sortant de cette réqion quand la surface piézomé- triaue "monte" ou "descend". 126 Le deuxième terme de 1'équation 6 se calcule à l'aide de [aJ qui la "matrice des divergences" I A J qui a la dimension des trans- missivités (m/s) et qui dépend uniquement du champ T (ou K) et de la géométrie des éléments. En multipliant les hauteurs no- dales h. par la matrice I A I nous obtenons la "somme des diver- gences" du flux q dans la réqion d'influence de chaque noeud, Il s'aqit du débit qui représenta la différence entre les flux q entrant par les limites d'influence d'un noeud et les flux qui sortent par ces limites. Ce débit, positif si les "sorties" sont supérieures aux "entrées", est déterminé par le champ q uniquement. Le troisième terme est la "fonction de source" q. qui a la 3 X dimension d'un débit (m /s), Il s1agit du débit total qui entre dans la région d'influence d'un noeud ou quitte cette région et ce débit total, déterminé à la fois par la variation du "niveau de la nappe" et par la divergence du champ q, doit être égal è la somme des deux premiers termes de 1'équation 6. Dans 1'écoule- ment permanent le niveau de la nappe ne varie pas en fonction du temps et les débits q. doivent être égaux aux "divergences" du champ q: 5.3. Les conditions imposées et les unités UHG-4 Puisque dans un système d'équations linéaires an ne peut déterminer qu'une seule valeur inconnue par équation» il s'en- suit que: - è chaque noeud d'un modèle il faut imposer soit la hau- teur piézomètrique h., soit le débit nodal q,. "Ne rien imposer' équivaut a imposer q- = D! - è un noeud donné du modèle on ne peut imposer à la fois les deux valeurs h. et q.. - en imposant la hauteur h. on calcule le débit q. et inversement. 127 En vertu de ces dernières remarques nous décidons de rempla- cer le terme "conditions aux limites" par le terme "conditions imposées" chaque fois que nous parlerons des modèles d'écoule- ment. En effet, dans les modèles on doit imposer des conditions non seulement aux limites de la région d 'écoulement, mais par- tout, à touB les noeuds, même à l'intérieur du modèle. On pourrait objecter que les noeuds avec q. £ 0 ne font pas partie de l'in- térieur du modèle, mais c'est une question de terminologie. Ce qui est, par contre, important, c'est 1'existenee d'un véritable champ des conditions imposées dans la région d'écoulement simu- lée que l1on peut partitionner en UHG-4. c'est-à-dire en "unités des conditions imposées". Les principaux types d'UHG-4 que 1'on doit identifier sont : UHG-4 • 1 ^ on impose q. = Q, le noeud est sur les limites de la région d'écoulement simulée (limite "imperméable"). UHG-4 . 2 ^ on impose q. = 0, le noeud est à l'intérieur du mo- dèle (pas d'alimentation, pas de prélèvement à l'intérieur de la nappe simulée). UHG-4 • 3 ^q, ^ CU c'est-à-dire on impose un débit connu (alimentation ou prélèvement dans la nappe simulée). UHG-4 • 4 ^ on impose les hauteurs h. (ou les potentiels 0.) à l'intérieur ou sur les limites de la nappe simulée, Il faut remarquer que ces UHG sont définies dans le modèle parce que nous en avons besoin pour faire fonctionner le modèle. Mais il va sans dire que si l'on veut "mettre en modèle" une nappe réelle, les mêmes types de conditions doivent être identi- fiés dans la nappe réelle, c'est-à-dire on doit délimiter les réqionB des "imperméables"; les régions où la nappe n'est pas alimentée, ni exploitée; les régions où la nappe est alimentée ou exploitée avec des débits connus et, finalement, les régions où la hauteur de la nappe est connue avec une exactitude suffi- sante. Toutes ces réqions doivent être considérées comme des "unités hydrDgéologiques" importantes, du type UHG-4, dont la cartographie est absolument indispensable pour la mise en mo- dèle des aquifères. 126 En revenant aux problèmes des modèles, il faut insister '.e fait qu'un déi parties différentes : sur Ib fait qu'un débit imposé q. peut être composé de trois a) un "débit concentré" QN, {m /s) peut être imposé ponctuelle- ment au noeud i simulant, par exemple, les pompages au in- j ections dans un puits. Dans ce cas les "puits" devraient être placés aux sommets des éléments et non sur les côtés. b) Sur les limites de la réqion d'écoulement on peut imposer un débit QL (m /s * m) sortant de la nappe ou entrant dans la nappe par unité de longueur de limite. A chaque noeud se situant sur ces limites on attribue alors un débit QL • LN. i où LN. est la "lonqueur d'influence" du noeud i. Cette "lonqueur d'influence" LN. varie passablement suivant que le noeud se trouve sur un côté "linéaire" ou sur un côté "quadratique" plus ou moins déformé. La figure 5 illustre la variation des LN. pour les cas les plus simples, mais les valeurs indiquées pour les côtés quadratiques ne sont vala- bles que si le noeud central se trouve exactement "eu milieu", c) Enfin, un débit distribué QD (m /s . m ) peut être imposé uniformément sur toute la surface d'un élément. On attribue alors è chaque noeud de 1'élément un débit QD ¦ 5N,, c'est- à-di re un débit proportionnel à la"surface d'influence" SN. du noeud i. Ces surfaces d'influence changent beaucoup avec le type et la distorsion des éléments, souvent de manière tout è fait contraire aux idées dictées par le "bon sens". Ainsi, sur la figure 5 où l'on a représenté la répartition des surfaces d'influence SN. pour quelques éléments simples et réguliers, on découvre avec étonnement 1'attribution de surfaces d'influence négatives aux noeuds situés sur les sommets du rectangle quadratique à B noeuds. Cela signifie que si l'on admet une infiltration uniforme dans la nappe, il faut imposer des débits sortant de la nappe aux sommets des éléments quadratiques è B noeuds. Ce n'est pas "logique" pour "l'honnête homme", mais c'est cohérent dans le cadre des méthodes employées et cela donne la meilleure approximation 129 ¦£!0 — ><5-J- 0' O^ Fig:5 130 des hauteurs h pour 1'infiItretion admise, du moins en l'ab- sence d'un réseau karstique très perméable. La distinction entre les trois parties du débit total im- posé JqJ = JQN1J+ QL .JLN j + QD -JsnJ est, par conséquent, très importante dans 1'interprétation résultats obtenus par le modèle, en particulier lors de la mulation des aquifères karstiques. 5.4. Les éléments du modèle karstique et les UHG-2 Dans le modèle d'un aquifère karstique nous devons dispo- ser d'éléments qui simulent convenablement le réseau connexe, très perméable, mais de faible volume des zones fortement kars- tifiées, Si le modèle est en deux dimensions, le réseau très per- méable est simulé à l'aide de segments unidimensionnels (voir figure 4 ) introduits en "sandwich" entre les éléments bi- dimensionnels qui représentent les "blocs" peu perméables; La transmissivité attribuée à chaque segment unidimensionnel repré- sente l'intégrale des perméabilités à la verticale du segment. On admet, donc, implicitement que les zones très fissurées et les zaneB où les conduits karstiques sont particulièrement fré- quents sont plutôt subverticales, surtout si l'extension laté- rale de l'aquifère est beaucoup plus grande que son épaisseur. Dans les modèles tridimensionnels les conduits karstiques peuvent être simulés è 1'aide de segments unidimensionnels et les zones très fissurées {failles, décrochements) à l'aide de "feuilles" bidimensionnelles placées en sandwich Entre les blocs peu perméables. H s'agit des mêmes éléments que sur le fig 4 mais cette fois ils sont plongés et déformés dans l'espace tri- dimensionnel, Mous admettons que dans les éléments "karstiques" le gradient hydraulique et le vecteur q sont, en tous les points, perpendiculaires à la normale des courbes et des surfaces qui 131 représentent ces éléments. Dans ce cas nous pouvons calculer les matrices [bJ , I AJ et I d] pour chaque élément "karstique" sans trop de difficultés . Remarquons, toutefois, que par la suite nous n 'utilisons que des modèles bidimensionnels. L'introduction du réseau karstique dans le modèle partition- ne, pratiquement, le champ des perméabilités en deux types d'UHG-2, c'est-à-dire en "unités de perméabilités": UHG-2 ¦ 1 ^"blocs ou volumes peu perméables" UHG-2 ¦ 2 ^ zones karstifiées, très perméables La perméabilité varie, bien sûr, à l'intérieur de chacun de ces deux types, mais les ordres de grandeur sont tellement différents d1un type à 1'autre (les rapports de perméabilités arrivent autour de ID à ID ! ) que cette distinction reste justifiée. La véritable difficulté réside, bien entendu, dans 1'identification de ces unités sur le terrain. Les énormes différences de perméabilité entre les deux types d'unité sont, d'ailleurs, la cause principale des diffi- cultés que l'on rencontre lors de la simulation d'un aquifère karstique. Tout d'abord, il sera pour ainsi dire impossible de simu- ler correctement les très forts gradients à la limite du réseau karstique et des blocs peu perméables. En effet, la simulation correcte de ces gradients (si la nappe réelle est continue à cette limite, ce qui n'est pas assuré du tout!) exigerait la subdivision de la région d'écoulement en de si nombreux éléments que le coût des calculs serait prohibitif. L'allure de la surface piézornétrique dans les blocs peu perméables sera donc approxi- mative déjà â cause de la grandeur deB éléments employés et il serait tout à fait illusoire de vouloir "ajuster" le modèle pour simuler "exactement" les hauteurs de nappe mesuréeB dans quel- ques piézomètres ou forages. Dans des aquifères aussi hétéro- gènes, où les variations des hauteurs de la nappe Bont aussi grandes, la compréhension (même approximative) du comportement dynamique doit primer sur la reconstitution "photographique", 132 La deuxième difficulté principale concerne 1'imposition des alimentations "distribuées" (infiltration) en fonction des élé- ments utilisés pour la simulation des régions peu perméables. En effet, ces dernières sont, généralement, simulées par l'as- semblage d'un petit nombre d'éléments bidimensionnels (4 ou 5 en principe), cet assemblage étant entouré, en partie ou entièrement, par des segments a forte transmissivité simulant le réseau kars- tique. Inévitablement, un certain nombre de noeuds des éléments à faible transmissiwité se situera aussi sur le réseau karstique et un débit proportionnel à la surface d'influence de ces noeuds sera "injecté" directement dans le réseau très perméable où il pourra transiter rapidement vers l'exutDire. Dr, ce débit qui entre directement dans le réseau n'influence pour ainsi dire pas le niveau de la nappe dans les blocs et nous avons là une Bouree d'erreur de plus pour fausser les hauteurs de nappe simulées. Le pourcentage de cette alimentation "directe" dans le réseau peut être très élevé pour certains assemblages comme 1'indique la figure 6 . Cette "difficulté" dans 1'imposition des alimentations dis- tribuées n'est, au fond» pas si négative que cela, car elle nous amène à poser une question très importante: est-ce qu'il est si faux»d'admettre qu'un certain pourcentage des infiltrations ar- rive directement ou très rapidement dans le réseau karstique ? Certainement pas! Au contraire, les théories génétiques sur le développement du karst, l'étude de " 1'hydrogramme des sources karstiques et, précisément, la simulation de la réaction des sources karstiques (KlRALY et MORCL, 1976a et 1976b) semblent indiquer qu'une partie importante de l'alimentation de la zone saturée arrive directement dans le réseau très perméable, c'est- à-dire une partie importante (plus de 20 %) des infiltrations do;it Être imposée directement aux noeuds du réseau. Nous pouvons. donc utiliser le modèle à réseau karstique sans faire trop d'er- reurs è condition de pouvoir contrôler, dans une certaine mesure, le pourcentage des infiltrations arrivant directement dans le réseau. Avec un assemblsqe de 5 éléments par "bloc" peu perméa- ble c'est presque toujours passible (voir figure 6 ) 133 o o i i A ' » 6 c Fig. 6 Pourcentage de !'alimentation distribuée arrivant directement dans le réseau karstique. 13« Regardons maintenant si toutes ces difficultés et approxi- mations permettent encore de donner un sens à la notion de système d'écoulement dans un aquifère karstique. 5.5. Les systèmes d'écoulement dans les aquifères sans réseau karstique 5ur la base de le figure 1 et de la définition du chapi- tre III, la notion du système d1écoulement est facilement compréhensible, du moins dans un profil bidimensionnel, A chaque exutoire, principal ou secondaire, ponctuel ou curviligne, appartient un ou plusieurs systèmes d'écoulement, chaque système d'écoulement ayant sa propre région alimentaire bien définie. Des régions alimentaires voisines peuvent appar- tenir à des systèmes d'écoulement différents, donc, les eaux s'infiltrant en des endroits rapprochés peuvent aboutir à des exutoires très éloignés les uns des autres (voir figure 1), D'après TOTH (1962), un système d'écoulement dont la région d'exutoire et la région alimentaire sont adjacentes est un sys- tème local. 5i la région d'exutoire est séparée de la région alimentaire par un ou plusieurs systèmes locaux, c'est-à-dire les lignes d'écoulement passent au-dessous des .systèmes locaux, on parle de système intermédiaire au de système régional. La figure 1 montre bien que le sens de l'écoulement peut être différent dans le système régional et dans les systèmes locaux superposés. Gn comprend intuitivement que la température, le chimisme et le temps de séjour (iqe, renouvellement) des eaux peuvent varier beaucoup suivant la position du point de mesure dans 1'aquifère structuré en systèmes d'écoulement locaux ou régionaux. Les eaux aboutissant à un exutoire donné sont des mélanges provenant, souvent, de plusieurs systèmes d'écoulement et suivant l'importance relative des débits fournis par les systèmes la "composition" des eaux peut varier passablement. Il faut souligner que les limites entre les systèmes d'écoulement sont des limites "hydrauliques" qui passent par 135 quelques points singuliers où ni le gradient ni la vitesse ne sont définis. La position de ces limites est influencée aussi bien par les conditions imposées au système que par le champ des propriétés physiques de 1'aquifère. Etant donné que les conditions imposées, en particulier l'alimentation, varient dans le temps, la position des frontières et l'extension des systèmes d'écoulement changent aussi dans le temps, "oscillent" autour d'une position "moyenne". Les fluctuations peuvent être parti- culièrement importantes pour les systèmes locaux tandis que les limites latérales des systèmes régionaux, déterminées par les exutoires principaux et par les lignes de partage des eaux souterraines majeures, sont relativement stables (HUBBEHT, 1940). Les limites latérales des systèmes régionaux, formées d'un "rideau" de lignes de courant subverticales, peuvent être consi- dérées comme des "surfaces imperméables" et nous les utilisons pour la délimitation latérale des "bassins hydrogéologiques" (UHG-5) en l'absence de couches géologiques "étanches". Ainsi, une région où l'on connaît les exutoires principaux et les lignes de partage des eaux souterraines majeures (souvent dé- ductibles de la géomorphologie) forme un cadre idéal et "naturel" pour la mise en modèle des aquifères. L'influence de 1'allure de la surface de la nappe et l'influence de la géométrie du bassin hydrogéologique sur les systèmes d'écoulement d'un aquifère "homogène" sont décrites en détail par TDTH (1963). Parmi ses nombreuses conclusions, il faut mentionner l'augmentation de la profondeur des systèmes locaux avec la différence d'altitude entre la région alimen- taire et l'exutoire des systèmes, ainsi que la disparition du système régional avec la diminution de l'épaisseur de l'aqui- fère. L'influence de 1'hétérogénéité et de l'anisotropie de la perméabilité sur la forme et 1'extension des systèmes d'écoulement est illustrée par les modèles de FREEZE et WITHERSPOON (1966, 1967, 196B) et par ceux de KIRALY (1970, 1971). Etant donné que la structure d'un champ de perméabilités hétérogène peut être tr*s variée, il est difficile de tirer des conclusions générales: chaque bassin hydrogéologique "hétérogène" 136 Fig. 7 (a,b,e,d,e,f ) Equipotentielles (—¦__J) et systèmes d'écoulement (-..*") en nappes libres sons réseau korstique (coupes verticales), 137 est un cas particulier. Si les séries "perméables" et les séries "peu perméables" sont continues et sont disposées selon la struc- ture géologique, l'intérêt principal de la détermination des systèmes d'écoulement réside dans la représentation des communi- cations entre les parties libres et les parties captives de la zone saturée (communications entre les "nappes superposées" à travers les "séries peu perméables"). Les figures 7a, b, c, d, e, f montrent, à titre d'exempie et en coupes verticales, la disposition des systèmes d'écoule- ment dans quelques nappes libres théoriques. Les régions d'écou- lement, longues de 2 km, sont subdivisées en éléments isopara- métriques quadratiques bidimensionnels, représentés sur chaque figure. Les conditions imposées sont simples: les limites infé- rieures et latérales des aquifères sont imperméables, les SUr- faces libres reçoivent une alimentation moyenne de 3,2 *1D m /s*m (environ 1'0DD mm d'infiltration par an). Les potentiels sont imposés aux exutoires marqués sur chaque dessin. Le modèle mathématique {programme TLOW 1 du Centre d'Hydrogéologie de Neuchâtel) permet de calculer la position de la surface libre, le champ des potentiels à l'intérieur de la nappe et les débits aux exutoires pour un champ de perméabilités donné. Sur la base de ces résultats, la séparation des systèmes d'écoulement se fait assez facilement. Les figures 7 a et b montrent l'influence de l'épaisseur de la nappe (150 m et 300 m) sur les systèmes d'écoulement dans . un aquifère à perméabilité homogène et isotrope da 5*ID m/s : on remarque l'apparition d'un "faible" système régional avec 1'augmentation de 1'épaisseur de la nappe. Il est à souligner que ce système régional occupera une place de plus en plus im- portante si l'on diminue l'alimentation de la surface libre: en effet, la pente générale de la surface libre et la profondeur des systèmes locaux diminuent dans ce cas. La figure 7 c illustre 1'influence de 1'anisotropie du milieu sur les systèmes d'écoulement. Les ccnditions imposées sont les mimes que pour la figure 7 b, mais la perméabilité (homogène) est anisotrope: Kl = 9-lCT m/s; K2 = 1-10~ m/s et la perméabilité principale Kl fait un angle de 45 avec l'axe +X. Les lignes d'écoulement ne sont, généralement, pas perpendiculaires aux équipotentielles et le système régional n'existe plus, malgré la grande épaisseur de la nappe. Les figures 7 d, e, f montrent la situation dans un aqui- fère à perméabilité hétérogène: la perméabilité des séries 1 e est de 10 m/s tandis que la perméabilité des couches 2 et d est de 10 m/s. L'alimentation de la surface libre reste 3,2-10 m /s*m et l'épaisseur de la nappe est de 4 00 m en- viron. Sur la fiqure 7 d, les deux exutoires se situent dans la série 1 et la "nappe" de la série 3 devient fortement cap- tive partout où elle est recouverte de la couche peu perméa- ble 2. Sur la fiqure 7 e, la couche peu perméable 2 est tra- versée par une faille (f) ayant pour effet une forte diminutio des potentiels dans la "nappe" 3. Dans les deux cas (figures 7d, e), la limite des systèmes d'écoulement se situe dsns la région axiale de l'anticlinal. Le fiqure 7 f montre l'effet d'un troisième exutaire à la limite des séries 2 et 3: trois systèmes locaux et un système régional faiblement développé apparaissent dans la nappe, avec une diminution des pertes de charge à travers la couche peu perméable. Soulignons que dans les trois cas, l'aquifère supérieur 1 est alimenté aussi par l'aquifère inférieur 3. Il apparaît, donc, clairement que dans un aquifère sans réseau de drainage organisé (sans réseau karstique), les sys- tèmes d'écoulement peuvent être délimités de façon univoque à l'aide des modèles de simulation. 5.6. Les systèmes d'écoulement dans les aquifères karstiques Une première série de modèles illustre la situation en coupes verticales e travers des nappes karstiques théoriques à surface libre (figures 8, 9 et 10). Dans chaque modèle, les 139 limites inférieures et latérales de la région d'écoulement sont imperméables et la surface libre de la nappe est alimentée par un débit moyen de 3,2-1O- m /s m . Il s'agit, au fond, de coupes transversales dans un "synclinal" karstique large de 2 km environ, Chaque coupe est perpendiculaire eux branches principales du réseau très perméable et ces dernières apparais- sent dans le modèle comme des "exutoires souterrains" ponctuels à potentiel imposé. La perméabilité du milieu fissuré entourant les branches du réseau karstique est de 5-10* m/s. Le programme FLOW I calcule, pour chaque variante, la position de la surface libre, la distribution des potentiels dans 1'aquifère et les débits qui arrivent dans le réseau karstique (par unité de longueur). Sur les figures B a, b, c, la nappe est "drainée" par trais chenaux perpendiculaires au plan du dessin. A chaque branche du réseau appartient un système d'écoulement et les figures montrent l'influence de l'épaisseur de la nappe sur les limites de ces systèmes locaux. Il faut noter les faibles gradients (donc, les faibles vitesses) dans la partie centrale de 1'aquifère, "protégée" par les drains latéraux. Les figures 9 a, b, c montrent la distribution des poten- tiels en fonction de l'épaisseur de la nappe dans un aquifère drainé par deux chenaux parallèles. Si les potentiels sont à peu près IeB mêmes dans les deux branches du réseau, la limite des systèmes locaux est subverticale. Sur la figure 9 d, le système d'écoulement du chenal B devient moins étendu car le potentiel en B est de 5 m supérieur au potentiel imposé en A. 5i, tout en gardant la même différence de potentiel de 5 m entre les branches A et B, on diminue 1'alimentation de la surface libre, 1'extension du système local de la branche B devient presque négligeable (voir figure 9 e). En comparant les figures 9 b, 9 d et 9 e, on imagine facilement les fluc- tuations passibles de la limite des systèmes d'écoulement dans le temps et un volume de l'aqui fere qui se situe entre les deux branches du réseau pourrait appartenir alternativement à l'un au à l'autre des systèmes locaux suivant la fluctuation UD Fig. 8 (o,b,c) : Equipotentietles {/ ') et systèmes d'écoulement (*'•¦•') en nappes libres karstiques (coupes vertîcoles). hmmmmmmœ Fig. 9 (a,b,c,d,e ) : Equipotentielles (/"""N) et systèmes d'écoulement ('¦.***) en nappes libres karstiques {coupes verticales) 142 de leur limite. 5i les chenaux aboutissent à des exutoires superficiels distincts, les fluctuations des systèmes locaux se traduisent par des variations plus ou moina rapides du débit et de la composition chimique de l'eau des sources. La figure 10 a montre l'effet d'une zone verticale très perméable sur la distribution des potentiels dans une nappe libre, épaisse de 150 m environ. En remplaçant cette zone ver- ticale par trois chenaux karstiques superposés, on ne modifie pas de façon radicale la configuration des äquipotentielles (voir figure 10 b) : elles restent subverticales dans la plus grande partie de l'aquifère. Les figures 10 c et 10 d montrent 1'influence de l1épaisseur de la nappe sur le distribution des potentiels et sur la limite des systèmes locaux qui se forment autour de chaque branche du réseau karstique. Les quelques modèles que nous venons de présenter, bien que théoriques et extrêmement simplifiés, laissent déjà entre- voir certaines modifications que l'on devrait apporter à la notion da "système d'écoulement" en milieu karstique : - Les "systèmes locaux" sa développent dans la roche fissurée relativement peu perméable, autour de chaque branche du réseau karstique organisé, exactement de la même façon que les systèmes locaux de TOTK se développent autour de cha- que branche du réseau hydrographique superficiel dans les ter- rains non karstiques. Le réseau karstique organisé et à forte conductivité hydraulique joue donc le rôle des zones d'exutoires par rapport aux volumes de roche peu perméables î - Les zones d'exutoires des systèmes locaux (c* est-è-dire les branches du réseau karstique) sont souterraines I donc, dif- ficilement localisables), de géométrie très compliquée (jonc- tions ou subdivisions des segments du réseau) et les variations du potentiel y sont, souvent, très importantes (plusieurs dizaines de mètres). Autrement dit, dans un aquifère karstique, il y aura beaucoup de petits systèmes locaux, de forine et de dimension très variables, branchés sur le réseau principal or- ganisé, mais qu'il sera impossible de simuler en détail sans la connaissance exacte du réseau karstique. 143 a ^c / / 'ILlUIU ;ffiSIÏÏY\ IA V V v J^ \ UuMëi :b:::ÏÏ: ::::ffflW ttw^ «MMgMfflffi^ Fig. 10 (a,b,c,cl) : Equipotentîelles ("*~\) et systèmes d'écoulement { "**•.) en nappes libres karstiques (coupes verticales). 144 - LeB sources ou les résurgences sont les exutoires "ponctuels" et superficiels du réseau karstique organisé et elles représentent les "vrais" exutoires ou exutoires principaux de 1'aquifère karstique. C'est ici qu'apparaît la principale différence entre un aquifère karstique et un aquifère non kars- tique. Dans ce dernier, l'exutoire principal est alimenté par un système local voisin et par un système régional plus ou moins développé (voir figures 1 et 7). Dans le caB de !'aquifère kars- tique, l'exutoire est alimenté par un "chapelet" de systèmes locaux branchés sur le réseau aboutissant h la source et la réaction de la source karstique sera, déterminée par les réactions intégrées de tous les systèmes locaux (voir aussi figures 11 b, c, d, e). Une source karstique ponctuelle n'a pas de système local ou de système régional au sens habituel du terme, car aucune ligne de courant n'arrive directement à la source : partant des régions alimentaires, les lignes de courant abou- tissent d'abord au réseau karstique (dont les branches repré- sentent les exutoires "locaux") et, à l'intérieur du réseau karstique, la séparation des lignes de courant n'a, évidemment, pas de sens. La détermination des systèmes d'écoulement, au sens habituel du terme, n'est possible que là où l'on peut partitionner les lignes de courant, c'est-à-dire entre les régions alimentaires et le réseau karstique (systèmes locaux). Quand nous parlerons du "système d'écoulement d'une source karstique", nous donnerons è ce terme le sens suivant : région d'influence de la source, "Einzugsgebiet" de la source, portion du réseau karstique qui est drainée par la source, 1'ensemble des "chapelets de systèmes locaux" ou 1'ensemble des "chaînes de systèmes locaux" aboutissant à la source. Si un aquifère karstique a plusieurs exutgires principaux (sources), à chaque exutoire devrait appartenir une région d1influence, une portion du réseau karstique plus ou moins connexe, un ensemble de "chaînes de systèmes locaux". Examinons, h l'aide de quelques modèles simples, s'il est possible de déter- miner les "systèmes d'écoulements" des sources karstiques de manière univoque. 145 La figure 11 a représente la surface de la nappe dans un aquifère bidimenaionnel horizontal à perméabilité homogène et d'épaisseur constante (150 m) pour une alimentation de 3,2*10 m /s*m . La région d'influence de chacun des deux exutoires "ponctuels" est déterminée de façon univoque. La figure 11 b représente la situation dans un aquifère karstique horizontal, de perméabilité K = 5*10 m/s, d'épais- seur e = 150 m, alimenté par une infiltration moyenne de — B 3 2 3,2 "10 m /s'm et drainé par un réseau karstique connexe de maille kilométrique. Le sens des écoulements est indiqué par des flèches aussi bien dans le réseau karstique que dans les "blocs" peu perméables. L'aquifère a deux exutoires principaux et la figure 11 b montre clairement que la région d'influence des sources karstiques ne peut pas être déterminée, dans une telle situation, de façon univoque : en effet, en injectant du colorant en un point quelconque des régions marquées en pointillés, toutes les deux sources seront colorées. Cela ' provient du fait que deux des points d'intersection des bran- ches du réseau sont des points de diffluence, donc les eaux drainées par les segments du réseau situés en amont des dif- fluences peuvent aboutir indifféremment à l'un ou a l'autre des exutoires, après s'être "mélangées" dans le réseau kars- tique. La situation est encore plus compliquée en régime tran- sitoire, car les points de diffluence peuvent se déplacer en fonction de 1'alimentation et en fonction de l1évolution des potentiels dans le réseau principal. Dans le modèle de la figure 11 c, nous avons introduit un réseau karstique plus dense, de maille x = 0,5 km, tout en gardant les mêmes conditions aux limites et la même per- méabilité des "blocs". Les points de diffluence se multiplient et la région d'alimentation "commune" aux deux sources augmente. La figure 11 d montre l'effet de deux réseaux karstiques "emboîtés", de mailles et de conductivités hydrauliques dif- férentes : le réseau principal a une maille de x = 1 km et une conductivité hydraulique 100 fois plus grande que celle us Fig. 11 : Equìpotentìelles i^s) ßt systèmes d'écoulement ( *«#) dans des nappes "bidîménsionnèlles" horizontales' sans réseau karstique ( a ) et avec réseau karstique { b,c,d,e,f ) ¦ En pointillés' : systèmes karstiques ' J' diffluents. . f ¦ ¦ 'il ^j • ?\~ " j. ¦ % *¦• ¦ s.-.......-< * - .... ¦„.,-,„ I, v: '\''^n •fi, '-jj ?j'. t|- >-:t»* ji'y^' 147 \^x ¦Î4B f 149 du réseau secondaire, de maille x = 0,5 km. La répartition des régions "communes" est conditionnée surtout par le réseau prin- cipal. Dans le modèle de la figure 11 e, nous gardons les deux réseaux emboîtés mais, cette fois, le réseau principal est discontinu : chaque exutoire possède son propre réseau prin- cipal. Les régions "communes", beaucoup plus petites que pré- cédemment, sont déterminées par les points de diffluence du réseau secondaire. Enfin, dans le modèle de la figure 11 f, on diminue encore la connexité des réseaux principaux : certains segments à forte conductivité hydraulique restent isolés dans un milieu moins perméable. La surface de la nappe devient plus irrêgulière et les régions "communes" sont légèrement déplacées par rapport au modèle précédent. Conclusions Les quelques modèles présentés montrent clairement qu'une source karstique ne possède pas de "système d'écoulement" au sens habituel de ce terme : elle est 1'exutoire "ponctuel" d'un ensemble de "chaînes de systèmes locaux" (figure 11 b) , cet ensemble étant structuré par le réseau karstique qui re- présente les zones d'exutoire souterraines pour les systèmes locaux. La structure des écoulements est, par conséquent, beau- coup plus complexe dans un aquifère karstique que dans un aquifère sans réseau de drainage organisé. Toutefois, dans des aquifères plus complexes nous pouvons définir des "unités" plus complexes, les "systèmes d'écoulement karstiques" : un système d'écoulement karstique ou système karstique est un ensemble dB systèmes locaux reliés par le réseau très perméa- ble de telle façon que les eaux souterraines des systèmes lo- caux aboutiront soit à une même source unique { système kars- tique simple ou confluent ), soit â un mSme groupe de plusieurs sources ( système karstique diffluent ) . 150 Les systèmes karstiques diffluents se situeront toujours "en amont" des systèmes simples Ou confluents. La région alimentaire d'un système confluent alimente une seule source (c'est le terme "Einzugsgebiet" qui convient), tandis que la réqion alimentaire d'un système diffluent alimentera toujours plusieurs sources à la foi>s (dans ce cas seul le terme "Ab- flussgebiet" peut avoir un sens, voir aussi ZOTL, 1961) . Théo- riquement nous pouvons donc toujours partitionner la réqion alimentairE globale en "Einzugsgebiete" et en "Abflussgebiete" et nous pouvons toujours partitionner l'aquifère en systèmes karstiques confluents et en systèmes karstiques diffluents. 5i les systèmes locaux sont des classes d'équivalence dans 1'ensemble des lignes de courant, c'est-à-dire qu'ils sont des UHG-I, alors les Systèmes karstiques seront des clas- ses d'équivalence dans 1'ensemble des systèmes locaux (dans 1'ensemble des UHG-I) et formeront des "unités" ayant un "ordre de grandeur" supérieur à celui des systèmes d'écoulement hydro- dynamiques définis par J. TOTH. Il s'agit, en fait, d'une clas- sification emboîtée selon la définition du paragraphe 2.1. (pages 92 et 93) : un premier ensemble-quotient est formé par les UHG-I (c'est-à-dire par des ensembles de lignes de courant) et un second ensemble-quotient est formé par les systèmes karsti- ques (c'est-à-dire par des ensembles de UHG-I) . Cela signifie que les systèmes d'écoulement hydrodynamiques (UHG-I) et les systèmes karstiques sont définis à des niveaux différents et la présence des uns n'exclut pas la présence des autres, bien au contraire, dans un aquifère karstique. CHAPITRE VI 151 LE CHAMP DES CARACTERES PHYSIQUE5 DANS LES AQUIFERES KARSTIQUE5 6.1. Introduction Ce chapitre est consacré aux problèmes liés à la détermi- nation des caractères physiques dans les aquifères karstiques, en particulier à la détermination indirecte du champ des perméa- bilités. Ainsi qu'il a été, déjà, s.ouligné dans les chapitres précédents, on peut schématiser la distribution des perméabilités dans le karst comme le résultat de la superposition d'un champ "continu" de faibles perméabilités et d'un réseau "discontinu" très perméable. La détermination du champ des perméabilités pour- rait être divisée, par conséquent, en deux problèmes distincts qu'il s'agirait d'aborder à des échelles différentes et avec des méthodes assez dissemblables ; a) Estimation du champ des "faibles perméabilités" dont Ib par- tit ion en classes d*équivalence donnera les UHG - 2*1 et qui sera découpé, par la superposition du réseau karstique, en "blocs peu perméables". Les volumes à faibles perméabilités occupent presque la totalité de l'aquifere. b) Estimation de la structure (densité, degré d'organisation) du réseau karstique "très perméable" qui formera les UHG - 2*2 et qui sera superposé au champ "continu" des faibles perméabilités. Le réseau occupe une petite partie seulement du volume total de 1'aquifère. Pour l'estimation du champ des faibles perméabilités on peut utiliser les méthodes habituellement employées dans les aqui- fères à perméabilités d'interstices: essais de pompage pour la détermination directe et "ponctuelle" des valeurs et, si les fo- rages ne sont pas très éloignés les uns des autres, interpolation statistique entre les valeurs ponctuelles (moyennes pondérées, 152 krigeage, régressions, etc). Il faut souligner cependant qu'à l'exception de certaines régions "privilégiées" (régions minières, sites de barrages, aquifères calcaires exploités pour les grandes villes), ou bien on ne possède pas de mesures de perméabilité dans les calcaires, ou bien la densité des mesures ponctuelles est tellement faible que 1'interpolation purement statistique est â peu près inutilisable. Par conséquent nous cherchons, dans ce chapitre, â préciser les relations représentées symbolique- ment dans la partie droite du schéma qénéral de la figure 3 et, en particulier, nous examinons si les relations analysées sont, oui ou non, utilisables comme fonctions d'interpolation pour 1'estimation indirecte du champ des perméabilités. Nous ana- lysons , plus spécialement, deux groupes de relations: a) Les relations entre perméabilité et certaines variables carac- térisant la distribution des vides {principalement 1'orientation, la fréquence, 1'ouverture et l'extension latérale des fractures). b) Les relations qui montrent 1'influence des écoulements et des facteurs géologiques (lithologie, structure, etc) sur 1'orienta- tion, la fréquence et, surtout, 1'ouverture des fractures. Il faut souligner tout de suite que les relations du pre- mier groupe (a) sont des relations quantitatives (SNQW, 1969; KlRALY, 1969), donc la connaissance de certains paramètres de la fissuration permet de calculer la valeur du tenseur de perméa- bilité dans le milieu anisotrope des roches fracturées. Dans le second groupe nous trouvons des relations "géologiques" par ex- cellence et si nos collèques géologues pouvaient les transformer en relations quantitatives, ils rendraient un signalé service aux hydrogéologues du karst. Il est, par conséquent, compréhen- sible que dans le présent chapitre nous nous intéressons tout particulièrement à ces deux groupes de relations et nous laissons de côté le problème de la détermination "directe" des perméabili- tés par des essais de pompage, sujet traité dans d'excellents ouvrages spécialisés (voir, par exemple: KRU5EMAN et OE RIDDEN, 153 6.2. Distribution des valeurs mesurées. Effet d'échelle. Les publications sur la distribution des valeurs de perméa- bilité dans les calcaires sont, malheureusement, très rares et le plus souvent on ne trouve que des valeurs moyennes, sans in- dications de la dispersion des mesures et sans renseignement sur le plan d'échantillonnage. Tout d'abord, les perméabilités peuvent être déterminées en laboratoire, sur des échantillons prélevés dans les forages. Le volume intéressé par 1'essai étant très faible (de 1'ordre -4 -3 3 de 10 à 10 m ), la probabilité de trouver de grandes perméa- bilités est, évidemment, minime. Ces essais ne sont, toutefois, pas dépourvus d'intérêt car les maerofissures étant absentes, on ne mesure que la perméabilité "primaire" (perméabilité d'in- terstices, de microfissures) et la comparaison avec des valeurs obtenues sur le terrain montrera la contribution des fissures plus ou moins karstifiées à la perméabilité moyenne globale (effet d'échelle). Les valeurs de perméabilité obtenues en labo- ratoire sont qénéralement très faibles: elles varient entre 10 m/s et 10 m/s fvDir par exemple NEWBERRY, 1968; MDRRI5 et JOHNSON, 1967) , la arande majorité des valeurs étant inférieure à 10 m/s. Les essais de pompaqe ou d'injection dans les puits ou dans lés forages nous renseiqnent sur la perméabilité moyenne de vo- lumes beaucoup plus importants (10 à 10 m ) où l'effet des macro- fissures plus ou moins karstifiées doit être prépondérant. Etant donné que la qrande majorité des forages ne recoupe pas le ré- seau karstique, les résultats des essais de pompage ou d'injec- tion caractérisent surtout les blocs relativement peu perméables. BDRELLl et PAULIN (1967) présentent la distribution de plu- sieurs centaines de valeurs de perméabilité (obtenues à partir d'essais d'injection dans des forages) de façon exemplaire, à l'aide de courbes de fréquence cumulatives. 80 à 95 % de toutes les valeurs, mesurées dans trois régions du karst yougoslave, se situent entre lu m/s et 5« 1G~ m/s (les valeurs originelles sont données en unités Lugeon). La distribution des valeurs est 154 loq-normale en première approximation et la"queue de distribu- tion" est dans la direction des grandes perméabilités indiquant par là que la probabilité de trouver des zones très perméables dans les forages est faible. Dans sa thèse TRIPET (1971) donne les valeurs de perméabi- lité mesurées par essais d'injection {passes de 3 à 5 m) dans quelques forages du Jura neuchâtelois. Environ BO % de toutes les valeurs sont inférieures à 2•10 m/s et les plus grandes valeurs ne dépassent guère 5-10 m/s. Les résultats obtenus par essais de pompage confirment 1'ordre de grandeur des perméabilités. SIMEDWI (1976) montre que les valeurs diminuent jusqu'à une pro- fondeur de <150 m environ, mais ensuite la perméabilité moyenne reste statistiquement "constante" (la variance autour de la moyen- ne étant» toutefois, assez importante}. La variabilité spatiale de la perméabilité dans le sens la- téral étant aussi un caractère important du champ K, nous présen- tons brièvement l'analyse de variance par "modèle emboîté" de 399 mesures de perméabilité effectuées dans B foraqes profonds, eux-mêmes répartis dans trois synclinaux du Jura neuchâtelois. Les forages étant éloignés de plusieurs kilomètres les uns des autres, le problème était de savoir à quel niveau se situait la plus grande hétérogénéité ou la plus grande variabilité des va- leurs K: entre les valeurs moyennes des synclinaux, entre les valeurs moyennes des forages à l'intérieur des synclinaux ou entre les mesures individuelles à l'intérieur des forages ? Le tableau 1, calculé par SIMEONI (1977) donne les résultats de l'analyse de variance: - la variabilité est maximale entre les moyennes des synclinaux ( (¾ = 75,8-10"-1 ), c'est-à-dire entre les bassins hydrogéolo- giques, ce qui semble démontrer le caractère régional du dé- veloppement de la perméabilité (à l'échelle des "blocs") - la variabilité est minimale entre les moyennes des forages -2 -14 { O8 = 10,3*10 ) indiquant per là que les perméabilités moyen- nes des "blocs" varient relativement peu à 1'intérieur d'un même bassin hydroqéoloqique 155 ¦* ¦**¦ ¦* *~ •7 '7 O O O « » »" *" "~ lance timée co .ics SO CM r- O Os m W CSI £ CJ CM -j ¦ • b cu a 0 CSI O Vt o — in T- £ C csi" UÌ O) « TT n u Z * TJ -» Tf ^J- 7 O *~ _l O O a ¦**». t— Os CSI Ü » CO 9039 Os Os CSI a> X) * _J -CJ *- ,— CO CvI in O) Ol co co O = Q «J- ¦** TT -d- (0 I I I I 0) O O O O "D en ^ *- T- W Eie t— I^ Os irs *i #fc E « •* tj- T" IfS Os t— ¦>* Os O 1— ti t»- CTs CO O CSl m TJ C w O OJ ^ O .2? D 'J) Q) Q) CTi 3 CtJ en O LL QJ 5 O) Q) « — C C O 156 - Ia variance des mesures à 1 'intérieur des forages ( Oc = 29,3*10" ) est presque trois fois plus grande que la variance entre les valeurs moyennes des forages, démontrant ainsi une hétéroqénéité locale importante du champ des per- méabilités . Par conséquent, si les perméabilités moyennes des "blocs" sont relativement homoqènes à 1'intérieur des bassins, il semble qu'il y ait deux niveaux d'hétérogénéité importants dans le champ des "faibles perméabilités": a) Une hétérogénéité "locale" qui rend très difficile 1'estima- tion statistique des valeurs "ponctuelles", car cette estimation se fera toujours avec une erreur considérable. b) Une hétérogénéité "régionale" qui peut entraîner d 'importan- tes erreurs lorsqu'on veut extrapoler les perméabilités moyennes d'un bassin hydrogéologique à un autre. Ces "conclusions" ne sont, évidemment, valables que pour les synclinaux et pour les forages étudiés, mais il est souhai- table que de telles analyses soient publiées pour d'autres ré- gions karstiques afin que l'on puisse les comparer sur des bases plus objectives que par le passé. La perméabilité moyenne peut être déterminée pour des vo- 1 urnes de terrains beaucoup plus grands (de l'ordre du km ) que le volume de blocs peu perméables et à cette échelle 1'effet du réseau karstique devient généralement prépondérant. Apres avoir effectué un essai de pompaqe directement sur la source de Verdier (France), PALDC 1196d) trouve une valeur de perméabilité de l'ordre de 10~ m/s. KIJATGVIC (1970) calcule la transmissivité d'après les courbes de décrue de deux sources du Liban et trouve des valeurs de T = 1 à d , 7 m"/s correspondant à des perméabilités moyennes de l'ordre de 10 m/s. Dans ces exemples, les valeurs moyennes obtenues concernent l'ensemble de l'aquifère (réseau karstique + blocs peu perméables) dont la source est l'exutoire principal, 157 L*emploi des modèles pour la détermination de la perméabilité moyenne équivalente à 1'échelle du bassin donne des résultats particulièrement intéressants. Par exemple, TRIPET (1971), en voulant simuler les écoulements dans le bassin karstique de la source de l'A reuse (Jura neuchStelois) à 1 * aide d' un modèle élec- trique, avait dû utiliser des perméabilités équivalentes de l'ordre de IQ m/s alors que les essais de pompage dans les fo- rsqeB profonds donnaient des valeurs de 3-1U m/s à 5-10 m/s '. En comparant ces dernières valeurs avec les résultats obtenus dans les forages ou en laboratoire, 1'effet d'échelle sur la perméabilité ne fait aucun doute. Cet effet d'échelle est dû, probablement, à l'existence d'un réseau karstique connexe, de densité et de volume très faibles, mais de perméabilité très grande qui draîne des "blocs" volumineux (hectométriques ou kilo- métriques), mais relativement peu perméables. La comparaison des perméabilités moyennes obtenues a partir d'échantillons d'ordres de grandeur différents (essais de laboratoire, essais dans les forages, valeur moyenne qlobale pour l'aquifère) nous renseigne à la fois sur le "degré de karstl'f ication" . sur 1 ' effet d'échelle et sur la structure du champ des perméabilités de l'aquifère. La comparaison des valeurs moyennes devrait, toujours, être ac- compagnée d'une analyse de variance par modèle emboîté pour pou- voir décider à quel "niveau de 1'échantillonnage" 1'hétérogénéité spatiale du champ K est la plus grande: à l'intérieur des forages, entre les forages ou,' éventuellement, entre les systèmes d'écou- lement. - ' Les moyennes des valeurs empiriques citées, dans ce paragraphe sont représentées d'une façon synthétique sur la figure 12 où la distinction entre "effet' de pores", "effet.: de fissures" et "effet du réseau karstique"*h'est, bien entendu, qu'une première approximation destinée à montrer les ordres de grandeur. Les /ones en pointillé sont entièrement hypothétiques: elles repré- sentent le distribution possible des valeurs autour des moyennes. En fait, de telles "courbes'd'effet d'échelle" pourraient (ou devraient) être construites pour chaque aquifère karstique parti- culier et ces courbes pourraient servir de base à une classifi- cation hydrogéologique objective des bassins karstiques. 1 SB 3 * n a :• i --— m Ul t_ 5 5 il O J _ 3 ^ ** n ui a E pli > \ l l 1 I t I t 1 I 1 I t £ l. o t 2 i 41 s \ "(A * U) * us \ m \ \ V 1 ........Ja -Mt UÌ Ä C r*E Forages ^—,— -----1------------(------------r— —--I- « Fig. 12 Effet d'échelle sur Io perméabilité dans le karst. 159 Dans les exemples cites on admet implicitement que la per- méabilité est isotrope. Ce n'est certainement pas le cas dans la réalité, mais la détermination pratique du tenseur de perméa- bilité (9 valeurs, dont 6 indépendantes, en trois dimensions ï est très difficile et les techniques ou méthodes existantes (PAPADOPOULDS, 1967; KRUSEMAN et DE RIDDEN, 197Oj LOUIS, 1974) ne sont pas employées en routine par les hydrogéologues. Les valeurs "mesurées" du coefficient d'emmaqasinement 5 (défini, en principe, dans les nappes captives) et de la porosité efficace m (dans le cas des nappes libres) montrent une varia- tion beaucoup moins grande que les valeurs de perméabilité. TRIPET (1971) donne, pour le synclinal de la Brévine (Jura neu- châtelois) un coefficient d'emmaqasinement de 0,0035 a 0,004 ; BORELLI et PAVLIN (1967) indiquent pour la région de Busko Blato (Yougoslavie) une porosité efficace de 0,01 et PALOC (1964) trouve, pour le bassin de la source de Verdier (France), un coefficient d'emmagasinement de 0,03B. Ces valeurs sont plus faibles que les valeurs de porosité totale obtenues en laboratoire (0,03 à 0,33) sur des petits échantillons. 6.3. La nature de la perméabilité anisotrope Etant donné que dans les roches fissurées la conductivité hydraulique est qénéralement anisotrope il faut préciser, en quelques mots, la nature du tenseur de perméabilité. En remplaçant le milieu poreux discontinu par un milieu continu fictif et en remplaçant la vitesse réelle du liquide dans les pores par le vecteur vitesse de filtration fictif q, HUBBERT (1940, 1957) généralise la loi de Darcy en montrant qu'elle relie le champ des vecteurs gradients g rad1!*= J au champ des vecteurs vitesses de filtration q par l'intermédiaire d'un opérateur linéaire appelé conductivité hydraulique ou per- méabilité. SCHEIDEGGER (1954, I960) montre que dans les terrains anisotropes cet opérateur linéaire est un tenseur symétrique de 160 second ordre et les champs J et q sont reliés par 1'intermèdiaire du champ tensoriel des perméabilités. Dans ce cas les lignes d'écoulement ne sont, qénéralement, pas perpendiculaires aux lignes équipotentielles. Pour visualiser 1 * effet du tenseur de perméabilité nous pre- nons un cas bidimensionnel simple représenté sur la figure 13 ô. Sur cette figure les perméabilités principales Kl = 0,5 et K2 = 2 sont parallèles aux axes des coordonnées (x,y) et le ten- seur Ik! prend la forme d'une matrice diagonale: [¦]¦- Kl 0 0 K2 0,5 0 0 2 La figure 13a montre clairement que si l'extrémité du vec- teur gradient J décrit un cercle (de rayon unité dans notre cas), alors 1'extrémité du vecteur vitesse de filtration q = I K [¦ J parcourt une ellipse, les deux axes de l'ellipse étant parallèles aux perméabilités principales. Si le Gradient J n'est pas pa- rai lele à une des perméabilités principales, alors la direction de 1 ' écoulement q sera oblique aux équipotentielles. Il faut souliqner que les perméabilités directionnelles ne forment pas une ellipse et l'utilisation de l'expression "ellipse de perméabilités" peut induire en erreur. La définition même du terme "perméabilité directionnelle" n'est, malheureuse- ment, pes univoque car il y a, au moins, deux définitions théo- riques possibles (voir aussi SCHEI DEGGEP , 196D). En effet, la perméabilité K dans une direction n (n est un vecteur unité) n est le rapport entre la composante de q dans la direction n et la composante de J dans la direction n. Théoriquement deux cas extrêmes peuvent se présenter : a) C'est le vecteur q qui est parallèle à n et nous devons pro- jeter le gradient J=JK • q dans la direction do q. Dans ce cas le perméabilité directionnelle K est calculée par T[K] ou ---- = n-l K J - n 161 Fig. 13 Relations entre gradient ( J ) et vecteur vitesse de filtration (q } en milieu anisotrope. 162 et la figure 13a montre que les valeurs K ne forment pas une ellipse. Par contre, en reportant les valeurs */K dans chaque direction on obtient une ellipse, de demi-axes JKp et *jK2? Cette ellipse (voir figure 13 b ) correspond à l'équation 2 2 T ' |"k] T= 1 ou —— + —ï- = 1 U J Kl K2 et peut être utilisée pour la détermination graphique de la di- rection de 1'écoulement q directement à partir des équipotentiel- les (KIRALY, 1970 et 1971). En effet, tout vecteur "ä*sa,tisf aisant à 1'équation ci-dessus a son extrémité (de coordonnées q-, q ) ^x y située sur l'ellipse et en ce point (q , q ) le vecteur gradient ~~~ r "i -i —" x y J=Kj .q est perpendiculaire à la tangente de l'ellipse. La figure 13b indique la marche a suivre pour obtenir la direction de 1 ' écoulement : - on translate 1'ellipse représentative jusqu'à ce que la ligne équipotentielle soit tangente à l'ellipse - le rayon vecteur aboutissant au point de contact indique la direction (mais non la grandeur!) de q. b) Dans le second cas théorique, c'est le gradient qui est paral- lèle à n et l'on doit projeter le vecteur vitesse q = [.Kj*J dans la direction de J. La perméabilité directionnelle K se cal- cule par KJ a W ou K = n [KJ n et la figure13a montre que les valeurs K. ne forment pas une ellipse, non plus. Par contre, en reportant les valeurs 1/*lK . dans chaque direction on obtient de nouveau une ellipse, mais cette fois avec des demi-axes l/JTPet l/JTF^fvoir fi gure 13 c] Elle correspond ô l'équation M J=I ou K1*J2 + K2* J2 x y 163 et LIAK0P0ULG5 (1965Ì l'a utilisée pour déterminer graphiquement la direction d'écoulement a partir de la direction du gradient (voir figure 13c) Dans le cas où les perméabilités principales ne sont pas parallèles aux axes (x, y) du système de référence, la matrice I K I n'est plus diaconale, mais peut être calculée facilement si l'on connaît la valeur des perméabilités principales Kl, K 2 et l'angle Cl que fait la direction de Kl avec l'axe + x: W = cos a -sin a sin a cos a Kl 0 a K2 cosa sina -sina cos (L KXX KYX KXY KYY avec KXX = Kl-cos OL + K2«sin a KXY = KYX = (Kl - K2) sintt • cosa KYY = Kl-sin2a + K2-cos2a Si K1=K2 alors la perméabilité est, évidemment, isotrope. Il arrive que l'on connaisse le tenseur de perméabilité W- KXX KXY KYX KTY dans une base (x,y) quelconque (calculé à partir de la fissura- tion, par exemple) et l'on aimerait connaître la grandeur et 1'orientation des perméabilités principales Kl, K2. Dans ce cas les valeurs propres de la matrice[K]donnent la valeur des per- méabilités principales et les vecteurs propres associés indiquent l'orientation des perméabilités principales (pour le calcul des valeurs propres et des vecteurs propres voir, par exemple, LIP5CHUTZ, 196B). Dans un cas tridimensionnel les opérations sont analogues (voir, par exemple, TEICHMANN, 1964) . 164 6.4. Relations entre caractères physiques et fissuration 5ur la base des valeurs empiriques citées au paragraphe 6.2, nous admettons que la perméabilité d'interstices de la "matrice rocheuse" est négligeable par rapport à la perméabilité de fis- sures et de chenaux. Nous pouvons calculer le tenseur de perméa- bilité è partir de la fissuration (SNDW, 1969; KIRALY, 1969) si nous admettons certaines hypothèses simplificatrices qui rempla- cent le milieu fissuré réel par un milieu fissuré idéalisé: - l'eau souterraine circule dans des fissures qui sont planes et continues à l'intérieur d'un "volume élémentaire représen- tatif" de 1'aquifère - la conductivité hydraulique est isotrope dans le plan des fissures - la vitesse moyenne V^ de l'écoulement varie linéairement avec la projei fissure. la projection J du vecteur gradient J dans le plan de la WITTKE et LOUIS (196R) ont confirmé que dans une seule fissure 1'écoulement se fait dans la direction J avec la vites- P se moyenne 12 V d2 .T p © g = accélération due à la gravité (9,81 m/s } d = ouverture de la fissure (m) V = viscosité cinématique de l'eau (10 m /s) J = proj ection du gradient général J dans le plan de la fissure ( "sBf>9 dimension" ) . 5i nous avons f fissures parallèles par mètre, la section 2 -•¦ d'écoulement effective par m de surface perpendiculaire a J est f • d et le vecteur vitesse de filtration fictif sera: . _^ q q = f . d • V =-- 12V . f d3 . J P La figure 14 montre que la composante J peut être expr mée en fonction du qradient général J et de la normale n des fissures parallèles (n est un vecteur unité, de composantes n n2, n3) : J=J- (J.n)TT= [l - nBTTj . T où J I J= matrice unité nBn = produit tensoriel de la normale par elle-même. Plus explicitement : Il - nBn J = 1 0 D 0 1 0 - 0 D 1 r"l nln2 nln3 Ln3nl n3n2 En substituant J nous avons pour q: T = M-T =----- ¦ f d3 . [i -~nBn"j ¦ "J* 12 v et le tenseur de perméabilité tridimensionnel devient, pour un système de fissures parallèles: ^kJ = ------- . f.d3 . ^1 _-n*87Tj 12 V Si nous avons N systèmes de fissures, le tenseur de permea bilité global se calcule par simple sommation des H matrices IK ------------- i=l____________________ © f = fréquence moyenne du groupe i (fissures/m) d. = ouverture moyenne du groupe i (m) n. = normale moyenne du groupe i ("sans dimensions") 166 \ t"ï K Fig. 14 Schéma pour le calcul du tenseur de perméabilité en milieu fissuré. 167 Dans la plupart des cas on peut admettre que g et V sont constantes, donc le tenseur de perméabilité ne dépend que des trois paramètres n (orientations). f (fréquence) et d {ouverture) de la fissuration. Parfois il est raisonnable de supposer que l'écoulement se fait surtout dans les intersections des fissures et dans les che- naux plus ou moins karstifiés. Nous pouvons calculer le t enseur de perméabilité anisotrope pour un tel milieu (KIRALY, 1969) si nous admettons les hypothèses simplificatrices suivantes: - l'écoulement se fait dans les intersections des fissures et dans les intersections des joints de stratification avec les fissures - les intersections de deux systèmes de fissures forment un fais- ceau de tubes cylindriques et parallèles, de diamètre D - les faisceaux d'intersections sont continus dans un volume élémentaire représentatif de l'aquifère - la vitesse moyenne V de l'écoulement varie linéairement avec m t la projection J du vecteur qradient J sur le faisceau d'in- tersections . N groupes de fissures, de normales n. et de fréquences f., déterminent M = N ¦ (N-I)/2 faisceaux d'intersections et ces M faisceaux représentent un réseau tridimensionnel de "tubes" qui communiquent entre eux. L'orientation moyenne m des fais- ceaux se calcule à l'aide des normales moyennes des fissures 2A5 * Ik B Bre«ine • Gorget du Seyon l>2[nq x—~~* • ~9~ t^^ ' -**^ /^ / \ """ / °*K \ W\ A-1" t I D ^1 v A / • ,1 l * ^J ¦¦• /A to. \ .:' $A J>\ <*• 1F a; » $ 1 • / if lì • •; Vol I / . J ; J \J Yi ^ \ ^V V / /T\ y \ i- i • \ / V''^ • XJ ,s ¦ */ Fîg. 21 Orientation moyenne des principaux groupes de fissures aux stations de mesure ( Jura suisse ). 185 Les directions des contraintes principales O« \ O9 ^ CJjsont reliées à ces axes par 1'hypothèse de travail suivante: direc- tion de O = X; direction de O1 = Y; direction de O- = Z. Dans cette base, l1orientation théorique des principaux groupes de fissures est montrée par la fiqure 20 : 50 = stratification 51 et 52 = fissures de tension 53-S4 et 5 5-56 = deux paires de fissures de cisaillement conjuguées (dextre et senestre). 57-58 et 59-SlO = deux paires de fiBsures de cisaillement produites expérimentalement par BELITSKI (voir ASHGIREl, 1963; page 55). S11-512 = une paire de fissures de cisaillement et S13 = fissure de tension, développées à 1 'échelle du banc seulement. Elles sont dues au "glissement de couche sur couche" lors du plissement, L*analyse de quelque 3'5OO mesures à 70 stations réparties dans le Jura neuchStelois et en Ajoie montre que, dans la base des couches, les qroupes de fissures è grande fréquence (plus ' de deux fissures par mètre) occupent une position conforme au schéma théorique: sur le diagramme de la figure 21 on dis- tingue nettement les accumulations de pBles moyens correspondant aux groupes 51, S2, 53-54 et 55-56. Sur la figure 22 nous représentons les trajectoires moyennes hypothétiques des con- traintes O. et O- dans le Jura d'après les grands cisaillements dextres et senestres dessinés par PAVONI (1961). Le groupe Sl est effectivement parallèle à O,; le groupe 52 est parallèle à O-î les groupes 53 et 55 sont subparallèles aux grands décro- chements senestres et les groupes 54 et S6 sont subparallèles aux grands cisaillements dextres. La cohérence entre l'image théorique et les mesures sur le terrain montre que la théorie génétique proposée par PAVONI "{1961) pour le Jura (voir figure23 ) représente une excellente fonction d'interpolation pour l'orien- tation n des principaux qroupes de fissures. Cela signifie que 186 \ \ Fig. 22 Trajectoires hypothétiques des contraintes principales dans le Jura d'oprès les grands cisaillements dextres et senestres. Plan de décrochements : d'après N. PAVONI, 1961. 187 Fig. 23 Schémas génétiques théoriques des décrochements dextres et senesrres dans le Jura : a) d'après N. PAVONI (1961) b) avec des zones de cisaillement majeures délimitant la région fortement déformée. 188 l'orientation moyenne des groupes- à grande fréquence est prévi- sible dans le Jura avec une probabilité assez grande. La fréquence f La fréquence f (50) des joints de stratification est, sans doute, liée au microfaciès des séries lithologiques, mais les relations sont, pour le moment, qualitatives et hypothétiques, donc eHeB représentent des fonctions d'interpolation assez mé- diocres . En première approximation nous admettons les relations qualitatives suivantes pour le Jura neuchâtelois (communication orale du professeur B. KUBLER): - marno-calcaires : f(S0) est grande (bancs minces) - calcaires micritiques : f (SO) est petite (bancs épais) - si le rapport "élément biodétritiques/micrites" augmente, f (SQ) Bugmente aussi. Les mesures effectuées à 67 stations situées dans les cal- caires du Portlandien, du Kimmeridgîen et du Séquanien indiquent une distribution "loq-normale" des f(50), avec une valeur moyenne de 2,5 joints/m. Des levés plus détaillés seraient, toutefois, nécessaires pour connaître 1•épaisseur des bancs à 1'intérieur de chaque série. La fréquence des fissures tectoniques dépend du matériel déformé et des processus de déformation. Les relations avec la lithologie sont qualitatives: la fréquence des fissures aug- mente généralement dans la série conglomérats —¦- grès —» calcaires —»dalomies—•- mernocalcaîres . L'influence de 1"épaisseur des bancs est assez controversée. Selon PR ICE (1966) le fréquence des fissures eBt inversement pro- portionnelle à l'épaisseur des bancs, mais RAT5 et CHERNYA5H0V (1967) montrent que la distance moyenne x des fissures est une fonction non linéaire de l'épaisseur E des strates: laq x = e - loa E + b (a et b sont des constantes) m 189 Les quelques mesures dont nous disposons confirment les résultats de PATS et CHERNYA5H0V, à condition de prendre en con- sideration le "deqré de tectonisation" des affleurements- La fiqure 24 montre le plan d'échantillonnage "emboîté" du MaIm supérieur dans la région des Ponts-de-Martel (KIRALY, 1969). Les quatre groupes de fissures les plus importantes sont Sl, 52, S3-S5, 54-S6 et la figure 25 montre la fréquence de ces groupes au niveau des "échantillons" B... Les stations de mesures sont ré- parties dans trois ensembles: I {peu tectonisées): D28, 032, D33, 034, 037, 03B, 040, 041, 042 II (moyennement tectonisées): 021, 022, 023, 035, 036 III (fortement tectonisées): 024, 025, 026, 027, 029, 030, 031 Sur la figure 26a nous avons reporté, pour chaque station, le logarithme des fréquences cumulées des quatre groupes de fis- sures ( log^_^ f) en fonction de log f(50), sans avoir séparé les ensembles I, II et III. La corrélation est très faible. Sur les figures 26 b ,26 c et 26 d les ensembles I, II et III sont sépa- rés et les tendances propres à chaque ensemble se marquent plus nettement: l'influence de l'épaisseur des bancs diminue si la "tectonisation" augmente et dans les volumes fortement tectonisés la fréquence cumulée des fissures est statistiquement indépen- dante de l'épaisseur des bancs. Finalement, la relation entre fréquence des fissures et épaisseur des bancs est une "fonction d'interpolation" assez médiocre et qualitative seulement: les bancs épais sont probablement moins fissurés que les strates plus minces. L'influence de la tectonique est intuitivement claire, mais difficilement quantifiable: sur la figure 25 les fissures de tension Sl sont particulièrement développées dans le flanc -ren- versé de l'anticlinal (région B 13) et près d'un petit décroche- ment (région B 22), mais rien ne permet de prévoir quantitative- ment la fréquence. Cette fréquence varie, d'ailleurs, très rapi- dement dans l'espace comme l'indique la figure 27 où nous repré- sentons les résultats d'un levé continu sur une distance de 200 m. Le diagramme montre la variation de l'espacement moyen x 190 #*f Fig. 24 191 O" fe S m = I— 3 U < < UJ CO CC > U. = NOd Uj _ co UJ M LU MJ S O HlN AOIAI LU < LU X LU P U ^ — UJ ^j CO < Z CO > UJ LU a a Fig. 25 192 Fig. 26 Relations entre Fréquence des bancs et fissuration dans la Vallée des Ponts {Jura neuchâtelois). a,b,C,d, If1 100- 80 60- 40- 20- 10- 8 6- 4- 2- ZU WO- 80- 60- 40- 20- 10- B- 6- 4- 2- 25 30 I ?136 35.-2¾1--¾ '40 23 28^' ..^.' .28 h Ml 42 •33 «37 r= 0.389 log (Iti ) = 1.165« 0.324 log f (S0) f (S, i~ i i î r i I-IiI 2 4 6 BlO 20 40 6090IQC 3? ^0 **?/¦ / S& *41 y / • *33 •37 / / / / r=0.817 log (Ifj)=0.475 «0,664 log f(Se) 4 6 B ' K) "S 40 60 80 ' f (S0) —^- 100 193 , I fi C 100- 80-60- 36___-— 2J _—- • — 35/—- 40-20- .---------- 23__¦------22 M - • 0- r ¦- 0,982 8-6- log (Ifi)= 1,340 ? 0,2437 log (S0) 4- 2-¦ !(S0) I- ---------r~ 2 I I I I 1------------1-------|—1—I--------------^" 4 6 8 io 20 40 60 80100 100- 80- 60- 40H 20H 10- 8- 6- 2H Ii 27 30 31 '25 2& *24 . *s 26 r = -0.01 log (0^1,8346-0,0018 log (S0) r (S, -i------1—r 4 6 8 10 ~i----------1------1—r 20 40 60 80 100 194 Fig. 27 Variation de la distance mutuelle des fissures de tension SI près d'un décrochement dexrre (Creux-du-Von, Jura neuchâtelois ), 195 des fissures de tension 51 dans le voisinage d'un décrochement dextre {Creux-du-Van, Neuchâtel), ainsi que la variation de 1'écart-type des distances mutuelles. Il apparaît clairement que la zone très fissurée (f(51) /ID fissure s /m être!) est rela- tivement étroite (une centaine de m) et la transition vers les fréquences normales {f(Sl)"., 2 fissures/m) se fait très rapidement (sur 100 m environ). Cela confirme les résultats de l'analyse de variance par modèle emboîté (voir KIRALY, 1969a) dont le plan d'échantillonnage est représenté sur la figure 24 . Le problème était de savoir où se situe la plus grande variabilité: entre les stations de mesure rapprochées, entre les moyennes des groupes de stations B.. ou entre les moyennes des réqions A.. La réponse était très nette: pour les quatre groupes de fissures déjà men- tionnés (Sl, 52, 53-55 et 54-56) la variabilité maximale de la fréquence se situe entre les stations de mesure voisines (dis- tantes de quelques dizaines de mètres les unes des autres) et la fissuration peut être considérée comme homogène au niveau des ré- qions A. d*extension kilométrique. Autrement dit, il est facile de prévoir la fréquence moyenne régionale, mais il est pour ainsi dire impossible de prédire la fréquence locale. La figure 28 montre la distribution des fréquences pour les quatre systèmes Sl, 52, 53-S5 et Sd-56 d'après les résultats de 70 stations de mesure. La distribution semble être loq-normale et, sauf pour les fissures de tension 51, les fréquences inférieures à 2 fissures/m dominent. La figure 29 indique la valeur moyenne régionale et l'écart-type de 1 ' espacement des fissures pour les quatre groupes mentionnés, ainsi que la distance mutuelle des fissures de tension 51 dans le voisinage du décrochement de Creux-du-Van. Les conclusions les plus importantes que l'on peut dégager de ces diagrammes sont: - la prédominance des fissures de tension 51, 52 qui forment, avec la stratification SO,' un réseau de fissures mutuellement perpendiculaires et régionalement développé - 1'augmentation très nette de la fréquence de fissures de tension Sl dans le voisinage des décrochements ss 90 i: !" 10 h ì a tension T û ten*ion L • citatile ment 5 a cisaillement 0 — ----- / k - i AA'' 4' h • *"•'''. T * t,''.''..'* ' ¦* ,¦', ' T o* ,' ~1-' y . ' , / 03 05 densité [tissures/me Ire] Fig. 28 Distribution de ta densité des fissures dans le Jura suisse pour quatre groupes de fissures ( 70 stations de mesure ). 200- 150 - S [cm] écart type (7)=valeurs régionales T : tension (transversale) L: tension (longitudinale) S: cisaillement (senestre) D: cisaillement ( deatre ) 100 - 60 - © © 2) = valeurs locales (sur tronçons de 20 m ) du groupe T près d'un décrochement —r 50 i 100 150 AX [cm] écartement moyen 29 Distance mutuelle des fissures dans le Jura suisse. des fissures 197 - une certaine similitude entre la structure du champ des fréquences et la structure du champ des "faibles perméabilités" une grande hétérogénéité locale des valeurs à 1'échelle mé- trique ou décamétrique et une relative homogénéité des valeurs moyennes des renions d'extension kilométrique - parmi tous les "facteurs" géologiques énumérés seule la struc- ture (décrochement, faille, flanc renversé) permet de prévoir les régions où les fissures (en particulier le groupe 51) ont une grande fréquence. En ce qui concerne la variation de la fréquence des fissures avec le profondeur, nous admettons une diminution de UQ % à 60 % pour tous les systèmes (voir KRAUSE, 1966; JAMIER, 1975), L'ouverture d D'après RATS et CHERNYA5H0V (1967) l'ouverture moyenne d m des fissures dans les qrès augmente avec 1'épaisseur E des bancs selon une équation du type log d a • log E + b a,b = constantes) Etant donné que les mêmes auteurs proposent une relation semblable pour la distance moyenne x des fissures: log x = a ' .log E + b' (a',b' = constantes) il est tout naturel de combiner les deux équations pour exprimer 1'ouverture moyenne d des fissures en fonction de leur distance J m moyenne x : loqd = a"•loq x +b" (a",b" = constantes) Cette dernière relation implique que l'ouverture moyenne des fissures augmente avec leur distance moyenne, c'est-à-dire dans un banc très fissuré 1'ouverture moyenne des fissures sera très faible. A notre avis, cette manière de voir les choses 196 n'est, peut-être, pas la meilleure car elle ne rend pas compte de 1'existenee des réseaux de fissures emboîtés, de maille x . et d'ouverture d. (voir paragraphe 6.d.). Pour 1*hydrogéologue il serait beaucoup plus important de connaître la distance mu- tuelle x. entre les fissures de "même ouverture" d. dans diffé- î î rentes situations géologiques (lithologie, épaisseur des bancs, structure, etc). BOCKER (1973) établit une telle relation statistique basée sur 1'étude des microfissures ayant des ouvertures inférieures à 0,G3 mm: loq x. = 0,92 ¦ loq d. + 3,1 En admettant la validité de cette relation pour les macro- fissures on trouve, par exemple, une distance moyenne de 20 m environ entre les ouvertures centimétriques et une distance moyenne de 1,26 km entre les ouvertures de 1 m. La formule ne dit rien sur l'extension des ouvertures, mais malgré cela, son importance paraît évidente et de telles relations statistiques devraient être établies pour plusieurs régions du Jura. Dans les roches carbonatées c'est certainement la dissolution qui influence le plus fortement 1'ouverture des fissures ou le diamètre des intersections et la meilleure "fonction d'interpo- lation" serait une théorie qénétique qui prévoit la distribu- tion spatiale de l'intensité de la. karstification. Il existe beau- coup de théories sur la karstification, mais la plupart d'entre elles se limitent à 1'explication de phénomènes particuliers déjà localisés et elles ont peu de valeur quand il s'agit de prévoir les endroits où 1'élargissement des fissures par la dissolution est particulièrement important. Notre préférence va aux idées développées dans le chapitre IV et nous admettons avec MAWDEL (1967) et BEDINGER (1966) que, en première approximation, 1'élargisse- ment des fissures par dissolution est proportionnel au débit qui s'écoule parallèlement à ces fissures, pour autant que le direc- tion qénérale des écoulements reste la même pendant assez long- temps. BEDINGER (1966) examine les conséquences d'une telle hypo- thèse à l'aide de modèles électriques et il arrive aux conclu- sions suivantes : 199 a) La zone de dissolution la plus active se situe à faible profon- deur sous la surface de la nappe (les vitesses de filtration q sont généralement grandes). b) Les chenaux de dissolution augmentent en dimension (diamètre) et diminuent en nombre depuis la région alimentaire vers l'exu- toire du système d'écoulement. c) L'orientation des chenaux de dissolution est déterminée par la direction réqionale du mouvement de l'eau souterraine. Le modèle de Bedinqer simule, en fait, l'évolution du sys- tème autoréqulateur de la figure 3 et les conclusions a, b, c décrivent le résultat des auto-réglages successifs entre le champ des vecteurs vitesses q et le champ des perméabilités K dans des conditions "géoloqiques", "morphologiques" et "climatiques" très simples : - les perméabilités sont homogènes et isotropes au départ de la simulation, donc on admet une "lithologie" uniforme et une "fissuration" homogène pour l'aquifère initial - il n'y a qu'un seul système d'écoulement, donc l'aquifère ne possède qu'un seul exutoire - 1'altitude de 1'exutoire, le débit d'alimentation et l'exten- sion de la région alimentaire restent constants dans le temps, donc il n'y a pas de changements dans les conditions "morpholo- giques" et "climatiques" pendant 1'évolution du système d'écou- lement - l'élargissement des fissures ou des chenaux est proportionnel au débit, donc "1'agressivité" des eaux souterraines est homo- gène dens l'aquifère. Si l'on peut admettre, en première approximation, des con- ditions aussi simples pour un aqûifère réel, les résultats de Bedinqer peuvent être utilisés comme des fonctions d'interpola- tion qualitatives et asse? qrossières, mais somme toute "utiles", pour la répartition spatiale des ouvertures (et des perméabilités): 2OD a) L'ouverture des vides devrait être plus grande près de la sur- face de la nappe que dans les zones profondes de l'aquifère. Si- gnalons que d 'autres considérations théoriques renforcent cette hypothèse, notamment 1'agressivité théoriquement plus élevée des eaux souterraines près de la surface de la nappe et une certaine augmentation de l'ouverture des fissures (même en l'absence de la dissolution) dans la zone de "décompression" des roches, par la suite dE la diminution des contraintes effectives vers la sur- face topDgraphique. La diminution rapide de la perméabilité avec la profondeur dans les "blocs peu perméables" (BORELLl et PAVLIN, 1967; LOUIS, 1974; 5IHEDNI1 197fi) semble confirmer ces hypothèses: b) La deuxième conclusion de Redinqer, c'est-à-dire la diminution du nombre et l'augmentation de l'ouverture des chenaux de dissolu- tion vers la région d'exutoire, doit it re appliquée à deux échel- les différentes, pour deux réseaux emboîtés: - à l'échelle du bassin hydroqéoloqique pour le réseau karstique très perméable (ouvertures décimétriqués ou métriques) avec la source karstique comme exutoire - à l'échelle des blocs à faibles perméabilités où l'exutoire direct des chenaux de dissolution (ouvertures millimétriques ou cent!métriques) n'est pas la source, mais le réseau karstique principal qui les draine. c) La troisième conclusion de Bedinger doit être nuancée. 5'ïl est vrai que 1'orientation générale des branches principales du réseau karstique est déterminée par la forme du bassin et par la position de la source (donc par la direction moyenne du gradient hydraulique), il est tout aussi vrai que, localement, l'orienta- tion des gradients et des lignes d'écoulement est profondément modifiée autour de chaque zone très perméable (voir chapitre V). Cela signifie que l'orientation des branches secondaires du ré- seau karstique et 1'orientation des chenaux de dissolution dens les blocs peu perméables pourraient être complètement différentes de la "direction régionale" du mouvement de l'eau souterraine. Autrement dit, tant que la théorie qénétique ne permet pas de prévoir 1'emplacement et 1'orientation des segments très per- 201 méables du réseau karstique, elle sera inutilisable comme "fonction d'interpolation" pour 1'orientation des chenaux de dissolution (et pour l'anisotropie de la perméabilité) dans les blocs peu karstifiés. Dans le cas des aquifères "réels" 1'évolution du sys- tème karstique autorégulateur est, généralement, plus complexe et seule une théorie génétique, basée sur 1 ' étude géologique, géo- morphologique et climatique de la région, permet de comprendre les changements des "conditions aux limites" (vitesses et durées) et d'évaluer, qualitativement et grossièrement, 1'effet de ces changements sur les systèmes d'écoulement et sur l'élargissement des fissures et des chenaux par dissolution. Résumé Le problème de la détermination indirecte des paramètres de la fissuration à l'aide des facteurs géologiques peut être résumé de la façon suivante: orientation n: Une théorie génétique vérifiée par des études sta- tistiques permet de prévoir, avec une exactitude suffisante, l'o- rientation des principaux groupes de fissures dans le Jura à par- tir de la structure géologique cartographiée (plis, décrochements, orientation des couches, etc.). fréquence f: Seule l'étude statistique permet de prévoir quanti- tativement les fréquences moyennes régionales (échelle kilométri- que) des principaux groupes de fissures dans le Jura. La prévi- sion quantitative des valeurs ponctuelles est empêchée par la qrande hétéroqénéité locale du champ des fréquences. Les facteurs qéologiques sont utilisables pour une estimation qualitative des fréquences (localement aussi!): augmentation de la fréquence dans la zone des décrochements; diminution de la fréquence avec l'aug- mentation de l'épaisseur des bancs dans les régions peu tectonisées Il arrive que la cartographie de la fissuration permette de mettre en évidence des zones particulièrement fracturées même la où le carte géologigue n'indique ni décrochement, ni faille, ni déformations importantes des strates. La figure 30 202 Nombre do fractures par 10m. 566.000 Fig. 30 Fissuration dans la région de la Grotte de Mîlandre (Ajoie, Jura tabulaire). 203 représente les résultats du lever de fissuration à 10 stations de mesure dans le jura tabulaire, en Ajoie (région de la Grotte de Milandre) d'après KIRALY - MATHEY - TRIPET (1971). La repré- sentation cartographique de la fréquence réelle des fissures fait apparaître une zone fracturée d'orientation NE - SW où la "prédiction" des grandes fréquences est relativement facile pour tous les points situés entre les stations de mesure 7, 5, A et 9. Il faut remarquer qu'aucun des groupes de fissureB représentés n'est strictement parallèle à la zone fracturée. ouverture d; Seuls les modèles statistiques permettent de prévoir la "densité moyenne probable" des différentes ouvertures dans une région (BOEKER. 1973). Les théories génétiques, basées sur 1'interaction des facteurs géologiques, morphologiques et climatiques avec les systèmes d'écoulement et les processus de dissolution, donnent des prévisions qualitatives sur la dis- tribution des chenaux karstifiés. Dans ces prévisions, la kars- tification est une "fonction" du lieu dans le système d'écou- lement général (BEDINGER, 1966). Il y a une très nette corré- lation entre 1'orientation des branches principales du réseau karstique et 1'orientation des intersections fissures-strates subparalièies au gradient hydraulique général (KIRALY. 1968 i KIRALY - MATHEY - TRIPET, 1971). 6.6. Les problèmes liés à la détermination indirecte du réseau très perméable La détermination directe de la densité et de 1'organisation du réseau karstique dans un Bquifère calcaire est impossible dans l'état actuel de nos connaissances. Certaines méthodes indirectes (géologiques, géomorpholo- giques, géophysiques) permettent parfais de localiser, sur une courte distance et è faible profondeur, quelques segments du réseau karstique, mais les renseignements ainsi obtenus sont beaucoup trop fragmentaires pour pouvoir les utiliser dans les modèles de simulation. 204 L'impossibilité de déterminer, directement ou indirecte- ment, la densité et le structure réelles du réseau très per- méable dans un aquifère karstique nous a amenés à proposer 1'u- tilisation des modèles de simulation pour la détermination indirecte d'un "réseau équivalent" (KIRALY, 1975; KIRALY - MOREL, 1976 a et 1976 b). Cette méthode donne des résultats satisfaisants dans un aquifère bien délimité par des séries peu perméables et ne pos- sédant qu'un seul exutoire principal. Dn doit connaître la perméabilité moyenne dans les blocs peu perméables (mesures dans les forages ou estimation indirecte "raisonnable"), ainsi que 1'hydrogramme de l'exutoire et le niveau de la nappe en quelques endroits de 1'aquifère. Avec les données disponibles, on construit un modèle mathématique de l'aquifère et l'on y introduit un réseau hypothétique très perméable, de maille plus ou moins régulière, la conductivité hydraulique variant entre 10 m/s et 10 m/s. En utiliBant un programme de calcul pour écoulement transitoire, on simule, dans la mesure du possible, 1'hydrogramme de Is source karstique pour une alimentation isolée (orage) de la nappe. La comparaison entre 1'hydrogramme réel de la source et 1'hydrogramme simulé (rapidité de la crue, débit maximal, partie non-exponentielle de la décrue et coefficient de taris- sement de la partie exponentielle) permet de corriger 1'hy- pothèse initiale sur la densité et la conductivité hydraulique du réseau. Dans la plupart des cas, on doit admettre qu'un certain pourcentage des infiltrations arrive de façon déjà "concentrée" à la zone saturée. En employant cette méthode (KIRALY - MOREL-, 1976 a), nous avons trouvé une "maille" variant de 0,8 à 1,5 km pour le réseau karstique "équivalent" de l'aquifère calcaire de la source de l'Areuse (Jura neuchâtelois). L'analyse de l'hydrogramme des sources karstiques et 1'utilisation des modèles de simulation sont apparemment les 205 seules méthodes indirectes qui permettent de contrôler, dans une certaine mesure, les hypothèses que l'on fait sur la structure du champ des propriétés physiques de 1'aquifère et sur la structure des "systèmes d'écoulement karstiques". 206 CHAPITRE VII CONCLUSIONS Etant donné que nous connaissons peu de chose sur 1 ' empla- cement et la forme exacte du réseau karstique réel, il serait tentant de conclure que "dans l'état actuel de nos connaissances il est impossible de déterminer les systèmes d'écoulement dans le Jura" et de nous consoler en disant que, finalement, 1 * impor- tant n'est pas tellement d'arriver au but, mais d'avoir fait un beau voyage à travers les problèmes hydroqéologiques du karst. Pourtant, il reste de ce "voyage" un enseignement précieux: si les restrictions imposées par les problèmes de calcul, par la connaissance approximative des conditions aux limites et par 1'ignorance de la géométrie réelle du réseau karstique ne permet- tent pas la détermination "exacte" des systèmes d'Écoulement, ces restrictions ne nous empêchent pas pour autant de mettre en mo- dèle un aquifère avec les données que nous possédons et de déter- miner les potentiels ou les débits au niveau et avec 1'exactitu- de de nos connaissances. Or, la mise en modèle d'une nappe exige, selon l'expression de A. Bürger, la synthèse de tous les rensei- gnements, empiriques au théoriques, que nous avons sur 1 'aquifère. Quel que soit le niveau ou "l'exactitude" de ces renseignements, ils doivent s'ordonner d'une façon cohérente dans le modèle de simulation et tester la cohérence des renseignements ou des hy- pothèses sur le champ des caractères physiques et sur les condi- tions aux limites est un des buts de 1'emploi des modèles. L'utilisation des modèles est possible dès que l'on formule les premières hypothèses sur la géométrie, sur le champ des perméa- bilités et sur les conditions aux limites de 1'aquifère. Due ces hypothèses reposent sur des théories génétiques expliquant la kBrstification, sur 1'intuition et 1'expérience d'un hydrogéolo- gue ou sur une carte des perméabilités obtenue par krigeage à partir d'un réseau de forage, le modèle donnera une "réponse" 207 qui représente les conséquences véri fiables des hypothèses uti- lisées. Finalement c'est la vérification de la réponse du mo- dèle, c'est-à-dire la comparaison de la réponse avec le comporte- ment réel de la nappe qui permettra de décider si les hypothèses sont acceptables ou non. Cela signifie que les modèles peuvent être employés pour la détermination approximative des systèmes d'écoulement déjà avant que tout soit prêt pour la "grande simu- lation finale et infaillible" et les premiers résultats hypothé- tiques pourraient suggérer la marche à suivre pour obtenir des renseignements plus précis sur le champ des caractères physiques, sur les conditions aux limites et sur la géométrie de l'aquifère. Ainsi, la détermination des systèmes d'écoulement et des autres unités hydroqéologiques définies (UHG-2, UHG-3, UHG-â, UHG-5) devrait être considérée comme un processus dans 1'acqui- sition des connaissances de plus en plus précises, de plus en plus "exactes" sur l'aquifere. Dans ce processus toutes les me- thod ea d'estimation décrites, directes ou indirectes, sont néces- saires car elles entrent en jeu, à un moment ou à un autre, au fur et è mesure que notre connaissance d'une région progresse. Le point de départ pour la détermination des unités hydro- qéoloqiques dans une région est fourni par les théories généti- ques (chapitre IV} et par la connaissance des facteurs géologi- ques» morphologiques et climatiques. A l'aide des théories gé- nétiques et des relations qualitatives entre facteurs géologi- ques et distribution des vides (paragraphe 6.5. et schéma géné- ral de la figure 3 ) on transforme 1'image géologique (morpho- logique, climatique} de la région en une répartition hypothéti- que (mais plus ou moins probable) et.qualitative des perméabili- tés, des porosités, des régions alimentaires et des régions d'exutoire. Les unités hydrogéalogiques les plus importantes sont esquissées et déjà è ce niveau de la connaissance on peut songer à estimer l'ordre de grandeur des alimentations, des per- méabilités ou des porosités efficaces et des débits aux exutoires en vue de tester la cohérence de toutes ces hypothèses à l'aide de modèles. 208 TrèB souvent Ib manque d'observations pertinentes sur le comportement réel de la nsppe {niveau piézométrique dans les calcaires, régime des exutoires) limite singulièrement 1'effi- cacité de la vérification des résultats du modèle car on ne - peut pas faire un choix strictement univoque parmi toutes les hypothèses compatibles avec les données observées. Malgré cette indétermination, le résultat général doit être ccnaidéré com- me positif. En effet, è ce stade "initial" la réqion est struc- turée par les cinq types d 'unites hydrogéologiques "fondamen- tales", hypothétiques certes, mais qui sont organisées en un système interactif cohérent par l'intermédiaire du modèle d'é- coulement et qui forment un cadre général où l'on peut intégrer chaque information complémentaire, chaque nouvelle donnée, cha- que nouvelle relation entre les variables et chaque amélioration des techniques de simulation. Avec 1*euqmentation progressive du nombre des données disponibles (niveau de la nappe; valeurs de perméabilité, de porosité, de débits; régime thermique et chimique des eaux souterraines) qui doivent être cohérentes entre elles, l'incertitude concernant la valeur et l'extension réelles ("exactes") des unités hydrogéologiques diminuera aussi- En définitive, dans une région donnée et à un moment donné, il peut y avoir plusieurs unités hydrogéologiques du mSme type (par exemple des "unités de perméabilité" UHG-2), juxtaposées ou super- posées, mais qui ont été déterminées è des niveaux différents, avec des méthodes différentes, avec des exactitudes différentes. L'important est de ne pas les confondre, de ne pas leur donner la même signification et de ne pas les utiliser pour la résolu- tion de problèmes qui exigeraient un autre niveau de précision. En guise d'exemple concret nous mentionnons la carte hy- drogéologique du Canton de Neuchatel (KIRALY, 1973) où nous re- présentons la géométrie des principaux aquifères du Jura neu- châteloiB, leurs régions.alimentaires, leurs exutoires, Is dis- tribution des perméabilités dans les calcaires fissurés, ainsi que l'effet d'échelle dû è la présence du réseau karstique. Cer- tains des aquifères représentés sur la carte hydrogéologique sont "mis en modèle"actuellement (KIRALY-MOREL, 1976). 209 REMERCIEMENTS Mes remerciements les plus aineères vont à Monsieur le professeur A. BURGER, directeur du Centre d'Hydrogeologie, pour m'avoir donné la possibilité de réaliser ce travail et pour 1'intérêt bienveillant qu'il témoignait à 1'égard de mes recherchas, ainsi qu'à Monsieur le professeur J.-P. SCHAER pour ses encouragements constants et patients. L'appui financier du Fonds national suisse de la recher- che scientifique était indispensable pour mener à bien les travaux de recherche. Due cette institution trouve ici 1'ex- pression de ma gratitude.. Monsieur le professeur BANDERET, directeur du Centre de calcul de Neuchâtel, m'a considérablement facilité 1'accès à l'ordinateur et le travail avec les modèles mathématiques ; qu'il en soit vivement remercié. Je tiens à remercier également, Madame J.-G. GREMAUD pour la dactylographie du texte, ainsi que M,-A. BETf)IX pour le des- sin de certaines fiqures et pour se contribution à la présenta- tion définitive de ce travail. Enfin, j'aimerais exprimer mes remerciements les plus sin- cères à tous mes collègues du Centre d'Hydrogéologie, en parti- culier à J.-P. TPIPET et è J,-P. SIMEGNI, pour les discussions tonifiantes qui étaient souvent, à l'origine des idées dévelop- pées dans les chapitres précédents. 210 BIBLIOGRAPHIE ASHGIFlEI, G.D. (1963) - 5trukturgeologie. VEB Deutscher Verlag der Wissenschaften, Berlin 1963, 572 s., 369 abb., 22 taf. AURDUZE, J. (1966) - Les conditions d'existence des nappes aquifères et la notion de piège aquifère. Sciences de la Terre, TXl (1966), No. 1, p. 19-40, 6 fig., I tabi., Nancy. BEAR, J-, ZASLAWSKY, D., IRMAY, S. (196B) - Physical principles of water percolation and seepage. UNESCD; 465 p. BEDINGER, M.S. (1966) - Electric analoq study of cave formation. Nat. speleol. Soc. Bull. Vol. 28, no 3, pp. 127- 132, 4 fiq. BERKALOFF, E. (1967) - Limite de validité des formules courantes de ta ri ss ein ent de débit. Chronique d 'Hydrogéologie, 10 : 31-41. SORELLI, M., PAVLIN1 H. (1967) - Approach to the problem of under- ground water lakage from the storages in karst regions. AIH pubi. 73 : 120-13B. BOCKER, T. (1973) - Theorical model for karstic rock. "Karst és Barlangkutatas" 7 : 93-104. BRAIT5CH, 0. (1956) - Quantitative Auswertung einfacher Gefü'ge- diagramme. Heidelberger Heìtr. Minerai. Petrogr. 5 : 210-226. BURGER, A. (1956) - Interprétation mathématique de la courbe de décroissance du débit de l'A reuse, Jura neuehâtelois (Suisse). Bull. Soc. neuch. Sc. nat. 79 : 49-54. - (1959) - Hydrogeologie du bassin de l'Areuse. Bull. de la Soc. neuch. de Géographie. Nouvelle série, no. Il/Tome LII, fase. 1-, 1959. pp. 5-304, 29 fig., 6 pi. CASTANY, G. (1963)— Traité pratique des eaux souterraines. Dunod, Paris, 1963, 658 p. - (196B) - Prospection et exploitation des eaux sou- terraines. Dunod, Paris, 717 p. 211 CLOOS, E. (1955} - Experimental analysis of fracture patterns. Bull. Geol. Soc. Amer. 66/3, 241-256, B pi., 11 fig., 21 réf. DELAR05IERE, D. (1960) - Contribution à l'étude du Bassin du Doubs. Thèse, Paris, Faculté des sciences. DESAI, 5.C., ABEL, J,F. (1972) - Introduction to the finite element method. Van Nostrand Reinhold Co., London, 477 p. DROGUE, C. (1967) - Essai de détermination des composantes de l'écoulement des sources karstiques. Chronique d'Hydrogéologie 10 : 42-47. - (1972) - Analyse statistique des hydrogrammes de décrue des sources karstiques. J, of Hydrology, 15 : 49-6B. ERGAT0UDI5, B. et alii. (1968) - Curved, isoparametric, "quadri- lateral" elements for finite element analysis. Int. J. Solids Structures 4 : 31-42. F0RKA5IEWICZ, J., PALOC, H. (1967) - Le régime de tarissement de la Foux de la Vis. Chronique d'Hydrogeologie 10 : 59-73. FREEZE, A.R., WITHERSPODN. P.A. (1966-1967) - Theorical analysis of Regional groundwater flow. Water Resources Research (1) ; 2/d, pp. 641-656, 9 fig.: (2) : 3/2, pp. 623-634, 6 fig. GIESSLER, A. (1957) - Das unterirdische Wasser. VEB Deutscher Verlag der Wissenschaften Berlin, 1957, 1B7 p., 82 fig. GILES, R.V. (1962) - Fluid mechanics and hydraulics. Schaum Pubi. Co., New York, 274 p. HITCHON, B. (1969) - Fluid Flow in the Western Canada Sedimentary Basin. Effect of Topography. Water Resources Research 5/1 : 186-195. HUBBERT, M.K. (1940) - Theory of groundwater motion. Jour. Geol. 48, pp. 785-944. - (1967) - Darcy's law and the field equations of the flow of underground fluids. Bull. AIHS, 1957/No. 5, pp. 24-59, 11 fig. 212 I.S.R.M. (1975) - Suggested methods for the description of rock masses, joints and discontinuities. Second draft of ISRH Working Party, Oslo. JAMIER, D. (1975) - Etude de la fissuration, de l*hydrogéologie et de la géochimie des eaux profondes des massifs de l'Arpille et du Mont Blanc. Thèse, Neuehatal, 153 p. JAVANDEL, I., WITHER5P00N, P.A. (1968) - Application de la méthode des éléments finis aux écoulements transitoires en milieu poreux. Revue I.F.P. Vol. XXIII No. 12, 1509-1529. KARL, F. (1964) - Anwendung der Gefügekunde in der Petrotektonik. Tektonische Hefte 5; 142 p; Clausthal-Zellerfeld. KIRALY, L. .(1968) - Eléments structuraux et alignement de phéno- mènes karstiques (Région du gouffre du Petit-Pré de St. Livres. Jura vaudois). Bull. Soc . Neueh. Sc. nat. 91 : 127-146; 10 fig. - (1969)a - Statistical analysis of fractures (orien- tation and density). Geoi. Rundschau, 59/1 : 125-151. - (I969)b - Anisotropie et hétérogénéité de le per- méabilité dans les calcaires fissurés. Eclogae Geoi. HeIv. - (1970) - L'influence de l'hétérogénéité et de 1'anisotropie de la perméabilité sur Ibs sys- tèmes d'écoulement. Bull. Ver. Schweiz. Petrol. Geol.-v. Ing., vol.' 37., Nr. 91 : 50-57. - (1971) - Groundwater flow in heterogeneous, anisotropic fractured media: a simple two- dimensional electric analog. Journal of Hydrology 12 : 255-261. - (1973) - Notice explicative de la courbe hydro- qéologique du canton de Neuchâtel. Suppl. Bull. soc. neuch. 5c. nat. 96; 15 p., 6 fig., I tabi. 1 carte. KlRALT, L., MATHEY, B., TRIPET, J.P. (1971) - Fissuration et orien- tation des cavités souterraines : région de la Grotte de Milandre (Jura tabulaire). Bull. Soc. neuch. 5c. nat., 94 : 99-114 213 KIRALY, L-, SIMEONI, G.-P. )1971) - Structure géologique et orien- tation des cavités karstiques : la grotte de Chez-le-Brandt (Jura neuchâtelois). Bull. Soc. neuen. 5c. nat., 94 : 91-97. KIRALY, L., MOREL, G. (1976a) - Etude de la régularisation de l'Areuse par-modèle mathématique. Bull. Centre d 'Hydrogeologie de Neuchâtel 1; 19-35. KIRALY, L., MOREL, G. (1976b) - Remarques sur 1'hydrogramme des sources karstiques simulées par modèles mathémat- ques. Bull. Centre d'Hydrogéologie de Neuchâtel 1: 37-60. KLINGBEIL, E. (1966) - Tensorrechnung für Ingenieure. Bibi. Instit. Mannheim, 194 p. KNORING, L.D. (1968) - Relations entre orientation des fractures et éléments de la structure géologique (en russe), In: "Problèmes de géoloqie mathématique", Leningrad. KRAU5E, H. (1966) - Oberflächennahe Auflocherungserscheinungen in Sedimentgesteinen Baden-Württembergs. JahreBhefte geol. Landesamt Baden-Württemberg, B, 269-323 S., Abb.. 56-69, Tafel 15-23, Tab. 51, Freiburg im Breis- gau 1966. KRUSEMAN, G.P., DE RIDDEN, N.A. (197Q) - Analysis and evaluation of pumping test data. Internat. Inst, for Land Recla- mation and Improvement, Wageningen; Bull. 11, 200 p. LATTMAN, L.H., FARIZEK, R.R. (1964) - Relationship between fracture traces and the occurence of ground water in carbo- nate racks. J. of Hydrology; 2/3 : 73-91. LeGRAND, H.E., STRINGFIELD, V.T. (1966) - Development of permea- bility and storage in the Tertiary limestones of the south eastern States, U.S.A. Bull. AIHS, 11/4 : ^ 61-73. LEINFELLNER, W. (1965) - Einführung in die Erkenntnis-und WiBsen- schaftstheorie. B. I. Mannheim, 207 p. LIAK0POUL05, A. (1965) - Variation of the permeability tensor ellipsoid in homogeneous anisotropic soils. Water Resources Research 1/1, pp. 135-151, 3 fig. LIPSCHUTZ, 5. (196B) - Linear algebra. Mc. Graw-Hill, 334 p. LOUIS, Cl. (1974) - Introduction a l'hydraulique des roches. Bull. B.R.G.M. (2), 111/4 : 2B3-356. MANDEL, 5. (1966) - A conceptual model of ksrstic erosion by ground water. Fluii. AIH5, XI/1 : 5-7. MANGIN, A., (1975) - Contribution à l'étude hydrodynamique des aquifères karstiques. Thèse, Institut des sciences de la Terre de l'Université de Dijon, MAXEY, G.B. (1964) - Hydrostratigraphic units. Jour, of Hydrology 2/2, pp. 127-129, 2 fig., 3 réf. MEYBDOM, P. (1966) - Current trends in Hydroqeology. Earth- Science Reviews 2 : 34 5-356. MIJATOVIC, S/ (1970) - A method of studying the hydrodynamic regime of karst aquifers by analysis of the discharge curve and level fluctuations during recession. Bull. 1nst. Geoi. and Geaphys. Res ., ser. B, No 8 : 41-74. MCRRI5, D.A., JOHNSON, A.I. (1967) - Summary of hydroloqic and physical properties of rock and soil material as analysed by the Hydroloqic Laboratory of the U.S. Geological Survey 19^0-60. Geol. Survey Water-Supply Paper 1B39-D; 42 p. MULLER, L. (1963) - Der felsbau. 1 Bd. Theoretischer Teil. 624 5., 307 Abb, 22 Taf. Stuttgart 1963. NEWBERRY, B (1966) - The perched water table in the Upper lime- stones aquifer of Malta. J. Inst. Water Engineers; 22/B : 551-570. PALOC, H. (1964) - Caractéristiques hydroqéologiques des dolomies de la région languedocienne. Mém. Centre d'Etudes et de Rech. Hydrogéol., Montpellier, 1 ; 123-127. FAPADCP0UL0S. 1.5. (1967) - Nonsteady flow to a well in an infinite anisotropic aquifer. AIHS pubi. No 73 ; 21-31. PAVONI, N. (1961) - Faltunq durch Harizontaiverschiebung. EGH 54, pp. 513-534, 9 fiq. PRICE, N.J. (1966) - Fault and joint development in brittle and semi-brittle rock. Pergamon Press, Dxford, etc., 1966, 176 p., 58 fig. nATS, M.V., CHERNYA5H0V, S.N, (1967) - Statistical aspect of the problem on the permeability Df the jointy rocks. Publication AIH5 No. 73 (Colloque de Dubrovnik 1965, Hq. des roches fissurées) pp. 227-236. 215 RHODES, R.t SINACORI, M.N. (1941) - Pattern of ground-water flow and solution. J. of Geology 49/8 : 7B5-794 5CHEIDEGGER, A.E. (I960} - The physics of flow through porous media. University Df Toronto Press, 313 p. - (1965) - On the statistics of the orientation of beddinq planes, grain axes and similar sedimento- logie al data. U.S, Geol. Survey Prof. Paper 525-C : 164-167. SCHOLLER, H, (1962) - Les eaux souterraines. Hydrologie dynamique et chimique. Recherche, exploitation et évaluation des ressources. Mass dp S. Cie. , Paris, 1962, 642 p., 1B7 fiq. - (1967) - Hydrodynamique dans le Warst. (Ecoulement et emmagasinement). Publications AIH5 No. 73, (Hydrologie des roches fissurées, Dubrovnik 1965) pp 3-20, 7 fiq., et: Chronique d'Hydrogéologie 10 : 7-21, 7 fig. (1967). 5IME0NI, G.P. (1976) - Etude de la perméabilité des formations calcaires du Jura neuchâtelois. Bull. Centre d'Hydrogeologie de Neuchâtel, 1 :.9-16. - (1977) - Analyse statistique des perméabilités ponctuelles des calcaires aquifères du Jura neuchâtelois. Ann. scient. Univ. Besançon, Fase . 25, 3ème série. SNOW, D.T. (196B) - Rock fracture spacinqs, openings, and poro- sities. Jour, of the Soil Mechanics and foundations division. Vol. 94, no. SMl, January 1968; pp. 73-91, 7 fig. - (1969) - Anisotropic permeability of fractured media. Water Resources Res., 5 : 1273-1289. 5WINNERTDN, A.C. (1949) - Hydroloqy of limestone terrains. (in: Meinzer, D.E. editor, Hydmloqy, Physics of the Earth; New York, Dover Pubi.: 656-C77. TCICHMANN, H. (1964) - Physikalische Anwendungen der Vektor- und Tensorrechnunq. Bibl. Inst. Mannheim, 231 p. THDMA5, H.E., LEOPOLD, L.D. (1964) - Groundwater in North America. Science, 143 : 1001-1006. 216 ULLMO, J. (1967) - Les concepts en physique. In: Logique et connaissance scientifique, Encyclopédie de la pléiade, p. 623-705. TRIPET, J.-P. (1972) - Etude hydrogéologique du bassin de la source de 1'Areuse. Université de Neuchâtel, thèse, 183 p. VISTELIUS, A.B. (1966) - Structurel Diagrams. Pergamon Press Ltd., Oxford, London 1966, XI + ITB p., àl fig. 16 tabi., réf. WEGMANN, E.C. (1961) - Anatomie comparée des hypothèses sur les plissements de couverture (le Jura plissé). Reprinted from the Hull. Geol. Inst. Univ. Uppsala Vol. XL: 196-182, 2 fig. WITTKE, W., LOUIS, Cl. (1968) - Modellversuche zur Durchsträmmung klüftiger Medien. Felsinee h. u. Ing. Geologie, Suppl. IV: 52-7B, 2d fig. Z1ENKIEWICZ, O.C. (1971) - The finite element méthode in engineering science. Mc. Graw-Hill, London; 521 p. 20TL, J. (1961) - Die Hydrographie des nordostalpinen Karstes. 5teierische Reitrage zur Hydrogeol., (1960/61 : Heft 2) .