Dynamic urban

bld/CLMBuildNamelist.pm

@@ -2455,6 +2455,7 @@ sub setup_logic_dynamic_subgrid {
    setup_logic_do_transient_pfts($opts, $nl_flags, $definition, $defaults, $nl);
    setup_logic_do_transient_crops($opts, $nl_flags, $definition, $defaults, $nl);
    setup_logic_do_transient_lakes($opts, $nl_flags, $definition, $defaults, $nl);
+   setup_logic_do_transient_urban($opts, $nl_flags, $definition, $defaults, $nl);
    setup_logic_do_harvest($opts, $nl_flags, $definition, $defaults, $nl);
 
    add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'reset_dynbal_baselines');
@@ -2683,6 +2684,69 @@ sub setup_logic_do_transient_lakes {
    }
 }
 
+sub setup_logic_do_transient_urban {
+   #
+   # Set do_transient_urban default value, and perform error checking on do_transient_urban
+   #
+   # Assumes the following are already set in the namelist (although it's okay
+   # for them to be unset if that will be their final state):
+   # - flanduse_timeseries
+   #
+   # NOTE(kwo, 2021-08-11) I based this function on setup_logic_do_transient_lakes. 
+   # As in NOTE(wjs, 2020-08-23) I'm not sure if all of the checks here are truly important 
+   # for transient urban (in particular, my guess is that collapse_urban could probably be done with transient
+   # urban - as well as transient pfts and transient crops for that matter), but some of
+   # the checks probably are needed, and it seems best to keep transient urban consistent
+   # with other transient areas in this respect.
+   my ($opts, $nl_flags, $definition, $defaults, $nl) = @_;
+
+   my $var = 'do_transient_urban';
+
+   # cannot_be_true will be set to a non-empty string in any case where
+   # do_transient_urban should not be true; if it turns out that
+   # do_transient_urban IS true in any of these cases, a fatal error will be
+   # generated
+   my $cannot_be_true = "";
+
+   my $n_dom_pfts = $nl->get_value( 'n_dom_pfts' );
+   my $n_dom_landunits = $nl->get_value( 'n_dom_landunits' );
+   my $toosmall_soil = $nl->get_value( 'toosmall_soil' );
+   my $toosmall_crop = $nl->get_value( 'toosmall_crop' );
+   my $toosmall_glacier = $nl->get_value( 'toosmall_glacier' );
+   my $toosmall_lake = $nl->get_value( 'toosmall_lake' );
+   my $toosmall_wetland = $nl->get_value( 'toosmall_wetland' );
+   my $toosmall_urban = $nl->get_value( 'toosmall_urban' );
+
+   if (string_is_undef_or_empty($nl->get_value('flanduse_timeseries'))) {
+      $cannot_be_true = "$var can only be set to true when running a transient case (flanduse_timeseries non-blank)";
+   }
+
+   if (!$cannot_be_true) {
+      # Note that, if the variable cannot be true, we don't call add_default
+      # - so that we don't clutter up the namelist with variables that don't
+      # matter for this case
+      add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, $var);
+   }
+
+   # Make sure the value is false when it needs to be false - i.e., that the
+   # user hasn't tried to set a true value at an inappropriate time.
+
+   if (&value_is_true($nl->get_value($var)) && $cannot_be_true) {
+      $log->fatal_error($cannot_be_true);
+   }
+
+   # if do_transient_urban is .true. and any of these (n_dom_* or toosmall_*)
+   # are > 0 or collapse_urban = .true., then give fatal error
+   if (&value_is_true($nl->get_value($var))) {
+      if (&value_is_true($nl->get_value('collapse_urban'))) {
+         $log->fatal_error("$var cannot be combined with collapse_urban");
+      }
+      if ($n_dom_pfts > 0 || $n_dom_landunits > 0 || $toosmall_soil > 0 || $toosmall_crop > 0 || $toosmall_glacier > 0 || $toosmall_lake > 0 || $toosmall_wetland > 0 || $toosmall_urban > 0) {
+         $log->fatal_error("$var cannot be combined with any of the of the following > 0: n_dom_pfts > 0, n_dom_landunit > 0, toosmall_soil > 0._r8, toosmall_crop > 0._r8, toosmall_glacier > 0._r8, toosmall_lake > 0._r8, toosmall_wetland > 0._r8, toosmall_urban > 0._r8");
+      }
+   }
+}
+
 sub setup_logic_do_harvest {
    #
    # Set do_harvest default value, and perform error checking on do_harvest

bld/namelist_files/namelist_defaults_ctsm.xml

@@ -2569,6 +2569,7 @@ lnd/clm2/surfdata_map/release-clm5.0.30/surfdata_ne0np4.CONUS.ne30x8_hist_78pfts
 <!-- are set in the BuildNamelist code)         -->
 <!-- =========================================  -->
 <do_transient_lakes>.false.</do_transient_lakes>
+<do_transient_urban>.false.</do_transient_urban>
 <reset_dynbal_baselines>.false.</reset_dynbal_baselines>
 
 <!-- =========================================  -->

bld/namelist_files/namelist_definition_ctsm.xml

@@ -2484,6 +2484,12 @@ If TRUE, apply transient lakes from flanduse_timeseries file.
 (Only valid for transient runs, where there is a flanduse_timeseries file.)
 </entry>
 
+<entry id="do_transient_urban" type="logical" category="physics"
+       group="dynamic_subgrid" valid_values="" >
+If TRUE, apply transient urban from flanduse_timeseries file.
+(Only valid for transient runs, where there is a flanduse_timeseries file.)
+</entry>
+
 <entry id="do_harvest" type="logical" category="physics"
        group="dynamic_subgrid" valid_values="" >
 If TRUE, apply harvest from flanduse_timeseries file.

cime_config/testdefs/testlist_clm.xml

@@ -1669,6 +1669,16 @@
       </machine>
     </machines>
   </test>
+  <test name="ERS_Lm25" grid="1x1_smallvilleIA" compset="IHistClm50BgcCropQianRs" testmods="clm/smallville_dynurban_monthly">
+    <machines>
+      <machine name="cheyenne" compiler="gnu" category="aux_clm">
+        <options>
+          <option name="wallclock">0:50:00</option>
+          <option name="comment">Include a test of transient urban</option>
+        </options>
+      </machine>
+    </machines>
+  </test>
   <test name="SMS_D_P48x1_Ld5" grid="f10_f10_mg37" compset="I2000Clm50BgcCrop" testmods="clm/irrig_spunup">
     <machines>
       <machine name="izumi" compiler="nag" category="aux_clm"/>

cime_config/testdefs/testmods_dirs/clm/smallville_dynurban_monthly/include_user_mods

@@ -0,0 +1 @@
+../monthly

cime_config/testdefs/testmods_dirs/clm/smallville_dynurban_monthly/modify_smallville_with_dynurban.ncl

@@ -0,0 +1,63 @@
+; NCL script
+; modify_smallville_with_dynurban.ncl
+; Keith Oleson, Dec 2021
+; Purpose is to create a dynamic urban file for the smallville grid for test
+; ERS_Lm25.1x1_smallvilleIA.IHistClm50BgcCropQianRs.cheyenne_gnu.clm-smallville_dynurban_monthly
+;**************************************
+
+load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
+load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
+load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
+load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
+
+begin
+
+ print ("=========================================")
+ print ("Start Time: "+systemfunc("date") )
+ print ("=========================================")
+
+ infile  = "/glade/p/cgd/tss/people/oleson/modify_surfdata/landuse.timeseries_1x1_smallvilleIA_hist_78pfts_simyr1850-1855_c160127.nc"
+ outfile = "/glade/p/cgd/tss/people/oleson/modify_surfdata/landuse.timeseries_1x1_smallvilleIA_hist_78pfts_simyr1850-1855_dynUrban_c211206.nc" 
+
+ system("cp " + infile + " " + outfile)
+
+ outf = addfile(outfile,"w")
+
+ numurbl = 3
+
+ pct_crop = outf->PCT_CROP 
+ printVarSummary(pct_crop)
+ pct_urban = new((/dimsizes(pct_crop(:,0,0)),numurbl,dimsizes(pct_crop(0,:,0)),dimsizes(pct_crop(0,0,:))/),double,"No_FillValue")
+ pct_urban!0 = "time"
+ pct_urban&time = pct_crop&time
+ pct_urban!1 = "numurbl"
+ pct_urban!2 = pct_crop!1
+ pct_urban!3 = pct_crop!2
+ pct_urban@long_name = "percent urban for each density type (tbd, hd, md)"
+ pct_urban@units = "unitless"
+ printVarSummary(pct_urban)
+
+ hasurban = new((/numurbl,dimsizes(pct_crop(0,:,0)),dimsizes(pct_crop(0,0,:))/),double,"No_FillValue")
+ hasurban!0 = pct_urban!1
+ hasurban!1 = pct_urban!2
+ hasurban!2 = pct_urban!3
+ hasurban = 1.d
+ printVarSummary(hasurban)
+
+ pct_urban(:,0,0,0) = (/0.d,20.d,10.d,10.d,10.d,10.d/)
+ pct_urban(:,1,0,0) = (/0.d,15.d, 8.d, 8.d, 8.d, 8.d/)
+ pct_urban(:,2,0,0) = (/0.d,10.d, 5.d, 5.d, 5.d, 5.d/)
+
+ pct_crop(:,0,0) = (/0.,25.,12.,12.,12.,12./)
+
+ outf->HASURBAN  = hasurban
+ outf->PCT_URBAN = pct_urban
+ outf->PCT_CROP = pct_crop
+
+ outf@history = "This file was created with the following NCL script: /glade/p/cgd/tss/people/oleson/modify_surfdata/modify_smallville_with_dynurban.ncl.  The file used as a template is: /glade/p/cesm/cseg/inputdata/lnd/clm2/surfdata_map/landuse.timeseries_1x1_smallvilleIA_hist_78pfts_simyr1850-1855_c160127.nc.  Key points are that urban area starts as 0, increases after the first year, then decreases after the second year.  PCT_CROP is also changed so that PCT_URBAN + PCT_CROP <= 100. (Here, PCT_CROP increases and decreases at the same time as PCT_URBAN in order to exercise the simultaneous increase or decrease of two landunits, but that isn't a critical part of this test.). Note that the use of this file means that this testmod can only be used with the 1x1_smallvilleIA grid."
+
+ print ("=========================================")
+ print ("Finish Time: "+systemfunc("date") )
+ print ("=========================================")
+
+end

cime_config/testdefs/testmods_dirs/clm/smallville_dynurban_monthly/user_nl_clm

@@ -0,0 +1,12 @@
+do_transient_urban = .true.
+!KO The following run_zero_weight_urban setting is temporary until the HASURBAN methdology is implemented.
+run_zero_weight_urban = .true.
+
+! This file was created with the following NCL script:
+! /glade/p/cgd/tss/people/oleson/modify_surfdata/modify_smallville_with_dynurban.ncl
+! The file used as a template is: 
+! /glade/p/cesm/cseg/inputdata/lnd/clm2/surfdata_map/landuse.timeseries_1x1_smallvilleIA_hist_78pfts_simyr1850-1855_c160127.nc
+! Key points are that urban area starts as 0, increases after the first year, then decreases after the second year.
+! PCT_CROP is also changed so that PCT_URBAN + PCT_CROP <= 100. (Here, PCT_CROP increases and decreases at the same time as PCT_URBAN in order to exercise the simultaneous increase or decrease of two landunits, but that isn't a critical part of this test.)
+! Note that the use of this file means that this testmod can only be used with the 1x1_smallvilleIA grid.
+flanduse_timeseries = '$DIN_LOC_ROOT/lnd/clm2/surfdata_map/landuse.timeseries_1x1_smallvilleIA_hist_78pfts_simyr1850-1855_dynUrban_c211206.nc'

doc/design/dynamic_urban.rst

@@ -0,0 +1,45 @@
+.. sectnum::
+
+.. contents::
+
+==================================
+ Overview of this design document
+==================================
+
+This documents some of the high-level design decisions made during implementation of 
+dynamic urban landunits.
+
+============================================================================
+ The use of dzsoi_decomp for urban landunits to calculate totcolch4 in ch4Mod.F90
+============================================================================
+During the first test simulation for dynamic urban, we encountered a methane conservation 
+error the first time PCT_URBAN changed. The dynamic adjustments for conc_ch4_sat_col and 
+conc_ch4_unsat_col (the column_state_updater in subroutine DynamicColumnAdjustments within
+ch4Mod.F90) were distributing non-zero values for roof and walls for layers 1,nlevsoi.
+When the total column ch4 is summed over the soil layers (or in this case, urban layers), the 
+summation is done over nlevsoi, not nlevurb, using dz. dz is 1.e36 for roof/wall layers 
+that are greater than nlevurb, thus creating an imbalance.
+
+Rather than trying to keep the BGC variables physically meaningful in urban landunits,
+we will just pack these variables in a way that should conserve these variables, even if
+the values in each of the urban columns is somewhat nonsensical. Specifically: we'll take
+col%wtgcell at face value in urban columns in dynColumnStateUpdaterMod - i.e., for the sake
+of storing / conserving these BGC variables, we'll act as if that gives the true column
+weight on the grid cell. This way we'll end up storing all of the C & N from the vegetated
+column in the urban columns, and there shouldn't be any that is lost from the system. If that
+urban landunit later shrinks, the stored C & N should be restored symmetrically. It shouldn't
+really matter that it was stored in a non-physical way (e.g., with some C & N stored in urban
+walls), since the BGC variables are irrelevant over the urban areas and we just want to be able
+to restore the amount that was originally stored if an urban landunit grows and then later shrinks.
+But for this to work right, we need to treat the relevant BGC variables as having the same dz over
+all urban columns as over the soil column. Note that there already seems to be an implicit assumption
+that dz is the same for all columns in the dynamic column state updates, in that dz doesn't enter
+into the conservation equations. In terms of what needs to change, we think that the only relevant
+code is the code that sums up total C / N / CH4 for the sake of balance checks: these balance checks
+need to be consistent with the assumptions made in the conservation code. The C and N summations
+already use dzsoi_decomp, which is the same for all columns, so this is already what we want.
+The only thing that needs to change is the use of dz in totcolch4 in ch4Mod.F90: we've changed that to now use
+dzsoi_decomp over urban columns. (This begs the question of why this isn't already using
+dzsoi_decomp for consistency with the C & N code; we're not sure about this.)
+
+See issue #1445 for the original discussion on this topic.

src/biogeochem/ch4Mod.F90

@@ -4371,7 +4371,9 @@ subroutine ch4_totcolch4(bounds, num_nolakec, filter_nolakec, num_lakec, filter_
     ! its original values
     !
     ! !USES:
-    use ch4varcon  , only : allowlakeprod
+    use ch4varcon       , only : allowlakeprod
+    use clm_varcon      , only : dzsoi_decomp
+    use landunit_varcon , only : isturb_tbd, isturb_hd, isturb_md
     !
     ! !ARGUMENTS:
     type(bounds_type) , intent(in)    :: bounds   
@@ -4384,7 +4386,7 @@ subroutine ch4_totcolch4(bounds, num_nolakec, filter_nolakec, num_lakec, filter_
     !
     ! !LOCAL VARIABLES:
     integer :: fc, c
-    integer :: j
+    integer :: j, l
 
     character(len=*), parameter       :: subname = 'ch4_totcolch4'
     !-----------------------------------------------------------------------
@@ -4411,9 +4413,18 @@ subroutine ch4_totcolch4(bounds, num_nolakec, filter_nolakec, num_lakec, filter_
     do j = 1, nlevsoi
        do fc = 1, num_nolakec
           c = filter_nolakec(fc)
-          totcolch4(c) = totcolch4(c) + &
-               (finundated(c)*conc_ch4_sat(c,j) + (1._r8-finundated(c))*conc_ch4_unsat(c,j)) * &
-               dz(c,j)*catomw
+          l = col%landunit(c)
+          ! See doc/design/dynamic_urban.rst for an explanation of why we use dzsoi_decomp instead of dz for 
+          ! urban landunits.
+          if (lun%itype(l) .eq. isturb_tbd .or. lun%itype(l) .eq. isturb_hd .or. lun%itype(l) .eq. isturb_md) then
+             totcolch4(c) = totcolch4(c) + &
+                  (finundated(c)*conc_ch4_sat(c,j) + (1._r8-finundated(c))*conc_ch4_unsat(c,j)) * &
+                  dzsoi_decomp(j)*catomw
+          else
+             totcolch4(c) = totcolch4(c) + &
+                  (finundated(c)*conc_ch4_sat(c,j) + (1._r8-finundated(c))*conc_ch4_unsat(c,j)) * &
+                  dz(c,j)*catomw
+          end if
           ! mol CH4 --> g C
        end do
 

src/biogeophys/UrbanParamsType.F90

@@ -357,9 +357,8 @@ subroutine Init(this, bounds)
        end if
     end do
 
-    ! Deallocate memory for urbinp datatype
-
-    call UrbanInput(bounds%begg, bounds%endg, mode='finalize')
+    ! Note that we don't deallocate memory for urbinp datatype (call UrbanInput with
+    ! mode='finalize') because the arrays are needed for dynamic urban landunits.
 
   end subroutine Init
 

src/dyn_subgrid/dynConsBiogeophysMod.F90

@@ -382,9 +382,7 @@ subroutine set_glacier_baselines(bounds, num_icec, filter_icec, &
     ! all we really need are averages for the natural veg landunit. But since this
     ! subroutine is only called in initialization, code reuse/simplicity is more important
     ! than performance.
-    !
-    ! No thought has been given to c2l_scale_type here. This may need to be fixed if we
-    ! ever start caring about urban averages in this routine.
+
     call c2l(bounds = bounds, &
          carr = vals_col(bounds%begc:bounds%endc), &
          larr = vals_lun(bounds%begl:bounds%endl), &