Merge pull request #539 from ImperialCollegeLondon/512-litter-phosphorus

Add phosphorus to litter model
ImperialCollegeLondon · Sep 10, 2024 · 1b9487d · 1b9487d
2 parents ca04b51 + 76f76f9
commit 1b9487d
Show file tree

Hide file tree

Showing 16 changed files with 749 additions and 100 deletions.
diff --git a/docs/source/refs.bib b/docs/source/refs.bib
@@ -839,4 +839,22 @@ @article{harfoot_madingley_2014
   publisher={Public Library of Science},
   doi={10.1371/journal.pbio.1001841},
   url={https://journals.plos.org/plosbiology/article?id=10.1371/journal.pbio.1001841}
-}
+}
+
+@article{cleveland_cnp_2007,
+	title = {C:{N}:{P} stoichiometry in soil: is there a “{Redfield} ratio” for the microbial biomass?},
+	volume = {85},
+	copyright = {http://www.springer.com/tdm},
+	issn = {0168-2563, 1573-515X},
+	shorttitle = {C},
+	url = {http://link.springer.com/10.1007/s10533-007-9132-0},
+	doi = {10.1007/s10533-007-9132-0},
+	language = {en},
+	number = {3},
+	urldate = {2024-08-22},
+	journal = {Biogeochemistry},
+	author = {Cleveland, Cory C. and Liptzin, Daniel},
+	month = aug,
+	year = {2007},
+	pages = {235--252},
+}
diff --git a/tests/models/litter/conftest.py b/tests/models/litter/conftest.py
@@ -17,7 +17,7 @@ def fixture_litter_model(dummy_litter_data, fixture_core_components):
 
     # Build the config object
     config = Config(
-        cfg_strings="[core]\n[core.timing]\nupdate_interval = '24 hours'\n[litter]\n"
+        cfg_strings="[core]\n[core.timing]\nupdate_interval = '48 hours'\n[litter]\n"
     )
     core_components = CoreComponents(config)
 
@@ -54,6 +54,11 @@ def dummy_litter_data(fixture_core_components):
         "c_n_ratio_woody": [55.5, 63.3, 47.3, 59.1],
         "c_n_ratio_below_metabolic": [10.7, 11.3, 15.2, 12.4],
         "c_n_ratio_below_structural": [50.5, 55.6, 73.1, 61.2],
+        "c_p_ratio_above_metabolic": [57.3, 68.7, 100.1, 95.8],
+        "c_p_ratio_above_structural": [337.5, 473.2, 415.8, 570.2],
+        "c_p_ratio_woody": [555.5, 763.3, 847.3, 599.1],
+        "c_p_ratio_below_metabolic": [310.7, 411.3, 315.2, 412.4],
+        "c_p_ratio_below_structural": [550.5, 595.6, 773.1, 651.2],
         "decomposed_excrement": [8e-07, 8.42857e-07, 3.28571e-05, 3.28571e-05],
         "decomposed_carcasses": [1.0714e-4, 4.8571e-4, 1.15714e-3, 1.15714e-3],
         "deadwood_production": [0.075, 0.099, 0.063, 0.033],
@@ -68,6 +73,10 @@ def dummy_litter_data(fixture_core_components):
         "leaf_turnover_c_n_ratio": [15.0, 25.5, 43.1, 57.4],
         "plant_reproductive_tissue_turnover_c_n_ratio": [12.5, 23.8, 15.7, 18.2],
         "root_turnover_c_n_ratio": [30.3, 45.6, 43.3, 37.1],
+        "deadwood_c_p_ratio": [856.5, 675.4, 933.2, 888.8],
+        "leaf_turnover_c_p_ratio": [415.0, 327.4, 554.5, 380.9],
+        "plant_reproductive_tissue_turnover_c_p_ratio": [125.5, 105.0, 145.0, 189.2],
+        "root_turnover_c_p_ratio": [656.7, 450.6, 437.3, 371.9],
     }
 
     for var, vals in pool_values.items():
@@ -140,6 +149,11 @@ def metabolic_splits(dummy_litter_data):
             "plant_reproductive_tissue_turnover_c_n_ratio"
         ].to_numpy(),
         root_turnover_c_n_ratio=dummy_litter_data["root_turnover_c_n_ratio"].to_numpy(),
+        leaf_turnover_c_p_ratio=dummy_litter_data["leaf_turnover_c_p_ratio"].to_numpy(),
+        reproduct_turnover_c_p_ratio=dummy_litter_data[
+            "plant_reproductive_tissue_turnover_c_p_ratio"
+        ].to_numpy(),
+        root_turnover_c_p_ratio=dummy_litter_data["root_turnover_c_p_ratio"].to_numpy(),
         constants=LitterConsts,
     )
 
@@ -197,3 +211,36 @@ def input_c_n_ratios(dummy_litter_data, metabolic_splits, litter_chemistry):
     )
 
     return input_c_n_ratios
+
+
+@pytest.fixture
+def input_c_p_ratios(dummy_litter_data, metabolic_splits, litter_chemistry):
+    """Carbon:nitrogen ratio of each input flow."""
+
+    input_c_p_ratios = litter_chemistry.calculate_litter_input_phosphorus_ratios(
+        metabolic_splits=metabolic_splits,
+        struct_to_meta_phosphorus_ratio=LitterConsts.structural_to_metabolic_p_ratio,
+    )
+
+    return input_c_p_ratios
+
+
+@pytest.fixture
+def updated_pools(dummy_litter_data, decay_rates, plant_inputs):
+    """Updated carbon mass of each pool."""
+    from virtual_ecosystem.models.litter.carbon import calculate_updated_pools
+
+    updated_pools = calculate_updated_pools(
+        above_metabolic=dummy_litter_data["litter_pool_above_metabolic"].to_numpy(),
+        above_structural=dummy_litter_data["litter_pool_above_structural"].to_numpy(),
+        woody=dummy_litter_data["litter_pool_woody"].to_numpy(),
+        below_metabolic=dummy_litter_data["litter_pool_below_metabolic"].to_numpy(),
+        below_structural=dummy_litter_data["litter_pool_below_structural"].to_numpy(),
+        decomposed_excrement=dummy_litter_data["decomposed_excrement"].to_numpy(),
+        decomposed_carcasses=dummy_litter_data["decomposed_carcasses"].to_numpy(),
+        decay_rates=decay_rates,
+        plant_inputs=plant_inputs,
+        update_interval=2.0,
+    )
+
+    return updated_pools
diff --git a/tests/models/litter/test_carbon.py b/tests/models/litter/test_carbon.py
@@ -105,11 +105,11 @@ def test_calculate_updated_pools(dummy_litter_data, decay_rates, plant_inputs):
     from virtual_ecosystem.models.litter.carbon import calculate_updated_pools
 
     expected_pools = {
-        "above_metabolic": [0.31632696, 0.15296346, 0.08537701, 0.08087947],
-        "above_structural": [0.50453639, 0.25006367, 0.09842669, 0.11162423],
+        "above_metabolic": [0.31567198, 0.1529074957, 0.0813030042, 0.0736771942],
+        "above_structural": [0.50519138, 0.25011962, 0.10250070, 0.11882651],
         "woody": [4.77403361, 11.89845863, 7.3598224, 7.3298224],
-        "below_metabolic": [0.40174907, 0.36687303, 0.06792061, 0.08224246],
-        "below_structural": [0.60638765, 0.31821335, 0.02010401, 0.03038216],
+        "below_metabolic": [0.3976309, 0.3630269, 0.06787947, 0.07794085],
+        "below_structural": [0.61050583, 0.32205947352, 0.02014514530, 0.03468376530],
     }
 
     actual_pools = calculate_updated_pools(

diff --git a/tests/models/litter/test_chemistry.py b/tests/models/litter/test_chemistry.py
@@ -27,27 +27,24 @@ def test_calculate_litter_chemistry_factor():
 
 
 def test_calculate_new_pool_chemistries(
-    dummy_litter_data, plant_inputs, metabolic_splits, litter_chemistry
+    dummy_litter_data, plant_inputs, metabolic_splits, updated_pools, litter_chemistry
 ):
     """Test that function to calculate updated pool chemistries works correctly."""
 
-    updated_pools = {
-        "above_metabolic": np.array([0.32072786, 0.15473132, 0.08523907, 0.08074153]),
-        "above_structural": np.array([0.5047038, 0.25068224, 0.09843778, 0.11163532]),
-        "woody": np.array([4.774517, 11.898729, 7.361411, 7.331411]),
-        "below_metabolic": np.array([0.4090768, 0.37287148, 0.06883228, 0.08315412]),
-        "below_structural": np.array([0.6066315, 0.31860251, 0.02010566, 0.03038382]),
-    }
-
     expected_chemistries = {
-        "lignin_above_structural": [0.49790843, 0.10067782, 0.70495536, 0.71045831],
-        "lignin_woody": [0.49580586, 0.79787834, 0.35224223, 0.35012603],
-        "lignin_below_structural": [0.50313604, 0.26586391, 0.7499951, 0.82142894],
-        "c_n_ratio_above_metabolic": [7.42828417, 8.93702902, 11.13974273, 10.28862942],
-        "c_n_ratio_above_structural": [37.5698310, 43.3465444, 49.0206010, 54.4471558],
-        "c_n_ratio_woody": [55.58168366, 63.25507083, 47.52080006, 59.08199528],
-        "c_n_ratio_below_metabolic": [10.9044015, 11.4675610, 15.2070612, 12.6623415],
-        "c_n_ratio_below_structural": [50.7755820, 56.387878, 73.1837156, 64.0424461],
+        "lignin_above_structural": [0.49726219, 0.10065698, 0.67693666, 0.6673972],
+        "lignin_woody": [0.49580543, 0.7978783, 0.35224272, 0.35012606],
+        "lignin_below_structural": [0.49974338, 0.26270806, 0.74846367, 0.71955592],
+        "c_n_ratio_above_metabolic": [7.39175978, 8.93054462, 10.40414542, 9.86041981],
+        "c_n_ratio_above_structural": [37.5547150, 43.3448492, 48.0974058, 52.0359678],
+        "c_n_ratio_woody": [55.5816919, 63.2550698, 47.5208477, 59.0819914],
+        "c_n_ratio_below_metabolic": [10.7299421, 11.3394567, 15.1984024, 12.2222413],
+        "c_n_ratio_below_structural": [50.6228215, 55.9998994, 73.0948342, 58.6661277],
+        "c_p_ratio_above_metabolic": [69.957176, 68.5502416, 107.1709994, 96.55826106],
+        "c_p_ratio_above_structural": [346.048307, 472.496124, 465.834123, 525.882608],
+        "c_p_ratio_woody": [560.22870571, 762.56863636, 848.03530307, 600.40427444],
+        "c_p_ratio_below_metabolic": [308.200782, 405.110726, 314.824814, 372.870229],
+        "c_p_ratio_below_structural": [563.06464, 597.68324, 772.78968, 609.82810],
     }
 
     actual_chemistries = litter_chemistry.calculate_new_pool_chemistries(
@@ -63,20 +60,14 @@ def test_calculate_new_pool_chemistries(
 
 
 def test_calculate_lignin_updates(
-    dummy_litter_data, plant_inputs, input_lignin, litter_chemistry
+    dummy_litter_data, plant_inputs, input_lignin, updated_pools, litter_chemistry
 ):
     """Test that the function to calculate the lignin updates works as expected."""
 
-    updated_pools = {
-        "above_structural": np.array([0.5047038, 0.25068224, 0.09843778, 0.11163532]),
-        "woody": np.array([4.774517, 11.898729, 7.361411, 7.331411]),
-        "below_structural": np.array([0.6066315, 0.31860251, 0.02010566, 0.03038382]),
-    }
-
     expected_lignin = {
-        "above_structural": [-0.00209157, 0.00067782, 0.00495532, 0.01045834],
-        "woody": [-0.00419414, -0.00212166, 0.00224223, 0.00012603],
-        "below_structural": [3.1360386e-3, 1.5863906e-2, -4.90160482e-6, 7.1428885e-2],
+        "above_structural": [-0.00273781, 0.00065698, -0.02306334, -0.03260280],
+        "woody": [-0.00419457, -0.0021217, 0.00224272, 0.00012606],
+        "below_structural": [-0.00025662, 0.01270806, -0.00153633, -0.03044408],
     }
 
     actual_lignin = litter_chemistry.calculate_lignin_updates(
@@ -113,24 +104,16 @@ def test_calculate_change_in_chemical_concentration(dummy_litter_data):
 
 
 def test_calculate_c_n_ratio_updates(
-    dummy_litter_data, plant_inputs, input_c_n_ratios, litter_chemistry
+    dummy_litter_data, plant_inputs, input_c_n_ratios, updated_pools, litter_chemistry
 ):
     """Test that calculation of C:N ratio updates works properly."""
 
-    updated_pools = {
-        "above_metabolic": np.array([0.32072786, 0.15473132, 0.08523907, 0.08074153]),
-        "above_structural": np.array([0.5047038, 0.25068224, 0.09843778, 0.11163532]),
-        "woody": np.array([4.774517, 11.898729, 7.361411, 7.331411]),
-        "below_metabolic": np.array([0.4090768, 0.37287148, 0.06883228, 0.08315412]),
-        "below_structural": np.array([0.6066315, 0.31860251, 0.02010566, 0.03038382]),
-    }
-
     expected_change = {
-        "above_metabolic": [0.12828416, 0.23702901, 1.03974239, 0.48862956],
-        "above_structural": [0.06983094, 0.14654437, 3.22060275, 4.24715499],
-        "woody": [0.081683655, -0.04492917, 0.220800061, -0.01800472],
-        "below_metabolic": [0.20440145, 0.16756069, 0.00706121, 0.26234147],
-        "below_structural": [0.27558203, 0.78787769, 0.08371555, 2.8424462],
+        "above_metabolic": [0.09175978, 0.23054462, 0.30414542, 0.06041981],
+        "above_structural": [0.05471499, 0.14484922, 2.29740576, 1.835967773],
+        "woody": [0.0816919, -0.0449302, 0.2208477, -0.0180086],
+        "below_metabolic": [0.02994209, 0.03945672, -0.00159759, -0.17775875],
+        "below_structural": [0.12282146, 0.39989943, -0.00516585, -2.53387232],
     }
 
     actual_change = litter_chemistry.calculate_c_n_ratio_updates(
@@ -145,6 +128,31 @@ def test_calculate_c_n_ratio_updates(
         assert np.allclose(actual_change[key], expected_change[key])
 
 
+def test_calculate_c_p_ratio_updates(
+    dummy_litter_data, plant_inputs, input_c_p_ratios, updated_pools, litter_chemistry
+):
+    """Test that calculation of C:P ratio updates works properly."""
+
+    expected_change = {
+        "above_metabolic": [12.657176, -0.14975840, 7.07099940, 0.75826106],
+        "above_structural": [8.5483073, -0.7038764, 50.034123, -44.317392],
+        "woody": [4.72870571, -0.73136364, 0.73530307, 1.30427444],
+        "below_metabolic": [-2.49921796, -6.18927446, -0.37518617, -39.52977135],
+        "below_structural": [12.56464272, 2.08324337, -0.31032454, -41.37190224],
+    }
+
+    actual_change = litter_chemistry.calculate_c_p_ratio_updates(
+        plant_inputs=plant_inputs,
+        input_c_p_ratios=input_c_p_ratios,
+        updated_pools=updated_pools,
+    )
+
+    assert set(expected_change.keys()) == set(actual_change.keys())
+
+    for key in actual_change.keys():
+        assert np.allclose(actual_change[key], expected_change[key])
+
+
 def test_calculate_N_mineralisation(dummy_litter_data, decay_rates, litter_chemistry):
     """Test that function to calculate nitrogen mineralisation rate works properly."""
 
@@ -158,14 +166,27 @@ def test_calculate_N_mineralisation(dummy_litter_data, decay_rates, litter_chemi
     assert np.allclose(actual_n_mineral, expected_n_mineral)
 
 
+def test_calculate_P_mineralisation(dummy_litter_data, decay_rates, litter_chemistry):
+    """Test that function to calculate phosphorus mineralisation rate works properly."""
+
+    expected_p_mineral = [4.39937479e-4, 2.13832149e-4, 6.40698004e-5, 6.56405873e-5]
+
+    actual_p_mineral = litter_chemistry.calculate_P_mineralisation(
+        decay_rates=decay_rates,
+        active_microbe_depth=CoreConsts.max_depth_of_microbial_activity,
+    )
+
+    assert np.allclose(actual_p_mineral, expected_p_mineral)
+
+
 def test_calculate_litter_input_lignin_concentrations(
     dummy_litter_data, plant_inputs, litter_chemistry
 ):
     """Check calculation of lignin concentrations of each plant flow to litter."""
 
     expected_woody = [0.233, 0.545, 0.612, 0.378]
-    expected_concs_above_struct = [0.28329484, 0.23062465, 0.75773447, 0.75393599]
-    expected_concs_below_struct = [0.77196233, 0.80040249, 0.74908861, 0.95895666]
+    expected_concs_above_struct = [0.24971768, 0.22111396, 0.51122474, 0.56571041]
+    expected_concs_below_struct = [0.48590258, 0.56412613, 0.54265483, 0.67810978]
 
     actual_concs = litter_chemistry.calculate_litter_input_lignin_concentrations(
         plant_input_below_struct=plant_inputs["below_ground_structural"],
@@ -184,10 +205,10 @@ def test_calculate_litter_input_nitrogen_ratios(
 
     expected_c_n_ratios = {
         "woody": [60.7, 57.9, 73.1, 55.1],
-        "below_metabolic": [14.879783, 16.587126, 17.733169, 13.903046],
-        "below_structural": [74.398916, 82.935630, 88.665843, 69.515230],
-        "above_metabolic": [8.9373399, 14.343140, 15.968877, 13.520689],
-        "above_structural": [44.735092, 71.440811, 83.323241, 72.103527],
+        "below_metabolic": [11.449427, 13.09700, 14.48056, 11.04331],
+        "below_structural": [57.24714, 65.48498, 72.40281, 55.21655],
+        "above_metabolic": [8.48355299, 14.17116914, 12.3424635, 11.10877484],
+        "above_structural": [42.5018709, 69.9028550, 64.6044513, 57.7622747],
     }
 
     actual_c_n_ratios = litter_chemistry.calculate_litter_input_nitrogen_ratios(
@@ -201,6 +222,30 @@ def test_calculate_litter_input_nitrogen_ratios(
         assert np.allclose(actual_c_n_ratios[key], expected_c_n_ratios[key])
 
 
+def test_calculate_litter_input_phosphorus_ratios(
+    dummy_litter_data, metabolic_splits, litter_chemistry
+):
+    """Check function to calculate the C:P ratios of input to each litter pool works."""
+
+    expected_c_p_ratios = {
+        "woody": [856.5, 675.4, 933.2, 888.8],
+        "below_metabolic": [248.1465, 129.418998, 146.243645, 110.700999],
+        "below_structural": [1240.73249721, 647.09498874, 731.2182237, 553.50499377],
+        "above_metabolic": [220.55713162, 65.14600889, 152.23446238, 112.22496062],
+        "above_structural": [1118.95921, 343.440873, 825.333331, 387.658509],
+    }
+
+    actual_c_p_ratios = litter_chemistry.calculate_litter_input_phosphorus_ratios(
+        metabolic_splits=metabolic_splits,
+        struct_to_meta_phosphorus_ratio=LitterConsts.structural_to_metabolic_p_ratio,
+    )
+
+    assert set(expected_c_p_ratios.keys()) == set(actual_c_p_ratios.keys())
+
+    for key in actual_c_p_ratios.keys():
+        assert np.allclose(actual_c_p_ratios[key], expected_c_p_ratios[key])
+
+
 def test_calculate_nutrient_split_between_litter_pools(
     dummy_litter_data, metabolic_splits
 ):
@@ -209,8 +254,8 @@ def test_calculate_nutrient_split_between_litter_pools(
         calculate_nutrient_split_between_litter_pools,
     )
 
-    expected_meta_c_n = np.array([14.879783, 16.587126, 17.733169, 13.903046])
-    expected_struct_c_n = np.array([74.398916, 82.935630, 88.665843, 69.515230])
+    expected_meta_c_n = np.array([11.449427, 13.09700, 14.48056, 11.04331])
+    expected_struct_c_n = np.array([57.24714, 65.48498, 72.40281, 55.21655])
 
     actual_meta_c_n, actual_struct_c_n = calculate_nutrient_split_between_litter_pools(
         input_c_nut_ratio=dummy_litter_data["root_turnover_c_n_ratio"],