\def\ie{i.~e.} \def\domega{\partial \Omega} \def\domegai{\partial \Omega_i} \def\domegao{\partial \Omega_o} \def\scal#1{\left< #1 \right>} \def\ts{t_S} \def\td{t_D} \def\ex#1#2#3#4{\chi(#1,#2,#3,#4)} \def\dep#1{\left(#1\right)} \def\depb#1{\left[#1\right]} \def\depc#1{\left\{#1\right\}} \def\dem{1/2} \def\eq#1{Eq.~\ref{eq:#1}} \def\sec#1{Section~\ref{sec:#1}} \def\dq{{\dep{\delta c}}} \def\grad{{\mbox{\bf grad}}} \def\div#1{{\mbox{div}}\dep{#1}} \def\kdiff#1{\frac{1}{\rho}\frac{\partial}{\partial z} \dep{\rho K_z \frac{\partial #1}{\partial z}}} \def\deriv#1#2{\frac{\partial #1}{\partial #2}} \def\x{{\bf x}} \def\V{{\bf v}} \def\n{{\bf n}} \def\dx{{\bf dx}} \def\dxdt{\dx dt} \def\dt#1{\frac{\partial #1}{\partial t}} \def\local{retro} \def\retrodir{\local/FIGURES} \chapter{Inversion du transport atmosphérique\label{ch:retro}} \section{Contexte} Le troisième volet de ce travail concerne l'inversion du transport atmosphérique. Il a été initié suite à une demande du CEA visant à évaluer l'efficacité du réseau global de mesure de la radioactivité atmosphérique déployé dans le cadre du Traité d'Interdiction Complète des Essais nucléaires (TICE). La surveillance des essais reposera sur quatre réseaux complémentaires. Les trois premiers mesureront les ondes sismiques, hydroacoustiques et infra-sonores et devraient permettre de localiser relativement bien les explosions pour des essais effectués respectivement sous terre, dans la mer ou dans l'atmosphère. En parallèle de ces trois technologies, un réseau global de 80 stations, en cours de déploiement, mesurera en permanence la concentration en radio-éléments dans l'atmosphère permettant de mieux caractériser le type d'évènement détecté. Toutes ces stations détecterons les aérosols radioactifs. Un sous-réseau de 40 stations détectera également les isotopes du xénon dont on pense qu'il sont relâchés en quantité significative dans l'atmosphère même lors d'essais souterrains ou sous-marins. La question qui nous était posée par le CEA était l'évaluation de la capacité de détection de ce réseau. A l'époque (été 1997), nous disposions d'une première version de LMDZT qui permettait d'effectuer des calcul de dispersion de polluants ou de radio-éléments atmosphériques. Pour répondre à la question, une approche directe aurait consisté à simuler des essais nucléaires en injectant un traceur en chaque point d'un maillage de la sphère (en chaque point de grille du modèle de transport par exemple) et de comptabiliser les fois où la concentration simulée aux différentes stations du réseau excédait le seuil de détection. On voit vite qu'une telle approche conduit à des coûts informatiques prohibitifs. Avec Robert Sadourny (LMD), nous nous sommes convaincus que ce problème pouvait être abordé en inversant le sens du temps dans notre modèle de transport Eulérien, le calcul du transport étant effectué alors en remontant à rebours dans le temps le long des trajectoires atmosphériques. L'idée d'utiliser des {\em modèles de récepteurs} pour ce type de problème n'est pas nouvelle. De façon générale, l'identification des sources pour un traceur atmosphérique (localisation spatiale et temporelle des sources ainsi que la quantification de la quantité de traceur émise) est une question très importante pour beaucoup d'aspects des sciences de l'environnement. La caractérisation des puits et sources naturels et anthropiques de CO$_2$ \cite[]{Kami:99b,Rayn:99,Bous:00,Gurn:02} est par exemple un sujet très sensible dans la perspective du contrôle des émissions de gaz à effet de serre. La surveillance de possibles rejets accidentels de pollution par des centrales nucléaires ou des installations chimiques est un autre bon exemple. Les inventaires d'émissions d'espèces chimiques qui alimentent les modèles de prévision de la pollution atmosphérique doivent aussi souvent être en partie construits en utilisant des méthodes inverses, à partir de la comparaison des prévisions du modèle et des observations \cite[voir e.~g.][]{Menu:00}. Autre exemple encore, l'interprétation fine en termes de paléoclimats des carottages effectués sur les calottes de glace en Antarctique ou au Gro\"enland, nécessite qu'on soit capable de remonter à l'origine de l'eau qui a précipité sur les calottes ou des molécules gazeuses et aérosols qui s'y sont sont accumulés. Nous nous intéressons donc dans ce chapitre à {\em l'inversion du transport atmosphérique}, c'est à dire aux approches permettant de remonter à l'identification des sources à partir de mesures de concentration dans l'atmosphère. A noter que l'inversion du transport, au-delà des nombreuses applications mentionnées ci-dessus, est un outil très intéressant pour analyser la dynamique et la chimie de l'atmosphère. Une première approche, fréquemment utilisée pour interpréter des mesures de concentration dans l'atmosphère, consiste à remonter à rebours dans le temps le long des trajectoires des masses d'air. On parle alors de {\em rétro-transport Lagrangien} ou de {\em rétro-trajectoires}. Cette technique permet d'obtenir facilement une description qualitative de l'origine de la masse d'air échantillonnée lors d'une mesure de concentration \cite[]{Hess:96,Merr:94,Ramo:96,Chia:97,Vero:92}. Comme dans les modèles lagrangiens directs, la composante turbulente du transport peut être prise en compte au moyen de perturbations aléatoires des trajectoires \cite[]{Fles:95,Vaut:01,Sieb:03}. Cependant, les calculs de rétro-trajectoires sont souvent limités à la composante grande échelle. Curieusement, alors que de nombreux modèles directs sont basés sur une description eulérienne du transport, le rétro-transport est généralement réservé à une description lagrangienne. On montre ici qu'on peut de façon générale, comme dans le monde direct, définir le {\em rétro-transport} le long des trajectoires d'air parcourues à rebours dans le temps, et lui appliquer les outils et approches utilisées habituellement pour le transport direct~: formalismes Lagrangien ou Eulérien, décomposition de Reynolds entre transport grande échelle et turbulent, schémas numériques sophistiqués pour l'advection (monotones, conservatifs, peu diffusifs, etc.). L'équation du rétro-transport se déduit de l'équation directe par des transformations simples. La diffusion turbulente, basée sur l'idée de mélange par des fluctuations symétriques du champs de vent, produit par exemple la même diffusion pour le transport direct et à rebours. Pour les paramétrisations en flux de masse de la convection, il faut en revanche inverser la direction des flux de masse dans les différents compartiments. C'est sur ces idées physiques que nous avons développé une version rétro-transport du modèle de transport LMDZT \cite[]{Hour:99CRAS} dans sa version débranchée (il s'agit d'inverser le transport de traceurs à météo connue). Le rétro-transport, ainsi défini, peut être utilisé pour des inversions quantitatives en tirant avantage d'une {\em symétrie temporelle du transport}. En absence d'autres sources ou puits, la concentration d'un traceur mesurée par un détecteur à un instant $\td$, faisant suite à l'injection d'une certaine quantité de traceur à un point source à un instant antérieur $\ts$, peut être calculée également comme la concentration au niveau de la source à l'instant $\ts$ d'un rétro-traceur émis en même quantité au lieu de détection à l'instant $\td$. Ce rétro-traceur correspond (à un facteur multiplicatif près) à la distribution de l'air qui va être échantillonné ultérieurement au niveau du détecteur. A noter qu'on ne viole pas ici l'irréversibilité du transport atmosphérique. Une plus forte dispersion ou diffusion dans le monde direct, correspondant à un gain d'entropie, se traduit également par une origine de l'air échantillonnée plus diffuse et donc également par un gain d'entropie pour le rétro-transport. Le gain d'entropie ou la perte d'information est en fait identique dans les deux cas. Une seconde approche utilisée classiquement pour l'inversion du transport atmosphérique consiste à partir d'un modèle existant du transport direct, aussi complexe soit-il, et à lui appliquer les techniques classiques de l'assimilation des observations telles que des variantes de l'interpolation optimale, ou l'assimilation variationnelle \cite[voir par exemple][]{Petr:02}, les filtres ou lissages de Kalman \cite[]{Haas:96,Zhan:99}. En particulier, les techniques adjointes, qui remontent aux travaux de Lions et Marchuk \cite[se reporter par exemple à][]{Lion:71,Marc:74,Marc:82}, fournissent une méthode systématique et puissante pour déterminer les sensibilités d'un modèle par rapport à ses variables d'état ou à des paramètres de contrôle. Ces méthodes sont généralement introduite sur des bases purement mathématiques sans référence à une quelconque signification physique. Les techniques adjointes sont utilisées pour de nombreuses applications tant météorologiques qu'océanographiques, en particulier pour l'assimilation variationnelle d'observations \cite[]{Pene:76,LeDi:86,Tala:87,Cour:87}. Les techniques adjointes ont également été utilisées dans de nombreuses études relatives à l'inversion du transport atmosphérique \cite[]{Ulia:91,Pudy:98,Kami:99a,Kami:99b,Vuki:00,Robe:91,Houw:99}. Si l'observable considérée est une mesure de concentration à un instant donné et à une station particulière, l'équation de transport adjointe est en fait un modèle "orienté récepteur" de la mesure, restreint aux processus qui vont effectivement influencer cette observation particulière. Le calcul adjoint détermine la sensibilité ou encore la {\em fonction d'influence}, qui, une fois combinée avec les sources et la concentration initiale du traceur, permet de calculer effectivement l'équivalent de la mesure. A partir de la propriété de symétrie du transport, on montre facilement que l'équation de rétro-transport, définie plus haut sur des considérations physiques, est en fait l'adjoint de l'équation du transport direct pour un produit scalaire particulier, le {\em produit scalaire pondéré par l'air} \begin{equation}\label{eq:scal1} \scal{\phi,\psi}=\int \rho \phi \psi \dxdt \end{equation} où $\phi$ et $\psi$ sont des concentrations massiques de traceurs. Dans ce cas particulier, l'adjoint d'un modèle peut donc être obtenu sur des considérations purement physiques, sans faire appel au techniques algébriques utilisées classiquement, et qui consistent en particulier à enchaîner des intégrations par partie. De plus, pour ce produit scalaire particulier, le modèle adjoint se déduit du modèle direct en changeant simplement le signe de certains termes ce qui évite en pratique certaines difficultés relatives au développement et à la maintenance des codes adjoints. Si on utilise un produit scalaire plus classique de la forme $[\phi,\psi]=\int\phi\psi \dxdt$, la symétrie est mise à mal et les équations directes et adjointes prennent des formes différentes \cite[cf. e.~g.][]{Pudy:98,Vuki:00}. \cite{Ulia:91} étaient en fait déjà arrivés sur un jeu d'équations symétriques en utilisant pourtant un produit scalaire non pondéré, mais ils utilisaient une approximation de Boussinesq pour le fluide, ce qui rend en fait les deux produits scalaires équivalents. \cite{Ulia:91} avaient également déjà remarqué l'essentiel à savoir que ``la fonction d'influence peut se calculer à partir de rétro-trajectoires particulaires quand on utilise un modèle Lagrangien. Les modèles Eulériens, gouvernés par des équations aux dérivées partielles, sont formulés dans un cadre variationnel et, dans ce cas, la fonction d'influence peut s'obtenir comme solution de l'équation adjointe à rebours dans le temps, en utilisant le récepteur comme source." Aux vues des considérations ci-dessus, on est tenté de récrire ce paragraphe sous une forme un peu plus symétrique~: la fonction d'influence peut être obtenue soit en suivant à rebours dans le temps la masse d'air échantillonnée au niveau du détecteur (rétro-transport) soit comme solution de l'équation adjointe, et ce indépendamment du cadre (Eulérien ou Lagrangien) choisi pour représenter le transport atmosphérique. La symétrie du transport à la base de cette équivalence peut cependant être perdue dans le monde numérique, c'est à dire que le modèle numérique intégré à rebours dans le temps peut être différent de l'adjoint du code numérique direct. En particulier, seuls les opérateurs linéaires par rapport à la concentration du traceur ont un adjoint et peuvent donc respecter la symétrie temporelle. Les schémas sophistiqués utilisés aujourd'hui pour l'advection des traceurs dans la plupart des modèles introduisent des non linéarités pour garantir un meilleur comportement physique (conservation, positivité, monotonie, faible diffusion). Ces schémas ne sont donc pas symétriques. C'est le cas en particulier du schéma de Van Leer utilisé dans LMDZ. On montre en revanche que les codes numériques utilisés dans LMDZ pour calculer le transport turbulent ou convectif respectent bien l'équivalence entre transport rétro et adjoint. Dans le cas où cette symétrie n'est pas respectée, la question se pose de savoir s'il est préférable d'utiliser le code adjoint du code direct ou d'intégrer à rebours dans le temps le modèle direct. Dans le cas des schémas d'advection, l'utilisation du rétro-transport garantira la positivité des rétro-panaches alors que le calcul adjoint pourra faire apparaître des valeurs négatives non physiques (une mesure de concentration pouvant augmenter à la suite d'une émission moindre). Dans certains algorithmes d'inversion, la positivité de la source peut être une contrainte importante, auquel cas on aura tendance à privilégier le rétro-transport. Pour des algorithmes de minimisation linéaires, utilisant des descentes de gradient, on peut penser a contrario que le code adjoint, qui fournit un calcul exact du gradient, doit être privilégié. On montre ici, sur des cas académiques, que cette conclusion est un peu trop rapide et que, même dans ce cas, l'utilisation d'un gradient approché mais garantissant certaines propriétés physiques peut permettre une minimisation moins poussée mathématiquement mais plus robuste. Dans ce chapitre, nous revenons en détail sur l'ensemble des aspects introduits ci-dessus. Nous commençons par introduire (\sec{transport}) le rétro-transport et la symétrie du transport atmosphérique pour des sources et détecteurs ponctuels et montrons comment le rétro-transport ainsi défini s'étend aux processus diffusifs. Nous expliquons dans la même section le lien avec l'équation adjointe et étendons la théorie à des sources diffuses. Nous explicitons en particulier ce lien entre rétro-transport et équations adjointes sur le cas des schémas de convection en flux de masse (\sec{param}). La discussion, développée d'abord dans le monde analytique, est ensuite étendue au monde numérique en présentant tout d'abord une illustration numérique du rétro-transport dans le cas de la campagne ETEX (\sec{numeric}). On démontre la symétrie des algorithmes utilisés dans LMDZ pour le transport turbulent et convectif. On étudie ensuite le cas des schémas non symétriques en se focalisant sur le schéma de Van Leer utilisé dans LMDZ (\sec{inversion}). Le chapitre se termine par quelques illustrations relatives à des applications menées autour de l'utilisation inverse du modèle LMDZ, et notamment des résultats relatifs à la surveillance des essais nucléaires. Ce travail doit beaucoup aux discussions initiales avec Robert Sadourny. La sollicitation et le financement des études par Jean-Pierre Issartel puis Philippe Heinrich (CEA/DAM) ont permis de mener le travail au-delà du simple cadre académique. Ce travail a aussi bénéficié de nombreuses et fructueuses discussions avec Bertrand Cabrit. Les tests et illustrations ETEX ont été réalisés par Abderrahmane Idelkadi pendant sa thèse. Le travail sur le TICE a vu passer plusieurs stagiaires comme Alexandre Maes et Elie Anselin. Enfin Olivier Talagrand a été d'une aide précieuse pour mener à bien la mise en forme de la partie adjointe et l'écriture d'un article en deux parties sur le lien entre transport rétro et adjoint \cite[]{Hour:05I,Hour:05II}, article sur lequel repose largement le présent chapitre.