module UrbanInitMod !----------------------------------------------------------------------- !BOP ! ! !MODULE: UrbanInitMod ! ! !DESCRIPTION: ! Initialize urban data ! ! !USES: use shr_kind_mod, only : r8 => shr_kind_r8 use abortutils , only : endrun use shr_sys_mod , only : shr_sys_flush use clm_varctl , only : iulog use UrbanMod, only : urban_traffic, urban_hac, urban_hac_off ! ! !PUBLIC TYPES: implicit none save private ! ! !PUBLIC MEMBER FUNCTIONS: public :: UrbanInitTimeVar ! Initialize urban time varying variables public :: UrbanInitTimeConst ! Initialize urban time constant variables public :: UrbanInitAero ! Calculate urban landunit aerodynamic constants ! !EOP !----------------------------------------------------------------------- contains !----------------------------------------------------------------------- !BOP ! ! !IROUTINE: UrbanInitAero ! ! !INTERFACE: subroutine UrbanInitAero( ) ! ! !DESCRIPTION: ! Calculate urban land unit aerodynamic constants using Macdonald (1998) as used in ! Grimmond and Oke (1999) ! ! !USES: use clmtype , only : clm3 use clm_varcon, only : isturb, vkc use decompMod , only : get_proc_bounds ! ! !ARGUMENTS: implicit none ! ! local pointers to original implicit in arguments (urban clump) ! real(r8), pointer :: ht_roof(:) ! height of urban roof (m) real(r8), pointer :: canyon_hwr(:) ! ratio of building height to street width (-) integer , pointer :: ltype(:) ! landunit type ! ! local pointers to original implicit out arguments ! real(r8), pointer :: z_0_town(:) ! urban landunit momentum roughness length (m) real(r8), pointer :: z_d_town(:) ! urban landunit displacement height (m) ! ! !CALLED FROM: ! subroutine initialize ! ! !REVISION HISTORY: ! Created by Keith Oleson January 2005 ! ! ! !LOCAL VARIABLES: !EOP real(r8), parameter :: alpha = 4.43_r8 ! coefficient used to calculate z_d_town real(r8), parameter :: beta = 1.0_r8 ! coefficient used to calculate z_d_town real(r8), parameter :: C_d = 1.2_r8 ! drag coefficient as used in Grimmond and Oke (1999) real(r8) :: plan_ai ! plan area index - ratio building area to plan area (-) real(r8) :: frontal_ai ! frontal area index of buildings (-) real(r8) :: build_lw_ratio ! building short/long side ratio (-) integer :: l,g ! indices integer :: begp, endp ! clump beginning and ending pft indices integer :: begc, endc ! clump beginning and ending column indices integer :: begl, endl ! clump beginning and ending landunit indices integer :: begg, endg ! clump beginning and ending gridcell indices !----------------------------------------------------------------------- ! Assign local pointers to derived type members (landunit level) ltype => clm3%g%l%itype z_0_town => clm3%g%l%z_0_town z_d_town => clm3%g%l%z_d_town ht_roof => clm3%g%l%ht_roof canyon_hwr => clm3%g%l%canyon_hwr call get_proc_bounds(begg, endg, begl, endl, begc, endc, begp, endp) do l = begl, endl if (ltype(l) == isturb) then ! Calculate plan area index plan_ai = canyon_hwr(l)/(canyon_hwr(l) + 1._r8) ! Building shape shortside/longside ratio (e.g. 1 = square ) ! This assumes the building occupies the entire canyon length build_lw_ratio = plan_ai ! Calculate frontal area index frontal_ai = (1._r8 - plan_ai) * canyon_hwr(l) ! Adjust frontal area index for different building configuration frontal_ai = frontal_ai * sqrt(1/build_lw_ratio) * sqrt(plan_ai) ! Calculate displacement height #if (defined VANCOUVER) z_d_town(l) = 3.5_r8 #elif (defined MEXICOCITY) z_d_town(l) = 10.9_r8 #else z_d_town(l) = (1._r8 + alpha**(-plan_ai) * (plan_ai - 1._r8)) * ht_roof(l) #endif ! Calculate the roughness length #if (defined VANCOUVER) z_0_town(l) = 0.35_r8 #elif (defined MEXICOCITY) z_0_town(l) = 2.2_r8 #else z_0_town(l) = ht_roof(l) * (1._r8 - z_d_town(l) / ht_roof(l)) * & exp(-1.0_r8 * (0.5_r8 * beta * C_d / vkc**2 * & (1 - z_d_town(l) / ht_roof(l)) * frontal_ai)**(-0.5_r8)) #endif end if end do end subroutine UrbanInitAero !----------------------------------------------------------------------- !BOP ! ! !IROUTINE: UrbanInitTimeConst ! ! !INTERFACE: subroutine UrbanInitTimeConst() ! ! !DESCRIPTION: ! Initialize urban time-constant variables ! ! !USES: use clmtype , only : clm3 use clm_varcon , only : isturb, icol_roof, icol_sunwall, icol_shadewall, & icol_road_perv, icol_road_imperv, spval use decompMod , only : get_proc_bounds, ldecomp use UrbanInputMod, only : urbinp ! ! !ARGUMENTS: implicit none ! ! !LOCAL VARIABLES: ! ! local pointers to original implicit in arguments ! integer , pointer :: gdc(:) ! grid index for landunit integer , pointer :: coli(:) ! beginning column index for landunit integer , pointer :: colf(:) ! ending column index for landunit integer , pointer :: ctype(:) ! column type integer , pointer :: ltype(:) ! landunit type index integer , pointer :: lgridcell(:) ! gridcell of corresponding landunit ! ! local pointers to original implicit out arguments ! real(r8), pointer :: canyon_hwr(:) ! urban canyon height to width ratio real(r8), pointer :: emg(:) ! ground emissivity real(r8), pointer :: wtroad_perv(:) ! weight of pervious column to total road real(r8), pointer :: ht_roof(:) ! height of urban roof (m) real(r8), pointer :: wtlunit_roof(:) ! weight of roof with respect to landunit real(r8), pointer :: wind_hgt_canyon(:) ! height above road at which wind in canyon is to be computed (m) real(r8), pointer :: eflx_traffic_factor(:) ! multiplicative factor for sensible heat flux from urban traffic real(r8), pointer :: t_building_max(:) ! maximum internal building temperature (K) real(r8), pointer :: t_building_min(:) ! minimum internal building temperature (K) real(r8), pointer :: tk_wall(:,:) ! thermal conductivity of urban wall (W/m/K) real(r8), pointer :: tk_roof(:,:) ! thermal conductivity of urban roof (W/m/K) real(r8), pointer :: tk_improad(:,:) ! thermal conductivity of urban impervious road (W/m/K) real(r8), pointer :: cv_wall(:,:) ! thermal conductivity of urban wall (J/m^3/K) real(r8), pointer :: cv_roof(:,:) ! thermal conductivity of urban roof (J/m^3/K) real(r8), pointer :: cv_improad(:,:) ! thermal conductivity of urban impervious road (J/m^3/K) real(r8), pointer :: thick_wall(:) ! thickness of urban wall (m) real(r8), pointer :: thick_roof(:) ! thickness of urban roof (m) integer, pointer :: nlev_improad(:) ! number of impervious road layers (-) ! ! ! !OTHER LOCAL VARIABLES !EOP integer :: nc,fl,ib,l,c,p,g ! indices integer :: ier ! error status integer :: begp, endp ! clump beginning and ending pft indices integer :: begc, endc ! clump beginning and ending column indices integer :: begl, endl ! clump beginning and ending landunit indices integer :: begg, endg ! clump beginning and ending gridcell indices ! Assign local pointers to derived type members (landunit-level) ltype => clm3%g%l%itype lgridcell => clm3%g%l%gridcell coli => clm3%g%l%coli colf => clm3%g%l%colf canyon_hwr => clm3%g%l%canyon_hwr wtroad_perv => clm3%g%l%wtroad_perv ht_roof => clm3%g%l%ht_roof wtlunit_roof => clm3%g%l%wtlunit_roof wind_hgt_canyon => clm3%g%l%wind_hgt_canyon eflx_traffic_factor => clm3%g%l%lef%eflx_traffic_factor t_building_max => clm3%g%l%lps%t_building_max t_building_min => clm3%g%l%lps%t_building_min canyon_hwr => clm3%g%l%canyon_hwr tk_wall => clm3%g%l%lps%tk_wall tk_roof => clm3%g%l%lps%tk_roof tk_improad => clm3%g%l%lps%tk_improad cv_wall => clm3%g%l%lps%cv_wall cv_roof => clm3%g%l%lps%cv_roof cv_improad => clm3%g%l%lps%cv_improad thick_wall => clm3%g%l%lps%thick_wall thick_roof => clm3%g%l%lps%thick_roof nlev_improad => clm3%g%l%lps%nlev_improad ! Assign local pointers to derived type members (column-level) ctype => clm3%g%l%c%itype emg => clm3%g%l%c%cps%emg ! Initialize time constant urban variables call get_proc_bounds(begg, endg, begl, endl, begc, endc, begp, endp) do l = begl, endl if (ltype(l) == isturb) then g = clm3%g%l%gridcell(l) canyon_hwr(l) = urbinp%canyon_hwr(g) wtroad_perv(l) = urbinp%wtroad_perv(g) ht_roof(l) = urbinp%ht_roof(g) wtlunit_roof(l) = urbinp%wtlunit_roof(g) wind_hgt_canyon(l) = urbinp%wind_hgt_canyon(g) tk_wall(l,:) = urbinp%tk_wall(g,:) tk_roof(l,:) = urbinp%tk_roof(g,:) tk_improad(l,:) = urbinp%tk_improad(g,:) cv_wall(l,:) = urbinp%cv_wall(g,:) cv_roof(l,:) = urbinp%cv_roof(g,:) cv_improad(l,:) = urbinp%cv_improad(g,:) thick_wall(l) = urbinp%thick_wall(g) thick_roof(l) = urbinp%thick_roof(g) nlev_improad(l) = urbinp%nlev_improad(g) t_building_min(l) = urbinp%t_building_min(g) t_building_max(l) = urbinp%t_building_max(g) do c = coli(l),colf(l) if (ctype(c) == icol_roof ) emg(c) = urbinp%em_roof(g) if (ctype(c) == icol_sunwall ) emg(c) = urbinp%em_wall(g) if (ctype(c) == icol_shadewall ) emg(c) = urbinp%em_wall(g) if (ctype(c) == icol_road_imperv) emg(c) = urbinp%em_improad(g) if (ctype(c) == icol_road_perv ) emg(c) = urbinp%em_perroad(g) end do ! Inferred from Sailor and Lu 2004 if (urban_traffic) then eflx_traffic_factor(l) = 3.6_r8 * (canyon_hwr(l)-0.5_r8) + 1.0_r8 else eflx_traffic_factor(l) = 0.0_r8 end if #if (defined VANCOUVER || defined MEXICOCITY) ! Freely evolving t_building_max(l) = 380.00_r8 t_building_min(l) = 200.00_r8 #else if (urban_hac == urban_hac_off) then ! Overwrite values read in from urbinp by freely evolving values t_building_max(l) = 380.00_r8 t_building_min(l) = 200.00_r8 end if #endif else eflx_traffic_factor(l) = spval t_building_max(l) = spval t_building_min(l) = spval end if end do end subroutine UrbanInitTimeConst !----------------------------------------------------------------------- !BOP ! ! !IROUTINE: UrbanInitTimeVar ! ! !INTERFACE: subroutine UrbanInitTimeVar( ) ! ! !DESCRIPTION: ! Initialize urban time-varying variables ! ! !USES: use clmtype , only : clm3 use clm_varcon, only : isturb, spval, icol_road_perv use decompMod , only : get_proc_bounds ! ! !ARGUMENTS: implicit none ! ! local pointers to original implicit in arguments (urban clump) ! integer , pointer :: ltype(:) ! landunit type integer , pointer :: lgridcell(:) ! gridcell of corresponding landunit integer , pointer :: clandunit(:) ! landunit index of corresponding column integer , pointer :: plandunit(:) ! landunit index of corresponding pft integer , pointer :: ctype(:) ! column type ! ! local pointers to original implicit out arguments ! real(r8), pointer :: taf(:) ! urban canopy air temperature (K) real(r8), pointer :: qaf(:) ! urban canopy air specific humidity (kg/kg) real(r8), pointer :: eflx_building_heat(:) ! heat flux from urban building interior to walls, roof (W/m**2) real(r8), pointer :: eflx_urban_ac(:) ! urban air conditioning flux (W/m**2) real(r8), pointer :: eflx_urban_heat(:) ! urban heating flux (W/m**2) real(r8), pointer :: fcov(:) ! fractional impermeable area real(r8), pointer :: fsat(:) ! fractional area with water table at surface real(r8), pointer :: qcharge(:) ! aquifer recharge rate (mm/s) real(r8), pointer :: t_building(:) ! internal building temperature (K) real(r8), pointer :: eflx_traffic(:) ! traffic sensible heat flux (W/m**2) real(r8), pointer :: eflx_wasteheat(:) ! sensible heat flux from urban heating/cooling sources of waste heat (W/m**2) real(r8), pointer :: eflx_wasteheat_pft(:) ! sensible heat flux from urban heating/cooling sources of waste heat at pft level (W/m**2) real(r8), pointer :: eflx_heat_from_ac_pft(:) ! sensible heat flux put back into canyon due to removal by AC (W/m**2) real(r8), pointer :: eflx_traffic_pft(:) ! sensible heat flux from traffic (W/m**2) real(r8), pointer :: eflx_anthro(:) ! total anthropogenic heat flux (W/m**2) real(r8), pointer :: t_ref2m_u(:) ! Urban 2 m height surface air temperature (Kelvin) real(r8), pointer :: t_ref2m_min_u(:) ! Urban daily minimum of average 2 m height surface air temperature (K) real(r8), pointer :: t_ref2m_max_u(:) ! Urban daily maximum of average 2 m height surface air temperature (K) real(r8), pointer :: rh_ref2m_u(:) ! Urban 2 m height surface relative humidity (%) real(r8), pointer :: t_grnd_u(:) ! Urban ground temperature (Kelvin) real(r8), pointer :: qflx_runoff_u(:) ! Urban total runoff (qflx_drain+qflx_surf) (mm H2O /s) real(r8), pointer :: fsa_u(:) ! Urban absorbed solar radiation (W/m**2) real(r8), pointer :: eflx_lwrad_net_u(:) ! Urban net longwave radiation (W/m**2) real(r8), pointer :: eflx_lh_tot_u(:) ! Urban latent heat flux (W/m**2) real(r8), pointer :: eflx_sh_tot_u(:) ! Urban sensible heat flux (W/m**2) real(r8), pointer :: eflx_soil_grnd_u(:) ! Urban ground heat flux (W/m**2) real(r8), pointer :: eflx_snomelt_u(:) ! Urban snow melt heat flux (W/m**2) ! ! !CALLED FROM: ! subroutine initialize ! ! !REVISION HISTORY: ! Created by Keith Oleson February 2005 ! ! ! !LOCAL VARIABLES: !EOP integer :: l,g,c,p ! indices integer :: begp, endp ! clump beginning and ending pft indices integer :: begc, endc ! clump beginning and ending column indices integer :: begl, endl ! clump beginning and ending landunit indices integer :: begg, endg ! clump beginning and ending gridcell indices !----------------------------------------------------------------------- ! Assign local pointers to derived type members (landunit level) taf => clm3%g%l%lps%taf qaf => clm3%g%l%lps%qaf ltype => clm3%g%l%itype lgridcell => clm3%g%l%gridcell t_building => clm3%g%l%lps%t_building eflx_traffic => clm3%g%l%lef%eflx_traffic eflx_wasteheat => clm3%g%l%lef%eflx_wasteheat ! Assign local pointers to derived type members (column level) clandunit => clm3%g%l%c%landunit eflx_building_heat => clm3%g%l%c%cef%eflx_building_heat eflx_urban_ac => clm3%g%l%c%cef%eflx_urban_ac eflx_urban_heat => clm3%g%l%c%cef%eflx_urban_heat fcov => clm3%g%l%c%cws%fcov fsat => clm3%g%l%c%cws%fsat qcharge => clm3%g%l%c%cws%qcharge ctype => clm3%g%l%c%itype t_grnd_u => clm3%g%l%c%ces%t_grnd_u qflx_runoff_u => clm3%g%l%c%cwf%qflx_runoff_u eflx_snomelt_u => clm3%g%l%c%cef%eflx_snomelt_u ! Assign local pointers to derived type members (pft level) t_ref2m_u => clm3%g%l%c%p%pes%t_ref2m_u t_ref2m_min_u => clm3%g%l%c%p%pes%t_ref2m_min_u t_ref2m_max_u => clm3%g%l%c%p%pes%t_ref2m_max_u rh_ref2m_u => clm3%g%l%c%p%pes%rh_ref2m_u plandunit => clm3%g%l%c%p%landunit eflx_wasteheat_pft => clm3%g%l%c%p%pef%eflx_wasteheat_pft eflx_heat_from_ac_pft => clm3%g%l%c%p%pef%eflx_heat_from_ac_pft eflx_traffic_pft => clm3%g%l%c%p%pef%eflx_traffic_pft eflx_anthro => clm3%g%l%c%p%pef%eflx_anthro fsa_u => clm3%g%l%c%p%pef%fsa_u eflx_lwrad_net_u => clm3%g%l%c%p%pef%eflx_lwrad_net_u eflx_lh_tot_u => clm3%g%l%c%p%pef%eflx_lh_tot_u eflx_sh_tot_u => clm3%g%l%c%p%pef%eflx_sh_tot_u eflx_soil_grnd_u => clm3%g%l%c%p%pef%eflx_soil_grnd_u call get_proc_bounds(begg, endg, begl, endl, begc, endc, begp, endp) do l = begl, endl g = lgridcell(l) if (ltype(l) == isturb) then #if (defined VANCOUVER) taf(l) = 297.56_r8 qaf(l) = 0.0111_r8 #elif (defined MEXICOCITY) taf(l) = 289.46_r8 qaf(l) = 0.00248_r8 #else taf(l) = 283._r8 ! Arbitrary set since forc_q is not yet available qaf(l) = 1.e-4_r8 #endif else t_building(l) = spval eflx_traffic(l) = spval eflx_wasteheat(l) = spval end if end do do c = begc, endc l = clandunit(c) if (ltype(l) == isturb) then eflx_building_heat(c) = 0._r8 eflx_urban_ac(c) = 0._r8 eflx_urban_heat(c) = 0._r8 ! ! Set hydrology variables for urban to spvalue -- as only valid for pervious road ! if (ctype(c) /= icol_road_perv )then fcov(c) = spval fsat(c) = spval qcharge(c) = spval end if else eflx_building_heat(c) = spval eflx_urban_ac(c) = spval eflx_urban_heat(c) = spval t_grnd_u(c) = spval qflx_runoff_u(c) = spval eflx_snomelt_u(c) = spval end if end do do p = begp, endp l = plandunit(p) if (ltype(l) /= isturb) then t_ref2m_u(p) = spval t_ref2m_min_u(p) = spval t_ref2m_max_u(p) = spval rh_ref2m_u(p) = spval eflx_wasteheat_pft(p) = spval eflx_heat_from_ac_pft(p) = spval eflx_traffic_pft(p) = spval eflx_anthro(p) = spval fsa_u(p) = spval eflx_lwrad_net_u(p) = spval eflx_lh_tot_u(p) = spval eflx_sh_tot_u(p) = spval eflx_soil_grnd_u(p) = spval end if end do end subroutine UrbanInitTimeVar end module UrbanInitMod