Changeset 6


Ignore:
Timestamp:
Oct 22, 2010, 12:47:12 PM (14 years ago)
Author:
slebonnois
Message:

cf commit_v6.log :

  • manipulation traceurs
  • homogeneisation .def
  • bilan_dyn
  • etats initiaux start.nc
  • appels specifiques pour physique
Location:
trunk
Files:
1 added
17 edited

Legend:

Unmodified
Added
Removed
  • trunk/deftanks/titan/gcm.def

    r5 r6  
    11## $Header$
     2#
     3## Planete:
     4planet_type=titan
     5#
    26## nombre de pas par jour (multiple de iperiod) ( ici pour  dt ~ 2 min )     
    37day_step=12000
    4 iapp_tracvl=5
    58## periode pour le pas Matsuno (en pas)
    69iperiod=5
     
    2528## choix du shema d'integration temporelle (Matsuno ou Matsuno-leapfrog)
    2629purmats=n
    27 ## avec ou sans physique                                                 
     30## avec ou sans physique
     31## 0: pas de physique (e.g. en mode Shallow Water)
     32## 1: avec physique (e.g. physique phylmd)
     33## 2: avec rappel newtonien dans la dynamique                                         
    2834iflag_phys=1
     35## avec ou sans fichiers de demarrage (start.nc, startphy.nc) ?
     36## (sans fichiers de demarrage, initialisation des champs par iniacademic
     37##  dans la dynamique) PAS AU POINT POUR TITAN
     38read_start=y
    2939## periode de la physique (en pas)                                       
    3040iphysiq=20
    3141## avec ou sans traceurs                                                 
    3242iflag_trac=1
    33 ##  Cycle diurne  ou non                 
    34 cycle_diurne=y
    35 ##  Soil Model  ou non               
    36 soil_model=y
    37 ##  Orodr  ou  non   pour l orographie             
    38 ok_orodr=n
    39 ##  Orolf  ou  non   pour l orographie             
    40 ok_orolf=n
    41 ## Nombre  d'appels des routines de rayonnements ( par jour)                 
    42 nbapp_rad=10
    43 ## Nombre  d'appels des routines de chimie ( par jour)                 
    44 nbapp_chim=1
     43#
    4544## longitude en degres du centre du zoom                                 
    4645clon=0.
  • trunk/deftanks/titan/physiq.def

    r5 r6  
     1## $Header: /home/cvsroot/LMDZ4/physiq.def,v 1.2 2004/06/22 11:45:18 lmdzadmin Exp $
    12#
    2 # $Header: /home/cvsroot/LMDZ4/physiq.def,v 1.2 2004/06/22 11:45:18 lmdzadmin Exp $
     3# PARAMETRES ANCIENNEMENT DANS gcm.def
     4##  Cycle diurne  ou non                 
     5cycle_diurne=y
     6##  Soil Model  ou non               
     7soil_model=y
     8##  Orodr  ou  non   pour l orographie             
     9ok_orodr=n
     10##  Orolf  ou  non   pour l orographie             
     11ok_orolf=n
     12##  Gravity Waves non-orographiques
     13ok_gw_nonoro=n
     14## Nombre  d'appels des routines de rayonnements ( par jour)                 
     15nbapp_rad=10
     16## Nombre  d'appels des routines de chimie ( par jour)                 
     17nbapp_chim=1
     18##  Flag  pour la convection : 1 pour LMD, 2 pour Tiedtke, 3 KE(nvlle version JYG), 30 KE(version IPCC AR4), 4 KE vect
     19iflag_con=0
    320#
    4 
     21#
     22# Parametres fichiers de sortie
     23#
     24### OK_mensuel= y sortir fichier mensuel histmth.nc, =n pas de fichier histmth.nc
    525OK_mensuel=y
     26### OK_journe= y sortir fichier journalier histday.nc, =n pas de fichier histday.nc
    627OK_journe=y
     28### OK_instan=y, ecrire sorties "instantanees" (chaque pas de temps de la  physique)
    729OK_instan=n
    8 if_ebil=0
    9 #
    10 # parametres outputs
     30#
     31# Parametres niveau de sorties differents fichiers
    1132#
    1233#niveau de sortie "mth" lev_histmth
  • trunk/deftanks/titan/run.def

    r5 r6  
    88##  Annee de l'etat  initial (   avec  4  chiffres   )
    99anneeref=1979
     10## Remise a zero de la date initiale
     11raz_date=0
    1012## Nombre de jours d'integration
    1113nday=11
    1214## periode de sortie des variables de controle (en pas)
    1315iconser=500
    14 ## periode d'ecriture du fichier histoire (en jour) -OLD-
    15 iecri=1000
    16 ## periode de stockage fichier dyn_hist_ave et dynzon (en jour)
    17 periodav=0.0016666666667
     16#
     17## sorties instantanees dans la dynamique (fichiers dyn_hist.nc and co.)
     18ok_dyn_ins=n
     19## periode d'ecriture des sorties instantanees dans la dynamique
     20## (en pas dynamiques)
     21iecri= A REPRENDRE !!
     22## sorties de valeurs moyennes dans la dynamique (fichiers dyn_hist_ave.nc and co.)
     23ok_dyn_ave=n
     24## periode de stockage des moyennes dans la dynamique et dans dynzon (en jour)
     25periodav= 10.
     26## flag de sortie dynzon
     27ok_dynzon=n
  • trunk/deftanks/venus/gcm.def

    r5 r6  
    11## $Header$
     2#
     3## Planete:
     4planet_type=venus
     5#
    26## nombre de pas par jour (multiple de iperiod) ( ici pour  dt = 1 min )     
    37day_step=120000
    4 iapp_tracvl=5
    58## periode pour le pas Matsuno (en pas)
    69iperiod=5
     
    2629purmats=n
    2730## avec ou sans physique                                                 
     31## 0: pas de physique (e.g. en mode Shallow Water)
     32## 1: avec physique (e.g. physique phylmd)
     33## 2: avec rappel newtonien dans la dynamique                                         
    2834iflag_phys=1
     35## avec ou sans fichiers de demarrage (start.nc, startphy.nc) ?
     36## (sans fichiers de demarrage, initialisation des champs par iniacademic
     37##  dans la dynamique) PAS AU POINT POUR VENUS
     38read_start=y
    2939## periode de la physique (en pas)                                       
    3040iphysiq=20
    31 ## frequence (en  jours ) de l'ecriture du fichier histmth               
    32 # regle la frequence d'ecriture dans histmth (mensuel)
    33 ecritphy=0.1
    3441## avec ou sans traceurs                                                 
    3542iflag_trac=0
    36 ##  Cycle diurne  ou non                 
    37 cycle_diurne=y
    38 ##  Soil Model  ou non               
    39 soil_model=y
    40 ##  Choix ou non  de  New oliq               
    41 new_oliq=y
    42 ##  Orodr  ou  non   pour l orographie             
    43 ok_orodr=n
    44 ##  Orolf  ou  non   pour l orographie             
    45 ok_orolf=n
    46 ##   Si = .T. ,  lecture du fichier limit avec la bonne annee             
    47 ok_limitvrai=n
    48 ## Nombre  d'appels des routines de rayonnements ( par jour)                 
    49 nbapp_rad=2000
    50 ##  Flag  pour la convection (1 pour LMD, 2 pour Tiedtke, 3 KE, 4 KE vect)
    51 iflag_con=0
     43#
    5244## longitude en degres du centre du zoom                                 
    5345clon=0.
  • trunk/deftanks/venus/physiq.def

    r5 r6  
     1## $Header: /home/cvsroot/LMDZ4/physiq.def,v 1.2 2004/06/22 11:45:18 lmdzadmin Exp $
    12#
    2 # $Header: /home/cvsroot/LMDZ4/physiq.def,v 1.2 2004/06/22 11:45:18 lmdzadmin Exp $
     3# PARAMETRES ANCIENNEMENT DANS gcm.def
     4##  Cycle diurne  ou non                 
     5cycle_diurne=y
     6##  Soil Model  ou non               
     7soil_model=y
     8##  Orodr  ou  non   pour l orographie             
     9ok_orodr=n
     10##  Orolf  ou  non   pour l orographie             
     11ok_orolf=n
     12##  Gravity Waves non-orographiques
     13ok_gw_nonoro=n
     14## Nombre  d'appels des routines de rayonnements ( par jour)                 
     15nbapp_rad=2000
     16## Nombre  d'appels des routines de chimie ( par jour)                 
     17nbapp_chim=1
     18##  Flag  pour la convection : 1 pour LMD, 2 pour Tiedtke, 3 KE(nvlle version JYG), 30 KE(version IPCC AR4), 4 KE vect
     19iflag_con=0
    320#
    4 
     21#
     22# Parametres fichiers de sortie
     23#
     24### OK_journe= y sortir fichier journalier histday.nc, =n pas de fichier histday.nc
     25### Pour Venus, ca n'a pas de sens
    526OK_journe=n
     27### OK_mensuel= y sortir fichier mensuel histmth.nc, =n pas de fichier histmth.nc
     28### Pour Venus, seul fichier de sorties moyennees
    629OK_mensuel=y
     30## frequence (en  jours ) de l'ecriture du fichier histmth pour Venus               
     31# regle la frequence d'ecriture dans histmth (mensuel)
     32ecritphy=0.1
     33### OK_instan=y, ecrire sorties "instantanees" (chaque pas de temps de la  physique)
    734OK_instan=n
    8 if_ebil=0
    935#
    10 # parametres outputs
     36# parametres niveau de sorties differents fichiers
    1137#
    12 #niveau de sortie "hf" lev_histhf avec
    13 # - lev_histhf=3 => defaut
    14 # - lev_histhf=4 => histhf3d.nc champs 3d niveaux modele
    15 lev_histhf=0
    1638#niveau de sortie "day" lev_histday
    1739# - lev_histday=1 => basiques, 2D
  • trunk/deftanks/venus/run.def

    r5 r6  
    44INCLUDEDEF=physiq.def
    55INCLUDEDEF=gcm.def
    6 INCLUDEDEF=orchidee.def
    76## Jour de l'etat initial ( = 350  si 20 Decembre ,par expl. ,comme ici )
    87dayref=1
    98##  Annee de l'etat  initial (   avec  4  chiffres   )
    109anneeref=1979
     10## Remise a zero de la date initiale
     11raz_date=0
    1112## Nombre de jours d'integration
    1213nday=2
    1314## periode de sortie des variables de controle (en pas)
    1415iconser=500
    15 ## periode d'ecriture du fichier histoire (en jour)
    16 iecri=1000
    17 ## periode de stockage fichier dyn_hist_ave et dynzon (en jour)
    18 periodav=5.
     16#
     17## sorties instantanees dans la dynamique (fichiers dyn_hist.nc and co.)
     18ok_dyn_ins=n
     19## periode d'ecriture des sorties instantanees dans la dynamique
     20## (en pas dynamiques)
     21iecri= A REPRENDRE !!
     22## sorties de valeurs moyennes dans la dynamique (fichiers dyn_hist_ave.nc and co.)
     23ok_dyn_ave=n
     24## periode de stockage des moyennes dans la dynamique et dans dynzon (en jour)
     25periodav= 2.
     26## flag de sortie dynzon
     27ok_dynzon=n
  • trunk/libf/dyn3d/addfi.F

    r1 r6  
    2727c      pdufi(ip1jmp1,llm)     |
    2828c      pdvfi(ip1jm,llm)       |   respective
    29 c      pdhfi(ip1jmp1)         |      tendencies
     29c      pdhfi(ip1jmp1)         |      tendencies  (unit/s)
    3030c      pdtsfi(ip1jmp1)        |
    3131c
     
    116116      ENDDO
    117117 
    118       DO iq = 1, 2
     118      DO iq = 1, nqtot
     119       IF ((planet_type.eq.'earth').and.(iq.le.2)) THEN
    119120         DO k = 1,llm
    120121            DO j = 1,ip1jmp1
     
    123124            ENDDO
    124125         ENDDO
    125       ENDDO
    126 
    127       DO iq = 3, nqtot
     126       ELSE
    128127         DO k = 1,llm
    129128            DO j = 1,ip1jmp1
     
    133132         ENDDO
    134133      ENDDO
    135 
    136134
    137135      DO  ij   = 1, iim
  • trunk/libf/dyn3d/bilan_dyn.F

    r1 r6  
    22! $Id: bilan_dyn.F 1403 2010-07-01 09:02:53Z fairhead $
    33!
    4       SUBROUTINE bilan_dyn (ntrac,dt_app,dt_cum,
    5      s  ps,masse,pk,flux_u,flux_v,teta,phi,ucov,vcov,trac)
     4      SUBROUTINE bilan_dyn (dt_app,dt_cum,
     5     s  ps,masse,pk,flux_u,flux_v,teta,phi,ucov,vcov,
     6     s  ducovdyn,ducovdis,ducovspg,ducovphy)
     7c si besoin des traceurs:
     8c      SUBROUTINE bilan_dyn (ntrac,dt_app,dt_cum,
     9c     s  ps,masse,pk,flux_u,flux_v,teta,phi,ucov,vcov,trac,
     10c     s  ducovdyn,ducovdis,ducovspg,ducovphy)
    611
    712c   AFAIRE
     
    3944c   ===========
    4045
    41       integer ntrac
     46c      integer ntrac
    4247      real dt_app,dt_cum
    4348      real ps(iip1,jjp1)
     
    4954      real ucov(iip1,jjp1,llm)
    5055      real vcov(iip1,jjm,llm)
    51       real trac(iip1,jjp1,llm,ntrac)
     56c      real trac(iip1,jjp1,llm,ntrac)
     57c Tendances en m/s2 :
     58      real ducovdyn(iip1,jjp1,llm)
     59      real ducovdis(iip1,jjp1,llm)
     60      real ducovspg(iip1,jjp1,llm)
     61      real ducovphy(iip1,jjp1,llm)
    5262
    5363c   Local :
     
    5565
    5666      integer icum,ncum
     67      save icum,ncum
     68
     69      integer i,j,l,iQ,num
     70      real zz,zqy,zfactv(jjm,llm),zfactw(jjm,llm)
     71      character*2 strd2
     72      real ww
     73
    5774      logical first
    58       real zz,zqy,zfactv(jjm,llm)
    59 
    60       integer nQ
    61       parameter (nQ=7)
    62 
    63 
    64 cym      character*6 nom(nQ)
    65 cym      character*6 unites(nQ)
    66       character*6,save :: nom(nQ)
    67       character*6,save :: unites(nQ)
    68 
    69       character*10 file
    70       integer ifile
    71       parameter (ifile=4)
    72 
    73       integer itemp,igeop,iecin,iang,iu,iovap,iun
     75      save first
     76      data first/.true./
     77
    7478      integer i_sortie
    75 
    76       save first,icum,ncum
    77       save itemp,igeop,iecin,iang,iu,iovap,iun
    7879      save i_sortie
     80      data i_sortie/1/
    7981
    8082      real time
     
    8385      data time,itau/0.,0/
    8486
    85       data first/.true./
    86       data itemp,igeop,iecin,iang,iu,iovap,iun/1,2,3,4,5,6,7/
    87       data i_sortie/1/
    88 
    89       real ww
     87! facteur = -1. pour Venus
     88      real    fact_geovenus
     89      save    fact_geovenus
    9090
    9191c   variables dynamiques intermédiaires
     92c -----------------------------------
    9293      REAL vcont(iip1,jjm,llm),ucont(iip1,jjp1,llm)
    9394      REAL ang(iip1,jjp1,llm),unat(iip1,jjp1,llm)
     
    9596      REAL vorpot(iip1,jjm,llm)
    9697      REAL w(iip1,jjp1,llm),ecin(iip1,jjp1,llm),convm(iip1,jjp1,llm)
    97       REAL bern(iip1,jjp1,llm)
     98      real temp(iip1,jjp1,llm)
     99      real dudyn(iip1,jjp1,llm)
     100      real dudis(iip1,jjp1,llm)
     101      real duspg(iip1,jjp1,llm)
     102      real duphy(iip1,jjp1,llm)
     103
     104c CHAMPS SCALAIRES Q ADVECTES
     105c ----------------------------
     106      integer nQ
     107c avec tous les composes, ca fait trop.... Je les enleve
     108c     parameter (nQ=6+nqmx)
     109      parameter (nQ=6)
     110
     111      character*6,save :: nom(nQ)
     112      character*6,save :: unites(nQ)
     113
     114      integer itemp,igeop,iecin,iang,iu,iun
     115      save itemp,igeop,iecin,iang,iu,iun
     116      data itemp,igeop,iecin,iang,iu,iun/1,2,3,4,5,6/
    98117
    99118c   champ contenant les scalaires advectés.
     
    105124      real flux_u_cum(iip1,jjp1,llm)
    106125      real flux_v_cum(iip1,jjm,llm)
     126      real flux_w_cum(iip1,jjp1,llm)
    107127      real Q_cum(iip1,jjp1,llm,nQ)
    108128      real flux_uQ_cum(iip1,jjp1,llm,nQ)
     
    113133      save ps_cum,masse_cum,flux_u_cum,flux_v_cum
    114134      save Q_cum,flux_uQ_cum,flux_vQ_cum
    115 
    116 c   champs de tansport en moyenne zonale
     135      save flux_w_cum,flux_wQ_cum
     136
     137c   champs de transport en moyenne zonale
    117138      integer ntr,itr
    118139      parameter (ntr=5)
    119140
    120 cym      character*10 znom(ntr,nQ)
    121 cym      character*20 znoml(ntr,nQ)
    122 cym      character*10 zunites(ntr,nQ)
    123141      character*10,save :: znom(ntr,nQ)
    124142      character*20,save :: znoml(ntr,nQ)
    125143      character*10,save :: zunites(ntr,nQ)
     144      character*10,save :: znom2(ntr,nQ)
     145      character*20,save :: znom2l(ntr,nQ)
     146      character*10,save :: zunites2(ntr,nQ)
     147      character*10,save :: znom3(ntr,nQ)
     148      character*20,save :: znom3l(ntr,nQ)
     149      character*10,save :: zunites3(ntr,nQ)
    126150
    127151      integer iave,itot,immc,itrs,istn
     
    131155
    132156      real zvQ(jjm,llm,ntr,nQ),zvQtmp(jjm,llm)
     157      real zwQ(jjm,llm,ntr,nQ),zwQtmp(jjm,llm)
    133158      real zavQ(jjm,ntr,nQ),psiQ(jjm,llm+1,nQ)
    134       real zmasse(jjm,llm),zamasse(jjm)
    135 
    136       real zv(jjm,llm),psi(jjm,llm+1)
    137 
    138       integer i,j,l,iQ
    139 
     159      real zawQ(jjm,llm,ntr,nQ)
     160      real zdQ(jjm,llm,nQ)
     161      real zmasse(jjm,llm),zavmasse(jjm),zawmasse(llm)
     162      real psbar(jjm)
     163
     164      real zv(jjm,llm),zw(jjp1,llm),psi(jjm,llm+1)
     165
     166c TENDANCES POUR MOMENT CINETIQUE
     167c -------------------------------
     168
     169      integer ntdc,itdc
     170      parameter (ntdc=4)
     171
     172      integer itdcdyn,itdcdis,itdcspg,itdcphy
     173      data    itdcdyn,itdcdis,itdcspg,itdcphy/1,2,3,4/
     174
     175      character*6,save :: nomtdc(ntdc)
     176
     177c   champ contenant les tendances du moment cinetique.
     178      real    tdc(iip1,jjp1,llm,ntdc)
     179      real    ztdc(jjm,llm,ntdc)   ! moyenne zonale
     180   
     181c   champs cumulés
     182      real tdc_cum(iip1,jjp1,llm,ntdc)
     183      save tdc_cum
     184
     185c   integrations completes
     186      real mtot,mctot,dmctot(ntdc)
    140187
    141188c   Initialisation du fichier contenant les moyennes zonales.
     
    149196
    150197      integer ndex3d(jjm*llm)
     198      real    ztmp3d(jjm,llm)
    151199
    152200C   Variables locales
     
    154202      integer tau0
    155203      real zjulian
    156       character*3 str
    157       character*10 ctrac
    158       integer ii,jj
    159204      integer zan, dayref
    160205C
     
    167212c=====================================================================
    168213
    169       time=time+dt_app
    170       itau=itau+1
    171 cIM
    172214      ndex3d=0
    173215
    174216      if (first) then
    175217
     218        if (planet_type.eq."venus") then
     219            fact_geovenus = -1.
     220        else
     221            fact_geovenus = 1.
     222        endif
    176223
    177224        icum=0
     
    180227c   ncum est la frequence de stokage en pas de temps
    181228        ncum=dt_cum/dt_app
    182         if (abs(ncum*dt_app-dt_cum).gt.1.e-5*dt_app) then
     229        if (abs(ncum*dt_app-dt_cum).gt.1.e-2*dt_app) then
     230         if (abs((ncum+1)*dt_app-dt_cum).lt.1.e-2*dt_app) then
     231           ncum = ncum+1
     232         elseif (abs((ncum-1)*dt_app-dt_cum).lt.1.e-2*dt_app) then
     233           ncum = ncum-1
     234         else
    183235           WRITE(lunout,*)
    184236     .            'Pb : le pas de cumule doit etre multiple du pas'
    185237           WRITE(lunout,*)'dt_app=',dt_app
    186238           WRITE(lunout,*)'dt_cum=',dt_cum
     239           WRITE(lunout,*)'ncum*dt_app=',ncum*dt_app
     240           WRITE(lunout,*)'ncum=',ncum
    187241           stop
     242         endif
    188243        endif
    189244
    190         if (i_sortie.eq.1) then
    191          file='dynzon'
    192          call inigrads(ifile,1
    193      s  ,0.,180./pi,0.,0.,jjm,rlatv,-90.,90.,180./pi
    194      s  ,llm,presnivs,1.
    195      s  ,dt_cum,file,'dyn_zon ')
    196         endif
    197 
    198         nom(itemp)='T'
     245c VARIABLES ADVECTEES:
     246
     247        nom(itemp)='temp'
    199248        nom(igeop)='gz'
    200         nom(iecin)='K'
     249        nom(iecin)='ecin'
    201250        nom(iang)='ang'
    202251        nom(iu)='u'
    203         nom(iovap)='ovap'
    204252        nom(iun)='un'
    205253
     
    209257        unites(iang)='ang'
    210258        unites(iu)='m/s'
    211         unites(iovap)='kg/kg'
    212259        unites(iun)='un'
    213260
     261c avec tous les composes, ca fait trop.... Je les enleve
     262c       do num=1,ntrac
     263c        write(strd2,'(i2.2)') num
     264c        nom(6+num)='trac'//strd2
     265c        unites(6+num)='kg/kg'
     266c       enddo
     267
     268c TENDANCES MOMENT CIN:
     269       
     270        nomtdc(itdcdyn) ='dmcdyn'
     271        nomtdc(itdcdis) ='dmcdis'
     272        nomtdc(itdcspg) ='dmcspg'
     273        nomtdc(itdcphy) ='dmcphy'
    214274
    215275c   Initialisation du fichier contenant les moyennes zonales.
     
    221281      dayref = day_ref
    222282      CALL ymds2ju(zan, 1, dayref, 0.0, zjulian)
    223       tau0 = itau_dyn
     283c     tau0 = itau_dyn
     284      tau0 = 0
     285      itau = tau0
    224286     
    225287      rlong=0.
    226       rlatg=rlatv*180./pi
     288      rlatg=rlatv*180./pi*fact_geovenus
    227289       
    228290      call histbeg(infile, 1, rlong, jjm, rlatg,
     
    237299C
    238300C  Appels a histdef pour la definition des variables a sauvegarder
     301
    239302      do iQ=1,nQ
    240303         do itr=1,ntr
    241304            if(itr.eq.1) then
    242                znom(itr,iQ)=nom(iQ)
    243                znoml(itr,iQ)=nom(iQ)
    244                zunites(itr,iQ)=unites(iQ)
     305               znom(itr,iQ)    =nom(iQ)
     306               znoml(itr,iQ)   =nom(iQ)
     307               zunites(itr,iQ) =unites(iQ)
    245308            else
    246                znom(itr,iQ)=ctrs(itr)//'v'//nom(iQ)
    247                znoml(itr,iQ)='transport : v * '//nom(iQ)//' '//ctrs(itr)
    248                zunites(itr,iQ)='m/s * '//unites(iQ)
     309           znom(itr,iQ)    =ctrs(itr)//'v'//nom(iQ)
     310           znoml(itr,iQ)   ='transport : v * '//nom(iQ)//' '//ctrs(itr)
     311           zunites(itr,iQ) ='m/s * '//unites(iQ)
     312           znom2(itr,iQ)   =ctrs(itr)//'w'//nom(iQ)
     313           znom2l(itr,iQ)  ='transport: w * '//nom(iQ)//' '//ctrs(itr)
     314           zunites2(itr,iQ)='Pa/s * '//unites(iQ)
    249315            endif
    250316         enddo
     317               znom3(iQ)='d'//nom(iQ)
     318               znom3l(iQ)='convergence: '//nom(iQ)
     319               zunites3(iQ)=unites(iQ)//' /s'
     320c          print*,'DEBUG:',znom3(iQ),znom3l(iQ),zunites3(iQ)
    251321      enddo
    252322
    253323c   Declarations des champs avec dimension verticale
    254 c      print*,'1HISTDEF'
    255       do iQ=1,nQ
     324
     325      if (1.eq.0) then  ! on les sort, ou pas...
     326
     327c     do iQ=1,nQ
     328c !!!! JE NE SORS ICI QUE temp et ang POUR CAUSE DE PLACE !
     329      do iQ=1,4,3
    256330         do itr=1,ntr
    257331      IF (prt_level > 5)
     
    262336     .        32,'ave(X)',dt_cum,dt_cum)
    263337         enddo
     338c transport vertical:
     339         do itr=2,ntr
     340      IF (prt_level > 5)
     341     . WRITE(lunout,*)'var ',itr,iQ
     342     .      ,znom2(itr,iQ),znom2l(itr,iQ),zunites2(itr,iQ)
     343            call histdef(fileid,znom2(itr,iQ),znom2l(itr,iQ),
     344     .        zunites2(itr,iQ),1,jjm,thoriid,llm,1,llm,zvertiid,
     345     .        32,'ins(X)',dt_cum,dt_cum)
     346         enddo
     347
     348c Declarations pour convergences
     349      IF (prt_level > 5)
     350     . WRITE(lunout,*)'var ',iQ
     351     .      ,znom3(iQ),znom3l(iQ),zunites3(iQ)
     352            call histdef(fileid,znom3(iQ),znom3l(iQ),
     353     .        zunites3(iQ),1,jjm,thoriid,llm,1,llm,zvertiid,
     354     .        32,'ins(X)',dt_cum,dt_cum)
     355
    264356c   Declarations pour les fonctions de courant
    265 c      print*,'2HISTDEF'
    266           call histdef(fileid,'psi'//nom(iQ)
    267      .      ,'stream fn. '//znoml(itot,iQ),
    268      .      zunites(itot,iQ),1,jjm,thoriid,llm,1,llm,zvertiid,
    269      .      32,'ave(X)',dt_cum,dt_cum)
    270       enddo
    271 
     357c   Non sorties ici...
     358c          call histdef(fileid,'psi'//nom(iQ)
     359c     .      ,'stream fn. '//znoml(itot,iQ),
     360c     .      zunites(itot,iQ),1,jjm,thoriid,llm,1,llm,zvertiid,
     361c     .      32,'ave(X)',dt_cum,dt_cum)
     362
     363      enddo ! iQ
     364
     365      endif ! 1=1 sortie ou non...
    272366
    273367c   Declarations pour les champs de transport d'air
    274 c      print*,'3HISTDEF'
    275368      call histdef(fileid, 'masse', 'masse',
    276369     .             'kg', 1, jjm, thoriid, llm, 1, llm, zvertiid,
     
    279372     .             'm/s', 1, jjm, thoriid, llm, 1, llm, zvertiid,
    280373     .             32, 'ave(X)', dt_cum, dt_cum)
    281 c   Declarations pour les fonctions de courant
    282 c      print*,'4HISTDEF'
     374      call histdef(fileid, 'w', 'w',
     375     .             'Pa/s', 1, jjm, thoriid, llm, 1, llm, zvertiid,
     376     .             32, 'ins(X)', dt_cum, dt_cum)
     377
     378c   Declarations pour la fonction de courant
    283379          call histdef(fileid,'psi','stream fn. MMC ','mega t/s',
    284380     .      1,jjm,thoriid,llm,1,llm,zvertiid,
     
    286382
    287383
     384c   Declarations pour les tendances de moment cinetique
     385      do itdc=1,ntdc
     386      call histdef(fileid, nomtdc(itdc), nomtdc(itdc),
     387     .             'ang/s', 1, jjm, thoriid, llm, 1, llm, zvertiid,
     388     .             32, 'ins(X)', dt_cum, dt_cum)
     389      enddo
     390
    288391c   Declaration des champs 1D de transport en latitude
    289 c      print*,'5HISTDEF'
    290       do iQ=1,nQ
     392c     do iQ=1,nQ
     393c !!!! JE NE SORS ICI QUE temp et ang POUR CAUSE DE PLACE !
     394      do iQ=1,4,3
    291395         do itr=2,ntr
    292396            call histdef(fileid,'a'//znom(itr,iQ),znoml(itr,iQ),
    293397     .        zunites(itr,iQ),1,jjm,thoriid,1,1,1,-99,
    294398     .        32,'ave(X)',dt_cum,dt_cum)
    295          enddo
    296       enddo
    297 
    298 
    299 c      print*,'8HISTDEF'
     399c JE VIRE LE VERTICAL POUR L'INSTANT
     400c           call histdef(fileid,'a'//znom2(itr,iQ),znom2l(itr,iQ),
     401c    .        zunites2(itr,iQ),1,jjm,thoriid,llm,1,llm,zvertiid,
     402c    .        32,'ins(X)',dt_cum,dt_cum)
     403         enddo
     404      enddo
     405
    300406               CALL histend(fileid)
    301407
    302408
    303       endif
     409      endif  ! first
    304410
    305411
     
    313419      CALL enercin(vcov,ucov,vcont,ucont,ecin)
    314420
    315 c   moment cinétique
     421c   moment cinétique et tendances
     422      dudyn = 0.
     423      dudis = 0.
     424      duspg = 0.
     425      duphy = 0.
    316426      do l=1,llm
    317          ang(:,:,l)=ucov(:,:,l)+constang(:,:)
    318427         unat(:,:,l)=ucont(:,:,l)*cu(:,:)
    319       enddo
    320 
    321       Q(:,:,:,itemp)=teta(:,:,:)*pk(:,:,:)/cpp
    322       Q(:,:,:,igeop)=phi(:,:,:)
    323       Q(:,:,:,iecin)=ecin(:,:,:)
    324       Q(:,:,:,iang)=ang(:,:,:)
    325       Q(:,:,:,iu)=unat(:,:,:)
    326       Q(:,:,:,iovap)=trac(:,:,:,1)
    327       Q(:,:,:,iun)=1.
    328 
     428         dudyn(:,2:jjm,l) = ducovdyn(:,2:jjm,l)/cu(:,2:jjm)
     429         dudis(:,2:jjm,l) = ducovdis(:,2:jjm,l)/cu(:,2:jjm)
     430         duspg(:,2:jjm,l) = ducovspg(:,2:jjm,l)/cu(:,2:jjm)
     431         duphy(:,2:jjm,l) = ducovphy(:,2:jjm,l)/cu(:,2:jjm)
     432         do j=1,jjp1
     433          ang(:,j,l)= rad*cos(rlatu(j))*
     434     .     ( unat(:,j,l) + rad*cos(rlatu(j))*omeg )
     435          tdc(:,j,l,1) = rad*cos(rlatu(j))*dudyn(:,j,l)
     436          tdc(:,j,l,2) = rad*cos(rlatu(j))*dudis(:,j,l)
     437          tdc(:,j,l,3) = rad*cos(rlatu(j))*duspg(:,j,l)
     438          tdc(:,j,l,4) = rad*cos(rlatu(j))*duphy(:,j,l)
     439         enddo
     440      enddo
     441c Normalisation:
     442      ang = ang / (2./3. *rad*rad*omeg)
     443      do itdc=1,ntdc
     444        tdc(:,:,:,itdc)=tdc(:,:,:,itdc) / (2./3. *rad*rad*omeg)
     445      enddo
     446
     447! ADAPTATION GCM POUR CP(T)
     448      call tpot2t(ip1jmp1*llm,teta,temp,pk)
     449      Q(:,:,:,itemp) = temp(:,:,:)
     450      Q(:,:,:,igeop) =phi(:,:,:)
     451      Q(:,:,:,iecin) =ecin(:,:,:)
     452      Q(:,:,:,iang)  =ang(:,:,:)
     453      Q(:,:,:,iu)    =unat(:,:,:)
     454      Q(:,:,:,iun)   =1.
     455c avec tous les composes, ca fait trop.... Je les enleve
     456c     do num=1,ntrac
     457c      Q(:,:,:,6+num)=trac(:,:,:,num)
     458c     enddo
     459
     460c   calcul du flux de masse vertical (+ vers le bas)
     461      call convmas(flux_u,flux_v,convm)
     462      CALL vitvert(convm,w)
    329463
    330464c=====================================================================
     
    333467c
    334468      if(icum.EQ.0) then
    335          ps_cum=0.
    336          masse_cum=0.
    337          flux_u_cum=0.
    338          flux_v_cum=0.
    339          Q_cum=0.
    340          flux_vQ_cum=0.
    341          flux_uQ_cum=0.
     469         ps_cum      = 0.
     470         masse_cum   = 0.
     471         flux_u_cum  = 0.
     472         flux_v_cum  = 0.
     473         flux_w_cum  = 0.
     474         Q_cum       = 0.
     475         flux_vQ_cum = 0.
     476         flux_uQ_cum = 0.
     477         flux_wQ_cum = 0.
     478         tdc_cum     = 0.
    342479      endif
    343480
     
    347484
    348485c   accumulation des flux de masse horizontaux
    349       ps_cum=ps_cum+ps
    350       masse_cum=masse_cum+masse
    351       flux_u_cum=flux_u_cum+flux_u
    352       flux_v_cum=flux_v_cum+flux_v
     486      ps_cum          = ps_cum     + ps
     487      masse_cum       = masse_cum  + masse
     488      flux_u_cum      = flux_u_cum + flux_u
     489      flux_v_cum      = flux_v_cum + flux_v
     490      flux_w_cum      = flux_w_cum + w
    353491      do iQ=1,nQ
    354       Q_cum(:,:,:,iQ)=Q_cum(:,:,:,iQ)+Q(:,:,:,iQ)*masse(:,:,:)
     492      Q_cum(:,:,:,iQ) = Q_cum(:,:,:,iQ) + Q(:,:,:,iQ)*masse(:,:,:)
     493      enddo
     494      do itdc=1,ntdc
     495      tdc_cum(:,:,:,itdc) =
     496     .       tdc_cum(:,:,:,itdc) + tdc(:,:,:,itdc)*masse(:,:,:)
    355497      enddo
    356498
     
    386528      enddo
    387529
     530c   Flux vertical
     531c   -------------
     532      do iQ=1,nQ
     533         do l=2,llm
     534            do j=1,jjp1
     535               do i=1,iip1
     536                  flux_wQ_cum(i,j,l,iQ)=flux_wQ_cum(i,j,l,iQ)
     537     s            +w(i,j,l)*0.5*(Q(i,j,l-1,iQ)+Q(i,j,l,iQ))
     538               enddo
     539            enddo
     540         enddo
     541         flux_wQ_cum(:,:,1,iQ)=0.0e0
     542      enddo
    388543
    389544c    tendances
     
    397552      CALL vitvert(convm,w)
    398553
     554c  ajustement tendances (vitesse verticale)
    399555      do iQ=1,nQ
    400556         do l=1,llm-1
     
    410566      IF (prt_level > 5)
    411567     . WRITE(lunout,*)'Apres les calculs fait a chaque pas'
     568
    412569c=====================================================================
    413570c   PAS DE TEMPS D'ECRITURE
     
    415572      if (icum.eq.ncum) then
    416573c=====================================================================
     574
     575      time=time+dt_cum
     576      itau=itau+1
    417577
    418578      IF (prt_level > 5)
     
    421581c   Normalisation
    422582      do iQ=1,nQ
    423          Q_cum(:,:,:,iQ)=Q_cum(:,:,:,iQ)/masse_cum(:,:,:)
    424       enddo
     583         Q_cum(:,:,:,iQ) = Q_cum(:,:,:,iQ)/masse_cum(:,:,:)
     584         dQ(:,:,:,iQ)    = dQ(:,:,:,iQ)   /masse_cum(:,:,:)
     585      enddo
     586      do itdc=1,ntdc
     587         tdc_cum(:,:,:,itdc) = tdc_cum(:,:,:,itdc)/masse_cum(:,:,:)
     588      enddo
     589
    425590      zz=1./REAL(ncum)
    426       ps_cum=ps_cum*zz
    427       masse_cum=masse_cum*zz
    428       flux_u_cum=flux_u_cum*zz
    429       flux_v_cum=flux_v_cum*zz
    430       flux_uQ_cum=flux_uQ_cum*zz
    431       flux_vQ_cum=flux_vQ_cum*zz
    432       dQ=dQ*zz
    433 
    434 
    435 c   A retravailler eventuellement
    436 c   division de dQ par la masse pour revenir aux bonnes grandeurs
    437       do iQ=1,nQ
    438          dQ(:,:,:,iQ)=dQ(:,:,:,iQ)/masse_cum(:,:,:)
    439       enddo
    440  
     591      ps_cum      = ps_cum      *zz
     592      masse_cum   = masse_cum   *zz
     593      flux_u_cum  = flux_u_cum  *zz
     594      flux_v_cum  = flux_v_cum  *zz
     595      flux_w_cum  = flux_w_cum  *zz
     596      flux_uQ_cum = flux_uQ_cum *zz
     597      flux_vQ_cum = flux_vQ_cum *zz
     598      flux_wQ_cum = flux_wQ_cum *zz
     599
     600c Integration complete
     601      mtot  = 0.
     602      mctot  = 0.
     603      dmctot = 0.
     604      do l=1,llm
     605       do j=2,jjm
     606        do i=1,iim
     607          mtot  = mtot  + masse_cum(i,j,l)
     608          mctot = mctot + Q_cum(i,j,l,iang)*masse_cum(i,j,l)
     609        enddo
     610       enddo
     611      enddo
     612      mctot = mctot/mtot
     613      do itdc=1,ntdc
     614      do l=1,llm
     615       do j=2,jjm
     616        do i=1,iim
     617          dmctot(itdc) = dmctot(itdc)
     618     .               + tdc_cum(i,j,l,itdc)*masse_cum(i,j,l)/mtot
     619        enddo
     620       enddo
     621      enddo
     622      enddo
     623
    441624c=====================================================================
    442625c   Transport méridien
     
    446629c   ----------------------------------
    447630      zv=0.
     631      zw=0.
    448632      zmasse=0.
    449633      call massbar(masse_cum,massebx,masseby)
     634
     635c moy zonale de la ps cumulee
     636         do j=1,jjm
     637            psbar(j)=0.
     638            do i=1,iim
     639               psbar(j)=psbar(j)+ps_cum(i,j)/iim
     640            enddo
     641         enddo
     642
    450643      do l=1,llm
    451644         do j=1,jjm
     
    453646               zmasse(j,l)=zmasse(j,l)+masseby(i,j,l)
    454647               zv(j,l)=zv(j,l)+flux_v_cum(i,j,l)
     648               zw(j,l)=zw(j,l)+flux_w_cum(i,j,l)
    455649            enddo
    456650            zfactv(j,l)=cv(1,j)/zmasse(j,l)
    457          enddo
     651            zfactw(j,l)=((ap(l)-ap(l+1))+psbar(j)*(bp(l)-bp(l+1)))
     652     .                    /zmasse(j,l)
     653         enddo
     654            do i=1,iim
     655               zw(jjp1,l)=zw(jjp1,l)+flux_w_cum(i,jjp1,l)
     656            enddo
    458657      enddo
    459658
     
    480679c    la variable zfactv transforme un transport meridien cumule
    481680c    en kg/s * unte-du-champ-transporte en m/s * unite-du-champ-transporte
     681c    la variable zfactw transforme un transport vertical cumule
     682c    en kg/s * unte-du-champ-transporte en Pa/s * unite-du-champ-transporte
    482683c
    483684c   --------------------------------------------------------------
     
    489690
    490691      zvQ=0.
     692      zwQ=0.
     693      zdQ=0.
    491694      psiQ=0.
     695
    492696      do iQ=1,nQ
     697
     698c   transport meridien
    493699         zvQtmp=0.
    494700         do l=1,llm
    495701            do j=1,jjm
    496702c              print*,'j,l,iQ=',j,l,iQ
    497 c   Calcul des moyennes zonales du transort total et de zvQtmp
     703c   Calcul des moyennes zonales du transport total et de zvQtmp
    498704               do i=1,iim
    499705                  zvQ(j,l,itot,iQ)=zvQ(j,l,itot,iQ)
     
    518724         do l=llm,1,-1
    519725            do j=1,jjm
    520                psiQ(j,l,iQ)=psiQ(j,l+1,iQ)+zvQ(j,l,itot,iQ)
     726             psiQ(j,l,iQ)=psiQ(j,l+1,iQ)+zvQ(j,l,itot,iQ)/zfactv(j,l)
     727            enddo
     728         enddo
     729      enddo
     730
     731c   transport vertical
     732         zwQtmp=0.
     733         do l=1,llm
     734            do j=1,jjm
     735c              print*,'j,l,iQ=',j,l,iQ
     736c   Calcul des moyennes zonales du transport vertical total et de zwQtmp
     737               do i=1,iim
     738                  zwQ(j,l,itot,iQ)=zwQ(j,l,itot,iQ)
     739     s                            +flux_wQ_cum(i,j,l,iQ)
     740                  zqy=      0.5*(Q_cum(i,j,l,iQ)*masse_cum(i,j,l)+
     741     s                           Q_cum(i,j+1,l,iQ)*masse_cum(i,j+1,l))
     742                  zwQtmp(j,l)=zwQtmp(j,l)+flux_w_cum(i,j,l)*zqy
     743     s             /(0.5*(masse_cum(i,j,l)+masse_cum(i,j+1,l)))
     744                  zwQ(j,l,iave,iQ)=zwQ(j,l,iave,iQ)+zqy
     745               enddo
     746c   Decomposition
     747               zwQ(j,l,iave,iQ)=zwQ(j,l,iave,iQ)/zmasse(j,l)
     748               zwQ(j,l,itot,iQ)=zwQ(j,l,itot,iQ)*zfactw(j,l)
     749               zwQtmp(j,l)=zwQtmp(j,l)*zfactw(j,l)
     750               zwQ(j,l,immc,iQ)=zw(j,l)*zwQ(j,l,iave,iQ)*zfactw(j,l)
     751               zwQ(j,l,itrs,iQ)=zwQ(j,l,itot,iQ)-zwQtmp(j,l)
     752               zwQ(j,l,istn,iQ)=zwQtmp(j,l)-zwQ(j,l,immc,iQ)
     753            enddo
     754         enddo
     755
     756c   convergence
     757c   Calcul moyenne zonale de la convergence totale
     758         do l=1,llm
     759            do j=1,jjm
     760c              print*,'j,l,iQ=',j,l,iQ
     761               do i=1,iim
     762                  zdQ(j,l,iQ)=zdQ(j,l,iQ) +
     763     .                   ( dQ(i,j,l,iQ)   * masse_cum(i,j,l)
     764     .                   + dQ(i,j+1,l,iQ) * masse_cum(i,j+1,l))
     765     .                 / ( masse_cum(i,j,l)+masse_cum(i,j+1,l))
     766               enddo
    521767            enddo
    522768         enddo
     
    527773      do l=llm,1,-1
    528774         do j=1,jjm
    529             psi(j,l)=psi(j,l+1)+zv(j,l)
    530             zv(j,l)=zv(j,l)*zfactv(j,l)
     775            psi(j,l)= psi(j,l+1)+zv(j,l)
     776            zv(j,l) = zv(j,l)*zfactv(j,l)
     777            zw(j,l) = 0.5*(zw(j,l)+zw(j+1,l))*zfactw(j,l)
     778         enddo
     779      enddo
     780
     781c   Calcul moyenne zonale des tendances moment cin.
     782      ztdc=0.
     783      do itdc=1,ntdc
     784         do l=1,llm
     785            do j=1,jjm
     786               do i=1,iim
     787                  ztdc(j,l,itdc)=ztdc(j,l,itdc) +
     788     .            ( tdc_cum(i,j,l,itdc)   * masse_cum(i,j,l)
     789     .            + tdc_cum(i,j+1,l,itdc) * masse_cum(i,j+1,l))
     790     .          / ( masse_cum(i,j,l)+masse_cum(i,j+1,l))
     791               enddo
     792            enddo
    531793         enddo
    532794      enddo
    533795
    534796c     print*,'4OK'
     797
     798c--------------------------------------
     799c--------------------------------------
    535800c   sorties proprement dites
     801c--------------------------------------
     802c--------------------------------------
     803
    536804      if (i_sortie.eq.1) then
    537       do iQ=1,nQ
    538          do itr=1,ntr
    539             call histwrite(fileid,znom(itr,iQ),itau,zvQ(:,:,itr,iQ)
     805
     806c sortie des integrations completes dans le listing
     807      write(*,'(A12,5(1PE11.4,X))') "BILANMCDYN  ",mctot,dmctot
     808
     809c sorties dans fichier dynzon
     810
     811      if (1.eq.0) then  ! on les sort, ou pas...
     812
     813c avec tous les composes, ca fait trop.... Je les enleve
     814c      do iQ=1,nQ
     815c      do iQ=1,6
     816c !!!! JE NE SORS ICI QUE temp et ang POUR CAUSE DE PLACE !
     817      do iQ=1,4,3
     818
     819         ztmp3d(:,:)= zvQ(:,:,1,iQ) ! valeur moyenne
     820            call histwrite(fileid,znom(1,iQ),itau,ztmp3d
    540821     s      ,jjm*llm,ndex3d)
    541          enddo
    542          call histwrite(fileid,'psi'//nom(iQ),itau,psiQ(:,1:llm,iQ)
     822       do itr=2,ntr
     823         ztmp3d(:,:)= zvQ(:,:,itr,iQ)*fact_geovenus ! transport horizontal
     824            call histwrite(fileid,znom(itr,iQ),itau,ztmp3d
    543825     s      ,jjm*llm,ndex3d)
    544       enddo
    545 
    546       call histwrite(fileid,'masse',itau,zmasse
     826       enddo
     827
     828       do itr=2,ntr
     829         ztmp3d(:,:)=zwQ(:,:,itr,iQ)
     830            call histwrite(fileid,znom2(itr,iQ),itau,ztmp3d
     831     s      ,jjm*llm,ndex3d)
     832       enddo
     833       
     834         ztmp3d(:,:)= zdQ(:,:,iQ)
     835            call histwrite(fileid,znom3(iQ),itau,ztmp3d
     836     s      ,jjm*llm,ndex3d)
     837
     838c        ztmp3d(:,:)= psiQ(:,1:llm,iQ)*fact_geovenus
     839c        call histwrite(fileid,'psi'//nom(iQ),itau,ztmp3d
     840c    s      ,jjm*llm,ndex3d)
     841      enddo
     842
     843      endif ! 1=1 sortie ou non...
     844
     845      ztmp3d=zmasse
     846      call histwrite(fileid,'masse',itau,ztmp3d
    547847     s   ,jjm*llm,ndex3d)
    548       call histwrite(fileid,'v',itau,zv
     848     
     849      ztmp3d= zv*fact_geovenus
     850      call histwrite(fileid,'v',itau,ztmp3d
    549851     s   ,jjm*llm,ndex3d)
    550       psi=psi*1.e-9
    551       call histwrite(fileid,'psi',itau,psi(:,1:llm),jjm*llm,ndex3d)
    552 
    553       endif
     852      ztmp3d(:,:)=zw(1:jjm,:)
     853      call histwrite(fileid,'w',itau,ztmp3d
     854     s   ,jjm*llm,ndex3d)
     855      ztmp3d= psi(:,1:llm)*1.e-9*fact_geovenus
     856      call histwrite(fileid,'psi',itau,ztmp3d,jjm*llm,ndex3d)
     857
     858      do itdc=1,ntdc
     859         ztmp3d(:,:)= ztdc(:,:,itdc)
     860         call histwrite(fileid,nomtdc(itdc),itau,ztmp3d
     861     s    ,jjm*llm,ndex3d)
     862      enddo
     863
     864      endif ! i_sortie
    554865
    555866
     
    558869c   -----------------
    559870
    560       zamasse=0.
     871      zavmasse=0.
    561872      do l=1,llm
    562          zamasse(:)=zamasse(:)+zmasse(:,l)
     873         zavmasse(:)=zavmasse(:)+zmasse(:,l)
    563874      enddo
    564875      zavQ=0.
    565       do iQ=1,nQ
     876
     877c avec tous les composes, ca fait trop.... Je les enleve
     878c      do iQ=1,nQ
     879c      do iQ=1,6
     880c !!!! JE NE SORS ICI QUE temp et ang POUR CAUSE DE PLACE !
     881      do iQ=1,4,3
    566882         do itr=2,ntr
    567883            do l=1,llm
    568884               zavQ(:,itr,iQ)=zavQ(:,itr,iQ)+zvQ(:,l,itr,iQ)*zmasse(:,l)
    569885            enddo
    570             zavQ(:,itr,iQ)=zavQ(:,itr,iQ)/zamasse(:)
    571             call histwrite(fileid,'a'//znom(itr,iQ),itau,zavQ(:,itr,iQ)
    572      s      ,jjm*llm,ndex3d)
    573          enddo
    574       enddo
    575 
    576 c     on doit pouvoir tracer systematiquement la fonction de courant.
     886            zavQ(:,itr,iQ)=zavQ(:,itr,iQ)/zavmasse(:)
     887      if (i_sortie.eq.1) then
     888         ztmp3d=0.0
     889         ztmp3d(:,1)= zavQ(:,itr,iQ)*fact_geovenus
     890         call histwrite(fileid,'a'//znom(itr,iQ),itau,ztmp3d
     891     .      ,jjm*llm,ndex3d)     
     892      endif
     893         enddo
     894      enddo
     895
     896c   ------------------
     897c   Moyenne meridienne
     898c   ------------------
     899
     900      zawmasse=0.
     901      do j=1,jjm
     902           do l=1,llm
     903         zawmasse(l)=zawmasse(l)+zmasse(j,l)
     904           enddo
     905      enddo
     906      zawQ=0.
     907
     908c avec tous les composes, ca fait trop.... Je les enleve
     909c      do iQ=1,nQ
     910c      do iQ=1,6
     911c !!!! JE NE SORS ICI QUE temp et ang POUR CAUSE DE PLACE !
     912      do iQ=1,4,3
     913         do itr=2,ntr
     914           do l=1,llm
     915            do j=1,jjp1
     916          zawQ(1,l,itr,iQ)=zawQ(1,l,itr,iQ)+zwQ(j,l,itr,iQ)*zmasse(j,l)
     917            enddo
     918            zawQ(1,l,itr,iQ)=zawQ(1,l,itr,iQ)/zawmasse(l)
     919           enddo
     920      if (i_sortie.eq.1) then
     921         ztmp3d=0.0
     922           do l=1,llm
     923         ztmp3d(1,l)=zawQ(1,l,itr,iQ)
     924           enddo
     925c JE VIRE LE VERTICAL POUR L'INSTANT
     926c        call histwrite(fileid,'a'//znom2(itr,iQ),itau,ztmp3d
     927c    .      ,jjm*llm,ndex3d)     
     928      endif
     929         enddo
     930      enddo
     931
     932      call histsync(fileid)
    577933
    578934c=====================================================================
  • trunk/libf/dyn3d/caladvtrac.F

    r1 r6  
    5252        dq = 0.
    5353
     54      IF (planet_type.eq."earth") THEN
     55! Earth-specific treatment of first 2 tracers (water)
     56
    5457        CALL SCOPY( 2 * ijp1llm, q, 1, dq, 1 )
    5558
     
    7881          ENDDO
    7982         
    80           if (planet_type.eq."earth") then
    81 ! Earth-specific treatment of first 2 tracers (water)
    82             CALL qminimum( q, 2, finmasse )
    83           endif
     83          CALL qminimum( q, 2, finmasse )
    8484
    8585          CALL SCOPY   ( ip1jmp1*llm, masse, 1, finmasse,       1 )
     
    100100           ENDDO
    101101c
    102          ELSE
     102         ELSE 
    103103           DO iq = 1 , 2
    104104           DO l  = 1, llm
     
    110110
    111111
    112          ENDIF
     112         ENDIF ! iapptrac VS iapp_tracvl
     113
     114      ELSE ! not Earth
     115
     116c   advection
     117
     118        CALL advtrac( pbaru,pbarv,
     119     *       p,  masse,q,iapptrac, teta,
     120     .       flxw, pk)
     121c
     122
     123      ENDIF ! planet_type
    113124
    114125c
  • trunk/libf/dyn3d/calfis.F

    r5 r6  
    116116      REAL pducov(iip1,jjp1,llm)
    117117      REAL pdteta(iip1,jjp1,llm)
     118! commentaire SL: pdq ne sert que pour le calcul de pcvgq,
     119! qui lui meme ne sert a rien dans la routine telle qu'elle est
     120! ecrite, et que j'ai donc commente....
    118121      REAL pdq(iip1,jjp1,llm,nqtot)
    119122c
     
    122125      REAL ppk(iip1,jjp1,llm)
    123126c
     127c TENDENCIES in */s
    124128      REAL pdvfi(iip1,jjm,llm)
    125129      REAL pdufi(iip1,jjp1,llm)
     
    476480c   ---------------------
    477481
    478 #ifdef CPP_EARTH
     482! Appel de la physique: pose probleme quand on tourne
     483! SANS physique, car physiq.F est dans le repertoire phy[]...
     484! Il faut une cle CPP_PHYS
     485
     486! Le fait que les arguments de physiq soient differents selon les planetes
     487! ne pose pas de probleme a priori.
     488
     489! #ifdef CPP_PHYS
    479490
    480491!      write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys
     
    491502         lafin_split=lafin.and.isplit==nsplit_phys
    492503
     504      if (planet_type.eq."earth") then
    493505         CALL physiq (ngridmx,
    494506     .             llm,
     
    518530     .             PVteta)
    519531
    520          zufi(:,:)=zufi(:,:)+zdufi(:,:)*zdt_split
    521          zvfi(:,:)=zvfi(:,:)+zdvfi(:,:)*zdt_split
    522          ztfi(:,:)=ztfi(:,:)+zdtfi(:,:)*zdt_split
    523          zqfi(:,:,:)=zqfi(:,:,:)+zdqfi(:,:,:)*zdt_split
    524 
    525          zdufic(:,:)=zdufic(:,:)+zdufi(:,:)
    526          zdvfic(:,:)=zdvfic(:,:)+zdvfi(:,:)
    527          zdtfic(:,:)=zdtfic(:,:)+zdtfi(:,:)
    528          zdqfic(:,:,:)=zdqfic(:,:,:)+zdqfi(:,:,:)
    529 
    530       enddo
    531       zdufi(:,:)=zdufic(:,:)/nsplit_phys
    532       zdvfi(:,:)=zdvfic(:,:)/nsplit_phys
    533       zdtfi(:,:)=zdtfic(:,:)/nsplit_phys
    534       zdqfi(:,:,:)=zdqfic(:,:,:)/nsplit_phys
    535 
    536 #else 
    537 ! si il n'y a pas la cle CPP_EARTH
    538 
    539 ! ADAPTATION GCM POUR CP(T)
    540       CALL physiq (ngridmx,
     532      else  ! a moduler pour Mars !!
     533
     534         CALL physiq (ngridmx,
    541535     .             llm,
    542536     .             nqtot,
    543      .             debut,
    544      .             lafin,
     537     .             debut_split,
     538     .             lafin_split,
    545539     .             jD_cur,
    546      .             jH_cur,
    547      .             dtphys,
     540     .             jH_cur_split,
     541     .             zdt_split,
    548542     .             zplev,
    549543     .             zplay,
     
    564558     .             zdpsrf)
    565559
    566 #endif
     560      endif ! planet_type
     561
     562         zufi(:,:)=zufi(:,:)+zdufi(:,:)*zdt_split
     563         zvfi(:,:)=zvfi(:,:)+zdvfi(:,:)*zdt_split
     564         ztfi(:,:)=ztfi(:,:)+zdtfi(:,:)*zdt_split
     565         zqfi(:,:,:)=zqfi(:,:,:)+zdqfi(:,:,:)*zdt_split
     566
     567         zdufic(:,:)=zdufic(:,:)+zdufi(:,:)
     568         zdvfic(:,:)=zdvfic(:,:)+zdvfi(:,:)
     569         zdtfic(:,:)=zdtfic(:,:)+zdtfi(:,:)
     570         zdqfic(:,:,:)=zdqfic(:,:,:)+zdqfi(:,:,:)
     571
     572      enddo
     573      zdufi(:,:)=zdufic(:,:)/nsplit_phys
     574      zdvfi(:,:)=zdvfic(:,:)/nsplit_phys
     575      zdtfi(:,:)=zdtfic(:,:)/nsplit_phys
     576      zdqfi(:,:,:)=zdqfic(:,:,:)/nsplit_phys
     577
     578! #endif ! CPP_PHYS
    567579
    568580500   CONTINUE
     
    622634      pdqfi(:,:,:,:)=0.
    623635C
    624       DO iq=1,nqtot
     636       DO iq=1,nqtot
    625637         iiq=niadv(iq)
    626638         DO l=1,llm
     
    637649            ENDDO
    638650         ENDDO
    639       ENDDO
     651       ENDDO
    640652
    641653c   65. champ u:
  • trunk/libf/dyn3d/clesph0.h

    r1 r6  
    22! $Header$
    33!
     4! NE SERT A RIEN !! A VIRER... PAS A JOUR !!!
     5
    46c..include clesph0.h
    57c
  • trunk/libf/dyn3d/conf_gcm.F

    r1 r6  
    164164
    165165!Config  Key  = nsplit_phys
    166 !Config  Desc = nombre de pas par jour
     166!Config  Desc = nombre de subdivisions par pas physique
    167167!Config  Def  = 1
    168 !Config  Help = nombre de pas par jour (multiple de iperiod) (
    169 !Config          ici pour  dt = 1 min )
     168!Config  Help = nombre de subdivisions par pas physique
    170169       nsplit_phys = 1
    171170       CALL getin('nsplit_phys',nsplit_phys)
     
    361360       ip_ebil_dyn = 0
    362361       CALL getin('ip_ebil_dyn',ip_ebil_dyn)
    363 !
    364362
    365363      DO i = 1, longcles
     
    815813!Config  Help = active la version stratosphérique de LMDZ de F. Lott
    816814
    817       ok_strato=.FALSE.
     815      ok_strato=.TRUE.
    818816      CALL getin('ok_strato',ok_strato)
    819817
  • trunk/libf/dyn3d/defrun.F

    r1 r6  
    66      SUBROUTINE defrun( tapedef, etatinit, clesphy0 )
    77c
     8! ========================== ATTENTION =============================
     9! COMMENTAIRE SL :
     10! NE SERT PLUS APPAREMMENT
     11! DONC PAS MIS A JOUR POUR L'UTILISATION AVEC LES PLANETES
     12! ==================================================================
     13
    814      USE control_mod
    915 
  • trunk/libf/dyn3d/dynredem.F

    r1 r6  
    136136c
    137137      ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 27,
    138      .                       "Fichier demmarage dynamique")
     138     .                       "Fichier demarrage dynamique")
    139139c
    140140c Definir les dimensions du fichiers:
     
    536536#include "iniprint.h"
    537537
    538 
    539538      INTEGER l
    540539      REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm)
  • trunk/libf/dyn3d/gcm.F

    r5 r6  
    253253        endif
    254254
    255         if (planet_type.eq."earth") then
    256 #ifdef CPP_EARTH
    257 ! Load an Earth-format start file
     255        if (planet_type.eq."mars") then
     256! POUR MARS, METTRE UNE FONCTION A PART, genre dynetat0_mars
     257         abort_message = 'dynetat0_mars A FAIRE'
     258         call abort_gcm(modname,abort_message,0)
     259        else
    258260         CALL dynetat0("start.nc",vcov,ucov,
    259261     &              teta,q,masse,ps,phis, time_0)
    260 #else
    261         ! SW model also has Earth-format start files
    262         ! (but can be used without the CPP_EARTH directive)
    263           if (iflag_phys.eq.0) then
    264             CALL dynetat0("start.nc",vcov,ucov,
    265      &              teta,q,masse,ps,phis, time_0)
    266           endif
    267 #endif
    268         endif ! of if (planet_type.eq."earth")
     262        endif ! of if (planet_type.eq."mars")
    269263       
    270264c       write(73,*) 'ucov',ucov
     
    441435         WRITE(lunout,*)
    442436     .       'GCM: WARNING!!! vitesse verticale nulle dans la physique'
    443 ! Earth:
    444          if (planet_type.eq."earth") then
    445 #ifdef CPP_EARTH
     437
     438! Initialisation de la physique: pose probleme quand on tourne
     439! SANS physique, car iniphysiq.F est dans le repertoire phy[]...
     440! Il faut une cle CPP_PHYS
     441! #ifdef CPP_PHYS
    446442         CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys ,
    447443     ,                latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp     )
    448 #endif
    449          endif ! of if (planet_type.eq."earth")
     444! #endif ! CPP_PHYS
    450445         call_iniphys=.false.
    451446      ENDIF ! of IF (call_iniphys.and.(iflag_phys.eq.1))
     
    472467#endif
    473468
    474       if (planet_type.eq."earth") then
     469      if (planet_type.eq."mars") then
     470! POUR MARS, METTRE UNE FONCTION A PART, genre dynredem0_mars
     471         abort_message = 'dynredem0_mars A FAIRE'
     472         call abort_gcm(modname,abort_message,0)
     473      else
    475474        CALL dynredem0("restart.nc", day_end, phis)
    476       endif
     475      endif ! of if (planet_type.eq."mars")
    477476
    478477      ecripar = .TRUE.
  • trunk/libf/dyn3d/infotrac.F90

    r1 r6  
    8686    descrq(30)='PRA'
    8787   
    88 
    89     IF (config_inca=='none') THEN
     88    IF (planet_type=='earth') THEN
     89     IF (config_inca=='none') THEN
    9090       type_trac='lmdz'
     91     ELSE
     92       type_trac='inca'
     93     END IF
    9194    ELSE
    92        type_trac='inca'
    93     END IF
     95     type_trac='plnt'  ! planets... May want to dissociate between each later.
     96    ENDIF
    9497
    9598!-----------------------------------------------------------------------
     
    99102!
    100103!-----------------------------------------------------------------------
    101     IF (type_trac == 'lmdz') THEN
     104    IF (planet_type=='earth') THEN
     105     IF (type_trac == 'lmdz') THEN
    102106       OPEN(90,file='traceur.def',form='formatted',status='old', iostat=ierr)
    103107       IF(ierr.EQ.0) THEN
     
    109113          nqtrue=4 ! Defaut value
    110114       END IF
    111        ! Attention! Only for planet_type=='earth'
    112115       nbtr=nqtrue-2
    113     ELSE
     116     ELSE
    114117       ! nbtr has been read from INCA by init_cont_lmdz() in gcm.F
    115118       nqtrue=nbtr+2
    116     END IF
    117 
    118     IF (nqtrue < 2) THEN
     119     END IF
     120
     121     IF (nqtrue < 2) THEN
    119122       WRITE(lunout,*) 'nqtrue=',nqtrue, ' is not allowded. 2 tracers is the minimum'
    120123       CALL abort_gcm('infotrac_init','Not enough tracers',1)
    121     END IF
     124     END IF
     125
     126    ELSE  ! not Earth
     127       OPEN(90,file='traceur.def',form='formatted',status='old', iostat=ierr)
     128       IF(ierr.EQ.0) THEN
     129          WRITE(lunout,*) 'Open traceur.def : ok'
     130          READ(90,*) nqtrue
     131       ELSE
     132          WRITE(lunout,*) 'Problem in opening traceur.def'
     133          WRITE(lunout,*) 'ATTENTION using defaut values: nqtrue=1'
     134          nqtrue=1 ! Defaut value
     135       END IF
     136       nbtr=nqtrue
     137     
     138    ENDIF  ! planet_type
    122139!
    123140! Allocate variables depending on nqtrue and nbtr
     
    154171!    Get choice of advection schema from file tracer.def or from INCA
    155172!---------------------------------------------------------------------
    156     IF (type_trac == 'lmdz') THEN
     173    IF (planet_type=='earth') THEN
     174     IF (type_trac == 'lmdz') THEN
    157175       IF(ierr.EQ.0) THEN
    158176          ! Continue to read tracer.def
     
    182200       END DO
    183201
    184     ELSE  ! type_trac=inca : config_inca='aero' ou 'chem'
     202     ELSE  ! type_trac=inca : config_inca='aero' ou 'chem'
    185203! le module de chimie fournit les noms des traceurs
    186204! et les schemas d'advection associes.
     
    201219       END DO
    202220
    203     END IF ! type_trac
     221     END IF ! type_trac
     222
     223    ELSE  ! not Earth
     224       IF(ierr.EQ.0) THEN
     225          ! Continue to read tracer.def
     226          DO iq=1,nqtrue
     227             READ(90,999) hadv(iq),vadv(iq),tnom_0(iq)
     228          END DO
     229          CLOSE(90) 
     230       ELSE ! Without tracer.def
     231          hadv(1) = 10
     232          vadv(1) = 10
     233          tnom_0(1) = 'dummy'
     234       END IF
     235       
     236       WRITE(lunout,*) 'Valeur de traceur.def :'
     237       WRITE(lunout,*) 'nombre de traceurs ',nqtrue
     238       DO iq=1,nqtrue
     239          WRITE(lunout,*) hadv(iq),vadv(iq),tnom_0(iq)
     240       END DO
     241
     242    ENDIF  ! planet_type
    204243
    205244!-----------------------------------------------------------------------
  • trunk/libf/dyn3d/leapfrog.F

    r5 r6  
    9999      REAL massem1(ip1jmp1,llm)
    100100
    101 c   tendances dynamiques
     101c   tendances dynamiques en */s
    102102      REAL dv(ip1jm,llm),du(ip1jmp1,llm)
    103103      REAL dteta(ip1jmp1,llm),dq(ip1jmp1,llm,nqtot),dp(ip1jmp1)
    104104
    105 c   tendances de la dissipation
     105c   tendances de la dissipation en */s
    106106      REAL dvdis(ip1jm,llm),dudis(ip1jmp1,llm)
    107107      REAL dtetadis(ip1jmp1,llm)
    108108
    109 c   tendances physiques
     109c   tendances de la couche superieure */s
     110      REAL dvtop(ip1jm,llm),dutop(ip1jmp1,llm)
     111      REAL dtetatop(ip1jmp1,llm)
     112      REAL dqtop(ip1jmp1,llm,nqtot),dptop(ip1jmp1)
     113
     114c   tendances physiques */s
    110115      REAL dvfi(ip1jm,llm),dufi(ip1jmp1,llm)
    111116      REAL dtetafi(ip1jmp1,llm),dqfi(ip1jmp1,llm,nqtot),dpfi(ip1jmp1)
     
    182187      data dissip_conservative/.true./
    183188
    184       LOGICAL prem
    185       save prem
    186       DATA prem/.true./
    187189      INTEGER testita
    188190      PARAMETER (testita = 9)
     
    190192      logical , parameter :: flag_verif = .false.
    191193     
    192 
    193       integer itau_w   ! pas de temps ecriture = itap + itau_phy
    194 
    195194
    196195      itaufin   = nday*day_step
     
    198197      modname="leapfrog"
    199198     
     199c INITIALISATIONS
     200        dudis(:,:)   =0.
     201        dvdis(:,:)   =0.
     202        dtetadis(:,:)=0.
     203        dutop(:,:)   =0.
     204        dvtop(:,:)   =0.
     205        dtetatop(:,:)=0.
     206        dqtop(:,:,:) =0.
     207        dptop(:)     =0.
     208        dufi(:,:)   =0.
     209        dvfi(:,:)   =0.
     210        dtetafi(:,:)=0.
     211        dqfi(:,:,:) =0.
     212        dpfi(:)     =0.
    200213
    201214      itau = 0
     
    382395
    383396
    384 c   Inbterface avec les routines de phylmd (phymars ... )
     397c   Interface avec les routines de phylmd (phymars ... )
    385398c   -----------------------------------------------------
    386399
     
    400413cIM : pour sortir les param. du modele dans un fis. netcdf 110106
    401414         IF (first) THEN
    402           first=.false.
    403415#include "ini_paramLMDZ_dyn.h"
    404416         ENDIF
     
    408420#endif
    409421! #endif of #ifdef CPP_IOIPSL
     422
    410423         CALL calfis( lafin , jD_cur, jH_cur,
    411424     $               ucov,vcov,teta,q,masse,ps,p,pk,phis,phi ,
     
    414427     $               clesphy0, dufi,dvfi,dtetafi,dqfi,dpfi  )
    415428
     429c      ajout des tendances physiques:
     430c      ------------------------------
     431          CALL addfi( dtphys, leapf, forward   ,
     432     $                  ucov, vcov, teta , q   ,ps ,
     433     $                 dufi, dvfi, dtetafi , dqfi ,dpfi  )
     434
     435c      Couche superieure :
     436c      -------------------
    416437         IF (ok_strato) THEN
    417            CALL top_bound( vcov,ucov,teta,masse,dufi,dvfi,dtetafi)
     438           CALL top_bound( vcov,ucov,teta,masse,dutop,dvtop,dtetatop)
    418439         ENDIF
     440c dqtop=0, dptop=0
    419441       
    420442c      ajout des tendances physiques:
     
    422444          CALL addfi( dtphys, leapf, forward   ,
    423445     $                  ucov, vcov, teta , q   ,ps ,
    424      $                 dufi, dvfi, dtetafi , dqfi ,dpfi  )
     446     $                 dutop, dvtop, dtetatop , dqtop ,dptop  )
    425447c
    426448c  Diagnostique de conservation de l'énergie : difference
     
    448470        ! Sponge layer (if any)
    449471        IF (ok_strato) THEN
    450           dufi(:,:)=0.
    451           dvfi(:,:)=0.
    452           dtetafi(:,:)=0.
    453           dqfi(:,:,:)=0.
    454           dpfi(:)=0.
    455           CALL top_bound(vcov,ucov,teta,masse,dufi,dvfi,dtetafi)
     472          CALL top_bound(vcov,ucov,teta,masse,dutop,dvtop,dtetatop)
    456473          CALL addfi( dtvr, leapf, forward   ,
    457474     $                  ucov, vcov, teta , q   ,ps ,
    458      $                 dufi, dvfi, dtetafi , dqfi ,dpfi  )
     475     $                 dutop, dvtop, dtetatop , dqtop ,dptop  )
    459476        ENDIF ! of IF (ok_strato)
    460477      ENDIF ! of IF (iflag_phys.EQ.2)
     
    485502        ucov=ucov+dudis
    486503        vcov=vcov+dvdis
    487 c       teta=teta+dtetadis
    488 
     504        dudis=dudis/dtdiss   ! passage en (m/s)/s
     505        dvdis=dvdis/dtdiss   ! passage en (m/s)/s
    489506
    490507c------------------------------------------------------------------------
     
    506523        endif
    507524        teta=teta+dtetadis
     525        dtetadis=dtetadis/dtdiss   ! passage en K/s
    508526c------------------------------------------------------------------------
    509527
     
    600618               IF (ok_dynzon) THEN
    601619#ifdef CPP_IOIPSL
    602                  CALL bilan_dyn(2,dtvr*iperiod,dtvr*day_step*periodav,
    603      &                 ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q)
     620c les traceurs ne sont pas sortis, trop lourd.
     621c Peut changer eventuellement si besoin.
     622                 CALL bilan_dyn(dtvr*iperiod,dtvr*day_step*periodav,
     623     &                 ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,
     624     &                 du,dudis,duspg,dufi)
    604625#endif
    605626               END IF
     
    651672
    652673
    653               if (planet_type.eq."earth") then
    654 ! Write an Earth-format restart file
     674              if (planet_type.eq."mars") then
     675! POUR MARS, METTRE UNE FONCTION A PART, genre dynredem1_mars
     676                abort_message = 'dynredem1_mars A FAIRE'
     677                call abort_gcm(modname,abort_message,0)
     678              else
    655679                CALL dynredem1("restart.nc",0.0,
    656680     &                         vcov,ucov,teta,q,masse,ps)
    657               endif ! of if (planet_type.eq."earth")
     681              endif ! of if (planet_type.eq."mars")
    658682
    659683              CLOSE(99)
     
    726750               IF (ok_dynzon) THEN
    727751#ifdef CPP_IOIPSL
    728                  CALL bilan_dyn(2,dtvr*iperiod,dtvr*day_step*periodav,
    729      &                 ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q)
     752c les traceurs ne sont pas sortis, trop lourd.
     753c Peut changer eventuellement si besoin.
     754                 CALL bilan_dyn(dtvr*iperiod,dtvr*day_step*periodav,
     755     &                 ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,
     756     &                 du,dudis,duspg,dufi)
    730757#endif
    731758               ENDIF
     
    766793
    767794              IF(itau.EQ.itaufin) THEN
    768                 if (planet_type.eq."earth") then
     795                if (planet_type.eq."mars") then
     796! POUR MARS, METTRE UNE FONCTION A PART, genre dynredem1_mars
     797                  abort_message = 'dynredem1_mars A FAIRE'
     798                  call abort_gcm(modname,abort_message,0)
     799                else
    769800                  CALL dynredem1("restart.nc",0.0,
    770      &                           vcov,ucov,teta,q,masse,ps)
    771                 endif ! of if (planet_type.eq."earth")
     801     &                         vcov,ucov,teta,q,masse,ps)
     802                endif ! of if (planet_type.eq."mars")
    772803              ENDIF ! of IF(itau.EQ.itaufin)
    773804
     
    779810      END IF ! of IF(.not.purmats)
    780811
     812      first=.false.
     813
    781814      STOP
    782815      END
Note: See TracChangeset for help on using the changeset viewer.