module simpleclouds_mod implicit none contains subroutine simpleclouds(ngrid,nlay,ptimestep, & pplay,pzlay,pt,pdt, & pq,pdq,pdqcloud,pdtcloud, & nq,tau,rice) USE updaterad, ONLY: updaterice_typ USE watersat_mod, ONLY: watersat use tracer_mod, only: igcm_h2o_vap, igcm_h2o_ice, & igcm_hdo_vap, igcm_hdo_ice, & qperemin use comcstfi_h, only: cpp use dimradmars_mod, only: naerkind use callkeys_mod, only: hdo, hdofrac implicit none c------------------------------------------------------------------ c This routine is used to form clouds when a parcel of the GCM is c saturated. It is a simplified scheme, and there is almost no c microphysics involved. When the air is saturated, water-ice c clouds form on a fraction of the dust particles, specified by c the constant called "ccn_factor". There is no supersaturation, c and no nucleation rates computed. A more accurate scheme can c be found in the routine called "improvedclouds.F". c Authors: Franck Montmessin (water ice scheme) c Francois Forget (changed nuclei density & outputs) c Ehouarn Millour (sept.2008, tracers are now handled c by name and not fixed index) c J.-B. Madeleine (developed a single routine called c simpleclouds.F, and corrected calculations c of the typical CCN profile, Oct. 2011) c------------------------------------------------------------------ c Arguments: c --------- c Inputs: integer,intent(in) :: ngrid ! number of atmospheric columns integer,intent(in) :: nlay ! number of atmospheric layers integer,intent(in) :: nq ! number of tracers real,intent(in) :: ptimestep ! physics time step (s) real,intent(in) :: pplay(ngrid,nlay) ! pressure at mid-layer (Pa) real,intent(in) :: pzlay(ngrid,nlay) ! altitude of the layers (m) real,intent(in) :: pt(ngrid,nlay) ! input temperature at mid-layer (K) real,intent(in) :: pdt(ngrid,nlay) ! tendency on temperature from previous paramatrizations (K/s) real,intent(in) :: pq(ngrid,nlay,nq) ! input tracer mixing ratio (kg/kg) real,intent(in) :: pdq(ngrid,nlay,nq) ! tendency on tracers from previous parametrizations (kg/kg.s-1) real,intent(in) :: tau(ngrid,naerkind) ! Column dust optical depth in each column c Output: real,intent(out) :: rice(ngrid,nlay) ! Water ice mass mean radius (m) ! (r_c in montmessin_2004) real,intent(out) :: pdqcloud(ngrid,nlay,nq) ! tendencies due to H2O condensation ! and sublimation (kg/kg.s-1) real,intent(out) :: pdtcloud(ngrid,nlay) ! tendency on temperature due ! to latent heat (K/s) c------------------------------------------------------------------ c Local variables: REAL rhocloud(ngrid,nlay) ! Cloud density (kg.m-3) INTEGER ig,l REAL zq(ngrid,nlay,nq) ! local value of tracers REAL zq0(ngrid,nlay,nq) ! local initial value of tracers REAL zt(ngrid,nlay) ! local value of temperature REAL zqsat(ngrid,nlay) ! saturation REAL*8 dzq ! masse de glace echangee (kg/kg) REAL lw !Latent heat of sublimation (J.kg-1) REAL,PARAMETER :: To=273.15 ! reference temperature, T=273.15 K real rdusttyp(ngrid,nlay) ! Typical dust geom. mean radius (m) REAL ccntyp(ngrid,nlay) ! Typical dust number density (#/kg) REAL alpha_c(ngrid,nlay) !!MARGAUX: alpha_c as a spatial variable c CCN reduction factor c REAL, PARAMETER :: ccn_factor = 4.5 !! comme TESTS_JB // 1. avant REAL DoH_vap(ngrid,nlay) REAL DoH_ice(ngrid,nlay) c----------------------------------------------------------------------- c 1. initialisation c ----------------- c Update values of water vapor and ice, and temperature. do l=1,nlay do ig=1,ngrid zq(ig,l,igcm_h2o_vap)= & pq(ig,l,igcm_h2o_vap)+pdq(ig,l,igcm_h2o_vap)*ptimestep zq(ig,l,igcm_h2o_vap)=max(zq(ig,l,igcm_h2o_vap),1.E-30) ! FF 12/2004 zq0(ig,l,igcm_h2o_vap)=zq(ig,l,igcm_h2o_vap) zt(ig,l)=pt(ig,l)+ pdt(ig,l)*ptimestep zq(ig,l,igcm_h2o_ice)= & pq(ig,l,igcm_h2o_ice)+pdq(ig,l,igcm_h2o_ice)*ptimestep zq(ig,l,igcm_h2o_ice)=max(zq(ig,l,igcm_h2o_ice),0.) ! FF 12/2004 zq0(ig,l,igcm_h2o_ice)=zq(ig,l,igcm_h2o_ice) if (hdo) then zq(ig,l,igcm_hdo_vap)= & pq(ig,l,igcm_hdo_vap)+pdq(ig,l,igcm_hdo_vap)*ptimestep zq(ig,l,igcm_hdo_vap)=max(zq(ig,l,igcm_hdo_vap),1e-30) ! FF 12/2004 zq0(ig,l,igcm_hdo_vap)=zq(ig,l,igcm_hdo_vap) zq(ig,l,igcm_hdo_ice)= & pq(ig,l,igcm_hdo_ice)+pdq(ig,l,igcm_hdo_ice)*ptimestep zq(ig,l,igcm_hdo_ice)=max(zq(ig,l,igcm_hdo_ice),1e-30) ! FF 12/2004 zq0(ig,l,igcm_hdo_ice)=zq(ig,l,igcm_hdo_ice) endif !hdo enddo enddo pdqcloud(1:ngrid,1:nlay,1:nq)=0 pdtcloud(1:ngrid,1:nlay)=0 alpha_c(1:ngrid,1:nlay)=0. c ---------------------------------------------- c c Compute saturation mixing ratio in each layer c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ call watersat(ngrid*nlay,zt,pplay,zqsat) c compute condensation rates (kg/kg/s-1) in each layer c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ do l=1,nlay do ig=1,ngrid if (zq(ig,l,igcm_h2o_vap).ge.zqsat(ig,l))then ! Condensation dzq=zq(ig,l,igcm_h2o_vap)-zqsat(ig,l) elseif(zq(ig,l,igcm_h2o_vap).lt.zqsat(ig,l))then ! Sublimation dzq=-min(zqsat(ig,l)-zq(ig,l,igcm_h2o_vap), & zq(ig,l,igcm_h2o_ice)) endif c Water Mass change c ~~~~~~~~~~~~~~~~~ zq(ig,l,igcm_h2o_ice)=zq(ig,l,igcm_h2o_ice)+dzq zq(ig,l,igcm_h2o_vap)=zq(ig,l,igcm_h2o_vap)-dzq enddo ! of do ig=1,ngrid enddo ! of do l=1,nlay c Final tendency c ~~~~~~~~~~~~~~~ do l=1, nlay do ig=1,ngrid pdqcloud(ig,l,igcm_h2o_vap)=(zq(ig,l,igcm_h2o_vap) & -zq0(ig,l,igcm_h2o_vap))/ptimestep pdqcloud(ig,l,igcm_h2o_ice) = & (zq(ig,l,igcm_h2o_ice) - zq0(ig,l,igcm_h2o_ice))/ptimestep lw=(2834.3-0.28*(zt(ig,l)-To)-0.004*(zt(ig,l)-To)**2)*1.e+3 pdtcloud(ig,l)=-pdqcloud(ig,l,igcm_h2o_vap)*lw/cpp if (hdo) then if (pdqcloud(ig,l,igcm_h2o_ice).gt.0.0) then !condens if (hdofrac) then ! do we use fractionation? c alpha_c(ig,l) = exp(16288./zt(ig,l)**2.-9.34d-2) alpha_c(ig,l) = exp(13525./zt(ig,l)**2.-5.59d-2) !Lamb else alpha_c(ig,l) = 1.d0 endif if (zq0(ig,l,igcm_h2o_vap).gt.qperemin) then pdqcloud(ig,l,igcm_hdo_ice)= & pdqcloud(ig,l,igcm_h2o_ice)*alpha_c(ig,l)* & ( zq0(ig,l,igcm_hdo_vap) & /zq0(ig,l,igcm_h2o_vap) ) else pdqcloud(ig,l,igcm_hdo_ice)= 0.0 endif pdqcloud(ig,l,igcm_hdo_ice) = & min(pdqcloud(ig,l,igcm_hdo_ice), & zq0(ig,l,igcm_hdo_vap)/ptimestep) pdqcloud(ig,l,igcm_hdo_vap)= & -pdqcloud(ig,l,igcm_hdo_ice) else ! sublimation 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) ) else pdqcloud(ig,l,igcm_hdo_ice)= 0. endif pdqcloud(ig,l,igcm_hdo_ice) = & max(pdqcloud(ig,l,igcm_hdo_ice), & -zq0(ig,l,igcm_hdo_ice)/ptimestep) pdqcloud(ig,l,igcm_hdo_vap)= & -pdqcloud(ig,l,igcm_hdo_ice) endif ! condensation/sublimation endif ! hdo enddo ! of do ig=1,ngrid enddo ! of do l=1,nlay c ice crystal radius do l=1, nlay do ig=1,ngrid call updaterice_typ(zq(ig,l,igcm_h2o_ice), & tau(ig,1),pzlay(ig,l),rice(ig,l)) end do end do c if (hdo) then c CALL WRITEDIAGFI(ngrid,'alpha_c', c & 'alpha_c', c & ' ',3,alpha_c) c endif !hdo c------------------------------------------------------------------ end subroutine simpleclouds end module simpleclouds_mod