diff --git a/tests/models/soil/test_carbon.py b/tests/models/soil/test_carbon.py
index 288b14ee1..135c4d600 100644
--- a/tests/models/soil/test_carbon.py
+++ b/tests/models/soil/test_carbon.py
@@ -15,9 +15,9 @@ def test_calculate_soil_carbon_updates(dummy_carbon_data, top_soil_layer_index):
from virtual_rainforest.models.soil.carbon import calculate_soil_carbon_updates
change_in_pools = {
- "soil_c_pool_lmwc": [-3.70969075e-2, -2.2791314e-2, -0.15100693, 7.2342438e-4],
- "soil_c_pool_maom": [1.6992235e-4, 4.45995156e-3, 6.252756e-3, 2.0712399e-5],
- "soil_c_pool_microbe": [-0.03828773, -0.01245439, -0.06446385, -0.00711458],
+ "soil_c_pool_lmwc": [-0.00371115, 0.00278502, -0.01849181, 0.00089995],
+ "soil_c_pool_maom": [-1.28996257e-3, 2.35822401e-3, 1.5570399e-3, 1.2082886e-5],
+ "soil_c_pool_microbe": [-0.04978105, -0.02020101, -0.10280967, -0.00719517],
"soil_c_pool_pom": [0.04809165, 0.01023544, 0.07853728, 0.01167564],
"soil_enzyme_pom": [1.18e-8, 1.67e-8, 1.8e-9, -1.12e-8],
"soil_enzyme_maom": [-0.00031009, -5.09593e-5, 0.0005990658, -3.72112e-5],
@@ -165,8 +165,8 @@ def test_calculate_microbial_carbon_uptake(
"""Check microbial carbon uptake calculates correctly."""
from virtual_rainforest.models.soil.carbon import calculate_microbial_carbon_uptake
- expected_uptake = [0.04484178, 0.03190812, 0.18552910, 0.00022563]
- expected_assimilation = [0.01614304, 0.01052968, 0.05565873, 1.08303e-4]
+ expected_uptake = [1.29159055e-2, 8.43352433e-3, 5.77096991e-2, 5.77363558e-5]
+ expected_assimilation = [4.64972597e-3, 2.78306303e-3, 1.73129097e-2, 2.77134508e-5]
actual_uptake, actual_assimilation = calculate_microbial_carbon_uptake(
soil_c_pool_lmwc=dummy_carbon_data["soil_c_pool_lmwc"],
@@ -200,7 +200,11 @@ def test_calculate_enzyme_mediated_decomposition(
pH_factor=environmental_factors["pH"],
clay_factor_saturation=environmental_factors["clay_saturation"],
soil_temp=dummy_carbon_data["soil_temperature"][top_soil_layer_index],
- constants=SoilConsts,
+ reference_temp=SoilConsts.arrhenius_reference_temp,
+ max_decomp_rate=SoilConsts.max_decomp_rate_pom,
+ activation_energy_rate=SoilConsts.activation_energy_pom_decomp_rate,
+ half_saturation=SoilConsts.half_sat_pom_decomposition,
+ activation_energy_sat=SoilConsts.activation_energy_pom_decomp_saturation,
)
assert np.allclose(actual_decomp, expected_decomp)
diff --git a/tests/models/soil/test_soil_model.py b/tests/models/soil/test_soil_model.py
index cb607367a..d07ef3ab4 100644
--- a/tests/models/soil/test_soil_model.py
+++ b/tests/models/soil/test_soil_model.py
@@ -294,22 +294,22 @@ def test_update(mocker, soil_model_fixture, dummy_carbon_data):
Dataset(
data_vars=dict(
lmwc=DataArray(
- [0.03460932, 0.01198232, 0.04607224, 0.00520078], dims="cell_id"
+ [0.04823845, 0.02119714, 0.0895863, 0.00528887], dims="cell_id"
),
maom=DataArray(
- [2.50009515, 1.70223563, 4.50321333, 0.50001044], dims="cell_id"
+ [2.49936689, 1.70118553, 4.50085129, 0.50000614], dims="cell_id"
),
microbe=DataArray(
- [5.779758, 2.2926573, 11.26080985, 0.99645009], dims="cell_id"
+ [5.77512315, 2.2899636, 11.24827514, 0.99640928], dims="cell_id"
),
pom=DataArray(
- [0.12398572, 1.00509262, 0.73902017, 0.35583213], dims="cell_id"
+ [0.12397575, 1.00508662, 0.7389913, 0.35583206], dims="cell_id"
),
enzyme_pom=DataArray(
- [0.02267854, 0.00957584, 0.05005016, 0.00300993], dims="cell_id"
+ [0.02267842, 0.00957576, 0.05004963, 0.00300993], dims="cell_id"
),
enzyme_maom=DataArray(
- [0.03544542, 0.0116745, 0.02538689, 0.00454144], dims="cell_id"
+ [0.0354453, 0.01167442, 0.02538637, 0.00454144], dims="cell_id"
),
)
),
@@ -424,18 +424,18 @@ def test_construct_full_soil_model(dummy_carbon_data, top_soil_layer_index):
from virtual_rainforest.models.soil.soil_model import construct_full_soil_model
delta_pools = [
- -3.70969075e-2,
- -2.2791314e-2,
- -0.15100693,
- 7.2342438e-4,
- 1.6992235e-4,
- 4.45995156e-3,
- 6.252756e-3,
- 2.0712399e-5,
- -3.82877348e-2,
- -1.24543917e-2,
- -6.44638512e-2,
- -7.11457978e-3,
+ -0.00371115,
+ 0.00278502,
+ -0.01849181,
+ 0.00089995,
+ -1.28996257e-3,
+ 2.35822401e-3,
+ 1.5570399e-3,
+ 1.2082886e-5,
+ -0.04978105,
+ -0.02020101,
+ -0.10280967,
+ -0.00719517,
4.80916464e-2,
1.02354410e-2,
7.85372753e-2,
diff --git a/virtual_rainforest/models/soil/carbon.py b/virtual_rainforest/models/soil/carbon.py
index 377bd8d76..5bea5394c 100644
--- a/virtual_rainforest/models/soil/carbon.py
+++ b/virtual_rainforest/models/soil/carbon.py
@@ -38,7 +38,9 @@ class MicrobialBiomassLoss:
"""Rate at which biomass is lost to the POM pool [kg C m^-3 day^-1]."""
-# TODO - This function should probably be shortened
+# TODO - This function should probably be shortened. First step is group the calculation
+# of environmental factors. Otherwise there's probably stuff that I'm missing, which I
+# should tackle later
def calculate_soil_carbon_updates(
soil_c_pool_lmwc: NDArray[np.float32],
soil_c_pool_maom: NDArray[np.float32],
@@ -150,26 +152,45 @@ def calculate_soil_carbon_updates(
pH_factor=pH_factor,
clay_factor_saturation=clay_factor_saturation,
soil_temp=soil_temp,
- constants=model_constants,
+ reference_temp=model_constants.arrhenius_reference_temp,
+ max_decomp_rate=model_constants.max_decomp_rate_pom,
+ activation_energy_rate=model_constants.activation_energy_pom_decomp_rate,
+ half_saturation=model_constants.half_sat_pom_decomposition,
+ activation_energy_sat=model_constants.activation_energy_pom_decomp_saturation,
)
-
+ # Calculate how pom decomposition is split between lmwc and maom pools
pom_decomposition_to_lmwc = (
pom_decomposition_rate * model_constants.pom_decomposition_fraction_lmwc
)
pom_decomposition_to_maom = pom_decomposition_rate * (
1 - model_constants.pom_decomposition_fraction_lmwc
)
- # TODO - Use enzyme decomp equation to calculate lmwc to maom transfer
+ maom_decomposition_to_lmwc = calculate_enzyme_mediated_decomposition(
+ soil_c_pool=soil_c_pool_maom,
+ soil_enzyme=soil_enzyme_maom,
+ water_factor=water_factor,
+ pH_factor=pH_factor,
+ clay_factor_saturation=clay_factor_saturation,
+ soil_temp=soil_temp,
+ reference_temp=model_constants.arrhenius_reference_temp,
+ max_decomp_rate=model_constants.max_decomp_rate_maom,
+ activation_energy_rate=model_constants.activation_energy_maom_decomp_rate,
+ half_saturation=model_constants.half_sat_maom_decomposition,
+ activation_energy_sat=model_constants.activation_energy_maom_decomp_saturation,
+ )
# Determine net changes to the pools
delta_pools_ordered["soil_c_pool_lmwc"] = (
pom_decomposition_to_lmwc
+ biomass_losses.necromass_decay_to_lmwc
+ pom_enzyme_turnover
+ + maom_decomposition_to_lmwc
- microbial_uptake
- labile_carbon_leaching
)
- delta_pools_ordered["soil_c_pool_maom"] = pom_decomposition_to_maom
+ delta_pools_ordered["soil_c_pool_maom"] = (
+ pom_decomposition_to_maom - maom_decomposition_to_lmwc
+ )
delta_pools_ordered["soil_c_pool_microbe"] = (
microbial_assimilation - biomass_losses.maintenance_synthesis
)
@@ -394,7 +415,11 @@ def calculate_enzyme_mediated_decomposition(
pH_factor: NDArray[np.float32],
clay_factor_saturation: NDArray[np.float32],
soil_temp: NDArray[np.float32],
- constants: SoilConsts,
+ reference_temp: float,
+ max_decomp_rate: float,
+ activation_energy_rate: float,
+ half_saturation: float,
+ activation_energy_sat: float,
) -> NDArray[np.float32]:
"""Calculate rate of a enzyme mediated decomposition process.
@@ -413,7 +438,15 @@ def calculate_enzyme_mediated_decomposition(
clay_factor_saturation: A factor capturing the impact of soil clay fraction on
enzyme saturation constants [unitless]
soil_temp: soil temperature for each soil grid cell [degrees C]
- constants: Set of constants for the soil model.
+ reference_temp: The reference temperature that enzyme rates were determined
+ relative to [degrees C]
+ max_decomp_rate: The maximum rate of substrate decomposition (at the reference
+ temperature) [day^-1]
+ activation_energy_rate: Activation energy for maximum decomposition rate
+ [J K^-1]
+ half_saturation: Half saturation constant for decomposition (at the reference
+ temperature) [kg C m^-3]
+ activation_energy_sat: Activation energy for decomposition saturation [J K^-1]
Returns:
The rate of decomposition of the organic matter pool in question [kg C m^-3
@@ -423,23 +456,19 @@ def calculate_enzyme_mediated_decomposition(
# Calculate the factors which impact the rate and saturation constants
temp_factor_rate = calculate_temperature_effect_on_microbes(
soil_temperature=soil_temp,
- activation_energy=constants.activation_energy_pom_decomp_rate,
- reference_temperature=constants.arrhenius_reference_temp,
+ activation_energy=activation_energy_rate,
+ reference_temperature=reference_temp,
)
temp_factor_saturation = calculate_temperature_effect_on_microbes(
soil_temperature=soil_temp,
- activation_energy=constants.activation_energy_pom_decomp_saturation,
- reference_temperature=constants.arrhenius_reference_temp,
+ activation_energy=activation_energy_sat,
+ reference_temperature=reference_temp,
)
# Calculate the adjusted rate and saturation constants
- rate_constant = (
- constants.max_decomp_rate_pom * temp_factor_rate * water_factor * pH_factor
- )
+ rate_constant = max_decomp_rate * temp_factor_rate * water_factor * pH_factor
saturation_constant = (
- constants.half_sat_pom_decomposition
- * temp_factor_saturation
- * clay_factor_saturation
+ half_saturation * temp_factor_saturation * clay_factor_saturation
)
return (
diff --git a/virtual_rainforest/models/soil/constants.py b/virtual_rainforest/models/soil/constants.py
index df7a60b1f..3b51823e0 100644
--- a/virtual_rainforest/models/soil/constants.py
+++ b/virtual_rainforest/models/soil/constants.py
@@ -7,7 +7,7 @@
from virtual_rainforest.core.constants_class import ConstantsDataclass
-# TODO - Need to figure out a sensible area to volume conversion
+# TODO - Once lignin is tracked a large number of constants will have to be duplicated
@dataclass(frozen=True)
@@ -77,14 +77,14 @@ class SoilConsts(ConstantsDataclass):
:attr:`max_uptake_rate_labile_C`.
"""
- half_sat_labile_C_uptake: float = 0.091
+ half_sat_labile_C_uptake: float = 0.364
"""Half saturation constant for microbial uptake of labile carbon (LMWC).
- [kg C m^-2]. This was calculated from the value provided in
- :cite:t:`wang_development_2013` assuming an average bulk density of 1400 [kg m^-3],
- and a topsoil depth of 0.25 [m]. The reference temperature is given by
- :attr:`arrhenius_reference_temp`, and the corresponding activation energy is given
- by :attr:`activation_energy_labile_C_saturation`.
+ [kg C m^-3]. This was calculated from the value provided in
+ :cite:t:`wang_development_2013` assuming an average bulk density of 1400 [kg m^-3].
+ The reference temperature is given by :attr:`arrhenius_reference_temp`, and the
+ corresponding activation energy is given by
+ :attr:`activation_energy_labile_C_saturation`.
"""
activation_energy_labile_C_saturation: float = 30000
@@ -93,15 +93,13 @@ class SoilConsts(ConstantsDataclass):
Taken from :cite:t:`wang_development_2013`.
"""
- # TODO - Add another set of constants once we start tracking lignin
half_sat_pom_decomposition: float = 70.0
- """Half saturation constant for POM decomposition to LMWC [kg C m^-2].
+ """Half saturation constant for POM decomposition to LMWC [kg C m^-3].
This was calculated from the value provided in :cite:t:`wang_development_2013`
- assuming an average bulk density of 1400 [kg m^-3], and a topsoil depth of 0.25 [m].
- The reference temperature is given by :attr:`arrhenius_reference_temp`, and the
- corresponding activation energy is given by
- :attr:`activation_energy_pom_decomp_saturation`.
+ assuming an average bulk density of 1400 [kg m^-3]. The reference temperature is
+ given by :attr:`arrhenius_reference_temp`, and the corresponding activation energy
+ is given by :attr:`activation_energy_pom_decomp_saturation`.
"""
activation_energy_pom_decomp_saturation: float = 30000
@@ -110,7 +108,6 @@ class SoilConsts(ConstantsDataclass):
Taken from :cite:t:`wang_development_2013`.
"""
- # TODO - Add another set of constants once we start tracking lignin
max_decomp_rate_pom: float = 60.0
"""Maximum rate for particulate organic matter break down (at reference temp).
@@ -128,6 +125,36 @@ class SoilConsts(ConstantsDataclass):
Taken from :cite:t:`wang_development_2013`.
"""
+ half_sat_maom_decomposition: float = 350.0
+ """Half saturation constant for MAOM decomposition to LMWC [kg C m^-3].
+
+ This was calculated from the value provided in :cite:t:`wang_development_2013`
+ assuming an average bulk density of 1400 [kg m^-3]. The reference temperature is
+ given by :attr:`arrhenius_reference_temp`, and the corresponding activation energy
+ is given by :attr:`activation_energy_maom_decomp_saturation`.
+ """
+
+ activation_energy_maom_decomp_saturation: float = 30000
+ """Activation energy for MAOM decomposition saturation constant [J K^-1].
+
+ Taken from :cite:t:`wang_development_2013`.
+ """
+
+ max_decomp_rate_maom: float = 24.0
+ """Maximum rate for mineral associated organic matter decomposition enzyme.
+
+ Units of [day^-1]. The rate is for a reference temperature which is given by
+ :attr:`arrhenius_reference_temp`, and the corresponding activation energy is given
+ by :attr:`activation_energy_maom_decomp_rate`. The value is taken from
+ :cite:t:`wang_development_2013`.
+ """
+
+ activation_energy_maom_decomp_rate: float = 47000
+ """Activation energy for decomposition of mineral associated organic matter.
+
+ Units of [J K^-1]. Taken from :cite:t:`wang_development_2013`.
+ """
+
# TODO - Split this and the following into 2 constants once fungi are introduced
microbial_turnover_rate: float = 0.005
"""Microbial turnover rate at reference temperature [day^-1].