diff --git a/components/elm/src/biogeochem/CH4Mod.F90 b/components/elm/src/biogeochem/CH4Mod.F90 index a5cbd2f041c2..f6203665b7ad 100644 --- a/components/elm/src/biogeochem/CH4Mod.F90 +++ b/components/elm/src/biogeochem/CH4Mod.F90 @@ -790,7 +790,7 @@ subroutine InitCold(this, bounds, cellorg_col) this%o2_decomp_depth_unsat_col(c,:)= spval l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%conc_ch4_sat_col (c,1:nlevsoi) = 0._r8 this%conc_ch4_unsat_col (c,1:nlevsoi) = 0._r8 @@ -806,7 +806,7 @@ subroutine InitCold(this, bounds, cellorg_col) ! Note that finundated will be overwritten with this%fsat_bef_col upon reading ! a restart file - either in a continuation, branch or startup spun-up case - else if (lun_pp%itype(l) == istdlak) then + else if (col_pp%is_lake(c)) then this%conc_ch4_sat_col(c,1:nlevsoi) = 0._r8 this%conc_o2_sat_col (c,1:nlevsoi) = 0._r8 @@ -851,7 +851,7 @@ subroutine InitCold(this, bounds, cellorg_col) this%ch4stress_unsat_col (c,nlevsoi+1:nlevgrnd) = 0._r8 this%ch4stress_sat_col (c,nlevsoi+1:nlevgrnd) = 0._r8 - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%conc_ch4_lake_col (c,:) = spval this%conc_o2_lake_col (c,:) = spval @@ -860,7 +860,7 @@ subroutine InitCold(this, bounds, cellorg_col) this%ch4_prod_depth_lake_col (c,:) = spval this%ch4_oxid_depth_lake_col (c,:) = spval - else if (lun_pp%itype(l) == istdlak .and. allowlakeprod) then + else if (col_pp%is_lake(c) .and. allowlakeprod) then this%ch4_prod_depth_unsat_col (c,:) = spval this%ch4_oxid_depth_unsat_col (c,:) = spval diff --git a/components/elm/src/biogeochem/CNCarbonFluxType.F90 b/components/elm/src/biogeochem/CNCarbonFluxType.F90 index a6c56a935f05..f2af432a869b 100644 --- a/components/elm/src/biogeochem/CNCarbonFluxType.F90 +++ b/components/elm/src/biogeochem/CNCarbonFluxType.F90 @@ -1035,7 +1035,7 @@ subroutine InitCold(this, bounds) this%xsmrpool_c13ratio_patch(p) = spval endif end if - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then this%tempsum_npp_patch(p) = 0._r8 this%annsum_npp_patch(p) = 0._r8 this%availc_patch(p) = 0._r8 @@ -1057,9 +1057,9 @@ subroutine InitCold(this, bounds) end if this%fphr_col(c,nlevdecomp+1:nlevgrnd) = 0._r8 !used to be in CH4Mod - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%fphr_col(c,nlevdecomp+1:nlevgrnd) = 0._r8 - else if (lun_pp%itype(l) == istdlak .and. allowlakeprod) then + else if (col_pp%is_lake(c) .and. allowlakeprod) then this%fphr_col(c,:) = spval else ! Inactive CH4 columns this%fphr_col(c,:) = spval @@ -1067,7 +1067,7 @@ subroutine InitCold(this, bounds) ! also initialize dynamic landcover fluxes so that they have ! real values on first timestep, prior to calling pftdyn_cnbal - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%dwt_conv_cflux_col(c) = 0._r8 this%dwt_prod10c_gain_col(c) = 0._r8 this%dwt_prod100c_gain_col(c) = 0._r8 diff --git a/components/elm/src/biogeochem/CNCarbonStateType.F90 b/components/elm/src/biogeochem/CNCarbonStateType.F90 index b733773fb20d..baa2b3d9b816 100644 --- a/components/elm/src/biogeochem/CNCarbonStateType.F90 +++ b/components/elm/src/biogeochem/CNCarbonStateType.F90 @@ -450,7 +450,7 @@ subroutine InitCold(this, bounds, ratio, c12_carbonstate_vars) this%leafcmax_patch(p) = 0._r8 l = veg_pp%landunit(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then if (veg_pp%itype(p) == noveg) then this%leafc_patch(p) = 0._r8 @@ -569,7 +569,7 @@ subroutine InitCold(this, bounds, ratio, c12_carbonstate_vars) ! initialize column-level variables do c = bounds%begc, bounds%endc l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%cwdc_col(c) = 0._r8 diff --git a/components/elm/src/biogeochem/CNNitrogenStateType.F90 b/components/elm/src/biogeochem/CNNitrogenStateType.F90 index a7f7315c6c83..93911c4f6390 100644 --- a/components/elm/src/biogeochem/CNNitrogenStateType.F90 +++ b/components/elm/src/biogeochem/CNNitrogenStateType.F90 @@ -562,7 +562,7 @@ subroutine InitCold(this, bounds, & do p = bounds%begp,bounds%endp l = veg_pp%landunit(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then if (veg_pp%itype(p) == noveg) then this%leafn_patch(p) = 0._r8 this%leafn_storage_patch(p) = 0._r8 @@ -635,7 +635,7 @@ subroutine InitCold(this, bounds, & do c = bounds%begc, bounds%endc l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then ! column nitrogen state variables this%ntrunc_col(c) = 0._r8 diff --git a/components/elm/src/biogeochem/CropType.F90 b/components/elm/src/biogeochem/CropType.F90 index 1c5e6b36ad44..ba4d337c3d8d 100644 --- a/components/elm/src/biogeochem/CropType.F90 +++ b/components/elm/src/biogeochem/CropType.F90 @@ -331,7 +331,7 @@ subroutine InitCold(this, bounds) t = veg_pp%topounit(p) topi = grc_pp%topi(g) - if (lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_crop_col(p)) then ti = t - topi + 1 m = veg_pp%itype(p) this%fertnitro_patch(p) = fert_cft(g,ti,m) diff --git a/components/elm/src/biogeochem/DUSTMod.F90 b/components/elm/src/biogeochem/DUSTMod.F90 index f53e970423be..16a8b4ef5c34 100644 --- a/components/elm/src/biogeochem/DUSTMod.F90 +++ b/components/elm/src/biogeochem/DUSTMod.F90 @@ -184,15 +184,13 @@ subroutine InitCold(this, bounds) type(bounds_type), intent(in) :: bounds ! ! !LOCAL VARIABLES: - integer :: c,l + integer :: c !----------------------------------------------------------------------- ! Set basin factor to 1 for now do c = bounds%begc, bounds%endc - l = col_pp%landunit(c) - - if (.not.lun_pp%lakpoi(l)) then + if (.not.col_pp%is_lake(c)) then this%mbl_bsn_fct_col(c) = 1.0_r8 end if end do @@ -350,7 +348,7 @@ subroutine DustEmission (bounds, & ! linearly from 1 to 0 as VAI(=tlai+tsai) increases from 0 to vai_mbl_thr ! if ice sheet, wetland, or lake, no dust allowed - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then if (tlai_lu(l) < vai_mbl_thr) then lnd_frc_mbl(p) = 1.0_r8 - (tlai_lu(l))/vai_mbl_thr else diff --git a/components/elm/src/biogeochem/FanMod.F90 b/components/elm/src/biogeochem/FanMod.F90 old mode 100755 new mode 100644 diff --git a/components/elm/src/biogeochem/FanUpdateMod.F90 b/components/elm/src/biogeochem/FanUpdateMod.F90 index fbb0e0467ab9..95981b052966 100644 --- a/components/elm/src/biogeochem/FanUpdateMod.F90 +++ b/components/elm/src/biogeochem/FanUpdateMod.F90 @@ -264,7 +264,7 @@ subroutine fan_eval(bounds, num_soilc, filter_soilc, & if (.not. col_pp%active(c) .or. col_pp%wtgcell(c) < 1e-6_r8) cycle g = col_pp%gridcell(c) - if (lun_pp%itype(l) == istsoil) then + if (col_pp%is_soil(c)) then col_nf%manure_n_grz(c) & = atm2lnd_vars%forc_ndep_past_grc(g) / col_pp%wtgcell(c) * 1e3_r8 ! kg to g if (col_nf%manure_n_grz(c) > 1e12 .or. isnan(col_nf%manure_n_grz(c))) then @@ -314,12 +314,12 @@ subroutine fan_eval(bounds, num_soilc, filter_soilc, & c = filter_soilc(fc) l = col_pp%landunit(c) g = col_pp%gridcell(c) - if (.not. (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop)) cycle + if (.not. (col_pp%is_soil(c) .or. col_pp%is_crop(c))) cycle if (.not. col_pp%active(c) .or. col_pp%wtgcell(c) < 1e-15_r8) cycle ! Find and average the atmospheric resistances Rb and Ra. ! - if (lun_pp%itype(col_pp%landunit(c)) == istcrop) then + if (veg_pp%is_on_crop_col(p)) then ! Crop column, only one patch p = col_pp%pfti(c) if (p /= col_pp%pftf(c)) call endrun(msg='Strange patch for crop') @@ -775,7 +775,7 @@ subroutine handle_storage(bounds, frictionvel_vars, dt, & col_grass = ispval l = grc_pp%landunit_indices(istsoil, g) do c = lun_pp%coli(l), lun_pp%colf(l) - if (col_pp%itype(c) == istsoil) then + if (col_pp%is_soil(c)) then col_grass = c exit end if @@ -885,11 +885,12 @@ subroutine update_summary(filter_soilc, num_soilc) integer, intent(in) :: num_soilc ! number of soil columns in filter integer, intent(in) :: filter_soilc(:) ! filter for soil columns - integer :: c, fc + integer :: c, l, fc real(r8) :: total, fluxout, fluxin, flux_loss do fc = 1, num_soilc c = filter_soilc(fc) + l = col_pp%landunit(c) if (.not. col_pp%active(c) .or. col_pp%wtgcell(c) < 1.e-15_r8) cycle total = col_ns%tan_g1(c) + col_ns%tan_g2(c) + col_ns%tan_g3(c) total = total + col_ns%manure_u_grz(c) + col_ns%manure_a_grz(c) + col_ns%manure_r_grz(c) @@ -900,7 +901,7 @@ subroutine update_summary(filter_soilc, num_soilc) total = total + col_ns%manure_n_stored(c) col_ns%fan_totn(c) = total - if (lun_pp%itype(col_pp%landunit(c)) == istcrop) then + if (col_pp%is_crop(c)) then ! no grazing, manure_n_appl is from the same column and not counted as input fluxin = col_nf%manure_n_mix(c) + col_nf%fert_n_appl(c) else @@ -939,7 +940,7 @@ subroutine fan_to_sminn(bounds, filter_soilc, num_soilc, nfertilization) ! patch level fertilizer application + manure production real(r8), intent(inout) :: nfertilization(bounds%begp:) - integer :: c, fc, p + integer :: c, l, fc, p real(r8) :: flux_manure, flux_fert, manure_prod logical :: included @@ -947,12 +948,13 @@ subroutine fan_to_sminn(bounds, filter_soilc, num_soilc, nfertilization) do fc = 1, num_soilc c = filter_soilc(fc) + l = col_pp%landunit(c) flux_manure = col_nf%manure_no3_to_soil(c) + col_nf%manure_nh4_to_soil(c) flux_fert = col_nf%fert_no3_to_soil(c) + col_nf%fert_nh4_to_soil(c) manure_prod = col_nf%manure_n_barns(c) + col_nf%manure_n_grz(c) - included = (lun_pp%itype(col_pp%landunit(c)) == istcrop .and. fan_to_bgc_crop) & - .or. (lun_pp%itype(col_pp%landunit(c)) == istsoil .and. fan_to_bgc_veg) + included = (col_pp%is_crop(c) .and. fan_to_bgc_crop) & + .or. (col_pp%is_soil(c) .and. fan_to_bgc_veg) if (included) then col_nf%fert_to_sminn(c) = flux_fert + flux_manure diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 old mode 100755 new mode 100644 index c39487e5e91d..d16e2eaef4a4 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -462,7 +462,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & snow_sinks(c) = qflx_sub_snow(c) + qflx_evap_grnd(c) + qflx_snow_melt(c) & + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) + qflx_sl_top_soil(c) - if (lun_pp%itype(l) == istdlak) then + if (col_pp%is_lake(c)) then if (.not. use_firn_percolation_and_compaction) then if ( do_capsnow(c)) then snow_sources(c) = qflx_snow_grnd_col(c) & @@ -490,7 +490,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & endif endif - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop .or. lun_pp%itype(l) == istwet ) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c) .or. lun_pp%itype(l) == istwet ) then if (.not. use_firn_percolation_and_compaction) then if (do_capsnow(c)) then snow_sources(c) = frac_sno_eff(c) * (qflx_dew_snow(c) + qflx_dew_grnd(c) ) & diff --git a/components/elm/src/biogeophys/BareGroundFluxesMod.F90 b/components/elm/src/biogeophys/BareGroundFluxesMod.F90 index 03ecff38f137..aff2447831b6 100644 --- a/components/elm/src/biogeophys/BareGroundFluxesMod.F90 +++ b/components/elm/src/biogeophys/BareGroundFluxesMod.F90 @@ -414,7 +414,7 @@ subroutine BareGroundFluxes(bounds, num_nolakeurbanp, filter_nolakeurbanp, & rh_ref2m(p) = min(100._r8, q_ref2m(p) / qsat_ref2m * 100._r8) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then rh_ref2m_r(p) = rh_ref2m(p) t_ref2m_r(p) = t_ref2m(p) end if diff --git a/components/elm/src/biogeophys/CanopyFluxesMod.F90 b/components/elm/src/biogeophys/CanopyFluxesMod.F90 old mode 100755 new mode 100644 diff --git a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 old mode 100755 new mode 100644 index 343380fa18fb..4a2cd4a5fd7a --- a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 +++ b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 @@ -278,8 +278,8 @@ subroutine CanopyHydrology(bounds, & ! Canopy interception and precipitation onto ground surface ! Add precipitation to leaf water - if (ltype(l)==istsoil .or. ltype(l)==istwet .or. urbpoi(l) .or. & - ltype(l)==istcrop) then + if (veg_pp%is_on_soil_col(p) .or. ltype(l) == istwet .or. urbpoi(l) .or. & + ltype(l) == istcrop) then qflx_candrip(p) = 0._r8 ! rate of canopy runoff qflx_through_snow(p) = 0._r8 ! snow precipitation direct through canopy @@ -334,7 +334,7 @@ subroutine CanopyHydrology(bounds, & end if end if - else if (ltype(l)==istice .or. ltype(l)==istice_mec) then + else if (ltype(l) == istice .or. ltype(l) == istice_mec) then h2ocan(p) = 0._r8 qflx_candrip(p) = 0._r8 @@ -642,7 +642,7 @@ subroutine CanopyHydrology(bounds, & end if !end of do_capsnow construct ! set frac_sno_eff variable - if (ltype(l) == istsoil .or. ltype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then if (subgridflag ==1) then frac_sno_eff(c) = frac_sno(c) else @@ -652,7 +652,7 @@ subroutine CanopyHydrology(bounds, & frac_sno_eff(c) = 1._r8 endif - if (ltype(l)==istwet .and. t_grnd(c)>tfrz) then + if (ltype(l) == istwet .and. t_grnd(c)>tfrz) then h2osno(c)=0._r8 snow_depth(c)=0._r8 end if @@ -826,7 +826,7 @@ subroutine FracH2OSfc(bounds, num_h2osfc, filter_h2osfc, & l = col_pp%landunit(c) qflx_h2osfc2topsoi(c) = 0._r8 ! h2osfc only calculated for soil vegetated land units - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then ! Use newton-raphson method to iteratively determine frac_h20sfc ! based on amount of surface water storage (h2osfc) and diff --git a/components/elm/src/biogeophys/CanopyStateType.F90 b/components/elm/src/biogeophys/CanopyStateType.F90 index 803366c05000..b07db4b7ad04 100644 --- a/components/elm/src/biogeophys/CanopyStateType.F90 +++ b/components/elm/src/biogeophys/CanopyStateType.F90 @@ -514,7 +514,7 @@ subroutine InitCold(this, bounds) this%dewmx_patch(p) = 0.1_r8 this%vegwp_patch(p,:) = -2.5e4_r8 - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then this%laisun_patch(p) = 0._r8 this%laisha_patch(p) = 0._r8 end if @@ -531,7 +531,7 @@ subroutine InitCold(this, bounds) do c = bounds%begc, bounds%endc l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%alt_col(c) = 0._r8 this%altmax_col(c) = 0._r8 this%altmax_lastyear_col(c) = 0._r8 diff --git a/components/elm/src/biogeophys/CanopyTemperatureMod.F90 b/components/elm/src/biogeophys/CanopyTemperatureMod.F90 index 017c3770a91d..63083feb42f9 100644 --- a/components/elm/src/biogeophys/CanopyTemperatureMod.F90 +++ b/components/elm/src/biogeophys/CanopyTemperatureMod.F90 @@ -244,7 +244,7 @@ subroutine CanopyTemperature(bounds, & if (lun_pp%itype(l)/=istwet .AND. lun_pp%itype(l)/=istice & .AND. lun_pp%itype(l)/=istice_mec) then - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then wx = (h2osoi_liq(c,1)/denh2o+h2osoi_ice(c,1)/denice)/dz(c,1) fac = min(1._r8, wx/watsat(c,1)) fac = max( fac, 0.01_r8 ) @@ -297,7 +297,7 @@ subroutine CanopyTemperature(bounds, & end if ! compute humidities individually for snow, soil, h2osfc for vegetated landunits - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then call QSat(t_soisno(c,snl(c)+1), forc_pbot(t), eg, degdT, qsatg, qsatgdT) if (qsatg > forc_q(t) .and. forc_q(t) > qsatg) then @@ -354,7 +354,7 @@ subroutine CanopyTemperature(bounds, & ! Urban emissivities are currently read in from data file if (.not. urbpoi(l)) then - if (lun_pp%itype(l)==istice .or. lun_pp%itype(l)==istice_mec) then + if (lun_pp%itype(l) == istice .or. lun_pp%itype(l) == istice_mec) then emg(c) = 0.97_r8 else emg(c) = (1._r8-frac_sno(c))*0.96_r8 + frac_sno(c)*0.97_r8 @@ -414,15 +414,16 @@ subroutine CanopyTemperature(bounds, & eflx_sh_tot(p) = 0._r8 l = veg_pp%landunit(p) + if (urbpoi(l)) then eflx_sh_tot_u(p) = 0._r8 - else if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + else if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then eflx_sh_tot_r(p) = 0._r8 end if eflx_lh_tot(p) = 0._r8 if (urbpoi(l)) then eflx_lh_tot_u(p) = 0._r8 - else if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + else if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then eflx_lh_tot_r(p) = 0._r8 end if eflx_sh_veg(p) = 0._r8 @@ -454,7 +455,7 @@ subroutine CanopyTemperature(bounds, & t = veg_pp%topounit(p) l = veg_pp%landunit(p) c = veg_pp%column(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then if (frac_veg_nosno(p) == 0) then forc_hgt_u_patch(p) = forc_hgt_u(t) + z0mg(c) + displa(p) forc_hgt_t_patch(p) = forc_hgt_t(t) + z0mg(c) + displa(p) diff --git a/components/elm/src/biogeophys/FrictionVelocityType.F90 b/components/elm/src/biogeophys/FrictionVelocityType.F90 index f00afe453f8c..73954cabfdb6 100644 --- a/components/elm/src/biogeophys/FrictionVelocityType.F90 +++ b/components/elm/src/biogeophys/FrictionVelocityType.F90 @@ -314,7 +314,7 @@ subroutine InitCold(this, bounds) type(bounds_type), intent(in) :: bounds ! ! !LOCAL VARIABLES: - integer :: p, c, l ! indices + integer :: p, c ! indices !----------------------------------------------------------------------- ! Added 5/4/04, PET: initialize forc_hgt_u (gridcell-level), @@ -331,8 +331,7 @@ subroutine InitCold(this, bounds) end if do c = bounds%begc, bounds%endc - l = col_pp%landunit(c) - if (lun_pp%lakpoi(l)) then !lake + if (col_pp%is_lake(c)) then !lake this%z0mg_col(c) = 0.0004_r8 end if end do diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 old mode 100755 new mode 100644 index e14c229160ce..eb8bb116f53e --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -222,7 +222,7 @@ subroutine HydrologyDrainage(bounds, & qflx_glcice_frz(c) = 0._r8 qflx_glcice_frz_diag(c) = 0._r8 - if (lun_pp%itype(l)==istice .and. qflx_snwcp_ice(c) > 0.0_r8) then + if (lun_pp%itype(l) == istice .and. qflx_snwcp_ice(c) > 0.0_r8) then qflx_glcice_frz_diag(c) = qflx_snwcp_ice(c) qflx_glcice_diag(c) = qflx_glcice_diag(c) + qflx_glcice_frz_diag(c) endif @@ -239,7 +239,7 @@ subroutine HydrologyDrainage(bounds, & if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 end if - !if (lun_pp%itype(l)==istice) then + !if (lun_pp%itype(l) == istice) then ! qflx_glcice_frz_diag(c) = qflx_snwcp_ice(c) ! qflx_glcice_diag(c) = qflx_glcice_diag(c) + qflx_glcice_frz_diag(c) !endif @@ -259,8 +259,8 @@ subroutine HydrologyDrainage(bounds, & tpu_ind = top_pp%topo_grc_ind(t) !Get topounit index on the grid g = col_pp%gridcell(c) - if (lun_pp%itype(l)==istwet .or. lun_pp%itype(l)==istice & - .or. lun_pp%itype(l)==istice_mec) then + if (lun_pp%itype(l) == istwet .or. lun_pp%itype(l) == istice & + .or. lun_pp%itype(l) == istice_mec) then qflx_drain(c) = 0._r8 qflx_drain_perched(c) = 0._r8 @@ -276,7 +276,7 @@ subroutine HydrologyDrainage(bounds, & ! glc_dyn_runoff_routing = true: in this case, melting ice runs off, and excess ! snow is sent to CISM, where it is converted to ice. These corrections are ! done here: - if (glc_dyn_runoff_routing(g) .and. lun_pp%itype(l)==istice_mec) then + if (glc_dyn_runoff_routing(g) .and. lun_pp%itype(l) == istice_mec) then ! this allows GLC melt to runoff to qflx_qrgwl! ! If glc_dyn_runoff_routing=T, add meltwater from istice_mec ice columns to the runoff. @@ -333,12 +333,12 @@ subroutine HydrologyDrainage(bounds, & qflx_runoff(c) = qflx_drain(c) + qflx_surf(c) + qflx_h2osfc_surf(c) + qflx_qrgwl(c) + qflx_drain_perched(c) - if ((lun_pp%itype(l)==istsoil .or. lun_pp%itype(l)==istcrop) .and. col_pp%active(c)) then + if ((col_pp%is_soil(c) .or. col_pp%is_crop(c)) .and. col_pp%active(c)) then qflx_irr_demand(c) = -1.0_r8 * f_surf(g,tpu_ind)*qflx_irrig(c) !surface water demand send to MOSART end if if (lun_pp%urbpoi(l)) then qflx_runoff_u(c) = qflx_runoff(c) - else if (lun_pp%itype(l)==istsoil .or. lun_pp%itype(l)==istcrop) then + else if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then qflx_runoff_r(c) = qflx_runoff(c) end if diff --git a/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 index fa6fe25c4ae7..33c27defeec2 100644 --- a/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 @@ -448,7 +448,7 @@ subroutine HydrologyNoDrainage(bounds, & t_soi_10cm(c) = t_soi_10cm(c)/0.1_r8 tsoi17(c) = tsoi17(c)/0.17_r8 ! F. Li and S. Levis end if - if (lun_pp%itype(l)==istsoil .or. lun_pp%itype(l)==istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then t_grnd_r(c) = t_soisno(c,snl(c)+1) end if diff --git a/components/elm/src/biogeophys/LakeStateType.F90 b/components/elm/src/biogeophys/LakeStateType.F90 index 405dd8807706..3ae7fbb44868 100644 --- a/components/elm/src/biogeophys/LakeStateType.F90 +++ b/components/elm/src/biogeophys/LakeStateType.F90 @@ -227,8 +227,7 @@ subroutine InitCold(this, bounds) !------------------------------------------------- do c = bounds%begc, bounds%endc - l = col_pp%landunit(c) - if (lun_pp%lakpoi(l)) then + if (col_pp%is_lake(c)) then ! Set lake ice fraction and top eddy conductivity from previous timestep ! Always initialize with no ice to prevent excessive ice sheets from forming when ! starting with old lake model that has unrealistically cold lake conseratures. diff --git a/components/elm/src/biogeophys/PhotosynthesisType.F90 b/components/elm/src/biogeophys/PhotosynthesisType.F90 index 5cb94628d8d2..c788726838de 100644 --- a/components/elm/src/biogeophys/PhotosynthesisType.F90 +++ b/components/elm/src/biogeophys/PhotosynthesisType.F90 @@ -10,6 +10,7 @@ module PhotosynthesisType use elm_varcon , only : spval use LandunitType , only : lun_pp use VegetationType , only : veg_pp + use ColumnDataType ,only : col_pp ! implicit none save @@ -420,12 +421,13 @@ subroutine TimeStepInit (this, bounds) type(bounds_type) , intent(in) :: bounds ! ! !LOCAL VARIABLES: - integer :: p,l ! indices + integer :: p,c,l ! indices !----------------------------------------------------------------------- do p = bounds%begp, bounds%endp l = veg_pp%landunit(p) - if (.not. lun_pp%lakpoi(l)) then + c = veg_pp%column(p) + if (.not. col_pp%is_lake(c)) then this%psnsun_patch(p) = 0._r8 this%psnsun_wc_patch(p) = 0._r8 this%psnsun_wj_patch(p) = 0._r8 @@ -454,7 +456,7 @@ subroutine TimeStepInit (this, bounds) this%c14_psnsha_patch(p) = 0._r8 endif end if - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop & + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p) & .or. lun_pp%itype(l) == istice .or. lun_pp%itype(l) == istice_mec & .or. lun_pp%itype(l) == istwet) then if (use_c13) then @@ -480,12 +482,13 @@ subroutine photosyns_vars_TimeStepInit (photosyns_vars, bounds) type(bounds_type) , intent(in) :: bounds ! ! !LOCAL VARIABLES: - integer :: p,l ! indices + integer :: p,c,l ! indices !----------------------------------------------------------------------- do p = bounds%begp, bounds%endp l = veg_pp%landunit(p) - if (.not. lun_pp%lakpoi(l)) then + c = veg_pp%column(p) + if (.not. col_pp%is_lake(c)) then photosyns_vars%psnsun_patch(p) = 0._r8 photosyns_vars%psnsun_wc_patch(p) = 0._r8 photosyns_vars%psnsun_wj_patch(p) = 0._r8 @@ -512,7 +515,7 @@ subroutine photosyns_vars_TimeStepInit (photosyns_vars, bounds) photosyns_vars%c14_psnsha_patch(p) = 0._r8 endif end if - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop & + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p) & .or. lun_pp%itype(l) == istice .or. lun_pp%itype(l) == istice_mec & .or. lun_pp%itype(l) == istwet) then if (use_c13) then diff --git a/components/elm/src/biogeophys/SedYieldMod.F90 b/components/elm/src/biogeophys/SedYieldMod.F90 index 250929f4a423..83c6ae220fa5 100644 --- a/components/elm/src/biogeophys/SedYieldMod.F90 +++ b/components/elm/src/biogeophys/SedYieldMod.F90 @@ -157,7 +157,7 @@ subroutine SoilErosion (bounds, num_soilc, filter_soilc, & flx_sed_yld(c) = 0._r8 ! check landunit type and ground covered by snow/ice - if ( lun_pp%itype(l)/=istsoil .and. lun_pp%itype(l)/=istcrop ) then + if ( .not.col_pp%is_soil(c) .and. .not.col_pp%is_crop(c) ) then cycle end if diff --git a/components/elm/src/biogeophys/SnowHydrologyMod.F90 b/components/elm/src/biogeophys/SnowHydrologyMod.F90 index 57e486d0ce90..8a0ba08b8d4b 100644 --- a/components/elm/src/biogeophys/SnowHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SnowHydrologyMod.F90 @@ -650,7 +650,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & if (.not. use_firn_percolation_and_compaction) then ! I don't think the next 5 lines are necessary (removed in CLMv5) l = col_pp%landunit(c) - if (ltype(l)==istice_mec .and. void < 0._r8) then + if (ltype(l) == istice_mec .and. void < 0._r8) then dz(c,j) = h2osoi_ice(c,j)/denice + h2osoi_liq(c,j)/denh2o void = 0._r8 endif @@ -705,7 +705,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & ! Compaction occurring during melt if (imelt(c,j) == 1) then - if(subgridflag==1 .and. (ltype(col_pp%landunit(c)) == istsoil .or. ltype(col_pp%landunit(c)) == istcrop)) then + if(subgridflag==1 .and. (col_pp%is_soil(c) .or. col_pp%is_crop(c))) then ! first term is delta mass over mass ddz3 = max(0._r8,min(1._r8,(swe_old(c,j) - wx)/wx)) @@ -858,7 +858,7 @@ subroutine CombineSnowLayers(bounds, num_snowc, filter_snowc, & if (num_snowc > 0) then c = filter_snowc(1) l = col_pp%landunit(c) - if (ltype(l) == istdlak) then ! Called from LakeHydrology + if (col_pp%is_lake(c)) then ! Called from LakeHydrology if (.not. use_extrasnowlayers) then dzminloc(:) = dzmin(:) + lsadz else @@ -884,7 +884,7 @@ subroutine CombineSnowLayers(bounds, num_snowc, filter_snowc, & do j = msn_old(c)+1,0 ! use 0.01 to avoid runaway ice buildup if (h2osoi_ice(c,j) <= .01_r8) then - if (ltype(l) == istsoil .or. urbpoi(l) .or. ltype(l) == istcrop) then + if (col_pp%is_soil(c) .or. urbpoi(l) .or. col_pp%is_crop(c)) then h2osoi_liq(c,j+1) = h2osoi_liq(c,j+1) + h2osoi_liq(c,j) h2osoi_ice(c,j+1) = h2osoi_ice(c,j+1) + h2osoi_ice(c,j) @@ -911,7 +911,7 @@ subroutine CombineSnowLayers(bounds, num_snowc, filter_snowc, & mss_dst4(c,j+1) = mss_dst4(c,j+1) + mss_dst4(c,j) end if - else if (ltype(l) /= istsoil .and. .not. urbpoi(l) .and. ltype(l) /= istcrop .and. j /= 0) then + else if (.not.col_pp%is_soil(c) .and. .not. urbpoi(l) .and. .not.col_pp%is_crop(c) .and. j /= 0) then h2osoi_liq(c,j+1) = h2osoi_liq(c,j+1) + h2osoi_liq(c,j) h2osoi_ice(c,j+1) = h2osoi_ice(c,j+1) + h2osoi_ice(c,j) @@ -935,7 +935,7 @@ subroutine CombineSnowLayers(bounds, num_snowc, filter_snowc, & ! urban, soil or crop, the h2osoi_liq and h2osoi_ice associated with this layer is sent ! to qflx_qrgwl later on in the code. To keep track of this for the snow balance ! error check, we add this to qflx_sl_top_soil here - if (ltype(l) /= istsoil .and. ltype(l) /= istcrop .and. .not. urbpoi(l) .and. i == 0) then + if (.not.col_pp%is_soil(c) .and. .not.col_pp%is_crop(c) .and. .not. urbpoi(l) .and. i == 0) then qflx_sl_top_soil(c) = (h2osoi_liq(c,i) + h2osoi_ice(c,i))/dtime end if @@ -987,7 +987,7 @@ subroutine CombineSnowLayers(bounds, num_snowc, filter_snowc, & c = filter_snowc(fc) l = col_pp%landunit(c) if (snow_depth(c) > 0._r8) then - if ((ltype(l) == istdlak .and. snow_depth(c) < 0.01_r8 + lsadz ) .or. & + if ((col_pp%is_lake(c) .and. snow_depth(c) < 0.01_r8 + lsadz ) .or. & ((ltype(l) /= istdlak) .and. ((frac_sno_eff(c)*snow_depth(c) < 0.01_r8) & .or. (h2osno(c)/(frac_sno_eff(c)*snow_depth(c)) < 50._r8)))) then @@ -1005,7 +1005,7 @@ subroutine CombineSnowLayers(bounds, num_snowc, filter_snowc, & if (h2osno(c) <= 0._r8) snow_depth(c) = 0._r8 ! this is where water is transfered from layer 0 (snow) to layer 1 (soil) - if (ltype(l) == istsoil .or. urbpoi(l) .or. ltype(l) == istcrop) then + if (col_pp%is_soil(c) .or. urbpoi(l) .or. col_pp%is_crop(c)) then h2osoi_liq(c,0) = 0.0_r8 h2osoi_liq(c,1) = h2osoi_liq(c,1) + zwliq(c) qflx_snow2topsoi(c) = zwliq(c)/dtime @@ -1014,7 +1014,7 @@ subroutine CombineSnowLayers(bounds, num_snowc, filter_snowc, & if (ltype(l) == istwet) then h2osoi_liq(c,0) = 0.0_r8 endif - if (ltype(l) == istice .or. ltype(l)==istice_mec) then + if (ltype(l) == istice .or. ltype(l) == istice_mec) then h2osoi_liq(c,0) = 0.0_r8 endif endif @@ -2148,7 +2148,7 @@ subroutine InitSnowLayers (bounds, snow_depth) real(r8) , intent(in) :: snow_depth(bounds%begc:) ! ! LOCAL VARAIBLES: - integer :: c,l,j ! indices + integer :: c,j ! indices real(r8) :: minbound, maxbound ! helper variables !SHR_ASSERT_ALL((ubound(snow_depth) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) @@ -2160,14 +2160,13 @@ subroutine InitSnowLayers (bounds, snow_depth) zi => col_pp%zi & ! Output: [real(r8) (:,:) ] interface level below a "z" level (m) (-nlevsno+0:nlevgrnd) ) loop_columns: do c = bounds%begc,bounds%endc - l = col_pp%landunit(c) dz(c,-nlevsno+1: 0) = spval z (c,-nlevsno+1: 0) = spval zi(c,-nlevsno :-1) = spval ! Special case: lake - if (lun_pp%lakpoi(l)) then + if (col_pp%is_lake(c)) then snl(c) = 0 dz(c,-nlevsno+1:0) = 0._r8 z(c,-nlevsno+1:0) = 0._r8 diff --git a/components/elm/src/biogeophys/SoilFluxesMod.F90 b/components/elm/src/biogeophys/SoilFluxesMod.F90 index 1b221a617a83..f18ead86f1e2 100644 --- a/components/elm/src/biogeophys/SoilFluxesMod.F90 +++ b/components/elm/src/biogeophys/SoilFluxesMod.F90 @@ -288,7 +288,7 @@ subroutine SoilFluxes (bounds, num_urbanl, filter_urbanl, & - emg(c)*sb*lw_grnd - emg(c)*sb*t_grnd0(c)**3*(4._r8*tinc(c)) & - (eflx_sh_grnd(p)+qflx_evap_soi(p)*htvp(c)) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then eflx_soil_grnd_r(p) = eflx_soil_grnd(p) end if else @@ -312,7 +312,7 @@ subroutine SoilFluxes (bounds, num_urbanl, filter_urbanl, & eflx_sh_tot(p) = eflx_sh_veg(p) + eflx_sh_grnd(p) qflx_evap_tot(p) = qflx_evap_veg(p) + qflx_evap_soi(p) eflx_lh_tot(p)= hvap*qflx_evap_veg(p) + htvp(c)*qflx_evap_soi(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then eflx_lh_tot_r(p)= eflx_lh_tot(p) eflx_sh_tot_r(p)= eflx_sh_tot(p) else if (lun_pp%urbpoi(l)) then @@ -377,7 +377,7 @@ subroutine SoilFluxes (bounds, num_urbanl, filter_urbanl, & errsoi_patch(p) = errsoi_patch(p)+eflx_h2osfc_to_snow_col(c) ! For urban sunwall, shadewall, and roof columns, the "soil" energy balance check ! must include the heat flux from the interior of the building. - if (col_pp%itype(c)==icol_sunwall .or. col_pp%itype(c)==icol_shadewall .or. col_pp%itype(c)==icol_roof) then + if (col_pp%itype(c) == icol_sunwall .or. col_pp%itype(c) == icol_shadewall .or. col_pp%itype(c) == icol_roof) then errsoi_patch(p) = errsoi_patch(p) + eflx_building_heat(c) end if end do @@ -432,7 +432,7 @@ subroutine SoilFluxes (bounds, num_urbanl, filter_urbanl, & + 4._r8*emg(c)*sb*t_grnd0(c)**3*tinc(c) eflx_lwrad_net(p) = eflx_lwrad_out(p) - forc_lwrad(t) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then eflx_lwrad_net_r(p) = eflx_lwrad_out(p) - forc_lwrad(t) eflx_lwrad_out_r(p) = eflx_lwrad_out(p) end if diff --git a/components/elm/src/biogeophys/SoilHydrologyMod.F90 b/components/elm/src/biogeophys/SoilHydrologyMod.F90 index e226adf96ebe..b43eafe7d39e 100644 --- a/components/elm/src/biogeophys/SoilHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SoilHydrologyMod.F90 @@ -440,10 +440,11 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f do fc = 1, num_hydrologyc c = filter_hydrologyc(fc) g = cgridcell(c) + l = col_pp%landunit(c) pc = pc_grid(g) ! partition moisture fluxes between soil and h2osfc - if (lun_pp%itype(col_pp%landunit(c)) == istsoil .or. lun_pp%itype(col_pp%landunit(c))==istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then ! explicitly use frac_sno=0 if snl=0 if (snl(c) >= 0) then @@ -1068,7 +1069,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte ! ! !LOCAL VARIABLES: !character(len=32) :: subname = 'Drainage' ! subroutine name - integer :: c,j,fc,i ! indices + integer :: c,l,j,fc,i ! indices integer :: nlevbed ! # layers to bedrock real(r8) :: xs(bounds%begc:bounds%endc) ! water needed to bring soil moisture to watmin (mm) real(r8) :: dzmm(bounds%begc:bounds%endc,1:nlevgrnd) ! layer thickness (mm) @@ -1571,6 +1572,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte do fc = 1, num_hydrologyc c = filter_hydrologyc(fc) + l = col_pp%landunit(c) !scs: watmin addition to fix water balance errors xs1(c) = max(max(h2osoi_liq(c,1)-watmin,0._r8)- & @@ -1578,7 +1580,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte if (use_vsfm) xs1(c) = 0._r8 h2osoi_liq(c,1) = h2osoi_liq(c,1) - xs1(c) - if (lun_pp%urbpoi(col_pp%landunit(c))) then + if (lun_pp%urbpoi(l)) then qflx_rsub_sat(c) = xs1(c) / dtime else if(h2osfcflag == 1) then @@ -1596,7 +1598,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte ! add in ice check xs1(c) = max(max(h2osoi_ice(c,1),0._r8)-max(0._r8,(pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1))),0._r8) h2osoi_ice(c,1) = min(max(0._r8,pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1)), h2osoi_ice(c,1)) - if ( (lun_pp%itype(col_pp%landunit(c)) == istice .or. lun_pp%itype(col_pp%landunit(c)) == istice_mec) .or. (.not. use_firn_percolation_and_compaction)) then + if ( (lun_pp%itype(l) == istice .or. lun_pp%itype(l) == istice_mec) .or. (.not. use_firn_percolation_and_compaction)) then qflx_snwcp_ice(c) = qflx_snwcp_ice(c) + xs1(c) / dtime else qflx_ice_runoff_xs(c) = qflx_ice_runoff_xs(c) + xs1(c) / dtime diff --git a/components/elm/src/biogeophys/SoilHydrologyType.F90 b/components/elm/src/biogeophys/SoilHydrologyType.F90 index 203e0a7d52b1..398ab68005fc 100644 --- a/components/elm/src/biogeophys/SoilHydrologyType.F90 +++ b/components/elm/src/biogeophys/SoilHydrologyType.F90 @@ -289,7 +289,7 @@ subroutine InitCold(this, bounds) do c = bounds%begc,bounds%endc l = col_pp%landunit(c) nlevbed = col_pp%nlevbed(c) - if (.not. lun_pp%lakpoi(l)) then !not lake + if (.not. col_pp%is_lake(c)) then !not lake if (lun_pp%urbpoi(l)) then if (col_pp%itype(c) == icol_road_perv) then this%wa_col(c) = 0._r8 @@ -314,7 +314,7 @@ subroutine InitCold(this, bounds) do c = bounds%begc,bounds%endc l = col_pp%landunit(c) nlevbed = col_pp%nlevbed(c) - if (.not. lun_pp%lakpoi(l)) then !not lake + if (.not. col_pp%is_lake(c)) then !not lake if (lun_pp%urbpoi(l)) then if (col_pp%itype(c) == icol_road_perv) then this%wa_col(c) = 4800._r8 @@ -462,7 +462,7 @@ subroutine InitCold(this, bounds) ti = t - topi + 1 if (lun_pp%itype(l) /= istdlak) then ! soil columns of both urban and non-urban types - if (lun_pp%itype(l)==istwet .or. lun_pp%itype(l)==istice .or. lun_pp%itype(l)==istice_mec) then + if (lun_pp%itype(l) == istwet .or. lun_pp%itype(l) == istice .or. lun_pp%itype(l) == istice_mec) then ! do nothing else if (lun_pp%urbpoi(l) .and. (col_pp%itype(c) /= icol_road_perv) .and. (col_pp%itype(c) /= icol_road_imperv) )then ! do nothing @@ -510,7 +510,7 @@ subroutine InitCold(this, bounds) if (lun_pp%itype(l) /= istdlak) then ! soil columns of both urban and non-urban types if (lun_pp%urbpoi(l)) then - if (col_pp%itype(c)==icol_sunwall .or. col_pp%itype(c)==icol_shadewall .or. col_pp%itype(c)==icol_roof) then + if (col_pp%itype(c) == icol_sunwall .or. col_pp%itype(c) == icol_shadewall .or. col_pp%itype(c) == icol_roof) then ! do nothing else this%depth_col(c, 1:nlayer) = dzvic diff --git a/components/elm/src/biogeophys/SoilStateType.F90 b/components/elm/src/biogeophys/SoilStateType.F90 index 24f7343db813..bea48018a550 100644 --- a/components/elm/src/biogeophys/SoilStateType.F90 +++ b/components/elm/src/biogeophys/SoilStateType.F90 @@ -407,9 +407,9 @@ subroutine InitCold(this, bounds) do c = bounds%begc,bounds%endc this%rootfr_col (c,nlevsoi+1:nlevgrnd) = 0._r8 - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%rootfr_col (c,nlevsoi+1:nlevgrnd) = 0._r8 - else if (lun_pp%itype(l) == istdlak .and. allowlakeprod) then + else if (col_pp%is_lake(c) .and. allowlakeprod) then this%rootfr_col (c,:) = spval else ! Inactive CH4 columns this%rootfr_col (c,:) = spval @@ -588,7 +588,7 @@ subroutine InitCold(this, bounds) topi = grc_pp%topi(g) ti = t - topi + 1 - if (lun_pp%itype(l)==istwet .or. lun_pp%itype(l)==istice .or. lun_pp%itype(l)==istice_mec) then + if (lun_pp%itype(l) == istwet .or. lun_pp%itype(l) == istice .or. lun_pp%itype(l) == istice_mec) then do lev = 1,nlevgrnd this%bsw_col(c,lev) = spval @@ -613,7 +613,7 @@ subroutine InitCold(this, bounds) this%tkmg_col(c,lev) = spval this%tksatu_col(c,lev) = spval this%tkdry_col(c,lev) = spval - if (lun_pp%itype(l)==istwet .and. lev > nlevbed) then + if (lun_pp%itype(l) == istwet .and. lev > nlevbed) then this%csol_col(c,lev) = csol_bedrock else this%csol_col(c,lev)= spval @@ -690,7 +690,7 @@ subroutine InitCold(this, bounds) endif end if - if (lun_pp%itype(l) == istdlak) then + if (col_pp%is_lake(c)) then if (lev <= nlevsoi) then this%cellsand_col(c,lev) = sand @@ -699,7 +699,7 @@ subroutine InitCold(this, bounds) this%cellorg_col(c,lev) = om_frac*organic_max end if - else if (lun_pp%itype(l) /= istdlak) then ! soil columns of both urban and non-urban types + else ! soil columns of both urban and non-urban types if (lun_pp%urbpoi(l)) then om_frac = 0._r8 ! No organic matter for urban @@ -810,7 +810,7 @@ subroutine InitCold(this, bounds) g = col_pp%gridcell(c) l = col_pp%landunit(c) - if (lun_pp%itype(l)==istdlak) then + if (col_pp%is_lake(c)) then do lev = 1,nlevgrnd if ( lev <= nlevsoi )then diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index 6826b3e22e63..37878d4e77ab 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -677,9 +677,9 @@ subroutine SoilTemperature(bounds, num_urbanl, filter_urbanl, num_nolakec, filte if (j == 1) then ! this only needs to be done once eflx_fgr12(c) = -cnfac*fn(c,1) - (1._r8-cnfac)*fn1(c,1) end if - if (j > 0 .and. j < nlevgrnd .and. (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop)) then + if (j > 0 .and. j < nlevgrnd .and. (col_pp%is_soil(c) .or. col_pp%is_crop(c))) then eflx_fgr(c,j) = -cnfac*fn(c,j) - (1._r8-cnfac)*fn1(c,j) - else if (j == nlevgrnd .and. (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop)) then + else if (j == nlevgrnd .and. (col_pp%is_soil(c) .or. col_pp%is_crop(c))) then eflx_fgr(c,j) = 0._r8 end if @@ -1465,7 +1465,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & ! from Zhao (1997) and Koren (1999) supercool(c,j) = 0.0_r8 - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop .or. col_pp%itype(c) == icol_road_perv) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c) .or. col_pp%itype(c) == icol_road_perv) then if(t_soisno(c,j) < tfrz) then smp = hfus*(tfrz-t_soisno(c,j))/(grav*t_soisno(c,j)) * 1000._r8 !(mm) supercool(c,j) = watsat(c,j)*(smp/sucsat(c,j))**(-1._r8/bsw(c,j)) @@ -1648,7 +1648,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & ! as computed in HydrologyDrainageMod.F90. l = col_pp%landunit(c) - if ( lun_pp%itype(l)==istice_mec) then + if ( lun_pp%itype(l) == istice_mec) then if (j>=1 .and. h2osoi_liq(c,j) > 0._r8) then ! ice layer with meltwater ! melting corresponds to a negative ice flux qflx_glcice_melt(c) = qflx_glcice_melt(c) + h2osoi_liq(c,j)/dtime @@ -1662,7 +1662,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & endif ! istice_mec ! for diagnostic QICE SMB output only - ! these are to calculate SMB even without MECs - if ( lun_pp%itype(l)==istice) then + if ( lun_pp%itype(l) == istice) then if (j>=1 .and. h2osoi_liq(c,j) > 0._r8) then ! ice layer with meltwater ! melting corresponds to a negative ice flux qflx_glcice_melt_diag(c) = qflx_glcice_melt_diag(c) + h2osoi_liq(c,j)/dtime @@ -1682,7 +1682,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & l = col_pp%landunit(c) if (lun_pp%urbpoi(l)) then eflx_snomelt_u(c) = eflx_snomelt(c) - else if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + else if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then eflx_snomelt_r(c) = eflx_snomelt(c) end if end do diff --git a/components/elm/src/biogeophys/SurfaceAlbedoMod.F90 b/components/elm/src/biogeophys/SurfaceAlbedoMod.F90 index 4aed983b0ae0..9bcf1e4df24e 100644 --- a/components/elm/src/biogeophys/SurfaceAlbedoMod.F90 +++ b/components/elm/src/biogeophys/SurfaceAlbedoMod.F90 @@ -116,7 +116,7 @@ subroutine SurfaceAlbedo(bounds, & integer :: i ! index for layers [idx] integer :: aer ! index for sno_nbr_aer real(r8) :: extkn ! nitrogen allocation coefficient - integer :: fp,fc,g,c,p,iv ! indices + integer :: fp,fc,g,c,l,p,iv ! indices integer :: ib ! band index integer :: ic ! 0=unit incoming direct; 1=unit incoming diffuse real(r8) :: dinc ! lai+sai increment for canopy layer @@ -652,6 +652,7 @@ subroutine SurfaceAlbedo(bounds, & do ib = 1, nband do fc = 1,num_nourbanc c = filter_nourbanc(fc) + l = col_pp%landunit(c) if (coszen_col(c) > 0._r8) then ! ground albedo was originally computed in SoilAlbedo, but is now computed here ! because the order of SoilAlbedo and SNICAR_RT/SNICAR_AD_RT was switched for SNICAR/SNICAR_AD_RT. @@ -683,7 +684,7 @@ subroutine SurfaceAlbedo(bounds, & ! weight snow layer radiative absorption factors based on snow fraction and soil albedo ! (NEEDED FOR ENERGY CONSERVATION) do i = -nlevsno+1,1,1 - if (subgridflag == 0 .or. lun_pp%itype(col_pp%landunit(c)) == istdlak) then + if (subgridflag == 0 .or. col_pp%is_lake(c)) then if (ib == 1) then flx_absdv(c,i) = flx_absd_snw(c,i,ib)*frac_sno(c) + & ((1.-frac_sno(c))*(1-albsod(c,ib))*(flx_absd_snw(c,i,ib)/(1.-albsnd(c,ib)))) @@ -733,8 +734,8 @@ subroutine SurfaceAlbedo(bounds, & do fp = 1,num_nourbanp p = filter_nourbanp(fp) if (coszen_patch(p) > 0._r8) then - if ((lun_pp%itype(veg_pp%landunit(p)) == istsoil .or. & - lun_pp%itype(veg_pp%landunit(p)) == istcrop ) & + if ((veg_pp%is_on_soil_col(p) .or. & + veg_pp%is_on_crop_col(p) ) & .and. (elai(p) + esai(p)) > 0._r8) then num_vegsol = num_vegsol + 1 filter_vegsol(num_vegsol) = p @@ -1051,7 +1052,7 @@ subroutine SoilAlbedo (bounds, & pi = ldomain%pftm(g) if (coszen(c) > 0._r8) then l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop .and. top_pp%active(t)) then ! soil + if (col_pp%is_soil(c) .or. col_pp%is_crop(c) .and. top_pp%active(t)) then ! soil inc = max(0.11_r8-0.40_r8*h2osoi_vol(c,1), 0._r8) soilcol = isoicol(c) ! changed from local variable to elm_type: @@ -1066,7 +1067,7 @@ subroutine SoilAlbedo (bounds, & albsod(c,ib) = albice(ib) albsoi(c,ib) = albsod(c,ib) ! unfrozen lake, wetland - else if (t_grnd(c) > tfrz .or. (lakepuddling .and. lun_pp%itype(l) == istdlak .and. t_grnd(c) == tfrz .and. & + else if (t_grnd(c) > tfrz .or. (lakepuddling .and. col_pp%is_lake(c) .and. t_grnd(c) == tfrz .and. & lake_icefrac(c,1) < 1._r8 .and. lake_icefrac(c,2) > 0._r8) .and. top_pp%active(t)) then albsod(c,ib) = 0.05_r8/(max(0.001_r8,coszen(c)) + 0.15_r8) @@ -1078,7 +1079,7 @@ subroutine SoilAlbedo (bounds, & ! ZMS: Attn EK, currently restoring this for wetlands even though it is wrong in order to try to get ! bfb baseline comparison when no lakes are present. I'm assuming wetlands will be phased out anyway. - if (lun_pp%itype(l) == istdlak) then + if (col_pp%is_lake(c)) then albsoi(c,ib) = 0.10_r8 else albsoi(c,ib) = albsod(c,ib) @@ -1090,7 +1091,7 @@ subroutine SoilAlbedo (bounds, & ! Tenatively I'm restricting this to lakes because I haven't tested it for wetlands. But if anything ! the albedo should be lower when melting over frozen ground than a solid frozen lake. ! - if (lun_pp%itype(l) == istdlak .and. .not. lakepuddling .and. snl(c) == 0 .and. top_pp%active(t)) then + if (col_pp%is_lake(c) .and. .not. lakepuddling .and. snl(c) == 0 .and. top_pp%active(t)) then ! Need to reference snow layers here because t_grnd could be over snow or ice ! but we really want the ice surface temperature with no snow sicefr = 1._r8 - exp(-calb * (tfrz - t_grnd(c))/tfrz) diff --git a/components/elm/src/biogeophys/SurfaceRadiationMod.F90 b/components/elm/src/biogeophys/SurfaceRadiationMod.F90 index 47fe244c530e..d4a25c54192c 100644 --- a/components/elm/src/biogeophys/SurfaceRadiationMod.F90 +++ b/components/elm/src/biogeophys/SurfaceRadiationMod.F90 @@ -484,7 +484,7 @@ subroutine SurfaceRadiation(bounds, num_nourbanp, filter_nourbanp, & sabg(p) = 0._r8 sabv(p) = 0._r8 fsa(p) = 0._r8 - if (lun_pp%itype(l)==istsoil .or. lun_pp%itype(l)==istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then fsa_r(p) = 0._r8 end if sabg_lyr(p,:) = 0._r8 @@ -521,7 +521,7 @@ subroutine SurfaceRadiation(bounds, num_nourbanp, filter_nourbanp, & if (ib == 1) then parveg(p) = cad(p,ib) + cai(p,ib) end if - if (lun_pp%itype(l)==istsoil .or. lun_pp%itype(l)==istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then fsa_r(p) = fsa_r(p) + cad(p,ib) + cai(p,ib) end if @@ -538,7 +538,7 @@ subroutine SurfaceRadiation(bounds, num_nourbanp, filter_nourbanp, & absrad = trd(p,ib)*(1._r8-albgrd(c,ib)) + tri(p,ib)*(1._r8-albgri(c,ib)) sabg(p) = sabg(p) + absrad fsa(p) = fsa(p) + absrad - if (lun_pp%itype(l)==istsoil .or. lun_pp%itype(l)==istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then fsa_r(p) = fsa_r(p) + absrad end if if (snl(c) == 0) then @@ -637,7 +637,7 @@ subroutine SurfaceRadiation(bounds, num_nourbanp, filter_nourbanp, & ! If shallow snow depth, all solar radiation absorbed in top or top two snow layers ! to prevent unrealistic timestep soil warming - if (subgridflag == 0 .or. lun_pp%itype(l) == istdlak) then + if (subgridflag == 0 .or. col_pp%is_lake(c)) then if (snow_depth(c) < 0.10_r8) then if (snl(c) == 0) then sabg_lyr(p,-nlevsno+1:0) = 0._r8 @@ -678,7 +678,7 @@ subroutine SurfaceRadiation(bounds, num_nourbanp, filter_nourbanp, & endif ! Diagnostic: shortwave penetrating ground (e.g. top layer) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then sabg_pen(p) = sabg(p) - sabg_lyr(p, snl(c)+1) end if diff --git a/components/elm/src/biogeophys/SurfaceResistanceMod.F90 b/components/elm/src/biogeophys/SurfaceResistanceMod.F90 index b4217890fca3..27b9dccd90fa 100644 --- a/components/elm/src/biogeophys/SurfaceResistanceMod.F90 +++ b/components/elm/src/biogeophys/SurfaceResistanceMod.F90 @@ -141,7 +141,7 @@ subroutine calc_beta_leepielke1992(bounds, num_nolakec, filter_nolakec, & l = col_pp%landunit(c) if (lun_pp%itype(l)/=istwet .AND. lun_pp%itype(l)/=istice & .AND. lun_pp%itype(l)/=istice_mec) then - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then wx = (h2osoi_liq(c,1)/denh2o+h2osoi_ice(c,1)/denice)/col_pp%dz(c,1) fac = min(1._r8, wx/watsat(c,1)) fac = max( fac, 0.01_r8 ) diff --git a/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 b/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 index f222dca037e2..5287b74c13c6 100644 --- a/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 +++ b/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 @@ -267,7 +267,7 @@ subroutine ComputeLiqIceMassNonLake(bounds, num_nolakec, filter_nolakec, & do fc = 1, num_nolakec c = filter_nolakec(fc) l = col_pp%landunit(c) - if ( (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop ) & + if ( (col_pp%is_soil(c) .or. col_pp%is_crop(c) ) & .or. (lun_pp%itype(l) == istwet ) & .or. (lun_pp%itype(l) == istice ) & .or. (lun_pp%itype(l) == istice_mec ) & @@ -276,7 +276,7 @@ subroutine ComputeLiqIceMassNonLake(bounds, num_nolakec, filter_nolakec, & end if l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then ! note: soil specified at LU level + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then ! note: soil specified at LU level do p = col_pp%pfti(c),col_pp%pftf(c) ! loop over patches if (veg_pp%active(p)) then liquid_mass(c) = liquid_mass(c) + h2ocan_patch(p) * veg_pp%wtcol(p) @@ -564,7 +564,7 @@ subroutine ComputeHeatNonLake(bounds, num_nolakec, filter_nolakec, & c = filter_nolakec(fc) l = col_pp%landunit(c) - if (col_pp%itype(c)==icol_sunwall .or. col_pp%itype(c)==icol_shadewall) then + if (col_pp%itype(c) == icol_sunwall .or. col_pp%itype(c) == icol_shadewall) then has_h2o = .false. if (j <= nlevurb) then heat_dry_mass(c) = heat_dry_mass(c) + & diff --git a/components/elm/src/biogeophys/WaterfluxType.F90 b/components/elm/src/biogeophys/WaterfluxType.F90 index 8dfb1fb00c80..d165c0750c61 100644 --- a/components/elm/src/biogeophys/WaterfluxType.F90 +++ b/components/elm/src/biogeophys/WaterfluxType.F90 @@ -373,7 +373,7 @@ subroutine InitCold(this, bounds) ! needed for NitrogenLeaching do c = bounds%begc, bounds%endc l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%qflx_drain_col(c) = 0._r8 this%qflx_surf_col(c) = 0._r8 this%qflx_irr_demand_col(c) = 0._r8 @@ -383,7 +383,7 @@ subroutine InitCold(this, bounds) do p = bounds%begp, bounds%endp l = veg_pp%landunit(p) - if (lun_pp%itype(l)==istsoil) then + if (veg_pp%is_on_soil_col(p)) then this%n_irrig_steps_left_patch(p) = 0 this%irrig_rate_patch(p) = 0.0_r8 end if diff --git a/components/elm/src/cpl/lnd_disagg_forc.F90 b/components/elm/src/cpl/lnd_disagg_forc.F90 old mode 100755 new mode 100644 diff --git a/components/elm/src/data_types/CNStateType.F90 b/components/elm/src/data_types/CNStateType.F90 index 2f880568f30d..660f3eac1e40 100644 --- a/components/elm/src/data_types/CNStateType.F90 +++ b/components/elm/src/data_types/CNStateType.F90 @@ -1020,7 +1020,7 @@ subroutine initCold(this, bounds) end do end if - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%annsum_counter_col(c) = 0._r8 this%annavg_t2m_col(c) = 280._r8 @@ -1098,7 +1098,7 @@ subroutine initCold(this, bounds) do p = bounds%begp,bounds%endp l = veg_pp%landunit(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then this%rc14_atm_patch(p) = c14ratio diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index 3f886d04c7de..dc35d59ffe51 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -1234,9 +1234,9 @@ subroutine col_es_init(this, begc, endc) end if ! Below snow temperatures - nonlake points (lake points are set below) - if (.not. lun_pp%lakpoi(l)) then + if (.not. col_pp%is_lake(c)) then - if (lun_pp%itype(l)==istice .or. lun_pp%itype(l)==istice_mec) then + if (lun_pp%itype(l) == istice .or. lun_pp%itype(l) == istice_mec) then this%t_soisno(c,1:nlevgrnd) = 250._r8 else if (lun_pp%itype(l) == istwet) then @@ -1292,7 +1292,7 @@ subroutine col_es_init(this, begc, endc) this%t_grnd(c) = this%t_soisno(c,snl(c)+1) endif - if (lun_pp%lakpoi(l)) then ! special handling for lake points + if (col_pp%is_lake(c)) then ! special handling for lake points this%t_soisno(c,1:nlevgrnd) = 277._r8 this%t_lake(c,1:nlevlak) = 277._r8 this%t_grnd(c) = 277._r8 @@ -1733,10 +1733,10 @@ subroutine col_ws_init(this, begc, endc, h2osno_input, snow_depth_input, watsat_ this%h2osoi_liq(c,-nlevsno+1:) = spval this%h2osoi_ice(c,-nlevsno+1:) = spval - if (.not. lun_pp%lakpoi(l)) then !not lake + if (.not. col_pp%is_lake(c)) then !not lake nlevbed = col_pp%nlevbed(c) ! volumetric water - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then nlevs = nlevgrnd do j = 1, nlevs if (j > nlevbed) then @@ -1823,7 +1823,7 @@ subroutine col_ws_init(this, begc, endc, h2osno_input, snow_depth_input, watsat_ !-------------------------------------------- ! Set Lake water !-------------------------------------------- - if (lun_pp%lakpoi(l)) then + if (col_pp%is_lake(c)) then do j = -nlevsno+1, 0 if (j > col_pp%snl(c)) then this%h2osoi_ice(c,j) = col_pp%dz(c,j)*bdsno @@ -2071,7 +2071,7 @@ subroutine col_ws_restart(this, bounds, ncid, flag, watsat_input) end if do j = 1,nlevs l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%h2osoi_liq(c,j) = max(0._r8,this%h2osoi_liq(c,j)) this%h2osoi_ice(c,j) = max(0._r8,this%h2osoi_ice(c,j)) this%h2osoi_vol(c,j) = this%h2osoi_liq(c,j)/(col_pp%dz(c,j)*denh2o) & @@ -2482,7 +2482,7 @@ subroutine col_cs_init(this, begc, endc, carbon_type, ratio, c12_carbonstate_var do c = begc, endc l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then if (.not. present(c12_carbonstate_vars)) then ! initializing a c12 type do j = 1, nlevdecomp do k = 1, ndecomp_pools @@ -3699,7 +3699,7 @@ subroutine col_ns_init(this, begc, endc, col_cs) do c = begc, endc l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then ! column nitrogen state variables this%ntrunc(c) = 0._r8 @@ -4850,7 +4850,7 @@ subroutine col_ps_init(this, begc, endc, col_cs) do c = begc, endc l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then ! column phosphorus state variables this%ptrunc(c) = 0._r8 @@ -6032,7 +6032,7 @@ subroutine col_wf_init(this, begc, endc) ! needed for CNNLeaching do c = begc, endc l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%qflx_drain(c) = 0._r8 this%qflx_surf(c) = 0._r8 end if @@ -7189,9 +7189,9 @@ subroutine col_cf_init(this, begc, endc, carbon_type) end if this%fphr(c,nlevdecomp+1:nlevgrnd) = 0._r8 !used to be in ch4Mod - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%fphr(c,nlevdecomp+1:nlevgrnd) = 0._r8 - else if (lun_pp%itype(l) == istdlak .and. allowlakeprod) then + else if (col_pp%is_lake(c) .and. allowlakeprod) then this%fphr(c,:) = spval else ! Inactive CH4 columns this%fphr(c,:) = spval @@ -7199,7 +7199,7 @@ subroutine col_cf_init(this, begc, endc, carbon_type) ! also initialize dynamic landcover fluxes so that they have ! real values on first timestep, prior to calling pftdyn_cnbal - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then this%dwt_conv_cflux(c) = 0._r8 this%dwt_prod10c_gain(c) = 0._r8 this%dwt_prod100c_gain(c) = 0._r8 diff --git a/components/elm/src/data_types/ColumnType.F90 b/components/elm/src/data_types/ColumnType.F90 index 21e9d24bd8f4..8edc53511d3e 100644 --- a/components/elm/src/data_types/ColumnType.F90 +++ b/components/elm/src/data_types/ColumnType.F90 @@ -76,6 +76,11 @@ module ColumnType logical, pointer :: is_fates(:) => null() ! True if this column is associated with a FATES active column ! False if otherwise. If fates is turned off, this array is ! all false + + logical, pointer :: is_soil(:) => null() ! True if the column is a soil + logical, pointer :: is_crop(:) => null() ! True if the column is a crop + logical, pointer :: is_lake(:) => null() ! True if the column is a lake + contains procedure, public :: Init => col_pp_init @@ -139,6 +144,10 @@ subroutine col_pp_init(this, begc, endc) ! Assume that columns are not fates columns until fates initialization begins allocate(this%is_fates(begc:endc)); this%is_fates(:) = .false. + allocate(this%is_soil (begc:endc)); this%is_soil (:) = .false. + allocate(this%is_crop (begc:endc)); this%is_crop (:) = .false. + allocate(this%is_lake (begc:endc)); this%is_lake (:) = .false. + end subroutine col_pp_init !------------------------------------------------------------------------ @@ -178,6 +187,9 @@ subroutine col_pp_clean(this) deallocate(this%meangradz ) deallocate(this%hydrologically_active) deallocate(this%is_fates) + deallocate(this%is_soil) + deallocate(this%is_crop) + deallocate(this%is_lake) end subroutine col_pp_clean diff --git a/components/elm/src/data_types/VegetationDataType.F90 b/components/elm/src/data_types/VegetationDataType.F90 index a1b90d925b78..814bd201859e 100644 --- a/components/elm/src/data_types/VegetationDataType.F90 +++ b/components/elm/src/data_types/VegetationDataType.F90 @@ -29,6 +29,7 @@ module VegetationDataType use CNStateType , only: cnstate_type use SpeciesMod , only : species_from_string use VegetationType , only : veg_pp + use ColumnType , only : col_pp use VegetationPropertiesType , only : veg_vp use LandunitType , only : lun_pp use GridcellType , only : grc_pp @@ -2460,7 +2461,7 @@ subroutine veg_cs_init(this, begp, endp, carbon_type, ratio) this%leafcmax(p) = 0._r8 l = veg_pp%landunit(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then if (veg_pp%itype(p) == noveg) then this%leafc(p) = 0._r8 @@ -3933,7 +3934,7 @@ subroutine veg_ns_init(this, begp, endp, veg_cs) do p = begp,endp l = veg_pp%landunit(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then if (veg_pp%itype(p) == noveg) then this%leafn(p) = 0._r8 this%leafn_storage(p) = 0._r8 @@ -4414,7 +4415,7 @@ subroutine veg_ps_init(this, begp, endp, veg_cs) type(vegetation_carbon_state), intent(in) :: veg_cs ! ! !LOCAL VARIABLES: - integer :: fp,l,p ! indices + integer :: fp,l,c,p ! indices integer :: num_special_patch ! number of good values in special_patch filter integer :: special_patch (endp-begp+1) ! special landunit filter - patches !------------------------------------------------------------------------ @@ -4616,7 +4617,7 @@ subroutine veg_ps_init(this, begp, endp, veg_cs) do p = begp,endp l = veg_pp%landunit(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then if (veg_pp%itype(p) == noveg) then this%leafp(p) = 0._r8 @@ -5583,7 +5584,7 @@ subroutine veg_wf_init(this, begp, endp) do p = begp, endp l = veg_pp%landunit(p) - if (lun_pp%itype(l)==istsoil) then + if (veg_pp%is_on_soil_col(p)) then this%n_irrig_steps_left(p) = 0 this%irrig_rate(p) = 0.0_r8 end if @@ -7949,7 +7950,7 @@ subroutine veg_cf_init(this, begp, endp, carbon_type) this%xsmrpool_c13ratio(p) = spval endif end if - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then this%tempsum_npp(p) = 0._r8 this%annsum_npp(p) = 0._r8 this%availc(p) = 0._r8 @@ -9497,7 +9498,7 @@ subroutine veg_nf_init(this, begp, endp) this%soyfixn(p) = 0._r8 end if - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then this%fert_counter(p) = 0._r8 end if @@ -9879,7 +9880,7 @@ subroutine veg_pf_init(this, begp, endp) integer, intent(in) :: begp,endp ! ! !LOCAL VARIABLES: - integer :: p,l ! indices + integer :: p,c,l ! indices integer :: fp ! filter indices integer :: num_special_patch ! number of good values in special_patch filter integer :: special_patch(endp-begp+1) ! special landunit filter - patches @@ -10558,7 +10559,7 @@ subroutine veg_pf_init(this, begp, endp) end if - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then this%fert_p_counter(p) = 0._r8 end if diff --git a/components/elm/src/data_types/VegetationType.F90 b/components/elm/src/data_types/VegetationType.F90 index dd8bd7d842d4..38d2babc75a1 100644 --- a/components/elm/src/data_types/VegetationType.F90 +++ b/components/elm/src/data_types/VegetationType.F90 @@ -67,6 +67,9 @@ module VegetationType integer , pointer :: mxy (:) => null() ! m index for laixy(i,j,m),etc. (undefined for special landunits) logical , pointer :: active (:) => null() ! true=>do computations on this patch + logical, pointer :: is_on_soil_col(:) => null() ! true, if it is on a column that is soil + logical, pointer :: is_on_crop_col(:) => null() ! ture, if it is on a column that is assocaited with crops + ! Fates relevant types logical , pointer :: is_veg (:) => null() ! This is an ACTIVE fates patch logical , pointer :: is_bareground (:) => null() ! ? @@ -119,6 +122,9 @@ subroutine veg_pp_init(this, begp, endp) allocate(this%mxy (begp:endp)); this%mxy (:) = ispval allocate(this%active (begp:endp)); this%active (:) = .false. + allocate(this%is_on_soil_col(begp:endp)); this%is_on_soil_col(:) = .false. + allocate(this%is_on_crop_col(begp:endp)); this%is_on_crop_col(:) = .false. + allocate(this%is_fates (begp:endp)); this%is_fates (:) = .false. if (use_fates) then allocate(this%is_veg (begp:endp)); this%is_veg (:) = .false. @@ -148,6 +154,9 @@ subroutine veg_pp_clean(this) deallocate(this%mxy ) deallocate(this%active ) + deallocate(this%is_on_soil_col) + deallocate(this%is_on_crop_col) + deallocate(this%is_fates) if (use_fates) then deallocate(this%is_veg) diff --git a/components/elm/src/dyn_subgrid/dynConsBiogeochemMod.F90 b/components/elm/src/dyn_subgrid/dynConsBiogeochemMod.F90 index b8db8144bb87..5f8fb02a4cc1 100644 --- a/components/elm/src/dyn_subgrid/dynConsBiogeochemMod.F90 +++ b/components/elm/src/dyn_subgrid/dynConsBiogeochemMod.F90 @@ -245,7 +245,7 @@ subroutine dyn_cnbal_patch(bounds, & p = filter_soilp_with_inactive(fp) c = veg_pp%column(p) l = veg_pp%landunit(p) - if (.not.(lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) ) Then + if (.not.(veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) ) Then print *, "istsoil,istcrop",istsoil, istcrop print *, lun_pp%itype(l) end if diff --git a/components/elm/src/dyn_subgrid/dynEDMod.F90 b/components/elm/src/dyn_subgrid/dynEDMod.F90 index 29882dedb4c7..1cfde410e718 100644 --- a/components/elm/src/dyn_subgrid/dynEDMod.F90 +++ b/components/elm/src/dyn_subgrid/dynEDMod.F90 @@ -30,7 +30,7 @@ subroutine dyn_ED( bounds ) do p = bounds%begp,bounds%endp c = veg_pp%column(p) - if (col_pp%itype(c) == istsoil .and. col_pp%active(c) ) then + if (col_pp%is_soil(c) .and. col_pp%active(c) ) then if ( veg_pp%is_veg(p) .or. veg_pp%is_bareground(p)) then veg_pp%wtcol(p) = veg_pp%wt_ed(p) else diff --git a/components/elm/src/dyn_subgrid/dyncropFileMod.F90 b/components/elm/src/dyn_subgrid/dyncropFileMod.F90 old mode 100755 new mode 100644 index b8f6cf44ff0e..02b815f6565f --- a/components/elm/src/dyn_subgrid/dyncropFileMod.F90 +++ b/components/elm/src/dyn_subgrid/dyncropFileMod.F90 @@ -194,7 +194,7 @@ subroutine dyncrop_interp(bounds,crop_vars) topi = grc_pp%topi(g) ti = t - topi + 1 - if (lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_crop_col(p)) then m = veg_pp%itype(p) ! The following assumes there is a single CFT on each crop column. The diff --git a/components/elm/src/dyn_subgrid/dynpftFileMod.F90 b/components/elm/src/dyn_subgrid/dynpftFileMod.F90 index 69a3535fd51a..b01c72515c30 100644 --- a/components/elm/src/dyn_subgrid/dynpftFileMod.F90 +++ b/components/elm/src/dyn_subgrid/dynpftFileMod.F90 @@ -283,7 +283,7 @@ subroutine dynpft_interp(bounds) ! (if there is one) ! (However, currently [as of 5-9-13] the code won't let you run with transient ! Patches combined with create_crop_landunit anyway, so it's a moot point.) - if (lun_pp%itype(l) == istsoil) then + if (veg_pp%is_on_soil_col(p)) then m = veg_pp%itype(p) ! Note that the following assignment assumes that all Patches share a single column diff --git a/components/elm/src/external_models/emi/src/emi/ExternalModelInterfaceMod.F90 b/components/elm/src/external_models/emi/src/emi/ExternalModelInterfaceMod.F90 index c2c5c7d20642..46476a8ed9d2 100644 --- a/components/elm/src/external_models/emi/src/emi/ExternalModelInterfaceMod.F90 +++ b/components/elm/src/external_models/emi/src/emi/ExternalModelInterfaceMod.F90 @@ -371,8 +371,8 @@ subroutine EMI_Init_EM(em_id) do c = bounds_clump%begc,bounds_clump%endc if (col_pp%active(c)) then l = col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. col_pp%itype(c) == icol_road_perv .or. & - lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%itype(c) == icol_road_perv .or. & + col_pp%is_crop(c)) then num_e2l_filter_col = num_e2l_filter_col + 1 tmp_col(c) = 1 end if diff --git a/components/elm/src/main/controlMod.F90 b/components/elm/src/main/controlMod.F90 old mode 100755 new mode 100644 diff --git a/components/elm/src/main/elm_initializeMod.F90 b/components/elm/src/main/elm_initializeMod.F90 old mode 100755 new mode 100644 diff --git a/components/elm/src/main/elm_instMod.F90 b/components/elm/src/main/elm_instMod.F90 index ca74dda6ccf4..bf621ad0f474 100644 --- a/components/elm/src/main/elm_instMod.F90 +++ b/components/elm/src/main/elm_instMod.F90 @@ -309,10 +309,10 @@ subroutine elm_inst_biogeophys(bounds_proc) g = col_pp%gridcell(c) if (.not. use_extrasnowlayers) then ! original 5 layer shallow snowpack model - if (lun_pp%itype(l)==istice) then + if (lun_pp%itype(l) == istice) then h2osno_col(c) = h2osno_max - elseif (lun_pp%itype(l)==istice_mec .or. & - (lun_pp%itype(l)==istsoil .and. ldomain%glcmask(g) > 0._r8)) then + elseif (lun_pp%itype(l) == istice_mec .or. & + (col_pp%is_soil(c) .and. ldomain%glcmask(g) > 0._r8)) then ! Initialize a non-zero snow thickness where the ice sheet can/potentially operate. ! Using glcmask to capture all potential vegetated points around GrIS (ideally ! we would use icemask from CISM, but that isn't available until after initialization.) @@ -341,11 +341,11 @@ subroutine elm_inst_biogeophys(bounds_proc) ! a small amount of snow in places that are likely to be snow-covered for much or ! all of the year. ! amschnei@uci.edu: Initializing "deep firn" for glacier columns - if (lun_pp%itype(l)==istice .or. lun_pp%itype(l)==istice_mec) then + if (lun_pp%itype(l) == istice .or. lun_pp%itype(l) == istice_mec) then ! land ice (including multiple elevation classes, i.e. glacier_mec) h2osno_col(c) = 0.5_r8*h2osno_max ! start with half full snow column, representing deep firn snow_depth_col(c) = h2osno_col(c) / bdfirn - else if (lun_pp%itype(l)==istsoil .and. grc_pp%latdeg(g) >= 44._r8) then + else if (col_pp%is_soil(c) .and. grc_pp%latdeg(g) >= 44._r8) then ! Northern hemisphere seasonal snow h2osno_col(c) = 50._r8 snow_depth_col(c) = h2osno_col(c) / bdsno diff --git a/components/elm/src/main/elm_interface_pflotranMod.F90 b/components/elm/src/main/elm_interface_pflotranMod.F90 index 4298107c0a06..f3303afb4f96 100644 --- a/components/elm/src/main/elm_interface_pflotranMod.F90 +++ b/components/elm/src/main/elm_interface_pflotranMod.F90 @@ -484,7 +484,7 @@ subroutine interface_init(bounds) !write (iulog,*) 'WARNING: SOIL/CROP column with wtgcell <= 0 or inactive... within the domain' !write (iulog,*) 'ELM-- PFLOTRAN does not include such a SOIL/CROP column, AND will skip it' - elseif ( .not.(ltype(l)==istsoil .or. ltype(l)==istcrop) ) then + elseif ( .not.(col_pp%is_soil(c) .or. col_pp%is_crop(c)) ) then !write (iulog,*) 'WARNING: non-SOIL/CROP column found in filter%num_soilc: nc, l, ltype', nc, l, ltype(l) !write (iulog,*) 'ELM-- PFLOTRAN does not include such a SOIL/CROP column, AND will skip it' @@ -548,7 +548,7 @@ subroutine interface_init(bounds) g = cgridcell(c) gcount = g - bounds%begg + 1 - if ((.not.(ltype(l)==istsoil)) .and. (.not.(ltype(l)==istcrop)) ) then + if ((.not.(col_pp%is_soil(c))) .and. (.not.(col_pp%is_crop(c))) ) then !write (iulog,*) 'WARNING: Land Unit type of Non-SOIL/CROP... within the domain' !write (iulog,*) 'ELM-- PFLOTRAN does not support this land unit at present, AND will skip it' @@ -598,7 +598,7 @@ subroutine interface_init(bounds) g = cgridcell(c) gcount = g-bounds%begg+1 - if( (ltype(l)==istsoil .or. ltype(l)==istcrop) .and. & + if( (col_pp%is_soil(c) .or. col_pp%is_crop(c)) .and. & (cactive(c) .and. cwtgcell(c)>0._r8) ) then mapped_gid(gcount) = grc_pp%gindex(g) ! this is the globally grid-index, i.e. 'an' in its original calculation @@ -1856,7 +1856,7 @@ subroutine get_elm_soil_properties(elm_interface_data, bounds, filters) g = cgridcell(c) l = clandunit(c) - if ( (ltype(l)==istsoil .or. ltype(l)==istcrop) .and. & + if ( (col_pp%is_soil(c) .or. col_pp%is_crop(c)) .and. & (cactive(c) .and. cwtgcell(c)>0._r8) ) then ! skip inactive or zero-weighted column (may be not needed, but in case) #ifdef COLUMN_MODE diff --git a/components/elm/src/main/elm_varsur.F90 b/components/elm/src/main/elm_varsur.F90 old mode 100755 new mode 100644 diff --git a/components/elm/src/main/elmfates_interfaceMod.F90 b/components/elm/src/main/elmfates_interfaceMod.F90 index 414a8d79f97c..bc151ec8e2d0 100644 --- a/components/elm/src/main/elmfates_interfaceMod.F90 +++ b/components/elm/src/main/elmfates_interfaceMod.F90 @@ -940,7 +940,7 @@ subroutine init(this, bounds_proc, flandusepftdat) ! INTERF-TODO: WE HAVE NOT FILTERED OUT FATES SITES ON INACTIVE COLUMNS.. YET ! NEED A RUN-TIME ROUTINE THAT CLEARS AND REWRITES THE SITE LIST - if ( (lun_pp%itype(l) == istsoil) .and. (col_pp%active(c)) ) then + if ( (col_pp%is_soil(c)) .and. (col_pp%active(c)) ) then s = s + 1 collist(s) = c this%f2hmap(nc)%hsites(c) = s diff --git a/components/elm/src/main/filterMod.F90 b/components/elm/src/main/filterMod.F90 index e1dde9445598..4909c82682fd 100644 --- a/components/elm/src/main/filterMod.F90 +++ b/components/elm/src/main/filterMod.F90 @@ -309,7 +309,7 @@ subroutine setFiltersOneGroup(bounds, this_filter, include_inactive, icemask_grc if (top_pp%active(t)) then if (col_pp%active(c) .or. include_inactive) then l =col_pp%landunit(c) - if (lun_pp%lakpoi(l)) then + if (col_pp%is_lake(c)) then fl = fl + 1 this_filter(nc)%lakec(fl) = c else @@ -329,10 +329,11 @@ subroutine setFiltersOneGroup(bounds, this_filter, include_inactive, icemask_grc fnlu = 0 do p = bounds%begp,bounds%endp t =veg_pp%topounit(p) + c = veg_pp%column(p) if (top_pp%active(t)) then if (veg_pp%active(p) .or. include_inactive) then l =veg_pp%landunit(p) - if (lun_pp%lakpoi(l) ) then + if (lun_pp%lakpoi(l) .and. col_pp%is_lake(c)) then fl = fl + 1 this_filter(nc)%lakep(fl) = p else @@ -358,7 +359,7 @@ subroutine setFiltersOneGroup(bounds, this_filter, include_inactive, icemask_grc if (top_pp%active(t)) then if (col_pp%active(c) .or. include_inactive) then l =col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%is_crop(c)) then fs = fs + 1 this_filter(nc)%soilc(fs) = c end if @@ -375,7 +376,7 @@ subroutine setFiltersOneGroup(bounds, this_filter, include_inactive, icemask_grc if (top_pp%active(t)) then if (veg_pp%active(p) .or. include_inactive) then l =veg_pp%landunit(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then fs = fs + 1 this_filter(nc)%soilp(fs) = p end if @@ -393,8 +394,8 @@ subroutine setFiltersOneGroup(bounds, this_filter, include_inactive, icemask_grc if (top_pp%active(t)) then if (col_pp%active(c) .or. include_inactive) then l =col_pp%landunit(c) - if (lun_pp%itype(l) == istsoil .or. col_pp%itype(c) == icol_road_perv .or. & - lun_pp%itype(l) == istcrop) then + if (col_pp%is_soil(c) .or. col_pp%itype(c) == icol_road_perv .or. & + col_pp%is_crop(c)) then f = f + 1 this_filter(nc)%hydrologyc(f) = c @@ -422,7 +423,7 @@ subroutine setFiltersOneGroup(bounds, this_filter, include_inactive, icemask_grc if (veg_pp%active(p) .or. include_inactive) then if (.not. iscft(veg_pp%itype(p))) then l =veg_pp%landunit(p) - if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then + if (veg_pp%is_on_soil_col(p) .or. veg_pp%is_on_crop_col(p)) then fnc = fnc + 1 this_filter(nc)%soilnopcropp(fnc) = p end if @@ -530,7 +531,7 @@ subroutine setFiltersOneGroup(bounds, this_filter, include_inactive, icemask_grc l = col_pp%landunit(c) g = col_pp%gridcell(c) if ( lun_pp%itype(l) == istice_mec .or. & - (lun_pp%itype(l) == istsoil .and. icemask_grc(g) > 0.)) then + (col_pp%is_soil(c) .and. icemask_grc(g) > 0.)) then f = f + 1 this_filter(nc)%do_smb_c(f) = c end if diff --git a/components/elm/src/main/initGridCellsMod.F90 b/components/elm/src/main/initGridCellsMod.F90 index 2d8547334eea..a0e99dd2102d 100644 --- a/components/elm/src/main/initGridCellsMod.F90 +++ b/components/elm/src/main/initGridCellsMod.F90 @@ -391,14 +391,14 @@ subroutine set_landunit_veg_compete (ltype, gi, ti,topo_ind, li, ci, pi, setdata call add_landunit(li=li, ti=ti, ltype=ltype, wttopounit=wtlunit2topounit) ! Assume one column on the landunit - call add_column(ci=ci, li=li, ctype=1, wtlunit=1.0_r8) + call add_column(ci=ci, li=li, ctype=1, wtlunit=1.0_r8, is_soil=.true.) do m = natpft_lb,natpft_ub if(use_fates .and. .not.use_fates_sp)then p_wt = 1.0_r8/real(natpft_size,r8) else p_wt = wt_nat_patch(gi,topo_ind,m) end if - call add_patch(pi=pi, ci=ci, ptype=m, wtcol=p_wt) + call add_patch(pi=pi, ci=ci, ptype=m, wtcol=p_wt, is_on_soil_col=.true.) end do ! add polygonal landunits and columns if feature turned on @@ -409,7 +409,7 @@ subroutine set_landunit_veg_compete (ltype, gi, ti,topo_ind, li, ci, pi, setdata ! get new weight for wttopounit: wtpoly2lndunit = wt_lunit(gi, topo_ind, z) call add_polygon_landunit(li=li, ti=ti, ltype=ltype, wttopounit=wtpoly2lndunit, polytype = z - max_non_poly_lunit) - call add_column(ci=ci, li=li, ctype=1, wtlunit=1.0_r8) + call add_column(ci=ci, li=li, ctype=1, wtlunit=1.0_r8, is_soil=.true.) ! add patch: do m = natpft_lb,natpft_ub if(use_fates .and. .not.use_fates_sp) then @@ -417,7 +417,7 @@ subroutine set_landunit_veg_compete (ltype, gi, ti,topo_ind, li, ci, pi, setdata else p_wt = wt_nat_patch(gi,topo_ind,m) end if - call add_patch(pi=pi, ci=ci, ptype=m, wtcol=p_wt) + call add_patch(pi=pi, ci=ci, ptype=m, wtcol=p_wt, is_on_soil_col=.true.) end do end do end if @@ -459,6 +459,7 @@ subroutine set_landunit_wet_ice_lake (ltype, gi, ti,topo_ind, li, ci, pi, setdat integer :: npfts ! number of pfts in landunit real(r8) :: wtlunit2topounit ! landunit weight in topounit real(r8) :: wtcol2lunit ! col weight in landunit + logical :: is_lake_col !------------------------------------------------------------------------ ! Set decomposition properties @@ -467,9 +468,11 @@ subroutine set_landunit_wet_ice_lake (ltype, gi, ti,topo_ind, li, ci, pi, setdat ! gridcell as the new landunit weights on each topounit. ! Later, this information will come from new surface datasat. + is_lake_col = .false. if (ltype == istwet) then call subgrid_get_topounitinfo(ti, gi,tgi=topo_ind, nwetland=npfts) else if (ltype == istdlak) then + is_lake_col = .true. call subgrid_get_topounitinfo(ti, gi,tgi=topo_ind, nlake=npfts) else if (ltype == istice) then call subgrid_get_topounitinfo(ti, gi,tgi=topo_ind, nglacier=npfts) @@ -522,7 +525,7 @@ subroutine set_landunit_wet_ice_lake (ltype, gi, ti,topo_ind, li, ci, pi, setdat ! and that each column has its own pft call add_landunit(li=li, ti=ti, ltype=ltype, wttopounit=wtlunit2topounit) - call add_column(ci=ci, li=li, ctype=ltype, wtlunit=1.0_r8) + call add_column(ci=ci, li=li, ctype=ltype, wtlunit=1.0_r8, is_lake=is_lake_col) call add_patch(pi=pi, ci=ci, ptype=noveg, wtcol=1.0_r8) end if ! ltype = istice_mec @@ -591,8 +594,8 @@ subroutine set_landunit_crop_noncompete (ltype, gi, ti,topo_ind, li, ci, pi, set if (create_crop_landunit) then do m = cft_lb, cft_ub - call add_column(ci=ci, li=li, ctype=((istcrop*100) + m), wtlunit=wt_cft(gi,topo_ind,m)) - call add_patch(pi=pi, ci=ci, ptype=m, wtcol=1.0_r8) + call add_column(ci=ci, li=li, ctype=((istcrop*100) + m), wtlunit=wt_cft(gi,topo_ind,m), is_crop=.true.) + call add_patch(pi=pi, ci=ci, ptype=m, wtcol=1.0_r8, is_on_crop_col=.true.) end do end if @@ -1165,10 +1168,11 @@ subroutine initGhostPatch() grid_count(:) = 0.d0 last_lun_type = -1 do c = bounds_proc%begc_all, bounds_proc%endc_all - g = col_pp%gridcell(c) - if (last_lun_type /= lun_pp%itype(col_pp%landunit(c))) then + g = col_pp%gridcell(c) + l = col_pp%landunit(c) + if (last_lun_type /= lun_pp%itype(l)) then grid_count(:) = 0.d0 - last_lun_type = lun_pp%itype(col_pp%landunit(c)) + last_lun_type = lun_pp%itype(l) endif grid_count(g) = grid_count(g) + 1.d0 col_rank(c) = grid_count(g) diff --git a/components/elm/src/main/initSubgridMod.F90 b/components/elm/src/main/initSubgridMod.F90 index 1dddae8fee97..e496bce442eb 100644 --- a/components/elm/src/main/initSubgridMod.F90 +++ b/components/elm/src/main/initSubgridMod.F90 @@ -537,7 +537,7 @@ end subroutine add_polygon_landunit !----------------------------------------------------------------------- - subroutine add_column(ci, li, ctype, wtlunit) + subroutine add_column(ci, li, ctype, wtlunit, is_soil, is_crop, is_lake) ! ! !DESCRIPTION: ! Add an entry in the column-level arrays. ci gives the index of the last column @@ -549,6 +549,9 @@ subroutine add_column(ci, li, ctype, wtlunit) integer , intent(in) :: li ! landunit index on which this column should be placed (assumes this landunit has already been created) integer , intent(in) :: ctype ! column type real(r8) , intent(in) :: wtlunit ! weight of the column relative to the landunit + logical , optional :: is_soil ! true for a soil column + logical , optional :: is_crop ! true for a crop column + logical , optional :: is_lake ! true for a lake column ! ! !LOCAL VARIABLES: character(len=*), parameter :: subname = 'add_column' @@ -562,11 +565,15 @@ subroutine add_column(ci, li, ctype, wtlunit) col_pp%wtlunit(ci) = wtlunit col_pp%itype(ci) = ctype - + + if (present(is_soil)) col_pp%is_soil(ci) = is_soil + if (present(is_crop)) col_pp%is_crop(ci) = is_crop + if (present(is_lake)) col_pp%is_lake(ci) = is_lake + end subroutine add_column !----------------------------------------------------------------------- - subroutine add_patch(pi, ci, ptype, wtcol) + subroutine add_patch(pi, ci, ptype, wtcol, is_on_soil_col, is_on_crop_col) ! ! !DESCRIPTION: ! Add an entry in the patch-level arrays. pi gives the index of the last patch added; the @@ -582,6 +589,8 @@ subroutine add_patch(pi, ci, ptype, wtcol) integer , intent(in) :: ci ! column index on which this patch should be placed (assumes this column has already been created) integer , intent(in) :: ptype ! patch type real(r8) , intent(in) :: wtcol ! weight of the patch relative to the column + logical, optional :: is_on_soil_col + logical, optional :: is_on_crop_col ! ! !LOCAL VARIABLES: integer :: li ! landunit index, for convenience @@ -608,6 +617,9 @@ subroutine add_patch(pi, ci, ptype, wtcol) veg_pp%mxy(pi) = ispval end if + if (present(is_on_soil_col)) veg_pp%is_on_soil_col(pi) = is_on_soil_col + if (present(is_on_crop_col)) veg_pp%is_on_crop_col(pi) = is_on_crop_col + end subroutine add_patch diff --git a/components/elm/src/main/initVerticalMod.F90 b/components/elm/src/main/initVerticalMod.F90 old mode 100755 new mode 100644 index b24aad2d4622..f0458d024c7f --- a/components/elm/src/main/initVerticalMod.F90 +++ b/components/elm/src/main/initVerticalMod.F90 @@ -312,7 +312,7 @@ subroutine initVertical(bounds, snow_depth, thick_wall, thick_roof) l = col_pp%landunit(c) if (lun_pp%urbpoi(l)) then - if (col_pp%itype(c)==icol_sunwall .or. col_pp%itype(c)==icol_shadewall) then + if (col_pp%itype(c) == icol_sunwall .or. col_pp%itype(c) == icol_shadewall) then col_pp%z(c,1:nlevurb) = zurb_wall(l,1:nlevurb) col_pp%zi(c,0:nlevurb) = ziurb_wall(l,0:nlevurb) col_pp%dz(c,1:nlevurb) = dzurb_wall(l,1:nlevurb) @@ -321,7 +321,7 @@ subroutine initVertical(bounds, snow_depth, thick_wall, thick_roof) col_pp%zi(c,nlevurb+1:nlevgrnd) = spval col_pp%dz(c,nlevurb+1:nlevgrnd) = spval end if - else if (col_pp%itype(c)==icol_roof) then + else if (col_pp%itype(c) == icol_roof) then col_pp%z(c,1:nlevurb) = zurb_roof(l,1:nlevurb) col_pp%zi(c,0:nlevurb) = ziurb_roof(l,0:nlevurb) col_pp%dz(c,1:nlevurb) = dzurb_roof(l,1:nlevurb) @@ -335,7 +335,7 @@ subroutine initVertical(bounds, snow_depth, thick_wall, thick_roof) col_pp%zi(c,0:nlevgrnd) = zisoi(0:nlevgrnd) col_pp%dz(c,1:nlevgrnd) = dzsoi(1:nlevgrnd) end if - else if (lun_pp%itype(l) /= istdlak) then + else if (.not.col_pp%is_lake(c)) then col_pp%z(c,1:nlevgrnd) = zsoi(1:nlevgrnd) col_pp%zi(c,0:nlevgrnd) = zisoi(0:nlevgrnd) col_pp%dz(c,1:nlevgrnd) = dzsoi(1:nlevgrnd) @@ -428,7 +428,7 @@ subroutine initVertical(bounds, snow_depth, thick_wall, thick_roof) do c = bounds%begc,bounds%endc l = col_pp%landunit(c) - if (lun_pp%itype(l) == istdlak) then + if (col_pp%is_lake(c)) then if (col_pp%lakedepth(c) == spval) then col_pp%lakedepth(c) = zlak(nlevlak) + 0.5_r8*dzlak(nlevlak) @@ -476,13 +476,12 @@ subroutine initVertical(bounds, snow_depth, thick_wall, thick_roof) associate(snl => col_pp%snl) ! Output: [integer (:) ] number of snow layers do c = bounds%begc,bounds%endc - l = col_pp%landunit(c) col_pp%dz(c,-nlevsno+1: 0) = spval col_pp%z (c,-nlevsno+1: 0) = spval col_pp%zi(c,-nlevsno :-1) = spval - if (.not. lun_pp%lakpoi(l)) then + if (.not. col_pp%is_lake(c)) then if (snow_depth(c) < 0.01_r8) then snl(c) = 0 col_pp%dz(c,-nlevsno+1:0) = 0._r8 @@ -646,7 +645,7 @@ subroutine initVertical(bounds, snow_depth, thick_wall, thick_roof) if (lun_pp%urbpoi(l) .and. col_pp%itype(c) /= icol_road_imperv .and. col_pp%itype(c) /= icol_road_perv) then col_pp%nlevbed(c) = nlevurb - else if (lun_pp%itype(l) == istdlak) then + else if (col_pp%is_lake(c)) then col_pp%nlevbed(c) = nlevlak else if (lun_pp%itype(l) == istice_mec) then col_pp%nlevbed(c) = 5 @@ -686,7 +685,7 @@ subroutine initVertical(bounds, snow_depth, thick_wall, thick_roof) do c = begc,endc l = col_pp%landunit(c) - if (lun_pp%itype(l)==istice_mec) then + if (lun_pp%itype(l) == istice_mec) then ! ice_mec columns already account for subgrid topographic variability through ! their use of multiple elevation classes; thus, to avoid double-accounting for ! topographic variability in these columns, we ignore topo_std and use a value diff --git a/components/elm/src/main/lnd2glcMod.F90 b/components/elm/src/main/lnd2glcMod.F90 index f938bb7abf23..9b9b4673f737 100644 --- a/components/elm/src/main/lnd2glcMod.F90 +++ b/components/elm/src/main/lnd2glcMod.F90 @@ -182,7 +182,7 @@ subroutine update_lnd2glc(this, bounds, num_do_smb_c, filter_do_smb_c, init) if (lun_pp%itype(l) == istice_mec) then n = col_itype_to_icemec_class(col_pp%itype(c)) flux_normalization = 1.0_r8 - else if (lun_pp%itype(l) == istsoil) then + else if (col_pp%is_soil(c)) then n = 0 !0-level index (bareland information) flux_normalization = bareland_normalization(c) else diff --git a/components/elm/src/main/subgridRestMod.F90 b/components/elm/src/main/subgridRestMod.F90 old mode 100755 new mode 100644 index 0dbd0888ead0..1bb19f8fdf44 --- a/components/elm/src/main/subgridRestMod.F90 +++ b/components/elm/src/main/subgridRestMod.F90 @@ -686,7 +686,7 @@ subroutine check_weights(bounds) do p = bounds%begp, bounds%endp l = veg_pp%landunit(p) - if (lun_pp%itype(l) == istsoil) then + if (veg_pp%is_on_soil_col(p)) then diff = abs(veg_pp%wtlunit(p) - pft_wtlunit_before_rest_read(p)) if (diff > tol) then write(iulog,*) 'ERROR: PFT weights are SIGNIFICANTLY different between the restart (finidat) file' diff --git a/components/elm/src/main/subgridWeightsMod.F90 b/components/elm/src/main/subgridWeightsMod.F90 old mode 100755 new mode 100644 index 0367efd71cbe..338d00789b2d --- a/components/elm/src/main/subgridWeightsMod.F90 +++ b/components/elm/src/main/subgridWeightsMod.F90 @@ -958,10 +958,10 @@ subroutine set_pct_pft_diagnostics(bounds) !topi = grc_pp%topi(g) !ti = t - topi + 1 ptype = veg_pp%itype(p) - if (lun_pp%itype(l) == istsoil) then + if (veg_pp%is_on_soil_col(p)) then ptype_1indexing = ptype + (1 - natpft_lb) subgrid_weights_diagnostics%pct_nat_pft(t, ptype_1indexing) = veg_pp%wtlunit(p) * 100._r8 - else if (lun_pp%itype(l) == istcrop) then + else if (veg_pp%is_on_crop_col(p)) then ptype_1indexing = ptype + (1 - cft_lb) subgrid_weights_diagnostics%pct_cft(t, ptype_1indexing) = veg_pp%wtlunit(p) * 100._r8 end if diff --git a/components/elm/src/main/surfrdMod.F90 b/components/elm/src/main/surfrdMod.F90 old mode 100755 new mode 100644 diff --git a/components/elm/src/utils/domainMod.F90 b/components/elm/src/utils/domainMod.F90 old mode 100755 new mode 100644