Skip to content

Energy-balance multi-layer snowpack model implementation#15

Open
jpnousu wants to merge 60 commits intomainfrom
snowmodel
Open

Energy-balance multi-layer snowpack model implementation#15
jpnousu wants to merge 60 commits intomainfrom
snowmodel

Conversation

@jpnousu
Copy link
Copy Markdown
Collaborator

@jpnousu jpnousu commented Mar 12, 2026

The original degreeday snow model approach has been complemented with an option to use an energy-balance, multi-layer snowpack model based on FSM2 (Essery et al., 2025). The original snow process modules written in Fortran were translated to Python, made consistent with pyAPES coding style, and coupled to pyAPES' forestfloor computations. Tests: forestfloor runs standalone with both FSM2 and degreeday snowpack models. Full multi-layer canopy model with prescribed canopy also works with both snowmodel options.

Essery, R., Mazzotti, G., Barr, S., Jonas, T., Quaife, T., and Rutter, N.: A Flexible Snow Model (FSM 2.1.1) including a forest canopy, Geosci. Model Dev., 18, 3583–3605, https://doi.org/10.5194/gmd-18-3583-2025, 2025.

Copy link
Copy Markdown
Contributor

@khaahti khaahti left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks good, only minor comments to check before merge!


return fluxes, states

def heat_and_water_exchange_under_snow(self, dt: float, forcing: Dict, parameters: Dict, sub_dt: float=60.0) -> Tuple:
Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe add here in documentation that this is not used do to instability.

If we at some point want to return to this, we should check iterative solution of temperature need when freezing/thawing is happening. Implemented in heat_and_water_exchange().


if self.snow_model == 'degreeday':
fluxes_snow['snow_heat_flux'] = 0.
snow_bottom_temperature = 0.
Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we get reasonable soil temperatures under degree-day snow if we keep these zero? Or are we fine with that.

kmelt melting coefficient [m degC-1 s-1]
kfreeze (float): freezing coefficient coefficient [m degC-1 s-1]
retention (float): max fraction of liquid water in snow [-]
Tmelt (float): melting temperature (~0.0 degC) [degC]
Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

convert units in documentation to Kelvin

}

fluxes = {'potential_infiltration': snow_fluxes['Roff'],
'snow_heat_flux': snow_fluxes['Gsoil'], # heat flux to organiclayer
Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

or actually directly to soil now, right?

# --- Model control flags
ctr = {'Eflow': True, # use ensemble flow statistics; i.e fixed ratio of Utop/ustar.
'WMA': False, # assume air-space scalar profiles well-mixed
'Ebal': True, # computes leaf temperature by solving energy balance
Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Did you check with Ebal=False? When energy balance is not solved the snow model should be forced to degree-day i think.

dat['Prec'] = dat['Prec'] / dt

cols = ['doy', 'Prec', 'P', 'Tair', 'Tdaily', 'U', 'Ustar', 'H2O', 'CO2', 'Zen',
cols = ['doy', 'Prec', 'P', 'Tair', 'Tdaily', 'U', 'Ustar', 'H2O', 'RH', 'CO2', 'Zen',
Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

can be calculated if we want to avoid adding variables to forcing.

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

See also my comment above to mlm_canopy. I think we want to avoid adding variables to forcing to maintain as much backwards compatibility as possible

Copy link
Copy Markdown
Collaborator

@Altesmi Altesmi left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks almost ready, great work!

Two important things

  • Add comments to new if statements what real-life state the if/elif/else block is modelling e.g., if swe <= 0 # no snow -> do not solve snowpack
  • Remove RH dependency from forcing. IMO we don't need both [H2O]_air and RH in forcing files

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Delete this data before merge

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Delete this data before merge

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Delete this result file before merge

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I've made changes to this notebook. Remember to fetch and merge main before merging

if controls['energy_balance']:
# calculate moss / litter energy and water balance
fluxes, states = self.heat_and_water_exchange(
if forcing['snow_water_equivalent'] > 0.:
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could you add comment here (and perhaps elsewhere) what we are doing and why when we are entering a certain if/elif/else block

@@ -180,6 +180,11 @@
['ffloor_nir_albedo', 'NIR albedo (forest floor) [-]', ('date', 'simulation')],
['ffloor_groundtypes', 'forestfloor groundtype names', ('groundtype')],

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Are these the only outputs from fsm2?

'air_pressure': forcing['air_pressure'],
'wind_speed': U[1],
'friction_velocity': ustar[1],
'relative_humidity': forcing['relative_humidity'],
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we need relative humidity from forcing or should we calculate it from [H2O]_air. At minimum relative humidity should be added to the run function documentation that it is needed. Does this brake old run scripts that do not have relative humidity in forcing?

'r10': 2.5, # [umol m-2 s-1]
'q10': 2.0, # [-]
'moisture_coeff': [3.83, 4.43, 1.25, 0.854] # Skopp moisture function param [a ,b, d, g]}
'moisture_coeff': [3.83, 4.43, 1.25, 0.854], # Skopp moisture function param [a ,b, d, g]}
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is skopp moisture function generally known or should we add a reference?

Args:
snowpara (dict):
'initial_conditions' (dict)
'Sice' (np.ndarray): # Snow ice content (kg/m^2)
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Here the style of units has been kg m-2

dat['Prec'] = dat['Prec'] / dt

cols = ['doy', 'Prec', 'P', 'Tair', 'Tdaily', 'U', 'Ustar', 'H2O', 'CO2', 'Zen',
cols = ['doy', 'Prec', 'P', 'Tair', 'Tdaily', 'U', 'Ustar', 'H2O', 'RH', 'CO2', 'Zen',
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

See also my comment above to mlm_canopy. I think we want to avoid adding variables to forcing to maintain as much backwards compatibility as possible

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