Dynamic lakes - without tools changes

bld/CLMBuildNamelist.pm

@@ -2382,6 +2382,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_harvest($opts, $nl_flags, $definition, $defaults, $nl);
 
    add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'reset_dynbal_baselines');
@@ -2547,6 +2548,69 @@ sub setup_logic_do_transient_crops {
    }
 }
 
+sub setup_logic_do_transient_lakes {
+   #
+   # Set do_transient_lakes default value, and perform error checking on do_transient_lakes
+   #
+   # 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(wjs, 2020-08-23) I based this function on setup_logic_do_transient_crops. I'm
+   # not sure if all of the checks here are truly important for transient lakes (in
+   # particular, my guess is that collapse_urban could probably be done with transient
+   # lakes - 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 lakes consistent
+   # with other transient areas in this respect.
+   my ($opts, $nl_flags, $definition, $defaults, $nl) = @_;
+
+   my $var = 'do_transient_lakes';
+
+   # cannot_be_true will be set to a non-empty string in any case where
+   # do_transient_lakes should not be true; if it turns out that
+   # do_transient_lakes 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_lakes 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

@@ -3599,6 +3599,7 @@ lnd/clm2/surfdata_map/release-clm5.0.30/surfdata_ne0np4.CONUS.ne30x8_hist_78pfts
 <!-- (some other defaults for dynamic subgrid   -->
 <!-- are set in the BuildNamelist code)         -->
 <!-- =========================================  -->
+<do_transient_lakes>.false.</do_transient_lakes>
 <reset_dynbal_baselines>.false.</reset_dynbal_baselines>
 
 <!-- =========================================  -->

bld/namelist_files/namelist_definition_ctsm.xml

@@ -2424,6 +2424,12 @@ If TRUE, apply transient crops from flanduse_timeseries file.
 (Only valid for transient runs, where there is a flanduse_timeseries file.)
 </entry>
 
+<entry id="do_transient_lakes" type="logical" category="physics"
+       group="dynamic_subgrid" valid_values="" >
+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_harvest" type="logical" category="physics"
        group="dynamic_subgrid" valid_values="" >
 If TRUE, apply harvest from flanduse_timeseries file.

cime_config/testdefs/testlist_clm.xml

@@ -1622,6 +1622,16 @@ for ERS test as otherwise it won't work for a sub-day test"</option>
       <option name="comment"  >Want a debug test that tests a number of aspects of transient crops, including a new crop landunit and shifting PCT_CFT; move to CLM50 once we can get it fast enough (see bug 2391)</option>
     </options>
   </test>
+  <test name="ERS_Lm25" grid="1x1_smallvilleIA" compset="IHistClm50BgcCropQianRsGs" testmods="clm/smallville_dynlakes_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 lakes</option>
+        </options>
+      </machine>
+    </machines>
+  </test>
   <test name="SMS_D_P48x1_Ld5" grid="f10_f10_musgs" compset="I2000Clm50BgcCrop" testmods="clm/irrig_spunup">
     <machines>
       <machine name="izumi" compiler="nag" category="aux_clm"/>

cime_config/testdefs/testmods_dirs/clm/smallville_dynlakes_monthly/include_user_mods

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

cime_config/testdefs/testmods_dirs/clm/smallville_dynlakes_monthly/user_nl_clm

@@ -0,0 +1,12 @@
+do_transient_lakes = .true.
+
+! This file was created with the following command:
+! ncap2 -s 'PCT_LAKE=array(0.0,0.0,PCT_CROP); PCT_LAKE={0.,50.,25.,25.,25.,25.}; HASLAKE=array(1.,1.,AREA); PCT_CROP=array(0.0,0.0,PCT_LAKE); PCT_CROP={0.,25.,12.,12.,12.,12.}' landuse.timeseries_1x1_smallvilleIA_hist_78pfts_simyr1850-1855_c160127.nc landuse.timeseries_1x1_smallvilleIA_hist_78pfts_simyr1850-1855_dynLakes_c200928.nc
+! Key points are that lake area starts as 0, increases after the first year, then decreases after the second year.
+! PCT_CROP is also changed so that PCT_LAKE + PCT_CROP <= 100. (Here, PCT_CROP increases and decreases at the same time as PCT_LAKE 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_dynLakes_c200928.nc'
+
+! BUG(wjs, 2020-09-25, ESCOMP/CTSM#43) Dynamic lakes don't work when methane is active,
+! so for now disable methane for this test.
+use_lch4 = .false.

src/biogeochem/ch4Mod.F90

@@ -1251,7 +1251,7 @@ subroutine DynamicColumnAdjustments(this, bounds, clump_index, column_state_upda
     character(len=*), parameter :: subname = 'DynamicColumnAdjustments'
     !-----------------------------------------------------------------------
 
-    ! BUG(wjs, 2016-02-16, bugz 2283) Need to do some special handling of finundated for
+    ! BUG(wjs, 2016-02-16, ESCOMP/CTSM#43) Need to do some special handling of finundated for
     ! increases in lake area, since lakes are assumed to be 100% inundated. Probably it's
     ! most appropriate for this special handling to happen elsewhere - i.e., within this
     ! routine, we do the standard adjustments as they are currently done, but then in the

src/biogeophys/BalanceCheckMod.F90

@@ -20,6 +20,7 @@ module BalanceCheckMod
   use SoilHydrologyType  , only : soilhydrology_type
   use SurfaceAlbedoType  , only : surfalb_type
   use WaterStateType     , only : waterstate_type
+  use LakestateType      , only : lakestate_type
   use WaterDiagnosticBulkType, only : waterdiagnosticbulk_type
   use WaterDiagnosticType, only : waterdiagnostic_type
   use Wateratm2lndType   , only : wateratm2lnd_type
@@ -122,7 +123,7 @@ end function GetBalanceCheckSkipSteps
   !-----------------------------------------------------------------------
   subroutine BeginWaterBalance(bounds, &
        num_nolakec, filter_nolakec, num_lakec, filter_lakec, &
-       water_inst, soilhydrology_inst, &
+       water_inst, soilhydrology_inst, lakestate_inst, &
        use_aquifer_layer)
     !
     ! !DESCRIPTION:
@@ -136,6 +137,7 @@ subroutine BeginWaterBalance(bounds, &
     integer                   , intent(in)    :: num_lakec            ! number of column lake points in column filter
     integer                   , intent(in)    :: filter_lakec(:)      ! column filter for lake points
     type(water_type)          , intent(inout) :: water_inst
+    type(lakestate_type)      , intent(in)    :: lakestate_inst
     type(soilhydrology_type)  , intent(in)    :: soilhydrology_inst
     logical                   , intent(in)    :: use_aquifer_layer    ! whether an aquifer layer is used in this run
     !
@@ -150,6 +152,7 @@ subroutine BeginWaterBalance(bounds, &
             num_nolakec, filter_nolakec, &
             num_lakec, filter_lakec, &
             soilhydrology_inst, &
+            lakestate_inst, &
             water_inst%bulk_and_tracers(i)%waterstate_inst, &
             water_inst%bulk_and_tracers(i)%waterdiagnostic_inst, &
             water_inst%bulk_and_tracers(i)%waterbalance_inst, &
@@ -161,7 +164,8 @@ end subroutine BeginWaterBalance
   !-----------------------------------------------------------------------
   subroutine BeginWaterBalanceSingle(bounds, &
        num_nolakec, filter_nolakec, num_lakec, filter_lakec, &
-       soilhydrology_inst, waterstate_inst, waterdiagnostic_inst, waterbalance_inst, &
+       soilhydrology_inst, lakestate_inst, waterstate_inst, & 
+       waterdiagnostic_inst, waterbalance_inst, &
        use_aquifer_layer)
     !
     ! !DESCRIPTION:
@@ -175,6 +179,7 @@ subroutine BeginWaterBalanceSingle(bounds, &
     integer                   , intent(in)    :: num_lakec            ! number of column lake points in column filter
     integer                   , intent(in)    :: filter_lakec(:)      ! column filter for lake points
     type(soilhydrology_type)  , intent(in)    :: soilhydrology_inst
+    type(lakestate_type)      , intent(in)    :: lakestate_inst
     class(waterstate_type)    , intent(inout) :: waterstate_inst
     class(waterdiagnostic_type), intent(in)   :: waterdiagnostic_inst
     class(waterbalance_type)  , intent(inout) :: waterbalance_inst
@@ -210,8 +215,8 @@ subroutine BeginWaterBalanceSingle(bounds, &
          water_mass = begwb(bounds%begc:bounds%endc))
 
     call ComputeWaterMassLake(bounds, num_lakec, filter_lakec, &
-         waterstate_inst, &
-         subtract_dynbal_baselines = .false., &
+         waterstate_inst, lakestate_inst, &
+         add_lake_water_and_subtract_dynbal_baselines = .false., &
          water_mass = begwb(bounds%begc:bounds%endc))
 
     call waterstate_inst%CalculateTotalH2osno(bounds, num_nolakec, filter_nolakec, &

src/biogeophys/CMakeLists.txt

@@ -12,6 +12,7 @@ list(APPEND clm_sources
   InfiltrationExcessRunoffMod.F90
   IrrigationMod.F90
   LakeCon.F90
+  LakeStateType.F90
   QSatMod.F90
   RootBiophysMod.F90
   SaturatedExcessRunoffMod.F90

src/biogeophys/LakeHydrologyMod.F90

@@ -164,7 +164,7 @@ subroutine LakeHydrology(bounds, &
          t_soisno             =>  temperature_inst%t_soisno_col         , & ! Output: [real(r8) (:,:) ]  snow temperature (Kelvin)             
          dTdz_top             =>  temperature_inst%dTdz_top_col         , & ! Output: [real(r8) (:)   ]  temperature gradient in top layer K m-1] !TOD 
          snot_top             =>  temperature_inst%snot_top_col         , & ! Output: [real(r8) (:)   ]  snow temperature in top layer [K]  !TODO
-         t_sno_mul_mss        => temperature_inst%t_sno_mul_mss_col     , & ! Output: [real(r8) (:)   ]  col snow temperature multiplied by layer mass, layer sum (K * kg/m2) 
+         t_sno_mul_mss        =>  temperature_inst%t_sno_mul_mss_col     , & ! Output: [real(r8) (:)   ]  col snow temperature multiplied by layer mass, layer sum (K * kg/m2) 
 
          begwb                =>  b_waterbalance_inst%begwb_col             , & ! Input:  [real(r8) (:)   ]  water mass begining of the time step    
          endwb                =>  b_waterbalance_inst%endwb_col             , & ! Output: [real(r8) (:)   ]  water mass end of the time step         
@@ -649,8 +649,8 @@ subroutine LakeHydrology(bounds, &
     ! Determine ending water balance and volumetric soil water
 
     call ComputeWaterMassLake(bounds, num_lakec, filter_lakec, &
-         b_waterstate_inst, &
-         subtract_dynbal_baselines = .false., &
+         b_waterstate_inst, lakestate_inst, &
+         add_lake_water_and_subtract_dynbal_baselines = .false., &
          water_mass = endwb(bounds%begc:bounds%endc))
 
     do j = 1, nlevgrnd

src/biogeophys/LakeTemperatureMod.F90

@@ -134,6 +134,7 @@ subroutine LakeTemperature(bounds, num_lakec, filter_lakec, num_lakep, filter_la
     type(energyflux_type)  , intent(inout) :: energyflux_inst
     type(temperature_type) , intent(inout) :: temperature_inst
     type(lakestate_type)   , intent(inout) :: lakestate_inst
+
     !
     ! !LOCAL VARIABLES:
     real(r8), parameter :: p0 = 1._r8                                      ! neutral value of turbulent prandtl number
@@ -247,7 +248,7 @@ subroutine LakeTemperature(bounds, num_lakec, filter_lakec, num_lakep, filter_la
          t_grnd          =>   temperature_inst%t_grnd_col          , & ! Input:  [real(r8) (:)   ]  ground temperature (Kelvin)             
          t_soisno        =>   temperature_inst%t_soisno_col        , & ! Output: [real(r8) (:,:) ]  soil (or snow) temperature (Kelvin)   
          t_lake          =>   temperature_inst%t_lake_col          , & ! Output: [real(r8) (:,:) ]  col lake temperature (Kelvin)             
-
+        
          beta            =>   lakestate_inst%betaprime_col         , & ! Output: [real(r8) (:)   ]  col effective beta: sabg_lyr(p,jtop) for snow layers, beta otherwise
          lake_icefrac    =>   lakestate_inst%lake_icefrac_col      , & ! Output: [real(r8) (:,:) ]  col mass fraction of lake layer that is frozen
          lake_icefracsurf => lakestate_inst%lake_icefracsurf_col   , & ! Output: [real(r8) (:,:) ]  col mass fraction of surface lake layer that is frozen
@@ -999,8 +1000,10 @@ subroutine LakeTemperature(bounds, num_lakec, filter_lakec, num_lakep, filter_la
           ncvts(c) = ncvts(c) + cv_lake(c,j)*(t_lake(c,j)-tfrz) &
                    + cfus*dz_lake(c,j)*(1._r8-lake_icefrac(c,j)) 
           fin(c) = fin(c) + phi(c,j)
+
        end do
     end do
+
 
     call waterstatebulk_inst%CalculateTotalH2osno(bounds, num_lakec, filter_lakec, &
          caller = 'LakeTemperature-2', &
@@ -1022,6 +1025,7 @@ subroutine LakeTemperature(bounds, num_lakec, filter_lakec, num_lakep, filter_la
     end do
 
 
+
     ! Check energy conservation.
     do fp = 1, num_lakep
        p = filter_lakep(fp)

src/biogeophys/TemperatureType.F90

@@ -116,7 +116,7 @@ module TemperatureType
      real(r8), pointer    :: xmf_h2osfc_col        (:)   ! latent heat of phase change of surface water
      real(r8), pointer    :: fact_col              (:,:) ! used in computing tridiagonal matrix
      real(r8), pointer    :: c_h2osfc_col          (:)   ! heat capacity of surface water
-
+     
    contains
 
      procedure, public  :: Init
@@ -620,7 +620,6 @@ subroutine InitHistory(this, bounds, is_simple_buildtemp, is_prog_buildtemp )
             ptr_patch=this%t_veg10_night_patch, default='inactive')
     endif
 
-
   end subroutine InitHistory
 
   !-----------------------------------------------------------------------