Mécanique statistique de modèles d'écoulements océaniques

Nicolas Sauvage

Description

Nicolas SAUVAGE
UNIVERSITÉ CLAUDE BERNARD LYON I 2e année de Master
ECOLE NORMALE SUPÉRIEURE DE LYON Sciences de la matière
Rapport de Stage
sujet :
Mécanique statistique de modèles d’écoulements océaniques
Laboratoire de Physique de l’ENS Lyon Juillet 2005
1Table des matières
Résumé 3
1 Modélisation des écoulements océaniques :le modèle Quasi Géostrophique 5
1.1 Description de l’équilibre géostrophique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.2 Modèle QG équivalent barotrope et vorticité potentielle . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.3 Contraintes et quantités conservées . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
2 Mécanique statistique des écoulements 7
2.1 Modèle à deux niveaux de vorticité potentielle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
2.2 Etats d’équilibre statistique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
2.2.1 Etats de Gibbs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
2.2.2 Application à plusieurs problèmes physiques . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
3 Résolution numérique des équations 9
3.1 Principe de la résolution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
3.2 Résolution du problème variationnel linéarisé . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
3.3 ...

Sujets

• Circulation atmosphérique
• PVA
• Topographie
• Ailier fort
• Champéry
• Igny

Informations

• Publié par: Chawyog
• Nombre de lectures: 226
• Langue: Français
• Poids de l'ouvrage: 3 Mo

Extrait

U INEVSRTIÉ C ALDUE B NAERRD L YON I E CELO N LEORMA S UPÉRIEURE DE L YON

Rapport de Stage

Nicolas SAUVAGE

2e année de Master Sciences de la matière

sujet : Mécanique statistique de modèles d’écoulements océaniques

Laboratoire de Physique de l’ENS Lyon

1

Juillet 2005

5

5.1 Ordres de grandeur et adimensionnalisation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23

5.1.1 Caractérisation du Gulf Stream . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 5.1.2 Choix du système d’unités . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 5.2 Structures des états de Gibbs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23

5.2.1 Sans topographie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

23

5.2.2 Avec effetβ. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

25

28

3 Résolution numérique des équations

16

4.1 Energie associée à une forme de jet donnée pourR→0. . . . . . . . . . . . . . . . . . . . . . . . . . . 16 4.2 Structure des écoulements maximisant l’énergie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18

4.2.1 Solutions pourβ= 0 18 . . . . . . . . . . . . . . . . . . . .:Transition de phase du premier ordre 4.2.2 Inuence de la topographie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

Ecoulements d’énergies maximales

4

22

Structures d’équilibre du modèle QG barotrope

2.2.1 Etats de Gibbs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

9

2.2.2 Application à plusieurs problèmes physiques . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9

3.2 Résolution du problème variationnel linéarisé . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10

3.1 Principe de la résolution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10

3.5 Calcul des états d’énergie maximale . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13

3.3 L’équation de Helmoltz . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 3.4 Extrémalisation numérique d’une fonction :méthode de Newton . . . . . . . . . . . . . . . . . . . . . . . 13

Modèle QG équivalent-barotrope et vorticité potentielle . . . . . . . . . . . . . . . . . . . . . . . . . . .

1.1

6

1.2

Contraintes et quantités conservées . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7

Mécanique statistique des écoulements 7 2.1 Modèle à deux niveaux de vorticité potentielle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 2.2 Etats d’équilibre statistique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

1.3

Reférences

Résumé

5

2

6 Conclusion

3

2

29

Table des matières

5

Description de l’équilibre géostrophique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1 Modélisation des écoulements océaniques :le modèle Quasi-Géostrophique

Remerciements

Ce stage s’est effectué au Laboratoire de Physique de l’Ecole Normale Supérieure de Lyon, d’avril à juillet 2005.

Je tiens à remercier tout d’abord mon maître de stage Freddy Bouchet, Chargé de Recherche à l’ENSL, qui a guidé mes investigations durant ces quatre mois en manifestant une attention et un soutient permanents à mon travail.

Je salue également Thierry Dumont, chercheur à l’UCBL, pour ses suggestions judicieuses en analyse numérique.

J’exprime enn toute ma reconnaissance à l’ensemble du personnel du laboratoire, enseignants, doctorants et techni- ciens, pour leur accueil sympathique.

Résumé

Ce travail est une application de la mécanique statistique à l’étude des ots géophysiques à grande échelle. Nous souhaitons expliquer l’existence à moyennes latitudes de jets intenses dirigés vers l’est dans les océans.

Nous calculons numériquement les états d’équilibre statistiques d’un modèle simplié d’océan (Quasi-Géostrophique barotrope). Ces états d’équilibre statistiques correspondent aux écoulements stationnaires d’un uide inertiel non forcé. Dans la limite des faibles rayons de Rossby nous présentons également des arguments analytiques pour comprendre la structure de ces écoulements stationnaires. Nous expliquons la forme de jets intérieurs par l’organisation de la vorticité potentielle à l’intérieur d’un bassin océanique rectangulaire. Nous mettons ainsi en évidence des transitions entre plusieurs structures d’équilibre, en fonction notamment de la circulation et de la topographie.

Nous montrons que ce modèle ne permet pas d’états d’équilibre avec un jet océanique intérieur dirigé vers l’est, contrairement aux modèles Quasi-Géostrophique barotropes avec faibles forcage et dissipation visqueuse. Ceux-ci sont donc fortement hors-équilibre. La forte sensibilité des états d’équilibre de notre modèle vis à vis de la circulation à la frontière du bassin suggère d’utiliser ultérieurement des modèles d’océans globaux intégrant l’équateur (type Shalow- Water), pour tenter de justier l’existence de jets vers l’est par des états d’équilibre statistique.

3

Introduction

Les écoulements géophysiques externes (océans, atmosphères) ont la caractéristique commune de former des struc- tures cohérentes à grande échelle. Cette organisation en cyclones-anticyclones, courants marins, jets, peut être comprise grâce à une théorie récente basée sur la mécanique statistique des structures à petites échelles, qui permet de prédire l’état le plus probable pour l’écoulement à grande échelle. Cette approche a permis de prédire avec succès l’organisation de l’atmosphère de Jupiter (par exemple la Grande Tache Rouge).

Il nous semble donc judicieux d’appliquer cette approche à l’étude des courants océaniques. L’objectif du stage est de déterminer les états d’équilibres statistiques d’un modèle simplié d’océan. Nous pourrons alors savoir si ces équilibres reproduisent bien les courants océaniques réels. En particulier, nous souhaiterions expliquer l’existence à moyennes lati- tudes de jets intenses dirigés vers l’est (comme le Gulf Stream dans l’océan atlantique) grâce à des arguments statistiques.

A plus long terme, cette théorie pourrait permettre d’étudier la dynamique hors équilibre des océans. Ceux -ci sont en effet soumis à des forçages incessants : vents, ensoleillement, échanges thermiques provoquent des uctuations qu’il serait intéressant de prédire.

A l’échelle des océans, les courants sont quasi-bidimensionnels et dominés par la force de Coriolis. Pour modéliser des écoulements de ce type, les océanographes utilisent des modèles de complexité très différentes, déduit des équations de Navier Stokes. On peut utiliser par exemple des équations primitives (météo), ou des modèles en couches (par exemple les modèles isopycnaux ou shallow water multicouches). Dans ce travail, nous choisirons d’utiliser le modèle le plus simple contenant la physique de la dynamique à grande échelle des océans aux latitudes moyennes : le modèle Quasi- Géostrophique équivalent barotrope. Certains travaux récents ([19],[18] et [20]) ont utilisé ce type de modèle d’océan avec un faible forçage et une faible dissipation visqueuse . Les différentes instabilités et leur bifurcation ont été étudiées numériquement.

Le forçage et la dissipation étant relativement faibles pour les océans réels, nous allons ici rechercher des solutions stationnaires en supposant le uide inertiel et non forcé. Certaines de ces solutions stationnaires ont été étudiées dans les travaux classiques de Fofonoff [2] . Ces études étaient limités à des solutions ayant une relation linéaire entre vorticité potentielle et fonction courant. Nous souhaitons utiliser la théorie statistique pour trouver les états stationnaires du cas général non-linéaire : ces états d’équilibre statistiques n’ont encore jamais été étudiés à notre connaissance.

Dans une première partie nous présentons le modèle retenu pour étudier les courants marins. Nous expliquerons le rôle fondamental de la force de Coriolis et les conséquences de la stratication sur les écoulements. Nous verrons que le modèle Quasi-Géostrophique équivalent barotrope est le plus simple qui prenne en compte ces phénomènes. Nous utiliserons ensuite la théorie statistique de Robert, Miller et Sommeria pour prédire l’écoulement géophysique le plus probable (état d’équilibre thermodynamique) (partie 2). La description probabiliste des écoulements à grande échelles permettra de réduire les degrés de liberté de ce système, tout en conservant la notion de mélange turbulent à petite échelle. Dans la partie 3 nous présenterons les méthodes utilisées pour calculer numériquement les états d’équilibre statistique. Nous chercherons à déterminer les écoulements stationnaires contenus dans un domaine rectangulaire, en fonction des quantités conservées que sont l’énergie et la circulation. Les techniques utilisées pour résoudre les différentes équations du problème sont en majorité inspirées de travaux précédents, mais nous avons dû les compléter par des méthodes origi- nales. Nous consacrons la partie 4 aux résultats obtenus concernant les états de plus hautes énergie. Nous confrontons des approches analytiques et numériques qui révèlent des phénomènes de transition de phase du 1er ordre entre différentes structures d’écoulement. Nous verrons dans la partie 5 que ces structures conditionnent fortement celles des états d’éner- gies inférieures. Nous discuterons de l’ensemble des écoulements envisageables, pour conclure sur la pertinence de notre approche.

4

1 Modélisation des écoulements océaniques : le modèle Quasi-Géostrophique

Les écoulements géophysiques sont dominés par la force de Coriolis. Elle est le moteur principal des structures rota- tionnelles à grande échelle bien connues de l’atmosphère (cyclones, anticyclones), ainsi que des océans. Nous présentons ici le modèle Quasi-Géostrophique, le plus simple qui prenne en compte ses effets, et nous en expliciterons les équations résultantes permettant d’étudier la première couche d’un océan stratié. Pour approfondir les concepts qui suivent, nous conseillons les ouvrages de Pedlosky [1] et de Gill [3].

1.1 Description de l’équilibre géostrophique Le Gulf Stream est un jet de vitesse typique1ms −1 allant du golfe du mexique aux côtes nord-européennes. Sa largeur est de l’ordre de50km. L’écoulement est donc associé à un nombre de Reynolds très élevé (typiquementRe≈10 11 ). Les effets de dissipation visqueuse seront donc dans une très bonne approximation négligés. En notantvla vitesse du uide,ρsa densité,Pla pression,gla gravité etΩla vitesse de rotation planétaire, les champs satisfont : ∂ ∂ v t + (v∙ r)v= − r ρ P −2Ω×v+g(1) r ∙v= 0(2)

Le nombre de Rossbymesure le rapport entre l’accélération convective(v∙ r)vet la force de Coriolis2Ω×v. Pour des courants océaniques de vitesse typiqueV= 1ms −1 organisés sur une échelle de longueur typiqueR= 50km, = 2 V ΩR ≈ 0 1 1 , et diminue encore pour des échelles plus importantes. Ainsi, dans une bonne approximation, l’équation d’Euler est dominée par l’équilibre géostrophique entre le gradient de la pression et la force de Coriolis :

rP + 2Ω×v=0 ρ

(3)

Dans ces conditions, la force de Coriolis importante rend le mouvement quasi-bidimensionnel (effet Proudman- Taylor). Cet effet est accentué par la stratication du uide : des variations brusques de température ou/et de salinité séparent des couches de uide dont l’épaisseur moyenne est faible devant l’extension horizontale. On décrit donc l’écou- lement moyen d’une couche dans un plan horizontal (x,y) d’altitude constante grâce à la fonction de courantψtelle que : v=r hor ×(ψe z ) =−e z × r hor ψ(4) oùe z désigne un vecteur unitaire suivant la verticale ascendante. Dans une région de latitudeλ, on déni la vorticité planétairef= 2Ωsinλ. On déduit alors la fonction de courant à partir de l’équilibre géostrophique : r hor P=−2ρΩsinλe z ×(−e z × r hor ψ) =−fr hor ψ⇒ψ=− P f +cste(5)

L’équilibre géostrophique impose donc aux courants de suivre des lignes de pression constante, ce qui rappelle un comportement analogue pour des mouvements atmosphériques bien connus : cyclones et anticyclones sont repérés par les lignes isobares. Grâce à cet équilibre, les uides soumis à de fortes rotations pourront atteindre un équilibre dynamique avec des isobares non horizontales. La force de gravitation lie le gradient de pression à la déformation des isobares (équi- libre hydrostatique). La dimension caractéristique de cette déformation peut être évaluée par analyse dimensionnelle : R= √ gH /2Ω.Rest le rayon de déformation de Rossby, dépendant de la vitesse de rotation de la planèteΩ, de l’épais- seurHde la couche de uide et de la gravitationg. Cette longueur est essentielle dans l’étude des ots géophysiques : on l’évalue de20à100kmeffectivement proche de la largeur du Gulf Stream.pour les océans, ce qui est Pour des bassins océaniques offrant des variations de latitude importante, il faut prendre en compte la variation de la vorticité planétairef= 2Ωsinλ. La première correction qui prenne en compte ce phénomène consiste à linéariser la variation autour une latitude moyenne. En orientant la coordonnéeyvers le nord : f(y)≈f(y 0 )− 2Ω R c T y osλ =f(y 0 )−βyavecR T  (6)le rayon de la planète On effectue ainsi une approximationβ-plan qui sera valide si l’extention méridionnale du bassin reste faible devantR T .

5

1.2 Modèle QG équivalent-barotrope et vorticité potentielle

La stratication de l’océan est la conséquence de plusieurs effets : l’équilibre géostrophique combiné à la faible extension verticale de l’océan (≈4km) conduit à la création de couches de température et de salinité différentes. La couche supérieure, offrant en plus des possibilités d’échange avec l’atmosphère, possède la dynamique la plus importante : on l’appèle Couche de Mélange Océanique (CMO). Elle est séparée des couches inférieures par la thermocline principale, qui consiste en une brusque variation de température au voisinage d’une profondeurH= 500m:

F IG . 1 –Modèle océanique en couches

Le modèle Quasi-Géostrophique barotrope ne prend en compte que l’écoulement moyen de la CMO de densité constante. Des rafnements sont possibles en prenant en compte l’écoulement de couches inférieures aux densités diffé- rentes, ainsi que le forçage des vents atmosphériques. Nous renvoyons le lecteur aux ouvrages de Pedlosky [1] et Dijkstra [5] pour plus de précisions sur ces modèles.

L’altitude de l’interface air-eau est notéeη(x,y), oùxetyposition selon les axes est-ouest et nord-sud.représentent la Les variations deηseront supposées faibles devant l’épaisseurHde la CMO. Sous certaines hypothèses (bidimension- nalité, incompressibilité, équilibre hydrostatique,→0), on peut montrer (voir Pedlosky [1]) que la dynamique de la couche de uide est régie par les équations Quasi-Géostrophiques :

∂ ∂ q t +v∙ rq= 0 q=−Δψ+ R ψ 2 −h(y) v=−e z × rψ

(7)

(8)

(9)

oùqest la vorticité potentielle (PV), advectée par le champ de vitesse non divergentv, etRest la rayon de déformation de Rossby. La relation (9) traduit l’équilibre géostrophique. Le termeh(y)est appelé topographie équivalente, car dans des modèles plus rafnés avec plusieurs couches, il combine la topographieh B de la couche profonde et l’effetβ. Pour un modèle QG barotrope,h B = 0, mais nous pouvons considérer l’effetβcomme une topographie qui décroit en allant vers le nord. Nous travaillerons uniquement avech(y) =−βy, en choisissant l’origine deyadaptée pour que la topographie soit de moyenne nulle sur le domaineDqu’occupe la CMO : < h(y)>= Z D h(y) d 2 r=−β Z D ydxdy= 0(10) Pour un écoulement stationnaire, l’équation (7) se résume àv∙ rq= 0, ou encore en utilisant la fonction de courant rq× rψ= 0. Cela signie que les lignes de niveau deψetqet donc que les deux champs sontsont identiques, directement liés par une relation du typeq=f(ψ). Le problème revient à exprimer la vorticité potentielle avec la fonction de courant pour résoudre : −Δψ+ ψ R 2 =f(ψ) +hsur un domaineD(11)

6

1.3 Contraintes et quantités conservées

Nous avons vu que les variations deψsont proportionnelles à celles de la pression géostrophique, et donc proportion- nelles au niveauη(x,y)en équilibre hydrostatique. Par conservation de la masse du uide incompressible, lede l’ocean niveau moyen de la surface doit rester constant, ce que nous traduirons sans perte de généralité en imposant : < ψ >= Z ψ(x,y) dxdy= 0(12) D D’autre part nous supposerons que sur les bords la vitesse sera toujours tangentielle à la frontière, ce qui se traduit (voir (9)) parψ=cste=ψ f r Les conditions aux limites de (8) seront donc. ψ=cste=ψ f r sur∂Det< ψ >= 0(13) En utilisant (8) et (13), nous pouvons estimer la circulation sur la frontière∂Ddu domaine en fonction deq: Γ = Z ∂D v∙dl= Z D −Δψdxdy=< q+h >− <ψ R 2 > =< q >(14) L’énergie de l’écoulement comporte des parties cinétique et potentielle : U=U cin +U po ρ 2 H Z (rψ) 2 dxdy+ ρ 2 H Z D ψ R 2 2 dxdy t = D U= ρ 2 H Z ∂D h ψrψ∙n i dl + ρ 2 H Z D h (−Δψ+ R ψ 2 )ψ i dxdy= ρ 2 H h −ψ f r Γ+<(q+h)ψ > i (15) Dans la suite, nous choisirons l’unité de masse adaptée pour que l’énergie adimensionnée s’exprime par U= 1 2 h <(q+h)ψ >−ψ f r Γ i (16) D’autres quantités conservées sont susceptibles d’exister, en fonction des symétries du domaine. Nous ne les utiliserons pas ici.

La vorticité potentielle est un scalaire, advectée par une vitesse à divergence nulle. En conséquence toutes les fonc- tionnelles de la forme : C f [q] = Z f(q) d 2 r(17) D sont conservées. La conservation de ces fonctionnelles est également une conséquence de la dégénérescence du modèle QG considéré comme un système hamiltonien de dimension innie. Ces fonctionnelles apparaissent alors liées au crochet de Poisson (Casimirs).

2 Mécanique statistique des écoulements

Les écoulements atmosphériques et océaniques sont souvents organisés en jets, et leur stabilité en milieu fortement turbulent est étonnante : la Grande tache rouge est en effet observée depuis plusieurs siècles. Nous expliquons ce phéno- mène grâce à une approche de mécanique statistique, où le jet apparaît comme l’état le plus probable de l’écoulement, étant donné des contraintes sur l’énergie et la circulation conservées par la dynamique. Une telle approche statistique a pour la première fois été proposée pour les équations d’Euler en 2D par Kuz-min (1982), Robert, Miller et Sommeria (1991). Les résultats conrmaient l’organisation de la turbulence en écoulement moyen stable superposé par des uctuations de vorticité potentielle aux échelles microscopiques.

2.1 Modèle à deux niveaux de vorticité potentielle

Comme tous les écoulements à hauts nombre de Reynolds, les équations QG conduisent à des structures de vorti- cité très complexes, avec l’apparition de laments fortement instationnnaires à petites échelles. Nous considérons ces structures à petites échelles comme microscopiques. Décrire correctement l’écoulement jusqu’à ces échelles nécessiterait un nombre gigantesque de degrés de liberté, ce qui est irréalisable et surtout peu approprié lorque l’on veut se focaliser

7

sur l’organisation à grande échelle. C’est pourquoi il est préférable d’adopter une vision probabiliste pour le champ de vorticité potentielle.

On considère que microscopiquement les niveaux accessibles par la PV sont mélangés par la turbulence. Nous cher- cherons donc à caractériser l’écoulement par la moyenne localeq¯de la vorticité potentielle (coarse-graining), qui constitue notre description macroscopique. L’état le plus probable dépend alors uniquement des contraintes sur l’énergieEet la circulationγrecherchées, ainsi que des niveaux de PV microscopiques accessibles. Par simplicité, nous choisissons de travailler avec deux niveaux microscopiques disponibles, hypothèse que nous justierons par la suite. Nous discuterons à partir de 2 niveaux adimensionnnésa + eta − en choisissant l’unité de temps adéquate : a + +a − =B 2 et a + 2 −a − = 1(18) oùBcaractérise la dissymétrie des deux niveaux. Dans une région de taille macroscopique repérée parr, la probabilité de trouver le niveau le plus haut estp(r), si bien que la valeur globale de PV pour cette région est donnée par :

1 ¯ q(r) =a + p(r) +a − (1−p(r)) = 2(p− +2 )B(19) Les uctuations à petite échelle de la PV ont un effet négligeable sur les valeurs de la fonction de courant (le calcul de ψportée), si bien que les champs macroscopiques satisfont toujours lesfait intervenir une moyenne spatiale à longue équation QG moyennées : ¯ q¯ =−Δψ ¯ + ψ h ; ¯ v=−e z × rψ ¯ (20) − R 2 Nous utiliserons alors la circulationΓ =< q¯>et l’énergieU= 1 2 h <(q¯ +h)ψ ¯ >−ψ f r Γ i en négligeant l’énergie des uctuations microscopiques. Cette approximation a été justiée par Robert et Someria (1991), pourvu que le "cutoff" soit effectué a des échelles très inférieures au rayon de RossbyR, qui correspond à la largeur typique des jets. Dorénavant, ¯ nous travaillerons uniquement avec les quantités macroscopiques et nous poseronsq=q¯etψ=ψ.

2.2 Etats d’équilibre statistique

Nous allons désormais chercher à caractériser les états d’équilibre statistique. Nous montrerons que les champs ma- croscopiquesq,vetψsont solutions des équations QG stationnaires.

2.2.1 Etats de Gibbs

D’un point de vue statistique, le nombre d’états macroscopiques correspondant à un champ de probabilitépest donné par l’entropie associée au système à deux niveaux : S=− Z D [p(r) lnp(r) + (1−p(r)) ln(1−p(r))] d 2 r(21) Nous nous plaçons dans l’ensemble microcanonique, c’est à dire que nous xons les quantités invariantes de la dyna- mique que sont la circulation et l’énergie. Pour un modèle à deux niveaux, la conservation de la circulation implique automatiquement la conservation des Casimirs, c’est pourquoi ils n’apparaitront pas dans nos contraintes. Pour obtenir des équilibres associés à une circulationγet une énergieEdonnées, nous devons donc calculer les points stationnaires de l’énergie libreFconstruite avec deux multiplicateurs de Lagrangeαet−C/R 2 : F[p(r)] =−S[p(r)]− R C 2 U[p(r)] +αΓ[p(r)](22) Nous recherchons doncp(r) satisfaisantδF=−δS− RC 2 δU+αδΓ = 0, avec : δS=− Z D [lnp−ln(1−p)]δpd 2 r;δΓ = Z D δqd 2 r= Z D 2δpd 2 r; δU= Z (ψ−ψ f r )δqd 2 r= Z D 2(ψ−ψ f r )δpd 2 r D

8

L’expression pourδUest obtenue après une intégration par partie et en utilisant (19). On doit donc annuler l’intégrale R D [ln(1−p)−lnp+ 2α−2 RC 2 (ψ−ψ f r )]δpd 2 rquelle que soit la différentielleδp, ce qui revient à annuler l’intégrant. On en déduit donc pour l’équilibre statistique : p = 1−tanh  α− 2 RC 2 (ψ−ψ f )  ce qui est équivalent àq=B−tanh  α− R C 2 (ψ−ψ f r )  (23) r Cette dernière expression offre un lien du typeq=f(ψ)montrant que la vorticité potentielle macroscopique peut varier continuement entrea + eta − . Il s’agit d’une solution stationnaire paramétrée parαetC. Leurs valeurs sont telles que les contraintes sur l’énergie et la circulation sont respectées. Ainsi, la mécanique statistique prédit un écoulement stationnaire que nous allons comparer aux courants réels. Si nous avions plus de 2 niveaux microscopiques accessibles, la forme def en serait légèrement affectée, mais ses caractériques principales conservées : monotonie et limites asymptotiques vers les niveaux extrêmes de PV. Ces modications n’apportent donc pas de changement physiques sufsants pour justier des développements techniques supplémentaires.

2.2.2 Application à plusieurs problèmes physiques

L’équation des états d’équilibres a été utilisé dans plusieurs travaux antérieurs, qui selon les auteurs sont focalisés sur les problèmes analytique ou numérique de résolution des équations Quasi-Géostrophiques :

?Aspects analytiques

– Dans le casR=∞(équation d’Euler sans vorticité planétaire), Chavanis etal[13] ont proposé des solutions théoriques pour une fonctionf(ψ)linéarisée, ce qui correspond à des états d’équilibre de faibles énergies, comme nous le verrons dans la suite. – Dans la limite des faibles rayons de déformation (R→0), les calculs de Bouchet [6] [7] ont expliqué l’existence de jets stables dans l’atmosphère jovienne fortement turbulente. Les prédictions théoriques sont en accord précis avec les observations des structures réelles de la Grande Tache Rouge.

?Aspects numériques

– Des méthodes de relaxation ont été utilisés par Sommeria etal: l’ajout d’une dissipation visqueuse articielle[11] fait converger les solutions vers les états stationnaires recherchés. L’étude montre la destabilisation d’un double jet en tourbillon. – Thess et Sommeria proposent une méthode alternative dite de continuation [12] : en partant des énergies faibles la procédure identie les éventuelles bifurcations donnant naissance à des états métastables. – Une autre méthode a été mise au point par Whitaker et Turkington [17] dans le casR=∞, et utilisée récemment par Tatekawa etal[8] pour étudier la thermodynamique de systèmes gravitationnels. A l’inverse de la prédédente, cette procédure nécessite le calcul des états d’énergies maximales.

Nous allons nous inspirer de ces travaux pour déterminer les équilibres du modèle QG avec des paramètres qui corres- pondent à un bassin océanique, ce qui n’a jamais été fait auparavant. Les écoulements stationnaires seront calculés numé- riquement selon une méthode proche de Whitaker et Turkington, pour un rayon de déformation de Rossby quelconque. Dans la limiteR→0réaliste pour les courants marins, nous effectuerons parallèlement des prévisions analytiques sur les états d’équilibres d’énergies maximales pour expliquer les résultats numériques. Nous discuterons de la forme des jets obtenus en fonction des paramètres du problème (circulationγ, dissymétrieBdes niveaux, effetβ) pour déterminer si les courant longitudinaux à moyennes latitudes correspondent à un état stationnaire de ce modèle.

3

Résolution numérique des équations

Pour calculer les états d’équilibre statistique, nous devons maximiser l’entropie S[q] associée au système à deux niveaux, avec des contraintes sur les invariants de la dynamique (énergieU[q]et circulationΓ[q]). Nous avons mis au point une méthode numérique pour trouver les solutions de ce problème variationnel. Elle est inspirée de celle utilisée par Whitaker et Turkington [17] pour résoudre un problème analogue de maximisation d’entropie concernant des écoulements bidimensionnels. Nous nous sommes également inspirés des idées de Tatekawa, Bouchet, Dauxois et Ruffo [8] dans le cadre d’un modèle de particules auto-gravitantes.

9

Dans un premier temps nous exposons le principe de la méthode itérative utilisée, qui calcule des champs de vorticité potentielleq k à chaque pas d’itération un problème varia-qui convergent vers le maximum d’entropie. Nous résolvons tionnel aux contraintes linéarisées, dont nous préciserons un algorithme de résolution. Nous nous attarderons plus par- ticulièrement sur les étapes numériquement délicates : l’équation de Helmoltz avec conditions aux limites non-standard ainsi que la maximisation d’une fonction concave. Nous présentons enn une méthode de calcul des états d’équilibre d’énergie maximale, dont les solutions permettent une bonne initialisation du calcul des maxima d’entropie.

3.1 Principe de la résolution

La fonctionnelle à maximiserS[q]strictement concave, et nous xons de plus une contrainte linéaireest Γ[q] =γ ainsi qu’une non-linéaireU[q] =E. C’est celle-ci qui rend le problème variationnel difcile à résoudre, aussi allons-nous considérer une contrainte sur l’énergie linéarisée autour d’une distributionq k provenant d’un pas précédent du calcul itératif. La distribution suivante est la solution du problème variationnel suivant : δU max ( S[q]|Γ[q] =γ , U[q k ] + Z δq  q k (q−q k ) d 2 r=E ) (24)

où la dérivée fonctionnelle vaut :

δ δq U  q k = (ψ k+1 −ψ k+1f r )avec−Δψ k+1 + R ψ k 2 1+ =q k +h(y)et< ψ k+1 >= 0(25) La maximisation d’une fonctionnelle strictement concave avec des contraintes linéaires possède ainsi une unique solution q k+1 . On peut également montrer que l’algorithme fait croître l’entropie (voir [8]) :S[q k+1 ]> S[q k ]. Comme l’entropie est majorée parAln 2,Aétant l’aire du bassin, la convergence du processus est théoriquement assurée. Signalons d’autre part que l’énergieU[q]est convexe, ce qui implique : U[q k+1 ] =U[q k + (q k+1 −q k )]≥U[q k ] + Z δ q U δ  q k (q k+1 −q k ) d 2 r⇒U[q k+1 ]≥E(26) L’énergie est donc toujours approchée par valeurs supérieures. Si nous connaissons les états d’énergie maximale, nous sommes alors certains de pouvoir initialiser avec une énergie supérieure à celle recherchée puisqueU[q 1 ] =E max ≥E. Ceci motive donc notre intérêt pour la recherche annexe des maxima d’énergie exposée dans le dernier paragraphe. Nous avons vérié numériquement (gure 2) queS[q k+1 ]> S[q k ]et queU[q k+1 ]≥Eà chaque pas du calcul. Nous allons maintenant détailler comment nous obtenons les champsq k .

F IG . 2 –EntropieS[q k ]et énergie relative U[q k E ]−E pour chaque pas d’itérationk≥1. Les niveaux de PV sont symé- triques (B= 0) et répartis sur une surface rectangulaire d’aireA= 1, pourγ= 0.5,E=E max /2etβ= 0.

3.2 Résolution du problème variationnel linéarisé

Comme nous l’avons expliqué précédemment, à chaque itération notre algorithme résoud un problème variationnel linéarisé. Nous détaillons ici les étapes de cette résolution. Nous supposons connu un champ de vorticité potentielleq k , et

10

nous en déduisons la fonction de courantψ k+1 en résolvant numériquement l’équation de Helmoltz (voir le paragraphe suivant): −Δψ k+1 + ψ R k1+ 2 =q k +h(y)avec< ψ k+1 >= 0etψ k+1 =cste=ψ k+1f r sur∂D(27) Dénissons le Lagrangien : L k [q](α,C) =−S[q]− R C 2
U[q k ] + Z δ q U δ (q−q k ) d 2 r−E ! +α(Γ[q]−γ)(28)  q k Par un calcul analogue à celui effectué dans la partie 2 sur la mécanique statistique, nous pouvons montrer que le champ de vorticité potentielleqlagrangien est de la forme :minimisant ce q=B−tanh  α− R C 2 (ψ k+1 −ψ k+1f r )  (29) A ce stade,αetCsont encore indéterminés. Intéressons-nous alors à : L k ? (α,C) = i q nf{L k [q](α,C)}=L k  q= tanh  α− C R 2 (ψ k+1 −ψ k+1f r )  (α,C);(30) On peut alors montrer [9] queL ? k est concave et queα k+1 etC k+1 correspondent à l’unique maxima deL ? k . En pratique nous utilisons l’expression suivante : L ? k (α,C) =− Z D ln  cosh  α− R C 2 (ψ k+1 −ψ k+1f r )  d 2 r+α Z D B+h(y) d 2 r  + R C Z D ψ k+1f r )}d 2 r−U[q k ]−E  (31) 2 {(B+h(y)) (ψ k+1 − An de maximiser cette fonction, nous utilisons un algorithme basé sur la méthode itérative de Newton qui consiste à annuler les dérivées premières. Une fois obtenus ces paramètres de Lagrange, nous calculons un nouveau champ de vorticité potentielle : q k+1 =B−tanh  α k+1 − C k R +1 2 (ψ k+1 −ψ k+1f r )  (32) et nous continuons le processus jusqu’à la convergence. En pratique nous arrêtons les calculs lorsque la quatité suivante est inférieure à la précisionrequise sur les champs : 2 s  α k+1 α k −α k  2 +  C k+1 C k −C k  +  U[q k+ U 1 [ ] q − k ] U[q k ]  2 + (Γ[q k+1 ]−Γ[q k ]) 2 < (33) Aux vues des expériences numériques, nous considérons que= 10 −4 est un critère d’arrêt sufsament contraignant (nous reprenons l’exemple numérique précédent sur la gure 3)

F IG . 3 –Quantité de convergence en fonction du nombre d’itérations, pourB= 0,γ= 0.5,E=E max /2etβ= 0

Les étapes numériquement difciles sont la résolution de l’équation de Helmoltz et la recherche du maximum de L k ? (α,C). Nous détaillons les techniques utilisés pour résoudre ces problèmes dans les deux prochains paragraphes.

11