MODULE watercloud_mod IMPLICIT NONE REAL,SAVE,ALLOCATABLE :: zdqcloud(:,:,:) ! tendencies on pq due to condensation of H2O(kg/kg.s-1) REAL,SAVE,ALLOCATABLE :: zdqscloud(:,:) ! tendencies on qsurf (calculated only by calchim but declared here) !$OMP THREADPRIVATE(zdqcloud,zdqscloud) CONTAINS SUBROUTINE watercloud(ngrid,nlay,ptimestep, & pplev,pplay,pdpsrf,pzlay,pt,pdt, & pq,pdq,pdqcloud,pdtcloud, & nq,tau,tauscaling,rdust,rice,nuice, & rsedcloud,rhocloud,totcloudfrac) USE ioipsl_getin_p_mod, ONLY : getin_p USE updaterad, ONLY: updaterdust, updaterice_micro, & updaterice_typ USE improvedclouds_mod, ONLY: improvedclouds USE watersat_mod, ONLY: watersat use tracer_mod, only: nqmx, igcm_h2o_vap, igcm_h2o_ice, & igcm_hdo_vap, igcm_hdo_ice, & igcm_dust_mass, igcm_dust_number, & igcm_ccn_mass, igcm_ccn_number, & rho_dust, nuice_sed, nuice_ref, & qperemin use dimradmars_mod, only: naerkind use conc_mod, only: mmean use write_output_mod, only: write_output IMPLICIT NONE c======================================================================= c Water-ice cloud formation c c Includes two different schemes: c - A simplified scheme (see simpleclouds.F) c - An improved microphysical scheme (see improvedclouds.F) c c There is a time loop specific to cloud formation c due to timescales smaller than the GCM integration timestep. c c Authors: Franck Montmessin, Francois Forget, Ehouarn Millour, c J.-B. Madeleine, Thomas Navarro c c 2004 - 2012 c c 2023: J. Naar, adding different subtimestep on each grid cell c plus not doing microphysics if no water present c plus simpleclouds no longer in the loop for microphysics c======================================================================= c----------------------------------------------------------------------- c declarations: c ------------- include "callkeys.h" c Inputs/outputs: c ------ INTEGER, INTENT(IN) :: ngrid,nlay INTEGER, INTENT(IN) :: nq ! nombre de traceurs REAL, INTENT(IN) :: ptimestep ! pas de temps physique (s) REAL, INTENT(IN) :: pplev(ngrid,nlay+1) ! pression aux inter-couches (Pa) REAL, INTENT(IN) :: pplay(ngrid,nlay) ! pression au milieu des couches (Pa) REAL, INTENT(IN) :: pdpsrf(ngrid) ! tendence surf pressure REAL, INTENT(IN) :: pzlay(ngrid,nlay) ! altitude at the middle of the layers REAL, INTENT(IN) :: pt(ngrid,nlay) ! temperature at the middle of the layers (K) REAL, INTENT(IN) :: pdt(ngrid,nlay) ! tendence temperature des autres param. REAL, INTENT(IN) :: pq(ngrid,nlay,nq) ! traceur (kg/kg) rEAL, INTENT(IN) :: pdq(ngrid,nlay,nq) ! tendence avant condensation (kg/kg.s-1) REAL, INTENT(IN) :: tau(ngrid,naerkind) ! Column dust optical depth at each point REAL, INTENT(IN) :: tauscaling(ngrid) ! Convertion factor for dust amount REAL, INTENT(INOUT) :: rdust(ngrid,nlay) ! Dust geometric mean radius (m) REAL, INTENT(OUT) :: pdqcloud(ngrid,nlay,nq) ! tendence de la condensation H2O(kg/kg.s-1) REAL, INTENT(OUT) :: pdtcloud(ngrid,nlay) ! tendence temperature due ! a la chaleur latente REAL, INTENT(INOUT) :: rice(ngrid,nlay) ! Ice mass mean radius (m) ! (r_c in montmessin_2004) REAL, INTENT(OUT) :: nuice(ngrid,nlay) ! Estimated effective variance ! of the size distribution REAL, INTENT(OUT) :: rsedcloud(ngrid,nlay) ! Cloud sedimentation radius REAL, INTENT(OUT) :: rhocloud(ngrid,nlay) ! Cloud density (kg.m-3) REAL, INTENT(INOUT):: totcloudfrac(ngrid) ! Cloud fraction (A. Pottier 2013) c Locals: c ------ ! for ice radius computation REAL Mo,No REAl ccntyp ! for time loop INTEGER microstep ! time subsampling step variable INTEGER,SAVE :: imicro ! time subsampling for coupled water microphysics & sedimentation REAL,SAVE :: microtimestep ! integration timestep for coupled water microphysics & sedimentation REAL,SAVE :: microtimestep_prev=-999 !$OMP THREADPRIVATE(imicro,microtimestep) !$OMP THREADPRIVATE(microtimestep_prev) ! tendency given by clouds (inside the micro loop) REAL subpdqcloud(ngrid,nlay,nq) ! cf. pdqcloud REAL subpdtcloud(ngrid,nlay) ! cf. pdtcloud ! global tendency (clouds+physics) ! JN : keeping this for simpleclouds scheme REAL sum_subpdq(ngrid,nlay,nq) ! cf. pdqcloud REAL sum_subpdt(ngrid,nlay) ! cf. pdtcloud ! no supersaturation when option supersat is false REAL zt(ngrid,nlay) ! local value of temperature REAL zqsat(ngrid,nlay) ! saturation INTEGER iq,ig,l LOGICAL,SAVE :: firstcall=.true. !$OMP THREADPRIVATE(firstcall) ! HDO cycle REAL :: alpha(ngrid,nlay) ! fractionation coefficient for HDO REAL :: zq0(ngrid,nlay,nq) ! Initial mixing ratio: intermediate variable for HDO ! Representation of sub-grid water ice clouds A. Pottier 2013 REAL :: ztclf(ngrid, nlay) REAL :: zqclf(ngrid, nlay,nq) REAL :: zdelt REAL :: norm REAL :: ponder REAL :: tcond(ngrid,nlay) REAL :: zqvap(ngrid,nlay) REAL :: zqice(ngrid,nlay) REAL :: spant ! delta T for the temperature distribution ! REAL :: zqsat(ngrid,nlay) ! saturation REAL :: pteff(ngrid, nlay)! effective temperature in the cloud,neb REAL :: pqeff(ngrid, nlay, nq)! effective tracers quantities in the cloud REAL :: cloudfrac(ngrid,nlay) ! cloud fraction REAL :: mincloud ! min cloud frac LOGICAL, save :: flagcloud=.true. !$OMP THREADPRIVATE(flagcloud) ! Scheme for adaptative timestep J. Naar 2023 c LOGICAL :: computed_micro(ngrid,nlay) ! Check if microphy was done in this cell REAL :: computed_micro(ngrid,nlay) ! Check if microphy was done inthis cell (logical) REAL :: zt_micro(ngrid,nlay) ! Temperature during microphysics (K) REAL :: zq_micro(ngrid,nlay,nq) ! Tracers during microphysics (kg/kg) REAL :: zqsat_micro(ngrid,nlay) ! Theoritical q(h2o_vap) at saturation (kg/kg) INTEGER :: zimicro(ngrid,nlay) ! Number of ptimestep divisions REAL :: zpotcond_inst(ngrid,nlay) ! condensable water at the beginning of physics (kg/kg) REAL :: zpotcond_full(ngrid,nlay) ! condensable water with integrated tendancies (kg/kg) REAL :: zpotcond(ngrid,nlay) ! maximal condensable water, used to compute adaptative subdivision of ptimestep c ** un petit test de coherence c -------------------------- IF (firstcall) THEN if (nq.gt.nqmx) then write(*,*) 'stop in watercloud (nq.gt.nqmx)!' write(*,*) 'nq=',nq,' nqmx=',nqmx call abort_physic("watercloud","nq.gt.nqmx",1) endif write(*,*) "watercloud: igcm_h2o_vap=",igcm_h2o_vap write(*,*) " igcm_h2o_ice=",igcm_h2o_ice write(*,*) "time subsampling for microphysic ?" #ifdef MESOSCALE imicro = 2 #else imicro = 30 #endif call getin_p("imicro",imicro) write(*,*)"watercloud: imicro = ",imicro firstcall=.false. ENDIF ! of IF (firstcall) !! AS: moved out of firstcall to allow nesting+evoluting timestep !! TBD: consider possible diff imicro with domains? microtimestep = ptimestep/real(imicro) if (microtimestep/=microtimestep_prev) then ! only tell the world if microtimestep has changed write(*,*)"watercloud: Physical timestep is ",ptimestep write(*,*)"watercloud: Microphysics timestep is ",microtimestep microtimestep_prev=microtimestep endif c-----Initialization sum_subpdq(1:ngrid,1:nlay,1:nq) = 0 sum_subpdt(1:ngrid,1:nlay) = 0 ! default value if no ice rhocloud(1:ngrid,1:nlay) = rho_dust c------------------------------------------------------------------- c 0. Representation of sub-grid water ice clouds c------------------ c-----Initialization pteff(1:ngrid,1:nlay) = 0 pqeff(1:ngrid,1:nlay,1:nq) = 0 DO l=1,nlay DO ig=1,ngrid pteff(ig,l)=pt(ig,l) END DO END DO DO l=1,nlay DO ig=1,ngrid DO iq=1,nq pqeff(ig,l,iq)=pq(ig,l,iq) ENDDO ENDDO ENDDO c-----Tendencies DO l=1,nlay DO ig=1,ngrid ztclf(ig,l)=pt(ig,l)+ pdt(ig,l)*ptimestep ENDDO ENDDO DO l=1,nlay DO ig=1,ngrid DO iq=1,nq zqclf(ig,l,iq)=pq(ig,l,iq)+pdq(ig,l,iq)*ptimestep ENDDO ENDDO ENDDO c-----Effective temperature calculation IF (CLFvarying) THEN spant=3.0 ! delta T for the temprature distribution mincloud=0.1 ! min cloudfrac when there is ice IF (flagcloud) THEN write(*,*) "Delta T", spant write(*,*) "mincloud", mincloud flagcloud=.false. END IF !CALL watersat(ngrid*nlay,ztclf,pplay,zqsat) !MV17: we dont need zqsat in the CLFvarying scheme zqvap=zqclf(:,:,igcm_h2o_vap) zqice=zqclf(:,:,igcm_h2o_ice) CALL tcondwater(ngrid*nlay,pplay,zqvap+zqice,tcond) DO l=1,nlay DO ig=1,ngrid zdelt=spant !MAX(spant*ztclf(ig,l),1.e-12), now totally in K. Fixed width IF (tcond(ig,l) .ge. (ztclf(ig,l)+zdelt)) THEN pteff(ig,l)=ztclf(ig,l) cloudfrac(ig,l)=1. ELSE IF (tcond(ig,l) .le. (ztclf(ig,l)-zdelt)) THEN pteff(ig,l)=ztclf(ig,l)-zdelt cloudfrac(ig,l)=mincloud ELSE cloudfrac(ig,l)=(tcond(ig,l)-ztclf(ig,l)+zdelt)/ & (2.0*zdelt) pteff(ig,l)=(tcond(ig,l)+ztclf(ig,l)-zdelt)/2. END IF pteff(ig,l)=pteff(ig,l)-pdt(ig,l)*ptimestep IF (cloudfrac(ig,l).le.mincloud) THEN !MV17: replaced .le.0 by .le.mincloud cloudfrac(ig,l)=mincloud ELSE IF (cloudfrac(ig,l).gt.1) THEN cloudfrac(ig,l)=1. END IF ENDDO ENDDO c-----Calculation of the total cloud coverage of the column DO ig=1,ngrid totcloudfrac(ig) = 0. norm=0. DO l=1,nlay ponder=zqice(ig,l)*(pplev(ig,l) - pplev(ig,l+1)) totcloudfrac(ig) = totcloudfrac(ig) & + cloudfrac(ig,l)*ponder norm=norm+ponder ENDDO totcloudfrac(ig)=MAX(totcloudfrac(ig)/norm,1.e-12) ! min value if NaNs ENDDO c-----Effective tracers quantities in the cloud fraction IF (microphys) THEN pqeff(:,:,igcm_ccn_mass)=pq(:,:,igcm_ccn_mass)/ & cloudfrac(:,:) pqeff(:,:,igcm_ccn_number)=pq(:,:,igcm_ccn_number)/ & cloudfrac(:,:) END IF ! end if (microphys) pqeff(:,:,igcm_h2o_ice)=pq(:,:,igcm_h2o_ice)/ & cloudfrac(:,:) !! CLFvarying outputs ! CALL write_output('pqeffice','pqeffice', ! & 'kg/kg',pqeff(:,:,igcm_h2o_ice)) ! CALL write_output('pteff','pteff', ! & 'K',pteff(:,:)) ! CALL write_output('tcond','tcond', ! & 'K',tcond(:,:)) ! CALL write_output('cloudfrac','cloudfrac', ! & 'K',cloudfrac(:,:)) END IF ! end if (CLFvarying) c------------------------------------------------------------------ c Time subsampling for microphysics c------------------------------------------------------------------ rhocloud(1:ngrid,1:nlay) = rho_dust c Initialisation of all the stuff JN2023 c computed_micro(1:ngrid,1:nlay)=.false. computed_micro(1:ngrid,1:nlay)=0. zt_micro(:,:)=pt(:,:) zq_micro(:,:,:)=pq(:,:,:) zq_micro(:,:,:)=pq(:,:,:) call watersat(ngrid*nlay,zt_micro,pplay,zqsat_micro) zpotcond_inst=zq_micro(:,:,igcm_h2o_vap) - zqsat_micro call watersat(ngrid*nlay,zt_micro+pdt*ptimestep,pplay,zqsat_micro) zpotcond_full=(zq_micro(:,:,igcm_h2o_vap)+ & pdq(:,:,igcm_h2o_vap)*ptimestep) - zqsat_micro zimicro(1:ngrid,1:nlay)=imicro if (cloud_adapt_ts) then call adapt_imicro(ptimestep,zpotcond(ig,l), $ zimicro(ig,l)) endif! (cloud_adapt_ts) then DO l=1,nlay DO ig=1,ngrid c Start by computing the condensable water vapor amount if (zpotcond_full(ig,l).gt.0.) then zpotcond(ig,l)=max(zpotcond_inst(ig,l),zpotcond_full(ig,l)) else if (zpotcond_full(ig,l).le.0.) then zpotcond(ig,l)=min(zpotcond_inst(ig,l),zpotcond_full(ig,l)) endif! (zpotcond_full.gt.0.) then microtimestep=ptimestep/real(zimicro(ig,l)) c Check if microphysics is even needed, that is if enough action c is happening water-wise if ((pq(ig,l,igcm_h2o_ice)+pdq(ig,l,igcm_h2o_ice)*ptimestep & .gt.1e-22) .or. (abs(zpotcond(ig,l)).gt.1e-22)) then c Eventuellement sortir simpleclouds de la boucle egalement computed_micro(ig,l)=1. DO microstep=1,zimicro(ig,l) c JN : incrementing after main microphysics scheme c Previously we were incrementing tendencies, we now c increment tracers and temperature directly c We are thus starting at the end of the first iteration c c------------------------------------------------------------------- c 1. Main call to the different cloud schemes: c------------------------------------------------ IF (microphys) THEN CALL improvedclouds(microtimestep, & pplay(ig,l),zt_micro(ig,l), & zq_micro(ig,l,:),subpdqcloud(ig,l,:), & subpdtcloud(ig,l),nq,tauscaling(ig),mmean(ig,l)) ELSE c Simpleclouds should maybe be taken out and put in a specific loop ? CALL simpleclouds(ngrid,nlay,microtimestep, & pplay,pzlay,pteff,sum_subpdt, & pqeff,sum_subpdq,subpdqcloud,subpdtcloud, & nq,tau,rice) ENDIF c------------------------------------------------------------------- c 2. Updating tracers and temperature after cloud scheme: c----------------------------------------------- IF (microphys) THEN zq_micro(ig,l,igcm_dust_mass) = & zq_micro(ig,l,igcm_dust_mass)+(pdq(ig,l,igcm_dust_mass) & +subpdqcloud(ig,l,igcm_dust_mass))*microtimestep zq_micro(ig,l,igcm_dust_number) = & zq_micro(ig,l,igcm_dust_number) & +(pdq(ig,l,igcm_dust_number) & + subpdqcloud(ig,l,igcm_dust_number))*microtimestep zq_micro(ig,l,igcm_ccn_mass) = & zq_micro(ig,l,igcm_ccn_mass) + & (pdq(ig,l,igcm_ccn_mass) & +subpdqcloud(ig,l,igcm_ccn_mass))*microtimestep zq_micro(ig,l,igcm_ccn_number) = & zq_micro(ig,l,igcm_ccn_number) + & (pdq(ig,l,igcm_ccn_number) & + subpdqcloud(ig,l,igcm_ccn_number))*microtimestep ENDIF zq_micro(ig,l,igcm_h2o_ice) = & zq_micro(ig,l,igcm_h2o_ice)+ & (pdq(ig,l,igcm_h2o_ice) & + subpdqcloud(ig,l,igcm_h2o_ice))*microtimestep zq_micro(ig,l,igcm_h2o_vap) = & zq_micro(ig,l,igcm_h2o_vap)+ & (pdq(ig,l,igcm_h2o_vap) & + subpdqcloud(ig,l,igcm_h2o_vap))*microtimestep IF (hdo) THEN zq_micro(ig,l,igcm_hdo_ice) = & zq_micro(ig,l,igcm_hdo_ice)+ & (pdq(ig,l,igcm_hdo_ice) & + subpdqcloud(ig,l,igcm_hdo_ice))*microtimestep zq_micro(ig,l,igcm_hdo_vap) = & zq_micro(ig,l,igcm_hdo_vap)+ & (pdq(ig,l,igcm_hdo_vap) & + subpdqcloud(ig,l,igcm_hdo_vap))*microtimestep ENDIF ! hdo c Could also set subpdtcloud to 0 if not activice to make it simpler zt_micro(ig,l) = zt_micro(ig,l)+ & pdt(ig,l)*microtimestep IF (activice) THEN zt_micro(ig,l) = zt_micro(ig,l)+ & subpdtcloud(ig,l)*microtimestep ENDIF ! !! Example of how to use writediagmicrofi useful to ! !! get outputs at each microphysical sub-timestep (better to be used in 1D) ! CALL WRITEDIAGMICROFI(ngrid,imicro,microstep, ! & microtimestep,'subpdtcloud', ! & 'subpdtcloud','K/s',1,subpdtcloud(:,:)) ENDDO ! of DO microstep=1,imicro endif! (zq(ig,l,igcm_h2o_ice)+pdq(ig,l,igcm_h2o_ice)*ptimestep ! & .gt.1e-22).or.(abs(zpotcond).gt.1e-22) then ENDDO ! ig=1,ngrid ENDDO ! l=1,nlay c------ Useful outputs to check how it went call write_output("computed_micro","computed_micro "// & "after microphysics","logical",computed_micro(:,:)) call write_output("zimicro","Used number of subtimestep "// & "in cloud microphysics","integer",real(zimicro(:,:))) c------------------------------------------------------------------- c 3. Compute final tendencies after time loop: c------------------------------------------------ c------ Temperature tendency pdtcloud DO l=1,nlay DO ig=1,ngrid pdtcloud(ig,l) = -pdt(ig,l)+ & (zt_micro(ig,l)-pt(ig,l)) / ptimestep ENDDO ENDDO c------ Tracers tendencies pdqcloud DO l=1,nlay DO ig=1,ngrid pdqcloud(ig,l,igcm_h2o_ice) = & -pdq(ig,l,igcm_h2o_ice)+ & (zq_micro(ig,l,igcm_h2o_ice) - & pq(ig,l,igcm_h2o_ice)) / ptimestep pdqcloud(ig,l,igcm_h2o_vap) = & -pdq(ig,l,igcm_h2o_vap)+ & (zq_micro(ig,l,igcm_h2o_vap) - & pq(ig,l,igcm_h2o_vap)) / ptimestep IF (hdo) THEN pdqcloud(ig,l,igcm_hdo_ice) = & -pdq(ig,l,igcm_hdo_ice)+ & (zq_micro(ig,l,igcm_hdo_ice) - & pq(ig,l,igcm_hdo_ice)) / ptimestep pdqcloud(ig,l,igcm_hdo_vap) = & -pdq(ig,l,igcm_hdo_vap)+ & (zq_micro(ig,l,igcm_hdo_vap) - & pq(ig,l,igcm_hdo_vap)) / ptimestep ENDIF !hdo ENDDO ENDDO IF(microphys) THEN DO l=1,nlay DO ig=1,ngrid pdqcloud(ig,l,igcm_ccn_mass) = & -pdq(ig,l,igcm_ccn_mass)+ & (zq_micro(ig,l,igcm_ccn_mass) - & pq(ig,l,igcm_ccn_mass)) / ptimestep pdqcloud(ig,l,igcm_ccn_number) = & -pdq(ig,l,igcm_ccn_number)+ & (zq_micro(ig,l,igcm_ccn_number) - & pq(ig,l,igcm_ccn_number)) / ptimestep ENDDO ENDDO ENDIF IF(scavenging) THEN DO l=1,nlay DO ig=1,ngrid pdqcloud(ig,l,igcm_dust_mass) = & -pdq(ig,l,igcm_dust_mass)+ & (zq_micro(ig,l,igcm_dust_mass) - & pq(ig,l,igcm_dust_mass)) / ptimestep pdqcloud(ig,l,igcm_dust_number) = & -pdq(ig,l,igcm_dust_number)+ & (zq_micro(ig,l,igcm_dust_number) - & pq(ig,l,igcm_dust_number)) / ptimestep ENDDO ENDDO ENDIF c------- Due to stepped entry, other processes tendencies can add up to negative values c------- Therefore, enforce positive values and conserve mass IF(microphys) THEN DO l=1,nlay DO ig=1,ngrid IF ((pq(ig,l,igcm_ccn_number) + & ptimestep* (pdq(ig,l,igcm_ccn_number) + & pdqcloud(ig,l,igcm_ccn_number)) .le. 1.) & .or. (pq(ig,l,igcm_ccn_mass) + & ptimestep* (pdq(ig,l,igcm_ccn_mass) + & pdqcloud(ig,l,igcm_ccn_mass)) .le. 1.e-20)) THEN pdqcloud(ig,l,igcm_ccn_number) = & - pq(ig,l,igcm_ccn_number)/ptimestep & - pdq(ig,l,igcm_ccn_number) + 1. pdqcloud(ig,l,igcm_dust_number) = & -pdqcloud(ig,l,igcm_ccn_number) pdqcloud(ig,l,igcm_ccn_mass) = & - pq(ig,l,igcm_ccn_mass)/ptimestep & - pdq(ig,l,igcm_ccn_mass) + 1.e-20 pdqcloud(ig,l,igcm_dust_mass) = & -pdqcloud(ig,l,igcm_ccn_mass) ENDIF ENDDO ENDDO ENDIF IF(scavenging) THEN DO l=1,nlay DO ig=1,ngrid IF ((pq(ig,l,igcm_dust_number) + & ptimestep* (pdq(ig,l,igcm_dust_number) + & pdqcloud(ig,l,igcm_dust_number)) .le. 1.) & .or. (pq(ig,l,igcm_dust_mass) + & ptimestep* (pdq(ig,l,igcm_dust_mass) + & pdqcloud(ig,l,igcm_dust_mass)) .le. 1.e-20)) THEN pdqcloud(ig,l,igcm_dust_number) = & - pq(ig,l,igcm_dust_number)/ptimestep & - pdq(ig,l,igcm_dust_number) + 1. pdqcloud(ig,l,igcm_ccn_number) = & -pdqcloud(ig,l,igcm_dust_number) pdqcloud(ig,l,igcm_dust_mass) = & - pq(ig,l,igcm_dust_mass)/ptimestep & - pdq(ig,l,igcm_dust_mass) + 1.e-20 pdqcloud(ig,l,igcm_ccn_mass) = & -pdqcloud(ig,l,igcm_dust_mass) ENDIF ENDDO ENDDO ENDIF DO l=1,nlay DO ig=1,ngrid IF (pq(ig,l,igcm_h2o_ice) + ptimestep* & (pdq(ig,l,igcm_h2o_ice) + pdqcloud(ig,l,igcm_h2o_ice)) & .le. qperemin) THEN pdqcloud(ig,l,igcm_h2o_ice) = & - pq(ig,l,igcm_h2o_ice)/ptimestep - pdq(ig,l,igcm_h2o_ice) pdqcloud(ig,l,igcm_h2o_vap) = -pdqcloud(ig,l,igcm_h2o_ice) if (hdo) then pdqcloud(ig,l,igcm_hdo_ice) = & - pq(ig,l,igcm_hdo_ice)/ptimestep - pdq(ig,l,igcm_hdo_ice) pdqcloud(ig,l,igcm_hdo_vap) = -pdqcloud(ig,l,igcm_hdo_ice) endif ENDIF IF (pq(ig,l,igcm_h2o_vap) + ptimestep* & (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap)) & .le. qperemin) THEN pdqcloud(ig,l,igcm_h2o_vap) = & - pq(ig,l,igcm_h2o_vap)/ptimestep - pdq(ig,l,igcm_h2o_vap) pdqcloud(ig,l,igcm_h2o_ice) = -pdqcloud(ig,l,igcm_h2o_vap) if (hdo) then pdqcloud(ig,l,igcm_hdo_vap) = & - pq(ig,l,igcm_hdo_vap)/ptimestep - pdq(ig,l,igcm_hdo_vap) pdqcloud(ig,l,igcm_hdo_ice) = -pdqcloud(ig,l,igcm_hdo_vap) endif ENDIF ENDDO ENDDO c------Update the ice and dust particle size "rice" for output or photochemistry c------Only rsedcloud is used for the water cycle IF(scavenging) THEN DO l=1, nlay DO ig=1,ngrid call updaterdust( & pq(ig,l,igcm_dust_mass) + ! dust mass & (pdq(ig,l,igcm_dust_mass) + ! dust mass & pdqcloud(ig,l,igcm_dust_mass))*ptimestep, ! dust mass & pq(ig,l,igcm_dust_number) + ! dust number & (pdq(ig,l,igcm_dust_number) + ! dust number & pdqcloud(ig,l,igcm_dust_number))*ptimestep, ! dust number & rdust(ig,l)) ENDDO ENDDO ENDIF IF(microphys) THEN ! In case one does not want to allow supersatured water when using microphysics. ! Not done by default. IF(.not.supersat) THEN ! !! initial mixing ratios for initial D/H ratio calculation zq0(:,:,:) = pq(:,:,:) + pdq(:,:,:)*ptimestep zt = pt + (pdt+pdtcloud)*ptimestep call watersat(ngrid*nlay,zt,pplay,zqsat) DO l=1, nlay DO ig=1,ngrid IF (pq(ig,l,igcm_h2o_vap) & + (pdq(ig,l,igcm_h2o_vap) + pdqcloud(ig,l,igcm_h2o_vap)) & * ptimestep .ge. zqsat(ig,l)) THEN pdqcloud(ig,l,igcm_h2o_vap) = & (zqsat(ig,l) - pq(ig,l,igcm_h2o_vap))/ptimestep & - pdq(ig,l,igcm_h2o_vap) pdqcloud(ig,l,igcm_h2o_ice) = & -pdqcloud(ig,l,igcm_h2o_vap) ! !! HDO ice flux has to be modified in consequence IF (hdo) THEN ! !! Logically only condensation can happen in this case IF (pdqcloud(ig,l,igcm_h2o_ice) .gt. 0.0) THEN IF ( zq0(ig,l,igcm_h2o_vap) .gt. qperemin ) THEN ! !! Lamb et al. 2017 alpha(ig,l) = exp(13525./zt(ig,l)**2.-5.59d-2) pdqcloud(ig,l,igcm_hdo_ice) = & pdqcloud(ig,l,igcm_h2o_ice)*alpha(ig,l)* & ( zq0(ig,l,igcm_hdo_vap) & /zq0(ig,l,igcm_h2o_vap) ) pdqcloud(ig,l,igcm_hdo_ice) = & min(pdqcloud(ig,l,igcm_hdo_ice), & zq0(ig,l,igcm_hdo_vap)/ptimestep) ELSE pdqcloud(ig,l,igcm_hdo_ice) = 0. ENDIF !! sublimation ELSE IF ( zq0(ig,l,igcm_h2o_ice) .gt. qperemin ) THEN pdqcloud(ig,l,igcm_hdo_ice) = & pdqcloud(ig,l,igcm_h2o_ice)* & ( zq0(ig,l,igcm_hdo_ice) & /zq0(ig,l,igcm_h2o_ice) ) pdqcloud(ig,l,igcm_hdo_ice) = & max(pdqcloud(ig,l,igcm_hdo_ice), & -zq0(ig,l,igcm_hdo_ice)/ptimestep) ELSE pdqcloud(ig,l,igcm_hdo_ice) = 0. ENDIF ENDIF !IF (pdqcloud(ig,l,igcm_h2o_ice).gt.0.) pdqcloud(ig,l,igcm_hdo_vap) = & -pdqcloud(ig,l,igcm_hdo_ice) ENDIF !IF (hdo) ! no need to correct ccn_number, updaterad can handle this properly. ENDIF ENDDO ENDDO ENDIF DO l=1, nlay DO ig=1,ngrid call updaterice_micro( & pq(ig,l,igcm_h2o_ice) + ! ice mass & (pdq(ig,l,igcm_h2o_ice) + ! ice mass & pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass & pq(ig,l,igcm_ccn_mass) + ! ccn mass & (pdq(ig,l,igcm_ccn_mass) + ! ccn mass & pdqcloud(ig,l,igcm_ccn_mass))*ptimestep, ! ccn mass & pq(ig,l,igcm_ccn_number) + ! ccn number & (pdq(ig,l,igcm_ccn_number) + ! ccn number & pdqcloud(ig,l,igcm_ccn_number))*ptimestep, ! ccn number & tauscaling(ig),rice(ig,l),rhocloud(ig,l)) ENDDO ENDDO ELSE ! no microphys DO l=1,nlay DO ig=1,ngrid call updaterice_typ( & pq(ig,l,igcm_h2o_ice) + ! ice mass & (pdq(ig,l,igcm_h2o_ice) + ! ice mass & pdqcloud(ig,l,igcm_h2o_ice))*ptimestep, ! ice mass & tau(ig,1),pzlay(ig,l),rice(ig,l)) ENDDO ENDDO ENDIF ! of IF(microphys) c A correction if a lot of subliming CO2 fills the 1st layer FF04/2005 c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ c Then that should not affect the ice particle radius do ig=1,ngrid if(pdpsrf(ig)*ptimestep.gt.0.9*(pplev(ig,1)-pplev(ig,2)))then if(pdpsrf(ig)*ptimestep.gt.0.9*(pplev(ig,1)-pplev(ig,3))) & rice(ig,2)=rice(ig,3) rice(ig,1)=rice(ig,2) end if end do DO l=1,nlay DO ig=1,ngrid rsedcloud(ig,l)=max(rice(ig,l)* & (1.+nuice_sed)*(1.+nuice_sed)*(1.+nuice_sed), & rdust(ig,l)) ! rsedcloud(ig,l)=min(rsedcloud(ig,l),1.e-4) ENDDO ENDDO ! used for rad. transfer calculations ! nuice is constant because a lognormal distribution is prescribed nuice(1:ngrid,1:nlay)=nuice_ref c------Update tendencies for sub-grid water ice clouds IF (CLFvarying) THEN DO ig=1,ngrid DO l=1,nlay pdqcloud(ig,l,igcm_dust_mass)=pdqcloud(ig,l,igcm_dust_mass) & *cloudfrac(ig,l) pdqcloud(ig,l,igcm_ccn_mass)=pdqcloud(ig,l,igcm_ccn_mass) & *cloudfrac(ig,l) pdqcloud(ig,l,igcm_dust_number)=pdqcloud(ig,l, & igcm_dust_number) *cloudfrac(ig,l) pdqcloud(ig,l,igcm_ccn_number)=pdqcloud(ig,l, & igcm_ccn_number) *cloudfrac(ig,l) pdqcloud(ig,l,igcm_h2o_vap)=pdqcloud(ig,l, & igcm_h2o_vap) *cloudfrac(ig,l) pdqcloud(ig,l,igcm_h2o_ice)=pdqcloud(ig,l, & igcm_h2o_ice) *cloudfrac(ig,l) ENDDO ENDDO pdtcloud(:,:)=pdtcloud(:,:)*cloudfrac(:,:) ENDIF #ifndef MESOSCALE c======================================================================= call write_output("watercloud_pdqh2oice","pdqcloud_h2o_ice "// & "after microphysics","kg/kg.s-1",pdqcloud(:,:,igcm_h2o_ice)) call write_output("watercloud_pdqh2ovap","pdqcloud_h2o_vap "// & "after microphysics","kg/kg.s-1",pdqcloud(:,:,igcm_h2o_vap)) if (hdo) then call write_output("watercloud_pdqhdoice","pdqcloud_hdo_ice "// & "after microphysics","kg/kg.s-1",pdqcloud(:,:,igcm_hdo_ice)) call write_output("watercloud_pdqhdovap","pdqcloud_hdo_vap "// & "after microphysics","kg/kg.s-1",pdqcloud(:,:,igcm_hdo_vap)) endif ! call write_output("pdqccn2","pdqcloudccn apres microphysique" ! & ,"kg/kg.s-1",pdqcloud(:,:, ! & igcm_ccn_mass)) ! call write_output("pdqccnN2","pdqcloudccnN apres "// ! & "microphysique","nb/kg.s-1",pdqcloud(:,:, ! & igcm_ccn_number)) ! call write_output("pdqdust2", "pdqclouddust apres "// ! & "microphysique","kg/kg.s-1",pdqcloud(:,:, ! & igcm_dust_mass)) ! call write_output("pdqdustN2", "pdqclouddustN apres "// ! & "microphysique","nb/kg.s-1",pdqcloud(:,:, ! & igcm_dust_number)) c======================================================================= #endif END SUBROUTINE watercloud subroutine ini_watercloud_mod(ngrid,nlayer,nq) implicit none integer,intent(in) :: ngrid ! number of atmospheric columns integer,intent(in) :: nlayer ! number of atmospheric layers integer,intent(in) :: nq ! number of tracers allocate(zdqcloud(ngrid,nlayer,nq)) zdqcloud(:,:,:)=0 allocate(zdqscloud(ngrid,nq)) zdqscloud(:,:)=0 end subroutine ini_watercloud_mod subroutine end_watercloud_mod implicit none if (allocated(zdqcloud)) deallocate(zdqcloud) if (allocated(zdqscloud)) deallocate(zdqscloud) end subroutine end_watercloud_mod SUBROUTINE adapt_imicro(ptimestep,potcond, $ zimicro) c Pas de temps adaptatif pour les nuages real,intent(in) :: ptimestep ! total duration of physics (sec) real,intent(in) :: potcond ! total duration of physics (sec) real :: alpha, beta ! total duration of physics (sec) integer,intent(out) :: zimicro ! number of ptimestep division c zimicro = ptimestep*alpha*potcond**beta zimicro = 30 END SUBROUTINE adapt_imicro END MODULE watercloud_mod