Small answer changes in preparation for adding option for bio-mass heat storage

bld/unit_testers/build-namelist_test.pl

@@ -138,7 +138,7 @@ sub make_config_cache {
 #
 # Figure out number of tests that will run
 #
-my $ntests = 1516;
+my $ntests = 1513;
 if ( defined($opts{'compare'}) ) {
    $ntests += 1017;
 }

cime_config/config_compsets.xml

@@ -36,6 +36,11 @@
 
 <!-- I single point forcing -->
 
+  <compset>
+    <alias>I1PtClm51SpRs</alias>
+    <lname>2000_DATM%1PT_CLM51%SP_SICE_SOCN_SROF_SGLC_SWAV</lname>
+  </compset>
+
   <compset>
     <alias>I1PtClm50SpRs</alias>
     <lname>2000_DATM%1PT_CLM50%SP_SICE_SOCN_SROF_SGLC_SWAV</lname>
@@ -59,6 +64,17 @@
     <lname>2000_DATM%GSWP3v1_CLM50%SP_SICE_SOCN_SROF_SGLC_SWAV</lname>
   </compset>
 
+  <compset>
+    <alias>I2000Clm51Sp</alias>
+    <lname>2000_DATM%GSWP3v1_CLM51%SP_SICE_SOCN_MOSART_SGLC_SWAV</lname>
+  </compset>
+
+  <!-- Primarily for testing -->
+  <compset>
+    <alias>I2000Clm51SpRs</alias>
+    <lname>2000_DATM%GSWP3v1_CLM51%SP_SICE_SOCN_SROF_SGLC_SWAV</lname>
+  </compset>
+
   <!-- Primarily for testing -->
   <compset>
     <alias>I2000Clm50SpRtm</alias>
@@ -128,6 +144,11 @@
     <science_support grid="f19_g17"/>
   </compset>
 
+  <compset>
+    <alias>I1850Clm51Sp</alias>
+    <lname>1850_DATM%GSWP3v1_CLM51%SP_SICE_SOCN_MOSART_SGLC_SWAV</lname>
+  </compset>
+
   <compset>
     <alias>I1850Clm51Bgc</alias>
     <lname>1850_DATM%GSWP3v1_CLM51%BGC_SICE_SOCN_MOSART_SGLC_SWAV</lname>
@@ -177,6 +198,10 @@
   </compset>
 
   <!---I FATES compsets -->
+  <compset>
+    <alias>I2000Clm51Fates</alias>
+    <lname>2000_DATM%GSWP3v1_CLM51%FATES_SICE_SOCN_MOSART_SGLC_SWAV</lname>
+  </compset>
   <compset>
     <alias>I2000Clm50Fates</alias>
     <lname>2000_DATM%GSWP3v1_CLM50%FATES_SICE_SOCN_MOSART_SGLC_SWAV</lname>
@@ -213,6 +238,16 @@
     <science_support grid="f19_g17"/>
   </compset>
 
+  <compset>
+     <alias>I1850Clm51SpNoAnthro</alias>
+     <lname>1850_DATM%GSWP3v1_CLM51%SP-NOANTHRO_SICE_SOCN_MOSART_SGLC_SWAV</lname>
+   </compset>
+
+  <compset>
+    <alias>IHistClm51Sp</alias>
+    <lname>HIST_DATM%GSWP3v1_CLM51%SP_SICE_SOCN_MOSART_SGLC_SWAV</lname>
+  </compset>
+
   <compset>
     <alias>IHistClm51Bgc</alias>
     <lname>HIST_DATM%GSWP3v1_CLM51%BGC_SICE_SOCN_MOSART_SGLC_SWAV</lname>

cime_config/testdefs/testlist_clm.xml

@@ -381,12 +381,13 @@
       <option name="wallclock">00:20:00</option>
     </options>
   </test>
-  <test name="ERP_D_Ld5" grid="f10_f10_musgs" compset="I2000Clm50Sp" testmods="clm/decStart">
+  <test name="ERP_D_Ld5" grid="f10_f10_musgs" compset="I2000Clm51Sp" testmods="clm/decStart">
     <machines>
       <machine name="cheyenne" compiler="intel" category="aux_clm"/>
     </machines>
     <options>
       <option name="wallclock">00:20:00</option>
+      <option name="comment">2000 Sp test for CLM51</option>
     </options>
   </test>
   <test name="ERP_D_Ld5" grid="f10_f10_musgs" compset="I2000Clm50Sp" testmods="clm/reduceOutput">
@@ -423,6 +424,15 @@
       <option name="wallclock">00:20:00</option>
     </options>
   </test>
+  <test name="ERP_D_Ld5" grid="f10_f10_musgs" compset="IHistClm51Sp" testmods="clm/default">
+    <machines>
+      <machine name="cheyenne" compiler="intel" category="aux_clm"/>
+    </machines>
+    <options>
+      <option name="wallclock">00:20:00</option>
+      <option name="comment"  >Test Hist compset with Sp for CLM5.1</option>
+    </options>
+  </test>
   <test name="SMS_Ld5" grid="f09_g17" compset="IHistClm50SpCru" testmods="clm/default">
     <machines>
       <machine name="cheyenne" compiler="intel" category="ctsm_sci"/>

doc/ChangeLog

@@ -1,4 +1,112 @@
 ===============================================================
+Tag name:  ctsm5.1.dev020
+Originator(s):  erik (Erik Kluzek,UCAR/TSS,303-497-1326)
+Date: Wed Dec 30 00:42:16 MST 2020
+One-line Summary: Potential roundoff changes in preparation for bio-mass heat storage option
+
+Purpose of changes
+------------------
+
+This adds in some changes that may be roundoff different in preparation for bringing in the changes
+that allow turning on the bio-mass heat storage option. There's a few zeroed terms added to some calculations
+in CanopyFluxes. The new terms are identical to what will come in for BHS except the zeroed terms are
+not arrays here but scalars. I showed the new terms are only off by roundoff with a previous commit that added
+explicit tests for the terms and verified it for two tests.
+
+
+Bugs fixed or introduced
+------------------------
+
+Significant changes to scientifically-supported configurations
+--------------------------------------------------------------
+
+Does this tag change answers significantly for any of the following physics configurations?
+(Details of any changes will be given in the "Answer changes" section below.)
+
+    [Put an [X] in the box for any configuration with significant answer changes.]
+
+[ ] clm5_1
+
+[ ] clm5_0
+
+[ ] ctsm5_0-nwp
+
+[ ] clm4_5
+
+Notes of particular relevance for users
+---------------------------------------
+
+Changes to CTSM's user interface (e.g., new/renamed XML or namelist variables): Added new Clm51 compsets
+     I1PtClm51SpRs
+     I2000Clm51Sp
+     I2000Clm51SpRs
+     I1850Clm51Sp
+     I2000Clm51Fates
+     I1850Clm51SpNoAnthro
+     IHistClm51Sp
+
+Notes of particular relevance for developers: (including Code reviews and testing)
+---------------------------------------------
+NOTE: Be sure to review the steps in README.CHECKLIST.master_tags as well as the coding style in the Developers Guide
+
+Caveats for developers (e.g., code that is duplicated that requires double maintenance):
+   Here I set some terms to zero that will be non-zero even when BHS is off
+   Earlier commits show the smaller set of changes that are needed.
+
+Changes to tests or testing: Switched one Clm50Sp test to Clm51Sp and added one HistClm51Sp test
+
+CTSM testing: regular
+
+ [PASS means all tests PASS and OK means tests PASS other than expected fails.]
+
+  build-namelist tests:
+
+    cheyenne - PASS
+
+  python testing (see instructions in python/README.md; document testing done):
+
+    cheyenne - PASS
+
+  regular tests (aux_clm):
+
+    cheyenne ---- OK
+    izumi ------- OK
+
+If the tag used for baseline comparisons was NOT the previous tag, note that here:
+
+
+Answer changes
+--------------
+
+Changes answers relative to baseline:
+
+  Summarize any changes to answers, i.e.,
+    - what code configurations: all
+    - what platforms/compilers: izumi-PGI
+    - nature of change: roundoff
+
+   If bitwise differences were observed, how did you show they were no worse
+   than roundoff?
+
+    Earlier commit tested the new terms versus the old and showed they were
+    order e-11 for absolute difference and e-14 relative difference.
+
+   Tests that change from baseline....
+
+   ERS_Ly5_P144x1.f10_f10_musgs.IHistClm51BgcCrop.cheyenne_intel.clm-cropMonthOutput
+   SMS.f10_f10_musgs.I2000Clm50BgcCrop.izumi_pgi.clm-crop
+   SMS_D.f10_f10_musgs.I2000Clm51BgcCrop.izumi_pgi.clm-crop
+
+Detailed list of changes
+------------------------
+
+Pull Requests that document the changes (include PR ids): #1241
+(https://github.com/ESCOMP/ctsm/pull)
+
+  #1241 -- Small answer changes in preparation for adding option for bio-mass heat storage
+
+===============================================================
+===============================================================
 Tag name:  ctsm5.1.dev019
 Originator(s):  sacks (Bill Sacks)
 Date:  Sat Dec 19 06:55:46 MST 2020

doc/ChangeSum

@@ -1,5 +1,6 @@
 Tag                   Who      Date  Summary
 ============================================================================================================================
+  ctsm5.1.dev020     erik 12/30/2020 Potential roundoff changes in preparation for bio-mass heat storage option
   ctsm5.1.dev019    sacks 12/19/2020 Fix ndep from coupler
   ctsm5.1.dev018   slevis 12/08/2020 Add ACTIVE (T/F) column to master hist fields table and alphabetize
   ctsm5.1.dev017   slevis 11/17/2020 Write history fields master list to separate optional file

src/biogeophys/CanopyFluxesMod.F90

@@ -402,6 +402,14 @@ subroutine CanopyFluxes(bounds,  num_exposedvegp, filter_exposedvegp,
     real(r8) :: dt_veg_temp(bounds%begp:bounds%endp)
     integer  :: iv
     logical  :: is_end_day                               ! is end of current day
+    real(r8) :: tl_ini                                   ! leaf temperature from beginning of time step [K]
+    real(r8) :: ts_ini                                   ! stem temperature from beginning of time step [K]
+    real(r8) :: cp_leaf                                  !heat capacity of leaves
+    real(r8) :: dt_stem                                  !change in stem temperature
+    real(r8) :: frac_rad_abs_by_stem                     !fraction of incoming radiation absorbed by stems
+    real(r8) :: lw_stem                                  !internal longwave stem
+    real(r8) :: lw_leaf                                  !internal longwave leaf
+    real(r8) :: sa_internal                              !min(sa_stem,sa_leaf)
 
     integer :: dummy_to_make_pgi_happy
     !------------------------------------------------------------------------------
@@ -627,6 +635,14 @@ subroutine CanopyFluxes(bounds,  num_exposedvegp, filter_exposedvegp,
          obuold(p) = 0._r8
          btran(p)  = btran0
       end do
+      frac_rad_abs_by_stem = 0._r8
+      lw_leaf              = 0._r8
+      cp_leaf              = 0._r8
+      lw_stem              = 0._r8
+      sa_internal          = 0._r8
+      dt_stem              = 0._r8
+      tl_ini               = 0._r8
+      ts_ini               = 0._r8
 
       ! calculate daylength control for Vcmax
       do f = 1, fn
@@ -1046,20 +1062,28 @@ subroutine CanopyFluxes(bounds,  num_exposedvegp, filter_exposedvegp,
                  +(1._r8-frac_sno(c)-frac_h2osfc(c))*t_soisno(c,1)**4 &
                  +frac_h2osfc(c)*t_h2osfc(c)**4)
 
-            dt_veg(p) = (sabv(p) + air(p) + bir(p)*t_veg(p)**4 + &
-                 cir(p)*lw_grnd - efsh - efe(p)) / &
-                 (- 4._r8*bir(p)*t_veg(p)**3 +dc1*wtga +dc2*wtgaq*qsatldT(p))
+            dt_veg(p) = ((1._r8-frac_rad_abs_by_stem)*(sabv(p) + air(p) &
+                  + bir(p)*t_veg(p)**4 + cir(p)*lw_grnd) &
+                  - efsh - efe(p) - lw_leaf + lw_stem &
+                  - (cp_leaf/dtime)*(t_veg(p) - tl_ini)) &
+                  / ((1._r8-frac_rad_abs_by_stem)*(- 4._r8*bir(p)*t_veg(p)**3 &
+                  + 4._r8*sa_internal*emv(p)*sb*t_veg(p)**3 &
+                  +dc1*wtga+dc2*wtgaq*qsatldT(p))+ cp_leaf/dtime)
             t_veg(p) = tlbef(p) + dt_veg(p)
             dels = dt_veg(p)
             del(p)  = abs(dels)
             err(p) = 0._r8
             if (del(p) > delmax) then
                dt_veg(p) = delmax*dels/del(p)
                t_veg(p) = tlbef(p) + dt_veg(p)
-               err(p) = sabv(p) + air(p) + bir(p)*tlbef(p)**3*(tlbef(p) + &
-                    4._r8*dt_veg(p)) + cir(p)*lw_grnd - &
-                    (efsh + dc1*wtga*dt_veg(p)) - (efe(p) + &
-                    dc2*wtgaq*qsatldT(p)*dt_veg(p))
+               err(p) = (1._r8-frac_rad_abs_by_stem)*(sabv(p) + air(p) &
+                    + bir(p)*tlbef(p)**3*(tlbef(p) + &
+                    4._r8*dt_veg(p)) + cir(p)*lw_grnd) &
+                    -sa_internal*emv(p)*sb*tlbef(p)**3*(tlbef(p) + 4._r8*dt_veg(p)) &
+                    + lw_stem &
+                    - (efsh + dc1*wtga*dt_veg(p)) - (efe(p) + &
+                    dc2*wtgaq*qsatldT(p)*dt_veg(p)) &
+                    - (cp_leaf/dtime)*(t_veg(p) - tl_ini)
             end if
 
             ! Fluxes from leaves to canopy space
@@ -1183,9 +1207,10 @@ subroutine CanopyFluxes(bounds,  num_exposedvegp, filter_exposedvegp,
               +(1._r8-frac_sno(c)-frac_h2osfc(c))*t_soisno(c,1)**4 &
               +frac_h2osfc(c)*t_h2osfc(c)**4)
 
-         err(p) = sabv(p) + air(p) + bir(p)*tlbef(p)**3*(tlbef(p) + 4._r8*dt_veg(p)) &
-                                !+ cir(p)*t_grnd(c)**4 - eflx_sh_veg(p) - hvap*qflx_evap_veg(p)
-              + cir(p)*lw_grnd - eflx_sh_veg(p) - hvap*qflx_evap_veg(p)
+         err(p) = (1.0_r8-frac_rad_abs_by_stem)*(sabv(p) + air(p) + bir(p)*tlbef(p)**3 &
+              *(tlbef(p) + 4._r8*dt_veg(p)) + cir(p)*lw_grnd) &
+                - lw_leaf + lw_stem - eflx_sh_veg(p) - hvap*qflx_evap_veg(p) &
+                - ((t_veg(p)-tl_ini)*cp_leaf/dtime)
 
          ! Fluxes from ground to canopy space
 
@@ -1269,14 +1294,20 @@ subroutine CanopyFluxes(bounds,  num_exposedvegp, filter_exposedvegp,
 
          ! Downward longwave radiation below the canopy
 
-         dlrad(p) = (1._r8-emv(p))*emg(c)*forc_lwrad(c) + &
-              emv(p)*emg(c)*sb*tlbef(p)**3*(tlbef(p) + 4._r8*dt_veg(p))
+         dlrad(p) = (1._r8-emv(p))*emg(c)*forc_lwrad(c) &
+              + emv(p)*emg(c)*sb*tlbef(p)**3*(tlbef(p) + 4._r8*dt_veg(p)) &
+              *(1.0_r8-frac_rad_abs_by_stem) &
+              + emv(p)*emg(c)*sb*ts_ini**3*(ts_ini + 4._r8*dt_stem) &
+              *frac_rad_abs_by_stem
 
          ! Upward longwave radiation above the canopy
 
          ulrad(p) = ((1._r8-emg(c))*(1._r8-emv(p))*(1._r8-emv(p))*forc_lwrad(c) &
-              + emv(p)*(1._r8+(1._r8-emg(c))*(1._r8-emv(p)))*sb*tlbef(p)**3*(tlbef(p) + &
-              4._r8*dt_veg(p)) + emg(c)*(1._r8-emv(p))*sb*lw_grnd)
+              + emv(p)*(1._r8+(1._r8-emg(c))*(1._r8-emv(p)))*sb &
+              *tlbef(p)**3*(tlbef(p) + 4._r8*dt_veg(p))*(1._r8-frac_rad_abs_by_stem) &
+              + emv(p)*(1._r8+(1._r8-emg(c))*(1._r8-emv(p)))*sb &
+              *ts_ini**3*(ts_ini+ 4._r8*dt_stem)*frac_rad_abs_by_stem &
+              + emg(c)*(1._r8-emv(p))*sb*lw_grnd)
 
          ! Calculate the skin temperature as a weighted sum of all the ground and vegetated fraction
          ! The weight is the so-called vegetation emissivity, but not that emv is actually an attentuation