Code Dynamique 3-D

Nul besoin de lire cette partie technique pour travailler avec le GCM!

Discrétisation des équations de la dynamique

Extrait de la note de Robert Sadourny, Phu Le Van et Frédéric Hourdin, Laboratoire de Météorologie Dynamique.

Le modèle climatique du LMD est bâti, comme tous les modèles de circulation générale atmosphérique, sur la résolution numérique des équations primitives de la météorologie décrites dans de nombreux ouvrages (6). L'analyse présentée ici a été menée sur la nouvelle version de la dynamique du LMD écrite par Phu Le Van (11) sur une formulation de Robert Sadourny. Cette formulation diffère de l'ancienne essentiellement par deux points: dans la nouvelle formulation, la répartition des points en longitude et en latitude peut être changée arbitrairement. L'autre modification porte sur la répartition des points aux pôles^3.1.

La coordonnée verticale du modèle est la pression normalisée par sa valeur à la surface: $\sigma = p/p_s$. On utilise en fait $\sigma$ aux niveaux inter-couches et $s=\sigma^\kappa$ au milieu des couches. On note $X$ et $Y$ les coordonnées horizontales:

Figure 3.1: Grille obtenue avec 96 points en longitude et 73 en latitude et un zoom d'un facteur 3 centré sur la méditérannée (grille utilisée au laboratoire par Ali Harzallah)

.

$X$ (resp. $Y$) est une fonction biunivoque de la longitude $\lambda$ (resp. de la latitude $\phi$). Ces deux fonctions peuvent être choisies de façon arbitraire dans le modèle LMDZ ce qui permet d'effectuer un zoom sur une région du globe particulière. Une grille de ce type est montrée sur la Figure 3.1. Les variables scalaires (température potentielle $\theta = c_p T/{p_s}^\kappa$, géopotentiel $\Phi$ et pression de surface $p_s$) sont évaluées aux points correspondant à des couples de valeurs entières $(X,Y)=(i,j)$. Les variables dynamiques sont décalées par rapport aux variables scalaires en utilisant une grille $C$ dans la définition de Arakawa (): le vent zonal est calculé aux points $(X,Y)=(i+1/2,j)$ et le vent méridien aux points $(X,Y)=(i,j+1/2)$. La disposition des variables sur la grille est illustrée sur la Figure 3.2.

Figure 3.2: Disposition des variables dans la grille du LMD

On utilise en fait les composantes covariantes ($\tilde{u}$ et $\tilde{v}$) et contravariantes ( $\tilde{\tilde{u}}$ et $\tilde{\tilde{v}}$) du vent définies par

$$\begin{array}{llllllllll} \tilde{u}= c_uu & \mbox{et} & \til... ...c_v& \mbox{avec} & c_v= a \left( d\phi / dY \right) \end{array}$$ (3.1)

où $u$ et $v$ sont les composantes physiques du vecteur vent horizontal. On introduit également:

la pression extensive:

$\tilde{p}_s$ (pression au sol multipliée par l'aire de la maille).

les trois composantes du flux de masse:

$$U={\overline{\tilde{p}_s}}^{ X } \tilde{\tilde{u}}, V= {\ov... ...}_s\dot{\sigma} \mbox{avec} \dot{\sigma}=\frac{d\sigma}{dt}$$ (3.2)

le facteur de Coriolis multiplié par l'aire de la maille:

$f=2\Omega \sin{\phi} c_uc_v$
où $\Omega$ est la vitesse de rotation de la planète.

la vorticité potentielle absolue:

$$Z=\frac{{\cal F}\left( \delta_X \tilde{v}- \delta_Y \tilde{u} \right) + f}{{\overline{\tilde{p}_s}}^{ X,Y }}$$ (3.3)

l'énergie cinétique

$$K=\frac{1}{2} \left( {\overline{\tilde{u}\tilde{\tilde{u}}}}^{ X } + {\overline{\tilde{v}\tilde{\tilde{v}}}}^{ Y } \right)$$ (3.4)

La notation $\delta X$ signifie simplement qu'on effectue la différence entre deux points consécutifs suivant la direction $X$. La notation ${\overline{a}}^{ X }$ signifie qu'on prend la moyenne arithmétique de la quantité $a$ suivant la direction $X$. ${\cal F}$ est un filtre longitudinale appliqué dans les régions polaires. Les équations discrétisées sont écrites sous la forme suivante:

équations du mouvement:

$$\frac{\partial \tilde{u}}{\partial t} - {\overline{Z}}^{ Y }... ...\overline{\tilde{p}_s}}^{ X } \delta_Z \sigma } =S_{\tilde{u}}$$ (3.5)

oú $\tilde{u}_{a}$ est la composante zonale covariante du vecteur vent absolu: $\tilde{u}_{a}=\tilde{u}+c_ua \Omega \cos{\phi}$ et

$$\frac{\partial \tilde{v}}{\partial t} + {\overline{Z}}^{ X }... ...\overline{\tilde{p}_s}}^{ X } \delta_Z \sigma } =S_{\tilde{v}}$$ (3.6)

équation thermodynamique:

$$\frac{\partial \left( \tilde{p}_s\theta \right) }{\partial t... ...verline{\theta }}^{ Z } W \right) }{\delta_Z \sigma}=S_\theta$$ (3.7)

équation hydrostatique:

$$\delta_Z \Phi=-p_s^\kappa {\overline{\theta }}^{ z } \delta_Z s$$ (3.8)

équations de continuité:

$$\frac{\partial p_s}{\partial t} = {\cal F}\left[ \sum_z{\delta_Z \sigma \left( \delta_X U+ \delta_Y V \right) } \right]$$ (3.9)

$$\delta_Z W = -\delta_Z \sigma \left[ {\cal F}\left( \delta_X U+ \delta_Y V \right) + \frac{\partial p_s}{\partial t} \right]$$ (3.10)

On a noté $S$ les termes sources dans les différentes équations. Dans ces termes sources, on distingue 1) d'une part les paramétrisations physiques mentionnées plus haut et qui font intervenir pour une maille donnée du modèle, tous les points situés sur une même verticale mais ceux-là seulement; 2) les opérateurs de dissipation horizontale, censés rendre compte des échanges entre échelles explicitement représentées dans le modèle et échelles sous-mailles. Ces opérateurs ont la structure de Laplaciens agissant sur des plans horizontaux c'est à dire qu'il font intervenir un voisin de chaque côté dans les deux directions horizontales. Cet opérateur est généralement itéré pour le rendre plus sélectif en échelle (plus on itère un laplacien et plus son effet sur les petites échelles devient important relativement).

Filtres aux hautes latitudes

Extrait adapté de Forget et al. [1999]

At high latitude a filter is applied near the singularity in the grid at the pole in order to satisfy the Courant-Friedrichs-Lewy numerical stability criterion without going to an excessively small timestep. In the original version of the dynamical code a classical Fourier filter was used, but we found that because the Martian polar atmosphere appears to be much more dynamically unstable than the Earth's polar atmosphere, a more efficient formulation (based on the grouping of adjacent gridpoints together) was necessary to avoid numerical instability.

En pratique, la technique suivante est utilisée dans la subroutine nommée groupeun.F :

• Les points sont regroupés par paquets de $2^{\mbox{ngroup}}$ aux poles (e.g. ngroup=3 $\rightarrow$ paquets de 8), puis $2^{\mbox{ngroup-1}}$, $2^{\mbox{ngroup-2}}$, etc... aux latitude plus basse en s'eloignant du pôle

• Plus ngroup est élevé, plus le lissage est efficace, et plus le modèle est stable.

• MAIS, il faut iim divisible par $2^{\mbox{ngroup}}$ !!!

Dissipation

Extrait adapté de Forget et al. [1999]

In the LMD grid point model, nonlinear interactions between explicitly resolved scales and subgrid-scale processes are parameterized by applying a scale-selective horizontal dissipation operator based on an $n$ time iterated Laplacian $\Delta^{n}$. For the grid point model, for instance, this can be written ${\partial q}/{\partial t} = ([-1]^{n}/ {\tau_{\mbox{\scriptsize diss}}}) (\delta x)^{2n} \Delta^{n} q$ where $\delta x$ is the smallest horizontal distance represented in the model and $\tau_{\mbox{\scriptsize diss}}$ is the dissipation timescale for a st ructure of scale $\delta x$. These operators are necessary to ensure the grid point model numerical stability. In practice, the operator is separately applied to (1) potential temperature, (2) the divergence of the flow, and (3) its vorticity. We respectively use $n=2$, $n=1$, and $n=2$ in the grid point model.

Note : En pratique, les valeurs des $n$ et de $\tau_{\mbox{\scriptsize diss}}$ sont réglables et prescrites au début de chaque run dans le fichier de définition du run ``run.def'' (cf. 6.2.2)

Sponge layer

Uniquement pour le modèle martien. Adapté de Forget et al. [1999]

In the upper levels a sponge layer is also used in both models in an attempt to reduce spurious reflections of vertically propagating waves from the model top. Unlike the traditional Rayleigh friction formulation, this operates as a linear drag solely on the eddy components of the vorticity and divergence fields and is not scale-selective. The timescales on which it operates are typically half a day, 1 day, and 2 days at the three uppermost levels, respectively.

Note : les valeurs des ``timescale'' de la sponge layer et son extension en altitude sont réglables et prescrites au début de chaque run dans le fichier de définition du run ``run.def'' (cf. 6.2.2)