Changeset 1673


Ignore:
Timestamp:
Oct 27, 2012, 4:23:07 PM (12 years ago)
Author:
Laurent Fairhead
Message:

Fin du phasage de la dynamique parallele localisee (petite memoire) avec le tronc LMDZ5 r1671
Il reste quelques routines a verifier (en particulier ce qui touche a l'etude des cas academiques)
et la validation a effectuer

End of the phasing of the localised (low memory) parallel dynamics package with the
LMDZ5 trunk (r1671)
Some routines still need some checking (in particular the academic cases) and some
validation is still required

Location:
LMDZ5/trunk/libf/dyn3dmem
Files:
9 added
8 deleted
46 edited

Legend:

Unmodified
Added
Removed

  • LMDZ5/trunk/libf/dyn3dmem/abort_gcm.F

    r1658 r1673  
    11!
    2 ! $Id: abort_gcm.F 1425 2010-09-02 13:45:23Z lguez $
     2! $Id$
    33!
    44c
    … …  
    4545      if (ierr .eq. 0) then
    4646        write(lunout,*) 'Everything is cool'
    47         stop
    4847      else
    4948        write(lunout,*) 'Houston, we have a problem ', ierr

  • LMDZ5/trunk/libf/dyn3dmem/academic.h

    r1632 r1673  
    11!
    2 ! $Header$
     2! $Id$
    33!
    4       real tetarappel(ip1jmp1,llm),taurappel
    5       common/academic/tetarappel,taurappel
     4      common/academic/tetarappel,knewt_t,kfrict,knewt_g,clat4
     5      real :: tetarappel(ip1jmp1,llm)
     6      real :: knewt_t(llm)
     7      real :: kfrict(llm)
     8      real :: knewt_g
     9      real :: clat4(ip1jmp1)

  • LMDZ5/trunk/libf/dyn3dmem/addfi_loc.F

    r1632 r1673  
    11!
    2 ! $Header$
     2! $Id$
    33!
    44      SUBROUTINE addfi_loc(pdt, leapf, forward,
    … …  
    77      USE parallel
    88      USE infotrac, ONLY : nqtot
     9      USE control_mod, ONLY : planet_type
    910      IMPLICIT NONE
    1011c
    … …  
    156157c$OMP END MASTER
    157158 
     159      if (planet_type=="earth") then
     160      ! earth case, special treatment for first 2 tracers (water)
    158161      DO iq = 1, 2
    159162c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    … …  
    177180c$OMP END DO NOWAIT
    178181      ENDDO
     182      else
     183      ! general case, treat all tracers equally)
     184       DO iq = 1, nqtot
     185c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     186         DO k = 1,llm
     187            DO j = ijb,ije
     188               pq(j,k,iq)= pq(j,k,iq) + pdqfi(j,k,iq) * pdt
     189               pq(j,k,iq)= AMAX1( pq(j,k,iq), qtestt )
     190            ENDDO
     191         ENDDO
     192c$OMP END DO NOWAIT
     193       ENDDO
     194      endif ! of if (planet_type=="earth")
    179195
    180196c$OMP MASTER

  • LMDZ5/trunk/libf/dyn3dmem/advtrac_loc.F

    r1632 r1673  
    11!
    2 ! $Id: advtrac_p.F 1299 2010-01-20 14:27:21Z fairhead $
     2! $Id$
    33!
    44c
    … …  
    342342c$OMP BARRIER
    343343
     344      if (planet_type=="earth") then
     345
    344346      ijb=ij_begin
    345347      ije=ij_end
    … …  
    355357       CALL qminimum_loc( q, 2, finmasse )
    356358
     359      endif ! of if (planet_type=="earth")
    357360
    358361       RETURN

  • LMDZ5/trunk/libf/dyn3dmem/bands.F90

    r1632 r1673  
    11!
    2 ! $Id: bands.F90 1279 2009-12-10 09:02:56Z fairhead $
     2! $Id$
    33!
    44  module Bands
    … …  
    105105   SUBROUTINE  Set_Bands
    106106     USE parallel
    107 #ifdef CPP_EARTH
    108 ! Ehouarn: what follows is only related to // physics; for now only for Earth
     107#ifdef CPP_PHYS
     108! Ehouarn: what follows is only related to // physics
    109109     USE mod_phys_lmdz_para, ONLY : jj_para_begin,jj_para_end
    110110#endif
    … …  
    118118      enddo
    119119         
    120 #ifdef CPP_EARTH
    121 ! Ehouarn: what follows is only related to // physics; for now only for Earth         
     120#ifdef CPP_PHYS
    122121      do i=0,MPI_Size-1
    123122        jj_Nb_physic(i)=jj_para_end(i)-jj_para_begin(i)+1
    … …  
    374373    subroutine AdjustBands_physic
    375374      use times
    376 #ifdef CPP_EARTH
    377 ! Ehouarn: what follows is only related to // physics; for now only for Earth
     375#ifdef CPP_PHYS
     376! Ehouarn: what follows is only related to // physics
    378377      USE mod_phys_lmdz_para, only : klon_mpi_para_nb
    379378#endif
    … …  
    401400      medium=medium/mpi_size     
    402401      NbTot=0
    403 #ifdef CPP_EARTH
    404 ! Ehouarn: what follows is only related to // physics; for now only for Earth
     402#ifdef CPP_PHYS
    405403      do i=0,mpi_size-1
    406404        Inc(i)=nint(klon_mpi_para_nb(i)*(medium-value(i))/value(i))

  • LMDZ5/trunk/libf/dyn3dmem/bilan_dyn_loc.F

    r1632 r1673  
    421421          Q_cum(:,jjb:jje,:,l)=0.
    422422          flux_uQ_cum(:,jjb:jje,l,:)=0.
    423           flux_v_cum(:,jjb:jje,l)=0.
    424423          if (pole_sud) jje=jj_end-1
    425424          flux_v_cum(:,jjb:jje,l)=0.

  • LMDZ5/trunk/libf/dyn3dmem/caladvtrac_loc.F

    r1632 r1673  
    88     *                   flxw, pk, iapptrac)
    99      USE parallel
    10       USE infotrac
    11       USE control_mod
     10      USE infotrac, ONLY : nqtot
     11      USE control_mod, ONLY : iapp_tracvl,planet_type
    1212      USE caladvtrac_mod
    1313      USE mod_hallo
    … …  
    3838      REAL :: masse(ijb_u:ije_u,llm)
    3939      REAL :: p( ijb_u:ije_u,llmp1)
    40       REAL :: q( ijb_u:ije_u,llm,nqtot),dq( ijb_u:ije_u,llm,2 )
     40      REAL :: q( ijb_u:ije_u,llm,nqtot),dq( ijb_u:ije_u,llm, nqtot )
    4141      REAL :: teta( ijb_u:ije_u,llm),pk( ijb_u:ije_u,llm)
    4242      REAL :: flxw(ijb_u:ije_u,llm)

  • LMDZ5/trunk/libf/dyn3dmem/calfis_loc.F

    r1658 r1673  
    2727     $                  pdqfi,
    2828     $                  pdpsfi)
    29 #ifdef CPP_EARTH
    30 ! Ehouarn: For now, calfis_p needs Earth physics
     29#ifdef CPP_PHYS
     30! If using physics
    3131c
    3232c    Auteur :  P. Le Van, F. Hourdin
    … …  
    3636      USE parallel, ONLY : omp_chunk, using_mpi,jjb_u,jje_u,jjb_v,jje_v
    3737      USE mod_interface_dyn_phys
     38      USE IOPHY
     39#endif
    3840      USE Write_Field
    3941      Use Write_field_p
    4042      USE Times
    41       USE IOPHY
    4243      USE infotrac
    4344      USE control_mod
    … …  
    145146
    146147
     148#ifdef CPP_PHYS
     149! Ehouarn: for now calfis_p needs some informations from physics to compile
    147150c    Local variables :
    148151c    -----------------
    … …  
    220223      PARAMETER(ntetaSTD=3)
    221224      REAL rtetaSTD(ntetaSTD)
    222       DATA rtetaSTD/350., 380., 405./
     225      DATA rtetaSTD/350., 380., 405./ ! Earth-specific values, beware !!
    223226      REAL PVteta(klon,ntetaSTD)
    224227     
    … …  
    243246      REAL,ALLOCATABLE,SAVE:: zdufi2(:,:),zdvfi2(:,:)
    244247      integer :: k,kstart,kend
    245       INTEGER :: offset 
     248      INTEGER :: offset
     249
     250      LOGICAL tracerdyn 
    246251c
    247252c-----------------------------------------------------------------------
    … …  
    512517
    513518
    514       IF (is_sequential) THEN
    515 c
     519      IF (is_sequential.and.(planet_type=="earth")) THEN
     520#ifdef CPP_PHYS
     521! PVtheta calls tetalevel, which is in the physics
    516522cIM calcul PV a teta=350, 380, 405K
    517523        CALL PVtheta(ngridmx,llm,pucov,pvcov,pteta,
    … …  
    519525     $           ntetaSTD,rtetaSTD,PVteta)
    520526c
     527#endif
    521528      ENDIF
    522529
    … …  
    662669c$OMP BARRIER
    663670     
    664       if (planet_type=="earth") then
    665 #ifdef CPP_EARTH
    666 
    667671
    668672!$OMP MASTER
    … …  
    675679      zdqfic_omp(:,:,:)=0.
    676680
     681#ifdef CPP_PHYS
    677682      do isplit=1,nsplit_phys
    678683
    … …  
    681686         lafin_split=lafin.and.isplit==nsplit_phys
    682687
     688      if (planet_type=="earth") then
    683689
    684690      CALL physiq (klon,
    … …  
    711717     .             PVteta)
    712718
     719      else if ( planet_type=="generic" ) then
     720
     721      CALL physiq (klon,     !! ngrid
     722     .             llm,            !! nlayer
     723     .             nqtot,          !! nq
     724     .             tname,          !! tracer names from dynamical core (given in infotrac)
     725     .             debut_split,    !! firstcall
     726     .             lafin_split,    !! lastcall
     727     .             float(day_ini), !! pday <-- day_ini (dans temps.h)
     728     .             jH_cur_split,   !! ptime "fraction of day"
     729     .             zdt_split,      !! ptimestep
     730     .             zplev_omp,  !! pplev
     731     .             zplay_omp,  !! pplay
     732     .             zphi_omp,   !! pphi
     733     .             zufi_omp,   !! pu
     734     .             zvfi_omp,   !! pv
     735     .             ztfi_omp,   !! pt
     736     .             zqfi_omp,   !! pq
     737     .             flxwfi_omp, !! pw !! or 0. anyway this is for diagnostic. not used in physiq.
     738     .             zdufi_omp,  !! pdu
     739     .             zdvfi_omp,  !! pdv
     740     .             zdtfi_omp,  !! pdt
     741     .             zdqfi_omp,  !! pdq
     742     .             zdpsrf_omp, !! pdpsrf
     743     .             tracerdyn)      !! tracerdyn <-- utilite ???
     744
     745      endif ! of if (planet_type=="earth")
     746
     747
    713748         zufi_omp(:,:)=zufi_omp(:,:)+zdufi_omp(:,:)*zdt_split
    714749         zvfi_omp(:,:)=zvfi_omp(:,:)+zdvfi_omp(:,:)*zdt_split
    … …  
    723758      enddo
    724759
     760#endif
     761! of #ifdef CPP_PHYS
     762
     763
    725764      zdufi_omp(:,:)=zdufic_omp(:,:)/nsplit_phys
    726765      zdvfi_omp(:,:)=zdvfic_omp(:,:)/nsplit_phys
    … …  
    728767      zdqfi_omp(:,:,:)=zdqfic_omp(:,:,:)/nsplit_phys
    729768
    730 #endif
    731       endif !of if (planet_type=="earth")
    732769c$OMP BARRIER
    733770
    … …  
    11791216      firstcal = .FALSE.
    11801217
    1181 #else
    1182       write(*,*) "calfis_p: for now can only work with parallel physics"
    1183       write(lunout,*)
    1184    & "calfis_p: for now can only work with parallel physics"
    1185       stop
    1186 #endif
    1187 ! of #ifdef CPP_EARTH
     1218#else
     1219      write(lunout,*)
     1220     & "calfis_p: for now can only work with parallel physics"
     1221      stop
     1222#endif
     1223! of #ifdef CPP_PHYS
    11881224      RETURN
    11891225      END

  • LMDZ5/trunk/libf/dyn3dmem/call_calfis_mod.F90

    r1632 r1673  
    136136  !$OMP END MASTER
    137137   
    138     jD_cur = jD_ref + day_ini - day_ref + int (itau * dtvr / daysec)
    139     jH_cur = jH_ref + (itau * dtvr / daysec - int(itau * dtvr / daysec))
     138           jD_cur = jD_ref + day_ini - day_ref                           &
     139     &        + itau/day_step
     140           jH_cur = jH_ref + start_time +                                &
     141     &              mod(itau,day_step)/float(day_step)
     142    if (jH_cur > 1.0 ) then
     143      jD_cur = jD_cur +1.
     144      jH_cur = jH_cur -1.
     145    endif
    140146
    141147!   Inbterface avec les routines de phylmd (phymars ... )

  • LMDZ5/trunk/libf/dyn3dmem/call_dissip_mod.F90

    r1632 r1673  
    240240    !$OMP END DO NOWAIT
    241241
     242         if (1 == 0) then
     243!!! Ehouarn: lines here 1) kill 1+1=2 in the dynamics
     244!!!                     2) should probably not be here anyway
     245!!! but are kept for those who would want to revert to previous behaviour
    242246    !$OMP MASTER               
    243247      DO ij =  1,iim
    … …  
    251255    !$OMP END MASTER
    252256   
    253     ENDIF
     257    ENDIF ! of if (1 == 0)
     258    endif ! of of (pole_nord)
    254259       
    255260    IF (pole_sud) THEN
    … …  
    269274    !$OMP END DO NOWAIT
    270275
     276    if (1 == 0) then
     277!!! Ehouarn: lines here 1) kill 1+1=2 in the dynamics
     278!!!                     2) should probably not be here anyway
     279!!! but are kept for those who would want to revert to previous behaviour
    271280    !$OMP MASTER               
    272281      DO ij =  1,iim
    … …  
    279288      ENDDO
    280289    !$OMP END MASTER
    281     ENDIF
     290    ENDIF ! of if (1 == 0)
     291    endif ! of if (pole_sud)
    282292
    283293

  • LMDZ5/trunk/libf/dyn3dmem/ce0l.F90

    r1658 r1673  
    11!
    2 ! $Id: ce0l.F90 1425 2010-09-02 13:45:23Z lguez $
     2! $Id$
    33!
    44!-------------------------------------------------------------------------------
    … …  
    1919  USE dimphy
    2020  USE comgeomphy
    21   USE mod_phys_lmdz_para
     21  USE mod_phys_lmdz_para, mpi_root_xx=>mpi_root
    2222  USE mod_const_mpi
    2323  USE infotrac
     24  USE parallel, ONLY: finalize_parallel
    2425
    2526#ifdef CPP_IOIPSL
    … …  
    3031  IMPLICIT NONE
    3132#ifndef CPP_EARTH
     33#include "iniprint.h"
    3234  WRITE(lunout,*)'limit_netcdf: Earth-specific routine, needs Earth physics'
    3335#else
    … …  
    4143#include "temps.h"
    4244#include "logic.h"
     45#ifdef CPP_MPI
     46      include 'mpif.h'
     47#endif
     48
    4349  INTEGER, PARAMETER            :: longcles=20
     50  INTEGER                       :: ierr
    4451  REAL,    DIMENSION(longcles)  :: clesphy0
    4552  REAL,    DIMENSION(iip1,jjp1) :: masque
    4653  CHARACTER(LEN=15)             :: calnd
     54  REAL,    DIMENSION(iip1,jjp1) :: phis ! geopotentiel au sol
    4755!-------------------------------------------------------------------------------
    4856  CALL conf_gcm( 99, .TRUE. , clesphy0 )
    4957
     58#ifdef CPP_MPI
    5059  CALL init_mpi
     60#endif
    5161
    5262  CALL Init_Phys_lmdz(iim,jjp1,llm,1,(/(jjm-1)*iim+2/))
    … …  
    5565       CALL abort_gcm('ce0l','In parallel mode,                         &
    5666 &                 ce0l must be called only                             &
    57  &                 for 1 process and 1 task')
     67 &                 for 1 process and 1 task',1)
    5868  ENDIF
    5969
    … …  
    7686#endif
    7787
    78   IF (config_inca /= 'none') THEN
     88  IF (type_trac == 'inca') THEN
    7989#ifdef INCA
    80     CALL init_const_lmdz(nbtr,anneeref,dayref,iphysiq,day_step,nday)
    81     CALL init_inca_para(iim,jjm+1,klon,1,klon_mpi_para_nb,0)
    82     WRITE(lunout,*)'nbtr =' , nbtr
     90      CALL init_const_lmdz( &
     91         nbtr,anneeref,dayref,&
     92         iphysiq,day_step,nday,&
     93         nbsrf, is_oce,is_sic,&
     94         is_ter,is_lic)
     95     
    8396#endif
    8497  END IF
    … …  
    90103  WRITE(lunout,'(//)')
    91104  WRITE(lunout,*) ' interbar = ',interbar
    92   CALL etat0_netcdf(interbar,masque,ok_etat0)
     105  CALL etat0_netcdf(interbar,masque,phis,ok_etat0)
    93106
    94107  IF(ok_limit) THEN
    … …  
    101114  END IF
    102115
     116  IF (grilles_gcm_netcdf) THEN
     117     WRITE(lunout,'(//)')
     118     WRITE(lunout,*) '  ***************************  '
     119     WRITE(lunout,*) '  ***  grilles_gcm_netcdf ***  '
     120     WRITE(lunout,*) '  ***************************  '
     121     WRITE(lunout,'(//)')
     122     CALL grilles_gcm_netcdf_sub(masque,phis)
     123  END IF
     124 
     125#ifdef CPP_MPI
     126!$OMP MASTER
     127  CALL MPI_FINALIZE(ierr)
     128  IF (ierr /= 0) CALL abort_gcm('ce0l','Error in MPI_FINALIZE',1)
     129!$OMP END MASTER
     130#endif
     131
    103132#endif
    104133! of #ifndef CPP_EARTH #else

  • LMDZ5/trunk/libf/dyn3dmem/comconst.h

    r1632 r1673  
    11!
    2 ! $Id: comconst.h 1279 2009-12-10 09:02:56Z fairhead $
     2! $Id$
    33!
    44!-----------------------------------------------------------------------
    55! INCLUDE comconst.h
    66
    7       COMMON/comconst/im,jm,lllm,imp1,jmp1,lllmm1,lllmp1,lcl,           &
    8      & dtvr,daysec,                                                     &
     7      COMMON/comconsti/im,jm,lllm,imp1,jmp1,lllmm1,lllmp1,lcl,          &
     8     &                 iflag_top_bound
     9      COMMON/comconstr/dtvr,daysec,                                     &
    910     & pi,dtphys,dtdiss,rad,r,cpp,kappa,cotot,unsim,g,omeg              &
    1011     &                   ,dissip_factz,dissip_deltaz,dissip_zref        &
    11      &                   ,iflag_top_bound,tau_top_bound
     12     &                   ,tau_top_bound,                                &
     13     & daylen,year_day,molmass, ihf
    1214
    1315
    1416      INTEGER im,jm,lllm,imp1,jmp1,lllmm1,lllmp1,lcl
    15       REAL dtvr,daysec
    16       REAL pi,dtphys,dtdiss,rad,r,cpp,kappa
    17       REAL cotot,unsim,g,omeg
     17      REAL dtvr ! dynamical time step (in s)
     18      REAL daysec !length (in s) of a standard day
     19      REAL pi    ! something like 3.14159....
     20      REAL dtphys ! (s) time step for the physics
     21      REAL dtdiss ! (s) time step for the dissipation
     22      REAL rad ! (m) radius of the planet
     23      REAL r ! Reduced Gas constant r=R/mu
     24             ! with R=8.31.. J.K-1.mol-1, mu: mol mass of atmosphere (kg/mol)
     25      REAL cpp   ! Specific heat Cp (J.kg-1.K-1)
     26      REAL kappa ! kappa=R/Cp
     27      REAL cotot
     28      REAL unsim ! = 1./iim
     29      REAL g ! (m/s2) gravity
     30      REAL omeg ! (rad/s) rotation rate of the planet
    1831      REAL dissip_factz,dissip_deltaz,dissip_zref
    1932      INTEGER iflag_top_bound
    2033      REAL tau_top_bound
     34      REAL daylen ! length of solar day, in 'standard' day length
     35      REAL year_day ! Number of standard days in a year
     36      REAL molmass ! (g/mol) molar mass of the atmosphere
    2137
     38      REAL ihf ! (W/m2) Intrinsic heat flux (for giant planets)
    2239
    2340!-----------------------------------------------------------------------

  • LMDZ5/trunk/libf/dyn3dmem/comdissipn.h

    r1632 r1673  
    11!
    2 ! $Header$
     2! $Id$
    33!
    4 c-----------------------------------------------------------------------
    5 c INCLUDE comdissipn.h
     4!  Attention : ce fichier include est compatible format fixe/format libre
     5!                 veillez à n'utiliser que des ! pour les commentaires
     6!                 et à bien positionner les & des lignes de continuation
     7!                 (les placer en colonne 6 et en colonne 73)
     8!-----------------------------------------------------------------------
     9! INCLUDE comdissipn.h
    610
    711      REAL  tetaudiv, tetaurot, tetah, cdivu, crot, cdivh
    8 c
    9       COMMON/comdissipn/ tetaudiv(llm),tetaurot(llm),tetah(llm)   ,
    10      1                        cdivu,      crot,         cdivh
     12!
     13      COMMON/comdissipn/ tetaudiv(llm),tetaurot(llm),tetah(llm)   ,     &
     14     &                        cdivu,      crot,         cdivh
    1115
    12 c
    13 c    Les parametres de ce common proviennent des calculs effectues dans
    14 c             Inidissip  .
    15 c
    16 c-----------------------------------------------------------------------
     16!
     17!    Les parametres de ce common proviennent des calculs effectues dans
     18!             Inidissip  .
     19!
     20!-----------------------------------------------------------------------

  • LMDZ5/trunk/libf/dyn3dmem/comvert.h

    r1632 r1673  
    11!
    2 ! $Header$
     2! $Id$
    33!
    44!-----------------------------------------------------------------------
    55!   INCLUDE 'comvert.h'
    66
    7       COMMON/comvert/ap(llm+1),bp(llm+1),presnivs(llm),dpres(llm),       &
    8      &               pa,preff,nivsigs(llm),nivsig(llm+1)
     7      COMMON/comvertr/ap(llm+1),bp(llm+1),presnivs(llm),dpres(llm),     &
     8     &               pa,preff,nivsigs(llm),nivsig(llm+1),               &
     9     &               aps(llm),bps(llm),scaleheight,pseudoalt(llm)
    910
    10       common/comverti/disvert_type
     11      common/comverti/disvert_type, pressure_exner
    1112
    12       REAL ap,bp,presnivs,dpres,pa,preff,nivsigs,nivsig
     13      real ap     ! hybrid pressure contribution at interlayers
     14      real bp     ! hybrid sigma contribution at interlayer
     15      real presnivs ! (reference) pressure at mid-layers
     16      real dpres
     17      real pa     ! reference pressure (Pa) at which hybrid coordinates
     18                  ! become purely pressure
     19      real preff  ! reference surface pressure (Pa)
     20      real nivsigs
     21      real nivsig
     22      real aps    ! hybrid pressure contribution at mid-layers
     23      real bps    ! hybrid sigma contribution at mid-layers
     24      real scaleheight ! atmospheric (reference) scale height (km)
     25      real pseudoalt ! for planets
    1326
    1427      integer disvert_type ! type of vertical discretization:
    … …  
    1730                           ! 2: Planets (default for planet_type!=earth),
    1831                           !     using 'z2sig.def' (or 'esasig.def) file
    19 !-----------------------------------------------------------------------
     32
     33      logical pressure_exner
     34!     compute pressure inside layers using Exner function, else use mean
     35!     of pressure values at interfaces
     36
     37 !-----------------------------------------------------------------------

  • LMDZ5/trunk/libf/dyn3dmem/conf_gcm.F

    r1658 r1673  
    11!
    2 ! $Id: conf_gcm.F 1403 2010-07-01 09:02:53Z fairhead $
     2! $Id$
    33!
    44c
    … …  
    1818      use parallel, ONLY : omp_chunk
    1919      USE control_mod
     20      USE infotrac, ONLY : type_trac
    2021      IMPLICIT NONE
    2122c-----------------------------------------------------------------------
    … …  
    103104      CALL getin('lunout', lunout)
    104105      IF (lunout /= 5 .and. lunout /= 6) THEN
    105         OPEN(lunout,FILE='lmdz.out')
     106        OPEN(UNIT=lunout,FILE='lmdz.out_0000',ACTION='write',
     107     &          STATUS='unknown',FORM='formatted')
     108
    106109      ENDIF
    107110
    … …  
    166169      CALL getin('nday',nday)
    167170
     171!Config  Key  = starttime
     172!Config  Desc = Heure de depart de la simulation
     173!Config  Def  = 0
     174!Config  Help = Heure de depart de la simulation
     175!Config         en jour
     176      starttime = 0
     177      CALL getin('starttime',starttime)
     178
    168179!Config  Key  = day_step
    169180!Config  Desc = nombre de pas par jour
    … …  
    226237       CALL getin('output_grads_dyn',output_grads_dyn)
    227238
    228 !Config  Key  = idissip
     239!Config  Key  = dissip_period
    229240!Config  Desc = periode de la dissipation
    230 !Config  Def  = 10
     241!Config  Def  = 0
    231242!Config  Help = periode de la dissipation
    232 !Config         (en pas) ... a completer !
    233        idissip = 10
    234        CALL getin('idissip',idissip)
     243!Config  dissip_period=0 => la valeur sera calcule dans inidissip       
     244!Config  dissip_period>0 => on prend cette valeur
     245       dissip_period = 0
     246       CALL getin('dissip_period',dissip_period)
    235247
    236248ccc  ....   P. Le Van , modif le 29/04/97 .pour la dissipation  ...
    … …  
    579591       offline = .FALSE.
    580592       CALL getin('offline',offline)
     593       IF (offline .AND. adjust) THEN
     594          WRITE(lunout,*)
     595     &         'WARNING : option offline does not work with adjust=y :'
     596          WRITE(lunout,*) 'the files defstoke.nc, fluxstoke.nc ',
     597     &         'and fluxstokev.nc will not be created'
     598          WRITE(lunout,*)
     599     &         'only the file phystoke.nc will still be created '
     600       END IF
     601       
     602!Config  Key  = type_trac
     603!Config  Desc = Choix de couplage avec model de chimie INCA ou REPROBUS
     604!Config  Def  = lmdz
     605!Config  Help =
     606!Config         'lmdz' = pas de couplage, pur LMDZ
     607!Config         'inca' = model de chime INCA
     608!Config         'repr' = model de chime REPROBUS
     609      type_trac = 'lmdz'
     610      CALL getin('type_trac',type_trac)
     611
    581612
    582613!Config  Key  = config_inca
    … …  
    628659      write(lunout,*)' periodav = ', periodav
    629660      write(lunout,*)' output_grads_dyn = ', output_grads_dyn
    630       write(lunout,*)' idissip = ', idissip
     661      write(lunout,*)' dissip_period = ', dissip_period
    631662      write(lunout,*)' lstardis = ', lstardis
    632663      write(lunout,*)' nitergdiv = ', nitergdiv
    … …  
    651682      write(lunout,*)' tauyy = ', tauyy
    652683      write(lunout,*)' offline = ', offline
     684      write(lunout,*)' type_trac = ', type_trac
    653685      write(lunout,*)' config_inca = ', config_inca
    654686      write(lunout,*)' ok_dynzon = ', ok_dynzon
    … …  
    769801       offline = .FALSE.
    770802       CALL getin('offline',offline)
     803       IF (offline .AND. adjust) THEN
     804          WRITE(lunout,*)
     805     &         'WARNING : option offline does not work with adjust=y :'
     806          WRITE(lunout,*) 'the files defstoke.nc, fluxstoke.nc ',
     807     &         'and fluxstokev.nc will not be created'
     808          WRITE(lunout,*)
     809     &         'only the file phystoke.nc will still be created '
     810       END IF
     811
     812!Config  Key  = type_trac
     813!Config  Desc = Choix de couplage avec model de chimie INCA ou REPROBUS
     814!Config  Def  = lmdz
     815!Config  Help =
     816!Config         'lmdz' = pas de couplage, pur LMDZ
     817!Config         'inca' = model de chime INCA
     818!Config         'repr' = model de chime REPROBUS
     819      type_trac = 'lmdz'
     820      CALL getin('type_trac',type_trac)
    771821
    772822!Config  Key  = config_inca
    … …  
    875925      CALL getin('ok_etat0',ok_etat0)
    876926
     927!Config  Key  = grilles_gcm_netcdf
     928!Config  Desc = creation de fichier grilles_gcm.nc dans create_etat0_limit
     929!Config  Def  = n
     930      grilles_gcm_netcdf = .FALSE.
     931      CALL getin('grilles_gcm_netcdf',grilles_gcm_netcdf)
     932
    877933      write(lunout,*)' #########################################'
    878934      write(lunout,*)' Configuration des parametres de cel0'
    … …  
    889945      write(lunout,*)' periodav = ', periodav
    890946      write(lunout,*)' output_grads_dyn = ', output_grads_dyn
    891       write(lunout,*)' idissip = ', idissip
     947      write(lunout,*)' dissip_period = ', dissip_period
    892948      write(lunout,*)' lstardis = ', lstardis
    893949      write(lunout,*)' nitergdiv = ', nitergdiv
    … …  
    912968      write(lunout,*)' tauy = ', tauy
    913969      write(lunout,*)' offline = ', offline
     970      write(lunout,*)' type_trac = ', type_trac
    914971      write(lunout,*)' config_inca = ', config_inca
    915972      write(lunout,*)' ok_dynzon = ', ok_dynzon
    … …  
    923980      write(lunout,*)' ok_limit = ', ok_limit
    924981      write(lunout,*)' ok_etat0 = ', ok_etat0
     982      write(lunout,*)' grilles_gcm_netcdf = ', grilles_gcm_netcdf
    925983c
    926984      RETURN

  • LMDZ5/trunk/libf/dyn3dmem/control_mod.F90

    r1657 r1673  
    1010  IMPLICIT NONE
    1111
    12   REAL    :: periodav
     12  REAL    :: periodav, starttime
    1313  INTEGER :: nday,day_step,iperiod,iapp_tracvl,nsplit_phys
    14   INTEGER :: iconser,iecri,idissip,iphysiq,iecrimoy
     14  INTEGER :: iconser,iecri,dissip_period,iphysiq,iecrimoy
    1515  INTEGER :: dayref,anneeref, raz_date, ip_ebil_dyn
    1616  LOGICAL :: offline

  • LMDZ5/trunk/libf/dyn3dmem/defrun.F

    r1658 r1673  
    11!
    2 ! $Id: defrun.F 1403 2010-07-01 09:02:53Z fairhead $
     2! $Id$
    33!
    44c
    … …  
    132132
    133133      READ (tapedef,9001) ch1,ch4
    134       READ (tapedef,*)    idissip
    135       WRITE(tapeout,9001) ch1,'idissip'
    136       WRITE(tapeout,*)    idissip
     134      READ (tapedef,*)    dissip_period
     135      WRITE(tapeout,9001) ch1,'dissip_period'
     136      WRITE(tapeout,*)    dissip_period
    137137
    138138ccc  ....   P. Le Van , modif le 29/04/97 .pour la dissipation  ...

  • LMDZ5/trunk/libf/dyn3dmem/dynetat0_loc.F

    r1658 r1673  
    55     .                    teta,q,masse,ps,phis,time)
    66      USE infotrac
     7      use control_mod, only : planet_type
    78      USE parallel
    89      IMPLICIT NONE
    … …  
    5758      REAL,ALLOCATABLE :: phis_glo(:)
    5859
     60      INTEGER idecal
     61
    5962c-----------------------------------------------------------------------
    6063c  Ouverture NetCDF du fichier etat initial
    … …  
    8487      ENDIF
    8588
     89      !!! AS: idecal is a hack to be able to read planeto starts...
     90      !!!     .... while keeping everything OK for LMDZ EARTH
     91      if (planet_type.eq."generic") then
     92          print*,'NOTE NOTE NOTE : Planeto-like start files'
     93          idecal = 4
     94          annee_ref  = 2000
     95      else
     96          print*,'NOTE NOTE NOTE : Earth-like start files'
     97          idecal = 5
     98          annee_ref  = tab_cntrl(5)
     99      endif
     100
     101
    86102      im         = tab_cntrl(1)
    87103      jm         = tab_cntrl(2)
    88104      lllm       = tab_cntrl(3)
    89105      day_ref    = tab_cntrl(4)
    90       annee_ref  = tab_cntrl(5)
    91       rad        = tab_cntrl(6)
    92       omeg       = tab_cntrl(7)
    93       g          = tab_cntrl(8)
    94       cpp        = tab_cntrl(9)
    95       kappa      = tab_cntrl(10)
    96       daysec     = tab_cntrl(11)
    97       dtvr       = tab_cntrl(12)
    98       etot0      = tab_cntrl(13)
    99       ptot0      = tab_cntrl(14)
    100       ztot0      = tab_cntrl(15)
    101       stot0      = tab_cntrl(16)
    102       ang0       = tab_cntrl(17)
    103       pa         = tab_cntrl(18)
    104       preff      = tab_cntrl(19)
    105 c
    106       clon       = tab_cntrl(20)
    107       clat       = tab_cntrl(21)
    108       grossismx  = tab_cntrl(22)
    109       grossismy  = tab_cntrl(23)
    110 c
    111       IF ( tab_cntrl(24).EQ.1. )  THEN
     106      rad        = tab_cntrl(idecal+1)
     107      omeg       = tab_cntrl(idecal+2)
     108      g          = tab_cntrl(idecal+3)
     109      cpp        = tab_cntrl(idecal+4)
     110      kappa      = tab_cntrl(idecal+5)
     111      daysec     = tab_cntrl(idecal+6)
     112      dtvr       = tab_cntrl(idecal+7)
     113      etot0      = tab_cntrl(idecal+8)
     114      ptot0      = tab_cntrl(idecal+9)
     115      ztot0      = tab_cntrl(idecal+10)
     116      stot0      = tab_cntrl(idecal+11)
     117      ang0       = tab_cntrl(idecal+12)
     118      pa         = tab_cntrl(idecal+13)
     119      preff      = tab_cntrl(idecal+14)
     120c
     121      clon       = tab_cntrl(idecal+15)
     122      clat       = tab_cntrl(idecal+16)
     123      grossismx  = tab_cntrl(idecal+17)
     124      grossismy  = tab_cntrl(idecal+18)
     125c
     126      IF ( tab_cntrl(idecal+19).EQ.1. )  THEN
    112127        fxyhypb  = . TRUE .
    113128c        dzoomx   = tab_cntrl(25)
    … …  
    118133        fxyhypb = . FALSE .
    119134        ysinus  = . FALSE .
    120         IF( tab_cntrl(27).EQ.1. ) ysinus = . TRUE.
     135        IF( tab_cntrl(idecal+22).EQ.1. ) ysinus = . TRUE.
    121136      ENDIF
    122137
    … …  
    266281      ierr = NF_INQ_VARID (nid, "temps", nvarid)
    267282      IF (ierr .NE. NF_NOERR) THEN
    268          write(lunout,*)"dynetat0_loc: Le champ <temps> est absent"
    269          CALL abort
     283         write(lunout,*)"dynetat0: Le champ <temps> est absent"
     284         write(lunout,*)"dynetat0: J essaie <Time>"
     285         ierr = NF_INQ_VARID (nid, "Time", nvarid)
     286         IF (ierr .NE. NF_NOERR) THEN
     287            write(lunout,*)"dynetat0: Le champ <Time> est absent"
     288            CALL abort
     289         ENDIF
    270290      ENDIF
    271291#ifdef NC_DOUBLE

  • LMDZ5/trunk/libf/dyn3dmem/dynredem_loc.F

    r1657 r1673  
    11!
    2 ! $Id: dynredem_p.F 1299 2010-01-20 14:27:21Z fairhead $
     2! $Id$
    33!
    44c
    … …  
    126126       tab_cntrl(30) =  REAL(iday_end)
    127127       tab_cntrl(31) =  REAL(itau_dyn + itaufin)
     128c start_time: start_time of simulation (not necessarily 0.)
     129       tab_cntrl(32) = start_time
    128130c
    129131c    .........................................................
    … …  
    635637      CALL dynredem_write_u(nid,"ps",ps,1)
    636638
    637       IF (config_inca == 'none') THEN
     639      IF (type_trac /= 'inca') THEN
    638640        DO iq=1,nqtot
    639641          CALL dynredem_write_u(nid,tname(iq),q(:,:,iq),llm)

  • LMDZ5/trunk/libf/dyn3dmem/ener.h

    r1632 r1673  
    11!
    2 ! $Header$
     2! $Id$
    33!
    4 !-----------------------------------------------------------------------
     4!  ATTENTION!!!!: ce fichier include est compatible format fixe/format libre
     5!                 veillez à n'utiliser que des ! pour les commentaires
     6!                 et à bien positionner les & des lignes de continuation
     7!                 (les placer en colonne 6 et en colonne 73)
     8!
    59! INCLUDE 'ener.h'
    610
    7       COMMON/ener/ang0,etot0,ptot0,ztot0,stot0,                          &
    8      &            ang,etot,ptot,ztot,stot,rmsdpdt ,                      &
     11      COMMON/ener/ang0,etot0,ptot0,ztot0,stot0,                         &
     12     &            ang,etot,ptot,ztot,stot,rmsdpdt ,                     &
    913     &            rmsv,gtot(llmm1)
    10 
    11       REAL ang0,etot0,ptot0,ztot0,stot0,                                 &
     14      REAL ang0,etot0,ptot0,ztot0,stot0,                                &
    1215     &     ang,etot,ptot,ztot,stot,rmsdpdt,rmsv,gtot
    1316

  • LMDZ5/trunk/libf/dyn3dmem/etat0_netcdf.F90

    r1658 r1673  
    11!
    2 ! $Id: etat0_netcdf.F90 1425 2010-09-02 13:45:23Z lguez $
     2! $Id$
    33!
    44!-------------------------------------------------------------------------------
    55!
    6 SUBROUTINE etat0_netcdf(ib, masque, letat0)
     6SUBROUTINE etat0_netcdf(ib, masque, phis, letat0)
    77!
    88!-------------------------------------------------------------------------------
    … …  
    3737  LOGICAL,                    INTENT(IN)    :: ib     ! barycentric interpolat.
    3838  REAL, DIMENSION(iip1,jjp1), INTENT(INOUT) :: masque ! land mask
     39  REAL, DIMENSION(iip1,jjp1), INTENT(OUT)   :: phis   ! geopotentiel au sol
    3940  LOGICAL,                    INTENT(IN)    :: letat0 ! F: masque only required
    4041#ifndef CPP_EARTH
    … …  
    5152  REAL,    DIMENSION(klon)                 :: tsol, qsol
    5253  REAL,    DIMENSION(klon)                 :: sn, rugmer, run_off_lic_0
    53   REAL,    DIMENSION(iip1,jjp1)            :: orog, rugo, psol, phis
     54  REAL,    DIMENSION(iip1,jjp1)            :: orog, rugo, psol
    5455  REAL,    DIMENSION(iip1,jjp1,llm+1)      :: p3d
    5556  REAL,    DIMENSION(iip1,jjp1,llm)        :: uvent, t3d, tpot, qsat, qd
    … …  
    9899  REAL    :: dummy
    99100  LOGICAL :: ok_newmicro, ok_journe, ok_mensuel, ok_instan, ok_hf
    100   LOGICAL :: ok_LES, ok_ade, ok_aie, aerosol_couple, new_aod
     101  LOGICAL :: ok_LES, ok_ade, ok_aie, aerosol_couple, new_aod, callstats
    101102  INTEGER :: iflag_radia, flag_aerosol
    102103  REAL    :: bl95_b0, bl95_b1, fact_cldcon, facttemps, ratqsbas, ratqshaut
    … …  
    130131!--- CONSTRUCT A GRID
    131132  CALL conf_phys(  ok_journe, ok_mensuel, ok_instan, ok_hf, ok_LES,     &
     133                   callstats,                                           &
    132134                   solarlong0,seuil_inversion,                          &
    133135                   fact_cldcon, facttemps,ok_newmicro,iflag_radia,      &
    … …  
    137139                   flag_aerosol, new_aod,                               &
    138140                   bl95_b0, bl95_b1,                                    &
    139                    iflag_thermals,nsplit_thermals,tau_thermals,         &
    140                    iflag_thermals_ed,iflag_thermals_optflux,            &
    141                    iflag_coupl,iflag_clos,iflag_wake, read_climoz,      &
     141                   read_climoz,                                         &
    142142                   alp_offset)
    143143
    … …  
    251251!*******************************************************************************
    252252  CALL pression(ip1jmp1, ap, bp, psol, p3d)
    253   CALL exner_hyb(ip1jmp1, psol, p3d, alpha, beta, pks, pk, y)
     253  if (pressure_exner) then
     254    CALL exner_hyb(ip1jmp1, psol, p3d, alpha, beta, pks, pk, y)
     255  else
     256    CALL exner_milieu(ip1jmp1,psol,p3d,beta,pks,pk,y)
     257  endif
    254258  pls(:,:,:)=preff*(pk(:,:,:)/cpp)**(1./kappa)
    255259!  WRITE(lunout,*) 'P3D :', p3d(10,20,:)

  • LMDZ5/trunk/libf/dyn3dmem/exner_hyb.F

    r1658 r1673  
    11!
    2 ! $Id: exner_hyb.F 1403 2010-07-01 09:02:53Z fairhead $
     2! $Id$
    33!
    44      SUBROUTINE  exner_hyb ( ngrid, ps, p,alpha,beta, pks, pk, pkf )
    … …  
    5151      REAL SSUM
    5252c
     53      logical,save :: firstcall=.true.
     54      character(len=*),parameter :: modname="exner_hyb"
     55     
     56      ! Sanity check
     57      if (firstcall) then
     58        ! sanity checks for Shallow Water case (1 vertical layer)
     59        if (llm.eq.1) then
     60          if (kappa.ne.1) then
     61            call abort_gcm(modname,
     62     &      "kappa!=1 , but running in Shallow Water mode!!",42)
     63          endif
     64          if (cpp.ne.r) then
     65            call abort_gcm(modname,
     66     &      "cpp!=r , but running in Shallow Water mode!!",42)
     67          endif
     68        endif ! of if (llm.eq.1)
     69
     70        firstcall=.false.
     71      endif ! of if (firstcall)
    5372
    5473      if (llm.eq.1) then
    55         ! Specific behaviour for Shallow Water (1 vertical layer) case
    56      
    57         ! Sanity checks
    58         if (kappa.ne.1) then
    59           call abort_gcm("exner_hyb",
    60      &    "kappa!=1 , but running in Shallow Water mode!!",42)
    61         endif
    62         if (cpp.ne.r) then
    63         call abort_gcm("exner_hyb",
    64      &    "cpp!=r , but running in Shallow Water mode!!",42)
    65         endif
    6674       
    6775        ! Compute pks(:),pk(:),pkf(:)
    … …  
    7785        ! our work is done, exit routine
    7886        return
     87
    7988      endif ! of if (llm.eq.1)
    8089
     90!!!! General case:
    8191     
    8292      unpl2k    = 1.+ 2.* kappa

  • LMDZ5/trunk/libf/dyn3dmem/exner_hyb_loc.F

    r1657 r1673  
    5757      EXTERNAL SSUM
    5858      INTEGER ije,ijb,jje,jjb
     59      logical,save :: firstcall=.true.
     60!$OMP THREADPRIVATE(firstcall)
     61      character(len=*),parameter :: modname="exner_hyb_loc"
    5962c
    6063c$OMP BARRIER           
    6164
     65      ! Sanity check
     66      if (firstcall) then
     67        ! sanity checks for Shallow Water case (1 vertical layer)
     68        if (llm.eq.1) then
     69          if (kappa.ne.1) then
     70            call abort_gcm(modname,
     71     &      "kappa!=1 , but running in Shallow Water mode!!",42)
     72          endif
     73          if (cpp.ne.r) then
     74            call abort_gcm(modname,
     75     &      "cpp!=r , but running in Shallow Water mode!!",42)
     76          endif
     77        endif ! of if (llm.eq.1)
     78
     79        firstcall=.false.
     80      endif ! of if (firstcall)
     81
     82c$OMP BARRIER
     83
     84! Specific behaviour for Shallow Water (1 vertical layer) case
    6285      if (llm.eq.1) then
    63         ! Specific behaviour for Shallow Water (1 vertical layer) case
    64      
    65         ! Sanity checks
    66         if (kappa.ne.1) then
    67           call abort_gcm("exner_hyb",
    68      &    "kappa!=1 , but running in Shallow Water mode!!",42)
    69         endif
    70         if (cpp.ne.r) then
    71         call abort_gcm("exner_hyb",
    72      &    "cpp!=r , but running in Shallow Water mode!!",42)
    73         endif
    7486       
    7587        ! Compute pks(:),pk(:),pkf(:)
    … …  
    111123      endif
    112124!$OMP END MASTER
    113 
     125!$OMP BARRIER
    114126        jjb=jj_begin
    115127        jje=jj_end

  • LMDZ5/trunk/libf/dyn3dmem/exner_milieu_loc.F

    r1632 r1673  
    11!
    2 ! $Id $
     2! $Id$
    33!
    44      SUBROUTINE  exner_milieu_loc ( ngrid, ps, p,beta, pks, pk, pkf )
    … …  
    5454      logical,save :: firstcall=.true.
    5555!$OMP THREADPRIVATE(firstcall)
    56       character(len=*),parameter :: modname="exner_milieu_p"
     56      character(len=*),parameter :: modname="exner_milieu_loc"
    5757
    5858      ! Sanity check
    5959      if (firstcall) then
    60         ! check that vertical discretization is compatible
    61         ! with this routine
    62         if (disvert_type.ne.2) then
    63           call abort_gcm(modname,
    64      &     "this routine should only be called if disvert_type==2",42)
    65         endif
    6660       
    6761        ! sanity checks for Shallow Water case (1 vertical layer)
    … …  
    123117      endif
    124118!$OMP END MASTER
    125 
     119!$OMP BARRIER
    126120        jjb=jj_begin
    127121        jje=jj_end
    … …  
    171165      endif
    172166c$OMP END MASTER
     167c$OMP BARRIER
    173168c
    174169c

  • LMDZ5/trunk/libf/dyn3dmem/filtreg_p.F

    r1632 r1673  
    208208               IF( ifiltre.EQ.-2 )   THEN
    209209                  DO j = jdfil,jffil
     210#ifdef BLAS
    210211                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    211212     &                    matrinvn(1,1,j), iim,
    212213     &                    champ_loc(1,j,1), iip1*nlat, 0.0,
    213214     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
     215#else
     216                     champ_fft(:,j-jdfil+1,:)
     217     &                    =matmul(matrinvn(:,:,j),champ_loc(:iim,j,:))
     218#endif
    214219                  ENDDO
    215220                 
    216221               ELSE IF ( griscal )     THEN
    217222                  DO j = jdfil,jffil
     223#ifdef BLAS
    218224                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    219225     &                    matriceun(1,1,j), iim,
    220226     &                    champ_loc(1,j,1), iip1*nlat, 0.0,
    221227     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
     228#else
     229                     champ_fft(:,j-jdfil+1,:)
     230     &                    =matmul(matriceun(:,:,j),champ_loc(:iim,j,:))
     231#endif
    222232                  ENDDO
    223233                 
    224234               ELSE
    225235                  DO j = jdfil,jffil
     236#ifdef BLAS
    226237                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    227238     &                    matricevn(1,1,j), iim,
    228239     &                    champ_loc(1,j,1), iip1*nlat, 0.0,
    229240     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
     241#else
     242                     champ_fft(:,j-jdfil+1,:)
     243     &                    =matmul(matricevn(:,:,j),champ_loc(:iim,j,:))
     244#endif
    230245                  ENDDO
    231246                 
    … …  
    236251               IF( ifiltre.EQ.-2 )   THEN
    237252                  DO j = jdfil,jffil
     253#ifdef BLAS
    238254                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    239255     &                    matrinvs(1,1,j-jfiltsu+1), iim,
    240256     &                    champ_loc(1,j,1), iip1*nlat, 0.0,
    241257     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
     258#else
     259                     champ_fft(:,j-jdfil+1,:)
     260     &                    =matmul(matrinvs(:,:,j-jfiltsu+1),
     261     &                            champ_loc(:iim,j,:))
     262#endif
    242263                  ENDDO
    243264                 
    … …  
    245266                 
    246267                  DO j = jdfil,jffil
     268#ifdef BLAS
    247269                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    248270     &                    matriceus(1,1,j-jfiltsu+1), iim,
    249271     &                    champ_loc(1,j,1), iip1*nlat, 0.0,
    250272     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
     273#else
     274                     champ_fft(:,j-jdfil+1,:)
     275     &                    =matmul(matriceus(:,:,j-jfiltsu+1),
     276     &                            champ_loc(:iim,j,:))
     277#endif
    251278                  ENDDO
    252279                 
    … …  
    254281                 
    255282                  DO j = jdfil,jffil
     283#ifdef BLAS
    256284                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    257285     &                    matricevs(1,1,j-jfiltsv+1), iim,
    258286     &                    champ_loc(1,j,1), iip1*nlat, 0.0,
    259287     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
     288#else
     289                     champ_fft(:,j-jdfil+1,:)
     290     &                    =matmul(matricevs(:,:,j-jfiltsv+1),
     291     &                            champ_loc(:iim,j,:))
     292#endif
    260293                  ENDDO
    261294                 

  • LMDZ5/trunk/libf/dyn3dmem/friction_loc.F

    r1632 r1673  
    66      USE parallel
    77      USE control_mod
     8#ifdef CPP_IOIPSL
     9      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
    814      IMPLICIT NONE
    915
    10 c=======================================================================
    11 c
    12 c
    13 c   Objet:
    14 c   ------
    15 c
    16 c  ***********
    17 c    Friction
    18 c  ***********
    19 c
    20 c=======================================================================
     16!=======================================================================
     17!
     18!   Friction for the Newtonian case:
     19!   --------------------------------
     20!    2 possibilities (depending on flag 'friction_type'
     21!     friction_type=0 : A friction that is only applied to the lowermost
     22!                       atmospheric layer
     23!     friction_type=1 : Friction applied on all atmospheric layer (but
     24!       (default)       with stronger magnitude near the surface; see
     25!                       iniacademic.F)
     26!=======================================================================
    2127
    2228#include "dimensions.h"
    … …  
    2430#include "comgeom2.h"
    2531#include "comconst.h"
    26 
    27       REAL pdt
     32#include "iniprint.h"
     33#include "academic.h"
     34
     35! arguments:
     36      REAL,INTENT(inout) :: ucov( iip1,jjb_u:jje_u,llm )
     37      REAL,INTENT(inout) :: vcov( iip1,jjb_v:jje_v,llm )
     38      REAL,INTENT(in) :: pdt ! time step
     39
     40! local variables:
     41
    2842      REAL modv(iip1,jjb_u:jje_u),zco,zsi
    2943      REAL vpn,vps,upoln,upols,vpols,vpoln
    3044      REAL u2(iip1,jjb_u:jje_u),v2(iip1,jjb_v:jje_v)
    31       REAL ucov( iip1,jjb_u:jje_u,llm ),vcov( iip1,jjb_v:jje_v,llm )
    32       INTEGER  i,j
    33       REAL cfric
    34       parameter (cfric=1.e-5)
     45      INTEGER  i,j,l
     46      REAL,PARAMETER :: cfric=1.e-5
     47      LOGICAL,SAVE :: firstcall=.true.
     48      INTEGER,SAVE :: friction_type=1
     49      CHARACTER(len=20) :: modname="friction_p"
     50      CHARACTER(len=80) :: abort_message
     51!$OMP THREADPRIVATE(firstcall,friction_type)
    3552      integer :: jjb,jje
    3653
    37 
     54!$OMP SINGLE
     55      IF (firstcall) THEN
     56        ! set friction type
     57        call getin("friction_type",friction_type)
     58        if ((friction_type.lt.0).or.(friction_type.gt.1)) then
     59          abort_message="wrong friction type"
     60          write(lunout,*)'Friction: wrong friction type',friction_type
     61          call abort_gcm(modname,abort_message,42)
     62        endif
     63        firstcall=.false.
     64      ENDIF
     65!$OMP END SINGLE COPYPRIVATE(friction_type,firstcall)
     66
     67      if (friction_type.eq.0) then ! friction on first layer only
     68!$OMP SINGLE
    3869c   calcul des composantes au carre du vent naturel
    3970      jjb=jj_begin
    … …  
    138169         vcov(iip1,j,1)=vcov(1,j,1)
    139170      enddo
     171!$OMP END SINGLE
     172      endif ! of if (friction_type.eq.0)
     173
     174      if (friction_type.eq.1) then
     175       ! for ucov()
     176        jjb=jj_begin
     177        jje=jj_end
     178        if (pole_nord) jjb=jj_begin+1
     179        if (pole_sud) jje=jj_end-1
     180
     181!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     182        do l=1,llm
     183          ucov(1:iip1,jjb:jje,l)=ucov(1:iip1,jjb:jje,l)*
     184     &                                  (1.-pdt*kfrict(l))
     185        enddo
     186!$OMP END DO NOWAIT
     187       
     188       ! for vcoc()
     189        jjb=jj_begin
     190        jje=jj_end
     191        if (pole_sud) jje=jj_end-1
     192       
     193!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     194        do l=1,llm
     195          vcov(1:iip1,jjb:jje,l)=vcov(1:iip1,jjb:jje,l)*
     196     &                                  (1.-pdt*kfrict(l))
     197        enddo
     198!$OMP END DO
     199      endif ! of if (friction_type.eq.1)
    140200
    141201      RETURN

  • LMDZ5/trunk/libf/dyn3dmem/gcm.F

    r1658 r1673  
    11!
    2 ! $Id: gcm.F 1403 2010-07-01 09:02:53Z fairhead $
     2! $Id$
    33!
    44c
    … …  
    2020      USE control_mod
    2121
    22 ! Ehouarn: for now these only apply to Earth:
    23 #ifdef CPP_EARTH
     22#ifdef CPP_PHYS
    2423      USE mod_grid_phy_lmdz
    2524      USE mod_phys_lmdz_para, ONLY : klon_mpi_para_nb
    … …  
    8786
    8887      REAL zdtvr
    89 c      INTEGER nbetatmoy, nbetatdem,nbetat
    90       INTEGER nbetatmoy, nbetatdem
    9188
    9289c   variables dynamiques
    … …  
    189186      call ini_getparam("out.def")
    190187      call Read_Distrib
    191 ! Ehouarn : temporarily (?) keep this only for Earth
    192       if (planet_type.eq."earth") then
    193 #ifdef CPP_EARTH
     188
     189#ifdef CPP_PHYS
    194190        CALL Init_Phys_lmdz(iim,jjp1,llm,mpi_size,distrib_phys)
    195191#endif
    196       endif ! of if (planet_type.eq."earth")
    197192      CALL set_bands
    198 #ifdef CPP_EARTH
    199 ! Ehouarn: For now only Earth physics is parallel
     193#ifdef CPP_PHYS
    200194      CALL Init_interface_dyn_phys
    201195#endif
    … …  
    209203c$OMP END PARALLEL
    210204
    211 ! Ehouarn : temporarily (?) keep this only for Earth
    212       if (planet_type.eq."earth") then
    213 #ifdef CPP_EARTH
     205#ifdef CPP_PHYS
    214206c$OMP PARALLEL
    215207      call InitComgeomphy
    216208c$OMP END PARALLEL
    217209#endif
    218       endif ! of if (planet_type.eq."earth")
    219210
    220211c-----------------------------------------------------------------------
    … …  
    240231#endif
    241232
    242       IF (config_inca /= 'none') THEN
     233      IF (type_trac == 'inca') THEN
    243234#ifdef INCA
    244235         call init_const_lmdz(
    … …  
    282273        endif
    283274
    284         if (planet_type.eq."earth") then
    285 #ifdef CPP_EARTH
     275!        if (planet_type.eq."earth") then
    286276! Load an Earth-format start file
    287277         CALL dynetat0_loc("start.nc",vcov,ucov,
    288278     &              teta,q,masse,ps,phis, time_0)
    289 #else
    290         ! SW model also has Earth-format start files
    291         ! (but can be used without the CPP_EARTH directive)
    292           if (iflag_phys.eq.0) then
    293             CALL dynetat0_loc("start.nc",vcov,ucov,
    294      &              teta,q,masse,ps,phis, time_0)
    295           endif
    296 #endif
    297         endif ! of if (planet_type.eq."earth")
     279!        endif ! of if (planet_type.eq."earth")
     280
    298281c       write(73,*) 'ucov',ucov
    299282c       write(74,*) 'vcov',vcov
    … …  
    337320C on remet le calendrier à zero si demande
    338321c
     322      IF (start_time /= starttime) then
     323        WRITE(lunout,*)' GCM: Attention l''heure de depart lue dans le'
     324     &,' fichier restart ne correspond pas à celle lue dans le run.def'
     325        IF (raz_date == 1) then
     326          WRITE(lunout,*)'Je prends l''heure lue dans run.def'
     327          start_time = starttime
     328        ELSE
     329          WRITE(lunout,*)'Je m''arrete'
     330          CALL abort
     331        ENDIF
     332      ENDIF
    339333      IF (raz_date == 1) THEN
    340334        annee_ref = anneeref
    … …  
    404398#endif
    405399
    406 c  nombre d'etats dans les fichiers demarrage et histoire
    407       nbetatdem = nday / iecri
    408       nbetatmoy = nday / periodav + 1
    409400
    410401c-----------------------------------------------------------------------
    411402c   Initialisation des constantes dynamiques :
    412403c   ------------------------------------------
    413       dtvr = zdtvr
    414       CALL iniconst
     404        dtvr = zdtvr
     405        CALL iniconst
    415406
    416407c-----------------------------------------------------------------------
    417408c   Initialisation de la geometrie :
    418409c   --------------------------------
    419       CALL inigeom
     410        CALL inigeom
    420411
    421412c-----------------------------------------------------------------------
    422413c   Initialisation du filtre :
    423414c   --------------------------
    424       CALL inifilr
     415        CALL inifilr
    425416c
    426417c-----------------------------------------------------------------------
    … …  
    434425c   Initialisation de la physique :
    435426c   -------------------------------
    436       IF (call_iniphys.and.iflag_phys.eq.1) THEN
     427      IF (call_iniphys.and.(iflag_phys==1.or.iflag_phys>=100)) THEN
    437428         latfi(1)=rlatu(1)
    438429         lonfi(1)=0.
    … …  
    455446         WRITE(lunout,*)
    456447     .       'GCM: WARNING!!! vitesse verticale nulle dans la physique'
    457 ! Earth:
    458          if (planet_type.eq."earth") then
    459 #ifdef CPP_EARTH
    460          CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys ,
    461      ,                latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp     )
    462 #endif
    463          endif ! of if (planet_type.eq."earth")
     448! Physics:
     449#ifdef CPP_PHYS
     450         CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys,
     451     &                latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp,
     452     &                iflag_phys)
     453#endif
    464454         call_iniphys=.false.
    465       ENDIF ! of IF (call_iniphys.and.(iflag_phys.eq.1))
     455      ENDIF ! of IF (call_iniphys.and.(iflag_phys==1.or.iflag_phys>=100))
    466456
    467457
    … …  
    469459c   Initialisation des dimensions d'INCA :
    470460c   --------------------------------------
    471       IF (config_inca /= 'none') THEN
     461      IF (type_trac == 'inca') THEN
    472462!$OMP PARALLEL
    473463#ifdef INCA
    … …  
    496486#endif
    497487
    498       if (planet_type.eq."earth") then
     488!      if (planet_type.eq."earth") then
     489! Write an Earth-format restart file
    499490        CALL dynredem0_loc("restart.nc", day_end, phis)
    500       endif
     491!      endif
    501492
    502493      ecripar = .TRUE.
    … …  
    544535c       write(78,*) 'q',q
    545536
    546 c$OMP PARALLEL DEFAULT(SHARED) COPYIN(/temps/,/logic/)
     537c$OMP PARALLEL DEFAULT(SHARED) COPYIN(/temps/,/logici/,/logicl/)
    547538      CALL leapfrog_loc(ucov,vcov,teta,ps,masse,phis,q,clesphy0,
    548539     .              time_0)

  • LMDZ5/trunk/libf/dyn3dmem/gr_dyn_fi_p.F

    r1632 r1673  
    11!
    2 ! $Id: gr_dyn_fi_p.F 1279 2009-12-10 09:02:56Z fairhead $
     2! $Id$
    33!
    44      SUBROUTINE gr_dyn_fi_p(nfield,im,jm,ngrid,pdyn,pfi)
    5 #ifdef CPP_EARTH
     5#ifdef CPP_PHYS
    66! Interface with parallel physics,
    7 ! for now this routine only works with Earth physics
    87      USE mod_interface_dyn_phys
    98      USE dimphy
    … …  
    4039      ENDDO
    4140c$OMP END DO NOWAIT
    42 #else
    43       write(lunout,*) "gr_fi_dyn_p : This routine should not be called",
    44      &   "without parallelized physics"
    45       stop
    4641#endif
    47 ! of #ifdef CPP_EARTH
     42! of #ifdef CPP_PHYS
    4843      RETURN
    4944      END

  • LMDZ5/trunk/libf/dyn3dmem/gr_fi_dyn_p.F

    r1632 r1673  
    11!
    2 ! $Id: gr_fi_dyn_p.F 1279 2009-12-10 09:02:56Z fairhead $
     2! $Id$
    33!
    44      SUBROUTINE gr_fi_dyn_p(nfield,ngrid,im,jm,pfi,pdyn)
    5 #ifdef CPP_EARTH
     5#ifdef CPP_PHYS
    66! Interface with parallel physics,
    7 ! for now this routine only works with Earth physics
    87      USE mod_interface_dyn_phys
    98      USE dimphy
    … …  
    5251      ENDDO
    5352c$OMP END DO NOWAIT
    54 #else
    55       write(lunout,*) "gr_fi_dyn_p : This routine should not be called",
    56      &   "without parallelized physics"
    57       stop
    5853#endif
    59 ! of #ifdef CPP_EARTH
     54! of #ifdef CPP_PHYS
    6055      RETURN
    6156      END

  • LMDZ5/trunk/libf/dyn3dmem/grid_noro.F

    r1658 r1673  
    11!
    2 ! $Id: grid_noro.F 1403 2010-07-01 09:02:53Z fairhead $
     2! $Id$
    33!
    44c
    … …  
    458458C MAKE A MOVING AVERAGE OVER 9 GRIDPOINTS OF THE X FIELDS
    459459
    460       PARAMETER (ISMo=300,JSMo=200)
    461       REAL X(IMAR,JMAR),XF(ISMo,JSMo)
     460      REAL X(IMAR,JMAR),XF(IMAR,JMAR)
    462461      real WEIGHTpb(-1:1,-1:1)
    463462
    464       if(imar.gt.ismo) stop'surdimensionner ismo dans mva9 (grid_noro)'
    465       if(jmar.gt.jsmo) stop'surdimensionner jsmo dans mva9 (grid_noro)'
    466463     
    467464      SUM=0.

  • LMDZ5/trunk/libf/dyn3dmem/guide_loc_mod.F90

    r1632 r1673  
    467467!       Calcul niveaux pression milieu de couches
    468468        CALL pression_loc( ijnb_u, ap, bp, ps, p )
    469         if (disvert_type==1) then
     469        if (pressure_exner) then
    470470          CALL exner_hyb_loc(ip1jmp1,ps,p,alpha,beta,pks,pk,pkf)
    471471        else

  • LMDZ5/trunk/libf/dyn3dmem/infotrac.F90

    r1632 r1673  
    3232  SUBROUTINE infotrac_init
    3333    USE control_mod
     34#ifdef REPROBUS
     35    USE CHEM_REP, ONLY : Init_chem_rep_trac
     36#endif
    3437    IMPLICIT NONE
    3538!=======================================================================
    … …  
    6164    CHARACTER(len=1), DIMENSION(3)  :: txts
    6265    CHARACTER(len=2), DIMENSION(9)  :: txtp
    63     CHARACTER(len=13)               :: str1,str2
     66    CHARACTER(len=23)               :: str1,str2
    6467 
    6568    INTEGER :: nqtrue  ! number of tracers read from tracer.def, without higer order of moment
    6669    INTEGER :: iq, new_iq, iiq, jq, ierr
    67     INTEGER, EXTERNAL :: lnblnk
    68  
     70
     71    character(len=*),parameter :: modname="infotrac_init"
    6972!-----------------------------------------------------------------------
    7073! Initialization :
    … …  
    8588   
    8689
    87     IF (config_inca=='none') THEN
    88        type_trac='lmdz'
     90    ! Coherence test between parameter type_trac, config_inca and preprocessing keys
     91    IF (type_trac=='inca') THEN
     92       WRITE(lunout,*) 'You have choosen to couple with INCA chemestry model : type_trac=', &
     93            type_trac,' config_inca=',config_inca
     94       IF (config_inca/='aero' .AND. config_inca/='chem') THEN
     95          WRITE(lunout,*) 'Incoherence between type_trac and config_inca. Model stops. Modify run.def'
     96          CALL abort_gcm('infotrac_init','Incoherence between type_trac and config_inca',1)
     97       END IF
     98#ifndef INCA
     99       WRITE(lunout,*) 'To run this option you must add cpp key INCA and compile with INCA code'
     100       CALL abort_gcm('infotrac_init','You must compile with cpp key INCA',1)
     101#endif
     102    ELSE IF (type_trac=='repr') THEN
     103       WRITE(lunout,*) 'You have choosen to couple with REPROBUS chemestry model : type_trac=', type_trac
     104#ifndef REPROBUS
     105       WRITE(lunout,*) 'To run this option you must add cpp key REPROBUS and compile with REPRPBUS code'
     106       CALL abort_gcm('infotrac_init','You must compile with cpp key REPROBUS',1)
     107#endif
     108    ELSE IF (type_trac == 'lmdz') THEN
     109       WRITE(lunout,*) 'Tracers are treated in LMDZ only : type_trac=', type_trac
    89110    ELSE
    90        type_trac='inca'
    91     END IF
     111       WRITE(lunout,*) 'type_trac=',type_trac,' not possible. Model stops'
     112       CALL abort_gcm('infotrac_init','bad parameter',1)
     113    END IF
     114
     115
     116    ! Test if config_inca is other then none for run without INCA
     117    IF (type_trac/='inca' .AND. config_inca/='none') THEN
     118       WRITE(lunout,*) 'config_inca will now be changed to none as you do not couple with INCA model'
     119       config_inca='none'
     120    END IF
     121
    92122
    93123!-----------------------------------------------------------------------
    … …  
    97127!
    98128!-----------------------------------------------------------------------
    99     IF (type_trac == 'lmdz') THEN
     129    IF (type_trac == 'lmdz' .OR. type_trac == 'repr') THEN
    100130       OPEN(90,file='traceur.def',form='formatted',status='old', iostat=ierr)
    101131       IF(ierr.EQ.0) THEN
    102           WRITE(lunout,*) 'Open traceur.def : ok'
     132          WRITE(lunout,*) trim(modname),': Open traceur.def : ok'
    103133          READ(90,*) nqtrue
    104134       ELSE
    105           WRITE(lunout,*) 'Problem in opening traceur.def'
    106           WRITE(lunout,*) 'ATTENTION using defaut values'
    107           nqtrue=4 ! Defaut value
    108        END IF
    109        ! Attention! Only for planet_type=='earth'
    110        nbtr=nqtrue-2
    111     ELSE
    112        ! nbtr has been read from INCA by init_cont_lmdz() in gcm.F
     135          WRITE(lunout,*) trim(modname),': Problem in opening traceur.def'
     136          WRITE(lunout,*) trim(modname),': WARNING using defaut values'
     137          if (planet_type=='earth') then
     138            nqtrue=4 ! Default value for Earth
     139          else
     140            nqtrue=1 ! Default value for other planets
     141          endif
     142       END IF
     143       if ( planet_type=='earth') then
     144         ! For Earth, water vapour & liquid tracers are not in the physics
     145         nbtr=nqtrue-2
     146       else
     147         ! Other planets (for now); we have the same number of tracers
     148         ! in the dynamics than in the physics
     149         nbtr=nqtrue
     150       endif
     151    ELSE ! type_trac=inca
     152       ! nbtr has been read from INCA by init_const_lmdz() in gcm.F
    113153       nqtrue=nbtr+2
    114154    END IF
    115155
    116     IF (nqtrue < 2) THEN
    117        WRITE(lunout,*) 'nqtrue=',nqtrue, ' is not allowded. 2 tracers is the minimum'
     156    IF ((planet_type=="earth").and.(nqtrue < 2)) THEN
     157       WRITE(lunout,*) trim(modname),': nqtrue=',nqtrue, ' is not allowded. 2 tracers is the minimum'
    118158       CALL abort_gcm('infotrac_init','Not enough tracers',1)
    119159    END IF
     160   
     161! Transfert number of tracers to Reprobus
     162    IF (type_trac == 'repr') THEN
     163#ifdef REPROBUS
     164       CALL Init_chem_rep_trac(nbtr)
     165#endif
     166    END IF
     167       
    120168!
    121169! Allocate variables depending on nqtrue and nbtr
    … …  
    152200!    Get choice of advection schema from file tracer.def or from INCA
    153201!---------------------------------------------------------------------
    154     IF (type_trac == 'lmdz') THEN
     202    IF (type_trac == 'lmdz' .OR. type_trac == 'repr') THEN
    155203       IF(ierr.EQ.0) THEN
    156204          ! Continue to read tracer.def
    157205          DO iq=1,nqtrue
    158              READ(90,999) hadv(iq),vadv(iq),tnom_0(iq)
     206             READ(90,*) hadv(iq),vadv(iq),tnom_0(iq)
    159207          END DO
    160208          CLOSE(90) 
    161        ELSE ! Without tracer.def
     209       ELSE ! Without tracer.def, set default values
     210         if (planet_type=="earth") then
     211          ! for Earth, default is to have 4 tracers
    162212          hadv(1) = 14
    163213          vadv(1) = 14
    … …  
    172222          vadv(4) = 10
    173223          tnom_0(4) = 'PB'
     224         else ! default for other planets
     225          hadv(1) = 10
     226          vadv(1) = 10
     227          tnom_0(1) = 'dummy'
     228         endif ! of if (planet_type=="earth")
    174229       END IF
    175230       
    176        WRITE(lunout,*) 'Valeur de traceur.def :'
    177        WRITE(lunout,*) 'nombre de traceurs ',nqtrue
     231       WRITE(lunout,*) trim(modname),': Valeur de traceur.def :'
     232       WRITE(lunout,*) trim(modname),': nombre de traceurs ',nqtrue
    178233       DO iq=1,nqtrue
    179234          WRITE(lunout,*) hadv(iq),vadv(iq),tnom_0(iq)
    … …  
    217272          new_iq=new_iq+10 ! 9 tracers added
    218273       ELSE
    219           WRITE(lunout,*) 'This choice of advection schema is not available'
     274          WRITE(lunout,*) trim(modname),': This choice of advection schema is not available',iq,hadv(iq),vadv(iq)
    220275          CALL abort_gcm('infotrac_init','Bad choice of advection schema - 1',1)
    221276       END IF
    … …  
    227282       nqtot = new_iq
    228283
    229        WRITE(lunout,*) 'The choice of advection schema for one or more tracers'
     284       WRITE(lunout,*) trim(modname),': The choice of advection schema for one or more tracers'
    230285       WRITE(lunout,*) 'makes it necessary to add tracers'
    231        WRITE(lunout,*) nqtrue,' is the number of true tracers'
    232        WRITE(lunout,*) nqtot, ' is the total number of tracers needed'
     286       WRITE(lunout,*) trim(modname)//': ',nqtrue,' is the number of true tracers'
     287       WRITE(lunout,*) trim(modname)//': ',nqtot, ' is the total number of tracers needed'
    233288
    234289    ELSE
    … …  
    258313          iadv(new_iq)=11
    259314       ELSE
    260           WRITE(lunout,*)'This choice of advection schema is not available'
     315          WRITE(lunout,*)trim(modname),': This choice of advection schema is not available',iq,hadv(iq),vadv(iq)
     316
    261317          CALL abort_gcm('infotrac_init','Bad choice of advection schema - 2',1)
    262318       END IF
    … …  
    265321       tname(new_iq)= tnom_0(iq)
    266322       IF (iadv(new_iq)==0) THEN
    267           ttext(new_iq)=str1(1:lnblnk(str1))
     323          ttext(new_iq)=trim(str1)
    268324       ELSE
    269           ttext(new_iq)=str1(1:lnblnk(str1))//descrq(iadv(new_iq))
     325          ttext(new_iq)=trim(tnom_0(iq))//descrq(iadv(new_iq))
    270326       END IF
    271327
    … …  
    276332             new_iq=new_iq+1
    277333             iadv(new_iq)=-20
    278              ttext(new_iq)=str2(1:lnblnk(str2))//txts(jq)
    279              tname(new_iq)=str1(1:lnblnk(str1))//txts(jq)
     334             ttext(new_iq)=trim(str2)//txts(jq)
     335             tname(new_iq)=trim(str1)//txts(jq)
    280336          END DO
    281337       ELSE IF (iadv(new_iq)==30) THEN
    … …  
    283339             new_iq=new_iq+1
    284340             iadv(new_iq)=-30
    285              ttext(new_iq)=str2(1:lnblnk(str2))//txtp(jq)
    286              tname(new_iq)=str1(1:lnblnk(str1))//txtp(jq)
     341             ttext(new_iq)=trim(str2)//txtp(jq)
     342             tname(new_iq)=trim(str1)//txtp(jq)
    287343          END DO
    288344       END IF
    … …  
    303359
    304360
    305     WRITE(lunout,*) 'Information stored in infotrac :'
    306     WRITE(lunout,*) 'iadv  niadv tname  ttext :'
     361    WRITE(lunout,*) trim(modname),': Information stored in infotrac :'
     362    WRITE(lunout,*) trim(modname),': iadv  niadv tname  ttext :'
    307363    DO iq=1,nqtot
    308        WRITE(lunout,*) iadv(iq),niadv(iq), tname(iq), ttext(iq)
     364       WRITE(lunout,*) iadv(iq),niadv(iq),&
     365       ' ',trim(tname(iq)),' ',trim(ttext(iq))
    309366    END DO
    310367
    … …  
    315372    DO iq=1,nqtot
    316373       IF (iadv(iq)/=10 .AND. iadv(iq)/=14 .AND. iadv(iq)/=0) THEN
    317           WRITE(lunout,*)'STOP : The option iadv=',iadv(iq),' is not tested in this version of LMDZ'
     374          WRITE(lunout,*)trim(modname),' STOP : The option iadv=',iadv(iq),' is not tested in this version of LMDZ'
    318375          CALL abort_gcm('infotrac_init','In this version only iadv=10 and iadv=14 is tested!',1)
    319376       ELSE IF (iadv(iq)==14 .AND. iq/=1) THEN
    320           WRITE(lunout,*)'STOP : The option iadv=',iadv(iq),' is not tested in this version of LMDZ'
     377          WRITE(lunout,*)trim(modname),'STOP : The option iadv=',iadv(iq),' is not tested in this version of LMDZ'
    321378          CALL abort_gcm('infotrac_init','In this version iadv=14 is only permitted for water vapour!',1)
    322379       END IF
    … …  
    329386    DEALLOCATE(tracnam)
    330387
    331 999 FORMAT (i2,1x,i2,1x,a15)
    332 
    333388  END SUBROUTINE infotrac_init
    334389

  • LMDZ5/trunk/libf/dyn3dmem/inigrads.F

    r1632 r1673  
    99      implicit none
    1010
    11       integer if,im,jm,lm,i,j,l,lnblnk
     11      integer if,im,jm,lm,i,j,l
    1212      real x(im),y(jm),z(lm),fx,fy,fz,dt
    1313      real xmin,xmax,ymin,ymax
    … …  
    4040      ivar(if)=0
    4141
    42       fichier(if)=file(1:lnblnk(file))
     42      fichier(if)=trim(file)
    4343
    4444      firsttime(if)=.true.
    … …  
    7070
    7171      print*,4*(ifd(if)-iid(if))*(jfd(if)-jid(if))
    72       print*,file(1:lnblnk(file))//'.dat'
     72      print*,trim(file)//'.dat'
    7373
    74       OPEN (unit(if)+1,FILE=file(1:lnblnk(file))//'.dat'
     74      OPEN (unit(if)+1,FILE=trim(file)//'.dat'
    7575     s   ,FORM='unformatted',
    7676     s   ACCESS='direct'

  • LMDZ5/trunk/libf/dyn3dmem/initfluxsto_p.F

    r1632 r1673  
    11!
    2 ! $Id: initfluxsto_p.F 1279 2009-12-10 09:02:56Z fairhead $
     2! $Id$
    33!
    44      subroutine initfluxsto_p
    … …  
    203203     .              llm, nivsigs, zvertiid)
    204204c pour le fichier def
    205       nivd(1) = 1
    206       call histvert(filedid, 'sig_s', 'Niveaux sigma',
    207      .  'sigma_level',
    208      .              1, nivd, dvertiid)
    209 
     205      if (mpi_rank==0) then
     206         nivd(1) = 1
     207         call histvert(filedid, 'sig_s', 'Niveaux sigma',
     208     .        'sigma_level',
     209     .        1, nivd, dvertiid)
     210      endif
    210211C
    211212C  Appels a histdef pour la definition des variables a sauvegarder
    … …  
    282283      call histend(fileid)
    283284      call histend(filevid)
    284       call histend(filedid)
     285      if (mpi_rank==0) call histend(filedid)
    285286      if (ok_sync) then
    286287        call histsync(fileid)
    287288        call histsync(filevid)
    288         call histsync(filedid)
     289        if (mpi_rank==0) call histsync(filedid)
    289290      endif
    290291       

  • LMDZ5/trunk/libf/dyn3dmem/integrd_loc.F

    r1632 r1673  
    44      SUBROUTINE integrd_loc
    55     $  (  nq,vcovm1,ucovm1,tetam1,psm1,massem1,
    6      $     dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps0,masse,phis,finvmaold)
     6     $     dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps0,masse,phis) !,finvmaold)
    77      USE parallel
    88      USE control_mod
    … …  
    3737#include "temps.h"
    3838#include "serre.h"
     39#include "iniprint.h"
    3940      include 'mpif.h'
    4041
    … …  
    4243c   ----------
    4344
    44       INTEGER nq
    45 
    46       REAL vcov(ijb_v:ije_v,llm),ucov(ijb_u:ije_u,llm)
    47       REAL teta(ijb_u:ije_u,llm)
    48       REAL q(ijb_u:ije_u,llm,nq)
    49       REAL ps0(ijb_u:ije_u),masse(ijb_u:ije_u,llm),phis(ijb_u:ije_u)
    50 
    51       REAL vcovm1(ijb_v:ije_v,llm),ucovm1(ijb_u:ije_u,llm)
    52       REAL tetam1(ijb_u:ije_u,llm),psm1(ijb_u:ije_u)
    53       REAL massem1(ijb_u:ije_u,llm)
    54 
    55       REAL dv(ijb_v:ije_v,llm),du(ijb_u:ije_u,llm)
    56       REAL dteta(ijb_u:ije_u,llm),dp(ijb_u:ije_u)
    57       REAL dq(ijb_u:ije_u,llm,nq), finvmaold(ijb_u:ije_u,llm)
     45      INTEGER,intent(in) :: nq ! number of tracers to handle in this routine
     46
     47      REAL,INTENT(INOUT) :: vcov(ijb_v:ije_v,llm) ! covariant meridional wind
     48      REAL,INTENT(INOUT) :: ucov(ijb_u:ije_u,llm) ! covariant zonal wind
     49      REAL,INTENT(INOUT) :: teta(ijb_u:ije_u,llm) ! potential temperature
     50      REAL,INTENT(INOUT) :: q(ijb_u:ije_u,llm,nq) ! advected tracers
     51      REAL,INTENT(INOUT) :: ps0(ijb_u:ije_u) ! surface pressure
     52      REAL,INTENT(INOUT) :: masse(ijb_u:ije_u,llm) ! atmospheric mass
     53      REAL,INTENT(INOUT) :: phis(ijb_u:ije_u) ! ground geopotential !!! unused
     54      ! values at previous time step
     55      REAL,INTENT(INOUT) :: vcovm1(ijb_v:ije_v,llm)
     56      REAL,INTENT(INOUT) :: ucovm1(ijb_u:ije_u,llm)
     57      REAL,INTENT(INOUT) :: tetam1(ijb_u:ije_u,llm)
     58      REAL,INTENT(INOUT) :: psm1(ijb_u:ije_u)
     59      REAL,INTENT(INOUT) :: massem1(ijb_u:ije_u,llm)
     60      ! the tendencies to add
     61      REAL,INTENT(INOUT) :: dv(ijb_v:ije_v,llm)
     62      REAL,INTENT(INOUT) :: du(ijb_u:ije_u,llm)
     63      REAL,INTENT(INOUT) :: dteta(ijb_u:ije_u,llm)
     64      REAL,INTENT(INOUT) :: dp(ijb_u:ije_u)
     65      REAL,INTENT(INOUT) :: dq(ijb_u:ije_u,llm,nq) !!! unused
     66!      REAL,INTENT(INOUT) ::finvmaold(ijb_u:ije_u,llm) !!! unused
    5867
    5968c   Local:
    … …  
    6271      REAL vscr( ijb_v:ije_v ),uscr( ijb_u:ije_u )
    6372      REAL hscr( ijb_u:ije_u ),pscr(ijb_u:ije_u)
    64       REAL massescr( ijb_u:ije_u,llm ), finvmasse(ijb_u:ije_u,llm)
     73      REAL massescr( ijb_u:ije_u,llm )
     74!      REAL finvmasse(ijb_u:ije_u,llm)
    6575      REAL tpn,tps,tppn(iim),tpps(iim)
    6676      REAL qpn,qps,qppn(iim),qpps(iim)
    6777
    68       INTEGER  l,ij,iq
     78      INTEGER  l,ij,iq,i,j
    6979
    7080      REAL SSUM
    … …  
    7484      LOGICAL,SAVE :: checksum_all=.TRUE.
    7585      INTEGER :: stop_it
    76       INTEGER :: ierr,j
     86      INTEGER :: ierr
    7787
    7888c-----------------------------------------------------------------------
    … …  
    137147!     &                   MPI_LOGICAL,MPI_LOR,COMM_LMDZ,ierr)
    138148      IF( .NOT. checksum ) THEN
    139          PRINT*,' Au point ij = ',stop_it, ' , pression sol neg. '
    140      &         , ps(stop_it)
    141          STOP' dans integrd'
    142       ENDIF
     149          write(lunout,*) "integrd: negative surface pressure ",
     150     &                                                ps(stop_it)
     151         write(lunout,*) " at node ij =", stop_it
     152         ! since ij=j+(i-1)*jjp1 , we have
     153!         j=modulo(stop_it,jjp1)
     154!         i=1+(stop_it-j)/jjp1
     155!         write(lunout,*) " lon = ",rlonv(i)*180./pi, " deg",
     156!     &                   " lat = ",rlatu(j)*180./pi, " deg"
     157      ENDIF
     158
    143159c$OMP END MASTER
    144160c$OMP BARRIER
    … …  
    160176        call WriteField_u('int_dteta',dteta)
    161177        call WriteField_u('int_dp',dp)
    162         call WriteField_u('int_finvmaold',finvmaold)
     178!        call WriteField_u('int_finvmaold',finvmaold)
    163179        do j=1,nq
    164180          call WriteField_u('int_q'//trim(int2str(j)),
    … …  
    206222      CALL massdair_loc (     p  , masse         )
    207223
    208 c      CALL   SCOPY( ijp1llm  , masse, 1, finvmasse,  1      )
    209       ijb=ij_begin
    210       ije=ij_end
    211      
    212 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    213       DO  l = 1,llm
    214         finvmasse(ijb:ije,l)=masse(ijb:ije,l)
    215       ENDDO
    216 c$OMP END DO NOWAIT
    217 
    218       jjb=jj_begin
    219       jje=jj_end
    220       CALL filtreg_p( finvmasse,jjb_u,jje_u,jjb,jje, jjp1, llm,
    221      &                -2, 2, .TRUE., 1  )
     224! Ehouarn : we don't use/need finvmaold and finvmasse,
     225!           so might as well not compute them
     226!c      CALL   SCOPY( ijp1llm  , masse, 1, finvmasse,  1      )
     227!      ijb=ij_begin
     228!      ije=ij_end
     229!     
     230!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     231!      DO  l = 1,llm
     232!        finvmasse(ijb:ije,l)=masse(ijb:ije,l)
     233!      ENDDO
     234!c$OMP END DO NOWAIT
     235
     236!      jjb=jj_begin
     237!      jje=jj_end
     238!      CALL filtreg_p( finvmasse,jjb_u,jje_u,jjb,jje, jjp1, llm,
     239!     &                -2, 2, .TRUE., 1  )
    222240c
    223241
    … …  
    320338
    321339          CALL qminimum_loc( q, nq, deltap )
    322          endif ! of if (planet_type.eq."earth")
    323340c
    324341c    .....  Calcul de la valeur moyenne, unique  aux poles pour  q .....
    … …  
    371388      ENDIF
    372389     
    373 c         CALL  SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
    374 
    375 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    376       DO l = 1, llm     
    377         finvmaold(ijb:ije,l)=finvmasse(ijb:ije,l)       
    378       ENDDO
    379 c$OMP END DO NOWAIT
     390! Ehouarn: forget about finvmaold
     391!c         CALL  SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
     392
     393!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     394!      DO l = 1, llm     
     395!        finvmaold(ijb:ije,l)=finvmasse(ijb:ije,l)       
     396!      ENDDO
     397!c$OMP END DO NOWAIT
     398
     399      endif ! of if (planet_type.eq."earth")
     400
    380401c
    381402c

  • LMDZ5/trunk/libf/dyn3dmem/leapfrog_loc.F

    r1659 r1673  
    11!
    2 ! $Id: leapfrog_p.F 1299 2010-01-20 14:27:21Z fairhead $
     2! $Id$
    33!
    44c
    … …  
    119119
    120120c   tendances physiques
    121 !      REAL,SAVE,ALLOCATABLE :: dvfi(:,:),dufi(:,:)
    122 !      REAL,SAVE,ALLOCATABLE :: dtetafi(:,:)
    123 !      REAL,SAVE,ALLOCATABLE :: dpfi(:)
    124 !      REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: dqfi
     121      REAL,SAVE,ALLOCATABLE :: dvfi(:,:),dufi(:,:)
     122      REAL,SAVE,ALLOCATABLE :: dtetafi(:,:)
     123      REAL,SAVE,ALLOCATABLE :: dpfi(:)
     124      REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: dqfi
    125125
    126126c   variables pour le fichier histoire
    … …  
    150150      REAL :: secondes
    151151
     152      logical :: physic
    152153      LOGICAL first,callinigrads
    153154
    … …  
    174175
    175176      character*80 dynhist_file, dynhistave_file
    176       character*20 modname
     177      character(len=*),parameter :: modname="leapfrog_loc"
    177178      character*80 abort_message
    178179
    … …  
    195196
    196197      INTEGER :: true_itau
    197       LOGICAL :: verbose=.true.
    198198      INTEGER :: iapptrac
    199199      INTEGER :: AdjustCount
    … …  
    215215      itaufin   = nday*day_step
    216216      itaufinp1 = itaufin +1
    217       modname="leapfrog_p"
    218217
    219218      itau = 0
     219      physic=.true.
     220      if (iflag_phys==0.or.iflag_phys==2) physic=.false.
    220221      CALL init_nan
    221222      CALL leapfrog_allocate
    … …  
    252253!      ALLOCATE(dvdis(ijb_v:ije_v,llm),dudis(ijb_u:ije_u,llm))
    253254!      ALLOCATE(dtetadis(ijb_u:ije_u,llm))
    254 !      ALLOCATE(dvfi(ijb_v:ije_v,llm),dufi(ijb_u:ije_u,llm))
    255 !      ALLOCATE(dtetafi(ijb_u:ije_u,llm))
    256 !      ALLOCATE(dpfi(ijb_u:ije_u))
     255      ALLOCATE(dvfi(ijb_v:ije_v,llm),dufi(ijb_u:ije_u,llm))
     256      ALLOCATE(dtetafi(ijb_u:ije_u,llm))
     257      ALLOCATE(dpfi(ijb_u:ije_u))
    257258!      ALLOCATE(dq(ijb_u:ije_u,llm,nqtot))
    258 !      ALLOCATE(dqfi(ijb_u:ije_u,llm,nqtot))
     259      ALLOCATE(dqfi(ijb_u:ije_u,llm,nqtot))
    259260!      ALLOCATE(dqfi_tmp(iip1,llm,nqtot))
    260261!      ALLOCATE(finvmaold(ijb_u:ije_u,llm))
    … …  
    277278
    278279c$OMP MASTER
    279       dq=0.
     280      dq(:,:,:)=0.
    280281      CALL pression ( ijnb_u, ap, bp, ps, p       )
    281282c$OMP END MASTER
     283      if (pressure_exner) then
    282284      CALL exner_hyb_loc( ijnb_u, ps, p,alpha,beta, pks, pk, pkf)
    283 
     285      else
     286        CALL exner_milieu_loc( ijnb_u, ps, p, beta, pks, pk, pkf )
     287      endif
    284288c-----------------------------------------------------------------------
    285289c   Debut de l'integration temporelle:
    … …  
    287291c et du parallelisme !!
    288292
    289    1  CONTINUE
    290 
    291       jD_cur = jD_ref + day_ini - day_ref + int (itau * dtvr / daysec)
    292       jH_cur = jH_ref +                                                 &
    293      &          (itau * dtvr / daysec - int(itau * dtvr / daysec))
     293   1  CONTINUE ! Matsuno Forward step begins here
     294
     295      jD_cur = jD_ref + day_ini - day_ref +                             &
     296     &          itau/day_step
     297      jH_cur = jH_ref + start_time +                                    &
     298     &         mod(itau,day_step)/float(day_step)
     299      if (jH_cur > 1.0 ) then
     300        jD_cur = jD_cur +1.
     301        jH_cur = jH_cur -1.
     302      endif
     303
    294304
    295305#ifdef CPP_IOIPSL
    … …  
    323333         psm1= ps
    324334         
    325          finvmaold = masse
    326 c$OMP END MASTER
    327 c$OMP BARRIER
    328          CALL filtreg_p ( finvmaold ,jjb_u,jje_u,jjb_u,jje_u,jjp1, llm,
    329      &                    -2,2, .TRUE., 1 )
     335! Ehouarn: finvmaold is actually not used       
     336!         finvmaold = masse
     337c$OMP END MASTER
     338c$OMP BARRIER
     339!         CALL filtreg_p ( finvmaold ,jjb_u,jje_u,jjb_u,jje_u,jjp1, llm,
     340!     &                    -2,2, .TRUE., 1 )
    330341       else
    331342! Save fields obtained at previous time step as '...m1'
    … …  
    343354           tetam1   (ijb:ije,l) = teta  (ijb:ije,l)
    344355           massem1  (ijb:ije,l) = masse (ijb:ije,l)
    345            finvmaold(ijb:ije,l)=masse(ijb:ije,l)
     356!           finvmaold(ijb:ije,l)=masse(ijb:ije,l)
    346357                 
    347358           if (pole_sud) ije=ij_end-iip1
    … …  
    353364
    354365
    355           CALL filtreg_p(finvmaold ,jjb_u,jje_u,jj_begin,jj_end,jjp1,
    356      .                    llm, -2,2, .TRUE., 1 )
     366! Ehouarn: finvmaold not used
     367!          CALL filtreg_p(finvmaold ,jjb_u,jje_u,jj_begin,jj_end,jjp1,
     368!     .                    llm, -2,2, .TRUE., 1 )
    357369
    358370       endif ! of if (FirstCaldyn)
    … …  
    370382cym      call minmax(ijp1llm,q(:,:,3),zqmin,zqmax)
    371383
    372    2  CONTINUE
     384   2  CONTINUE ! Matsuno backward or leapfrog step begins here
    373385
    374386c$OMP MASTER
    … …  
    399411      ! Purely Matsuno time stepping
    400412         IF( MOD(itau,iconser) .EQ.0.AND.  forward    ) conser = .TRUE.
    401          IF( MOD(itau,idissip ).EQ.0.AND..NOT.forward ) apdiss = .TRUE.
     413         IF( MOD(itau,dissip_period ).EQ.0.AND..NOT.forward )
     414     s        apdiss = .TRUE.
    402415         IF( MOD(itau,iphysiq ).EQ.0.AND..NOT.forward
    403      s          .and. iflag_phys.EQ.1                 ) apphys = .TRUE.
     416     s          .and. physic                        ) apphys = .TRUE.
    404417      ELSE
    405418      ! Leapfrog/Matsuno time stepping
    406419         IF( MOD(itau   ,iconser) .EQ. 0              ) conser = .TRUE.
    407          IF( MOD(itau+1,idissip)  .EQ. 0              ) apdiss = .TRUE.
    408          IF( MOD(itau+1,iphysiq).EQ.0.AND.iflag_phys.EQ.1) apphys=.TRUE.
     420         IF( MOD(itau+1,dissip_period).EQ.0 .AND. .NOT. forward )
     421     s        apdiss = .TRUE.
     422         IF( MOD(itau+1,iphysiq).EQ.0.AND.physic) apphys=.TRUE.
    409423      END IF
    410424
    … …  
    450464c$OMP MASTER
    451465           call allgather_timer_average
    452         verbose=.TRUE.
    453         if (Verbose) then
     466
     467        if (prt_level > 9) then
    454468       
    455469        print *,'*********************************'
    … …  
    622636      call start_timer(timer_caldyn)
    623637
     638      ! compute geopotential phi()
    624639      CALL geopot_loc  ( ip1jmp1, teta  , pk , pks,  phis  , phi   )
    625640
    … …  
    697712
    698713       CALL integrd_loc ( 2,vcovm1,ucovm1,tetam1,psm1,massem1 ,
    699      $         dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis ,
    700      $              finvmaold                                    )
     714     $         dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis)
     715!     $              finvmaold                                    )
    701716
    702717!       CALL FTRACE_REGION_END("integrd")
    … …  
    10811096!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    10821097       do l=1,llm
    1083        teta(ijb:ije,l)=teta(ijb:ije,l)
    1084      &  -iphysiq*dtvr*(teta(ijb:ije,l)-tetarappel(ijb:ije,l))/taurappel
     1098       teta(ijb:ije,l)=teta(ijb:ije,l) -dtvr*
     1099     &        (teta(ijb:ije,l)-tetarappel(ijb:ije,l))*
     1100     &                 (knewt_g+knewt_t(l)*clat4(ijb:ije))       
    10851101       enddo
    10861102!$OMP END DO
     1103
     1104!$OMP MASTER
     1105       if (planet_type.eq."giant") then
     1106         ! add an intrinsic heat flux at the base of the atmosphere
     1107         teta(ijb:ije,1) = teta(ijb:ije,1)
     1108     &        + dtvr * aire(ijb:ije) * ihf / cpp / masse(ijb:ije,1)
     1109       endif
     1110!$OMP END MASTER
     1111!$OMP BARRIER
     1112
    10871113
    10881114       call Register_Hallo_u(ucov,llm,0,1,1,0,Request_Physic)
    … …  
    10921118       call WaitRequest(Request_Physic)     
    10931119c$OMP BARRIER
    1094        call friction_loc(ucov,vcov,iphysiq*dtvr)
     1120       call friction_loc(ucov,vcov,dtvr)
    10951121!$OMP BARRIER
     1122
     1123        ! Sponge layer (if any)
     1124        IF (ok_strato) THEN
     1125          ! set dufi,dvfi,... to zero
     1126          ijb=ij_begin
     1127          ije=ij_end
     1128!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     1129          do l=1,llm
     1130            dufi(ijb:ije,l)=0
     1131            dtetafi(ijb:ije,l)=0
     1132            dqfi(ijb:ije,l,1:nqtot)=0
     1133          enddo
     1134!$OMP END DO
     1135!$OMP MASTER
     1136          dpfi(ijb:ije)=0
     1137!$OMP END MASTER
     1138          ijb=ij_begin
     1139          ije=ij_end
     1140          if (pole_sud) ije=ije-iip1
     1141!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     1142          do l=1,llm
     1143            dvfi(ijb:ije,l)=0
     1144          enddo
     1145!$OMP END DO
     1146
     1147          CALL top_bound_loc(vcov,ucov,teta,masse,dufi,dvfi,dtetafi)
     1148          CALL addfi_loc( dtvr, leapf, forward   ,
     1149     $                  ucov, vcov, teta , q   ,ps ,
     1150     $                 dufi, dvfi, dtetafi , dqfi ,dpfi  )
     1151!$OMP BARRIER
     1152        ENDIF ! of IF (ok_strato)
    10961153      ENDIF ! of IF(iflag_phys.EQ.2)
    10971154
    … …  
    10991156        CALL pression_loc ( ip1jmp1, ap, bp, ps, p                  )
    11001157c$OMP BARRIER
    1101         CALL exner_hyb_loc( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
     1158        if (pressure_exner) then
     1159        CALL exner_hyb_loc( ijnb_u, ps, p,alpha,beta, pks, pk, pkf )
     1160        else
     1161          CALL exner_milieu_loc( ijnb_u, ps, p, beta, pks, pk, pkf )
     1162        endif
    11021163c$OMP BARRIER
    11031164
    … …  
    14961557c$OMP BARRIER
    14971558
    1498               if (planet_type.eq."earth") then
     1559!              if (planet_type.eq."earth") then
    14991560! Write an Earth-format restart file
    15001561                CALL dynredem1_loc("restart.nc",0.0,
    15011562     &                           vcov,ucov,teta,q,masse,ps)
    1502               endif ! of if (planet_type.eq."earth")
     1563!              endif ! of if (planet_type.eq."earth")
    15031564
    15041565!              CLOSE(99)
    … …  
    16081669
    16091670              IF(itau.EQ.itaufin) THEN
    1610                 if (planet_type.eq."earth") then
     1671!                if (planet_type.eq."earth") then
    16111672                   CALL dynredem1_loc("restart.nc",0.0,
    16121673     .                               vcov,ucov,teta,q,masse,ps)
    1613                 endif ! of if (planet_type.eq."earth")
     1674!               endif ! of if (planet_type.eq."earth")
    16141675              ENDIF ! of IF(itau.EQ.itaufin)
    16151676

  • LMDZ5/trunk/libf/dyn3dmem/limit_netcdf.F90

    r1658 r1673  
    11!
    2 ! $Id: limit_netcdf.F90 1425 2010-09-02 13:45:23Z lguez $
     2! $Id$
    33!-------------------------------------------------------------------------------
    44!
    … …  
    4242  REAL, DIMENSION(iip1,jjp1), INTENT(IN) :: masque   ! land mask
    4343#ifndef CPP_EARTH
    44   WRITE(lunout,*)'limit_netcdf: Earth-specific routine, needs Earth physics'
     44  CALL abort_gcm('limit_netcdf','Earth-specific routine, needs Earth physics',1)
    4545#else
    4646!-------------------------------------------------------------------------------
    … …  
    5252#include "indicesol.h"
    5353
    54 !--- For fractionary sub-cell use (old coding used soil index: 0,1,2,3) --------
    55   LOGICAL, PARAMETER :: fracterre=.TRUE.
    56 
    5754!--- INPUT NETCDF FILES NAMES --------------------------------------------------
    5855  CHARACTER(LEN=25) :: icefile, sstfile, dumstr
    5956  CHARACTER(LEN=25), PARAMETER :: famipsst='amipbc_sst_1x1.nc        ',        &
    6057                                  famipsic='amipbc_sic_1x1.nc        ',        &
    61                                   fclimsst='amipbc_sst_1x1_clim.nc   ',        &
    62                                   fclimsic='amipbc_sic_1x1_clim.nc   ',        &
    6358                                  fcpldsst='cpl_atm_sst.nc           ',        &
    6459                                  fcpldsic='cpl_atm_sic.nc           ',        &
     60                                  fhistsst='histmth_sst.nc           ',        &
     61                                  fhistsic='histmth_sic.nc           ',        &
    6562                                  frugo   ='Rugos.nc                 ',        &
    6663                                  falbe   ='Albedo.nc                '
    67 
     64  CHARACTER(LEN=10) :: varname
    6865!--- OUTPUT VARIABLES FOR NETCDF FILE ------------------------------------------
    6966  REAL,   DIMENSION(klon)                :: fi_ice, verif
    … …  
    8077  INTEGER :: id_FOCE, id_FSIC, id_FTER, id_FLIC
    8178  INTEGER :: NF90_FORMAT
    82   LOGICAL :: lCPL                    !--- T: IPCC-IPSL cpl model output files
    8379  INTEGER :: ndays                   !--- Depending on the output calendar
    8480
    … …  
    9793  kappa = 0.2857143
    9894  cpp   = 1004.70885
    99   dtvr  = daysec/FLOAT(day_step)
     95  dtvr  = daysec/REAL(day_step)
    10096  CALL inigeom
    10197
    … …  
    104100
    105101!--- RUGOSITY TREATMENT --------------------------------------------------------
    106   WRITE(lunout,*) 'Traitement de la rugosite'
    107   CALL get_2Dfield(frugo,'RUG',interbar,ndays,phy_rug,mask=masque(1:iim,:))
     102  IF (prt_level>1) WRITE(lunout,*) 'Traitement de la rugosite'
     103  varname='RUGOS'
     104  CALL get_2Dfield(frugo,varname,'RUG',interbar,ndays,phy_rug,mask=masque(1:iim,:))
    108105
    109106!--- OCEAN TREATMENT -----------------------------------------------------------
    110   PRINT*, 'Traitement de la glace oceanique' ; icefile=''; lCPL=.FALSE.
     107  IF (prt_level>1) WRITE(lunout,*) 'Traitement de la glace oceanique'
    111108
    112109! Input SIC file selection
    113   icefile='fake'
    114   IF(NF90_OPEN(famipsic,NF90_NOWRITE,nid)==NF90_NOERR) icefile=TRIM(famipsic)
    115   IF(NF90_OPEN(fclimsic,NF90_NOWRITE,nid)==NF90_NOERR) icefile=TRIM(fclimsic)
    116   IF(NF90_OPEN(fcpldsic,NF90_NOWRITE,nid)==NF90_NOERR) icefile=TRIM(fcpldsic)
    117   SELECT CASE(icefile)
    118     CASE(famipsic); dumstr='Amip.'
    119     CASE(fclimsic); dumstr='Amip climatologique.'
    120     CASE(fcpldsic); dumstr='de sortie du modele couplé IPSL/IPCC.';lCPL=.TRUE.
    121     CASE('fake');   CALL abort_gcm('limit_netcdf','Fichier SIC non reconnu.',1)
    122   END SELECT
     110! Open file only to test if available
     111  IF ( NF90_OPEN(TRIM(famipsic),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
     112     icefile=TRIM(famipsic)
     113     varname='sicbcs'
     114  ELSE IF( NF90_OPEN(TRIM(fcpldsic),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
     115     icefile=TRIM(fcpldsic)
     116     varname='SIICECOV'
     117  ELSE IF ( NF90_OPEN(TRIM(fhistsic),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
     118     icefile=TRIM(fhistsic)
     119     varname='pourc_sic'
     120  ELSE
     121     WRITE(lunout,*) 'ERROR! No sea-ice input file was found.'
     122     WRITE(lunout,*) 'One of following files must be availible : ',trim(famipsic),', ',trim(fcpldsic),', ',trim(fhistsic)
     123     CALL abort_gcm('limit_netcdf','No sea-ice file was found',1)
     124  END IF
    123125  ierr=NF90_CLOSE(nid)
    124   WRITE(lunout,*)'Pour la glace de mer a ete choisi un fichier '//TRIM(dumstr)
    125 
    126   CALL get_2Dfield(icefile,'SIC',interbar,ndays,phy_ice,flag=oldice,lCPL=lCPL)
     126  IF (prt_level>=0) WRITE(lunout,*)'Pour la glace de mer a ete choisi le fichier '//TRIM(icefile)
     127
     128  CALL get_2Dfield(icefile,varname, 'SIC',interbar,ndays,phy_ice,flag=oldice)
    127129
    128130  ALLOCATE(pctsrf_t(klon,nbsrf,ndays))
    129131  DO k=1,ndays
    130     fi_ice=phy_ice(:,k)
    131     WHERE(fi_ice>=1.0  ) fi_ice=1.0
    132     WHERE(fi_ice<EPSFRA) fi_ice=0.0
    133     IF(fracterre) THEN
    134       pctsrf_t(:,is_ter,k)=pctsrf(:,is_ter)       ! land soil
    135       pctsrf_t(:,is_lic,k)=pctsrf(:,is_lic)       ! land ice
    136       IF(lCPL) THEN                               ! SIC=pICE*(1-LIC-TER)
    137         pctsrf_t(:,is_sic,k)=fi_ice*(1-pctsrf(:,is_lic)-pctsrf(:,is_ter))
    138       ELSE                                        ! SIC=pICE-LIC
     132     fi_ice=phy_ice(:,k)
     133     WHERE(fi_ice>=1.0  ) fi_ice=1.0
     134     WHERE(fi_ice<EPSFRA) fi_ice=0.0
     135     pctsrf_t(:,is_ter,k)=pctsrf(:,is_ter)       ! land soil
     136     pctsrf_t(:,is_lic,k)=pctsrf(:,is_lic)       ! land ice
     137     IF (icefile==trim(fcpldsic)) THEN           ! SIC=pICE*(1-LIC-TER)
     138        pctsrf_t(:,is_sic,k)=fi_ice(:)*(1.-pctsrf(:,is_lic)-pctsrf(:,is_ter))
     139     ELSE IF (icefile==trim(fhistsic)) THEN      ! SIC=pICE
     140        pctsrf_t(:,is_sic,k)=fi_ice(:)
     141     ELSE ! icefile==famipsic                    ! SIC=pICE-LIC
    139142        pctsrf_t(:,is_sic,k)=fi_ice-pctsrf_t(:,is_lic,k)
    140       END IF
    141       WHERE(pctsrf_t(:,is_sic,k)<=0) pctsrf_t(:,is_sic,k)=0.
    142       WHERE(1.0-zmasq<EPSFRA)
     143     END IF
     144     WHERE(pctsrf_t(:,is_sic,k)<=0) pctsrf_t(:,is_sic,k)=0.
     145     WHERE(1.0-zmasq<EPSFRA)
    143146        pctsrf_t(:,is_sic,k)=0.0
    144147        pctsrf_t(:,is_oce,k)=0.0
    145       ELSEWHERE
     148     ELSEWHERE
    146149        WHERE(pctsrf_t(:,is_sic,k)>=1.0-zmasq)
    147           pctsrf_t(:,is_sic,k)=1.0-zmasq
    148           pctsrf_t(:,is_oce,k)=0.0
     150           pctsrf_t(:,is_sic,k)=1.0-zmasq
     151           pctsrf_t(:,is_oce,k)=0.0
    149152        ELSEWHERE
    150           pctsrf_t(:,is_oce,k)=1.0-zmasq-pctsrf_t(:,is_sic,k)
    151           WHERE(pctsrf_t(:,is_oce,k)<EPSFRA)
    152             pctsrf_t(:,is_oce,k)=0.0
    153             pctsrf_t(:,is_sic,k)=1.0-zmasq
    154           END WHERE
     153           pctsrf_t(:,is_oce,k)=1.0-zmasq-pctsrf_t(:,is_sic,k)
     154           WHERE(pctsrf_t(:,is_oce,k)<EPSFRA)
     155              pctsrf_t(:,is_oce,k)=0.0
     156              pctsrf_t(:,is_sic,k)=1.0-zmasq
     157           END WHERE
    155158        END WHERE
    156       END WHERE
    157       nbad=COUNT(pctsrf_t(:,is_oce,k)<0.0)
    158       IF(nbad>0) WRITE(lunout,*) 'pb sous maille pour nb point = ',nbad
    159       nbad=COUNT(abs(sum(pctsrf_t(:,:,k),dim=2)-1.0)>EPSFRA)
    160       IF(nbad>0) WRITE(lunout,*) 'pb sous surface pour nb points = ',nbad
    161     ELSE
    162       pctsrf_t(:,is_ter,k)=pctsrf(:,is_ter)
    163       WHERE(NINT(pctsrf(:,is_ter))==1)
    164         pctsrf_t(:,is_sic,k)=0.
    165         pctsrf_t(:,is_oce,k)=0.                 
    166         WHERE(fi_ice>=1.e-5)
    167           pctsrf_t(:,is_lic,k)=fi_ice
    168           pctsrf_t(:,is_ter,k)=pctsrf_t(:,is_ter,k)-pctsrf_t(:,is_lic,k)
    169         ELSEWHERE
    170           pctsrf_t(:,is_lic,k)=0.0
    171         END WHERE
    172       ELSEWHERE
    173         pctsrf_t(:,is_lic,k) = 0.0
    174         WHERE(fi_ice>=1.e-5)
    175           pctsrf_t(:,is_ter,k)=0.0
    176           pctsrf_t(:,is_sic,k)=fi_ice
    177           pctsrf_t(:,is_oce,k)=1.0-pctsrf_t(:,is_sic,k)
    178         ELSEWHERE
    179           pctsrf_t(:,is_sic,k)=0.0
    180           pctsrf_t(:,is_oce,k)=1.0
    181         END WHERE
    182       END WHERE
    183       verif=sum(pctsrf_t(:,:,k),dim=2)
    184       nbad=COUNT(verif<1.0-1.e-5.OR.verif>1.0+1.e-5)
    185       IF(nbad>0) WRITE(lunout,*) 'pb sous maille pour nb point = ',nbad
    186     END IF
     159     END WHERE
     160     nbad=COUNT(pctsrf_t(:,is_oce,k)<0.0)
     161     IF(nbad>0) WRITE(lunout,*) 'pb sous maille pour nb point = ',nbad
     162     nbad=COUNT(abs(sum(pctsrf_t(:,:,k),dim=2)-1.0)>EPSFRA)
     163     IF(nbad>0) WRITE(lunout,*) 'pb sous surface pour nb points = ',nbad
    187164  END DO
    188165  DEALLOCATE(phy_ice)
    189166
    190167!--- SST TREATMENT -------------------------------------------------------------
    191   WRITE(lunout,*) 'Traitement de la sst' ; sstfile=''; lCPL=.FALSE.
     168  IF (prt_level>1) WRITE(lunout,*) 'Traitement de la sst'
    192169
    193170! Input SST file selection
    194   sstfile='fake'
    195   IF(NF90_OPEN(famipsst,NF90_NOWRITE,nid)==NF90_NOERR) sstfile=TRIM(famipsst)
    196   IF(NF90_OPEN(fclimsst,NF90_NOWRITE,nid)==NF90_NOERR) sstfile=TRIM(fclimsst)
    197   IF(NF90_OPEN(fcpldsst,NF90_NOWRITE,nid)==NF90_NOERR) sstfile=TRIM(fcpldsst)
    198   SELECT CASE(icefile)
    199     CASE(famipsic); dumstr='Amip.'
    200     CASE(fclimsic); dumstr='Amip climatologique.'
    201     CASE(fcpldsic); dumstr='de sortie du modele couplé IPSL/IPCC.';lCPL=.TRUE.
    202     CASE('fake');   CALL abort_gcm('limit_netcdf','Fichier SST non reconnu',1)
    203   END SELECT
     171! Open file only to test if available
     172  IF ( NF90_OPEN(TRIM(famipsst),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
     173     sstfile=TRIM(famipsst)
     174     varname='tosbcs'
     175  ELSE IF ( NF90_OPEN(TRIM(fcpldsst),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
     176     sstfile=TRIM(fcpldsst)
     177     varname='SISUTESW'
     178  ELSE IF ( NF90_OPEN(TRIM(fhistsst),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
     179     sstfile=TRIM(fhistsst)
     180     varname='tsol_oce'
     181  ELSE
     182     WRITE(lunout,*) 'ERROR! No sst input file was found.'
     183     WRITE(lunout,*) 'One of following files must be availible : ',trim(famipsst),trim(fcpldsst),trim(fhistsst)
     184     CALL abort_gcm('limit_netcdf','No sst file was found',1)
     185  END IF
    204186  ierr=NF90_CLOSE(nid)
    205   WRITE(lunout,*)'Pour la temperature de mer a ete choisi un fichier '//TRIM(dumstr)
    206 
    207   CALL get_2Dfield(trim(sstfile),'SST',interbar,ndays,phy_sst,flag=extrap)
     187  IF (prt_level>=0) WRITE(lunout,*)'Pour la temperature de mer a ete choisi un fichier '//TRIM(sstfile)
     188
     189  CALL get_2Dfield(sstfile,varname,'SST',interbar,ndays,phy_sst,flag=extrap)
    208190
    209191!--- ALBEDO TREATMENT ----------------------------------------------------------
    210   WRITE(lunout,*) 'Traitement de l albedo'
    211   CALL get_2Dfield(falbe,'ALB',interbar,ndays,phy_alb)
     192  IF (prt_level>1) WRITE(lunout,*) 'Traitement de l albedo'
     193  varname='ALBEDO'
     194  CALL get_2Dfield(falbe,varname,'ALB',interbar,ndays,phy_alb)
    212195
    213196!--- REFERENCE GROUND HEAT FLUX TREATMENT --------------------------------------
    … …  
    215198
    216199!--- OUTPUT FILE WRITING -------------------------------------------------------
    217   WRITE(lunout,*) 'Ecriture du fichier limit'
     200  IF (prt_level>5) WRITE(lunout,*) 'Ecriture du fichier limit : debut'
    218201
    219202  !--- File creation
    … …  
    264247  ierr=NF90_CLOSE(nid)
    265248
     249  IF (prt_level>5) WRITE(lunout,*) 'Ecriture du fichier limit : fin'
     250
    266251  DEALLOCATE(pctsrf_t,phy_sst,phy_bil,phy_alb,phy_rug)
    267 #endif
    268 ! of #ifdef CPP_EARTH
    269252
    270253
    … …  
    278261!-------------------------------------------------------------------------------
    279262!
    280 SUBROUTINE get_2Dfield(fnam, mode, ibar, ndays, champo, flag, mask, lCPL)
     263SUBROUTINE get_2Dfield(fnam, varname, mode, ibar, ndays, champo, flag, mask)
    281264!
    282265!-----------------------------------------------------------------------------
    … …  
    306289! Arguments:
    307290  CHARACTER(LEN=*),  INTENT(IN)     :: fnam     ! NetCDF file name
     291  CHARACTER(LEN=10), INTENT(IN)     :: varname  ! NetCDF variable name
    308292  CHARACTER(LEN=3),  INTENT(IN)     :: mode     ! RUG, SIC, SST or ALB
    309293  LOGICAL,           INTENT(IN)     :: ibar     ! interp on pressure levels
    310294  INTEGER,           INTENT(IN)     :: ndays    ! current year number of days
    311   REAL,    POINTER,  DIMENSION(:, :) :: champo   ! output field = f(t)
     295  REAL,    POINTER,  DIMENSION(:, :) :: champo  ! output field = f(t)
    312296  LOGICAL, OPTIONAL, INTENT(IN)     :: flag     ! extrapol. (SST) old ice (SIC)
    313297  REAL,    OPTIONAL, DIMENSION(iim, jjp1), INTENT(IN) :: mask
    314   LOGICAL, OPTIONAL, INTENT(IN)     :: lCPL     ! Coupled model flag (for ICE)
    315298!------------------------------------------------------------------------------
    316299! Local variables:
    … …  
    318301  INTEGER :: ncid, varid                  ! NetCDF identifiers
    319302  CHARACTER(LEN=30)               :: dnam       ! dimension name
    320   CHARACTER(LEN=80)               :: varname    ! NetCDF variable name
    321303!--- dimensions
    322304  INTEGER,           DIMENSION(4) :: dids       ! NetCDF dimensions identifiers
    … …  
    333315!--- input files
    334316  CHARACTER(LEN=20)                 :: cal_in   ! calendar
     317  CHARACTER(LEN=20)                 :: unit_sic ! attribute unit in sea-ice file
    335318  INTEGER                           :: ndays_in ! number of days
    336319!--- misc
    … …  
    339322  CHARACTER(LEN=25)                 :: title    ! for messages
    340323  LOGICAL                           :: extrp    ! flag for extrapolation
     324  LOGICAL                           :: oldice   ! flag for old way ice computation
    341325  REAL                              :: chmin, chmax
    342326  INTEGER ierr
    343327  integer n_extrap ! number of extrapolated points
    344328  logical skip
     329
    345330!------------------------------------------------------------------------------
    346331!---Variables depending on keyword 'mode' -------------------------------------
    347332  NULLIFY(champo)
     333
    348334  SELECT CASE(mode)
    349     CASE('RUG'); varname='RUGOS'; title='Rugosite'
    350     CASE('SIC'); varname='sicbcs'; title='Sea-ice'
    351     CASE('SST'); varname='tosbcs'; title='SST'
    352     CASE('ALB'); varname='ALBEDO'; title='Albedo'
     335  CASE('RUG'); title='Rugosite'
     336  CASE('SIC'); title='Sea-ice'
     337  CASE('SST'); title='SST'
     338  CASE('ALB'); title='Albedo'
    353339  END SELECT
     340 
     341
    354342  extrp=.FALSE.
     343  oldice=.FALSE.
    355344  IF ( PRESENT(flag) ) THEN
    356345    IF ( flag .AND. mode=='SST' ) extrp=.TRUE.
     346    IF ( flag .AND. mode=='SIC' ) oldice=.TRUE.
    357347  END IF
    358348
    359349!--- GETTING SOME DIMENSIONAL VARIABLES FROM FILE -----------------------------
     350  IF (prt_level>5) WRITE(lunout,*) ' Now reading file : ',fnam
    360351  ierr=NF90_OPEN(fnam, NF90_NOWRITE, ncid);             CALL ncerr(ierr, fnam)
    361   ierr=NF90_INQ_VARID(ncid, varname, varid);            CALL ncerr(ierr, fnam)
     352  ierr=NF90_INQ_VARID(ncid, trim(varname), varid);            CALL ncerr(ierr, fnam)
    362353  ierr=NF90_INQUIRE_VARIABLE(ncid, varid, dimids=dids); CALL ncerr(ierr, fnam)
     354
     355!--- Read unit for sea-ice variable only
     356  IF (mode=='SIC') THEN
     357     ierr=NF90_GET_ATT(ncid, varid, 'units', unit_sic)
     358     IF(ierr/=NF90_NOERR) THEN
     359        IF (prt_level>5) WRITE(lunout,*) 'No unit was given in sea-ice file. Take percentage as default value'
     360        unit_sic='X'
     361     ELSE
     362        IF (prt_level>5) WRITE(lunout,*) ' Sea-ice cover has unit=',unit_sic
     363     END IF
     364  END IF
    363365
    364366!--- Longitude
    … …  
    367369  ierr=NF90_INQ_VARID(ncid, dnam, varid);                CALL ncerr(ierr, fnam)
    368370  ierr=NF90_GET_VAR(ncid, varid, dlon_ini);              CALL ncerr(ierr, fnam)
    369   WRITE(lunout, *) 'variable ', dnam, 'dimension ', imdep
     371  IF (prt_level>5) WRITE(lunout, *) 'variable ', dnam, 'dimension ', imdep
    370372
    371373!--- Latitude
    … …  
    374376  ierr=NF90_INQ_VARID(ncid, dnam, varid);                CALL ncerr(ierr, fnam)
    375377  ierr=NF90_GET_VAR(ncid, varid, dlat_ini);              CALL ncerr(ierr, fnam)
    376   WRITE(lunout, *) 'variable ', dnam, 'dimension ', jmdep
     378  IF (prt_level>5) WRITE(lunout, *) 'variable ', dnam, 'dimension ', jmdep
    377379
    378380!--- Time (variable is not needed - it is rebuilt - but calendar is)
    … …  
    387389      CASE('SIC', 'SST'); cal_in='gregorian'
    388390    END SELECT
    389     WRITE(lunout, *)'ATTENTION: variable "time" sans attribut "calendrier" ' &
     391    IF (prt_level>5) WRITE(lunout, *)'ATTENTION: variable "time" sans attribut "calendrier" ' &
    390392         // 'dans '//TRIM(fnam)//'. On choisit la valeur par defaut.'
    391393  END IF
    392   WRITE(lunout, *) 'variable ', dnam, 'dimension ', lmdep, 'calendrier ', &
     394  IF (prt_level>5) WRITE(lunout, *) 'variable ', dnam, 'dimension ', lmdep, 'calendrier ', &
    393395       cal_in
    394396
     397 
    395398!--- CONSTRUCTING THE INPUT TIME VECTOR FOR INTERPOLATION --------------------
    396399  !--- Determining input file number of days, depending on calendar
    … …  
    400403!--- If input records are not monthly, time sampling has to be constant !
    401404  timeyear=mid_months(anneeref, cal_in, lmdep)
    402   IF (lmdep /= 12) WRITE(lunout, '(a, i3, a)') 'Note : les fichiers de ' &
    403        // TRIM(mode) // ' ne comportent pas 12, mais ', lmdep, &
    404        ' enregistrements.'
     405  IF (lmdep /= 12) WRITE(lunout,*) 'Note : les fichiers de ', TRIM(mode), &
     406       ' ne comportent pas 12, mais ', lmdep, ' enregistrements.'
    405407
    406408!--- GETTING THE FIELD AND INTERPOLATING IT ----------------------------------
    … …  
    408410  IF(extrp) ALLOCATE(work(imdep, jmdep))
    409411
    410   WRITE(lunout, *)
    411   WRITE(lunout, '(a, i3, a)')'LECTURE ET INTERPOLATION HORIZ. DE ', lmdep, &
    412        ' CHAMPS.'
     412  IF (prt_level>5) WRITE(lunout, *)
     413  IF (prt_level>5) WRITE(lunout,*)'LECTURE ET INTERPOLATION HORIZ. DE ', lmdep, ' CHAMPS.'
    413414  ierr=NF90_INQ_VARID(ncid, varname, varid);             CALL ncerr(ierr, fnam)
    414415  DO l=1, lmdep
    … …  
    421422         work)
    422423
    423     IF(ibar.AND..NOT.(mode=='SIC'.AND.flag)) THEN
    424       IF(l==1) THEN
    425         WRITE(lunout, *)                                                      &
    426   '-------------------------------------------------------------------------'
    427         WRITE(lunout, *)                                                     &
    428   'Utilisation de l''interpolation barycentrique pour '//TRIM(title)//' $$$'
    429         WRITE(lunout, *)                                                      &
    430   '-------------------------------------------------------------------------'
     424    IF(ibar .AND. .NOT.oldice) THEN
     425      IF(l==1 .AND. prt_level>5) THEN
     426        WRITE(lunout, *) '-------------------------------------------------------------------------'
     427        WRITE(lunout, *) 'Utilisation de l''interpolation barycentrique pour ',TRIM(title),' $$$'
     428        WRITE(lunout, *) '-------------------------------------------------------------------------'
    431429      END IF
    432430      IF(mode=='RUG') champ=LOG(champ)
    … …  
    455453
    456454!--- TIME INTERPOLATION ------------------------------------------------------
    457   WRITE(lunout, *)
    458   WRITE(lunout, *)'INTERPOLATION TEMPORELLE.'
    459   WRITE(lunout, "(2x, ' Vecteur temps en entree: ', 10f6.1)") timeyear
    460   WRITE(lunout, "(2x, ' Vecteur temps en sortie de 0 a ', i3)") ndays
     455  IF (prt_level>5) THEN
     456     WRITE(lunout, *)
     457     WRITE(lunout, *)'INTERPOLATION TEMPORELLE.'
     458     WRITE(lunout, *)' Vecteur temps en entree: ', timeyear
     459     WRITE(lunout, *)' Vecteur temps en sortie de 0 a ', ndays
     460  END IF
     461
    461462  ALLOCATE(yder(lmdep), champan(iip1, jjp1, ndays))
    462463  skip = .false.
    … …  
    473474  END DO
    474475  if (n_extrap /= 0) then
    475      print *, "get_2Dfield pchfe_95: n_extrap = ", n_extrap
     476     WRITE(lunout,*) "get_2Dfield pchfe_95: n_extrap = ", n_extrap
    476477  end if
    477478  champan(iip1, :, :)=champan(1, :, :)
    … …  
    481482  DO j=1, jjp1
    482483    CALL minmax(iip1, champan(1, j, 10), chmin, chmax)
    483     WRITE(lunout, *)' '//TRIM(title)//' au temps 10 ', chmin, chmax, j
     484    IF (prt_level>5) WRITE(lunout, *)' ',TRIM(title),' au temps 10 ', chmin, chmax, j
    484485  END DO
    485486
    486487!--- SPECIAL FILTER FOR SST: SST>271.38 --------------------------------------
    487488  IF(mode=='SST') THEN
    488     WRITE(lunout, *) 'Filtrage de la SST: SST >= 271.38'
     489    IF (prt_level>5) WRITE(lunout, *) 'Filtrage de la SST: SST >= 271.38'
    489490    WHERE(champan<271.38) champan=271.38
    490491  END IF
    … …  
    492493!--- SPECIAL FILTER FOR SIC: 0.0<SIC<1.0 -------------------------------------
    493494  IF(mode=='SIC') THEN
    494     WRITE(lunout, *) 'Filtrage de la SIC: 0.0 < Sea-ice < 1.0'
    495     IF(.NOT.lCPL) champan(:, :, :)=champan(:, :, :)/100.
     495    IF (prt_level>5) WRITE(lunout, *) 'Filtrage de la SIC: 0.0 < Sea-ice < 1.0'
     496
     497    IF (unit_sic=='1') THEN
     498       ! Nothing to be done for sea-ice field is already in fraction of 1
     499       ! This is the case for sea-ice in file cpl_atm_sic.nc
     500       IF (prt_level>5) WRITE(lunout,*) 'Sea-ice field already in fraction of 1'
     501    ELSE
     502       ! Convert sea ice from percentage to fraction of 1
     503       IF (prt_level>5) WRITE(lunout,*) 'Transformt sea-ice field from percentage to fraction of 1.'
     504       champan(:, :, :)=champan(:, :, :)/100.
     505    END IF
     506
    496507    champan(iip1, :, :)=champan(1, :, :)
    497508    WHERE(champan>1.0) champan=1.0
    498509    WHERE(champan<0.0) champan=0.0
    499   END IF
     510 END IF
    500511
    501512!--- DYNAMICAL TO PHYSICAL GRID ----------------------------------------------
    … …  
    592603
    593604!--- Mid-months times
    594   mid_months(1)=0.5*FLOAT(mnth(1))
     605  mid_months(1)=0.5*REAL(mnth(1))
    595606  DO k=2,nm
    596     mid_months(k)=mid_months(k-1)+0.5*FLOAT(mnth(k-1)+mnth(k))
     607    mid_months(k)=mid_months(k-1)+0.5*REAL(mnth(k-1)+mnth(k))
    597608  END DO
    598609
    … …  
    626637!-------------------------------------------------------------------------------
    627638
     639#endif
     640! of #ifdef CPP_EARTH
    628641
    629642END SUBROUTINE limit_netcdf

  • LMDZ5/trunk/libf/dyn3dmem/limy.F

    r1632 r1673  
    1 !
    2 ! $Header$
    3 !
     1c
     2c $Id$
     3c
    44      SUBROUTINE limy(s0,sy,sm,pente_max)
    55c
    … …  
    4040      REAL qbyv(ip1jm,llm)
    4141
    42       REAL qpns,qpsn,apn,aps,dyn1,dys1,dyn2,dys2
     42      REAL qpns,qpsn,appn,apps,dyn1,dys1,dyn2,dys2
    4343      Logical extremum,first
    4444      save first
    … …  
    5252      REAL      SSUM
    5353      integer ismax,ismin
    54       EXTERNAL  SSUM, ismin,ismax
     54      EXTERNAL  SSUM, convflu,ismin,ismax
     55      EXTERNAL filtreg
    5556
    5657      data first/.true./
    … …  
    116117
    117118c     print*,dyqv(iip1+1)
    118 c     apn=abs(dyq(1)/dyqv(iip1+1))
     119c     appn=abs(dyq(1)/dyqv(iip1+1))
    119120c     print*,dyq(ip1jm+1)
    120121c     print*,dyqv(ip1jm-iip1+1)
    121 c     aps=abs(dyq(ip1jm+1)/dyqv(ip1jm-iip1+1))
     122c     apps=abs(dyq(ip1jm+1)/dyqv(ip1jm-iip1+1))
    122123c     do ij=2,iim
    123 c        apn=amax1(abs(dyq(ij)/dyqv(ij)),apn)
    124 c        aps=amax1(abs(dyq(ip1jm+ij)/dyqv(ip1jm-iip1+ij)),aps)
     124c        appn=amax1(abs(dyq(ij)/dyqv(ij)),appn)
     125c        apps=amax1(abs(dyq(ip1jm+ij)/dyqv(ip1jm-iip1+ij)),apps)
    125126c     enddo
    126 c     apn=min(pente_max/apn,1.)
    127 c     aps=min(pente_max/aps,1.)
     127c     appn=min(pente_max/appn,1.)
     128c     apps=min(pente_max/apps,1.)
    128129
    129130
    … …  
    131132
    132133c     if(dyqv(ismin(iim,dyqv,1))*dyqv(ismax(iim,dyqv,1)).le.0.)
    133 c    &   apn=0.
     134c    &   appn=0.
    134135c     if(dyqv(ismax(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1)*
    135136c    &   dyqv(ismin(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1).le.0.)
    136 c    &   aps=0.
     137c    &   apps=0.
    137138
    138139c   limitation des pentes aux poles
    139140c     do ij=1,iip1
    140 c        dyq(ij)=apn*dyq(ij)
    141 c        dyq(ip1jm+ij)=aps*dyq(ip1jm+ij)
     141c        dyq(ij)=appn*dyq(ij)
     142c        dyq(ip1jm+ij)=apps*dyq(ip1jm+ij)
    142143c     enddo
    143144

  • LMDZ5/trunk/libf/dyn3dmem/logic.h

    r1658 r1673  
    11!
    2 ! $Id: logic.h 1319 2010-02-23 21:29:54Z fairhead $
     2! $Id$
    33!
    44!
    5 !
     5! NB: keep items of different kinds in seperate common blocs to avoid
     6!     "misaligned commons" issues
    67!-----------------------------------------------------------------------
    78! INCLUDE 'logic.h'
    89
    9       COMMON/logic/ purmats,iflag_phys,forward,leapf,apphys,            &
     10      COMMON/logicl/ purmats,forward,leapf,apphys,                      &
    1011     &  statcl,conser,apdiss,apdelq,saison,ecripar,fxyhypb,ysinus       &
    1112     &  ,read_start,ok_guide,ok_strato,ok_gradsfile                     &
    12      &  ,ok_limit,ok_etat0
     13     &  ,ok_limit,ok_etat0,grilles_gcm_netcdf,hybrid
    1314
     15      COMMON/logici/ iflag_phys,iflag_trac
     16     
    1417      LOGICAL purmats,forward,leapf,apphys,statcl,conser,               &
    1518     & apdiss,apdelq,saison,ecripar,fxyhypb,ysinus                      &
    1619     &  ,read_start,ok_guide,ok_strato,ok_gradsfile                     &
    17      &  ,ok_limit,ok_etat0
     20     &  ,ok_limit,ok_etat0,grilles_gcm_netcdf
     21      logical hybrid ! vertical coordinate is hybrid if true (sigma otherwise)
     22                     ! (only used if disvert_type==2)
    1823
    19       INTEGER iflag_phys
    20 !$OMP THREADPRIVATE(/logic/)
     24      integer iflag_phys,iflag_trac
     25!$OMP THREADPRIVATE(/logicl/)
     26!$OMP THREADPRIVATE(/logici/)
    2127!-----------------------------------------------------------------------

  • LMDZ5/trunk/libf/dyn3dmem/mod_interface_dyn_phys.F90

    r1632 r1673  
    11!
    2 ! $Id: mod_interface_dyn_phys.F90 1279 2009-12-10 09:02:56Z fairhead $
     2! $Id$
    33!
    44MODULE mod_interface_dyn_phys
    … …  
    77 
    88 
    9 #ifdef CPP_EARTH
     9#ifdef CPP_PHYS
    1010! Interface with parallel physics,
    11 ! for now this routine only works with Earth physics
    1211CONTAINS
    1312 
    … …  
    5655  END SUBROUTINE Init_interface_dyn_phys
    5756#endif
    58 ! of #ifdef CPP_EARTH
     57! of #ifdef CPP_PHYS
    5958END MODULE mod_interface_dyn_phys

  • LMDZ5/trunk/libf/dyn3dmem/parallel.F90

    r1632 r1673  
    11!
    2 ! $Id: parallel.F90 1279 2009-12-10 09:02:56Z fairhead $
     2! $Id$
    33!
    44  module parallel
    … …  
    9494      integer, dimension(3) :: blocklen,type
    9595      integer :: comp_id
    96 
     96      character(len=4)  :: num
     97      character(len=20) :: filename
     98 
    9799#ifdef CPP_OMP   
    98100      INTEGER :: OMP_GET_NUM_THREADS
    … …  
    126128        mpi_rank=0
    127129      ENDIF
    128  
     130
     131
     132! Open text output file with mpi_rank in suffix of file name
     133      IF (lunout /= 5 .and. lunout /= 6) THEN
     134         WRITE(num,'(I4.4)') mpi_rank
     135         filename='lmdz.out_'//num
     136         IF (mpi_rank .NE. 0) THEN
     137            OPEN(UNIT=lunout,FILE=TRIM(filename),ACTION='write', &
     138               STATUS='unknown',FORM='formatted',IOSTAT=ierr)
     139         ENDIF
     140      ENDIF
     141
    129142     
    130143      allocate(jj_begin_para(0:mpi_size-1))
    … …  
    376389      integer :: ierr
    377390      integer :: i
    378       deallocate(jj_begin_para)
    379       deallocate(jj_end_para)
    380       deallocate(jj_nb_para)
     391
     392      if (allocated(jj_begin_para)) deallocate(jj_begin_para)
     393      if (allocated(jj_end_para))   deallocate(jj_end_para)
     394      if (allocated(jj_nb_para))    deallocate(jj_nb_para)
    381395
    382396      if (type_ocean == 'couple') then
    … …  
    717731       
    718732   
    719     /* 
    720   Subroutine verif_hallo(Field,ij,ll,up,down)
    721     implicit none
    722 #include "dimensions.h"
    723 #include "paramet.h"   
    724     include 'mpif.h'
    725    
    726       INTEGER :: ij,ll
    727       REAL, dimension(ij,ll) :: Field
    728       INTEGER :: up,down
    729      
    730       REAL,dimension(ij,ll): NewField
    731      
    732       NewField=0
    733      
    734       ijb=ij_begin
    735       ije=ij_end
    736       if (pole_nord)
    737       NewField(ij_be       
    738 */
     733!  Subroutine verif_hallo(Field,ij,ll,up,down)
     734!    implicit none
     735!#include "dimensions.h"
     736!#include "paramet.h"   
     737!    include 'mpif.h'
     738!   
     739!      INTEGER :: ij,ll
     740!      REAL, dimension(ij,ll) :: Field
     741!      INTEGER :: up,down
     742!     
     743!      REAL,dimension(ij,ll): NewField
     744!     
     745!      NewField=0
     746!     
     747!      ijb=ij_begin
     748!      ije=ij_end
     749!      if (pole_nord)
     750!      NewField(ij_be       
     751
    739752  end module parallel

  • LMDZ5/trunk/libf/dyn3dmem/paramet.h

    r1632 r1673  
    11!
    2 ! $Header$
     2! $Id$
    33!
    44!

  • LMDZ5/trunk/libf/dyn3dmem/temps.h

    r1632 r1673  
    11!
    2 ! $Id: temps.h 1279 2009-12-10 09:02:56Z fairhead $
     2! $Id$
    33!
    44!  ATTENTION!!!!: ce fichier include est compatible format fixe/format libre
    … …  
    1414
    1515      COMMON/temps/itaufin, dt, day_ini, day_end, annee_ref, day_ref,   &
    16      &             itau_dyn, itau_phy, jD_ref, jH_ref, calend
     16     &             itau_dyn, itau_phy, jD_ref, jH_ref, calend,          &
     17     &             start_time
     18
    1719
    1820      INTEGER   itaufin
    1921      INTEGER itau_dyn, itau_phy
    2022      INTEGER day_ini, day_end, annee_ref, day_ref
    21       REAL      dt, jD_ref, jH_ref
     23      REAL      dt, jD_ref, jH_ref, start_time
    2224      CHARACTER (len=10) :: calend
    2325

  • LMDZ5/trunk/libf/dyn3dmem/vlsplt_loc.F

    r1632 r1673  
     1!
     2! $Id$
     3!
    14      SUBROUTINE vlx_loc(q,pente_max,masse,u_m,ijb_x,ije_x)
    25
    … …  
    372375      REAL qbyv(ijb_v:ije_v,llm)
    373376
    374       REAL qpns,qpsn,apn,aps,dyn1,dys1,dyn2,dys2,newmasse,fn,fs
     377      REAL qpns,qpsn,appn,apps,dyn1,dys1,dyn2,dys2,newmasse,fn,fs
    375378c     REAL newq,oldmasse
    376379      Logical extremum,first,testcpu
    … …  
    543546C     PRINT*,dyq(1)
    544547C     PRINT*,dyqv(iip1+1)
    545 C     apn=abs(dyq(1)/dyqv(iip1+1))
     548C     appn=abs(dyq(1)/dyqv(iip1+1))
    546549C     PRINT*,dyq(ip1jm+1)
    547550C     PRINT*,dyqv(ip1jm-iip1+1)
    548 C     aps=abs(dyq(ip1jm+1)/dyqv(ip1jm-iip1+1))
     551C     apps=abs(dyq(ip1jm+1)/dyqv(ip1jm-iip1+1))
    549552C     DO ij=2,iim
    550 C        apn=amax1(abs(dyq(ij)/dyqv(ij)),apn)
    551 C        aps=amax1(abs(dyq(ip1jm+ij)/dyqv(ip1jm-iip1+ij)),aps)
     553C        appn=amax1(abs(dyq(ij)/dyqv(ij)),appn)
     554C        apps=amax1(abs(dyq(ip1jm+ij)/dyqv(ip1jm-iip1+ij)),apps)
    552555C     ENDDO
    553 C     apn=min(pente_max/apn,1.)
    554 C     aps=min(pente_max/aps,1.)
     556C     appn=min(pente_max/appn,1.)
     557C     apps=min(pente_max/apps,1.)
    555558C
    556559C
    … …  
    558561C
    559562C     IF(dyqv(ismin(iim,dyqv,1))*dyqv(ismax(iim,dyqv,1)).le.0.)
    560 C    &   apn=0.
     563C    &   appn=0.
    561564C     IF(dyqv(ismax(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1)*
    562565C    &   dyqv(ismin(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1).le.0.)
    563 C    &   aps=0.
     566C    &   apps=0.
    564567C
    565568C   limitation des pentes aux poles
    566569C     DO ij=1,iip1
    567 C        dyq(ij)=apn*dyq(ij)
    568 C        dyq(ip1jm+ij)=aps*dyq(ip1jm+ij)
     570C        dyq(ij)=appn*dyq(ij)
     571C        dyq(ip1jm+ij)=apps*dyq(ip1jm+ij)
    569572C     ENDDO
    570573C

  • LMDZ5/trunk/libf/dyn3dmem/vlspltqs_loc.F

    r1632 r1673  
    549549C     PRINT*,dyq(1)
    550550C     PRINT*,dyqv(iip1+1)
    551 C     apn=abs(dyq(1)/dyqv(iip1+1))
     551C     appn=abs(dyq(1)/dyqv(iip1+1))
    552552C     PRINT*,dyq(ip1jm+1)
    553553C     PRINT*,dyqv(ip1jm-iip1+1)
    554 C     aps=abs(dyq(ip1jm+1)/dyqv(ip1jm-iip1+1))
     554C     apps=abs(dyq(ip1jm+1)/dyqv(ip1jm-iip1+1))
    555555C     DO ij=2,iim
    556 C        apn=amax1(abs(dyq(ij)/dyqv(ij)),apn)
    557 C        aps=amax1(abs(dyq(ip1jm+ij)/dyqv(ip1jm-iip1+ij)),aps)
     556C        appn=amax1(abs(dyq(ij)/dyqv(ij)),appn)
     557C        apps=amax1(abs(dyq(ip1jm+ij)/dyqv(ip1jm-iip1+ij)),apps)
    558558C     ENDDO
    559 C     apn=min(pente_max/apn,1.)
    560 C     aps=min(pente_max/aps,1.)
     559C     appn=min(pente_max/appn,1.)
     560C     apps=min(pente_max/apps,1.)
    561561C
    562562C
    … …  
    564564C
    565565C     IF(dyqv(ismin(iim,dyqv,1))*dyqv(ismax(iim,dyqv,1)).le.0.)
    566 C    &   apn=0.
     566C    &   appn=0.
    567567C     IF(dyqv(ismax(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1)*
    568568C    &   dyqv(ismin(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1).le.0.)
    569 C    &   aps=0.
     569C    &   apps=0.
    570570C
    571571C   limitation des pentes aux poles
    572572C     DO ij=1,iip1
    573 C        dyq(ij)=apn*dyq(ij)
    574 C        dyq(ip1jm+ij)=aps*dyq(ip1jm+ij)
     573C        dyq(ij)=appn*dyq(ij)
     574C        dyq(ip1jm+ij)=apps*dyq(ip1jm+ij)
    575575C     ENDDO
    576576C

  • LMDZ5/trunk/libf/dyn3dmem/wrgrads.F

    r1632 r1673  
    2222c   local
    2323
    24       integer im,jm,lm,i,j,l,lnblnk,iv,iii,iji,iif,ijf
     24      integer im,jm,lm,i,j,l,iv,iii,iji,iif,ijf
    2525
    2626      logical writectl
    … …  
    5555            nvar(if)=ivar(if)
    5656            var(ivar(if),if)=name
    57             tvar(ivar(if),if)=titlevar(1:lnblnk(titlevar))
     57            tvar(ivar(if),if)=trim(titlevar)
    5858            nld(ivar(if),if)=nl
    5959            print*,'initialisation ecriture de ',var(ivar(if),if)
    … …  
    9696      file=fichier(if)
    9797c   WARNING! on reecrase le fichier .ctl a chaque ecriture
    98       open(unit(if),file=file(1:lnblnk(file))//'.ctl'
     98      open(unit(if),file=trim(file)//'.ctl'
    9999     &         ,form='formatted',status='unknown')
    100100      write(unit(if),'(a5,1x,a40)')
    101      &       'DSET ','^'//file(1:lnblnk(file))//'.dat'
     101     &       'DSET ','^'//trim(file)//'.dat'
    102102
    103103      write(unit(if),'(a12)') 'UNDEF 1.0E30'
Note: See TracChangeset for help on using the changeset viewer.