GROUPE DE PROJET MESONH
Compte-rendu de la réunion du 11 mai 2010
Rédacteur C.Lac, CNRM

Présents : B.Barret, S.Bieilli, A.Boilley, C.Bovalo, J.-P. Chaboureau, M.Chong, V.Ducrocq, J. Escobar, C. Lac, D.Lambert, M. Leriche, C. Mari, V. Masson, J.Payart, J.-P. Pinty, D.Ricard, E. Richard, S.Riette, Y.Seity, G. Tanguy, O. Thouron

1. Informations générales

Licences

Deux nouvelles licences ont été signées, la première avec l'Université Antilles-Guyane sur les thématiques pollution atmosphérique et dispersion des boosters d'Ariane, la seconde avec Numtech, pour des études comparatives avec WRF.

Workshop ALADIN

Du 17 au 20 mai a lieu en Croatie un atelier ALADIN/HIRLAM sur la stratégie de la Prévision Numérique à haute résolution (3km-500m). Le DWD, UKMO et ECMWF sont également représentés. Les principaux sujets sont la précision et l'efficacité des schémas numériques (approches spectral/semi-lagrangien/semi-implicite vs. Point de grille/eulérien/time-splitting), la scalabilité des codes sur machines massivement parallèles et l'évolution de la physique à ces résolutions.

AROME

Un atelier organisé par PREVI/LABO de retour d'expérience des prévisionnistes sur le modèle AROME a eu lieu du 3 au 5 mai. Le principal problème soulevé concerne la prévision de nuages bas en situation hivernale, avec une sous-estimation de la nébulosité pour des humidités proches de 100%. En l'absence de convection et en couche limite stable, le schéma statistique ne produit pas de nuage. Une solution possible va être examinée par S.Riette, consistant à introduire une variabilité sous-maille liée à l'hétérogénéité de surface, en plus de la turbulence et de la convection. Les brouillards de rayonnement ont également été sous-représentés cet hiver. En revanche, les prévisions d'AROME en situation convective apportent une information supplémentaire, dans les zones orographiques plus que dans les zones de plaine, les décalages dans l'espace et le temps étant principalement liés aux limites de prévisibilité de la convection.

Ergonomie de Méso-NH

S.Bielli présente un état des lieux de l'ergonomie de Méso-NH, à partir d'une comparaison avec WRF (20100511_Bielli_Ergonomie.ppt). Concernant le téléchargement du code, le A-Install n'est pas évident à trouver lors de la 1ère connexion. Il manque un fichier « Getting started » qui guide l'utilisateur pas à pas du téléchargement jusqu'au post-traitement. Avec WRF, la compilation peut être réalisée sur une partie du code, ce qui permet de limiter le temps de compilation. En terme d'utilisation, il serait utile d'associer la signification des variables pour les namelists des cas tests. Les noms des variables sont plus standardisés dans WRF. Du point de vue des fichiers de sortie, WRF différencie les fichiers restart des fichiers output. Le format de sortie Netcdf est un format direct des fichiers output de WRF, alors qu'il nécessite des outils de transformation dans Méso-NH. Ceci constitue un avantage pour l'utilisation de logiciel graphique comme NCL ou VAPOR. Une discussion s'est également engagée sur les avantages respectifs de DIAPROG et de NCL, sans qu'aucun consensus ne se dégage : NCL fonctionne également sous forme de lignes de commande et n'est pas parallélisé (« il faut penser parallèle », dixit Juan !).

SURFEX

Actuellement, SURFEX est parallélisé d'un point de vue MPI, depuis le modèle coupleur (piloté par le découpage en processeurs du modèle appelant). Chaque Surfex sur chaque processeur découpé en MPI est indépendant des autres Surfex des autres processeurs. Au sein de chaque découpage MPI, on ne peut donc appeler qu'une seule partie de la surface par module (successivement). Mais le code SURFEX n'est pas parallélisé dans son fonctionnement off-line et présente une mauvaise scalabilité lorsqu'il est couplé à un modèle atmosphérique. Une action de parallélisation par OPEN-MP est en cours au CERFACS.

Une autre action concernant SURFEX est également en cours au CERFACS dans le cadre du projet ACCLIMAT : il s'agit d'introduire le coupleur PALM entre SURFEX et Méso-NH afin de pouvoir faire tourner les modèles de manière indépendante. SURFEX pourra ainsi tourner à une résolution plus forte que Méso-NH (échéance fin été).

Valentina-Méso-NH

Suite au couplage par PALM de Méso-NH et SURFEX, il est envisagé d'utiliser le système d'assimilation VALENTINA avec Méso-NH (autre action CERFACS). Les principales caractéristiques du système Valentina sont :

L'intérêt du couplage Valentina-Méso-NH est fort pour la chimie, moins immédiat pour les variables météorologiques (1 seule variable à la fois). Pour plus d'information, contacter Brice Barret.

Stages

Le prochain stage Méso-NH aura lieu du 4 au 7 octobre 2010. Le prochain stage SURFEX suivra du 11 au 15 octobre.

2. Développements en cours 

Schémas d'advection

J.-P. Pinty présente l'historique des innovations des schémas numériques dans Méso-NH (20100511_Pinty_WENO.pdf). Dans la version courante coexistent actuellement le schéma PPM pour le transport des variables scalaires, associé à un schéma temporel FIT, et le schéma centré CEN4TH pour le transport de la quantité de mouvement, associé au Leap-Frog (LF). Cette situation n'est pas complètement satisfaisante car hybride entre 2 schémas temporels, et nécessitant un filtre d'Asselin pour le LF et de la diffusion numérique à régler pour le schéma de transport centré. Le développement en cours dans Méso-NH, réalisé par F.Visentin, consiste à introduire les schémas de type WENO (Weighted Essentially Non Oscillating) pour l'advection du vent (PPM étant instable pour les variables flux en grille C) : WENO3 (3ème ordre) et WENO5 (5ème ordre). WENO3 associé au schéma temporel FIT montre d'excellents résultats : sur des cas idéalisés (bulle froide de Straka, ondes de gravité piégées), il permet d'atteindre des nombres de Courant de 1 en s'affranchissant de la diffusion explicite, avec une très bonne précision. Il fonctionne en 3D pour tout type de conditions aux limites latérales, est parallélisé et documenté. WENO5 est en revanche instable en FIT : il doit donc être associé à un schéma temporel spécifique. L'introduction d'un schéma d'intégration temporelle de type RK3 est en cours, puis WENO5 devra être parallélisé (introduction de HALO3). Les tests vont se poursuivre sur des cas réels. Actuellement l'association PPM+WENO3+FIT est satisfaisante. Reste à poursuivre l'effort pour intégrer WENO5+RKX.

Spectres d'énergie dans Méso-NH et AROME

Un outil diagnostic de calcul de spectres d'énergie cinétique a été développé dans Méso-NH (à partir du DIAG ; 20100511_Ricard_Spectre.pdf). Le calcul, basé sur une décomposition en transformée de cosinus (DCT), est classiquement appliqué aux composantes résolues du vent, mais peut être appliqué aussi aux autres variables. Il a également été utilisé sur les sorties d'AROME. Il s'avère un outil pertinent pour visualiser la distribution d'énergie du modèle : comparaison aux pentes théoriques (k-3 dans les grandes échelles, k-5/3 à méso-échelle, k étant le nombre d'onde), comportement dans la queue du spectre (petites échelles), où est artificiellement dissipée l'énergie (impact de la diffusion numérique). Il permet également de définir une résolution effective, correspondant à la longueur d'onde à partir de laquelle la distribution d'énergie s'écarte de la pente théorique : les mouvements ne sont plus totalement résolus et la variabilité sous-maille est réduite. Cette résolution effective, estimée autour de 7dx avec WRF[1] (Skamarock, 2004), est voisine de 6dx avec Méso-NH, et de 9dx avec AROME, et est indépendante de la résolution. La différence entre les deux modèles est, selon toute vraisemblance, liée à la diffusion implicite du semi-lagrangien. La descente en résolution avec Méso-NH (1km-500m) tend à montrer un excès d'énergie dans les plus courtes longueurs d'onde, qui semble être lié à la surestimation de la part résolue dans la zone grise de la turbulence. Le calcul de spectres sera disponible dans la prochaine version de Méso-NH.

Tests sur la grêle dans AROME

Un certain nombre de cas très convectifs donnant lieu à des observations de grêle au sol, ont été testés avec AROME et le schéma microphysique ICE4 représentant explicitement la grêle (20100511_Seity_Grele.ppt). Si le modèle représente en général assez bien la présence de grêle, il surestime quasi-systématiquement l'occurrence de grêle sur l'ensemble du domaine. Des réglages s'imposent donc sur la limitation de la formation de la grêle et sur la possibilité de reconvertir de la grêle en graupel. Un autre défaut observé est la dépendance du processus de croissance humide au pas de temps du modèle (dépendance qui n'existe pas avec ICE3, opérationnel dans AROME). Cette dépendance disparaît si l'on introduit un time-splitting sur l'ensemble de la microphysique (dt=60s -> dt=15s, et dt=15s équivalent à dt=5s). Les réglages et l'évaluation du schéma ICE4 vont se poursuivre avec AROME sur l'année 2009 grâce à des données de grêlimètres et radars.

Impact de RRTM-SW

Le code RRTM pour le SW (RRTM-SW), développé au CEPMMT, a été implanté dans Méso-NH (20100511_Thouron_Rayonnement.ppt) en plus du schéma de Morcrette (MORC-SW), nécessitant au préalable une mise à jour de l'ensemble du code radiatif du Centre Européen. Les principales différences entre les 2 schémas sont que, d'une part MORC-SW considère un angle zénithal effectif (AZE) alors que RRTM-SW fonctionne en Ordonnées Discrètes et sépare le rayonnement direct du rayonnement diffus, et d'autre part MORC-SW dispose de 6 bandes spectrales pour l'absorption des gaz, contre 14 pour RRTM-SW. Le transfert radiatif dans les nuages étant maintenant traité suivant la méthode McICA (Monte Carlo Independent Column Approximation) dans la nouvelle version du code, l'hypothèse du recouvrement nuageux n'était plus présente. Le recouvrement maximum aléatoire a donc été réintroduit dans l'interface de Méso-NH, en attendant de tester dans le futur la nouvelle méthode McICA. L'impact du RRTM-SW a été évalué sur 3 situations idéalisées, en LES et en 1D : le cas de Sc marins de FIRE, le cas de Cu continentaux d'ARM et un cas de CL convective sèche à Niamey. Les simulations ont montré un impact négligeable de RRTM-SW comparativement à MORC-SW, pour un coût numérique double du schéma de rayonnement. Ce résultat n'est pas en contradiction avec Morcrette et al. (MWR, 2008), l'amélioration sur des simulations climatiques étant principalement liée au schéma de recouvrement McICA et aux paramétrisations des propriétés optiques. Le nouveau code RRTM-SW sera disponible dans la prochaine version de Méso-NH.

Simulation Xynthia-Pyrénées

Les modèles AROME et Méso-NH à 2.5km de maille horizontale se sont bien comportés lors de la tempête Xynthia sur les Pyrénées (20100511_Lac_Vent.ppt). Ils ont simulé des diagnostics de rafales de vent réalistes, jusqu'à 210km/h sur certains sommets (209km/h mesuré au Pic du Midi vers 21h). Une caractéristique exceptionnelle de cette situation est la bande de vents forts le long de la barrière des Pyrénées, au nord de la chaîne dans le flux Sud-Nord. Les simulations montrent clairement des ondes de gravité piégées, qui génèrent des vents extrêmement forts dans les vallées, expliquant les violentes rafales relevées par exemple à Luchon. La simulation en hydrostatique dans AROME (au lieu du NH par défaut) montre une structure différente avec l'absence de piégeage d'ondes dans les basses couches, et des déferlements d'ondes en haute altitude beaucoup plus présents : les vents dans les vallées sont alors plus faibles. La simulation NH se révèle plus réaliste (confirmé par les scores). Les simulations Méso-NH et AROME à 500m sont en cours d'évaluation, et nécessitent au préalable de disposer d'un fichier relief à 250m au minimum sur la zone, pour représenter correctement la longueur de rugosité liée au relief.

Développement grandes grilles

Actuellement, grâce au projet PRACE (http://www.prace-project.eu/news/prace-granted-over-4-3-million-core-hours-to-prace-prototype-systems), un million d'heures de calcul sont disponibles sur JUGENE/IBM-BlueGene/P à 290 000 coeurs pour 4 mois (jusqu'au 31 mai 2010), afin d'optimiser la scalabilité de Méso-NH. L'objectif du projet est d'optimiser le code au-delà des 100 000 coeurs de l'architecture BlueGene et de mener des tests initiaux sur la nouvelle architecture de Cray XT5.
Les tests préliminaires réalisés uniquement sur des FFT-3D parallèles pour le solveur de pression (code P3DFFT) sur une grille 4096x4096x128 montrent que la scalabilité des FFT-3D est au mieux en SQRT[nombre de processeurs] sur la BG/P. Cela est dû à l'algorithme[2], mais surtout à la topologie du réseau des machines BG/P = TORE 3D[3]. Pour une grille 4096x4096x128, Juan a mesuré environ un facteur 3 de speed-up en passant de 16384 processeurs (= 16Kprocs) à 131072 processeurs ( = 128Kprocs), ce qui n'est pas satisfaisant. Pour autant, il reste des raisons d'espérer :

  1. D'une part, la scalabilité est liée à la topologie du réseau : ainsi, sur les machines de type SGI/ICE = JADE machine du CINES (plus grosse machine de France pour la recherche actuellement), le réseau est un HYPERCUBE . Sur ce type de machine/topologie , on trouve dans la littérature une scalabilité linéaire en fonction du nombre de processeurs (Scalabilité = NB_processeurs). La contre-partie de ce type de réseau est qu'il faut connecter chaque processeur avec Log(Nb_processeurs) : cela coute plus cher qu'un TORE 3D (où chaque processeur n'est lié qu'à ses 6 plus proches voisins). Par exemple, un hypercube de dimension 7 (2^7 = 128 processeurs) a déjà 7 liens par processeur, alors qu'un TORE 3D (4 * 4 * 8 = 128 processeurs) n'en que 6. Techniquement, il faut de longs câbles capables de relier des noeuds distants. Actuellement, la machine SGI/ICE est annoncée supporter jusqu'à 32Kprocesseurs sans rajouter de matériel autre que les câbles pour relier les noeuds. Par ailleurs, la prochaine machine de IBM (IBM BlueGene/Q & 1.6 Million de coeurs) aura un TORE 5D, soit en théorie une scalabilité de (NB_processeurs) ^ 0.8, ce qui est plus près d'une scalabilité idéale linéaire que celle du TORE-3D.
  2. D'autre part, on peut envisager de remplacer la FFT3D, qui sert de pré-conditionneur/accélérateur aux méthodes itératives de gradient conjugué (GC), par des solveurs Multi-Grilles, beaucoup plus scalables. Mais le coût numérique du solveur Multi-Grille+GC est d'au moins un facteur 5 par rapport à la FFT+GC. Des tests supplémentaires pourraient permettre d'améliorer la méthode. Enfin, les solveurs Multi-Grilles sont aussi fortement recommandés dans l'optique calcul sur Carte-Graphique car il faut sur ce type d'architecture éviter au maximum de déplacer sans cesse les données de la mémoire du processeur vers celle de la carte graphique (goulot d'étranglement de la connexion CPU<=>PCI-express<=>Carte-Graphique).

[1] La résolution effective est un multiple de dx car le modèle filtre les échelles en fonction de la résolution
[2] FFT3D = transposition de toutes les données au départ suivant l'axe des Z vers l'axe de X +FFT1D(X) + transposition de toutes les données au départ suivant l'axe des X vers l'axe de Y + FFT1D(Y) + transposition de toutes les données au départ suivant l'axe des Y vers l'axe de Z + FFT1D(Z)
[3] Si on note D ( = 3 ) la dimension de TORE Scalabilité(TORE D-dimension) = (NB_processeurs) ^ [( D-1)/D]= (NB_processeurs)^(0.6666)