Ignore:
Timestamp:
May 19, 2017, 11:19:17 AM (8 years ago)
Author:
mvals
Message:

Mars GCM

  • implementation of a sub-grid water cloud fraction scheme (by A. Pottier): If CLFvarying is set to true in callphys.def, then the sub-grid cloud fraction routine is applied in watercloud.F and aeropacity.F.
  • accordingly modified files are: watercloud.F, aeropacity.F, callradite.F, conf_phys.F, phyetat0.F, phyredem.F90, physiq_mod.F, testphys1d.F, callkeys.h, newtsart.F, start2archive.F, lect_start_archive.F
  • added file: tcondwater.F90, used by watercloud.F to calculate the condensation temperature of water
  • watercloud.F, aeropacity.F, callradite.F are converted to module files as watercloud_mod.F, aeropacity_mod.F, callradite_mod.F

MV

Location:
trunk/LMDZ.MARS/libf/phymars
Files:
1 added
6 edited
3 moved

Legend:

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

    r1710 r1711  
     1      MODULE aeropacity_mod
     2
     3      IMPLICIT NONE
     4
     5      CONTAINS
     6
    17      SUBROUTINE aeropacity(ngrid,nlayer,nq,zday,pplay,pplev,ls,
    28     &    pq,tauscaling,tauref,tau,taucloudtes,aerosol,dsodust,reffrad,
    3      &    nueffrad,QREFvis3d,QREFir3d,omegaREFvis3d,omegaREFir3d)
     9     &    nueffrad,QREFvis3d,QREFir3d,omegaREFvis3d,omegaREFir3d,
     10     &    clearsky,totcloudfrac)
    411                                                   
    512! to use  'getin'
     
    8491      REAL omegaREFvis3d(ngrid,nlayer,naerkind)
    8592      REAL omegaREFir3d(ngrid,nlayer,naerkind)
     93      real,intent(in) :: totcloudfrac(ngrid) ! total cloud fraction
     94      logical,intent(in) :: clearsky ! true for part without clouds,
     95                                     ! false for part with clouds (total or sub-grid clouds)
     96      real :: CLFtot ! total cloud fraction
    8697c
    8798c    Local variables :
     
    417428        taucloudtes(1:ngrid) = 0.
    418429c       2. Opacity calculation
    419         DO ig=1, ngrid
    420           DO l=1,nlayer
    421             aerosol(ig,l,iaer) = max(1E-20,
    422      &        (  0.75 * QREFvis3d(ig,l,iaer) /
    423      &        ( rho_ice * reffrad(ig,l,iaer) )  ) *
    424      &        pq(ig,l,i_ice) *
    425      &        ( pplev(ig,l) - pplev(ig,l+1) ) / g
     430        ! NO CLOUDS
     431        IF (clearsky) THEN
     432            aerosol(1:ngrid,1:nlayer,iaer) =1.e-9
     433        ! CLOUDSs
     434        ELSE ! else (clearsky)
     435          DO ig=1, ngrid
     436            DO l=1,nlayer
     437              aerosol(ig,l,iaer) = max(1E-20,
     438     &          (  0.75 * QREFvis3d(ig,l,iaer) /
     439     &          ( rho_ice * reffrad(ig,l,iaer) )  ) *
     440     &          pq(ig,l,i_ice) *
     441     &          ( pplev(ig,l) - pplev(ig,l+1) ) / g
    426442     &                              )
    427             taucloudvis(ig) = taucloudvis(ig) + aerosol(ig,l,iaer)
    428             taucloudtes(ig) = taucloudtes(ig) + aerosol(ig,l,iaer)*
    429      &        QREFir3d(ig,l,iaer) / QREFvis3d(ig,l,iaer) *
    430      &        ( 1.E0 - omegaREFir3d(ig,l,iaer) )
    431           ENDDO
    432         ENDDO
     443              taucloudvis(ig) = taucloudvis(ig) + aerosol(ig,l,iaer)
     444              taucloudtes(ig) = taucloudtes(ig) + aerosol(ig,l,iaer)*
     445     &          QREFir3d(ig,l,iaer) / QREFvis3d(ig,l,iaer) *
     446     &          ( 1.E0 - omegaREFir3d(ig,l,iaer) )
     447            ENDDO
     448          ENDDO
     449          ! SUB-GRID SCALE CLOUDS
     450          IF (CLFvarying) THEN
     451             DO ig=1, ngrid
     452                DO l=1,nlayer-1
     453                   CLFtot  = max(totcloudfrac(ig),0.01)
     454                   aerosol(ig,l,iaer)=
     455     &                    aerosol(ig,l,iaer)/CLFtot
     456                   aerosol(ig,l,iaer) =
     457     &                    max(aerosol(ig,l,iaer),1.e-9)
     458                ENDDO
     459             ENDDO
     460!          ELSE ! else (CLFvarying)
     461!             DO ig=1, ngrid
     462!                DO l=1,nlayer-1 ! to stop the rad tran bug
     463!                   CLFtot  = CLFfixval
     464!                   aerosol(ig,l,iaer)=
     465!     &                    aerosol(ig,l,iaer)/CLFtot
     466!                   aerosol(ig,l,iaer) =
     467!     &                    max(aerosol(ig,l,iaer),1.e-9)
     468!                ENDDO
     469!             ENDDO
     470          ENDIF ! end (CLFvarying)             
     471        ENDIF ! end (clearsky)
     472
    433473c       3. Outputs -- Now done in physiq.F
    434474!        IF (ngrid.NE.1) THEN
     
    726766c      ENDIF ! of IF (ngrid.NE.1)
    727767c -----------------------------------------------------------------
    728       return
    729       end
     768
     769      END SUBROUTINE aeropacity
     770     
     771      END MODULE aeropacity_mod
  • trunk/LMDZ.MARS/libf/phymars/callkeys.h

    r1684 r1711  
    1414     &   ,activice,water,tifeedback,microphys,supersat,caps,photochem   &
    1515     &   ,calltherm,callrichsl,callslope,tituscap,callyamada4,co2clouds,&
    16      &   microphysco2
     16     &   microphysco2,CLFvarying
    1717     
    1818      COMMON/callkeys_i/iradia,iaervar,iddist,ilwd,ilwb,ilwn,ncouche    &
     
    2020     
    2121      COMMON/callkeys_r/topdustref,semi,alphan,euveff,                  &
    22      &   tke_heat_flux,dustrefir,fixed_euv_value
     22     &   tke_heat_flux,dustrefir,fixed_euv_value,CLFfixval
    2323     
    2424      LOGICAL callrad,calldifv,calladj,callcond,callsoil,               &
     
    4141      real euveff
    4242      real tke_heat_flux
     43      real CLFfixval
    4344
    4445      integer iddist
     
    5960      logical sedimentation
    6061      logical activice,tifeedback,supersat,caps
     62      logical CLFvarying
    6163      logical co2clouds
    6264      logical water
  • trunk/LMDZ.MARS/libf/phymars/callradite_mod.F

    r1710 r1711  
     1      MODULE callradite_mod
     2
     3      IMPLICIT NONE
     4
     5      CONTAINS
     6
    17      SUBROUTINE callradite(icount,ngrid,nlayer,nq,zday,ls,pq,albedo,
    28     $     emis,mu0,pplev,pplay,pt,tsurf,fract,dist_sol,igout,
    39     $     dtlw,dtsw,fluxsurf_lw,fluxsurf_sw,fluxtop_lw,fluxtop_sw,
    410     &     tauref,tau,aerosol,dsodust,tauscaling,taucloudtes,rdust,rice,
    5      &     nuice,co2ice)
    6 
     11     &     nuice,co2ice,clearsky,totcloudfrac)
     12
     13      use aeropacity_mod
    714      use dimradmars_mod, only: ndomainsz, nflev, nsun, nir
    815      use dimradmars_mod, only: naerkind, name_iaer,
     
    186193      REAL,INTENT(OUT) :: nuice(ngrid,nlayer)  ! Estimated effective variance
    187194      REAL,INTENT(IN) :: co2ice(ngrid)           ! co2 ice surface layer (kg.m-2)
     195c     sub-grid scale water ice clouds
     196      real,intent(out) :: totcloudfrac(ngrid)
     197      logical,intent(in) :: clearsky
    188198
    189199c
     
    376386      CALL aeropacity(ngrid,nlayer,nq,zday,pplay,pplev,ls,
    377387     &     pq,tauscaling,tauref,tau,taucloudtes,aerosol,dsodust,reffrad,
    378      &     nueffrad,QREFvis3d,QREFir3d,omegaREFvis3d,omegaREFir3d)
     388     &     nueffrad,QREFvis3d,QREFir3d,omegaREFvis3d,omegaREFir3d,
     389     &     clearsky,totcloudfrac)
    379390
    380391c     Starting loop on sub-domain
     
    545556
    546557
    547       return
    548       end
     558      END SUBROUTINE callradite
     559
     560      END MODULE callradite_mod
  • trunk/LMDZ.MARS/libf/phymars/conf_phys.F

    r1684 r1711  
    438438         write(*,*) " water = ",water
    439439
     440! sub-grid cloud fraction: fixed
     441         write(*,*) "Fixed cloud fraction?"
     442         CLFfixval=1.0 ! default value
     443         call getin("CLFfixval",CLFfixval)
     444         write(*,*) "CLFfixval=",CLFfixval
     445! sub-grid cloud fraction: varying
     446         write(*,*) "Use partial nebulosity?"
     447         CLFvarying=.false. ! default value
     448         call getin("CLFvarying",CLFvarying)
     449         write(*,*)"CLFvarying=",CLFvarying
     450
    440451!CO2 clouds scheme?
    441452         write(*,*) "Compute CO2 clouds ?"
  • trunk/LMDZ.MARS/libf/phymars/dyn1d/testphys1d.F

    r1621 r1711  
    118118      character(len=44) :: txt
    119119
     120!MV17 nuages alizee
     121!     included by AP14 for cloud fraction computations- see BC generic
     122      !real,allocatable,dimension(:,:),save :: cloudfrac
     123      real,allocatable,dimension(:),save :: totcloudfrac
    120124c=======================================================================
    121125
     
    225229        allocate(dqdyn(nlayer,nq))
    226230        allocate(mqtot(nq))
     231!MV17 nuages alizee
     232        !allocate(cloudfrac(ngrid,nlayermx))!essai allocate AP14
     233        allocate(totcloudfrac(ngrid))
    227234       
    228235        ! read tracer names from file traceur.def
     
    711718      call physdem1("startfi.nc",nsoilmx,ngrid,llm,nq,
    712719     .              dtphys,time,
    713      .              tsurf,tsoil,co2ice,emis,q2,qsurf,tauscaling)
     720     .              tsurf,tsoil,co2ice,emis,q2,qsurf,tauscaling,
     721     .              totcloudfrac) !MV17 nuages alizee
    714722
    715723c=======================================================================
  • trunk/LMDZ.MARS/libf/phymars/phyetat0.F90

    r1621 r1711  
    11subroutine phyetat0 (fichnom,tab0,Lmodif,nsoil,ngrid,nlay,nq, &
    22                     day_ini,time0,tsurf,tsoil,emis,q2,qsurf,co2ice, &
    3                      tauscaling)
     3                     tauscaling,totcloudfrac)
    44!  use netcdf
    55  use tracer_mod, only: noms ! tracer names
     
    4646  real,intent(out) :: co2ice(ngrid) ! co2 ice cover
    4747  real,intent(out) :: tauscaling(ngrid) ! dust conversion factor
    48 
     48  real,intent(out) :: totcloudfrac(ngrid) ! total cloud fraction
    4949!======================================================================
    5050!  Local variables:
     
    241241endif
    242242
     243! Sub-grid cloud fraction
     244call get_field("totcloudfrac",totcloudfrac,found,indextime)
     245if (.not.found) then
     246  write(*,*) "phyetat0: <totcloudfrac> not in file WARNING put to 1"
     247  totcloudfrac(:) = 1.0 !valeur par defaut (CLFfixval par defaut)
     248else
     249  write(*,*) "phyetat0: total cloud fraction <totcloudfrac> range:", &
     250             minval(totcloudfrac), maxval(totcloudfrac)
     251endif
    243252
    244253! Surface temperature :
  • trunk/LMDZ.MARS/libf/phymars/phyredem.F90

    r1621 r1711  
    145145subroutine physdem1(filename,nsoil,ngrid,nlay,nq, &
    146146                    phystep,time,tsurf,tsoil,co2ice,emis,q2,qsurf,&
    147                     tauscaling)
     147                    tauscaling,totcloudfrac)
    148148  ! write time-dependent variable to restart file
    149149  use iostart, only : open_restartphy, close_restartphy, &
     
    165165  real,intent(in) :: qsurf(ngrid,nq)
    166166  real,intent(in) :: tauscaling(ngrid)
     167  real, intent(in) :: totcloudfrac(ngrid)
    167168 
    168169  integer :: iq
     
    194195  ! Planetary Boundary Layer
    195196  call put_field("q2","pbl wind variance",q2,time)
    196  
     197
     198  ! Sub-grid cloud fraction
     199  call put_field("totcloudfrac","Total cloud fraction",totcloudfrac,time)
     200 
    197201  ! Dust conversion factor
    198202  ! Only to be read by newstart to convert to actual dust values
  • trunk/LMDZ.MARS/libf/phymars/physiq_mod.F

    r1708 r1711  
    1414     $            ,pdu,pdv,pdt,pdq,pdpsrf)
    1515
     16      use watercloud_mod
     17      use aeropacity_mod
     18      use callradite_mod
    1619      use tracer_mod, only: noms, mmol, igcm_co2, igcm_n2, igcm_co2_ice,
    1720     &                      igcm_co, igcm_o, igcm_h2o_vap, igcm_h2o_ice,
     
    386389!      REAL zu2(ngrid)
    387390
     391c sub-grid scale water ice clouds (A. Pottier 2013)
     392      logical clearsky
     393      ! flux for the part without clouds
     394      real zdtsw1(ngrid,nlayer)
     395      real zdtlw1(ngrid,nlayer)
     396      real fluxsurf_lw1(ngrid)     
     397      real fluxsurf_sw1(ngrid,2) 
     398      real fluxtop_lw1(ngrid)
     399      real fluxtop_sw1(ngrid,2)
     400      REAL taucloudtes1(ngrid)
     401      ! tf: fraction of clouds, ntf: fraction without clouds
     402      real tf, ntf
     403      real,allocatable,dimension(:),save :: totcloudfrac ! total cloud fraction over the column
     404      REAL rave2(ngrid), totrave2(ngrid) ! Mean water ice mean radius (m)
     405
    388406c=======================================================================
    389407
     
    394412c  ---------------------------------------
    395413      IF (firstcall) THEN
     414        ALLOCATE(totcloudfrac(ngrid)) !A. Pottier: souci avec ca a priori, tentative
     415                                      ! de les allouer dans dimradmars_mod.F90
     416        !ouverture de nebuldata.out et nettoyage du precedent
     417         open(144,file="nebuldata.out")
     418         close(144,status="delete")
    396419
    397420c        variables set to 0
     
    412435     &         nsoilmx,ngrid,nlayer,nq,
    413436     &         day_ini,time_phys,
    414      &         tsurf,tsoil,emis,q2,qsurf,co2ice,tauscaling)
     437     &         tsurf,tsoil,emis,q2,qsurf,co2ice,tauscaling,
     438     &         totcloudfrac)
    415439
    416440         if (pday.ne.day_ini) then
     
    552576         call solarlong(float(day_ini),zls)
    553577      end if
     578      !A. Pottier: checking nebuldata.out
     579      open(14,file="nebuldata.out",position="append")
     580      write(14,*) "LsEtJour ",zls*180/pi,zday
     581      close(14)
    554582
    555583c     Initialize pressure levels
     
    699727c          Radiative transfer
    700728c          ------------------
    701 
     729           ! callradite for the part with clouds
     730           clearsky=.false. ! part with clouds for both cases CLFvarying true/false
    702731           CALL callradite(icount,ngrid,nlayer,nq,zday,zls,pq,albedo,
    703732     $     emis,mu0,zplev,zplay,pt,tsurf,fract,dist_sol,igout,
    704733     $     zdtlw,zdtsw,fluxsurf_lw,fluxsurf_sw,fluxtop_lw,fluxtop_sw,
    705734     $     tauref,tau,aerosol,dsodust,tauscaling,taucloudtes,rdust,rice,
    706      $     nuice,co2ice)
    707 
     735     $     nuice,co2ice,clearsky,totcloudfrac)
     736           ! case of sub-grid water ice clouds: callradite for the clear case
     737            IF (CLFvarying) THEN
     738               ! ---> PROBLEMS WITH ALLOCATED ARRAYS
     739               ! (temporary solution in callcorrk: do not deallocate
     740               ! if
     741               ! CLFvarying ...) ?? AP ??
     742               clearsky=.true. !
     743               CALL callradite(icount,ngrid,nlayer,nq,zday,zls,pq,
     744     &              albedo,emis,mu0,zplev,zplay,pt,tsurf,fract,
     745     &              dist_sol,igout,zdtlw1,zdtsw1,fluxsurf_lw1,
     746     &              fluxsurf_sw1,fluxtop_lw1,fluxtop_sw1,tauref,tau,
     747     &              aerosol,dsodust,tauscaling,taucloudtes1,rdust,
     748     &              rice,nuice,co2ice, clearsky, totcloudfrac)
     749               clearsky = .false.  ! just in case.
     750               ! Sum the fluxes and heating rates from cloudy/clear
     751               ! cases
     752               DO ig=1,ngrid
     753                  tf=totcloudfrac(ig)
     754                  ntf=1.-tf
     755                  fluxsurf_lw(ig) = ntf*fluxsurf_lw1(ig)
     756     &                                          + tf*fluxsurf_lw(ig)
     757                  fluxsurf_sw(ig,1) = ntf*fluxsurf_sw1(ig,1)
     758     &                                          + tf*fluxsurf_sw(ig,1)
     759                  fluxsurf_sw(ig,2) = ntf*fluxsurf_sw1(ig,2)
     760     &                                          + tf*fluxsurf_sw(ig,2)
     761                  fluxtop_lw(ig)  = ntf*fluxtop_lw1(ig)
     762     &                                          + tf*fluxtop_lw(ig)
     763                  fluxtop_sw(ig,1)  = ntf*fluxtop_sw1(ig,1) 
     764     &                                          + tf*fluxtop_sw(ig,1)
     765                  fluxtop_sw(ig,2)  = ntf*fluxtop_sw1(ig,2) 
     766     &                                          + tf*fluxtop_sw(ig,2)
     767                  taucloudtes(ig) = ntf*taucloudtes1(ig)
     768     &                                          + tf*taucloudtes(ig)
     769                  zdtlw(ig,1:nlayer) = ntf*zdtlw1(ig,1:nlayer)
     770     &                                         + tf*zdtlw(ig,1:nlayer)
     771                  zdtsw(ig,1:nlayer) = ntf*zdtsw1(ig,1:nlayer)
     772     &                                         + tf*zdtsw(ig,1:nlayer)
     773               ENDDO
     774
     775            ENDIF ! (CLFvarying)
    708776
    709777#ifdef DUSTSTORM
     
    10821150     &                pq,pdq,zdqcloud,zdtcloud,
    10831151     &                nq,tau,tauscaling,rdust,rice,nuice,
    1084      &                rsedcloud,rhocloud)
     1152     &                rsedcloud,rhocloud,totcloudfrac)
    10851153     
    10861154c Temperature variation due to latent heat release
     
    17401808          call physdem1("restartfi.nc",nsoilmx,ngrid,nlayer,nq,
    17411809     .                ptimestep,ztime_fin,
    1742      .                tsurf,tsoil,co2ice,emis,q2,qsurf,tauscaling)
     1810     .                tsurf,tsoil,co2ice,emis,q2,qsurf,tauscaling,
     1811     .               totcloudfrac)
    17431812         
    17441813         ENDIF
     
    18441913               enddo
    18451914             endif ! of if (scavenging)
     1915
     1916           !Alternative A. Pottier weighting
     1917           rave2(:) = 0.
     1918           totrave2(:) = 0.
     1919           do ig=1,ngrid
     1920              do l=1,nlayer
     1921              rave2(ig) =rave2(ig)+ zq(ig,l,igcm_h2o_ice)*rice(ig,l)
     1922              totrave2(ig) = totrave2(ig) + zq(ig,l,igcm_h2o_ice)
     1923              end do
     1924              rave2(ig)=max(rave2(ig)/max(totrave2(ig),1.e-30),1.e-30)
     1925           end do
    18461926
    18471927           endif ! of if (water)
     
    23882468     &                       'Mean reff',
    23892469     &                       'm',2,rave)
     2470!A. Pottier
     2471             CALL WRITEDIAGFI(ngrid,'rmoym',
     2472     &                      'alternative reffice',
     2473     &                      'm',2,rave2)
    23902474            call WRITEDIAGFI(ngrid,'saturation',
    23912475     &           'h2o vap saturation ratio','dimless',3,satu)
     
    24142498     &                       'J.s-1/2.m-2.K-1',3,inertiesoil)
    24152499              endif
     2500!A. Pottier
     2501              if (CLFvarying) then !AP14 nebulosity
     2502                 call WRITEDIAGFI(ngrid,'totcloudfrac',
     2503     &                'Total cloud fraction',
     2504     &                       ' ',2,totcloudfrac)
     2505              endif !clf varying
     2506
    24162507           endif !(water)
    24172508
     
    28852976     &                      'm',0,rave)
    28862977
     2978          !Alternative A. Pottier weighting
     2979           rave2 = 0.
     2980           totrave2 = 0.
     2981           do l=1,nlayer
     2982              rave2 =rave2+ zq(1,l,igcm_h2o_ice)*rice(1,l)
     2983              totrave2 = totrave2 + zq(1,l,igcm_h2o_ice)
     2984           end do
     2985           rave2=max(rave2/max(totrave2,1.e-30),1.e-30)
     2986          CALL WRITEDIAGFI(ngrid,'rmoym',
     2987     &                      'reffice',
     2988     &                      'm',0,rave2)
     2989
    28872990           do iq=1,nq
    28882991             call WRITEDIAGFI(ngrid,trim(noms(iq))//'_s',
     
    28983001     &                          +zdqcloud(1,:,igcm_h2o_vap))
    28993002     
    2900           call WRITEDIAGFI(ngrid,"rice","ice radius","m",1,
     3003        call WRITEDIAGFI(ngrid,"rice","ice radius","m",1,
    29013004     &                    rice)
     3005             
     3006              if (CLFvarying) then
     3007                 call WRITEDIAGFI(ngrid,'totcloudfrac',
     3008     &                'Total cloud fraction',
     3009     &                       ' ',0,totcloudfrac)
     3010              endif !clfvarying
     3011
    29023012         ENDIF ! of IF (water)
    29033013         
  • trunk/LMDZ.MARS/libf/phymars/watercloud_mod.F

    r1710 r1711  
     1       MODULE watercloud_mod
     2
     3       IMPLICIT NONE
     4
     5       CONTAINS
     6
    17       SUBROUTINE watercloud(ngrid,nlay,ptimestep,
    28     &                pplev,pplay,pdpsrf,pzlay,pt,pdt,
    39     &                pq,pdq,pdqcloud,pdtcloud,
    410     &                nq,tau,tauscaling,rdust,rice,nuice,
    5      &                rsedcloud,rhocloud)
     11     &                rsedcloud,rhocloud,totcloudfrac)
    612! to use  'getin'
    713      USE ioipsl_getincom
     
    7278      real rhocloud(ngrid,nlay)  ! Cloud density (kg.m-3)
    7379
     80      REAL, INTENT(INOUT):: totcloudfrac(ngrid) ! Cloud fraction (A. Pottier 2013)
    7481c   local:
    7582c   ------
     
    101108      LOGICAL,SAVE :: firstcall=.true.
    102109
     110! Representation of sub-grid water ice clouds A. Pottier 2013
     111!      REAL :: zt(ngrid, nlay)
     112      REAL :: zq(ngrid, nlay,nq)
     113      REAL :: zdelt 
     114      REAL :: norm
     115      REAL :: ponder
     116      REAL :: tcond(ngrid,nlay)
     117      REAL ::  zqvap(ngrid,nlay)
     118      REAL :: zqice(ngrid,nlay)
     119      REAL ::  spant ! delta T for the temperature distribution
     120!      REAL :: zqsat(ngrid,nlay) ! saturation
     121      REAL :: zteff(ngrid, nlay)! effective temperature in the cloud,neb
     122      REAL :: pqeff(ngrid, nlay, nq)! effective tracers quantities in the cloud
     123      REAL :: cloudfrac(ngrid,nlay) ! cloud fraction
     124      REAL :: mincloud ! min cloud frac
     125      LOGICAL, save :: flagcloud=.true.
    103126c    ** un petit test de coherence
    104127c       --------------------------
     
    138161      rhocloud(1:ngrid,1:nlay) = rho_dust
    139162
    140 
     163c-------------------------------------------------------------------
     164c   0.  Representation of sub-grid water ice clouds
     165c------------------
     166c-----Tendencies
     167      DO l=1,nlay
     168       DO ig=1,ngrid
     169          zt(ig,l)=pt(ig,l)+ pdt(ig,l)*ptimestep
     170       ENDDO
     171      ENDDO
     172      DO l=1,nlay
     173        DO ig=1,ngrid
     174          DO iq=1,nq
     175             zq(ig,l,iq)=pq(ig,l,iq)+pdq(ig,l,iq)*ptimestep
     176          ENDDO
     177        ENDDO
     178      ENDDO
     179c-----Effective temperature calculation
     180      IF (CLFvarying) THEN
     181         spant=3.0 ! delta T for the temprature distribution
     182         mincloud=0.1 ! min cloudfrac when there is ice
     183         IF (flagcloud) THEN
     184             write(*,*) "Delta T", spant
     185             write(*,*) "mincloud", mincloud
     186             flagcloud=.false.
     187         END IF
     188         CALL watersat(ngrid*nlay,zt,pplay,zqsat)
     189         zqvap=zq(:,:,igcm_h2o_vap)
     190         zqice=zq(:,:,igcm_h2o_ice)
     191         CALL tcondwater(ngrid*nlay,pplay,zqvap+zqice,tcond)
     192         DO l=1,nlay
     193           DO ig=1,ngrid
     194              zdelt=spant !MAX(spant*zt(ig,l),1.e-12), now totally in K. Fixed width
     195              IF (tcond(ig,l) .ge. (zt(ig,l)+zdelt)) THEN
     196                 zteff(ig,l)=zt(ig,l)
     197                 cloudfrac(ig,l)=1.
     198              ELSE IF (tcond(ig,l) .le. (zt(ig,l)-zdelt)) THEN
     199                 zteff(ig,l)=zt(ig,l)-zdelt
     200                 cloudfrac(ig,l)=mincloud
     201              ELSE
     202                 cloudfrac(ig,l)=(tcond(ig,l)-zt(ig,l)+zdelt)/
     203     &                           (2.0*zdelt)
     204                 zteff(ig,l)=(tcond(ig,l)+zt(ig,l)-zdelt)/2.
     205              END IF
     206              zteff(ig,l)=zteff(ig,l)-pdt(ig,l)*ptimestep
     207              IF (cloudfrac(ig,l).le.0) THEN
     208                 cloudfrac(ig,l)=mincloud
     209              ELSE IF (cloudfrac(ig,l).gt.1) THEN
     210                 cloudfrac(ig,l)=1.
     211              END IF
     212           ENDDO
     213         ENDDO
     214c-----Calculation of the total cloud coverage of the column
     215         DO ig=1,ngrid
     216            totcloudfrac(ig) = 0.
     217            norm=0.
     218            DO l=1,nlay
     219               ponder=zqice(ig,l)*(pplev(ig,l) - pplev(ig,l+1))
     220               totcloudfrac(ig) = totcloudfrac(ig)
     221     &                   + cloudfrac(ig,l)*ponder
     222               norm=norm+ponder
     223            ENDDO
     224            totcloudfrac(ig)=MAX(totcloudfrac(ig)/norm,1.e-12) ! min value if NaNs
     225         ENDDO
     226c-----No sub-grid cloud representation (CLFvarying=false)
     227      ELSE
     228         DO l=1,nlay
     229            DO ig=1,ngrid
     230               zteff(ig,l)=pt(ig,l)
     231            END DO
     232         END DO
     233      END IF ! end if (CLFvarying)
    141234
    142235c------------------------------------------------------------------
    143236c Time subsampling for microphysics
    144 c------------------------------------------------------------------
     237      rhocloud(1:ngrid,1:nlay) = rho_dust
    145238      DO microstep=1,imicro
    146239     
     
    189282          ENDDO
    190283        ENDDO
    191        
     284c------ Effective tracers quantities in the cloud fraction
     285        IF (CLFvarying) THEN     
     286           pqeff(:,:,:)=pq(:,:,:) ! prevent from buggs (A. Pottier)
     287           pqeff(:,:,igcm_ccn_mass) =pq(:,:,igcm_ccn_mass)/
     288     &                                cloudfrac(:,:)
     289           pqeff(:,:,igcm_ccn_number)=
     290     &     pq(:,:,igcm_ccn_number)/cloudfrac(:,:)
     291           pqeff(:,:,igcm_h2o_ice)= pq(:,:,igcm_h2o_ice)/
     292     &                               cloudfrac(:,:)
     293        ELSE
     294           pqeff(:,:,:)=pq(:,:,:)
     295           pqeff(:,:,igcm_ccn_mass)= pq(:,:,igcm_ccn_mass)
     296           pqeff(:,:,igcm_ccn_number)= pq(:,:,igcm_ccn_number)
     297           pqeff(:,:,igcm_h2o_ice)= pq(:,:,igcm_h2o_ice)
     298        END IF     
    192299       
    193300c-------------------------------------------------------------------
     
    196303        IF (microphys) THEN
    197304           CALL improvedclouds(ngrid,nlay,microtimestep,
    198      &             pplay,pt,subpdt,
    199      &             pq,subpdq,subpdqcloud,subpdtcloud,
     305     &             pplay,zteff,subpdt,
     306     &             pqeff,subpdq,subpdqcloud,subpdtcloud,
    200307     &             nq,tauscaling)
    201308
    202309        ELSE
    203310           CALL simpleclouds(ngrid,nlay,microtimestep,
    204      &             pplay,pzlay,pt,subpdt,
    205      &             pq,subpdq,subpdqcloud,subpdtcloud,
     311     &             pplay,pzlay,zteff,subpdt,
     312     &             pqeff,subpdq,subpdqcloud,subpdtcloud,
    206313     &             nq,tau,rice)
    207314        ENDIF
    208 
    209315
    210316c-------------------------------------------------------------------
     
    305411c------- Due to stepped entry, other processes tendencies can add up to negative values
    306412c------- Therefore, enforce positive values and conserve mass
    307 
    308 
    309413       IF(microphys) THEN
    310414        DO l=1,nlay
     
    479583      nuice(1:ngrid,1:nlay)=nuice_ref
    480584
    481 
    482 
     585c------Update tendencies for sub-grid water ice clouds
     586      IF (CLFvarying) THEN
     587        DO ig=1,ngrid
     588          DO l=1,nlay
     589            pdqcloud(ig,l,igcm_dust_mass)=pdqcloud(ig,l,igcm_dust_mass)
     590     &          *cloudfrac(ig,l)
     591            pdqcloud(ig,l,igcm_ccn_mass)=pdqcloud(ig,l,igcm_ccn_mass)
     592     &          *cloudfrac(ig,l)
     593            pdqcloud(ig,l,igcm_dust_number)=pdqcloud(ig,l,
     594     &           igcm_dust_number) *cloudfrac(ig,l)
     595            pdqcloud(ig,l,igcm_ccn_number)=pdqcloud(ig,l,
     596     &           igcm_ccn_number) *cloudfrac(ig,l)
     597            pdqcloud(ig,l,igcm_h2o_vap)=pdqcloud(ig,l,
     598     &           igcm_h2o_vap) *cloudfrac(ig,l)
     599            pdqcloud(ig,l,igcm_h2o_ice)=pdqcloud(ig,l,
     600     &           igcm_h2o_ice) *cloudfrac(ig,l)
     601          ENDDO
     602        ENDDO   
     603        pdtcloud(:,:)=pdtcloud(:,:)*cloudfrac(:,:)
     604      ENDIF
    483605c=======================================================================
    484 
    485       END
    486  
     606      call WRITEDIAGFI(ngrid,"pdqice2","pdqcloudice apres microphysique"
     607     &      ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,igcm_h2o_ice))
     608      call WRITEDIAGFI(ngrid,"pdqvap2","pdqcloudvap apres microphysique"
     609     &      ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
     610     &      igcm_h2o_vap))
     611      call WRITEDIAGFI(ngrid,"pdqccn2","pdqcloudccn apres microphysique"
     612     &      ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
     613     &      igcm_ccn_mass))
     614      call WRITEDIAGFI(ngrid,"pdqccnN2","pdqcloudccnN apres
     615     &      microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
     616     &      igcm_ccn_number))
     617      call WRITEDIAGFI(ngrid,"pdqdust2", "pdqclouddust apres
     618     &      microphysique","kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
     619     &      igcm_dust_mass))
     620      call WRITEDIAGFI(ngrid,"pdqdustN2", "pdqclouddustN apres
     621     &      microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
     622     &      igcm_dust_number))
     623
     624
     625c=======================================================================
     626
     627      END SUBROUTINE watercloud
     628     
     629      END MODULE watercloud_mod
Note: See TracChangeset for help on using the changeset viewer.