Formulation incrémentale du 4D-Var

La formulation incrémentale a déjà été abordée rapidement à propos du 3D-Var dans la section 5.5.1, mais elle prend avec le 4D-Var tout son sens.

En effet, l'introduction de l'approche incrémentale en météorologie a été motivée par la réduction de coût qu'elle propose. Dans le cadre du 4D-Var classique décrit précédemment, à chaque itération de la minimisation de la fonction coût (Eqs. 5.27 et 5.28), l'intégration du modèle direct non-linéaire et de l'adjoint du modèle linéarisé est très coûteuse. Les non-linéarités des modèles numériques atmosphériques peuvent conduire à des fonctions coûts complexes. Les minimiseurs utilisés sur ces fonctions n'aboutissent pas forcément à une minimisation fiable et ce, de surcroît, à un coût élevé. Ces non-linéarités compliquent de plus, lourdement la tâche de l'écriture de l'adjoint (Xu, 1996). Dans l'approche incrémentale, la fonction coût est rendue quadratique, ce qui garantit l'identification d'un minimum unique par une méthode de descente pour un coût de calcul notablement inférieur à celui du problème non-linéaire. Le modèle linéaire-tangent est écrit avec une physique simplifiée, ce qui facilite grandement l'écriture de l'adjoint du modèle. De plus, une approche communément choisie en météorologie est d'utiliser un modèle linéaire-tangent à une résolution plus basse que celle du modèle non-linéaire. Le coût de la minimisation du 4D-Var en est significativement réduit. C'est formulation incrémentale de l'approche variationnelle du 4D-Var qui a permis de le rendre applicable de façon opérationnelle pour la prévision météorologique (Courtier , 1994 et Rabier , 2000).

Dans la formulation incrémentale de l'assimilation variationnelle, l'objectif est de minimiser la fonction coût, non plus par rapport à la variable d'état ${\mathbf x}$, mais par rapport à un incrément $\delta {\mathbf x}$ tel que ${\mathbf x}={\mathbf x}^b+\delta{\mathbf x}$. L'hypothèse principale de la formulation incrémentale est d'utiliser un modèle d'évolution et un opérateur d'observation linéarisé pour propager l'incrément mais de conserver le modèle non-linéaire pour la propagation de l'ébauche ${\mathbf x}^b$. La solution de la minimisation est l'incrément d'analyse $\delta{\mathbf x}^a$ à $t_0$ tel que le vecteur d'analyse ${\mathbf x}^a$ soit

$${\mathbf x}^a = {\mathbf x}^b+\delta{\mathbf x}^a.$$ (5.35)

En supposant que l'ébauche ${\mathbf x}^b$ est une bonne approximation a priori de l'état optimal du système au sens des moindres carrés, l'incrément $\delta {\mathbf x}$ devrait être petit. Les opérateurs non-linéaires d'observation $H$ et du modèle $M$ sont linéarisés au voisinage de l'ébauche de sorte que pour tout état du modèle ${\mathbf x}$, à chaque instant $t_i$ de la fenêtre temporelle d'assimilation,

$${\mathbf x}(t_{i+1}) = M_{0\to i}{\mathbf x}$$

$$= M_{0\to i}({\mathbf x}^b+\delta{\mathbf x})$$

$$\approx M_{0\to i}{\mathbf x}^b+{\mathbf M}_{0\to i}\delta{\mathbf x}.$$ (5.36)

Par conséquent

$$H_i{\mathbf x}(t_i) = H_i{\mathbf x}^b(t_i)+{\mathbf H}_i \delta{\mathbf x}+\parallel\delta{\mathbf x}\parallel^2+\cdots,$$ (5.37)

$$G_i{\mathbf x}(t_i) = G_i{\mathbf x}^b(t_i)+{\mathbf G}_i \delta{\mathbf x}+\parallel\delta{\mathbf x}\parallel^2+\cdots,$$ (5.38)

où ${\mathbf M}_{0\to i}$, ${\mathbf H}_i$ et ${\mathbf G}_i$ sont les opérateurs linéarisés à $t_0$ autour de l'ébauche $x^b$ de $M_{0 \to i}$, $H_i$ et $G_i$.

$${\mathbf M}_{0\to i} = \left.\dpdp{M_{0\to i}}{{\mathbf x}}\right\vert _{{\mathbf x}={\mathbf x}^b},$$

$${\mathbf H}_i = \left.\dpdp{H_i}{{\mathbf x}}\right\vert _{{\mathbf x}={\mathbf x}^b},$$

$${\mathbf G}_i = \left.\dpdp{G_i}{{\mathbf x}}\right\vert _{{\mathbf x}={\mathbf x}^b}.$$ (5.39)

En insérant les Eqs. 5.36, 5.37 et 5.38 dans la formulation de la fonction coût donnée par les Eqs. 5.27 et 5.28, elle peut être reformulée de manière incrémentale

$$J^b({\mathbf x}) = \frac{1}{2}\delta{\mathbf x}^T{\mathbf B}^{-1}\delta{\mathbf x},$$ (5.40)

$$J^o({\mathbf x}) = \frac{1}{2}\sum_{i=0}^N \left({\mathbf d}_i-{\mathbf G}_i \delta ... ...{\mathbf R}_i^{-1} \left({\mathbf d}_i-{\mathbf G}_i \delta {\mathbf x}\right),$$ (5.41)

où ${\mathbf d}_i={\mathbf y}^o_i-G_i{\mathbf x}^b={\mathbf y}^o_i-H_i{\mathbf x}^b(t_i)$ représente l'innovation au temps $t_i$, c'est-à-dire l'écart entre les observations et l'équivalent de l'ébauche donné par l'opérateur d'observation généralisé $G_i$ dans l'espace des observations à chaque temps $t_i$.

La fonction coût $J=J^b+J^o$ du 4D-Var incrémental à contrainte forte est quadratique et la minimisation possède une solution unique. Si la linéarisation des opérateurs $H_i$, $M_{0 \to i}$ et $G_i$ est exacte, alors la solution est identique à celle obtenue par le filtre de Kalman étendu.

L'incrément d'analyse qui minimise la fonction coût donné par les Eqs. 5.40 et 5.41 est

$$\delta{\mathbf x}^a = {\mathbf B}{\mathbf G}^T\left({\mathbf G}{\mathbf B}{\mathbf G}^T+{\mathbf R}\right)^{-1}{\mathbf d}.$$ (5.42)

La fonction coût incrémentale (Eqs. 5.40 et 5.41) est minimisée par une méthode itérative de descente. Cette minimisation nécessite le calcul de la fonction coût et de son gradient à chaque itération de la minimisation comme pour l'algorithme classique du 4D-Var. Avant le début de la minimisation, l'état d'ébauche à $t_0$ noté ${\mathbf x}^b$ est propagé par le modèle non-linéaire permettant le calcul des innovations ${\mathbf d}_i$ à chaque instant d'observation $t_i$.

La figure 5.8 permet de représenter simplement l'utilisation du 4D-Var incrémental. Elle est à comparer aux Figs. 5.6 et 5.7 qui représentent, selon les mêmes codes, les 3D-Var FGAT incrémental avec ou sans IAU.

Figure: Illustration de la procédure pour cycler le 4D-Var incrémental. Pour chaque cycle $c$, le modèle d'évolution est intégré de $t_0$ à $t_N$ à partir de l'état initial d'ébauche ${\mathbf x}^b_c(t_0)$ (courbe noir pleine) et le vecteur d'innovation $\mathbf{d}_i$ est calculé pour les différentes observations ${\mathbf y}^o_i$ avec $i=1,\cdots,N$ (ligne fine verticale). L'analyse est effectuée à l'instant $t_0$ en ramenant les innovations à l'instant $t_0$ à l'aide du modèle adjoint. Après la minimisation, un incrément est obtenu qui est rajouté à l'état de l'ébauche initial pour obtenir l'état analysé ${\mathbf x}^a_c(t_0)={\mathbf x}^b_c(t_0)+\delta{\mathbf x}^a$. Cet état analysé tient compte de la dynamique du modèle de sorte que la trajectoire analysée (courbe grise pointillée) minimise au mieux l'écart aux observations tout au long du cycle d'assimilation (de $t_0$ à $t_N$). L'état analysé ${\mathbf x}^a_c(t_N)$ est ensuite utilisé comme état initial d'ébauche pour le cycle suivant.