Ignore:
Timestamp:
Mar 30, 2009, 4:46:54 PM (16 years ago)
Author:
Ehouarn Millour
Message:

Premiere vaque de modifications pour l'unification des dynamiques (planetes-Terre) et un peu de netoyage ...

  • modified 'makegcm' and 'makegcm_fcm' to remove 'CPP_PHYS' key and add 'CPP_EARTH' preprocessing key instead
  • updated 'diagedyn.F' (in dyn3d and dyn3dpar) to use 'CPP_EARTH' key
  • added 'ioipsl_getincom.F90' and 'ioipsl_stringop.F90' to 'dyn3d' and 'dyn3dpar' for future possibility of running without IOIPSL library
  • modified conf_gcm.F ( in d'yn3d' and 'dyn3dpar') to read in flag 'planet_type' (default=='earth') (flag added in 'control.h')
  • modified 'gcm.F' (in 'dyn3d' and 'dyn3dpar') so that flags so that 'read_start' and 'iflag_phys' (known from conf_gcm.F) are used
  • added flag 'output_grads_dyn' (read by conf_gcm.F, stored in 'control.h') to write grads outputs from 'leapfrog.F' and 'leapfrog_p.F'
  • removed 'comdiss.h' from 'dyn3d' and 'dyn3dpar' (it is not used)
  • removed variable 'lstardis' from 'comdissip.h' (it is also in

'comdissnew.h'), in dyn3d as well as in dyn3dpar

  • adapted 'dyn3d/iniacademic.F' to not use 'inicons0.F' but 'iniconst.F'
  • updated 'dyn3d/etat0_netcdf.F' to not use 'inicons0' but 'iniconst' (added prerequisite pa=50000 instruction) and added #ifdef CPP_EARTH keys
  • removed 'inicons0.F' and 'disvert0.F' (not used any more)
Location:
LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar
Files:
2 added
3 deleted
8 edited

Legend:

Unmodified
Added
Removed

  • LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/comdissip.h

    r774 r1140  
    22! $Header$
    33!
    4 c-----------------------------------------------------------------------
    5 c INCLUDE dissip.h
     4!-----------------------------------------------------------------------
     5! INCLUDE comdissip.h
    66
    7       COMMON/comdissip/
    8      $    lstardis,niterdis,coefdis,tetavel,tetatemp,gamdissip
     7      COMMON/comdissip/                                                 &
     8     &    niterdis,coefdis,tetavel,tetatemp,gamdissip
    99
    1010
    11       LOGICAL lstardis
    1211      INTEGER niterdis
    1312
    1413      REAL tetavel,tetatemp,coefdis,gamdissip
    1514
    16 c-----------------------------------------------------------------------
     15!-----------------------------------------------------------------------

  • LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/conf_gcm.F

    r1111 r1140  
    66      SUBROUTINE conf_gcm( tapedef, etatinit, clesphy0 )
    77c
     8#ifdef CPP_IOIPSL
    89      use IOIPSL
     10#else
     11! if not using IOIPSL, we still need to use (a local version of) getin
     12      use ioipsl_getincom
     13#endif
    914      use misc_mod
    1015      use mod_filtre_fft, ONLY : use_filtre_fft
    … …  
    99104
    100105!Config  Key  = prt_level
    101 !Config  Desc = niveau d'impressions de d�bogage
    102 !Config  Def  = 0
    103 !Config  Help = Niveau d'impression pour le d�bogage
     106!Config  Desc = niveau d'impressions de débogage
     107!Config  Def  = 0
     108!Config  Help = Niveau d'impression pour le débogage
    104109!Config         (0 = minimum d'impression)
    105110      prt_level = 0
    … …  
    109114c  Parametres de controle du run:
    110115c-----------------------------------------------------------------------
     116!Config  Key  = planet_type
     117!Config  Desc = planet type ("earth", "mars", "venus", ...)
     118!Config  Def  = earth
     119!Config  Help = this flag sets the type of atymosphere that is considered
     120      planet_type="earth"
     121      CALL getin('planet_type',planet_type)
    111122
    112123!Config  Key  = dayref
    … …  
    189200       CALL getin('periodav',periodav)
    190201
     202!Config  Key  = output_grads_dyn
     203!Config  Desc = output dynamics diagnostics in 'dyn.dat' file
     204!Config  Def  = n
     205!Config  Help = output dynamics diagnostics in Grads-readable 'dyn.dat' file
     206       output_grads_dyn=.false.
     207       CALL getin('output_grads_dyn',output_grads_dyn)
     208
    191209!Config  Key  = idissip
    192210!Config  Desc = periode de la dissipation
    … …  
    284302c    ...............................................................
    285303
     304!Config  Key  =  read_start
     305!Config  Desc = Initialize model using a 'start.nc' file
     306!Config  Def  = y
     307!Config  Help = y: intialize dynamical fields using a 'start.nc' file
     308!               n: fields are initialized by 'iniacademic' routine
     309       read_start= .true.
     310       CALL getin('read_start',read_start)
     311
    286312!Config  Key  = iflag_phys
    287313!Config  Desc = Avec ls physique
    … …  
    341367c
    342368      IF( ABS(clat - clatt).GE. 0.001 )  THEN
    343         PRINT *,' La valeur de clat passee par run.def est differente de
    344      * celle lue sur le fichier  start '
     369        write(lunout,*)'conf_gcm: La valeur de clat passee par run.def',
     370     &    ' est differente de celle lue sur le fichier  start '
    345371        STOP
    346372      ENDIF
    … …  
    356382
    357383      IF( ABS(grossismx - grossismxx).GE. 0.001 )  THEN
    358         PRINT *,' La valeur de grossismx passee par run.def est differente 
    359      * de celle lue sur le fichier  start '
     384        write(lunout,*)'conf_gcm: La valeur de grossismx passee par ',
     385     &  'run.def est differente de celle lue sur le fichier  start '
    360386        STOP
    361387      ENDIF
    … …  
    370396
    371397      IF( ABS(grossismy - grossismyy).GE. 0.001 )  THEN
    372         PRINT *,' La valeur de grossismy passee par run.def est differen
    373      * te de celle lue sur le fichier  start '
     398        write(lunout,*)'conf_gcm: La valeur de grossismy passee par ',
     399     & 'run.def est differente de celle lue sur le fichier  start '
    374400        STOP
    375401      ENDIF
    376402     
    377403      IF( grossismx.LT.1. )  THEN
    378         PRINT *,' ***  ATTENTION !! grossismx < 1 .   *** '
     404        write(lunout,*)
     405     &       'conf_gcm: ***  ATTENTION !! grossismx < 1 .   *** '
    379406         STOP
    380407      ELSE
    … …  
    384411
    385412      IF( grossismy.LT.1. )  THEN
    386         PRINT *,' ***  ATTENTION !! grossismy < 1 .   *** '
     413        write(lunout,*)
     414     &       'conf_gcm: ***  ATTENTION !! grossismy < 1 .   *** '
    387415         STOP
    388416      ELSE
    … …  
    390418      ENDIF
    391419
    392       PRINT *,' alphax alphay defrun ',alphax,alphay
     420      write(lunout,*)'conf_gcm: alphax alphay',alphax,alphay
    393421c
    394422c    alphax et alphay sont les anciennes formulat. des grossissements
    … …  
    405433
    406434      IF( .NOT.fxyhypb )  THEN
    407            IF( fxyhypbb )     THEN
    408               PRINT *,' ********  PBS DANS  DEFRUN  ******** '
    409               PRINT *,' *** fxyhypb lu sur le fichier start est F ',
    410      *       'alors  qu il est  T  sur  run.def  ***'
     435         IF( fxyhypbb )     THEN
     436            write(lunout,*)' ********  PBS DANS  CONF_GCM  ******** '
     437            write(lunout,*)' *** fxyhypb lu sur le fichier start est ',
     438     *       'F alors  qu il est  T  sur  run.def  ***'
    411439              STOP
    412            ENDIF
     440         ENDIF
    413441      ELSE
    414            IF( .NOT.fxyhypbb )   THEN
    415               PRINT *,' ********  PBS DANS  DEFRUN  ******** '
    416               PRINT *,' ***  fxyhypb lu sur le fichier start est T ',
    417      *        'alors  qu il est  F  sur  run.def  ****  '
     442         IF( .NOT.fxyhypbb )   THEN
     443            write(lunout,*)' ********  PBS DANS  CONF_GCM  ******** '
     444            write(lunout,*)' ***  fxyhypb lu sur le fichier start est ',
     445     *        'T alors  qu il est  F  sur  run.def  ****  '
    418446              STOP
    419            ENDIF
     447         ENDIF
    420448      ENDIF
    421449c
    … …  
    430458      IF( fxyhypb )  THEN
    431459       IF( ABS(dzoomx - dzoomxx).GE. 0.001 )  THEN
    432         PRINT *,' La valeur de dzoomx passee par run.def est differente
    433      *  de celle lue sur le fichier  start '
     460        write(lunout,*)'conf_gcm: La valeur de dzoomx passee par ',
     461     *  'run.def est differente de celle lue sur le fichier  start '
    434462        STOP
    435463       ENDIF
    … …  
    446474      IF( fxyhypb )  THEN
    447475       IF( ABS(dzoomy - dzoomyy).GE. 0.001 )  THEN
    448         PRINT *,' La valeur de dzoomy passee par run.def est differente
    449      * de celle lue sur le fichier  start '
     476        write(lunout,*)'conf_gcm: La valeur de dzoomy passee par ',
     477     * 'run.def est differente de celle lue sur le fichier  start '
    450478        STOP
    451479       ENDIF
    … …  
    461489      IF( fxyhypb )  THEN
    462490       IF( ABS(taux - tauxx).GE. 0.001 )  THEN
    463         PRINT *,' La valeur de taux passee par run.def est differente
    464      * de celle lue sur le fichier  start '
     491        write(lunout,*)'conf_gcm: La valeur de taux passee par ',
     492     * 'run.def est differente de celle lue sur le fichier  start '
    465493        STOP
    466494       ENDIF
    … …  
    476504      IF( fxyhypb )  THEN
    477505       IF( ABS(tauy - tauyy).GE. 0.001 )  THEN
    478         PRINT *,' La valeur de tauy passee par run.def est differente
    479      * de celle lue sur le fichier  start '
     506        write(lunout,*)'conf_gcm: La valeur de tauy passee par ',
     507     * 'run.def est differente de celle lue sur le fichier  start '
    480508        STOP
    481509       ENDIF
    … …  
    495523
    496524        IF( .NOT.ysinus )  THEN
    497            IF( ysinuss )     THEN
    498               PRINT *,' ********  PBS DANS  DEFRUN  ******** '
    499               PRINT *,' *** ysinus lu sur le fichier start est F ',
    500      *       'alors  qu il est  T  sur  run.def  ***'
     525          IF( ysinuss )     THEN
     526            write(lunout,*)' ********  PBS DANS  CONF_GCM  ******** '
     527            write(lunout,*)' *** ysinus lu sur le fichier start est F',
     528     *       ' alors  qu il est  T  sur  run.def  ***'
     529            STOP
     530          ENDIF
     531        ELSE
     532          IF( .NOT.ysinuss )   THEN
     533            write(lunout,*)' ********  PBS DANS  CONF_GCM  ******** '
     534            write(lunout,*)' *** ysinus lu sur le fichier start est T',
     535     *        ' alors  qu il est  F  sur  run.def  ****  '
    501536              STOP
    502            ENDIF
    503         ELSE
    504            IF( .NOT.ysinuss )   THEN
    505               PRINT *,' ********  PBS DANS  DEFRUN  ******** '
    506               PRINT *,' ***  ysinus lu sur le fichier start est T ',
    507      *        'alors  qu il est  F  sur  run.def  ****  '
    508               STOP
    509            ENDIF
     537          ENDIF
    510538        ENDIF
    511       ENDIF
     539      ENDIF ! of IF( .NOT.fxyhypb  )
    512540c
    513541!Config  Key  = offline
    … …  
    532560      write(lunout,*)' #########################################'
    533561      write(lunout,*)' Configuration des parametres du gcm: '
     562      write(lunout,*)' planet_type = ', planet_type
    534563      write(lunout,*)' dayref = ', dayref
    535564      write(lunout,*)' anneeref = ', anneeref
    … …  
    540569      write(lunout,*)' iecri = ', iecri
    541570      write(lunout,*)' periodav = ', periodav
     571      write(lunout,*)' output_grads_dyn = ', output_grads_dyn
    542572      write(lunout,*)' idissip = ', idissip
    543573      write(lunout,*)' lstardis = ', lstardis
    … …  
    550580      write(lunout,*)' coefdis = ', coefdis
    551581      write(lunout,*)' purmats = ', purmats
     582      write(lunout,*)' read_start = ', read_start
    552583      write(lunout,*)' iflag_phys = ', iflag_phys
    553584      write(lunout,*)' clonn = ', clonn
    … …  
    600631
    601632      IF( grossismx.LT.1. )  THEN
    602         PRINT *,' ***  ATTENTION !! grossismx < 1 .   *** '
     633        write(lunout,*)
     634     &   'conf_gcm: ***  ATTENTION !! grossismx < 1 .   *** '
    603635         STOP
    604636      ELSE
    … …  
    608640
    609641      IF( grossismy.LT.1. )  THEN
    610         PRINT *,' ***  ATTENTION !! grossismy < 1 .   *** '
     642        write(lunout,*)
     643     &  'conf_gcm: ***  ATTENTION !! grossismy < 1 .   *** '
    611644         STOP
    612645      ELSE
    … …  
    614647      ENDIF
    615648
    616       PRINT *,' alphax alphay defrun ',alphax,alphay
     649      write(lunout,*)'conf_gcm: alphax alphay ',alphax,alphay
    617650c
    618651c    alphax et alphay sont les anciennes formulat. des grossissements
    … …  
    697730        write(lunout,*)"Le zoom en longitude est incompatible",
    698731     &                 " avec l'utilisation du filtre FFT ",
    699      &                 "---> filtre FFT désactivé "
     732     &                 "---> filtre FFT désactivé "
    700733       use_filtre_fft=.FALSE.
    701734      ENDIF
    … …  
    704737     
    705738!Config  Key  = use_mpi_alloc
    706 !Config  Desc = Utilise un buffer MPI en m�moire globale
     739!Config  Desc = Utilise un buffer MPI en m�moire globale
    707740!Config  Def  = false
    708741!Config  Help = permet d'activer l'utilisation d'un buffer MPI
    709 !Config         en m�moire globale a l'aide de la fonction MPI_ALLOC.
    710 !Config         Cela peut am�liorer la bande passante des transferts MPI
     742!Config         en m�moire globale a l'aide de la fonction MPI_ALLOC.
     743!Config         Cela peut am�liorer la bande passante des transferts MPI
    711744!Config         d'un facteur 2 
    712745      use_mpi_alloc=.FALSE.
    … …  
    716749!Config  Desc = taille des blocs openmp
    717750!Config  Def  = 1
    718 !Config  Help = defini la taille des packets d'it�ration openmp
    719 !Config         distribu�e � chaque t�che lors de l'entr�e dans une
    720 !Config         boucle parall�lis�e
     751!Config  Help = defini la taille des packets d'it�ration openmp
     752!Config         distribu�e � chaque t�che lors de l'entr�e dans une
     753!Config         boucle parall�lis�e
    721754 
    722755      omp_chunk=1
    … …  
    726759!Config  Desc = activation de la version strato
    727760!Config  Def  = .FALSE.
    728 !Config  Help = active la version stratosph�rique de LMDZ de F. Lott
     761!Config  Help = active la version stratosph�rique de LMDZ de F. Lott
    729762
    730763      ok_strato=.FALSE.
    … …  
    741774      write(lunout,*)' #########################################'
    742775      write(lunout,*)' Configuration des parametres du gcm: '
     776      write(lunout,*)' planet_type = ', planet_type
    743777      write(lunout,*)' dayref = ', dayref
    744778      write(lunout,*)' anneeref = ', anneeref
    … …  
    749783      write(lunout,*)' iecri = ', iecri
    750784      write(lunout,*)' periodav = ', periodav
     785      write(lunout,*)' output_grads_dyn = ', output_grads_dyn
    751786      write(lunout,*)' idissip = ', idissip
    752787      write(lunout,*)' lstardis = ', lstardis
    … …  
    759794      write(lunout,*)' coefdis = ', coefdis
    760795      write(lunout,*)' purmats = ', purmats
     796      write(lunout,*)' read_start = ', read_start
    761797      write(lunout,*)' iflag_phys = ', iflag_phys
    762798      write(lunout,*)' clon = ', clon
    … …  
    764800      write(lunout,*)' grossismx = ', grossismx
    765801      write(lunout,*)' grossismy = ', grossismy
    766       write(lunout,*)' fxyhypbb = ', fxyhypbb
     802      write(lunout,*)' fxyhypb = ', fxyhypb
    767803      write(lunout,*)' dzoomx = ', dzoomx
    768804      write(lunout,*)' dzoomy = ', dzoomy

  • LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/control.h

    r985 r1140  
    1414     &              iperiod,iapp_tracvl,iconser,iecri,idissip,iphysiq , &
    1515     &              periodav,iecrimoy,dayref,anneeref,                  &
    16      &              raz_date,offline,ip_ebil_dyn,config_inca
     16     &              raz_date,offline,ip_ebil_dyn,config_inca,           &
     17     &              planet_type,output_grads_dyn
    1718
    1819      INTEGER   nday,day_step,iperiod,iapp_tracvl,iconser,iecri,        &
    … …  
    2122      REAL periodav
    2223      logical offline
    23       CHARACTER*4 config_inca
     24      CHARACTER (len=4) :: config_inca
     25      CHARACTER(len=10) :: planet_type ! planet type ('earth','mars',...)
     26      LOGICAL :: output_grads_dyn ! output dynamics diagnostics in
     27                                  ! binary grads file 'dyn.dat' (y/n)
    2428!-----------------------------------------------------------------------

  • LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/diagedyn.F

    r774 r1140  
    5858#include "paramet.h"
    5959#include "comgeom.h"
    60 
    61 #ifdef CPP_PHYS
     60#include "iniprint.h"
     61
     62#ifdef CPP_EARTH
    6263#include "../phylmd/YOMCST.h"
    6364#include "../phylmd/YOETHF.h"
    … …  
    139140
    140141
    141 #ifdef CPP_PHYS
     142#ifdef CPP_EARTH
    142143c======================================================================
    143144C     Compute Kinetic enrgy
    … …  
    314315C
    315316#else
    316       print*,'Pour l instant diagedyn a besoin de la physique'
     317      write(lunout,*),'diagedyn: Needs Earth physics to function'
    317318#endif
     319! #endif of #ifdef CPP_EARTH
    318320      RETURN
    319321      END

  • LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/etat0_netcdf.F

    r1117 r1140  
    55c
    66      SUBROUTINE etat0_netcdf (interbar, masque)
    7    
     7#ifdef CPP_EARTH       
    88      USE startvar
    99      USE ioipsl
    … …  
    1414      USE filtreg_mod
    1515      USE infotrac
     16#endif
     17!#endif of #ifdef CPP_EARTH
    1618      !
    1719      IMPLICIT NONE
    … …  
    2527!     .KLON=KFDIA-KIDIA+1,KLEV=llm
    2628      !
     29#ifdef CPP_EARTH   
    2730#include "comgeom2.h"
    2831#include "comvert.h"
    … …  
    3134#include "dimsoil.h"
    3235#include "temps.h"
    33       !
     36#endif
     37!#endif of #ifdef CPP_EARTH
     38      ! arguments:
    3439      LOGICAL interbar
     40      REAL :: masque(iip1,jjp1)
     41
     42#ifdef CPP_EARTH
     43      ! local variables:
    3544      REAL :: latfi(klon), lonfi(klon)
    36       REAL :: orog(iip1,jjp1), rugo(iip1,jjp1), masque(iip1,jjp1),
     45      REAL :: orog(iip1,jjp1), rugo(iip1,jjp1),
    3746     . psol(iip1, jjp1), phis(iip1, jjp1)
    3847      REAL :: p3d(iip1, jjp1, llm+1)
    … …  
    144153      !
    145154      preff     = 101325.
     155      pa        =  50000.
    146156      unskap = 1./kappa
    147157      !
    … …  
    167177      print*,'dtvr',dtvr
    168178
    169       CALL inicons0()
     179      CALL iniconst()
    170180      CALL inigeom()
    171181      !
    … …  
    754764      print*,'entree histclo'
    755765      CALL histclo
     766
     767#endif
     768!#endif of #ifdef CPP_EARTH
    756769      RETURN
    757770      !

  • LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/gcm.F

    r1129 r1140  
    99      USE IOIPSL
    1010#endif
     11
    1112      USE mod_const_mpi, ONLY: init_const_mpi
    1213      USE parallel
    1314      USE mod_phys_lmdz_para, ONLY : klon_mpi_para_nb
     15      USE infotrac
     16      USE mod_interface_dyn_phys
     17      USE mod_hallo
     18      USE Bands
     19
     20      USE filtreg_mod
     21
     22! Ehouarn: for now these only apply to Earth:
     23#ifdef CPP_EARTH
    1424      USE mod_grid_phy_lmdz
    1525      USE mod_phys_lmdz_omp_data, ONLY: klon_omp
    1626      USE dimphy
    17       USE infotrac
    18       USE mod_interface_dyn_phys
    1927      USE comgeomphy
    20       USE mod_hallo
    21       USE Bands
    22 
    23       USE filtreg_mod
    24 
     28#endif
    2529      IMPLICIT NONE
    2630
    … …  
    160164      dynhistave_file = 'dyn_hist_ave'
    161165
    162 c--------------------------------------------------------------------------
    163 c   Iflag_phys controle l'appel a la physique :
    164 c   -------------------------------------------
    165 c      0 : pas de physique
    166 c      1 : Normale (appel a phylmd, phymars ...)
    167 c      2 : rappel Newtonien pour la temperature + friction au sol
    168       iflag_phys=1
    169 
    170 c--------------------------------------------------------------------------
    171 c   Lecture de l'etat initial :
    172 c   ---------------------------
    173 c     T : on lit start.nc
    174 c     F : le modele s'autoinitialise avec un cas academique (iniacademic)
    175 #ifdef CPP_IOIPSL
    176       read_start=.true.
    177 #else
    178       read_start=.false.
    179 #endif
    180166
    181167c-----------------------------------------------------------------------
    … …  
    194180c  ---------------------------------------
    195181c
    196 #ifdef CPP_IOIPSL
     182! Ehouarn: dump possibility of using defrun
     183!#ifdef CPP_IOIPSL
    197184      CALL conf_gcm( 99, .TRUE. , clesphy0 )
    198 #else
    199       CALL defrun( 99, .TRUE. , clesphy0 )
    200 #endif
     185!#else
     186!      CALL defrun( 99, .TRUE. , clesphy0 )
     187!#endif
    201188c
    202189c
    … …  
    208195      call init_parallel
    209196      call Read_Distrib
    210       CALL Init_Phys_lmdz(iim,jjp1,llm,mpi_size,distrib_phys)
     197! Ehouarn : temporarily (?) keep this only for Earth
     198      if (planet_type.eq."earth") then
     199#ifdef CPP_EARTH
     200        CALL Init_Phys_lmdz(iim,jjp1,llm,mpi_size,distrib_phys)
     201#endif
     202      endif ! of if (planet_type.eq."earth")
    211203      CALL set_bands
    212204      CALL Init_interface_dyn_phys
    … …  
    220212c$OMP END PARALLEL
    221213
     214! Ehouarn : temporarily (?) keep this only for Earth
     215      if (planet_type.eq."earth") then
     216#ifdef CPP_EARTH
    222217c$OMP PARALLEL
    223218      call InitComgeomphy
    224219c$OMP END PARALLEL
     220#endif
     221      endif ! of if (planet_type.eq."earth")
    225222
    226223      IF (config_inca /= 'none') THEN
    … …  
    252249c  lecture du fichier start.nc
    253250      if (read_start) then
    254 #ifdef CPP_IOIPSL
     251      ! we still need to run iniacademic to initialize some
     252      ! constants & fields, if we run the 'newtonian' case:
     253        if (iflag_phys.eq.2) then
     254          CALL iniacademic(vcov,ucov,teta,q,masse,ps,phis,time_0)
     255        endif
     256!#ifdef CPP_IOIPSL
     257        if (planet_type.eq."earth") then
     258#ifdef CPP_EARTH
     259! Load an Earth-format start file
    255260         CALL dynetat0("start.nc",vcov,ucov,
    256261     .              teta,q,masse,ps,phis, time_0)
     262#endif
     263        endif ! of if (planet_type.eq."earth")
    257264c       write(73,*) 'ucov',ucov
    258265c       write(74,*) 'vcov',vcov
    … …  
    261268c       write(77,*) 'q',q
    262269
    263 #endif
    264       endif
     270      endif ! of if (read_start)
    265271
    266272c le cas echeant, creation d un etat initial
    267273      IF (prt_level > 9) WRITE(lunout,*)
    268      .                 'AVANT iniacademic AVANT AVANT AVANT AVANT'
     274     .              'GCM: AVANT iniacademic AVANT AVANT AVANT AVANT'
    269275      if (.not.read_start) then
    270276         CALL iniacademic(vcov,ucov,teta,q,masse,ps,phis,time_0)
    … …  
    351357c   Initialisation de la physique :
    352358c   -------------------------------
    353 #ifdef CPP_PHYS
    354359      IF (call_iniphys.and.iflag_phys.eq.1) THEN
    355360         latfi(1)=rlatu(1)
    … …  
    372377
    373378         WRITE(lunout,*)
    374      .           'WARNING!!! vitesse verticale nulle dans la physique'
    375 
     379     .       'GCM: WARNING!!! vitesse verticale nulle dans la physique'
     380! Earth:
     381         if (planet_type.eq."earth") then
     382#ifdef CPP_EARTH
    376383         CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys ,
    377384     ,                latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp     )
    378 
     385#endif
     386         endif ! of if (planet_type.eq."earth")
    379387         call_iniphys=.false.
    380 
    381       ENDIF
    382 #endif
     388      ENDIF ! of IF (call_iniphys.and.(iflag_phys.eq.1))
    383389
    384390
    … …  
    402408      day_end = day_ini + nday
    403409      WRITE(lunout,300)day_ini,day_end
     410 300  FORMAT('1'/,15x,'run du jour',i7,2x,'au jour',i7//)
     411
     412!#ifdef CPP_IOIPSL
     413      if (planet_type.eq."earth") then
     414#ifdef CPP_EARTH
     415      CALL dynredem0_p("restart.nc", day_end, phis)
     416#endif
     417      endif
     418
     419      ecripar = .TRUE.
    404420
    405421#ifdef CPP_IOIPSL
    406       CALL dynredem0_p("restart.nc", day_end, phis)
    407 
    408       ecripar = .TRUE.
    409 
    410422      if ( 1.eq.1) then
    411423      time_step = zdtvr
    … …  
    425437
    426438#endif
     439! #endif of #ifdef CPP_IOIPSL
    427440
    428441c  Choix des frequences de stokage pour le offline
    … …  
    450463
    451464
    452  300  FORMAT('1'/,15x,'run du pas',i7,2x,'au pas',i7,2x,
    453      . 'c''est a dire du jour',i7,3x,'au jour',i7//)
    454465      END
    455466

  • LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/iniacademic.F

    r1114 r1140  
    4545#include "temps.h"
    4646#include "control.h"
     47#include "iniprint.h"
    4748
    4849c   Arguments:
    … …  
    5758      REAL ps(ip1jmp1)                       ! pression  au sol
    5859      REAL masse(ip1jmp1,llm)                ! masse d'air
     60      REAL phis(ip1jmp1)                     ! geopotentiel au sol
     61
     62c   Local:
     63c   ------
     64
    5965      REAL p (ip1jmp1,llmp1  )               ! pression aux interfac.des couches
    6066      REAL pks(ip1jmp1)                      ! exner au  sol
    6167      REAL pk(ip1jmp1,llm)                   ! exner au milieu des couches
    6268      REAL pkf(ip1jmp1,llm)                  ! exner filt.au milieu des couches
    63       REAL phis(ip1jmp1)                     ! geopotentiel au sol
    6469      REAL phi(ip1jmp1,llm)                  ! geopotentiel
    65 
    66 
    67 
    68 
    69 
    70 c   Local:
    71 c   ------
    72 
    7370      REAL ddsin,tetarappelj,tetarappell,zsig
    7471      real tetajl(jjp1,llm)
    … …  
    8178
    8279c-----------------------------------------------------------------------
     80! 1. Initializations for Earth-like case
     81! --------------------------------------
     82      if (planet_type=="earth") then
     83c
     84        time_0=0.
    8385
    84 c
    85       time_0=0.
     86        im         = iim
     87        jm         = jjm
     88        day_ini    = 0
     89        omeg       = 4.*asin(1.)/86400.
     90        rad    = 6371229.
     91        g      = 9.8
     92        daysec = 86400.
     93        dtvr    = daysec/FLOAT(day_step)
     94        zdtvr=dtvr
     95        kappa  = 0.2857143
     96        cpp    = 1004.70885
     97        preff     = 101325.
     98        pa        =  50000.
     99        etot0      = 0.
     100        ptot0      = 0.
     101        ztot0      = 0.
     102        stot0      = 0.
     103        ang0       = 0.
    86104
    87       im         = iim
    88       jm         = jjm
    89       day_ini    = 0
    90       omeg       = 4.*asin(1.)/86400.
    91       rad    = 6371229.
    92       g      = 9.8
    93       daysec = 86400.
    94       dtvr    = daysec/FLOAT(day_step)
    95       zdtvr=dtvr
    96       kappa  = 0.2857143
    97       cpp    = 1004.70885
    98       preff     = 101325.
    99       pa        =  50 000.
    100       etot0      = 0.
    101       ptot0      = 0.
    102       ztot0      = 0.
    103       stot0      = 0.
    104       ang0       = 0.
    105       pa         = 0.
     105        CALL iniconst
     106        CALL inigeom
     107        CALL inifilr
    106108
    107       CALL inicons0
    108       CALL inigeom
    109       CALL inifilr
    110 
    111       ps=0.
    112       phis=0.
     109        ps=0.
     110        phis=0.
    113111c---------------------------------------------------------------------
    114112
    115       taurappel=10.*daysec
     113        taurappel=10.*daysec
    116114
    117115c---------------------------------------------------------------------
    … …  
    119117c   --------------------------------------
    120118
    121       DO l=1,llm
    122        zsig=ap(l)/preff+bp(l)
    123        if (zsig.gt.0.3) then
    124          lsup=l
    125          tetarappell=1./8.*(-log(zsig)-.5)
    126          DO j=1,jjp1
    127             ddsin=sin(rlatu(j))-sin(pi/20.)
    128             tetajl(j,l)=300.*(1+1./18.*(1.-3.*ddsin*ddsin)+tetarappell)
    129          ENDDO
    130         else
     119        DO l=1,llm
     120         zsig=ap(l)/preff+bp(l)
     121         if (zsig.gt.0.3) then
     122           lsup=l
     123           tetarappell=1./8.*(-log(zsig)-.5)
     124           DO j=1,jjp1
     125             ddsin=sin(rlatu(j))-sin(pi/20.)
     126             tetajl(j,l)=300.*(1+1./18.*(1.-3.*ddsin*ddsin)+tetarappell)
     127           ENDDO
     128          else
    131129c   Choix isotherme au-dessus de 300 mbar
    132          do j=1,jjp1
    133             tetajl(j,l)=tetajl(j,lsup)*(0.3/zsig)**kappa
    134          enddo
    135         endif
    136       ENDDO
     130           do j=1,jjp1
     131             tetajl(j,l)=tetajl(j,lsup)*(0.3/zsig)**kappa
     132           enddo
     133          endif ! of if (zsig.gt.0.3)
     134        ENDDO ! of DO l=1,llm
    137135
    138       do l=1,llm
    139          do j=1,jjp1
    140             do i=1,iip1
    141                ij=(j-1)*iip1+i
    142                tetarappel(ij,l)=tetajl(j,l)
    143             enddo
    144          enddo
    145       enddo
     136        do l=1,llm
     137           do j=1,jjp1
     138              do i=1,iip1
     139                 ij=(j-1)*iip1+i
     140                 tetarappel(ij,l)=tetajl(j,l)
     141              enddo
     142           enddo
     143        enddo
    146144
    147 c     call dump2d(jjp1,llm,tetajl,'TEQ   ')
     145c       call dump2d(jjp1,llm,tetajl,'TEQ   ')
    148146
    149       ps=1.e5
    150       phis=0.
    151       CALL pression ( ip1jmp1, ap, bp, ps, p       )
    152       CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
    153       CALL massdair(p,masse)
     147        ps=1.e5
     148        phis=0.
     149        CALL pression ( ip1jmp1, ap, bp, ps, p       )
     150        CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
     151        CALL massdair(p,masse)
    154152
    155153c  intialisation du vent et de la temperature
    156       teta(:,:)=tetarappel(:,:)
    157       CALL geopot(ip1jmp1,teta,pk,pks,phis,phi)
    158       call ugeostr(phi,ucov)
    159       vcov=0.
    160       q(:,:,1   )=1.e-10
    161       q(:,:,2   )=1.e-15
    162       q(:,:,3:nqtot)=0.
     154        teta(:,:)=tetarappel(:,:)
     155        CALL geopot(ip1jmp1,teta,pk,pks,phis,phi)
     156        call ugeostr(phi,ucov)
     157        vcov=0.
     158        q(:,:,1   )=1.e-10
     159        q(:,:,2   )=1.e-15
     160        q(:,:,3:nqtot)=0.
    163161
    164162
    165 c   perturbation al\351atoire sur la temp\351rature
    166       idum  = -1
    167       zz = ran1(idum)
    168       idum  = 0
    169       do l=1,llm
    170          do ij=iip2,ip1jm
    171             teta(ij,l)=teta(ij,l)*(1.+0.005*ran1(idum))
    172          enddo
    173       enddo
     163c   perturbation aleatoire sur la temperature
     164        idum  = -1
     165        zz = ran1(idum)
     166        idum  = 0
     167        do l=1,llm
     168           do ij=iip2,ip1jm
     169              teta(ij,l)=teta(ij,l)*(1.+0.005*ran1(idum))
     170           enddo
     171        enddo
    174172
    175       do l=1,llm
    176          do ij=1,ip1jmp1,iip1
    177             teta(ij+iim,l)=teta(ij,l)
    178          enddo
    179       enddo
     173        do l=1,llm
     174           do ij=1,ip1jmp1,iip1
     175              teta(ij+iim,l)=teta(ij,l)
     176           enddo
     177        enddo
    180178
    181179
    … …  
    187185
    188186c   initialisation d'un traceur sur une colonne
    189       j=jjp1*3/4
    190       i=iip1/2
    191       ij=(j-1)*iip1+i
    192       q(ij,:,3)=1.
    193 
     187        j=jjp1*3/4
     188        i=iip1/2
     189        ij=(j-1)*iip1+i
     190        q(ij,:,3)=1.
     191     
     192      else
     193        write(lunout,*)"iniacademic: planet types other than earth",
     194     &                 " not implemented (yet)."
     195        stop
     196      endif ! of if (planet_type=="earth")
    194197      return
    195198      END

  • LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/leapfrog_p.F

    r1134 r1140  
    44c
    55c
    6 #define IO_DEBUG
    7 
    8 #undef CPP_IOIPSL
    9 #define CPP_IOIPSL
    106
    117      SUBROUTINE leapfrog_p(ucov,vcov,teta,ps,masse,phis,q,clesphy0,
    … …  
    285281c$OMP BARRIER
    286282       else
    287          
     283! Save fields obtained at previous time step as '...m1'
    288284         ijb=ij_begin
    289285         ije=ij_end
    … …  
    312308     .                    llm, -2,2, .TRUE., 1 )
    313309
    314        endif
     310       endif ! of if (FirstCaldyn)
    315311       
    316312      forward = .TRUE.
    … …  
    356352         IF( MOD(itau,idissip ).EQ.0.AND..NOT.forward ) apdiss = .TRUE.
    357353         IF( MOD(itau,iphysiq ).EQ.0.AND..NOT.forward
    358      s          .and. iflag_phys.NE.0                 ) apphys = .TRUE.
     354     s          .and. iflag_phys.EQ.1                 ) apphys = .TRUE.
    359355      ELSE
    360356         IF( MOD(itau   ,iconser) .EQ. 0              ) conser = .TRUE.
    361357         IF( MOD(itau+1,idissip)  .EQ. 0              ) apdiss = .TRUE.
    362          IF( MOD(itau+1,iphysiq).EQ.0.AND.iflag_phys.NE.0) apphys=.TRUE.
     358         IF( MOD(itau+1,iphysiq).EQ.0.AND.iflag_phys.EQ.1) apphys=.TRUE.
    363359      END IF
    364360
    … …  
    537533
    538534c$OMP MASTER
    539       IF (prt_level>9)WRITE(lunout,*)"Iteration No",True_itau
     535      IF (prt_level>9) THEN
     536        WRITE(lunout,*)"leapfrog_p: Iteration No",True_itau
     537      ENDIF
    540538
    541539
    … …  
    594592
    595593
    596          ENDIF
    597 c
    598       ENDIF
     594         ENDIF ! of IF (offline)
     595c
     596      ENDIF ! of IF( forward. OR . leapf )
    599597cc$OMP END PARALLEL
    600598
    … …  
    673671c$OMP MASTER
    674672         call suspend_timer(timer_caldyn)
    675          WRITE(lunout,*)'Entree dans la physique : Iteration No ',      &
    676      &                   true_itau
     673
     674         write(lunout,*)
     675     &   'leapfrog_p: Entree dans la physique : Iteration No ',true_itau
    677676c$OMP END MASTER
    678677
    … …  
    693692c   -----------------------------------------------------
    694693
    695 #ifdef CPP_PHYS
    696694c+jld
    697695
    … …  
    699697      IF (ip_ebil_dyn.ge.1 ) THEN
    700698          ztit='bil dyn'
    701           CALL diagedyn(ztit,2,1,1,dtphys
    702      e  , ucov    , vcov , ps, p ,pk , teta , q(:,:,1), q(:,:,2))
     699! Ehouarn: be careful, diagedyn is Earth-specific (includes ../phylmd/..)!
     700           IF (planet_type.eq."earth") THEN
     701            CALL diagedyn(ztit,2,1,1,dtphys
     702     &    , ucov    , vcov , ps, p ,pk , teta , q(:,:,1), q(:,:,2))
     703           ENDIF
    703704      ENDIF
    704705c-jld
    … …  
    787788          dpfi_tmp(1:iip1)     = dpfi(ijb:ijb+iim) 
    788789c$OMP END MASTER
    789         endif
     790        endif ! of if ( .not. pole_nord)
    790791
    791792c$OMP BARRIER
    … …  
    843844c$OMP END MASTER
    844845         
    845         endif
     846        endif ! of if (.not. pole_nord)
    846847c$OMP BARRIER
    847848cc$OMP MASTER       
    … …  
    932933cc$OMP END MASTER
    933934
    934 #else
    935 
     935
     936c-jld
     937c$OMP MASTER
     938         call resume_timer(timer_caldyn)
     939         if (FirstPhysic) then
     940           ok_start_timer=.TRUE.
     941           FirstPhysic=.false.
     942         endif
     943c$OMP END MASTER
     944       ENDIF ! of IF( apphys )
     945
     946      IF(iflag_phys.EQ.2) THEN ! "Newtonian" case
    936947c   Calcul academique de la physique = Rappel Newtonien + fritcion
    937948c   --------------------------------------------------------------
    … …  
    949960
    950961       call friction_p(ucov,vcov,iphysiq*dtvr)
    951 
    952 #endif
    953 
    954 c-jld
    955 c$OMP MASTER
    956          call resume_timer(timer_caldyn)
    957          if (FirstPhysic) then
    958            ok_start_timer=.TRUE.
    959            FirstPhysic=.false.
    960          endif
    961 c$OMP END MASTER
    962        ENDIF
     962      ENDIF ! of IF(iflag_phys.EQ.2)
     963
    963964
    964965        CALL pression_p ( ip1jmp1, ap, bp, ps, p                  )
    … …  
    13121313c      IF( MOD(itau,iecri         ).EQ.0) THEN
    13131314
    1314            IF( MOD(itau,iecri*day_step).EQ.0) THEN
    1315 c$OMP BARRIER
    1316 c$OMP MASTER
    1317                nbetat = nbetatdem
    1318         CALL geopot_p  ( ip1jmp1, teta  , pk , pks,  phis  , phi )
     1315            IF( MOD(itau,iecri*day_step).EQ.0) THEN
     1316c$OMP BARRIER
     1317c$OMP MASTER
     1318              nbetat = nbetatdem
     1319              CALL geopot_p(ip1jmp1,teta,pk,pks,phis,phi)
    13191320       
    13201321cym        unat=0.
    13211322       
    1322         ijb=ij_begin
    1323         ije=ij_end
    1324        
    1325         if (pole_nord) then
    1326           ijb=ij_begin+iip1
    1327           unat(1:iip1,:)=0.
    1328         endif
    1329        
    1330         if (pole_sud) then
    1331           ije=ij_end-iip1
    1332           unat(ij_end-iip1+1:ij_end,:)=0.
    1333         endif
     1323              ijb=ij_begin
     1324              ije=ij_end
     1325       
     1326              if (pole_nord) then
     1327                ijb=ij_begin+iip1
     1328                unat(1:iip1,:)=0.
     1329              endif
     1330       
     1331              if (pole_sud) then
     1332                ije=ij_end-iip1
     1333                unat(ij_end-iip1+1:ij_end,:)=0.
     1334              endif
    13341335           
    1335         do l=1,llm
    1336            unat(ijb:ije,l)=ucov(ijb:ije,l)/cu(ijb:ije)
    1337         enddo
    1338 
    1339         ijb=ij_begin
    1340         ije=ij_end
    1341         if (pole_sud) ije=ij_end-iip1
    1342        
    1343         do l=1,llm
    1344            vnat(ijb:ije,l)=vcov(ijb:ije,l)/cv(ijb:ije)
    1345         enddo
     1336              do l=1,llm
     1337                unat(ijb:ije,l)=ucov(ijb:ije,l)/cu(ijb:ije)
     1338              enddo
     1339
     1340              ijb=ij_begin
     1341              ije=ij_end
     1342              if (pole_sud) ije=ij_end-iip1
     1343       
     1344              do l=1,llm
     1345                vnat(ijb:ije,l)=vcov(ijb:ije,l)/cv(ijb:ije)
     1346              enddo
    13461347       
    13471348#ifdef CPP_IOIPSL
    13481349 
    1349         CALL writehist_p(histid,histvid, itau,vcov,
    1350      s                       ucov,teta,phi,q,masse,ps,phis)
     1350              CALL writehist_p(histid,histvid, itau,vcov,
     1351     &                         ucov,teta,phi,q,masse,ps,phis)
    13511352
    13521353#endif
    1353 c$OMP END MASTER
    1354            ENDIF
     1354! For some Grads outputs of fields
     1355              if (output_grads_dyn) then
     1356! Ehouarn: hope this works the way I think it does:
     1357                  call Gather_Field(unat,ip1jmp1,llm,0)
     1358                  call Gather_Field(vnat,ip1jm,llm,0)
     1359                  call Gather_Field(teta,ip1jmp1,llm,0)
     1360                  call Gather_Field(ps,ip1jmp1,1,0)
     1361                  do iq=1,nqtot
     1362                    call Gather_Field(q(1,1,iq),ip1jmp1,llm,0)
     1363                  enddo
     1364                  if (mpi_rank==0) then
     1365#include "write_grads_dyn.h"
     1366                  endif
     1367              endif ! of if (output_grads_dyn)
     1368c$OMP END MASTER
     1369            ENDIF ! of IF(MOD(itau,iecri).EQ.0)
    13551370
    13561371            IF(itau.EQ.itaufin) THEN
    … …  
    13591374c$OMP MASTER
    13601375
    1361 c#ifdef CPP_IOIPSL
    1362 
    1363        CALL dynredem1_p("restart.nc",0.0,
    1364      ,                     vcov,ucov,teta,q,masse,ps)
    1365 c#endif
     1376              if (planet_type.eq."earth") then
     1377#ifdef CPP_EARTH
     1378! Write an Earth-format restart file
     1379                CALL dynredem1_p("restart.nc",0.0,
     1380     &                           vcov,ucov,teta,q,masse,ps)
     1381
     1382#endif
     1383              endif ! of if (planet_type.eq."earth")
    13661384
    13671385              CLOSE(99)
    13681386c$OMP END MASTER
    1369             ENDIF
     1387            ENDIF ! of IF (itau.EQ.itaufin)
    13701388
    13711389c-----------------------------------------------------------------------
    … …  
    13981416                 dt  = 2.*dtvr
    13991417                 GO TO 2
    1400             END IF
    1401 
    1402       ELSE
     1418            END IF ! of IF (MOD(itau,iperiod).EQ.0)
     1419                   !    ELSEIF (MOD(itau-1,iperiod).EQ.0)
     1420
     1421
     1422      ELSE ! of IF (.not.purmats)
    14031423
    14041424c       ........................................................
    … …  
    14271447               GO TO 2
    14281448
    1429             ELSE
    1430 
    1431             IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN
     1449            ELSE ! of IF(forward)
     1450
     1451              IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN
    14321452               IF(itau.EQ.itaufin) THEN
    14331453                  iav=1
    … …  
    14381458c$OMP BARRIER
    14391459
    1440               call Register_Hallo(vcov,ip1jm,llm,1,0,0,1,TestRequest)
    1441               call SendRequest(TestRequest)
    1442 c$OMP BARRIER
    1443               call WaitRequest(TestRequest)
    1444 
    1445 c$OMP BARRIER
    1446 c$OMP MASTER
    1447               CALL writedynav_p(histaveid, itau,vcov ,
     1460               call Register_Hallo(vcov,ip1jm,llm,1,0,0,1,TestRequest)
     1461               call SendRequest(TestRequest)
     1462c$OMP BARRIER
     1463               call WaitRequest(TestRequest)
     1464
     1465c$OMP BARRIER
     1466c$OMP MASTER
     1467               CALL writedynav_p(histaveid, itau,vcov ,
    14481468     ,                          ucov,teta,pk,phi,q,masse,ps,phis)
    14491469               call bilan_dyn_p (2,dtvr*iperiod,dtvr*day_step*periodav,
    … …  
    14511471c$OMP END MASTER
    14521472#endif
    1453             ENDIF
     1473              ENDIF ! of IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin)
     1474
    14541475
    14551476c               IF(MOD(itau,iecri         ).EQ.0) THEN
    … …  
    14571478c$OMP BARRIER
    14581479c$OMP MASTER
    1459                   nbetat = nbetatdem
    1460        CALL geopot_p( ip1jmp1, teta  , pk , pks,  phis  , phi   )
     1480                nbetat = nbetatdem
     1481                CALL geopot_p(ip1jmp1,teta,pk,pks,phis,phi)
    14611482
    14621483cym        unat=0.
    1463         ijb=ij_begin
    1464         ije=ij_end
    1465        
    1466         if (pole_nord) then
    1467           ijb=ij_begin+iip1
    1468           unat(1:iip1,:)=0.
    1469         endif
    1470        
    1471         if (pole_sud) then
    1472           ije=ij_end-iip1
    1473           unat(ij_end-iip1+1:ij_end,:)=0.
    1474         endif
     1484                ijb=ij_begin
     1485                ije=ij_end
     1486       
     1487                if (pole_nord) then
     1488                  ijb=ij_begin+iip1
     1489                  unat(1:iip1,:)=0.
     1490                endif
     1491       
     1492                if (pole_sud) then
     1493                  ije=ij_end-iip1
     1494                  unat(ij_end-iip1+1:ij_end,:)=0.
     1495                endif
    14751496           
    1476         do l=1,llm
    1477            unat(ijb:ije,l)=ucov(ijb:ije,l)/cu(ijb:ije)
    1478         enddo
    1479 
    1480         ijb=ij_begin
    1481         ije=ij_end
    1482         if (pole_sud) ije=ij_end-iip1
    1483        
    1484         do l=1,llm
    1485            vnat(ijb:ije,l)=vcov(ijb:ije,l)/cv(ijb:ije)
    1486         enddo
     1497                do l=1,llm
     1498                  unat(ijb:ije,l)=ucov(ijb:ije,l)/cu(ijb:ije)
     1499                enddo
     1500
     1501                ijb=ij_begin
     1502                ije=ij_end
     1503                if (pole_sud) ije=ij_end-iip1
     1504       
     1505                do l=1,llm
     1506                  vnat(ijb:ije,l)=vcov(ijb:ije,l)/cv(ijb:ije)
     1507                enddo
    14871508
    14881509#ifdef CPP_IOIPSL
    14891510
    1490        CALL writehist_p( histid, histvid, itau,vcov ,
    1491      ,                           ucov,teta,phi,q,masse,ps,phis)
    1492 c#else
    1493 c      call Gather_Field(unat,ip1jmp1,llm,0)
    1494 c      call Gather_Field(vnat,ip1jm,llm,0)
    1495 c      call Gather_Field(teta,ip1jmp1,llm,0)
    1496 c      call Gather_Field(ps,ip1jmp1,1,0)
    1497 c      do iq=1,nqtot
    1498 c        call Gather_Field(q(1,1,iq),ip1jmp1,llm,0)
    1499 c      enddo
     1511                CALL writehist_p(histid, histvid, itau,vcov ,
     1512     &                           ucov,teta,phi,q,masse,ps,phis)
     1513#endif
     1514! For some Grads output (but does it work?)
     1515                if (output_grads_dyn) then
     1516                  call Gather_Field(unat,ip1jmp1,llm,0)
     1517                  call Gather_Field(vnat,ip1jm,llm,0)
     1518                  call Gather_Field(teta,ip1jmp1,llm,0)
     1519                  call Gather_Field(ps,ip1jmp1,1,0)
     1520                  do iq=1,nqtot
     1521                    call Gather_Field(q(1,1,iq),ip1jmp1,llm,0)
     1522                  enddo
    15001523c     
    1501 c      if (mpi_rank==0) then
    1502 c#include "write_grads_dyn.h"
    1503 c      endif
    1504 #endif
    1505 
    1506 c$OMP END MASTER
    1507                ENDIF
    1508 
    1509                  IF(itau.EQ.itaufin) THEN
     1524                  if (mpi_rank==0) then
     1525#include "write_grads_dyn.h"
     1526                  endif
     1527                endif ! of if (output_grads_dyn)
     1528
     1529c$OMP END MASTER
     1530              ENDIF ! of IF(MOD(itau,iecri*day_step).EQ.0)
     1531
     1532              IF(itau.EQ.itaufin) THEN
     1533                if (planet_type.eq."earth") then
     1534#ifdef CPP_EARTH
    15101535c$OMP MASTER
    15111536                   CALL dynredem1_p("restart.nc",0.0,
    15121537     .                               vcov,ucov,teta,q,masse,ps)
    15131538c$OMP END MASTER
    1514                  ENDIF
    1515                  forward = .TRUE.
    1516                  GO TO  1
    1517 
    1518             ENDIF
    1519 
    1520       END IF
    1521 c$OMP MASTER
    1522         call finalize_parallel
    1523 c$OMP END MASTER
    1524         RETURN
     1539#endif
     1540                endif ! of if (planet_type.eq."earth")
     1541              ENDIF ! of IF(itau.EQ.itaufin)
     1542
     1543              forward = .TRUE.
     1544              GO TO  1
     1545
     1546            ENDIF ! of IF (forward)
     1547
     1548      END IF ! of IF(.not.purmats)
     1549c$OMP MASTER
     1550      call finalize_parallel
     1551c$OMP END MASTER
     1552      RETURN
    15251553      END
Note: See TracChangeset for help on using the changeset viewer.