Modélisation des bilans de surface et des débits sur la France, application à la prévision d’ensemble des débits

La modélisation hydro-météorologique, présentation de la chaˆıne SIM 

La modélisation hydro-météorologique est engagée depuis plusieurs années à Météo-France, avec le développement et la validation de la chaˆıne SIM (SAFRAN-ISBA-MODCOU). Le schéma de surface ISBA (Interactions Surface Biosphère Atmosphère, [NOILHAN and PLANTON 1988]) a été développé au CNRM. Il est relativement ”simple”, car il avait pour objectif le couplage avec des modèles de circulation générale, pour l’étude du climat ou pour la prévision numérique opérationnelle, avec de fortes contraintes numériques. ISBA est le schéma de surface du modèle de climat de Météo-France ARPEGE-CLIMAT, et il est opérationnel depuis 1998 dans les modèle de prévision de Météo-France ARPEGE et ALADIN. ISBA, dans la version que j’ai utilisée dans ma thèse et qui est présentée en détails dans le chapitre 2, prend en compte trois couches dans le sol, un seul bilan d’énergie pour le sol et la végétation, et le calcul des transferts d’eau et de chaleur par la méthode Force-Restore. ISBA est très modulaire, et dispose de nombreuses autres options. Ainsi, il existe notamment la version ISBA-DF qui calcule par diffusion les flux d’eau et de chaleur dans le sol. Dans la version ISBA-A-gs ([CALVET et al. 1998], [CALVET et al. 2001]), la représentation de la végétation a été améliorée, en tenant compte de l’assimilation nette du CO2 et du stress hydrique pour décrire l’ouverture et la fermeture des stomates et pour simuler l’évolution de la biomasse et de l’indice foliaire. De plus, le schéma de surface ISBA est également couplé au modèle hydro-géologique MODCOU (MODèle COUplé, [LEDOUX et al. 1984], [LEDOUX 2003], [HABETS 1998]), développé par le CG/ENSMP (Centre de Géosciences de l’Ecole Nationale Supérieure des Mines de Paris). Le fonctionnement de MODCOU est détaillé au chapitre 3. Un tel couplage, dont le principe est décrit au chapitre 4, permet de valider à l’échelle régionale les différents processus d’écoulement paramétrisés par le modèle (superficiels et souterrains). Les débits intègrent en effet cette information et sont facilement comparables à des observations grˆace à un réseau assez dense de stations de mesure. Le couplage entre ISBA et MODCOU a été initié en 1998 sur les bassins de l’Adour et du Rhône ([HABETS 1998]). Ces deux modèles sont complétés par le système SAFRAN (chapitre 1) pour l’analyse des paramètres météorologiques de surface, pour constituer la chaˆıne SIM. Après avoir été validée sur trois grands bassins versants fran¸cais, le Rhône ([ETCHEVERS 2000]), l’Adour-Garonne ([MOREL 2002]) et la Seine (partie III de ce manuscrit et [ROUSSET et al. 2004]), la chaˆıne SIM a été étendue à l’ensemble de la France. De plus, elle est en partie exploitée de fa¸con opérationnelle à Météo-France, afin de fournir un suivi en temps réel des composantes des bilans d’eau et d’énergie, en particulier des réserves en eau du sol. La chaˆıne SIM France et son exploitation opérationnelle sont décrites au chapitre 5. 9 10 Deuxième partie : La modélisation hydro-météorologique, présentation de la chaˆıne SIM.  Le système SAFRAN (Système d’Analyse Fournissant des Renseignements Atmosphériques à la Neige) a été développé à l’origine dans le cadre de la prévision d’avalanches sur le massif alpin ([DURAND et al. 1993], [DURAND 1995]). Il avait pour but d’analyser sur les zones de relief des paramètres météorologiques de surface : température de l’air au voisinage de la surface, vent à 10m, humidité spécifique de l’air, précipitations solides et liquides, nébulosité totale, rayonnements solaire (direct et diffus) et infrarouge, au pas de temps horaire, et avec une discrétisation verticale variable (généralement fixée à 300m). Ensuite, ce système a été étendu sur toute la France, afin de fournir une analyse du for¸cage atmosphérique sur l’ensemble du territoire ([ETCHEVERS 2000], [LEMOIGNE 2002]). L’analyse des paramètres météorologiques effectuée par SAFRAN sur la France ainsi que sa validation sont décrites en détails dans [QUINTANA-SEGUI et al. 2007]. SAFRAN ne travaille pas sur un maillage régulier, mais sur des ”massifs”, des zones géographiquement et climatiquement homogènes. C’est le zonage SYMPOSIUM II de Météo-France qui est utilisé. Sur l’ensemble de la France, il y a 615 zones SAFRAN, de forme irrégulière, d’une surface en général inférieure à 1 000 km2 , et sur lesquelles les gradients horizontaux climatologiques (en particulier de précipitations) sont faibles. La carte 1.1 montre le zonage SAFRAN/France.  Le système SAFRAN fait appel à la fois à des observations et à des sorties de modèles (analyses ARPEGE ou CEPMMT). A partir de ces données, il fournit pour chaque massif des profils verticaux des paramètres météorologiques au pas de temps horaire. L’analyse SAFRAN se fait en plusieurs étapes : – Analyse des profils verticaux de température, vent, humidité, rayonnement et nébulosité Pour chaque massif, une ébauche du modèle est comparée aux observations, aux niveaux o`u ont été faites ces mesures, afin de vérifier la cohérence des observations et d’éliminer les observations douteuses. Ensuite, l’analyse des paramètres est effectuée à l’aide des observations valides et de l’ébauche en utilisant la technique de l’interpolation optimale, au pas de temps de 6 heures, et avec une résolution verticale de 300 m. Cas particulier du rayonnement : en raison de la trop faible densité de postes météorologiques fournissant des observations de rayonnement (atmosphérique et solaire), ces derniers sont calculés à l’aide d’un modèle de transfert radiatif ([RITTER and GELEYN 1992]). Celui-ci utilise les profils verticaux analysés de température et d’humidité, permettant l’estimation d’un profil de nébulosité. – Analyse des précipitations : L’analyse des précipitations est faite au pas de temps journalier, afin de pouvoir prendre en compte les points de mesure o`u seules des observations journalières sont disponibles. Ici, l’ébauche n’est pas fournie par des modèles, SAFRAN utilise un gradient climatologique moyen (comme décrit dans [ETCHEVERS 2000]). – Interpolation horaire : L’ensemble des paramètres sont interpolés au pas de temps horaire. Cette interpolation se fait de manière linéaire, sauf pour la température et les précipitations. Pour la température, les variations diurnes sont estimées à l’aide d’un modèle de transfert radiatif et d’un terme de rappel à l’équilibre ([MARTIN 1988]). La température à 2 m à 12 UTC est corrigée en fonction de la température maximale observée. La répartition horaire des précipitations est déterminée à partir de l’humidité spécifique. La limite pluie/neige est estimée à partir de l’altitude de l’isotherme 0.5 ˚C issue de l’analyse de la température, altitude ajustée si besoin par les observations de pluie et de neige dans la zone concernée. A l’issue de ces étapes, SAFRAN fournit des analyses horaires des paramètres météorologiques sur l’ensemble des massifs, sous forme de profils verticaux de résolution 300 m. Ces données sont ensuite interpolées sur la grille régulière 8 km sur laquelle travaille ISBA. Chaque maille appartient à un massif SAFRAN, pour lequel on connait le profil vertical des paramètres météorologiques. Une maille ISBA ayant une altitude connue, l’interpolation se fait entre les deux niveaux SAFRAN entourant cette altitude (distants de 300 m). Chapitre 2 Le schéma de surface ISBA

 Présentation 

Le schéma de surface ISBA a été developpé afin de simuler les échanges d’eau et d’énergie entre le sol, la végétation et les basses couches de l’atmosphère. Il est couplé avec les modèles opérationnels de prévision numérique utilisés à Météo-France, ARPEGE et ALADIN, mais aussi avec des modèles comme ARPEGE-climat (modèle de circulation générale), Méso-NH (modèle de méso-échelle non-hydrostatique), ou MODCOU (modèle hydro-géologique) (chapitre 4). 

 Les paramètres d’ISBA 

Les paramètres d’entrée d’ISBA sont classés en deux catégories : Les paramètres primaires : ce sont le type de sol (via les pourcentages de sable et d’argile) et le type de végétation ; Les paramètres secondaires : à partir de l’information texturale du sol décrite ci-dessus, on peut estimer d’autres paramètres à l’aide de relations continues ([NOILHAN and LACARRERE 1995]). 

 les paramètres liés au sol 

on les détermine à partir des relations de [CLAPP and HORNBERGER 1978] et du type de sol. – le contenu en eau du sol à la saturation wsat (en m3/m3 ), c’est le contenu en eau maximal du sol, en dessous de ce seuil il y a drainage gravitationnel, jusqu’à ce que le contenu en eau diminue et atteigne la capacité au champ ; – le contenu en eau du sol à la capacité au champ wfc (en m3/m3 ), à partir duquel le drainage gravitationnel devient très faible. Cette eau peut ˆetre utilisée par la plante par absorption par les racines mais il y a stress hydrique, jusqu’à ce que le contenu en eau atteigne le point de flétrissement ; – le contenu en eau du sol au point de flétrissement wwilt (en m3/m3 ), en dessous de ce seuil l’eau contenue dans le sol est inutilisable par les plantes, car elle forme des films très mince autour des particules du sol et les forces capillaires deviennent très fortes. On définit alors la réserve utile RU par le produit de la différence entre capacité au champ et point de flétrissement par la profondeur racinaire d2 : RU = d2×(wfc−wwilt). Elle représente la quantité d’eau disponible pour la plante.

les paramètres liés à la végétation :

– la fraction de végétation dans la maille veg (0 ≤ veg ≤ 1) ; – l’indice foliaire LAI (Leaf Area Index) est le rapport entre la surface de feuilles vertes et la surface du sol ; 

 – la résistance stomatique minimale Rsmin est définie comme la résistance au transfert de l’eau des racines vers les feuilles, et ce dans un environnement optimal (rayonnement solaire important, sol très alimenté en eau, air saturé et température optimale) ; – les longueurs de rugosité dynamique z0m et thermique z0h interviennent dans les équations de flux de la couche limite, et représentent les influences liées à la hauteur de la végétation sur les échanges turbulents avec l’atmosphère ; – les profondeurs des trois couches de sol, d1, d2, et d3. La couche de surface a une profondeur d1 au maximum d’1 cm et contrôle l’évaporation du sol nu. Elle fait partie de la couche racinaire de profondeur d2 qui permet de définir la quantité d’eau utilisable par la plante (RU). Enfin, la dernière couche a une épaisseur d3 − d2, d3 étant la profondeur totale du sol. 3. Enfin, à cette liste s’ajoutent l’albédo α et l’émissivité ǫ de la surface. 

 Les variables pronostiques d’ISBA

 ISBA calcule l’évolution temporelle de 11 variables pronostiques : – Ts, T2 : températures respectivement de surface et profonde ; – wg, w2, w3, wr : contenus en eau respectivement de la couche de surface, de la zone racinaire, du réservoir profond et du réservoir d’interception de la pluie par la végétation ; – wsf et w2f les équivalents en eau du réservoir de glace respectivement à la surface et dans la zone racinaire ; – Ws, Ds et Hs : le manteau neigeux est discrétisé en trois couches pour lesquelles sont définis respectivement l’équivalent en eau, l’épaisseur et la quantité de chaleur stockée par la neige. Pour une meilleure compréhension du fonctionnement d’ISBA, la figure 2.1 synthétise les composantes des bilans hydrique et énergétique. Fig. 2.1 – Principe du schéma de surface ISBA La version d’ISBA utilisée possède trois réservoirs d’eau dans le sol. Pour établir le bilan hydrique, ISBA calcule la quantité d’eau interceptée par la végétation (Wr), l’équivalent en eau du manteau neigeux total (Ws), la quantité d’eau présente dans le sol (wg,w2 et w3), ainsi que les flux  d’eau ruisselée (Qr), drainée (D) et évaporée, par la surface du sol nu d’une part (Eg, et sublimation de la neige Es), et par la végétation d’autre part (transpiration Etr et évaporation de l’eau interceptée Er) Le bilan d’énergie d’ISBA fait intervenir le rayonnement net (donnée d’entrée), le flux de chaleur sensible (H), le flux de chaleur latente (lié aux différents termes d’évaporation) et le flux de conduction dans le sol (G). C’est le terme de chaleur latente (ie d’évaporation) qui permet de coupler ces deux bilans. Le détail de ce fonctionnement est donné ci-après. 2.3.1 Aspect énergétique 1. Les équations d’évolution de Ts et T2 ∂Ts ∂t = CT (Rn − H − LE) − 2π τ (Ts − T2) (2.1) ∂T2 ∂t = 1 τ (Ts − T2) (2.2) La température de surface (équation 2.1) évolue à la fois sous l’influence du flux de chaleur dans le sol Rn − H −LE = G (avec Rn bilan radiatif à la surface, H et LE les flux de chaleur sensible et latent), et d’un terme de rappel à la température du sol profond T2. Le coefficient CT est la capacité thermique de l’ensemble sol-végétation-couverture neigeuse. La température profonde T2 (équation 2.2) est la moyenne journalière de la température de surface ; son évolution dépend d’un terme de rappel vers Ts avec une constante de temps τ = 1 jour. 2. Le bilan énergétique dans ISBA Il s’exprime par la relation suivante : Rn = RG(1 − αt) + ǫ(RAT − σT4 s ) = H + LE + G Ts est la température unique du milieu sol-végétation-neige. avec :    σ constante de Stefan Boltzmann αt albédo total (sol, végétation et neige) RG(1 − αt) rayonnement solaire descendant absorbé ǫRAT rayonnement atmosphérique absorbé ǫσT4 s radiations infrarouges émises par la surface H représente le flux de chaleur sensible, LE le flux de chaleur latente (évaporation) et G le flux de conduction dans le sol (stockage de chaleur dans le sol). Ces flux sont déterminés par des lois aérodynamiques : H = ρaCpCH Va(Ts − Ta) o`u Cp est la chaleur spécifique de l’air, ρa, Va et Ta sont la densité de l’air, la vitesse du vent et la température de l’air. CH est un coefficient d’échange qui dépend de la stabilité thermique de l’atmosphère et des rugosités z0m et z0h. Le flux de chaleur latente est décrit par la suite, car en effet il assure le lien entre bilan énergétique et bilan hydrique. 2.3.2 La gestion de l’eau 1. Les équations d’évolution de wg, w2, w3 et wr 16 Deuxième partie : La modélisation hydro-météorologique, présentation de la chaˆıne SIM. ∂wg ∂t = C1 ρwd1 (Pg − Eg − Fgw) − C2 τ (wg − wgeq) 0 ≤ wg ≤ wsat (2.3) ∂w2 ∂t = 1 ρwd2 (Pg − Eg − Etr − F2w) − Dr1 − Df1 0 ≤ w2 ≤ wsat (2.4) ∂w3 ∂t = d2 d3 − d2 (Dr1 + Df1) − Dr2 0 ≤ w3 ≤ wsat (2.5) ∂wr ∂t = vegP − Er − R − r 0 ≤ wr ≤ wrmax (2.6) o`u ρw est la masse volumique de l’eau. Les termes de drainage Dr1 et Dr2 sont définis par : Dr1 = C3 d2τ max [wdrain,(w2 − wfc)] (2.7) Dr2 = C3 τ (d3 − d2) max [wdrain,(w3 − wfc)] (2.8) (voir paragraphe Le drainage gravitationnel page 17) et la diffusion entre la zone profonde et la zone racinaire est exprimée par : Df1 = C4 τ (w2 − w3) (2.9) L’évolution du contenu en eau du réservoir d’interception Wr (équation 2.6) dépend du taux de précipitations P et de l’évaporation du réservoir Er. Lorsque wr excède le contenu maximal wrmax, une quantité Rr est ruisselée vers la surface. wrmax dépend uniquement de la densité du couvert végétal : wrmax = 0.2 × LAI × veg. La partie des précipitations qui n’est pas interceptée ((1− veg)P) atteint le sol nu, une partie de cette eau va ruisseler en surface (Qr), et seule la quantité Pg = (1 − veg)P + Rr − Qr s’infiltre dans le sol. Le contenu en eau de surface wg (équation 2.3) évolue en fonction du flux d’eau de surface (les précipitations Pg qui s’infiltrent, l’évaporation du sol nu Eg, et le flux Fgw d’eau se transformant en glace), et de la diffusion vers la couche inférieure (rappel vers le contenu en eau wgeq à l’équilibre entre les forces de gravité et de capillarité). Les coefficients C1 et wgeq ont été calibrés, ils dépendent de la texture du sol et du contenu en eau. Le contenu en eau de la zone racinaire w2 (équation 2.4) (incluant le réservoir de surface) évolue en fonction d’un terme source, le taux de précipitations infiltrées Pg, et de plusieurs termes puits : l’évaporation du sol nu Eg, la transpiration Etr, l’eau se transformant en glace F2w, le drainage gravitationnel vers la couche profonde et les échanges d’eau par capillarité entre ces deux couches. La zone profonde (équation 2.5) est alimentée en eau par le drainage depuis la zone racinaire (Dr1), en perd par sa base par drainage gravitationnel (Dr2) et peut en échanger par diffusion avec la couche racinaire (Df1). Les coefficients C2, C3 et C4 ont été calibrés et dépendent uniquement de la nature du sol. 

 Le bilan hydrique dans ISBA 

L’équation du bilan hydrique sur une durée t est : XPrécipitations = XEvaporation + XRuissellement + XDrainage + ∆w ∆w représente la variation du stock d’eau dans le réservoir total, il a été décomposé réservoir par réservoir dans le paragraphe précédent. Détaillons les autres termes