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/dyn3d
Files:
2 added
3 deleted
8 edited

Legend:

Unmodified
Added
Removed

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

    r524 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/dyn3d/conf_gcm.F

    r1046 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      IMPLICIT NONE
    1015c-----------------------------------------------------------------------
    … …  
    99104c  Parametres de controle du run:
    100105c-----------------------------------------------------------------------
     106!Config  Key  = planet_type
     107!Config  Desc = planet type ("earth", "mars", "venus", ...)
     108!Config  Def  = earth
     109!Config  Help = this flag sets the type of atymosphere that is considered
     110      planet_type="earth"
     111      CALL getin('planet_type',planet_type)
    101112
    102113!Config  Key  = dayref
    … …  
    179190       CALL getin('periodav',periodav)
    180191
     192!Config  Key  = output_grads_dyn
     193!Config  Desc = output dynamics diagnostics in 'dyn.dat' file
     194!Config  Def  = n
     195!Config  Help = output dynamics diagnostics in Grads-readable 'dyn.dat' file
     196       output_grads_dyn=.false.
     197       CALL getin('output_grads_dyn',output_grads_dyn)
     198
    181199!Config  Key  = idissip
    182200!Config  Desc = periode de la dissipation
    … …  
    274292c    ...............................................................
    275293
     294!Config  Key  =  read_start
     295!Config  Desc = Initialize model using a 'start.nc' file
     296!Config  Def  = y
     297!Config  Help = y: intialize dynamical fields using a 'start.nc' file
     298!               n: fields are initialized by 'iniacademic' routine
     299       read_start= .true.
     300       CALL getin('read_start',read_start)
     301
    276302!Config  Key  = iflag_phys
    277303!Config  Desc = Avec ls physique
    … …  
    330356c
    331357      IF( ABS(clat - clatt).GE. 0.001 )  THEN
    332         PRINT *,' La valeur de clat passee par run.def est differente de
    333      * celle lue sur le fichier  start '
     358        write(lunout,*)'conf_gcm: La valeur de clat passee par run.def',
     359     &    ' est differente de celle lue sur le fichier  start '
    334360        STOP
    335361      ENDIF
    … …  
    345371
    346372      IF( ABS(grossismx - grossismxx).GE. 0.001 )  THEN
    347         PRINT *,' La valeur de grossismx passee par run.def est differente 
    348      * de celle lue sur le fichier  start '
     373        write(lunout,*)'conf_gcm: La valeur de grossismx passee par ',
     374     &  'run.def est differente de celle lue sur le fichier  start '
    349375        STOP
    350376      ENDIF
    … …  
    359385
    360386      IF( ABS(grossismy - grossismyy).GE. 0.001 )  THEN
    361         PRINT *,' La valeur de grossismy passee par run.def est differen
    362      * te de celle lue sur le fichier  start '
     387        write(lunout,*)'conf_gcm: La valeur de grossismy passee par ',
     388     & 'run.def est differente de celle lue sur le fichier  start '
    363389        STOP
    364390      ENDIF
    365391     
    366392      IF( grossismx.LT.1. )  THEN
    367         PRINT *,' ***  ATTENTION !! grossismx < 1 .   *** '
     393        write(lunout,*)
     394     &       'conf_gcm: ***  ATTENTION !! grossismx < 1 .   *** '
    368395         STOP
    369396      ELSE
    … …  
    373400
    374401      IF( grossismy.LT.1. )  THEN
    375         PRINT *,' ***  ATTENTION !! grossismy < 1 .   *** '
     402        write(lunout,*)
     403     &       'conf_gcm: ***  ATTENTION !! grossismy < 1 .   *** '
    376404         STOP
    377405      ELSE
    … …  
    379407      ENDIF
    380408
    381       PRINT *,' alphax alphay defrun ',alphax,alphay
     409      write(lunout,*)'conf_gcm: alphax alphay',alphax,alphay
    382410c
    383411c    alphax et alphay sont les anciennes formulat. des grossissements
    … …  
    394422
    395423      IF( .NOT.fxyhypb )  THEN
    396            IF( fxyhypbb )     THEN
    397               PRINT *,' ********  PBS DANS  DEFRUN  ******** '
    398               PRINT *,' *** fxyhypb lu sur le fichier start est F ',
    399      *       'alors  qu il est  T  sur  run.def  ***'
     424         IF( fxyhypbb )     THEN
     425            write(lunout,*)' ********  PBS DANS  CONF_GCM  ******** '
     426            write(lunout,*)' *** fxyhypb lu sur le fichier start est ',
     427     *       'F alors  qu il est  T  sur  run.def  ***'
    400428              STOP
    401            ENDIF
     429         ENDIF
    402430      ELSE
    403            IF( .NOT.fxyhypbb )   THEN
    404               PRINT *,' ********  PBS DANS  DEFRUN  ******** '
    405               PRINT *,' ***  fxyhypb lu sur le fichier start est T ',
    406      *        'alors  qu il est  F  sur  run.def  ****  '
     431         IF( .NOT.fxyhypbb )   THEN
     432            write(lunout,*)' ********  PBS DANS  CONF_GCM  ******** '
     433            write(lunout,*)' ***  fxyhypb lu sur le fichier start est ',
     434     *        'T alors  qu il est  F  sur  run.def  ****  '
    407435              STOP
    408            ENDIF
     436         ENDIF
    409437      ENDIF
    410438c
    … …  
    419447      IF( fxyhypb )  THEN
    420448       IF( ABS(dzoomx - dzoomxx).GE. 0.001 )  THEN
    421         PRINT *,' La valeur de dzoomx passee par run.def est differente
    422      *  de celle lue sur le fichier  start '
     449        write(lunout,*)'conf_gcm: La valeur de dzoomx passee par ',
     450     *  'run.def est differente de celle lue sur le fichier  start '
    423451        STOP
    424452       ENDIF
    … …  
    435463      IF( fxyhypb )  THEN
    436464       IF( ABS(dzoomy - dzoomyy).GE. 0.001 )  THEN
    437         PRINT *,' La valeur de dzoomy passee par run.def est differente
    438      * de celle lue sur le fichier  start '
     465        write(lunout,*)'conf_gcm: La valeur de dzoomy passee par ',
     466     * 'run.def est differente de celle lue sur le fichier  start '
    439467        STOP
    440468       ENDIF
    … …  
    450478      IF( fxyhypb )  THEN
    451479       IF( ABS(taux - tauxx).GE. 0.001 )  THEN
    452         PRINT *,' La valeur de taux passee par run.def est differente
    453      * de celle lue sur le fichier  start '
     480        write(lunout,*)'conf_gcm: La valeur de taux passee par ',
     481     * 'run.def est differente de celle lue sur le fichier  start '
    454482        STOP
    455483       ENDIF
    … …  
    465493      IF( fxyhypb )  THEN
    466494       IF( ABS(tauy - tauyy).GE. 0.001 )  THEN
    467         PRINT *,' La valeur de tauy passee par run.def est differente
    468      * de celle lue sur le fichier  start '
     495        write(lunout,*)'conf_gcm: La valeur de tauy passee par ',
     496     * 'run.def est differente de celle lue sur le fichier  start '
    469497        STOP
    470498       ENDIF
    … …  
    484512
    485513        IF( .NOT.ysinus )  THEN
    486            IF( ysinuss )     THEN
    487               PRINT *,' ********  PBS DANS  DEFRUN  ******** '
    488               PRINT *,' *** ysinus lu sur le fichier start est F ',
    489      *       'alors  qu il est  T  sur  run.def  ***'
     514          IF( ysinuss )     THEN
     515            write(lunout,*)' ********  PBS DANS  CONF_GCM  ******** '
     516            write(lunout,*)' *** ysinus lu sur le fichier start est F',
     517     *       ' alors  qu il est  T  sur  run.def  ***'
     518            STOP
     519          ENDIF
     520        ELSE
     521          IF( .NOT.ysinuss )   THEN
     522            write(lunout,*)' ********  PBS DANS  CONF_GCM  ******** '
     523            write(lunout,*)' *** ysinus lu sur le fichier start est T',
     524     *        ' alors  qu il est  F  sur  run.def  ****  '
    490525              STOP
    491            ENDIF
    492         ELSE
    493            IF( .NOT.ysinuss )   THEN
    494               PRINT *,' ********  PBS DANS  DEFRUN  ******** '
    495               PRINT *,' ***  ysinus lu sur le fichier start est T ',
    496      *        'alors  qu il est  F  sur  run.def  ****  '
    497               STOP
    498            ENDIF
     526          ENDIF
    499527        ENDIF
    500       ENDIF
     528      ENDIF ! of IF( .NOT.fxyhypb  )
    501529c
    502530!Config  Key  = offline
    … …  
    521549      write(lunout,*)' #########################################'
    522550      write(lunout,*)' Configuration des parametres du gcm: '
     551      write(lunout,*)' planet_type = ', planet_type
    523552      write(lunout,*)' dayref = ', dayref
    524553      write(lunout,*)' anneeref = ', anneeref
    … …  
    529558      write(lunout,*)' iecri = ', iecri
    530559      write(lunout,*)' periodav = ', periodav
     560      write(lunout,*)' output_grads_dyn = ', output_grads_dyn
    531561      write(lunout,*)' idissip = ', idissip
    532562      write(lunout,*)' lstardis = ', lstardis
    … …  
    539569      write(lunout,*)' coefdis = ', coefdis
    540570      write(lunout,*)' purmats = ', purmats
     571      write(lunout,*)' read_start = ', read_start
    541572      write(lunout,*)' iflag_phys = ', iflag_phys
    542573      write(lunout,*)' iphysiq = ', iphysiq
    … …  
    590621
    591622      IF( grossismx.LT.1. )  THEN
    592         PRINT *,' ***  ATTENTION !! grossismx < 1 .   *** '
     623        write(lunout,*)
     624     &   'conf_gcm: ***  ATTENTION !! grossismx < 1 .   *** '
    593625         STOP
    594626      ELSE
    … …  
    598630
    599631      IF( grossismy.LT.1. )  THEN
    600         PRINT *,' ***  ATTENTION !! grossismy < 1 .   *** '
     632        write(lunout,*)
     633     &  'conf_gcm: ***  ATTENTION !! grossismy < 1 .   *** '
    601634         STOP
    602635      ELSE
    … …  
    604637      ENDIF
    605638
    606       PRINT *,' alphax alphay defrun ',alphax,alphay
     639      write(lunout,*)'conf_gcm: alphax alphay ',alphax,alphay
    607640c
    608641c    alphax et alphay sont les anciennes formulat. des grossissements
    … …  
    693726      write(lunout,*)' #########################################'
    694727      write(lunout,*)' Configuration des parametres du gcm: '
     728      write(lunout,*)' planet_type = ', planet_type
    695729      write(lunout,*)' dayref = ', dayref
    696730      write(lunout,*)' anneeref = ', anneeref
    … …  
    701735      write(lunout,*)' iecri = ', iecri
    702736      write(lunout,*)' periodav = ', periodav
     737      write(lunout,*)' output_grads_dyn = ', output_grads_dyn
    703738      write(lunout,*)' idissip = ', idissip
    704739      write(lunout,*)' lstardis = ', lstardis
    … …  
    711746      write(lunout,*)' coefdis = ', coefdis
    712747      write(lunout,*)' purmats = ', purmats
     748      write(lunout,*)' read_start = ', read_start
    713749      write(lunout,*)' iflag_phys = ', iflag_phys
    714750      write(lunout,*)' iphysiq = ', iphysiq
    715       write(lunout,*)' clonn = ', clonn
    716       write(lunout,*)' clatt = ', clatt
     751      write(lunout,*)' clon = ', clon
     752      write(lunout,*)' clat = ', clat
    717753      write(lunout,*)' grossismx = ', grossismx
    718754      write(lunout,*)' grossismy = ', grossismy
    719       write(lunout,*)' fxyhypbb = ', fxyhypbb
     755      write(lunout,*)' fxyhypb = ', fxyhypb
    720756      write(lunout,*)' dzoomx = ', dzoomx
    721757      write(lunout,*)' dzoomy = ', dzoomy

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

    r962 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,        &
    … …  
    2223      logical offline
    2324      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/dyn3d/diagedyn.F

    r524 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/dyn3d/etat0_netcdf.F

    r1114 r1140  
    55c
    66      SUBROUTINE etat0_netcdf (interbar, masque)
    7    
     7#ifdef CPP_EARTH       
    88      USE startvar
    99      USE ioipsl
    … …  
    1414      USE phys_state_var_mod
    1515      USE filtreg_mod
     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      REAL :: psol(iip1, jjp1), phis(iip1, jjp1)
    3847      REAL :: p3d(iip1, jjp1, llm+1)
    … …  
    143152      !
    144153      preff     = 101325.
     154      pa        =  50000.
    145155      unskap = 1./kappa
    146156      !
    … …  
    168178
    169179
    170       CALL inicons0()
     180      CALL iniconst()
    171181      CALL inigeom()
    172182
    … …  
    759769      DEALLOCATE(q3d)
    760770
     771#endif
     772!#endif of #ifdef CPP_EARTH
    761773      RETURN
    762774      !

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

    r1114 r1140  
    88#ifdef CPP_IOIPSL
    99      USE IOIPSL
     10#else
     11! if not using IOIPSL, we still need to use (a local version of) getin
     12      USE ioipsl_getincom
    1013#endif
    1114
    … …  
    1720! A nettoyer. On ne veut qu'une ou deux routines d'interface
    1821! dynamique -> physique pour l'initialisation
    19 #ifdef CPP_PHYS
     22! Ehouarn: for now these only apply to Earth:
     23#ifdef CPP_EARTH
    2024      USE dimphy
    2125      USE comgeomphy
    … …  
    158162
    159163
    160 c--------------------------------------------------------------------------
    161 c   Iflag_phys controle l'appel a la physique :
    162 c   -------------------------------------------
    163 c      0 : pas de physique
    164 c      1 : Normale (appel a phylmd, phymars ...)
    165 c      2 : rappel Newtonien pour la temperature + friction au sol
    166       iflag_phys=1
    167 
    168 c--------------------------------------------------------------------------
    169 c   Lecture de l'etat initial :
    170 c   ---------------------------
    171 c     T : on lit start.nc
    172 c     F : le modele s'autoinitialise avec un cas academique (iniacademic)
    173       read_start=.true.
    174 #ifdef CPP_IOIPSL
    175 #else
    176       read_start=.false.
    177 #endif
    178 #ifdef CPP_PHYS
    179 #else
    180       read_start=.false.
    181 #endif
    182164
    183165c-----------------------------------------------------------------------
    … …  
    196178c  ---------------------------------------
    197179c
    198 #ifdef CPP_IOIPSL
     180! Ehouarn: dump possibility of using defrun
     181!#ifdef CPP_IOIPSL
    199182      CALL conf_gcm( 99, .TRUE. , clesphy0 )
    200 #else
    201       CALL defrun( 99, .TRUE. , clesphy0 )
    202 #endif
     183!#else
     184!      CALL defrun( 99, .TRUE. , clesphy0 )
     185!#endif
    203186
    204187!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    … …  
    206189! A nettoyer. On ne veut qu'une ou deux routines d'interface
    207190! dynamique -> physique pour l'initialisation
    208 #ifdef CPP_PHYS
     191! Ehouarn : temporarily (?) keep this only for Earth
     192      if (planet_type.eq."earth") then
     193#ifdef CPP_EARTH
    209194      CALL Init_Phys_lmdz(iim,jjp1,llm,1,(jjm-1)*iim+2)
    210195      call InitComgeomphy
    211196#endif
     197      endif
    212198!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    213199
    … …  
    242228c  lecture du fichier start.nc
    243229      if (read_start) then
    244 #ifdef CPP_IOIPSL
     230      ! we still need to run iniacademic to initialize some
     231      ! constants & fields, if we run the 'newtonian' case:
     232        if (iflag_phys.eq.2) then
     233          CALL iniacademic(vcov,ucov,teta,q,masse,ps,phis,time_0)
     234        endif
     235!#ifdef CPP_IOIPSL
     236        if (planet_type.eq."earth") then
     237#ifdef CPP_EARTH
     238! Load an Earth-format start file
    245239         CALL dynetat0("start.nc",vcov,ucov,
    246240     .              teta,q,masse,ps,phis, time_0)
     241#endif
     242        endif ! of if (planet_type.eq."earth")
    247243c       write(73,*) 'ucov',ucov
    248244c       write(74,*) 'vcov',vcov
    … …  
    251247c       write(77,*) 'q',q
    252248
    253 #endif
    254       endif
     249      endif ! of if (read_start)
    255250
    256251      IF (config_inca /= 'none') THEN
    … …  
    264259c le cas echeant, creation d un etat initial
    265260      IF (prt_level > 9) WRITE(lunout,*)
    266      .                 'AVANT iniacademic AVANT AVANT AVANT AVANT'
     261     .              'GCM: AVANT iniacademic AVANT AVANT AVANT AVANT'
    267262      if (.not.read_start) then
    268263         CALL iniacademic(vcov,ucov,teta,q,masse,ps,phis,time_0)
    … …  
    298293      if (annee_ref .ne. anneeref .or. day_ref .ne. dayref) then
    299294        write(lunout,*)
    300      .  ' Attention les dates initiales lues dans le fichier'
     295     .  'GCM: Attention les dates initiales lues dans le fichier'
    301296        write(lunout,*)
    302297     .  ' restart ne correspondent pas a celles lues dans '
    … …  
    304299        if (raz_date .ne. 1) then
    305300          write(lunout,*)
    306      .    ' On garde les dates du fichier restart'
     301     .    'GCM: On garde les dates du fichier restart'
    307302        else
    308303          annee_ref = anneeref
    … …  
    313308          time_0 = 0.
    314309          write(lunout,*)
    315      .   ' On reinitialise a la date lue dans gcm.def'
     310     .   'GCM: On reinitialise a la date lue dans gcm.def'
    316311        endif
    317312      ELSE
    … …  
    350345c   Initialisation de la physique :
    351346c   -------------------------------
    352 #ifdef CPP_PHYS
    353       IF (call_iniphys.and.iflag_phys.eq.1) THEN
     347
     348      IF (call_iniphys.and.(iflag_phys.eq.1)) THEN
    354349         latfi(1)=rlatu(1)
    355350         lonfi(1)=0.
    … …  
    370365         CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,aire,airefi)
    371366         WRITE(lunout,*)
    372      .           'WARNING!!! vitesse verticale nulle dans la physique'
     367     .       'GCM: WARNING!!! vitesse verticale nulle dans la physique'
     368! Earth:
     369         if (planet_type.eq."earth") then
     370#ifdef CPP_EARTH
    373371         CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys ,
    374372     ,                latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp     )
     373#endif
     374         endif ! of if (planet_type.eq."earth")
    375375         call_iniphys=.false.
    376       ENDIF
    377 #endif
     376      ENDIF ! of IF (call_iniphys.and.(iflag_phys.eq.1))
     377!#endif
    378378
    379379c  numero de stockage pour les fichiers de redemarrage:
    … …  
    386386      day_end = day_ini + nday
    387387      WRITE(lunout,300)day_ini,day_end
     388 300  FORMAT('1'/,15x,'run du jour',i7,2x,'au jour',i7//)
     389
     390      if (planet_type.eq."earth") then
     391#ifdef CPP_EARTH
     392      CALL dynredem0("restart.nc", day_end, phis)
     393#endif
     394      endif
     395
     396      ecripar = .TRUE.
    388397
    389398#ifdef CPP_IOIPSL
    390       CALL dynredem0("restart.nc", day_end, phis)
    391 
    392       ecripar = .TRUE.
    393 
    394399      if ( 1.eq.1) then
    395400      time_step = zdtvr
    … …  
    409414
    410415#endif
     416! #endif of #ifdef CPP_IOIPSL
    411417
    412418c  Choix des frequences de stokage pour le offline
    … …  
    432438     .              time_0)
    433439
    434 
    435 
    436  300  FORMAT('1'/,15x,'run du pas',i7,2x,'au pas',i7,2x,
    437      . 'c''est a dire du jour',i7,3x,'au jour',i7//)
    438440      END
    439441

  • LMDZ4/branches/LMDZ4-dev/libf/dyn3d/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/dyn3d/leapfrog.F

    r1114 r1140  
    77
    88cIM : pour sortir les param. du modele dans un fis. netcdf 110106
    9       USE IOIPSL
     9#ifdef CPP_IOIPSL
     10      use IOIPSL
     11#endif
    1012      USE infotrac
    1113
    … …  
    161163
    162164      character*80 dynhist_file, dynhistave_file
    163       character*20 modname
     165      character(len=20) :: modname
    164166      character*80 abort_message
    165167
    … …  
    217219        call guide(itau,ucov,vcov,teta,q,masse,ps)
    218220      else
    219         IF(prt_level>9)WRITE(*,*)'attention on ne guide pas les ',
    220      .    '6 dernieres heures'
     221        IF(prt_level>9)WRITE(lunout,*)'leapfrog: attention on ne ',
     222     .    'guide pas les 6 dernieres heures'
    221223      endif
    222224#endif
    … …  
    227229c     ENDIF
    228230c
     231
     232! Save fields obtained at previous time step as '...m1'
    229233      CALL SCOPY( ijmllm ,vcov , 1, vcovm1 , 1 )
    230234      CALL SCOPY( ijp1llm,ucov , 1, ucovm1 , 1 )
    … …  
    242246      CALL filtreg ( finvmaold ,jjp1, llm, -2,2, .TRUE., 1 )
    243247
    244       call minmax(ijp1llm,q(:,:,3),zqmin,zqmax)
     248! Ehouarn: what is this for? zqmin & zqmax are not used anyway ...
     249!      call minmax(ijp1llm,q(:,:,3),zqmin,zqmax)
    245250
    246251   2  CONTINUE
    … …  
    302307
    303308
    304          ENDIF
    305 c
    306       ENDIF
     309         ENDIF ! of IF (offline)
     310c
     311      ENDIF ! of IF( forward. OR . leapf )
    307312
    308313
    … …  
    350355c   -----------------------------------------------------
    351356
    352 #ifdef CPP_PHYS
    353357c+jld
    354358
    355359c  Diagnostique de conservation de l'énergie : initialisation
    356       IF (ip_ebil_dyn.ge.1 ) THEN
     360         IF (ip_ebil_dyn.ge.1 ) THEN
    357361          ztit='bil dyn'
    358           CALL diagedyn(ztit,2,1,1,dtphys
    359      e  , ucov    , vcov , ps, p ,pk , teta , q(:,:,1), q(:,:,2))
    360       ENDIF
     362! Ehouarn: be careful, diagedyn is Earth-specific (includes ../phylmd/..)!
     363           IF (planet_type.eq."earth") THEN
     364            CALL diagedyn(ztit,2,1,1,dtphys
     365     &    , ucov    , vcov , ps, p ,pk , teta , q(:,:,1), q(:,:,2))
     366           ENDIF
     367         ENDIF ! of IF (ip_ebil_dyn.ge.1 )
    361368c-jld
     369#ifdef CPP_IOIPSL
    362370cIM : pour sortir les param. du modele dans un fis. netcdf 110106
    363       IF (first) THEN
    364        first=.false.
     371         IF (first) THEN
     372          first=.false.
    365373#include "ini_paramLMDZ_dyn.h"
    366       ENDIF
     374         ENDIF
    367375c
    368376#include "write_paramLMDZ_dyn.h"
    369377c
    370 
    371         CALL calfis( lafin ,rdayvrai,time  ,
     378#endif
     379! #endif of #ifdef CPP_IOIPSL
     380         CALL calfis( lafin ,rdayvrai,time  ,
    372381     $               ucov,vcov,teta,q,masse,ps,p,pk,phis,phi ,
    373382     $               du,dv,dteta,dq,
    … …  
    375384     $               clesphy0, dufi,dvfi,dtetafi,dqfi,dpfi  )
    376385
    377        IF (ok_strato) THEN
    378          CALL top_bound( vcov,ucov,teta, dufi,dvfi,dtetafi)
    379        ENDIF
     386         IF (ok_strato) THEN
     387           CALL top_bound( vcov,ucov,teta, dufi,dvfi,dtetafi)
     388         ENDIF
    380389       
    381390c      ajout des tendances physiques:
    … …  
    386395c
    387396c  Diagnostique de conservation de l'énergie : difference
    388       IF (ip_ebil_dyn.ge.1 ) THEN
     397         IF (ip_ebil_dyn.ge.1 ) THEN
    389398          ztit='bil phys'
    390           CALL diagedyn(ztit,2,1,1,dtphys
    391      e  , ucov    , vcov , ps, p ,pk , teta , q(:,:,1), q(:,:,2))
    392       ENDIF
    393 #endif
     399          IF (planet_type.eq."earth") THEN
     400           CALL diagedyn(ztit,2,1,1,dtphys
     401     &     , ucov    , vcov , ps, p ,pk , teta , q(:,:,1), q(:,:,2))
     402          ENDIF
     403         ENDIF ! of IF (ip_ebil_dyn.ge.1 )
     404
    394405       ENDIF ! of IF( apphys )
    395406
    396       IF(iflag_phys.EQ.2) THEN ! "Newtonian physics" case
     407      IF(iflag_phys.EQ.2) THEN ! "Newtonian" case
    397408c   Calcul academique de la physique = Rappel Newtonien + friction
    398409c   --------------------------------------------------------------
    … …  
    472483
    473484
    474       END IF
     485      END IF ! of IF(apdiss)
    475486
    476487c ajout debug
    … …  
    545556c           IF( MOD(itau,iecri*day_step).EQ.0) THEN
    546557
    547                nbetat = nbetatdem
    548        CALL geopot  ( ip1jmp1, teta  , pk , pks,  phis  , phi        )
    549         unat=0.
    550         do l=1,llm
    551            unat(iip2:ip1jm,l)=ucov(iip2:ip1jm,l)/cu(iip2:ip1jm)
    552            vnat(:,l)=vcov(:,l)/cv(:)
    553         enddo
    554 #ifdef CPP_IOIPSL
    555 c        CALL writehist(histid,histvid,itau,vcov,
    556 c     s                       ucov,teta,phi,q,masse,ps,phis)
    557 #else
     558              nbetat = nbetatdem
     559              CALL geopot(ip1jmp1,teta,pk,pks,phis,phi)
     560              unat=0.
     561              do l=1,llm
     562                unat(iip2:ip1jm,l)=ucov(iip2:ip1jm,l)/cu(iip2:ip1jm)
     563                vnat(:,l)=vcov(:,l)/cv(:)
     564              enddo
     565#ifdef CPP_IOIPSL
     566c             CALL writehist(histid,histvid,itau,vcov,
     567c     &                      ucov,teta,phi,q,masse,ps,phis)
     568#endif
     569! For some Grads outputs of fields
     570             if (output_grads_dyn) then
    558571#include "write_grads_dyn.h"
    559 #endif
    560 
    561 
    562             ENDIF
     572             endif
     573
     574            ENDIF ! of IF(MOD(itau,iecri).EQ.0)
    563575
    564576            IF(itau.EQ.itaufin) THEN
    565577
    566578
    567 #ifdef CPP_IOIPSL
    568        CALL dynredem1("restart.nc",0.0,
    569      ,                     vcov,ucov,teta,q,masse,ps)
    570 #endif
     579              if (planet_type.eq."earth") then
     580#ifdef CPP_EARTH
     581! Write an Earth-format restart file
     582                CALL dynredem1("restart.nc",0.0,
     583     &                         vcov,ucov,teta,q,masse,ps)
     584#endif
     585              endif ! of if (planet_type.eq."earth")
    571586
    572587              CLOSE(99)
    573             ENDIF
     588            ENDIF ! of IF (itau.EQ.itaufin)
    574589
    575590c-----------------------------------------------------------------------
    … …  
    593608                        leapf =  .TRUE.
    594609                        dt  =  2.*dtvr
    595                         GO TO 2
    596                    END IF
     610                        GO TO 2 
     611                   END IF ! of IF (forward)
    597612            ELSE
    598613
    … …  
    602617                 dt  = 2.*dtvr
    603618                 GO TO 2
    604             END IF
    605 
    606       ELSE
     619            END IF ! of IF (MOD(itau,iperiod).EQ.0)
     620                   !    ELSEIF (MOD(itau-1,iperiod).EQ.0)
     621
     622      ELSE ! of IF (.not.purmats)
    607623
    608624c       ........................................................
    … …  
    627643               GO TO 2
    628644
    629             ELSE
    630 
    631             IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN
     645            ELSE ! of IF(forward)
     646
     647              IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN
    632648               IF(itau.EQ.itaufin) THEN
    633649                  iav=1
    … …  
    636652               ENDIF
    637653#ifdef CPP_IOIPSL
    638               CALL writedynav(histaveid, itau,vcov ,
     654               CALL writedynav(histaveid, itau,vcov ,
    639655     ,                          ucov,teta,pk,phi,q,masse,ps,phis)
    640656               call bilan_dyn (2,dtvr*iperiod,dtvr*day_step*periodav,
    … …  
    642658#endif
    643659
    644             ENDIF
    645 
    646                IF(MOD(itau,iecri         ).EQ.0) THEN
     660              ENDIF ! of IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin)
     661
     662              IF(MOD(itau,iecri         ).EQ.0) THEN
    647663c              IF(MOD(itau,iecri*day_step).EQ.0) THEN
    648                   nbetat = nbetatdem
    649        CALL geopot  ( ip1jmp1, teta  , pk , pks,  phis  , phi       )
    650         unat=0.
    651         do l=1,llm
    652            unat(iip2:ip1jm,l)=ucov(iip2:ip1jm,l)/cu(iip2:ip1jm)
    653            vnat(:,l)=vcov(:,l)/cv(:)
    654         enddo
    655 #ifdef CPP_IOIPSL
    656 c       CALL writehist( histid, histvid, itau,vcov ,
    657 c    ,                           ucov,teta,phi,q,masse,ps,phis)
    658 #else
     664                nbetat = nbetatdem
     665                CALL geopot(ip1jmp1,teta,pk,pks,phis,phi)
     666                unat=0.
     667                do l=1,llm
     668                  unat(iip2:ip1jm,l)=ucov(iip2:ip1jm,l)/cu(iip2:ip1jm)
     669                  vnat(:,l)=vcov(:,l)/cv(:)
     670                enddo
     671#ifdef CPP_IOIPSL
     672c               CALL writehist( histid, histvid, itau,vcov ,
     673c    &                           ucov,teta,phi,q,masse,ps,phis)
     674#endif
     675! For some Grads outputs
     676                if (output_grads_dyn) then
    659677#include "write_grads_dyn.h"
    660 #endif
    661 
    662 
    663                ENDIF
    664 
    665 #ifdef CPP_IOIPSL
    666                  IF(itau.EQ.itaufin)
    667      . CALL dynredem1("restart.nc",0.0,
    668      .                     vcov,ucov,teta,q,masse,ps)
    669 #endif
    670 
    671                  forward = .TRUE.
    672                  GO TO  1
    673 
    674             ENDIF
    675 
    676       END IF
     678                endif
     679
     680              ENDIF ! of IF(MOD(itau,iecri         ).EQ.0)
     681
     682              IF(itau.EQ.itaufin) THEN
     683                if (planet_type.eq."earth") then
     684#ifdef CPP_EARTH
     685                  CALL dynredem1("restart.nc",0.0,
     686     &                           vcov,ucov,teta,q,masse,ps)
     687#endif
     688                endif ! of if (planet_type.eq."earth")
     689              ENDIF ! of IF(itau.EQ.itaufin)
     690
     691              forward = .TRUE.
     692              GO TO  1
     693
     694            ENDIF ! of IF (forward)
     695
     696      END IF ! of IF(.not.purmats)
    677697
    678698      STOP
Note: See TracChangeset for help on using the changeset viewer.