Index: trunk/LMDZ.MARS/libf/phymars/conf_phys.F =================================================================== --- trunk/LMDZ.MARS/libf/phymars/conf_phys.F (revision 3316) +++ trunk/LMDZ.MARS/libf/phymars/conf_phys.F (revision 3325) @@ -38,5 +38,5 @@ USE mod_phys_lmdz_transfert_para, ONLY: bcast USE paleoclimate_mod,ONLY: paleoclimate,albedo_perennialco2, - & lag_layer + & lag_layer, include_waterbuoyancy use microphys_h, only: mteta use comsoil_h, only: adsorption_soil, choice_ads,ads_const_D, @@ -364,11 +364,10 @@ lag_layer=.false. call getin_p("lag_layer",lag_layer) - write(*,*) " lag_layer = ", lag_layer - - write(*,*)"Is it paleoclimate run?" + write(*,*) "lag_layer = ", lag_layer + + write(*,*)"Is it paleoclimate run?" paleoclimate=.false. ! default value call getin_p("paleoclimate",paleoclimate) - write(*,*)" paleoclimate = ",paleoclimate - + write(*,*)"paleoclimate = ",paleoclimate write(*,*)"Albedo for perennial CO2 ice?" @@ -376,4 +375,16 @@ call getin_p("albedo_perennialco2",albedo_perennialco2) write(*,*)"albedo_perennialco2 = ",albedo_perennialco2 + + write(*,*)"Using lag layer??" + lag_layer=.false. + call getin_p("include_waterbuoyancy",include_waterbuoyancy) + write(*,*) "include_waterbuoyancy = ", include_waterbuoyancy + + if (include_waterbuoyancy.and.(.not.paleoclimate)) then + print*,'You should not include the water buoyancy effect + & on current climate (model not tunned with this effect' + call abort_physic(modname, + & "include_waterbuoyancy must be used with paleoclimate",1) + endif ! TRACERS: Index: trunk/LMDZ.MARS/libf/phymars/paleoclimate_mod.F90 =================================================================== --- trunk/LMDZ.MARS/libf/phymars/paleoclimate_mod.F90 (revision 3316) +++ trunk/LMDZ.MARS/libf/phymars/paleoclimate_mod.F90 (revision 3325) @@ -15,10 +15,11 @@ !$OMP THREADPRIVATE(paleoclimate) - real, save, allocatable, dimension(:,:) :: h2o_ice_depth ! Thickness of the lag before H2O ice [m] - real, save, allocatable, dimension(:,:) :: lag_co2_ice ! Thickness of the lag before CO2 ice [m] - real, save, allocatable, dimension(:,:) :: d_coef ! Diffusion coeficent - real, save :: albedo_perennialco2 ! Albedo for perennial co2 ice [1] - logical, save :: lag_layer ! Does lag layer is present? -!$OMP THREADPRIVATE(h2o_ice_depth,d_coef,lag_co2_ice,albedo_perennialco2) + real, save, allocatable, dimension(:,:) :: h2o_ice_depth ! Thickness of the lag before H2O ice [m] + real, save, allocatable, dimension(:,:) :: lag_co2_ice ! Thickness of the lag before CO2 ice [m] + real, save, allocatable, dimension(:,:) :: d_coef ! Diffusion coeficent + real, save :: albedo_perennialco2 ! Albedo for perennial co2 ice [1] + logical, save :: lag_layer ! Does lag layer is present? + logical, save :: include_waterbuoyancy ! Include the effect of water buoyancy when computing the sublimation of water ice ? +!$OMP THREADPRIVATE(h2o_ice_depth,d_coef,lag_co2_ice,albedo_perennialco2,include_waterbuoyancy) !======================================================================= Index: trunk/LMDZ.MARS/libf/phymars/pbl_parameters_mod.F90 =================================================================== --- trunk/LMDZ.MARS/libf/phymars/pbl_parameters_mod.F90 (revision 3316) +++ trunk/LMDZ.MARS/libf/phymars/pbl_parameters_mod.F90 (revision 3325) @@ -19,5 +19,5 @@ CONTAINS - SUBROUTINE pbl_parameters(ngrid,nlay,ps,pplay,pz0,pg,zzlay,zzlev,pu,pv,wstar_in,hfmax,zmax,tke,pts,ph,z_out,n_out, & + SUBROUTINE pbl_parameters(ngrid,nlay,ps,pplay,pz0,pg,zzlay,zzlev,pu,pv,wstar_in,hfmax,zmax,tke,pts,ph,pqvap,pqsurf,mumean,z_out,n_out, & T_out,u_out,ustar,tstar,vhf,vvv) @@ -26,4 +26,6 @@ use turb_mod, only: turb_resolved use lmdz_atke_turbulence_ini, only : smmin, ric, cinf, cepsilon, pr_slope, pr_asym, pr_neut, ri0, ri1, cn, rpi + use watersat_mod, only: watersat + use paleoclimate_mod, only: include_waterbuoyancy IMPLICIT NONE @@ -92,4 +94,7 @@ REAL, INTENT(IN) :: pts(ngrid),ph(ngrid,nlay) ! surface temperature, potentiel temperature REAL, INTENT(IN) :: z_out(n_out) ! altitudes of the interpolation in the surface layer + REAL, INTENT(IN) :: mumean(ngrid) ! Molecular mass of the atmosphere (kg/mol) + REAL, INTENT(IN) :: pqvap(ngrid,nlay) ! Water vapor mixing ratio (kg/kg) to account for vapor flottability + REAL, INTENT(IN) :: pqsurf(ngrid) ! Water ice frost on the surface (kg/m^2) to account for vapor flottability ! Outputs: @@ -142,4 +147,13 @@ REAL f_ri_cd_min ! minimum of the stability function with the ATKE scheme REAL ric_4interp ! critical richardson number which is either the one from Colaitis et al. (2013) or ATKE (1) + + + REAL tsurf_v(ngrid) ! Virtual surface temperature, accounting for vapor flottability + REAL temp_v(ngrid) ! Potential virtual air temperature in the frist layer, accounting for vapor flottability + REAL mu_h2o ! Molecular mass of water vapor (kg/mol) + REAL tol_frost ! Tolerance to consider the effect of frost on the surface + REAL qsat(ngrid) ! saturated mixing ratio of water vapor + + ! Code: ! -------- @@ -166,4 +180,8 @@ tol_iter = 0.01 f_ri_cd_min = 0.01 + mu_h2o = 18e-3 + tol_frost = 1e-4 + tsurf_v(:) = 0. + temp_v(:) = 0. ! this formulation assumes alphah=1., implying betah=betam @@ -183,5 +201,21 @@ endif - + + + IF(include_waterbuoyancy) then + temp_v(:) = ph(:,1)*(1.+pqvap(:,1)/(mu_h2o/mumean(:)))/(1.+pqvap(:,1)) + DO ig = 1,ngrid + IF(pqsurf(ig).gt.tol_frost) then + call watersat(1,pts(ig),pplay(ig,1),qsat(ig)) + tsurf_v(ig) = pts(ig)*(1.+qsat(ig)/(mu_h2o/mumean(ig)))/(1.+qsat(ig)) + ELSE + tsurf_v(ig) = pts(ig)*(1.+pqvap(ig,1)/(mu_h2o/mumean(ig)))/(1.+pqvap(ig,1)) + ENDIF + ENDDO + ELSE + tsurf_v(:) = pts(:) + temp_v(:) = ph(:,1) + ENDIF + DO ig=1,ngrid ite = 0. @@ -199,5 +233,5 @@ IF ((pu(ig,1) .ne. 0.) .or. (pv(ig,1) .ne. 0.)) then ! Richardson number formulation proposed by D.E. England et al. (1995) - rib(ig) = (pg/pts(ig))*sqrt(zzlay(ig,1)*pz0(ig))*(((log(zzlay(ig,1)/pz0(ig)))**2)/(log(zzlay(ig,1)/pz0t)))*(ph(ig,1)-pts(ig))/zu2(ig) + rib(ig) = (pg/tsurf_v(ig))*sqrt(zzlay(ig,1)*pz0(ig))*(((log(zzlay(ig,1)/pz0(ig)))**2)/(log(zzlay(ig,1)/pz0t)))*(temp_v(ig)-tsurf_v(ig))/zu2(ig) ELSE print*,'warning, infinite Richardson at surface' Index: trunk/LMDZ.MARS/libf/phymars/physiq_mod.F =================================================================== --- trunk/LMDZ.MARS/libf/phymars/physiq_mod.F (revision 3316) +++ trunk/LMDZ.MARS/libf/phymars/physiq_mod.F (revision 3325) @@ -2584,5 +2584,6 @@ call pbl_parameters(ngrid,nlayer,ps,zplay,z0, & g,zzlay,zzlev,zu,zv,wstar,hfmax_th,zmax_th,q2,tsurf(:,iflat), - & zh,z_out,n_out,T_out,u_out,ustar,tstar,vhf,vvv) + & zh,zq(:,1,igcm_h2o_vap),qsurf(:,igcm_h2o_ice,iflat),mmean(:,1), + & z_out,n_out,T_out,u_out,ustar,tstar,vhf,vvv) ! pourquoi ustar recalcule ici? fait dans vdifc. Index: trunk/LMDZ.MARS/libf/phymars/vdif_cd_mod.F90 =================================================================== --- trunk/LMDZ.MARS/libf/phymars/vdif_cd_mod.F90 (revision 3316) +++ trunk/LMDZ.MARS/libf/phymars/vdif_cd_mod.F90 (revision 3325) @@ -16,5 +16,5 @@ CONTAINS - SUBROUTINE vdif_cd(ngrid,nlay,nslope,pz0,pg,pz,pp,pu,pv,wstar,pts,ph,virtual,mumean,pqvap,pqsurf,pcdv,pcdh) + SUBROUTINE vdif_cd(ngrid,nlay,nslope,pz0,pg,pz,pp,pu,pv,wstar,pts,ph,mumean,pqvap,pqsurf,write_outputs,pcdv,pcdh) use comcstfi_h @@ -22,7 +22,11 @@ use watersat_mod, only: watersat use lmdz_atke_turbulence_ini, only : smmin, ric, cinf, cepsilon, pr_slope, pr_asym, pr_neut, ri0, ri1, cn, rpi - + use paleoclimate_mod, only: include_waterbuoyancy + use write_output_mod, only: write_output + use comslope_mod, ONLY: iflat IMPLICIT NONE include "callkeys.h" + + !======================================================================= ! @@ -47,5 +51,4 @@ ! pts(ngrid) surface temperature ! ph(ngrid) potential temperature T*(p/ps)^kappa -! virtual Boolean to account for vapor flottability ! mumean Molecular mass of the atmosphere (kg/mol) ! pqvap(ngrid,nlay) Water vapor mixing ratio (kg/kg) to account for vapor flottability @@ -76,8 +79,8 @@ REAL, INTENT(IN) :: pts(ngrid,nslope),ph(ngrid,nlay) ! Surface Temperature and atmospheric temperature (K) REAL, INTENT(IN) :: wstar(ngrid) ! Vertical velocity due to turbulence (m/s) - LOGICAL, INTENT(IN) :: virtual ! Boolean to account for vapor flottability REAL, INTENT(IN) :: mumean(ngrid) ! Molecular mass of the atmosphere (kg/mol) REAL, INTENT(IN) :: pqvap(ngrid,nlay) ! Water vapor mixing ratio (kg/kg) to account for vapor flottability REAL, INTENT(IN) :: pqsurf(ngrid,nslope) ! Water ice frost on the surface (kg/m^2) to account for vapor flottability + LOGICAL, INTENT(IN) :: write_outputs ! write_output in xios or not. ! Outputs: @@ -100,7 +103,8 @@ !$OMP THREADPRIVATE(karman,nu) - REAL rib(ngrid) ! Bulk Richardson number (1) - REAL fm(ngrid) ! stability function for momentum (1) - REAL fh(ngrid) ! stability function for heat (1) + REAL rib(ngrid,nslope) ! Bulk Richardson number [virtual or dry] (1) + REAL rib_dry(ngrid,nslope) ! Dry bulk Richardson number [virtual or dry] (1) + REAL fm(ngrid,nslope) ! stability function for momentum (1) + REAL fh(ngrid,nslope) ! stability function for heat (1) ! Formalism for stability functions fm= 1/phim^2; fh = 1/(phim*phih) where @@ -111,28 +115,29 @@ REAL betam, betah, alphah, bm, bh, lambda - REAL pz0t ! initial thermal roughness length z0t for the iterative scheme (m) - REAL ric_colaitis ! critical richardson number (1) - REAL reynolds(ngrid) ! Reynolds number (1) - REAL prandtl(ngrid) ! Prandtl number (1) - REAL pz0tcomp(ngrid) ! computed z0t during the iterative scheme (m) - REAL ite ! Number of iteration (1) - REAL itemax ! Maximum number of iteration (1) - REAL residual ! Residual between pz0t and pz0tcomp (m) - REAL tol_iter ! Tolerance for the residual to ensure convergence (1= - - REAL zu2(ngrid) ! Near surface winds (m/s) - - REAL cdn(ngrid) ! neutral momentum drag coefficient (1) - REAL chn(ngrid) ! neutral heat drag coefficient (1) - REAL z1z0,z1z0t ! ratios z1/z0 and z1/z0T - REAL z1,zcd0 ! Neutral roughness (m) and Cd/Ch coefficient when call richls is not called - REAL tsurf_v(ngrid,nslope) ! Virtual surface temperature, accounting for vapor flottability - REAL temp_v(ngrid) ! Potential virtual air temperature in the frist layer, accounting for vapor flottability - REAL mu_h2o ! Molecular mass of water vapor (kg/mol) - REAL tol_frost ! Tolerance to consider the effect of frost on the surface - REAL qsat(ngrid) ! saturated mixing ratio of water vapor + REAL pz0t ! initial thermal roughness length z0t for the iterative scheme (m) + REAL ric_colaitis ! critical richardson number (1) + REAL reynolds(ngrid,nslope) ! Reynolds number (1) + REAL prandtl(ngrid) ! Prandtl number (1) + REAL pz0tcomp(ngrid) ! computed z0t during the iterative scheme (m) + REAL z0t(ngrid,nslope) ! computed z0t at the last step the iterative scheme (m) + REAL ite ! Number of iteration (1) + REAL itemax ! Maximum number of iteration (1) + REAL residual ! Residual between pz0t and pz0tcomp (m) + REAL tol_iter ! Tolerance for the residual to ensure convergence (1= + + REAL zu2(ngrid) ! Near surface winds (m/s) + + REAL cdn(ngrid) ! neutral momentum drag coefficient (1) + REAL chn(ngrid) ! neutral heat drag coefficient (1) + REAL z1z0,z1z0t ! ratios z1/z0 and z1/z0T + REAL z1,zcd0 ! Neutral roughness (m) and Cd/Ch coefficient when call richls is not called + REAL tsurf_v(ngrid,nslope) ! Virtual surface temperature, accounting for vapor flottability + REAL temp_v(ngrid) ! Potential virtual air temperature in the frist layer, accounting for vapor flottability + REAL mu_h2o ! Molecular mass of water vapor (kg/mol) + REAL tol_frost ! Tolerance to consider the effect of frost on the surface + REAL qsat(ngrid) ! saturated mixing ratio of water vapor - REAL sm ! Stability function computed with the ATKE scheme - REAL f_ri_cd_min ! minimum of the stability function with the ATKE scheme + REAL sm ! Stability function computed with the ATKE scheme + REAL f_ri_cd_min ! minimum of the stability function with the ATKE scheme ! Code: @@ -144,10 +149,12 @@ mu_h2o = 18e-3 tol_frost = 1e-4 - reynolds(:) = 0. + reynolds(:,:) = 0. pz0t = 0. pz0tcomp(:) = 0.1*pz0(:) - rib(:) = 0. + rib(:,:) = 0. + rib_dry(:,:) = 0. pcdv(:,:) = 0. pcdh(:,:) = 0. + z0t(:,:) = 0. f_ri_cd_min = 0.01 ! this formulation assumes alphah=1., implying betah=betam @@ -162,5 +169,5 @@ ric_colaitis = betah/(betam**2) - IF(virtual) then + IF(include_waterbuoyancy) then temp_v(:) = ph(:,1)*(1.+pqvap(:,1)/(mu_h2o/mumean(:)))/(1.+pqvap(:,1)) DO islope = 1,nslope @@ -210,9 +217,11 @@ IF(turb_resolved) zu2(ig)=MAX(zu2(ig),1.) ! Richardson number formulation proposed by D.E. England et al. (1995) - rib(ig) = (pg/tsurf_v(ig,islope))*sqrt(pz(ig,1)*pz0(ig))*(((log(pz(ig,1)/pz0(ig)))**2)/(log(pz(ig,1)/pz0t)))*(temp_v(ig)-tsurf_v(ig,islope))/zu2(ig) + rib(ig,islope) = (pg/tsurf_v(ig,islope))*sqrt(pz(ig,1)*pz0(ig))*(((log(pz(ig,1)/pz0(ig)))**2)/(log(pz(ig,1)/pz0t)))*(temp_v(ig)-tsurf_v(ig,islope))/zu2(ig) + rib_dry(ig,islope) = (pg/pts(ig,islope))*sqrt(pz(ig,1)*pz0(ig))*(((log(pz(ig,1)/pz0(ig)))**2)/(log(pz(ig,1)/pz0t)))*(ph(ig,1)-pts(ig,islope))/zu2(ig) ELSE print*,'warning, infinite Richardson at surface' print*,pu(ig,1),pv(ig,1) - rib(ig) = ric_colaitis + rib(ig,islope) = ric_colaitis + rib_dry(ig,islope) = ric_colaitis ENDIF @@ -221,42 +230,42 @@ IF(callatke) THEN ! Computation following parameterizaiton from ATKE - IF(rib(ig) .gt. 0) THEN - ! STABLE BOUNDARY LAYER : - sm = MAX(smmin,cn*(1.-rib(ig)/ric)) + IF(rib(ig,islope) .gt. 0) THEN + ! STABLE BOUNDARY LAYER :include_waterbuoyancy + sm = MAX(smmin,cn*(1.-rib(ig,islope)/ric)) ! prandlt expression from venayagamoorthy and stretch 2010, Li et al 2019 - prandtl(ig) = pr_neut*exp(-pr_slope/pr_neut*rib(ig)+rib(ig)/pr_neut) + rib(ig) * pr_slope - fm(ig) = max((sm**(3./2.) * sqrt(cepsilon) * (1 - rib(ig) / prandtl(ig))**(1./2.)),f_ri_cd_min) - fh(ig) = max((fm(ig) / prandtl(ig)), f_ri_cd_min) + prandtl(ig) = pr_neut*exp(-pr_slope/pr_neut*rib(ig,islope)+rib(ig,islope)/pr_neut) + rib(ig,islope) * pr_slope + fm(ig,islope) = max((sm**(3./2.) * sqrt(cepsilon) * (1 - rib(ig,islope) / prandtl(ig))**(1./2.)),f_ri_cd_min) + fh(ig,islope) = max((fm(ig,islope) / prandtl(ig)), f_ri_cd_min) ELSE ! UNSTABLE BOUNDARY LAYER : - sm = 2./rpi * (cinf - cn) * atan(-rib(ig)/ri0) + cn - prandtl(ig) = -2./rpi * (pr_asym - pr_neut) * atan(rib(ig)/ri1) + pr_neut - fm(ig) = MAX((sm**(3./2.) * sqrt(cepsilon) * (1 - rib(ig) / prandtl(ig))**(1./2.)),f_ri_cd_min) - fh(ig) = MAX((fm(ig) / prandtl(ig)), f_ri_cd_min) + sm = 2./rpi * (cinf - cn) * atan(-rib(ig,islope)/ri0) + cn + prandtl(ig) = -2./rpi * (pr_asym - pr_neut) * atan(rib(ig,islope)/ri1) + pr_neut + fm(ig,islope) = MAX((sm**(3./2.) * sqrt(cepsilon) * (1 - rib(ig,islope) / prandtl(ig))**(1./2.)),f_ri_cd_min) + fh(ig,islope) = MAX((fm(ig,islope) / prandtl(ig)), f_ri_cd_min) ENDIF ! Rib < 0 ELSE ! Computation following parameterizaiton from from D.E. England et al. (95) - IF (rib(ig) .gt. 0.) THEN + IF (rib(ig,islope) .gt. 0.) THEN ! STABLE BOUNDARY LAYER : prandtl(ig) = 1. - IF(rib(ig) .lt. ric_colaitis) THEN + IF(rib(ig,islope) .lt. ric_colaitis) THEN ! Assuming alphah=1. and bh=bm for stable conditions : - fm(ig) = ((ric_colaitis-rib(ig))/ric_colaitis)**2 - fh(ig) = ((ric_colaitis-rib(ig))/ric_colaitis)**2 + fm(ig,islope) = ((ric_colaitis-rib(ig,islope))/ric_colaitis)**2 + fh(ig,islope) = ((ric_colaitis-rib(ig,islope))/ric_colaitis)**2 ELSE ! For Ri>Ric, we consider Ri->Infinity => no turbulent mixing at surface - fm(ig) = 1. - fh(ig) = 1. + fm(ig,islope) = 1. + fh(ig,islope) = 1. ENDIF ELSE ! UNSTABLE BOUNDARY LAYER : - fm(ig) = sqrt(1.-lambda*bm*rib(ig)) - fh(ig) = (1./alphah)*((1.-lambda*bh*rib(ig))**0.5)*(1.-lambda*bm*rib(ig))**0.25 - prandtl(ig) = alphah*((1.-lambda*bm*rib(ig))**0.25)/((1.-lambda*bh*rib(ig))**0.5) + fm(ig,islope) = sqrt(1.-lambda*bm*rib(ig,islope)) + fh(ig,islope) = (1./alphah)*((1.-lambda*bh*rib(ig,islope))**0.5)*(1.-lambda*bm*rib(ig,islope))**0.25 + prandtl(ig) = alphah*((1.-lambda*bm*rib(ig,islope))**0.25)/((1.-lambda*bh*rib(ig,islope))**0.5) ENDIF ! Rib < 0 ENDIF ! atke ! Recompute the reynolds and z0t - reynolds(ig) = karman*sqrt(fm(ig))*sqrt(zu2(ig))*pz0(ig)/(log(z1z0)*nu) - pz0tcomp(ig) = pz0(ig)*exp(-karman*7.3*(reynolds(ig)**0.25)*(prandtl(ig)**0.5)+5*karman) + reynolds(ig,islope) = karman*sqrt(fm(ig,islope))*sqrt(zu2(ig))*pz0(ig)/(log(z1z0)*nu) + pz0tcomp(ig) = pz0(ig)*exp(-karman*7.3*(reynolds(ig,islope)**0.25)*(prandtl(ig)**0.5)+5*karman) residual = abs(pz0t-pz0tcomp(ig)) ite = ite+1 @@ -264,10 +273,29 @@ ! Compute the coefficients Cdv; Cdh : neutral coefficient x stability functions - pcdv(ig,islope)=cdn(ig)*fm(ig) - pcdh(ig,islope)=chn(ig)*fh(ig) + pcdv(ig,islope)=cdn(ig)*fm(ig,islope) + pcdh(ig,islope)=chn(ig)*fh(ig,islope) + z0t(ig,islope) = pz0tcomp(ig) pz0t = 0. ! for next grid point ENDDO ! of ngrid enddo ENDIF !of if call richardson surface layer + + IF(write_outputs) then + + call write_output("rib_dry_vdif_cd","Dry Richardson number in vdif_cd","1",rib_dry(:,iflat)) + + call write_output("rib_vdif_cd","Richardson number in vdif_cd","1",rib(:,iflat)) + + call write_output("fm_vdif_cd","Momentum stability function in vdif_cd","1",fm(:,iflat)) + + call write_output("fh_vdif_cd","Heat stability function in vdif_cd","1",fh(:,iflat)) + + call write_output("z0t_vdif_cd","Thermal roughness length in vdif_cd","m",z0t(:,iflat)) + + call write_output("z0_vdif_cd","Momentum roughness length in vdif_cd","m",pz0(:)) + + call write_output("Reynolds_vdif_cd","Reynolds number in vdif_cd","1",reynolds(:,iflat)) + + ENDIF RETURN Index: trunk/LMDZ.MARS/libf/phymars/vdifc_mod.F =================================================================== --- trunk/LMDZ.MARS/libf/phymars/vdifc_mod.F (revision 3316) +++ trunk/LMDZ.MARS/libf/phymars/vdifc_mod.F (revision 3325) @@ -235,6 +235,4 @@ LOGICAL :: writeoutput ! boolean to say to soilexchange.F if we are at the last iteration and thus if he can write in the diagsoil REAL, INTENT(out) :: pore_icefraction(ngrid,nsoil,nslope) ! ice filling fraction in the pores -!! Water buyoncy - LOGICAL :: virtual !! Subsurface ice interactions @@ -309,5 +307,4 @@ zq1temp_regolith(1:ngrid)=0 zdqsdif_tot(1:ngrid)=0 - virtual = .false. pore_icefraction(:,:,:) = 0. zdqsdif_ssi_atm(:,:) = 0. @@ -434,6 +431,6 @@ CALL vdif_cd(ngrid,nlay,nslope,pz0,g,pzlay,pplay,pu,pv,wstar, - & ptsrf,ph,virtual,mmean(:,1),zq(:,:,igcm_h2o_vap), - & pqsurf(:,igcm_h2o_ice,:), + & ptsrf,ph,mmean(:,1),zq(:,:,igcm_h2o_vap), + & pqsurf(:,igcm_h2o_ice,:), .false., & zcdv_true,zcdh_true) @@ -612,6 +609,6 @@ CALL vdif_cd(ngrid,nlay,nslope,pz0,g,pzlay,pplay,zu,zv,wstar, - & ptsrf,ph,virtual,mmean(:,1),zq(:,:,igcm_h2o_vap), - & pqsurf(:,igcm_h2o_ice,:), + & ptsrf,ph,mmean(:,1),zq(:,:,igcm_h2o_vap), + & pqsurf(:,igcm_h2o_ice,:), .true., & zcdv_true,zcdh_true)