Environment">
RP 60222 FR
RP 60222 FR
RP 60222 FR
Vérificateur : Approbateur :
Nom : C. Lamotte Nom : M. Audibert
I
Mots clés : karst, modélisation numérique, hydrogéologie, crue, pompage, modèle réservoir,
modèle distribué, hydrodynamique
© BRGM, 2012, ce document ne peut être reproduit en totalité ou en partie sans l’autorisation expresse du BRGM.
Modèles numériques du Lez
Synthèse
- modèle spatialisé,
Sommaire
1. Introduction .......................................................................................................... 21
3.4. CONCLUSION................................................................................................. 71
et la courbe orange en pointillés est celle simulée par le modèle VENSIM . Les barres
inversées représentent le cumul de pluie journalier (pluie efficace). ......................................... 125
Figure 62 : Piézométrie simulée à partir de la pluie pondérée et de l’évapotranspiration
a priori sur l’année test 2003. La courbe noire correspond à la piézométrie mesurée
dans le drain principal de la source du Lez, la courbe verte en tirets est la piézométrie
simulée par le réseau de neurones, et la courbe orange en pointillés est celle simulée
par le modèle VENSIM. Les barres inversées représentent le cumul de pluie journalier
(pluie pondérée). ........................................................................................................................ 126
Figure 63 : Schéma de fonctionnement du GR4J , Perrin (2002).............................................. 130
Figure 64 : Tableau de calage du GR4J .................................................................................... 131
Figure 65 : Résultats des simulations du débit du Lez à Lavalette avec le GR4J
(Observations en bleu, simulation en rouge) ............................................................................. 132
Figure 66 : Débit simulé par le GR4J calé en fonction du débit observé à Lavalette ................ 133
Figure 67 : Débit simulé par le GR4J calé (rouge) et débit observé à Lavalette (bleu) au
cours de l'automne - hiver 2002 ................................................................................................. 134
Figure 68 : Schématisation de la prise en compte de l'état de remplissage du souterrain ....... 137
Figure 69 : Schéma de fonctionnement du modèle SMA-SCS.................................................. 138
Figure 70 : Schéma représentant la fonction de production SCS aménagée pour le
bassin du Lez ............................................................................................................................. 140
Figure 71 : schéma de la fonction de transfert « Lag and Route Simple » ................................ 142
Figure 72 : Hydrogrammes simulés et observés pour la crue de décembre 2002. La
zone entourée repère une production inapropriée de ruissellement direct alors que ce
dernier n’a pas commencé en réalité. (en bleu, les débits observés; en gris les débits
simulés avec SMA SCS) ............................................................................................................ 144
Figure 73 : Exemple de l’épisode de décembre 2002 illustrant le problème lié à la
vidange du réservoir de production pour le modèle SMA. Le graphique du haut
compare les hydrogrammes observé (tirets bleus) et simulés avec une vidange de ds =
-1 -1
0.28 j (courbe grise avec les cercles) ou de ds = 0.01 j (courbe grise avec les croix).
Le graphique du bas compare le niveau piézométrique normalisé à Gour Noir (tirets
bleus) aux niveaux normalisés dans le réservoir « sol » pour une vidange de ds = 0.28
-1 -1
j (courbe grise avec les cercles) ou ds = 0.01 j (courbe grise avec les croix). ...................... 146
Figure 74 : Hydrogrammes observé (en tirets bleus) et simulé avec la fonction de
production SMA (en trait continu gris) ou SCS aménagée (en trait continu noir) pour
une des crues testées. La zone entourée repère les périodes en début d’épisode où le
SMA-SCS produisait du ruissellement direct alors que ce dernier n’a pas commencé en
réalité.......................................................................................................................................... 148
Figure 75 : Graphiques permettant de comparer, pour chacun des 6 épisodes de crue,
l’évolution de la piézométrie à Gour Noir avec celle du niveau dans le réservoir « sol »
des fonctions de production SMA-SCS (en gris) et SCS aménagée (en noir). ......................... 148
Figure 76 : Comparaison de l’hydrogramme observé (en tirets bleu) avec celui de la
simulation « de référence » (en noir) pour l’épisode d’octobre 2001 ......................................... 149
Figure 77 : Distribution de l’occupation des sols du bassin versant du Lez. Le karst est
en rose et la couverture en bordeaux. ....................................................................................... 151
Figure 78 : Evolution des critères de Nash obtenus en fonction de la distribution de S ........... 153
b
Figure 105 : Hydrogrammes observé (courbe bleue), simulé à partir de l’ébauche S =
a
160 mm (courbe noire) et simulé à partir de l’analyse S = 131 mm (courbe rouge)
après assimilation d’une donnée (croix bleue) pour l’épisode d’octobre 2001. Le trait
vertical noir sépare la période d’assimilation de la période de prévision ................................... 211
Figure 106 : Graphiques représentant la fonction coût J(x) (courbe rouge) et son
minimum (rond rouge) ainsi que l’évolution de la fonction coût incrémentale J(xl ; x)
a
(courbe bleue), de son minimum x (rond bleu) et de sa tangente (droite noire) au point
de linéarisation xl (croix noire) pour les 2 premières itérations de la boucle externe pour
l’épisode d’octobre 2001 ............................................................................................................ 212
Figure 107 : Résultats de l’assimilation de données pour l’épisode de décembre 1997.
La courbe bleue correspond aux observations ; la courbe noire est l’hydrogramme
simulé à partir de l’ébauche (avant assimilation) et la courbe rouge est l’hydrogramme
simulé à partir de l’analyse (après assimilation). Les croix bleues sont les données
assimilées................................................................................................................................... 217
Figure 108 : Hydrogrammes issus du cas fictif réalisé en « assimilation groupée » (a) et
en « assimilation séparée » (b). La courbe rouge représente l’hydrogramme obtenu à
partir de l’analyse, la courbe bleue représente l’état vrai, les croix bleues correspondent
aux observations, les ronds bleus sont les observations assimilées. Le trait noir vertical
représente l’instant de prévision. ............................................................................................... 219
Figure 109 : Evolution des charges hydrauliques mesurées par le réseau de suivi du
Lez. Les évolutions piézométriques de Fontbonne (données CAM) et de la source des
Fontanilles (Onema-Brgm) sont également reportées. .............................................................. 235
Figure 110 : Courbe de tarage théorique du déversoir de la source du lez .............................. 237
Figure 111 : Réponses impulsionnelles des composantes lentes des modèles de
transferts (ModHE_Ginger Intermédiaire_128j ; ModHE_ Ginger Intermédiaire _256j ;
ModHE_ Ginger Intermédiaire _366j) utilisés pour reconstituer les débits naturels de
Lez en période de basses eaux ................................................................................................. 238
Figure 112 : Réponses impulsionnelles des composantes rapides des modèles de
transferts (ModHE_Ginger Intermédiaire_128j ; ModHE_ Ginger Intermédiaire _256j ;
ModHE_ Ginger Intermédiaire _366j) utilisés pour reconstituer les débits naturels de
Lez en période de basses eaux ................................................................................................. 238
Figure 113 : Exemples de reconstitutions de débit naturels obtenus à l’aide des
modèles de transfert qualifiés par differents temps de régulation (compris entre 128 j et
366j). .......................................................................................................................................... 239
Figure 114 : Cartes de localisation des postes pluviométriques utilisés pour le modèle
hydrologique ATHYS. ................................................................................................................. 249
Figure 115 : Carte de localisation des postes pluviométriques utilisés pour les modèles
hydrogéologiques TEMPO, VENSIM et Réseau de Neurones. ................................................. 251
1. Introduction
L’hydrosystème du Lez, fleuve côtier dont la source d’origine karstique est exploitée
par une gestion active depuis 1981 par la ville de Montpellier puis l’Agglomération de
Montpellier, constitue un site de référence pour étudier, dans le cadre d’un projet
d’étude pluridisciplinaire (géologie, hydrogéologie, hydrologie, hydraulique, biologie et
économie), les possibilités d’une gestion multi-usages.
Le projet de recherche « Gestion multi-usages des aquifères karstiques
méditerranéens – Le Lez, son bassin versant et son bassin d’alimentation » regroupe
un partenariat scientifique et technique constitué par le BRGM, les UMR
Hydrosciences, G-EAU, TETIS, le CERFACS et BIOTOPE, ainsi q'une collaboration
supplémentaire avec l'Ecole des Mines d'Alès (EMA). Les questions d'intérêt qui se
posent dans ce cadre concernent :
- l'amélioration de la connaissance de l'aquifère karstique du Lez,
- l'évaluation de la vulnérabilité de l’aquifère karstique et les contraintes associées en
termes d’aménagement du territoire,
- l'évaluation de la ressource en eau exploitable au sein de l’aquifère,
- la caractérisation des impacts éventuels de nouveaux scénarios d’exploitation,
- la mitigation des risques d’inondation,
- la caractérisation des interactions du pompage de la source du Lez sur d’autres
systèmes connexes,
- la connaissance de la biodiversité souterraine et la définition de bio-indicateurs des
eaux,
- la gestion concertée de la ressource.
La mise en place des différents modèles hydrogéologiques, hydrologiques et
hydrauliques reproduisant le fonctionnement de l'hydrosystème karstique du Lez et
permettant d'apporter des éléments de réponse à ces questions font l'objet de l'atelier
AT4.
Les résultats attendus pour cet atelier sont :
i) reproduire et expliquer les fonctionnements et évolutions observés concernant la
gestion de la ressource en eau souterraine et les crues de surface ;
ii) évaluer et proposer des scénarios de gestion alternatifs des eaux souterraines et de
surface, en contexte ou non de changement climatique.
Dans le présent rapport L4.1 sont présentés : une description de ces modèles, leurs
hypothèses, leur sensibilité, leur validation et les résultats d'application obtenus sur
une période rétrospective.
Les résultats des scénarios prospectifs et les mises en cascade des modèles feront
l'objet d'un prochain rapport (L4.2).
1
Temps de réponse du bassin versant estimé à 2 heures pour les crues des 17/12/1996 et 9/12/2002 par
exemple
2
Module annuel sur la période 1975 - 2011 de 2.14 m3/s source : banquehydro et Débit maximum de crue
de 467 m3/s pour la crue de septembre 2005
Les tâches 2, 3, 4, 6 et 7 font l’objet du présent rapport. Les résultats de la tâche 1 ont
été présentés dans le livrable de l’AT1 (Rapport BRGM/RP-60041-FR). Les chroniques
de pluie utilisées sont décrites dans l’annexe 2 du présent rapport.
Figure 2 : Evolution de la piézométrie (en m NGF) du drain karstique à la source du Lez (suivi
véolia), au piézomètre de Claret (suivi véolia) et au forage de secours du Suquet (suivi CG34)
« État hydrologique »
Condition de débit
réservé
Piézométrie
Pompage
« Systeme karstique »
Demande AEP
En étiage, lorsque la source est tarie, le débit de pompage est supérieur au débit qui
se serait naturellement écoulé en absence de pompage. Le débit pompé lors des
périodes de tarissement de la source peut s’interpréter comme la somme du débit qui
se serait naturellement écoulé en absence de pompage et du débit d’eau issue des
réserves (déstockage des drains karstiques et mobilisation de la matrice carbonatée
dans laquelle s’est développé le réseau de drainage karstique). Les pompages en été
mobilisent les réserves en eau du système naturellement inaccessibles en étiage, ce
qui se traduit, par une diminution de la charge piézométrique dans le drain karstique.
Une fois le débit naturel calculé, un deuxième modèle de transfert est utilisé :
« piézométrie du Lez ». Ce modèle comporte 4 entrées :
- la pluie efficace. Dans le modèle piézométrique du Lez, la contribution de la pluie
efficace a été fixée à 1% (contribution négligeable) car l’information recharge par les
pluies efficaces est déjà portée par les chroniques piézométriques utilisées en
entrée pour décrire le fonctionnement du drain karstique ;
- la fonction de sollicitation des réserves du karst calculée à l’aide des débits naturels
modélisés et des pompages (qui engendre une diminution du niveau piézométrique
lorsque le débit de pompage devient supérieur au débit naturel de la source) ;
- la contribution en provenance du compartiment Jurassique situé à l’Ouest de la faille
de Corconne. Cette contribution est estimée à l’aide de la variable piézométrique
Suquet qui est générée par la méthode de transfert « piézométrie Suquet » (sous
module de transfert) ;
- la contribution en provenance du compartiment Crétacé à l’aide de la variable
piézométrique Laudou qui est générée par la méthode de transfert « Piézométrie
Laudou » (sous module de transfert). Cette entrée prend en compte les effets liés
au pompage par l’intermédiaire de prise en compte de la fonction de sollicitation.
Les variables hydro-climatiques utilisées dans les simulations sont les suivantes :
- variable climatique : température moyenne mensuelle de la station de St-Martin de
Londres (1960-2008). L’éloignement du site considéré n’est pas préjudiciable dans
l’approche de modélisation utilisée car seules les variations relatives sont prises en
compte dans les calculs (les variables d’entrées sont dites « centrées») ;
- variable pluie : pluie à St-Martin-de-Londres (données Météo France, 1962-2008),
Montpellier Fréjorgues (données Météo France, 1962-2008); Valfaunes plaine
(données Météo France, 1960-2008). La fonction de pondération utilisée (cf. §4.3.2.
rapport BRGM/RP-60041-FR, Jourde et al., 2011) est la suivante : Pluie_BV =
0.32*St-Martin-de-Londres + 0.56*Valfaunes + 0.13*Montpellier
3
This work was funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under
grants no. LI 727/11-1 and SA 501/24-1.
La méthodologie utilisée pour estimer les débits naturels est présentée en détail dans
le rapport BRGM/RP-60041-FR (Jourde et al., 2011). Une nouvelle courbe de tarage
du seuil du Lez a récemment été produite par Ginger Environnement. Cette courbe de
tarage dite intermédiaire (tableau en annexe) a été établie en considérant des
variations progressives du coefficient de débit (Kv) compris entre 0.34 et 0.42 entre les
basses et hautes hauteurs d’eau à la vasque. Cette courbe de tarage a été jugée plus
adaptée à la quantification des débits de débordement à la vasque.
Mois 7 8 9 10 11 12
Nb valeurs journalières 1085 1085 1050 1085 1020 1054
minimum 103 81 83 98 83 83
Décile 0.1 403 215 214 272 338 469
Decile 0.2 507 300 283 385 587 821
1er Quartile 533 330 327 439 860 985
Médiane 704 485 462 1408 2071 2066
3ème Quartile 883 609 730 2508 3388 3609
Décile 0.8 934 636 1243 3014 3704 3916
Decile 0.9 1101 746 2215 4092 5070 5041
maximum 3934 2662 10043 9963 10459 16358
Moyenne 771 514 882 1807 2419 2554
ecart-type 441 303 1170 1709 1885 2054
Tableau 2 : Statistiques descriptives des débits naturalisés de la source du Lez (valeurs en l/s,
période 1976-2009)
b) Débit de pompage
Figure 6 : Evolution de la piézométrie à la source du Lez depuis 1970 et des débits naturels
estimés depuis 1975. Les débits de pompage ainsi que les débits de sollicitation des réserves
sont également reportés (depuis 1974).
Prélèvement
Mois 1 2 3 4 5 6
Nb valeurs journalières 1085 989 1085 1050 1085 1050
minimum 466 555 571 318 530 414
Décile 0.1 835 852 886 887 919 887
Decile 0.2 896 902 927 932 960 963
1er Quartile 915 920 949 952 984 998
Médiane 998 1013 1027 1033 1065 1105
3ème Quartile 1060 1078 1105 1117 1168 1222
Décile 0.8 1078 1109 1137 1149 1195 1245
Decile 0.9 1168 1178 1201 1217 1254 1307
maximum 1419 1387 1385 1529 1498 1525
Moyenne 992 1006 1035 1041 1077 1100
ecart-type 133 132 121 133 131 167
Mois 7 8 9 10 11 12
Nb valeurs journalières 1085 1085 1050 1085 1020 1054
minimum 388 435 245 140 316 459
Décile 0.1 766 645 596 649 722 763
Decile 0.2 910 731 709 812 869 864
1er Quartile 966 781 771 849 894 886
Médiane 1097 1005 964 945 972 988
3ème Quartile 1184 1094 1041 1035 1040 1060
Décile 0.8 1205 1118 1054 1054 1058 1077
Decile 0.9 1260 1169 1109 1117 1114 1117
maximum 1453 1369 1343 1337 1577 1335
Moyenne 1055 947 907 914 949 959
ecart-type 193 198 191 192 159 147
7 8 9 10 11 12
Nb valeurs journalières 1085 1085 1050 1085 1020 1054
minimum -1001 -1018 -1037 -952 -905 -934
Décile 0.1 -737 -808 -784 -639 -447 -298
Decile 0.2 -613 -714 -669 -406 -180 -112
1er Quartile -565 -671 -640 -282 -103 -15
Médiane -350 -513 -416 0 0 0
3ème Quartile -107 -234 -47 0 0 0
Décile 0.8 -54 -162 0 0 0 0
Decile 0.9 0 -57 0 0 0 0
maximum 0 0 0 0 0 0
Moyenne -356 -466 -383 -168 -103 -73
ecart-type 269 271 305 255 208 164
Tableau 4 : Statistiques descriptives des débits de sollicitations des réserves du karst (valeurs
en l/s, periode 1976-2009)
Le terme modélisation inverse recouvre ici les techniques visant à reproduire un débit,
un niveau piézométrique ou un flux (transfert de masse) à partir d’une chronique de
pluie et d’ETP (ou de température moyenne journalière), éventuellement de
concentrations ou de paramètres physico-chimiques. La modélisation inverse fait appel
à des techniques numériques complexes s'appuyant sur différentes méthodes de
régularisation (Pinault et al 2001a, b)
Dans un aquifère karstique, le transfert rapide est souvent associé à l'effet de chasse
qui peut se produire lorsqu'une continuité hydraulique s'établit entre la zone non
saturée et la zone noyée, suite à une recharge importante (mise en connexion d’une
nappe temporaire dans l’épikarst avec la zone noyée au travers de fractures et/ou
conduits subverticaux (Pinault et al, 2001).
Le transfert lent met en jeu les différents processus d'infiltration induisant un retard plus
ou moins important entre la pluie et le débit (ou la variation du niveau piézométrique)
observés en sortie.
Le seuil Ω (t) est calculé à partir de la pluie et de l’ETP (ou de la température de l’air
moyenne) des jours précédents de sorte que :
(Figure 10). Ce seuil Oméga est estimé à partir de la pluie et d’une variable climatique
(ETP, Montpellier Fréjorgues) de telle façon que seule la hauteur de pluie située au-
dessus de ce seuil génère un débit à l’exutoire. Le seuil est calculé à partir d’une
réponse impulsionnelle à la pluie et à la température de l’air (Figure 7), également
obtenue par inversion.
L’option « modèle paramétrique » de Tempo a été utilisée afin de bien contraindre les
débits d’étiage du Lez. Dans l’approche paramétrique utilisée, la recharge par les
précipitations est décrite par une réponse impulsionnelle qui convolue une fonction
analytique complexe constituée de deux réservoirs (Modèle de Samani et Ebrahimi,
1996), la vidange de la zone noyée étant décrite par une fonction exponentielle
décroissante.
4
Pas d'échantillonnage: 1.00 j ; Nombre d'années (calage) : 12.01342 [21/07/1995-04/09/2006]
TRANSFERT; Système non linéaire ; Minimisation de la distance S ; Aire(modèle)=Aire(Observations)
VARIABLE de SORTIE :Q_Nat_HE_BE
VARIABLE d’ENTREE :Pluie :pl_Lez, ETP :TMM_SML (Température moyenne mensuelle à St-Martin de
Londres
Surface en km2 : 115.5622
$$$ Q_Lez_Final Q_Nat_HE_BE ; écart-type: 0.64955E+03 l/s R2= 0.926 Nash=0.8595
sont donc inférieurs à 1 an. L’effet pluriannuel d’une succession d’années déficitaires
en précipitation ne pourra donc pas être pris en compte de manière précise.
Figure 10 : Calage. Comparaison du débit modélisé au débit naturel estimé de la source du Lez
Nous avons supposé, sur la base des résultats de la chronique de l’ouvrage Suquet
(suivi CG34) que le compartiment jurassique situé à l’Ouest de la Faille Corconne n’est
pas influencé par les pompages effectués à la source du Lez.
N 2
1
H (t ) * Reff
Le niveau piézométrique en sortie est exprimé par rapport au niveau de basses eaux
(H0) qui est estimé à partir de la moyenne m et l’écart-type du niveau piézométrique
mesuré (Pinault et al, 2005) :
H0 m
Dans le cas d’un aquifère présentant des temps de régulation inférieur à 1 an, le
coefficient est égal à 2.5 ce qui correspond à une période de retour de 160 ans de H0
à condition que H(t) suit une distribution gaussienne. Cette formulation est bien
adaptée pour les aquifères poreux (Pinault et Schomburgk (2006). Dans le cas des
karsts, les variations piézométriques suivent rarement une distribution gaussienne, le
niveau de base est donc fixé en fonction de la connaissance de la structure et du
fonctionnement du milieu. Pour le cas de l’ouvrage Suquet (09903X0111/F), le niveau
de base est celui mesuré en basses eaux lorsque la charge devient constante (H0= 60
m NGF), suite à la déconnexion piézométrique du compartiment jurassique situé à l’est
de la faille de Corconne.
Le modèle a été utilisé pour estimer la piézométrie au sein du karst pour la période
antérieure à 2002 et postérieure à 2004 (Figure 14, Figure 20).
2.4. DISCUSSION
Dans le but d’évaluer la sensibilité du modèle au débit naturel estimé, nous avons
simulé la piézométrie du Lez en considérant que le débit de la source est connu à
±10% près. Les résultats des nouvelles simulations piézométriques indiquent que le
modèle est globalement peu sensible à l’incertitude sur les débits estimés du Lez
(Figure 22). Par contre, les résultats indiquent que les biais observés dans la
piézométrie sont causés soit par une mauvaise estimation de la fonction de production
de pluie efficace (seuil oméga trop faible) ou par une surestimation du cumul de pluie
considérée en entrée du modèle. Pour conclure, le modèle de transfert « débit du
Lez » apparaît optimiste lorsque des épisodes pluvieux significatifs (H> 25 mm) se
produisent au cours de l’étiage. La piézométrie simulée lors de ces périodes apparaît
alors surestimée.
2.5. PERSPECTIVES
Les modèles développés dans cette étude seront ultérieurement utilisés dans le
cadre du projet pour aborder l’effet du changement climatique et des scénarios
de gestion sur les ressources en eau du système karstique du Lez.
3.1. RESUME
Le modèle est constitué d’un réservoir principal assimilé à la zone noyée (ou zone des
réserves) de l’aquifère karstique. Sa structure schématique est fournie dans la Figure
23 ci-après. Les différents modules présentés dans le schéma seront expliqués dans
les paragraphes suivants. Le modèle simule ainsi les variations de niveau dans la zone
noyée ainsi que les débits à la source mesurés à la station de la DIREN.
Pluie ETP
Peff
Hsol min
Sol
Hpompé
Débordement
xdébordement * Peff Hseuil
α
Hsortant Hréservoir
0
Débit sortant réserves
3.2.1. Alimentation
Les pluies sont définies à partir d’une pondération de la pluie mesurée sur les stations
de Valflaunès, Prades le Lez et Saint Martin de Londres (Fleury et al., 2009). Cette
dernière est obtenue par le logiciel TEMPO© du BRGM en se basant sur le calcul des
corrélations croisées entre les séries chronologiques fournies (débit mesuré, mesures
de pluie issues de plusieurs pluviomètres, évapotranspiration). La pondération optimale
définie sur ces trois postes est : 42,7 % Valflaunès, 18,7 % Prades le Lez et 38,6% St
Martin de Londres. Cette pondération est différente de celle utilisée par le modèle
Tempo au chapitre précédent car la période de simulation du modèle Vensim est plus
courte et l’utilisation de la station de Prades le Lez est donc possible.
Selon les modèles, l’ETP sera soit une ETP artificielle issue d’une fonction sinusoïdale,
soit la donnée de la station Météo France de Montpellier Fréjorgues. L’ETP sinusoïdale
est estimée à partir de la formule de Thornthwaite en utilisant une interpolation grâce à
une fonction sinusoïdale des données d’ETP de la station Météo France (Tritz et al.,
2011).
Le niveau d’eau dans le réservoir sol est ainsi calculé pour chaque pas de temps selon
l’équation suivante :
Hsol (t+1) = Hsol (t) + Pluie (t+1) – ETP (t+1) – Peff (t+1) ;
et Peff = Hsol si Hsol > 0 sinon Peff = 0.
Aussi le réservoir sol est caractérisé par un niveau minimum (Hsol min) défini par calage.
Lorsque ce réservoir est saturé (niveau positif) alors la hauteur d’eau est évacuée vers
le réservoir zone noyée, cette hauteur correspond à la Pluie efficace (Peff).
3.2.2. Prélèvements
Les prélèvements réalisés dans le système sont également pris en compte. Le débit
prélevé est transformé en hauteur d’eau unitaire (débit prélevé divisé par la surface du
bassin d’alimentation) et est soustrait à chaque pas de temps à la hauteur dans le
réservoir principal. Il s’agit de la composante Hpompé présentée dans l’illustration
précédente.
Le débit est obtenu en intégrant la hauteur sortant (Hsortant) sur la surface du bassin
d’alimentation (110 km²). Aussi, un seuil de débordement a dû être ajouté afin de
représenter le comportement du système lors des fortes crues. En effet, l’étude
hydrogéologique, conformément aux observations de terrain, révèle lors des fortes
crues, une augmentation moins forte des niveaux piézométriques et des débits à la
source. Parallèlement de nombreux trop-pleins se mettent en place sur le système
(source du Lirou, source de la Fleurette). Ainsi pour représenter ce phénomène un
seuil avec évacuation d’une partie du flux a été mis en place dans le modèle. Une fois
le seuil atteint, une partie de la recharge par les pluies est évacuée. La valeur de ce
seuil (Hseuil) est déduite d’un calage tout comme la part du flux évacuée du système
(xdébordement).
Trois paramètres calés sont nécessaires pour simuler les débits, il s’agit des valeurs de
Hseuil, xdébordement et α.
Un même calage est effectué à la fois pour le débit et la piézométrie, de sorte que le
modèle permette de simuler en même temps ces deux grandeurs sans modifier les
paramètres.
Trois simulations vont être testées : la première sera réalisée à partir de l’ETP
sinusoïdale, la seconde à partir de l’ETP Mauguio, enfin pour la dernière simulation, la
structure du modèle va être simplifiée : le réservoir sol sera supprimé, les flux
Les modèles vont être calés sur l’ensemble de la chronique à savoir du 01/09/1996 au
31/05/2010. La qualité du modèle sera approchée à l’aide du calcul du critère de Nash.
Ce critère sera calculé sur 3 périodes distinctes à savoir : la totalité de la chronique, la
période 1 du 01/09/2002 à 31/08/2003 et la période 2 du 01/09/2003 au 31/08/2004.
Pour chacune des simulations les valeurs des paramètres calés ainsi que les critères
de qualité seront fournis.
Une simulation de bonne qualité est obtenue avec les paramètres calés suivants :
- Hsol min = - 30 mm ;
- α = 0,04 ;
- Hseuil = 0,09 m ;
- Xdébordement = 0,8 ;
- P1 = 0,1 ;
- P2 = 0,05.
Les résultats graphiques de cette simulation sont présentés Figure 25 et Figure 26.
On vérifie graphiquement que les débits sont relativement bien reproduits. Les étiages
simulés sont conformes aux étiages mesurés. Dans l’ensemble, les premières crues
automnales simulées (faisant suite à la période d’étiage) sont synchrones avec celles
observées. Cette bonne représentation est à mettre en relation avec la simulation de la
piézométrie. En effet lors du calage du modèle pour la représentation de la
piézométrie, un grand intérêt a porté sur le calage temporel des périodes de
sollicitation (début de prélèvement sur les réserves et désaturation du karst) et de
recharge du karst (remontée des niveaux piézométriques aboutissant à la recharge du
karst).
Aussi, il est nécessaire de remarquer que l’intensité des quatre plus fortes crues sur la
période est sous-estimée. Ceci est lié à la mise en place du seuil. Plusieurs essais ont
été réalisés afin de diminuer son impact sur ces quatre crues (diminution de la fraction
de débordement ou augmentation du niveau du seuil), ils ont tous conduit à une forte
augmentation des crues significatives si bien qu’à l’inverse le modèle surestimait
globalement l’ensemble des débits des crues. Ainsi, bien qu’imparfait pour la
simulation des quatre plus fortes crues, ce calage a été conservé car jugé satisfaisant
sur l’ensemble de la chronique.
La qualité du modèle est estimée à partir du critère de Nash (exprimé dans le tableau
ci-dessous). Ces résultats valident le modèle mis en place. En effet les critères sont
proches de 80 % pour les débits comme pour la piézométrie quelle que soit la
chronique utilisée. On remarque toutefois des écarts significatifs selon la période de
calcul, montrant l’intérêt de travailler sur l’ensemble de la chronique pour caler le
modèle et non pas sur une période réduite.
Piézométrie 78 % 90 % 71 %
Débit 78 % 77 % 83 %
Figure 25 : Piézométrie simulée et mesurée sur la période septembre 1996 – mai 2010 (ETP sinusoïdale)
Figure 26 : Débit simulé et mesuré sur la période septembre 1996 – mai 2010 (ETP sinusoïdale)
Une simulation de bonne qualité est obtenue avec les paramètres calés suivant :
- Hsol min = - 10 mm ;
- α = 0,04 ;
- Hseuil = 0,09 m ;
- Xdébordement = 0,8 ;
- P1 = 0,1 ;
- P2 = 0,05.
Cette simulation présente des résultats acceptables pour les mêmes paramètres que la
simulation précédente, à l’exception du réservoir sol dont la valeur minimale est égale
à -10 mm.
ce décalage induit une moins bonne représentation des débits sur les périodes de
transition (petites crues en période de basses eaux notamment).
Piézométrie 67 % 89 % 59 %
Débit 75 % 73 % 79 %
Figure 28 : Piézométrie simulée et mesurée sur la période septembre 1996 – mai 2010
Figure 29 : Débit simulé et mesuré sur la période septembre 1996 – mai 2010 (ETP Météo France)
Pluie ETP
Hpompé
Débordement
0.8 * Peff Hseuil = 0.09
α = 0.04
Hsortant Hréservoir
0
Débit sortant réserves
Piézométrie 75 % 96 % 83 %
Débit 79 % 81 % 75 %
Aussi, bien que le modèle soit de bonne qualité, sa structure est difficilement justifiable.
En effet, il n’est pas concevable d’appliquer directement l’ETP à la zone noyée et
permettre à une partie du flux de quitter le système. Ce modèle ainsi que ses résultats
sont toutefois fournis à titre informatif.
Figure 31 : Piézométrie simulée et mesurée sur la période septembre 1996 – mai 2010 – modèle sans réservoir sol
Figure 32 : Débit simulé et mesuré sur la période septembre 1996 – mai 2010 – modèle sans réservoir sol
3.4. CONCLUSION
Deux autres paramètres (porosité équivalente) sont utilisés afin de convertir les
niveaux dans le réservoir en niveau piézométrique.
Ainsi ce modèle pourra être utilisé à la fois pour la gestion des crues mais aussi
des étiages.
4.1. INTRODUCTION
L’objectif de ce modèle est de simuler au mieux l’impact des cycles de pompages sur
la piézométrie au sein du conduit karstique du Lez. La solution analytique classique de
type radial (solution de Theis) pour l’interprétation d’un pompage d’essai peut
représenter le rabattement mesuré dans le conduit karstique à la source du Lez lors
des cycles de pompages (référence au rapport AT1). Cependant, l’ajustement est
effectué au prix de l’utilisation d’un coefficient d’emmagasinement dont la valeur
(supérieure à 10000) n’a aucun sens physique. Cette solution n’est donc pas en
mesure de représenter de façon réaliste et physiquement basée la réponse
hydrodynamique globale du système karstique du Lez sollicité, notamment le
stockage/déstockage d’eau dans le réseau de drains karstiques et leur échange avec
la matrice rocheuse.
Lors du pompage d’essai, on considère que le débit Q pompé dans le drain karstique
provient (Figure 33) :
- de la contribution des blocs de calcaires (et/ou dolomies) poreux et/ou fissurés,
et/ou des systèmes annexes au drainage (SAD) situés entre les drains karstiques et
alimentant ceux-ci ;
- du déstockage de l’eau contenue dans le drain karstique.
5
Le terme « bloc » ne préjuge pas de la nature karstique des vides alimentant le drain : fissures ou
porosité de la matrice ou systèmes annexes au drainage
Fonctionnalités
Il possède une fonction capacitive liée au stock d‘eau situé dans les cavités karstiques
connectées entre elles, constituant le drain karstique.
Description mathématique
Sdrain
Qb Qpomp
Qsource
Hseuil
Hd Drain
karstique
0
Si Hd>Hseuil alors Qsource est calculé par la courbe de tarage à la vasque, sinon Qsource =0
H
(éq. 1)
dH d Qb Q pomp Qsource
(éq. 2)
dt S drain
Le débit de contribution des blocs (Qb) est calculé dans le réservoir « blocs » décrit ci-
dessous. Il comprend le débit naturalisé (Qnat) des blocs (débit de tarissement naturel
de la source, en l’absence de pompage, calculé par l’approche Tempo) ainsi que le
débit résultant de leur sollicitation du fait des rabattements induits par le pompage.
Seuls trois paramètres caractérisent ce réservoir « drain karstique » (Tableau 5). Ils
sont décrits ci-dessous :
- la fonction capacitive du drain karstique correspond à son emmagasinement,
exprimé en m2. Il s’agit du volume d’eau libéré par le réseau karstique lorsque le
niveau d’eau dans le drain baisse d’un mètre, exprimé en m3/m, soit en m2. Le drain
karstique étant constitué de vides de porosité égale à 100 %, son emmagasinement
correspond donc à la surface totale des vides : c’est la surface dénoyée du drain
Sdrain ;
Sur la base des trois paramètres décrits ci-dessus et des débits d’échange entre les
blocs et le drain calculés dans l’autre réservoir, le réservoir « drain karstique » calcule
le niveau d’eau dans le drain (Hd) au moyen de la relation (éq. 2) et le débit de la
source (Qsource) au moyen de la relation (éq. 1).
Les niveaux d’eau calculés dans le drain karstique pourront alors être comparés à ceux
mesurés à la vasque.
Fonctionnalités
Ce réservoir correspond aux blocs poreux (et/ou fissurés) et/ou aux systèmes annexes
au drainage (SAD) situés entre les drains du réseau karstique. Ceux-ci alimentent le
drain karstique selon une relation dépendant de la différence de charge hydraulique qui
règne entre blocs et drain.
6
Au sein de ce chapitre, les termes « hauteur d’eau » et « niveau d’eau » sont fréquemment utilisés en
tant que synonymes de « charge hydraulique » ou de « niveau piézométrique ».
Description mathématique
Simpactée
Recharge=0
Hb
Blocs Qb
0 Blocs
dH b R Qb
(eq. 3)
dt Simpactée
Le débit des blocs vers le drain évoluant en fonction des variations de différence de
charge hydraulique entre blocs et drain, la solution analytique utilisée pour calculer les
échanges entre les deux réservoirs est une solution de type potentiel imposé. Par
ailleurs, il a été vu précédemment (Rapport AT1) que les courbes de rabattements
pouvaient être simulées assez correctement par le modèle de Theis, à savoir avec une
hypothèse d’écoulements radiaux autour des drains karstiques. Dès lors, la solution
d’écoulements radiaux vers un puits artésien, initialement publiée par Jacon et Lohman
(1952) est la plus adaptée7. Etant donnée sa complexité numérique, une version
simplifiée par Perrochet (2005) a été implémentée dans le modèle.
2 T
Q(t ) H0 h0 (éq. 4)
Tt
ln 1
Sr02
Avec :
r0 (m) : rayon du puits, ou, dans le cas présent, le rayon moyen des drains karstiques,
7
A noter qu’une solution considérant une géométrie d’écoulements linéaires vers une structure drainante
de type tranchée a également été testée. Les résultats obtenus, décevants, sont brièvement décrits dans
la section validation.
H0 (m) : charge hydraulique initiale dans l’aquifère, dans les blocs par exemple,
h0 (m) : charge hydraulique constante dans le puits, dans le drain dans le cas présent.
Cette solution a été obtenue pour une nappe captive par intégration de l’équation de
Darcy avec les conditions initiales et aux limites suivantes (schématisées au sein de la
Figure 35) :
H (r ,0) H0 H ( , t) H0
et H (r0 , t ) h0
Figure 35 : Schéma d’un écoulement transitoire vers un puits de rayon r0, dont le potentiel,
initialement à H0, est instantanément imposé à h0 (le niveau piézométrique dans l’aquifère
évolue au cours du temps t1, t2, t3)
Cependant, cette solution analytique ne peut être appliquée aussi simplement au cas
d’un drain karstique sollicité par un pompage. Dans tout système karstique, les
conditions initiales avant pompage ne sont pas des conditions hydrostatiques de type
H(r,0) = H0 partout avec un débit nul des blocs vers le drain. Au contraire, il existe,
initialement en t = 0, une différence de potentiel entre les blocs et le drain, à l’origine
du débit s’écoulant à la source.
Dès lors, la relation (éq. 4) doit été appliquée « en perturbation » pour calculer le débit
lié au rabattement supplémentaire H b (t ) H d (t ) , où les indices « b » et « d » font
respectivement référence aux « blocs » et au « drain », par rapport à un rabattement
initial H b 0 H d 0 générant le débit naturalisé de la source Qnat. De cette façon, le
débit s’écoulant des blocs vers le drain peut s’écrire :
2 T
Qb (t ) Qnat (t ) H b (t ) H d (t ) H b0 Hd0 (éq. 5)
Tt
ln 1
Sr02
En supposant que, durant un cycle de pompage, la variation de niveau dans les blocs
est négligeable en face du rabattement total dans le drain :
H b (t ) H b 0 0
2 T
Qb (t ) Qnat sd (t ) (éq. 6)
Tt
ln 1
Sr02
Où :
- s d (t ) Hd0 H d (t ) désigne le rabattement dans le drain, égal aussi au
rabattement mesuré à la vasque ;
Dans la mesure où les rabattements imposés dans le drain par le pompage varient au
cours du temps, la relation (éq. 6) a été mise en œuvre au moyen du principe de
superposition. On peut se référer à Maréchal et al. (2008) pour plus de précisions.
Le modèle a été testé sur plusieurs cycles de pompages observés ces dernières
années et non perturbés par des pluies. Il a été initialisé avec le niveau d’eau observé
dans la vasque à l’instant t = 0 choisi comme démarrage de la modélisation. Le débit
naturalisé Qnat a été calculé par l’approche Tempo avec la nouvelle courbe de tarage
dite « intermédiaire » (§ 2.2.3. du présent rapport).
Les années 2006 et 2007 ont été choisies comme périodes de calibration car des
données de débits horaires étaient disponibles en 2007, permettant un meilleur calage
de certains paramètres comme nous le verrons dans la suite. Ensuite, le modèle a été
validé sur les périodes 2001 et 2002.
Les résultats de la phase de calibration sont présentés aux Figure 36 et Figure 37.
La simulation de 2007 débute le 15 juin 2007 (Figure 37a). Elle est similaire à la
précédente avec un décrochage des niveaux simulés qui apparaît dès les premières
heures.
On observe, sur les deux graphes, que le niveau d’eau simulé dans le drain karstique
tend sur le long terme vers les niveaux mesurés dans la vasque de la source du Lez.
Ceci est satisfaisant. Cependant, au début des deux modélisations, on observe une
baisse plus précoce des niveaux d’eau simulés. Ceci est dû à un décrochage des
niveaux au moment où le débit naturalisé Qnat devient inférieur au débit de pompage,
c’est à ce moment précis que les niveaux simulés baissent.
Cet écart entre le modèle et les observations est lié à une représentation incertaine du
débit naturalisé. Cependant, dans un tel hydrosystème exploité en permanence, il nous
est impossible d’avoir précisément accès au débit naturel du système. Dès lors, il est
nécessaire d’opérer un choix sur la période qu’il nous semble préférable de privilégier
lors du calage. Nous avons choisi de calibrer le modèle sur le long terme de façon à
mieux représenter les niveaux les plus bas. De plus, il nous semble logique de
privilégier la calibration de ce modèle de simulation de la sollicitation des réserves sur
la période où cette sollicitation est la plus forte, à savoir en fin de cycle de pompage au
moment où la sollicitation des réserves devient nettement supérieure au débit
naturalisé sur lequel pèse une incertitude.
Figure 36 : Evolution des niveaux d’eau observés et simulés pendant le cycle de pompage de
2001. (a) Les débits naturalisés, de pompages et de sollicitation des réserves sont également
représentés. (b) L’erreur relative sur le rabattement simulé est également représentée.
Figure 37 : Evolution des niveaux d’eau observés et simulés pendant le cycle de pompage de
2001. (a) Les débits naturalisés, de pompages et de sollicitation des réserves sont également
représentés. (b) L’erreur relative sur le rabattement simulé est également représentée.
Les paramètres ajustés avec ce modèle sont les suivants (Tableau 7).
La transmissivité obtenue pour les blocs, Tb , oscille entre 7 et 8 x 10-2 m2/s, valeur très
élevée. Cependant, dans ce modèle à écoulement radial vers les drains, cette valeur
de transmissivité élevée intègre la longueur du réseau de drains sollicité. Dès lors, la
perméabilité correspondante est sans doute 2 à 3 ordres de grandeur inférieure (10-3 à
10-4 m/s).
Le rayon moyen des conduits constituant le réseau de drains karstiques est inconnu.
Nous l’avons arbitrairement fixé à r0 = 1 m. Dans l’équation (6), ce terme intervient
sous un logarithme, tout comme l’emmagasinement, son influence sur le résultat est
donc faible.
Figure 39 : Modèle conceptuel de l’aquifère karstique du Lez et paramètres utilisés pour le calage du modèle
Les résultats de la phase de validation durant les périodes 2001 et 2002 sont
présentés respectivement aux Figure 40 et Figure 41. On constate que les niveaux
calculés évoluent très similairement aux niveaux mesurés, avec respectivement des
erreurs relatives moyennes de 9 et 15 % pour les rabattements significatifs. Ceci
permet de valider le modèle mathématique dans la gamme des niveaux d’eau
observés, c’est-à-dire entre 65 et 45 m NGF. L’utilisation du modèle pour des niveaux
d’eau inférieurs à 45 m NGF implique l’hypothèse que le milieu se comporte
similairement d’un point de vue hydrodynamique dans ces conditions de basses eaux.
Pour mémoire, un modèle à écoulements linaires vers une structure drainante de type
tranchée a également été testé. Contrairement à la solution radiale, il ne permet pas de
simuler de façon aussi satisfaisante les rabattements observés en maintenant des
paramètres de calage dans une faible gamme de variation d’un cycle hydrologique à
un autre. La solution radiale a donc été retenue.
Figure 40 : Evolution des niveaux d’eau observés et simulés pendant le cycle de pompage de
2001. (a) Les débits naturalisés, de pompages et de sollicitation des réserves sont également
représentés. (b) L’erreur relative sur le rabattement simulé est également représentée.
Figure 41 : Evolution des niveaux d’eau observés et simulés pendant le cycle de pompage de
2001. (a) Les débits naturalisés, de pompages et de sollicitation des réserves sont également
représentés. (b) L’erreur relative sur le rabattement simulé est également représentée.
Il faut noter qu’à cette période, les débits de pompages horaires n’étaient pas
disponibles, ce qui explique que les fluctuations journalières de niveau d’eau ne
peuvent être simulées par le modèle.
Afin d’observer le comportement du modèle dans des conditions de très basses eaux
(<45 m NGF), nous avons exploré les années durant lesquelles les niveaux d’eau
mesurés à la vasque ont baissé en-dessous de ce niveau. Malheureusement, chacune
des récessions observées jusqu’à ce niveau a été perturbée à un moment ou un autre
par un épisode de pluie (pour rappel, le modèle POKA n’est pas adapté à la prise en
compte des phases de recharge de l’aquifère) ou par une baisse du débit naturalisé
calculé non synchronisée avec la baisse de niveau observé à la vasque. Cependant,
l’année 1998 est sans aucun doute la moins perturbée par des pluies et comporte une
bonne synchronisation des signaux. Même si ces conditions diffèrent des conditions
optimales d’utilisation du modèle, celui-ci a été testé sur cette période (Figure 42).
L’erreur relative moyenne de cette simulation est de 7 % pour les rabattements
significatifs (niveau d’eau inférieur à 50 m NGF). Sur la fin de la simulation, pour des
niveaux proches de 40 m NGF, l’erreur relative est de l’ordre de 10 %. Sur l’ensemble
des années testées (1998-2001, 2002, 2006 et 2007), l’erreur relative en fin de
simulation pour des rabattements significatifs est systématiquement inférieure ou égale
à 10 %. Nous considérons qu’il est raisonnable de retenir cette valeur de 10 % d’erreur
relative sur les rabattements simulés (soit 2,5 m pour un rabattement de 25 m, c’est-à-
dire un niveau dans la vasque à 40 m NGF) pour des niveaux d’eau compris entre 40
et 65 m NGF. Au-delà, c’est-à-dire pour des niveaux en-dessous de 40 m NGF, les
résultats fournis par le modèle devront être considérés avec prudence car en-dehors
de la gamme de calibration et validation du modèle.
Figure 42 : Evolution des niveaux d’eau observés et simulés pendant le cycle de pompage de
1998. (a) Les débits naturalisés, de pompages et de sollicitation des réserves sont également
représentés. (b) L’erreur relative sur le rabattement simulé est également représentée.
4.4. SYNTHESE
Le modèle hybride distribué est réalisé avec le logiciel FEFLOW (Finite Element
subsurface FLOW) [Diersch, 1998a]. Le réseau de drainage karstique est représenté
par des éléments discrets monodimensionnels qui sont couplés à une matrice
tridimensionnelle représentative des blocs fracturés en supposant que les charges
hydrauliques sont continues. Seul le réseau principal de drainage (et d’exutoires) a été
pris en compte.
Ce modèle a été développé dans le cadre de la thèse de Naomi Mazzilli (2011) portant
sur la sensibilité et l’incertitude de modélisation sur les bassins versants
méditerranéens à forte composante karstique.
Le modèle couvre une étendue de 1110 km2 qui s'étend entre la faille des Cévennes
au Nord, le pli de Montpellier au Sud, et les rivières de l’Hérault et du Vidourle à
l'Ouest et à l'Est respectivement. La modélisation géologique 3-D a été réalisée à l’aide
du logiciel GMS [Owen et al.1996]. Du sommet à la base des séries, les unités
géologiques suivantes sont représentées dans le modèle géologique 3-D (Figure 43):
- l'aquitard superficiel qui comprend l'Hauterivien et les unités du Quaternaire. Les
unités de l’Hauterivien supérieur et du Lutétien qui peuvent être considérées comme
déconnectées de l’aquifère du Lez ne sont pas prises en compte dans cette
analyse ;
- l'aquifère perché qui comprend les unités du Valanginien supérieur jusqu’aux
formations de l’Hauterivien du synclinal de Sauteyrargues. Les aquifères du causse
de l’Hortus et du bassin de Liouc ne sont pas considérés ;
- le toit de l’aquifère qui correspond à l'unité marneuse du Valanginien inférieur ;
- l'aquifère du Lez comprend les unités du Jurassique supérieur et du Berriasien
inférieur. Les unités du Jurassique moyen et inférieur (partiellement séparées de
l'unité Jurassique supérieur par la série de l’oxfordo-Callovien) sont négligées. A
noter que le synclinal de Sauteyrargues, qui correspond à l’aquifère perché, est relié
à l’aquifère du Lez par l’intermédiaire de la faille de Corconne.
Figure 43 : Géométrie du modèle: a) vue en plan, b-c) vue en coupe des unités
hydrogéologiques.
Ecoulement dans les conduits. L'écoulement dans les conduits karstiques est
modélisé à l'aide de l’onde diffusive, issue d’une approximation des équations de
Saint-Venant. Les équations de Saint-Venant sont des équations différentielles
partielles, fondées sur la conservation de la masse et du moment ; elles sont
couramment utilisées pour décrire les écoulements à surface libre dans des canaux
découverts [Chow, 1959]. Pour des régimes d’écoulement faibles à modérés, les
termes inertiels peuvent être négligés par rapport aux termes gravitationnels, à la
friction et la pression. Par ailleurs, la viscosité et le cisaillement à l'interface air/fluide
peuvent être négligés par rapport aux contraintes de cisaillement à l'interface
roche/fluide.
Dans le modèle, les éléments discrets représentatifs des conduits karstiques sont
considérés comme totalement saturés ; ils sont désactivés dès qu’un nœud est
asséché (charge hydraulique plus faible que la cote du nœud).
Figure 44 : Configuration du Modèle, étendue du modèle, conditions aux limites (voir aussi
Tableau 9) et réseau de conduits karstiques modélisés.
.
Les conduits karstiques modélisés, ainsi que les différents exutoires de l’hydrosystème
sont présentés dans la Figure 44. Les éléments représentatifs des conduits karstiques
sont situés en profondeur dans l’unité "aquifère du Lez", à l'interface entre la couche
"épikarst" et la couche "aquifère profond".
Le drainage au travers du pli de Montpellier est modélisé par une condition de flux (voir
Figure 44). Les flux sortants associés à cette condition à la limite de transfert sont
calibrés durant la phase de simulations en régime permanent.
Figure 45 : Géométrie du modèle: aire de recharge (jaune) et frontières imperméables (en gris).
Tableau 8 : Conditions aux limites et contraintes de flux associées aux exutoires karstiques. La
charge imposée est liée au niveau de débordement de la source. La contrainte sur le flux
dépend des pompages effectués à l’exutoire ou du débit moyen interannuel.
Figure 46 : Raffinement de maillage autour du réseau de conduits. Les éléments ont une taille
comprise entre 30 mètres au voisisnage du conduit karstique jusqu à une taille moyenne de
l’ordre de 300 mètres.
Conditions aux limites en régime permanent : les conditions aux limites pour la
modélisation de l'hydrosystème karstique sont résumées dans le Tableau 8 et la Figure
44. Les paramètres (flux sortant) associés à la condition de transfert au niveau du pli
de Montpellier sont déduits de la phase de calibration.
1 - Système Lez 4
2 - Système Lirou 3
3 - Système Fontbonne 0.5
4 - Système Montlobre 1
5 - Système Fontanilles 1
6 - Système Vernade 1
7 - Système Sauve 1
1. Le fait que le modèle sous-estime les débits de pointe à la source du Lez, du Lirou
et probablement aussi aux sources de Sauve et des Fontanilles suggère que :
b. soit la section des conduits karstiques est trop faible ce qui a pour
conséquence de limiter les flux. Néanmoins, le comportement hydrodynamique
cohérent de la source de Fontbonne suggère que la section des conduits
karstiques contrôle convenablement les écoulements vers la source.
2. Le fait que les périodes de débordement soient trop longues et que le débit moyen
de débordement à la source du Lez soit trop élevé peut être expliqué par :
b. une limitation des débits de pointe liée à une trop faible section des conduits
karstiques,
5.4. CONCLUSION
2. bien que tous les exutoires et conduits karstiques connus aient été considérés pour
cette modélisation hybride distribuée, le réseau de drainage modélisé ne comporte
que les exutoires et conduits karstiques principaux sans les différents réseaux de
drainage associés ; on s'attend par conséquent à des flux sortants au niveau de
ces exutoires principaux, supérieurs à la réalité. Ce biais sera notamment important
durant les crues, lorsque nombre de sources temporaires deviennent actives. Par
exemple, comme le modèle ne tient pas compte des différents exutoires
temporaires le long du fleuve Hérault, cela devrait se traduire par des débits
supérieurs à la réalité au niveau des exutoires considérés par le modèle (sources
de Vernède, Uglas, et Fontanilles) ;
Une étude de sensibilité menée dans le cadre de la thèse de Naomi Mazzilli (2011) a
montré que l’influence des propriétés hydrodynamiques de la couche "épikarst" sur les
débits simulés aux exutoires est notablement moins importante que la section
transversale des conduits karstiques. Les différences de comportement entre
Ainsi, alors que le pli de Montpellier est traditionnellement considéré comme étant une
limite imperméable ou à flux très limité, ce modèle montre qu'il pourrait être associé à
une condition de transfert, notamment dans la partie ouest du pli. Cette limite pourrait
donc être le siège de circulations d’eau conséquentes vers les formations aquifères
avales (Causse d’Aumelas, Issanka). Les flux considérés dans cette modélisation (0.4
m3/s environ) pourraient être sous-estimés par rapport à la réalité.
6.1. RESUME
Cette partie fait l’objet d’une réflexion menée au-delà du cadre prévu par le projet Lez.
Le travail sur les Réseaux de Neurones a été développé par Line Kong A Siou dans le
cadre de sa thèse de doctorat (école doctorale SIBAGHE, financement Ecole des
Mines d'Alès et encadrement de thèse assuré par l’EMA et HSM).
Le lecteur se réfèrera au manuscrit de thèse de Line Kong A Siou pour de plus amples
descriptions et expérimentations, qui vont bien au-delà de l’application traitée.
Dans cette partie, il est proposé de mener à bien une comparaison de deux approches
pour la modélisation du système karstique du Lez : le modèle de type modèle à
réservoir VENSIM et un modèle par réseau de neurones.
L’objectif de ces simulations est de comparer, d’une part, le fonctionnement des deux
modèles, et d’autre part, leur sensibilité aux variables d’entrée. Pour cela, trois
chroniques de pluie sont considérées successivement : la pluie brute, la pluie efficace,
et la pluie brute à laquelle s’ajoute une donnée « d’évapotranspiration présupposée »
qui se présente sous la forme d’une courbe gaussienne ou sinusoïdale, maximale
lorsque l’évapotranspiration est la plus intense. Nous pourrons donc évaluer la
capacité des modèles à exploiter la pluie brute, la pluie efficace et une information
« qualitative » sur l’évapotranspiration. Pour rendre les simulations comparables pour
les deux modèles, les données d’entrée et de sortie, les ensembles de test et de
calage/apprentissage sont identiques.
Le modèle par réseaux de neurones est présenté dans un premier paragraphe, ainsi
que ses propriétés principales. Le second paragraphe est consacré à la simulation du
débit et le troisième à la simulation de la piézométrie.
a) Neurone formel
Tout d’abord, le neurone calcule son potentiel : la somme pondérée de ses entrées.
Les paramètres, coefficients de pondération, ou coefficients synaptiques, sont calculés
durant la phase d’apprentissage.
Le potentiel du neurone subit ensuite une transformation, opérée par une fonction f,
appelée fonction d’évaluation ou fonction d’activation du neurone. Plusieurs fonctions
d’activation existent : la fonction d’activation de Heaviside est un seuil binaire, par
exemple à valeurs 0 ou 1, selon que le potentiel du neurone est positif ou négatif. La
fonction d’activation sigmoïde est en général une tangente hyperbolique, présentant
l’avantage d’être continue et dérivable. Enfin la fonction d’activation peut être la
fonction identité ; dans ce cas le neurone est linéaire.
b) Réseau de neurones
Le réseau de neurones utilisé dans cette étude est un perceptron multicouches, choisi
pour ses propriétés d’approximation universelle et de parcimonie vis-à-vis d’autres
modèles statistiques non-linéaire.
Ce réseau particulier est présenté en Figure 54, les neurones cachés sont non-
linéaires alors que le neurone de sortie est linéaire. Tous les neurones d’une couche
sont reliés à tous les neurones de la couche suivante.
estimés du jour k-1 pour estimer les débits du jour k. Les débits estimés q(k-1)sont tout
simplement ceux qui sont estimés lors du calcul précédent du réseau.
Pour le débit calculé en sortie du modèle, le débit restitué de 160 L/s à la source est
systématiquement ajouté au débit simulé par le réseau lorsque celui-ci est inférieur à
160 L/s.
Le graphique suivant présente les résultats des simulations effectuées avec VENSIM
et le réseau de neurones, prenant tous deux comme variable d’entrée la pluie
pondérée. Cette variable est également représentée sur le graphique en histogramme
inversé.
Figure 56 : Hydrogrammes simulés à partir de la pluie pondérée calculée par TEMPO sur
l’ensemble du test. La courbe noire correspond au débit de débordement mesuré à la source du
Lez, la courbe verte en tirets est l'hydrogramme simulé par le réseau de neurones, et la courbe
orange en pointillés est l'hydrogramme simulé par le modèle VENSIM. Les barres inversées
représentent le cumul de pluie journalier (pluie pondérée).
Pour la crue principale de décembre 2003, les deux modèles ont un décalage de un
jour par rapport au pic de crue observé. Le modèle VENSIM a une amplitude plus
proche du débit observé (le pic simulé atteint 97% du pic mesuré) que le modèle par
réseau de neurones qui sous-estime la crue (81% du pic observé). Enfin, les crues de
printemps sont sous-estimées par les deux modèles.
Les coefficients de Nash sont de 0,75 pour le modèle à réservoir et de 0,79 pour le
modèle par réseaux de neurones. On peut remarquer que les valeurs de ces
coefficients traduisent assez mal la comparaison entre les deux modèles que nous
venons de proposer et qui était plutôt en faveur d’une meilleure modélisation par le
modèle VENSIM.
Le graphique permet de constater que les données du pompage sont bien prises en
compte par les deux modèles. En effet, une pluie au mois de Juin fait réagir les deux
modèles mais de façon bien moins importante qu’une pluie similaire survenant au mois
de septembre. Il se peut donc que cette réaction soit uniquement due au manque
d’information sur l’évapotranspiration.
A la fin de l’été, on constate une différence de comportement pour les deux modèles :
le modèle VENSIM ne réagit pas aux pluies du mois d’août alors que le réseau de
neurones continue de réagir. Cela pourrait suggérer que le réseau de neurones serait
plus sensible à l’absence d’information d’évapotranspiration que VENSIM.
Concernant les périodes pluvieuses, les deux modèles ont un comportement proche
sauf pour la première crue d’automne et pour la plus forte crue qui survient au mois de
décembre. La première crue d’automne est particulière car le karst est très déprimé à
cette époque par le pompage, l’eau doit donc d’abord remplir le réservoir avant que le
débit ne réagisse en conséquence. De ce fait, le débit observé forme d’abord un petit
pic de crue (< 200 000 m3/jour) avant de passer à une augmentation plus importante
(environ 400 000 m3/jour). Ce comportement est mal représenté par VENSIM qui ne
réagit pas du tout à la première pluie et atténue fortement le pic suivant. Le modèle par
réseaux de neurones réagit de façon plus satisfaisante même s’il surestime le premier
pic et sous-estime le second.
Les modèles utilisés sont les mêmes que précédemment à savoir VENSIM et le
modèle neuronal, mais les calages et apprentissages ont bien entendu été effectués
avec la chronique de pluie efficace. La sélection de la complexité du réseau de
neurones avec la chronique reff en entrée a donné un modèle plus parcimonieux que
pour la pluie brute : seulement quatre jours de pluie et deux neurones cachés. Les
performances sont pourtant équivalentes en termes de critère de Nash qui est égal à
0,78. Ce critère ne traduit pas l’amélioration apportée sur la simulation du pic de crue
principal ainsi qu’à l’étiage.
Le pic de crue principal est toujours décalé d’un jour pour les deux modèles, mais
l’amplitude est mieux représentée (95% du pic pour VENSIM et 101% pour le réseau
de neurones qui estime donc le niveau parfaitement). C’est un point d’amélioration
pour le modèle neuronal qui avait une amplitude trop faible dans la simulation
précédente. De plus, les crues de printemps sont mieux simulées que précédemment.
Sur ce point les deux modèles ont été améliorés. Cependant le modèle VENSIM
simule des décrues un peu trop lentes, tandis que le modèle neuronal au contraire les
décrues sont un peu rapides par rapport au débit observé.
Figure 57 : Hydrogrammes simulés à partir de la pluie efficace calculée par TEMPO sur
l’ensemble de test. La courbe noire correspond au débit de débordement mesuré à la source du
Lez, la courbe verte en tirets est l'hydrogramme simulé par le réseau de neurones, et la courbe
orange en pointillés est l'hydrogramme simulé par le modèle VENSIM. Les barres inversées
représentent le cumul de pluie journalier (pluie efficace).
Concernant la période estivale, on constate que les pics de crue du mois d’août ont
quasiment disparu, l’évapotranspiration ayant enlevé la plupart des jours de pluie à
cette période, en revanche le pic du mois de Juin est encore présent pour les deux
modèles, bien que l’intensité de la pluie correspondante ait été diminuée par
l’évapotranspiration. Ceci infirme donc l’hypothèse formulée au paragraphe précédent
qui supposait que la prise en compte de l’évapotranspiration permettrait de supprimer
le pic de crue indésirable de la fin du mois de Juin. Il reste à vérifier si
l’évapotranspiration « artificielle » n’est pas plus efficace dans ce cas, ce qui sera
discuté au paragraphe suivant.
pour que le modèle ne réagisse. Encore une fois, l’évapotranspiration artificielle pourra
éventuellement apporter une amélioration.
Le coefficient de Nash pour le modèle VENSIM est de 0,55 alors que celui pour le
réseau de neurone est de 0,78. La valeur particulièrement basse du coefficient de
Nash pour VENSIM est essentiellement due au manque de réaction à l’automne.
Les simulations pour le modèle à réservoir ont cette fois été effectuées avec le modèle
VENSIM qui utilise en entrée une chronique de pluie efficace «artificielle» décrite au
chapitre 3. Le réseau de neurones possède une entrée ETP artificielle distincte de
l’entrée de la pluie ce qui lui confère plus de souplesse quand à l’utilisation de l’ETP
artificielle. Ces simulations permettront d’évaluer si l’ajout d’une ETP « artificielle »,
donc disponible pour n’importe quel système hydrologique améliore le résultat
comparativement aux simulations faites avec la seule pluie brute pondérée.
Dans un deuxième temps, seront observées les différences entre ces simulations et
celles effectuées avec la pluie efficace « réelle », issue de l’estimation de
l’évapotranspiration sur le bassin d’alimentation de la source.
Concernant l’étiage, on peut observer que le modèle VENSIM ne réagit plus du tout
aux pluies d’été ce qui constitue une amélioration par rapport aux modèles utilisant reff
et rpond. En revanche, le réseau de neurones réagit à chaque pluie, comme pour la
simulation avec rpond.
Les épisodes d’automne ne sont pas mieux simulés qu’avec rpond et reff pour le modèle
par réseau de neurones qui donne une simulation très proche de celle qu’il a fournie
avec rpond. Le modèle VENSIM réagit un peu plus rapidement à l’automne que le
modèle VENSIM, utilisant la reff , mais cette réaction est encore trop tardive puisque le
premier grand pic d’automne n’est pas du tout représenté.
Enfin, concernant les épisodes de printemps, ils sont correctement simulés par
VENSIM alors que le modèle par réseaux de neurones sous-estime tous les pics.
VENSIM 97 % 95 % 50 %
RNF 81 % 101 % 81 %
Tableau 16 : Pourcentage du pic principal de débit obtenu pour la crue de décembre 2003 à la
source du Lez par le modèle VENSIM et par le modèle par réseaux de neurones, en fonction de
la chronique de pluie utilisée.
Ces chiffres, comme les observations notées dans les paragraphes précédents,
traduisent bien la sensibilité des modèles aux variables d’entrée : le modèle par
réseaux de neurones est très peu affecté par les changements de variable d’entrée.
Son coefficient de Nash est quasiment constant alors que le pourcentage du pic estimé
varie de 20%, avec une excellente simulation pour le modèle utilisant la pluie efficace.
Il semble donc que le modèle neuronal ne trouve pas l’information qui lui manque sur
l’évapotranspiration quelles que soient les variables proposées. Ceci a l’avantage de
fournir un modèle robuste, peu dépendant des données d’entrées. En contrepartie il
semble difficile dans ces conditions d’obtenir de meilleures simulations de l’étiage et
des premières crues d’automne. Il faut cependant remarquer que ces simulations sont
très satisfaisantes. On doit noter que ce résultat est tout à fait remarquable pour un
modèle récurrent n’utilisant pas le débit précédent observé (avec un modèle dirigé,
c'est-à-dire utilisant en entrée les débits précédents observés, le réseau de neurones
effectue une simulation avec un critère de Nash de 0,95 et 87 % du pic retrouvé). Par
conséquent, le modèle neuronal mis en œuvre apparaît comme un outil fiable, et peu
sensible aux variables d’entrées ce qui est très attractif dans un domaine où les
données sont souvent très coûteuses. Une fois l’apprentissage effectué, le modèle
peut simuler efficacement le débit à la source du Lez en utilisant uniquement la pluie
brute.
Dans tous les cas, l’information sur le remplissage du karst semble indispensable pour
simuler correctement les premières crues d’automne. Au vu des résultats qui viennent
d’être présentés, cette information ne peut venir des estimations actuelles de
l’évapotranspiration.
L’objectif de ces simulations est d’évaluer la capacité des modèles à représenter l’état
de remplissage du karst, ici représenté par la piézométrie. Les variables d’entrée sont
les mêmes que précédemment : pluviométrie, débit pompé et éventuellement
évapotranspiration. Si les modèles parviennent à faire des simulations correctes de la
piézométrie, il sera alors envisageable d’inclure la sortie de ce modèle au modèle de
simulation du débit, permettant à ce dernier d’avoir une information sur le remplissage
du karst. Dans le cas contraire, on déduira que ces modèles n’ont pas été capables de
représenter l’état de remplissage du karst et ne pourront donc pas être utilisés pour
améliorer les modèles de simulation du débit. Dans ce paragraphe, une attention
particulière sera donnée à la qualité des simulations lors des premières pluies
d’automne, période qui semble la plus problématique d’après les résultats des
simulations du débit.
La démarche suivie est la même que pour le débit : en premier lieu c’est la pluie brute
pondérée qui est utilisée comme variable d’entrée, puis la pluie efficace calculée grâce
à l’ETP estimée à Mauguio et enfin la pluie brute pondérée à laquelle on ajoute une
évapotranspiration artificielle.
Figure 59 : Débit pompé et piézométrie à la source du Lez pour l'année civile 2003. La courbe
verte correspond à la piézométrie dans le drain principal et la courbe bleue correspond au débit
journalier pompé.
Les simulations ont été effectuées avec VENSIM et avec le modèle par réseaux de
neurones, en effectuant un nouveau calage/apprentissage de façon à simuler la
piézométrie. Le nombre d’entrées et de neurones cachés pour chaque réseau est
reporté dans le Tableau 18 et sera discuté par la suite.
Figure 60 : Piézométrie simulée à partir de la pluie pondérée sur l’année test 2003. La courbe
noire correspond à la piézométrie mesurée dans le drain principal de la source du Lez, la
courbe verte en tirets est la piézométrie simulée par le réseau de neurones, et la courbe orange
en pointillés est celle simulée par le modèle VENSIM. Les barres inversées représentent le
cumul de pluie journalier (pluie pondérée).
La piézométrie constante à 65 m NGF est mal représentée par les deux modèles qui
continuent de réagir aux pluies, surtout pour le modèle par réseaux de neurones qui
simule de grandes variations de piézométrie lors de l’épisode pluvieux du mois de
décembre. La variation par rapport au niveau 65m NGF du pic principal est surestimée
par le réseau de neurones et représente plus de 500% de la variation observée (contre
150% pour le modèle VENSIM). Ce problème pourra éventuellement être réglé par
l’utilisation de la variable pluie efficace en entrée ou par l’ajout d’une entrée
évapotranspiration, au risque de compromettre la qualité du modèle en étiage qui est
bien simulé. A l’automne, une brutale remontée du niveau piézométrique est
observée : la niveau passe de 40 m NGF à 60 m NGF en une seule journée et retrouve
le niveau de 65m NGF en seulement deux jours. VENSIM simule bien cette remontée,
alors que le réseau de neurone simule une remontée beaucoup plus lente, plus de 30
jours sont nécessaires pour retrouver le niveau de 65m NGF. Par contre le modèle
VENSIM surestime le déficit de l’été sur le niveau piézométrique : le 8 septembre,
VENSIM prévoit un déficit pratiquement deux fois supérieur au déficit observé. Le
modèle par réseaux de neurones est meilleur sur ce point, prévoyant en moyenne un
déficit égal à celui observé. Les valeurs des critères de Nash sont de 0,59 pour
VENSIM et 0,81 pour le modèle par réseaux de neurones. Encore une fois ceci traduit
assez mal les observations effectuées : le modèle VENSIM paraît très mauvais à
cause de l’importance du déficit qu’il a simulé mais reste meilleur que le réseau de
neurones sur les points les plus importants de la simulation, à savoir le début de
l’automne et la grande crue du mois de décembre.
L’introduction de la pluie efficace en entrée des modèles est effectuée dans l’espoir de
simuler moins de réactions aux pluies que précédemment, la piézométrie étant souvent
constante ou connaissant de très faibles variations lors d’épisodes pluvieux importants.
Figure 61 : Piézométrie simulée à partir de la pluie efficace sur l’année test 2003. La courbe
noire correspond à la piézométrie mesurée dans le drain principal de la source du Lez, la
courbe verte en tirets est la piézométrie simulée par le réseau de neurones, et la courbe orange
en pointillés est celle simulée par le modèle VENSIM . Les barres inversées représentent le
cumul de pluie journalier (pluie efficace).
Le réseau de neurones donne une simulation peu différente de celle effectuée avec la
pluie brute pondérée. La remontée en automne est plus longue qu’avec la pluie brute,
et similaire à celle simulée par VENSIM. Une forte réaction est toujours observée au
moment de l’épisode de décembre 2003.
L’introduction de la pluie efficace n’a pas donné de bons résultats, pour le réseau de
neurones comme pour VENSIM, les résultats sont moins bons qu’avec la pluie
pondérée. Les coefficients de Nash (0,59 pour VENSIM et 0,81 pour le modèle
neuronal) sont d’ailleurs les plus bas obtenus de toutes les simulations (voir Tableau
17).
Sur ces dernières simulations, le réseau de neurones ne donne pas de résultats très
différents de ceux présentés précédemment. Le modèle VENSIM (simplifié de façon à
ne comporter qu'un seul coefficient de tarissement) a été amélioré par rapport aux
deux autres simulations au niveau de la baisse du niveau piézométrique en été et du
niveau le plus bas atteint qui était très surestimé lors de l’utilisation de la chronique
Ppond. Cependant, le point le plus important de l’ensemble du test n’a pas été amélioré :
la remontée rapide est meilleure que celle effectuée avec Peff mais moins bonne que
celle simulée avec Ppond.
L’ajout d’un signal d’évapotranspiration, qu’il soit estimé ou artificiel n’a pas amélioré
significativement les simulations, il paraît donc préférable d’utiliser uniquement la pluie
brute pondérée en entrée, qui donne de meilleurs résultats particulièrement à
l’automne.
Comme cela a été constaté pour le débit, le modèle par réseaux de neurones est peu
sensible aux variables d’entrée, ceci est confirmé par les valeurs des critères de Nash
données dans le Tableau 17. De même, le modèle VENSIM est plus sensible aux
données d’entrées que le modèle neuronal. Par contre il a été en mesure de donner de
meilleurs résultats sur deux points importants : la remontée brutale du niveau
piézométrique en automne, et des variations modérées au-dessus de 65 m NGF.
Tableau 17 : Coefficients de Nash obtenus par les deux modèles pour la simulation de la
piézométrie en fonction de la chronique de pluie utilisée.
Nombre d’entrées 14 4 9 12 8 6
Nombre de neurones 10 2 2 1 3 3
cachés
7. Modèle Hydrologique
7.1. RESUME
Un modèle Pluie-Débit mono réservoir (appelé 'ATHYS' car développé sous cette
plateforme de modélisation) a été développé pour reproduire la genèse et la
propagation des crues de surface du Lez à Lavalette.
- la base de données AT1 a été complétée et interprétée (sondes d’humidité, pluies
radar…) ;
- le modèle développé est robuste (Nash très satisfaisants et pics de crues très bien
reproduits), parcimonieux (4 paramètres de calage invariants pour le Lez et une
condition initiale qu’il convient de définir avant chaque crue), distribué (par la prise
en compte de la lame d’eau précipitée spatialisée) ;
- la condition initiale peut être pré-déterminée (avant le début de la crue) par le biais
d’indicateurs de l’état de saturation antérieur du système. Les indicateurs testés les
plus pertinents sont le niveau piézométrique et l’humidité des sols.
7.2.1. Principe
La feuille "GR4J" permet de faire des simulations de débit au pas de temps journalier à
l'aide du modèle GR4J (cf schéma du modèle ci-dessous). La version utilisée ici est
celle présentée par Perrin (2002) et Perrin et al. (2003) et développée sous Excel. On
se référera à ces articles pour une explication des équations représentées.8
E P
interception
En Pn
Es Ps Pn-Ps
Réservoir de
X1
production S
Perc Pr
0.9 0.1
HU1 HU2
X4 2.X4
Q9 Q1
Réservoir
X3 F(X2) F(X2)
de routage R
Qr Qd
Quatre critères d'efficacité sont calculés pour juger de la qualité des simulations :
- critère de Nash-Sutcliffe calculé sur les débits ;
- critère de Nash-Sutcliffe calculé sur les racines carrées des débits ;
- critère de Nash-Sutcliffe calculé sur les logarithmes des débits ;
8
Cemagref, UR Hydrosystèmes et Bioprocédés, Antony. Contacts:
charles.perrin@cemagref.fr, vazken.andreassian@cemagref.fr, version 1.11 – 11/2009
- critère de bilan.
Les paramètres de calage sont X1, X2, X3, X4 représentant respectivement la taille du
réservoir de production (au sein duquel se développe l’ETP), les interactions entre la
surface et le souterrain, la taille du réservoir de routage (représentant la composante
souterraine des écoulements), la dynamique temporelle des transferts.
Valeurs initiales
Taux de remplissage initial S0/x1 0,60
Taux de remplissage initial R0/x3 0,70
Période
Longueur de la période de mise en route (j) 365
Durée de la période test (j) 5479
Date de départ 01/09/1995
Date de fin 01/09/2010
200 0
180
160 50
140
120 100
pluie
100 Qobs
Qsim
80 150
60
40 200
20
0 250
1/9/94 0:00
10/12/94 0:00
20/3/95 0:00
28/6/95 0:00
6/10/95 0:00
14/1/96 0:00
23/4/96 0:00
1/8/96 0:00
9/11/96 0:00
17/2/97 0:00
28/5/97 0:00
5/9/97 0:00
14/12/97 0:00
24/3/98 0:00
2/7/98 0:00
10/10/98 0:00
18/1/99 0:00
28/4/99 0:00
6/8/99 0:00
14/11/99 0:00
22/2/00 0:00
1/6/00 0:00
9/9/00 0:00
18/12/00 0:00
28/3/01 0:00
6/7/01 0:00
14/10/01 0:00
22/1/02 0:00
2/5/02 0:00
10/8/02 0:00
18/11/02 0:00
26/2/03 0:00
6/6/03 0:00
14/9/03 0:00
23/12/03 0:00
1/4/04 0:00
10/7/04 0:00
18/10/04 0:00
26/1/05 0:00
6/5/05 0:00
14/8/05 0:00
22/11/05 0:00
2/3/06 0:00
10/6/06 0:00
18/9/06 0:00
27/12/06 0:00
6/4/07 0:00
15/7/07 0:00
23/10/07 0:00
31/1/08 0:00
10/5/08 0:00
18/8/08 0:00
26/11/08 0:00
6/3/09 0:00
14/6/09 0:00
22/9/09 0:00
31/12/09 0:00
10/4/10 0:00
Figure 65 : Résultats des simulations du débit du Lez à Lavalette avec le GR4J (Observations 19/7/10 0:00
Le critère de Nash associé à cette chronique est de 0.75 (pour une valeur maximale de
1), ce qui est un résultat correct.
Toutefois les périodes de fortes crues sont mal simulées (observer les maximums qui
ne sont systématiquement reproduits qu’à hauteur de la moitié de leur amplitude), ainsi
que de nombreuses petites crues pour lesquelles des écoulements sont observés à
Lavalette, mais ne sont pas reproduits par le modèle. Enfin, la précision sur les basses
eaux n’est pas fine, les prélèvements et restitutions ne sont pas pris en compte.
La dispersion entre les débits observés et simulés est d’autant plus forte que les débits
sont importants, comme en témoigne la figure ci-dessous.
25
y = 0,7435x + 0,2572
R2 = 0,756
20
15
10
0
-1 1 3 5 7 9 11 13 15 17 19 21 23 25
Débit observé (mm/j)
Figure 66 : Débit simulé par le GR4J calé en fonction du débit observé à Lavalette
Enfin le pas de temps journalier n’est pas adapté à la simulation des crues, comme en
témoigne la chronique de l’automne – hiver 2002, au cours de laquelle les débits de
pointe observés en crue ne sont pas expliqués par le modèle (Figure 67).
250 0
50
200
100
150
pluie
Qobs 150
Qsim
100
200
50
250
0 300
1/9/02 0:00
8/9/02 0:00
15/9/02 0:00
22/9/02 0:00
29/9/02 0:00
6/10/02 0:00
13/10/02 0:00
20/10/02 0:00
27/10/02 0:00
3/11/02 0:00
10/11/02 0:00
17/11/02 0:00
24/11/02 0:00
1/12/02 0:00
8/12/02 0:00
15/12/02 0:00
22/12/02 0:00
29/12/02 0:00
Figure 67 : Débit simulé par le GR4J calé (rouge) et débit observé à Lavalette (bleu) au cours
de l'automne - hiver 2002
7.2.3. Conclusion
Dans notre cas d’étude, le Lez est un bassin à forte composante karstique, qui
complexifie les échanges entre la surface et le souterrain, ainsi que les circulations
souterraines. La réduction de la représentation de son fonctionnement, souterrain et de
surface, à un modèle global avec 2 réservoirs et 4 paramètres n’est pas pleinement
satisfaisante. De plus, les conditions météorologiques fortement variables qui y sont
rencontrées (été chaud et sec, automne orageux) engendrent des réactions du bassin
elles aussi fortement variables et non linéaires. L’extrapolation d’un fonctionnement
« normal » à un fonctionnement de crue est alors hasardeuse. Enfin, la dynamique
rapide du bassin du Lez (temps de réponse inférieur à 6 heures en cas de forte crue)
n’est pas adaptée à une simulation journalière, les cumuls de pluie n’étant pas les
variables explicatives des débits (il faudrait aussi considérer les intensités).
Pour la suite de l’étude, nous nous concentrerons sur les épisodes de crues
importantes (> 40 m3/s) nécessitant un début de vigilance pour la sécurité des
personnes et des biens. Pour ce faire, nous allons nous concentrer sur les processus
dominants lors de la genèse et la propagation de ces crues et négliger les autres, on
parle alors de modèle perceptuel. Ce modèle sera donc nécessairement événementiel,
c'est-à-dire qu’il ne simulera que les événements de crue et fonctionnera avec des
données horaires.
Les caractéristiques des épisodes de crue à Lavalette étudiés sont présentées dans le
Tableau 19.
Début Fin QHp tr P Vr Cr
3 6 3 6 3
(m /s) (h) (10 m ) (10 m ) (ad.)
18/10/1994 06:00 26/10/1994 20:00 124 3 24.1 19.5 0.81
27/10/1994 06:00 11/11/1994 11:00 99.8 4 19.4 17.5 0.90
15/12/1995 06:00 25/12/1995 00:00 53.2 5 11.3 8.5 0.75
13/03/1996 06:00 22/03/1996 03:00 40.3 4 11.3 7.9 0.70
17/12/1996 06:00 28/12/1996 17:00 139 2 21.6 22.5 1.04
24/11/1997 06:00 04/12/1997 12:00 53.9 5 14.5 7.9 0.55
16/12/1997 06:00 27/12/1997 00:00 122 5 21.0 22.0 1.05
11/11/1999 06:00 23/11/1999 00:00 42.8 19 14.5 11.3 0.78
28/09/2000 06:00 03/10/2000 19:00 51.5 8 9.1 3.8 0.41
23/12/2000 06:00 31/12/2000 23:00 48.3 21 11.2 10.4 0.93
16/01/2001 06:00 24/01/2001 14:00 93.1 8 10.8 11.1 1.03
09/10/2001 06:00 14/10/2001 14:00 238 4 11.6 6.9 0.60
08/09/2002 06:00 12/09/2002 10:00 103 6 15.1 6.3 0.41
08/10/2002 06:00 14/10/2002 17:00 43.0 18 14.3 8.1 0.57
09/12/2002 06:00 21/12/2002 01:00 376 2 36.7 44.4 1.21
22/09/2003 06:00 25/09/2003 15:00 91.5 3 13.3 3.4 0.26
15/11/2003 06:00 20/11/2003 04:00 64.1 4 9.6 8.8 0.91
29/11/2003 06:00 10/12/2003 04:00 424 3 31.1 31.6 1.01
05/09/2005 06:00 07/09/2005 15:00 467 4 40.7 19.8 0.49
27/01/2006 06:00 07/02/2006 17:00 52.5 16 13.5 12.1 0.90
19/10/2008 06:00 25/10/2008 06:00 109 4 23.4 5.7 0.24
Tableau 19 : Caractéristiques des événements de crue étudiées : avec QHp: le débit du pic de
crue, tr le temps de réponse du bassin, P la lame d'eau moyenne précipitée sur le bassin, Vr le
volume ruissellé et Cr le coefficient de ruissellement.
7.3.2. Concept
Le bassin topographique du Lez a une petite taille (114 km2) et une réponse rapide
(env. 6h). L’exutoire retenu est la station de Lavalette en amont de Montpellier. Les
processus de production des écoulements y sont dominants dans la réponse du bassin
à une sollicitation pluvieuse [Refsgaard, 1997]. Un effort particulier a ainsi été fait pour
déterminer la fonction de production la plus adaptée au Lez. Dans cette partie, la
fonction de transfert se limitera donc a une fonction « Lag and Route Simple ».
Toutefois, le modèle de propagation en rivière sera étudié précisément par la suite
7.3.3. Equations
Une première analyse des données laisse suggérer une influence importante de l’état
hydrique sur la réponse du bassin, notamment en début d’épisode. Les modèles à
réservoir permettent d’introduire de façon explicite cette influence de l’état hydrique via
un ou plusieurs réservoirs. Certains modèles vont chercher à représenter la capacité
de stockage totale de l’hydrosystème et l’évolution de la teneur en eau absolue de
l’hydrosystème, alors que d’autres vont chercher à reproduire l’évolution du déficit
hydrique de l’hydrosystème. La figure ci-dessous illustre la différence entre ces deux
approches :
Déficit hydrique en
début d’épisode
Capacité de
stockage totale de
l’hydrosystème
Teneur en eau en
début d’épisode
S
S
h(t0)
ib(t)
ie(t)
id(t)
ib(t) – ie(t)
V(t)
ds . V(t)
Sa = S/3
(1-w) . ds . V(t)
V (t ) S a V Sa
ie (t ) ib (t ) 2 si V (t ) Sa
S S
ie (t ) 0 sinon
où ie(t) est l’intensité de pluie efficace (ruissellement direct), ib(t) est l’intensité de pluie
brute, V(t) est le niveau dans le réservoir « sol », Sa=S/3 est le seuil du réservoir
« sol » au-delà duquel le ruissellement se déclenche et S est la taille du réservoir sol
au-dessus du seuil de ruissellement. On s’aperçoit ici que le coefficient de
V (t ) S a V Sa
C (t ) 2
S S
L’évolution du niveau V(t) est régie par l’équation différentielle ordinaire (EDO)
suivante :
dV (t )
ib (t ) ie (t ) ds . V (t )
dt
où V(t) est le niveau dans le réservoir « sol » au temps t, ib(t) est l’intensité de pluie
brute à t, ie(t) est l’intensité de pluie efficace (ruissellement direct). Une vidange a été
ajoutée afin de permettre la décroissance du niveau dans le réservoir et donc celle du
coefficient de ruissellement au cours d’un épisode lorsqu’il ne pleut pas. Cette vidange
est proportionnelle au niveau V(t) dans le réservoir « sol ». Le coefficient de vidange ds
est supposé constant pour un bassin donné.
Afin de représenter les écoulements retardés id(t) dus au drainage des sols et à la
vidange de l’aquifère karstique, une partie w de la vidange du réservoir « sol » participe
à la crue de surface.
id (t ) w . ds .V (t )
it (t ) ie (t ) id (t )
ib(t)
ib(t)
0.2 S ds . P(t)
P(t)
ie(t) = C(t) . ib(t)
id(t)
ib(t) – ie(t)
stoc(t) ds . stoc(t)
(1-min(1,w/S)) . ds . stoc(t)
La fonction de production décrite ici (Figure 70) est dérivée de la fonction SCS
d’origine. Elle est représentée par un réservoir de cumul de pluie de capacité infinie qui
détermine la pluie efficace (ou ruissellement direct) [Gaume et al., 2004] :
ie (t ) C (t ) . ib (t )
où ib(t) est l’intensité de pluie brute, ie(t) est l’intensité de pluie efficace et C(t) est le
coefficient de ruissellement :
0, si P(t ) 0.2 S
C (t ) P(t ) 0.2 S P (t ) 0.2 S
2 , si P (t ) 0.2 S
P (t ) 0.8S P(t ) 0.8S
dP(t )
ib (t ) ds.P(t )
dt
où ds est le coefficient de vidange appliqué au cumul de pluie P(t). Pour tous les
événements, ce coefficient de vidange est constant et le cumul de pluie au début de
chaque événement est nul : P(0) = 0.
dstoc(t )
ib (t ) ie (t ) ds.stoc(t )
dt
w
id (t ) min 1, .ds.stoc(t )
S
Ainsi, plus S est grand (plus le bassin est sec), plus le rapport w/S est faible (moins la
récession de l’hydrogramme est soutenue).
Ce ruissellement retardé vient ensuite s’ajouter au ruissellement ie(t) déjà produit pour
représenter le ruissellement total it(t).
it (t ) ie (t ) id (t )
Qm
t0
t
T
m
t0 t0 + t
Tm
Figure 71 : schéma de la fonction de transfert « Lag and Route Simple »
L'écoulement total it(t) produit par chaque maille est ensuite acheminé à l’exutoire du
bassin par une fonction de transfert de type « lag and route » [Maidment 1992 [26-9];
Bouvier and Delclaux, 1996; Lhomme et al., 2004]. La contribution d’une maille
parvient à l'exutoire après un temps de propagation Tm, calculé à partir de la vitesse
de transfert V et de la longueur de la trajectoire séparant la maille de l’exutoire :
lk
Tm
k V
où lk est la distance parcourue par l’eau dans chacune des k mailles de la trajectoire.
Cette contribution est également amortie par application d’un temps de diffusion Km
proportionnel au temps de propagation Tm.
Km K 0 .Tm
La fonction de transfert est donc pilotée par deux paramètres : V pour ajuster le temps
de propagation Tm et K0 pour ajuster le temps de diffusion Km. L’hydrogramme
élémentaire produit par chaque maille à l’exutoire du bassin peut alors être
calculé par :
it (t 0 ) t (t 0 Tm )
Qm (t ) . exp .A
Km Km
où A désigne l’aire de la maille. La
contribution Qm(t) de chaque maille k est ensuite additionnée pour donner
l’hydrogramme complet de crue.
a) Modèle SMA-SCS
On constate que le niveau initial V0 dans le réservoir est faible pour les épisodes de
début d’automne alors qu’il est élevé pour les épisodes de fin d’automne. La valeur de
w, qui varie d’un épisode sur l’autre, semble en rapport avec l’état hydrique initial du
bassin. Lorsque V0 est élevé (V0 = 400 mm) la valeur de w est forte (>0.5). Ce
comportement du modèle semble dû à la saisonnalité observée lors de l’analyse des
courbes de récession (après les premières pluies d’automne, un bassin humide
(représenté par un V0 élevé dans le modèle) conduit à des décrues plus soutenues
Du ruissellement direct est alors généré dès le début de l’épisode et s’ajoute aux
écoulements retardés générés par la vidange du réservoir. Sur l’hydrogramme
observé, le ruissellement direct ne commence apparemment pas dès le début de
l’épisode. Celui-ci se déclenche après une certaine quantité de pluie. Cette particularité
se reproduit sur d’autres crues observées. C'est la raison de la performance
particulièrement mauvaise du modèle pour la crue de décembre 2003 (Nash=0.02),
puisque durant toute la période précédant la montée de la crue, du ruissellement est
simulé (car le Vo est très important et supérieur au seuil d'apparition du ruissellement)
alors qu'en réalité, le ruissellement de surface n'a pas commencé.
Toutefois, alors que la valeur de ds = 0.28 j-1 calée à partir des courbes de récession
des hydrogrammes observés permet de bien soutenir les écoulements retardés, elle
entraîne pour les épisodes de fin d’automne une vidange trop importante du réservoir
« sol » de la fonction de production. Le second calage (ds= 0.01j-1) permet au modèle
de simuler un niveau dans le réservoir de production plus proche de la dynamique de
la piézométrie. Néanmoins, la vidange étant bien plus faible, les écoulements retardés
sont aussi fortement diminués. Le modèle ne reproduit plus de façon satisfaisante les
écoulements retardés observés.
On retrouve dans ce double calage possible, les 2 dynamiques que l’on retrouvait dans
les modèles de l’aquifère karstique et de la production des débits à la source.
Toutefois, l’objet ici, n’étant pas de reproduire correctement l’intégralité de
l’écoulement, mais bien sa part en crue, un autre modèle mono-réservoir va être testé.
Ceci permettra de conserver le caractère parcimonieux de ce modèle de genèse des
crues de surface et donc d’envisager des applications pour la prévision des crues. On
suppose que les écoulements de crue sont expliqués par un seul réservoir et se
La taille S du réservoir « sol » représentant le déficit hydrique initial est importante (S >
200 mm) pour les épisodes de septembre suivant un étiage prolongé (lorsque le bassin
est sec) et elle est plus petite (S < 200 mm) pour les épisodes intervenant après les
premières pluies de l’automne comme octobre 2001, décembre 2002 ou décembre
2003 (lorsque le bassin est plus humide). L’aménagement min(1,w/S) réalisé sur la
participation de la vidange du réservoir « sol » aux écoulements de surface permet à la
vidange de varier d’un épisode sur l’autre. Ainsi, la vidange participe à 40 % pour les
épisodes de septembre pour lesquels les écoulements retardés observés sont peu
soutenus et elle participe à 100 % pour les épisodes de décembre pour lesquels les
écoulements retardés sont plus soutenus. Le modèle reflète ainsi la variation
saisonnière observée sur le soutien des écoulements retardés. Les critères de Nash de
ce modèle sont globalement meilleurs que ceux du modèle précédent.
Figure 74 : Hydrogrammes observé (en tirets bleus) et simulé avec la fonction de production
SMA (en trait continu gris) ou SCS aménagée (en trait continu noir) pour une des crues testées.
La zone entourée repère les périodes en début d’épisode où le SMA-SCS produisait du
ruissellement direct alors que ce dernier n’a pas commencé en réalité.
Cette fonction de production semble donc résoudre les problèmes rencontrés avec la
fonction de production SMA-SCS. Un réservoir de production initialement vide autorise
des pertes systématiques au début de tous les épisodes quel que soit l’état hydrique
initial et empêche le déclenchement prématuré du ruissellement direct. Avec cette
fonction de production, une valeur du coefficient de vidange ds permet à la fois un bon
soutien des écoulements retardés et une bonne représentation de l’évolution de l’état
7.3.5. Sensibilité
Une étude de sensibilité des débits simulés aux différents paramètres du modèle a été
effectuée. Elle permet de qualifier et de quantifier l’effet de chacun des paramètres sur
l’hydrogramme de crue simulé à l’exutoire du bassin.
Cette étude de sensibilité a été effectuée sur l’épisode d’octobre 2001, représentatif
des épisodes de crue importants se produisant sur le bassin.
Figure 76 : Comparaison de l’hydrogramme observé (en tirets bleu) avec celui de la simulation
« de référence » (en noir) pour l’épisode d’octobre 2001
Une première simulation dite « de référence » (courbe noire) est effectuée avec les
valeurs calibrées des paramètres w, ds, V et K pour l’ensemble des épisodes) et la
valeur de S calibrée pour l’épisode d’octobre 2001. Malgré un démarrage tardif de la
montée de crue simulée, cette simulation représente de façon satisfaisante les débits
observés à l’exutoire avec un Nash de 0.94.
Les débits du pic de crue sont surtout sensibles aux paramètres S et V. Une
calibration sur les débits élevés (> 50 m3/s) permettra de s’affranchir de
l’influence des paramètres w, ds et K0 et d’avoir une première estimation de S et
V.
Le bassin versant a été découpé en deux parties : les causses en altitude pour
lesquels la partie karstique est à une altitude de 201 à 623 m et les garrigues pour
lesquelles la partie disposant d’une couverture végétale est à une altitude de -2 à 200
m.
Figure 77 : Distribution de l’occupation des sols du bassin versant du Lez. Le karst est en rose
et la couverture en bordeaux.
Méthode 1 : Le fonctionnement du karst est déterminé par l’état de l’aquifère. Cet état
est exprimé par la piézomètre à la source du Lez qui déborde à 65 m NGF. Le karst est
supposé pouvoir stocker 95mm de pluie entre 46 et 65m NGF (cf. étude abaques dans
le livrable L4.2). Le S karst est décrit donc par :
et la valeur de S couv est ajustée pour respecter S avg = S opt crue par crue.
Méthode 2 : La valeur de S karst est fixée au maximum de tous les S opt de l’étude
précédente, et la valeur de S couv est ajustée pour respecter S avg = S opt crue par
crue.
Méthode 3 : Enfin, la seule contrainte retenue est S avg = S opt et les S couv et S karst
sont optimisés.
S karst fixé par débit source S karst fixé à 405 mm Savg fixé à S opt
date de l'épisode S optimal Nash opt piezo source S karst S couv Nash S couv Nash S couv (mm) S karst (mm) Nash
1 oct-94 200 0,65 xx xx xx xx 199 0,67 250 94 0,66
2 nov-94 121 0,59 xx xx xx xx 104 0,55 131 99 0,59
3 déc-95 135 0,81 xx xx xx xx 120 0,79 156 91 0,82
4 mars-96 154 0,84 xx xx xx xx 151 0,81 179 100 0,85
5 déc-96 146 0,81 xx xx xx xx 135 0,82 173 88 0,81
6 nov-97 266 0,82 xx xx xx xx 262 0,85 266 267 0,83
7 déc-97 150 0,67 xx xx xx xx 135 0,66 174 99 0,68
8 nov-99 168 0,74 xx xx xx xx 151 0,7 186 129 0,75
9 sept-00 267 0,91 52,28 63,6 516 0,35 262 0,89 293 211 0,91
10 déc-00 117 0,6 65,25 0 167 0,16 104 0,53 125 99 0,61
11 janv-01 101 0,83 65,58 0 135 0,69 104 0,78 108 86 0,83
12 sept-02 238 0,9 61,18 19,1 453 0,79 231 0,85 291 125 0,91
13 oct-02 196 0,33 63,75 6,25 389 -0,3 231 0,57 220 145 0,38
14 déc-02 95 0,87 65,62 0 104 0,86 72 0,83 72 143 0,87
15 sept-03 254 0,9 39,55 127,25 294 0,8 231 0,94 218 332 0,92
16 nov-03 112 0,88 65,38 0 151 0,59 104 0,83 112 112 0,88
17 déc-03 101 0,89 65,85 0 135 0,87 88 0,87 101 101 0,89
18 sept-05 246 0,78 45,88 95,6 246 0,79 183 0,8 185 378 0,8
19 janv-06 132 0,85 65,38 0 199 0,5 120 0,81 149 96 0,85
20 oct-08 386 0,81 46,03 94,875 754 0,67 389 0,81 377 405 0,82
A partir de ces valeurs de S, les simulations de crues à Lavalette ont été menées. Leur
qualité a été estimée à travers le critère de Nash :
Comparaison des Nash obtenus pour différents calage de S (en noir, un seul S pour tout le bassin
optimisé - les autres couleurs correspondent à une spatialisation du S suivant 2 valeurs : S karst et
S couverture)
1,2
0,8
0,6
Nash opt
Nash (Skarst=f(piezo source)
0,4
Nash (Skarst=405)
Nash (Savg)
0,2
0
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
-0,2
-0,4
les performances du modèle, et les améliore même légèrement dans certains cas.
La recherche du jeu optimal de paramètres est toutefois plus longue, pour un gain
moyen de 1.6% sur le Nash, ce qui reste faible.
R2 = 0,4805
300
R2 = 0,6945
200
R2 = 0,6309
100
R2 = 0,7696
0
0 10 20 30 40 50 60 70 80 R2 = 90
0,7894 100
-100
Les relations entre la piézométrie mesurée à Bois Saint Mathieu, Claret et la source
ainsi que l’humidité fournie par SIM avant la crue ont été tracées en fonction de la
valeur des différents S distribués calés (cf. exemple graphe ci-dessus). La qualité de la
régression entre indicateur et condition initiale est estimée avec le coefficient de
détermination R2 (cf. tableau ci-dessous) :
Il en ressort que :
- l’utilisation de la piézométrie à la source pour pré-déterminer la valeur de S karst
(étude abaque), permet d’améliorer très fortement la corrélation entre indicateur (S
karst) et condition initiale du modèle. Ce résultat s’explique par la construction
même de la relation. Toutefois, le résultat sur les simulations de débits à Lavalette
étaient détériorées. Cette solution ne permet pas d’améliorer simultanément la
simulation des débits de surface et la corrélation entre l’initialisation du modèle et le
niveau piézométrique ;
- L’utilisation d’un S karst constant et maximal détériore les relations (S karst et S
couv). Or précédemment, nous avions montré que cette distribution de S détériorait
aussi les simulations de débits à Lavalette. Elle est donc sans intérêt pour la
présente étude ;
- L’utilisation d’un couple de S optimaux (S karst, S couv) conduit à une amélioration
de la relation entre la condition initiale S karst et la piézométrie par rapport à
l'utilisation d'un unique S uniforme sur le bassin. En particulier pour l'initialisation du
modèle à partir de la connaissance de la piézométrie à la source (R2 = 0.67 pour S
uniforme --> R2 = 0.88 pour S karts optimal). Seule la piézomètrie de Bois Saint
Mathieu s'éloigne de la condition initiale du modèle avec S distribué.
De plus, nous avions vu que les simulations en débit étaient améliorées par la prise
en compte de cette distribution de S. Cette distribution du S pourrait donc être une
piste à explorer. Toutefois dans la suite de l’étude, nous avons préféré ne plus la
considérer, car elle nécessite l’ajout d’un paramètre de plus au modèle pour un gain
limité.
300
250
R2 = 0,878
200
150
R2 = 0,669
100
50
R2 = 0,9994
0
0,00 10,00 20,00 30,00 40,00 50,00 60,00 70,00
-50
7.3.6. Calage
Parmi les 27 crues sélectionnées, 6 crues présentent un débit de pointe très faible (<
30 m3/s) comparé au reste de l’échantillon ce qui laisse suggérer un comportement
différent du bassin lors de ces crues. La calibration est alors menée sur les 21 crues
restantes. Les données de pluies utilisées en entrée du modèle sont soit les pluies
radar corrigées par le « Mean Field Bias » (MFB) soit les données de pluie au sol
fournies par les pluviographes. Ces pluies sont interpolées par la méthode des
polygones de Thiessen. Le modèle s’applique à pas de temps horaire sur un Modèle
Numérique de Terrain (MNT) de 75 m de résolution. Le critère optimisé est celui de
Nash.
(Qsim (t ) Qobs (t ) 2
- critère de Nash : Nash 1
(Qobs (t ) Qobs ) 2
- ds = 0.28 j-1
- K0 = 0.3
- V = 1.3 m/s
- w = 101 mm
- La condition initiale est variable selon l’événement.
Tableau 24 : Valeurs de la condition initiale et du critère de Nash obtenues avec les données
de pluie au sol et de pluie radar (traitement HYDRAM ou CALAMAR) pour le jeu de paramètres
-1
suivant : ds = 0.28 j ; w = 101 mm ; V = 1.3 m/s ; K0 = 0.3
7.3.7. Initialisation
Figure 81 : régression linéaire établie entre la valeur calibrée de la condition initiale du modèle
S et la valeur de l’indice Hu2 pris en début d’épisode.
Une régression linéaire a été établie entre l’indicateur Hu2 et la condition initiale du
modèle S. Le coefficient de détermination R2 est de 0.69.
Figure 82 : régression linéaire établie entre la valeur calibrée de la condition initiale du modèle
S et la piézométrie à Claret prise en début d’épisode.
Plusieurs incertitudes sont supposées réduire la qualité des régressions linéaires. Ces
incertitudes concernent l’estimation des indicateurs et/ou celle de la condition initiale
optimale du modèle :
- l’indicateur Hu2, en tant que sortie de ISBA, dépend des forçages météorologiques
(pluies, température …) et géographiques (hypothèses concernant la structure des
sols ou le type de végétation). Des incertitudes sur ces forçages ou sur la structure
du modèle ISBA peuvent affecter l’estimation de cet indicateur ;
- la piézométrie est sans doute l’indicateur connu localement avec le plus de
précision. Sa capacité à rendre compte du remplissage de l’ensemble de l’aquifère
peut être remise en cause par la nature hétérogène de l’aquifère karstique et les
particularités du piézomètre. Néanmoins, les bonnes corrélations obtenues entre la
source du Lez (représentative de la dynamique de l’ensemble de l’aquifère) et la
plupart des piézomètres semblent indiquer que ces piézomètres sont représentatifs
du remplissage de l’ensemble de l’aquifère ;
- l’estimation de la condition initiale S est aussi sujette à de nombreux biais.
Premièrement, les incertitudes concernant l’estimation de la pluie sur
7.3.8. Validation
Pour contrôler la robustesse du modèle, trois tests de validation croisée (« split sample
test ») sont réalisés. L’échantillon de 21 épisodes est séparé en deux parties : 14
épisodes font servir à la calibration et 7 vont servir à la validation. Trois
échantillonnages différents sont réalisés. La calibration s’effectue sur les épisodes 1 à
14, pour le 1er échantillonnage, sur les épisodes 1 à 7 et 15 à 21, pour le 2nd
échantillonnage et sur les épisodes 8 à 21 pour le 3ème échantillonnage. Pour chaque
échantillonnage la validation se fait sur les 7 épisodes restants.
Le karst peut aussi contribuer à la crue de surface en se vidangeant par ses sources
qui constituent des exutoires pérennes ou temporaires de l’aquifère. On peut alors se
demander si la vidange par ces sources permet au karst de contribuer à la pointe de la
crue de surface. A la différence du mécanisme de participation par saturation de la
capacité de stockage, la contribution par les sources peut permettre un transfert
souterrain rapide de la crue par effet piston (transfert de pression). Ce phénomène est
évoqué dans les études géochimiques récentes menées à la source du Lez. En
supposant qu’un tel phénomène existe pour l’aquifère du Lez, le transfert souterrain
(transfert de pression) peut être plus rapide que le transfert en surface (transfert de
masse) intervenant dans le mécanisme de saturation de la capacité de stockage. La
vidange du karst par les sources pourrait donc permettre une contribution du karst plus
rapide que la saturation de la capacité de stockage et permettrait au karst de
contribuer à la montée ou à la pointe de crue.
Pour savoir si cette vidange des sources peut permettre une contribution du karst à la
pointe de la crue de surface, l’hydrogramme observé de la source du Lez est injecté
dans le modèle distribué et transféré jusqu’à l’exutoire. Les temps d’arrivée des pics de
crue sont ensuite comparés. Dans l’hypothèse où la fonction de transfert utilisée est
représentative du transfert en surface du bassin, si le pic de crue (de l’hydrogramme
transféré) de la source du Lez arrive plus tôt que celui de Lavalette, alors le karst
contribue davantage à la montée ou à la pointe de crue qu’à la décrue. En revanche, si
le pic de la source du Lez arrive plus tard, alors le karst contribue davantage à la
décrue qu’à la pointe ou à la montée. De plus, la source du Lez étant l’exutoire du karst
situé le plus en aval de l’aquifère et du bassin topographique, c’est elle qui peut
contribuer le plus rapidement à la crue de surface. En effet, par rapport aux autres
sources situées plus en amont, son débordement est plus rapide et le transfert de son
hydrogramme en surface (supposé plus lent que le transfert de pression souterrain) est
limité.
Tableau 26 : Valeurs et temps d’arrivée des débits de pointe de chaque pic de crue. tQHXL et
QHXL sont respectivement le temps d’arrivée et la magnitude du pic de crue à Lavalette ;
tQHXtS est le temps d’arrivée du pic de crue de la source du Lez transféré à Lavalette ;
t(QHX) est la différence entre le temps d’arrivée du pic de crue de la source du Lez (transféré
à Lavalette) et celui du pic de crue à Lavalette.
Remarque concernant le rôle des eaux karstiques à la genèse des crues du fleuve
Lez : Il conviendrait également de conduire une analyse sur la source du Lirou. La
dynamique de réponse de la source de trop plein du compartiment jurassique situé à
l’ouest de la faille de Corconne pourrait être différente de celle du Lez (plus rapide). Il
est probable que les débits de pointes de crue de la source du Lirou soient supérieurs
à ceux du Lez (Qlirou>20 m3/s en crue, le débit max de la source du Lez étant de
l’ordre de 16 m3/s). Dans ce contexte, la contribution des eaux des sources karstiques
au débit du Lez à Garrigliano pourrait être supérieure à la petite dizaine de %.
L'analyse des récentes instrumentations mises en place sur le Lirou pourront permettre
d'apporter des éléments de réponse à cette interrogation dans un futur proche
(nécessité de disposer de chroniques suffisamment longues). Cette discussion pourra
être complétée par les données de l’ouvrage Suquet qui permettent de qualifier
l’évolution de la charge hydraulique de ce compartiment.
7.5. CONCLUSION
8. Modélisation hydraulique
8.1. RESUME
Le modèle hydraulique de propagation d'une crue dans le Lez est ici présenté. Il s'agit
d'un modèle basé sur les équations de Barré de Saint Venant appliquées pour décrire
les écoulements à surface libre dans le cours d'eau du Lez. Le principe du modèle y
est présenté. Tous les ouvrages y ont été également introduits à partir des mesures
effectuées dans le cadre de l'AT1. L'influence sur la ligne d'eau de ces ouvrages est
prise en compte et analysée.
8.2. INTRODUCTION
Figure 85 : Interface de SIC pour les données topologiques du Lez de Lavalette à la 3ème
écluse
Nous avons également indiqué dans le rapport de l'AT1 l'importance des 25 ouvrages
en travers situés le long du Lez entre la confluence Lez-Lirou et la 3ème écluse, dont
21 sur le tronçon retenu pour la modélisation dans SIC (Station de la Valette à la
Station de la 3ème écluse). Ces seuils ont été entièrement relevés au GPS différentiel
et au Tachéomètre Laser. Des données complémentaires, en particulier sur le seuil de
la 3ème écluse ont été fournies par la DREAL.
(date début, date fin, stations disponibles avec les codes suivants : L=Lavalette,
G=Garigliano, E=3ème Ecluse) :
> Du 18/10/1994 06:00 au 26/10/1994 20:00 L
> Du 27/10/1994 06:00 au 11/11/1994 11:00 L
Mais parmi ces 21 crues, 1 seule (Du 19/10/2008 06:00 au 25/10/2008 06:00) présente
des enregistrements aux 3 stations (Lavalette, Garigliano, 3ème Ecluse), et 12 crues
aux 2 stations (Lavalette, Garigliano). Nous utiliserons ces 12 crues pour le calage-
validation du modèle hydraulique. Ce travail sera bien entendu, pour être cohérent,
réalisé sur le modèle hydraulique correspondant aux géométries des lits antérieurs à
2008, donc ne comprenant pas les aménagements récents de recalibrage et
d'endiguement vers Lattes. Ces aménagements seront saisis dans une autre version
des données géométriques (à partir des données BRL 2009 et des plans de
recollement des travaux BEC), et utilisée probablement pour les scénarios prospectifs.
3
Tableau 29 : Amplitudes des crues retenues (débits en m /s et cotes en m NGF)
Pour valider un modèle il est nécessaire d’avoir des données hydrauliques auxquelles
les données du modèle pourront être comparées. Sur la partie du fleuve Lez étudié,
trois stations hydrométriques sont répertoriées dans la Banque Hydro : la station dite
de « Lavalette », la station de « Garigliano » et la station de « La troisième écluse » (en
souligné dans le Tableau 27). Les données cote/débit à ces stations ont pu être
récupérées auprès de la Banque Hydro. Les cotes sont mesurées par des capteurs de
niveaux (indiquées en général en mm par rapport à un zéro échelle, cf. Tableau 30), et
les débits sont déduits de ces niveaux via une courbe de tarage mise au point à partir
de jaugeages (au moulinet pour les faibles débits et au Saumon ou à l'ADCP9 pour les
forts débits). L'aide de la DREAL a été très importante et utile dans la collecte et
l'interprétation de ces données.
ZeroEchelle_Lavalette=34,431 m NGF;
ZeroEchelle_Garigliano=12,881 m NGF;
ZeroEchelle_3emeEcluse=-0,0275 m NGF;
On peut remarquer que la courbe de tarage est la même pour la crue et la décrue.
Dans certains cas les courbes de tarage sont définies avec des hystérésis, ce qui n'est
pas le cas sur les stations étudiées ici. Voici le résultat au niveau des trois stations
étudiées (croix en bleu des points de la Q(Z) de la DREAL et points en rouge des
points de toutes les chroniques des 12 crues retenues ci-dessus) :
9
ADCP = Acoustic Doppler Current Profiler, il s'agit de mesurer le profil vertical de vitesses de
l’eau à partir des particules en suspension dans l’eau par effet Doppler
Courbe de tarage - L
600
Point Hydro
Q(Z) Dreal
500
400
Débit (m3/s)
300
200
100
0
34.5 35 35.5 36 36.5 37 37.5 38 38.5 39 39.5
Hauteur (m NGF)
Courbe de tarage - G
800
Point Hydro
Q(Z) Dreal
700
600
500
Débit (m3/s)
400
300
200
100
0
13 13.5 14 14.5 15 15.5 16 16.5 17 17.5
Hauteur (m NGF)
Courbe de tarage - E
800
Point Hydro
Q(Z) Dreal
700
600
500
Débit (m3/s)
400
300
200
100
0
2.5 3 3.5 4 4.5 5 5.5 6 6.5 7
Hauteur (m NGF)
Les écoulements du Lez sont modélisés grâce au logiciel SIC basé sur les équations
1D de Saint-Venant, avec gestion des lits multiples (lits mineur, moyen et majeur), ainsi
que des casiers aux nœuds.
¶S ¶Q
Continuité : ¶t + ¶x = q [1]
¶Q ¶Q²/S ¶z
Dynamique : ¶t + ¶x + g.S ¶x + g.S.J = k.q.V [2]
où :
Q le débit (m3/s),
q > 0 : apports
q le débit latéral par unité de longueur (m2/s) : q < 0 : pertes
J la pente de frottement,
Q² n² [3]
J = S² R4/3
1
où R est le rayon hydraulique (m) et n le coefficient de Manning (n = K, où K est le
coefficient de Strickler).
Les ouvrages qui peuvent être introduits dans le modèle hydraulique sont de type,
seuils, déversoirs, orifices...
Les seuils en travers sont modélisés dans SIC grâce aux formulations décrites ci-
après. Ces équations coïncident avec les équations classiques rencontrées dans la
littérature, mais garantissent également les continuités entre toutes les conditions
hydrauliques possibles (surface libre / en charge et noyé / dénoyé). Nous explicitons ci-
après les équations utilisées pour les seuils rectangulaires ainsi que pour les seuils
triangulaires. Un seuil trapézoïdal sera considéré comme la somme d'un seuil
rectangulaire et d'un seuil triangulaire. Un seuil oblique sera considéré comme un
demi-seuil triangulaire. Un seuil oblique tronqué à une cote de berge maximum sera
considéré comme la différence de deux seuils obliques, avec les cotes de radier
correspondant aux cotes min et max de ce seuil.
Q = µF L 2g h13/2 [4]
2
Le passage noyé-dénoyé s'effectue pour h2 = 3 h1, on a alors :
3 3
µS = 2 µF pour µF = 0.4 => µS = 1.04
Q
µF =
L 2g h13/2
2
Orifice - régime dénoyé (h1 w et h2 3 h1)
h1
La continuité vers le fonctionnement à surface libre est assuré quand W = 1, on a alors
µ = µF
Il existe deux formulations suivant que l'on est partiellement noyé ou totalement noyé.
2 2 W
Régime partiellement noyé (h1 w et 3 h1 < h2 < 3 h1 + 3 )
3 3
Q = µF L 2g [ ((h1 - h2)1/2 h2) - (h1 - W)3/2] [7]
2
2 W
Régime totalement noyé (h1 w et 3 h1 + 3 < h2)
Le fonctionnement déversoir orifice est représenté par les équations ci-dessus. Quel
que soit le type d'écoulement en charge, on calcule un coefficient de débit dénoyé
équivalent correspondant à la formulation classique de l'orifice dénoyé :
Q
CF =
L 2g W (h1 - 0.5 W)1/2
Pour les seuils triangulaires on peut choisir, soit une formule de Gourley :
α 2.47
Q = Cd tg(2) h1 [9]
Soit une formule gérant toutes les conditions d'écoulement (et en particulier noyé-
dénoyé), par différence de deux seuils triangulaires sur le même principe que pour les
seuils rectangulaires comme présenté ci-dessus. Pour cette approche nous
choisissons les équations suivantes :
α 2.5
Q = Cd tg(2) h1 (dénoyé) [10]
1 α 2
Q = 0.286 Cd tg(2) h1-h2 h2 (noyé) [11]
Et une transition noyé-dénoyée obtenue pour h2 = 0,8 h1. On vérifie qu'il y a continuité
entre les 2 équations précédentes lors de la transition.
Pour construire un modèle dans SIC, la première étape consiste à rentrer la géométrie
du fleuve ou du cours d’eau étudié. Il est donc nécessaire d’avoir des relevés
topologiques de terrain afin de créer un modèle le plus fiable possible.
Dans le cas de l’étude du fleuve Lez, un premier modèle avait déjà été effectué avec
des données sous forme Abscisse/Cote de BCEOM-EGIS datant de 1988 recoupées
avec des données plus récentes (cf. Rapport de l'AT1). Pour le calage du modèle
réalisé avec des crues datant des années 1999 à 2008, nous avons donc dans un
premier temps travaillé sur ce modèle. Puis, de nouvelles données ont pu être
récupérées auprès de BRL datant de 2009 ainsi que les plans de recollement des
travaux effectués (CAM) au niveau des digues des écluses. Nous avons donc créé
ensuite une nouvelle géométrie plus récente comportant les nouveaux profils en
travers.
a) Les structures
La connaissance des profils en travers n’est pas suffisante pour une modélisation
correcte. En effet, des aménagements disposés le long du fleuve influencent fortement
son écoulement. Sur le Lez, de nombreux ouvrages sont présents (cf. Tableau 27) :
des seuils, des passerelles, des écluses. Que ce soit en période de crue, en étiage ou
en fonctionnement normal de la rivière, ces ouvrages jouent un rôle essentiel dans les
équilibres en régime permanent, ainsi que dans les dynamiques transitoires. Il est donc
important de les modéliser. Les profils de ces ouvrages ayant été relevés au GPS
différentiel et tachéomètres laser en 2009 et 2010 nous avons pu les rentrer dans le
logiciel SIC. La prise en compte des seuils trapézoïdaux dans la version 5.22 de SIC
permet de modéliser avec une plus grande précision ces ouvrages, en général
constitués de structures obliques sur les berges.
seuil triangulaire pour les écoulements latéraux dont la valeur par défaut et de 1.32
(formule de Gourley et Grimp : Q=1,32*tan(α/2)*h2.7).
Nous allons réaliser un calage de ces coefficients par optimisation sur les seuils de
Garigliano, de Lavalette et de la 3ème écluse.
Il a été décomposé en trois parties différentes, chacune ayant son propre coefficient de
débit.
La première peut être modélisée par un seuil rectangulaire dont les équations de débit
sont celles utilisés dans SIC pour T=1 (cf. Notice théorique du logiciel SIC). Les
caractéristiques de ce seuil sont L=2m, hr=13,23.
La deuxième partie peut aussi être modélisée par un seuil rectangulaire dont les
équations de débit sont celles de SIC pour T=1. Les caractéristiques de ce seuils sont
L=22,2 et hr=13,54.
La dernière partie (les 2 berges latérales obliques) peut se modéliser par deux demi-
seuils triangulaires dont les équations de débit sont celles de SIC pour T=11 (avec
puissance 2.5), voire avec l'équation proche de Gourley (T=10, avec puissance 2.47).
Les caractéristiques de ces seuils sont de demi-angle atan(3/2) et atan(7/(16.40-
13,54)) et hr=13,54. Pour simplifier, on supposera que leurs coefficients de débits sont
les mêmes m3=m4.
Le débit passant sur le seuil à Garigliano est donc la somme de ces trois débits.
Calage de m1
Pour le calage du coefficient de débit m1 sur le seuil central mince paroi, nous pouvons
utiliser les données de la Q(Z) pour les cotes situées entre 13,23 m et 13,54 m (8
couples de points Q-Z), ou également les données de la crue 28/09/2000 (3 couples de
points Q-Z), 05/09/2005 (9 couples de points Q-Z), 19/10/2008 (14 couples de points
Q-Z) qui comportent des points dans cette zone de cotes. Ce choix n'est pas neutre car
suivant le nombre et la répartition des points, le critère quadratique minimisé donnera
plus de poids aux forts débits et aux zones où les points sont plus nombreux.
Avec la loi de seuil classique de SIC (T=1) on obtient les résultats de calage suivants :
Nb de points 8 3 9 14
0.6
Erreur relative en Débit (%)
60
0.5
Débit (m3/s)
40 0.4
0.3
20
0.2
0
0.1
-20 0
13.2 13.25 13.3 13.35 13.4 13.45 13.5 13.55 13.2 13.25 13.3 13.35 13.4 13.45 13.5 13.55
Cote NGF (m) Cote NGF (m)
On observe que les résultats les plus satisfaisants sont obtenus sur les données de la
loi Q(Z) qui comprend 8 points dans cette gamme de cote. On a un coefficient de seuil
de 0,5137 à l'optimum et un Nash de 0,99357. Les écarts relatifs sont importants aux
très faibles débits mais deviennent très faibles aux débits plus importants. Le
coefficient de 0,5137 peut paraitre raisonnable, surtout compte tenu du fait que la
nappe n'est pas aérée (ce qui conduit à des coefficients plus forts).
Calage de m2 et m3
On fixe donc ici m1=0,5137 calé précédemment. Pour les cotes > 13,54 m NGF on
dispose de beaucoup plus de couples de données, soit issues des données de crues
(687 points au total sur les 12 crues), soit directement de la loi Q(Z) (23 points) qui fait
référence (et dont d'ailleurs les données de crues sont issues, mais avec des nombres
et répartitions de points différents).
Si l'on prend tous les couples de points avec des cotes > 13,54, sans limite supérieure,
donc y compris pour les forts débits, on obtient :
Nb de points 23 23 40 130
Erreur max (%) 37% puis 5% 65% puis 3% 16% puis 4% 16% puis 2%
On peut également caler les coefficients m2 et m3 sur la totalité des points des 12
crues retenues (687 points) :
Nb de points 687
m2 0,405
m3 9,29
Nash 0,99849
500
400
Débit (m3/s)
300
200
100
0
13 13.5 14 14.5 15 15.5 16 16.5 17
Cote NGF (m)
On voit que les résultats obtenus pour la crue 05/09/2005 ou pour la totalité des 12
crues donnent un bon Nash (0.99902 ou 0.99849), un coefficient de débit m2 de l'ordre
de 0,41 ou 0,44 et un coefficient de débit m3 de 8,4 ou 9,3. Mais on remarque
également (Figure 92) que ce calage se fait en dégradant les résultats pour les débits
dans les gammes 20-100 m3/s et 180-420 m3/s afin de favoriser un comportement
moyen satisfaisant. La forme de la courbe laisse pourtant apparaitre une inflexion au-
delà des valeurs 200 m3/s qui laisse supposer un ennoiement progressif du seuil aux
forts débits, ce qui parait possible. Le calage effectué à ce stade suppose en effet une
équation dénoyée, car la cote aval n'est pas mesurée. Ce point sera vérifié lors des
simulation-calage-validation sur le modèle complet SIC qui lui tient compte, bien
entendu, de l'ennoiement éventuel du seuil à partir de la cote aval calculée par le
logiciel.
Nb de points 17 23 22 129
Erreur max (%) 17% puis 3% 65% puis 3% 14% puis 2% 15% puis 3%
Sur la totalité des points des 12 crues compris entre Zmin=13.54 et Zmax=15.04, on
trouve :
Nb de points 597
m2 0,318
m3 12,21
Nash 0,9989
Resultats pour Cd=0.318169 12.2106 - Nash=0.9989 Resultats pour Cd=0.318169 12.2106 - Nash=0.99279
600 700
Q données totales (713) points Q données totales (713) points
Q données utilisées (597) points Q Equation SIC
600
500 Q Equation SIC
500
400
Débit (m3/s)
Débit (m3/s)
400
300
300
200
200
100
100
0 0
13 13.5 14 14.5 15 15.5 16 16.5 17 13 13.5 14 14.5 15 15.5 16 16.5 17
Cote NGF (m) Cote NGF (m)
Figure 93 : Résultats du calage de m2=0.318 et m3=12.21 sur données 12 crues limitées entre
Zmin=13.54 et Zmax=15.04 et vérification de la loi ainsi obtenue sur tous les points
Il reste néanmoins vrai que le coefficient m2=0.318 peut paraitre faible, et que le
coefficient m3=12.21 parait très fort et bien éloigné du coefficient par défaut de la loi de
Gourley (m3=1.32).
Une explication possible est que les écoulements lors des forts débits génèrent des
débordements importants en rive droite qui sont bien pris en compte lors des
jaugeages réalisés à partir du pont (Figure 94), mais qui sortent de l’ébauche du seuil
dont la géométrie retenue est représentée Figure 90. Mais cette raison ne nous semble
pas entièrement satisfaisante car nous avons tenté de caler l’angle du seuil oblique
rive droite et obtenu des valeurs incompatibles avec les observations terrain (173°). Il
est également vrai que les formules de seuils triangulaires (type Gourley) sont
adaptées aux seuils triangulaires avec contraction et pelle aval, ce qui n’est pas
vraiment le cas ici.
Le capteur de Lavalette est placé en amont du seuil de Gasconnet. Nous allons donc
utiliser les données hydrauliques de la station de la banque Hydro de Lavalette pour
caler les coefficients de débit du seuil de Gasconnet.
C'est un seuil qui est beaucoup moins régulier que celui de Garigliano, et nous allons
devoir faire un certain nombre de simplifications géométriques pour pouvoir
représenter son fonctionnement hydraulique. Ce seuil sera donc composé de six
parties différentes :
- un petit seuil central avec L=3 m et hr=34,76 m (coefficient de débit m1) ;
- un seuil bas avec L=15,5 (7+8,5) m et hr=34,96 m (coefficient de débit m2) ;
- un seuil triangulaire avec α=atan(3,8/0,7) (3,8 m étant la longueur provenant de 40-
4,2-7-3-8,5-13,5) et hr=34,96 m (coefficient de débit m3) ;
- un seuil haut avec L=13,5 m et hr=35,66 m (coefficient de débit m4) ;
- un seuil muret avec L=4,2 m et hr=35,96 m (coefficient de débit m5) ;
- un seuil triangulaire mimant les rives obliques, d'angle approximatif 70 (coefficient
de débit m6).
Le débit passant sur le seuil à Lavalette est donc la somme des débits passant sur ces
six seuils. Le capteur au niveau de Lavalette est placé à la cote NGF : 34,431
Si on cale dans un premier temps le seuil central seul (pour les points entre 34,76 et
34,96), on trouve m1=0,659. Le Nash est de mauvaise qualité (0,098), mais peu
significatif car correspondant à des débits très faibles :
15
10
Débit (m3/s)
-5
-10
-15
34.65 34.7 34.75 34.8 34.85 34.9 34.95 35 35.05 35.1
Cote NGF (m)
Dans un second temps on optimise les coefficients des débits supposés les mêmes
pour les seuils rectangulaire 2 et 5, et un autre coefficient pour le seuil 4. Pour le seuil
oblique, on lui donne la valeur de 1,32, valeur par défaut de la formule de Gourley.
Pour cette étape on se restreint aux cotes inférieures à 36,75 m NGF pour éviter
d'utiliser des points correspondant potentiellement à des écoulements noyés (cf.
discussion ci-après).
Les résultats obtenus par optimisation (toujours Simplex Nelder-Mead) sont les
suivant :
- m2 = m5 = 0,1871 pour les seuils 2 et 5,
- m3 = m6 = 1,32 pour les seuils obliques,
On obtient un Nash de 0,9984 avec ces valeurs sur les (716) points de calage, et de
0,9978 sur la totalité des (798) points disponibles.
400
350
300
Débit (m3/s)
250
200
150
100
50
0
34.5 35 35.5 36 36.5 37 37.5 38 38.5 39
Cote NGF (m)
500
400
Débit (m3/s)
300
200
100
0
34.5 35 35.5 36 36.5 37 37.5 38 38.5 39
Cote NGF (m)
Cette dernière courbe, obtenue par extrapolation (aux forts débits, correspondant à des
cotes > 36,75 m c'est-à-dire environ > 125 m3/s) des lois d'ouvrages calées comme
indiqué précédemment, montre un léger décrochage pour les cotes supérieures à
environ 38 m NGF, laissant supposer un ennoiement du seuil pour ces valeurs. Ce
point pourra être étudié lors de la modélisation complète dans le logiciel SIC.
L'ouverture en rive droite représente l’écluse qui est utilisée (peu fréquemment, et
encore moins en crues) pour le passage des bateaux. La modélisation se fera par :
- un seuil rectangulaire bas : L=57,1 et hr=2,70 (m1),
- un seuil rectangulaire : L=11+2,8=13,8 et hr=3,99 (m2),
- une vanne de caractéristique : cote de radier =1, Largeur=6.1, cote de
surverse =3,99 qui sera modélisé dans le programme Matlab par un seuil
rectangulaire : L=6,1 hr=3.99 (également m2),
- un seuil rectangulaire : L=2,,5 et hr=6,66 (m3),
- un seuil rectangulaire : L=4,7 et hr=2,99 (m4),
- un seuil rectangulaire : L=3 et hr=6,14 (m5),
- un demi-seuil oblique avec L=0, α=2*(90-tan-1(1,29/2,5))=125,41 et hr=2,70 (m6),
- un demi-seuil oblique avec L=0 α=2*(90-tan-1(2,6/6,8)) et hr=4,06 (m7),
- un demi-seuil oblique avec L=0, α=2*(90-tan-1(0,29/0,7))=135 et hr=2,70 (m8),
- un demi-seuil oblique avec L=0 α=2*(90-tan-1(3,15/5,3))=118,55 et hr=2,99 (m9).
Pour cet ouvrage, la station Hydro est récente (2008) et nous ne disposons que d'une
seule crue (19/10/2008), avec 18 points. Nous disposons également de la courbe de
tarage de la DREAL, mais qui est largement extrapolée pour les débits importants et
ne disposent que de peu de jaugeages (information discutée avec la DREAL, qui est
intéressée pour améliorer cette courbe de tarage).
Comme précédemment nous calons le premier coefficient m1, pour les cotes faibles
(entre 2,7 et 2,99 m NGF). Nous ne sommes pas arrivés à obtenir un calage
satisfaisant avec la cote hr=2,70. La forme de loi obtenue semblait indiquer que la cote
de seuil 2,70 m NGF n'était pas possible. N'étant pas sûrs de cette valeur, car ce seuil
est très grand (long et large) et diverses cotes ont été relevées sur ce seuil (entre 2.69
et 2,73 avec nos relevés de 2009-2010, et des points aux extrémités de 2,86 et 2,73 en
rive gauche et 2,78 et 2,75 en rive droite, données DREAL relevées en décembre
2008), nous avons tenté de faire une optimisation simultanée de cette cote hr ainsi que
du coefficient de débit m1. Nous obtenons un calage satisfaisant avec m1=0,34754 et
hr=2,8036, avec un Nash de 0,9975. Par la suite nous prendrons donc ces 2 valeurs
pour m1 ainsi que pour hr.
Nous effectuons ensuite le calage des autres coefficients pour les cotes comprises en
2,99 et 4,06 (4 points sur la Q(Z)). Pour les coefficients des seuils obliques, nous
supposons les valeurs par défaut Cd=1,32. Pour les coefficients concernés m2 et m4
nous trouvons Cd=-0,734 (coefficient négatif, impossible) pour un Nash de 0,9946 !
Ces résultats peuvent probablement être expliqués par le fait que les données utilisées
datent de 2008, alors que la géométrie du seuil, et en particulier les banquettes
latérales proviennent du recalibrage des rives du Lez dans cette zone, travaux réalisés
en 2009-2010. Faute de données récentes nous allons donc utiliser des coefficients
par défaut réalistes pour ces autres seuils, et comparer les résultats obtenus à la
courbe de tarage actuelle, qu'il faudra certainement améliorer avec de nouveaux
jaugeages. Par ailleurs, les écoulements sur ces banquettes latérales sont très noyés
(il n'y a pas de chute à ce niveau), ce qui correspond à des coefficients faibles. La
littérature manque en fait cruellement d'études sur ce type de seuils, qu'ils soient
horizontaux ou obliques. De manière très heuristique et sans fondement théorique
nous prenons un coefficient m2=m4=0,1. Avec cette valeur nous obtenons les résultats
présentés à la Figure 99.
15
10
Débit (m3/s)
-5
-10
-15
-20
Figure 100 : Résultats du calage de m1=0,3475 et hr=2,804, sur données limitées entre
Zmin=2,70 et Zmax=2,99
600
500
Débit (m3/s)
400
300
200
100
0
2.5 3 3.5 4 4.5 5 5.5 6 6.5 7
Cote NGF (m)
Figure 101 : Comparaison loi "calée" de SIC par rapport à la Q(Z) DREAL
Une fois ces données calées, il faut pouvoir pour chaque tronçon rentrer des
coefficients de Strickler et des infiltrations. Ce sont les paramètres du modèle sur
lesquels on interviendra pour caler au mieux le modèle à la réalité. Pour commencer,
les infiltrations sont calées à 0. Les Strickler dépendent des caractéristiques du lit du
fleuve. Les valeurs provenant des études et calages précédents sont comprises entre
15 et 25, et ponctuellement 75 dans les zones bétonnées (Cf. Mauris 2010, Monnier
2011).
Le pont de l’Evêque est situé dans la ville de Montpellier au niveau de l’hôtel de région.
Il est composé de trois vannes mobiles gérées par la Mairie de Montpellier. Le but de
cet ouvrage est de conserver une hauteur constante à l’amont de ce dernier.
1m 1m 1m 1m
Berges: minimum
3m de delta z sur
Sens d’écoulement 10m environ
(obliques)
Gauche
Droite
Les données hydrologiques que nous possédons sur cet ouvrage nous proviennent de
la mairie de Montpellier. Nous avons les débits, les cotes et les mouvements de seuils
sur la période 2000-2011. En regardant les données, nous pouvons vérifier que le but
de cet ouvrage est d’imposer une hauteur d’eau constante Z=2,5m au niveau de
l’amont de l’ouvrage.
Pour ce faire, en période de crue les seuils (vannes) sont baissés et en période
d’étiage, ils sont levés. En période d’étiage, lorsque le débit prend les valeurs 0,488,
0,732, 0,977, 1,221 m3/s les positions angulaires des seuils varient entre α=20° et 21°.
En période de crue, lorsqu’on regarde la crue de décembre 03, pour des débits de 100,
373, 500 m3/s on a des valeurs de 53, 76, 83,3°.
Nous recherchons la loi qui lie la hauteur d’eau sur le seuil (h), le niveau absolu de
l'eau (NGF) et l’angle du clapet. D’après les données de plan de l’ouvrage on a la
relation suivante :
h=(Niveau d’eau)-(L*cosα+9,85)
avec 9,85 étant la hauteur en mètres de l’axe de rotation du clapet et L étant la hauteur
des vannes. L=2,77 d’après les plans fournit. Pour trouver le débit qui passe sur le
seuil on applique la formule du déversoir en prenant en compte les piliers sur le côté
(0.2h) :
Q=m√(2g)*(ld-0.2h)*h^(3/2)
Avec ici m=0,41 et ld=18,35 m (largeur d'un seuil, sachant qu'il y en a 3 d'identiques).
Pour rentrer les données du fonctionnement de ce seuil dans SIC, nous avons besoin
de la hauteur du seuil en fonction du temps. Sachant que nous avons l’angle
d’ouverture de la vanne, pour avoir accès à la hauteur du seuil, il faut calculer
hr=9,85+2,77cosα.
8.6.1. La Modélisation
Une première modélisation sur les différentes crues répertoriées a été effectuée. Dans
un premier temps, nous avons effectué ces modélisations sur l’ancien modèle
(données géométriques antérieures à 2009, c'est-à-dire avant le recalibrage du lit à
l'aval vers Lattes).
Nous présentons ici quelques résultats (pour les crues de novembre 1999, janvier
2001, octobre 2002, décembre 2002, décembre 2003, et octobre 2008). Les autres
sont disponibles dans le document Malaterre & Monnier 2011.
50,000
45,000
40,000
35,000
30,000
Débit
Sic
25,000
Banque-hydro
20,000
15,000
10,000
5,000
0,000
0 2 4 6 8 10 12 14
Temps (jours)
Gargliano janv01
140,000
120,000
100,000
80,000
sic-T
Q
banque-hydro
60,000
40,000
20,000
0,000
0 1 2 3 4 5 6 7 8 9 10
temps(calibrage 6h30)
Gargliano oct02
60,000
50,000
40,000
sic-T
30,000
Q
banque Hydro
20,000
10,000
0,000
0 1 2 3 4 5 6 7 8
temps(calibrage 5h41)
Gargliano dec02
500,000
450,000
400,000
350,000
300,000
sic-T
250,000
Q
Banque-Hydro
200,000
150,000
100,000
50,000
0,000
0 2 4 6 8 10 12 14
temps(calibrage 6h13)
Gargliano-dec03
600,000
500,000
400,000
sic
300,000
Q
banque hydro
200,000
100,000
0,000
0 2 4 6 8 10 12 14
temps(calibrage06h01)
160,000
140,000
120,000
100,000
Sic-NR
80,000
Hydro
60,000
40,000
20,000
0,000
0 1 2 3 4 5 6
120
100
80
60
Q(m3/s)
banque hydro
sic
40
20
0
0 1 2 3 4 5 6 7
-20
temps (jours)
Nous observons que dans l’ensemble, les courbes obtenues grâce aux valeurs de SIC
sont de la même forme que celles obtenues par la banque Hydro.
Cependant, nous pouvons remarquer que généralement les premiers pics sont plus
grands avec la Banque hydro. Cela s’explique par le fait que l’on n’a pas encore
modélisé les apports au niveau de la Lironde et du Verdanson, faute de données
disponibles. Ces apports seront reconstitués dans la suite de l'étude (fin 2011, début
2012). Diverses approches seront mises en œuvre et comparées : par simple
différence, par assimilation de donnée et par modélisation hydrologique simplifiée des
sous-bassins correspondants.
Afin d’optimiser l’utilisation de SIC pour les différentes crues étudiées, un script a été
créé à l’aide de Matlab pour lancer les calculs de SIC en mode batch et un module de
régulation BOLST a été utilisé, pour choisir la crue souhaitée. Ce module est
disponible dans la bibliothèque des modules de régulation de SIC.
8.7. CONCLUSION
9.1. RESUME
Les crues « éclair » parfois dévastatrices qui touchent les bassins versants
méditerranéens du Sud de la France sont difficiles à anticiper. L’imprécision des
modèles pluie-débit servant à les prévoir conduit les hydrologues à mettre à jour les
résultats de leurs simulations par des techniques d’assimilation de données. Celles-ci
cherchent à réduire l’incertitude sur l’une des composantes du modèle hydrologique
(forçages, paramètres ou variables d’état) pour diminuer l’erreur commise sur les
sorties (les débits). Ce rapport présente une procédure d’assimilation de données, le
BLUE, utilisée pour améliorer la prévision du pic de débit réalisée par un modèle
hydrologique événementiel. À une date de prévision t0, certains paramètres du modèle
sont recalibrés grâce aux données de débits observés à l'exutoire du bassin. L'étude
est menée sur le bassin méditerranéen du Lez situé au Sud de la France. Recalibrer la
condition initiale du modèle événementiel s’avère être la correction la plus efficace.
Dans 75 % des cas, celle-ci permet d’améliorer de 12 % en moyenne la magnitude du
pic de crue. Une analyse détaillée a permis de mettre en évidence des cas
d’application favorables à l’approche choisie ainsi que des cas pour lesquels il est
nécessaire d’étendre l’espace de contrôle ou de corriger le modèle.
9.2. INTRODUCTION
modèle ATHYS
comme maquette démonstrative et pédagogique dans le cadre du présent projet
méthode
MERCEDES qui combine une fonction de production de type SCS et une fonction de
transfert de type translation- stockage. La fonction de production SCS est paramétrée
- être
êt
tres sont :
-
2008),
-
Cette deuxième stratégie est celle mise en œuvre dans le cadre du présent projet et
appliquée au bassin versant du Lez, notamment dans le cadre de la thèse de M.
Coustau. Les résultats présentés ici font l’objet d’un chapitre de la thèse de M.
Coustau.
Les bassins méditerranéens du Sud de la France sont soumis à des pluies intenses
responsables de crues « éclair » parfois dévastatrices comme la crue du Gard en 2002
ou celle de l’Aude en 1999 [Gaume 2009]. Afin de mieux prévoir ce type d’événements,
des modèles pluie-débit sont utilisés. Ces modèles sont soumis à un certain nombre
d’incertitudes liées à leur structure qui implique une simplification des processus
physiques réels, aux données utilisées pour leur forçage ou la calibration de leurs
paramètres.
L’assimilation de données peut être utilisée pour réduire ces incertitudes. S’inscrivant
dans le cadre général des problèmes inverses [Tarantola, 1987], cette discipline a pour
principe de combiner deux types d’informations, l’une venant d’observations et l’autre
d’un modèle numérique, pour proposer une meilleure estimation de l’état d’un système
[Bouttier et Courtier, 1999]. D’abord utilisées en météorologie [Daley, 1991] et
océanographie [Ghil et Malanotte-Rizzoli, 1991], les techniques d’assimilation de
données sont de plus en plus utilisées en géosciences notamment en raison de la
disponibilité croissante de données de télédétection [McLaughlin, 2002 ; Reichle, 2008]
mais aussi dans les sciences de l’eau comme l’hydrogéologie [De Marsily et al., 1999]
ou l’hydraulique [Malaterre et al., 2010 ; Nelly et al., 2010 ; Ricci et al., 2011]. Ces
Dans l’approche présentée ici, à une date de prévision t0, on cherche à recalibrer les
paramètres du modèle en fonction des observations déjà disponibles, notamment les
débits à l'exutoire du bassin, et ainsi améliorer les prévisions de débit réalisées par le
modèle hydrologique événementiel. Cette correction se fait par la technique du BLUE
(Best Linear Unbiased Estimator).
Dans cette étude, on utilise une technique simple d’assimilation de données permettant
d’effectuer une analyse à un instant donné d’observation. L’analyse consiste à estimer
la valeur optimale des paramètres du modèle étant donné les observations, les valeurs
a priori des paramètres et les erreurs dont ces informations sont respectivement
entachées. Pour ce faire, on utilise l’algorithme du BLUE dont on donne ici une
description en partant de l’approche classique variationnelle 3D-Var.
1 1 o
J( x) (x x b ) T B 1 (x x b ) (y H (x)) T R 1 (y o H (x))
2 2
1 T 1 1 o
J inc (x l ; δx) δx B δx (y H (x l ) H xl (δx x b x l )) T R 1 (y o H (x l ) H xl (δx x b x l ))
2 2
dH (x) H (x l dx ) H (x l )
H xl
dx dx
J inc(x l ; x a ) 0
Pour identifier xa, on peut utiliser un minimiseur qui nécessite à chaque itération le
calcul du gradient, ou bien, on peut expliciter la solution par un calcul matriciel, et ainsi
formuler le Best Linear Unbiased Estimator (BLUE). Dans notre cas, étant donné la
petite taille des matrices d’erreurs de covariance B et R, le calcul direct du BLUE
[Gelb, 1974 ; Talagrand, 1997] est utilisé pour trouver le minimum de la fonction coût
Jinc(xl ; x). L’analyse est donnée par :
xa xb Kd
T 1 T
K H xl R(B H xl R 1H xl ) 1
d yo [H (xl) H xl (x xl)]
L’analyse obtenue avec l’algorithme du BLUE est équivalente à celle obtenue si l’on
minimise Jinc(xl ; x). Cette solution est une approximation de l’estimateur optimal, qui
lui minimise J(x). Afin de se rapprocher du minimum de J(x), l’algorithme du BLUE est
muni d’une boucle externe qui met périodiquement à jour la valeur du point de
linéarisation xl. Ainsi lors de la première itération de la boucle externe, la linéarisation
se fait autour de l’ébauche xb alors que lors des itérations suivantes, la linéarisation se
fait autour de l’analyse de l’itération précédente. Ceci permet de prendre en compte les
non-linéarités de l’opérateur d’observation H. La mise en œuvre de cette technique est
présentée ci-dessous.
Cette étude de sensibilité fait déjà l’objet d’une section de ce rapport, au niveau du
modèle hydrologique (partie n°7 ) . Le lecteur s’y réfèrera pour la suite.
L’objectif est ici d’améliorer la prévision du pic de crue. L’assimilation des débits
observés débute au début de l’épisode jusqu’à une date t0, située 3h avant le pic de
crue. A partir de cette date, on effectue une prévision. La figure n°4 illustre le mode de
fonctionnement de l’assimilation pour l’épisode d’octobre 2001. On distingue la
« période d’assimilation » qui s’étend du début de l’épisode jusqu’à l’observation située
3h avant le pic de crue, et la « période de prévision ». Les débits simulés sur cette
période sont calculés en utilisant l'ensemble des pluies observées sur la totalité de
l'épisode (hypothèse de pluie future parfaitement connue).
moyenne des valeurs calibrées de S. L’écart-type de l’erreur commise sur Vb est fixée
à 0.2 m.s-1, valeur de l’écart-type obtenu après la calibration événement par
événement. On suppose que les erreurs sur Sb et Vb sont décorrélées.
Les observations utilisées, stockées dans le vecteur yo, correspondent aux n premiers
débits observés depuis le début d’un épisode à l’exutoire du bassin. Les erreurs sur les
débits observés sont supposées être décorrélées. L’erreur d’observation étant la
somme de l’erreur de mesure et de l’erreur de représentativité (liée à l’opérateur
d’observation c’est-à-dire dans notre cas, au modèle hydrologique), une erreur
d’observation importante a été attribuée aux débits pour lesquels l’une ou l’autre de
ces 2 erreurs étaient importantes. Nous avons retenus les seuils suivants :
- pour un débit observé supérieur à 300 m3/s, l’écart-type associé est infinie, car
l’erreur commise sur la courbe de tarage par extrapolation est importante ;
- pour un débit observé inférieur à 20 m3/s, l’écart-type associé est infinie, car le
modèle hydrologique n’a été calé que sur les débits de crue ;
- enfin, pour un débit observé compris entre 20 et 300 m3/s, l’écart-type associé est
de 20 m3/s. Cette valeur a été retenue suite à une étude de sensibilité effectuée sur
R (section 3.4).
L’opérateur d’observation H est le modèle hydrologique forcé par les pluies. Appliqué à
x, il permet d’obtenir une chronique de débits simulés y = H(x). La formulation de la
matrice jacobienne Hxl associée à l’opérateur d’observation non linéaire H requiert
plusieurs intégrations du modèle hydrologique :
- une intégration avec les valeurs des paramètres stockés dans xl
- une intégration supplémentaire par paramètre perturbé : xl + [dS , 0]T pour la
perturbation sur S et xl + [0 , dV]T pour la perturbation sur V.
Période Période de
d’assimilation prévision
b
Figure 105 : Hydrogrammes observé (courbe bleue), simulé à partir de l’ébauche S = 160 mm
a
(courbe noire) et simulé à partir de l’analyse S = 131 mm (courbe rouge) après assimilation
d’une donnée (croix bleue) pour l’épisode d’octobre 2001. Le trait vertical noir sépare la période
d’assimilation de la période de prévision
Lorsque la donnée de débit assimilée (croix bleue) est supérieure au débit simulé avec
l’ébauche (en noir), l’algorithme d’assimilation diminue la valeur du contrôle pour
augmenter la valeur du débit simulé au temps d’assimilation. Cette correction du
contrôle affecte les autres débits de crue qui augmentent eux aussi : les débits obtenus
avec l’analyse (en rouge) sont supérieurs à ceux obtenus avec l’ébauche (en noir). La
correction de S entraîne donc une correction monotone de la magnitude des débits de
crue. Dans le cas de cet épisode, l’ensemble des débits de crue étant sous-estimé,
l’assimilation d’un débit observé en début de crue permet d’améliorer l’ensemble de la
crue et notamment le débit de pointe. Avant assimilation, la simulation sous-estime le
débit de pointe de 16% alors qu’après assimilation elle le surestime de 7%. L’erreur sur
l’estimation du pic de crue est réduite de 9% par rapport à la valeur du débit de pointe
observé.
Itération 1
xl1 = xb
xa1
Itération 2
xl2 = xa1
Figure 106 : Graphiques représentant la fonction coût J(x) (courbe rouge) et son minimum (rond
rouge) ainsi que l’évolution de la fonction coût incrémentale J(xl ; x) (courbe bleue), de son
a
minimum x (rond bleu) et de sa tangente (droite noire) au point de linéarisation xl (croix noire)
pour les 2 premières itérations de la boucle externe pour l’épisode d’octobre 2001
Cette approche permet, au fil des itérations, de se déplacer sur la courbe décrite par
J(x) pour en proposer une meilleure approximation parabolique J(xl(2) ; x). On se
rapproche ainsi du minimum de J(x) (rond rouge). Dans les tests effectués, on ne
réalise que cinq itérations de la boucle externe pour des raisons de coût de calcul. Ce
nombre d’itération est toutefois suffisant pour que le minimum de la fonction coût
incrémentale Jinc(xl ; x) converge vers celui de la fonction coût J(x).
Il a été montré que, sur le bassin du Lez, les pluies radar de début d’automne étaient
de meilleure qualité que celles de fin d’automne et d’hiver. En effet, le bassin étant à
une distance importante du radar (~ 60 km), la faible extension verticale des nuages et
la faible altitude de l’isotherme 0°C en hiver détériorent la qualité de la mesure
[Coustau et al., 2011 soumis].
Les données de pluies utilisées sont les pluies radar pour les épisodes d’octobre 2001,
septembre 2003, septembre 2005 et octobre 2008 et les données de pluie au sol pour
les autres épisodes. Ces pluies au sol sont fournies à pas de temps horaire et
proviennent de 4 pluviographes. L’un d’entre eux est situé sur le bassin topographique,
les autres sont situés dans un rayon de 5 à 10 km à l’extérieur du bassin.
Le paramètre S jouant sur l’intensité du pic de crue, le gain de l'assimilation est évalué
avec un critère portant sur le débit de pointe de la crue.
Qp si m Qp o b s
PDIFF
Qp o b s
où Qpobs est le débit de pointe observé Qpsim le débit de pointe simulé. Qpsim peut
correspondre au débit de pointe simulé avant assimilation (Qpb) ou à celui simulé après
assimilation (Qpa). PDIFF est l’écart relatif au débit de pointe observé. Il peut être
calculé avant assimilation (PDIFFb) ou après assimilation (PDIFFa).
Ce critère exprime l’écart entre le pic simulé et le pic observé et permet de quantifier
l’amélioration apportée par assimilation de données.
Le paramètre V jouant sur le temps d’arrivée du pic de crue, le gain apporté par
l’assimilation de données est aussi évalué avec un critère portant sur le temps
d’arrivée du débit de pointe.
tp = |tpsim – tpobs|
où tpobs est le temps d’arrivée du débit de pointe observé et tpsim est le temps d’arrivée
du débit de pointe simulé. tpsim peut correspondre au temps d’arrivée du débit de pointe
simulé avant assimilation (tpb) ou après assimilation (tpa). tp est le décalage en temps
existant entre le débit de pointe calculé et le débit de pointe observé. Ce décalage peut
être calculé avant assimilation ( tpb) ou après assimilation ( tpa).
tpb et tpa quantifient l’avance ou le retard qu’a le pic simulé sur le pic observé. Leur
comparaison permet de quantifier l’amélioration apportée par assimilation de données.
9.4.3. Résultats
Nombre de b a
S S PDIFFb tpa tp)
Pics données PDIFFa (PDIFF) tpb(h)
(mm) (mm) (h) (h)
assimilées
Oct 94 Pic 1 9 184 178 0.04 0.08 0.04 0 0 0
Oct 94 Pic 2 24 184 186 0.61 0.60 -0.01 1 1 0
Nov 94 Pic 1 0 101 101 0.95 0.95 0 3 3 0
Nov 94 Pic 2 20 101 131 0.84 0.31 -0.53 -1 -1 0
Nov 94 Pic 3 34 101 130 0.36 0.09 -0.27 0 0 0
Déc 96 12 141 145 0.16 0.14 -0.02 -1 -1 0
Déc 97 20 184 172 0.09 0.15 0.06 0 0 0
Jan 01 7 107 109 0.19 0.18 -0.01 0 0 0
Oct 01 1 160 131 -0.16 0.07 -0.09 0 0 0
Sep 02 5 211 243 0.2 -0.04 -0.16 -2 -1 -1
Déc 02 Pic 1 5 119 121 0.13 0.11 -0.02 0 0 0
Déc 02 Pic 2 26 119 113 -0.38 -0.36 -0.02 1 1 0
Déc 02 Pic 3 37 119 91 -0.16 -0.11 -0.05 -2 -2 0
Déc 02 Pic 91 119 88 0.89 1.06 0.17 0 0 0
4
Sep 03 1 273 266 -0.23 -0.19 -0.04 1 1 0
Déc 03 41 64 111 0.22 0.07 -0.15 1 1 0
Sep 05 Pic 1 5 302 210 -0.37 -0.08 -0.29 -1 -1 0
Sep 05 Pic 10 302 205 0.58 0.99 0.41 0 0 0
2
Oct 08 Pic 1 0 304 304 0.45 0.45 0 0 0 0
Oct 08 Pic 2 8 304 346 0.65 0.39 -0.26 0 0 0
b a
Tableau 32 : Résultats de la correction de S ; S , valeur d’ébauche ; S , valeur de l’analyse ;
PDIFFb et PDIFFa, écarts relatifs au débit de pointe observé respectivement avant et après
assimilation ; (PDIFF) gain apporté par l’assimilation de données sur l’estimation du débit de
pointe ; tpb et tpa décalage en temps entre les débits de pointe simulé et observé
respectivement avant et après assimilation ; ( tp) gain apporté par l’assimilation sur
l’estimation du temps d’arrivée du pic de crue
( tp| = 0). Pour 14 des 20 pics testés, l’assimilation améliore l’estimation du pic de
crue de 12 % en moyenne. Pour 2 d’entre eux (pic 1 de novembre 1994 et d’octobre
2008), l’algorithme n’a aucun impact : l’assimilation n’a pas commencé 3h avant le pic
de crue car il n’y a pas d’observations supérieures à 20 m3/s. Enfin, pour 4 des 20 pics
testés, l’assimilation détériore l’estimation de départ. Les raisons de cette détérioration
seront exposées dans la partie 9.4.4.
Pour les 4 pics pour lesquels l’assimilation est prise en défaut, le modèle ne reproduit
pas la montée de crue à la même vitesse que la crue observée. Comme on le voit pour
décembre 1997 (Figure 107), le modèle sans assimilation (courbe noire) sous-estime
les débits de début de montée alors qu’il surestime les débits de fin de montée. On
peut schématiser cette situation en disant que, dans la phase de montée de crue, les
hydrogrammes observé et simulé (avant assimilation des débits) se croisent. La
période d’assimilation comprenant un grand nombre de débits sous-estimés par le
modèle, la technique d’assimilation a tendance à compenser cette sous-estimation en
diminuant la valeur de S. Cette correction monotone sur l’événement a pour effet
d’aggraver la surestimation des débits de crue (courbe rouge : hydrogramme après
assimilation). Une assimilation trop précoce (avant croisement des courbes de
montées de crue) des débits trop faibles aggrave la surestimation du pic de crue dans
ce cas.
Pour améliorer les performances de l'assimilation dans de tels cas, il est envisageable
d’augmenter la valeur seuil sur les observations assimilées pour n’assimiler que des
débits situé après le croisement des courbes de montée de crue.
Parmi les quatre pics pour lesquels l’assimilation des premiers débits détériore la
simulation du pic, deux d’entre eux surviennent après un premier pic de débits. La
détérioration observée après assimilation de données peut alors s’expliquer non
seulement par un problème de représentation de la montée de crue mais aussi par une
erreur sur les pluies différentes entre les 2 pics. Pour rendre compte des
conséquences de ces erreurs après assimilation de données, un cas fictif a été mis en
place (Figure 108).
Un épisode de pluie a été créé à partir du premier pic de pluie de septembre 2005,
répété 2 fois. Les observations de débits (croix bleues) sont créées à partir de débits «
vrais » simulés avec les paramètres suivants : S = 250 mm ; w = 101 mm ; ds = 0.28 j-
1 ; V0 = 1.3 m/s ; K0 = 0.3. Tous les débits sont perturbés d’un bruit blanc dont l’écart-
type vaut 3% de leur valeur moyenne et constitue l’écart-type des erreurs
d’observation. L’ébauche xb vaut 200 mm et l’écart-type de son erreur est fixé à 10 %
de sa valeur. La modélisation et l’assimilation se font sur le modèle hydrologique
simplifié à une maille de 114 km2 située à 10 km de l’exutoire. Dans cette expérience,
le premier pic de pluie n’est pas perturbé alors que le second est sous-estimé de 20 %,
10 % puis surestimé de 10 %, 20 %. La prévision se fait sur le 2ème pic de débit. Les
courbes de montées étant semblables, seule une erreur différente sur l’estimation des
pics de pluie est donc testée. Les données de débits sont ensuite assimilées de 2
manières différentes : i) l’assimilation porte sur les données du 1er pic de débit et
celles du début du 2nd (assimilation « groupée ») et ii) seules les données de débits du
début du 2nd pic sont assimilées (assimilation « séparée »). Dans le premier cas, les
données de débits de tous les pics sont assimilées ensemble, de façon « groupée »,
alors que dans le deuxième cas, les données de chacun des pics sont assimilées de
façon « séparée ».
Figure 108 : Hydrogrammes issus du cas fictif réalisé en « assimilation groupée » (a) et en
« assimilation séparée » (b). La courbe rouge représente l’hydrogramme obtenu à partir de
l’analyse, la courbe bleue représente l’état vrai, les croix bleues correspondent aux
observations, les ronds bleus sont les observations assimilées. Le trait noir vertical représente
l’instant de prévision.
Ce cas test permet de conclure que le traitement des pics de manière séparée permet
d’améliorer les résultats de l’assimilation. Ainsi, l’assimilation « séparée » des pics de
crue a été testée sur le pic 4 de décembre 2002 et le pic 2 de septembre 2005. Pour le
pic 2 de septembre 2005, l’assimilation « séparée » détériore toujours l’estimation du
pic de crue ( |PDIFF| = +0.17) mais de façon moins importante que l’assimilation
9.5. CONCLUSION
10. Conclusion
Le modèle par réseaux de neurones (RN) permet de simuler les débits journaliers à la
source du Lez à partir d’observations de pluies sur le bassin d’alimentation du karst. La
structuration du modèle proposée a été conçue de façon à capter une partie du
fonctionnement hydrogéologique du bassin karstique du Lez. Il en ressort la possibilité
de simuler d’autres grandeurs que les débits de la source, comme la piézométrie du
drain principal par exemple. Après les phases de calage et de validation croisée
relativement lourdes, les simulations à l’échelle d’un événement sont menées
rapidement. Les performances de ce modèle sont particulièrement bonnes et son
application pour la prévision des crues à la source sont satisfaisantes comme en
témoignent les valeurs du critère de persistance. Cette approche a été comparée avec
celle adoptée par TEMPO-VENSIM, l’analyse menée montre une bonne
complémentarité des deux modèles.
Le modèle hydrologique distribué (ATHYS – SCS LR) permet, à partir d’un formalisme
parcimonieux et robuste, de simuler les crues de surface rapides et importantes
survenant sur le bassin du Lez en amont de Montpellier. En effet, les performances, en
termes de modélisation, sont particulièrement bonnes avec des critères de Nash assez
élevés et une bonne anticipation du pic de crue. L’unicité du jeu de paramètres de
calage, quelque soit la crue considérée, et la bonne corrélation entre l’état hydrique
initial du système observé (par la piézométrie ou l’humidité des sols) avec la condition
initiale du modèle rendent possible son utilisation en mode opérationnel pour la
prévision des crues. En contrepartie, les simplifications opérées sur le fonctionnement
de ce bassin sont à l’origine d’incertitudes sur les gammes de crues absentes de
l’échantillon de calage. Afin de palier à ces difficultés, une procédure d’assimilation des
données de débits observés au début de la crue a été mise en place pour corriger les
défauts du modèle hydrologique. Cette chaîne d’assimilation a elle aussi été
développée avec un souci de robustesse et de rapidité des temps de calcul afin d’en
permettre son utilisation en mode opérationnel. La correction effectuée, basée sur la
technique du BLUE, a été adaptée au regard de la confiance que l’on peut avoir dans
les mesures des débits. Une interprétation au cas par cas de la correction proposée
par l’assimilation a été menée afin d’optimiser le mode d’utilisation de cette chaîne
d’assimilation en mode opérationnel. Elle permet une amélioration des performances
du système de prévision des crues.
Ces modèles seront utilisés par la suite pour simuler des scénarios qui seront
présentés dans le cadre du prochain rapport (livrable L4.2).
11. Bibliographie
Tikhonov, A., and V. Arsenine, (1976) Méthodes de résolution de problèmes mal posés, MIR,
Moscow, 1976.
Tikhonov,A . N., and A. V. Goncharsky(E ds.), Ill-Posed Problems in the Natural Sciences, MIR,
Moscow, 1987.
Fleury, P., Ladouche, B., Conroux, Y., Jourde, H., Dörfliger, N. 2009. Modelling the
hydrologic functions of a karst aquifer under active water management – the Lez
spring. Journal of Hydrology. 365, 235-243, doi:10.1016/j.jhydrol.2008.11.037.
Tritz S., Guinot V. & Jourde H. 2011. Modelling the behaviour of a karst system
catchment using non linear hysteretic conceptual model. Journal of Hydrology, 397,
250-262.
Jacob CE, Lohman SW (1952) Nonsteady flow to a well of constant drawdown in an extensive
aquifer. Trans AGU 33(4):559–56
Maréchal, J.C., Ladouche, B., Dörfliger, N., Lachassagne, P. (2008) - Interpretation of pumping
tests in a mixed flow karst system, Water Resources Research, 2008, 44, W05401,
doi:10.1029/2007WR006288, 2008.
Perrochet, P. 2005. A simple solution to tunnel or well discharge under constant drawdown.
Hydrogeology Journal, 13, 886-888. DOI:10.1007/s10040-004-0355-z.
Kiraly, L. (1998). _ Modeling karst aquifers by the combined discrete channel and continuum
approach. In: Bulletin du Centre d'Hydrogéologie de Neuchâtel 16, pp. 77_98. See pp. 4, 143.
Mazzilli, N (2011) Sensibilité et incertitude de modélisation dans les bassins versants
méditerranéens à forte composante karstique. thèse de doctorat UM2 - HSM.
Neuman, S. P. (1977). _ Theoretical derivation of Darcy's law _. In: Acta Mechanica 25.3, pp.
153_170. doi: 10.1007/BF01376989. See p. 142.
Owen, S., N. Jones and J. Holland (1996). _ A comprehensive modeling environment for the
simulation of groundwater _ow and transport _. In: Engineering with Computers 12.3-4, pp.
235_242. url: http://www.andrew.cmu.edu/user/sowen/papers/gms/gms.html. See p. 141.
Shewchuk, J. R. (1996). _ Triangle: Engineering a 2D quality mesh generator and Delaunay
triangulator _. In: Selected Papers from the Applied Computational Geometry: Towards
Geometric Engineering FCRC'96, WACG'96 workshop held in Philadelphia, USA (May 27-28,
1996). Ed. by M. C. Lin and D. Manocha. Vol. 1148. Lecture Notes in Computer Sciences.
Springer-Verlag, pp. 203_222. doi: 10.1007/BFb0014497. url:
http://www.cs.cmu.edu/~quake/tripaper/ triangle0.html. See p. 144.
Stüben, K., P. Delaney and S. Chmakov (2003). _ Algebraic multigrid (AMG) for ground water
_ow and oil reservoir simulation _. In: Proceedings of the MODFLOW and More 2003
Conference held in Golden, Colorado, USA (September 16-19, 2003). Ed. by J. Dougherty and
S. Seo. url: http://www.scai.fraunhofer.de/fileadmin/download/samg/paper/Modflow_Paper.pdf.
See p. 144.
Kong A Siu, L. (2011) : Modélisation des crues de bassins karstiques par réseaux de
neurones. Cas du bassin du Lez (France). Thèse de Doctorat, Ecole Doctorale
SIBAGHE - UM2. HydroSciences Montpellier - Ecole des Mines d'Alès. Soutenue le 21
octobre 2011, 228 pages.
Aubert, D., Loumagne, C., and Oudin, L. (2003): Sequential assimilation of soil
moisture and streamflow data in a conceptual rainfall-runoff model. Journal of
Hydrology, 280, 145–161.
Gaume E., Livet M., Desbordes M., Villeneuve J.P. (2004) : Hydrological analysis of
the river Aude, France, flash flood on 12 and 13 November 1999. Journal of Hydrology,
286, 135–154.
Lhomme J., Bouvier C., Perrin J.L. (2004) : Applying a GIS-based geomorphological
routing model in urban catchment. Journal of Hydrology, 299, 203-216.
Perrin, C. (2002). Vers une amélioration d'un modèle global pluie-débit au travers d'une
approche comparative. La Houille Blanche, n°6/7, 84-91.
Vigicrue : http://www.vigicrues.gouv.fr/
Géoportail : http://www.geoportail.fr/5069711/visu2D/afficher-en-2d.htm
Google Map
http://www.agroparistech.fr/coursenligne/hydraulique/degoutte1.pdf
http://engees.unistra.fr/site/fileadmin/user_upload/pdf/shu/Guide_technique.pdf
Florian Mauris (2010) Rapport de stage « Modélisation hydraulique d’un fleuve côtier
traversant la ville de Montpellier et ses alentours : Le Lez ». ENSE3.
Bailly-Comte, V., Borrell-Estupina, V., Jourde, H., and Pistre, S. (2011) Influence of
karst aquifers on flood genesis and propagation in an ephemeral Mediterranean River:
A semi-distributed conceptual model of the Coulazou River (Southern France).Water
Resour. Res., submitted.
Berthet L., Andréassian V., Perrin C., and Javelle P. (2009): How crucial is it to account
for the antecedent moisture conditions in flood forecasting? Comparison of event-
based and continuous approaches on 178 catchments. Hydrol. Earth Syst. Sci., 13,
819–831.
Bessière, H., Roux, H., and Dartus, D. (2007): Data assimilation and distributed flash
flood modeling. Proceedings of the second Space of Hydrology Workshop, Surface
Water Storage and Runoff: Modeling, In-situ data and Remote Sensing.
Borrell-Estupina, V., Chorda, J., and Dartus, D. (2005): Prévision des crues éclair. C.
R. Geoscience, 337, 1109–1119.
Bouttier, F. and Courtier, P. (1999): Data assimilation concepts and methods, ECMWF
Lecture Note.
Bouvier C. and Delclaux F.: ATHYS : a hydrological environment for spatial modelling
and coupling with a GIS. Proceedings HydroGIS96, Vienna, Austria, 19-28, IAHS
publication n°235, 1996.
Buis, S., Piacentini, A., and Déclat, D. (2006): PALM: A computational framework for
assembling high performance computing applications. Concurrency Computat.:Pract.
Exper., 18(2), 247–262.
Coustau, M., Bouvier, C., Borrell-Estupina, V. , and Jourde, H. (2011): Flood modelling
with a distributed event-based parsimonious rainfall-runoff model: case of the karstic
Lez river catchment. Nat. Hazards Earth Syst. Sci., submitted.
Da Ros D., and Borga, M. (1997): Adaptative use of a conceptual model for real time
flood forecasting. Nordic Hydrology, 28(3), 169–188.
De Marsily, G., Delhomme, J.-P., Delay, F., and Buoro, A. (1999): Regards sur 40 ans
de problèmes inverses en hydrogéologie. Comptes Renus de l’Académie des Sciences
– Series IIA – Earth and Planetary Science, 329(2), 73–87.
Fleury, P., Ladouche, B., and Courtois, N. (2007): Aléas inondations de la villes de
Nîmes par contribution des eaux souterraines. Rapport final, BRGM/RP55558-FR, 152
p., 112 ill., 5 ann.
Fouilloux, A., and Piacentini, A. (1999): The PALM Project: MPMD Paradigm for an
Oceanic Data Assimilation Software. Lecture Notes In Computer Science, 1685, 1423–
1430.
Gaume E., Bain, V., Bernardara P., Newinger, O., Barbuc, M., Bateman A.,
Blaskovicová, L., Blöschl, G., Borga, M., Dumitrescu, A., Daliakopoulos, I., Garcia, J.,
Irimescu, A., Kohnova, S., Koutroulis, A., Marchi, L., Matreata, S., Medina, V., Preciso,
E., Sempere-Torresm, D., Stancalie, G., Szolgay,J., Tsanis, I., Velascom, D., Viglione,
A. (2009): A compilation of data on European flash floods. Journal of Hydrology, 367,
70–78.
Gaume E., Livet M., Desbordes M., Villeneuve J.P.(2004): Hydrological analysis of the
river Aude, France, flash flood on 12 and 13 November 1999. Journal of Hydrology,
286, 135–154.
Gelb, A. (1974): Applied optimal estimation, Cambridge Mass.: MIT Press. 1974
Habets, F., Boone, A., Champeaux, J. L., Etchevers, P., Franchisteguy, L., Leblois, E.,
Ledoux, E., Le Moigne, P., Martin, E., Morel, S., Noilhan, J., Segui, P. Q., Rousset-
Regimbeau, F., Viennot, P. (2008):The SAFRAN-ISBA-MODCOU hydrometeorological
model applied over France, J. Geophys. Res. 113.
Khal, B., and Nachtnebel, H., P. (2008): Online updating procedure for a real-time
hydrological forecasting system. XXIVth Conference of the Danube Countries. IOP
Conf. Series: Earth and Environmental Science, 4, 012001.
Kobold, M., and Kay, S. (2004): Flood forecasting and uncertainty of precipitation
forecasts. BALWOIS, in Ohrid, FY Republic of Macedonia, 25–29 May 2004.
Lagarde, T., Piacentini, A., and Thual, O. (2001): A new representation of data
assimilation methods: the PALM flow charting approach, Q. J. R. M. S., 127, 189–207.
Lhomme J., Bouvier C., Perrin J.L.: Applying a GIS-based geomorphological routing
model in urban catchment. Journal of Hydrology, 299, 203-216, 2004.
Malaterre, P.-O., Baume J.-P., Nelly, J.-B., and Sau, J.(2010): Calibration of open
channel flow models: a system analysis and control engineering approach. SimHydro
2010: Hydraulic modeling and uncertainty, 2–4 June 2010, Sophia Antipolis.
Moore, R., J., Bell, V., A., and Jones, D., A.(2005): Forecasting for flood warning. C.R.
Geoscience, 337, 203–217.
Nelly, J.-B., Dorée, C., Sau, J., and Malaterre, P.-O. (2010): Data assimilation for the
real-time update of a 1D hydrodynamic model, fault detection – application to the
automatic control of hydropower plants on the Rhône river. SinHydro 2010: Hydraulic
modeling and uncertainty, 2–4 June 2010, Sophia Antipolis.
Piotrowski, A., Napiorkowski, J., J., and Rowinski, P., M. (2006): Flash-flood
forecasting by means of neural networks and nearest neighbour approach – a
comparative study. Non linear Processes in Geophysics, 13, 443–448.
Quintana Seguí, P., Martin, E., Habets, F., Noilhan, J. (2009): Improvement, calibration
and validation of a distributed hydrological model over France. Hydrol. Earth Syst. Sci.,
13, 163-181.
Ricci, S., Piacentini, A., Thual, O., Le Pape, E., and Jonville, G.(2010): Correction of
upstream flow with data assimilation in the context of flood forecasting. Hydrol. Earth
Syst. Sci., submitted.
Tarantola, A.(1987): Inverse problem theory and methods for model parameter
estimation, SIAM: Society for Industrial and Applied Mathematics.
Thirel, G., Martin, E., Mahfouf, J.-F., Massart, S., Ricci, S., and Habets, F. (2010): A
past discharge assimilation system for ensemble streamflow forecasts over France –
Part 1: Description and validation of the assimilation system, Hydrol. Earth Syst. Sci.
14, 1623–1637.
Toukourou, M. S., Johannet, A., Dreyfus, G., and Ayral, P.-A.(2009): Flash floof
forecasting by statistical learning in the absence of rainfall forecast: a case study.
Engeneering Applications of Neural Networks, 98–107.
Tramblay Y., Bouvier C., Martin C., Didon-Lescot J.F., Todorovik D., Domergue
J.M.(2010):Assessment of initial soil moisture conditions for event-based rainfall–runoff
modelling. Journal of Hydrology, Vol 387, Issues 3-4, 176-187.
Weiss, A., Oudin, L., and Loumagne, C. (2003): Assimilation of soil moisture into
hydrological models for flood forecasting: comparison of a conceptual rainfall-runoff
model and a model with an explicit counterpart for soil moisture. Revue des Sciences
de l’Eau, 16(2), 173–197.
Xiong, L., and O’Connor, K., M.(2002): Comparison of four updating models for real-
time river fow forecasting. Hydrological Sciences Journal, 47(4), 621–639.
Yang, X., and Michel, C. (2000): Flood forecasting with a watershed model: a new
method of parameter updating. Hydrological Sciences Journal, 45(4), 537–546.
Zehe, E., and Blöschl, G. (2004): Predictability of hydrologic response at the plot and
catchment scales : Role of initial conditions, Water Resour. Res., 40, W10202,
doi:10.1029/2003WR002869, 2004.
Annexe 1
Figure 109 : Evolution des charges hydrauliques mesurées par le réseau de suivi du Lez. Les
évolutions piézométriques de Fontbonne (données CAM) et de la source des Fontanilles
(Onema-Brgm) sont également reportées.
Producteur donnée Libellé Nature XL2E (m) YL2E (m) Altitude (m NGF)
ONEMA-BRGM Aven de la Sœur (Conquérac) Cavité Naturelle 728125 1883621 130.0
ONEMA-BRGM Claret (ONEMA) Piézomètre 726884 1873815 139.0
CG34 Sce Fontanilles Sce Karstique 703766 1862097 82.2
CG34 F. Fontbonne Forage AEP 733669 1866205 57.0
CG34 F. Suquet Forage 717798 1860042 165.0
CG34 F. Crouzette Forage AEP 726171 1849108 51.0
Veolia Sce du Lez Forage AEP 721649 1858552 64.0
Veolia St Clément de Rivière Piézomètre 722273 1855018 64.4
Veolia Ste Croix de Quintillargues F3 Piézomètre 727052 1863043 137.9
Veolia St Bauzille Bois des Rosiers Piézomètre 731495 1866382 87.8
Veolia Mas de Martin F84 Piézomètre 730553 1867132 99.1
Veolia Les Matelles 1 Piézomètre 718063 1859851 111.5
Veolia Laudou Piézomètre 727647 1870043 114.6
Veolia Gour Noir Piézomètre 723491 1858884 83.2
Veolia Fontanès F2 Piézomètre 728171 1866289 94.0
Veolia Bois des Avants Piézomètre 723728 1862299 125.9
Veolia Carnas stade Piézomètre 731572 1869676 156.8
Veolia Claret Brissac Piézomètre 727553 1873518 139.0
Veolia Coutach Piézomètre 731395 1877076 138.2
Veolia St Gély les Vautes Piézomètre 720350 1854906 95.0
Veolia Bois de St Mathieu Piézomètre 721903 1865652 131.7
Figure 111 : Réponses impulsionnelles des composantes lentes des modèles de transferts
(ModHE_Ginger Intermédiaire_128j ; ModHE_ Ginger Intermédiaire _256j ; ModHE_ Ginger
Intermédiaire _366j) utilisés pour reconstituer les débits naturels de Lez en période de basses
eaux
Figure 112 : Réponses impulsionnelles des composantes rapides des modèles de transferts
(ModHE_Ginger Intermédiaire_128j ; ModHE_ Ginger Intermédiaire _256j ; ModHE_ Ginger
Intermédiaire _366j) utilisés pour reconstituer les débits naturels de Lez en période de basses
eaux
Figure 113 : Exemples de reconstitutions de débit naturels obtenus à l’aide des modèles de
transfert qualifiés par differents temps de régulation (compris entre 128 j et 366j).
(Qdbt+Qp) : débit naturel reconstitué en périodes de hautes eaux [Qdbt calculé à partir de la
hauteur mesurée dans la vasque et de la nouvelle courbe de tarage (courbe intermédiaire(01-
08-2001), établie par Ginger Environnement) ; Qp= débit prélévé]
1981 2 118 646 1 819 933 4 995 793 6 113 465 3 712 179 3 256 476 4 516 370 1 820 680 1 314 444 2 347 063 1 089 865 4 652 869 37 757 785
1982 6 253 217 8 303 697 5 368 490 6 493 986 3 802 582 2 175 377 1 398 437 1 323 018 974 581 1 149 828 5 302 719 5 318 031 47 863 963
1983 2 620 698 3 408 972 5 952 303 3 643 544 3 885 754 2 095 547 1 305 095 847 726 1 572 897 919 810 827 793 3 205 321 30 285 461
1984 2 919 338 2 812 457 4 387 630 4 390 217 7 059 535 5 350 873 2 531 945 1 627 004 1 454 200 2 797 397 11 239 938 11 472 564 58 043 098
1985 6 324 247 4 652 476 4 494 416 4 907 434 6 145 874 4 254 215 2 113 900 1 249 153 670 996 717 630 621 889 3 481 061 39 633 290
1986 3 121 069 9 489 796 7 697 143 7 260 399 5 474 900 2 712 656 1 701 122 981 935 1 061 121 4 896 173 4 245 786 4 028 068 52 670 167
1987 4 107 491 9 560 946 6 815 318 6 332 636 4 139 035 2 400 532 2 455 201 1 756 158 1 325 032 5 784 457 10 041 313 13 907 894 68 626 013
1988 14 556 413 7 810 819 5 113 947 7 924 781 9 017 654 5 054 763 2 046 341 1 602 645 1 358 821 2 634 072 7 446 184 2 924 582 67 491 021
1989 3 549 594 2 507 713 2 907 039 5 363 182 4 958 013 1 880 279 1 125 039 894 867 1 919 881 2 413 873 6 042 755 5 161 242 38 723 478
1990 2 151 237 6 158 298 2 978 548 5 606 578 4 236 581 2 467 490 1 443 055 1 004 275 566 630 4 859 927 4 748 959 6 807 007 43 028 585
1991 5 825 641 6 121 751 10 553 535 7 212 370 6 715 812 2 566 820 1 552 478 1 149 969 684 938 1 368 895 2 558 629 1 621 469 47 932 306
1992 2 354 255 2 623 346 3 414 999 2 058 623 5 164 730 7 586 162 3 841 822 1 726 722 3 246 561 4 484 269 2 187 213 3 515 494 42 204 197
1993 2 094 959 2 420 114 3 825 124 5 762 502 8 812 066 4 725 933 2 426 042 1 488 272 2 540 203 5 144 449 12 731 636 8 219 318 60 190 619
1994 6 979 288 11 270 206 6 597 409 4 265 248 4 606 563 2 266 816 1 345 049 586 277 4 213 927 11 614 968 11 890 221 7 975 549 73 611 523
1995 5 682 061 3 964 516 3 157 783 2 448 185 2 296 887 996 210 497 396 433 325 4 196 122 9 641 022 7 999 456 12 417 538 53 730 501
1996 18 310 406 11 785 548 11 203 894 7 668 164 8 080 085 4 829 665 3 473 359 2 053 579 3 595 451 5 675 486 7 086 960 15 988 968 99 751 566
1997 15 095 652 6 686 721 4 097 512 2 309 289 1 912 939 3 878 819 2 172 284 1 484 703 688 324 1 772 216 8 652 031 12 480 817 61 231 308
1998 10 176 097 6 164 078 3 351 376 6 923 219 9 182 532 5 813 078 2 224 885 1 385 541 962 974 730 360 435 018 1 012 279 48 361 437
1999 6 131 886 2 496 052 3 256 458 3 837 478 9 358 831 3 608 912 2 041 875 1 372 201 1 183 486 5 251 257 9 153 562 4 439 007 52 131 006
2000 2 325 252 1 794 166 1 490 899 3 484 596 5 816 232 2 796 627 1 393 634 722 640 2 209 614 5 213 477 5 451 113 10 185 100 42 883 351
2001 15 013 797 7 396 972 10 103 858 5 399 888 3 838 933 2 387 331 2 100 080 1 106 339 1 437 010 10 115 824 4 219 106 1 990 482 65 109 621
2002 2 234 128 6 272 989 7 253 323 6 112 178 5 115 096 4 487 618 1 839 275 1 943 467 6 819 099 8 499 496 7 470 369 15 012 734 73 059 773
2003 7 457 149 5 526 662 7 571 724 5 800 931 3 789 556 1 799 384 1 105 349 504 251 2 481 170 8 281 207 12 653 703 14 683 239 71 654 325
2004 6 600 796 7 053 299 9 239 028 9 230 708 7 783 093 4 767 379 2 351 014 2 427 560 3 377 237 6 736 650 6 702 756 5 992 065 72 261 586
2005 1 985 983 1 791 256 1 352 664 2 452 773 3 022 659 2 838 600 1 409 084 790 975 8 767 962 4 876 731 8 282 209 3 701 619 41 272 514
2006 9 985 038 8 207 870 4 746 358 2 983 627 1 628 268 1 071 981 441 189 278 717 3 602 353 4 488 929 3 432 074 4 053 957 44 920 361
2007 3 140 134 3 129 068 2 955 702 3 509 777 5 273 562 3 496 807 1 806 541 771 894 754 373 1 329 989 2 380 933 1 619 115 30 167 895
2008 6 076 134 4 186 443 2 793 150 4 057 252 8 392 007 5 982 569 2 626 651 1 558 953 1 185 479 5 595 482 12 217 452 9 219 951 63 891 524
2009 10 874 045 12 471 132 7 438 167 10 256 466 6 685 563 3 759 893 2 169 700 1 426 010 857 256 4 720 706 2 049 874 2 437 992 65 146 804
2010 7 102 624 9 018 536 8 703 564 6 662 442 6 641 853 2 766 913 1 828 242 1 209 334 1 346 335 4 189 987 49 469 831
Volumes mensuels et annuels (m3) des débits naturels estimés de la source du Lez
1982 2 775 168 2 461 018 2 876 429 2 793 312 2 752 186 1 882 656 1 605 139 1 816 042 1 700 006 1 965 341 2 361 226 2 314 742 27 303 264
1983 2 416 349 2 141 338 2 608 675 2 590 190 2 661 732 2 641 200 2 351 224 1 805 047 2 445 761 2 013 612 1 780 691 2 097 256 27 553 076
1984 2 712 666 2 662 148 3 053 428 3 193 983 3 087 876 3 346 635 2 729 030 2 091 571 1 588 378 2 487 888 2 594 627 2 840 581 32 388 811
2010 2 338 794 2 072 298 2 323 877 2 314 240 2 475 716 2 719 659 2 973 091 2 871 429 2 604 652 2 475 733 25 169 488
241
Cum ul du Volum e de sollicitation des réserves (m 3)
242
Année\Mois 1 2 3 4 5 6 7 8 9 10 11 12 Total général
1976 0 0 0 0 0 -103 261 -152 758 -119 688 -50 543 0 0 0 -426 250
1977 0 0 0 0 0 0 -47 799 -218 940 -158 738 -23 646 0 0 -449 124
1978 0 0 0 0 0 -26 434 -88 291 -146 914 -279 374 -272 632 -421 870 -98 122 -1 333 637
1979 0 0 0 0 -36 625 -175 387 -256 161 -374 644 -441 961 -69 627 0 0 -1 354 406
1980 0 0 0 0 -308 871 -5 849 -151 044 -164 866 -272 252 -375 940 -64 257 -259 864 -1 602 945
1981 -182 010 -401 299 0 0 0 0 0 -183 597 -291 552 -128 411 -484 256 -179 699 -1 850 824
1982 0 0 0 0 0 -122 199 -289 167 -498 543 -725 426 -815 512 -130 645 0 -2 581 492
1983 -258 261 -508 376 0 0 -21 436 -656 886 -1 046 129 -957 322 -872 864 -1 093 802 -954 558 -506 285 -6 875 918
1984 -289 990 -145 429 -102 483 -90 876 -100 252 0 -247 041 -464 568 -452 752 -300 282 -18 826 0 -2 212 498
Modèles numériques du Lez
1985 0 0 -3 793 -5 158 0 -96 129 -1 006 436 -1 302 325 -1 548 361 -1 776 306 -1 737 349 -226 622 -7 702 479
1986 -547 959 0 0 0 0 -556 271 -1 591 755 -1 964 564 -1 671 366 -725 590 -100 316 -59 870 -7 217 690
1987 -226 921 0 0 0 -75 660 -583 625 -724 103 -1 315 708 -1 586 389 -133 814 0 0 -4 646 220
1988 0 0 0 0 0 -13 962 -1 154 088 -1 425 589 -1 537 135 -807 205 -249 136 -309 560 -5 496 675
1989 -475 316 -321 524 -94 147 -155 920 -77 344 -1 253 276 -1 957 367 -1 373 392 -647 567 -604 386 -31 286 -6 860 -6 998 386
1990 -439 155 0 -277 125 -38 452 0 -406 934 -1 355 000 -1 567 490 -1 670 009 -646 428 0 0 -6 400 595
1991 0 0 0 0 0 -433 039 -1 554 164 -1 901 851 -2 006 716 -1 626 913 -258 581 -760 657 -8 541 922
1992 -383 752 -112 946 -47 229 -441 983 -87 716 0 -373 651 -1 062 451 -1 053 601 -4 427 -836 406 -312 086 -4 716 249
1993 -731 111 -354 605 -6 275 -47 077 0 0 -485 137 -1 363 922 -1 158 684 -18 306 0 0 -4 165 116
1994 0 0 0 0 0 -494 833 -1 762 284 -2 372 897 -459 501 0 0 0 -5 089 514
1995 0 0 -54 388 -387 209 -385 620 -1 783 142 -2 274 843 -2 140 255 -717 712 0 0 0 -7 743 169
1996 0 0 0 0 0 0 -188 333 -532 455 -404 590 0 0 0 -1 125 377
1997 0 0 -21 534 -282 157 -680 327 -102 527 -567 651 -1 130 624 -1 915 664 -924 545 -71 715 0 -5 696 743
1998 0 0 -116 763 -9 344 0 -70 805 -1 065 919 -1 614 612 -1 641 527 -1 783 583 -1 943 130 -1 651 652 -9 897 333
1999 0 -141 978 -646 030 -50 456 0 -39 587 -1 083 170 -1 449 113 -1 449 121 -690 582 0 0 -5 550 036
2000 -108 256 -460 593 -1 064 587 -608 862 0 -522 183 -1 587 411 -2 049 507 -1 065 224 -7 208 -96 396 0 -7 570 226
2001 0 0 0 0 -20 676 -645 711 -1 109 136 -2 039 669 -1 621 458 0 -21 305 -456 228 -5 914 183
2002 -379 375 -25 524 0 0 0 -234 980 -1 327 112 -1 568 944 -57 572 -59 789 0 0 -3 653 295
2003 0 0 0 0 -31 484 -1 513 331 -2 213 301 -2 203 892 -1 208 603 0 0 0 -7 170 611
2004 0 0 0 0 0 -49 654 -1 042 666 -836 635 -564 020 -137 818 -78 088 -2 160 -2 711 041
2005 -838 602 -665 890 -1 393 232 -585 426 -382 719 -542 882 -807 452 -1 070 040 -163 348 -73 056 0 -326 165 -6 848 811
2006 -49 289 0 0 -195 756 -1 137 342 -1 759 196 -2 025 626 -1 744 979 -721 526 0 -72 096 0 -7 705 808
2007 -197 239 -77 701 -153 591 -175 198 -27 441 -103 740 -1 295 729 -1 999 734 -1 941 157 -1 233 024 -1 067 630 -643 073 -8 915 257
2008 -6 535 0 -188 539 -407 085 0 0 -525 867 -1 333 305 -1 488 807 -568 497 0 0 -4 518 635
2009 0 0 0 0 0 -62 555 -857 358 -1 560 017 -1 604 219 -443 057 -473 974 -873 253 -5 874 434
2010 0 0 0 0 0 -177 340 -1 144 848 -1 662 095 -1 258 317 -409 116 -4 651 716
Volumes mensuels et annuels (m3) de sollicitation des réserves du karst du Lez
1985 3 005 364 1 663 641 1 158 503 1 654 966 2 921 080 936 680 0 0 0 0 0 1 063 238 12 403 472
1986 861 719 6 568 353 4 464 919 4 172 204 2 075 147 78 952 0 0 0 2 548 860 1 380 767 1 068 517 23 219 438
1987 1 165 951 6 690 393 3 596 486 3 234 332 998 110 13 293 65 362 0 0 2 885 198 7 098 443 10 859 616 36 607 183
2008 3 637 383 1 648 374 421 513 2 016 280 5 791 024 3 427 510 6 286 0 0 3 686 334 9 835 154 6 887 475 37 357 332
2009 8 480 843 10 326 762 4 886 585 7 890 004 3 989 491 962 870 0 0 0 2 736 704 115 137 938 923 40 327 318
2010 4 763 831 6 946 238 6 379 687 4 348 202 4 166 137 224 595 0 0 0 2 123 370 28 952 059
243
Modèles numériques du Lez
Annexe 2
postes selon
Saint-Martin (2002 – 2008) disponibilité ATHYS
Tableau 33 : Tableau récapitulatif des données pluviométriques utilisées pour les modèles
TEMPO, VENSIM, Réseau de Neurones, ATHYS.
NB : le poste pluviométrique noté "MAUGUIO" est celui de l'Aéroport de Montpellier-Fréjorgues
Figure 114 : Cartes de localisation des postes pluviométriques utilisés pour le modèle
hydrologique ATHYS.
A gauche, les postes utilisés pour recalibrer les images radar. A droite, ceux utilisés
pour simuler les crues du Lez avec le modèle hydrologique ATHYS.
Les deux cartes de localisation des postes pluviométriques utilisés pour le modèle
hydrologique ATHYS sont représentées ci-dessus (Figure 114). Deux cas sont
différenciés selon que l'on travaille sur le bassin versant, et que les pluies sont
directement utilisées en entrée du modèle (carte de droite de la Figure 114) ou image
radar redécoupée selon cette même emprise (contour vert du bassin versant), ou bien
selon que les données de pluies sont utilisées pour calibrer les images radar,
initialement traitées par CALAMAR et HYDRAM. En effet, l'emprise géographique de
ces images radar est très vaste et donc les postes pluviométriques retenus doivent
coïncider avec ce secteur. Ensuite, est extrait de cette image recalibrée, le secteur du
Lez du Lez seulement (contour vert sur la figure de droite) qui servira à son tour
d'entrée de pluie dans le modèle hydrologique ATHYS. Bien que les images radar
présentent un échantillonnage infra-horaire, le travail de recalibration sur l'ensemble de
l'image ne se fait qu'avec des valeurs journalières de la pluviométrie (méthode du
MFB).
Concernant les postes pluviométriques utilisés pour simuler les crues du Lez avec le
modèle hydrologique, ils ont été retenus en fonction de deux critères : ceux qui
proposent une intensité horaire de pluie et ceux qui sont situés à l'intérieur du bassin
topographique ou à proximité immédiate de celui-ci. La lame d'eau spatialisée sur
l'ensemble du bassin a été calculée selon la méthode des polygones de Thiessen par
le logiciel ATHYS.
Enfin, les simulations globales du cycle hydrologique réalisées sous GR4J ont utilisé le
seul poste pluviométrique présent au sein du bassin versant du Lez (Prades le Lez)
avec une résolution temporelle journalière.
Les postes pluviométriques utilisés en entrée des modèles hydrogéologiques
recouvrent l'emprise, non plus du bassin topographique, mais du bassin d'alimentation
de la source du Lez (en rouge sur la Figure 114 carte de droite). Les pluviomètres
retenus proposent des valeurs de cumul de pluie journalier.
Figure 115 : Carte de localisation des postes pluviométriques utilisés pour les modèles
hydrogéologiques TEMPO, VENSIM et Réseau de Neurones.