Changeset 2984 for trunk/LMDZ.MARS/libf


Ignore:
Timestamp:
Jun 26, 2023, 5:44:59 PM (19 months ago)
Author:
jnaar
Message:

Adaptative timestep working and computed directly in improvedclouds. Simpleclouds working again. JN

Location:
trunk/LMDZ.MARS/libf/phymars
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.MARS/libf/phymars/improvedclouds_mod.F

    r2966 r2984  
    55      CONTAINS
    66
    7       subroutine improvedclouds(microtimestep,
    8      &             pplay,pt,pq,
    9      &             subpdqcloud,subpdtcloud,
    10      &             nq,tauscaling,mmean)
     7      subroutine improvedclouds(ngrid,nlay,ptimestep,
     8     &             pplay,pt,pdt,pq,pdq,nq,tauscaling,
     9     &             imicro,zt,zq)
    1110      USE updaterad, ONLY: updaterice_micro, updaterccn
    1211      USE watersat_mod, ONLY: watersat
     
    1716     &                      igcm_hdo_vap,igcm_hdo_ice,
    1817     &                      qperemin
     18      use conc_mod, only: mmean
    1919      use comcstfi_h, only: pi, cpp
    2020      use write_output_mod, only: write_output
     
    4141c           T. Navarro, debug,correction, new scheme (October-April 2011)
    4242c           A. Spiga, optimization (February 2012)
    43 c           J. Naar, from global to local (no more ngrid, nlay) to allow
    44 c           different microphysical timesteps (May 2023)
     43c           J. Naar, adaptative subtimestep now done here (June 2023)
    4544c------------------------------------------------------------------
    4645#include "callkeys.h"
     
    4948c     Inputs/outputs:
    5049
     50      INTEGER, INTENT(IN) :: ngrid,nlay
    5151      INTEGER, INTENT(IN) :: nq                 ! nombre de traceurs
    52       REAL, INTENT(IN) :: microtimestep         ! pas de temps physique (s)
    53       REAL, INTENT(IN) :: pplay     ! pression au milieu des couches (Pa)           
    54       REAL, INTENT(IN) :: pt ! temperature at the middle of the
     52      REAL, INTENT(IN) :: ptimestep             ! pas de temps physique (s)
     53      REAL, INTENT(IN) :: pplay(ngrid,nlay)     ! pression au milieu des couches (Pa)           
     54      REAL, INTENT(IN) :: pt(ngrid,nlay)        ! temperature at the middle of the
    5555                                                ! layers (K)
    56                                                 !   param.
    57       REAL, INTENT(IN) :: pq(nq)  ! traceur (kg/kg)
    58                                                        !   (kg/kg.s-1)
    59       REAL, INTENT(IN) :: tauscaling     ! Convertion factor for qdust and Ndust
    60       REAL, INTENT(IN) :: mmean ! Mean atmospheric mass
    61      
    62       REAL, INTENT(OUT) :: subpdqcloud(nq)  ! tendance de la condensation
    63                                                        !   H2O(kg/kg.s-1)
    64       REAL, INTENT(OUT) :: subpdtcloud     ! tendance temperature due
    65                                                        !   a la chaleur latente
     56      REAL, INTENT(IN) :: pdt(ngrid,nlay)       ! temperature tendency (K/s)
     57      REAL, INTENT(IN) :: pq(ngrid,nlay,nq)     ! tracer (kg/kg)
     58      REAL, INTENT(IN) :: pdq(ngrid,nlay,nq)    ! tracer tendency (kg/kg/s)
     59      REAL, INTENT(IN) :: tauscaling(ngrid)     ! Convertion factor for qdust and Ndust
     60      INTEGER, INTENT(IN) :: imicro             ! nb. microphy calls(retrocompatibility)
     61     
     62      REAL, INTENT(OUT) :: zq(ngrid,nlay,nq)  ! tracers post microphy (kg/kg)
     63      REAL, INTENT(OUT) :: zt(ngrid,nlay)     ! temperature post microphy (K)
    6664
    6765c------------------------------------------------------------------
     
    7977      INTEGER ig,l,i
    8078     
    81       REAL zq(nq)  ! local value of tracers
    82       REAL zq0(nq) ! local initial value of tracers
    83       REAL zt       ! local value of temperature
    84       REAL zqsat_tmp(1)    ! saturation
    85       REAL zqsat    ! saturation
    86       REAL lw                         !Latent heat of sublimation (J.kg-1)
     79      REAL zqsat(ngrid,nlay)    ! saturation
     80      REAL lw                   !Latent heat of sublimation (J.kg-1)
    8781      REAL cste
    8882      REAL dMice           ! mass of condensed ice
    8983      REAL dMice_hdo       ! mass of condensed HDO ice
    90       REAL alpha ! HDO equilibrium fractionation coefficient (Saturation=1)
    91       REAL alpha_c ! HDO real fractionation coefficient
     84      REAL alpha(ngrid,nlay) ! HDO equilibrium fractionation coefficient (Saturation=1)
     85      REAL alpha_c(ngrid,nlay) ! HDO real fractionation coefficient
    9286!      REAL sumcheck
    9387      REAL*8 ph2o          ! Water vapor partial pressure (Pa)
     
    10599      REAL seq
    106100
    107       REAL rice      ! Ice mass mean radius (m)
     101      REAL rice(ngrid,nlay)      ! Ice mass mean radius (m)
    108102                                 ! (r_c in montmessin_2004)
    109       REAL rhocloud  ! Cloud density (kg.m-3)
    110       REAL rdust ! Dust geometric mean radius (m)
     103      REAL rhocloud(ngrid,nlay)  ! Cloud density (kg.m-3)
     104      REAL rdust(ngrid,nlay) ! Dust geometric mean radius (m)
    111105
    112106      REAL res      ! Resistance growth
     
    136130
    137131
     132c----------------------------------     
     133c JN : used in subtimestep
     134      REAL :: microtimestep! subdivision of ptimestep (s)
     135      REAL :: subpdtcloud(ngrid,nlay)  ! Temperature variation due to latent heat
     136      REAL :: zq0(ngrid,nlay,nq) ! local initial value of tracers
     137c      REAL zt0(ngrid,nlay,nq) ! local initial value of temperature
     138      INTEGER :: zimicro(ngrid,nlay) ! Subdivision of ptimestep
     139      INTEGER :: count_micro(ngrid,nlay) ! Number of microphys calls
     140      REAL :: zpotcond_inst(ngrid,nlay) ! condensable water at the beginning of physics (kg/kg)
     141      REAL :: zpotcond_full(ngrid,nlay) ! condensable water with integrated tendancies (kg/kg)
     142      REAL :: zpotcond(ngrid,nlay) ! maximal condensable water (previous two)
     143      REAL :: zqsat_tmp(1) ! maximal condensable water (previous two)
     144      REAL :: spenttime ! time spent in while loop
     145      REAL :: zdq ! used to compute adaptative timestep
     146      LOGICAL :: ending_ts ! Condition to end while loop
     147
    138148     
    139149c----------------------------------     
     
    148158
    149159
    150       REAL satubf,satuaf
    151       REAL res_out
     160      REAL satubf(ngrid,nlay),satuaf(ngrid,nlay)
     161      REAL res_out(ngrid,nlay)
    152162 
    153163
     
    234244! 1. Initialisation
    235245!=============================================================
    236       cste = 4*pi*rho_ice*microtimestep
    237 
    238       res_out = 0
    239       rice = 1.e-8
    240 
    241 c     Initialize the tendencies
    242       subpdqcloud(1:nq)=0
    243       subpdtcloud=0
    244    
    245 c     temperature and tracers previously incremented here     
    246 c     now done outside, in watercloud_mod.F     
    247 c      zt =  pteff + sum_subpdt * microtimestep
    248 c      zq(1:nq) = pqeff(1:nq) + sum_subpdq(1:nq) * microtimestep
    249      
    250       zq(1:nq)=pq(1:nq)
    251       zt=pt
    252      
    253       WHERE( zq(1:nq) < 1.e-30 )
    254      &       zq(1:nq) = 1.e-30
    255 
    256       zq0(1:nq) = zq(1:nq)
     246
     247      res_out(:,:) = 0
     248      rice(:,:) = 1.e-8
     249
     250c     Initialize the temperature, tracers
     251      zt(1:ngrid,1:nlay)=pt(1:ngrid,1:nlay)
     252      zq(1:ngrid,1:nlay,1:nq)=pq(1:ngrid,1:nlay,1:nq)
     253      subpdtcloud(1:ngrid,1:nlay)=0
     254     
     255      WHERE( zq(1:ngrid,1:nlay,1:nq) < 1.e-30 )
     256     &       zq(1:ngrid,1:nlay,1:nq) = 1.e-30
     257
    257258     
    258259!=============================================================
     
    263264      dev2 = 1. / ( sqrt(2.) * sigma_ice )
    264265
    265       call watersat(1,(/zt/),(/pplay/),zqsat_tmp)
    266       zqsat=zqsat_tmp(1)     
     266c     Compute the condensable amount of water vapor
     267c     or the sublimable water ice (if negative)
     268      call watersat(ngrid*nlay,zt+pdt*ptimestep,pplay,zqsat)
     269      zpotcond_full=(zq(:,:,igcm_h2o_vap)+
     270     &             pdq(:,:,igcm_h2o_vap)*ptimestep) - zqsat
     271      zpotcond_full=max(zpotcond_full,-zq(:,:,igcm_h2o_ice))
     272      call watersat(ngrid*nlay,zt,pplay,zqsat)
     273      zpotcond_inst=zq(:,:,igcm_h2o_vap) - zqsat
     274      zpotcond_inst=max(zpotcond_inst,-zq(:,:,igcm_h2o_ice))
     275c     zpotcond is the most strict criterion between the two
     276      DO l=1,nlay
     277        DO ig=1,ngrid
     278          if (zpotcond_full(ig,l).gt.0.) then
     279            zpotcond(ig,l)=max(zpotcond_inst(ig,l),zpotcond_full(ig,l))
     280          else if (zpotcond_full(ig,l).le.0.) then
     281            zpotcond(ig,l)=min(zpotcond_inst(ig,l),zpotcond_full(ig,l))
     282          endif! (zpotcond_full.gt.0.) then
     283        ENDDO !ig=1,ngrid
     284      ENDDO !l=1,nlay
     285           
    267286      countcells = 0
    268 
     287      zimicro(1:ngrid,1:nlay)=imicro
     288      count_micro(1:ngrid,1:nlay)=0
     289
     290c     Main loop over the GCM's grid
     291      DO l=1,nlay
     292        DO ig=1,ngrid
     293c       Subtimestep : here we go
     294        IF (cloud_adapt_ts) then
     295                call adapt_imicro(ptimestep,zpotcond(ig,l),
     296     &             zimicro(ig,l))
     297        ENDIF! (cloud_adapt_ts) then
     298        spenttime = 0.
     299        ending_ts=.false.
     300        print*, "ig,l :", ig, l
     301        DO while (.not.ending_ts)
     302          microtimestep=ptimestep/real(zimicro(ig,l))
     303          zq0(1:ngrid,1:nlay,1:nq) = zq
     304
     305          ! Check if we are integrating over ptimestep
     306          if (spenttime+microtimestep.ge.ptimestep) then
     307                  microtimestep=ptimestep-spenttime
     308          !       If so : last call !
     309                  ending_ts=.true.
     310          endif! (spenttime+microtimestep.ge.ptimestep) then
     311          call watersat(1,(/zt(ig,l)/),(/pplay(ig,l)/),zqsat_tmp)
     312          zqsat(ig,l)=zqsat_tmp(1)
    269313c       Get the partial pressure of water vapor and its saturation ratio
    270       ph2o = zq(igcm_h2o_vap) * (mmean/18.) * pplay
    271       satu = zq(igcm_h2o_vap) / zqsat
     314        ph2o = zq(ig,l,igcm_h2o_vap) * (mmean(ig,l)/18.) * pplay(ig,l)
     315        satu = zq(ig,l,igcm_h2o_vap) / zqsat(ig,l)
    272316
    273317!=============================================================
     
    275319!=============================================================
    276320
    277       IF ( satu .ge. 1. ) THEN         ! if there is condensation
    278 
    279        call updaterccn(zq(igcm_dust_mass),
    280      &         zq(igcm_dust_number),rdust,tauscaling)
     321       IF ( satu .ge. 1. ) THEN         ! if there is condensation
     322
     323        call updaterccn(zq(ig,l,igcm_dust_mass),
     324     &          zq(ig,l,igcm_dust_number),rdust(ig,l),tauscaling(ig))
    281325
    282326
    283327c       Expand the dust moments into a binned distribution
    284         Mo = zq(igcm_dust_mass)* tauscaling   + 1.e-30
    285         No = zq(igcm_dust_number)* tauscaling + 1.e-30
    286         Rn = rdust
     328        Mo = zq(ig,l,igcm_dust_mass)* tauscaling(ig)   + 1.e-30
     329        No = zq(ig,l,igcm_dust_number)* tauscaling(ig) + 1.e-30
     330        Rn = rdust(ig,l)
    287331        Rn = -log(Rn)
    288332        Rm = Rn - 3. * sigma_ice*sigma_ice 
     
    324368 
    325369c       Get the rates of nucleation
    326         call nuclea(ph2o,zt,satu,n_aer,rate)
     370        call nuclea(ph2o,zt(ig,l),satu,n_aer,rate)
    327371       
    328372        dN = 0.
     
    335379
    336380c       Update Dust particles
    337         zq(igcm_dust_mass)   =
    338      &  zq(igcm_dust_mass)   + dM/ tauscaling !max(tauscaling,1.e-10)
    339         zq(igcm_dust_number) =
    340      &  zq(igcm_dust_number) + dN/ tauscaling !max(tauscaling,1.e-10)
     381        zq(ig,l,igcm_dust_mass)   =
     382     &  zq(ig,l,igcm_dust_mass)   + dM/ tauscaling(ig) !max(tauscaling(ig),1.e-10)
     383        zq(ig,l,igcm_dust_number) =
     384     &  zq(ig,l,igcm_dust_number) + dN/ tauscaling(ig) !max(tauscaling(ig),1.e-10)
    341385c       Update CCNs
    342         zq(igcm_ccn_mass)   =
    343      &  zq(igcm_ccn_mass)   - dM/ tauscaling !max(tauscaling,1.e-10)
    344         zq(igcm_ccn_number) =
    345      &  zq(igcm_ccn_number) - dN/ tauscaling !max(tauscaling,1.e-10)
     386        zq(ig,l,igcm_ccn_mass)   =
     387     &  zq(ig,l,igcm_ccn_mass)   - dM/ tauscaling(ig) !max(tauscaling(ig),1.e-10)
     388        zq(ig,l,igcm_ccn_number) =
     389     &  zq(ig,l,igcm_ccn_number) - dN/ tauscaling(ig) !max(tauscaling(ig),1.e-10)
    346390       
    347        ENDIF ! of is satu >1
     391        ENDIF ! of is satu >1
    348392
    349393!=============================================================
     
    355399c to avoid unrealistic value for nuclei radius and so on for cases that remain negligible.
    356400
    357        IF ( zq(igcm_ccn_number)*tauscaling.ge. 1.) THEN ! we trigger crystal growth
     401       IF ( zq(ig,l,igcm_ccn_number)*tauscaling(ig).ge. 1.) THEN ! we trigger crystal growth
    358402
    359403     
    360         call updaterice_micro(zq(igcm_h2o_ice),
    361      &    zq(igcm_ccn_mass),zq(igcm_ccn_number),
    362      &    tauscaling,rice,rhocloud)
    363 
    364         No   = zq(igcm_ccn_number)* tauscaling + 1.e-30
     404        call updaterice_micro(zq(ig,l,igcm_h2o_ice),
     405     &    zq(ig,l,igcm_ccn_mass),zq(ig,l,igcm_ccn_number),
     406     &    tauscaling(ig),rice(ig,l),rhocloud(ig,l))
     407
     408        No   = zq(ig,l,igcm_ccn_number)* tauscaling(ig) + 1.e-30
    365409
    366410c       saturation at equilibrium
    367411c       rice should not be too small, otherwise seq value is not valid
    368         seq  = exp(2.*sig(zt)*mh2o / (rho_ice*rgp*zt*
    369      &             max(rice,1.e-7)))
     412        seq  = exp(2.*sig(zt(ig,l))*mh2o / (rho_ice*rgp*zt(ig,l)*
     413     &             max(rice(ig,l),1.e-7)))
    370414       
    371415c       get resistance growth
    372         call growthrate(zt,pplay,
    373      &             real(ph2o/satu),rice,res,Dv)
    374 
    375         res_out = res
     416        call growthrate(zt(ig,l),pplay(ig,l),
     417     &             real(ph2o/satu),rice(ig,l),res,Dv)
     418
     419        res_out(ig,l) = res
    376420
    377421ccccccc  implicit scheme of mass growth
     422c         cste here must be computed at each step
     423          cste = 4*pi*rho_ice*microtimestep
     424c          print*, 'ig',ig
     425c          print*, 'l',l
     426c          print*, 'cste',cste
     427c          print*, 'seq',seq
     428c          print*, 'zqsat',zqsat(ig,l)
     429c          print*, 'zq vap',zq(ig,l,igcm_h2o_vap)
     430c          print*, 'res',res
     431c          print*, 'No',No
     432c          print*, 'rice',rice(ig,l)
    378433
    379434        dMice =
    380      & (zq(igcm_h2o_vap)-seq*zqsat)
    381      &   /(res*zqsat/(cste*No*rice) + 1.)
     435     & (zq(ig,l,igcm_h2o_vap)-seq*zqsat(ig,l))
     436     &   /(res*zqsat(ig,l)/(cste*No*rice(ig,l)) + 1.)
    382437
    383438
    384439! With the above scheme, dMice cannot be bigger than vapor,
    385440! but can be bigger than all available ice.
    386        dMice = max(dMice,-zq(igcm_h2o_ice))
    387        dMice = min(dMice,zq(igcm_h2o_vap)) ! this should be useless...
    388 
    389        zq(igcm_h2o_ice) = zq(igcm_h2o_ice)+dMice
    390        zq(igcm_h2o_vap) = zq(igcm_h2o_vap)-dMice
     441       dMice = max(dMice,-zq(ig,l,igcm_h2o_ice))
     442       dMice = min(dMice,zq(ig,l,igcm_h2o_vap)) ! this should be useless...
     443
     444       zq(ig,l,igcm_h2o_ice) = zq(ig,l,igcm_h2o_ice)+dMice
     445       zq(ig,l,igcm_h2o_vap) = zq(ig,l,igcm_h2o_vap)-dMice
    391446
    392447
     
    394449       
    395450       ! latent heat release
    396        lw=(2834.3-0.28*(zt-To)-
    397      &     0.004*(zt-To)*(zt-To))*1.e+3
    398        subpdtcloud= dMice*lw/cpp/microtimestep
     451       lw=(2834.3-0.28*(zt(ig,l)-To)-
     452     &     0.004*(zt(ig,l)-To)*(zt(ig,l)-To))*1.e+3
     453       subpdtcloud(ig,l)= dMice*lw/cpp/microtimestep
    399454         
    400455c Special case of the isotope of water HDO   
     
    405460           if (hdofrac) then
    406461             !! Calculation of the HDO vapor coefficient
    407              Dv_hdo = 1./3. * sqrt( 8*kbz*zt/(pi*mhdo/nav) )
    408      &          * kbz * zt /
    409      &          ( pi * pplay * (molco2+molhdo)*(molco2+molhdo)
     462             Dv_hdo = 1./3. * sqrt( 8*kbz*zt(ig,l)/(pi*mhdo/nav) )
     463     &          * kbz * zt(ig,l) /
     464     &          ( pi * pplay(ig,l) * (molco2+molhdo)*(molco2+molhdo)
    410465     &          * sqrt(1.+mhdo/mco2) )
    411466             !! Calculation of the fractionnation coefficient at equilibrium
    412 c             alpha = exp(16288./zt**2.-9.34d-2)  ! Merlivat and Nief et al. 1967
    413              alpha = exp(13525./zt**2.-5.59d-2)  ! Lamb et al. 2017
     467             alpha(ig,l) = exp(16288./zt(ig,l)**2.-9.34d-2)
     468c             alpha = exp(13525./zt(ig,l)**2.-5.59d-2)  !Lamb
    414469             !! Calculation of the 'real' fractionnation coefficient
    415              alpha_c = (alpha*satu)/
    416      &          ( (alpha*(Dv/Dv_hdo)*(satu-1.)) + 1.)
    417 c             alpha_c = alpha ! to test without the effect of cinetics
     470             alpha_c(ig,l) = (alpha(ig,l)*satu)/
     471     &          ( (alpha(ig,l)*(Dv/Dv_hdo)*(satu-1.)) + 1.)
     472c             alpha_c(ig,l) = alpha(ig,l) ! to test without the effect of cinetics
    418473           else
    419              alpha_c = 1.d0
     474             alpha_c(ig,l) = 1.d0
    420475           endif
    421            if (zq0(igcm_h2o_vap).gt.qperemin) then
     476           if (zq0(ig,l,igcm_h2o_vap).gt.qperemin) then
    422477              dMice_hdo=
    423      &          dMice*alpha_c*
    424      &     ( zq0(igcm_hdo_vap)
    425      &             /zq0(igcm_h2o_vap) )
     478     &          dMice*alpha_c(ig,l)*
     479     &     ( zq0(ig,l,igcm_hdo_vap)
     480     &             /zq0(ig,l,igcm_h2o_vap) )
    426481           else
    427482             dMice_hdo=0.
     
    429484         !! sublimation
    430485         else
    431            if (zq0(igcm_h2o_ice).gt.qperemin) then
     486           if (zq0(ig,l,igcm_h2o_ice).gt.qperemin) then
    432487             dMice_hdo=
    433488     &            dMice*
    434      &     ( zq0(igcm_hdo_ice)
    435      &             /zq0(igcm_h2o_ice) )
     489     &     ( zq0(ig,l,igcm_hdo_ice)
     490     &             /zq0(ig,l,igcm_h2o_ice) )
    436491           else
    437492             dMice_hdo=0.
     
    439494         endif !if (dMice.gt.0.0)
    440495
    441        dMice_hdo = max(dMice_hdo,-zq(igcm_hdo_ice))
    442        dMice_hdo = min(dMice_hdo,zq(igcm_hdo_vap))
    443 
    444        zq(igcm_hdo_ice) = zq(igcm_hdo_ice)+dMice_hdo
    445        zq(igcm_hdo_vap) = zq(igcm_hdo_vap)-dMice_hdo
     496       dMice_hdo = max(dMice_hdo,-zq(ig,l,igcm_hdo_ice))
     497       dMice_hdo = min(dMice_hdo,zq(ig,l,igcm_hdo_vap))
     498
     499       zq(ig,l,igcm_hdo_ice) = zq(ig,l,igcm_hdo_ice)+dMice_hdo
     500       zq(ig,l,igcm_hdo_vap) = zq(ig,l,igcm_hdo_vap)-dMice_hdo
    446501
    447502       endif ! if (hdo)       
     
    453508c         If all the ice particles sublimate, all the condensation
    454509c         nuclei are released:
    455           if (zq(igcm_h2o_ice).le.1.e-28) then
     510          if (zq(ig,l,igcm_h2o_ice).le.1.e-28) then
    456511         
    457512c           Water
    458             zq(igcm_h2o_vap) = zq(igcm_h2o_vap)
    459      &                            + zq(igcm_h2o_ice)
    460             zq(igcm_h2o_ice) = 0.
     513            zq(ig,l,igcm_h2o_vap) = zq(ig,l,igcm_h2o_vap)
     514     &                            + zq(ig,l,igcm_h2o_ice)
     515            zq(ig,l,igcm_h2o_ice) = 0.
    461516            if (hdo) then
    462               zq(igcm_hdo_vap) = zq(igcm_hdo_vap)
    463      &                            + zq(igcm_hdo_ice)
    464               zq(igcm_hdo_ice) = 0.
     517              zq(ig,l,igcm_hdo_vap) = zq(ig,l,igcm_hdo_vap)
     518     &                            + zq(ig,l,igcm_hdo_ice)
     519              zq(ig,l,igcm_hdo_ice) = 0.
    465520            endif
    466521c           Dust particles
    467             zq(igcm_dust_mass) = zq(igcm_dust_mass)
    468      &                              + zq(igcm_ccn_mass)
    469             zq(igcm_dust_number) = zq(igcm_dust_number)
    470      &                                + zq(igcm_ccn_number)
     522            zq(ig,l,igcm_dust_mass) = zq(ig,l,igcm_dust_mass)
     523     &                              + zq(ig,l,igcm_ccn_mass)
     524            zq(ig,l,igcm_dust_number) = zq(ig,l,igcm_dust_number)
     525     &                                + zq(ig,l,igcm_ccn_number)
    471526c           CCNs
    472             zq(igcm_ccn_mass) = 0.
    473             zq(igcm_ccn_number) = 0.
     527            zq(ig,l,igcm_ccn_mass) = 0.
     528            zq(ig,l,igcm_ccn_number) = 0.
    474529
    475530          endif
    476531         
    477532          ENDIF !of if Nccn>1
     533
    478534         
    479      
    480       ! Get cloud tendencies
    481         subpdqcloud(igcm_h2o_vap) =
    482      &   (zq(igcm_h2o_vap) -
    483      &    zq0(igcm_h2o_vap))/microtimestep
    484         subpdqcloud(igcm_h2o_ice) =
    485      &   (zq(igcm_h2o_ice) -
    486      &    zq0(igcm_h2o_ice))/microtimestep
    487         if (hdo) then
    488           subpdqcloud(igcm_hdo_vap) =
    489      &     (zq(igcm_hdo_vap) -
    490      &      zq0(igcm_hdo_vap))/microtimestep
    491           subpdqcloud(igcm_hdo_ice) =
    492      &     (zq(igcm_hdo_ice) -
    493      &      zq0(igcm_hdo_ice))/microtimestep
    494         endif
    495         subpdqcloud(igcm_ccn_mass) =
    496      &   (zq(igcm_ccn_mass) -
    497      &    zq0(igcm_ccn_mass))/microtimestep
    498         subpdqcloud(igcm_ccn_number) =
    499      &   (zq(igcm_ccn_number) -
    500      &    zq0(igcm_ccn_number))/microtimestep
    501      
    502       if (scavenging) then
    503      
    504         subpdqcloud(igcm_dust_mass) =
    505      &   (zq(igcm_dust_mass) -
    506      &    zq0(igcm_dust_mass))/microtimestep
    507         subpdqcloud(igcm_dust_number) =
    508      &   (zq(igcm_dust_number) -
    509      &    zq0(igcm_dust_number))/microtimestep
    510          
    511       endif
     535      ! No more getting tendency : we increment tracers & temp directly
     536
     537      ! Increment tracers & temp for following microtimestep
     538      ! Dust :
     539      ! Special treatment of dust_mass / number for scavenging ?
     540            IF (scavenging) THEN
     541              zq(ig,l,igcm_dust_mass) =zq(ig,l,igcm_dust_mass)+
     542     &         pdq(ig,l,igcm_dust_mass)*microtimestep
     543              zq(ig,l,igcm_dust_number) =zq(ig,l,igcm_dust_number)+
     544     &         pdq(ig,l,igcm_dust_number)*microtimestep
     545            ELSE
     546              zq(ig,l,igcm_dust_mass) =zq0(ig,l,igcm_dust_mass)
     547              zq(ig,l,igcm_dust_number) =zq0(ig,l,igcm_dust_number)
     548            ENDIF !(scavenging) THEN
     549              zq(ig,l,igcm_ccn_mass) = zq(ig,l,igcm_ccn_mass) +
     550     &         pdq(ig,l,igcm_ccn_mass)*microtimestep
     551              zq(ig,l,igcm_ccn_number) = zq(ig,l,igcm_ccn_number) +
     552     &          pdq(ig,l,igcm_ccn_number)*microtimestep
     553
     554      ! Water :
     555            zq(ig,l,igcm_h2o_ice) = zq(ig,l,igcm_h2o_ice)+
     556     &         pdq(ig,l,igcm_h2o_ice)*microtimestep
     557            zq(ig,l,igcm_h2o_vap) = zq(ig,l,igcm_h2o_vap)+
     558     &         pdq(ig,l,igcm_h2o_vap)*microtimestep
     559
     560            ! HDO (if computed) :
     561            IF (hdo) THEN
     562            zq(ig,l,igcm_hdo_ice) =
     563     &       zq(ig,l,igcm_hdo_ice)+
     564     &         pdq(ig,l,igcm_hdo_ice)*microtimestep
     565            zq(ig,l,igcm_hdo_vap) =
     566     &       zq(ig,l,igcm_hdo_vap)+
     567     &         pdq(ig,l,igcm_hdo_vap)*microtimestep
     568            ENDIF ! hdo
     569
     570c  Could also set subpdtcloud to 0 if not activice to make it simpler
     571c  or change name of the flag
     572            IF (.not.activice) THEN
     573                    subpdtcloud(ig,l)=0.
     574            ENDIF
     575      ! Temperature :
     576            zt(ig,l) = zt(ig,l)+(pdt(ig,l)+
     577     &          subpdtcloud(ig,l))*microtimestep
     578
     579c         Prevent negative tracers
     580          WHERE( zq(1:ngrid,1:nlay,1:nq) < 1.e-30 )
     581     &       zq(1:ngrid,1:nlay,1:nq) = 1.e-30
     582
     583c         Increment time spent in here
     584          spenttime=spenttime+microtimestep
     585          count_micro(ig,l)=count_micro(ig,l)+1
     586          IF (cloud_adapt_ts) then
     587c           Estimation of how much is actually condensing/subliming
     588c                print*, 'zq',zq(ig,l,igcm_h2o_vap)
     589c                print*, 'zq0',zq0(ig,l,igcm_h2o_vap)
     590c                print*, 'ptimestep',ptimestep
     591                zdq=(zq(ig,l,igcm_h2o_vap)-zq0(ig,l,igcm_h2o_vap))
     592     &                   *ptimestep/microtimestep
     593c           Estimation of how much is actually condensing/subliming
     594c                zdq=min(abs(dMice)*ptimestep,zpotcond(ig,l))
     595                call adapt_imicro(ptimestep,zdq,
     596     &             zimicro(ig,l))
     597          ENDIF! (cloud_adapt_ts) then
     598          ENDDO ! while (.not. ending_ts)
     599        ENDDO ! of ig loop
     600      ENDDO ! of nlayer loop
     601     
     602     
     603c------ Useful outputs to check how it went
     604      call write_output("zpotcond_inst","zpotcond_inst "//
     605     &   "microphysics","(kg/kg)",zpotcond_inst(:,:))
     606      call write_output("zpotcond_full","zpotcond_full "//
     607     &   "microphysics","(kg/kg)",zpotcond_full(:,:))
     608      call write_output("zpotcond","zpotcond "//
     609     &   "microphysics","(kg/kg)",zpotcond(:,:))
     610      call write_output("count_micro","count_micro "//
     611     &   "after microphysics","integer",count_micro(:,:))
    512612     
    513613     
     
    519619     
    520620!       error2d(:) = 0.
    521 !         error2d(ig) = max(abs(error_out),error2d(ig))
    522         satubf = zq0(igcm_h2o_vap)/zqsat
    523         satuaf = zq(igcm_h2o_vap)/zqsat
    524 
    525 c      print*, 'count is ',countcells, ' i.e. ',
    526 c     &     countcells*100/(nlay*ngrid), '% for microphys computation'
     621       DO l=1,nlay
     622       DO ig=1,ngrid
     623!         error2d(ig) = max(abs(error_out(ig,l)),error2d(ig))
     624          satubf(ig,l) = zq0(ig,l,igcm_h2o_vap)/zqsat(ig,l)
     625          satuaf(ig,l) = zq(ig,l,igcm_h2o_vap)/zqsat(ig,l)
     626       ENDDO
     627       ENDDO
     628
     629      print*, 'count is ',countcells, ' i.e. ',
     630     &     countcells*100/(nlay*ngrid), '% for microphys computation'
    527631
    528632#ifndef MESOSCALE
    529633!      IF (ngrid.ne.1) THEN ! 3D
    530 !         call write_output("satu","ratio saturation","",
     634!         call WRITEDIAGFI(ngrid,"satu","ratio saturation","",3,
    531635!     &                    satu_out)
    532 !         call write_output("dM","ccn variation","kg/kg",
     636!         call WRITEDIAGFI(ngrid,"dM","ccn variation","kg/kg",3,
    533637!     &                    dM_out)
    534 !         call write_output("dN","ccn variation","#",
     638!         call WRITEDIAGFI(ngrid,"dN","ccn variation","#",3,
    535639!     &                    dN_out)
    536 !         call write_output("error","dichotomy max error","%",
     640!         call WRITEDIAGFI(ngrid,"error","dichotomy max error","%",2,
    537641!     &                    error2d)
    538 !         call write_output("zqsat","zqsat","kg",
     642!         call WRITEDIAGFI(ngrid,"zqsat","zqsat","kg",3,
    539643!     &                    zqsat)
    540644!      ENDIF
    541645
    542646!      IF (ngrid.eq.1) THEN ! 1D
    543 !         call write_output("error","incertitude sur glace","%",
     647!         call WRITEDIAGFI(ngrid,"error","incertitude sur glace","%",1,
    544648!     &                    error_out)
    545          call write_output("resist","resistance","s/m2",
     649         call WRITEdiagfi(ngrid,"resist","resistance","s/m2",1,
    546650     &                    res_out)
    547          call write_output("satu_bf","satu before","kg/kg",
     651         call WRITEdiagfi(ngrid,"satu_bf","satu before","kg/kg",1,
    548652     &                    satubf)
    549          call write_output("satu_af","satu after","kg/kg",
     653         call WRITEdiagfi(ngrid,"satu_af","satu after","kg/kg",1,
    550654     &                    satuaf)
    551          call write_output("vapbf","h2ovap before","kg/kg",
    552      &                    zq0(igcm_h2o_vap))
    553          call write_output("vapaf","h2ovap after","kg/kg",
    554      &                    zq(igcm_h2o_vap))
    555          call write_output("icebf","h2oice before","kg/kg",
    556      &                    zq0(igcm_h2o_ice))
    557          call write_output("iceaf","h2oice after","kg/kg",
    558      &                    zq(igcm_h2o_ice))
    559          call write_output("ccnbf","ccn before","/kg",
    560      &                    zq0(igcm_ccn_number))
    561          call write_output("ccnaf","ccn after","/kg",
    562      &                    zq(igcm_ccn_number))
    563 c         call write_output("growthrate","growth rate","m^2/s",
     655         call WRITEdiagfi(ngrid,"vapbf","h2ovap before","kg/kg",1,
     656     &                    zq0(1,1,igcm_h2o_vap))
     657         call WRITEdiagfi(ngrid,"vapaf","h2ovap after","kg/kg",1,
     658     &                    zq(1,1,igcm_h2o_vap))
     659         call WRITEdiagfi(ngrid,"icebf","h2oice before","kg/kg",1,
     660     &                    zq0(1,1,igcm_h2o_ice))
     661         call WRITEdiagfi(ngrid,"iceaf","h2oice after","kg/kg",1,
     662     &                    zq(1,1,igcm_h2o_ice))
     663         call WRITEdiagfi(ngrid,"ccnbf","ccn before","/kg",1,
     664     &                    zq0(1,1,igcm_ccn_number))
     665         call WRITEdiagfi(ngrid,"ccnaf","ccn after","/kg",1,
     666     &                    zq(1,1,igcm_ccn_number))
     667c         call WRITEDIAGFI(ngrid,"growthrate","growth rate","m^2/s",1,
    564668c     &                    gr_out)
    565 c         call write_output("nuclearate","nucleation rate","",
     669c         call WRITEDIAGFI(ngrid,"nuclearate","nucleation rate","",1,
    566670c     &                    rate_out)
    567 c         call write_output("dM","ccn variation","kg",
     671c         call WRITEDIAGFI(ngrid,"dM","ccn variation","kg",1,
    568672c     &                    dM_out)
    569 c         call write_output("dN","ccn variation","#",
     673c         call WRITEDIAGFI(ngrid,"dN","ccn variation","#",1,
    570674c     &                    dN_out)
    571          call write_output("zqsat","p vap sat","kg/kg",
     675         call WRITEdiagfi(ngrid,"zqsat","p vap sat","kg/kg",1,
    572676     &                    zqsat)
    573 !         call write_output("satu","ratio saturation","",
    574 !     &                    satu_out(:,:))
    575          call write_output("rice","ice radius","m",
     677!         call WRITEDIAGFI(ngrid,"satu","ratio saturation","",1,
     678!     &                    satu_out)
     679         call WRITEdiagfi(ngrid,"rice","ice radius","m",1,
    576680     &                    rice)
    577 !         call write_output("rdust_sca","rdust","m",
     681!         call WRITEDIAGFI(ngrid,"rdust_sca","rdust","m",1,
    578682!     &                    rdust)
    579 !         call write_output("rsedcloud","rsedcloud","m",
     683!         call WRITEDIAGFI(ngrid,"rsedcloud","rsedcloud","m",1,
    580684!     &                    rsedcloud)
    581 !         call write_output("rhocloud","rhocloud","kg.m-3",
     685!         call WRITEDIAGFI(ngrid,"rhocloud","rhocloud","kg.m-3",1,
    582686!     &                    rhocloud)
    583687!      ENDIF
     
    636740     
    637741      END SUBROUTINE improvedclouds
     742
     743      SUBROUTINE adapt_imicro(ptimestep,potcond,
     744     $                     zimicro)
     745
     746c Adaptative timestep for water ice clouds.
     747c Works using a powerlaw to compute the minimal duration of subtimestep
     748c (in s) should all the avalaible vapor (resp. ice) condenses (resp.sublimates)
     749c Then, we use the instantaneous vap (ice) gradient extrapolated to the
     750c rest of duration to increase subtimestep duration, for computing
     751c efficiency.
     752
     753      real,intent(in) :: ptimestep ! total duration of physics (sec)
     754      real,intent(in) :: potcond ! condensible vapor / sublimable ice(kg/kg)
     755      real :: alpha, beta ! Coefficients for power law
     756      integer,intent(out) :: zimicro ! number of ptimestep division
     757
     758c       zimicro = 30
     759c      Coefficients good enough for present-day Mars :
     760       alpha = 1.87485684e+09
     761       beta = 1.45655856e+00
     762c      Coefficients covering high obliquity scenarios :
     763c       alpha=3.36711332e+15
     764c       beta=1.98636414e+00
     765c       nimicro=min(max(alpha*abs(potcond)**beta,5.),7000.)
     766c       zimicro=ceiling((ptimestep/timeleft)*nimicro)
     767c       zimicro=ceiling((timeleft/ptimestep)*nimicro)
     768       zimicro=ceiling(min(max(alpha*abs(potcond)**beta,5.),7000.))
     769
     770      END SUBROUTINE adapt_imicro
    638771     
    639772      END MODULE improvedclouds_mod
  • trunk/LMDZ.MARS/libf/phymars/watercloud_mod.F

    r2966 r2984  
    4747c  2004 - 2012
    4848c
    49 c     2023: J. Naar, adding different subtimestep on each grid cell
    50 c          plus not doing microphysics if no water present
    51 c          plus simpleclouds no longer in the loop for microphysics
     49c     2023: J. Naar, now with adaptative timestep for improvedclouds
     50c          (done in improvedclouds_mod).
    5251c=======================================================================
    5352
     
    146145! Scheme for adaptative timestep J. Naar 2023
    147146c      LOGICAL :: computed_micro(ngrid,nlay) ! Check if microphy was done in this cell
    148       REAL :: computed_micro(ngrid,nlay) ! Check if microphy was done inthis cell (logical)
     147      REAL :: count_micro(ngrid,nlay) ! Initially computed microtimestep
    149148      REAL :: zt_micro(ngrid,nlay) ! Temperature during microphysics (K)
    150149      REAL :: zq_micro(ngrid,nlay,nq) ! Tracers during microphysics (kg/kg)
     
    155154      REAL :: zpotcond(ngrid,nlay) ! maximal condensable water, used to
    156155compute adaptative subdivision of ptimestep
     156      REAL :: spenttime ! timespent
     157      REAL :: zdq ! used to compute adaptative timestep
    157158
    158159
     
    298299      END IF ! end if (CLFvarying)
    299300c------------------------------------------------------------------
    300 c Time subsampling for microphysics
     301c Cloud physics (nucleation, condensation / sublimation)
    301302c------------------------------------------------------------------
    302303      rhocloud(1:ngrid,1:nlay) = rho_dust
    303304
    304 
    305 c     Initialisation of all the stuff JN2023
    306 c      computed_micro(1:ngrid,1:nlay)=.false.
    307       computed_micro(1:ngrid,1:nlay)=0.
     305c     Initialisation of all the stuff (JN,2023)
    308306      zt_micro(:,:)=pt(:,:)
    309307      zq_micro(:,:,:)=pq(:,:,:)
    310       zq_micro(:,:,:)=pq(:,:,:)
    311       call watersat(ngrid*nlay,zt_micro,pplay,zqsat_micro)
    312       zpotcond_inst=zq_micro(:,:,igcm_h2o_vap) - zqsat_micro
    313       call watersat(ngrid*nlay,zt_micro+pdt*ptimestep,pplay,zqsat_micro)
    314       zpotcond_full=(zq_micro(:,:,igcm_h2o_vap)+
    315      &             pdq(:,:,igcm_h2o_vap)*ptimestep) - zqsat_micro
    316       zimicro(1:ngrid,1:nlay)=imicro
    317       if (cloud_adapt_ts) then
    318               call adapt_imicro(ptimestep,zpotcond(ig,l),
    319      $                   zimicro(ig,l))
    320       endif! (cloud_adapt_ts) then
    321       DO l=1,nlay
    322         DO ig=1,ngrid
    323 c         Start by computing the condensable water vapor amount
    324           if (zpotcond_full(ig,l).gt.0.) then
    325             zpotcond(ig,l)=max(zpotcond_inst(ig,l),zpotcond_full(ig,l))
    326           else if (zpotcond_full(ig,l).le.0.) then
    327             zpotcond(ig,l)=min(zpotcond_inst(ig,l),zpotcond_full(ig,l))
    328           endif! (zpotcond_full.gt.0.) then
    329           microtimestep=ptimestep/real(zimicro(ig,l))
    330 c         Check if microphysics is even needed, that is if enough action
    331 c         is happening water-wise
    332           if ((pq(ig,l,igcm_h2o_ice)+pdq(ig,l,igcm_h2o_ice)*ptimestep
    333      &      .gt.1e-22) .or. (abs(zpotcond(ig,l)).gt.1e-22)) then
    334 c         Eventuellement sortir simpleclouds de la boucle egalement
    335           computed_micro(ig,l)=1.
    336           DO microstep=1,zimicro(ig,l)
    337      
    338 
    339 c JN : incrementing after main microphysics scheme
    340 c Previously we were incrementing tendencies, we now
    341 c increment tracers and temperature directly
    342 c We are thus starting at the end of the first iteration
    343 c
     308
    344309c-------------------------------------------------------------------
    345310c   1.  Main call to the different cloud schemes:
    346311c------------------------------------------------
    347         IF (microphys) THEN
    348            CALL improvedclouds(microtimestep,
    349      &          pplay(ig,l),zt_micro(ig,l),
    350      &          zq_micro(ig,l,:),subpdqcloud(ig,l,:),
    351      &          subpdtcloud(ig,l),nq,tauscaling(ig),mmean(ig,l))
    352 
    353         ELSE
    354 c Simpleclouds should maybe be taken out and put in a specific loop ?
    355            CALL simpleclouds(ngrid,nlay,microtimestep,
     312c ds.
     313      IF (microphys) THEN
     314           CALL improvedclouds(ngrid,nlay,ptimestep,
     315     &          pplay,pt,pdt,pq,pdq,nq,tauscaling,imicro,
     316     &          zt_micro,zq_micro)
     317
     318      ELSE
     319
     320c Specific loop for simpleclouds.
     321       DO l=1,nlay
     322         DO ig=1,ngrid
     323           CALL simpleclouds(ngrid,nlay,ptimestep,
    356324     &             pplay,pzlay,pteff,sum_subpdt,
    357325     &             pqeff,sum_subpdq,subpdqcloud,subpdtcloud,
    358326     &             nq,tau,rice)
    359         ENDIF
    360 
    361327c-------------------------------------------------------------------
    362328c   2.  Updating tracers and temperature after cloud scheme:
     329c   For improved clouds (with microphysics) this is done directly
     330c   in the microphysics, during the subtimestep
     331c   I put it like that to be retrocompatible (JN)
    363332c-----------------------------------------------
    364333
    365         IF (microphys) THEN
    366               zq_micro(ig,l,igcm_dust_mass) =
    367      &         zq_micro(ig,l,igcm_dust_mass)+(pdq(ig,l,igcm_dust_mass)
    368      &         +subpdqcloud(ig,l,igcm_dust_mass))*microtimestep
    369               zq_micro(ig,l,igcm_dust_number) =
    370      &         zq_micro(ig,l,igcm_dust_number)
    371      &         +(pdq(ig,l,igcm_dust_number)
    372      &         + subpdqcloud(ig,l,igcm_dust_number))*microtimestep
    373               zq_micro(ig,l,igcm_ccn_mass) =
    374      &         zq_micro(ig,l,igcm_ccn_mass) +
    375      &         (pdq(ig,l,igcm_ccn_mass)
    376      &         +subpdqcloud(ig,l,igcm_ccn_mass))*microtimestep
    377               zq_micro(ig,l,igcm_ccn_number) =
    378      &          zq_micro(ig,l,igcm_ccn_number) +
    379      &         (pdq(ig,l,igcm_ccn_number)
    380      &          + subpdqcloud(ig,l,igcm_ccn_number))*microtimestep
    381         ENDIF
    382334            zq_micro(ig,l,igcm_h2o_ice) =
    383335     &       zq_micro(ig,l,igcm_h2o_ice)+
    384336     &         (pdq(ig,l,igcm_h2o_ice)
    385      &        + subpdqcloud(ig,l,igcm_h2o_ice))*microtimestep
     337     &        + subpdqcloud(ig,l,igcm_h2o_ice))*ptimestep
    386338            zq_micro(ig,l,igcm_h2o_vap) =
    387339     &       zq_micro(ig,l,igcm_h2o_vap)+
    388340     &         (pdq(ig,l,igcm_h2o_vap)
    389      &        + subpdqcloud(ig,l,igcm_h2o_vap))*microtimestep
     341     &        + subpdqcloud(ig,l,igcm_h2o_vap))*ptimestep
    390342
    391343            IF (hdo) THEN
     
    393345     &       zq_micro(ig,l,igcm_hdo_ice)+
    394346     &         (pdq(ig,l,igcm_hdo_ice)
    395      &        + subpdqcloud(ig,l,igcm_hdo_ice))*microtimestep
     347     &        + subpdqcloud(ig,l,igcm_hdo_ice))*ptimestep
    396348            zq_micro(ig,l,igcm_hdo_vap) =
    397349     &       zq_micro(ig,l,igcm_hdo_vap)+
    398350     &         (pdq(ig,l,igcm_hdo_vap)
    399      &        + subpdqcloud(ig,l,igcm_hdo_vap))*microtimestep
     351     &        + subpdqcloud(ig,l,igcm_hdo_vap))*ptimestep
    400352            ENDIF ! hdo
    401353
    402354c  Could also set subpdtcloud to 0 if not activice to make it simpler
    403             zt_micro(ig,l) = zt_micro(ig,l)+
    404      &           pdt(ig,l)*microtimestep
    405         IF (activice) THEN
    406               zt_micro(ig,l) = zt_micro(ig,l)+
    407      &           subpdtcloud(ig,l)*microtimestep
    408         ENDIF
    409 !      !! Example of how to use writediagmicrofi useful to
    410 !      !! get outputs at each microphysical sub-timestep (better to be used in 1D)
    411 !            CALL WRITEDIAGMICROFI(ngrid,imicro,microstep,
    412 !     &       microtimestep,'subpdtcloud',
    413 !     &      'subpdtcloud','K/s',1,subpdtcloud(:,:))     
    414  
    415           ENDDO ! of DO microstep=1,imicro
    416           endif! (zq(ig,l,igcm_h2o_ice)+pdq(ig,l,igcm_h2o_ice)*ptimestep
    417 !     &      .gt.1e-22).or.(abs(zpotcond).gt.1e-22) then
    418         ENDDO ! ig=1,ngrid
    419       ENDDO ! l=1,nlay
    420 
    421      
    422 c------ Useful outputs to check how it went
    423       call write_output("computed_micro","computed_micro "//
    424      &   "after microphysics","logical",computed_micro(:,:))
    425       call write_output("zimicro","Used number of subtimestep "//
    426      &   "in cloud microphysics","integer",real(zimicro(:,:)))
     355c  or change name of the flag
     356            IF (activice) THEN
     357                zt_micro(ig,l) = zt_micro(ig,l)+
     358     &             subpdtcloud(ig,l)*ptimestep
     359            ENDIF
     360
     361         ENDDO !ig=1,ngrid
     362       ENDDO !l=1,nlay
     363      ENDIF
     364
     365
     366     
    427367c-------------------------------------------------------------------
    428368c   3.  Compute final tendencies after time loop:
     
    786726
    787727
    788        SUBROUTINE adapt_imicro(ptimestep,potcond,
    789      $                     zimicro)
    790 
    791 c Pas de temps adaptatif pour les nuages
    792 
    793       real,intent(in) :: ptimestep ! total duration of physics (sec)
    794       real,intent(in) :: potcond ! total duration of physics (sec)
    795       real :: alpha, beta ! total duration of physics (sec)
    796       integer,intent(out) :: zimicro ! number of ptimestep division
    797 
    798 c       zimicro = ptimestep*alpha*potcond**beta
    799        zimicro = 30
    800 
    801        END SUBROUTINE adapt_imicro
    802 
    803 
    804728      END MODULE watercloud_mod
Note: See TracChangeset for help on using the changeset viewer.