Merge branch 'develop' into feature/update_example_data_code

docs/source/api/soil/env_factors.md

@@ -0,0 +1,23 @@
+---
+jupytext:
+  cell_metadata_filter: -all
+  formats: md:myst
+  main_language: python
+  text_representation:
+    extension: .md
+    format_name: myst
+    format_version: 0.13
+    jupytext_version: 1.13.8
+kernelspec:
+  display_name: vr_python3
+  language: python
+  name: vr_python3
+---
+
+# API documentation for the {mod}`~virtual_rainforest.models.soil.env_factors` module
+
+```{eval-rst}
+.. automodule:: virtual_rainforest.models.soil.env_factors
+    :autosummary:
+    :members:
+```

docs/source/index.md

@@ -100,6 +100,7 @@ team.
   Soil Overview <api/soil.md>
   Soil Model <api/soil/soil_model.md>
   Soil Carbon <api/soil/carbon.md>
+  Soil Environmental Factors <api/soil/env_factors.md>
   Soil Constants <api/soil/constants.md>
   Abiotic Simple Overview <api/abiotic_simple.md>
   Abiotic Simple Model <api/abiotic_simple/abiotic_simple_model.md>

docs/source/refs.bib

@@ -1,3 +1,36 @@
+@article{porporato_hydrologic_2003,
+	title = {Hydrologic controls on soil carbon and nitrogen cycles. {I}. {Modeling} scheme},
+	volume = {26},
+	issn = {03091708},
+	url = {https://linkinghub.elsevier.com/retrieve/pii/S0309170802000945},
+	doi = {10.1016/S0309-1708(02)00094-5},
+	language = {en},
+	number = {1},
+	urldate = {2023-11-14},
+	journal = {Advances in Water Resources},
+	author = {Porporato, A and D’Odorico, P and Laio, F and Rodriguez-Iturbe, I},
+	month = jan,
+	year = {2003},
+	pages = {45--58},
+}
+
+@article{orwin_organic_2011,
+	title = {Organic nutrient uptake by mycorrhizal fungi enhances ecosystem carbon storage: a model-based assessment},
+	volume = {14},
+	issn = {1461023X},
+	shorttitle = {Organic nutrient uptake by mycorrhizal fungi enhances ecosystem carbon storage},
+	url = {https://onlinelibrary.wiley.com/doi/10.1111/j.1461-0248.2011.01611.x},
+	doi = {10.1111/j.1461-0248.2011.01611.x},
+	language = {en},
+	number = {5},
+	urldate = {2023-10-18},
+	journal = {Ecology Letters},
+	author = {Orwin, Kate H. and Kirschbaum, Miko U. F. and St John, Mark G. and Dickie, Ian A.},
+	month = may,
+	year = {2011},
+	pages = {493--502},
+}
+
 
 @book{monteith_light_1969,
 	address = {Madison, Wisconsin, U.S.A.},
@@ -333,22 +366,6 @@ @book{campbell_introduction_1998
 	year = {1998},
 }
 
-@article{mayes_relation_2012,
-	title = {Relation between {Soil} {Order} and {Sorption} of {Dissolved} {Organic} {Carbon} in {Temperate} {Subsoils}},
-	volume = {76},
-	issn = {0361-5995, 1435-0661},
-	url = {https://onlinelibrary.wiley.com/doi/10.2136/sssaj2011.0340},
-	doi = {10.2136/sssaj2011.0340},
-	language = {en},
-	number = {3},
-	urldate = {2023-02-22},
-	journal = {Soil Science Society of America Journal},
-	author = {Mayes, Melanie A. and Heal, Katherine R. and Brandt, Craig C. and Phillips, Jana R. and Jardine, Philip M.},
-	month = may,
-	year = {2012},
-	pages = {1027--1037},
-}
-
 @article{davis_simple_2017,
 	title = {Simple process-led algorithms for simulating habitats ({SPLASH} v.1.0): robust indices of radiation, evapotranspiration and plant-available moisture},
 	volume = {10},

tests/conftest.py

@@ -172,12 +172,22 @@ def dummy_carbon_data(layer_roles_fixture):
     """Microbial biomass (carbon) pool (kg C m^-3)"""
     data["soil_c_pool_pom"] = DataArray([0.1, 1.0, 0.7, 0.35], dims=["cell_id"])
     """Particulate organic matter pool (kg C m^-3)"""
+    data["soil_enzyme_pom"] = DataArray(
+        [0.022679, 0.009576, 0.050051, 0.003010], dims=["cell_id"]
+    )
+    """Soil enzyme that breaks down particulate organic matter (kg C m^-3)"""
+    data["soil_enzyme_maom"] = DataArray(
+        [0.0356, 0.0117, 0.02509, 0.00456], dims=["cell_id"]
+    )
+    """Soil enzyme that breaks down mineral associated organic matter (kg C m^-3)"""
     data["pH"] = DataArray([3.0, 7.5, 9.0, 5.7], dims=["cell_id"])
     data["bulk_density"] = DataArray([1350.0, 1800.0, 1000.0, 1500.0], dims=["cell_id"])
-    data["percent_clay"] = DataArray([80.0, 30.0, 10.0, 90.0], dims=["cell_id"])
+    data["clay_fraction"] = DataArray([0.8, 0.3, 0.1, 0.9], dims=["cell_id"])
     data["litter_C_mineralisation_rate"] = DataArray(
         [0.00212106, 0.00106053, 0.00049000, 0.0055], dims=["cell_id"]
     )
+    # Data for combined vertical flow (for entire timestep)
+    data["vertical_flow"] = DataArray([3.0, 15.0, 75.0, 47.7], dims=["cell_id"])
 
     # The layer dependant data has to be handled separately
     data["soil_moisture"] = xr.concat(
@@ -186,7 +196,7 @@ def dummy_carbon_data(layer_roles_fixture):
             # At present the soil model only uses the top soil layer, so this is the
             # only one with real test values in
             DataArray(
-                [[0.472467929, 0.399900047, 0.256053640, 0.153616897]],
+                [[0.9304620050, 0.787549327, 0.504263188, 0.302527807]],
                 dims=["layers", "cell_id"],
             ),
             DataArray(np.full((1, 4), np.nan), dims=["layers", "cell_id"]),
@@ -200,7 +210,6 @@ def dummy_carbon_data(layer_roles_fixture):
             "cell_id": data.grid.cell_id,
         }
     )
-    # TODO - Eventually this should replace the dummy soil moisture entirely
     data["matric_potential"] = xr.concat(
         [
             DataArray(np.full((13, 4), np.nan), dims=["layers", "cell_id"]),

tests/core/test_constants_class.py

@@ -29,14 +29,14 @@ class Test(ConstantsDataclass):  # type: ignore [misc]
         pytest.param(
             {},
             does_not_raise(),
-            123.4,
+            0.25,
             (),
             id="defaults_with_no_config",
         ),
         pytest.param(
-            {"placeholder": 432.1},
+            {"depth_of_active_soil_layer": 1.55},
             does_not_raise(),
-            432.1,
+            1.55,
             (),
             id="configured",
         ),
@@ -70,6 +70,6 @@ def test_ConstantsDataclass_from_config(caplog, config, raises, exp_val, exp_log
         constants_instance = CoreConsts.from_config(config)
 
         if isinstance(raises, does_not_raise):
-            assert constants_instance.placeholder == exp_val
+            assert constants_instance.depth_of_active_soil_layer == exp_val
 
         log_check(caplog=caplog, expected_log=exp_log)

tests/core/test_constants_loader.py

@@ -17,16 +17,16 @@
             "core",
             "CoreConsts",
             does_not_raise(),
-            123.4,
+            0.25,
             ((INFO, "Initialised core.CoreConsts from config"),),
             id="default_values",
         ),
         pytest.param(
-            "[core.constants.CoreConsts]\nplaceholder=432.1",
+            "[core.constants.CoreConsts]\ndepth_of_active_soil_layer=1.5",
             "core",
             "CoreConsts",
             does_not_raise(),
-            432.1,
+            1.5,
             ((INFO, "Initialised core.CoreConsts from config"),),
             id="configured_value",
         ),
@@ -96,7 +96,7 @@ def test_load_constants(
             assert isinstance(constants_instance, CoreConsts)
             # The unconfigurable zero_Celsius should take the default value
             assert constants_instance.zero_Celsius == constants.zero_Celsius
-            # Check the placeholder constant has been configured
-            assert constants_instance.placeholder == exp_val
+            # Check the depth_of_active_soil_layer constant has been configured
+            assert constants_instance.depth_of_active_soil_layer == exp_val
 
         log_check(caplog=caplog, expected_log=exp_log)

tests/core/test_data.py

@@ -986,6 +986,8 @@ def test_output_current_state(mocker, dummy_carbon_data, time_index):
             "soil_c_pool_lmwc",
             "soil_c_pool_microbe",
             "soil_c_pool_pom",
+            "soil_enzyme_pom",
+            "soil_enzyme_maom",
         ],
         time_index,
     )

tests/models/litter/test_litter_model.py

@@ -33,8 +33,8 @@ def litter_model_fixture(dummy_litter_data):
 
 def test_litter_model_initialization(caplog, dummy_litter_data):
     """Test `LitterModel` initialization."""
-
     from virtual_rainforest.core.base_model import BaseModel
+    from virtual_rainforest.core.constants import CoreConsts
     from virtual_rainforest.models.litter.constants import LitterConsts
     from virtual_rainforest.models.litter.litter_model import LitterModel
 
@@ -43,7 +43,8 @@ def test_litter_model_initialization(caplog, dummy_litter_data):
         update_interval=pint.Quantity("1 week"),
         soil_layers=[-0.5, -1.0],
         canopy_layers=10,
-        constants=LitterConsts,
+        model_constants=LitterConsts,
+        core_constants=CoreConsts,
     )
 
     # In cases where it passes then checks that the object has the right properties
@@ -76,7 +77,7 @@ def test_litter_model_initialization(caplog, dummy_litter_data):
 
 def test_litter_model_initialization_no_data(caplog):
     """Test `LitterModel` initialization fails when all data is missing."""
-
+    from virtual_rainforest.core.constants import CoreConsts
     from virtual_rainforest.core.data import Data
     from virtual_rainforest.core.grid import Grid
     from virtual_rainforest.models.litter.constants import LitterConsts
@@ -94,7 +95,8 @@ def test_litter_model_initialization_no_data(caplog):
             update_interval=pint.Quantity("1 week"),
             soil_layers=2,  # FIXME - incorrect soil layer spec in model
             canopy_layers=10,
-            constants=LitterConsts,
+            model_constants=LitterConsts,
+            core_constants=CoreConsts,
         )
 
     # Final check that expected logging entries are produced
@@ -149,7 +151,7 @@ def test_litter_model_initialization_no_data(caplog):
 
 def test_litter_model_initialization_bad_pool_bounds(caplog, dummy_litter_data):
     """Test `LitterModel` initialization fails when litter pools are out of bounds."""
-
+    from virtual_rainforest.core.constants import CoreConsts
     from virtual_rainforest.models.litter.constants import LitterConsts
     from virtual_rainforest.models.litter.litter_model import LitterModel
 
@@ -164,7 +166,8 @@ def test_litter_model_initialization_bad_pool_bounds(caplog, dummy_litter_data):
             update_interval=pint.Quantity("1 week"),
             soil_layers=2,
             canopy_layers=10,
-            constants=LitterConsts,
+            model_constants=LitterConsts,
+            core_constants=CoreConsts,
         )
 
     # Final check that the last log entry is as expected
@@ -177,7 +180,7 @@ def test_litter_model_initialization_bad_pool_bounds(caplog, dummy_litter_data):
 
 def test_litter_model_initialization_bad_lignin_bounds(caplog, dummy_litter_data):
     """Test `LitterModel` initialization fails for lignin proportions not in bounds."""
-
+    from virtual_rainforest.core.constants import CoreConsts
     from virtual_rainforest.models.litter.constants import LitterConsts
     from virtual_rainforest.models.litter.litter_model import LitterModel
 
@@ -187,7 +190,14 @@ def test_litter_model_initialization_bad_lignin_bounds(caplog, dummy_litter_data
         # Put incorrect data in for woody lignin
         litter_data["lignin_woody"] = DataArray([0.5, 0.4, 1.1], dims=["cell_id"])
 
-        LitterModel(litter_data, pint.Quantity("1 week"), 2, 10, constants=LitterConsts)
+        LitterModel(
+            litter_data,
+            pint.Quantity("1 week"),
+            2,
+            10,
+            model_constants=LitterConsts,
+            core_constants=CoreConsts,
+        )
 
     # Final check that expected logging entries are produced
     log_check(
@@ -207,6 +217,7 @@ def test_litter_model_initialization_bad_lignin_bounds(caplog, dummy_litter_data
             does_not_raise(),
             (
                 (INFO, "Initialised litter.LitterConsts from config"),
+                (INFO, "Initialised core.CoreConsts from config"),
                 (
                     INFO,
                     "Information required to initialise the litter model successfully "
@@ -255,6 +266,7 @@ def test_litter_model_initialization_bad_lignin_bounds(caplog, dummy_litter_data
             does_not_raise(),
             (
                 (INFO, "Initialised litter.LitterConsts from config"),
+                (INFO, "Initialised core.CoreConsts from config"),
                 (
                     INFO,
                     "Information required to initialise the litter model successfully "
@@ -327,7 +339,7 @@ def test_generate_litter_model(
             pint.Quantity(config["core"]["timing"]["update_interval"]),
         )
         assert model.update_interval == time_interval
-        assert model.constants.litter_decomp_temp_response == temp_response
+        assert model.model_constants.litter_decomp_temp_response == temp_response
 
     # Final check that expected logging entries are produced
     log_check(caplog, expected_log_entries)

tests/models/litter/test_litter_pools.py

@@ -131,6 +131,7 @@ def test_calculate_change_in_litter_variables(
     dummy_litter_data, surface_layer_index, top_soil_layer_index
 ):
     """Test that litter pool update calculation is correct."""
+    from virtual_rainforest.core.constants import CoreConsts
     from virtual_rainforest.models.litter.litter_pools import (
         calculate_change_in_litter_variables,
     )
@@ -168,7 +169,8 @@ def test_calculate_change_in_litter_variables(
         decomposed_excrement=dummy_litter_data["decomposed_excrement"].to_numpy(),
         decomposed_carcasses=dummy_litter_data["decomposed_carcasses"].to_numpy(),
         update_interval=1.0,
-        constants=LitterConsts,
+        model_constants=LitterConsts,
+        core_constants=CoreConsts,
     )
 
     for name in expected_pools.keys():
@@ -206,15 +208,15 @@ def test_calculate_decay_rates(dummy_litter_data, temp_and_water_factors):
 
 def test_calculate_total_C_mineralised(decay_rates):
     """Test that calculation of total C mineralised is as expected."""
+    from virtual_rainforest.core.constants import CoreConsts
     from virtual_rainforest.models.litter.litter_pools import (
         calculate_total_C_mineralised,
     )
 
     expected_mineralisation = [0.0212182, 0.0274272, 0.00617274]
 
     actual_mineralisation = calculate_total_C_mineralised(
-        decay_rates=decay_rates,
-        constants=LitterConsts,
+        decay_rates=decay_rates, model_constants=LitterConsts, core_constants=CoreConsts
     )
 
     assert np.allclose(actual_mineralisation, expected_mineralisation)

tests/models/soil/conftest.py

@@ -0,0 +1,51 @@
+"""Collection of fixtures to assist the testing of the soil model."""
+
+import pytest
+
+from virtual_rainforest.models.soil.constants import SoilConsts
+
+
+@pytest.fixture
+def environmental_factors(dummy_carbon_data, top_soil_layer_index):
+    """Environmental factors based on dummy carbon data."""
+    from virtual_rainforest.models.soil.env_factors import (
+        calculate_clay_impact_on_enzyme_saturation,
+        calculate_clay_impact_on_necromass_decay,
+        calculate_pH_suitability,
+        calculate_water_potential_impact_on_microbes,
+    )
+
+    water_factors = calculate_water_potential_impact_on_microbes(
+        water_potential=dummy_carbon_data["matric_potential"][
+            top_soil_layer_index
+        ].to_numpy(),
+        water_potential_halt=SoilConsts.soil_microbe_water_potential_halt,
+        water_potential_opt=SoilConsts.soil_microbe_water_potential_optimum,
+        response_curvature=SoilConsts.microbial_water_response_curvature,
+    )
+
+    pH_factors = calculate_pH_suitability(
+        soil_pH=dummy_carbon_data["pH"].to_numpy(),
+        maximum_pH=SoilConsts.max_pH_microbes,
+        minimum_pH=SoilConsts.min_pH_microbes,
+        lower_optimum_pH=SoilConsts.lowest_optimal_pH_microbes,
+        upper_optimum_pH=SoilConsts.highest_optimal_pH_microbes,
+    )
+
+    clay_saturation_factors = calculate_clay_impact_on_enzyme_saturation(
+        clay_fraction=dummy_carbon_data["clay_fraction"].to_numpy(),
+        base_protection=SoilConsts.base_soil_protection,
+        protection_with_clay=SoilConsts.soil_protection_with_clay,
+    )
+
+    clay_decay_factors = calculate_clay_impact_on_necromass_decay(
+        clay_fraction=dummy_carbon_data["clay_fraction"].to_numpy(),
+        decay_exponent=SoilConsts.clay_necromass_decay_exponent,
+    )
+
+    return {
+        "water": water_factors,
+        "pH": pH_factors,
+        "clay_saturation": clay_saturation_factors,
+        "clay_decay": clay_decay_factors,
+    }