Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Feature/snow layer creation threshold #520

Merged

Conversation

wknoben
Copy link
Collaborator

@wknoben wknoben commented Mar 7, 2023

Make sure all the relevant boxes are checked (and only check the box if you actually completed the step):

  • Closes #xxx (identify the issue associated with this PR)
  • Code passes standard test cases (results are either bit-for-bit identical, or differences are explained in the PR comment)
  • New tests added (describe which tests were performed to test the changes)
  • Science test figures (add figures to PR comment and describe the tests)
  • Checked that the new code conforms to the SUMMA coding conventions
  • Describe the change in the release notes (use either ./summa/docs/whats-new.md or ./summa/docs/minor-changes.md depending on what changed)

Test cases

Synthetic test cases complete bit-for-bit identical.

WRR test cases do not, because changes to the snow layering scheme translate into changes in simulated snow (and related) variables. See example below. Using the new layering scheme (bottom plot, dashed lines in both plots) can lead to lower, identical and higher snow depths than the original scheme gives us. Differences are in the order of centimeters over a multi-year simulation (second figure).

WIP: add an explanation about why the snow sims are slightly different.

summa_snowLayer_PR

summa_snowLayer_PR_depthDiff

@andywood
Copy link
Collaborator

andywood commented Mar 8, 2023 via email

@wknoben
Copy link
Collaborator Author

wknoben commented Mar 8, 2023

It's a combination of hard-wired logic and user-defined parameters in localParamInfo.txt.

Logic lines - layerDivide():

zmax_lower = (/zmaxLayer1_lower, zmaxLayer2_lower, zmaxLayer3_lower, zmaxLayer4_lower/)
zmax_upper = (/zmaxLayer1_upper, zmaxLayer2_upper, zmaxLayer3_upper, zmaxLayer4_upper/)
! initialize the number of snow layers
nSnow = indx_data%var(iLookINDEX%nSnow)%dat(1)
nSoil = indx_data%var(iLookINDEX%nSoil)%dat(1)
nLayers = indx_data%var(iLookINDEX%nLayers)%dat(1)
! ***** special case of no snow layers
if(nSnow==0)then
! check if create the first snow layer
select case(ix_snowLayers)
case(sameRulesAllLayers); createLayer = (scalarSnowDepth > zmax)
case(rulesDependLayerIndex); createLayer = (scalarSnowDepth > zmaxLayer1_lower)
case default; err=20; message=trim(message)//'unable to identify option to combine/sub-divide snow layers'; return
end select ! (option to combine/sub-divide snow layers)

[...]

else ! if nSnow>0
! identify the number of layers to check for need for sub-division
nCheck = min(nSnow, maxSnowLayers-1) ! the depth of the last layer, if it exists, does not have a maximum value
! loop through all layers, and sub-divide a given layer, if necessary
do iLayer=1,nCheck
divideLayer=.false.
! identify the maximum depth of the layer
select case(ix_snowLayers)
case(sameRulesAllLayers)
if (nCheck >= maxSnowLayers-1) then
! make sure we don't divide so make very big
zmaxCheck = veryBig
else
zmaxCheck = zmax
end if
case(rulesDependLayerIndex)
if(iLayer == nSnow)then
zmaxCheck = zmax_lower(iLayer)
else
zmaxCheck = zmax_upper(iLayer)
end if
case default; err=20; message=trim(message)//'unable to identify option to combine/sub-divide snow layers'; return
end select ! (option to combine/sub-divide snow layers)

**Logic lines - ** layerMerge():

! identify algorithmic control parameters to syb-divide and combine snow layers
zminLayer = (/zminLayer1, zminLayer2, zminLayer3, zminLayer4, zminLayer5/)

[...]

do iSnow=kSnow+1,nCheck
! associate local variables with the information in the data structures
! NOTE: do this here, since the layer variables are re-defined
associate(&
mLayerDepth => prog_data%var(iLookPROG%mLayerDepth)%dat , & ! depth of each layer (m)
mLayerVolFracIce => prog_data%var(iLookPROG%mLayerVolFracIce)%dat , & ! volumetric fraction of ice in each layer (-)
mLayerVolFracLiq => prog_data%var(iLookPROG%mLayerVolFracLiq)%dat & ! volumetric fraction of liquid water in each layer (-)
) ! (associating local variables with the information in the data structures)
! check if the layer depth is less than the depth threshold
select case(ix_snowLayers)
case(sameRulesAllLayers); removeLayer = (mLayerDepth(iSnow) < zmin)
case(rulesDependLayerIndex); removeLayer = (mLayerDepth(iSnow) < zminLayer(iSnow))
case default; err=20; message=trim(message)//'unable to identify option to combine/sub-divide snow layers'; return
end select ! (option to combine/sub-divide snow layers)

User-specified parameters:

zmin                      |       0.0100 |       0.0050 |       0.1000
zmax                      |       0.0500 |       0.0100 |       0.5000
! ---
zminLayer1                |       0.0075 |       0.0075 |       0.0075
zminLayer2                |       0.0100 |       0.0100 |       0.0100
zminLayer3                |       0.0500 |       0.0500 |       0.0500
zminLayer4                |       0.1000 |       0.1000 |       0.1000
zminLayer5                |       0.2500 |       0.2500 |       0.2500
! ---
zmaxLayer1_lower          |       0.0500 |       0.0500 |       0.0500
zmaxLayer2_lower          |       0.2000 |       0.2000 |       0.2000
zmaxLayer3_lower          |       0.5000 |       0.5000 |       0.5000
zmaxLayer4_lower          |       1.0000 |       1.0000 |       1.0000
! ---
zmaxLayer1_upper          |       0.0300 |       0.0300 |       0.0300
zmaxLayer2_upper          |       0.1500 |       0.1500 |       0.1500
zmaxLayer3_upper          |       0.3000 |       0.3000 |       0.3000
zmaxLayer4_upper          |       0.7500 |       0.7500 |       0.7500

I'm not a 100% sure but I expect it may be difficult (if not impossible) to mimic the current behavior using the new code and custom parameter values. What are your thoughts on this? Is requiring recalibration a reason to not make this change?

@andywood
Copy link
Collaborator

andywood commented Mar 8, 2023 via email

@wknoben
Copy link
Collaborator Author

wknoben commented Mar 8, 2023

This change should lead to more intuitive logic, but you're right that the PR doesn't mention this yet. I'm still working on it.

What we currently have is the following:

  • If there are no snow layers (nSnow == 0),
  • And "snow-without-a-layer" (stored in scalarSnowDepth) is larger than the maximum depth we want the first layer to have (zmaxLayer1_lower),
  • Make a new layer.

If we then look at the case where there are snow layers, the logic is:

  • If there is exactly one snow layer,
  • Set the threshold for sub-dividing this layer to zmaxLayer1_lower
  • Then, if mLayerDepth for that 1 layer is > zmaxLayer1_lower,
  • Subdivide the layer.

So essentially with the current logic, the conditions that lead to creating the first layer would lead to subdividing that layer on the next time step if nothing changes. In other words, the same conditions lead to a different number of snow layers depending on if there are 0 or 1 snow layers when the check is made.

This change should lead to more consistent logic in how many snow layers you get for a given snow depth.

@andywood
Copy link
Collaborator

andywood commented Mar 8, 2023 via email

@martynpclark martynpclark marked this pull request as ready for review July 5, 2023 20:43
@martynpclark martynpclark merged commit ecfa1a1 into CH-Earth:develop Jul 5, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants