source: LMDZ.3.3/trunk/libf/dyn3d/LISEZMOI

 ..  Version  du 28/12/1998  , avec  nouveau fxyhyper appele par inigeom ...
c
c
 
   La dynamique  de LMDZ.3.2  est pratiquement identique  a celle de 
                    LMDZ.3   .

  Elle ne differe , en gros , que par   :

 a ) etat0_lmd   ( on doit introduire , en plus de dayref qui est le jour
     -----------
      de l'etat initial cree ,  anneeref  , qui est l'annee correspondante
      et qui servira eventuellement  a lire la bonne annee pour le fichier
      limit.dat , dic  ) . Il y a aussi le  stockage des coefficients pour
      l'orographie  de Francois Lott , sur startphy .  

 b ) calfis  qui prepare l'interface avec la physique  qui a ete restruc-
     -------
     turee  par Laurent Li , dans ses grandes lignes , et qui est maintenant

     adaptee pour les traceurs avec l'ajout du shema de convection de Tiedtke. 
   
     Pour plus de details concernant le cote pratique des simulations avec
     cette bib.  , priere de se reporter a la fin  de ce  LISEZMOI .


 **************************************************************************

 1) Le gcm LMDZ.3.2  a des coordonnees  hybrides  verticalement ,c.a.d  des
          ********
   coordonnees equivalentes a sigma pres de la surface et a  P plus haut .     

   Les pressions  p(ij,l) aux differents interfaces des  llm  couches verti- 
   cales sont donnees  par   :

    p( ij,l ) = ap(l) + bp(l) * ps(ij)    ,    1 << l << llm + 1    . 

                  avec   bp(l) =  EXP ( 1. - 1./ sigma(l) ** 2 )
                         ap(l) =  pa * (  sigma(l) -  bp(l)  )

                    pa   est actuellement mis a  50 000.  Pascals .
            
   La fonction d'Exner  est calculee aux milieux des llm couches vertic.

     pk(ij,l) = cp * ( p(ij,l)/preff ) ** kappa   ,  1 << l <<  llm

   Les pressions  PLS(ij,l)  definies aux milieux des couches sont egales a :

     pls(ij,l) = preff * ( pk(ij,l) /cp ) ** 1./kappa

   (  preff = 101 325. Pascals  dans  etat0_lmd.F  actuellement  )
   

 2) Une des variables prognostiques de LMDZ.3  est   teta  , temperature 
     potentielle  :    

                 teta (ij,l) = temperature(ij,l) * cp/ pk(ij,l)
 
N.b : Les parametres  PA et PREFF  figurent dans l'entete des fichiers
***   start ( et restart) et histphy ,  ainsi que  AP(l) et BP(l)  .

    Y figurent aussi  PRESNIVS(l) qui donnent les valeurs approximatives
    des pressions aux milieux des couches , ainsi que  SIG_S(l) qui donnent
    simplement les indices des llm couches et dont  la presence  permet 
    d'utiliser tel quel  le logiciel ' wrctl '  de LMDZ.2.2 qui permet de  
    fabriquer a partir du fichier ' histphy' ( dat et dic ) le  histphy.ctl  
    pret a etre exploite  par  GrADS .

   Les pressions PLS(ij,l)  sont stockees sur le fichier ' histphy ' car
   elles sont indispensables quand on veut faire une interpolation des varia-
   bles d'analyse comme  u, v ,teta  a des valeurs determinees de la pression.
   Elles figurent aussi eventuellement sur le fichier histins.gdat (.ctl ) 
   qui renferme des etats instantanes .
   Voir le parametre ok_inst dans physiq.f  ( L.Li ) .

   Un exemple d'interpolation de u,v, et teta  se trouve  dans le fichier 
   'analyse_lmdz_hyb' du repertoire de la dynamique , analyse_hyb etant 
   adapte d'un fichier d'analyse de histphy  de  Laurent Li .
    
 3) Dans  le programme principal gcm.F , il faut definir , au debut ,les para-
    metres  IADV ,pour choisir le shema d'advection pour l'eau et les traceurs
               :    IADV = 1    , humidite specifique
                    IADV = 2    , shema    amont
                    IADV = 3    , Van-Leer

 4) Il y a 2 shemas de dissipation   ( diffusion laterale )  :

 a) La dissipation star ( qui ne fait pas intervenir la geometrie de la     
    grille et qui est utilisee quand il y a un zoom ) 

    Par expl . la diffusion d'un scalaire  comme  teta est formulee ( en 
     appelant l'operateur  DIVGRAD2 )   comme :

    d( ro * teta ) = (-1)**(it-1) *cdiv * Lap * ro**1/2 * Lap_gam** (it-2 ) *
                           ro**1/2 * Lap * teta
 avec  :
 Lap     = 1/aire * Filt[ dx ( coefx*dx ) + dy( coefy*dy ) ] Filt
       
 Lap_gam = (1./aire)** -gam { dx[ coefx**(-gam) dx] + dy[ coefy**(-gam) dy ] } 


  Filt est l'operateur  filtre 
  ro   est la densite dp  de chaque couche
  it   = nombre d'iterations de l'operateur de dissipation ( ici divgrad2 )
  gam  = coefdis /( it -2 ) , coefdis = 0. si it = 2 , sinon = 1. habituellem.
  dx   = operateur difference finie :
    dx(A) = A( X +1/2 , Y, z ) - A( X -1/2, Y, z ) , dans le systeme de 
  travail  (X,Y) , avec l'increment DX = 1.

  Idem pour dy .

  coefx = beta(ij) / alpha(ij)  , coefy = alpha(ij)/beta(ij)
     avec    u_covariant(ij) = alpha(ij)* u_nat(ij) 
             v_covariant(ij) = beta(ij) * v_nat(ij) 

 Cette dissipation star est choisie en mettant le parametre lstardis a .TRUE.
  dans le fichier run.def  lu  par le gcm en entree .

 b) La dissipation ' non star '  fait intervenir les operateurs elementaires
   diverg,grad,rotat .   Elle est choisie  avec lstardis = .FALSE.

N.B : Dans les 2 shemas de dissipation , le gcm est plus stable  quand on
     n'itere pas l'operateur ' GRADIV ' ou ' GRADIV2 " .            
     
 5) Les constantes thermodynamiques Cp et kappa = R/Cp   utilisees dans la
   partie dynamique ont ete mises , ce qui est normal, aux memes valeurs
   qu'elles ont dans la partie ' physique ' , en changeant INICONST et en 
   les mettant dans l'entete de controle de start et startphy ,en changeant
   etat0_lmd , par rapport  a LMDZ.2.2 
       
 6) Pour finir , quelques commentaires sur les parametres qui figurent sur
    le fichier  run.def qui commandent le gcm .

  Dans le run.def donne en exemple dans /users/SOURCES/LMDZ.3/libf/dyn3d  ,

 a) day_step  correspond a un  dt  de  2 minutes .

 b) idissip correspond a a un appel de la dissipation toutes les 30 minutes.

N.B : Maintenant la routine inidissip calcule la meilleure valeur pour
  idissip , en fonction des valeurs de tetagdiv , tetagrot et tetatemp qui
  representent les temps de dissipation desirees pour dissiper les plus peti-
  tes longueurs d'ondes pour u ,v et teta  dans le cas, respectivement des 
  operateurs gradiv (gradiv2) , nxgrarot (nxgraro2) et divgrad ( divgrad2) .

 c) coefdis , coefficient qui intervient dans le calcul de gam  , utilise 
   dans les routines de dissipation star , vus plus haut .

 d) iphysiq correspond ici a un appel de la physique toutes les 30 minutes.

 e) alphax et alphay correspondent au grossissement pour le zoom dans la
   direction X et Y .
   En fait , les facteurs de grossissement sont egaux a  1./ ( 1.- alphax )
     et 1./( 1.- alphay )  .

 f) fxyhypb  est egal a .TRUE. si on choisit  pour la fonction f(x,y) qui 
   definit la grille du modele  une fonction a tangente hyperbolique , un peu
   adaptee pour le zoom , que la fonction sinusoidale , definie par fxy_new.h 
   dans /..../gcm/libf/grid  et lu par  makegcm  lors de la compilation du 
   code .


N.B: Pour l'instant , cette fonction hyperbolique n'est activee eventuellement
****            que pour  un zoom  latitudinal . 
   Pour un zoom longitudinal , la fonction sinusoidale precisee dans fxy_new.h
       est  choisie automatiquement .

 g) psi0 et delta sont des parametres  pour cette fonction hyberbolique .
   Ces valeurs de 0.3 et 0.5  ont ete adaptees pour 96x72  .
  


 7) La partie 'physique'  est quasiment  identique a celle de  LMDZ.2.2 
             ***********
    ( /users/lmdz/SOURCES/LMDZ.2.2/libf/phylmd ) geree par  Laurent Li

    sauf  PHYSIQ.F ( ou il y a l'ecriture en plus de pls  sur 'histphy' et 
    histins  ) et   PHYETAT0.F et PHYREDEM.F  ou il y a l'argument  RLON
    en plus  dans le tableau des aurguments  et WRITEFIA  ou nivsigs rempla-
    sig_s  pour les niveaux verticaux .   


 ****************************************************************************

 
 AA) Choisir les options du run  , avec   le fichier  run.def   .

     Il y a entre autres choses , le nombre de jours qu'on veut integrer lors 
       d'un run .
     Il y a le nombre de pas a faire pour 1 jour , ce qui calcule donc le dt
       temporel ( en secondes ) pour un pas dynamique .
     Iperiod  est conseille  d'etre a 5 pas ( mais non obligatoiremnet = 5 ) . 
     Idissip  peut etre mis a n'importe quelle valeur . Il est maintenant
     calcule dans inidissip de telle facon a optimiser  sa valeur .

     Les parametres( la plupart des cles ) pour la physique commencent avec
      la var. logique         **  physic  **

           entiere( en pas )  ** iphysiq  ** qui est calculee de telle facon a
                              appeler la physique tous les  30 minutes .

           ( en jours )       ** ecritphy ** qui indique la frequence en jours
                                du stockage   sur   histphy .  

      les variables logiques  ** cycle_diurne **  pour le cycle diurne

                              ** soil_model  **   pour un modele pour le sol

                              ** new_oliq   **  pour une nouvelle formulation
            pour l'eau liquide nuageuse , utilisee dans les rayonnements .

                              ** ok_orodr **   si on appelle l'orographie
                              ** ok_orolf **        IDEM

                              ** ok_limitvrai ** si on veut lire la bonne
                          annee  du fichier limit .


 Si ok_limitvrai = .TRUE.  , la routine  condsurf  de la physique  va lire

 la bonne annee  des conditions limites et il faut , dans le script qui lance

 le run , faire figurer  les fichiers limites attendus par le gcm .


 PLus clairement , si on lance un run  de 30 jours ,  du 25 Decembre 1986 au

 25 Janvier  1987  ,  ( actuellement dans le Gcm , il n'y a que des mois de

 30 jours ) ,   des que le  modele  arrive au jour 360 (qui correspond au

 30 Decembre) pour aller au jour 361  ( qui correspond au 1 Janvier 1987 ) ,

 il s'attend a lire le fichier  limit1987 ( avec 4 chiffres  pour l'an 2000).

 Donc , dans le script du run , on devra avoir, grosso modo  :

***********************************************************************
cd $TMPDIR
#
cp $workdir/mod/limit86.dat limit1986.dat
cp $workdir/mod/limit86.dic limit1986.dic
cp $workdir/mod/limit87.dat limit1987.dat
cp $workdir/mod/limit87.dic limit1987.dic
#
gcm.e < run.def
**********************************************************************

      la variable entiere   ** nbapp_rad ** qui est le nombre d'appels /jour
                  des routines de rayonnements .

      la variable entiere   ** iflag_con ** qui indique le shema de convection
        choisi  :     1   >  LMD ,   2  > Tiedtke  ,   3  >  CCM  .
    
                       

 BB) Creer l'etat initial , avec les donnees de l'ECMWF , en compilant 
     etat0_lmd  , en precisant auparavant  dayref  et anneeref .


 CC) Creer le fichier des conditions aux limites ( SST, Albedo, Rugosite,
     Nature du sol  ) avec  si possible la concordance entre l'annee des
     donnees des sst, albedo, etc ..  avec  les donnees pour l'etat initial .

 DD)  Compiler  le Gcm , lancer le run , prendre un bon verre et attendre
      la sortie  des resultats .



                 Phu  Le Van  .

                 E.N.S , 24 Rue Lhomond  , PARIS  5 ieme .
                 Tel :  01 44 32 22 67

                 4 Place Jussieu         , PARIS  5 ieme .
                 TEL :  01 44 27 47 72

                 Email :  levan@lmd.ens.fr