module UrbanParamsType !------------------------------------------------------------------------------ ! !DESCRIPTION: ! Urban Constants ! ! !USES: use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use abortutils , only : endrun use decompMod , only : bounds_type use clm_varctl , only : iulog, fsurdat use clm_varcon , only : namel, grlnd, spval use LandunitType , only : lun ! implicit none save private ! ! !PUBLIC MEMBER FUNCTIONS: public :: UrbanReadNML ! Read in the urban namelist items public :: UrbanInput ! Read in urban input data public :: CheckUrban ! Check validity of urban points public :: IsSimpleBuildTemp ! If using the simple building temperature method public :: IsProgBuildTemp ! If using the prognostic building temperature method ! ! !PRIVATE TYPE type urbinp_type real(r8), pointer :: canyon_hwr (:,:) real(r8), pointer :: wtlunit_roof (:,:) real(r8), pointer :: wtroad_perv (:,:) real(r8), pointer :: em_roof (:,:) real(r8), pointer :: em_improad (:,:) real(r8), pointer :: em_perroad (:,:) real(r8), pointer :: em_wall (:,:) real(r8), pointer :: alb_roof_dir (:,:,:) real(r8), pointer :: alb_roof_dif (:,:,:) real(r8), pointer :: alb_improad_dir (:,:,:) real(r8), pointer :: alb_improad_dif (:,:,:) real(r8), pointer :: alb_perroad_dir (:,:,:) real(r8), pointer :: alb_perroad_dif (:,:,:) real(r8), pointer :: alb_wall_dir (:,:,:) real(r8), pointer :: alb_wall_dif (:,:,:) real(r8), pointer :: ht_roof (:,:) real(r8), pointer :: wind_hgt_canyon (:,:) real(r8), pointer :: tk_wall (:,:,:) real(r8), pointer :: tk_roof (:,:,:) real(r8), pointer :: tk_improad (:,:,:) real(r8), pointer :: cv_wall (:,:,:) real(r8), pointer :: cv_roof (:,:,:) real(r8), pointer :: cv_improad (:,:,:) real(r8), pointer :: thick_wall (:,:) real(r8), pointer :: thick_roof (:,:) integer, pointer :: nlev_improad (:,:) real(r8), pointer :: t_building_min (:,:) end type urbinp_type type (urbinp_type), public :: urbinp ! urban input derived type ! !PUBLIC TYPE type, public :: urbanparams_type real(r8), allocatable :: wind_hgt_canyon (:) ! lun height above road at which wind in canyon is to be computed (m) real(r8), allocatable :: em_roof (:) ! lun roof emissivity real(r8), allocatable :: em_improad (:) ! lun impervious road emissivity real(r8), allocatable :: em_perroad (:) ! lun pervious road emissivity real(r8), allocatable :: em_wall (:) ! lun wall emissivity real(r8), allocatable :: alb_roof_dir (:,:) ! lun direct roof albedo real(r8), allocatable :: alb_roof_dif (:,:) ! lun diffuse roof albedo real(r8), allocatable :: alb_improad_dir (:,:) ! lun direct impervious road albedo real(r8), allocatable :: alb_improad_dif (:,:) ! lun diffuse impervious road albedo real(r8), allocatable :: alb_perroad_dir (:,:) ! lun direct pervious road albedo real(r8), allocatable :: alb_perroad_dif (:,:) ! lun diffuse pervious road albedo real(r8), allocatable :: alb_wall_dir (:,:) ! lun direct wall albedo real(r8), allocatable :: alb_wall_dif (:,:) ! lun diffuse wall albedo integer , pointer :: nlev_improad (:) ! lun number of impervious road layers (-) real(r8), pointer :: tk_wall (:,:) ! lun thermal conductivity of urban wall (W/m/K) real(r8), pointer :: tk_roof (:,:) ! lun thermal conductivity of urban roof (W/m/K) real(r8), pointer :: tk_improad (:,:) ! lun thermal conductivity of urban impervious road (W/m/K) real(r8), pointer :: cv_wall (:,:) ! lun heat capacity of urban wall (J/m^3/K) real(r8), pointer :: cv_roof (:,:) ! lun heat capacity of urban roof (J/m^3/K) real(r8), pointer :: cv_improad (:,:) ! lun heat capacity of urban impervious road (J/m^3/K) real(r8), pointer :: thick_wall (:) ! lun total thickness of urban wall (m) real(r8), pointer :: thick_roof (:) ! lun total thickness of urban roof (m) real(r8), pointer :: vf_sr (:) ! lun view factor of sky for road real(r8), pointer :: vf_wr (:) ! lun view factor of one wall for road real(r8), pointer :: vf_sw (:) ! lun view factor of sky for one wall real(r8), pointer :: vf_rw (:) ! lun view factor of road for one wall real(r8), pointer :: vf_ww (:) ! lun view factor of opposing wall for one wall real(r8), pointer :: t_building_min (:) ! lun minimum internal building air temperature (K) real(r8), pointer :: eflx_traffic_factor (:) ! lun multiplicative traffic factor for sensible heat flux from urban traffic (-) contains procedure, public :: Init end type urbanparams_type ! ! !Urban control variables character(len= *), parameter, public :: urban_hac_off = 'OFF' character(len= *), parameter, public :: urban_hac_on = 'ON' character(len= *), parameter, public :: urban_wasteheat_on = 'ON_WASTEHEAT' character(len= 16), public :: urban_hac = urban_hac_off logical, public :: urban_traffic = .false. ! urban traffic fluxes ! !PRIVATE MEMBER DATA: logical, private :: ReadNamelist = .false. ! If namelist was read yet or not integer, parameter, private :: BUILDING_TEMP_METHOD_SIMPLE = 0 ! Simple method introduced in CLM4.5 integer, parameter, private :: BUILDING_TEMP_METHOD_PROG = 1 ! Prognostic method introduced in CLM5.0 integer, private :: building_temp_method = BUILDING_TEMP_METHOD_PROG ! Method to calculate the building temperature character(len=*), parameter, private :: sourcefile = & __FILE__ !----------------------------------------------------------------------- contains !----------------------------------------------------------------------- subroutine Init(this, bounds) ! ! Allocate module variables and data structures ! ! !USES: use shr_infnan_mod , only : nan => shr_infnan_nan, assignment(=) use clm_varpar , only : nlevcan, nlevcan, numrad, nlevgrnd, nlevurb use clm_varpar , only : nlevsoi, nlevgrnd use clm_varctl , only : use_vancouver, use_mexicocity use clm_varcon , only : vkc use column_varcon , only : icol_roof, icol_sunwall, icol_shadewall use column_varcon , only : icol_road_perv, icol_road_imperv, icol_road_perv use landunit_varcon , only : isturb_MIN ! ! !ARGUMENTS: class(urbanparams_type) :: this type(bounds_type) , intent(in) :: bounds ! ! !LOCAL VARIABLES: integer :: j,l,c,p,g ! indices integer :: nc,fl,ib ! indices integer :: dindx ! urban density type index integer :: ier ! error status real(r8) :: sumvf ! sum of view factors for wall or road 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 :: begl, endl integer :: begc, endc integer :: begp, endp integer :: begg, endg !--------------------------------------------------------------------- begp = bounds%begp; endp = bounds%endp begc = bounds%begc; endc = bounds%endc begl = bounds%begl; endl = bounds%endl begg = bounds%begg; endg = bounds%endg ! Allocate urbanparams data structure if ( nlevurb > 0 )then allocate(this%tk_wall (begl:endl,nlevurb)) ; this%tk_wall (:,:) = nan allocate(this%tk_roof (begl:endl,nlevurb)) ; this%tk_roof (:,:) = nan allocate(this%cv_wall (begl:endl,nlevurb)) ; this%cv_wall (:,:) = nan allocate(this%cv_roof (begl:endl,nlevurb)) ; this%cv_roof (:,:) = nan end if allocate(this%t_building_min (begl:endl)) ; this%t_building_min (:) = nan allocate(this%tk_improad (begl:endl,nlevurb)) ; this%tk_improad (:,:) = nan allocate(this%cv_improad (begl:endl,nlevurb)) ; this%cv_improad (:,:) = nan allocate(this%thick_wall (begl:endl)) ; this%thick_wall (:) = nan allocate(this%thick_roof (begl:endl)) ; this%thick_roof (:) = nan allocate(this%nlev_improad (begl:endl)) ; this%nlev_improad (:) = huge(1) allocate(this%vf_sr (begl:endl)) ; this%vf_sr (:) = nan allocate(this%vf_wr (begl:endl)) ; this%vf_wr (:) = nan allocate(this%vf_sw (begl:endl)) ; this%vf_sw (:) = nan allocate(this%vf_rw (begl:endl)) ; this%vf_rw (:) = nan allocate(this%vf_ww (begl:endl)) ; this%vf_ww (:) = nan allocate(this%wind_hgt_canyon (begl:endl)) ; this%wind_hgt_canyon (:) = nan allocate(this%em_roof (begl:endl)) ; this%em_roof (:) = nan allocate(this%em_improad (begl:endl)) ; this%em_improad (:) = nan allocate(this%em_perroad (begl:endl)) ; this%em_perroad (:) = nan allocate(this%em_wall (begl:endl)) ; this%em_wall (:) = nan allocate(this%alb_roof_dir (begl:endl,numrad)) ; this%alb_roof_dir (:,:) = nan allocate(this%alb_roof_dif (begl:endl,numrad)) ; this%alb_roof_dif (:,:) = nan allocate(this%alb_improad_dir (begl:endl,numrad)) ; this%alb_improad_dir (:,:) = nan allocate(this%alb_perroad_dir (begl:endl,numrad)) ; this%alb_perroad_dir (:,:) = nan allocate(this%alb_improad_dif (begl:endl,numrad)) ; this%alb_improad_dif (:,:) = nan allocate(this%alb_perroad_dif (begl:endl,numrad)) ; this%alb_perroad_dif (:,:) = nan allocate(this%alb_wall_dir (begl:endl,numrad)) ; this%alb_wall_dir (:,:) = nan allocate(this%alb_wall_dif (begl:endl,numrad)) ; this%alb_wall_dif (:,:) = nan allocate(this%eflx_traffic_factor (begl:endl)) ; this%eflx_traffic_factor (:) = nan ! Initialize time constant urban variables do l = bounds%begl,bounds%endl ! "0" refers to urban wall/roof surface and "nlevsoi" refers to urban wall/roof bottom if (lun%urbpoi(l)) then g = lun%gridcell(l) dindx = lun%itype(l) - isturb_MIN + 1 this%wind_hgt_canyon(l) = urbinp%wind_hgt_canyon(g,dindx) do ib = 1,numrad this%alb_roof_dir (l,ib) = urbinp%alb_roof_dir (g,dindx,ib) this%alb_roof_dif (l,ib) = urbinp%alb_roof_dif (g,dindx,ib) this%alb_improad_dir(l,ib) = urbinp%alb_improad_dir(g,dindx,ib) this%alb_perroad_dir(l,ib) = urbinp%alb_perroad_dir(g,dindx,ib) this%alb_improad_dif(l,ib) = urbinp%alb_improad_dif(g,dindx,ib) this%alb_perroad_dif(l,ib) = urbinp%alb_perroad_dif(g,dindx,ib) this%alb_wall_dir (l,ib) = urbinp%alb_wall_dir (g,dindx,ib) this%alb_wall_dif (l,ib) = urbinp%alb_wall_dif (g,dindx,ib) end do this%em_roof (l) = urbinp%em_roof (g,dindx) this%em_improad(l) = urbinp%em_improad(g,dindx) this%em_perroad(l) = urbinp%em_perroad(g,dindx) this%em_wall (l) = urbinp%em_wall (g,dindx) ! Landunit level initialization for urban wall and roof layers and interfaces lun%canyon_hwr(l) = urbinp%canyon_hwr(g,dindx) lun%wtroad_perv(l) = urbinp%wtroad_perv(g,dindx) lun%ht_roof(l) = urbinp%ht_roof(g,dindx) lun%wtlunit_roof(l) = urbinp%wtlunit_roof(g,dindx) this%tk_wall(l,:) = urbinp%tk_wall(g,dindx,:) this%tk_roof(l,:) = urbinp%tk_roof(g,dindx,:) this%tk_improad(l,:) = urbinp%tk_improad(g,dindx,:) this%cv_wall(l,:) = urbinp%cv_wall(g,dindx,:) this%cv_roof(l,:) = urbinp%cv_roof(g,dindx,:) this%cv_improad(l,:) = urbinp%cv_improad(g,dindx,:) this%thick_wall(l) = urbinp%thick_wall(g,dindx) this%thick_roof(l) = urbinp%thick_roof(g,dindx) this%nlev_improad(l) = urbinp%nlev_improad(g,dindx) this%t_building_min(l) = urbinp%t_building_min(g,dindx) ! Inferred from Sailor and Lu 2004 if (urban_traffic) then this%eflx_traffic_factor(l) = 3.6_r8 * (lun%canyon_hwr(l)-0.5_r8) + 1.0_r8 else this%eflx_traffic_factor(l) = 0.0_r8 end if if (use_vancouver .or. use_mexicocity) then ! Freely evolving this%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 this%t_building_min(l) = 200.00_r8 end if end if !---------------------------------------------------------------------------------- ! View factors for road and one wall in urban canyon (depends only on canyon_hwr) ! --------------------------------------------------------------------------------------- ! WALL | ! ROAD | ! wall | ! -----\ /----- - - |\----------/ ! | \ vsr / | | r | | \ vww / s ! | \ / | h o w | \ / k ! wall | \ / | wall | a | | \ / y ! |vwr \ / vwr| | d | |vrw \ / vsw ! ------\/------ - - |-----\/----- ! road wall | ! <----- w ----> | ! <---- h --->| ! ! vsr = view factor of sky for road vrw = view factor of road for wall ! vwr = view factor of one wall for road vww = view factor of opposing wall for wall ! vsw = view factor of sky for wall ! vsr + vwr + vwr = 1 vrw + vww + vsw = 1 ! ! Source: Masson, V. (2000) A physically-based scheme for the urban energy budget in ! atmospheric models. Boundary-Layer Meteorology 94:357-397 ! ! - Calculate urban land unit aerodynamic constants using Macdonald (1998) as used in ! Grimmond and Oke (1999) ! --------------------------------------------------------------------------------------- ! road -- sky view factor -> 1 as building height -> 0 ! and -> 0 as building height -> infinity this%vf_sr(l) = sqrt(lun%canyon_hwr(l)**2 + 1._r8) - lun%canyon_hwr(l) this%vf_wr(l) = 0.5_r8 * (1._r8 - this%vf_sr(l)) ! one wall -- sky view factor -> 0.5 as building height -> 0 ! and -> 0 as building height -> infinity this%vf_sw(l) = 0.5_r8 * (lun%canyon_hwr(l) + 1._r8 - sqrt(lun%canyon_hwr(l)**2+1._r8)) / lun%canyon_hwr(l) this%vf_rw(l) = this%vf_sw(l) this%vf_ww(l) = 1._r8 - this%vf_sw(l) - this%vf_rw(l) ! error check -- make sure view factor sums to one for road and wall sumvf = this%vf_sr(l) + 2._r8*this%vf_wr(l) if (abs(sumvf-1._r8) > 1.e-06_r8 ) then write (iulog,*) 'urban road view factor error',sumvf write (iulog,*) 'clm model is stopping' call endrun(decomp_index=l, clmlevel=namel, msg=errmsg(sourcefile, __LINE__)) endif sumvf = this%vf_sw(l) + this%vf_rw(l) + this%vf_ww(l) if (abs(sumvf-1._r8) > 1.e-06_r8 ) then write (iulog,*) 'urban wall view factor error',sumvf write (iulog,*) 'clm model is stopping' call endrun(decomp_index=l, clmlevel=namel, msg=errmsg(sourcefile, __LINE__)) endif !---------------------------------------------------------------------------------- ! Calculate urban land unit aerodynamic constants using Macdonald (1998) as used in ! Grimmond and Oke (1999) !---------------------------------------------------------------------------------- ! Calculate plan area index plan_ai = lun%canyon_hwr(l)/(lun%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) * lun%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 (use_vancouver) then lun%z_d_town(l) = 3.5_r8 else if (use_mexicocity) then lun%z_d_town(l) = 10.9_r8 else lun%z_d_town(l) = (1._r8 + alpha**(-plan_ai) * (plan_ai - 1._r8)) * lun%ht_roof(l) end if ! Calculate the roughness length if (use_vancouver) then lun%z_0_town(l) = 0.35_r8 else if (use_mexicocity) then lun%z_0_town(l) = 2.2_r8 else lun%z_0_town(l) = lun%ht_roof(l) * (1._r8 - lun%z_d_town(l) / lun%ht_roof(l)) * & exp(-1.0_r8 * (0.5_r8 * beta * C_d / vkc**2 * & (1 - lun%z_d_town(l) / lun%ht_roof(l)) * frontal_ai)**(-0.5_r8)) end if else ! Not urban point this%eflx_traffic_factor(l) = spval this%t_building_min(l) = spval this%vf_sr(l) = spval this%vf_wr(l) = spval this%vf_sw(l) = spval this%vf_rw(l) = spval this%vf_ww(l) = spval end if end do ! Deallocate memory for urbinp datatype call UrbanInput(bounds%begg, bounds%endg, mode='finalize') end subroutine Init !----------------------------------------------------------------------- subroutine UrbanInput(begg, endg, mode) ! ! !DESCRIPTION: ! Allocate memory and read in urban input data ! ! !USES: use clm_varpar , only : numrad, nlevurb use landunit_varcon , only : numurbl use fileutils , only : getavu, relavu, getfil, opnfil use spmdMod , only : masterproc use domainMod , only : ldomain use ncdio_pio , only : file_desc_t, ncd_io, ncd_inqvdlen, ncd_inqfdims use ncdio_pio , only : ncd_pio_openfile, ncd_pio_closefile, ncd_inqdid, ncd_inqdlen ! ! !ARGUMENTS: implicit none integer, intent(in) :: begg, endg character(len=*), intent(in) :: mode ! ! !LOCAL VARIABLES: character(len=256) :: locfn ! local file name type(file_desc_t) :: ncid ! netcdf id integer :: dimid ! netCDF id integer :: nw,n,k,i,j,ni,nj,ns ! indices integer :: nlevurb_i ! input grid: number of urban vertical levels integer :: numrad_i ! input grid: number of solar bands (VIS/NIR) integer :: numurbl_i ! input grid: number of urban landunits integer :: ier,ret ! error status logical :: isgrid2d ! true => file is 2d logical :: readvar ! true => variable is on dataset logical :: has_numurbl ! true => numurbl dimension is on dataset character(len=32) :: subname = 'UrbanInput' ! subroutine name !----------------------------------------------------------------------- if ( nlevurb == 0 ) return if (mode == 'initialize') then ! Read urban data if (masterproc) then write(iulog,*)' Reading in urban input data from fsurdat file ...' end if call getfil (fsurdat, locfn, 0) call ncd_pio_openfile (ncid, locfn, 0) if (masterproc) then write(iulog,*) subname,trim(fsurdat) end if ! Check whether this file has new-format urban data call ncd_inqdid(ncid, 'numurbl', dimid, dimexist=has_numurbl) ! If file doesn't have numurbl, then it is old-format urban; ! in this case, set nlevurb to zero if (.not. has_numurbl) then nlevurb = 0 if (masterproc) write(iulog,*)'PCT_URBAN is not multi-density, nlevurb set to 0' end if if ( nlevurb == 0 ) return ! Allocate dynamic memory allocate(urbinp%canyon_hwr(begg:endg, numurbl), & urbinp%wtlunit_roof(begg:endg, numurbl), & urbinp%wtroad_perv(begg:endg, numurbl), & urbinp%em_roof(begg:endg, numurbl), & urbinp%em_improad(begg:endg, numurbl), & urbinp%em_perroad(begg:endg, numurbl), & urbinp%em_wall(begg:endg, numurbl), & urbinp%alb_roof_dir(begg:endg, numurbl, numrad), & urbinp%alb_roof_dif(begg:endg, numurbl, numrad), & urbinp%alb_improad_dir(begg:endg, numurbl, numrad), & urbinp%alb_perroad_dir(begg:endg, numurbl, numrad), & urbinp%alb_improad_dif(begg:endg, numurbl, numrad), & urbinp%alb_perroad_dif(begg:endg, numurbl, numrad), & urbinp%alb_wall_dir(begg:endg, numurbl, numrad), & urbinp%alb_wall_dif(begg:endg, numurbl, numrad), & urbinp%ht_roof(begg:endg, numurbl), & urbinp%wind_hgt_canyon(begg:endg, numurbl), & urbinp%tk_wall(begg:endg, numurbl,nlevurb), & urbinp%tk_roof(begg:endg, numurbl,nlevurb), & urbinp%tk_improad(begg:endg, numurbl,nlevurb), & urbinp%cv_wall(begg:endg, numurbl,nlevurb), & urbinp%cv_roof(begg:endg, numurbl,nlevurb), & urbinp%cv_improad(begg:endg, numurbl,nlevurb), & urbinp%thick_wall(begg:endg, numurbl), & urbinp%thick_roof(begg:endg, numurbl), & urbinp%nlev_improad(begg:endg, numurbl), & urbinp%t_building_min(begg:endg, numurbl), & stat=ier) if (ier /= 0) then call endrun(msg="Allocation error "//errmsg(sourcefile, __LINE__)) endif call ncd_inqfdims (ncid, isgrid2d, ni, nj, ns) if (ldomain%ns /= ns .or. ldomain%ni /= ni .or. ldomain%nj /= nj) then write(iulog,*)trim(subname), 'ldomain and input file do not match dims ' write(iulog,*)trim(subname), 'ldomain%ni,ni,= ',ldomain%ni,ni write(iulog,*)trim(subname), 'ldomain%nj,nj,= ',ldomain%nj,nj write(iulog,*)trim(subname), 'ldomain%ns,ns,= ',ldomain%ns,ns call endrun(msg=errmsg(sourcefile, __LINE__)) end if call ncd_inqdid(ncid, 'nlevurb', dimid) call ncd_inqdlen(ncid, dimid, nlevurb_i) if (nlevurb_i /= nlevurb) then write(iulog,*)trim(subname)// ': parameter nlevurb= ',nlevurb, & 'does not equal input dataset nlevurb= ',nlevurb_i call endrun(msg=errmsg(sourcefile, __LINE__)) endif call ncd_inqdid(ncid, 'numrad', dimid) call ncd_inqdlen(ncid, dimid, numrad_i) if (numrad_i /= numrad) then write(iulog,*)trim(subname)// ': parameter numrad= ',numrad, & 'does not equal input dataset numrad= ',numrad_i call endrun(msg=errmsg(sourcefile, __LINE__)) endif call ncd_inqdid(ncid, 'numurbl', dimid) call ncd_inqdlen(ncid, dimid, numurbl_i) if (numurbl_i /= numurbl) then write(iulog,*)trim(subname)// ': parameter numurbl= ',numurbl, & 'does not equal input dataset numurbl= ',numurbl_i call endrun(msg=errmsg(sourcefile, __LINE__)) endif call ncd_io(ncid=ncid, varname='CANYON_HWR', flag='read', data=urbinp%canyon_hwr,& dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg='ERROR: CANYON_HWR NOT on fsurdat file '//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='WTLUNIT_ROOF', flag='read', data=urbinp%wtlunit_roof, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: WTLUNIT_ROOF NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='WTROAD_PERV', flag='read', data=urbinp%wtroad_perv, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: WTROAD_PERV NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='EM_ROOF', flag='read', data=urbinp%em_roof, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: EM_ROOF NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='EM_IMPROAD', flag='read', data=urbinp%em_improad, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: EM_IMPROAD NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='EM_PERROAD', flag='read', data=urbinp%em_perroad, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: EM_PERROAD NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='EM_WALL', flag='read', data=urbinp%em_wall, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: EM_WALL NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='HT_ROOF', flag='read', data=urbinp%ht_roof, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: HT_ROOF NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='WIND_HGT_CANYON', flag='read', data=urbinp%wind_hgt_canyon, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: WIND_HGT_CANYON NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='THICK_WALL', flag='read', data=urbinp%thick_wall, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: THICK_WALL NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='THICK_ROOF', flag='read', data=urbinp%thick_roof, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: THICK_ROOF NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='NLEV_IMPROAD', flag='read', data=urbinp%nlev_improad, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: NLEV_IMPROAD NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='T_BUILDING_MIN', flag='read', data=urbinp%t_building_min, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: T_BUILDING_MIN NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='ALB_IMPROAD_DIR', flag='read', data=urbinp%alb_improad_dir, & dim1name=grlnd, readvar=readvar) if (.not.readvar) then call endrun( msg=' ERROR: ALB_IMPROAD_DIR NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='ALB_IMPROAD_DIF', flag='read', data=urbinp%alb_improad_dif, & dim1name=grlnd, readvar=readvar) if (.not.readvar) then call endrun( msg=' ERROR: ALB_IMPROAD_DIF NOT on fsurdat file'//errmsg(sourcefile, __LINE__) ) end if call ncd_io(ncid=ncid, varname='ALB_PERROAD_DIR', flag='read',data=urbinp%alb_perroad_dir, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: ALB_PERROAD_DIR NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='ALB_PERROAD_DIF', flag='read',data=urbinp%alb_perroad_dif, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: ALB_PERROAD_DIF NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='ALB_ROOF_DIR', flag='read', data=urbinp%alb_roof_dir, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: ALB_ROOF_DIR NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='ALB_ROOF_DIF', flag='read', data=urbinp%alb_roof_dif, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: ALB_ROOF_DIF NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='ALB_WALL_DIR', flag='read', data=urbinp%alb_wall_dir, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: ALB_WALL_DIR NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='ALB_WALL_DIF', flag='read', data=urbinp%alb_wall_dif, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: ALB_WALL_DIF NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='TK_IMPROAD', flag='read', data=urbinp%tk_improad, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: TK_IMPROAD NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='TK_ROOF', flag='read', data=urbinp%tk_roof, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: TK_ROOF NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='TK_WALL', flag='read', data=urbinp%tk_wall, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: TK_WALL NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='CV_IMPROAD', flag='read', data=urbinp%cv_improad, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: CV_IMPROAD NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='CV_ROOF', flag='read', data=urbinp%cv_roof, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: CV_ROOF NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_io(ncid=ncid, varname='CV_WALL', flag='read', data=urbinp%cv_wall, & dim1name=grlnd, readvar=readvar) if (.not. readvar) then call endrun( msg=' ERROR: CV_WALL NOT on fsurdat file'//errmsg(sourcefile, __LINE__)) end if call ncd_pio_closefile(ncid) if (masterproc) then write(iulog,*)' Sucessfully read urban input data' write(iulog,*) end if else if (mode == 'finalize') then if ( nlevurb == 0 ) return deallocate(urbinp%canyon_hwr, & urbinp%wtlunit_roof, & urbinp%wtroad_perv, & urbinp%em_roof, & urbinp%em_improad, & urbinp%em_perroad, & urbinp%em_wall, & urbinp%alb_roof_dir, & urbinp%alb_roof_dif, & urbinp%alb_improad_dir, & urbinp%alb_perroad_dir, & urbinp%alb_improad_dif, & urbinp%alb_perroad_dif, & urbinp%alb_wall_dir, & urbinp%alb_wall_dif, & urbinp%ht_roof, & urbinp%wind_hgt_canyon, & urbinp%tk_wall, & urbinp%tk_roof, & urbinp%tk_improad, & urbinp%cv_wall, & urbinp%cv_roof, & urbinp%cv_improad, & urbinp%thick_wall, & urbinp%thick_roof, & urbinp%nlev_improad, & urbinp%t_building_min, & stat=ier) if (ier /= 0) then call endrun(msg='initUrbanInput: deallocation error '//errmsg(sourcefile, __LINE__)) end if else write(iulog,*)'initUrbanInput error: mode ',trim(mode),' not supported ' call endrun(msg=errmsg(sourcefile, __LINE__)) end if end subroutine UrbanInput !----------------------------------------------------------------------- subroutine CheckUrban(begg, endg, pcturb, caller) !----------------------------------------------------------------------- ! !DESCRIPTION: ! Confirm that we have valid urban data for all points with pct urban > 0. If this isn't ! true, abort with a message. ! ! !USES: use clm_instur , only : urban_valid use landunit_varcon , only : numurbl ! ! !ARGUMENTS: implicit none integer , intent(in) :: begg, endg ! beg & end grid cell indices real(r8) , intent(in) :: pcturb(begg:,:) ! % urban character(len=*), intent(in) :: caller ! identifier of caller, for more meaningful error messages ! ! !REVISION HISTORY: ! Created by Bill Sacks 7/2013, mostly by moving code from surfrd_special ! ! !LOCAL VARIABLES: logical :: found integer :: nl, n integer :: nindx, dindx integer :: nlev !----------------------------------------------------------------------- found = .false. do nl = begg,endg do n = 1, numurbl if ( pcturb(nl,n) > 0.0_r8 ) then if ( .not. urban_valid(nl) .or. & urbinp%canyon_hwr(nl,n) <= 0._r8 .or. & urbinp%em_improad(nl,n) <= 0._r8 .or. & urbinp%em_perroad(nl,n) <= 0._r8 .or. & urbinp%em_roof(nl,n) <= 0._r8 .or. & urbinp%em_wall(nl,n) <= 0._r8 .or. & urbinp%ht_roof(nl,n) <= 0._r8 .or. & urbinp%thick_roof(nl,n) <= 0._r8 .or. & urbinp%thick_wall(nl,n) <= 0._r8 .or. & urbinp%t_building_min(nl,n) <= 0._r8 .or. & urbinp%wind_hgt_canyon(nl,n) <= 0._r8 .or. & urbinp%wtlunit_roof(nl,n) <= 0._r8 .or. & urbinp%wtroad_perv(nl,n) <= 0._r8 .or. & any(urbinp%alb_improad_dir(nl,n,:) <= 0._r8) .or. & any(urbinp%alb_improad_dif(nl,n,:) <= 0._r8) .or. & any(urbinp%alb_perroad_dir(nl,n,:) <= 0._r8) .or. & any(urbinp%alb_perroad_dif(nl,n,:) <= 0._r8) .or. & any(urbinp%alb_roof_dir(nl,n,:) <= 0._r8) .or. & any(urbinp%alb_roof_dif(nl,n,:) <= 0._r8) .or. & any(urbinp%alb_wall_dir(nl,n,:) <= 0._r8) .or. & any(urbinp%alb_wall_dif(nl,n,:) <= 0._r8) .or. & any(urbinp%tk_roof(nl,n,:) <= 0._r8) .or. & any(urbinp%tk_wall(nl,n,:) <= 0._r8) .or. & any(urbinp%cv_roof(nl,n,:) <= 0._r8) .or. & any(urbinp%cv_wall(nl,n,:) <= 0._r8)) then found = .true. nindx = nl dindx = n exit else if (urbinp%nlev_improad(nl,n) > 0) then nlev = urbinp%nlev_improad(nl,n) if ( any(urbinp%tk_improad(nl,n,1:nlev) <= 0._r8) .or. & any(urbinp%cv_improad(nl,n,1:nlev) <= 0._r8)) then found = .true. nindx = nl dindx = n exit end if end if end if if (found) exit end if end do end do if ( found ) then write(iulog,*) trim(caller), ' ERROR: no valid urban data for nl=',nindx write(iulog,*)'density type: ',dindx write(iulog,*)'urban_valid: ',urban_valid(nindx) write(iulog,*)'canyon_hwr: ',urbinp%canyon_hwr(nindx,dindx) write(iulog,*)'em_improad: ',urbinp%em_improad(nindx,dindx) write(iulog,*)'em_perroad: ',urbinp%em_perroad(nindx,dindx) write(iulog,*)'em_roof: ',urbinp%em_roof(nindx,dindx) write(iulog,*)'em_wall: ',urbinp%em_wall(nindx,dindx) write(iulog,*)'ht_roof: ',urbinp%ht_roof(nindx,dindx) write(iulog,*)'thick_roof: ',urbinp%thick_roof(nindx,dindx) write(iulog,*)'thick_wall: ',urbinp%thick_wall(nindx,dindx) write(iulog,*)'t_building_min: ',urbinp%t_building_min(nindx,dindx) write(iulog,*)'wind_hgt_canyon: ',urbinp%wind_hgt_canyon(nindx,dindx) write(iulog,*)'wtlunit_roof: ',urbinp%wtlunit_roof(nindx,dindx) write(iulog,*)'wtroad_perv: ',urbinp%wtroad_perv(nindx,dindx) write(iulog,*)'alb_improad_dir: ',urbinp%alb_improad_dir(nindx,dindx,:) write(iulog,*)'alb_improad_dif: ',urbinp%alb_improad_dif(nindx,dindx,:) write(iulog,*)'alb_perroad_dir: ',urbinp%alb_perroad_dir(nindx,dindx,:) write(iulog,*)'alb_perroad_dif: ',urbinp%alb_perroad_dif(nindx,dindx,:) write(iulog,*)'alb_roof_dir: ',urbinp%alb_roof_dir(nindx,dindx,:) write(iulog,*)'alb_roof_dif: ',urbinp%alb_roof_dif(nindx,dindx,:) write(iulog,*)'alb_wall_dir: ',urbinp%alb_wall_dir(nindx,dindx,:) write(iulog,*)'alb_wall_dif: ',urbinp%alb_wall_dif(nindx,dindx,:) write(iulog,*)'tk_roof: ',urbinp%tk_roof(nindx,dindx,:) write(iulog,*)'tk_wall: ',urbinp%tk_wall(nindx,dindx,:) write(iulog,*)'cv_roof: ',urbinp%cv_roof(nindx,dindx,:) write(iulog,*)'cv_wall: ',urbinp%cv_wall(nindx,dindx,:) if (urbinp%nlev_improad(nindx,dindx) > 0) then nlev = urbinp%nlev_improad(nindx,dindx) write(iulog,*)'tk_improad: ',urbinp%tk_improad(nindx,dindx,1:nlev) write(iulog,*)'cv_improad: ',urbinp%cv_improad(nindx,dindx,1:nlev) end if call endrun(msg=errmsg(sourcefile, __LINE__)) end if end subroutine CheckUrban !----------------------------------------------------------------------- !----------------------------------------------------------------------- !BOP ! ! !IROUTINE: UrbanReadNML ! ! !INTERFACE: ! subroutine UrbanReadNML ( NLFilename ) ! ! !DESCRIPTION: ! ! Read in the urban namelist ! ! !USES: use shr_mpi_mod , only : shr_mpi_bcast use abortutils , only : endrun use spmdMod , only : masterproc, mpicom use fileutils , only : getavu, relavu, opnfil use shr_nl_mod , only : shr_nl_find_group_name use shr_mpi_mod , only : shr_mpi_bcast implicit none ! ! !ARGUMENTS: character(len=*), intent(IN) :: NLFilename ! Namelist filename ! ! !LOCAL VARIABLES: integer :: ierr ! error code integer :: unitn ! unit for namelist file character(len=32) :: subname = 'UrbanReadNML' ! subroutine name namelist / clmu_inparm / urban_hac, urban_traffic, building_temp_method !EOP !----------------------------------------------------------------------- ! ---------------------------------------------------------------------- ! Read namelist from input namelist filename ! ---------------------------------------------------------------------- if ( masterproc )then unitn = getavu() write(iulog,*) 'Read in clmu_inparm namelist' call opnfil (NLFilename, unitn, 'F') call shr_nl_find_group_name(unitn, 'clmu_inparm', status=ierr) if (ierr == 0) then read(unitn, clmu_inparm, iostat=ierr) if (ierr /= 0) then call endrun(msg="ERROR reading clmu_inparm namelist"//errmsg(sourcefile, __LINE__)) end if else call endrun(msg="ERROR finding clmu_inparm namelist"//errmsg(sourcefile, __LINE__)) end if call relavu( unitn ) end if ! Broadcast namelist variables read in call shr_mpi_bcast(urban_hac, mpicom) call shr_mpi_bcast(urban_traffic, mpicom) call shr_mpi_bcast(building_temp_method, mpicom) ! if (urban_traffic) then write(iulog,*)'Urban traffic fluxes are not implemented currently' call endrun(msg=errMsg(sourcefile, __LINE__)) end if ! if ( masterproc )then write(iulog,*) ' urban air conditioning/heating and wasteheat = ', urban_hac write(iulog,*) ' urban traffic flux = ', urban_traffic end if ReadNamelist = .true. end subroutine UrbanReadNML !----------------------------------------------------------------------- !----------------------------------------------------------------------- !BOP ! ! !IROUTINE: IsSimpleBuildTemp ! ! !INTERFACE: ! logical function IsSimpleBuildTemp( ) ! ! !DESCRIPTION: ! ! If the simple building temperature method is being used ! ! !USES: implicit none !EOP !----------------------------------------------------------------------- if ( .not. ReadNamelist )then write(iulog,*)'Testing on building_temp_method before urban namelist was read in' call endrun(msg=errMsg(sourcefile, __LINE__)) end if IsSimpleBuildTemp = building_temp_method == BUILDING_TEMP_METHOD_SIMPLE end function IsSimpleBuildTemp !----------------------------------------------------------------------- !----------------------------------------------------------------------- !BOP ! ! !IROUTINE: IsProgBuildTemp ! ! !INTERFACE: ! logical function IsProgBuildTemp( ) ! ! !DESCRIPTION: ! ! If the prognostic building temperature method is being used ! ! !USES: implicit none !EOP !----------------------------------------------------------------------- if ( .not. ReadNamelist )then write(iulog,*)'Testing on building_temp_method before urban namelist was read in' call endrun(msg=errMsg(sourcefile, __LINE__)) end if IsProgBuildTemp = building_temp_method == BUILDING_TEMP_METHOD_PROG end function IsProgBuildTemp !----------------------------------------------------------------------- end module UrbanParamsType