[WIP] Changes the way columns and patches are identified

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

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,17 +1057,17 @@ 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
        end if
 
        ! 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

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

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

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)

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

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,20 +940,21 @@ 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
 
     if (.not. (fan_to_bgc_veg .or. fan_to_bgc_crop)) return
 
     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

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) ) &

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

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

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

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)

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