ESCOMP / CTSM

Community Terrestrial Systems Model (includes the Community Land Model of CESM)
http://www.cesm.ucar.edu/models/cesm2.0/land/
Other
309 stars 313 forks source link

SnowCapping should probably use updated h2osno #736

Closed billsacks closed 5 years ago

billsacks commented 5 years ago

SnowCapping decides whether to remove excess snow based on h2osno. However, at the point where this is called, it looks like h2osno is out of sync with (h2osoi_liq + h2osoi_ice): Those variables can be changed slightly in SnowWater (which is called before SnowCapping), yet h2osno isn't recalculated based on those changes.

I don't think this causes a big problem in practice, but it seems like it would improve code understandability if we had the general rule that any code that references h2osno is using a version of h2osno that is in sync with (h2osoi_liq + h2osoi_ice). So in this particular case, I propose updating h2osno at the end of subroutine SnowWater.

I have confirmed that doing so will change answers, by running a one-year test with and without the following change (SMS_Ly1.f09_g17.I2000Clm50SpGs.cheyenne_intel.clm-monthly):

diff --git a/src/biogeophys/SnowHydrologyMod.F90 b/src/biogeophys/SnowHydrologyMod.F90
index 45edb7bb..56f80366 100644
--- a/src/biogeophys/SnowHydrologyMod.F90
+++ b/src/biogeophys/SnowHydrologyMod.F90
@@ -1887,6 +1887,8 @@ subroutine SnowCapping(bounds, num_initc, filter_initc, num_snowc, filter_snowc,
     real(r8)   :: frac_adjust                      ! fraction of mass remaining after capping
     real(r8)   :: rho                              ! partial density of ice (not scaled with frac_sno) [kg/m3]
     integer    :: fc, c                            ! counters
+    integer    :: j
+    real(r8)   :: h2osno_recalculated(bounds%begc:bounds%endc)
     real(r8)   :: h2osno_excess(bounds%begc:bounds%endc) ! excess snow that needs to be capped [mm H2O]
     logical    :: apply_runoff(bounds%begc:bounds%endc)  ! for columns with capping, whether the capping flux should be sent to runoff
     ! Always keep at least this fraction of the bottom snow layer when doing snow capping
@@ -1926,8 +1928,22 @@ subroutine SnowCapping(bounds, num_initc, filter_initc, num_snowc, filter_snowc,
        qflx_snwcp_discarded_liq(c) = 0.0_r8
     end do

+    do fc = 1, num_snowc
+       c = filter_snowc(fc)
+       h2osno_recalculated(c) = 0._r8
+    end do
+
+    do j = -nlevsno+1,0
+       do fc = 1, num_snowc
+          c = filter_snowc(fc)
+          if (j >= col%snl(c)+1) then
+             h2osno_recalculated(c) = h2osno_recalculated(c) + h2osoi_ice(c,j) + h2osoi_liq(c,j)
+          end if
+       end do
+    end do
+
     call SnowCappingExcess(bounds, num_snowc, filter_snowc, &
-         h2osno = h2osno(bounds%begc:bounds%endc), &
+         h2osno = h2osno_recalculated(bounds%begc:bounds%endc), &
          topo = topo(bounds%begc:bounds%endc), &
          h2osno_excess = h2osno_excess(bounds%begc:bounds%endc), &
          apply_runoff = apply_runoff(bounds%begc:bounds%endc))

@swensosc and @lvankampenhout can you please let me know if you see any problems with doing this extra update of h2osno before calling SnowCapping?

billsacks commented 5 years ago

@swensosc supports doing the extra update that I mention here. In fact, based on his input, I'll probably go a step further, and ensure that any routine that uses h2osno_total updates this value locally before using it.

lvankampenhout commented 5 years ago

Sorry for being late to the party, but I agree with these changes and with the general overhaul of the h2osno variable in ctsm1.0.dev046. I've always found the double bookkeeping unnecessary and confusing