diff --git a/Registry/Registry.EM_COMMON b/Registry/Registry.EM_COMMON index 889a92854b..8ff2204ac7 100644 --- a/Registry/Registry.EM_COMMON +++ b/Registry/Registry.EM_COMMON @@ -895,13 +895,23 @@ state real dfi_SMFR3D ilj misc 1 Z r "SMFR3D_df state real dfi_KEEPFR3DFLAG ilj misc 1 Z r "KEEPFR3DFLAG_dfi" "FLAG - 1. FROZEN SOIL YES, 0 - NO" "" # urban state variables -state real TSK_RURAL ij misc 1 - rd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "TSK_RURAL" "TSK for rural fraction" "K" state real TR_URB2D ij misc 1 - rd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "TR_URB" "URBAN ROOF SKIN TEMPERATURE" "K" state real TGR_URB2D ij misc 1 - rd=(interp_mask_land_field:lu_index)u=(copy_fcnm) "TGR_URB" "URBAN GREEN ROOF SKIN TEMPERATURE" "K" state real TB_URB2D ij misc 1 - rd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "TB_URB" "URBAN WALL SKIN TEMPERATURE" "K" state real TG_URB2D ij misc 1 - rd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "TG_URB" "URBAN ROAD SKIN TEMPERATURE" "K" state real TC_URB2D ij misc 1 - rd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "TC_URB" "URBAN CANOPY TEMPERATURE" "K" state real QC_URB2D ij misc 1 - rd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "QC_URB" "URBAN CANOPY HUMIDITY" "kg kg{-1}" +#CORDEX +state real TSK_RURAL ij misc 1 - rh "TSK_RURAL" "TSK for rural fraction" "K" +state real T2R_URB2D ij misc 1 - rh "T2R_URB2D" "2M ROOF TEMPERATURE" "K" +state real T2G_URB2D ij misc 1 - rh "T2G_URB2D" "2M ROAD TEMPERATURE" "K" +state real TGROUND_URB2D ij misc 1 - rh "TGROUND_URB2D" "ROAD SURFACE TEMPERATURE" "K" +state real TROOF_URB2D ij misc 1 - rh "TROOF_URB2D" "ROOF SURFACE TEMPERATURE" "K" +state real T2VEG_URB2D ij misc 1 - rh "T2VEG_URB2D" "2M VEG TEMPERATURE" "K" +state real TCAN_URB2D ij misc 1 - rh "TCAN_URB2D" "URBAN CANOPY TEMP" "K" +state real UCAN_URB2D ij misc 1 - rh "UCAN_URB2D" "URBAN CANOPY WIND" "M/S" +state real QCAN_URB2D ij misc 1 - rh "QCAN_URB2D" "URBAN CANOPY HUMIDITY" "Kg kg{-1}" +#CORDEX state real UC_URB2D ij misc 1 - rd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "UC_URB" "URBAN CANOPY WIND" "m s{-1}" state real XXXR_URB2D ij misc 1 - rd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "XXXR_URB" "M-O LENGTH ABOVE URBAN ROOF" "dimensionless" state real XXXB_URB2D ij misc 1 - rd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "XXXB_URB" "M-O LENGTH ABOVE URBAN WALL" "dimensionless" @@ -3094,7 +3104,9 @@ package sfclayscheme sf_sfclay_physics==91 - - package noahucmscheme sf_urban_physics==1 - state:trb_urb4d,tw1_urb4d,tw2_urb4d,tgb_urb4d,sfw1_urb3d,sfw2_urb3d,sfr_urb3d,sfg_urb3d,a_u_bep,a_v_bep,a_t_bep,a_q_bep,a_e_bep,b_u_bep,b_v_bep,b_t_bep,b_q_bep,b_e_bep,dlg_bep,dl_u_bep,sf_bep,vl_bep,mh_urb2d,stdh_urb2d,lf_urb2d,lp_urb2d,hgt_urb2d,lb_urb2d,tgr_urb2d,cmcr_urb2d,drelr_urb2d,drelb_urb2d,drelg_urb2d,flxhumr_urb2d,flxhumb_urb2d,flxhumg_urb2d,tgrl_urb3d,smr_urb3d,cmgr_sfcdif,chgr_sfcdif,trl_urb3d,tgl_urb3d,tbl_urb3d package bepscheme sf_urban_physics==2 - state:a_u_bep,a_v_bep,a_t_bep,a_q_bep,a_e_bep,b_u_bep,b_v_bep,b_t_bep,b_q_bep,b_e_bep,dlg_bep,dl_u_bep,sf_bep,vl_bep,trb_urb4d,tw1_urb4d,tw2_urb4d,tgb_urb4d,sfw1_urb3d,sfw2_urb3d,sfr_urb3d,sfg_urb3d,hi_urb2d,lp_urb2d,hgt_urb2d,lb_urb2d,trl_urb3d,tgl_urb3d,tbl_urb3d,tsk_rural -package bep_bemscheme sf_urban_physics==3 - state:a_u_bep,a_v_bep,a_t_bep,a_q_bep,a_e_bep,b_u_bep,b_v_bep,b_t_bep,b_q_bep,b_e_bep,dlg_bep,dl_u_bep,sf_bep,vl_bep,trb_urb4d,tw1_urb4d,tw2_urb4d,tgb_urb4d,tlev_urb3d,qlev_urb3d,tw1lev_urb3d,tw2lev_urb3d,tglev_urb3d,tflev_urb3d,sf_ac_urb3d,lf_ac_urb3d,cm_ac_urb3d,sfvent_urb3d,lfvent_urb3d,sfwin1_urb3d,sfwin2_urb3d,sfw1_urb3d,sfw2_urb3d,sfr_urb3d,sfg_urb3d,hi_urb2d,lp_urb2d,hgt_urb2d,lb_urb2d,trl_urb3d,tgl_urb3d,tbl_urb3d,tsk_rural,ep_pv_urb3d,t_pv_urb3d,trv_urb4d,qr_urb4d,qgr_urb3d,tgr_urb3d,drain_urb4d,draingr_urb3d,sfrv_urb3d,lfrv_urb3d,dgr_urb3d,dg_urb3d,lfr_urb3d,lfg_urb3d +#CORDEX +package bep_bemscheme sf_urban_physics==3 - state:a_u_bep,a_v_bep,a_t_bep,a_q_bep,a_e_bep,b_u_bep,b_v_bep,b_t_bep,b_q_bep,b_e_bep,dlg_bep,dl_u_bep,sf_bep,vl_bep,trb_urb4d,tw1_urb4d,tw2_urb4d,tgb_urb4d,tlev_urb3d,qlev_urb3d,tw1lev_urb3d,tw2lev_urb3d,tglev_urb3d,tflev_urb3d,sf_ac_urb3d,lf_ac_urb3d,cm_ac_urb3d,sfvent_urb3d,lfvent_urb3d,sfwin1_urb3d,sfwin2_urb3d,sfw1_urb3d,sfw2_urb3d,sfr_urb3d,sfg_urb3d,hi_urb2d,lp_urb2d,hgt_urb2d,lb_urb2d,trl_urb3d,tgl_urb3d,tbl_urb3d,tsk_rural,ep_pv_urb3d,t_pv_urb3d,trv_urb4d,qr_urb4d,qgr_urb3d,tgr_urb3d,drain_urb4d,draingr_urb3d,sfrv_urb3d,lfrv_urb3d,dgr_urb3d,dg_urb3d,lfr_urb3d,lfg_urb3d,t2r_urb2d,t2g_urb2d,t2veg_urb2d,tcan_urb2d,ucan_urb2d,qcan_urb2d,troof_urb2d,tground_urb2d +#CORDEX package nolsmscheme sf_surface_physics==0 - - diff --git a/dyn_em/module_first_rk_step_part1.F b/dyn_em/module_first_rk_step_part1.F index f5eb26734d..d3db03ab50 100644 --- a/dyn_em/module_first_rk_step_part1.F +++ b/dyn_em/module_first_rk_step_part1.F @@ -775,6 +775,15 @@ SUBROUTINE first_rk_step_part1 ( grid , config_flags & & ,TG_URB2D=grid%tg_urb2d & !H urban & ,TC_URB2D=grid%tc_urb2d ,QC_URB2D=grid%qc_urb2d & !H urban & ,UC_URB2D=grid%uc_urb2d & !H urban +!CORDEX + & ,T2G_URB2D=grid%t2g_urb2d,T2R_URB2D=grid%t2r_urb2d & + & ,T2VEG_URB2D=grid%t2veg_urb2d & + & ,TCAN_URB2D=grid%tcan_urb2d & + & ,UCAN_URB2D=grid%ucan_urb2d & + & ,QCAN_URB2D=grid%qcan_urb2d & + & ,TGROUND_URB2D=grid%tground_urb2d & + & ,TROOF_URB2D=grid%troof_urb2d & +!CORDEX & ,XXXR_URB2D=grid%xxxr_urb2d & & ,XXXB_URB2D=grid%xxxb_urb2d & !H urban & ,XXXG_URB2D=grid%xxxg_urb2d & diff --git a/dyn_em/module_initialize_real.F b/dyn_em/module_initialize_real.F index 96629232bf..3c7a8d557b 100644 --- a/dyn_em/module_initialize_real.F +++ b/dyn_em/module_initialize_real.F @@ -3117,9 +3117,20 @@ SUBROUTINE init_domain_rk ( grid & DO j = jts , MIN(jde-1,jte) DO i = its , MIN(ide-1,ite) IF ( skip_middle_points_t ( ids , ide , jds , jde , i , j , em_width , hold_ups ) ) CYCLE - DO k = 1, 15 - grid%HI_URB2D(i,k,j) = grid%URB_PARAM(i,k+117,j) +!CORDEX +! DO k = 1, 15 +! grid%HI_URB2D(i,k,j) = grid%URB_PARAM(i,k+117,j) +! END DO +! force to one level + DO k = 1, 6!15 + grid%HI_URB2D(i,k,j) =0.! grid%URB_PARAM(i,k+117,j) END DO + k=int(grid%HGT_URB2D(i,j)/10.)+1 + if(k.lt.1)k=1 + if(k.gt.4)k=4 + grid%HI_URB2D(i,k,j) =1. +! force to one level + END DO END DO ENDIF diff --git a/phys/module_sf_bep_bem.F b/phys/module_sf_bep_bem.F index f249700a08..55835b2da7 100644 --- a/phys/module_sf_bep_bem.F +++ b/phys/module_sf_bep_bem.F @@ -42,12 +42,12 @@ MODULE module_sf_bep_bem integer ngb_u !Number of grid levels in the ground below building (BEM) parameter (ngb_u=10) - +! real dz_u ! Urban grid resolution - parameter (dz_u=5.) + parameter (dz_u=10.) !CORDEX integer nbui_max !maximum number of types of buildings in an urban class - parameter (nbui_max=15) !must be less or equal than nz_um + parameter (nbui_max=1) !CORDEX!must be less or equal than nz_um real h_water @@ -117,6 +117,10 @@ subroutine BEP_BEM(FRC_URB2D,UTYPE_URB2D,itimestep,dz8w,dt,u_phy,v_phy, & a_u,a_v,a_t,a_e,b_u,b_v, & b_t,b_e,b_q,dlg,dl_u,sf,vl, & rl_up,rs_abs,emiss,grdflx_urb,qv_phy, & +!CORDEX + tground_urb2d,troof_urb2d,tcan_urb2d, & + t2g_urb2d,t2r_urb2d,ucan_urb2d,qcan_urb2d, & +!CORDEX ids,ide, jds,jde, kds,kde, & ims,ime, jms,jme, kms,kme, & its,ite, jts,jte, kts,kte) @@ -178,6 +182,15 @@ subroutine BEP_BEM(FRC_URB2D,UTYPE_URB2D,itimestep,dz8w,dt,u_phy,v_phy, & REAL, DIMENSION( ims:ime, 1:urban_map_zdf, jms:jme ), INTENT(INOUT) :: drain_urb4d REAL, DIMENSION( ims:ime, jms:jme ), INTENT(IN) :: rainbl REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: draingr_urb3d +!CORDEX + REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: troof_urb2d + REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: tground_urb2d + REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: tcan_urb2d + REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: t2g_urb2d + REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: t2r_urb2d + REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: ucan_urb2d + REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: qcan_urb2d +!CORDEX !New variables used for BEM REAL, DIMENSION( ims:ime, kms:kme, jms:jme ):: qv_phy REAL, DIMENSION( ims:ime, 1:urban_map_bd, jms:jme ), INTENT(INOUT) :: tlev_urb3d @@ -945,7 +958,10 @@ subroutine BEP_BEM(FRC_URB2D,UTYPE_URB2D,itimestep,dz8w,dt,u_phy,v_phy, & dlg1D,dl_u1D,sf1D,vl1D,rl_up(ix,iy), & rs_abs(ix,iy),emiss(ix,iy),grdflx_urb(ix,iy), & qv1D,tlev1D,qlev1D,sflev1D,lflev1D,consumlev1D, & - eppvlev1D,tpvlev1D,sfvlev1D,lfvlev1D,twlev1D,tglev1D,tflev1D,sfwin1D,tair1D,sfr_indoor1D,sfrpv1D,gfr1D) + eppvlev1D,tpvlev1D,sfvlev1D,lfvlev1D,twlev1D,tglev1D,tflev1D,sfwin1D,tair1D,sfr_indoor1D,sfrpv1D,gfr1D, & +!CORDEX + t2g_urb2d(ix,iy),t2r_urb2d(ix,iy),tcan_urb2d(ix,iy),ucan_urb2d(ix,iy),qcan_urb2d(ix,iy)) +!CORDEX do ibui=1,nbui_max !type of building do iz=1,nz_um !vertical levels @@ -987,22 +1003,31 @@ subroutine BEP_BEM(FRC_URB2D,UTYPE_URB2D,itimestep,dz8w,dt,u_phy,v_phy, & enddo enddo - +!CORDEX - initalize to zero the ground and roof average surface temperature + tground_urb2d(ix,iy)=0. + troof_urb2d(ix,iy)=0. +!CORDEX -end do id=1,ndm do ig=1,ng_u - tgb_urb4d(ix,ind_gd(ig,id),iy)=tg1D(id,ig) enddo +!CORDEX - averaging it over the different street directions + tground_urb2d(ix,iy)=tground_urb2d(ix,iy)+tg1D(id,ng_u)/ndm +!CORDEX -end do iz_u=1,nz_um do ir=1,nwr_u trb_urb4d(ix,ind_zrd(iz_u,ir,id),iy)=tr1D(id,iz_u,ir) enddo - if(gr_flag_u.eq.1)then +!CORDEX - averaging it over the different street directions, and different heights + troof_urb2d(ix,iy)=troof_urb2d(ix,iy)+ss(iz_u)*tr1D(id,iz_u,nwr_u)/ndm + +!CORDEX -end + if(gr_flag_u.eq.1)then do ir=1,ngr_u trv_urb4d(ix,ind_grd(iz_u,ir,id),iy)=trv1D(id,iz_u,ir) qr_urb4d(ix,ind_grd(iz_u,ir,id),iy)=qr1D(id,iz_u,ir) enddo - endif + endif enddo enddo ! @@ -1182,7 +1207,10 @@ subroutine BEP1D(itimestep,ix,iy,iurb,kms,kme,kts,kte,z,dt,ua,va,pt,da,pr,pt0, b_u,b_v,b_t,b_ac,b_e,b_q, & dlg,dl_u,sf,vl,rl_up,rs_abs,emiss,grdflx_urb, & qv,tlev,qlev,sflev,lflev,consumlev, & - eppvlev,tpvlev,sfvlev,lfvlev,twlev,tglev,tflev,sfwin,tmp_u,sfr_indoor,sfrpv,gfr) + eppvlev,tpvlev,sfvlev,lfvlev,twlev,tglev,tflev,sfwin,tmp_u,sfr_indoor,sfrpv,gfr, & +!CORDEX + t2g_urb2d,t2r_urb2d,tcan_urb2d,ucan_urb2d,qcan_urb2d) +!CORDEX ! print*,'SFR_AFT',sfr(id,iz) @@ -1228,7 +1256,9 @@ subroutine BEP1D(itimestep,ix,iy,iurb,kms,kme,kts,kte,z,dt,ua,va,pt,da,pr,pt0, real emg ! Emissivity of ground real emw ! Emissivity of wall real emr ! Emissivity of roof - +!CORDEX + real tcan_urb2d,t2g_urb2d,t2r_urb2d,ucan_urb2d,qcan_urb2d +!CORDEX ! fww,fwg,fgw,fsw,fsg are the view factors used to compute the long and ! short wave radation. @@ -1479,6 +1509,9 @@ subroutine BEP1D(itimestep,ix,iy,iurb,kms,kme,kts,kte,z,dt,ua,va,pt,da,pr,pt0, real tr_av(ndm,nz_um) real tr_avb(ndm,nbui_max) real sfr_avb(ndm,nbui_max) +!CORDEX + real zcenter, zcenter0,tc0,tc1,uc0,uc1,qc0,qc1 +!CORDEX ! ---------------------------------------------------------------------- ! END VARIABLES DEFINITIONS ! ---------------------------------------------------------------------- @@ -1515,6 +1548,38 @@ subroutine BEP1D(itimestep,ix,iy,iurb,kms,kme,kts,kte,z,dt,ua,va,pt,da,pr,pt0, call interpol(kms,kme,kts,kte,nzu,z,z_u,pr,pr_u) call interpol(kms,kme,kts,kte,nzu,z,z_u,da,da_u) call interpol(kms,kme,kts,kte,nzu,z,z_u,qv,qv_u) +!CORDEX - here I am computing the value of temperature at 2m and wind speed at 10m + tcan_urb2d=0. + ucan_urb2d=0. + qcan_urb2d=0. + do iz=1,nz_um-1 + if(z_u(iz).le.2.and.z_u(iz+1).gt.2)then + if(iz.eq.1)then + tcan_urb2d=pt_u(1)*(pr_u(1)/1.e5)**rcp_u + qcan_urb2d=qv_u(1) + else + zcenter0=(z_u(iz-1)+z_u(iz))/2. + tc0=pt_u(iz-1)*(pr_u(iz-1)/1.e5)**rcp_u + tc1=pt_u(iz)*(pr_u(iz)/1.e5)**rcp_u + qc0=qv_u(iz-1) + qc1=qv_u(iz) + tcan_urb2d=((zcenter-2.)*tc0+(2.-zcenter0)*tc1)/(zcenter-zcenter0) + qcan_urb2d=((zcenter-2.)*qc0+(2.-zcenter0)*qc1)/(zcenter-zcenter0) + endif + endif + if(z_u(iz).le.10.and.z_u(iz+1).gt.10)then + if(iz.eq.1)then + ucan_urb2d=sqrt(ua_u(iz)**2.+va_u(iz)**2.) + else + zcenter0=(z_u(iz-1)+z_u(iz))/2. + uc0=sqrt(ua_u(iz-1)**2.+va_u(iz-1)**2.) + uc1=sqrt(ua_u(iz)**2.+va_u(iz)**2.) + ucan_urb2d=((zcenter-10.)*uc0+(10.-zcenter0)*uc1)/(zcenter-zcenter0) + endif + endif + enddo +!CORDEX + ! Compute the modification of the radiation due to the buildings @@ -1585,6 +1650,8 @@ subroutine BEP1D(itimestep,ix,iy,iurb,kms,kme,kts,kte,z,dt,ua,va,pt,da,pr,pt0, call surf_temp(ndu,pr_u,dt, & rld,rsg,rlg, & tg,alag,csg,emg,albg,ptg,sfg,lfg,gfg) + if(ix.eq.50.and.iy.eq.42)write(98,'(10(2x,f15.8))')sfg(1),gfg(1),rsg(1),rlg(1),ptg(1) + if(ix.eq.50.and.iy.eq.42)write(99,'(10(2x,f15.8))')sfg(2),gfg(2),rsg(2),rlg(2),ptg(2) if(gr_flag.eq.1)then if(gr_frac_roof.gt.0.)then hfgr=0. @@ -1729,6 +1796,7 @@ subroutine BEP1D(itimestep,ix,iy,iurb,kms,kme,kts,kte,z,dt,ua,va,pt,da,pr,pt0, sfvlev1D,lfvlev1D,hfgrb(1,ibui),tr_avb(1,ibui), & tpv(ibui),sfpv(ibui),sfr_indoor(ibui)) + if(ix.eq.50.and.iy.eq.42)write(88,'(10(2x,f15.8))')sfr_avb(1,1),gfrb(1,1),rs,rld,ptr(1,1) ! !Temporal modifications @@ -1845,7 +1913,10 @@ subroutine BEP1D(itimestep,ix,iy,iurb,kms,kme,kts,kte,z,dt,ua,va,pt,da,pr,pt0, uhb_u,vhb_u,thb_u,ehb_u,ss,dt,sfw,sfg,sfr,sfrpv,sfrv,lfrv, & dgr,dg,lfr,lfg, & sfwin,pb,bs_u,dz_u,sflev,lflev,sfvlev,lfvlev,tvb_ac,ix,iy,rsg,rs,qr,gr_frac_roof, & - pv_frac_roof,gr_flag,gr_type) + pv_frac_roof,gr_flag,gr_type, & +!CORDEX + t2g_urb2d,t2r_urb2d) +!CORDEX @@ -1955,6 +2026,7 @@ subroutine param(iurb,nzu,nzurb,nzurban,ndu, & ! LOCAL: ! ---------------------------------------------------------------------- integer id,ig,ir,iw,iz,iflo,ihu + real sstot ! ---------------------------------------------------------------------- ! END VARIABLES DEFINITIONS ! ---------------------------------------------------------------------- @@ -2601,7 +2673,10 @@ subroutine buildings(iurb,nd,nz,z0,cdrag,ua_u,va_u,pt_u,pt0_u, & uhb_u,vhb_u,thb_u,ehb_u,ss,dt,sfw,sfg,sfr,sfrpv,sfrv,lfrv, & dgr,dg,lfr,lfg, & sfwin,pb,bs_u,dz_u,sflev,lflev,sfvlev,lfvlev,tvb_ac,ix,iy,rsg,rs,qr,gr_frac_roof, & - pv_frac_roof,gr_flag,gr_type) + pv_frac_roof,gr_flag,gr_type, & +!CORDEX + t2g_urb2d,t2r_urb2d) +!CORDEX ! ---------------------------------------------------------------------- ! This routine computes the sources or sinks of the different quantities @@ -2644,7 +2719,10 @@ subroutine buildings(iurb,nd,nz,z0,cdrag,ua_u,va_u,pt_u,pt0_u, & real trv(ndm,nz_um,ngr_u) ! Ground Soil Moisture real roof_frac real road_frac -! +!CORDEX + real t2g_urb2d,t2r_urb2d + real t2m +!CORDEX !New variables (BEM) ! real bs_u(ndm,nurbm) ! Building width [m] @@ -2722,6 +2800,8 @@ subroutine buildings(iurb,nd,nz,z0,cdrag,ua_u,va_u,pt_u,pt0_u, & real qr_tmp(ngr_u) data rsv /0,1/ real fh,ric,utot + !local + real sstot !------------------------------------------------------------------ ! END VARIABLES DEFINITIONS ! ---------------------------------------------------------------------- @@ -2763,11 +2843,19 @@ subroutine buildings(iurb,nd,nz,z0,cdrag,ua_u,va_u,pt_u,pt0_u, & enddo ! Calculation at the ground surfaces +!CORDEX + t2g_urb2d=0. + t2r_urb2d=0. +!CORDEX do id=1,nd call flux_flat(dz,z0(id,1),ua_u(1),va_u(1),pt_u(1),pt0_u(1), & ptg(id),qv_u(1),uhb(id,1), & - vhb(id,1),sfg(id),lfg(id),ehb(id,1),da_u(1),pr_u(1)) + vhb(id,1),sfg(id),lfg(id),ehb(id,1),da_u(1),pr_u(1), & +!CORDEX + t2m) + t2g_urb2d=t2g_urb2d+t2m/ndm +!CORDEX if(dg(id).gt.0)then wfg=dg(id)/dgmax lfg(id)=-da_u(1)*latent*(-(wfg*lfg(id))/(da_u(1)*latent)) @@ -2780,49 +2868,53 @@ subroutine buildings(iurb,nd,nz,z0,cdrag,ua_u,va_u,pt_u,pt0_u, & uhb_u(id,1)=uhb(id,1) ehb_u(id,1)=ehb(id,1) qhb_u(id,1)=-lfg(id)/(da_u(1)*latent) + do iz=2,nz - if(ss(iz).gt.0)then - + if(ss(iz).gt.0)then call flux_flat(dz,z0(id,iz),ua_u(iz),& va_u(iz),pt_u(iz),pt0_u(iz), & ptr(id,iz),qv_u(iz),uhb(id,iz), & - vhb(id,iz),sfr(id,iz),lfr(id,iz),ehb(id,iz),da_u(iz),pr_u(iz)) - if(dgr(id,iz).gt.0)then - wfr=dgr(id,iz)/drmax - lfr(id,iz)=-da_u(iz)*latent*(-(wfr*lfr(id,iz))/(da_u(iz)*latent)) - else - lfr(id,iz)=0. - endif - if(gr_flag.eq.1.and.gr_frac_roof.gt.0.)then - do il=1,ngr_u - qr_tmp(il)=qr(id,iz,il) - enddo + vhb(id,iz),sfr(id,iz),lfr(id,iz),ehb(id,iz),da_u(iz),pr_u(iz), & +!CORDEX + t2m) + t2r_urb2d=t2r_urb2d+t2m/ndm*ss(iz) +!CORDEX + if(dgr(id,iz).gt.0)then + wfr=dgr(id,iz)/drmax + lfr(id,iz)=-da_u(iz)*latent*(-(wfr*lfr(id,iz))/(da_u(iz)*latent)) + else + lfr(id,iz)=0. + endif + if(gr_flag.eq.1.and.gr_frac_roof.gt.0.)then + do il=1,ngr_u + qr_tmp(il)=qr(id,iz,il) + enddo call flux_flat_roof(dz,z0v,ua_u(iz),va_u(iz),pt_u(iz),pt0_u(iz), & ptrv(id,iz),uhbv(id,iz), & vhbv(id,iz),sfrv(id,iz),lfrv(id,iz),ehbv(id,iz),da_u(iz),qv_u(iz),pr_u(iz),rs,qr_tmp,resg,rsveg,f1,f2,f3,f4,gr_type,pv_frac_roof) - sfr(id,iz)=sfr(id,iz)+pv_frac_roof*sfrpv(id,iz) - thb_u(id,iz)=-((1.-gr_frac_roof)*sfr(id,iz)+gr_frac_roof*sfrv(id,iz))/(da_u(iz)*cp_u) - vhb_u(id,iz)=(1.-gr_frac_roof)*vhb(id,iz)+gr_frac_roof*vhbv(id,iz) - uhb_u(id,iz)=(1.-gr_frac_roof)*uhb(id,iz)+gr_frac_roof*uhbv(id,iz) - ehb_u(id,iz)=(1.-gr_frac_roof)*ehb(id,iz)+gr_frac_roof*ehbv(id,iz) - qhb_u(id,iz)=-(gr_frac_roof*lfrv(id,iz)+(1.-gr_frac_roof)*lfr(id,iz))/(da_u(iz)*latent) - sfr(id,iz)=sfr(id,iz)-pv_frac_roof*sfrpv(id,iz) - else - sfr(id,iz)=sfr(id,iz)+pv_frac_roof*sfrpv(id,iz) - thb_u(id,iz)=-sfr(id,iz)/(da_u(iz)*cp_u) - vhb_u(id,iz)=vhb(id,iz) - uhb_u(id,iz)=uhb(id,iz) - ehb_u(id,iz)=ehb(id,iz) - qhb_u(id,iz)=-lfr(id,iz)/(da_u(iz)*latent) - sfr(id,iz)=sfr(id,iz)-pv_frac_roof*sfrpv(id,iz) - endif - else + sfr(id,iz)=sfr(id,iz)+pv_frac_roof*sfrpv(id,iz) + thb_u(id,iz)=-((1.-gr_frac_roof)*sfr(id,iz)+gr_frac_roof*sfrv(id,iz))/(da_u(iz)*cp_u) + vhb_u(id,iz)=(1.-gr_frac_roof)*vhb(id,iz)+gr_frac_roof*vhbv(id,iz) + uhb_u(id,iz)=(1.-gr_frac_roof)*uhb(id,iz)+gr_frac_roof*uhbv(id,iz) + ehb_u(id,iz)=(1.-gr_frac_roof)*ehb(id,iz)+gr_frac_roof*ehbv(id,iz) + qhb_u(id,iz)=-(gr_frac_roof*lfrv(id,iz)+(1.-gr_frac_roof)*lfr(id,iz))/(da_u(iz)*latent) + sfr(id,iz)=sfr(id,iz)-pv_frac_roof*sfrpv(id,iz) + else + sfr(id,iz)=sfr(id,iz)+pv_frac_roof*sfrpv(id,iz) + thb_u(id,iz)=-sfr(id,iz)/(da_u(iz)*cp_u) + vhb_u(id,iz)=vhb(id,iz) + uhb_u(id,iz)=uhb(id,iz) + ehb_u(id,iz)=ehb(id,iz) + qhb_u(id,iz)=-lfr(id,iz)/(da_u(iz)*latent) + sfr(id,iz)=sfr(id,iz)-pv_frac_roof*sfrpv(id,iz) + endif + else uhb_u(id,iz) = 0.0 vhb_u(id,iz) = 0.0 thb_u(id,iz) = 0.0 ehb_u(id,iz) = 0.0 qhb_u(id,iz) = 0.0 - endif + endif enddo @@ -4594,7 +4686,10 @@ end subroutine flux_flat_roof ! ===6=8===============================================================72 subroutine flux_flat(dz,z0,ua,va,pt,pt0,ptg,qv, & - uhb,vhb,sf,lf,ehb,da,pr) + uhb,vhb,sf,lf,ehb,da,pr, & +!CORDEX + t2m) +!CORDEX ! ---------------------------------------------------------------------- ! Calculation of the flux at the ground @@ -4610,8 +4705,13 @@ subroutine flux_flat(dz,z0,ua,va,pt,pt0,ptg,qv, & real ua ! wind speed real va ! wind speed real z0 ! Roughness length + real z0t ! Roughness length for heat real da ! air density real qv +!CORDEX + real t2m + real a2 +!CORDEX ! ---------------------------------------------------------------------- ! OUTPUT: ! ---------------------------------------------------------------------- @@ -4632,6 +4732,7 @@ subroutine flux_flat(dz,z0,ua,va,pt,pt0,ptg,qv, & ! LOCAL: ! ---------------------------------------------------------------------- real aa + real ah real al real buu real c @@ -4648,11 +4749,18 @@ subroutine flux_flat(dz,z0,ua,va,pt,pt0,ptg,qv, & real zz real qvsg,qvs,es,esa,fbqq,tmp,resg real b,cm,ch,rr,tol + real akanda,sqvisc,ust_g + integer ikanda parameter(b=9.4,cm=7.4,ch=5.3,rr=0.74,tol=.001) ! ---------------------------------------------------------------------- ! END VARIABLES DEFINITIONS ! ---------------------------------------------------------------------- + akanda=1.29 + sqvisc = 258.2 + + ! ikanda=0 !ikanda=0 => same roughness length for momentum and heat - original BEP-BEM formulation + ikanda=1 !ikanda=1 => Kanda et al. (2007) fromula to estimate roughness length for heat - same as for SLUCM ! computation of the ground temperature @@ -4673,6 +4781,14 @@ subroutine flux_flat(dz,z0,ua,va,pt,pt0,ptg,qv, & aa=vk/log(zz/z0) + if(ikanda.eq.0)then + z0t=z0 + elseif(ikanda.eq.1)then + ust_g=utot*vk/alog(zz/z0) + z0t = exp (2.0-akanda*(sqvisc**2 * ust_g * z0)**0.25)* z0 + endif + + ah=vk/log(zz/z0t) tmp=ptg*(pr/(1.e+5))**(rcp_u)-273.15 @@ -4689,13 +4805,14 @@ subroutine flux_flat(dz,z0,ua,va,pt,pt0,ptg,qv, & else c=b*cm*aa*aa*(zz/z0)**.5 fm=1.-b*ric/(1.+c*(-ric)**.5) + c=b*cm*aa*ah*(zz/z0t)**.5 c=c*ch/cm - fh=1.-b*ric/(1.+c*(-ric)**.5) + fh=1-b*ric/(1+c*(-ric)**.5) endif resg= rr/(aa*aa*utot*fh) fbuw=-aa*aa*utot*utot*fm - fbpt=-aa*aa*utot*(pt-ptg)*fh/rr + fbpt=-aa*ah*utot*(pt-ptg)*fh/rr fbqq=-(qv-qvsg)/(resg) ustar=(-fbuw)**.5 @@ -4711,6 +4828,12 @@ subroutine flux_flat(dz,z0,ua,va,pt,pt0,ptg,qv, & lf= ustar*qstar*da*latent ehb=buu !!!!!!!!!!!!!!! +!CORDEX + a2=vk/log(2./z0t) + t2m=ah/a2*pt+(1.-ah/a2)*ptg + t2m=t2m*(pr/1.e5)**rcp_u +!CORDEX + return end subroutine flux_flat diff --git a/phys/module_sf_noahdrv.F b/phys/module_sf_noahdrv.F index 5c7df673a2..632305f0c2 100644 --- a/phys/module_sf_noahdrv.F +++ b/phys/module_sf_noahdrv.F @@ -68,6 +68,11 @@ SUBROUTINE lsm(DZ8W,QV3D,P8W3D,T3D,TSK, & !Optional Urban TR_URB2D,TB_URB2D,TG_URB2D,TC_URB2D,QC_URB2D, & !H urban UC_URB2D, & !H urban +!CORDEX + T2G_URB2D,T2R_URB2D,T2VEG_URB2D, & + TCAN_URB2D,UCAN_URB2D,QCAN_URB2D, & + TGROUND_URB2D,TROOF_URB2D, & +!CORDEX XXXR_URB2D,XXXB_URB2D,XXXG_URB2D,XXXC_URB2D, & !H urban TRL_URB3D,TBL_URB3D,TGL_URB3D, & !H urban SH_URB2D,LH_URB2D,G_URB2D,RN_URB2D,TS_URB2D, & !H urban @@ -535,6 +540,16 @@ SUBROUTINE lsm(DZ8W,QV3D,P8W3D,T3D,TSK, & REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: TB_URB2D REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: TG_URB2D REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: TC_URB2D +!CORDEX + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: T2G_URB2D + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: T2R_URB2D + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: T2VEG_URB2D + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: TCAN_URB2D + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: UCAN_URB2D + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: QCAN_URB2D + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: TGROUND_URB2D + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: TROOF_URB2D +!CORDEX REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: QC_URB2D REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: UC_URB2D REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: XXXR_URB2D @@ -686,6 +701,9 @@ SUBROUTINE lsm(DZ8W,QV3D,P8W3D,T3D,TSK, & REAL :: CMR_URB, CHR_URB, CMC_URB, CHC_URB, CMGR_URB, CHGR_URB REAL :: frc_urb,lb_urb REAL :: check +!CORDEX + real aa_veg,a2_veg +!CORDEX ! ---------------------------------------------------------------------- ! DECLARATIONS END - urban ! ---------------------------------------------------------------------- @@ -1658,6 +1676,10 @@ SUBROUTINE lsm(DZ8W,QV3D,P8W3D,T3D,TSK, & b_t_bep,b_e_bep,b_q_bep,dlg_bep, & dl_u_bep,sf_bep,vl_bep, & rl_up_urb,rs_abs_urb,emiss_urb,grdflx_urb,qv3d, & +!CORDEX + tground_urb2d,troof_urb2d,tcan_urb2d, & + t2g_urb2d,t2r_urb2d,ucan_urb2d,qcan_urb2d, & +!CORDEX ids,ide, jds,jde, kds,kde, & ims,ime, jms,jme, kms,kme, & its,ite, jts,jte, kts,kte ) @@ -1752,6 +1774,20 @@ SUBROUTINE lsm(DZ8W,QV3D,P8W3D,T3D,TSK, & ! print*,'ust',ust(i,j) ! print*,'swdown,glw',swdown(i,j),glw(i,j) ! endif +!CORDEX here I do a simple log interpolation between the surface temperature of the vegetated frction, tsk_rural, and the lowest +!model level value of temperature + aa_veg=alog(dz8w(i,1,j)/2./znt(i,j)) + a2_veg=alog(2/znt(i,j)) + t2veg_urb2d(i,j)=(1.-aa_veg/a2_veg)*tsk_rural(i,j)+aa_veg/a2_veg*th_phy(i,1,j) +! write(*,*)'znt=',znt(i,j) +! write(*,*)'aa_veg=',aa_veg +! write(*,*)'a2_veg=',a2_veg +! write(*,*)'tsk_rural=',tsk_rural(i,j) +! write(*,*)'th_phy=',th_phy(i,1,j) +! write(*,*)'p_phy=',p_phy(i,1,j) + + t2veg_urb2d(i,j)=t2veg_urb2d(i,j)*(p_phy(i,1,j)/1.e5)**0.286 +!CORDEX else SH_URB2D(I,J) = 0. LH_URB2D(I,J) = 0. diff --git a/phys/module_surface_driver.F b/phys/module_surface_driver.F index f1592a1f00..899acc8b6b 100644 --- a/phys/module_surface_driver.F +++ b/phys/module_surface_driver.F @@ -187,6 +187,11 @@ SUBROUTINE surface_driver( & & ,num_road_layers, dzr, dzb, dzg & !I urban & ,tr_urb2d,tb_urb2d,tg_urb2d,tc_urb2d,qc_urb2d & !H urban & ,uc_urb2d & !H urban +!CORDEX + & ,t2g_urb2d,t2r_urb2d,t2veg_urb2d & + & ,tcan_urb2d,ucan_urb2d,qcan_urb2d & + & ,tground_urb2d,troof_urb2d & +!CORDEX & ,xxxr_urb2d,xxxb_urb2d,xxxg_urb2d,xxxc_urb2d & !H urban & ,cmcr_urb2d,tgr_urb2d,tgrl_urb3d,smr_urb3d & !H urban & ,julian,julyr,drelr_urb2d,drelb_urb2d,drelg_urb2d & !H urban @@ -360,7 +365,12 @@ SUBROUTINE surface_driver( & USE module_sf_gfs USE module_sf_noahdrv ! danli mosaic, the " ,only : lsm " needs to be deleted USE module_sf_noahlsm, only : LCZ_1,LCZ_2,LCZ_3,LCZ_4,LCZ_5,LCZ_6,LCZ_7,LCZ_8,LCZ_9,LCZ_10,LCZ_11 - USE module_sf_noahmpdrv, only : noahmplsm, noahmp_urban + USE module_sf_noahmpdrv, only : noahmplsm,noahmp_urban +!CORDEX- BUG + USE NOAHMP_TABLES, ONLY: ISURBAN_TABLE,LCZ_1_TABLE,LCZ_2_TABLE,LCZ_3_TABLE,LCZ_4_TABLE, & + LCZ_5_TABLE,LCZ_6_TABLE,LCZ_7_TABLE,LCZ_8_TABLE, & + LCZ_9_TABLE,LCZ_10_TABLE,LCZ_11_TABLE +!CORDEX- BUG USE module_sf_noahmp_groundwater USE module_sf_noah_seaice_drv #ifdef WRF_USE_CLM @@ -1273,6 +1283,16 @@ SUBROUTINE surface_driver( & REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: TC_URB2D !urban REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: QC_URB2D !urban REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: UC_URB2D !urban +!CORDEX + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: T2G_URB2D !urban + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: T2R_URB2D !urban + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: T2VEG_URB2D !urban + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: TCAN_URB2D !urban + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: UCAN_URB2D !urban + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: QCAN_URB2D !urban + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: TGROUND_URB2D !urban + REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: TROOF_URB2D !urban +!CORDEX REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: XXXR_URB2D !urban REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: XXXB_URB2D !urban REAL, OPTIONAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT):: XXXG_URB2D !urban @@ -2838,6 +2858,11 @@ SUBROUTINE surface_driver( & ,cmgr_sfcdif,chgr_sfcdif & ,tr_urb2d,tb_urb2d,tg_urb2d,tc_urb2d,qc_urb2d, & !H urban uc_urb2d, & !H urban +!CORDEX + t2g_urb2d,t2r_urb2d,t2veg_urb2d, & + tcan_urb2d,ucan_urb2d,qcan_urb2d, & + tground_urb2d,troof_urb2d, & +!CORDEX xxxr_urb2d,xxxb_urb2d,xxxg_urb2d,xxxc_urb2d, & !H urban trl_urb3d,tbl_urb3d,tgl_urb3d, & !H urban sh_urb2d,lh_urb2d,g_urb2d,rn_urb2d,ts_urb2d, & !H urban @@ -3008,15 +3033,28 @@ SUBROUTINE surface_driver( & IF((SF_URBAN_PHYSICS.eq.2).OR.(SF_URBAN_PHYSICS.eq.3)) THEN DO j=j_start(ij),j_end(ij) !urban DO i=i_start(ij),i_end(ij) !urban +! write(*,*)'ivgtyp',IVGTYP(I,J),LCZ_2,LCZ_6,LCZ_8 IF( IVGTYP(I,J) == ISURBAN .or. IVGTYP(I,J) == LCZ_1 .or. IVGTYP(I,J) == LCZ_2 .or. & IVGTYP(I,J) == LCZ_3 .or. IVGTYP(I,J) == LCZ_4 .or. IVGTYP(I,J) == LCZ_5 .or. & IVGTYP(I,J) == LCZ_6 .or. IVGTYP(I,J) == LCZ_7 .or. IVGTYP(I,J) == LCZ_8 .or. & IVGTYP(I,J) == LCZ_9 .or. IVGTYP(I,J) == LCZ_10 .or. IVGTYP(I,J) == LCZ_11 )THEN - T2(I,J) = TH_PHY(i,1,j)/((1.E5/PSFC(I,J))**RCP) !urban - TH2(I,J) = TH_PHY(i,1,j) !urban - Q2(I,J) = qv_curr(i,1,j) !urban - U10(I,J) = U_phy(I,1,J) !urban - V10(I,J) = V_phy(I,1,J) !urban + T2(I,J) = TH_PHY(i,1,j)/((1.E5/PSFC(I,J))**RCP) !urban + TH2(I,J) = TH_PHY(i,1,j) !urban + Q2(I,J) = qv_curr(i,1,j) !urban + U10(I,J) = U_phy(I,1,J) !urban + V10(I,J) = V_phy(I,1,J) !urban +!CORDEX ! for non urban points put 1.e20 for urban diagnostic variables + else + tsk_rural(i,j)=1.e20 + tcan_urb2d(i,j)=1.e20 + ucan_urb2d(i,j)=1.e20 + qcan_urb2d(i,j)=1.e20 + t2g_urb2d(i,j)=1.e20 + t2r_urb2d(i,j)=1.e20 + t2veg_urb2d(i,j)=1.e20 + troof_urb2d(i,j)=1.e20 + tground_urb2d(i,j)=1.e20 +!CORDEX END IF !urban ENDDO !urban ENDDO !urban @@ -3170,6 +3208,16 @@ SUBROUTINE surface_driver( & IF(SF_URBAN_PHYSICS > 0 ) THEN !urban + DO j=j_start(ij),j_end(ij) + DO i=i_start(ij),i_end(ij) + IF( IVGTYP(I,J) == ISURBAN .or. IVGTYP(I,J) == LCZ_1_TABLE .or. IVGTYP(I,J) == LCZ_2_TABLE .or. & + IVGTYP(I,J) == LCZ_3_TABLE .or. IVGTYP(I,J) == LCZ_4_TABLE .or. IVGTYP(I,J) == LCZ_5_TABLE .or. & + IVGTYP(I,J) == LCZ_6_TABLE .or. IVGTYP(I,J) == LCZ_7_TABLE .or. IVGTYP(I,J) == LCZ_8_TABLE .or. & + IVGTYP(I,J) == LCZ_9_TABLE .or. IVGTYP(I,J) == LCZ_10_TABLE .or. IVGTYP(I,J) == LCZ_11_TABLE )THEN + tsk_rural(i,j)=tsk(i,j) + ENDIF + ENDDO + ENDDO call noahmp_urban (sf_urban_physics, NUM_SOIL_LAYERS, IVGTYP,ITIMESTEP, & ! IN : Model configuration DTBL, COSZEN, XLAT_URB2D, & ! IN : Time/Space-related T_PHY, QV_CURR, U_PHY, V_PHY, SWDOWN, & ! IN : Forcing @@ -3184,6 +3232,12 @@ SUBROUTINE surface_driver( & chc_sfcdif, cmgr_sfcdif, chgr_sfcdif, & tr_urb2d, tb_urb2d, tg_urb2d, & !H urban tc_urb2d, qc_urb2d, uc_urb2d, & !H urban + tsk_rural, & +!CORDEX + t2g_urb2d, t2r_urb2d, & + tcan_urb2d, ucan_urb2d, qcan_urb2d, & + tground_urb2d, troof_urb2d, & +!CORDEX xxxr_urb2d, xxxb_urb2d, xxxg_urb2d, xxxc_urb2d, & !H urban trl_urb3d, tbl_urb3d, tgl_urb3d, & !H urban sh_urb2d, lh_urb2d, g_urb2d, rn_urb2d, ts_urb2d, & !H urban @@ -3223,7 +3277,6 @@ SUBROUTINE surface_driver( & a_e_bep, b_u_bep, b_v_bep, & !O multi-layer urban b_t_bep, b_q_bep, b_e_bep, dlg_bep, & !O multi-layer urban dl_u_bep, sf_bep, vl_bep) !O multi-layer urban - ENDIF IF ( iopt_run .EQ. 5 ) THEN @@ -3341,10 +3394,10 @@ SUBROUTINE surface_driver( & ENDIF TH2(I,J) = T2(I,J)*(1.E5/PSFC(I,J))**RCP ! ELSEIF (IVGTYP(I,J) == ISURBAN .OR. IVGTYP(I,J) == ISICE .OR. FVEGXY(I,J) == 0.0 ) THEN - ELSEIF( IVGTYP(I,J) == ISURBAN .or. IVGTYP(I,J) == LCZ_1 .or. IVGTYP(I,J) == LCZ_2 .or. & - IVGTYP(I,J) == LCZ_3 .or. IVGTYP(I,J) == LCZ_4 .or. IVGTYP(I,J) == LCZ_5 .or. & - IVGTYP(I,J) == LCZ_6 .or. IVGTYP(I,J) == LCZ_7 .or. IVGTYP(I,J) == LCZ_8 .or. & - IVGTYP(I,J) == LCZ_9 .or. IVGTYP(I,J) == LCZ_10 .or. IVGTYP(I,J) == LCZ_11 .or. & + ELSEIF( IVGTYP(I,J) == ISURBAN .or. IVGTYP(I,J) == LCZ_1_TABLE .or. IVGTYP(I,J) == LCZ_2_TABLE .or. & + IVGTYP(I,J) == LCZ_3_TABLE .or. IVGTYP(I,J) == LCZ_4_TABLE .or. IVGTYP(I,J) == LCZ_5_TABLE .or. & + IVGTYP(I,J) == LCZ_6_TABLE .or. IVGTYP(I,J) == LCZ_7_TABLE .or. IVGTYP(I,J) == LCZ_8_TABLE .or. & + IVGTYP(I,J) == LCZ_9_TABLE .or. IVGTYP(I,J) == LCZ_10_TABLE .or. IVGTYP(I,J) == LCZ_11_TABLE .or. & (IVGTYP(I,J) == ISICE .AND. XICE(I,J) .LT. XICE_THRESHOLD)) THEN Q2(I,J) = Q2MBXY(I,J) @@ -3369,10 +3422,10 @@ SUBROUTINE surface_driver( & IF(SF_URBAN_PHYSICS.eq.1) THEN DO j=j_start(ij),j_end(ij) !urban DO i=i_start(ij),i_end(ij) !urban - IF( IVGTYP(I,J) == ISURBAN .or. IVGTYP(I,J) == LCZ_1 .or. IVGTYP(I,J) == LCZ_2 .or. & - IVGTYP(I,J) == LCZ_3 .or. IVGTYP(I,J) == LCZ_4 .or. IVGTYP(I,J) == LCZ_5 .or. & - IVGTYP(I,J) == LCZ_6 .or. IVGTYP(I,J) == LCZ_7 .or. IVGTYP(I,J) == LCZ_8 .or. & - IVGTYP(I,J) == LCZ_9 .or. IVGTYP(I,J) == LCZ_10 .or. IVGTYP(I,J) == LCZ_11 )THEN + IF( IVGTYP(I,J) == ISURBAN .or. IVGTYP(I,J) == LCZ_1_TABLE .or. IVGTYP(I,J) == LCZ_2_TABLE .or. & + IVGTYP(I,J) == LCZ_3_TABLE .or. IVGTYP(I,J) == LCZ_4_TABLE .or. IVGTYP(I,J) == LCZ_5_TABLE .or. & + IVGTYP(I,J) == LCZ_6_TABLE .or. IVGTYP(I,J) == LCZ_7_TABLE .or. IVGTYP(I,J) == LCZ_8_TABLE .or. & + IVGTYP(I,J) == LCZ_9_TABLE .or. IVGTYP(I,J) == LCZ_10_TABLE .or. IVGTYP(I,J) == LCZ_11_TABLE )THEN Q2(I,J) = (FVEGXY(I,J)*Q2MVXY(I,J) + (1.-FVEGXY(I,J))*Q2MBXY(I,J))*(1.-FRC_URB2D(I,J)) + & Q2_URB2D(I,J)*FRC_URB2D(I,J) T2(I,J) = (FVEGXY(I,J)*T2MVXY(I,J) + (1.-FVEGXY(I,J))*T2MBXY(I,J))*(1.-FRC_URB2D(I,J)) + & @@ -3391,17 +3444,35 @@ SUBROUTINE surface_driver( & ENDIF IF((SF_URBAN_PHYSICS.eq.2).OR.(SF_URBAN_PHYSICS.eq.3)) THEN + ! write(*,*)'entering section',j_start(ij),j_end(ij),i_start(ij),i_end(ij) DO j=j_start(ij),j_end(ij) !urban DO i=i_start(ij),i_end(ij) !urban - IF( IVGTYP(I,J) == ISURBAN .or. IVGTYP(I,J) == LCZ_1 .or. IVGTYP(I,J) == LCZ_2 .or. & - IVGTYP(I,J) == LCZ_3 .or. IVGTYP(I,J) == LCZ_4 .or. IVGTYP(I,J) == LCZ_5 .or. & - IVGTYP(I,J) == LCZ_6 .or. IVGTYP(I,J) == LCZ_7 .or. IVGTYP(I,J) == LCZ_8 .or. & - IVGTYP(I,J) == LCZ_9 .or. IVGTYP(I,J) == LCZ_10 .or. IVGTYP(I,J) == LCZ_11 )THEN +!CORDEX - BUG + IF( IVGTYP(I,J) == ISURBAN .or. IVGTYP(I,J) == LCZ_1_TABLE .or. IVGTYP(I,J) == LCZ_2_TABLE .or. & + IVGTYP(I,J) == LCZ_3_TABLE .or. IVGTYP(I,J) == LCZ_4_TABLE .or. IVGTYP(I,J) == LCZ_5_TABLE .or. & + IVGTYP(I,J) == LCZ_6_TABLE .or. IVGTYP(I,J) == LCZ_7_TABLE .or. IVGTYP(I,J) == LCZ_8_TABLE .or. & + IVGTYP(I,J) == LCZ_9_TABLE .or. IVGTYP(I,J) == LCZ_10_TABLE .or. IVGTYP(I,J) == LCZ_11_TABLE )THEN +!CORDEX - BUG T2(I,J) = TH_PHY(i,1,j)/((1.E5/PSFC(I,J))**RCP) !urban TH2(I,J) = TH_PHY(i,1,j) !urban Q2(I,J) = qv_curr(i,1,j) !urban U10(I,J) = U_phy(I,1,J) !urban V10(I,J) = V_phy(I,1,J) !urban +!CORDEX + T2VEG_URB2D(I,J) = FVEGXY(I,J)*T2MVXY(I,J) + (1.-FVEGXY(I,J))*T2MBXY(I,J) +!CORDEX +!CORDEX ! for non urban points put 1.e20 for urban diagnostic variables + else + tsk_rural(i,j)=1.e20 + tcan_urb2d(i,j)=1.e20 + ucan_urb2d(i,j)=1.e20 + qcan_urb2d(i,j)=1.e20 + t2g_urb2d(i,j)=1.e20 + t2r_urb2d(i,j)=1.e20 + t2veg_urb2d(i,j)=1.e20 + troof_urb2d(i,j)=1.e20 + tground_urb2d(i,j)=1.e20 +!CORDEX END IF !urban ENDDO !urban ENDDO !urban diff --git a/phys/noahmp b/phys/noahmp index e348cac8f1..2e328ee6a1 160000 --- a/phys/noahmp +++ b/phys/noahmp @@ -1 +1 @@ -Subproject commit e348cac8f1949f88cd122dcf9c3a84bd3f3b0451 +Subproject commit 2e328ee6a1ed8a6c1fda8acba92553a8c2deddd8