\documentclass[a4paper,12pt]{article}
\usepackage{latexsym}
\usepackage{amsmath}
% amsfonts permet le \mathbb(C) attention, peu de lettres disponibles
\usepackage{amsfonts}
\usepackage{epsfig}
\usepackage{eepic,psfig,pstricks,multirow}
\usepackage{floatflt}
\usepackage{fancyhdr}
\usepackage{footnpag}
\usepackage{rotating}
\usepackage{euscript}
\usepackage{tabularx,xspace}
%\usepackage{minitoc}
\usepackage[frenchb]{babel}
\usepackage[latin1]{inputenc}

\topmargin -2cm
\addtolength{\textheight}{4cm}
\tolerance=8000
\pretolerance=4000
\textwidth 16cm
\oddsidemargin 0.5cm

%

\newfont{\mfman}{logo10 at 10pt}
\newcommand{\MF}{{\mfman MINI}\-{\mfman KER}}
\newcommand{\norm}[1]{\parallel{#1}\parallel_2}
\newcommand{\bra}[1]{<{#1}\mid}
\newcommand{\ket}[1]{\mid{#1}>}
\newcommand{\mat}[1]{\underline{\underline{\mathcal{#1}}}}
\newcommand{\diag}[1]{{#1}^{d}}
\newcommand{\mean}[1]{\overline{{#1}}}
\newcommand{\scal}[2]{{<\, #1 \mid #2>}}
\newcommand{\ext}[2]{{\mid\, #1 >< #2\mid}}
\newcommand{\scalm}[3]{{<\, #1 {#2}^{\dagger} \, #2 #3>}}
\newcommand{\qrayl}[3]{\frac{\scalm{#1}{#2}{#3}}{\scal{#1}{#3}}}
\newcommand{\trac}[1]{{\mathcal{T}r (#1)}}
\newcommand{\cad}{c'est-\`a-dire }
\newcommand{\dgdeta}{\overline{ \overline{ \frac{\partial G}{\partial \eta}}}}
\newcommand{\dgdphi}{\overline{ \overline{ \frac{\partial G}{\partial \varphi}}}}
\newcommand{\dgdt}{\frac{\partial \vec{G}}{\partial t}}
\newcommand{\dfdeta}{\overline{ \overline{ \frac{\partial f}{\partial \eta}}}}
\newcommand{\dfdphi}{\overline{ \overline{ \frac{\partial f}{\partial \varphi}}}}
\newcommand{\dfdt}{\vec{ \frac{\partial f}{\partial t}}}
\newcommand{\vecdeltaeta}{\overrightarrow{ \delta \eta}}
\newcommand{\vecdeltaetadec}{\overrightarrow{ \delta \eta}_{dec}}
\newcommand{\vecdeltaphi}{\overrightarrow{ \delta \varphi}} 
\newcommand{\vecdeltaphiins}{\overrightarrow{ \delta \varphi}_{ins}}
\newcommand{\vecdeltaz}{\overrightarrow{ \delta z}} 
\newcommand{\vecdeltazins}{\overrightarrow{ \delta z_{ins}}} 
\newcommand{\vecdeltazinsst}{\overrightarrow{ \delta z^*_{ins}}} 
\newcommand{\vecdeltax}{\overrightarrow{ \delta x}}  
\newcommand{\vecdeltaxins}{\overrightarrow{ \delta x_{ins}}} 
\newcommand{\matra}{\overline{ \overline{ A}}} 
\newcommand{\matramun}{\overline{ \overline{ A^{-1}}}}
\newcommand{\matrb}{\overline{ \overline{ B}}} 
\newcommand{\matrcplus}{\overline{ \overline{ C^{+}}}} 
\newcommand{\matrcbardt}{\overline{ \overline{\mathbb{C}_{(\delta t)}}}}
\newcommand{\matrcbardtst}{\overline{ \overline{\mathbb{C}^*_{(\delta t)}}}}
\newcommand{\matrd}{\overline{ \overline{ D}}} 
\newcommand{\matrf}{\overline{ \overline{ F}}}
\newcommand{\matri}{\mathbb{I}}
\newcommand{\flo}{\Psi^{F}}


\usepackage[round]{natbib}
%%\bibliographystyle{jas}

%
\def\todec{{\overline {\overline {\Upsilon  }}}_{,dec} }
\def\tox{{\overline {\overline {\Upsilon  }}}_{,x} }
\def\dt{{\delta t}}
\def\edtaux{{exp (-{\delta t  \over \taux } ) }}
\def\deta{{\delta \eta } }
\def\Deta{{\Delta \eta } }
\def\Dt{{\delta t } }
\def\dfi{{\delta \varphi } }
\def\Dfi{{\Delta \varphi } }
\def\dts{{\delta t\over 2 } }
\def\dtaf{{\overline {\overline {\partial _\eta f }}} }
\def\d(#1)G{{\overline {\overline {\partial_#1 G}}} }
\def\dA{{\overline {\overline {\partial A}}} }
\def\dtA{\partial_t A}
\def\VD{{\mathcal V \mathcal D}}
%à ù
\newcommand{\Minik}{Mini$\_$ker\xspace}
% essai pour subsection a la JAM, cad italique a. titre
\makeatletter
% d'abord, numero = alphabetique
\renewcommand{\thesubsubsection}{\alph{subsubsection}.}
% style de la forme subsection
%\renewcommand{\subsubsection}{\@startsubsection
% {subsubsection}%                      %the name
% {2}%                               % the level
% {0mm}%                             % the indent
% {-\baselineskip}%                  % the before skip
% {0.5\baselineskip}%                % the after skip
% {\normalfont\large\itshape}}% % the style
%%  {\normalfont\normalsize\itshape}}% % the style
\makeatother 

%

\title{Analyse de la rétroaction vapeur d´eau dans un modèle de climat soumis
à un forçage périodique}\date{Juin 2007}
\author{Pablo Garcia Ramos}
%
\begin{document}
%
\maketitle
%
{\it Un modèle de circulation générale atmosphérique (AGCM) est un modèle difficile à interpréter
du fait de la quantité de variables qui sont définies et la difficulté des couplages complexes
se produisant dans une gamme d'échelles spatio-temporelles très large.
On se retrouve dans la nécessité de faire des modèles plus simplifiés pour comprendre les phénomènes
physiques les plus importants vis à vis d'un processus comme celui du changement climatique.
Le but de mon PPL est d´essayer de comprendre l'influence de la vapeur d´eau dans la dynamique de
l´atmosphère en la supposant soumise à un forçage périodique saisonnier. Cette étude est réalisée
avec un modèle mono-dimensionnel du climat établi par par Stéphane Hallegate pendant sa thèse~(2003-2005).
Son modèle original ayant été implémenté dans un logiciel différent que celui que j'ai utilisé,
une partie de l'exercice a consisté à établir la cohérence entre les deux séries de résultats.

Stéphane a construit un modèle destiné à l'étude de la rétroaction vapeur d'eau en réponse
à un doublement de CO2. Le modèle comprend trois couches, cinq variables d´état et 12 variables
de transfert en utilisant
le formalisme du TEF\footnote{Transfer Evolution Formalism, voire annexe.} pour la définition et résolution du
problème dynamique. Le TEF est un formalisme qui nous permet de bien voir l'influence des variables
de flux qui couplent les éléments du modèle, à savoir ici l'océan, la troposphère et la stratosphère.
Cela est rendu possible en supprimant l'influence directe entre les variables d´état grâce à la
définition des variables de flux qui rétablissent les couplages entre les différentes variables
d´état.

De plus, pour voir l'influence de la vapeur d´eau dans ce modèle d´atmosphère simplifié, on a utilisé
le concept de sensibilité, et de gain de rétroaction. Ce dernier concept introduit
par~\cite{Bode45} pour analyser les circuits électriques linéaires a été
largement utilisé en climatologie depuis les
années~70~(\cite{Bony06,Hansen84,Manabe67}). Mais ce concept jusqu'à présent a été utilisé uniquement
dans le cas d´une atmosphère en
équilibre~(voir par exemple \cite{HLG05}) et on sait que l´atmosphère n'est jamais en équilibre.
On a donc profité de la nouvelle implémentation du modèle pour analyser les effets dùs aux
trajectoires non stationnaires du climat.

S'agissant d'un système périodique, on a ensuite
utilisé la théorie de Floquet pour décomposer l'influence de la vapeur d´eau dans une
atmosphère soumise au cycle saisonnier. Cette étude exploite ainsi la définition
du Système Linéaire Tangent aux trajectoires en profitant de la théorie de Floquet, aujourd'hui
abondament utilisée dans l'analyse de stabilité de nombreux systèmes physiques oscillants.
On a alors mis \oe uvre les algorithmes décrits dans les deux rapports de stage de Céline
Conrardy\footnote{cf pages web de ZOOM.}.} 

\newpage
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Le modèle de Stéphane Hallegate}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%
\subsection{Description et hypothèses}
%

\medskip
On a un modèle à une seule dimension, première hypothèse simplificatrice qui supprime les aspects
de circulation (perturbations synoptiques, vents) considérés ici comme phénomènes à temps
caractéristiques rapides (typiquement la semaine).
On considère ensuite trois couches~: Océan, Troposphère, et Stratosphère.
Le couplage de ces trois couches est fait par des équations de bilan de énergie en définissant les
transferts d'énergie et de masse entre la Troposphère et l'Océan, entre les trois couches et
l´espace (flux radiatifs courtes et grandes longueurs d'onde).
 Ces transferts sont définis
de façon très simplifiée sauf pour les échanges radiatifs, car ils présentent des processus cruciaux
pour la rétroaction vapeur d'eau. Pour simplifier les équations des transferts, les hypothèses
simplificatrices on a été faites en utilisant des résultats de GCM et de faits connus sur le
comportement moyen de l'atmosphère.

Surtout ces hypothèses prennent en compte la quasi impossibilité de modéliser la convection profonde
et la couche limite océanique de manière simplifiée. Ainsi il y a une grande difficulté pour
définir les échanges de chaleur et vapeur d´eau entre la Troposphère et l´Océan. Le modèle
a retenu comme principales hypothèses~:
%
\begin{itemize}
\item	les effets de la convection sont implicitement définis par les gradients verticaux de
 température moyens; en effet, il a été constaté dans les GCM que le rôle des phénomènes convectifs
revient a tendre vers ces gradients ``cibles'';
%
\item	la couche de mélange océanique (CMO) a une épaisseur fixée; son épaisseur est choisie en rapport
avec une épaisseur moyenne sur quelques mois;
%
\item	l'absorption du flux radiatif d´ondes courtes est fixé pour le Stratosphère~(présence
d'ozone); ses échanges avec la troposphère sont négligés~(rôle de la tropopause);
%
\item	il y a trois couches nuageuses constantes dans la troposphère, fixées à 2000, 6000 et
14000~m qui interviennent à la
fois dans le bilan de ``shortwave'' (SW) et de ``longwave''~(LW). On a la liberté de choisir les
paramètres optiques de ces trois couches, mais on les choisit pour établir des flux moyens réalistes.
Ces flux sont à l'origine du ``moteur thermique'' climatique;
%
\item	les équations concernant la précipitation et l'évaporation sont les plus difficiles à
construire, car elles doivent prendre en compte à la fois l'apport de vapeur d'eau dans la
troposphère et l'effet moyen des précipitations sur le réchauffement troposphérique. Nous donnons
à présent plus de détail sur ces deux modèles.
\end{itemize}
%
La précipitation est définie comme égale à l'évaporation plus un terme qui prends en compte
l'excés d'humidité relative par rapport à une valeur moyenne d'équilibre égale a 0.7,
constatée dans les résultats des GCM. En effet, avec les GCM, on peut voir que l'humidité relative
reste, dans un délai de l'ordre de la semaine, en équilibre, toujours proche de 0.7~.

L'évaporation est la source de vapeur d'eau de la troposphère. Avec cette définition on ne
prend pas en compte les fluctuations sur les temps courts et on garde une tendance globale. Cette
évaporation est définie proportionnellement à une différence entre la quantité d´eau dans la
couche saturante collée à l'océan et celle de la couche limite atmosphérique, partie basse de la
troposphère où s'effectuent ces échanges jour après jour, avec des temps courts.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{Analyse aux sensibilités}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\medskip
Dans notre étude on s´intéresse à la dynamique et surtout à l'évolution de notre système lorsque il
est perturbé. De cette façon on peut nous étudier l'évolution du climat lorsque on perturbe
la concentration de CO2. On a un système en équilibre et on perturbe soit les conditions initiales,
soit via un forçage le long la trajectoire. On suppose que l'amplitude des perturbations est dans
un ``rayon de linéarité'' et que SLTC nous donne alors cette réponse.
On a donc un système linéaire de cette façon (hors TEF)~:
%
\begin{eqnarray}
\partial_t x(t) &=&  A(t)\, x(t) + B(t)\, u(t)
\label{sltc}\end{eqnarray}
%
Pour étudier les sensibilités de ce système on va utiliser le concept de propagateur. Le concept
du propagateur nous permet d´introduire avec élégance le système des sensibilités. La définition
du propagateur est donnée par~:
%
\begin{equation}
x(T) = \Phi (T,0)x(0) + \int_0^{T} {\Phi (T,t)B(t)u(t)} dt
\label{sltphi}\end{equation}
%
dans laquelle $x(t)=\Delta \eta(t)$ est l'écart à la trajectoire d'une variable considérée,
$u(t)$ la perturbation d'un forçage ayant été appliqué pour trouver la
trajectoire de référence ($\eta(t), t\in [0,T]$).
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{Sensibilité aux états}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
\medskip
On considère l'effet d'une perturbation des seules Conditions Initiales (C.I.)
$x(\tau)$, et même d'une seule composante $x_j(\tau)$. Alors
%
\begin{equation}
x(t) = \Phi (t,\tau)x(\tau)
\end{equation}
%
où le vecteur $x(\tau)$ ne comporte qu'une composante non nulle.
On voit donc que $x(t)$ est l'effet à $t$ de cette perturbation de l'état
$x_j(\tau)$, et ainsi que la j-ième colonne de $\Phi$ est le vecteur :
%
\begin{equation}
\frac{\Delta \eta(t)}{\Delta \eta_j(\tau)} = \,\ket{\Phi (t,\tau)_j}
\end{equation}
%
dans laquelle il faut interpréter le terme de gauche comme rapport de variations au sens de Gâteaux.
On définit ainsi la sensibilité des états aux conditions initiales en $\tau=0$ du j-ième état.

Si on suppose que $\Phi$ doit verifier:
%
\begin{equation}
\partial_t \Phi (t,\tau) = A(t) \Phi (t,\tau);
\qquad  \Phi (\tau,\tau) = I ;
\label{epro}
\end{equation}
%
on voit que ce système,en prennant,(\ref{sltphi}), peut s'écrire colonne par colonne et s'identifie
au système~(\ref{sltc}) lorsqu'on prend comme état initial
des perturbations la j-ième colonne de la matrice identité (et $B$=0).
On a donc construit l'équation d'évolution des sensibilités à un état initial,
donné par une colonne $s_{\eta}(t,0)$ de $\Phi(t,0)$.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{Sensibilité à la commande, considérée comme transfert exogène}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\medskip
La formule~(\ref{sltphi}) suggère que la sensibilité au transfert scalaire
 ($u(t)$) est dépendante de la fonction, et pour le cas $u(t)=\varpi\Upsilon(t,0^+)$
(step initial), est une colonne $s$ de l'intégrale :
%
\begin{equation}
s_{\varphi}(t,0) = \int_0^{t} {\Phi_i (t,\tau)B(\tau)}d\tau
\label{sltiphi}\end{equation}
%
Le système différentiel correspondant s'obtient par application de $\partial t$
à~(\ref{sltiphi}) :
%
\begin{equation}
\partial t s_{\varphi}(t,0) = A(t)\int_0 ^{t} {\Phi_i (t,\tau)B(\tau)} d\tau\, +B(t)
 = A(t)s_{\varphi}(t,0)\, + B(t)
\label{sensff}\end{equation}
%
On obtient ainsi un système différent de celui des sensibilité aux états ($s_{\eta}$),
mais qui s'intègre de manière semblable, à la source ($B$) près et aux C.I. nulles.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{Sensibilités dans le TEF}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
\medskip
Dans le TEF notre système a en plus des variables de transfert, ou
variables de ``Flux''. On a vu que notre système linéarisé du TEF est:
%
\begin{equation}
\begin{bmatrix} \partial_t \!\cdot\, -A  & -B \\ -C^{\dagger} & 1-
D\end{bmatrix} \quad 
\begin{bmatrix} \Deta  \\ \Dfi \end{bmatrix} \quad = \quad
\begin{bmatrix} B^u u^{\eta}\\ D^u u^{\varphi} \end{bmatrix}
\label{sltc1} 
\end{equation}
%
Suivant le même schéma que précédemment, reprenons l'équation~(\ref{sltc1})
en remplaçant $\Deta$ et $\Dfi$ par $\Phi(t,\tau)$ et $\Psi(t,\tau)$, une matrice
que l'on introduit {\it a priori} arbitrairement :
%
\begin{equation}
\begin{bmatrix} \partial_t \Phi\\ \Psi \end{bmatrix} \quad = \quad
\begin{bmatrix}  A  & B \\ C^{\dagger} & D\end{bmatrix} \quad 
\begin{bmatrix} \Phi  \\ \Psi \end{bmatrix} (t,\tau) 
\label{sltctef} 
\end{equation}
%
On constate que l'élimination de $\Psi$ fournit l'équation d'évolution
classique du propagateur :
\begin{equation}
\partial_t \Phi (t,\tau) = [A + B(I-D)^{-1}C^{\dagger}](t)\, \Phi (t,\tau)
\label{tefpro}
\end{equation}
%
Le système~(\ref{sltctef}) est donc, colonne par colonne, celui des
sensibilités aux états. La matrice $\Psi$ est donc définie par le système
linéaire~(\ref{sltctef})~:
%
\begin{equation}
\left \lbrace
\begin{array}{l}
\partial_t \Phi(t,0) =  A\Phi(t,0) + B\Psi(t,0) \\
\Psi(t,0) = C^{\dagger}\Phi(t,0) + D\Psi(t,0)
\end{array}
\right .
\label{propatef}
\end{equation}
%
Dans~(\ref{sltctef},\ref{propatef}), une j-ième colonne de $\Phi$ représente l'impact sur
$\eta(t)$ d'une perturbation unité des C.I. en $\eta_j$.
%% $v$ est une ligne (j) de $\Phi(t,\tau)$ donnant le gradient de $\eta_j(t)$ par rapport
%% à un step sur chaque composante de $\eta$ appliqué en $\tau$.
Comment interpréter $\Psi$ ? Partant du résultat\footnote{on reviendra sur la forme explicite
de $B^u$.} :
\begin{equation}
\Deta(t) = \Phi (t,\tau)\Deta(\tau) + \int_{\tau}^{t} {\Phi (t,s)B^u (s)u^{\eta}(s)} ds
\label{resol}
\end{equation}
on veut aboutir au système~(\ref{sltc1}) que l'on réécrit ainsi
%
\begin{equation}
\begin{bmatrix}  \partial_t\Deta  \\ \Dfi \end{bmatrix} \quad =
\begin{bmatrix}   A  & B \\ C^{\dagger} &D \end{bmatrix} \quad 
\begin{bmatrix} \Deta  \\ \Dfi \end{bmatrix} \quad + \quad
\begin{bmatrix} B^u u^{\eta} \\ 0 \end{bmatrix}
\label{sltc3} 
\end{equation}
%
en choisissant de ne garder qu'une perturbation $u^{\eta}$ indépendante (on la suppose
de même dimension que $\varphi$).
Dérivons~(\ref{resol}) par rapport à $t$ en appliquant~(\ref{sltctef})
%
\begin{eqnarray}
\partial_t \Deta(t) =
   A(t)\left[\,\Phi (t,0)\Deta(0) + \int_{0}^{t} {\Phi (t,\tau)B^u (\tau)u^{\eta}(\tau)} d\tau\,\right]
+ B^u (t)u^{\eta}(t)\nonumber\\
 + B(t)\left[\,\Psi (t,0)\Deta(0) + \int_{0}^{t} {\Psi (t,\tau)B^u (\tau)u^{\eta}(\tau)} d\tau\,\right]
\label{resol2}
\end{eqnarray}
%
On reconnaît $\Deta(t)$ dans la première somme, ce qui permet d'identifier avec
le système linéaire~(\ref{sltc3}) en posant 
%
\begin{equation}
\begin{boxed}
{\Dfi(t) = \Psi (t,0)\Deta(0) + \int_{0}^{t} {\Psi (t,\tau)B^u (\tau)u^{\eta}(\tau)} d\tau}
\end{boxed}
\label{resolfi}
\end{equation}
permettant de réécrire :
\begin{eqnarray}
\partial_t \Deta(t) =
   A(t)\Deta(t) + B^u (t)u^{\eta}(t) + B(t)\Dfi(t)
\label{resol3}
\end{eqnarray}
%
On a ainsi l'interprétation de $\Psi(t,0)$ donnant la sensibilité de $\varphi(t)$
à la perturbation des C.I. en $\eta$.
Dans l'égalité précédente, on peut exprimer $\Psi$ grâce à~(\ref{sltctef})
%
\begin{equation}
\Dfi(t) = [C^{\dagger}\Phi + D\Psi]\Deta(0) + \int_{0}^{t}
                                {[C^{\dagger}\Phi + D\Psi]B^u (\tau)u(\tau)} d\tau
\label{resolfi2}
\end{equation}
%
et en identifiant les termes en $\Phi$ et $\Psi$ :
%
\begin{equation}
\Dfi(t) = C^{\dagger}\Deta + D\Dfi
\label{resolfi3}
\end{equation}
%
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{Problèmes avec le point d´équilibre}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
Au début de cette étude, on a notamment rencontré des problèmes pour trouver un point d´équilibre
réaliste.
On avait des différences avec les valeurs standards~(cf \ref{shortfix}), et ayant dans notre modèle
beaucoup plus de flux détaillés, il n´était pas évident de trouver une bonne répartition
tout en assurant le réalisme en terme global. Après quelques jours on a pensé que notre problème
était dans le flux d´ondes courtes. On trouve des valeurs de référence du rapport de GIEC,
Fig.{\ref{shortfix}, et ce sont ces valeurs que l'on a essayé de reproduire.
%
\begin{figure}[h!]
\centerline{\psfig{figure=PS/Figure7.eps,%
bbllx=23pt,bblly=27pt,bburx=577pt,bbury=355pt,height=8.5cm,clip=}}
\caption{bilan des flux de référence avec lesquels on a calibré le modèle.}
\label{shortfix}
\end{figure}

Mais le modèle de
Stéphane était alors très rigide, les trois couches nuageuses étant couplées sous la forme d'un nuage
équivalent, on n´avait pas beaucoup de degrés de liberté.
Le modèle de Stéphane de flux d´onde courte était~:
%
\footnotesize
\begin{verbatim}
! ................................................
  F_S = flux_solaire_entrant; "W/m2" "solar flux"
  a_grnd = albedo_ground;
!         [transm]:r=xtransm;             >>>
!         [absorption]:r=xabsorption;>    >>>
;
! nouveauté : nuages hauts transparents
 Tr_NH=0.3;
;
! On fixe comme parametre les Cld Cover   
                             "les parametres rn_B,rn_M,rn_H sont egal a 
                             les cobertures nuageuses bass,
                                                 moyenne et haut"
 rN_tot=1.-(1.-rN_B)*(1.-rN_M)*(1.-Tr_NH*rN_H);
                            "nuage equivalent qui a une coberture egal a 
                             rN_tot"
;
z_pr: rN_TOT ;
! absorption Strato Ozone == 40 W (=> abs_o3=012
! Albedo nuages A_cld = 0.88 (=> 0.12 !)
;
! ----------------------------------------------------------------------
  eqn: flux_SW_tropo  =  rN_tot*0.12*F_S;  "Flux SW -> Tropo"
  eqn: flux_SW_strato =  0.12*F_S;          "Flux SW -> Strato" "0.12
                                              c'est la valeur 
                                             de la absortion de
                                             la strato"
  eqn: flux_SW_ocean  = (1.-rN_tot)*(1.-a_grnd)*F_S;    "Flux SW -> Sol"
! ----------------------------------------------------------------------

\end{verbatim}
\normalsize
%
On a alors changé de modèle en détaillant les propriétés des nuages, tout en conservant la
simplicité des équations. On a découplé les trois couches nuageuses, ce qui nous donné liberté
pour fixer les paramètres optiques de ces trois couches~:
%
\footnotesize
\begin{verbatim}
!-----------------------------------------------
!NUAGE HAUT
albedo_h=rN_H*(1-AL_NH-Tr_NH)*(1.-0.12)*F_S;  
                               "AL_NH,Tr_NH sont les ctes de absorption et 
                                               transmision des nuages hautes"
abs_h=rN_H*AL_NH*(1.-0.12)*F_S;
                                    "ce terme est egal a le flux shortwave qui 
                                            est absorve par les nuages hautes"
trans_h=rN_H*Tr_NH*(1.-0.12)*F_S;
                                    "ce terme est egal a le flux shortwave qui 
                                         passe pour la transparence des nuages"
F_SM=trans_h+(1.-rN_H)*(1.-0.12)*F_S; 
                                    "ce terme est egal a le flux shortwave total 
                                                qui arrive a les nuages moyenes"
!-----------------------------------------------------------
!NUAGE MOYENNE                                    "meme structure"
albedo_m=rN_M*(1.-AL_NM-Tr_NM)*F_SM;
abs_m=rN_M*AL_NM*F_SM;
trans_m=rN_M*Tr_NM*F_SM;
F_SB=trans_m+(1.-rN_M)*F_SM;
!______________________________________________________________
!NUAGE BASSE                                      "meme structure"
albedo_b=rN_B*(1.-AL_NB-Tr_NB)*F_SB;
abs_b=rN_B*AL_NB*F_SB;
trans_b=rN_B*Tr_NB*F_SB;
F_SOCEAN=trans_b+(1.-rN_B)*F_SB;

 z_pr: albedo_h,albedo_m,albedo_b;

! absorption Strato Ozone == 40 W (=> abs_o3=012
! Albedo nuages A_cld = 0.88 (=> 0.12 !)
! *** revoir ecriture SW (mais les flux doivent etre les memes)
;
! ----------------------------------------------------------------------
  eqn: flux_SW_tropo  =  abs_h+abs_m+abs_b;  "Flux SW -> Tropo"
  eqn: flux_SW_strato =  0.12*F_S;                   "Flux SW -> Strato"
  eqn: flux_SW_ocean  = (1.-a_grnd)*F_SOCEAN;
!  "Flux SW -> Sol"
! ----------------------------------------------------------------------

\end{verbatim}
\normalsize
%
Comme on peut le constater, on a distingué les propriétés LW des nuages, données par
les nébulosités, et les propriétés SW, à savoir une transparence et une absorption,
ce qui détermine aussi leur albédo.

Avec ces nouvelles équations, il était plus simple de trouver les valeurs de flux appropriées.
On a trouvé ces valeurs en prenant finalement les paramètres suivants~:
%
\footnotesize
\begin{verbatim}
! Low, Medium and High Cloud cover fixes en parametres
 rN_B  = 0.35;
 rN_M  = 0.25;
 rN_H  = 0.30;
 Tr_NH=0.53;
 Tr_NM=0.28;
 Tr_NB=0.3;
 AL_NM=0.11;
 AL_NB=0.13;
 AL_NH=0.07;
\end{verbatim}
\normalsize
%
La comparaison avec les données admises par le GIEC sont rassemblées dans le tableau suivant~:
%
\begin{table}[ht]
\begin{tabular}{|c||c|c|c|}
\hline Mode d'échanges                                   & Salby      &  modèle      \\
% & &\\
%%--------------------------------------------------------------------------------------------
\hline\hline {\bf Flux SW} \hspace{7cm}                      &    &       \\
\hline Océan                                 & 168    &  150     \\
% & &\\
\hline Troposphère                         &  (87)    &  20 (61)    \\
% & &\\
\hline Stratosphère                               &  (87)    &  41 (61)    \\
% & &\\
\hline\hline {\bf Flux LW} \hspace{7cm}                                &    &       \\
\hline Océan vers troposphère                                & 350    &  360     \\
% & &\\
\hline Océan vers stratosphère                               &  ---    & 0.4     \\
% & &\\
\hline Stratosphère vers océan      &  --- &   0.07    \\
% & &\\
\hline Troposphère vers océan                             &  324    &   326    \\
% & &\\
\hline Troposphère vers espace        &   ---    &  102     \\
% & &\\
\hline \hline Atmosphère vers espace            &    195    &   143   \\
% & &\\
\hline Océan vers espace                      &  40    &  66    \\
% & &\\
% & &\\
\hline\hline {\bf Autres flux} \hspace{6.5cm}    & & \\
 \hline Chaleur sensible            &    24    &   16   \\
% & &\\
\hline Chaleur latente                      &  78    &  74    \\
% & &\\
\hline       
% & &\\
%%----------------------------------------------------------------------------------------------
\hline 
\end{tabular}
\caption{Comparaiison avec les données moyennes climatiques admises par le GIEC}
\label{idxerr}
\end{table}
Les différences résiduelles ne sont pas gènantes; elles tiennent à la prise en compte des effets
LW de CO2 et H2O, sans autres gaz à effet de serre, et à l'absence d'aérosols dans le modèle.
Par ailleurs, on verra que ces valeurs moyennes sont sensiblement différentes
pour un modèle non stationnaire soumis au cycle saisonnier.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{Analyse et compréhension du modèle.}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{Stratosphère}  
%%%%%%%%%%%%%%%%%%%%%%%%%%%%        

La dynamique de la stratosphère ne dépend que du flux LW parce que on a fixé l'absorption de flux
d´ondes courtes. Ce flux dépend de la quantité d'énergie qui arrive depuis l'océan et la
troposphère, et de sa température qui fixe son émission. La planète est à l'équilibre dès
que le LW (vers l'espace) égale le SW absorbé.

Pour fixer les idées, nous donnons dès à présent quelques exemples de sensibilité fournie
par le modèle; le flux
de LW dépend de la température et de l'absorption par le CO2 et le vapeur d'eau~:
%
\begin{figure}[h!]
\begin{minipage}[h]{.5\linewidth}
\centerline{\psfig{figure=PS/sensstratopertco2.eps,%
bbllx=50pt,bblly=50pt,bburx=409pt,bbury=300pt,width=7.5cm,height=8.5cm,clip=}}
\end{minipage}\hfill
\begin{minipage}[h]{.5\linewidth}
\centerline{\psfig{figure=PS/sensstratopertvapeau.eps,%
bbllx=50pt,bblly=50pt,bburx=409pt,bbury=300pt,width=7.5cm,height=8.5cm,clip=}}
\end{minipage}\hfill
\caption{On peut voir que la sensibilité a la perturbation de CO2 et H2O est assez forte parce que
moins de chaleur qui viens de l' océan arrive à la stratosphère.}
\label{sensstratoco2}
\end{figure}
%

Par contre les effets de la stratosphère sur la troposphère et l'océan sont assez faibles.
On peut le voir en regardant les sensibilités à la température de la troposphère et de l'océan,
a priori les variables les plus influentes, parce que cette sensibilité est régie par le bilan
de LW~:
%
\begin{figure}[h!]
\begin{minipage}[h]{.45\linewidth}
\centerline{\psfig{figure=PS/sensoceanperttempstrato.eps,%
bbllx=50pt,bblly=50pt,bburx=409pt,bbury=300pt,width=7cm,height=8.5cm,clip=}}
\end{minipage}\hfill
\begin{minipage}[h]{.45\linewidth}
\centerline{\psfig{figure=PS/senstropoperttempstrato.eps,%
bbllx=50pt,bblly=50pt,bburx=409pt,bbury=300pt,width=7cm,height=8.5cm,clip=}}
\end{minipage}\hfill
\caption{On a des perturbations de l'ordre de 0.002 pour un changement de un Kelvin, valeur
négligeable, sachant en plus que dans notre problème la variation
 de la stratosphère est plus petite que le degré.}
\label{sensoceatropastrato}
\end{figure}

En plus le moteur du déséquilibre est dans les autres deux couches; on peut dire que dans ce
modèle la stratosphère n´a pas d'influence sur la dynamique du modèle et qui seulement varie
en fonction de ce qui se passe plus bas pour équilibrer le bilan total de flux de grandes
longueurs d'onde. Cet équlibre dépend des températures, mais aussi du contenu troposphérique
en CO2 et surtout en vapeur d'eau, ces deux gaz formant ce qu'on appelle couramment
{\it l'effet de serre}.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{Troposphère}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     
La dynamique de la troposphère est beaucoup plus complexe que celle de la stratosphère.
Pour l´analyser il faut voir quelques propriétés des phénomènes comme la vitesse d´actuation, la puissance
et son effet (stabilisateur or moteur). Pour expliquer ça je vais prendre comme exemple la perturbation
de la concentration du CO2 parce que le transitoire est assez complexe.

Premièrement on fait l´analyse avec la température de l´océan constante pour voir les effets moins
forts et leur influence\footnote{par fixation d'une très grande épaisseur
de la couche de mélange océanique.}.
%
\begin{figure}[h!]
\begin{minipage}[h]{.45\linewidth}
\centerline{\psfig{figure=PS/temperaturetropodoubleco2.eps,%
bbllx=50pt,bblly=50pt,bburx=409pt,bbury=300pt,width=7cm,height=8cm,clip=}}
\end{minipage}\hfill
\begin{minipage}[h]{.55\linewidth}
\centerline{\psfig{figure=PS/clapeyron.eps,%
bbllx=2pt,bblly=2pt,bburx=316pt,bbury=286pt,width=8cm,height=8cm,clip=}}
\end{minipage}\hfill
\caption{L'évolution de la température de la troposphère et la courbe de Clausius-Clapeyron.}
\label{clapeyron}
\end{figure}
%

Le doublement initial de la concentration est une perturbation qui va être seulement vue au début
par le flux LW, donc le première chose qu´on va voir c´est une croissance de la température, les effets
dans le flux LW étant très rapides parce que ces échanges sont fait à la vitesse de la lumière.
Par contre une augmentation de la température a une sensibilité négative parce que cette perturbation
produit une diminution de la précipitation,  par le fait que l'humidité relative a baissé. En effet
la relation de Clausius-Clapeyron montre que la pression saturante croit rapidement avec la température;
alors, une croissance de la température fait diminuer l'humidité relative, et ainsi les précipitations.
La Fig.{\ref{clapeyron} nous montre que $P_{sat}$ croit avec la température.

 La précipitation a un
effet puissant parce qu'elle libère dans l'atmosphère la chaleur latente élevée de l'eau,
 de l´ordre de $10^{6} J kg^{-1}K{-1}$. Dans le modèle, son effet est appliqué avec une constante de temps,
le but de ce coefficient est de prendre en compte le temps de formation des phénomèmes précipitants comme
orages, tempêtes, cumulo-nimbus etc, qui bine sùr ne sont pas explicitement représentés
dans un modèle simplifié. On pense que ce décalage peut être mieux realisé en multipliant
par une fonction aléatoire --- comme par un tirage par Monte-Carlo --- mais ce n´est pas évident a régler.

Ainsi, revenant au transitoire lié à la perturbation de CO2, on peut voir que  la croissance de la
température provoque une diminution des précipitations et donc à une augmentation de la quantité de
vapeur d´eau dans l'atmosphère, l'évaporation restant initialement inchangée.
Maintenant on a l´effet de la précipitation qui tend à baisser la température, l´effet d'accumulation
de vapeur d´eau qui accroît l'écran de grandes longueurs d´ondes (effet de serre) et ainsi la dynamique
de croissance qu´on avait auparavant~:
%
\begin{figure}[h!]
	\centering
		\epsfig{figure=PS/preevapdoubleco2.eps,width=.8\textwidth}
		\caption{ Evaporation en rouge et précips en vert. La Fig.~\ref{clapeyron} montre les effets
de la précipitation sur la température de la troposphère.}
	\label{precipevap}
\end{figure}
%

Il y aura un moment où les trois phénomènes s´équilibrent et à ce moment le terme d'évaporation s'impose
et la température commence à baisser. Donc le  maximum local qu´on voit sur la figure est produit
par la précipitation. La température baisse par diminution de la précipitation, un effet qui va perdre
force parce que la température décroît et parce que le quantité d´eau a monté beaucoup. Au moment
où la précipitation recommencera à monter et a s'équilibrer avec l'évaporation. Pendant la durée
de décroissance de la température de la troposphère, celle-ci gagne néanmoins de l'énergie sous forme
de chaleur latente. Avec la croissance de la précipitation la troposphère tends a s´équilibrer.
On a vu que la précipitation à un gros effet sur la dynamique de temps courts parce que sa valeur
c´est liée à la chaleur latente relâchée sur quelques jours dans le modèle.

%%%%%%%%%%%%%%%%%%%%%
\subsubsection{Océan}
%%%%%%%%%%%%%%%%%%%%%%

\medskip

L´effet de l´océan est l'effet le plus important sur le long terme parce que l´océan a une grande
inertie~: au début il varie peu mais dès que sa température prend un profil de croissance, il est
difficile à arrêter comme on va l'expliquer (rétroaction positive). L'océan a trois modes
 d'interaction avec la troposphère,
deux directs et un indirect dans le modèle comme dans la réalité.

 Les deux modes directs sont le
flux de chaleur sensible et le flux LW. La façon indirecte est l'évaporation parce que la précipitation
dépend dans le modèle de l'évaporation, ce qui est une manière de prendre en compte le fait
qu'aux temps de réponse supérieurs à quelques jours, l'eau est recyclée entre la troposphère et l'océan.
Les flux LW et l'évaporation ont une valeur plus puissante que celle de le chaleur sensible, 
celle-ci tendant a égaliser continûment les températures de la troposphère et de l´océan par les
échanges turbulents au sein de la couche limite atmosphérique; sans transport de vapeur d'eau et
intervention de la chaleur latente, ces flux `secs'' sont nettement moins efficaces.

En regardant les sensibilités à une perturbation de la température océanique de ces flux
(Fig.~{\ref{sensco2}), on peut dire que le flux LW augmente avec la perturbation, donc c´est un flux
qui tends a faire croître la perturbation initiale --- on dit alors que la rétroaction est positive.

Par contre l'évaporation augmente également avec la température, mais comme elle a un signe négatif
dans le bilan d´énergie de l'océan, elle a un effet stabilisateur --- on dit alors qu'il s'agit
d'une rétroaction négative. Donc on peut voir que lorsque la température de l´océan commence a croître,
le rapport entre les deux flux tend à augmenter la température jusqu´a ce que l'évaporation, aidée
par la chaleur sensible, réussit à équilibrer la température océanique.

On peut à présent examiner l'interaction complète océan-troposphère. On va reprendre la perturbation de CO2
mais cette fois on va libérer la température de l´océan. La perturbation est produite dans la troposphère
au tout début, donc la température de l'océan va varier par réaction aux effets produits dans la
troposphère, et en premier par les flux LW reçus par océan de la troposphère.
%
\begin{figure}[h!]
\begin{minipage}[h]{.5\linewidth}
\centerline{\psfig{figure=PS/senslwoceanpertco2.eps,%
bbllx=50pt,bblly=50pt,bburx=409pt,bbury=300pt,width=7.5cm,height=8.5cm,clip=}}
\end{minipage}\hfill
\begin{minipage}[h]{.5\linewidth}
\centerline{\psfig{figure=PS/senpertco2precipevapoceanlibr.eps,%
bbllx=50pt,bblly=50pt,bburx=409pt,bbury=300pt,width=7.5cm,height=8.5cm,clip=}}
\end{minipage}\hfill
\begin{minipage}[h]{.5\linewidth}
\centerline{\psfig{figure=PS/senschaleursensibpertco2.eps,%
bbllx=50pt,bblly=50pt,bburx=409pt,bbury=300pt,width=7.5cm,height=8.5cm,clip=}}
\end{minipage}\hfill
\caption{On peut voir les différences de valeurs entre les sensibilités aux perturbations de CO2
et du longwave et l'effet de la précipitation plus grande que celle de la chaleur sensible.}
\label{sensco2}
\end{figure}

 Une fois que la température
de l´océan est perturbée, il est très difficile à arrêter parce que la sensibilité du flux LW presque
tout de suite commence à renforcer cet effet. La précipitation a encore un effet aux temps courts mais
cette fois se superpose la tendance croissante de la température induite par l´océan. A partir de
quelques quarante jours dans notre modèle et notre choix de la constante de temps des précipitations,
toute la dynamique est contrôlée par la croissance de température de l´océan difficile à arrêter de
par son inertie. La précipitation augmente parce que l'évaporation augmente et cette augmentation
fait un apport positif au le bilan de la troposphère, donc la température de la troposphère augmente aussi.

On arrivera à l´équilibre de la température de l´océan (moteur du changement climatique de long terme),
mais il faut toujours dire que très lentement, grâce a l'augmentation de l'évaporation qui équilibre
le flux de grandes longueurs d´ondes. Mais cet effet a aussi la conséquence de donner
de l'énergie à la troposphère sous forme de chaleur latente (augmentation de $Q$) et par augmentation
de sa température. Cette augmentation du volume de précipitation et d'évaporation et l'augmentation
de la température de la troposphère dans le modèle correspond à la production des orages et cyclones
avec plus de fréquence et de force dans le climat réel.

%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{résumé}
%%%%%%%%%%%%%%%%%%%%%%%%

\medskip
Le véritable agent du changement climatique est l´océan, qui induit les augmentations irréversibles de
température et de la vapeur d´eau. L´océan est un moteur qui tarde un peu en démarrer mais lorsque
il a commencé il est vraiment difficile d´arrêter. La précipitation est le terme principal qui régit
la dynamique de la troposphère parce que liée à la forte valeur de la chaleur latente de l'eau.
Il réagit on peut dire avec indépendance de la température de l´océan à court terme, mais lorsque
la température de l´océan monte, la précipitation est contrôlée par l'évaporation croissante. L'évaporation
est ainsi le terme qui tend à équilibrer la température de l´océan. Il est un terme puissant qui a
trois conséquences~: i) il peut équilibrer plus vite, mais par contre, ii) il donne plus d'énergie
à la troposphère, et cette énergie peut produire des phénomènes météorologiques adverses, et
enfin~: iii) il renforce l'effet de serre.
%               

%%%%%%%%%%%%%%%%%%
\section{Analyse du système périodique}
%%%%%%%%%%%%%%%%%%
On le sait, le système climatique n'est pas à l'équilibre. On a donc soumis le modèle à un forçage
saisonnier comme on va le décrire plus loin. On voit alors que la sensibilité du système
au doublement de CO2 est amplifiée par le cyclage~:
%
\begin{figure}[h!]
\centerline{\psfig{figure=PS/ClCyclStaReEff4b.eps,%
bbllx=50pt,bblly=50pt,bburx=554pt,bbury=775pt,width=11.5cm,angle=-90,clip=}}
\caption{On peut voir que la sensibilité a la perturbation de CO2 est amplifiée
par le cycle (température de la tropo).}
\label{senscystaco2}
\end{figure}
%
Après une dizaine d'années, la troposphère se stabilise à environ 0.3 K de plus qu'en perturbation statique,
et reste alors une oscillation résiduelle autour de $T_0$+1.5 K. La théorie de Floquet a pour
but d'analyser la stabilité des systèmes périodiques; c'est ce que nous abordons à présent.
%
%%%%%%%%%%%%%%%%%%
\subsection{Théorie de Floquet}
%%%%%%%%%%%%%%%%%%

\medskip
Dans cette partie notre objectif est d´analyser notre système avec la théorie de Floquet. 
L'analyse de Floquet s'applique aux systèmes dynamiques $\partial_t x = A(t)\, x$ lorsque
la matrice d'avance de phase est périodique : $A(t+T)=A(t) \quad\forall t$.
Les propriétés du propagateur associé sont classiquement obtenues en considérant la
matrice $X(t)$ des solutions fondamentales du système --- qui forme une base de toute les
solutions au système. On pose d'abord que $X(0)=I$ pour simplifier l'écriture.
Alors, on a 
\begin{eqnarray}
 \partial_t X(t+T) &=& A(t)\,X(t+T) 
\end{eqnarray}
%
qui prouve que $X(t+T)$ est également solution, et peut donc s'exprimer dans la base
 $X(t)$~:
%
\begin{eqnarray}
 X(t+T)&=&X(t)\, M 
\label{xpm}\end{eqnarray}
%
où $M$ est une matrice constante (et régulière, puisque les $X$ le sont par définition).
La théorie a été établie pour les systèmes linéaires.
Dans notre problème ce qu'on apelle $X(t)$ est ici $\delta\eta(t)$ parce que c'est la perturbation
qui est régie par un ssytème linéaire.
Cette perturbation n'est pas périodique alors que $\eta(t)$ l'est, et ainsi $A(t)$. Cela nous dit
qu' en connaissant une perturbation initiale $\delta\eta(0)$ et son évolution jusqu'à $\delta\eta(T)$,
on peut connaitre toute l'évolution du système par rapport a cette matrice $M$ cte.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{Factorisation de Floquet}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\medskip
La théorie de Floquet nous dit que lorsque on a un système de flot periodique, on peux décomposer
la matrice $X(t)$ en une matrice périodique et une autre exponentielle. On va le voir en définissant
la matrice constante $C$ (sur le corps des complexes) telle que $\displaystyle X(T)=e^{CT}$,
et on a en reprenant l'équation (\ref{xpm}) pour $t=0$ 
%
\begin{eqnarray}
 X(T)&=&M=e^{CT}     
\end{eqnarray}
%
De cette facon on peux maintenant définir la nouvelle matrice $B$ qui est périodique~:
%
\begin{eqnarray}
 X(t)&=&B(t)\, e^{Ct}  \nonumber\\
 &\Rightarrow& B(t+T)=B(t) 
\label{bct}\end{eqnarray}
%
On a en effet
%
\begin{eqnarray}
X(t+T)&=&B(t+T)e^{Ct}e^{CT}\nonumber\\
     &=&X(t)\, e^{CT}  \nonumber\\
     &=&B(t)e^{Ct} e^{CT}  \nonumber
\end{eqnarray}
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{Vecteurs de Floquet}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Lorsque on a maintenant la propriété~(\ref{bct}) et comme on sait que le
propagateur est lié à $X$ on va voir ce qu'on peut sortir du propagateur.
Par analogie avec la partie de $X$, il est logique de démarrer l'étude à
l'instant $t=T$. On diagonalise $\Phi(T,0)$ et on change la forme des valeurs propres
en cherchant quelque chose de similaire au $M$.

Donc on a 
des valeurs propres $\sigma_i =e^{\lambda_i T}$ du propagateur en $t=T$ et
on définit $Y$:
%
\begin{eqnarray}
 \Phi(T,0)Z&=& Z\,e^{\Lambda T} \nonumber\\
 Y(t)&\stackrel{def}{=}&\Phi(t,0)\,Z  \nonumber\\
 Y(0)&=&Z  \nonumber\\
\end{eqnarray}
%
Avec cette définition on a la relation de $Y(t)$ avec sa réplique à $t+T$ est alors
%
\begin{eqnarray}
 Y(t+T)&=&\Phi(t+T,0)\,Z = \Phi(t+T,T)\,\Phi(T,0)\,Z   \nonumber\\
       &=&\Phi(t,0)\,Z\, e^{\Lambda T}    \nonumber\\
 Y(t+T)&=& Y(t)\, e^{\Lambda T}
\label{floq1}\end{eqnarray}
%
Par ailleurs, $Y$ s'intègre comme le propagateur, car
%
\begin{eqnarray}
 \partial_t Y &=& A(t)\,\Phi(t,0)\, Z=A(t)\, Y(t)
\label{dty}\end{eqnarray}
%
Ainsi~(\ref{floq1}), $Y(t)$ n'est pas périodique mais est ``amorti'' à chaque période
par $\displaystyle e^{\Lambda T}$.
On a simplement amélioré la décomposition donnée par~(\ref{bct}) en découplant les modes d'amortissement.

Maintenant on peut construire à partir des $Y$ un autre vecteur non amorti pour capter
la structure du champ associée aux phénomènes fluctuants; ce sont les {\bf vecteurs de Floquet} $\Psi^F$
que l'on définit en ``dés-amortissant'' $Y$ :
%
\begin{eqnarray}
 \Psi^F (t) &\stackrel{def}{=}& Y(t)e^{-\Lambda t} 
\label{vfloq}\end{eqnarray}
%
on obtient bien un vecteur périodique, car :
%
\begin{eqnarray}
 \Psi^F (t+T) &=& Y(t+T)\, e^{-\Lambda T}e^{-\Lambda t}   \nonumber\\
             &=& Y(t)e^{-\Lambda t}=\Psi^F (t)   \nonumber
\end{eqnarray}
%                                                                          
Et on a arrive a décomposer, plus au moins, le propagateur en des orbites périodiques amplifiées
ou amorties; donc on a gagné pas mal d'information~:
%
\begin{equation}
\Phi (t) = \Psi^F (t)e^{\Lambda t}Z
\end{equation} 
%
C'est l'un au moins de ces vecteurs que l'on tente de déterminer en avant, mais notre développement considère
l'ensemble de ces vecteurs par souci de généralité~:
%
\begin{eqnarray}
 \partial_t \Psi^F &=& AY \,e^{-\Lambda t}- Y(t)\Lambda e^{-\Lambda t} \nonumber\\
                   &=& A\Psi^F (t) - \Psi^F e^{\Lambda t}\Lambda e^{-\Lambda t}  \nonumber\\
\label{floq2}\end{eqnarray}
%
et du fait de la diagonalité des trois matrices de droite, le résultat final est
%
\begin{equation}
\boxed{
\partial_t \Psi^F = A\Psi^F (t) - \Psi^F \Lambda \quad {\mbox avec}\,\Psi^F(0)=Z   
}
\label{floqfin}\end{equation}
%
Ainsi, connaissant la période $T$ de la matrice d'avance de phase d'un système dynamique linéaire
(ou du SLTC en perturbation), et ayant calculé un ou plusieurs vecteurs propres du propagateur
à $t=T$, on est capable d'intégrer en avant un ou plusieurs vecteurs de Floquet.

\paragraph{passage au TEF} On peut remarquer que l'on a un autre ``vecteur de Floquet'' lié au TEF,
à savoir $\zeta^F(t)=\Psi(t,0)Zexp(-\Lambda t)$. Ces transferts de Floquet correspondent aux
transferts perturbés des vecteurs de Floquet, et sont donc très intéressants à analyser pour
interprèter physiquement les Floquets.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{résumé}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
1.	La théorie de Floquet est applicable au systèmes ou la matrice d'avance de phase
  est périodique~: $A(t+T)=A(t) \quad\forall t$ . C´est l'unique condition qu´on a pour l´appliquer;

2.	Lorsque on a un le flot périodique la théorie nous dit qu´on peut décomposer le propagateur
 en vecteurs périodiques qui sont exponentiellement amplifiés ou réduits;

3.	A vecteurs périodiques, trajectoires périodiques en lesquelles on peut décomposer le propagateur,
 on les appelle Vecteurs de Floquet;

4.	les coefficients qui déterminent si une trajectoire périodique est amplifiée ou réduite sont liés
 aux valeurs propres du propagateur pour la période du système.

5.	Le propagateur nous donne les sensibilités du système aux perturbations des variables d´état.
 Maintenant ces sensibilités sont décomposées en trajectoires qui vont monter ou décroître.
 Donc comme premier résultat de l´analyse on peut voir l´évolution de la dimension du système avec le temps,
et connaitre ainsi sa stabilité;

6.	Par contre donner un sens physique aux vecteurs de Floquet n´est pas un problème facile.
 Et on ne sait pas (a priori) quelle information on peut sortir de notre analyse.
%
%%%%%%%%%%%%%%%%%%%%%%
\subsection{Problèmes du démarrage}
%%%%%%%%%%%%%%%%%%%%%%%%
Dans ce moment on est en train de commencer à étudier l'information que l´analyse de Floquet nous donne.
 Les résultats plus complets de notre étude, je les commenterai dans l´exposition du 9 Juillet.

Avant de procéder à l'analyse, on est obligé de faire que le système devienne périodique.
On est arrivé a le faire
en mettant un flux solaire périodique. De cette façon on essayé de reproduire le caractère saisonnier
du climat. Donc au lieu d´avoir~:
%
\begin{eqnarray}
FS & = & 343  \quad (W m^{-2})  \nonumber\\
FS & = & 343\,(1+\frac{1}{3}\sin(\frac{2 \pi t}{365}))  \quad (W m^{-2})  \nonumber
\end{eqnarray}
%

On a pensé que cette modification suffirait parce que le flux longwave, très important dans
l´évolution du système, est lié au shortwave. Pour tester si le système est périodique on exige
que toutes les variables d´état soient périodiques. De cette façon on est sùr que la
 matrice A est périodique, condition du Floquet. On a vu que les variables étaient décalées en
phase avec le forçage, ce qui provoque des transitoires assez longs à stabiliser parce que
la constante de temps de notre système est trés grande. La façon d'évoluer nous a donné
l'impression que le système avait un comportement pas très linéaire. Pour tester ça on a fait
un fit avec la fonction($a\sin{(\frac{2\pi}{365}x+c)}+d$) de notre signal pour avoir la valeur moyenne.
On a alors constaté une variation de valeur moyenne par rapport a celle du problème non périodique et en
plus que le signal n´était pas purement sinusoïdal. Donc notre système n´a pas un comportement linéaire
pour l'amplitude du forçage choisi.

En effet lorsque on a un système linéaire on peut découpler les équations de la perturbation
et la valeur moyenne~:
%
\begin{figure}[h!]
	\centering
		\epsfig{figure=PS/fitstrato.eps,width=.9\textwidth}
\caption{Le fit nous montre que l'évolution n'est pas vraiment sinusoïdale et aussi des
différences en valeurs moyennes.}
\label{fit}
\end{figure}

\begin{eqnarray}
x & = & \overline{ x} + x\prime \nonumber\\
\partial_t (\overline{x + x\prime}) & = & F (\overline{x} + x\prime) \nonumber\\
\partial_t \overline{x} + \partial{x\prime} & = & F(\overline{x}) + F(x\prime)  \nonumber
\end{eqnarray}
%
On peut faire ça si $F$ est un application linéaire. Et on a aussi si $F$ est linéaire~:
%
\begin{eqnarray}
\overline{\partial_t x} & = & \overline{ F ( x )} \nonumber\\
\partial_t \overline{x} & = & F (\overline{x}) \nonumber\\
\partial_t x\prime & = & F (x\prime) \nonumber
\end{eqnarray}
% 

On a vu que l´effet de la constante de temps est très forte; de plus notre système changeait
de tendance lorsque on le redémarrait avec le dernier point du calcul. Ça nous a fait réfléchir
sur notre choix de pas de temps. En plus quand on a démarré l´analyse de Floquet avec le système
finalement équilibré on a vu que les valeurs propres du propagateur avaient une sortie avec beaucoup
de bruit qui nous a fait penser à un pas de temps trop fort (un jour). L'évolution des valeurs propres
était très mauvaise pour notre étude parce que toutes les valeurs propres sauf une décroissaient
très rapidement et dans la première période il nous en restait seulement une. Donc de cette
façon notre analyse de Floquet n´est pas trop intéressante parce que il nous reste qu'une dimension
quand on commence l´étude. Ainsi on a fait des modifications pour gagner une deuxième valeur propre
pour notre étude. Le deuxième valeur propre perdait son influence à bout de 60 jours. En regardant
la sortie des sensibilités on a vu que le effet precip-evap étaient aussi valable dans cette période.
Par conséquent peut être possible que ce deuxième vecteur de Floquet garde l'information
sur l´effet de ces
deux flux et que le premier vecteur garde l'information sur l'évolution régie par la température
de l´océan. Pour résoudre ces problèmes on a :

1.	changé le pas de temps à 0.1 jour;

2.	changé la cte de temps de la précipitation à 180 jours pour ralentir son effet;

3.	augmenté l'amplitude du forçage pour perturber plus le système~:
%
\begin{equation}
%
FS  =  343\,(1+\frac{1}{2}\sin(\frac{2 \pi t}{365}))  \quad (W m^{-2})  \nonumber
\end{equation}

%
\begin{figure}[h!]
	\centering
		\epsfig{figure=PS/valeurpropreevolut.eps,width=.9\textwidth}
		\caption{La deuxième valeur propre en vert a un ordre de magnitude acceptable pour $t=365$.
Toutes les valeurs propres ont perdu le bruit.}
	\label{evolutvaleurpropre}
\end{figure}
%

De cette façon on a réussi à gagner le deuxième vecteur de Floquet mais comme notre système est toujours
très lent ça nous a pris quelques jours de calculs pour aboutir enfin au cycle limite.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsubsection{Analyse des vecteurs de  Floquet}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
\paragraph{Exposants de Floquet}\footnote{Al1~: cette analyse est ajoutée au rapport officiel
et expose les arguments du stagiaire; nous en discuterons ensuite.}
%
\medskip

Lorsque notre espace de variable d´état est de quatre dimensions on a quatre exposants de Floquet.
Notre analyse nous donne qu´en effet on garde que deux exposants de Floquet parce que les autres
deux sont rapidement amortis. De ce résultat on peut toute suite dire que le couplage entre les
variables d´état fait que notre système a deux degrés de liberté.
Il nous reste deux exposants de Floquet réels, résultat logique parce que notre système est
intrinsèquement non oscillatoire. Ce résultat est bien lié au fait qu´on trouve que le système
a une réponse avec le période du forçage. Les deux exposants de Floquet ont des constantes de temps égale à :
%
\begin{eqnarray}
\tau_1 &=& 2.3 \,\textit{ans} \\
\tau_2 &=& 80 \,\textit{jours}
\end{eqnarray}
La décomposition en vecteur de Floquet nous laisse qu´un seul vecteur pour les temps
supérieurs à 800 jours. Donc le couplage dans notre système fait qu´on a un système unidimensionnel
pour les temps très grands. De cette façon on peut dire que l'information de la trajectoire
asymptotique est gardée par le premier vecteur.
%
\paragraph {Premier vecteur compréhension}
%
\medskip

La simplicité de ce système nous peux aider a savoir la information qu´on peut gagner
en faisant une analyse de Floquet. Dans notre système on a un vecteur de Floquet pour
le trajectoire asymptotique. De toutes les perturbations avec lesquelles on peut
perturber le système Fl1 nous sélectionne les plus lentes. Les reste des directions,
si elles ne sont pas perpendiculaires à F1, vont avoir une évolution lente parce qu'elles
vont finir pour s'aligner avec la direction de F1,  évolution plus au moins lente  selon leur
produit escalaire avec F1.

Maintenant on peut regarder les différentes valeurs de F1, et comment  F1 est une perturbation
des variables d´état. D´abord c´est logique que F1 perturbe l'océan parce que c´est la variable
avec la plus grande inertie, celui qui donne l´effet lent, les autre vecteurs de Floquet
étant fortement amortis ne le perturbe pas. Mais pour que la perturbation soit la plus lente
possible, l´océan doit perturber les autres variables d´une façon optimale qui mette la machine
en route sans créer de perturbations rapides.

Cette façon on peut la voir en regardant quelques composantes de F1.
 La composante RH de F1 est presque nulle, ça veux dire que RH ne joue pas dans la trajectoire
asymptotique. Donc de toutes les façons possibles d´être couplées les perturbations T tropo
et Q restent couplées pour l'équation de Clausius Clapeyron.
%
\begin{figure}[h!]
\begin{minipage}[h]{.5\linewidth}
\centerline{\psfig{figure=PS/ClFloHumRel.eps,%
bbllx=50pt,bblly=50pt,bburx=405pt,bbury=295pt,width=8cm,height=8cm,clip=}}
\end{minipage}\hfill
\begin{minipage}[h]{.5\linewidth}
\centerline{\psfig{figure=PS/ClFloLWOc.eps,%
bbllx=50pt,bblly=50pt,bburx=405pt,bbury=295pt,width=8cm,height=8cm,clip=}}
\end{minipage}\hfill
\begin{minipage}[h]{.5\linewidth}
\centerline{\psfig{figure=PS/ClFloLWTropo.eps,%
bbllx=50pt,bblly=50pt,bburx=405pt,bbury=295pt,width=8cm,height=8cm,clip=}}
\end{minipage}\hfill
\begin{minipage}[h]{.5\linewidth}
\centerline{\psfig{figure=PS/ClFloLWStrato.eps,%
bbllx=50pt,bblly=50pt,bburx=405pt,bbury=295pt,width=8cm,height=8cm,clip=}}
\end{minipage}\hfill
\caption{Composantes des deux vecteurs de Floquet (océan bloqué en référence.}
\label{FloVecs}
\end{figure}

Les perturbations de T tropo et T océan sont couplées par le Bilan de LW et par l´effet que
l'evaporation est presque égale à la precipitation. D´où que notre dimension soit égale à un.

Comme résultat important, on peut dire que plus l'inertie de l´océan est forte,
plus il est puissant. On a pu voir que le premier vecteur de Floquet nous donne
la façon espécifique (les couplages entre variables d´état) de la tendance asymptotique
de notre système vers l´équilibre.
Et Floquet nous peux donner l´information parce que notre système a la particularité
d´avoir une variable beaucoup plus dificile à estabiliser.
%
\paragraph {Deuxième vecteur de Floquet compréhension~:} 
%
\medskip

Comme première information on peut dire que ce vecteur nous donne une direction de
perturbation rapidement amortie. Donc son effet n'est présent que de manière transitoire.
Et alors que le premier vecteur garde la façon que le système tombe asymptotiquement,
le deuxième représente une perturbation qui,  si on l'ajoute à la tendance asymptotique,
nous donne le transitoire le plus significatif.
Don c´est pas le vrai transitoire mais il caractèrise le transitoire. Par contre sa
direction, à différence du F1, on ne peut pas la lier à une seule variable, et c'est
le couplage entre les variables qui nous donne la cartéristique du transitoire.

Par rapport à F1, F2 perturbe la chaleur sensible, le flux LW de la stratosphère et le RH.
Donc interpréter ce mode est difficile lorsque toute est perturbé, comment savoir quelle est
la combinaison qui dirige et pourquoi.
On peut donc dire que l'analyse de Floquet n´est pas très utile pour analyser les transitoires
parce que on trouve une direction de perturbation difficilement assimilable à une ou deux
variables spécifiques.

A mon regard, ce F2 doit être lié au couplage de la température de Tropo et Q, après que
le temps d´action des perturbations induites par T-strato s´est amorti beaucoup plus rapide que ce F2.
%
\medskip

Maintenant on se pose la question, lorsque on remet la concentration de CO2 en variable
d'état, que doit-on trouver ?
Qu´est ce que Floquet peut nous dire sur l´effet serre ?

D´abord dans notre modèle le concentration de CO2 est fixée~: on peut la perturber,  elle gardera sa
valeur. C´est cet effet qui fait qu'on a des sensibilités asymptotiques différentes de zéro pour
la perturbation de CO2 dans notre propagateur.
Quand on a système comme ça, Floquet répond en nous donnant un multiplicateur de Floquet constant
et égal à 1. Et de cette façon quand on définit un système dynamique avec un variable bloquée,
Floquet nous va donner un vecteur de perturbations qui ne bouge pas ($\mu =1$),
donc on aura tout de suite la perturbation de notre système pour $t=\infty$.

On se retrouve aussi avec un deuxième vecteur avec multiplicateur égal à 0.6.
Donc comme on a pas changé une autre fois la variable océan plus dure à ce stabiliser
on trouvera maintenant avec ce vecteur (à présent F2), la façon que le système tombe
sur le valeur asymptotique : soit zéro si on ne perturbe pas le CO2 , \cad une perturbation
de direction parallèle a F1 si on perturbe le CO2. Lorsque on ne perturbe pas la concentration
de CO2 l'évolution de notre système est la même que le système sans avoir [CO2] comme variable d´état.
Donc on doit trouver les mêmes vecteurs de Floquet. Et donc comme les vecteurs de Floquet ne dépendent
pas de la direction de perturbation et seulement du système, on peut dire que dans notre nouveau
système on va trouver le même vecteur F1 d´avant cette fois comme F2.
De cette façon on peut dire qu'on avait bien l´évolution asymptotique du système et en
particulier aussi pour la tendance à cette valeur asymptotique qui donne l´effet serre.

Donc on a regardant F1 de notre premier système, l'information de comment définit le modèle
de Stephane la façon que le climat évolue avec l´effet serre.
 
%%%%

\newpage

%%%%%%%%%%%%%%%%%%
\section{annexe}
%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%
\subsection{Transfer Evolution Formalism}
%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%
\subsubsection{Introduction}
%%%%%%%%%%%%%%%%%

\medskip
En physique on partitionne toujours (en molécules, en particules fluides, en atomes, en éléments
différentiels) parce que c´est plus facile comprendre quand on a des éléments``isolés'' qui
sont liés les uns aux autres par des interactions. Donc pour partitionner on doit définir
les interactions entre les parties, et c´est elles qui nous aident a trouver toute l'information.

Le TEF est une formalisme de modélisation de systèmes qui suit ce principe mais dans
un espace mathématique.

On a un modèle défini par ses équations d´état.
%
\begin{equation}
\frac{d \vec{\eta} }{dt}  = \vec{ G} ( \vec{\eta},\vec{\psi}, t )
\label{deta}
\end{equation}
%   
En donnant des conditions initiales $\vec{\eta} (t_0)$, et avec
des conditions limites et forçages
connus $\vec{\psi} (t)$, la trajectoire du système est déterminée\footnote{entre
les éventuelles singularités} : on dit que le problème est ``bien posé''.

Ici on a un couplage parmi les variables d´état qui peut nous cacher des informations.
Donc on partitionne dedans notre espace de variables d´état, et on a l'espace vectoriel $\{\vec{\eta}\}$,
ce qui aboutit à la définition d'équations d'état pour chaque sous-espace $\{\vec{\eta_{\alpha}}\}$ :
%
\begin{equation}
\frac{d \vec{\eta}_{\alpha} }{dt}  = \vec{ G} ( \vec{\eta}_{\alpha},\vec{\varphi},\vec{\psi},t )
\label{deta2}
\end{equation}
%
où on remarque l'adjonction de conditions limites supplémentaires $\vec{\varphi}$. Le TEF exige
en effet que chaque élément $\alpha$ du système partitionné corresponde à un problème bien posé.
Cela a deux conséquences :
%
\begin{itemize}
\item a) chaque élément, dénommé ``cellule'', ne doit explicitement pas
dépendre de l'état d'un autre élément;
\item b) l'interface introduit par le découpage avec les autres
cellules représente de nouvelles conditions aux limites (C.L.) dépendantes.
\end{itemize}
%
Si ces nouvelles C.L. sont
connues, alors chaque cellule va suivre une trajectoire dite ``découplée''. Pour que les cellules
suivent leur trajectoires d'ensemble, déterminée par le système global, il faut donc que ces
conditions d'interface soient contraintes par le reste du système, ce qui se fait par l'introduction
de modèles d'interface, dénommés ``transferts'' :
%
\begin{equation}
\vec{\varphi}  = \vec{ f} ( \vec{\eta} , \vec{\varphi}, t )
\label{vecphi}
\end{equation}
%

L'équation~(\ref{vecphi}) exprime que les transferts sont à tout instant déterminés par les états,
et n'ajoutent donc pas de conditions initiales supplémentaires. La définition de ces transferts
correspond ainsi à la notion de {\bf raccordement}.
Ces nouvelles variables, ``les transferts'', induites par le découpage,
introduisent des modèles de couplage entre des éléments-cellules. Ils re\-présentent
couramment des contraintes de continuité, ou des flux et des forces dans les modèles physiques,
des échanges, des répartitions et des informations dans les modèles macro-économiques.

Mathématiquement, la solution du système ($\ref{deta2},\ref{vecphi}$) fournit la même trajectoire
que~(\ref{deta}). C'est de cette manière que le TEF met en \oe uvre la notion de
découpage-raccordement. On aboutit en résumé à deux classes distinctes de modèles, et aux deux
types de variables associés. 
Les {\bf cellules} sont représentées par leur état et définies par des
équations d'état, \cad donnant leur évolution, alors que les {\bf transferts}
sont des équations de contrainte statique\footnote{on peut préférer
les qualifier d'instantannées, le terme de ``statique'' étant utilisé
en théorie de la commande}.

%%%%%%%%%%%%%%%
\subsubsection{Linéarisation et résolution}
%%%%%%%%%%%%%%%

\medskip
L'étude des relations entre sous-systèmes est excessivement complexe
si on regarde le système complet non-linéaire. Toute l'étude va donc 
se faire sur le système linéarisé sur une portion de la trajectoire.
Le système linéaire obtenu dépend encore continuement du temps, et
dans un deuxième temps ce système linéaire est discrétisé. 


On considère le système $(\ref{deta2},\ref{vecphi})$, non linéaire, sur une petite
section de sa trajectoire (soit entre $t_0$ et $t_0+\delta t$), et on détermine l'accroissement $\vecdeltaeta$
de $\vec{\eta}$ et $\vecdeltaphi$ de $\vec{\varphi}$ par développement de Taylor :
%
\begin{equation}
\frac{\partial \vecdeltaeta }{\partial t}(\tau)  = \left. \vec{G} \right|_{t_0} +
                                 \left. \dgdeta \right|_{t_0}  \vecdeltaeta (\tau)+
                                 \left. \dgdphi \right|_{t_0}  \vecdeltaphi (\tau)+
                                 \left. \dgdt \right|_{t_0}  \tau
\label{pardeltaeta}
\end{equation}               

\begin{equation}
\vecdeltaphi (\tau) =                  \left. \dfdeta \right | _{t_0} \vecdeltaeta (\tau)+
                                 \left. \dfdphi \right | _{t_0} \vecdeltaphi (\tau)+
                                 \left. \dfdt \right | _{t_0} \tau 
\label{deltaphi}
\end{equation}
%
On obtient un système d'Equations aux Dérivées Partielles (EDP) baptisé Système Linéarisé Tangent
(SLT) en une portion de la trajectoire, avec $\tau \in [t,t+\delta t [$ paramétrisant
la trajectoire.


La résolution de ce système linéaire est classique et peut
être effectuée, soit en introduisant un schéma temporel, soit par transformation symbolique
(Laplace, Fourier, Borel ...). 
On considère l'application d'un schéma temporel, qui permettra de définir la trajectoire
par morceaux. Nous verrons plus loin qu'une transformation symbolique 
donnerait le même résultat.

L'application du schéma de Crank Nicolson au système~(\ref{pardeltaeta},\ref{deltaphi})
fournit le système algébrique suivant :
%                          
\begin{equation}
\left\lbrace 
\begin{array}{l}
\vspace{0.5cm}
\matra \vecdeltaeta  + \matrb \vecdeltaphi = \vec{\Gamma} \delta t  \\
\matrcplus \vecdeltaeta - (\matri  +\matrd) \vecdeltaphi = \vec{\Omega} \delta t
\end{array} 
\right.
\label{eqtef}
\end{equation}
%
avec les matrices et vecteurs introduits définis par :                                                  
% 
\begin{equation}
\left\lbrace 
\begin{array}{l}
\matra = \matri - \frac{\delta t}{2}  \left. \dgdeta \right|_{t_0}  \\
\matrb =  - \frac{\delta t}{2}  \left. \dgdphi \right|_{t_0} \\
\vec{\Gamma} = \left. \vec{G} \right|_{t_0}
                   + \frac{\delta t}{2} \left. \dgdt \right|_{_t0} \\
\matrcplus =  \left. \dfdeta \right | _{t_0}   \\
\matrd = - \left. \dfdphi \right | _{t_0}  \\
\vec{\Omega} = - \left. \dfdt \right | _{t_0}  
\end{array}              
\right.
\label{coefmat}
\end{equation}
%
La matrice $(\matri  + \matrd)$
provient de la partition effectuée sur l'ensemble des transferts, qui introduit des contraintes
statiques entre les ``objets transferts'' du TEF.

Cette formulation répond à notre programme, comme nous allons maintenant le montrer. 


Symboliquement la structure matricielle du système~(\ref{eqtef}) est la 
suivante:
%
\begin{equation*}
\left[
\begin{array}{cc}
\matra & \matrb  \\
-\matrcplus & \matri +\matrd
\end{array}              
\right]
\left[
\begin{array}{c}
\vecdeltaeta  \\
\vecdeltaphi
\end{array}
\right]
=
\left[
\begin{array}{c}
\Gamma \delta t  \\
-\Omega \delta t
\end{array}
\right]
\label{matef}
\end{equation*}
%
La résolution est conduite en éliminant d'abord les variables d'état, en 
utilisant la procédure qui permet également d'obtenir~:

\begin{equation*}
\vecdeltaeta = -\matramun \, \matrb \, \vecdeltaphi + \matramun \,  \vec{\Gamma} \, \delta t 
\end{equation*} 

$\vecdeltaeta$ est introduite dans la seconde équation:
%
\begin{equation*}        
\matrcplus \left( -\matramun \, \matrb \, \vecdeltaphi 
+ \matramun \vec{\Gamma} \delta t \right)
- \left( \matri + \matrd \right) \vecdeltaphi = \vec{\Omega} \, \delta t
\end{equation*}       

\begin{equation*}
\left(\matri + \matrd + \matrcplus \, \matramun \, \matrb \right) \vecdeltaphi =
\matrcplus \, \matramun \, \vec{\Gamma} \delta t - \vec{\Omega} \, \delta t
\end{equation*}
         
Cette équation linéaire est classiquement résolue et donne $\vecdeltaphi$. On peut alors résoudre chaque
équation d'état où tout est connu sauf $\vecdeltaeta$.
%

\clearpage
\newpage
%%%%%%%%%%%%%%%%%%%%%
\subsection{Les equations du modele}
%%%%%%%%%%%%%%%%%%%%%
\label{modelequat}
%
%* 
%* {\bf nomenclature}
%* 
%* \begin{tabular}{||ll||ll||}
%* \hline
%* \multicolumn{2}{||c||}{Subscripts \& superscripts} &
%* \multicolumn{2}{c||}{Physical quantities}     \\
%* \hline
%* $p$ & precipitation                           &  $T$ & temperature \\
%* $v$ & vapor or vaporisation                   &  $\theta$ & sea surface temperature \\
%* $e$ & evaporation                             &  $\psi$ & radiation budget \\
%* $s$ & sensible                                &  $\phi$ & Sensible and latent fluxes \\
%* $l$ & latent (for fluxes) or low (for clouds) &  $\alpha$ & Absorption coefficient \\
%* $m$ & medium                                  &  $\mathcal{T}$ & Transmission coefficient \\
%* $h$ & high                                    &  $A$ & Albedo \\
%* $C$ & clouds                                  &  $N$ & cloud cover \\
%* $S$ & short wave                              &  $C$ & heat capacity \\
%* $L$ & long wave                               &  $\tau$ & characteristic time \\
%* $\sigma$ & stratosphere                       &  & \\
%* $\tau$   & troposphere                        &  & \\
%* $\theta$ & ocean                              &  & \\
%* $b$ & bottom (e.g. of troposphere)            &  & \\
%* \hline
%* \end{tabular}
%* 


\tiny
\begin{table}[h!]
\begin{center}
%\begin{tabular}{|c|c|c|}\hline
\begin{tabular}{|l|p{3.6cm}|p{9cm}|}\hline
&& \\
Sub-system&State variable & Model \\ && \\ \hline && \\
{\it Stratosphere} & Average temperature & Capacitance with fixed and uniform lapse rate  \\ && \\
&  ($T_{\sigma}$) ($K$)  & $ C_{\sigma} \cdot \frac{dT_{\sigma}}{dt} = \psi_{S \sigma} + \psi_{L \sigma}  $ \\ && \\ \hline
&&\\
{\it Troposphere} & Average temperature & Capacitance with fixed and uniform lapse rate \\&&\\
& at 5000~$m$ ($T_{\tau}$) ($K$) &
$ C_{\tau} \cdot \frac{dT_{\tau}}{dt} = \psi_{S \tau} + \psi_{L \tau} + \phi_{s \tau} + \phi_{l \tau}  $  \\&&\\\hline
&&\\
{\it Ocean mixed layer} & Temperature & Capacitance with uniform temperature \\&&\\
& ($\theta$) ($K$) &
$C_{\theta} \cdot \frac{d\theta}{dt} = \psi_{S \theta} + \psi_{L \theta} - \phi_{s \tau} - \phi_{l \theta}  $\\&&\\\hline
&&\\
{\it Water vapor} & Water vapor & Capacitance with fixed relative humidity gradient ($\gamma _{r})$ \\&&\\
& content (Q) ($kg\;m^{-2}$)&
$ \frac{dQ}{dt} = \frac{1}{L_v} \cdot ( \phi_{l \theta} - \phi_{l \tau} ) $ \\&&\\\hline
&&\\
{\it Carbon dioxide} & CO$_2$ concentration & Fixed \\&&\\
& (CO$_2$) (ppmv)& \\&&\\ \hline
\end{tabular}
\caption{Model sub-systems with state variables and model types}
\label{tab_cell}
\end{center}
\end{table}
%
\begin{table}[h]
\begin{center}
\begin{tabular}{|l|p{3cm}|p{10cm}|}\hline &&\\
Transfer & Transfer      & Model \\ &&\\
         &  variable(s)  &       \\&&\\\hline&&\\

{\it Solar flux} & SW budget of & constant solar flux shared between ocean, \\ &&\\
&  troposphere, & troposphere, stratosphere and reflected flux to space \\&&\\
&  stratosphere, & $ \psi_{S \sigma}= \alpha^{\sigma}_{S} \cdot F_S $ \\&&\\
&  ocean & $ \psi_{S \tau}=  N \cdot \alpha^{C}_{S} (1 - \alpha^{\sigma}_{S}) F_S $ \\&&\\
&  ($\psi_{S \tau}$, $\psi_{S \sigma}$, $\psi_{S \theta}$) & $ \psi_{S \theta}= (1-A) \left( (1-N) + N \cdot \mathcal{T}^{C}_{S} \right)
(1 - \alpha^{\sigma}_{S}) F_S$ \\&&\\
&  ($W\;m^{-2}$) & $N=1-(1-N_l) (1-N_m) (1-N_h)$ \\&&\\ \hline&&\\



{\it LW flux}    & LW budget of & Malkmus narrow-band radiative module \\&&\\
& troposphere, &  \\&&\\
& stratosphere, & \\&&\\
& ocean  & \\&&\\
& ($\psi_{L \tau} $, $\psi_{L \sigma} $, $ \psi_{L \theta} $) & \\&&\\
& ($W\;m^{-2}$) & \\&&\\ \hline&&\\



{\it Sensible heat flux} & Surface heat flux & Diffusive term proportional to surface temperature \\&&\\
 & ($\phi_{s \tau}$) & difference \\&&\\
 & ($W\;m^{-2}$) & $ \phi_{s \tau} = h_b \cdot ( T_{\theta} - T_{b\tau} )$ \\&&\\
 & & $ T_{b\tau} = T_{\tau} - \gamma_{\tau} \cdot H_m $ \\ &&\\\hline

\end{tabular}
\caption{Transfer models, with transfer variables and model types (Part I)}
\label{tab_tran1}
\end{center}
\end{table}
%
\begin{table}[h!]
\begin{center}
\begin{tabular}{|l|p{3cm}|p{10cm}|}\hline &&\\
Transfer & Transfer      & Model \\ &&\\
         &  variable(s)  &       \\&&\\\hline&&\\

{\it Precipitation} & Surface & Newtonian term nudging tropospheric relative humidity \\&&\\
 & precipitation and & toward target relative humidity \\&&\\
 & corresponding & $\phi_{l \tau} = -\phi_{l \theta} + L_v \cdot \frac{1}{\tau_p} \cdot (Q - Q^{\infty}(T_{\tau})) $ \\&&\\
 & latent heat flux & \\&&\\
 & ($\phi_{l \tau}$) & $Q^{\infty}(T_{\tau})$ is the water vapor content corresponding \\&&\\
 & ($W\;m^{-2}$)  & to the target relative humidity. \\&&\\ \hline&&\\


{\it Evaporation} & Surface latent & Diffusive term proportional to surface humidity difference \\&&\\
 & heat flux & $\phi_{l \theta} = k_e \cdot L_v \cdot ( P_{sat} (\theta) - r_H \cdot P_{sat} (T_{surf}) ) $ \\&&\\
 & ($\phi_{l \theta}$) & \\&&\\
 & ($W\;m^{-2}$) & $P_{sat}(T)$ is the water vapor saturation pressure. \\&&\\\hline&&\\


{\it WV temperature} & $T_{WV}$ ($K$) & $T_{WV} = T_{\tau}$ \\&&\\ \hline&&\\
{\it Cloud cover} & High, medium& fixed:\\&&\\
 & and low level & $N_h=N_{h0}$, \\&&\\
  & cloud cover & $N_m=N_{m0}$ \\&&\\
 & ($N_h$, $N_m$, $N_l$) & $N_l=N_{l0}$ \\&&\\
 &  (no unit)& \\&&\\ \hline
%{\it CC temperature} & $T_{CC}$ ($K$) & $T_{CC} = T_{\tau}$ \\ &&\\\hline
\end{tabular}
\caption{Transfer models, with transfer variables and model types (Part II)}
\label{tab_tran2}
\end{center}
\end{table}
%**Al1 rev III et II demande de justifier les valeurs numeriques

\begin{table}[h]
\begin{center}
\begin{tabular}{|c|c|}\hline
Variable & Value \\\hline
Surface relative humidity at equilibrium   & 0.7    \\&\\
Top of the atmosphere altitude & 65,000~$m$ \\&\\
Tropopause altitude   & 15,000~$m$    \\&\\
Mean troposphere altitude ($H_m$) & 5000~$m$ \\&\\
Troposphere lapse rate ($\gamma_{\tau}$)   & $-6.5 \cdot 10^{-3} \; K\;m^{-1}$    \\&\\
Troposphere heat capacity ($ C_{\tau}$) & $1.\cdot 10^7 \; J\;K^{-1}\;m^{-2}$\\&\\
Stratosphere lapse rate ($\gamma_{\sigma}$)  & $3.5\cdot 10^{-3} \; K\;m^{-1}$    \\&\\
Stratosphere heat capacity ($ C_{\sigma}$) & $1.2\cdot 10^6 \; J\;K^{-1}\;m^{-2}$ \\&\\
Stratosphere SW absorption ($ \alpha^{\sigma}_{S}$) & 12\% \\&\\ \hline&\\
High level cloud cover ($N_{h0}$)  & 22\%    \\&\\
High level cloud altitude   & 14,000~$m$    \\&\\
Midlevel cloud cover ($N_{m0}$)  &  18\%   \\&\\
Midlevel cloud altitude   & 6,000~$m$    \\&\\
Low level cloud cover ($N_{l0}$)   & 26\%   \\&\\
Low level cloud altitude   &  2,000~$m$   \\&\\
Cloud LW absorption ($\alpha^C_{L}$)  &  100\%   \\&\\
Cloud SW transmitivity ($\mathcal{T}^{C}_{S}$) &  20\%   \\&\\
Cloud SW absorption ($\alpha^{C}_{S}$) &  10\%   \\&\\ \hline
\end{tabular}
\caption{Atmosphere parameter values. Sources on atmosphere parameters and cloud properties: \citet{Salby96};
Source on cloud covers: http://isccp.giss.nasa.gov. }
\label{para1}
\end{center}
\end{table}
%

\begin{table}[h!]
\begin{center}
\begin{tabular}{|c|c|}\hline&\\
Variable & Value \\&\\\hline&\\
Ocean albedo ($A$)  & 0.05    \\&\\
Ocean heat capacity ($ C_{\theta}$) & $2.\cdot 10^8 \; J\;K^{-1}\;m^{-2}$ \\&\\
Oceanic mixed layer thickness   & 50~$m$    \\&\\ \hline&\\
Atmosphere-Ocean sensible heat transfer coefficient ($h_b$)   &  $5 \; W\;m^{-2}\;K^{-1}$  \\&\\
Evaporation coefficient ($k_e$)  & $0.35 \cdot 10^{-2} \; s\;m^{-1} $    \\&\\
Precipitation characteristic time ($\tau_p$)  & 5 days      \\  &\\
Latent heat of evaporation ($L_v$) & $2.5 \cdot 10^6 \, J\;kg^{-1} $ \\&\\
Incoming solar flux ($F_S$) & $340~W\;m^{-2}$ \\&\\ \hline

\end{tabular}
\caption{Surface parameter values. Data from \citet{Salby96}, except: the atmosphere-Ocean sensible heat
transfer coefficient ($h_b$) and the evaporation coefficient ($k_e$), that are calibrated in
order to have realistic fluxes; and the precipitation characteristic time ($\tau_p$), that is
\textit{ad hoc}, equal to the mean synoptic time scale.}
\label{para2}
\end{center}
\end{table}
\normalsize


%\begin{figure}
%\begin{center}
%\includegraphics[height=13cm]{climsi.eps}
%\caption{Model Structure }
%\label{Struct_ClimSI}
%\end{center}
%\end{figure}



\clearpage
%
\newpage

\bibliographystyle{ametsoc}
%\bibliographystyle{jas}
%\input{BiblioD}
\bibliography{biblio}
%
\end{document}
%\bibliography{/home/al1/Cahiers/Gains/Articles/biblio}
\begin{thebibliography}{16}
\expandafter\ifx\csname natexlab\endcsname\relax\def\natexlab#1{#1}\fi
\expandafter\ifx\csname bibnamefont\endcsname\relax
  \def\bibnamefont#1{#1}\fi
\expandafter\ifx\csname bibfnamefont\endcsname\relax
  \def\bibfnamefont#1{#1}\fi
\expandafter\ifx\csname citenamefont\endcsname\relax
  \def\citenamefont#1{#1}\fi
\expandafter\ifx\csname url\endcsname\relax
  \def\url#1{\texttt{#1}}\fi
\expandafter\ifx\csname urlprefix\endcsname\relax\def\urlprefix{URL }\fi
\providecommand{\bibinfo}[2]{#2}
\providecommand{\eprint}[2][]{\url{#2}}

\bibitem[{\citenamefont{Bode}(1945)}]{Bode45}
\bibinfo{author}{\bibfnamefont{H.}~\bibnamefont{Bode}},
  \emph{\bibinfo{title}{Network Analysis and Feedback Amplifier Design}}
  (\bibinfo{publisher}{Van Nostrand, New York}, \bibinfo{year}{1945}).

\bibitem[{\citenamefont{Manabe and Wetherald}(1967)}]{manabe67}
\bibinfo{author}{\bibfnamefont{S.}~\bibnamefont{Manabe}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{R.}~\bibnamefont{Wetherald}},
  \bibinfo{journal}{J. Atmos. Sci.} \textbf{\bibinfo{volume}{24}},
  \bibinfo{pages}{241} (\bibinfo{year}{1967}).

\bibitem[{\citenamefont{Bony et~al.}(2006)\citenamefont{Bony, Colman, Kattsov,
  Allan, Bretherton, Dufresne, Hall, Hallegatte, Holland, Ingram
  et~al.}}]{Bony06}
\bibinfo{author}{\bibfnamefont{S.}~\bibnamefont{Bony}},
  \bibnamefont{et~al.}, \bibinfo{journal}{Journal of Climate}
  \textbf{\bibinfo{volume}{19}}, \bibinfo{pages}{3445} (\bibinfo{year}{2006}).

\bibitem[{\citenamefont{Cess et~al.}(1989)\citenamefont{Cess, Potter, Blanchet,
  Boer, Ghan, Kiehl, Le, Li, Liang, Mitchell et~al.}}]{Cess89}
\bibinfo{author}{\bibfnamefont{R.}~\bibnamefont{Cess}},
  \bibnamefont{et~al.}, \bibinfo{journal}{Science}
  \textbf{\bibinfo{volume}{245}}, \bibinfo{pages}{513} (\bibinfo{year}{1989}).

\bibitem[{\citenamefont{Hansen et~al.}(1984)\citenamefont{Hansen, Lacis, Rind,
  Rusell, Stone, Fung, Ruedy, and Lerner}}]{Hansen84}
\bibinfo{author}{\bibfnamefont{J.}~\bibnamefont{Hansen}},
 \bibnamefont{et~al.}, in
  \emph{\bibinfo{booktitle}{Climate Processes and Climate Sensitivity}}, edited
  by \bibinfo{editor}{\bibfnamefont{J.}~\bibnamefont{Hansen}} \bibnamefont{and}
  \bibinfo{editor}{\bibfnamefont{T.}~\bibnamefont{Takahashi}}
  (\bibinfo{publisher}{Geophysical monograph {\bf 29}, AGU},
  \bibinfo{year}{1984}), vol.~\bibinfo{volume}{5}, pp.
  \bibinfo{pages}{130--162}.

\bibitem[{\citenamefont{Hallegatte et~al.}(2006)\citenamefont{Hallegatte,
  Lahellec, and Grandpeix}}]{HLG05}
\bibinfo{author}{\bibfnamefont{S.}~\bibnamefont{Hallegatte}},
%  \bibinfo{author}{\bibfnamefont{A.}~\bibnamefont{Lahellec}}, \bibnamefont{and}
%  \bibinfo{author}{\bibfnamefont{J.-Y.} \bibnamefont{Grandpeix}},
\bibnamefont{et~al.},
  \bibinfo{journal}{J. Atmos. Sci.} \textbf{\bibinfo{volume}{63}},
  \bibinfo{pages}{1878} (\bibinfo{year}{2006}).


\bibitem[{\citenamefont{Luenberger}(1979)}]{Luenberger}
\bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Luenberger}},
  \emph{\bibinfo{title}{Introduction to Dynamic Systems}}
  (\bibinfo{publisher}{John Wiley \& Sons}, \bibinfo{year}{1979}).


\end{thebibliography}
\large

\bigskip
\end{document}

