Ignore:
Timestamp:
Jan 11, 2013, 10:19:19 AM (12 years ago)
Author:
Laurent Fairhead
Message:

Version testing basée sur la r1706


Testing release based on r1706

Location:
LMDZ5/branches/testing
Files:
26 deleted
95 edited
30 copied

Legend:

Unmodified
Added
Removed
  • LMDZ5/branches/testing

  • LMDZ5/branches/testing/libf/dyn3d/calfis.F

    r1669 r1707  
    507507     .             debut_split,    !! firstcall
    508508     .             lafin_split,    !! lastcall
    509      .             float(day_ini), !! pday <-- day_ini (dans temps.h)
     509     .             jD_cur,         !! pday. see leapfrog
    510510     .             jH_cur_split,   !! ptime "fraction of day"
    511511     .             zdt_split,      !! ptimestep
  • LMDZ5/branches/testing/libf/dyn3d/comconst.h

    r1505 r1707  
    2121      REAL dtdiss ! (s) time step for the dissipation
    2222      REAL rad ! (m) radius of the planet
    23       REAL r ! Gas constant R=8.31 J.K-1.mol-1
    24       REAL cpp   ! Cp
     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)
    2526      REAL kappa ! kappa=R/Cp
    2627      REAL cotot
  • LMDZ5/branches/testing/libf/dyn3d/comdissnew.h

    r1319 r1707  
    1212
    1313      COMMON/comdissnew/ lstardis,nitergdiv,nitergrot,niterh,tetagdiv,  &
    14      &                   tetagrot,tetatemp,coefdis 
     14     &                   tetagrot,tetatemp,coefdis, vert_prof_dissip
    1515
    1616      LOGICAL lstardis
    1717      INTEGER nitergdiv, nitergrot, niterh
     18
     19      integer vert_prof_dissip ! vertical profile of horizontal dissipation
     20!     Allowed values:
     21!     0: rational fraction, function of pressure
     22!     1: tanh of altitude
     23
    1824      REAL     tetagdiv, tetagrot,  tetatemp, coefdis
    1925
  • LMDZ5/branches/testing/libf/dyn3d/conf_gcm.F

    r1665 r1707  
    1414#endif
    1515      USE infotrac, ONLY : type_trac
     16      use assert_m, only: assert
     17
    1618      IMPLICIT NONE
    1719c-----------------------------------------------------------------------
     
    9395      CALL getin('lunout', lunout)
    9496      IF (lunout /= 5 .and. lunout /= 6) THEN
    95         OPEN(lunout,FILE='lmdz.out')
     97        OPEN(UNIT=lunout,FILE='lmdz.out',ACTION='write',
     98     &          STATUS='unknown',FORM='formatted')
    9699      ENDIF
    97100
     
    173176
    174177!Config  Key  = nsplit_phys
    175 !Config  Desc = nombre de pas par jour
    176 !Config  Def  = 1
    177 !Config  Help = nombre de pas par jour (multiple de iperiod) (
    178 !Config          ici pour  dt = 1 min )
    179178       nsplit_phys = 1
    180179       CALL getin('nsplit_phys',nsplit_phys)
     
    625624      CALL getin('ok_dyn_ave',ok_dyn_ave)
    626625
    627 
    628626      write(lunout,*)' #########################################'
    629627      write(lunout,*)' Configuration des parametres du gcm: '
     
    635633      write(lunout,*)' day_step = ', day_step
    636634      write(lunout,*)' iperiod = ', iperiod
     635      write(lunout,*)' nsplit_phys = ', nsplit_phys
    637636      write(lunout,*)' iconser = ', iconser
    638637      write(lunout,*)' iecri = ', iecri
     
    805804!Config  Desc = sortie des transports zonaux dans la dynamique
    806805!Config  Def  = n
    807 !Config  Help =
     806!Config  Help = Permet de mettre en route le calcul des transports
    808807!Config         
    809        ok_dynzon = .FALSE.
    810        CALL getin('ok_dynzon',ok_dynzon)
     808      ok_dynzon = .FALSE.
     809      CALL getin('ok_dynzon',ok_dynzon)
    811810
    812811!Config  Key  = ok_dyn_ins
     
    838837        write(lunout,*)'STOP !!!'
    839838        write(lunout,*)'use_filtre_fft n est pas implemente dans dyn3d'
    840         STOP
     839        STOP 1
    841840      ENDIF
    842841     
     
    848847      ok_strato=.FALSE.
    849848      CALL getin('ok_strato',ok_strato)
     849
     850      vert_prof_dissip = merge(1, 0, ok_strato .and. llm==39)
     851      CALL getin('vert_prof_dissip', vert_prof_dissip)
     852      call assert(vert_prof_dissip == 0 .or. vert_prof_dissip ==  1,
     853     $     "bad value for vert_prof_dissip")
    850854
    851855!Config  Key  = ok_gradsfile
  • LMDZ5/branches/testing/libf/dyn3d/fxhyp.F

    r1403 r1707  
    6868       xzoom    = xzoomdeg * pi/180.
    6969c
     70       if (iim==1) then
     71
     72          rlonm025(1)=-pi/2.
     73          rlonv(1)=0.
     74          rlonu(1)=pi
     75          rlonp025(1)=pi/2.
     76          rlonm025(2)=rlonm025(1)+depi
     77          rlonv(2)=rlonv(1)+depi
     78          rlonu(2)=rlonu(1)+depi
     79          rlonp025(2)=rlonp025(1)+depi
     80
     81          xprimm025(:)=1.
     82          xprimv(:)=1.
     83          xprimu(:)=1.
     84          xprimp025(:)=1.
     85          champmin=depi
     86          champmax=depi
     87          return
     88
     89       endif
     90
    7091           decalx   = .75
    7192       IF( grossism.EQ.1..AND.scal180 )  THEN
     
    286307
    287308
     309
    288310       IF(ik.EQ.1.and.grossism.EQ.1.)  THEN
    289311         xvrai(1)    = xvrai(iip1)-depi
    290312         xxprim(1)   = xxprim(iip1)
    291313       ENDIF
     314
    292315       DO i = 1 , iim
    293316        xlon(i)     = xvrai(i)
  • LMDZ5/branches/testing/libf/dyn3d/gcm.F

    r1665 r1707  
    405405
    406406      CALL inidissip( lstardis, nitergdiv, nitergrot, niterh   ,
    407      *                tetagdiv, tetagrot , tetatemp              )
     407     *                tetagdiv, tetagrot , tetatemp, vert_prof_dissip)
    408408
    409409c-----------------------------------------------------------------------
     
    433433! Physics:
    434434#ifdef CPP_PHYS
    435          CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys ,
    436      ,                latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp     )
     435         CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys,
     436     &                latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp,
     437     &                iflag_phys)
    437438#endif
    438439         call_iniphys=.false.
     
    457458 301  FORMAT('1'/,15x,'run du ', i2,'/',i2,'/',i4)
    458459 302  FORMAT('1'/,15x,'    au ', i2,'/',i2,'/',i4)
    459 #endif
    460 
    461 #ifdef CPP_PHYS
    462 ! Create start file (startphy.nc) and boundary conditions (limit.nc)
    463 ! for the Earth verstion
    464        if (iflag_phys>=100) then
    465           call iniaqua(ngridmx,latfi,lonfi,iflag_phys)
    466        endif
    467460#endif
    468461
  • LMDZ5/branches/testing/libf/dyn3d/groupe.F

    r524 r1707  
    3838      integer i,j,l
    3939
    40       logical firstcall
    41       save firstcall
     40      logical firstcall,groupe_ok
     41      save firstcall,groupe_ok
    4242
    4343      data firstcall/.true./
     44      data groupe_ok/.true./
     45
     46      if (iim==1) then
     47         groupe_ok=.false.
     48      endif
    4449
    4550      if (firstcall) then
    46          if(mod(iim,2**ngroup).ne.0) stop'probleme du nombre ede point'
     51         if (groupe_ok) then
     52           if(mod(iim,2**ngroup).ne.0) stop'probleme du nombre de point'
     53         endif
    4754         firstcall=.false.
    4855      endif
     56
    4957
    5058c   Champs 1D
     
    5260      call convflu(pbaru,pbarv,llm,zconvm)
    5361
    54 c
    5562      call scopy(ijp1llm,zconvm,1,zconvmm,1)
    5663      call scopy(ijmllm,pbarv,1,pbarvm,1)
    5764
    58 c
     65      if (groupe_ok) then
    5966      call groupeun(jjp1,llm,zconvmm)
    6067      call groupeun(jjm,llm,pbarvm)
    6168
    6269c   Champs 3D
    63 
    6470      do l=1,llm
    6571         do j=2,jjm
     
    7480         enddo
    7581      enddo
     82
     83      else
     84         pbarum(:,:,:)=pbaru(:,:,:)
     85         pbarvm(:,:,:)=pbarv(:,:,:)
     86      endif
    7687
    7788c    integration de la convergence de masse de haut  en bas ......
  • LMDZ5/branches/testing/libf/dyn3d/inidissip.F90

    r1665 r1707  
    33!
    44SUBROUTINE inidissip ( lstardis,nitergdiv,nitergrot,niterh  , &
    5      tetagdiv,tetagrot,tetatemp             )
     5     tetagdiv,tetagrot,tetatemp, vert_prof_dissip)
    66  !=======================================================================
    77  !   initialisation de la dissipation horizontale
     
    2525  INTEGER,INTENT(in) :: nitergdiv,nitergrot,niterh
    2626  REAL,INTENT(in) :: tetagdiv,tetagrot,tetatemp
     27
     28  integer, INTENT(in):: vert_prof_dissip
     29  ! Vertical profile of horizontal dissipation
     30  ! Allowed values:
     31  ! 0: rational fraction, function of pressure
     32  ! 1: tanh of altitude
    2733
    2834! Local variables:
     
    167173  !   --------------------------------------------------
    168174
    169   if (ok_strato .and. llm==39) then
     175  if (vert_prof_dissip == 1) then
    170176     do l=1,llm
    171177        pseudoz=8.*log(preff/presnivs(l))
  • LMDZ5/branches/testing/libf/dyn3d/leapfrog.F

    r1669 r1707  
    383383           jD_cur = jD_ref + day_ini - day_ref +                        &
    384384     &          itau/day_step
     385
     386           IF (planet_type .eq."generic") THEN
     387              ! AS: we make jD_cur to be pday
     388              jD_cur = int(day_ini + itau/day_step)
     389           ENDIF
     390
    385391           jH_cur = jH_ref + start_time +                               &
    386392     &              mod(itau,day_step)/float(day_step)
  • LMDZ5/branches/testing/libf/dyn3d/paramet.h

    r792 r1707  
    1717      INTEGER jcfil,jcfllm
    1818
    19       PARAMETER( iip1= iim+1-1/iim,iip2=iim+2,iip3=iim+3                &
     19      PARAMETER( iip1= iim+1,iip2=iim+2,iip3=iim+3                       &
    2020     &    ,jjp1=jjm+1-1/jjm)
    2121      PARAMETER( llmp1 = llm+1,  llmp2 = llm+2, llmm1 = llm-1 )
  • LMDZ5/branches/testing/libf/dyn3dmem/abort_gcm.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/academic.h

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/addfi_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/advtrac_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/bands.F90

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/bilan_dyn_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/caladvtrac_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/calfis_loc.F

    r1669 r1707  
    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     .             jD_cur,         !! pday. see leapfrog_p
     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/branches/testing/libf/dyn3dmem/call_calfis_mod.F90

    r1669 r1707  
    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
     141           IF (planet_type .eq."generic") THEN
     142              ! AS: we make jD_cur to be pday
     143              jD_cur = int(day_ini + itau/day_step)
     144           ENDIF
     145
     146           jH_cur = jH_ref + start_time +                                &
     147     &              mod(itau,day_step)/float(day_step)
     148    if (jH_cur > 1.0 ) then
     149      jD_cur = jD_cur +1.
     150      jH_cur = jH_cur -1.
     151    endif
    140152
    141153!   Inbterface avec les routines de phylmd (phymars ... )
  • LMDZ5/branches/testing/libf/dyn3dmem/call_dissip_mod.F90

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/ce0l.F90

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/comconst.h

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/comdissipn.h

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/comdissnew.h

    r1669 r1707  
    1212
    1313      COMMON/comdissnew/ lstardis,nitergdiv,nitergrot,niterh,tetagdiv,  &
    14      &                   tetagrot,tetatemp,coefdis 
     14     &                   tetagrot,tetatemp,coefdis, vert_prof_dissip
    1515
    1616      LOGICAL lstardis
    1717      INTEGER nitergdiv, nitergrot, niterh
     18
     19      integer vert_prof_dissip ! vertical profile of horizontal dissipation
     20!     Allowed values:
     21!     0: rational fraction, function of pressure
     22!     1: tanh of altitude
     23
    1824      REAL     tetagdiv, tetagrot,  tetatemp, coefdis
    1925
  • LMDZ5/branches/testing/libf/dyn3dmem/comvert.h

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/conf_gcm.F

    r1669 r1707  
    11!
    2 ! $Id: conf_gcm.F 1403 2010-07-01 09:02:53Z fairhead $
     2! $Id$
    33!
    44c
     
    66      SUBROUTINE conf_gcm( tapedef, etatinit, clesphy0 )
    77c
     8      USE control_mod
    89#ifdef CPP_IOIPSL
    910      use IOIPSL
     
    1718      use mod_hallo, ONLY : use_mpi_alloc
    1819      use parallel, ONLY : omp_chunk
    19       USE control_mod
     20      USE infotrac, ONLY : type_trac
     21      use assert_m, only: assert
     22
    2023      IMPLICIT NONE
    2124c-----------------------------------------------------------------------
     
    4346#include "serre.h"
    4447#include "comdissnew.h"
    45 !#include "clesphys.h"
    46 #include "iniprint.h"
    4748#include "temps.h"
    4849#include "comconst.h"
    4950
    5051! FH 2008/05/09 On elimine toutes les clefs physiques dans la dynamique
     52! #include "clesphys.h"
     53#include "iniprint.h"
    5154c
    5255c
     
    103106      CALL getin('lunout', lunout)
    104107      IF (lunout /= 5 .and. lunout /= 6) THEN
    105         OPEN(lunout,FILE='lmdz.out')
     108        OPEN(UNIT=lunout,FILE='lmdz.out_0000',ACTION='write',
     109     &          STATUS='unknown',FORM='formatted')
    106110      ENDIF
    107111
     
    166170      CALL getin('nday',nday)
    167171
     172!Config  Key  = starttime
     173!Config  Desc = Heure de depart de la simulation
     174!Config  Def  = 0
     175!Config  Help = Heure de depart de la simulation
     176!Config         en jour
     177      starttime = 0
     178      CALL getin('starttime',starttime)
     179
    168180!Config  Key  = day_step
    169181!Config  Desc = nombre de pas par jour
     
    175187
    176188!Config  Key  = nsplit_phys
    177 !Config  Desc = nombre d'iteration de la physique
    178 !Config  Def  = 240
    179 !Config  Help = nombre d'itration de la physique
    180 !
    181189       nsplit_phys = 1
    182190       CALL getin('nsplit_phys',nsplit_phys)
     
    226234       CALL getin('output_grads_dyn',output_grads_dyn)
    227235
    228 !Config  Key  = idissip
     236!Config  Key  = dissip_period
    229237!Config  Desc = periode de la dissipation
    230 !Config  Def  = 10
     238!Config  Def  = 0
    231239!Config  Help = periode de la dissipation
    232 !Config         (en pas) ... a completer !
    233        idissip = 10
    234        CALL getin('idissip',idissip)
     240!Config  dissip_period=0 => la valeur sera calcule dans inidissip       
     241!Config  dissip_period>0 => on prend cette valeur
     242       dissip_period = 0
     243       CALL getin('dissip_period',dissip_period)
    235244
    236245ccc  ....   P. Le Van , modif le 29/04/97 .pour la dissipation  ...
     
    314323       CALL getin('tau_top_bound',tau_top_bound)
    315324
    316 !
    317325!Config  Key  = coefdis
    318326!Config  Desc = coefficient pour gamdissip
     
    579587       offline = .FALSE.
    580588       CALL getin('offline',offline)
     589       IF (offline .AND. adjust) THEN
     590          WRITE(lunout,*)
     591     &         'WARNING : option offline does not work with adjust=y :'
     592          WRITE(lunout,*) 'the files defstoke.nc, fluxstoke.nc ',
     593     &         'and fluxstokev.nc will not be created'
     594          WRITE(lunout,*)
     595     &         'only the file phystoke.nc will still be created '
     596       END IF
     597       
     598!Config  Key  = type_trac
     599!Config  Desc = Choix de couplage avec model de chimie INCA ou REPROBUS
     600!Config  Def  = lmdz
     601!Config  Help =
     602!Config         'lmdz' = pas de couplage, pur LMDZ
     603!Config         'inca' = model de chime INCA
     604!Config         'repr' = model de chime REPROBUS
     605      type_trac = 'lmdz'
     606      CALL getin('type_trac',type_trac)
    581607
    582608!Config  Key  = config_inca
     
    628654      write(lunout,*)' periodav = ', periodav
    629655      write(lunout,*)' output_grads_dyn = ', output_grads_dyn
    630       write(lunout,*)' idissip = ', idissip
     656      write(lunout,*)' dissip_period = ', dissip_period
    631657      write(lunout,*)' lstardis = ', lstardis
    632658      write(lunout,*)' nitergdiv = ', nitergdiv
     
    651677      write(lunout,*)' tauyy = ', tauyy
    652678      write(lunout,*)' offline = ', offline
     679      write(lunout,*)' type_trac = ', type_trac
    653680      write(lunout,*)' config_inca = ', config_inca
    654681      write(lunout,*)' ok_dynzon = ', ok_dynzon
     
    769796       offline = .FALSE.
    770797       CALL getin('offline',offline)
     798       IF (offline .AND. adjust) THEN
     799          WRITE(lunout,*)
     800     &         'WARNING : option offline does not work with adjust=y :'
     801          WRITE(lunout,*) 'the files defstoke.nc, fluxstoke.nc ',
     802     &         'and fluxstokev.nc will not be created'
     803          WRITE(lunout,*)
     804     &         'only the file phystoke.nc will still be created '
     805       END IF
     806
     807!Config  Key  = type_trac
     808!Config  Desc = Choix de couplage avec model de chimie INCA ou REPROBUS
     809!Config  Def  = lmdz
     810!Config  Help =
     811!Config         'lmdz' = pas de couplage, pur LMDZ
     812!Config         'inca' = model de chime INCA
     813!Config         'repr' = model de chime REPROBUS
     814      type_trac = 'lmdz'
     815      CALL getin('type_trac',type_trac)
    771816
    772817!Config  Key  = config_inca
     
    781826
    782827!Config  Key  = ok_dynzon
    783 !Config  Desc = calcul et sortie des transports
     828!Config  Desc = sortie des transports zonaux dans la dynamique
    784829!Config  Def  = n
    785830!Config  Help = Permet de mettre en route le calcul des transports
     
    817862        write(lunout,*)"Le zoom en longitude est incompatible",
    818863     &                 " avec l'utilisation du filtre FFT ",
    819      &                 "---> filtre FFT désactivé "
     864     &                 "---> FFT filter not active"
    820865       use_filtre_fft=.FALSE.
    821866      ENDIF
     
    851896      CALL getin('ok_strato',ok_strato)
    852897
     898      vert_prof_dissip = merge(1, 0, ok_strato .and. llm==39)
     899      CALL getin('vert_prof_dissip', vert_prof_dissip)
     900      call assert(vert_prof_dissip == 0 .or. vert_prof_dissip ==  1,
     901     $     "bad value for vert_prof_dissip")
     902
    853903!Config  Key  = ok_gradsfile
    854904!Config  Desc = activation des sorties grads du guidage
     
    874924      ok_etat0 = .TRUE.
    875925      CALL getin('ok_etat0',ok_etat0)
     926
     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)
    876932
    877933      write(lunout,*)' #########################################'
     
    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
    915       write(lunout,*)' ok_dynzon = ', ok_dynzon 
     972      write(lunout,*)' ok_dynzon = ', ok_dynzon
    916973      write(lunout,*)' ok_dyn_ins = ', ok_dyn_ins
    917974      write(lunout,*)' ok_dyn_ave = ', ok_dyn_ave
     
    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/branches/testing/libf/dyn3dmem/control_mod.F90

    r1669 r1707  
    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
     
    2525                     ! in NetCDF files dyn_hist*ave.nc
    2626
    27   LOGICAL ok_dynzon  ! output zonal transports in dynzon.nc file
    28   LOGICAL ok_dyn_ins ! output instantaneous values of fields
    29                      ! in the dynamics in NetCDF files dyn_hist*nc
    30   LOGICAL ok_dyn_ave ! output averaged values of fields in the dynamics
    31                      ! in NetCDF files dyn_hist*ave.nc
    32 
    3327END MODULE
  • LMDZ5/branches/testing/libf/dyn3dmem/defrun.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/dynetat0_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/dynredem_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/dynredem_mod.F90

    r1669 r1707  
     1!
     2! $Id$
     3!
    14MODULE dynredem_mod
    2 
    3 
    4 
    5 
    65
    76CONTAINS
     
    3736   ENDIF
    3837!$OMP END MASTER
    39  
    40     ll=size(var,2)
    4138
    4239!$OMP MASTER
     
    105102!$OMP END MASTER
    106103 
    107     ll=size(var,2)
    108 
    109104!$OMP MASTER
    110105    ALLOCATE(var_tmp(ijb_v:ije_v,ll))
     
    172167!$OMP END MASTER
    173168 
    174     ll=size(var,2)
    175 
    176169!$OMP MASTER
    177170    ALLOCATE(var_tmp(ijb_u:ije_u,ll))
  • LMDZ5/branches/testing/libf/dyn3dmem/ener.h

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/exner_hyb.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/exner_hyb_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/exner_milieu_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/filtreg_p.F

    r1669 r1707  
    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(:iim,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(:iim,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(:iim,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(:iim,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(:iim,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(:iim,j-jdfil+1,:)
     290     &                    =matmul(matricevs(:,:,j-jfiltsv+1),
     291     &                            champ_loc(:iim,j,:))
     292#endif
    260293                  ENDDO
    261294                 
  • LMDZ5/branches/testing/libf/dyn3dmem/friction_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/gcm.F

    r1669 r1707  
    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-----------------------------------------------------------------------
     
    429420
    430421      CALL inidissip( lstardis, nitergdiv, nitergrot, niterh   ,
    431      *                tetagdiv, tetagrot , tetatemp              )
     422     *                tetagdiv, tetagrot , tetatemp, vert_prof_dissip)
    432423
    433424c-----------------------------------------------------------------------
    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/branches/testing/libf/dyn3dmem/gr_dyn_fi_p.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/gr_fi_dyn_p.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/grid_noro.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/guide_loc_mod.F90

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/infotrac.F90

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/inigrads.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/initfluxsto_p.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/integrd_loc.F

    r1669 r1707  
    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 'mpif.h'
     39#include "iniprint.h"
     40!      include 'mpif.h'
    4041
    4142c   Arguments:
    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/branches/testing/libf/dyn3dmem/leapfrog_loc.F

    r1669 r1707  
    11!
    2 ! $Id: leapfrog_p.F 1299 2010-01-20 14:27:21Z fairhead $
     2! $Id$
    33!
    44c
     
    7878#include "iniprint.h"
    7979#include "academic.h"
    80       include "mpif.h"
     80!      include "mpif.h"
    8181     
    8282      INTEGER         longcles
     
    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
     
    236237
    237238c Allocate variables depending on dynamic variable nqtot
    238 !c$OMP MASTER
     239!$OMP MASTER
     240      if (firstcall) then
    239241!     
    240242!      ALLOCATE(p(ijb_u:ije_u,llmp1))
     
    252254!      ALLOCATE(dvdis(ijb_v:ije_v,llm),dudis(ijb_u:ije_u,llm))
    253255!      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))
     256      ALLOCATE(dvfi(ijb_v:ije_v,llm),dufi(ijb_u:ije_u,llm))
     257      ALLOCATE(dtetafi(ijb_u:ije_u,llm))
     258      ALLOCATE(dpfi(ijb_u:ije_u))
    257259!      ALLOCATE(dq(ijb_u:ije_u,llm,nqtot))
    258 !      ALLOCATE(dqfi(ijb_u:ije_u,llm,nqtot))
     260      ALLOCATE(dqfi(ijb_u:ije_u,llm,nqtot))
    259261!      ALLOCATE(dqfi_tmp(iip1,llm,nqtot))
    260262!      ALLOCATE(finvmaold(ijb_u:ije_u,llm))
     
    265267!      ALLOCATE(vcont(ijb_v:ije_v,llm),ucont(ijb_u:ije_u,llm))
    266268!      ALLOCATE(vnat(ijb_v:ije_v,llm),unat(ijb_u:ije_u,llm))
    267 !c$OMP END MASTER     
    268 !c$OMP BARRIER
     269      endif
     270!$OMP END MASTER     
     271!$OMP BARRIER
    269272
    270273!                CALL dynredem1_loc("restart.nc",0.0,
     
    277280
    278281c$OMP MASTER
    279       dq=0.
     282      dq(:,:,:)=0.
    280283      CALL pression ( ijnb_u, ap, bp, ps, p       )
    281284c$OMP END MASTER
     285      if (pressure_exner) then
    282286      CALL exner_hyb_loc( ijnb_u, ps, p,alpha,beta, pks, pk, pkf)
    283 
     287      else
     288        CALL exner_milieu_loc( ijnb_u, ps, p, beta, pks, pk, pkf )
     289      endif
    284290c-----------------------------------------------------------------------
    285291c   Debut de l'integration temporelle:
     
    287293c et du parallelisme !!
    288294
    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))
     295   1  CONTINUE ! Matsuno Forward step begins here
     296
     297      jD_cur = jD_ref + day_ini - day_ref +                             &
     298     &          itau/day_step
     299      jH_cur = jH_ref + start_time +                                    &
     300     &         mod(itau,day_step)/float(day_step)
     301      if (jH_cur > 1.0 ) then
     302        jD_cur = jD_cur +1.
     303        jH_cur = jH_cur -1.
     304      endif
     305
    294306
    295307#ifdef CPP_IOIPSL
     
    323335         psm1= ps
    324336         
    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 )
     337! Ehouarn: finvmaold is actually not used       
     338!         finvmaold = masse
     339c$OMP END MASTER
     340c$OMP BARRIER
     341!         CALL filtreg_p ( finvmaold ,jjb_u,jje_u,jjb_u,jje_u,jjp1, llm,
     342!     &                    -2,2, .TRUE., 1 )
    330343       else
    331344! Save fields obtained at previous time step as '...m1'
     
    343356           tetam1   (ijb:ije,l) = teta  (ijb:ije,l)
    344357           massem1  (ijb:ije,l) = masse (ijb:ije,l)
    345            finvmaold(ijb:ije,l)=masse(ijb:ije,l)
     358!           finvmaold(ijb:ije,l)=masse(ijb:ije,l)
    346359                 
    347360           if (pole_sud) ije=ij_end-iip1
     
    353366
    354367
    355           CALL filtreg_p(finvmaold ,jjb_u,jje_u,jj_begin,jj_end,jjp1,
    356      .                    llm, -2,2, .TRUE., 1 )
     368! Ehouarn: finvmaold not used
     369!          CALL filtreg_p(finvmaold ,jjb_u,jje_u,jj_begin,jj_end,jjp1,
     370!     .                    llm, -2,2, .TRUE., 1 )
    357371
    358372       endif ! of if (FirstCaldyn)
     
    370384cym      call minmax(ijp1llm,q(:,:,3),zqmin,zqmax)
    371385
    372    2  CONTINUE
     386   2  CONTINUE ! Matsuno backward or leapfrog step begins here
    373387
    374388c$OMP MASTER
    375389      ItCount=ItCount+1
    376       if (MOD(ItCount,1)==0) then
     390      if (MOD(ItCount,1)==1) then
    377391        debug=.true.
    378392      else
     
    399413      ! Purely Matsuno time stepping
    400414         IF( MOD(itau,iconser) .EQ.0.AND.  forward    ) conser = .TRUE.
    401          IF( MOD(itau,idissip ).EQ.0.AND..NOT.forward ) apdiss = .TRUE.
     415         IF( MOD(itau,dissip_period ).EQ.0.AND..NOT.forward )
     416     s        apdiss = .TRUE.
    402417         IF( MOD(itau,iphysiq ).EQ.0.AND..NOT.forward
    403      s          .and. iflag_phys.EQ.1                 ) apphys = .TRUE.
     418     s          .and. physic                        ) apphys = .TRUE.
    404419      ELSE
    405420      ! Leapfrog/Matsuno time stepping
    406421         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.
     422         IF( MOD(itau+1,dissip_period).EQ.0 .AND. .NOT. forward )
     423     s        apdiss = .TRUE.
     424         IF( MOD(itau+1,iphysiq).EQ.0.AND.physic) apphys=.TRUE.
    409425      END IF
    410426
     
    450466c$OMP MASTER
    451467           call allgather_timer_average
    452         verbose=.TRUE.
    453         if (Verbose) then
     468
     469        if (prt_level > 9) then
    454470       
    455471        print *,'*********************************'
     
    622638      call start_timer(timer_caldyn)
    623639
     640      ! compute geopotential phi()
    624641      CALL geopot_loc  ( ip1jmp1, teta  , pk , pks,  phis  , phi   )
    625642
     
    697714
    698715       CALL integrd_loc ( 2,vcovm1,ucovm1,tetam1,psm1,massem1 ,
    699      $         dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis ,
    700      $              finvmaold                                    )
     716     $         dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis)
     717!     $              finvmaold                                    )
    701718
    702719!       CALL FTRACE_REGION_END("integrd")
     
    10811098!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    10821099       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
     1100       teta(ijb:ije,l)=teta(ijb:ije,l) -dtvr*
     1101     &        (teta(ijb:ije,l)-tetarappel(ijb:ije,l))*
     1102     &                 (knewt_g+knewt_t(l)*clat4(ijb:ije))       
    10851103       enddo
    10861104!$OMP END DO
     1105
     1106!$OMP MASTER
     1107       if (planet_type.eq."giant") then
     1108         ! add an intrinsic heat flux at the base of the atmosphere
     1109         teta(ijb:ije,1) = teta(ijb:ije,1)
     1110     &        + dtvr * aire(ijb:ije) * ihf / cpp / masse(ijb:ije,1)
     1111       endif
     1112!$OMP END MASTER
     1113!$OMP BARRIER
     1114
    10871115
    10881116       call Register_Hallo_u(ucov,llm,0,1,1,0,Request_Physic)
     
    10921120       call WaitRequest(Request_Physic)     
    10931121c$OMP BARRIER
    1094        call friction_loc(ucov,vcov,iphysiq*dtvr)
     1122       call friction_loc(ucov,vcov,dtvr)
    10951123!$OMP BARRIER
     1124
     1125        ! Sponge layer (if any)
     1126        IF (ok_strato) THEN
     1127          ! set dufi,dvfi,... to zero
     1128          ijb=ij_begin
     1129          ije=ij_end
     1130!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     1131          do l=1,llm
     1132            dufi(ijb:ije,l)=0
     1133            dtetafi(ijb:ije,l)=0
     1134            dqfi(ijb:ije,l,1:nqtot)=0
     1135          enddo
     1136!$OMP END DO
     1137!$OMP MASTER
     1138          dpfi(ijb:ije)=0
     1139!$OMP END MASTER
     1140          ijb=ij_begin
     1141          ije=ij_end
     1142          if (pole_sud) ije=ije-iip1
     1143!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     1144          do l=1,llm
     1145            dvfi(ijb:ije,l)=0
     1146          enddo
     1147!$OMP END DO
     1148
     1149          CALL top_bound_loc(vcov,ucov,teta,masse,dufi,dvfi,dtetafi)
     1150          CALL addfi_loc( dtvr, leapf, forward   ,
     1151     $                  ucov, vcov, teta , q   ,ps ,
     1152     $                 dufi, dvfi, dtetafi , dqfi ,dpfi  )
     1153!$OMP BARRIER
     1154        ENDIF ! of IF (ok_strato)
    10961155      ENDIF ! of IF(iflag_phys.EQ.2)
    10971156
     
    10991158        CALL pression_loc ( ip1jmp1, ap, bp, ps, p                  )
    11001159c$OMP BARRIER
    1101         CALL exner_hyb_loc( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
     1160        if (pressure_exner) then
     1161        CALL exner_hyb_loc( ijnb_u, ps, p,alpha,beta, pks, pk, pkf )
     1162        else
     1163          CALL exner_milieu_loc( ijnb_u, ps, p, beta, pks, pk, pkf )
     1164        endif
    11021165c$OMP BARRIER
    11031166
     
    14961559c$OMP BARRIER
    14971560
    1498               if (planet_type.eq."earth") then
     1561!              if (planet_type.eq."earth") then
    14991562! Write an Earth-format restart file
    15001563                CALL dynredem1_loc("restart.nc",0.0,
    15011564     &                           vcov,ucov,teta,q,masse,ps)
    1502               endif ! of if (planet_type.eq."earth")
     1565!              endif ! of if (planet_type.eq."earth")
    15031566
    15041567!              CLOSE(99)
     
    16081671
    16091672              IF(itau.EQ.itaufin) THEN
    1610                 if (planet_type.eq."earth") then
     1673!                if (planet_type.eq."earth") then
    16111674                   CALL dynredem1_loc("restart.nc",0.0,
    16121675     .                               vcov,ucov,teta,q,masse,ps)
    1613                 endif ! of if (planet_type.eq."earth")
     1676!               endif ! of if (planet_type.eq."earth")
    16141677              ENDIF ! of IF(itau.EQ.itaufin)
    16151678
  • LMDZ5/branches/testing/libf/dyn3dmem/limy.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/logic.h

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/mod_filtreg_p.F

    r1669 r1707  
    210210               IF( ifiltre.EQ.-2 )   THEN
    211211                  DO j = jdfil,jffil
    212                      CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
     212#ifdef BLAS
     213                     CALL SGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    213214     &                    matrinvn(1,1,j), iim,
    214215     &                    champ_loc(1,j,1), iip1*(jje-jjb+1), 0.0,
    215216     &                    champ_fft(1,j,1), iip1*(jje-jjb+1))
     217#else
     218                     champ_fft(:,j,:)=
     219     &                    matmul(matrinvn(:,:,j),champ_loc(:iim,j,:))
     220#endif
    216221                  ENDDO
    217222                 
    218223               ELSE IF ( griscal )     THEN
    219224                  DO j = jdfil,jffil
    220                      CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
     225#ifdef BLAS
     226                     CALL SGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    221227     &                    matriceun(1,1,j), iim,
    222228     &                    champ_loc(1,j,1), iip1*(jje-jjb+1), 0.0,
    223229     &                    champ_fft(1,j,1), iip1*(jje-jjb+1))
     230#else
     231                     champ_fft(:,j,:)=
     232     &                    matmul(matriceun(:,:,j),champ_loc(:iim,j,:))
     233#endif
    224234                  ENDDO
    225235                 
    226236               ELSE
    227237                  DO j = jdfil,jffil
    228                      CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
     238#ifdef BLAS
     239                     CALL SGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    229240     &                    matricevn(1,1,j), iim,
    230241     &                    champ_loc(1,j,1), iip1*(jje-jjb+1), 0.0,
    231242     &                    champ_fft(1,j,1), iip1*(jje-jjb+1))
     243#else
     244                     champ_fft(:,j,:)=
     245     &                    matmul(matricevn(:,:,j),champ_loc(:iim,j,:))
     246#endif
    232247                  ENDDO
    233248                 
     
    238253               IF( ifiltre.EQ.-2 )   THEN
    239254                  DO j = jdfil,jffil
    240                      CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
     255#ifdef BLAS
     256                     CALL SGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    241257     &                    matrinvs(1,1,j-jfiltsu+1), iim,
    242258     &                    champ_loc(1,j,1), iip1*(jje-jjb+1), 0.0,
    243259     &                    champ_fft(1,j,1), iip1*(jje-jjb+1))
     260#else
     261                     champ_fft(:,j,:)=
     262     &                    matmul(matrinvs(:,:,j-jfiltsu+1),
     263     &                           champ_loc(:iim,j,:))
     264#endif
    244265                  ENDDO
    245266                 
     
    247268                 
    248269                  DO j = jdfil,jffil
    249                      CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
     270#ifdef BLAS
     271                     CALL SGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    250272     &                    matriceus(1,1,j-jfiltsu+1), iim,
    251273     &                    champ_loc(1,j,1), iip1*(jje-jjb+1), 0.0,
    252274     &                    champ_fft(1,j,1), iip1*(jje-jjb+1))
     275#else
     276                     champ_fft(:,j,:)=
     277     &                    matmul(matriceus(:,:,j-jfiltsu+1),
     278     &                           champ_loc(:iim,j,:))
     279#endif
    253280                  ENDDO
    254281                 
     
    256283                 
    257284                  DO j = jdfil,jffil
    258                      CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
     285#ifdef BLAS
     286                     CALL SGEMM("N", "N", iim, nbniv_loc, iim, 1.0,
    259287     &                    matricevs(1,1,j-jfiltsv+1), iim,
    260288     &                    champ_loc(1,j,1), iip1*(jje-jjb+1), 0.0,
    261289     &                    champ_fft(1,j,1), iip1*(jje-jjb+1))
     290#else
     291                     champ_fft(:,j,:)=
     292     &                    matmul(matricevs(:,:,j-jfiltsv+1),
     293     &                           champ_loc(:iim,j,:))
     294#endif
    262295                  ENDDO
    263296                 
     
    269302               
    270303c     !-------------------------------------!
    271 c     ! Dés-agregation des niveau verticaux !
     304c     ! Dés-agregation des niveau verticaux !
    272305c     ! uniquement necessaire pour une      !
    273306c     ! execution OpenMP                    !
     
    402435      END SUBROUTINE filtreg_p
    403436      END MODULE mod_filtreg_p
     437
  • LMDZ5/branches/testing/libf/dyn3dmem/mod_interface_dyn_phys.F90

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/parallel.F90

    r1669 r1707  
    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
     
    643657          enddo
    644658         
    645         endif
     659        else
     660          ! Ehouarn: When in debug mode, ifort complains (for call MPI_GATHERV
     661          !          below) about Buffer_Recv() being not allocated.
     662          !          So make a dummy allocation.
     663          allocate(Buffer_Recv(1))
     664        endif ! of if (MPI_Rank==rank)
    646665 
    647666!$OMP CRITICAL (MPI)
     
    717736       
    718737   
    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 */
     738!  Subroutine verif_hallo(Field,ij,ll,up,down)
     739!    implicit none
     740!#include "dimensions.h"
     741!#include "paramet.h"   
     742!    include 'mpif.h'
     743!   
     744!      INTEGER :: ij,ll
     745!      REAL, dimension(ij,ll) :: Field
     746!      INTEGER :: up,down
     747!     
     748!      REAL,dimension(ij,ll): NewField
     749!     
     750!      NewField=0
     751!     
     752!      ijb=ij_begin
     753!      ije=ij_end
     754!      if (pole_nord)
     755!      NewField(ij_be       
     756
    739757  end module parallel
  • LMDZ5/branches/testing/libf/dyn3dmem/paramet.h

    r1669 r1707  
    11!
    2 ! $Header$
     2! $Id$
    33!
    44!
  • LMDZ5/branches/testing/libf/dyn3dmem/temps.h

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/vlsplt_loc.F

    r1669 r1707  
     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/branches/testing/libf/dyn3dmem/vlspltqs_loc.F

    r1669 r1707  
    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/branches/testing/libf/dyn3dmem/wrgrads.F

    r1669 r1707  
    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'
  • LMDZ5/branches/testing/libf/dyn3dpar/calfis_p.F

    r1669 r1707  
    684684     .             debut_split,    !! firstcall
    685685     .             lafin_split,    !! lastcall
    686      .             float(day_ini), !! pday <-- day_ini (dans temps.h)
     686     .             jD_cur,         !! pday. see leapfrog_p
    687687     .             jH_cur_split,   !! ptime "fraction of day"
    688688     .             zdt_split,      !! ptimestep
  • LMDZ5/branches/testing/libf/dyn3dpar/comconst.h

    r1505 r1707  
    2121      REAL dtdiss ! (s) time step for the dissipation
    2222      REAL rad ! (m) radius of the planet
    23       REAL r ! Gas constant R=8.31 J.K-1.mol-1
    24       REAL cpp   ! Cp
     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)
    2526      REAL kappa ! kappa=R/Cp
    2627      REAL cotot
  • LMDZ5/branches/testing/libf/dyn3dpar/comdissnew.h

    r1319 r1707  
    1212
    1313      COMMON/comdissnew/ lstardis,nitergdiv,nitergrot,niterh,tetagdiv,  &
    14      &                   tetagrot,tetatemp,coefdis 
     14     &                   tetagrot,tetatemp,coefdis, vert_prof_dissip
    1515
    1616      LOGICAL lstardis
    1717      INTEGER nitergdiv, nitergrot, niterh
     18
     19      integer vert_prof_dissip ! vertical profile of horizontal dissipation
     20!     Allowed values:
     21!     0: rational fraction, function of pressure
     22!     1: tanh of altitude
     23
    1824      REAL     tetagdiv, tetagrot,  tetatemp, coefdis
    1925
  • LMDZ5/branches/testing/libf/dyn3dpar/conf_gcm.F

    r1665 r1707  
    66      SUBROUTINE conf_gcm( tapedef, etatinit, clesphy0 )
    77c
     8      USE control_mod
    89#ifdef CPP_IOIPSL
    910      use IOIPSL
     
    1617      use mod_hallo, ONLY : use_mpi_alloc
    1718      use parallel, ONLY : omp_chunk
    18       USE control_mod
    1919      USE infotrac, ONLY : type_trac
     20      use assert_m, only: assert
     21
    2022      IMPLICIT NONE
    2123c-----------------------------------------------------------------------
     
    4345#include "serre.h"
    4446#include "comdissnew.h"
    45 !#include "clesphys.h"
    46 #include "iniprint.h"
    4747#include "temps.h"
    4848#include "comconst.h"
    4949
    5050! FH 2008/05/09 On elimine toutes les clefs physiques dans la dynamique
     51! #include "clesphys.h"
     52#include "iniprint.h"
    5153c
    5254c
     
    105107        OPEN(UNIT=lunout,FILE='lmdz.out_0000',ACTION='write',
    106108     &          STATUS='unknown',FORM='formatted')
    107 
    108109      ENDIF
    109110
     
    185186
    186187!Config  Key  = nsplit_phys
    187 !Config  Desc = nombre d'iteration de la physique
    188 !Config  Def  = 240
    189 !Config  Help = nombre d'itration de la physique
    190 !
    191188       nsplit_phys = 1
    192189       CALL getin('nsplit_phys',nsplit_phys)
     
    325322       CALL getin('tau_top_bound',tau_top_bound)
    326323
    327 !
    328324!Config  Key  = coefdis
    329325!Config  Desc = coefficient pour gamdissip
     
    608604      type_trac = 'lmdz'
    609605      CALL getin('type_trac',type_trac)
    610 
    611606
    612607!Config  Key  = config_inca
     
    830825
    831826!Config  Key  = ok_dynzon
    832 !Config  Desc = calcul et sortie des transports
     827!Config  Desc = sortie des transports zonaux dans la dynamique
    833828!Config  Def  = n
    834829!Config  Help = Permet de mettre en route le calcul des transports
     
    865860        write(lunout,*)"Le zoom en longitude est incompatible",
    866861     &                 " avec l'utilisation du filtre FFT ",
    867      &                 "---> filtre FFT désactivé "
     862     &                 "---> FFT filter not active"
    868863       use_filtre_fft=.FALSE.
    869864      ENDIF
     
    898893      ok_strato=.FALSE.
    899894      CALL getin('ok_strato',ok_strato)
     895
     896      vert_prof_dissip = merge(1, 0, ok_strato .and. llm==39)
     897      CALL getin('vert_prof_dissip', vert_prof_dissip)
     898      call assert(vert_prof_dissip == 0 .or. vert_prof_dissip ==  1,
     899     $     "bad value for vert_prof_dissip")
    900900
    901901!Config  Key  = ok_gradsfile
     
    968968      write(lunout,*)' type_trac = ', type_trac
    969969      write(lunout,*)' config_inca = ', config_inca
    970       write(lunout,*)' ok_dynzon = ', ok_dynzon 
     970      write(lunout,*)' ok_dynzon = ', ok_dynzon
    971971      write(lunout,*)' ok_dyn_ins = ', ok_dyn_ins
    972972      write(lunout,*)' ok_dyn_ave = ', ok_dyn_ave
  • LMDZ5/branches/testing/libf/dyn3dpar/filtreg_p.F

    r1665 r1707  
    214214     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
    215215#else
    216                      champ_fft(:,j-jdfil+1,:)
     216                     champ_fft(:iim,j-jdfil+1,:)
    217217     &                    =matmul(matrinvn(:,:,j),champ_loc(:iim,j,:))
    218218#endif
     
    227227     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
    228228#else
    229                      champ_fft(:,j-jdfil+1,:)
     229                     champ_fft(:iim,j-jdfil+1,:)
    230230     &                    =matmul(matriceun(:,:,j),champ_loc(:iim,j,:))
    231231#endif
     
    240240     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
    241241#else
    242                      champ_fft(:,j-jdfil+1,:)
     242                     champ_fft(:iim,j-jdfil+1,:)
    243243     &                    =matmul(matricevn(:,:,j),champ_loc(:iim,j,:))
    244244#endif
     
    257257     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
    258258#else
    259                      champ_fft(:,j-jdfil+1,:)
     259                     champ_fft(:iim,j-jdfil+1,:)
    260260     &                    =matmul(matrinvs(:,:,j-jfiltsu+1),
    261261     &                            champ_loc(:iim,j,:))
     
    272272     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
    273273#else
    274                      champ_fft(:,j-jdfil+1,:)
     274                     champ_fft(:iim,j-jdfil+1,:)
    275275     &                    =matmul(matriceus(:,:,j-jfiltsu+1),
    276276     &                            champ_loc(:iim,j,:))
     
    287287     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
    288288#else
    289                      champ_fft(:,j-jdfil+1,:)
     289                     champ_fft(:iim,j-jdfil+1,:)
    290290     &                    =matmul(matricevs(:,:,j-jfiltsv+1),
    291291     &                            champ_loc(:iim,j,:))
  • LMDZ5/branches/testing/libf/dyn3dpar/fxhyp.F

    r1403 r1707  
    4848c
    4949       REAL   dzoom
    50        REAL*8 xlon(iip1),xprimm(iip1),xuv
    51        REAL*8 xtild(0:nmax2)
    52        REAL*8 fhyp(0:nmax2),ffdx,beta,Xprimt(0:nmax2)
    53        REAL*8 Xf(0:nmax2),xxpr(0:nmax2)
    54        REAL*8 xvrai(iip1),xxprim(iip1)
    55        REAL*8 pi,depi,epsilon,xzoom,fa,fb
    56        REAL*8 Xf1, Xfi , a0,a1,a2,a3,xi2
     50       REAL(KIND=8) xlon(iip1),xprimm(iip1),xuv
     51       REAL(KIND=8) xtild(0:nmax2)
     52       REAL(KIND=8) fhyp(0:nmax2),ffdx,beta,Xprimt(0:nmax2)
     53       REAL(KIND=8) Xf(0:nmax2),xxpr(0:nmax2)
     54       REAL(KIND=8) xvrai(iip1),xxprim(iip1)
     55       REAL(KIND=8) pi,depi,epsilon,xzoom,fa,fb
     56       REAL(KIND=8) Xf1, Xfi , a0,a1,a2,a3,xi2
    5757       INTEGER i,it,ik,iter,ii,idif,ii1,ii2
    58        REAL*8 xi,xo1,xmoy,xlon2,fxm,Xprimin
    59        REAL*8 champmin,champmax,decalx
     58       REAL(KIND=8) xi,xo1,xmoy,xlon2,fxm,Xprimin
     59       REAL(KIND=8) champmin,champmax,decalx
    6060       INTEGER is2
    6161       SAVE is2
    6262
    63        REAL*8 heavyside
     63       REAL(KIND=8) heavyside
    6464
    6565       pi       = 2. * ASIN(1.)
     
    6868       xzoom    = xzoomdeg * pi/180.
    6969c
     70       if (iim==1) then
     71
     72          print*,'Longitudes calculees a la main pour iim=1'
     73
     74          rlonm025(1)=-pi/2.
     75          rlonv(1)=0.
     76          rlonu(1)=pi
     77          rlonp025(1)=pi/2.
     78          rlonm025(2)=rlonm025(1)+depi
     79          rlonv(2)=rlonv(1)+depi
     80          rlonu(2)=rlonu(1)+depi
     81          rlonp025(2)=rlonp025(1)+depi
     82
     83          xprimm025(:)=1.
     84          xprimv(:)=1.
     85          xprimu(:)=1.
     86          xprimp025(:)=1.
     87          champmin=depi
     88          champmax=depi
     89          return
     90
     91       endif
     92
    7093           decalx   = .75
    7194       IF( grossism.EQ.1..AND.scal180 )  THEN
     
    286309
    287310
     311
    288312       IF(ik.EQ.1.and.grossism.EQ.1.)  THEN
    289313         xvrai(1)    = xvrai(iip1)-depi
    290314         xxprim(1)   = xxprim(iip1)
    291315       ENDIF
     316
    292317       DO i = 1 , iim
    293318        xlon(i)     = xvrai(i)
  • LMDZ5/branches/testing/libf/dyn3dpar/gcm.F

    r1665 r1707  
    418418
    419419      CALL inidissip( lstardis, nitergdiv, nitergrot, niterh   ,
    420      *                tetagdiv, tetagrot , tetatemp              )
     420     *                tetagdiv, tetagrot , tetatemp, vert_prof_dissip)
    421421
    422422c-----------------------------------------------------------------------
    423423c   Initialisation de la physique :
    424424c   -------------------------------
    425       IF (call_iniphys.and.iflag_phys.eq.1) THEN
     425      IF (call_iniphys.and.(iflag_phys==1.or.iflag_phys>=100)) THEN
    426426         latfi(1)=rlatu(1)
    427427         lonfi(1)=0.
     
    446446! Physics:
    447447#ifdef CPP_PHYS
    448          CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys ,
    449      ,                latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp     )
     448         CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys,
     449     &                latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp,
     450     &                iflag_phys)
    450451#endif
    451452         call_iniphys=.false.
    452       ENDIF ! of IF (call_iniphys.and.(iflag_phys.eq.1))
     453      ENDIF ! of IF (call_iniphys.and.(iflag_phys==1.or.iflag_phys>=100))
    453454
    454455
     
    481482 301  FORMAT('1'/,15x,'run du ', i2,'/',i2,'/',i4)
    482483 302  FORMAT('1'/,15x,'    au ', i2,'/',i2,'/',i4)
    483 #endif
    484 
    485 #ifdef CPP_PHYS
    486 ! Create start file (startphy.nc) and boundary conditions (limit.nc)
    487 ! for the Earth verstion
    488        if (iflag_phys>=100) then
    489           call iniaqua(ngridmx,latfi,lonfi,iflag_phys)
    490        endif
    491484#endif
    492485
  • LMDZ5/branches/testing/libf/dyn3dpar/groupe_p.F

    r764 r1707  
    3737      integer i,j,l
    3838
    39       logical firstcall
    40       save firstcall
    41 c$OMP THREADPRIVATE(firstcall)
     39      logical firstcall,groupe_ok
     40      save firstcall,groupe_ok
     41c$OMP THREADPRIVATE(firstcall,groupe_ok)
    4242
    4343      data firstcall/.true./
     44      data groupe_ok/.true./
     45
    4446      integer ijb,ije,jjb,jje
    4547     
     48      if (iim==1) then
     49         groupe_ok=.false.
     50      endif
     51
    4652      if (firstcall) then
    47          if(mod(iim,2**ngroup).ne.0) stop'probleme du nombre ede point'
     53         if (groupe_ok) then
     54           if(mod(iim,2**ngroup).ne.0) stop'probleme du nombre de point'
     55         endif
    4856         firstcall=.false.
    4957      endif
     
    6674c$OMP END DO NOWAIT
    6775
    68       call groupeun_p(jjp1,llm,jjb,jje,zconvmm)
     76      if (groupe_ok) then
     77         call groupeun_p(jjp1,llm,jjb,jje,zconvmm)
     78      endif
    6979     
    7080      jjb=jj_begin-1
     
    7888c$OMP END DO NOWAIT
    7989
    80       call groupeun_p(jjm,llm,jjb,jje,pbarvm)
     90      if (groupe_ok) then
     91         call groupeun_p(jjm,llm,jjb,jje,pbarvm)
     92      endif
    8193
    8294c   Champs 3D
     
    101113      enddo
    102114c$OMP END DO NOWAIT
     115
    103116c    integration de la convergence de masse de haut  en bas ......
    104117   
  • LMDZ5/branches/testing/libf/dyn3dpar/inidissip.F90

    r1665 r1707  
    33!
    44SUBROUTINE inidissip ( lstardis,nitergdiv,nitergrot,niterh  , &
    5      tetagdiv,tetagrot,tetatemp             )
     5     tetagdiv,tetagrot,tetatemp, vert_prof_dissip)
    66  !=======================================================================
    77  !   initialisation de la dissipation horizontale
     
    2525  INTEGER,INTENT(in) :: nitergdiv,nitergrot,niterh
    2626  REAL,INTENT(in) :: tetagdiv,tetagrot,tetatemp
     27
     28  integer, INTENT(in):: vert_prof_dissip
     29  ! Vertical profile of horizontal dissipation
     30  ! Allowed values:
     31  ! 0: rational fraction, function of pressure
     32  ! 1: tanh of altitude
    2733
    2834! Local variables:
     
    167173  !   --------------------------------------------------
    168174
    169   if (ok_strato .and. llm==39) then
     175  if (vert_prof_dissip == 1) then
    170176     do l=1,llm
    171177        pseudoz=8.*log(preff/presnivs(l))
  • LMDZ5/branches/testing/libf/dyn3dpar/leapfrog_p.F

    r1665 r1707  
    139139      REAL :: secondes
    140140
     141      logical :: physic
    141142      LOGICAL first,callinigrads
    142143
     
    208209
    209210      itau = 0
     211      physic=.true.
     212      if (iflag_phys==0.or.iflag_phys==2) physic=.false.
    210213!      iday = day_ini+itau/day_step
    211214!      time = REAL(itau-(iday-day_ini)*day_step)/day_step+time_0
     
    364367     s        apdiss = .TRUE.
    365368         IF( MOD(itau,iphysiq ).EQ.0.AND..NOT.forward
    366      s          .and. iflag_phys.EQ.1                 ) apphys = .TRUE.
     369     s          .and. physic                        ) apphys = .TRUE.
    367370      ELSE
    368371      ! Leapfrog/Matsuno time stepping
     
    370373         IF( MOD(itau+1,dissip_period).EQ.0 .AND. .NOT. forward )
    371374     s        apdiss = .TRUE.
    372          IF( MOD(itau+1,iphysiq).EQ.0.AND.iflag_phys.EQ.1) apphys=.TRUE.
     375         IF( MOD(itau+1,iphysiq).EQ.0.AND.physic) apphys=.TRUE.
    373376      END IF
    374377
     
    707710           jD_cur = jD_ref + day_ini - day_ref
    708711     $        + itau/day_step
     712
     713           IF (planet_type .eq."generic") THEN
     714              ! AS: we make jD_cur to be pday
     715              jD_cur = int(day_ini + itau/day_step)
     716           ENDIF
     717
    709718           jH_cur = jH_ref + start_time +                                &
    710719     &              mod(itau,day_step)/float(day_step)
  • LMDZ5/branches/testing/libf/dyn3dpar/parallel.F90

    r1664 r1707  
    489489          enddo
    490490         
    491         endif
     491        else
     492          ! Ehouarn: When in debug mode, ifort complains (for call MPI_GATHERV
     493          !          below) about Buffer_Recv() being not allocated.
     494          !          So make a dummy allocation.
     495          allocate(Buffer_Recv(1))
     496        endif ! of if (MPI_Rank==rank)
    492497 
    493498!$OMP CRITICAL (MPI)
  • LMDZ5/branches/testing/libf/dyn3dpar/paramet.h

    r792 r1707  
    1717      INTEGER jcfil,jcfllm
    1818
    19       PARAMETER( iip1= iim+1-1/iim,iip2=iim+2,iip3=iim+3                &
     19      PARAMETER( iip1= iim+1,iip2=iim+2,iip3=iim+3                      &
    2020     &    ,jjp1=jjm+1-1/jjm)
    2121      PARAMETER( llmp1 = llm+1,  llmp2 = llm+2, llmm1 = llm-1 )
  • LMDZ5/branches/testing/libf/filtrez/filtreg.F

    r1279 r1707  
    185185               DO j = jdfil,jffil
    186186#ifdef BLAS
    187                   CALL DGEMM("N", "N", iim, nbniv, iim, 1.0,
     187                  CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
    188188     &                 matrinvn(1,1,j),
    189189     &                 iim, champ(1,j,1), iip1*nlat, 0.0,
     
    199199               DO j = jdfil,jffil
    200200#ifdef BLAS
    201                   CALL DGEMM("N", "N", iim, nbniv, iim, 1.0,
     201                  CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
    202202     &                 matriceun(1,1,j),
    203203     &                 iim, champ(1,j,1), iip1*nlat, 0.0,
     
    213213               DO j = jdfil,jffil
    214214#ifdef BLAS
    215                   CALL DGEMM("N", "N", iim, nbniv, iim, 1.0,
     215                  CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
    216216     &                 matricevn(1,1,j),
    217217     &                 iim, champ(1,j,1), iip1*nlat, 0.0,
     
    231231               DO j = jdfil,jffil
    232232#ifdef BLAS
    233                   CALL DGEMM("N", "N", iim, nbniv, iim, 1.0,
     233                  CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
    234234     &                 matrinvs(1,1,j-jfiltsu+1),
    235235     &                 iim, champ(1,j,1), iip1*nlat, 0.0,
     
    247247               DO j = jdfil,jffil
    248248#ifdef BLAS
    249                   CALL DGEMM("N", "N", iim, nbniv, iim, 1.0,
     249                  CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
    250250     &                 matriceus(1,1,j-jfiltsu+1),
    251251     &                 iim, champ(1,j,1), iip1*nlat, 0.0,
     
    262262               DO j = jdfil,jffil
    263263#ifdef BLAS
    264                   CALL DGEMM("N", "N", iim, nbniv, iim, 1.0,
     264                  CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
    265265     &                 matricevs(1,1,j-jfiltsv+1),
    266266     &                 iim, champ(1,j,1), iip1*nlat, 0.0,
  • LMDZ5/branches/testing/libf/filtrez/filtreg_mod.F90

    r1665 r1707  
    11!
    2 ! $Id $
     2! $Id$
    33!
    44MODULE filtreg_mod
     
    1010
    1111  SUBROUTINE inifilr
    12   USE mod_filtre_fft
    13     !
     12  USE mod_filtre_fft, ONLY : use_filtre_fft,Init_filtre_fft
     13  USE mod_filtre_fft_loc, ONLY : Init_filtre_fft_loc=>Init_filtre_fft    !
    1414    !    ... H. Upadhyaya, O.Sharma   ...
    1515    !
     
    541541       CALL Init_filtre_fft(coefilu,modfrstu,jfiltnu,jfiltsu,  &
    542542                           coefilv,modfrstv,jfiltnv,jfiltsv)
     543       CALL Init_filtre_fft_loc(coefilu,modfrstu,jfiltnu,jfiltsu,  &
     544                           coefilv,modfrstv,jfiltnv,jfiltsv)
    543545    ENDIF
    544546
  • LMDZ5/branches/testing/libf/filtrez/mod_fft.F90

    r986 r1707  
    33#ifdef FFT_MATHKEISAN
    44  USE mod_fft_mathkeisan
    5 #elif FFT_FFTW
     5#else
     6#ifdef FFT_FFTW
    67  USE mod_fft_fftw
    7 #elif FFT_MKL
     8#else
     9#ifdef FFT_MKL
    810  USE mod_fft_mkl
    911#else
    1012  USE mod_fft_wrapper
    1113#endif
     14#endif
     15#endif
    1216
    1317END MODULE mod_fft
  • LMDZ5/branches/testing/libf/phy1d/1DUTILS.h_no_writelim

    r1665 r1707  
    5252!Config         (defaut sortie standard = 6)
    5353      lunout=6
    54       CALL getin('lunout', lunout)
     54!      CALL getin('lunout', lunout)
    5555      IF (lunout /= 5 .and. lunout /= 6) THEN
    5656        OPEN(lunout,FILE='lmdz.out')
  • LMDZ5/branches/testing/libf/phy1d/1DUTILS.h_with_writelim

    r1665 r1707  
    5252!Config         (defaut sortie standard = 6)
    5353      lunout=6
    54       CALL getin('lunout', lunout)
     54!      CALL getin('lunout', lunout)
    5555      IF (lunout /= 5 .and. lunout /= 6) THEN
    5656        OPEN(lunout,FILE='lmdz.out')
  • LMDZ5/branches/testing/libf/phy1d/1DUTILS.h_with_writelim_old

    r1665 r1707  
    5252!Config         (defaut sortie standard = 6)
    5353      lunout=6
    54       CALL getin('lunout', lunout)
     54!      CALL getin('lunout', lunout)
    5555      IF (lunout /= 5 .and. lunout /= 6) THEN
    5656        OPEN(lunout,FILE='lmdz.out')
  • LMDZ5/branches/testing/libf/phy1d/lmdz1d.F

    r1669 r1707  
    711711!---------------------------------------------------------------------
    712712
    713       fcoriolis=2.*sin(rpi*rlat(1)/180.)*romega
     713      fcoriolis=2.*sin(rpi*xlat/180.)*romega
    714714
    715715       if (forcing_radconv) then
  • LMDZ5/branches/testing/libf/phydev/iniphysiq.F

    r1665 r1707  
    22! $Id: iniphysiq.F 1403 2010-07-01 09:02:53Z fairhead $
    33!
    4 c
    5 c
    64      SUBROUTINE iniphysiq(ngrid,nlayer,
    75     $           punjours,
    86     $           pdayref,ptimestep,
    97     $           plat,plon,parea,pcu,pcv,
    10      $           prad,pg,pr,pcpp)
    11       USE dimphy
    12       USE mod_grid_phy_lmdz
    13       USE mod_phys_lmdz_para
    14       USE comgeomphy
     8     $           prad,pg,pr,pcpp,iflag_phys)
     9      USE dimphy, only : klev
     10      USE mod_grid_phy_lmdz, only : klon_glo
     11      USE mod_phys_lmdz_para, only : klon_omp,klon_omp_begin,
     12     &                               klon_omp_end,klon_mpi_begin
     13      USE comgeomphy, only : airephy,cuphy,cvphy,rlond,rlatd
     14      USE comcstphy, only : rradius,rg,rr,rcpp
    1515
    1616      IMPLICIT NONE
     
    1818c=======================================================================
    1919c
    20 c   subject:
    21 c   --------
     20c   Initialisation of the physical constants and some positional and
     21c   geometrical arrays for the physics
    2222c
    23 c   Initialisation for the physical parametrisations of the LMD
    24 c   martian atmospheric general circulation modele.
    25 c
    26 c   author: Frederic Hourdin 15 / 10 /93
    27 c   -------
    28 c
    29 c   arguments:
    30 c   ----------
    31 c
    32 c   input:
    33 c   ------
    3423c
    3524c    ngrid                 Size of the horizontal grid.
     
    3726c    nlayer                Number of vertical layers.
    3827c    pdayref               Day of reference for the simulation
    39 c    firstcall             True at the first call
    40 c    lastcall              True at the last call
    41 c    pday                  Number of days counted from the North. Spring
    42 c                          equinoxe.
    4328c
    4429c=======================================================================
    45 c
    46 c-----------------------------------------------------------------------
    47 c   declarations:
    48 c   -------------
     30 
    4931 
    5032cym#include "dimensions.h"
    5133cym#include "dimphy.h"
    5234cym#include "comgeomphy.h"
    53 #include "comcstphy.h"
    54       REAL prad,pg,pr,pcpp,punjours
    55  
    56       INTEGER ngrid,nlayer
    57       REAL plat(ngrid),plon(ngrid),parea(klon_glo)
    58       REAL pcu(klon_glo),pcv(klon_glo)
    59       INTEGER pdayref
     35#include "iniprint.h"
     36
     37      REAL,INTENT(IN) :: prad ! radius of the planet (m)
     38      REAL,INTENT(IN) :: pg ! gravitational acceleration (m/s2)
     39      REAL,INTENT(IN) :: pr ! ! reduced gas constant R/mu
     40      REAL,INTENT(IN) :: pcpp ! specific heat Cp
     41      REAL,INTENT(IN) :: punjours ! length (in s) of a standard day
     42      INTEGER,INTENT(IN) :: ngrid ! number of horizontal grid points in the physics
     43      INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers
     44      REAL,INTENT(IN) :: plat(ngrid) ! latitudes of the physics grid
     45      REAL,INTENT(IN) :: plon(ngrid) ! longitudes of the physics grid
     46      REAL,INTENT(IN) :: parea(klon_glo) ! area (m2)
     47      REAL,INTENT(IN) :: pcu(klon_glo) ! cu coeff. (u_covariant = cu * u)
     48      REAL,INTENT(IN) :: pcv(klon_glo) ! cv coeff. (v_covariant = cv * v)
     49      INTEGER,INTENT(IN) :: pdayref ! reference day of for the simulation
     50      REAL,INTENT(IN) :: ptimestep !physics time step (s)
     51      INTEGER,INTENT(IN) :: iflag_phys ! type of physics to be called
     52
    6053      INTEGER :: ibegin,iend,offset
    61  
    62       REAL ptimestep
    6354      CHARACTER (LEN=20) :: modname='iniphysiq'
    6455      CHARACTER (LEN=80) :: abort_message
    6556 
    6657      IF (nlayer.NE.klev) THEN
    67          PRINT*,'STOP in inifis'
    68          PRINT*,'Probleme de dimensions :'
    69          PRINT*,'nlayer     = ',nlayer
    70          PRINT*,'klev   = ',klev
     58         write(lunout,*) 'STOP in ',trim(modname)
     59         write(lunout,*) 'Problem with dimensions :'
     60         write(lunout,*) 'nlayer     = ',nlayer
     61         write(lunout,*) 'klev   = ',klev
    7162         abort_message = ''
    7263         CALL abort_gcm (modname,abort_message,1)
     
    7465
    7566      IF (ngrid.NE.klon_glo) THEN
    76          PRINT*,'STOP in inifis'
    77          PRINT*,'Probleme de dimensions :'
    78          PRINT*,'ngrid     = ',ngrid
    79          PRINT*,'klon   = ',klon_glo
     67         write(lunout,*) 'STOP in ',trim(modname)
     68         write(lunout,*) 'Problem with dimensions :'
     69         write(lunout,*) 'ngrid     = ',ngrid
     70         write(lunout,*) 'klon   = ',klon_glo
    8071         abort_message = ''
    8172         CALL abort_gcm (modname,abort_message,1)
    8273      ENDIF
    83 c$OMP PARALLEL PRIVATE(ibegin,iend)
    84 c$OMP+         SHARED(parea,pcu,pcv,plon,plat)
     74
     75!$OMP PARALLEL PRIVATE(ibegin,iend)
     76!$OMP+         SHARED(parea,pcu,pcv,plon,plat)
    8577     
    8678      offset=klon_mpi_begin-1
     
    9284      rlatd(1:klon_omp)=plat(offset+klon_omp_begin:offset+klon_omp_end)
    9385
    94 !     call suphel
    95 !     prad,pg,pr,pcpp
     86! copy some fundamental parameters to physics
    9687      rradius=prad
    9788      rg=pg
     
    9990      rcpp=pcpp
    10091
    101      
     92!$OMP END PARALLEL
    10293
    103 c$OMP END PARALLEL
     94!      print*,'ATTENTION !!! TRAVAILLER SUR INIPHYSIQ'
     95!      print*,'CONTROLE DES LATITUDES, LONGITUDES, PARAMETRES ...'
    10496
    105       print*,'ATTENTION !!! TRAVAILLER SUR INIPHYSIQ'
    106       print*,'CONTROLE DES LATITUDES, LONGITUDES, PARAMETRES ...'
     97! Additional initializations for aquaplanets
     98!$OMP PARALLEL
     99      if (iflag_phys>=100) then
     100        call iniaqua(klon_omp,rlatd,rlond,iflag_phys)
     101      endif
     102!$OMP END PARALLEL
    107103
    108       RETURN
    109 9999  CONTINUE
    110       abort_message ='Cette version demande les fichier rnatur.dat
    111      & et surf.def'
    112       CALL abort_gcm (modname,abort_message,1)
     104!      RETURN
     105!9999  CONTINUE
     106!      abort_message ='Cette version demande les fichier rnatur.dat
     107!     & et surf.def'
     108!      CALL abort_gcm (modname,abort_message,1)
    113109
    114110      END
  • LMDZ5/branches/testing/libf/phydev/phyaqua.F

    r1665 r1707  
    1 ! Routines complementaires pour la physique planetaire.
    2 
     1!
     2! $Id: $
     3!
    34
    45      subroutine iniaqua(nlon,latfi,lonfi,iflag_phys)
    56
    67!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    7 !  Creation d'un etat initial et de conditions aux limites
    8 !  (resp startphy.nc et limit.nc) pour des configurations idealisees
    9 ! du modele LMDZ dans sa version terrestre.
    10 !  iflag_phys est un parametre qui controle
    11 !  iflag_phys = N 
    12 !    de 100 a 199 : aqua planetes avec SST forcees
    13 !                 N-100 determine le type de SSTs
    14 !    de 200 a 299 : terra planetes avec Ts calcule
    15 !       
     8!  Create an initial state (startphy.nc) for the physics
     9!  Usefull for idealised cases (e.g. aquaplanets or testcases)
    1610!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    1711
     12      use phys_state_var_mod, only : rlat,rlon,
     13     &                               phys_state_var_init
     14      use mod_phys_lmdz_para, only : klon_omp
     15      use comgeomphy, only : rlond,rlatd
     16      implicit none
     17     
     18      integer,intent(in) :: nlon,iflag_phys
     19      real,intent(in) :: lonfi(nlon),latfi(nlon)
    1820
    19       integer nlon,iflag_phys
    20 cIM ajout latfi, lonfi
    21       REAL, DIMENSION (nlon) :: lonfi, latfi
     21! local variables
     22      real :: pi
     23
     24! initializations:
     25      pi=2.*asin(1.)
     26
     27      call phys_state_var_init()
     28
     29      rlat(1:klon_omp)=rlatd(1:klon_omp)*180./pi
     30      rlon(1:klon_omp)=rlond(1:klon_omp)*180./pi
    2231
    2332
    24       return
     33! Here you could create an initial condition for the physics
     34! ...
     35! ... fill in the fields...
     36! ...
     37! ... and create a "startphy.nc" file
     38      CALL phyredem ("startphy.nc")
     39
    2540      end
    2641
  • LMDZ5/branches/testing/libf/phydev/physiq.F90

    r1665 r1707  
    1111     &            , PVteta)
    1212
    13       USE dimphy
    14       USE infotrac
    15       USE comgeomphy
     13      USE dimphy, only : klon,klev
     14      USE infotrac, only : nqtot
     15      USE comgeomphy, only : rlatd
     16      USE comcstphy, only : rg
     17      USE iophy, only : histbeg_phy,histwrite_phy
     18      USE ioipsl, only : getin,histvert,histdef,histend,ymds2ju
     19      USE mod_phys_lmdz_para, only : jj_nb
     20      USE phys_state_var_mod, only : phys_state_var_init
    1621
    1722      IMPLICIT none
    18 !======================================================================
    19 ! Objet: Moniteur general de la physique du modele
    20 !======================================================================
     23#include "dimensions.h"
     24
     25      integer,parameter :: jjmp1=jjm+1-1/jjm
     26      integer,parameter :: iip1=iim+1
    2127!
    22 !  Arguments:
     28! Routine argument:
    2329!
    24 ! nlon----input-I-nombre de points horizontaux
    25 ! nlev----input-I-nombre de couches verticales, doit etre egale a klev
    26 ! debut---input-L-variable logique indiquant le premier passage
    27 ! lafin---input-L-variable logique indiquant le dernier passage
    28 ! jD_cur       -R-jour courant a l'appel de la physique (jour julien)
    29 ! jH_cur       -R-heure courante a l'appel de la physique (jour julien)
    30 ! pdtphys-input-R-pas d'integration pour la physique (seconde)
    31 ! paprs---input-R-pression pour chaque inter-couche (en Pa)
    32 ! pplay---input-R-pression pour le mileu de chaque couche (en Pa)
    33 ! pphi----input-R-geopotentiel de chaque couche (g z) (reference sol)
    34 ! pphis---input-R-geopotentiel du sol
    35 ! presnivs-input_R_pressions approximat. des milieux couches ( en PA)
    36 ! u-------input-R-vitesse dans la direction X (de O a E) en m/s
    37 ! v-------input-R-vitesse Y (de S a N) en m/s
    38 ! t-------input-R-temperature (K)
    39 ! qx------input-R-humidite specifique (kg/kg) et d'autres traceurs
    40 ! d_t_dyn-input-R-tendance dynamique pour "t" (K/s)
    41 ! d_q_dyn-input-R-tendance dynamique pour "q" (kg/kg/s)
    42 ! flxmass_w -input-R- flux de masse verticale
    43 ! d_u-----output-R-tendance physique de "u" (m/s/s)
    44 ! d_v-----output-R-tendance physique de "v" (m/s/s)
    45 ! d_t-----output-R-tendance physique de "t" (K/s)
    46 ! d_qx----output-R-tendance physique de "qx" (kg/kg/s)
    47 ! d_ps----output-R-tendance physique de la pression au sol
    48 !IM
    49 ! PVteta--output-R-vorticite potentielle a des thetas constantes
    50 !======================================================================
    51 #include "dimensions.h"
    52 #include "comcstphy.h"
     30      integer,intent(in) :: nlon ! number of atmospheric colums
     31      integer,intent(in) :: nlev ! number of vertical levels (should be =klev)
     32      real,intent(in) :: jD_cur ! current day number (Julian day)
     33      real,intent(in) :: jH_cur ! current time of day (as fraction of day)
     34      logical,intent(in) :: debut ! signals first call to physics
     35      logical,intent(in) :: lafin ! signals last call to physics
     36      real,intent(in) :: pdtphys ! physics time step (s)
     37      real,intent(in) :: paprs(klon,klev+1) ! interlayer pressure (Pa)
     38      real,intent(in) :: pplay(klon,klev) ! mid-layer pressure (Pa)
     39      real,intent(in) :: pphi(klon,klev) ! geopotential at mid-layer
     40      real,intent(in) :: pphis(klon) ! surface geopotential
     41      real,intent(in) :: presnivs(klev) ! pseudo-pressure (Pa) of mid-layers
     42      integer,parameter :: longcles=20
     43      real,intent(in) :: clesphy0(longcles) ! Not used
     44      real,intent(in) :: u(klon,klev) ! eastward zonal wind (m/s)
     45      real,intent(in) :: v(klon,klev) ! northward meridional wind (m/s)
     46      real,intent(in) :: t(klon,klev) ! temperature (K)
     47      real,intent(in) :: qx(klon,klev,nqtot) ! tracers (.../kg_air)
     48      real,intent(in) :: flxmass_w(klon,klev) ! vertical mass flux
     49      real,intent(out) :: d_u(klon,klev) ! physics tendency on u (m/s/s)
     50      real,intent(out) :: d_v(klon,klev) ! physics tendency on v (m/s/s)
     51      real,intent(out) :: d_t(klon,klev) ! physics tendency on t (K/s)
     52      real,intent(out) :: d_qx(klon,klev,nqtot) ! physics tendency on tracers
     53      real,intent(out) :: d_ps(klon) ! physics tendency on surface pressure
     54      real,intent(in) :: dudyn(iim+1,jjmp1,klev) ! Not used
     55!FH! REAL PVteta(klon,nbteta)
     56!      REAL PVteta(klon,1)
     57      real,intent(in) :: PVteta(klon,3) ! Not used ; should match definition
     58                                        ! in calfis.F
    5359
    54       integer jjmp1
    55       parameter (jjmp1=jjm+1-1/jjm)
    56       integer iip1
    57       parameter (iip1=iim+1)
     60integer,save :: itau=0 ! counter to count number of calls to physics
     61!$OMP THREADPRIVATE(itau)
     62real :: temp_newton(klon,klev)
     63integer :: k
     64logical, save :: first=.true.
     65!$OMP THREADPRIVATE(first)
    5866
    59       INTEGER ivap          ! indice de traceurs pour vapeur d'eau
    60       PARAMETER (ivap=1)
    61       INTEGER iliq          ! indice de traceurs pour eau liquide
    62       PARAMETER (iliq=2)
     67! For I/Os
     68integer :: itau0
     69real :: zjulian
     70real :: dtime
     71integer :: nhori ! horizontal coordinate ID
     72integer,save :: nid_hist ! output file ID
     73!$OMP THREADPRIVATE(nid_hist)
     74integer :: zvertid ! vertical coordinate ID
     75integer,save :: iwrite_phys ! output every iwrite_phys physics step
     76!$OMP THREADPRIVATE(iwrite_phys)
     77real :: t_ops ! frequency of the IOIPSL operations (eg average over...)
     78real :: t_wrt ! frequency of the IOIPSL outputs
    6379
    64 !
    65 !
    66 ! Variables argument:
    67 !
    68       INTEGER nlon
    69       INTEGER nlev
    70       REAL, intent(in):: jD_cur, jH_cur
     80! initializations
     81if (debut) then ! Things to do only for the first call to physics
     82! load initial conditions for physics (including the grid)
     83  call phys_state_var_init() ! some initializations, required before calling phyetat0
     84  call phyetat0("startphy.nc")
    7185
    72       REAL pdtphys
    73       LOGICAL debut, lafin
    74       REAL paprs(klon,klev+1)
    75       REAL pplay(klon,klev)
    76       REAL pphi(klon,klev)
    77       REAL pphis(klon)
    78       REAL presnivs(klev)
    79       REAL znivsig(klev)
    80       real pir
     86! Initialize outputs:
     87  itau0=0
     88  iwrite_phys=1 !default: output every physics timestep
     89  call getin("iwrite_phys",iwrite_phys)
     90  t_ops=pdtphys*iwrite_phys ! frequency of the IOIPSL operation
     91  t_wrt=pdtphys*iwrite_phys ! frequency of the outputs in the file
     92  ! compute zjulian for annee0=1979 and month=1 dayref=1 and hour=0.0
     93  !CALL ymds2ju(annee0, month, dayref, hour, zjulian)
     94  call ymds2ju(1979, 1, 1, 0.0, zjulian)
     95  dtime=pdtphys
     96  call histbeg_phy("histins.nc",itau0,zjulian,dtime,nhori,nid_hist)
     97!$OMP MASTER
     98  ! define vertical coordinate
     99  call histvert(nid_hist,"presnivs","Vertical levels","Pa",klev, &
     100                presnivs,zvertid,'down')
     101  ! define variables which will be written in "histins.nc" file
     102  call histdef(nid_hist,'temperature','Atmospheric temperature','K', &
     103               iim,jj_nb,nhori,klev,1,klev,zvertid,32, &
     104               'inst(X)',t_ops,t_wrt)
     105  call histdef(nid_hist,'u','Eastward Zonal Wind','m/s', &
     106               iim,jj_nb,nhori,klev,1,klev,zvertid,32, &
     107               'inst(X)',t_ops,t_wrt)
     108  call histdef(nid_hist,'v','Northward Meridional Wind','m/s', &
     109               iim,jj_nb,nhori,klev,1,klev,zvertid,32, &
     110               'inst(X)',t_ops,t_wrt)
     111  call histdef(nid_hist,'ps','Surface Pressure','Pa', &
     112               iim,jj_nb,nhori,1,1,1,zvertid,32, &
     113               'inst(X)',t_ops,t_wrt)
     114  ! end definition sequence
     115  call histend(nid_hist)
     116!$OMP END MASTER
     117endif ! of if (debut)
    81118
    82       REAL u(klon,klev)
    83       REAL v(klon,klev)
    84       REAL t(klon,klev),theta(klon,klev)
    85       REAL qx(klon,klev,nqtot)
    86       REAL flxmass_w(klon,klev)
    87       REAL omega(klon,klev) ! vitesse verticale en Pa/s
    88       REAL d_u(klon,klev)
    89       REAL d_v(klon,klev)
    90       REAL d_t(klon,klev)
    91       REAL d_qx(klon,klev,nqtot)
    92       REAL d_ps(klon)
    93       real da(klon,klev),phi(klon,klev,klev),mp(klon,klev)
    94 !IM definition dynamique o_trac dans phys_output_open
    95 !      type(ctrl_out) :: o_trac(nqtot)
    96 !FH! REAL PVteta(klon,nbteta)
    97       REAL PVteta(klon,1)
    98       REAL dudyn(iim+1,jjmp1,klev)
     119! increment counter itau
     120itau=itau+1
    99121
    100     INTEGER        longcles
    101     PARAMETER    ( longcles = 20 )
     122! set all tendencies to zero
     123d_u(1:klon,1:klev)=0.
     124d_v(1:klon,1:klev)=0.
     125d_t(1:klon,1:klev)=0.
     126d_qx(1:klon,1:klev,1:nqtot)=0.
     127d_ps(1:klon)=0.
    102128
    103 real temp_newton(klon,klev)
    104 integer k
    105 logical, save :: first=.true.
    106 
    107       REAL clesphy0( longcles      )
    108 
    109 d_u=0.
    110 d_v=0.
    111 d_t=0.
    112 d_qx=0.
    113 d_ps=0.
    114 
    115      d_u(:,1)=-u(:,1)/86400.
    116      do k=1,klev
    117         temp_newton(:,k)=280.+cos(rlatd(:))*40.-pphi(:,k)/rg*6.e-3
    118         d_t(:,k)=(temp_newton(:,k)-t(:,k))/1.e5
    119      enddo
     129! compute tendencies to return to the dynamics:
     130! "friction" on the first layer
     131d_u(1:klon,1)=-u(1:klon,1)/86400.
     132d_v(1:klon,1)=-v(1:klon,1)/86400.
     133! newtonian relaxation towards temp_newton()
     134do k=1,klev
     135  temp_newton(1:klon,k)=280.+cos(rlatd(1:klon))*40.-pphi(1:klon,k)/rg*6.e-3
     136  d_t(1:klon,k)=(temp_newton(1:klon,k)-t(1:klon,k))/1.e5
     137enddo
    120138
    121139
    122       print*,'COUCOU PHYDEV'
    123       return
    124       end
     140!print*,'PHYDEV: itau=',itau
     141
     142! write some outputs:
     143if (modulo(itau,iwrite_phys)==0) then
     144  call histwrite_phy(nid_hist,.false.,"temperature",itau,t)
     145  call histwrite_phy(nid_hist,.false.,"u",itau,u)
     146  call histwrite_phy(nid_hist,.false.,"v",itau,v)
     147  call histwrite_phy(nid_hist,.false.,"ps",itau,paprs(:,1))
     148endif
     149
     150! if lastcall, then it is time to write "restartphy.nc" file
     151if (lafin) then
     152  call phyredem("restartphy.nc")
     153endif
     154
     155end
  • LMDZ5/branches/testing/libf/phylmd/calcul_STDlev.h

    r1418 r1707  
    22c $Header$
    33c
    4 c
    54cIM on initialise les variables
     5c
     6        missing_val=nf90_fill_real
     7c
     8cIM freq_moyNMC = frequences auxquelles on moyenne les champs accumules
     9cIM               sur les niveaux de pression standard du NMC
     10      DO n=1, nout
     11       freq_moyNMC(n)=freq_outNMC(n)/freq_calNMC(n)
     12      ENDDO
    613c
    714        CALL ini_undefSTD(itap,freq_outNMC)
     
    157164     $     lwup,LWup200)
    158165c
     166      twriteSTD(:,:,1)=tsumSTD(:,:,1)
     167      qwriteSTD(:,:,1)=qsumSTD(:,:,1)
     168      rhwriteSTD(:,:,1)=rhsumSTD(:,:,1)
     169      phiwriteSTD(:,:,1)=phisumSTD(:,:,1)
     170      uwriteSTD(:,:,1)=usumSTD(:,:,1)
     171      vwriteSTD(:,:,1)=vsumSTD(:,:,1)
     172      wwriteSTD(:,:,1)=wsumSTD(:,:,1)
     173
     174      twriteSTD(:,:,2)=tsumSTD(:,:,2)
     175      qwriteSTD(:,:,2)=qsumSTD(:,:,2)
     176      rhwriteSTD(:,:,2)=rhsumSTD(:,:,2)
     177      phiwriteSTD(:,:,2)=phisumSTD(:,:,2)
     178      uwriteSTD(:,:,2)=usumSTD(:,:,2)
     179      vwriteSTD(:,:,2)=vsumSTD(:,:,2)
     180      wwriteSTD(:,:,2)=wsumSTD(:,:,2)
     181
     182      twriteSTD(:,:,3)=tlevSTD(:,:)
     183      qwriteSTD(:,:,3)=qlevSTD(:,:)
     184      rhwriteSTD(:,:,3)=rhlevSTD(:,:)
     185      phiwriteSTD(:,:,3)=philevSTD(:,:)
     186      uwriteSTD(:,:,3)=ulevSTD(:,:)
     187      vwriteSTD(:,:,3)=vlevSTD(:,:)
     188      wwriteSTD(:,:,3)=wlevSTD(:,:)
     189
     190      twriteSTD(:,:,4)=tlevSTD(:,:)
     191      qwriteSTD(:,:,4)=qlevSTD(:,:)
     192      rhwriteSTD(:,:,4)=rhlevSTD(:,:)
     193      phiwriteSTD(:,:,4)=philevSTD(:,:)
     194      uwriteSTD(:,:,4)=ulevSTD(:,:)
     195      vwriteSTD(:,:,4)=vlevSTD(:,:)
     196      wwriteSTD(:,:,4)=wlevSTD(:,:)
     197c
     198cIM initialisation 5eme fichier de sortie
     199      twriteSTD(:,:,5)=tlevSTD(:,:)
     200      qwriteSTD(:,:,5)=qlevSTD(:,:)
     201      rhwriteSTD(:,:,5)=rhlevSTD(:,:)
     202      phiwriteSTD(:,:,5)=philevSTD(:,:)
     203      uwriteSTD(:,:,5)=ulevSTD(:,:)
     204      vwriteSTD(:,:,5)=vlevSTD(:,:)
     205      wwriteSTD(:,:,5)=wlevSTD(:,:)
     206c
     207cIM initialisation 6eme fichier de sortie
     208      twriteSTD(:,:,6)=tlevSTD(:,:)
     209      qwriteSTD(:,:,6)=qlevSTD(:,:)
     210      rhwriteSTD(:,:,6)=rhlevSTD(:,:)
     211      phiwriteSTD(:,:,6)=philevSTD(:,:)
     212      uwriteSTD(:,:,6)=ulevSTD(:,:)
     213      vwriteSTD(:,:,6)=vlevSTD(:,:)
     214      wwriteSTD(:,:,6)=wlevSTD(:,:)
     215cIM for NMC files
     216      DO n=1, nlevSTD3
     217       DO k=1, nlevSTD
     218        if(rlevSTD3(n).EQ.rlevSTD(k)) THEN
     219         twriteSTD3(:,n)=tlevSTD(:,k)
     220         qwriteSTD3(:,n)=qlevSTD(:,k)
     221         rhwriteSTD3(:,n)=rhlevSTD(:,k)
     222         phiwriteSTD3(:,n)=philevSTD(:,k)
     223         uwriteSTD3(:,n)=ulevSTD(:,k)
     224         vwriteSTD3(:,n)=vlevSTD(:,k)
     225         wwriteSTD3(:,n)=wlevSTD(:,k)
     226        endif !rlevSTD3(n).EQ.rlevSTD(k)
     227       ENDDO
     228      ENDDO
     229c
     230      DO n=1, nlevSTD8
     231       DO k=1, nlevSTD
     232        if(rlevSTD8(n).EQ.rlevSTD(k)) THEN
     233         tnondefSTD8(:,n)=tnondef(:,k,2)
     234         twriteSTD8(:,n)=tsumSTD(:,k,2)
     235         qwriteSTD8(:,n)=qsumSTD(:,k,2)
     236         rhwriteSTD8(:,n)=rhsumSTD(:,k,2)
     237         phiwriteSTD8(:,n)=phisumSTD(:,k,2)
     238         uwriteSTD8(:,n)=usumSTD(:,k,2)
     239         vwriteSTD8(:,n)=vsumSTD(:,k,2)
     240         wwriteSTD8(:,n)=wsumSTD(:,k,2)
     241        endif !rlevSTD8(n).EQ.rlevSTD(k)
     242       ENDDO
     243      ENDDO
  • LMDZ5/branches/testing/libf/phylmd/change_srf_frac_mod.F90

    r1454 r1707  
    1212
    1313  SUBROUTINE change_srf_frac(itime, dtime, jour, &
    14        pctsrf, alb1, alb2, tsurf, u10m, v10m, pbl_tke)
     14       pctsrf, alb1, alb2, tsurf, ustar, u10m, v10m, pbl_tke)
    1515!
    1616! This subroutine is called from physiq.F at each timestep.
     
    4646    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: alb2   ! albedo second interval in SW spektrum
    4747    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: tsurf
     48    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: ustar
    4849    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: u10m
    4950    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: v10m
     
    150151!
    151152!****************************************************************************************
    152        CALL pbl_surface_newfrac(itime, pctsrf, pctsrf_old, tsurf, alb1, alb2, u10m, v10m, pbl_tke)
     153       CALL pbl_surface_newfrac(itime, pctsrf, pctsrf_old, tsurf, alb1, alb2, ustar, u10m, v10m, pbl_tke)
    153154
    154155    ELSE
  • LMDZ5/branches/testing/libf/phylmd/iniphysiq.F

    r1403 r1707  
    88     $           pdayref,ptimestep,
    99     $           plat,plon,parea,pcu,pcv,
    10      $           prad,pg,pr,pcpp)
    11       USE dimphy
    12       USE mod_grid_phy_lmdz
    13       USE mod_phys_lmdz_para
    14       USE comgeomphy
     10     $           prad,pg,pr,pcpp,iflag_phys)
     11      USE dimphy, only : klev
     12      USE mod_grid_phy_lmdz, only : klon_glo
     13      USE mod_phys_lmdz_para, only : klon_omp,klon_omp_begin,
     14     &                               klon_omp_end,klon_mpi_begin
     15      USE comgeomphy, only : airephy,cuphy,cvphy,rlond,rlatd
    1516
    1617      IMPLICIT NONE
     
    1819c=======================================================================
    1920c
    20 c   subject:
    21 c   --------
     21c   Initialisation of the physical constants and some positional and
     22c   geometrical arrays for the physics
    2223c
    23 c   Initialisation for the physical parametrisations of the LMD
    24 c   martian atmospheric general circulation modele.
    25 c
    26 c   author: Frederic Hourdin 15 / 10 /93
    27 c   -------
    28 c
    29 c   arguments:
    30 c   ----------
    31 c
    32 c   input:
    33 c   ------
    3424c
    3525c    ngrid                 Size of the horizontal grid.
     
    3727c    nlayer                Number of vertical layers.
    3828c    pdayref               Day of reference for the simulation
    39 c    firstcall             True at the first call
    40 c    lastcall              True at the last call
    41 c    pday                  Number of days counted from the North. Spring
    42 c                          equinoxe.
    4329c
    4430c=======================================================================
    45 c
    46 c-----------------------------------------------------------------------
    47 c   declarations:
    48 c   -------------
    4931 
    5032cym#include "dimensions.h"
     
    5234cym#include "comgeomphy.h"
    5335#include "YOMCST.h"
    54       REAL prad,pg,pr,pcpp,punjours
    55  
    56       INTEGER ngrid,nlayer
    57       REAL plat(ngrid),plon(ngrid),parea(klon_glo)
    58       REAL pcu(klon_glo),pcv(klon_glo)
    59       INTEGER pdayref
    60       INTEGER :: ibegin,iend,offset
    61  
    62       REAL ptimestep
     36#include "iniprint.h"
     37
     38      REAL,INTENT(IN) :: prad ! radius of the planet (m)
     39      REAL,INTENT(IN) :: pg ! gravitational acceleration (m/s2)
     40      REAL,INTENT(IN) :: pr ! ! reduced gas constant R/mu
     41      REAL,INTENT(IN) :: pcpp ! specific heat Cp
     42      REAL,INTENT(IN) :: punjours ! length (in s) of a standard day
     43      INTEGER,INTENT(IN) :: ngrid ! number of horizontal grid points in the physics
     44      INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers
     45      REAL,INTENT(IN) :: plat(ngrid) ! latitudes of the physics grid
     46      REAL,INTENT(IN) :: plon(ngrid) ! longitudes of the physics grid
     47      REAL,INTENT(IN) :: parea(klon_glo) ! area (m2)
     48      REAL,INTENT(IN) :: pcu(klon_glo) ! cu coeff. (u_covariant = cu * u)
     49      REAL,INTENT(IN) :: pcv(klon_glo) ! cv coeff. (v_covariant = cv * v)
     50      INTEGER,INTENT(IN) :: pdayref ! reference day of for the simulation
     51      REAL,INTENT(IN) :: ptimestep !physics time step (s)
     52      INTEGER,INTENT(IN) :: iflag_phys ! type of physics to be called
     53
     54      INTEGER :: ibegin,iend,offset
    6355      CHARACTER (LEN=20) :: modname='iniphysiq'
    6456      CHARACTER (LEN=80) :: abort_message
    6557 
    6658      IF (nlayer.NE.klev) THEN
    67          PRINT*,'STOP in inifis'
    68          PRINT*,'Probleme de dimensions :'
    69          PRINT*,'nlayer     = ',nlayer
    70          PRINT*,'klev   = ',klev
     59         write(lunout,*) 'STOP in ',trim(modname)
     60         write(lunout,*) 'Problem with dimensions :'
     61         write(lunout,*) 'nlayer     = ',nlayer
     62         write(lunout,*) 'klev   = ',klev
    7163         abort_message = ''
    7264         CALL abort_gcm (modname,abort_message,1)
     
    7466
    7567      IF (ngrid.NE.klon_glo) THEN
    76          PRINT*,'STOP in inifis'
    77          PRINT*,'Probleme de dimensions :'
    78          PRINT*,'ngrid     = ',ngrid
    79          PRINT*,'klon   = ',klon_glo
     68         write(lunout,*) 'STOP in ',trim(modname)
     69         write(lunout,*) 'Problem with dimensions :'
     70         write(lunout,*) 'ngrid     = ',ngrid
     71         write(lunout,*) 'klon   = ',klon_glo
    8072         abort_message = ''
    8173         CALL abort_gcm (modname,abort_message,1)
    8274      ENDIF
    83 c$OMP PARALLEL PRIVATE(ibegin,iend)
    84 c$OMP+         SHARED(parea,pcu,pcv,plon,plat)
     75
     76!$OMP PARALLEL PRIVATE(ibegin,iend)
     77!$OMP+         SHARED(parea,pcu,pcv,plon,plat)
    8578     
    8679      offset=klon_mpi_begin-1
     
    9285      rlatd(1:klon_omp)=plat(offset+klon_omp_begin:offset+klon_omp_end)
    9386
     87      ! suphel => initialize some physical constants (orbital parameters,
     88      !           geoid, gravity, thermodynamical constants, etc.) in the
     89      !           physics
    9490      call suphel
     91     
     92!$OMP END PARALLEL
    9593
    96 c$OMP END PARALLEL
     94      ! check that physical constants set in 'suphel' are coherent
     95      ! with values set in the dynamics:
     96      if (RDAY.ne.punjours) then
     97        write(lunout,*) "iniphysiq: length of day discrepancy!!!"
     98        write(lunout,*) "  in the dynamics punjours=",punjours
     99        write(lunout,*) "   but in the physics RDAY=",RDAY
     100        if (abs(RDAY-punjours).gt.0.01) then
     101          ! stop here if the relative difference is more than 1%
     102          abort_message = 'length of day discrepancy'
     103          CALL abort_gcm (modname,abort_message,1)
     104        endif
     105      endif
     106      if (RG.ne.pg) then
     107        write(lunout,*) "iniphysiq: gravity discrepancy !!!"
     108        write(lunout,*) "     in the dynamics pg=",pg
     109        write(lunout,*) "  but in the physics RG=",RG
     110        if (abs(RG-pg).gt.0.01) then
     111          ! stop here if the relative difference is more than 1%
     112          abort_message = 'gravity discrepancy'
     113          CALL abort_gcm (modname,abort_message,1)
     114        endif
     115      endif
     116      if (RA.ne.prad) then
     117        write(lunout,*) "iniphysiq: planet radius discrepancy !!!"
     118        write(lunout,*) "   in the dynamics prad=",prad
     119        write(lunout,*) "  but in the physics RA=",RA
     120        if (abs(RA-prad).gt.0.01) then
     121          ! stop here if the relative difference is more than 1%
     122          abort_message = 'planet radius discrepancy'
     123          CALL abort_gcm (modname,abort_message,1)
     124        endif
     125      endif
     126      if (RD.ne.pr) then
     127        write(lunout,*)"iniphysiq: reduced gas constant discrepancy !!!"
     128        write(lunout,*)"     in the dynamics pr=",pr
     129        write(lunout,*)"  but in the physics RD=",RD
     130        if (abs(RD-pr).gt.0.01) then
     131          ! stop here if the relative difference is more than 1%
     132          abort_message = 'reduced gas constant discrepancy'
     133          CALL abort_gcm (modname,abort_message,1)
     134        endif
     135      endif
     136      if (RCPD.ne.pcpp) then
     137        write(lunout,*)"iniphysiq: specific heat discrepancy !!!"
     138        write(lunout,*)"     in the dynamics pcpp=",pcpp
     139        write(lunout,*)"  but in the physics RCPD=",RCPD
     140        if (abs(RCPD-pcpp).gt.0.01) then
     141          ! stop here if the relative difference is more than 1%
     142          abort_message = 'specific heat discrepancy'
     143          CALL abort_gcm (modname,abort_message,1)
     144        endif
     145      endif
    97146
    98       print*,'ATTENTION !!! TRAVAILLER SUR INIPHYSIQ'
    99       print*,'CONTROLE DES LATITUDES, LONGITUDES, PARAMETRES ...'
     147! Additional initializations for aquaplanets
     148!$OMP PARALLEL
     149      if (iflag_phys>=100) then
     150        call iniaqua(klon_omp,rlatd,rlond,iflag_phys)
     151      endif
     152!$OMP END PARALLEL
    100153
    101       RETURN
    102 9999  CONTINUE
    103       abort_message ='Cette version demande les fichier rnatur.dat
    104      & et surf.def'
    105       CALL abort_gcm (modname,abort_message,1)
     154!      RETURN
     155!9999  CONTINUE
     156!      abort_message ='Cette version demande les fichier rnatur.dat
     157!     & et surf.def'
     158!      CALL abort_gcm (modname,abort_message,1)
    106159
    107160      END
  • LMDZ5/branches/testing/libf/phylmd/iophy.F90

    r1539 r1707  
    5151   
    5252!$OMP MASTER 
    53     ALLOCATE(io_lat(jjm+1-1/iim))
     53    ALLOCATE(io_lat(jjm+1-1/(iim*jjm)))
    5454    io_lat(1)=rlat_glo(1)
    55     io_lat(jjm+1-1/iim)=rlat_glo(klon_glo)
    56     IF (iim > 1) then
     55    io_lat(jjm+1-1/(iim*jjm))=rlat_glo(klon_glo)
     56    IF ((iim*jjm) > 1) then
    5757      DO i=2,jjm
    5858        io_lat(i)=rlat_glo(2+(i-2)*iim)
     
    6161
    6262    ALLOCATE(io_lon(iim))
    63     io_lon(:)=rlon_glo(2-1/iim:iim+1-1/iim)
     63    io_lon(:)=rlon_glo(2-1/(iim*jjm):iim+1-1/(iim*jjm))
    6464
    6565    ddid=(/ 1,2 /)
    66     dsg=(/ iim, jjm+1-1/iim /)
     66    dsg=(/ iim, jjm+1-1/(iim*jjm) /)
    6767    dsl=(/ iim, jj_nb /)
    6868    dpf=(/ 1,jj_begin /)
     
    8989  include 'dimensions.h'   
    9090    real,dimension(iim),intent(in) :: lon
    91     real,dimension(jjm+1-1/iim),intent(in) :: lat
     91    real,dimension(jjm+1-1/(iim*jjm)),intent(in) :: lat
    9292
    9393    INTEGER,DIMENSION(2) :: ddid
     
    100100
    101101!$OMP MASTER 
    102     allocate(io_lat(jjm+1-1/iim))
     102    allocate(io_lat(jjm+1-1/(iim*jjm)))
    103103    io_lat(:)=lat(:)
    104104    allocate(io_lon(iim))
     
    106106   
    107107    ddid=(/ 1,2 /)
    108     dsg=(/ iim, jjm+1-1/iim /)
     108    dsg=(/ iim, jjm+1-1/(iim*jjm) /)
    109109    dsl=(/ iim, jj_nb /)
    110110    dpf=(/ 1,jj_begin /)
     
    234234
    235235       CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon_glo,zx_lon)
    236        if (iim.gt.1) then
     236       if ((iim*jjm).gt.1) then
    237237       DO i = 1, iim
    238238         zx_lon(i,1) = rlon_glo(i+1)
  • LMDZ5/branches/testing/libf/phylmd/mod_grid_phy_lmdz.F90

    r1001 r1707  
    11!
    2 !$Header$
     2!$Id $
    33!
    44MODULE mod_grid_phy_lmdz
     5
     6  PUBLIC
     7  PRIVATE :: grid1dTo2d_glo_igen, grid1dTo2d_glo_rgen, grid1dTo2d_glo_lgen, &
     8             grid2dTo1d_glo_igen, grid2dTo1d_glo_rgen, grid2dTo1d_glo_lgen
     9 
    510  INTEGER,SAVE :: nbp_lon  ! == iim
    611  INTEGER,SAVE :: nbp_lat  ! == jjmp1
     
    271276  END SUBROUTINE grid2dTo1d_glo_l3
    272277
    273 END MODULE mod_grid_phy_lmdz
    274 
    275 
     278!----------------------------------------------------------------
     279!  Generic (private) fonctions
     280!----------------------------------------------------------------
    276281 
    277282  SUBROUTINE grid1dTo2d_glo_igen(VarIn,VarOut,dimsize)
    278     USE mod_grid_phy_lmdz
     283
    279284    IMPLICIT NONE
    280285
     
    311316
    312317  SUBROUTINE grid1dTo2d_glo_rgen(VarIn,VarOut,dimsize)
    313     USE mod_grid_phy_lmdz
     318
    314319    IMPLICIT NONE
    315320
     
    345350
    346351  SUBROUTINE grid1dTo2d_glo_lgen(VarIn,VarOut,dimsize)
    347     USE mod_grid_phy_lmdz
     352
    348353    IMPLICIT NONE
    349354   
     
    379384 
    380385  SUBROUTINE grid2dTo1d_glo_igen(VarIn,VarOut,dimsize)
    381     USE mod_grid_phy_lmdz
     386
    382387    IMPLICIT NONE
    383388
     
    402407 
    403408  SUBROUTINE grid2dTo1d_glo_rgen(VarIn,VarOut,dimsize)
    404     USE mod_grid_phy_lmdz
     409
    405410    IMPLICIT NONE
    406411
     
    425430   
    426431  SUBROUTINE grid2dTo1d_glo_lgen(VarIn,VarOut,dimsize)
    427     USE mod_grid_phy_lmdz
     432
    428433    IMPLICIT NONE
    429434
     
    446451   
    447452  END SUBROUTINE grid2dTo1d_glo_lgen   
     453
     454END MODULE mod_grid_phy_lmdz
  • LMDZ5/branches/testing/libf/phylmd/pbl_surface_mod.F90

    r1664 r1707  
    172172       t,         q,         u,        v,             &
    173173       pplay,     paprs,     pctsrf,                  &
    174        ts,        alb1,      alb2,     u10m,  v10m,  &
     174       ts,        alb1,      alb2,ustar, u10m, v10m,  &
    175175       lwdown_m,  cdragh,    cdragm,   zu1,    zv1,   &
    176176       alb1_m,    alb2_m,    zxsens,   zxevap,        &
     
    181181       s_capCL,   s_oliqCL,  s_cteiCL, s_pblT,        &
    182182       s_therm,   s_trmb1,   s_trmb2,  s_trmb3,       &
    183        zxrugs,    zu10m,     zv10m,    fder_print,    &
     183       zxrugs,zustar,zu10m,  zv10m,    fder_print,    &
    184184       zxqsurf,   rh2m,      zxfluxu,  zxfluxv,       &
    185185       rugos_d,   agesno_d,  sollw,    solsw,         &
     
    288288    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: alb1    ! albedo in visible SW interval
    289289    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: alb2    ! albedo in near infra-red SW interval
     290    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: ustar   ! u* (m/s)
    290291    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: u10m    ! u speed at 10m
    291292    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: v10m    ! v speed at 10m
     
    330331    REAL, DIMENSION(klon),        INTENT(OUT)       :: s_trmb3    ! point Omega, mean for each grid point
    331332    REAL, DIMENSION(klon),        INTENT(OUT)       :: zxrugs     ! rugosity at surface (m), mean for each grid point
     333    REAL, DIMENSION(klon),        INTENT(OUT)       :: zustar     ! u*
    332334    REAL, DIMENSION(klon),        INTENT(OUT)       :: zu10m      ! u speed at 10m, mean for each grid point
    333335    REAL, DIMENSION(klon),        INTENT(OUT)       :: zv10m      ! v speed at 10m, mean for each grid point
     
    10191021       t2m(:,nsrf)    = 0.
    10201022       q2m(:,nsrf)    = 0.
     1023       ustar(:,nsrf)   = 0.
    10211024       u10m(:,nsrf)   = 0.
    10221025       v10m(:,nsrf)   = 0.
    1023 
    10241026       pblh(:,nsrf)   = 0.        ! Hauteur de couche limite
    10251027       plcl(:,nsrf)   = 0.        ! Niveau de condensation de la CLA
     
    10691071         
    10701072          ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman
     1073          ustar(i,nsrf)=yustar(j)
    10711074          u10m(i,nsrf)=(yu10m(j) * uzon(j))/SQRT(uzon(j)**2+vmer(j)**2)
    10721075          v10m(i,nsrf)=(yu10m(j) * vmer(j))/SQRT(uzon(j)**2+vmer(j)**2)
     1076
    10731077       END DO
    10741078
     
    11501154    zxtsol(:) = 0.0  ; zxfluxlat(:) = 0.0
    11511155    zt2m(:) = 0.0    ; zq2m(:) = 0.0
    1152     zu10m(:) = 0.0   ; zv10m(:) = 0.0
     1156    zustar(:)=0.0 ; zu10m(:) = 0.0   ; zv10m(:) = 0.0
    11531157    s_pblh(:) = 0.0  ; s_plcl(:) = 0.0
    11541158    s_capCL(:) = 0.0 ; s_oliqCL(:) = 0.0
     
    11721176          zt2m(i)  = zt2m(i)  + t2m(i,nsrf)  * pctsrf(i,nsrf)
    11731177          zq2m(i)  = zq2m(i)  + q2m(i,nsrf)  * pctsrf(i,nsrf)
     1178          zustar(i) = zustar(i) + ustar(i,nsrf) * pctsrf(i,nsrf)
    11741179          zu10m(i) = zu10m(i) + u10m(i,nsrf) * pctsrf(i,nsrf)
    11751180          zv10m(i) = zv10m(i) + v10m(i,nsrf) * pctsrf(i,nsrf)
     
    13051310!****************************************************************************************
    13061311!
    1307   SUBROUTINE pbl_surface_newfrac(itime, pctsrf_new, pctsrf_old, tsurf, alb1, alb2, u10m, v10m, tke)
     1312  SUBROUTINE pbl_surface_newfrac(itime, pctsrf_new, pctsrf_old, tsurf, alb1, alb2, ustar, u10m, v10m, tke)
    13081313
    13091314    ! Give default values where new fraction has appread
     
    13231328    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT)        :: tsurf
    13241329    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT)        :: alb1, alb2
    1325     REAL, DIMENSION(klon,nbsrf), INTENT(INOUT)        :: u10m, v10m
     1330    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT)        :: ustar,u10m, v10m
    13261331    REAL, DIMENSION(klon,klev+1,nbsrf), INTENT(INOUT) :: tke
    13271332
     
    13691374                alb1(i,nsrf)  = alb1(i,nsrf_comp1)
    13701375                alb2(i,nsrf)  = alb2(i,nsrf_comp1)
     1376                ustar(i,nsrf)  = ustar(i,nsrf_comp1)
    13711377                u10m(i,nsrf)  = u10m(i,nsrf_comp1)
    13721378                v10m(i,nsrf)  = v10m(i,nsrf_comp1)
     
    13831389                alb1(i,nsrf)  = alb1(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + alb1(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
    13841390                alb2(i,nsrf)  = alb2(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + alb2(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
     1391                ustar(i,nsrf)  = ustar(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + ustar(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
    13851392                u10m(i,nsrf)  = u10m(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + u10m(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
    13861393                v10m(i,nsrf)  = v10m(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + v10m(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
  • LMDZ5/branches/testing/libf/phylmd/phyaqua.F

    r1530 r1707  
    1616!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    1717
    18       use comgeomphy
    19       use dimphy
     18      use comgeomphy, only : rlatd,rlond
     19      use dimphy, only : klon
    2020      use surface_data, only : type_ocean,ok_veget
    2121      use pbl_surface_mod, only : pbl_surface_init
    2222      USE fonte_neige_mod, only : fonte_neige_init
    2323      use phys_state_var_mod
    24       use control_mod
    25 
     24      use control_mod, only : dayref,nday,iphysiq
    2625
    2726      USE IOIPSL
     
    3534#include "dimsoil.h"
    3635#include "indicesol.h"
    37 
    38       integer nlon,iflag_phys
     36#include "temps.h"
     37
     38      integer,intent(in) :: nlon,iflag_phys
    3939cIM ajout latfi, lonfi
    40       REAL, DIMENSION (nlon) :: lonfi, latfi
     40      real,intent(in) :: lonfi(nlon),latfi(nlon)
     41
    4142      INTEGER type_profil,type_aqua
    4243
     
    7172!      integer demih_pas
    7273
    73       integer day_ini
    74 
    7574      CHARACTER*80 ans,file_forctl, file_fordat, file_start
    7675      character*100 file,var
     
    8887      REAL phy_flic(nlon,360)
    8988
    90       integer, save::  read_climoz ! read ozone climatology
     89      integer, save::  read_climoz=0 ! read ozone climatology
    9190
    9291
     
    131130      type_aqua=iflag_phys/100
    132131      type_profil=iflag_phys-type_aqua*100
    133       print*,'type_aqua, type_profil',type_aqua, type_profil
    134 
    135       if (klon.ne.nlon) stop'probleme de dimensions dans iniaqua'
     132      print*,'iniaqua:type_aqua, type_profil',type_aqua, type_profil
     133
     134      if (klon.ne.nlon) then
     135        write(*,*)"iniaqua: klon=",klon," nlon=",nlon
     136        stop'probleme de dimensions dans iniaqua'
     137      endif
    136138      call phys_state_var_init(read_climoz)
    137139
     
    154156
    155157         day_ini=dayref
     158         day_end=day_ini+nday
    156159         airefi=1.
    157160         zcufi=1.
     
    171174      radsol=0.
    172175      qsol_f=10.
    173       CALL getin('albedo',albedo)
     176!      CALL getin('albedo',albedo) ! albedo is set below, depending on type_aqua
    174177      alb_ocean=.true.
    175178      CALL getin('alb_ocean',alb_ocean)
     
    180183      qsol(:)    = qsol_f
    181184      rugsrel = 0.0    ! (rugsrel = rugoro)
     185      rugoro = 0.0
     186      u_ancien = 0.0
     187      v_ancien = 0.0
    182188      agesno  = 50.0
    183189! Relief plat
     
    308314     .     evap, frugs, agesno, tsoil)
    309315
    310         print*,'avant phyredem dans iniaqua'
     316        print*,'iniaqua: before phyredem'
    311317
    312318      falb1=albedo
     
    329335      CALL phyredem ("startphy.nc")
    330336
    331         print*,'apres phyredem'
     337        print*,'iniaqua: after phyredem'
    332338      call phys_state_var_end
    333339
     
    450456      RETURN
    451457      END
     458
     459!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     460
    452461      subroutine writelim
    453462     s   (klon,phy_nat,phy_alb,phy_sst,phy_bil,phy_rug,phy_ice,
    454463     s    phy_fter,phy_foce,phy_flic,phy_fsic)
    455464c
     465      use mod_phys_lmdz_para, only: is_mpi_root,is_omp_root
     466      use mod_grid_phy_lmdz, only : klon_glo
     467      use mod_phys_lmdz_transfert_para, only : gather
    456468!#include "dimensions.h"
    457469!#include "dimphy.h"
    458470#include "netcdf.inc"
    459471 
    460       integer klon
    461       REAL phy_nat(klon,360)
    462       REAL phy_alb(klon,360)
    463       REAL phy_sst(klon,360)
    464       REAL phy_bil(klon,360)
    465       REAL phy_rug(klon,360)
    466       REAL phy_ice(klon,360)
    467       REAL phy_fter(klon,360)
    468       REAL phy_foce(klon,360)
    469       REAL phy_flic(klon,360)
    470       REAL phy_fsic(klon,360)
    471  
     472      integer,intent(in) :: klon
     473      real,intent(in) :: phy_nat(klon,360)
     474      real,intent(in) :: phy_alb(klon,360)
     475      real,intent(in) :: phy_sst(klon,360)
     476      real,intent(in) :: phy_bil(klon,360)
     477      real,intent(in) :: phy_rug(klon,360)
     478      real,intent(in) :: phy_ice(klon,360)
     479      real,intent(in) :: phy_fter(klon,360)
     480      real,intent(in) :: phy_foce(klon,360)
     481      real,intent(in) :: phy_flic(klon,360)
     482      real,intent(in) :: phy_fsic(klon,360)
     483
     484      real :: phy_glo(klon_glo,360) ! temporary variable, to store phy_***(:)
     485                                    ! on the whole physics grid
    472486      INTEGER ierr
    473487      INTEGER dimfirst(3)
     
    480494      INTEGER id_FTER,id_FOCE,id_FSIC,id_FLIC
    481495 
    482       PRINT*, 'Ecriture du fichier limit'
    483 c
    484       ierr = NF_CREATE ("limit.nc", NF_CLOBBER, nid)
    485 c
    486       ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 30,
     496      if (is_mpi_root.and.is_omp_root) then
     497     
     498        PRINT*, 'writelim: Ecriture du fichier limit'
     499c
     500        ierr = NF_CREATE ("limit.nc", NF_CLOBBER, nid)
     501c
     502        ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 30,
    487503     .                       "Fichier conditions aux limites")
    488       ierr = NF_DEF_DIM (nid, "points_physiques", klon, ndim)
    489       ierr = NF_DEF_DIM (nid, "time", NF_UNLIMITED, ntim)
    490 c
    491       dims(1) = ndim
    492       dims(2) = ntim
     504!!        ierr = NF_DEF_DIM (nid, "points_physiques", klon, ndim)
     505        ierr = NF_DEF_DIM (nid, "points_physiques", klon_glo, ndim)
     506        ierr = NF_DEF_DIM (nid, "time", NF_UNLIMITED, ntim)
     507c
     508        dims(1) = ndim
     509        dims(2) = ntim
    493510c
    494511ccc      ierr = NF_DEF_VAR (nid, "TEMPS", NF_DOUBLE, 1,ntim, id_tim)
    495       ierr = NF_DEF_VAR (nid, "TEMPS", NF_FLOAT, 1,ntim, id_tim)
    496       ierr = NF_PUT_ATT_TEXT (nid, id_tim, "title", 17,
     512        ierr = NF_DEF_VAR (nid, "TEMPS", NF_FLOAT, 1,ntim, id_tim)
     513        ierr = NF_PUT_ATT_TEXT (nid, id_tim, "title", 17,
    497514     .                        "Jour dans l annee")
    498515ccc      ierr = NF_DEF_VAR (nid, "NAT", NF_DOUBLE, 2,dims, id_NAT)
    499       ierr = NF_DEF_VAR (nid, "NAT", NF_FLOAT, 2,dims, id_NAT)
    500       ierr = NF_PUT_ATT_TEXT (nid, id_NAT, "title", 23,
     516        ierr = NF_DEF_VAR (nid, "NAT", NF_FLOAT, 2,dims, id_NAT)
     517        ierr = NF_PUT_ATT_TEXT (nid, id_NAT, "title", 23,
    501518     .                        "Nature du sol (0,1,2,3)")
    502519ccc      ierr = NF_DEF_VAR (nid, "SST", NF_DOUBLE, 2,dims, id_SST)
    503       ierr = NF_DEF_VAR (nid, "SST", NF_FLOAT, 2,dims, id_SST)
    504       ierr = NF_PUT_ATT_TEXT (nid, id_SST, "title", 35,
     520        ierr = NF_DEF_VAR (nid, "SST", NF_FLOAT, 2,dims, id_SST)
     521        ierr = NF_PUT_ATT_TEXT (nid, id_SST, "title", 35,
    505522     .                        "Temperature superficielle de la mer")
    506523ccc      ierr = NF_DEF_VAR (nid, "BILS", NF_DOUBLE, 2,dims, id_BILS)
    507       ierr = NF_DEF_VAR (nid, "BILS", NF_FLOAT, 2,dims, id_BILS)
    508       ierr = NF_PUT_ATT_TEXT (nid, id_BILS, "title", 32,
     524        ierr = NF_DEF_VAR (nid, "BILS", NF_FLOAT, 2,dims, id_BILS)
     525        ierr = NF_PUT_ATT_TEXT (nid, id_BILS, "title", 32,
    509526     .                        "Reference flux de chaleur au sol")
    510527ccc      ierr = NF_DEF_VAR (nid, "ALB", NF_DOUBLE, 2,dims, id_ALB)
    511       ierr = NF_DEF_VAR (nid, "ALB", NF_FLOAT, 2,dims, id_ALB)
    512       ierr = NF_PUT_ATT_TEXT (nid, id_ALB, "title", 19,
     528        ierr = NF_DEF_VAR (nid, "ALB", NF_FLOAT, 2,dims, id_ALB)
     529        ierr = NF_PUT_ATT_TEXT (nid, id_ALB, "title", 19,
    513530     .                        "Albedo a la surface")
    514531ccc      ierr = NF_DEF_VAR (nid, "RUG", NF_DOUBLE, 2,dims, id_RUG)
    515       ierr = NF_DEF_VAR (nid, "RUG", NF_FLOAT, 2,dims, id_RUG)
    516       ierr = NF_PUT_ATT_TEXT (nid, id_RUG, "title", 8,
     532        ierr = NF_DEF_VAR (nid, "RUG", NF_FLOAT, 2,dims, id_RUG)
     533        ierr = NF_PUT_ATT_TEXT (nid, id_RUG, "title", 8,
    517534     .                        "Rugosite")
    518535
    519       ierr = NF_DEF_VAR (nid, "FTER", NF_FLOAT, 2,dims, id_FTER)
    520       ierr = NF_PUT_ATT_TEXT (nid, id_FTER, "title", 8,"Frac. Terre")
    521       ierr = NF_DEF_VAR (nid, "FOCE", NF_FLOAT, 2,dims, id_FOCE)
    522       ierr = NF_PUT_ATT_TEXT (nid, id_FOCE, "title", 8,"Frac. Terre")
    523       ierr = NF_DEF_VAR (nid, "FSIC", NF_FLOAT, 2,dims, id_FSIC)
    524       ierr = NF_PUT_ATT_TEXT (nid, id_FSIC, "title", 8,"Frac. Terre")
    525       ierr = NF_DEF_VAR (nid, "FLIC", NF_FLOAT, 2,dims, id_FLIC)
    526       ierr = NF_PUT_ATT_TEXT (nid, id_FLIC, "title", 8,"Frac. Terre")
    527 c
    528       ierr = NF_ENDDEF(nid)
    529 c
    530       DO k = 1, 360
    531 c
    532       debut(1) = 1
    533       debut(2) = k
    534       epais(1) = klon
    535       epais(2) = 1
    536 c
    537       print*,'Instant ',k
    538 #ifdef NC_DOUBLE
    539       print*,'NC DOUBLE'
    540       ierr = NF_PUT_VAR1_DOUBLE (nid,id_tim,k,DBLE(k))
    541       ierr = NF_PUT_VARA_DOUBLE (nid,id_NAT,debut,epais,phy_nat(1,k))
    542       ierr = NF_PUT_VARA_DOUBLE (nid,id_SST,debut,epais,phy_sst(1,k))
    543       ierr = NF_PUT_VARA_DOUBLE (nid,id_BILS,debut,epais,phy_bil(1,k))
    544       ierr = NF_PUT_VARA_DOUBLE (nid,id_ALB,debut,epais,phy_alb(1,k))
    545       ierr = NF_PUT_VARA_DOUBLE (nid,id_RUG,debut,epais,phy_rug(1,k))
    546       ierr = NF_PUT_VARA_DOUBLE (nid,id_FTER,debut,epais,phy_fter(1,k))
    547       ierr = NF_PUT_VARA_DOUBLE (nid,id_FOCE,debut,epais,phy_foce(1,k))
    548       ierr = NF_PUT_VARA_DOUBLE (nid,id_FSIC,debut,epais,phy_fsic(1,k))
    549       ierr = NF_PUT_VARA_DOUBLE (nid,id_FLIC,debut,epais,phy_flic(1,k))
    550 #else
    551       print*,'NC PAS DOUBLE'
    552       ierr = NF_PUT_VAR1_REAL (nid,id_tim,k,FLOAT(k))
    553       ierr = NF_PUT_VARA_REAL (nid,id_NAT,debut,epais,phy_nat(1,k))
    554       ierr = NF_PUT_VARA_REAL (nid,id_SST,debut,epais,phy_sst(1,k))
    555       ierr = NF_PUT_VARA_REAL (nid,id_BILS,debut,epais,phy_bil(1,k))
    556       ierr = NF_PUT_VARA_REAL (nid,id_ALB,debut,epais,phy_alb(1,k))
    557       ierr = NF_PUT_VARA_REAL (nid,id_RUG,debut,epais,phy_rug(1,k))
    558       ierr = NF_PUT_VARA_REAL (nid,id_FTER,debut,epais,phy_fter(1,k))
    559       ierr = NF_PUT_VARA_REAL (nid,id_FOCE,debut,epais,phy_foce(1,k))
    560       ierr = NF_PUT_VARA_REAL (nid,id_FSIC,debut,epais,phy_fsic(1,k))
    561       ierr = NF_PUT_VARA_REAL (nid,id_FLIC,debut,epais,phy_flic(1,k))
    562 
    563 #endif
    564 c
    565       ENDDO
    566 c
    567       ierr = NF_CLOSE(nid)
    568 c
    569       return
     536        ierr = NF_DEF_VAR (nid, "FTER", NF_FLOAT, 2,dims, id_FTER)
     537        ierr = NF_PUT_ATT_TEXT (nid, id_FTER, "title", 8,"Frac. Terre")
     538        ierr = NF_DEF_VAR (nid, "FOCE", NF_FLOAT, 2,dims, id_FOCE)
     539        ierr = NF_PUT_ATT_TEXT (nid, id_FOCE, "title", 8,"Frac. Terre")
     540        ierr = NF_DEF_VAR (nid, "FSIC", NF_FLOAT, 2,dims, id_FSIC)
     541        ierr = NF_PUT_ATT_TEXT (nid, id_FSIC, "title", 8,"Frac. Terre")
     542        ierr = NF_DEF_VAR (nid, "FLIC", NF_FLOAT, 2,dims, id_FLIC)
     543        ierr = NF_PUT_ATT_TEXT (nid, id_FLIC, "title", 8,"Frac. Terre")
     544c
     545        ierr = NF_ENDDEF(nid)
     546c
     547
     548! write the 'times'
     549        do k=1,360
     550#ifdef NC_DOUBLE
     551          ierr = NF_PUT_VAR1_DOUBLE (nid,id_tim,k,DBLE(k))
     552#else
     553          ierr = NF_PUT_VAR1_REAL (nid,id_tim,k,FLOAT(k))
     554#endif
     555        enddo
     556
     557      endif ! of if (is_mpi_root.and.is_omp_root)
     558
     559! write the fields, after having collected them on master
     560
     561      call gather(phy_nat,phy_glo)
     562      if (is_mpi_root.and.is_omp_root) then
     563#ifdef NC_DOUBLE
     564        ierr=NF_PUT_VAR_DOUBLE(nid,id_NAT,phy_glo)
     565#else
     566        ierr=NF_PUT_VAR_REAL(nid,id_NAT,phy_glo)
     567#endif
     568        if(ierr.ne.NF_NOERR) then
     569          write(*,*) "writelim error with phy_nat"
     570          write(*,*) NF_STRERROR(ierr)
     571        endif
     572      endif
     573
     574      call gather(phy_sst,phy_glo)
     575      if (is_mpi_root.and.is_omp_root) then
     576#ifdef NC_DOUBLE
     577        ierr=NF_PUT_VAR_DOUBLE(nid,id_SST,phy_glo)
     578#else
     579        ierr=NF_PUT_VAR_REAL(nid,id_SST,phy_glo)
     580#endif
     581        if(ierr.ne.NF_NOERR) then
     582          write(*,*) "writelim error with phy_sst"
     583          write(*,*) NF_STRERROR(ierr)
     584        endif
     585      endif
     586
     587      call gather(phy_bil,phy_glo)
     588      if (is_mpi_root.and.is_omp_root) then
     589#ifdef NC_DOUBLE
     590        ierr=NF_PUT_VAR_DOUBLE(nid,id_BILS,phy_glo)
     591#else
     592        ierr=NF_PUT_VAR_REAL(nid,id_BILS,phy_glo)
     593#endif
     594        if(ierr.ne.NF_NOERR) then
     595          write(*,*) "writelim error with phy_bil"
     596          write(*,*) NF_STRERROR(ierr)
     597        endif
     598      endif
     599
     600      call gather(phy_alb,phy_glo)
     601      if (is_mpi_root.and.is_omp_root) then
     602#ifdef NC_DOUBLE
     603        ierr=NF_PUT_VAR_DOUBLE(nid,id_ALB,phy_glo)
     604#else
     605        ierr=NF_PUT_VAR_REAL(nid,id_ALB,phy_glo)
     606#endif
     607        if(ierr.ne.NF_NOERR) then
     608          write(*,*) "writelim error with phy_alb"
     609          write(*,*) NF_STRERROR(ierr)
     610        endif
     611      endif
     612
     613      call gather(phy_rug,phy_glo)
     614      if (is_mpi_root.and.is_omp_root) then
     615#ifdef NC_DOUBLE
     616        ierr=NF_PUT_VAR_DOUBLE(nid,id_RUG,phy_glo)
     617#else
     618        ierr=NF_PUT_VAR_REAL(nid,id_RUG,phy_glo)
     619#endif
     620        if(ierr.ne.NF_NOERR) then
     621          write(*,*) "writelim error with phy_rug"
     622          write(*,*) NF_STRERROR(ierr)
     623        endif
     624      endif
     625
     626      call gather(phy_fter,phy_glo)
     627      if (is_mpi_root.and.is_omp_root) then
     628#ifdef NC_DOUBLE
     629        ierr=NF_PUT_VAR_DOUBLE(nid,id_FTER,phy_glo)
     630#else
     631        ierr=NF_PUT_VAR_REAL(nid,id_FTER,phy_glo)
     632#endif
     633        if(ierr.ne.NF_NOERR) then
     634          write(*,*) "writelim error with phy_fter"
     635          write(*,*) NF_STRERROR(ierr)
     636        endif
     637      endif
     638
     639      call gather(phy_foce,phy_glo)
     640      if (is_mpi_root.and.is_omp_root) then
     641#ifdef NC_DOUBLE
     642        ierr=NF_PUT_VAR_DOUBLE(nid,id_FOCE,phy_glo)
     643#else
     644        ierr=NF_PUT_VAR_REAL(nid,id_FOCE,phy_glo)
     645#endif
     646        if(ierr.ne.NF_NOERR) then
     647          write(*,*) "writelim error with phy_foce"
     648          write(*,*) NF_STRERROR(ierr)
     649        endif
     650      endif
     651
     652      call gather(phy_fsic,phy_glo)
     653      if (is_mpi_root.and.is_omp_root) then
     654#ifdef NC_DOUBLE
     655        ierr=NF_PUT_VAR_DOUBLE(nid,id_FSIC,phy_glo)
     656#else
     657        ierr=NF_PUT_VAR_REAL(nid,id_FSIC,phy_glo)
     658#endif
     659        if(ierr.ne.NF_NOERR) then
     660          write(*,*) "writelim error with phy_fsic"
     661          write(*,*) NF_STRERROR(ierr)
     662        endif
     663      endif
     664
     665      call gather(phy_flic,phy_glo)
     666      if (is_mpi_root.and.is_omp_root) then
     667#ifdef NC_DOUBLE
     668        ierr=NF_PUT_VAR_DOUBLE(nid,id_FLIC,phy_glo)
     669#else
     670        ierr=NF_PUT_VAR_REAL(nid,id_FLIC,phy_glo)
     671#endif
     672        if(ierr.ne.NF_NOERR) then
     673          write(*,*) "writelim error with phy_flic"
     674          write(*,*) NF_STRERROR(ierr)
     675        endif
     676      endif
     677
     678!  close file:
     679      if (is_mpi_root.and.is_omp_root) then
     680        ierr = NF_CLOSE(nid)
     681      endif
     682
    570683      end
     684
     685!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    571686
    572687      SUBROUTINE profil_sst(nlon, rlatd, type_profil, phy_sst)
  • LMDZ5/branches/testing/libf/phylmd/phyetat0.F

    r1665 r1707  
    7676c FH1D
    7777c     real iolat(jjm+1)
    78       real iolat(jjm+1-1/iim)
     78      real iolat(jjm+1-1/(iim*jjm))
    7979c
    8080c Ouvrir le fichier contenant l'etat initial:
  • LMDZ5/branches/testing/libf/phylmd/phys_output_mod.F90

    r1669 r1707  
    8181  type(ctrl_out),save :: o_sicf         = ctrl_out((/ 1, 1, 10, 10, 10, 10 /),'sicf')
    8282  type(ctrl_out),save :: o_q2m          = ctrl_out((/ 1, 1, 1, 5, 10, 10 /),'q2m')
     83  type(ctrl_out),save :: o_ustar        = ctrl_out((/ 1, 1, 1, 5, 10, 10 /),'ustar')
    8384  type(ctrl_out),save :: o_u10m         = ctrl_out((/ 1, 1, 1, 5, 10, 10 /),'u10m')
    8485  type(ctrl_out),save :: o_v10m         = ctrl_out((/ 1, 1, 1, 5, 10, 10 /),'v10m')
     
    8687  type(ctrl_out),save :: o_qsurf        = ctrl_out((/ 1, 10, 10, 10, 10, 10 /),'qsurf')
    8788
     89  type(ctrl_out),save,dimension(4) :: o_ustar_srf     = (/ ctrl_out((/ 10, 6, 10, 10, 10, 10 /),'ustar_ter'), &
     90       ctrl_out((/ 10, 6, 10, 10, 10, 10 /),'ustar_lic'), &
     91       ctrl_out((/ 10, 6, 10, 10, 10, 10 /),'ustar_oce'), &
     92       ctrl_out((/ 10, 6, 10, 10, 10, 10 /),'ustar_sic') /)
    8893  type(ctrl_out),save,dimension(4) :: o_u10m_srf     = (/ ctrl_out((/ 10, 6, 10, 10, 10, 10 /),'u10m_ter'), &
    8994       ctrl_out((/ 10, 6, 10, 10, 10, 10 /),'u10m_lic'), &
     
    585590
    586591  type(ctrl_out),save,allocatable :: o_trac(:)
     592  type(ctrl_out),save,allocatable :: o_trac_cum(:)
    587593
    588594  type(ctrl_out),save :: o_rsu        = ctrl_out((/ 4, 10, 10, 10, 10, 10 /),'rsu')
     
    719725
    720726    if (.not. allocated(o_trac)) ALLOCATE(o_trac(nqtot))
     727    if (.not. allocated(o_trac_cum)) ALLOCATE(o_trac_cum(nqtot))
    721728
    722729    levmax = (/ klev, klev, klev, klev, klev, klev /)
     
    960967          CALL histdef2d(iff,clef_stations(iff),o_sicf%flag,o_sicf%name, "Sea-ice fraction", "-" )
    961968          CALL histdef2d(iff,clef_stations(iff),o_q2m%flag,o_q2m%name, "Specific humidity 2m", "kg/kg")
     969          CALL histdef2d(iff,clef_stations(iff),o_ustar%flag,o_ustar%name, "Friction velocity", "m/s" )
    962970          CALL histdef2d(iff,clef_stations(iff),o_u10m%flag,o_u10m%name, "Vent zonal 10m", "m/s" )
    963971          CALL histdef2d(iff,clef_stations(iff),o_v10m%flag,o_v10m%name, "Vent meridien 10m", "m/s")
     
    969977          endif
    970978
     979             type_ecri(1) = 'inst(X)'
     980             type_ecri(2) = 'inst(X)'
     981             type_ecri(3) = 'inst(X)'
     982             type_ecri(4) = 'inst(X)'
     983             type_ecri(5) = 'inst(X)'
     984             type_ecri(6) = 'inst(X)'
    971985          CALL histdef2d(iff,clef_stations(iff),o_ndayrain%flag,o_ndayrain%name, "Number of dayrain(liq+sol)", "-")
     986             type_ecri(:) = type_ecri_files(:)
    972987          CALL histdef2d(iff,clef_stations(iff),o_precip%flag,o_precip%name, "Precip Totale liq+sol", "kg/(s*m2)" )
    973988          CALL histdef2d(iff,clef_stations(iff),o_plul%flag,o_plul%name, "Large-scale Precip.", "kg/(s*m2)")
     
    10271042                  o_tsol_srf(nsrf)%flag,o_tsol_srf(nsrf)%name,"Temperature "//clnsurf(nsrf),"K")
    10281043             CALL histdef2d(iff,clef_stations(iff), &
     1044                  o_ustar_srf(nsrf)%flag,o_ustar_srf(nsrf)%name,"Friction velocity "//clnsurf(nsrf),"m/s")
     1045             CALL histdef2d(iff,clef_stations(iff), &
    10291046                  o_u10m_srf(nsrf)%flag,o_u10m_srf(nsrf)%name,"Vent Zonal 10m "//clnsurf(nsrf),"m/s")
    10301047             CALL histdef2d(iff,clef_stations(iff), &
     
    17561773                o_trac(iq-2) = ctrl_out((/ 4, 5, 1, 1, 1, 10 /),tname(iiq))
    17571774                CALL histdef3d (iff,clef_stations(iff), &
    1758                      o_trac(iq-2)%flag,o_trac(iq-2)%name,'Tracer '//ttext(iiq), "-" )
     1775                o_trac(iq-2)%flag,o_trac(iq-2)%name,'Tracer '//ttext(iiq), "-" )
     1776                o_trac_cum(iq-2) = ctrl_out((/ 3, 4, 10, 10, 10, 10 /),'cum'//tname(iiq))
     1777                CALL histdef2d (iff,clef_stations(iff), &
     1778                o_trac_cum(iq-2)%flag,o_trac_cum(iq-2)%name,'Cumulated tracer '//ttext(iiq), "-" )
    17591779             ENDDO
    17601780          ENDIF
  • LMDZ5/branches/testing/libf/phylmd/phys_output_write.h

    r1669 r1707  
    101101      CALL histwrite_phy(nid_files(iff),clef_stations(iff),
    102102     $o_q2m%name,itau_w,zq2m)
     103       ENDIF
     104
     105       IF (o_ustar%flag(iff)<=lev_files(iff)) THEN
     106      CALL histwrite_phy(nid_files(iff),clef_stations(iff),
     107     $o_ustar%name,itau_w,zustar)
    103108       ENDIF
    104109
     
    437442     $      zx_tmp_fi2d)
    438443        ENDIF
     444
     445      IF (o_ustar_srf(nsrf)%flag(iff)<=lev_files(iff)) THEN
     446      zx_tmp_fi2d(1 : klon) = ustar(1 : klon, nsrf)
     447      CALL histwrite_phy(nid_files(iff),clef_stations(iff),
     448     $o_ustar_srf(nsrf)%name,
     449     $                 itau_w,zx_tmp_fi2d)
     450      ENDIF
    439451
    440452      IF (o_u10m_srf(nsrf)%flag(iff)<=lev_files(iff)) THEN
     
    22482260       ENDIF
    22492261         ENDDO
     2262         DO iq=3,nqtot
     2263       IF (o_trac_cum(iq-2)%flag(iff)<=lev_files(iff)) THEN
     2264         zx_tmp_fi2d=0.
     2265         do k=1,klev
     2266            zx_tmp_fi2d(:)=zx_tmp_fi2d(:)+zmasse(:,k)*qx(:,k,iq)
     2267         enddo
     2268         CALL histwrite_phy(nid_files(iff),clef_stations(iff),
     2269     s                  o_trac_cum(iq-2)%name,itau_w,zx_tmp_fi2d)
     2270
     2271       ENDIF
     2272         ENDDO
    22502273        endif
    22512274
  • LMDZ5/branches/testing/libf/phylmd/phys_state_var_mod.F90

    r1669 r1707  
    326326      REAL,SAVE,ALLOCATABLE :: newsst(:)
    327327!$OMP THREADPRIVATE(newsst)
    328       REAL,SAVE,ALLOCATABLE :: u10m(:,:), v10m(:,:)
    329 !$OMP THREADPRIVATE(u10m,v10m)
     328      REAL,SAVE,ALLOCATABLE :: ustar(:,:),u10m(:,:), v10m(:,:)
     329!$OMP THREADPRIVATE(ustar,u10m,v10m)
    330330!
    331331! ok_ade=T -ADE=topswad-topsw
     
    496496      ALLOCATE(rlonPOS(klon))
    497497      ALLOCATE(newsst(klon))
    498       ALLOCATE(u10m(klon,nbsrf), v10m(klon,nbsrf))
     498      ALLOCATE(ustar(klon,nbsrf),u10m(klon,nbsrf), v10m(klon,nbsrf))
    499499      ALLOCATE(topswad(klon), solswad(klon))
    500500      ALLOCATE(topswai(klon), solswai(klon))
     
    606606      deallocate(rlonPOS)
    607607      deallocate(newsst)
    608       deallocate(u10m, v10m)
     608      deallocate(ustar,u10m, v10m)
    609609      deallocate(topswad, solswad)
    610610      deallocate(topswai, solswai)
  • LMDZ5/branches/testing/libf/phylmd/physiq.F

    r1669 r1707  
    178178      save iflag_ratqs
    179179c$OMP THREADPRIVATE(iflag_ratqs)
    180       real facteur,zfratqs1,zfratqs2
     180      real facteur
    181181
    182182      REAL zz,znum,zden
     
    257257c variables a une pression donnee
    258258c
    259       real rlevSTD(nlevSTD)
    260       DATA rlevSTD/100000., 92500., 85000., 70000.,
    261      .60000., 50000., 40000., 30000., 25000., 20000.,
    262      .15000., 10000., 7000., 5000., 3000., 2000., 1000./
    263       SAVE rlevstd
    264 c$OMP THREADPRIVATE(rlevstd)
    265       CHARACTER*4 clevSTD(nlevSTD)
    266       DATA clevSTD/'1000','925 ','850 ','700 ','600 ',
    267      .'500 ','400 ','300 ','250 ','200 ','150 ','100 ',
    268      .'70  ','50  ','30  ','20  ','10  '/
    269       SAVE clevSTD
    270 c$OMP THREADPRIVATE(clevSTD)
     259#include "declare_STDlev.h"
    271260c
    272261      CHARACTER*4 bb2
    273262      CHARACTER*2 bb3
    274 
    275       real twriteSTD(klon,nlevSTD,nfiles)
    276       real qwriteSTD(klon,nlevSTD,nfiles)
    277       real rhwriteSTD(klon,nlevSTD,nfiles)
    278       real phiwriteSTD(klon,nlevSTD,nfiles)
    279       real uwriteSTD(klon,nlevSTD,nfiles)
    280       real vwriteSTD(klon,nlevSTD,nfiles)
    281       real wwriteSTD(klon,nlevSTD,nfiles)
    282 cIM for NMC files
    283       REAL geo500(klon)
    284       real :: rlevSTD3(nlevSTD3)
    285       DATA rlevSTD3/85000., 50000., 25000./
    286       SAVE rlevSTD3
    287 c$OMP THREADPRIVATE(rlevSTD3)
    288       real :: rlevSTD8(nlevSTD8)
    289       DATA rlevSTD8/100000., 85000., 70000., 50000., 25000., 10000.,
    290      $     5000., 1000./
    291       SAVE rlevSTD8
    292 c$OMP THREADPRIVATE(rlevSTD8)
    293       real twriteSTD3(klon,nlevSTD3)
    294       real qwriteSTD3(klon,nlevSTD3)
    295       real rhwriteSTD3(klon,nlevSTD3)
    296       real phiwriteSTD3(klon,nlevSTD3)
    297       real uwriteSTD3(klon,nlevSTD3)
    298       real vwriteSTD3(klon,nlevSTD3)
    299       real wwriteSTD3(klon,nlevSTD3)
    300 c
    301       real tnondefSTD8(klon,nlevSTD8)
    302       real twriteSTD8(klon,nlevSTD8)
    303       real qwriteSTD8(klon,nlevSTD8)
    304       real rhwriteSTD8(klon,nlevSTD8)
    305       real phiwriteSTD8(klon,nlevSTD8)
    306       real uwriteSTD8(klon,nlevSTD8)
    307       real vwriteSTD8(klon,nlevSTD8)
    308       real wwriteSTD8(klon,nlevSTD8)
    309 c
    310 c plevSTD3 END
    311 c
    312 c nout : niveau de output des variables a une pression donnee
    313       logical oknondef(klon,nlevSTD,nout)
    314 c
    315 c les produits uvSTD, vqSTD, .., T2STD sont calcules
    316 c a partir des valeurs instantannees toutes les 6 h
    317 c qui sont moyennees sur le mois
    318263c
    319264#include "radopt.h"
     
    958903      REAL snow_lsc(klon)
    959904c
    960       REAL ratqss(klon,klev),ratqsc(klon,klev)
     905      REAL ratqsc(klon,klev)
    961906      real ratqsbas,ratqshaut,tau_ratqs
    962907      save ratqsbas,ratqshaut,tau_ratqs
     
    1050995      REAL zx_tmp_fi3d(klon,klev) ! variable temporaire pour champs 3D
    1051996      REAL zx_tmp_fi3d1(klon,klev+1) !variable temporaire pour champs 3D (kelvp1)
    1052 c#ifdef histNMC
    1053 cym   A voir plus tard !!!!
    1054 cym      REAL zx_tmp_NC(iim,jjmp1,nlevSTD)
    1055       REAL zx_tmp_fiNC(klon,nlevSTD)
    1056 c#endif
    1057997      REAL(KIND=8) zx_tmp2_fi3d(klon,klev) ! variable temporaire pour champs 3D
    1058998      REAL zx_tmp_2d(iim,jjmp1), zx_tmp_3d(iim,jjmp1,klev)
    1059999      REAL zx_lon(iim,jjmp1), zx_lat(iim,jjmp1)
    1060 cIM for NMC files
    1061       REAL missing_val
    1062       REAL, SAVE :: freq_moyNMC(nout)
    1063 c$OMP THREADPRIVATE(freq_moyNMC)
    10641000c
    10651001      INTEGER nid_day, nid_mth, nid_ins, nid_mthnmc, nid_daynmc
     
    11371073      REAL q2m(klon,nbsrf)  ! humidite a 2m
    11381074
    1139 cIM: t2m, q2m, u10m, v10m et t2mincels, t2maxcels
     1075cIM: t2m, q2m, ustar, u10m, v10m et t2mincels, t2maxcels
    11401076      REAL zt2m(klon), zq2m(klon)             !temp., hum. 2m moyenne s/ 1 maille
    1141       REAL zu10m(klon), zv10m(klon)           !vents a 10m moyennes s/1 maille
     1077      REAL zustar(klon),zu10m(klon), zv10m(klon)  ! u* et vents a 10m moyennes s/1 maille
    11421078      CHARACTER*40 t2mincels, t2maxcels       !t2m min., t2m max
    11431079      CHARACTER*40 tinst, tave, typeval
     
    12551191      integer iostat
    12561192
    1257 cIM for NMC files
    1258       missing_val=nf90_fill_real
    12591193c======================================================================
    12601194! Gestion calendrier : mise a jour du module phys_cal_mod
     
    13261260      call phys_output_var_init
    13271261      print*, '================================================='
    1328 cIM for NMC files
    1329 cIM freq_moyNMC = frequences auxquelles on moyenne les champs accumules
    1330 cIM               sur les niveaux de pression standard du NMC
    1331       DO n=1, nout
    1332        freq_moyNMC(n)=freq_outNMC(n)/freq_calNMC(n)
    1333       ENDDO
    1334 c
    1335 cIM beg
     1262c
    13361263          dnwd0=0.0
    13371264          ftd=0.0
     
    13811308         lalim_conv(:)=1
    13821309cRC
     1310         ustar(:,:)=0.
    13831311         u10m(:,:)=0.
    13841312         v10m(:,:)=0.
     
    17681696!
    17691697      CALL change_srf_frac(itap, dtime, days_elapsed+1,
    1770      *     pctsrf, falb1, falb2, ftsol, u10m, v10m, pbl_tke)
     1698     *     pctsrf, falb1, falb2, ftsol, ustar, u10m, v10m, pbl_tke)
    17711699
    17721700
     
    20782006     e     t_seri,    q_seri,    u_seri,  v_seri,   
    20792007     e     pplay,     paprs,     pctsrf,           
    2080      +     ftsol,     falb1,     falb2,   u10m,   v10m,
     2008     +     ftsol,     falb1,     falb2,   ustar, u10m,   v10m,
    20812009     s     sollwdown, cdragh,    cdragm,  u1,    v1,
    20822010     s     albsol1,   albsol2,   sens,    evap, 
     
    20872015     d     s_capCL,   s_oliqCL,  s_cteiCL,s_pblT,
    20882016     d     s_therm,   s_trmb1,   s_trmb2, s_trmb3,
    2089      d     zxrugs,    zu10m,     zv10m,   fder,
     2017     d     zxrugs,    zustar, zu10m,     zv10m,   fder,
    20902018     d     zxqsurf,   rh2m,      zxfluxu, zxfluxv,
    20912019     d     frugs,     agesno,    fsollw,  fsolsw,
     
    28162744
    28172745c-------------------------------------------------------------------------
    2818 c  Caclul des ratqs
    2819 c-------------------------------------------------------------------------
    2820 
    2821 c      print*,'calcul des ratqs'
    2822 c   ratqs convectifs a l'ancienne en fonction de q(z=0)-q / q
    2823 c   ----------------
    2824 c   on ecrase le tableau ratqsc calcule par clouds_gno
    2825       if (iflag_cldcon.eq.1) then
    2826          do k=1,klev
    2827          do i=1,klon
    2828             if(ptconv(i,k)) then
    2829               ratqsc(i,k)=ratqsbas
    2830      s        +fact_cldcon*(q_seri(i,1)-q_seri(i,k))/q_seri(i,k)
    2831             else
    2832                ratqsc(i,k)=0.
    2833             endif
    2834          enddo
    2835          enddo
    2836 
    2837 c-----------------------------------------------------------------------
    2838 c  par nversion de la fonction log normale
    2839 c-----------------------------------------------------------------------
    2840       else if (iflag_cldcon.eq.4) then
    2841          ptconvth(:,:)=.false.
    2842          ratqsc(:,:)=0.
    2843          if(prt_level.ge.9) print*,'avant clouds_gno thermique'
    2844          call clouds_gno
    2845      s   (klon,klev,q_seri,zqsat,clwcon0th,ptconvth,ratqsc,rnebcon0th)
    2846          if(prt_level.ge.9) print*,' CLOUDS_GNO OK'
    2847        
    2848        endif
    2849 
    2850 c   ratqs stables
    2851 c   -------------
    2852 
    2853       if (iflag_ratqs.eq.0) then
    2854 
    2855 ! Le cas iflag_ratqs=0 correspond a la version IPCC 2005 du modele.
    2856          do k=1,klev
    2857             do i=1, klon
    2858                ratqss(i,k)=ratqsbas+(ratqshaut-ratqsbas)*
    2859      s         min((paprs(i,1)-pplay(i,k))/(paprs(i,1)-30000.),1.)
    2860             enddo
    2861          enddo
    2862 
    2863 ! Pour iflag_ratqs=1 ou 2, le ratqs est constant au dessus de
    2864 ! 300 hPa (ratqshaut), varie lineariement en fonction de la pression
    2865 ! entre 600 et 300 hPa et est soit constant (ratqsbas) pour iflag_ratqs=1
    2866 ! soit lineaire (entre 0 a la surface et ratqsbas) pour iflag_ratqs=2
    2867 ! Il s'agit de differents tests dans la phase de reglage du modele
    2868 ! avec thermiques.
    2869 
    2870       else if (iflag_ratqs.eq.1) then
    2871 
    2872          do k=1,klev
    2873             do i=1, klon
    2874                if (pplay(i,k).ge.60000.) then
    2875                   ratqss(i,k)=ratqsbas
    2876                else if ((pplay(i,k).ge.30000.).and.
    2877      s            (pplay(i,k).lt.60000.)) then
    2878                   ratqss(i,k)=ratqsbas+(ratqshaut-ratqsbas)*
    2879      s            (60000.-pplay(i,k))/(60000.-30000.)
    2880                else
    2881                   ratqss(i,k)=ratqshaut
    2882                endif
    2883             enddo
    2884          enddo
    2885 
    2886       else if (iflag_ratqs.eq.2) then
    2887 
    2888          do k=1,klev
    2889             do i=1, klon
    2890                if (pplay(i,k).ge.60000.) then
    2891                   ratqss(i,k)=ratqsbas
    2892      s            *(paprs(i,1)-pplay(i,k))/(paprs(i,1)-60000.)
    2893                else if ((pplay(i,k).ge.30000.).and.
    2894      s             (pplay(i,k).lt.60000.)) then
    2895                     ratqss(i,k)=ratqsbas+(ratqshaut-ratqsbas)*
    2896      s              (60000.-pplay(i,k))/(60000.-30000.)
    2897                else
    2898                     ratqss(i,k)=ratqshaut
    2899                endif
    2900             enddo
    2901          enddo
    2902 
    2903       else if (iflag_ratqs==3) then
    2904          do k=1,klev
    2905            ratqss(:,k)=ratqsbas+(ratqshaut-ratqsbas)
    2906      s     *min( ((paprs(:,1)-pplay(:,k))/70000.)**2 , 1. )
    2907          enddo
    2908 
    2909       else if (iflag_ratqs==4) then
    2910          do k=1,klev
    2911            ratqss(:,k)=ratqsbas+0.5*(ratqshaut-ratqsbas)
    2912      s     *( tanh( (50000.-pplay(:,k))/20000.) + 1.)
    2913          enddo
    2914 
    2915       endif
    2916 
    2917 
    2918 
    2919 
    2920 c  ratqs final
    2921 c  -----------
    2922 
    2923       if (iflag_cldcon.eq.1 .or.iflag_cldcon.eq.2
    2924      s    .or.iflag_cldcon.eq.4) then
    2925 
    2926 ! On ajoute une constante au ratqsc*2 pour tenir compte de
    2927 ! fluctuations turbulentes de petite echelle
    2928 
    2929          do k=1,klev
    2930             do i=1,klon
    2931                if ((fm_therm(i,k).gt.1.e-10)) then
    2932                   ratqsc(i,k)=sqrt(ratqsc(i,k)**2+0.05**2)
    2933                endif
    2934             enddo
    2935          enddo
    2936 
    2937 !   les ratqs sont une combinaison de ratqss et ratqsc
    2938        if(prt_level.ge.9)
    2939      $       write(lunout,*)'PHYLMD NOUVEAU TAU_RATQS ',tau_ratqs
    2940 
    2941          if (tau_ratqs>1.e-10) then
    2942             facteur=exp(-pdtphys/tau_ratqs)
    2943          else
    2944             facteur=0.
    2945          endif
    2946          ratqs(:,:)=ratqsc(:,:)*(1.-facteur)+ratqs(:,:)*facteur
    2947 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    2948 ! FH 22/09/2009
    2949 ! La ligne ci-dessous faisait osciller le modele et donnait une solution
    2950 ! assymptotique bidon et dépendant fortement du pas de temps.
    2951 !        ratqs(:,:)=sqrt(ratqs(:,:)**2+ratqss(:,:)**2)
    2952 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    2953          ratqs(:,:)=max(ratqs(:,:),ratqss(:,:))
    2954       else if (iflag_cldcon<=6) then
    2955 !   on ne prend que le ratqs stable pour fisrtilp
    2956          ratqs(:,:)=ratqss(:,:)
    2957       else
    2958           zfratqs1=exp(-pdtphys/10800.)
    2959           zfratqs2=exp(-pdtphys/10800.)
    2960 !         print*,'RAPPEL RATQS 1 ',zfratqs1,zfratqs2
    2961 !    s    ,ratqss(1,14),ratqs(1,14),ratqsc(1,14)
    2962           do k=1,klev
    2963              do i=1,klon
    2964                 if (ratqsc(i,k).gt.1.e-10) then
    2965                    ratqs(i,k)=ratqs(i,k)*zfratqs2
    2966      s             +(iflag_cldcon/100.)*ratqsc(i,k)*(1.-zfratqs2)
    2967                 endif
    2968                 ratqs(i,k)=min(ratqs(i,k)*zfratqs1
    2969      s          +ratqss(i,k)*(1.-zfratqs1),0.5)
    2970              enddo
    2971           enddo
    2972       endif
     2746! Computation of ratqs, the width (normalized) of the subrid scale
     2747! water distribution
     2748      CALL  calcratqs(klon,klev,prt_level,lunout,       
     2749     s     iflag_ratqs,iflag_con,iflag_cldcon,pdtphys,
     2750     s     ratqsbas,ratqshaut,tau_ratqs,fact_cldcon, 
     2751     s     ptconv,ptconvth,clwcon0th, rnebcon0th,   
     2752     s     paprs,pplay,q_seri,zqsat,fm_therm,
     2753     s     ratqs,ratqsc)
    29732754
    29742755
     
    38433624     I     cdragh,   coefh,     fm_therm, entr_therm,
    38443625     I     u1,       v1,        ftsol,    pctsrf,
     3626     I     ustar,     u10m,      v10m,
    38453627     I     rlat,     frac_impa, frac_nucl,rlon,
    38463628     I     presnivs, pphis,     pphi,     albsol1,
     
    39333715c
    39343716#include "calcul_STDlev.h"
    3935       twriteSTD(:,:,1)=tsumSTD(:,:,1)
    3936       qwriteSTD(:,:,1)=qsumSTD(:,:,1)
    3937       rhwriteSTD(:,:,1)=rhsumSTD(:,:,1)
    3938       phiwriteSTD(:,:,1)=phisumSTD(:,:,1)
    3939       uwriteSTD(:,:,1)=usumSTD(:,:,1)
    3940       vwriteSTD(:,:,1)=vsumSTD(:,:,1)
    3941       wwriteSTD(:,:,1)=wsumSTD(:,:,1)
    3942 
    3943       twriteSTD(:,:,2)=tsumSTD(:,:,2)
    3944       qwriteSTD(:,:,2)=qsumSTD(:,:,2)
    3945       rhwriteSTD(:,:,2)=rhsumSTD(:,:,2)
    3946       phiwriteSTD(:,:,2)=phisumSTD(:,:,2)
    3947       uwriteSTD(:,:,2)=usumSTD(:,:,2)
    3948       vwriteSTD(:,:,2)=vsumSTD(:,:,2)
    3949       wwriteSTD(:,:,2)=wsumSTD(:,:,2)
    3950 
    3951       twriteSTD(:,:,3)=tlevSTD(:,:)
    3952       qwriteSTD(:,:,3)=qlevSTD(:,:)
    3953       rhwriteSTD(:,:,3)=rhlevSTD(:,:)
    3954       phiwriteSTD(:,:,3)=philevSTD(:,:)
    3955       uwriteSTD(:,:,3)=ulevSTD(:,:)
    3956       vwriteSTD(:,:,3)=vlevSTD(:,:)
    3957       wwriteSTD(:,:,3)=wlevSTD(:,:)
    3958 
    3959       twriteSTD(:,:,4)=tlevSTD(:,:)
    3960       qwriteSTD(:,:,4)=qlevSTD(:,:)
    3961       rhwriteSTD(:,:,4)=rhlevSTD(:,:)
    3962       phiwriteSTD(:,:,4)=philevSTD(:,:)
    3963       uwriteSTD(:,:,4)=ulevSTD(:,:)
    3964       vwriteSTD(:,:,4)=vlevSTD(:,:)
    3965       wwriteSTD(:,:,4)=wlevSTD(:,:)
    3966 c
    3967 cIM initialisation 5eme fichier de sortie
    3968       twriteSTD(:,:,5)=tlevSTD(:,:)
    3969       qwriteSTD(:,:,5)=qlevSTD(:,:)
    3970       rhwriteSTD(:,:,5)=rhlevSTD(:,:)
    3971       phiwriteSTD(:,:,5)=philevSTD(:,:)
    3972       uwriteSTD(:,:,5)=ulevSTD(:,:)
    3973       vwriteSTD(:,:,5)=vlevSTD(:,:)
    3974       wwriteSTD(:,:,5)=wlevSTD(:,:)
    3975 c
    3976 cIM initialisation 6eme fichier de sortie
    3977       twriteSTD(:,:,6)=tlevSTD(:,:)
    3978       qwriteSTD(:,:,6)=qlevSTD(:,:)
    3979       rhwriteSTD(:,:,6)=rhlevSTD(:,:)
    3980       phiwriteSTD(:,:,6)=philevSTD(:,:)
    3981       uwriteSTD(:,:,6)=ulevSTD(:,:)
    3982       vwriteSTD(:,:,6)=vlevSTD(:,:)
    3983       wwriteSTD(:,:,6)=wlevSTD(:,:)
    3984 cIM for NMC files
    3985       DO n=1, nlevSTD3
    3986        DO k=1, nlevSTD
    3987         if(rlevSTD3(n).EQ.rlevSTD(k)) THEN
    3988          twriteSTD3(:,n)=tlevSTD(:,k)
    3989          qwriteSTD3(:,n)=qlevSTD(:,k)
    3990          rhwriteSTD3(:,n)=rhlevSTD(:,k)
    3991          phiwriteSTD3(:,n)=philevSTD(:,k)
    3992          uwriteSTD3(:,n)=ulevSTD(:,k)
    3993          vwriteSTD3(:,n)=vlevSTD(:,k)
    3994          wwriteSTD3(:,n)=wlevSTD(:,k)
    3995         endif !rlevSTD3(n).EQ.rlevSTD(k)
    3996        ENDDO
    3997       ENDDO
    3998 c
    3999       DO n=1, nlevSTD8
    4000        DO k=1, nlevSTD
    4001         if(rlevSTD8(n).EQ.rlevSTD(k)) THEN
    4002          tnondefSTD8(:,n)=tnondef(:,k,2)
    4003          twriteSTD8(:,n)=tsumSTD(:,k,2)
    4004          qwriteSTD8(:,n)=qsumSTD(:,k,2)
    4005          rhwriteSTD8(:,n)=rhsumSTD(:,k,2)
    4006          phiwriteSTD8(:,n)=phisumSTD(:,k,2)
    4007          uwriteSTD8(:,n)=usumSTD(:,k,2)
    4008          vwriteSTD8(:,n)=vsumSTD(:,k,2)
    4009          wwriteSTD8(:,n)=wsumSTD(:,k,2)
    4010         endif !rlevSTD8(n).EQ.rlevSTD(k)
    4011        ENDDO
    4012       ENDDO
    40133717c
    40143718c slp sea level pressure
  • LMDZ5/branches/testing/libf/phylmd/phytrac.F90

    r1665 r1707  
    88     cdragh,    coefh,    fm_therm, entr_therm,&
    99     yu1,       yv1,      ftsol,    pctsrf,    &
     10     ustar,     u10m,      v10m,               &
    1011     xlat,      frac_impa,frac_nucl,xlon,      &
    1112     presnivs,  pphis,    pphi,     albsol,    &
     
    119120!--------------
    120121!
    121   REAL,DIMENSION(klon),INTENT(IN)      :: cdragh ! coeff drag pour T et Q
    122   REAL,DIMENSION(klon,klev),INTENT(IN) :: coefh  ! coeff melange CL (m**2/s)
    123   REAL,DIMENSION(klon),INTENT(IN)      :: yu1    ! vents au premier niveau
    124   REAL,DIMENSION(klon),INTENT(IN)      :: yv1    ! vents au premier niveau
     122  REAL,DIMENSION(klon),INTENT(IN)     :: cdragh ! coeff drag pour T et Q
     123  REAL,DIMENSION(klon,klev),INTENT(IN):: coefh  ! coeff melange CL (m**2/s)
     124  REAL,DIMENSION(klon),INTENT(IN)     :: ustar,u10m,v10m ! u* & vent a 10m (m/s)
     125  REAL,DIMENSION(klon),INTENT(IN)     :: yu1    ! vents au premier niveau
     126  REAL,DIMENSION(klon),INTENT(IN)     :: yv1    ! vents au premier niveau
    125127!
    126128!Lessivage:
     
    244246     !    -- Traitement des traceurs avec traclmdz
    245247     CALL traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
    246           cdragh,  coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon, couchelimite, &
    247           sh, tr_seri, source, solsym, d_tr_cl, zmasse)
     248          cdragh,  coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon,couchelimite,sh,&
     249          rh, pphi, ustar, u10m, v10m, &
     250          tr_seri, source, solsym, d_tr_cl, zmasse)
    248251  CASE('inca')
    249252     !    -- CHIMIE INCA  config_inca = aero or chem --
  • LMDZ5/branches/testing/libf/phylmd/printflag.F

    r1403 r1707  
    8787!        radpas0  = NINT( 86400./tabcntr0(1)/INT( tabcntr0(6) ) )
    8888        PRINT 100
    89         PRINT 22, radpas0, radpas
     89!        PRINT 22, radpas0, radpas
    9090        PRINT 100
    9191       ENDIF
  • LMDZ5/branches/testing/libf/phylmd/traclmdz_mod.F90

    r1665 r1707  
    279279  SUBROUTINE traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
    280280       cdragh,  coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon, couchelimite, sh, &
     281       rh, pphi, ustar, zu10m, zv10m, &
    281282       tr_seri, source, solsym, d_tr_cl, zmasse)
    282283   
     
    315316!--------------
    316317!
    317     REAL,DIMENSION(klon),INTENT(IN)      :: cdragh     ! coeff drag pour T et Q
    318     REAL,DIMENSION(klon,klev),INTENT(IN) :: coefh      ! coeff melange CL (m**2/s)
    319     REAL,DIMENSION(klon),INTENT(IN)      :: yu1        ! vents au premier niveau
    320     REAL,DIMENSION(klon),INTENT(IN)      :: yv1        ! vents au premier niveau
     318    REAL,DIMENSION(klon),INTENT(IN)      :: cdragh  ! coeff drag pour T et Q
     319    REAL,DIMENSION(klon,klev),INTENT(IN) :: coefh   ! diffusivite turb (m**2/s)
     320    REAL,DIMENSION(klon),INTENT(IN)      :: yu1     ! vents au premier niveau
     321    REAL,DIMENSION(klon),INTENT(IN)      :: yv1     ! vents au premier niveau
    321322    LOGICAL,INTENT(IN)                   :: couchelimite
    322     REAL,DIMENSION(klon,klev),INTENT(IN) :: sh         ! humidite specifique
     323    REAL,DIMENSION(klon,klev),INTENT(IN) :: sh      ! humidite specifique
     324    REAL,DIMENSION(klon,klev),INTENT(IN) :: rh      ! Humidite relative
     325    REAL,DIMENSION(klon,klev),INTENT(IN) :: pphi    ! geopotentie
     326    REAL,DIMENSION(klon),INTENT(IN)      :: ustar   ! ustar (m/s)
     327    REAL,DIMENSION(klon),INTENT(IN)      :: zu10m   ! vent zonal 10m (m/s)
     328    REAL,DIMENSION(klon),INTENT(IN)      :: zv10m   ! vent zonal 10m (m/s)
    323329
    324330! Arguments necessaires pour les sources et puits de traceur:
Note: See TracChangeset for help on using the changeset viewer.