Je remercie mon encadrant Didier Paillard pour ses discussions très enrichissantes relatives au sujet du stage mais aussi à de nombreux autres domaines. Je remercie l'équipe CLIM et les membres du LSCE pour leur accueil chaleureux.

Table des matières 1 Introduction 1.1 Français

Dans les modèles usuels d'atmosphère et d'océan, les équations de la dynamique sont intégrées sur une grille avec le plus de points possibles. Pourtant la taille des mailles n'est pas suffisamment petite pour prendre en considération des phénomènes de petite échelle comme la turbulence. Les modélisateurs sont alors obligés d'utiliser des paramétrisations dont le choix est toujours un peu subjectif. Une autre approche consiste plutôt à utiliser une hypothèse thermodynamique, la maximisation de la production d'entropie, qui permet de s'affranchir de certaines paramétrisations, et qui pourrait prédire certains phénomènes simples de grande échelle. L'objectif du stage est d'incorporer le cycle de l'eau dans un modèle simple radiatif-convectif basé sur la maximisation de la production d'entropie. Avec la seule contrainte de conservation de l'énergie (Herbert, Paillard et Dubrulle 2013), il est assez facile d'établir un profil de température, mais les flux d'énergie associées sont opposés aux gradients d'énergie et ne peuvent donc pas être interprétés en terme de convection. Il est possible d'exprimer explicitement la convection ce qui impose une nouvelle contrainte sur les flux d'énergie (Labarre, Paillard et al. 2019a). Les flux de masse sont alors infinis ce qui empêche l'eau d'être conservée. En me basant sur le code Labarre, Paillard et al. 2019a amélioré par Karine Watrin, j'ai ajouté au modèle la conservation de l'eau, ce qui permet de calculer des précipitations en zone tropicale. J'ai aussi exploré des schémas de convection différents de Labarre, Paillard et al. 2019a.

English

In atmospheric or oceanic models, equations of motion are usually integrated into a grid with the most possible number of points. However the mesh size isn't small enough to take into account small scale processes such as turbulence. Modelers have to use parametrisations whose choice is always a bit biased. Another approach is to use a thermodynamic hypothesis, the maximisation of entropy production, which allows to get rid of some parametrisations, and could predict some large scale phenomena. The purpose of the internship is to incorporate the water cycle into a simple radiativ-convectiv model based on the maximisation of entropy production. With the only constraint of energy conservation (Herbert, Paillard et Dubrulle 2013), it is quite easy to create a temperature profile, but the energy fluxes are opposed to the energy gradient and cannot be interpreted as a consequence of convection. It is possible to explicitely express convection, leading to a new constraint on energy fluxes (

Les modèles de climat et l'hypothèse MEP 2.1 Les modèles d'atmosphère et les paramétrisations

Une manière classique de modéliser l'atmosphère est d'utiliser les équations de base de la thermodynamique et de la mécanique des fluides, c'est à dire essentiellement deux équations de bilan d'impulsion (seconde loi de Newton, dans ce cas précis équations de Navier-Stokes sur l'horizontale), deux bilans de masse pour l'air et l'eau, un bilan d'énergie (premier principe de la thermodynamique), et une équation d'état (comme la loi des gaz parfaits). Ces 6 équations dépendent de 6 variables essentielles qui sont les vitesses horizontales des vents ⃗ v, la masse volumique de l'air ρ, la température T , la pression p, et le rapport de mélange d'eau q. Ces équations aux dérivées temporelles peuvent être résolues numériquement, ce qui nécessite de les discrétiser sur une grille en 3 dimensions, où chaque point de grille possède une valeur pour chaque variable du problème. Comme les ordinateurs ne possèdent qu'une capacité finie de calcul et de mémoire, le nombre de points de grille est toujours limité donc la taille des boîtes est limitée. Pour donner un ordre de grandeur, la taille (horizontale) des boîtes utilisée dans les modèles de climat des rapports du GIEC était de 500 km en 1991, et devrait être de 30 km pour le prochain rapport AR6. Cela peut descendre jusqu'à quelques km actuellement pour un modèle de météo très résolu. Quelque soit la résolution du modèle, tous les processus physiques qui s'effectuent à une échelle inférieure à la taille des boîtes ne sont pas pris en compte lors de l'intégration des équations du mouvement. Or ces processus de petite échelle contribuent à modifier le système et la méthode généralement utilisée pour les décrire est d'utiliser des paramétrisations, qui sont justifiées à la fois de façon physique et empirique. Par exemple, à cause de la turbulence, les flux de masse ou d'énergie sont mal déterminés et les modèles se contentent en général de les décrire localement en fonction des variables d'intégration. La fonction qui les relie est appelée « hypothèse de fermeture » et nécessite l'utilisation de paramètres semi-empiriques. Ces paramètres sont calibrés sur les observations, ce qui fait que leur domaine de validité n'est pas forcément très clair. De plus ces paramétrisations peuvent posséder un côté « subjectif » et ne sont d'ailleurs pas exactement les mêmes selon les modèles. L'avantage de cette méthode est qu'il est toujours possible d'affiner ou de rajouter des paramètres pour décrire de mieux en mieux chaque processus, qui peut être traité indépendemment des autres. Si l'on croit qu'une description précise des choses doit se faire par la représentation de tous les processus existants, l'accumulation des paramétrisations paraît être une bonne voie. Le désavantage est que les modèles deviennent de plus en plus complexes et que de moins en moins de gens en ont une compréhension totale. Par conséquent l'adaptation de ces modèles à des situations différentes de celles pour lesquelles ils ont été créés est très compliquée. En effet les modèles d'atmosphère des General Circulation Model ou des Earth System Model sont historiquement d'abord des modèles de météo. Ils ont plus tard été adaptés comme modèles de climat pour prédire le climat des 100 prochaines années. La pertinence de cette adaptation est déjà questionnable. En effet, a-t-on réellement besoin d'intégrer des équations à partir de conditions initiales qui seront complètement oubliées au bout de quelques dizaines d'années ? A-ton besoin de connaître toutes les variables en tout point de l'espace et du temps lorsqu'on ne s'intéresse par exemple au final qu'à la moyenne des températures de surface ? Ensuite, ces modèles ont été adaptés pour étudier les climats du passé, parfois très anciens. Or les calibrations effectuées à partir du climat présent ne sont plus forcément valides pour des climats anciens très différents de l'actuel (ni même pour le climat futur !). Cela complexifie beaucoup leur adaptation.

Une des solutions envisagées pour régler le problème des paramétrisations serait de réduire la taille des mailles. Cela pourrait venir de l'amélioration des techniques numériques pour réduire les temps de calcul ou les besoins en stockage, ou bien de l'augmentation de la puissance des super-calculateurs. Cette dernière solution est peu écologique et peu économique, et bien qu'elle puisse améliorer les résultats, elle ne règle fondamentalement pas le problème des processus sous-maille, puisqu'on est loin de pouvoir atteindre des mailles d'une taille d'un mm, l'ordre de grandeur de la dissipation visqueuse. Pour toutes ces raisons, nous avons envie de créer un modèle de climat plus simple, plus facile à comprendre et à utiliser, plus adaptable, qui se passerait de paramétrisations. Il existe déjà une grande variété de modèles :

-

Le modèle de Paltridge

Le premier modèle de climat à utiliser l'hypothèse de maximisation de la production d'entropie (hypothèse MEP) est le modèle de Paltridge [START_REF] Paltridge | Global dynamics and climate -a system of minimum entropy exchange[END_REF]). C'est un modèle méridional 1D en boîtes, radiatif-convectif, qui prédit un profil de température et un profil de flux d'énergie méridionale. Pour comprendre ce dont il s'agit, on propose d'étudier un cas plus simple à deux boîtes (figure 1a). La boîte 1 représente les basses latitudes, et la boîte 2 représente les hautes latitudes. Les basses latitudes reçoivent plus d'énergie du Soleil et possèdent un albédo plus faible. Par conséquent elles sont chauffées plus fort ce qui crée un déséquilibre qui est compensé par un flux d'énergie des basses vers les hautes latitudes. La boîte i possède une température

T i et reçoit un flux d'énergie radiative R i (T i ) = S i (1 -α i ) -σT 4
i , où S i est la constante solaire et α i est l'albédo. Entre les deux boîtes, il y a un flux d'énergie inconnu F . Ce flux serait typiquement paramétrisé dans un modèle de climat mais ici il est trouvé en maximisant la production d'entropie sous contrainte. La première contrainte est la conservation de l'énergie qui implique que la somme des flux d'énergie radiative est nulle R 1 (T 1 ) + R 2 (T 2 ) = 0. On suppose qu'on est en régime stationnaire, donc la deuxième contrainte est le bilan d'énergie dans la boîte 1 qui impose que ce qui rentre égal ce qui sort R 1 (T 1 ) = F . La production d'entropie s'écrit : observations. Ce succès triomphal a fait l'objet d'un certain nombre de critiques. D'abord, un certain nombre de paramétrisations demeurent dans le modèle de Paltridge, comme par exemple l'effet des nuages sur le bilan radiatif, ou les flux entre l'océan et l'atmosphère. Ces paramétrisations très contraignantes pourraient expliquer la grande justesse des résultats. De plus, certains phénomènes ne sont pas du tout pris en compte dans le modèle. En effet, la viscosité du fluide ou la rotation de la Terre ne jouent absolument aucun rôle. Pour un physicien de l'atmosphère, il peut paraître aberrant d'obtenir les mêmes résultats pour une Terre qui tournerait 4 fois plus vite et dont l'atmosphère serait composée de miel... Enfin une dernière critique porte sur l'hypothèse MEP elle-même, qui jusqu'à ce jour n'a jamais été rigoureusement démontrée.

σ = dS dt = F 1 T 2 - 1 T 1 = F 1 R -1 2 (-F ) - 1 R -1 1 (F ) (1)
Le problème des paramétrisations additionnelles a été corrigé dans Herbert, Paillard, Kageyama et al. 2011, en arrangeant quelques éléments et surtout en supprimant les nuages dont l'effet sur le bilan radiatif est sans toute trop compliqué pour se passer de paramétrisation. Ce modèle de Paltridge amélioré ne possède pas de paramétrisations arbitraires et peut donc être utilisé pour prédire les climats du passé. Les résultats pour le dernier maximum glaciaire sont proches de ceux trouvé avec un GCM si l'on enlève l'effet des nuages.

Le modèle de Platridge est horizontal mais il est tout à fait possible de créer des modèles verticaux. Cela a été fait dans Herbert, Paillard et Dubrulle 2013. Une critique de ce modèle est l'absence de représentation de l'énergie interne, qui fait que les flux d'énergie ne peuvent pas être interprétés comme liés à de la convection. Cela a été corrigé dans Labarre, Paillard et al. 2019a, mais ce dernier modèle n'a pas de contrainte sur la conservation de l'eau. Dans mon stage, je suis parti de son modèle, dont les méthodes de résolution numériques ont été retravaillées par Karine Watrin, et j'y ai essentiellement rajouté une contrainte de conservation de l'eau, qui permet de prédire des précipitations.

L'hypothèse MEP

Le second principe de la thermodynamique énoncé par Carnot dit seulement que l'entropie d'un système fermé isolé ne peut que croître. Réinterprété par Boltzmann et Gibbs, l'entropie s'énonce en terme de probabilité ou de fréquence. En théorie de l'information (Shannon), l'entropie devient une mesure de notre manque d'information sur un système. Si le système explore tous les états possibles de façon équiprobable, alors l'état du système réel est celui qui maximise l'entropie (MaxEnt) sous certaines contraintes, qui représentent notre connaissance du système. Ce principe de maximisation de l'entropie n'est pas un principe physique, et il peut être appliqué à n'importe quel domaine sous couvert que le système étudié explore bien toutes les possibilités de façon équiprobable, et que les contraintes sont justes. Le lien entre la théorie de l'information et la physique statistique est attribué à Jaynes. Effectivement en physique statistique pour retrouver toutes les lois fondamentales de la thermodynamique, on utilise le principe de maximisation de l'entropie (MaxEnt).

L'hypothèse de maximisation de la production d'entropie (MEP) est finalement assez similaire. En effet on peut voir la production d'entropie comme étant une entropie associée aux trajectoires et non plus aux états. L'hypothèse MEP serait donc une hypothèse d'équiprobabilité des trajectoires possibles suivies par le système. Bien que cette hypothèse n'ait jamais été parfaitement démontrée, nous l'utilisons pour voir jusqu'où il est possible d'aller avec et s'il est possible d'obtenir des résultats qui ont du sens physique. Tout comme pour MaxEnt, à priori la réussite de MEP dépend très fortement du choix des bonnes contraintes qui doivent représenter la physique essentielle du problème.

Quoiqu'il en soit, essayer de prédire le climat avec des arguments thermodynamiques est une idée intéressante. De plus il est très courant que la nature maximise quelque chose. Peut-être que MEP n'est pas la bonne méthode, mais l'idée qu'un système optimise une certaine fonction est une idée très répandue dans de nombreux domaines. Dans le cas de l'atmosphère on peut penser par exemple à la maximisation de l'énergie potentielle de convection disponible (CAPE).

Description du modèle 3.1 Le code radiatif

La circulation atmosphérique est permise par un apport constant d'énergie électromagnétique dû au Soleil. Ce rayonnement solaire est en partie absorbé par l'atmosphère et la surface de la Terre avec un facteur (1 -α), où α est l'albédo, puis ré-émis vers l'espace dans le domaine infrarouge. Ces rayonnements solaires et infrarouges sont transformés en énergie thermique dans l'atmosphère par absorption par certains composés chimiques. Le code radiatif fournit une quantification de cet apport d'énergie dans la colonne verticale atmosphérique. Il prend en entrée des mesures expérimentales de CO

R i (T, q, O 3 , CO 2 , α) = SW i (q, O 3 , α) + LW i (T, q, CO 2 ) (2)
où T, q, O 3 et CO 2 sont les profils de température, du rapport de mélange de vapeur d'eau (q = m eau /m air ), d'ozone et de dioxyde de carbone dans toute la colonne atmosphérique, i.e T = {T j } j=0,...,n ), et SW i et LW i sont les composantes solaires et infrarouges du bilan radiatif. Le code radiatif ne prend pas en compte la présence d'eau liquide, c'est à dire les nuages. De plus l'humidité absolue ou bien l'humidité relative (h = q/q s (T )) est un paramètre fixé en entrée. Si l'humidité absolue (et donc q) est fixée, il n'y a pas de rétroaction avec la vapeur d'eau. Si l'humidité relative est fixée, la seule rétroaction est celle donnée par la loi de q S (T ). Dans les deux cas, une fois les différents paramètres fixés en entrée, R = R(T ) est une fonctionnelle de la température seulement. La vapeur d'eau comme variable de l'énergie [START_REF] Labarre | Supplementary Materials to : A Radiative Convective Model based on constrained Maximum Entropy Production[END_REF] Enfin cette formulation satisfait automatiquement les principes physiques de base : le principe de retour inverse de la lumière est respecté, ainsi que le premier et le second principes de la thermodynamique. Or la force à priori de l'hypothèse MEP repose sur la bonne qualité des contraintes imposées, c'est à dire sur leur justesse physique, particulièrement importante.

La maximisation de la production d'entropie

Dans l'atmosphère divisée en N + 1 boîtes sur la verticale, chaque boîte i possède une température T i . À l'interface entre la boîte i -1 et la boîte i, on écrit F S i le rayonnement solaire net (positif vers le bas par convention), F L i le rayonnement infrarouge net (positif vers haut), et F i le flux d'énergie convective (positif vers le haut). À l'état stationnaire la conservation de l'énergie impose que la somme des flux est nulle :

F S i = F L i + F i (3)
Dans chaque boîte on peut compter l'énergie gagnée par absorption en prenant le rayonnement entrant moins le rayonnement sortant. Décomposé entre la partie visible et infrarouge le bilan radiatif s'écrit :

SW i = F S i+1 -F S i (4) LW i = F L i -F L i+1 ( 5
)
Les équations 3,4,5 permettent alors de relier le bilan d'énergie radiative R et les flux d'énergie convective F :

R i = SW i + LW i = F i+1 -F i (6)
avec la notation F 0 = F n+1 = 0, notre colonne d'atmosphère étant un système fermé (il n'y a pas d'océan ou de flux géothermique).

L'entropie associée au transfert de chaleur Q i entre une boîte i -1 et i s'écrit :

S i = Q i 1 T i - 1 T i-1 (7)
Par conséquent la production d'entropie associée aux flux convectifs, s'écrit (c'est une entropie par unité de temps) :

σ = n i=1 F i 1 T i - 1 T i-1 = n i=0 F i -F i+1 T i (6) = n i=0 - R i (T ) T i (8)
La production d'entropie peut donc s'écrire comme une fonction de la température seulement. De plus la conservation de l'énergie impose que la somme des bilans radiatifs de chaque boîte est nulle :

n i=0 R i (T ) = 0 (9)
ce qui permet d'écrire le problème variationnel sous la forme :

max n i=0 - R i (T ) T i n i=0 R i (T ) = 0 (problème ENERGIE)
Ce modèle a déjà été étudié dans la thèse de Corentin Herbert (Herbert 2012, Herbert, Paillard et Dubrulle 2013).

Le déplacement d'air -schéma en ligne

Dans le modèle précédent, on a parlé de flux d'énergie convectifs sans expliciter la façon dont les masses d'air se déplacent. Dans cette formulation, comme aucune contrainte sur F n'est donnée, ce transfert d'énergie peut avoir n'importe quelle origine, y compris une « énergie magique » dont le processus est inconnu et dont la seule raison d'être serait de maximiser la production d'entropie.

Est-il possible d'expliciter un schéma de convection ? Il faut d'abord se demander comment la convection permet de transférer de l'énergie d'une boîte à une autre. Tout d'abord chaque boîte i de température T i , d'altitude z i , et de rapport de mélange d'eau q i possède une densité massique d'énergie :

e i = C p T i + gz i + Lq i ( 10
)
où C p = 1005 J/kg/K est la capacité calorifique de l'air, g = 9.81 m/s 2 est l'accélération de la pesanteur, L = 2.5 x 10 6 J/kg est la chaleur latente de vaporisation de l'eau. Pour simplifier, le rapport de mélange d'eau q i est fixé à saturation q i = q S (T i ) (ce qui est indépendant des valeurs prises dans le code radiatif). De plus en supposant que l'air est un gaz parfait, l'altitude z i peut s'exprimer en fonction des températures des niveaux inférieurs T j,j≤i (annexe A de Labarre, Paillard et al. 2019a). Finalement l'énergie est une fonction de la température uniquement. Tout échange de masses d'air entre 2 boîtes d'énergie différente est alors associé à un transfert d'énergie, qui est toujours noté F . Selon le schéma de convection de la figure 2, entre la boîte i -1 et la boîte i, il y a un débit massique m i vers le haut, et le même débit massique m i vers le bas, de sorte que la masse totale des boîtes reste constante. Si le transport est adiabatique, c'est à dire qu'il n'y a pas de perte d'énergie pendant le déplacement, alors toute l'énergie de la parcelle d'air dans la boîte i -1 est transférée à la boîte i suivante. Le débit massique m i entre la boîte i -1 et i est donc associé à un flux d'énergie vers le haut m i e i-1 . De même, le débit massique m i entre la boîte i et i -1 est associé à un flux d'énergie vers le bas m i e i . Le flux d'énergie net entre les boîtes i -1 et i s'écrit donc :

F i = m i (e i-1 -e i ) (11)

Le déplacement d'air -schéma en étoile

Dans le schéma de convection précédent, le transport s'effectuait de proche en proche, mais dans la nature il n'y a à priori rien qui oblige la convection à s'effectuer comme cela. Par exemple lors des orages dans les cumulonimbus il existe de forts courants ascendants et descendants qui peuvent traverser la totalité de la troposphère. Pour autoriser ce genre de phénomènes et voir s'ils sont en faveur d'une plus grande production d'entropie, on fabrique un schéma en étoile dans lequel les échanges (réciproques) de masses d'air sont possibles entre toutes les boîtes de l'atmosphère. Seule la boîte i = 0 ne peut échanger qu'avec la boîte i = 1. En effet dans le modèle la boîte i = 0 d'albédo α représente une couche limite de surface, et bien qu'elle se comporte comme une boîte d'air (elle possède une énergie [START_REF] Labarre | Supplementary Materials to : A Radiative Convective Model based on constrained Maximum Entropy Production[END_REF]

F = F 0 + Cl F (14) m = m 0 + Cl m ( 15
)
maximisent aussi la production d'entropie, quelques que soient l F , l m , vecteurs de dimension µ qui appartiennent à l'espace des cycles. C est la matrice des cycles élémentaires (33) définie dans l'annexe A. Finalement, les schémas de convection en ligne ou en étoile satisfont aux mêmes équations et se résolvent de manière identique. La différence fondamentale entre les deux est l'expression de la matrice ∇ (LIGNE ou ETOILE). Les deux problèmes variationnels s'écrivent :

max n i=0 - R i (T ) T i n i=0 R i (T ) = 0 ∃ m i ≥ 0, avec R = -∇ T F et F = -m • ∇e

La conservation de l'eau

L'objectif principal du stage du stage était de représenter le cycle de l'eau. Or dans le problème CONV l'eau peut être créée et disparaître à l'infini. Donc dans un premier temps, une contrainte de conservation de l'eau a été ajoutée. La quantité de vapeur d'eau présente dans une boîte modifie son énergie [START_REF] Labarre | Supplementary Materials to : A Radiative Convective Model based on constrained Maximum Entropy Production[END_REF] de façon non négligeable. C'est donc une donnée importante du problème. Dans un premier temps, pour simplifier, le rapport de mélange est fixé à saturation q i = q S (T i ) partout. Il dépend donc uniquement de la température. Cependant, dans ce cadre, la conservation de l'eau n'est pas imposée, c'est à dire que de la vapeur d'eau peut apparaître et disparaître n'importe où. La disparition de vapeur d'eau peut s'expliquer par le phénomène de précipitations : elle se transforme en eau liquide et tombe. En revanche, la création de vapeur d'eau ne devrait se faire qu'au niveau du sol, par le phénomène d'évaporation. Cela doit se traduire par une contrainte supplémentaire. Les flux d'eau F w s'écrivent :

F w = -m∇q (16)
Alors la convergence des flux d'eau est positive dans les boîtes où il ne peut que pleuvoir :

P i (q, m) = -∇ T (m • ∇q) i ≥ 0 ∀i ∈ [1, n] (17)
Et la nature de l'opérateur ∇ T implique que les précipitations sont égales à l'évaporation :

P 0 = - n i=1 P i (18)
Le problème variationnel avec la convection et la conservation de l'eau s'écrit :

max n i=0 - R i (T ) T i n i=0 R i (T ) = 0 ∃ m i , P i≥1 ≥ 0, avec R = -∇ T F , F = -m • ∇e
(problème PRECIP) Ce problème peut se résoudre aussi bien avec le schéma de convection en ligne qu'en étoile.

Les boucles de convection

Dans les schémas de convection précédents, un flux de masse dans un sens imposait un flux de masse dans l'autre sens. Cela permettait de conserver la masse très simplement. Cependant, avec un graphe en étoile il est possible d'imaginer des flux de masse qui forment un cycle. Soit a le vecteur de dimension M des débits massiques (orientés). Alors l'équation 13 devient :

F (e, m, a) = -m • ∇e - |a| • ∇e -a • (T + S)e 2 (19)
De même les précipitations qui sont égales à la convergence des flux d'eau (17) deviennent :

P (q, m, a) = -∇ T m • ∇q + |a| • ∇q -a • (T + S)q 2 (20)
La matrice ∇ (ETOILE) et les matrices T et S sont définies dans l'annexe A. La conservation de la masse pour les débits massiques a s'écrit ∃ l a , a = Cl a , où C est la matrice des cycles élémentaires (33), et l a est un vecteur de dimension µ. Le problème variationnel s'écrit comme pour le problème CONV et le problème PRECIP, avec le changement de l'expression de F et P , et la nouvelle contrainte a = Cl a . Bien-sûr ce problème ne peut se résoudre que lorsque la matrice C est définie, c'est à dire pour le schéma en étoile.

Le rapport de mélange q

L'enjeu principal du stage était d'incorporer le cycle de l'eau dans le modèle. Pour que cela soit plus réaliste, il ne faudrait plus laisser la vapeur d'eau être à saturation, mais laisser le rapport de mélange d'eau q varier :

0 ≤ q ≤ q S (T ) (21)
Les précipitations ne peuvent alors avoir lieu qu'à saturation, ce qui impose une nouvelle contrainte, qui est que soit q = q S (T ), soit P = 0, ce qui peut s'écrire :

(q -q S (T ))P (q, m, (a)) = 0

On remarque que cette contrainte est fortement non linéaire et non convexe.

La résolution numérique 4.1 L'optimisation sous contrainte

Le problème ENERGIE, le problème CONV et le problème PRECIP sont appelés problèmes variationnels ou d'optimisation sous contrainte. Ce type de problèmes est présent dans de nombreux domaines de la physique, de la biologie, de l'ingénierie, ou de l'économie..., et il existe de nombreuses méthodes pour les résoudre numériquement. Ce type de problème se résout bien lorsque la fonction à minimiser (ou à maximiser), appelée fonction objectif, est convexe (ou concave), et lorsque les contraintes sont linéaires ou bien convexes. Lorsque l'on s'éloigne de ces cas idéaux, la fiabilité des algorithmes est en général moins bonne et les résultats sont à analyser avec précaution, car il est possible de diverger ou bien de converger vers des maxima locaux. L'algorithme converge facilement et toujours au même endroit dans le cas le plus simple (problème ENERGIE). Cependant, à partir du moment où un schéma de convection est explicité et qu'il est donc nécessaire d'imposer une contrainte sur les flux de masse m ou a, le résultat devient moins fiable. À priori les sources d'erreurs peuvent venir de tous les extrema locaux dans lesquels il est possible de tomber selon le choix des conditions initiales, ou bien de la présence de vallées dans lesquelles il est difficile de savoir quelle direction explorer et où s'arrêter. Il est en général conseillé de choisir ses variables de sorte qu'elles soient toutes proches de 1, ce qui permet d'éviter certains problèmes numériques, et aussi d'empêcher certaines vallées d'exister.

La résolution avec un algorithme « boîte noire ».

Par conséquent, au lieu de travailler avec des températures T , il est plus judicieux d'utiliser la variable x définie par (notations vectorielles, la division est terme à terme) :

x = T 0 T , où T 0 = S 4σ 1/4 (23)
où S = 1368 W/m 2 est la constante solaire, σ = 5.67 10 -8 W/m 2 est la constante de Stefan-Boltzmann. Une possibilité offerte est de linéariser R autour d'un x 0 : R(x) ≃ R x 0 (x) = rx + r 0 (24) où r est une matrice de dimension (N + 1) × (N + 1), et r 0 est un vecteur de dimension N + 1. Bien-sûr r et r 0 dépendent de x 0 . L'avantage est que la production d'entropie devient alors quadratique :

σ(x) = -x T R(x) ≃ σ x 0 (x) = -x T (rx + r 0 ) (25)
Le but n'est pas de résoudre le problème avec une production d'entropie fausse mais plutôt de faciliter l'optimisation. Pour s'approcher du résultat optimal, l'algorithme d'optimisation est à nouveau lancé mais cette fois en linéarisant autour d'un nouveau x 0 égal au x f final trouvé précédemment. Ainsi de suite, la boucle linéarisation-optimisation tourne jusqu'à ce que x converge. Ce n'est pas obligatoire de linéariser mais cela peut parfois améliorer le résultat. La variable du problème variationnel est x, ce qui fonctionne dans le cas simple (problème ENERGIE), mais lorsque m ≥ 0 (problème CONV et problème PRECIP) il faut exprimer cette contrainte en fonction de x. Pour cela, on calcule déjà F (x) en inversant la relation 6 :

F (x) = -∇ T + R(x) (26)
où ∇ T + est la matrice pseudo-inverse de ∇ T , mais F n'est défini de façon unique que dans le cas du schéma en ligne (R de dimension N + 1 plus la contrainte 9 est équivalent à F de dimension N ; dans les autres cas F est défini à des boucles près (14).). Dans ce cas, F (x) peut se calculer avec la formule :

F i (x) = i-1 k=0 R k (x) (27)
On peut alors montrer que la contrainte m ≥ 0 est équivalente à :

-F (x) • ∇e(x) ≥ 0 (28) ce qui permet de résoudre le problème CONV pour le schéma en ligne.

Pour chacun des problèmes posés, il est possible de les résoudre de différentes manières en les formulant différemment. Par exemple un changement de variable permet de supprimer la contrainte de conservation de l'énergie [START_REF] Labarre | A radiative-convective model based on constrained maximum entropy production[END_REF]. En effet après avoir linéarisé R (24), l'équation 6 peut s'inverser :

x(F ) ≃ R -1 x 0 (-∇ T F ) = -r -1 (r 0 + ∇ T F ) (29)
ce qui permet d'exprimer la production d'entropie en fonction de F uniquement :

σ(F ) = x(F ) T ∇ T F (30)
La contrainte de conservation de l'énergie disparaît, car :

n i=0 (∇ T F ) i = F 0 -F n+1 = 0 (31)
Autrement dit le système est fermé, il n'y a pas d'énergie convective qui disparaît dans l'espace ou qui provient de la surface de la Terre. L'équivalent continu est le théorème de Green-Ostrogradski, qui dit que l'intégrale sur un volume d'une divergence est égale aux termes de bords, qui sont nuls ici. La contrainte de positivité des flux de masse s'écrit alors :

-F • ∇e(x(F )) ≥ 0 (32)
Encore plus de changements de variables ont été codés par Karine Watrin et permettent de résoudre chacun des problèmes de façon différente. Par exemple la variable xm est composée de la concaténation de x et de m, ce qui permet d'exprimer directement la contrainte m ≥ 0, mais cela nécessite de rajouter la contrainte qui correspond à la relation entre x et m. De même on peut concaténer m et F pour créer la variable mF , ou concaténer x, m et F pour créer la variable xmF . Les méthodes de résolution et les contraintes associées à chacune d'entre elles sont résumées dans la table 1, dans le cas du problème CONV. On constate que lorsque les variables sont xm ou xmF , la production d'entropie peut être exprimée de deux façons différentes. Cela peut sembler un détail mais cela modifie le chemin suivi par l'algorithme d'optimisation et donc le résultat. La contrainte de conservation de l'eau (problème PRECIP) n'est pas écrite dans le tableau, mais le principe de changement de variable est le même que pour les autres contraintes. Les précipitations P dépendent de m et de q s (x), et on obtient m(e(x), F ) directement de [START_REF] Lacis | A Parameterization for the Absorption of Solar Radiation in the Earth's Atmosphere[END_REF], puis F (x) (26) ou x(F ) (29) selon la variable souhaitée. On constate aussi que le schéma en étoile ne peut pas se résoudre avec x car il faudrait calculer F (x). De plus lorsque m et a sont présents, ce n'est pas possible de résoudre avec F car il dépend alors de m et de a et la contrainte m ≥ 0 ne peut plus s'écrire en fonction de F . variable production d'entropie contrainte m ≥ 0 autres contraintes

x σ = -x T R(x) -F (x) • ∇e(x) ≥ 0 n i=0 R i (x) = 0 F σ = x T (F )∇ T F -F • ∇e(x(F )) ≥ 0 xm σ = -x T R(x) m ≥ 0 R(x) = ∇ T (m • ∇e(x)) σ = -x T ∇ T (m • ∇e(x)) mF σ = x T (F )∇ T F m ≥ 0 F = -m • ∇e(x(F )) xmF σ = -x T R(x) m ≥ 0 R(x) = -∇ T F F = -m • ∇e(x) σ = x T ∇ T F Table 1
D'autres possibilités sont offertes par le code de Karine Watrin, comme par exemple multiplier m par un coefficient α de sorte à utiliser la variable z = αm au lieu de m, ou bien imposer un flux de masse maximal m ≤ m max pour restreindre la taille de l'espace des phases dans le cas où m est une variable explicite.

Toutes ces méthodes de résolution ont été codées afin d'effectuer des tests pour essayer de tendre vers de meilleurs résultats. En effet il n'est pas possible d'utiliser aveuglément l'algorithme d'optimisation pour n'importe quelle variable et n'importe quelle condition initiale et espérer trouver le résultat optimal à chaque fois. En fait l'algorithme garantit de trouver un maximum local (ou rien du tout), mais il n'y a aucun moyen de savoir s'il est global ou non. L'expérience montre que selon les options choisies, les résultats diffèrent légèrement, ou peuvent même être aberrants dans certains cas. Une façon simple de vérifier si un résultat est « bon » est de regarder la production d'entropie : si elle est plus grande le résultat est meilleur. Un exemple est donné figure 3, dans lequel le problème CONV_LIGNE est résolu avec 3 variables différentes : x, F , et xmF . La linéarisation de R est utilisée pour F et xmF mais pas pour x. Pour le profil de température, les résolutions avec F et xmF sont très proches, mais la résolution avec x diffère des autres. Dans ce cas c'est lié au choix de linéariser ou non. En terme de production d'entropie, la différence est de l'ordre de 0.1%, mais pour la température la différence peut aller jusqu'à 1%, ce qui correspond à une différence de température de 3 K. Même entre les 2 résolutions avec F et avec xmF pour lesquelles la différence de production d'entropie n'est que de 0.01%, et où les courbes de température semblent confondues, en réalité les différences entre les deux courbes vont jusqu'à ±1 K. Il est intéressant de constater que dans l'article de Une bonne partie de mon stage a consisté à tester pour chaque problème posé différentes options de résolution et à chercher la meilleure.

Résultats

J'ai essayé de choisir dans chaque cas les figures pour lesquelles j'avais trouvé la plus grande production d'entropie.

Il est possible d'utiliser dans le code radiatif différents profils atmosphériques issus des mesures de concentrations de CO 2 , d'H 2 O et d'O 3 dans des colonnes atmosphériques de différents endroits du globe [START_REF] Mcclatchey | Optical Properties of the Atmosphere (Third Edition)[END_REF]. En pratique il a été montré dans Labarre, Paillard et al. 2019a que pour les hautes latitudes, les résultats sont très éloignés de la réalité à cause de l'importance des flux d'énergie horizontaux, qui n'existent pas (encore) dans notre modèle 1D. C'est pourquoi tous les résultats présentés ici le sont pour une colonne tropicale.

Le schéma en ligne

Sur la figure 4 sont superposés les résultats du problème ENERGIE, du problème CONV et du problème PRECIP pour le schéma en ligne. De façon très logique, le problème ENERGIE possède la plus grande production d'entropie σ = 53.9 mW/m 2 /K, puisque c'est le problème le moins contraint. Ensuite vient le problème CONV avec σ = 44.3 mW/m 2 /K, puisqu'une contrainte sur les flux de masse m a été rajoutée, et enfin vient le problème PRECIP avec σ = 41.1 mW/m 2 /K, où une contrainte sur la conservation de l'eau est rajoutée. Aux courbes de températures ont été rajoutées les « observations » de McClatchey 1972. Déjà le problème CONV et le problème PRECIP sont beaucoup plus proches des observations que le problème ENERGIE. Notamment les valeurs dans la haute atmosphère et l'inversion de température liée à la couche d'ozone sont beaucoup plus réalistes. Globalement l'atmosphère est chauffée par le bas, ce qui crée un déséquilibre qui est compensé par des flux d'énergie vers le haut. C'est ce qui est observé figure 4b. Dans le problème ENERGIE, bien que cela n'ait pas de sens puisque cela conduit à m < 0, on peut définir une énergie [START_REF] Labarre | Supplementary Materials to : A Radiative Convective Model based on constrained Maximum Entropy Production[END_REF] et des « flux de masse » associés. On observe alors figure 4c que l'énergie décroît jusqu'à 700 mB, à cause de la baisse de température, puis croît, essentiellement à cause du potentiel de pesanteur. Lors de cette inversion du gradient d'énergie, comme les flux d'énergie ne changent pas de signe, les « flux de masse » deviennent négatifs. Par conséquent les flux d'énergie du problème ENERGIE ne peuvent pas être interprétés en terme de convection (telle que décrite dans le schéma en ligne), et d'ailleurs ils sont indépendants de l'énergie interne. Dans les deux autres problèmes, puisque m ≥ 0, pour permettre aux flux d'énergie de continuer d'être positifs non nuls, le gradient d'énergie ne s'inverse que beaucoup plus haut, vers 250 mB. C'est permis par la modification du gradient de température qui est plus important et qui contrecarre l'effet de la pesanteur. Au-dessus de 250 mB, les flux d'énergie, obligés d'être opposés aux gradients d'énergie, sont alors obligés d'être négatifs ou nuls. Lorsque l'atmosphère est chauffée par le bas, ils n'ont aucune raison d'être négatifs, et sont donc nuls. Tout en haut au dessus de 100 mB elle est chauffée par le haut à cause de l'ozone O 3 , et les flux d'énergie sont donc négatifs. C'est ce qu'on observe pour le problème problème CONV figure 4b. À cet endroit les flux de masse m ne sont pas nuls, mais justes très petits et compensés par le fort gradient d'énergie. En revanche pour le problème PRECIP, la conservation de l'eau force la nullité des flux d'énergie et de masse à partir 22) est fortement non linéaire et c'est probablement la raison de la difficulté à converger vers un maximum. Peut-être que nous manquons d'une astuce mathématique ou d'une « idée de génie » pour poser le problème d'une façon différente qui le rendrait bien plus facile à optimiser. À ce stade il y aurait deux directions à explorer, l'une serait mathématique et consisterait à cerner de façon plus précise ce qu'il se passe pendant l'optimisation, l'autre serait plus physique et consisterait à réfléchir à ce qu'est la physique essentielle du problème, quelles sont les contraintes qui manquent ou quelles sont celles en trop, s'il manque ou s'il y a trop de degrés de liberté. L'ajout dans le schéma en étoile de boucles de convection a été testé (section 3.6), mais d'une part les résultats de l'optimisation sont moins fiables, d'autre part ils ne semblent pas modifier radicalement la production d'entropie. À priori les boucles de convection seraient négligeables, mais il faut rester prudent car on a vu qu'un détail au niveau de la production d'entropie pouvait ne pas en être un pour certaines variables comme les précipitations.

Le schéma en étoile

Le CO 2

Des tests de sensibilité climatique (changement de température lié à un doublement de CO 2 ) ont été effectués pour le problème ENERGIE, le problème CONV et le problème PRECIP avec un schéma de convection en ligne (figure 6). Le niveau de CO 2 par défaut est celui pré-industriel de 280 ppm. Deux changements de CO 2 ont été testés, celui de 180 ppm du dernier maximum glaciaire et celui qui correspond au test classique de doublement de CO 2 (560 ppm). Pour interpréter les résultats et les comparer à ceux des GCMs, il faut garder en tête que de nombreuses rétroactions qui sont contenues dans les GCMs ne le sont pas ici, comme celle de l'albédo, des nuages, et toutes celles liées aux couplages avec les autre composantes que l'atmosphère. Dans le problème CONV et le problème PRECIP, la rétroaction principale est celle due à l'augmentation de la quantité de vapeur d'eau avec la température, qui est un gaz à effet de serre (car c'est l'humidité relative qui est fixée dans le code radiatif). Si l'on s'intéresse uniquement à la température à 1m50 du sol comme c'est le cas pour la sensibilité climatique, la modification du profil de température devient une rétroaction importante. Si l'on prend la boîte 1 pour évaluer la sensibilité climatique (puisque la boîte 0 n'est qu'une couche limite), alors la sensibilité climatique est de 1.1 K pour le problème ENERGIE, 0.7 K pour le problème CONV et 1.0 K pour le problème PRECIP. Elles sont toutes assez proches de 1 K. Par comparaison la sensibilité climatique des GCMs est en général comprise entre 1.5 et 4.5 K, certainement à cause des nombreuses rétroactions représentées en plus. En terme de profil de température, celui qui ressemble le plus qualitativement aux résultats des GCMs est le problème CONV (voir RCP4.5, CMIP5 du GIEC). Le problème PRECIP_LIGNE est certainement trop contraint et il faudrait lui rajouter des degrés de liberté comme un schéma en étoile ou une humidité non fixée à saturation.

Conclusion

Dans le monde de la recherche sur les modèles de climat, la majorité de l'effort actuel s'effectue sur l'amélioration des GCMs (General Circulation Model) ou des ESMs (Earth System Model), dont la partie atmosphère ou océan est basée sur la discrétisation des équations de la dynamique. Ces modèles ne peuvent se passer de paramétrisations, qui sont la source d'un certain nombre d'incertitudes, et dont le domaine de validité n'est pas toujours très clair. Ce sont ces modèles qui sont utilisés pour rédiger les rapports du GIEC. Une autre approche consisterait à chercher directement juste la physique nécessaire pour prédire le climat. Des principes variationnels peuvent être utilisés pour se passer de paramétrisations. Parmi eux, la maximisation de la production d'entropie semble prometteuse. Utilisée pour la première fois dans un modèle de climat dans Paltridge 1975, elle a ensuite été relativement abandonnée à la suite des critiques du modèle. Plus récemment, de nouveaux modèles qui corrigeaient certains problèmes ont été créés. Pour une atmosphère verticale, l'hypothèse MEP a été testée dans un modèle radiatif-convectif, dont le code radiatif est plus précis (Herbert, Paillard et Dubrulle 2013, problème ENERGIE). De la convection a été rajoutée dans Labarre, Paillard et al. 2019a (problème CONV). Dans mon stage, l'ajout d'une contrainte sur les précipitations (problème PRECIP) permet de prédire des précipitations cohérentes avec ce qui est observé en milieu tropical. Autoriser la convection à être plus libre et à effectuer de grands mouvements verticaux favorise le transfert d'eau directement du sol vers le haut la troposphère. Malheureusement, cela conduit à des précipitations qui tendent vers l'infini, ce qui indique qu'il manque certainement une contrainte, peut-être sur la mécanique ou sur la stratification (il est discutable d'échanger des masses d'air dont le profil est stable (e haut > e bas ), comme cela se passe dans le schéma en étoile ou dans la stratosphère du problème CONV_LIGNE. Peut-être qu'il faudrait rajouter une contrainte du type F ≥ 0).

Les possibilités futures pourraient être d'ajouter un cycle saisonnier, ou bien de faire un modèle en 2 dimensions. Rendre le modèle plus précis et pouvoir le comparer à des modèles plus compliqués serait intéressant. Une idée serait aussi de l'utiliser comme hypothèse de fermeture pour « tuner » un GCM. Enfin, une meilleure compréhension des mathématiques sous-jacentes au problème d'optimisation pourrait conduire à une amélioration des résultats et une plus grande facilité à modifier le modèle ultérieurement.

A Théorie des graphes

Un graphe est un ensemble de sommets reliés entre eux par des côtés. Dans le cas présent les sommets représentent les boîtes d'atmosphère et les côtés représentent les flux entre les boîtes. Deux graphes ont été étudiés et sont représentés figure 7. Si N + 1 est le nombre de sommets, alors le nombre de côtés est égal à M = N pour le graphe en ligne, et M = 1 + N (N -1)/2 pour le graphe en étoile. Un graphe peut être vu comme l'équivalent discret d'un espace continu, et par analogie l'ensemble des sommets est appelé champ scalaire, et l'ensemble des côtés est appelé champ vectoriel. Tout comme dans un espace continu, il peut être nécessaire de définir des relations (gradient, divergence, rotationnel) entre les champs scalaires et vectoriels. Pour cela des matrices sont utilisées. Les matrices S (start) et T (target) de tailles M × N transforment un champ scalaire en un champ vectoriel. La matrice S relie les sommets aux côtés qui en partent, et la matrice T relie les sommets aux côtés qui en arrivent. Mathématiquement, S i,j = 1 si le côté i a pour point de départ le sommet j, 0 sinon ; et T i,j = 1 si le côté i a pour point d'arrivée le sommet j, 0 sinon.

La matrice gradient est alors définie par ∇ = T -S. Soit un champ scalaire x, alors ∇x est un champ vectoriel qui vaut (x j -x i ) {i,j}∈{côtés} , ce qui correspond bien à ce qui est attendu d'un gradient. Avec un graphe en ligne, le gradient ∇ s'écrit : La transposée de la matrice gradient est appelée matrice de convergence ∇ T , c'est l'opposé de la divergence. Soit un champ vectoriel v = (v ij ) {i,j}∈{côtés} , alors ∇ T est un champ scalaire défini par (∇ T ) i = j u ij ϵ i , où ϵ i vaut 1 si le sommet i est un point d'arrivée de u ij , et vaut -1 si c'est un point de départ. Cela correspond bien à ce qui est attendu d'une convergence. Dans certains graphes la structure des côtés permet de définir des cycles. C'est le cas du graphe en étoile. Tout comme pour les sommets et les côtés, les cycles peuvent être numérotés. Ils sont décomposés en cycles élémentaires qui forment une base de l'ensemble des cycles. La matrice des cycles élémentaires C donne les côtés qui appartiennent à chaque cycle élémentaire :

∇ =        
C p,q =    0
si le côté p n'appartient pas au cycle q 1 si le côté p appartient au cycle q avec la même orientation -1 si le côté p appartient au cycle q avec l'orientation opposée

C =           0 0 0 1 0 0 -1 1 0 0 -1 0 1 0 1 0 0 -1 0 1 1           (34)

  Labarre, Paillard et al. 2019a). Mass fluxes are then infinite, avoiding water to be conserved. With the code of Labarre, Paillard et al. 2019a upgraded by Karine Watrin, I added the water conservation, allowing one to calculate precipitations in the tropical zone. I also explored other convection patterns, different than in Labarre, Paillard et al. 2019a.

  Assez logiquement, on voit figure 1b que la production d'entropie est nulle lorsque le flux d'énergie F est nul, et qu'elle est nulle lorsque les 2 boîtes sont isothermes. Elle est bien positive lorsque le flux F est dirigé du chaud vers le froid. Entre les deux valeurs nulles, il existe un maximum de production d'entropie. En plus d'obtenir un flux d'énergie F , on obtient très simplement les températures T 1 et T 2 en inversant R 1 et R 2 .

( a )

 a Schéma du modèle de Paltridge simplifié à 2 boîtes. La boîte 1 représente les basses latitudes, et la boîtes 2 les hautes latitudes. Les R i (T i ) sont les flux d'énergie radiative, et F est l'énergie d'échange enter les 2 boîtes. (b) Production d'entropie σ en fonction du flux d'énergie F , dans le modèle de Paltridge simplifié à deux boîtes.

Figure 1 Le

 1 Figure 1

  2 , d'O 3 et de vapeur d'eau dans la colonne atmosphérique (McClatchey 1972), l'albédo α de la surface, la constante solaire S c , et la température T dans la colonne. L'ozone O 3 absorbe dans le spectre solaire, le CO 2 dans l'infrarouge, et la vapeur d'eau H 2 O dans les deux spectres. Plus clairement, dans le code l'atmosphère est divisée en N + 1 boîtes entre le sol et l'espace, et le bilan radiatif R i (en W/m 2 ) dans la boîte i s'écrit :

Figure 2 -

 2 Figure 2 -Figure issue de Labarre, Paillard et al. 2019a. Atmosphère verticale divisée en N + 1 boîtes. La couche 0 est infiniment fine et possède un albédo α. Les R sont les flux d'énergie radiative, les F sont les flux d'énergie convective entre les boîtes. Les m sont les flux de masse réciproques entre chaque boîte. Chaque boîte i possède une température T i , une élévation z i et énergie massique interne e i (10). Les vecteurs R, F , m et e sont reliés par les équations 12 et 13.

  liée à sa température et sa quantité d'eau), on considère qu'elle n'est pas la base de grands courants d'airs. Les schémas en ligne et en étoile sont représentés figure 7 en annexe, sous forme de graphe. Pour utiliser le vocabulaire de la théorie des graphes, les boîtes d'atmosphère sont des sommets, et les flux sont des côtés. Les températures T et les énergies e sont des champs scalaires, et les flux d'énergie F et les débits massiques m sont des champs vectoriels. Les équations 6 et 11 peuvent se réécrire de façon compacte : R = -∇ T F (12) F = -m • ∇e (13) où • est un produit matriciel terme à terme (produit matriciel de Hadamar), et ∇ est la matrice gradient, dont la définition dépend du graphe. Elle est donnée dans l'annexe A pour le graphe en ligne (LIGNE) et pour le graphe en étoile (ETOILE). La matrice ∇ T est la transposée de ∇, et elle est appelée matrice de convergence (c'est l'opposé de la divergence). La matrice ∇ est de taille M × (N + 1), et M = N dans le schéma en ligne, alors que M = 1 + N (N -1)/2 dans le schéma en étoile. Cette différence introduit une « dégénérescence » des vecteurs F et m. En effet, lorsque des flux d'énergie forment une cycle, ils ne contribuent en rien à la production d'entropie (c'est lié aux équations 12 et 35). Par conséquent, si F 0 et m 0 maximisent la production d'entropie, alors

(

  problème CONV) Ce modèle a été étudié avec un schéma en ligne dans la thèse de Vincent Labarre (Labarre, Paillard et al. 2019a, Labarre 2020). Une partie du travail de stage a consisté à construire et étudier ce modèle avec le nouveau schéma en étoile.

  Labarre, Paillard et al. 2019a, la variable de résolution était F et que par conséquent, la figure 3 montre un résultat meilleur grâce la résolution en x sans linéariser.

Figure 3 -

 3 Figure 3 -Exemple des différences de résultats lorsque des méthodes de résolution différentes sont utilisées, pour le problème CONV_LIGNE avec N = 20. La courbe pleine en rouge correspond aux observations. 1) Résolution avec x, x ini = 1, pas de linéarisation. 2) Résolution avec F , x ini = 1, F ini = 1, avec linéarisation (obligatoire avec F ). 3) Résolution avec xmF , (x, m, F ) ini = (1, 0.15, 0.09), m max = ∞, z = 0.1m, sans linéarisation. Dans l'ordre σ = 44.304/44.258/44.255 mW/m 2 /K.

Figure 4 -

 4 Figure 4 -Schéma en ligne. N = 20. 3 cas différents : 1) problème ENERGIE, 2) problème CONV et 3) problème PRECIP. 1) x ini = 1. 2) Résolution avec x, x ini = 1, sans linéarisation. 3) Résolution avec xmF , (x, m, F ) ini = (1, 0.15, 0.09), avec linéarisation, m max = 0.33, z = 0.1m. Dans l'ordre σ = 53.917/44.304/41.074 mW/m 2 /K.

  Le problème CONV_ETOILE et le problème PRECIP_ETOILE permettent de rajouter des degrés de liberté au système. Bien-sûr cela rajoute des difficultés de résolution puisque le nombre de flux devient plus grand, M = 1 + N (N -1)/2 au lieu de M = N . C'est pourquoi on se contente de résoudre avec N = 15 (Note : La production d'entropie est plus grande lorsque le nombre de boîtes augmente, car le nombre de flux augmente. Des tests pour des cas simples, jusqu'à N = 80, semblent indiquer que σ converge quand N → ∞. Cela ne semble pas idiot puisque les différences de température entre les boîtes se réduisent. Les comparaisons sont donc toujours faites à nombre de boîtes égales.). Sans surprise les productions d'entropie trouvées sont plus grandes pour le schéma en étoile que pour le schéma en ligne : σ = 44.449 mW/m 2 /K contre σ = 41.289 mW/m 2 /K pour le problème CONV, et σ = 44.447 mW/m 2 /K contre σ = 38.934 mW/m 2 /K pour le problème PRECIP. Pour comparer, la production d'entropie la plus grande possible est

Figure 5 -

 5 Figure 5 -Schéma en étoile. N = 15. 1) problème ENERGIE pour comparer aux autres. 2) problème CONV, résolution en xmF , avec linéarisation, (x, m, F ) ini = (1, 0.15, 0.09), m max = 0.33, z = 0.1m. 3) problème PRECIP, résolution en xmF , avec linéarisation, (x, m, F ) ini = (1, 0.15, 1),m max = 0.33,z = 1m. Pour obtenir une plus grande production d'entropie (problème PRECIP), le flux de masse de 0 à 1 a été mis à la main vers l'infini. Dans les deux cas l'échelle pour les flux d'énergie a été tronquée car F 0->1 est trop grand. Dans les deux cas les précipitations divergent. Dans l'ordre σ = 52.4/44.449/44.447 mW/m 2 /K.

  280/560 = 43.0/44.3/46.0 mW/m 2 /K. (c) problème PRECIP_LIGNE. σ 180/280/560 = 39.6/41.1/43.1 mW/m 2 /K.

Figure 6 -

 6 Figure 6 -Différence de profil de température entre ceux avec un CO 2 à 180 ppm ou 560 ppm, et celui standard avec un CO 2 de 280 ppm.

  Graphe en ligne. M = 4.

  Graphe en étoile. M = 6.

Figure 7 -

 7 Figure 7 -Les deux graphes avec N = 4.

C

  est une matrice de taille M × µ, où µ = 1 + M -(N + 1) est le nombre de cycles élémentaires. Un exemple d'écriture de C pour le graphe en étoile et N = 4 est :

  Les ESMs (Earth System Model) ont pour but d'être le plus précis possible dans la représentation d'un maximum de phénomènes. Ils peuvent contenir un modèle d'atmosphère, d'océan, de surface continentale, de glace de mer, de calotte, de chimie, de biogéochimie marine et/ou de cycle du carbone. Le modèle étudié est un modèle radiatif-convectif, basé sur la maximisation de la production d'entropie, mais contrairement aux EBMs classiques et aux autres modèles, celui-ci se passe de paramétrisations.

Les modèles conceptuels sont composés de seulement quelques équations simples et sont utilisés pour rendre compte d'un phénomène particulier, par exemple les cycles glaciaires-interglaciaires. -Les EBMs (Energy Balance Model) sont uniquement thermodynamiques et sont basés sur la loi de conservation de l'énergie (premier principe). Ils peuvent être à 0D, 1D ou 2D. -Les modèles radiatif-convectifs sont des EBMs qui prennent en compte l'énergie radiative et une énergie due à la convection. Ils peuvent prendre en compte le cycle de l'eau. -Les modèles statistique-dynamiques à 2D sont une combinaison d'un EBM et de relations statistiques pour connaître la direction et la vitesse des vents. -Les GCMs (General Circulation Model) avec une atmosphère et/ou un océan, sont basés sur la discrétisation des équations de la mécanique des fluides, comme décrit dans le paragraphe précédent. Parfois l'atmosphère et l'océan sont couplés. -Les EMICs (Earth Model of Intermediate Complexity) sont des GCMs dont la résolution atmosphérique et/ou océanique a été réduite pour pouvoir effectuer des simulations sur plusieurs centaines ou milliers d'années, ou effectuer de nombreuses expériences de sensibilité. Ils peuvent être couplés à des modèles de végétation, de biogéochimie, de calotte... pour prendre en compte les processus à long terme. Ils sont très utilisés pour étudier les climats du passé. -

  Dufresne et al. 2005). Un modèle NEF fonctionne à partir des mêmes lois physiques de base qu'un modèle classique avec description locale des flux d'énergie, mais leur formulation est globale, c'est à dire que chaque apport d'énergie dans une couche d'atmosphère i est exprimé en fonction de toutes les contributions individuelles des différentes couches atmosphériques. Les composantes visibles et infrarouges sont calculées de 2 façons différentes. Dans le domaine infrarouge, le spectre est décomposé en 22 bandes. Le coefficient d'absorption est calculé avec le modèle de Goody 1952, à partir des données de Rodgers et Walshaw 1966. Dans ce cadre, l'approximation diffusive est réalisée avec un coefficient de diffusion standard µ = 1.66. Dans le domaine visible, la paramétrisation utilisée est celle de Lacis et Hansen 1974. L'intérêt du modèle NEF est multiple. Il est pour commencer plus précis qu'un modèle simple avec atmosphère grise dont les résultats avec MEP montraient une instabilité verticale (Herbert, Paillard et Dubrulle 2013). Ensuite, le temps de calcul (de R) est suffisamment rapide, et le résultat est dérivable et assez lisse pour que la recherche du maximum de production d'entropie soit numériquement efficace, contrairement aux codes radiatifs standards qui fournissent des résultats non dérivables et avec beaucoup de maxima locaux rendant donc inefficace l'algorithme d'optimisation.

dans le problème MEP est découplée de sa valeur dans le code radiatif. En pratique on fixera toujours l'humidité relative. Ces approximations sont faites pour simplifier l'interprétation des résultats et éviter davantage de problèmes de résolution, sans pour autant enlever au modèle l'essentiel de la physique.

Plus précisément le code radiatif utilisé (Herbert, Paillard et Dubrulle 2013) est basé sur une formulation en puissances nettes échangées (NEF)

[START_REF] Green | Division of radiative streams into internal transfer and cooling to space[END_REF]

  Naturellement, les débits massiques sont positifs, m i ≥ 0. Cela fournit une nouvelle contrainte sur F , et par extension, sur T et sur la production d'entropie, au travers des équations 6 et 8. Dans les modèles classiques de climat, les valeurs de F et m sont typiquement paramétrisées après résolution des équations du mouvement. Ici au contraire ce sont des variables du problème qui s'ajustent toutes seules pour maximiser la production d'entropie. Une conséquence d'un tel schéma de convection est que les flux d'énergie se font de l'endroit le plus énergétique vers le moins énergétique. D'autres schémas de convection sont envisageables, ce qui modifie les contraintes (voir sections 3.4 et 3.6). Le problème ENERGIE n'est basée sur aucun schéma de convection, c'est donc le problème le moins contraint, et donc celui pour lequel la production d'entropie est la plus grande possible.

  Il existe beaucoup d'algorithmes d'optimisation sous contrainte déjà codés et implantés dans de nombreuses librairies de programmation. Une grande partie d'entre eux fonctionne sur la méthode des multiplicateurs de Lagrange, et effectuent ensuite une descente de gradient à partir d'un point initial (d'où la nécessité d'un code radiatif lisse !). L'algorithme que nous avons utilisé est de ceux-ci, c'est un algorithme d'optimisation quadratique successive, décrit précisément dans Kraft 1988. Concrètement, il est utilisé comme une « boîte noire » fournie par la bibliothèque python scipy.optimize.minimize (méthode SLSQP, table 1 de Virtanen et al. 2020). L'algorithme prend en entrée d'une part la fonction objectif et des contraintes d'égalité et d'inégalité, qui doivent toutes être écrites en fonction d'une même variable, et d'autre part des conditions initiales. Un code python a été crée par Karine Watrin autour de cet algorithme, qui permet de résoudre chaque problème avec différentes méthodes de résolution, essentiellement par des changements de variable. Son code est lui-même inspiré de Labarre, Paillard et al. 2019b.

  de 250 mB à se propager jusqu'en haut de l'atmosphère. Dans le problème CONV, les gradients d'énergie sont nuls entre 700 et 250 mB, ce qui correspond à une évolution de température qui suit le gradient adiabatique humide. Cette zone correspond à des flux de masse m infinis, ce qui signifie que l'air est parfaitement mélangé (les valeurs infiniment négatives figure 4d sont justes des problèmes numériques de division par 0). Par suite, les précipitations sont infiniment négatives ou positives (figure4e). En revanche, dans le problème PRECIP, la conservation de l'eau empêche le mélange infini des masses d'air et donc le gradient d'énergie cesse d'être nul. Les flux de masse sont alors finis. Ils permettent le transport de vapeur d'eau du sol jusqu'à 300 mB, et comme le gradient de vapeur d'eau devient plus faible avec l'altitude, les flux de masse augmentent pour conserver l'eau. À partir de 300 mB, il n'y a plus assez d'énergie pour monter plus haut, et l'eau ne peut plus que précipiter. C'est ce qu'on observe figure4e, où les précipitations prédites sont de 1.3 m/an, ce qui est dans le bon ordre de grandeur des précipitations en région tropicale. Ce résultat est le résultat majeur du stage : finalement, avec seulement un code radiatif et une hypothèse physique assez simple, il est possible de prédire des précipitations qui ont du sens. À notre connaissance c'est la première fois que des précipitations sont prédites avec un modèle MEP.L'étape suivante aurait été d'ajouter des degrés de liberté au système en autorisant le rapport de mélange d'eau à varier (équations 21 et 22). Mais pour l'instant, l'algorithme d'optimisation n'arrive pas à converger, il se retrouve souvent bloqué dans des zones non faisables où les contraintes sont incompatibles, alors même qu'il part d'une zone faisable. Mathématiquement, la nouvelle contrainte (

  4/44.449/44.447 mW/m 2 /K. celle du problème ENERGIE, où σ = 52.4 mW/m 2 /K. De même les profils de température (figure5a) et d'énergie (figure5b) sont plus proches de ceux du problème ENERGIE. Comme il y a plus de flux, les gradients de e, ∇e, sont plus nombreux, il est plus facile de les avoir négatifs : il suffit que le point de départ soit plus énergétique que le point d'arrivée, ce qui est permis en de nombreux endroits par le profil en U de l'énergie. Cela permet l'existence de flux d'énergie positifs qui traversent l'atmosphère. C'est ce qu'on observe figures 5c et 5d, des flux d'énergie partent des boîtes 1 à 2 (1000-900 mB) et vont vers les boîtes 10 à 12 (400-200 mB). Mais le flux d'énergie le plus important est celui qui va de la boîte 0 (couche limite de surface) à la boîte 1. Il est tellement plus grand que les autres qu'il a été tronqué sur les figures 5c et 5d, il dépasse en réalité l'échelle de couleur. Il vaut F 0→1 = 82 W/m 2 pour le problème CONV et F 0→1 = 102 W/m 2 pour le problème PRECIP, ce qui est 4 à 5 fois plus grand que le maximum des autres flux d'énergie. Mais c'est dans le même ordre de grandeur que pour le schéma en ligne. Finalement, le schéma en étoile en autorisant beaucoup plus de flux dilue la force de chacun, et leur permet aussi d'être plus efficaces. Les flux de masse ne sont pas affichés mais peuvent être décrits qualitativement. Déjà les flux de masse ont été imposés inférieurs à 0.33 kg/m 2 /s parce que l'optimisation est meilleure dans ce cas, mais la raison est plus numérique que physique. Le problème CONV ne possède qu'un seul flux de masse vraiment significatif, c'est celui qui échange l'air entre la boîte 1 et la boîte 12, m 1→12 = 0.33 kg/m 2 /s. Les deux boîtes ont alors une énergie quasiment identique. Les autres flux de masse sont en général plus faibles et ne sont significatifs qu'aux endroits où ∇e est nul. Le problème PRECIP est un peu plus subtil. Les flux de masse ont été imposés inférieurs à 0.33 kg/m 2 /s dans un premier temps, puis dans un deuxième temps le flux de masse de la boîte 0 vers 1 qui valait 0.33 kg/m 2 /s a été manuellement positionné vers l'infini (ici m = 10 suffit). L'algorithme d'optimisation est alors à nouveau lancé avec ces nouvelles conditions initiales (et m max = ∞) et non seulement ce flux de masse reste quasiment le même, mais la production d'entropie devient légèrement plus grande ! Faire tendre m 0→1 vers l'infini produit donc plus d'entropie. Les autres flux de masse deviennent alors similaires au cas précédents, ils sont tous très faibles sauf un, m 1→13 = 0.36 kg/m 2 /s, qui était pourtant 10 fois plus faible lorsqu'il y avait m 0→1 = 0.33 kg/m 2 /s. Les deux boîtes 1 et 13 ont une énergie presque identique parce qu'elles sont très bien mélangées. En revanche leur rapport de mélange d'eau q S (T ) est très différent, il est beaucoup plus grand pour la boîte près du sol, ce qui fait tendre les précipitations dans la boîte 13 vers l'infini lorsque m 0→1 diverge. En fait ce flux de masse permet d'augmenter m 1→13 tout en conservant l'eau. Cette apparition de précipitations infinies montre qu'il manque certainement des contraintes dans le problème. Peut-être qu'il faut une condition pour que la couche limite de surface ne puisse pas trop se mélanger, ou bien il faudrait rajouter une contrainte générale sur la mécanique.

Remerciements