User defined top snow layer #792
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -24,6 +24,7 @@ module SnowHydrologyMod | |
| use clm_varpar , only : nlevsno, nlevsoi, nlevgrnd | ||
| use clm_varctl , only : iulog, use_subgrid_fluxes | ||
| use clm_varcon , only : namec, h2osno_max, hfus, denice, rpi, spval, tfrz | ||
| use clm_varcon, only: dzmin, dzmax_l, dzmax_u | ||
| use atm2lndType , only : atm2lnd_type | ||
| use AerosolMod , only : aerosol_type | ||
| use TemperatureType , only : temperature_type | ||
|
|
@@ -100,20 +101,6 @@ module SnowHydrologyMod | |
| integer, parameter, public :: LoTmpDnsSlater2017 = 2 ! For temperature below -15C use equation from Slater 2017 | ||
| integer, parameter, public :: LoTmpDnsTruncatedAnderson1976 = 1 ! Truncate low temp. snow density from the Anderson-1976 version at -15C | ||
|
|
||
| ! | ||
| ! !PRIVATE DATA MEMBERS: | ||
|
|
||
|
|
@@ -123,6 +110,7 @@ module SnowHydrologyMod | |
| ! If true, the density of new snow depends on wind speed, and there is also | ||
| ! wind-dependent snow compaction | ||
| logical :: wind_dependent_snow_density ! If snow density depends on wind or not | ||
| real(r8) :: dzmin_1, dzmax_l_1, dzmax_u_1 ! namelist-defined top snow layer information | ||
| integer :: overburden_compaction_method = -1 | ||
| integer :: new_snow_density = LoTmpDnsSlater2017 ! Snow density type | ||
| real(r8) :: upplim_destruct_metamorph = 100.0_r8 ! Upper Limit on Destructive Metamorphism Compaction [kg/m3] | ||
|
|
@@ -187,11 +175,15 @@ subroutine SnowHydrology_readnl( NLFilename) | |
| wind_dependent_snow_density, snow_overburden_compaction_method, & | ||
| lotmp_snowdensity_method, upplim_destruct_metamorph, & | ||
| overburden_compress_Tfactor, min_wind_snowcompact, & | ||
| reset_snow, reset_snow_glc, reset_snow_glc_ela, dzmin_1, dzmax_l_1, & | ||
| dzmax_u_1 | ||
|
|
||
| ! Initialize options to default values, in case they are not specified in the namelist | ||
| wind_dependent_snow_density = .false. | ||
| snow_overburden_compaction_method = ' ' | ||
| dzmin_1 = 0.01_r8 | ||
billsacks marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| dzmax_l_1 = 0.03_r8 | ||
| dzmax_u_1 = 0.02_r8 | ||
|
|
||
| if (masterproc) then | ||
| unitn = getavu() | ||
|
|
@@ -218,6 +210,9 @@ subroutine SnowHydrology_readnl( NLFilename) | |
| call shr_mpi_bcast (reset_snow , mpicom) | ||
| call shr_mpi_bcast (reset_snow_glc , mpicom) | ||
| call shr_mpi_bcast (reset_snow_glc_ela , mpicom) | ||
| call shr_mpi_bcast (dzmin_1, mpicom) | ||
| call shr_mpi_bcast (dzmax_l_1, mpicom) | ||
| call shr_mpi_bcast (dzmax_u_1, mpicom) | ||
|
|
||
| if (masterproc) then | ||
| write(iulog,*) ' ' | ||
|
|
@@ -243,6 +238,19 @@ subroutine SnowHydrology_readnl( NLFilename) | |
| errMsg(sourcefile, __LINE__)) | ||
| end if | ||
|
|
||
| if (dzmin_1 < 0.01_r8) then | ||
billsacks marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| call endrun(msg="ERROR dzmin_1 cannot be less than 0.01"// & | ||
| errMsg(sourcefile, __LINE__)) | ||
| end if | ||
| if (dzmax_l_1 < 0.03_r8) then | ||
| call endrun(msg="ERROR dzmax_l_1 cannot be less than 0.03"// & | ||
| errMsg(sourcefile, __LINE__)) | ||
| end if | ||
| if (dzmax_u_1 < 0.02_r8) then | ||
| call endrun(msg="ERROR dzmax_2_1 cannot be less than 0.02"// & | ||
| errMsg(sourcefile, __LINE__)) | ||
| end if | ||
|
|
||
| end subroutine SnowHydrology_readnl | ||
|
|
||
| !----------------------------------------------------------------------- | ||
|
|
@@ -800,12 +808,12 @@ subroutine BuildFilter_SnowpackInitialized(bounds, num_c, filter_c, & | |
| ! maintain answers as before. | ||
| snowpack_initialized(c) = ( & | ||
| snl(c) == 0 .and. & | ||
|
There was a problem hiding this comment. @billsacks I think the next four 0.01 values are what you wanted me to look for, right? More generally speaking, I think this code is ready for review. Also I would like to know your preference in this case between new Unit Test and new standard test-suite test. I will be working on that next. There was a problem hiding this comment. Yes, these are the four places where I know this value is hard-coded. If you haven't already, thenI also thought it would be worth grepping the code for As for testing: With some upcoming suggestions that I'm going to make (in which you don't have hard-coded defaults for these values, but rather have them always explicitly set in the namelist), your new code will be sufficiently covered by existing system tests, so no need to add additional tests. There was a problem hiding this comment. I looked for all the instances of 0.01 at the time, and there were some more that looked unrelated, so I left them alone. |
||
| frac_sno_eff(c)*snow_depth(c) >= (dzmin(1) + lsadz) .and. & | ||
| qflx_snow_grnd(c) > 0.0_r8) | ||
| else | ||
| snowpack_initialized(c) = ( & | ||
| snl(c) == 0 .and. & | ||
| frac_sno_eff(c)*snow_depth(c) >= dzmin(1)) | ||
| end if | ||
|
|
||
| end do | ||
|
|
@@ -1496,7 +1504,7 @@ subroutine CombineSnowLayers(bounds, num_snowc, filter_snowc, & | |
| real(r8):: h2osno_total(bounds%begc:bounds%endc) ! total snow water (mm H2O) | ||
| real(r8):: zwice(bounds%begc:bounds%endc) ! total ice mass in snow | ||
| real(r8):: zwliq (bounds%begc:bounds%endc) ! total liquid water in snow | ||
| real(r8):: dzminloc(nlevsno) ! minimum of top snow layer (local) | ||
| real(r8):: dtime !land model time step (sec) | ||
|
|
||
| !----------------------------------------------------------------------- | ||
|
|
@@ -1680,8 +1688,8 @@ subroutine CombineSnowLayers(bounds, num_snowc, filter_snowc, & | |
| c = filter_snowc(fc) | ||
| l = col%landunit(c) | ||
| if (snow_depth(c) > 0._r8) then | ||
| if ((ltype(l) == istdlak .and. snow_depth(c) < dzmin(1) + lsadz ) .or. & | ||
| ((ltype(l) /= istdlak) .and. ((frac_sno_eff(c)*snow_depth(c) < dzmin(1)) & | ||
| .or. (h2osno_total(c)/(frac_sno_eff(c)*snow_depth(c)) < 50._r8)))) then | ||
|
|
||
| snl(c) = 0 | ||
|
|
@@ -2234,6 +2242,26 @@ subroutine InitSnowLayers (bounds, snow_depth) | |
| zi => col%zi & ! Output: [real(r8) (:,:) ] interface level below a "z" level (m) (-nlevsno+0:nlevgrnd) | ||
| ) | ||
|
|
||
| ! These three variables determine the vertical structure of the snow pack: | ||
| ! dzmin: minimum of top snow layer | ||
| ! dzmax_l: maximum thickness of layer when no layers beneath | ||
| ! dzmax_u: maximum thickness of layer when layers beneath | ||
| dzmin(1) = dzmin_1 ! set in the namelist with default or user-defined value | ||
| dzmax_l(1) = dzmax_l_1 ! same comment | ||
| dzmax_u(1) = dzmax_u_1 ! same comment | ||
| do j = 2, nlevsno | ||
| dzmin(j) = dzmax_u(j-1) * 0.5_r8 | ||
| if (dzmin(j) <= dzmin(j-1)) then | ||
| dzmin(j) = (dzmin(j-1) + dzmax_u(j-1)) * 0.5_r8 | ||
| end if | ||
| dzmax_u(j) = 2._r8 * dzmax_u(j-1) + 0.01_r8 | ||
| if (j == 2) then | ||
| dzmax_l(j) = dzmax_u(j) + dzmax_l(j-1) - 0.01_r8 | ||
| else | ||
| dzmax_l(j) = dzmax_u(j) + dzmax_l(j-1) | ||
| end if | ||
| end do | ||
|
|
||
billsacks marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| loop_columns: do c = bounds%begc,bounds%endc | ||
| l = col%landunit(c) | ||
|
|
||
|
|
||
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
In reference to this comment, the original code assigned huge(1._r8) to dzmax_l(12) and dzmax_u(12). So to match the original behavior, I think I need to do the following:
if (nlevsno >= 12) thendzmax_l = huge(1._r8)dzmax_u = huge(1._r8)end ifIs this what you were suggesting @billsacks ?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
More precisely...
if (j >= 12) thendzmax_l(j) = ...dzmax_u(j) = ...There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Changing mind:
I prefer the first if statement with these corrections:
dzmax_l(nlevsno) = ...dzmax_u(nlevsno) = ...This replicates the current behavior for nlevsno = 12. If some day we allow nlevsno > 12, this just sets the bottom layer to huge.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes, the last comment is what I was thinking (where the ... should be
huge(1._r8)).There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
i.e., the bottom layer should be set to huge, regardless of the number of snow layers. I now realize that that isn't exactly what the original code did, but I think that was the intent.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I guess it would be good to verify that that does the right thing. Either (a) make that change and verify that it doesn't change answers for, say, a 5-layer case (a CLM45 case), or (b) don't make that change, and ensure that having dzmax_l and dzmax_u NOT set to huge for nlevsno gives the same answers as setting it to huge for a 12-layer case. I prefer (a), but don't feel strongly about it.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I also prefer (a).
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The hobart test-suite includes a few clm45 cases and they all PASS.