From 489780813971b0cd4fa7c1a2a0534cd963de3b94 Mon Sep 17 00:00:00 2001 From: Sally Date: Wed, 9 Apr 2025 13:58:05 +0100 Subject: [PATCH 01/22] StemStochiometry basic structure. --- virtual_ecosystem/models/plants/constants.py | 26 +++-- .../models/plants/plants_model.py | 13 +++ .../models/plants/stochiometry.py | 94 +++++++++++++++++++ 3 files changed, 123 insertions(+), 10 deletions(-) create mode 100644 virtual_ecosystem/models/plants/stochiometry.py diff --git a/virtual_ecosystem/models/plants/constants.py b/virtual_ecosystem/models/plants/constants.py index a282b3039..d914ee3e7 100644 --- a/virtual_ecosystem/models/plants/constants.py +++ b/virtual_ecosystem/models/plants/constants.py @@ -39,20 +39,17 @@ class PlantsConsts(ConstantsDataclass): deadwood_c_n_ratio: float = 56.5 """Carbon to Nitrogen ratio of deadwood.""" - leaf_turnover_c_n_ratio: float = 25.5 - """Carbon to Nitrogen ratio of leaf turnover.""" - - plant_reproductive_tissue_turnover_c_n_ratio: float = 12.5 - """Carbon to Nitrogen ratio of plant reproductive tissue turnover.""" - - root_turnover_c_n_ratio: float = 45.6 - """Carbon to Nitrogen ratio of root turnover.""" - deadwood_c_p_ratio: float = 856.5 """Carbon to Phosphorous ratio of deadwood.""" + leaf_turnover_c_n_ratio: float = 25.5 + """Carbon to Nitrogen ratio of leaf turnover (senesced leaves).""" + leaf_turnover_c_p_ratio: float = 415.0 - """Carbon to Phosphorous ratio of leaf turnover.""" + """Carbon to Phosphorous ratio of leaf turnover (senesced leaves).""" + + plant_reproductive_tissue_turnover_c_n_ratio: float = 12.5 + """Carbon to Nitrogen ratio of plant reproductive tissue turnover.""" plant_reproductive_tissue_turnover_c_p_ratio: float = 125.5 """Carbon to Phosphorous ratio of plant reproductive tissue turnover.""" @@ -60,6 +57,15 @@ class PlantsConsts(ConstantsDataclass): root_turnover_c_p_ratio: float = 656.7 """Carbon to Phosphorous ratio of root turnover.""" + root_turnover_c_n_ratio: float = 45.6 + """Carbon to Nitrogen ratio of root turnover.""" + + foliage_c_n_ratio: float = 15 + """Carbon to Nitrogen ratio of foliage on trees.""" + + foliage_c_p_ratio: float = 300 + """Carbon to Phosphorous ratio of foliage on trees.""" + root_exudates: float = 0.5 """Fraction of GPP topslice allocated to root exudates.""" diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index aa684001a..90129101f 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -29,6 +29,7 @@ from virtual_ecosystem.models.plants.communities import PlantCommunities from virtual_ecosystem.models.plants.constants import PlantsConsts from virtual_ecosystem.models.plants.functional_types import get_flora_from_config +from virtual_ecosystem.models.plants.stochiometry import StemStochiometry class PlantsModel( @@ -187,6 +188,8 @@ def __init__( self.communities: PlantCommunities """An instance of PlantCommunities providing dictionary access keyed by cell id to PlantCommunity instances for each cell.""" + self.stochiometries: dict[int, StemStochiometry] + """A dictionary keyed by cell id giving the stochiometry of each community.""" self._canopy_layer_indices: NDArray[np.bool_] """The indices of the canopy layers within wider vertical profile. This is a shorter reference to self.layer_structure.index_canopy.""" @@ -291,6 +294,16 @@ def _setup( self.communities = PlantCommunities( data=self.data, flora=self.flora, grid=self.grid ) + + # Initialize the stochiometries of each cohort + self.stochiometries = { + cell_id: StemStochiometry( + plant_constants=self.model_constants, + n_cohorts=self.communities[cell_id].number_of_cohorts, + ) + for cell_id in self.communities.keys() + } + # This is widely used internally so store it as an attribute. self._canopy_layer_indices = self.layer_structure.index_canopy diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py new file mode 100644 index 000000000..ba480f220 --- /dev/null +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -0,0 +1,94 @@ +"""The :mod:`~virtual_ecosystem.models.plants.stochiometry` module contains the class +for managing plant cohort stochiometry ratios. The carbon mass is stored in plant +alloemetry or allocation, so this class uses thoses as the anchor weights and stores +CN and CP ratios. + +The class holds current CN and CP ratios for foliage, reproductive tissue, wood, and +roots on the cohort level. Each tissue also has an idea CN and CP ratio, which is used +as a comparison in the case of any nutrient deficit. Senesced leaves also have fixed CN +and CP ratios, which is used for leaf turnover. + +In the future, the ideal CN and CP ratios will be PFT traits. +""" # noqa: D205 + +from dataclasses import dataclass, field + +import numpy as np +from numpy.typing import NDArray +from pyrealm.demography.core import CohortMethods, PandasExporter + +from virtual_ecosystem.models.plants.constants import PlantsConsts + + +@dataclass +class StemStochiometry(CohortMethods, PandasExporter): + """A class holding the ratios of Carbon to Nitrogen and Phosphorous for stems. + + This class holds the current ratios across tissue type for a community object, which + in essence is a series of cohorts. It acts in parallel with SteAllometry, a class + attribute of Community. + """ + + # Init vars + plant_constants: PlantsConsts + """The plant constants used in the model.""" + n_cohorts: int + """The number of cohorts in the community.""" + + # Post-init vars + cn_reproductive_tissue: NDArray[np.float64] = field(init=False) + """An array of carbon to nitrogen ratios for reproductive tissue.""" + cp_reproductive_tissue: NDArray[np.float64] = field(init=False) + """An array of carbon to phosphorous ratios for reproductive tissue.""" + + cn_foliage: NDArray[np.float64] = field(init=False) + """An array of carbon to nitrogen ratios for foliage.""" + cp_foliage: NDArray[np.float64] = field(init=False) + """An array of carbon to phosphorous ratios for foliage.""" + + cn_wood: NDArray[np.float64] = field(init=False) + """An array of carbon to nitrogen ratios for wood.""" + cp_wood: NDArray[np.float64] = field(init=False) + """An array of carbon to phosphorous ratios for wood.""" + + cn_roots: NDArray[np.float64] = field(init=False) + """An array of carbon to nitrogen ratios for roots.""" + cp_roots: NDArray[np.float64] = field(init=False) + """An array of carbon to phosphorous ratios for roots.""" + + def __post_init__( + self, + ) -> None: + """Initialise the stochiometry class. + + Args: + plant_constants: The plant constants used in the model. + len_cohort: The number of cohorts in the community. + """ + + # Initialise the arrays + self.cn_reproductive_tissue = np.full( + self.n_cohorts, + self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio, + ) + self.cp_reproductive_tissue = np.full( + self.n_cohorts, + self.plant_constants.plant_reproductive_tissue_turnover_c_p_ratio, + ) + + self.cn_foliage = np.full( + self.n_cohorts, self.plant_constants.foliage_c_n_ratio + ) + self.cp_foliage = np.full( + self.n_cohorts, self.plant_constants.foliage_c_p_ratio + ) + + self.cn_wood = np.full(self.n_cohorts, self.plant_constants.deadwood_c_n_ratio) + self.cp_wood = np.full(self.n_cohorts, self.plant_constants.deadwood_c_p_ratio) + + self.cn_roots = np.full( + self.n_cohorts, self.plant_constants.root_turnover_c_n_ratio + ) + self.cp_roots = np.full( + self.n_cohorts, self.plant_constants.root_turnover_c_p_ratio + ) From 6dc22b89e5199a36e1d7be2910e5b3bf3967e027 Mon Sep 17 00:00:00 2001 From: Sally Date: Tue, 29 Apr 2025 14:34:10 +0100 Subject: [PATCH 02/22] Moving stochiometry logic throuhgout the setup, estimate gpp, and allocate steps - DRAFT version --- .../models/plants/plants_model.py | 243 ++++++++++++++---- .../models/plants/stochiometry.py | 191 +++++++++++--- 2 files changed, 343 insertions(+), 91 deletions(-) diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index 90129101f..cb8f10b96 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -190,6 +190,8 @@ def __init__( to PlantCommunity instances for each cell.""" self.stochiometries: dict[int, StemStochiometry] """A dictionary keyed by cell id giving the stochiometry of each community.""" + self.allocations: dict[int, StemAllocation] + """A dictionary keyed by cell id giving the allocation of each community.""" self._canopy_layer_indices: NDArray[np.bool_] """The indices of the canopy layers within wider vertical profile. This is a shorter reference to self.layer_structure.index_canopy.""" @@ -269,6 +271,8 @@ def _setup( self.flora = flora self.model_constants = model_constants + self.temp_populate_symbiote_vars() + # Adjust flora turnover rates to timestep # TODO: Pyrealm provides annual turnover rates. Dividing by the number of # updates_per_year to get monthly turnover values is naive and will @@ -299,7 +303,9 @@ def _setup( self.stochiometries = { cell_id: StemStochiometry( plant_constants=self.model_constants, - n_cohorts=self.communities[cell_id].number_of_cohorts, + num_cohorts=self.communities[cell_id].number_of_cohorts, + stem_allometry=self.communities[cell_id].stem_allometry, + community=self.communities[cell_id], ) for cell_id in self.communities.keys() } @@ -613,6 +619,21 @@ def estimate_gpp(self, time_index: int) -> None: * n_seconds ).sum(axis=0) + # Calculate N uptake (g N per stem) due to transpiration. Conversion units: + # - Per stem transpiration (µmol H2O per stem) + # - Conversion factor from µmol H2O to m^3 (1.08015*10^-11) + # - Concentration of dissolved N (kg m^-3) + # - Kg to g (1000) + self.stochiometries[cell_id].n_surplus += ( + self.per_stem_transpiration[cell_id] + * (1.8015 * pow(10.0, -11)) + * ( + self.data["dissolved_ammonium"][cell_id] + + self.data["dissolved_nitrate"][cell_id] + ).item() + * 1000 + ) + # Calculate the total transpiration per layer in mm m2 in mm, converted from # an initial value is in µmol m2 s1.abs transpiration[active_layers, cell_id] = ( @@ -638,41 +659,34 @@ def allocate_gpp(self) -> None: turnover values. """ - # Allocate leaf and root turnover to per cell pools, merging across PFTs and - # cohorts. + # First, initialize all turnover variables to 0 with the proper shapes. + # Most variables are merged across PFTs and cohorts - one pool per cell. self.data["leaf_turnover"] = xr.full_like(self.data["elevation"], 0) self.data["root_turnover"] = xr.full_like(self.data["elevation"], 0) + self.data["root_carbohydrate_exudation"] = xr.full_like( + self.data["elevation"], 0 + ) + self.data["plant_symbiote_carbon_supply"] = xr.full_like( + self.data["elevation"], 0 + ) + self.data["fallen_non_propagule_c_mass"] = xr.full_like( + self.data["elevation"], 0 + ) - # Allocate reproductive tissue mass turnover - fallen propagules are stored per - # cell and per PFT, but fallen non-propagule reproductive tissue mass is merged - # into a single pool. + # Fallen propagules and canopy RT are stored per cell and per PFT. pft_cell_template = xr.DataArray( data=np.zeros((self.grid.n_cells, self.flora.n_pfts)), coords={"cell_id": self.data["cell_id"], "pft": self.flora.name}, ) - self.data["fallen_n_propagules"] = pft_cell_template.copy() - self.data["fallen_non_propagule_c_mass"] = xr.full_like( - self.data["elevation"], 0 - ) - - # Allocate canopy reproductive tissue mass. This is deliberately not - # partitioning tissue across canopy vertical layers. + # Canopy RT mass is deliberately not partitioned across canopy vertical layers. self.data["canopy_n_propagules"] = pft_cell_template.copy() self.data["canopy_non_propagule_c_mass"] = pft_cell_template.copy() - # Carbon supply to soil - self.data["root_carbohydrate_exudation"] = xr.full_like( - self.data["elevation"], 0 - ) - self.data["plant_symbiote_carbon_supply"] = xr.full_like( - self.data["elevation"], 0 - ) - - # Loop over each grid cell for cell_id in self.communities.keys(): community = self.communities[cell_id] cohorts = community.cohorts + stochiometry = self.stochiometries[cell_id] # Calculate the allocation of GPP per stem stem_allocation = StemAllocation( @@ -681,14 +695,12 @@ def allocate_gpp(self) -> None: at_potential_gpp=self.per_stem_gpp[cell_id], ) - # Grow the plants by increasing the stem dbh - # TODO: dimension mismatch (1d vs 2d array) - check in pyrealm - cohorts.dbh_values = cohorts.dbh_values + stem_allocation.delta_dbh - + # FIRST, ALLOCATE TO TURNOVER: # Sum of turnover from all cohorts in a grid cell self.data["leaf_turnover"][cell_id] = np.sum( stem_allocation.foliage_turnover * cohorts.n_individuals ) + self.data["root_turnover"][cell_id] = np.sum( stem_allocation.fine_root_turnover * cohorts.n_individuals ) @@ -746,6 +758,41 @@ def allocate_gpp(self) -> None: canopy_non_propagule_mass * cohort_n_stems ) + # SECOND: ALLOCATE N TO REGROW WHAT WAS LOST TO TURNOVER + # Nitrogen is lost from the tree in the form of turnover, and so an + # equivalent allocation of N is required to replace what was lost. To + # represent this process, N is allocated from the surplus store in the same + # quantities as turnover. + # TODO: is it correct to use the current ratios here? Should prob be a bit + # more complicated - need to replace N at the "ideal" ratio rate? + stochiometry.n_surplus = ( + stochiometry.n_surplus + - ( + stem_allocation.foliage_turnover + * (1 / self.model_constants.leaf_turnover_c_n_ratio) + ) + - ( + stem_allocation.fine_root_turnover + * (1 / stochiometry.cn_ratio_roots) + ) + - ( + stem_allocation.reproductive_tissue_turnover + * (1 / stochiometry.cn_ratio_reproductive_tissue) + ) + ) + + # N is lost in the form of senseced leaves, but some of the N emobodied in + # living leaves is returned to the surplus store before the leaves are shed. + stochiometry.n_surplus = ( + stochiometry.n_surplus + + stem_allocation.foliage_turnover + * ( + 1 / stochiometry.cn_ratio_foliage + - 1 / self.model_constants.leaf_turnover_c_n_ratio + ) + ) + + # THIRD, ALLOCATE GPP TO ACTIVE NUTRIENT PATHWAYS: # Allocate the topsliced GPP to root exudates with remainder as active # nutrient pathways self.data["root_carbohydrate_exudation"][cell_id] = np.sum( @@ -759,11 +806,90 @@ def allocate_gpp(self) -> None: * cohorts.n_individuals ) + # FOURTH, ALLOCATE TO GROWTH: + # Grow the plants by increasing the stem dbh + # TODO: dimension mismatch (1d vs 2d array) - check in pyrealm + cohorts.dbh_values = cohorts.dbh_values + stem_allocation.delta_dbh + + # Allocate N to growth and update stochiometry values + n_for_foliage_growth = ( + 1 / self.model_constants.foliage_c_n_ratio + ) * stem_allocation.delta_foliage_mass + n_for_stem_growth = ( + 1 / self.model_constants.deadwood_c_n_ratio + ) * stem_allocation.delta_stem_mass + n_for_rt_growth = ( + (1 / self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio) + * stem_allocation.delta_foliage_mass + * community.stem_traits.p_foliage_for_reproductive_tissue + ) + n_for_roots_growth = ( + (1 / self.model_constants.root_turnover_c_n_ratio) + * stem_allocation.delta_foliage_mass + * community.stem_traits.zeta + ) + + # Update the stochiometry values + stochiometry.n_surplus -= ( + n_for_foliage_growth + + n_for_stem_growth + + n_for_rt_growth + + n_for_roots_growth + ) + stochiometry.n_foliage = stochiometry.n_foliage + n_for_foliage_growth + stochiometry.n_wood = stochiometry.n_wood + n_for_stem_growth + stochiometry.n_reproductive_tissue = ( + stochiometry.n_reproductive_tissue + n_for_rt_growth + ) + stochiometry.n_roots = stochiometry.n_roots + n_for_roots_growth + + # Balance the N surplus/deficit with the symbiote carbon supply + self.data["plant_n_uptake_arbuscular"] = xr.full_like( + self.data["elevation"], 0 + ) + self.data["plant_n_uptake_ecto"] = xr.full_like(self.data["elevation"], 0) + n_available_from_symbiotes = ( + self.data["ecto_supply_limit_n"][cell_id] + + self.data["arbuscular_supply_limit_n"][cell_id] + ) + for cohort in range(0, len(stochiometry.n_surplus)): + if stochiometry.n_surplus[cohort] >= 0: + # Take no N from the symbiotes and keep the surplus + self.data["plant_n_uptake_arbuscular"][cell_id] += 0 + self.data["plant_n_uptake_ecto"][cell_id] += 0 + elif stochiometry.n_surplus[cohort] + n_available_from_symbiotes >= 0: + # Take N from arbuscular and ecto symbiotes in proportion to their + # supply, and reset the surplus to 0 + self.data["plant_n_uptake_arbuscular"] = stochiometry.n_surplus * ( + self.data["arbs_supply_limit_n"] + / ( + self.data["ecto_supply_limit_n"] + + self.data["arbs_supply_limit_n"] + ) + ) + self.data["plant_n_uptake_ecto"] = ( + stochiometry.n_surplus - self.data["plant_n_uptake_arbuscular"] + ) + stochiometry.n_surplus = 0 + else: + # NITROGEN DEFICIT!! + # Take all the N from the symbiotes, and decrease the deficit + self.data["plant_n_uptake_arbuscular"] = self.data[ + "arbuscular_supply_limit_n" + ] + self.data["plant_n_uptake_ecto"] = self.data["ecto_supply_limit_n"] + stochiometry.n_surplus = ( + stochiometry.n_surplus + n_available_from_symbiotes + ) + # TODO: what to do about the remaining deficit???? + # Update community allometry with new dbh values community.stem_allometry = StemAllometry( stem_traits=community.stem_traits, at_dbh=cohorts.dbh_values ) + self.allocations[cell_id] = stem_allocation + def apply_mortality(self) -> None: """Apply mortality to plant cohorts. @@ -795,6 +921,35 @@ def apply_mortality(self) -> None: mortality * community.stem_allometry.stem_mass ) + def update_cn_ratios(self) -> None: + # C:N and C:P ratios + self.data["deadwood_c_n_ratio"] = xr.full_like( + self.data["elevation"], self.model_constants.deadwood_c_n_ratio + ) + self.data["leaf_turnover_c_n_ratio"] = xr.full_like( + self.data["elevation"], self.model_constants.leaf_turnover_c_n_ratio + ) + self.data["plant_reproductive_tissue_turnover_c_n_ratio"] = xr.full_like( + self.data["elevation"], + self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio, + ) + self.data["root_turnover_c_n_ratio"] = xr.full_like( + self.data["elevation"], self.model_constants.root_turnover_c_n_ratio + ) + self.data["deadwood_c_p_ratio"] = xr.full_like( + self.data["elevation"], self.model_constants.deadwood_c_p_ratio + ) + self.data["leaf_turnover_c_p_ratio"] = xr.full_like( + self.data["elevation"], self.model_constants.leaf_turnover_c_p_ratio + ) + self.data["plant_reproductive_tissue_turnover_c_p_ratio"] = xr.full_like( + self.data["elevation"], + self.model_constants.plant_reproductive_tissue_turnover_c_p_ratio, + ) + self.data["root_turnover_c_p_ratio"] = xr.full_like( + self.data["elevation"], self.model_constants.root_turnover_c_p_ratio + ) + def calculate_turnover(self) -> None: """Calculate turnover of each plant biomass pool. @@ -828,33 +983,6 @@ def calculate_turnover(self) -> None: self.data["elevation"], self.model_constants.root_lignin ) - # C:N and C:P ratios - self.data["deadwood_c_n_ratio"] = xr.full_like( - self.data["elevation"], self.model_constants.deadwood_c_n_ratio - ) - self.data["leaf_turnover_c_n_ratio"] = xr.full_like( - self.data["elevation"], self.model_constants.leaf_turnover_c_n_ratio - ) - self.data["plant_reproductive_tissue_turnover_c_n_ratio"] = xr.full_like( - self.data["elevation"], - self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio, - ) - self.data["root_turnover_c_n_ratio"] = xr.full_like( - self.data["elevation"], self.model_constants.root_turnover_c_n_ratio - ) - self.data["deadwood_c_p_ratio"] = xr.full_like( - self.data["elevation"], self.model_constants.deadwood_c_p_ratio - ) - self.data["leaf_turnover_c_p_ratio"] = xr.full_like( - self.data["elevation"], self.model_constants.leaf_turnover_c_p_ratio - ) - self.data["plant_reproductive_tissue_turnover_c_p_ratio"] = xr.full_like( - self.data["elevation"], - self.model_constants.plant_reproductive_tissue_turnover_c_p_ratio, - ) - self.data["root_turnover_c_p_ratio"] = xr.full_like( - self.data["elevation"], self.model_constants.root_turnover_c_p_ratio - ) self.data["nitrogen_fixation_carbon_supply"] = xr.full_like( self.data["elevation"], 0.01 ) @@ -871,6 +999,7 @@ def calculate_nutrient_uptake(self) -> None: """ # Assume plants can take 0.1% of the available nutrient per day + # TODO: This is now set by transpiration? self.data["plant_ammonium_uptake"] = self.data["dissolved_ammonium"] * 0.01 self.data["plant_nitrate_uptake"] = self.data["dissolved_nitrate"] * 0.01 self.data["plant_phosphorus_uptake"] = self.data["dissolved_phosphorus"] * 0.01 @@ -897,3 +1026,11 @@ def partition_reproductive_tissue( ) return n_propagules, non_propagule_mass + + def temp_populate_symbiote_vars(self) -> None: + """Jacob has these vars in an open PR. Will delete before merging with develop.""" + + self.data["ecto_supply_limit_n"] = xr.full_like(self.data["elevation"], 1) + self.data["arbuscular_supply_limit_n"] = xr.full_like(self.data["elevation"], 1) + self.data["ecto_supply_limit_p"] = xr.full_like(self.data["elevation"], 1) + self.data["arbuscular_supply_limit_p"] = xr.full_like(self.data["elevation"], 1) diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index ba480f220..2340b0c18 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -15,7 +15,10 @@ import numpy as np from numpy.typing import NDArray +from pyrealm.demography.community import Community from pyrealm.demography.core import CohortMethods, PandasExporter +from pyrealm.demography.tmodel import StemAllometry, StemAllocation +from pyrealm.demography.community import Community from virtual_ecosystem.models.plants.constants import PlantsConsts @@ -32,63 +35,175 @@ class StemStochiometry(CohortMethods, PandasExporter): # Init vars plant_constants: PlantsConsts """The plant constants used in the model.""" - n_cohorts: int + num_cohorts: int """The number of cohorts in the community.""" + community: Community + """The Community object for the stochiometry.""" + stem_allometry: StemAllometry + """The StemAllometry object for the community.""" # Post-init vars - cn_reproductive_tissue: NDArray[np.float64] = field(init=False) - """An array of carbon to nitrogen ratios for reproductive tissue.""" - cp_reproductive_tissue: NDArray[np.float64] = field(init=False) - """An array of carbon to phosphorous ratios for reproductive tissue.""" - - cn_foliage: NDArray[np.float64] = field(init=False) - """An array of carbon to nitrogen ratios for foliage.""" - cp_foliage: NDArray[np.float64] = field(init=False) - """An array of carbon to phosphorous ratios for foliage.""" - - cn_wood: NDArray[np.float64] = field(init=False) - """An array of carbon to nitrogen ratios for wood.""" - cp_wood: NDArray[np.float64] = field(init=False) - """An array of carbon to phosphorous ratios for wood.""" - - cn_roots: NDArray[np.float64] = field(init=False) - """An array of carbon to nitrogen ratios for roots.""" - cp_roots: NDArray[np.float64] = field(init=False) - """An array of carbon to phosphorous ratios for roots.""" + n_reproductive_tissue: NDArray[np.float64] = field(init=False) + """Per stem Nitrogen mass of reproductive tissue for each cohort.""" + p_reproductive_tissue: NDArray[np.float64] = field(init=False) + """Per stem Phosphorous mass of reproductive tissue for each cohort.""" + + n_foliage: NDArray[np.float64] = field(init=False) + """Per stem Nitrogen mass of foliage for each cohort.""" + p_foliage: NDArray[np.float64] = field(init=False) + """Per stem Phosphorous mass of foliage for each cohort.""" + + n_wood: NDArray[np.float64] = field(init=False) + """Per stem Nitrogen mass of wood for each cohort.""" + p_wood: NDArray[np.float64] = field(init=False) + """Per stem Phosphorous mass of wood for each cohort.""" + + n_roots: NDArray[np.float64] = field(init=False) + """Per stem Nitrogen mass of roots for each cohort.""" + p_roots: NDArray[np.float64] = field(init=False) + """Per stem Phosphorous mass of roots for each cohort.""" + + n_surplus: NDArray[np.float64] = field(init=False) + """Per stem Nitrogen surplus (or deficit) for each cohort.""" + p_surplus: NDArray[np.float64] = field(init=False) + """Per stem Phosphorous surplus (or deficit) for each cohort.""" def __post_init__( self, ) -> None: """Initialise the stochiometry class. + TODO: Where do these nutrients come from? + Args: - plant_constants: The plant constants used in the model. - len_cohort: The number of cohorts in the community. + community: The Community object that parallels the Stochiometry. """ # Initialise the arrays - self.cn_reproductive_tissue = np.full( - self.n_cohorts, - self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio, + self.n_reproductive_tissue = np.full( + self.num_cohorts, + self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio + * self.stem_allometry.reproductive_tissue_mass, ) - self.cp_reproductive_tissue = np.full( - self.n_cohorts, - self.plant_constants.plant_reproductive_tissue_turnover_c_p_ratio, + self.p_reproductive_tissue = np.full( + self.num_cohorts, + self.plant_constants.plant_reproductive_tissue_turnover_c_p_ratio + * self.stem_allometry.reproductive_tissue_mass, ) - self.cn_foliage = np.full( - self.n_cohorts, self.plant_constants.foliage_c_n_ratio + self.n_foliage = np.full( + self.num_cohorts, + self.plant_constants.foliage_c_n_ratio * self.stem_allometry.foliage_mass, ) - self.cp_foliage = np.full( - self.n_cohorts, self.plant_constants.foliage_c_p_ratio + self.p_foliage = np.full( + self.num_cohorts, + self.plant_constants.foliage_c_p_ratio * self.stem_allometry.foliage_mass, ) - self.cn_wood = np.full(self.n_cohorts, self.plant_constants.deadwood_c_n_ratio) - self.cp_wood = np.full(self.n_cohorts, self.plant_constants.deadwood_c_p_ratio) + self.n_wood = np.full( + self.num_cohorts, + self.plant_constants.deadwood_c_n_ratio * self.stem_allometry.stem_mass, + ) + self.p_wood = np.full( + self.num_cohorts, + self.plant_constants.deadwood_c_p_ratio * self.stem_allometry.stem_mass, + ) - self.cn_roots = np.full( - self.n_cohorts, self.plant_constants.root_turnover_c_n_ratio + self.n_roots = np.full( + self.num_cohorts, + self.plant_constants.root_turnover_c_n_ratio + * self.community.stem_traits.zeta + * self.stem_allometry.foliage_mass, + ) + self.p_roots = np.full( + self.num_cohorts, + self.plant_constants.root_turnover_c_p_ratio + * self.community.stem_traits.zeta + * self.stem_allometry.foliage_mass, ) - self.cp_roots = np.full( - self.n_cohorts, self.plant_constants.root_turnover_c_p_ratio + + self.n_surplus = np.full(self.num_cohorts, 0.0) + self.p_surplus = np.full(self.num_cohorts, 0.0) + + def total_n(self) -> NDArray[np.float64]: + """Calculate the total nitrogen mass for each cohort. + + Returns: + The total nitrogen mass for each cohort. + """ + return self.n_foliage + self.n_wood + self.n_roots + self.n_reproductive_tissue + + def total_p(self) -> NDArray[np.float64]: + """Calculate the total phosphorous mass for each cohort. + + Returns: + The total phosphorous mass for each cohort. + """ + return self.p_foliage + self.p_wood + self.p_roots + self.p_reproductive_tissue + + def n_deficit(self) -> NDArray[np.float64]: + """Calculate the nitrogen deficit for each cohort. + + Returns: + The nitrogen deficit for each cohort. + """ + + ideal_n = ( + self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio + * self.stem_allometry.reproductive_tissue_mass + + self.plant_constants.foliage_c_n_ratio * self.stem_allometry.foliage_mass + + self.plant_constants.deadwood_c_n_ratio * self.stem_allometry.stem_mass + + ( + self.plant_constants.root_turnover_c_n_ratio + * self.community.stem_traits.zeta + * self.stem_allometry.foliage_mass + ) ) + + return ideal_n - self.total_n() + + def n_for_growth(self, allocation: StemAllocation) -> NDArray[np.float64]: + """Calculate the nitrogen required for growth for each cohort. + + Returns: + The nitrogen available for growth for each cohort. + """ + + n_needed_for_foliage_growth = allocation.delta_foliage_mass * ( + 1 / self.plant_consts["foliage_c_n_ratio"] + ) + + n_needed_for_stem_growth = allocation.delta_stem_mass * ( + 1 / self.plant_consts["deadwood_c_n_ratio"] + ) + + n_needed_for_rt_growth = ( + allocation.delta_foliage_mass + * (1 / self.plant_consts["plant_reproductive_tissue_turnover_c_n_ratio"]) + * self.stem_allometry.p_foliage_for_reproductive_tissue + ) + + n_needed_for_growth = ( + n_needed_for_foliage_growth + + n_needed_for_stem_growth + + n_needed_for_rt_growth + ) + + return n_needed_for_growth + + @property + def cn_ratio_foliage(self) -> NDArray[np.float64]: + """Get the carbon to nitrogen ratio for foliage.""" + return self.stem_allometry.foliage_mass / self.n_foliage + + @property + def cn_ratio_roots(self) -> NDArray[np.float64]: + """Get the carbon to nitrogen ratio for roots.""" + return ( + self.community.stem_traits.zeta * self.stem_allometry.foliage_mass + ) / self.n_roots + + @property + def cn_ratio_reproductive_tissue(self) -> NDArray[np.float64]: + """Get the carbon to nitrogen ratio for reproductive tissue.""" + return self.stem_allometry.reproductive_tissue_mass / self.n_reproductive_tissue From 9b8df282202f10cdfda7d6fc68bebfcead68c9a5 Mon Sep 17 00:00:00 2001 From: Sally Date: Thu, 8 May 2025 11:21:22 +0100 Subject: [PATCH 03/22] Draft of (N) stochiometry with tests passing. --- tests/models/plants/test_plants_model.py | 8 + virtual_ecosystem/main.py | 11 +- .../models/plants/plants_model.py | 161 +++++++------ .../models/plants/stochiometry.py | 222 ++++++++++++++---- 4 files changed, 281 insertions(+), 121 deletions(-) diff --git a/tests/models/plants/test_plants_model.py b/tests/models/plants/test_plants_model.py index 7a2691afc..12168032d 100644 --- a/tests/models/plants/test_plants_model.py +++ b/tests/models/plants/test_plants_model.py @@ -246,6 +246,14 @@ def test_PlantsModel_calculate_turnover(fxt_plants_model, fixture_config): ) assert np.allclose(fxt_plants_model.data["root_lignin"], consts.root_lignin) assert np.allclose(fxt_plants_model.data["leaf_lignin"], consts.leaf_lignin) + + +def test_PlantsModel_update_cn_ratios(fxt_plants_model, fixture_config): + """Test the update_cn_ratios method of the plants model.""" + + fxt_plants_model.update_cn_ratios() + consts = fxt_plants_model.model_constants + assert np.allclose( fxt_plants_model.data["deadwood_c_n_ratio"], consts.deadwood_c_n_ratio ) diff --git a/virtual_ecosystem/main.py b/virtual_ecosystem/main.py index f11309b2d..3bd238831 100644 --- a/virtual_ecosystem/main.py +++ b/virtual_ecosystem/main.py @@ -161,7 +161,7 @@ def ve_run( ) if progress: print("* Saved model inital state") - + print("checkpoint3") # If no path for saving continuous data is specified, fall back on using out_path if "out_folder_continuous" not in config["core"]["data_output_options"]: config["core"]["data_output_options"]["out_folder_continuous"] = str(out_path) @@ -188,32 +188,37 @@ def ve_run( time_index = 0 current_time = core_components.model_timing.start_time while current_time < core_components.model_timing.end_time: + print("starting update ", time_index) LOGGER.info(f"Starting update {time_index}: {current_time}") current_time += core_components.model_timing.update_interval # Run update() method for every model for model in models_update.values(): + print(model) LOGGER.info(f"Updating model {model.model_name}") model.update(time_index) + print("done with ", model.model_name) # With updates complete increment the time_index time_index += 1 + print("yo yo yo") # Append updated data to the continuous data file if config["core"]["data_output_options"]["save_continuous_data"]: + print("next checkpoint") outfile_path = data.output_current_state( variables_to_save, config["core"]["data_output_options"], time_index ) continuous_data_files.append(outfile_path) - + print("almost to end") pbar.update(n=1) pbar.close() if progress: print("* Simulation completed") - + print("checkpoint4") # Merge all files together based on a list if config["core"]["data_output_options"]["save_continuous_data"]: merge_continuous_data_files( diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index cb8f10b96..f9c478443 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -303,8 +303,6 @@ def _setup( self.stochiometries = { cell_id: StemStochiometry( plant_constants=self.model_constants, - num_cohorts=self.communities[cell_id].number_of_cohorts, - stem_allometry=self.communities[cell_id].stem_allometry, community=self.communities[cell_id], ) for cell_id in self.communities.keys() @@ -659,7 +657,7 @@ def allocate_gpp(self) -> None: turnover values. """ - # First, initialize all turnover variables to 0 with the proper shapes. + # First, initialize all turnover variables to 0 with the proper dimensions. # Most variables are merged across PFTs and cohorts - one pool per cell. self.data["leaf_turnover"] = xr.full_like(self.data["elevation"], 0) self.data["root_turnover"] = xr.full_like(self.data["elevation"], 0) @@ -674,12 +672,12 @@ def allocate_gpp(self) -> None: ) # Fallen propagules and canopy RT are stored per cell and per PFT. + # Canopy RT mass is deliberately not partitioned across canopy vertical layers. pft_cell_template = xr.DataArray( data=np.zeros((self.grid.n_cells, self.flora.n_pfts)), coords={"cell_id": self.data["cell_id"], "pft": self.flora.name}, ) self.data["fallen_n_propagules"] = pft_cell_template.copy() - # Canopy RT mass is deliberately not partitioned across canopy vertical layers. self.data["canopy_n_propagules"] = pft_cell_template.copy() self.data["canopy_non_propagule_c_mass"] = pft_cell_template.copy() @@ -762,34 +760,35 @@ def allocate_gpp(self) -> None: # Nitrogen is lost from the tree in the form of turnover, and so an # equivalent allocation of N is required to replace what was lost. To # represent this process, N is allocated from the surplus store in the same - # quantities as turnover. - # TODO: is it correct to use the current ratios here? Should prob be a bit - # more complicated - need to replace N at the "ideal" ratio rate? + # quantities as turnover. This uses the current ratios so that the tissue + # C:N ratios are maintained. stochiometry.n_surplus = ( stochiometry.n_surplus - ( - stem_allocation.foliage_turnover - * (1 / self.model_constants.leaf_turnover_c_n_ratio) - ) - - ( - stem_allocation.fine_root_turnover - * (1 / stochiometry.cn_ratio_roots) - ) - - ( - stem_allocation.reproductive_tissue_turnover - * (1 / stochiometry.cn_ratio_reproductive_tissue) - ) + ( + stem_allocation.foliage_turnover + * (1 / stochiometry.cn_ratio_foliage) + ) + - ( + stem_allocation.fine_root_turnover + * (1 / stochiometry.cn_ratio_roots) + ) + - ( + stem_allocation.reproductive_tissue_turnover + * (1 / stochiometry.cn_ratio_reproductive_tissue) + ) + ).squeeze() ) - # N is lost in the form of senseced leaves, but some of the N emobodied in # living leaves is returned to the surplus store before the leaves are shed. - stochiometry.n_surplus = ( - stochiometry.n_surplus - + stem_allocation.foliage_turnover - * ( - 1 / stochiometry.cn_ratio_foliage - - 1 / self.model_constants.leaf_turnover_c_n_ratio - ) + stochiometry.n_surplus = stochiometry.n_surplus + ( + ( + stem_allocation.foliage_turnover + * ( + 1 / stochiometry.cn_ratio_foliage + - 1 / self.model_constants.leaf_turnover_c_n_ratio + ) + ).squeeze() ) # THIRD, ALLOCATE GPP TO ACTIVE NUTRIENT PATHWAYS: @@ -800,11 +799,14 @@ def allocate_gpp(self) -> None: * self.model_constants.root_exudates * cohorts.n_individuals ) - self.data["plant_symbiote_carbon_supply"][cell_id] = np.sum( + carbon_to_symbiotes = ( stem_allocation.gpp_topslice * (1 - self.model_constants.root_exudates) * cohorts.n_individuals ) + self.data["plant_symbiote_carbon_supply"][cell_id] = np.sum( + carbon_to_symbiotes + ) # FOURTH, ALLOCATE TO GROWTH: # Grow the plants by increasing the stem dbh @@ -812,31 +814,37 @@ def allocate_gpp(self) -> None: cohorts.dbh_values = cohorts.dbh_values + stem_allocation.delta_dbh # Allocate N to growth and update stochiometry values - n_for_foliage_growth = ( - 1 / self.model_constants.foliage_c_n_ratio - ) * stem_allocation.delta_foliage_mass + n_for_stem_growth = ( 1 / self.model_constants.deadwood_c_n_ratio - ) * stem_allocation.delta_stem_mass + ) * stem_allocation.delta_stem_mass.squeeze() n_for_rt_growth = ( (1 / self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio) * stem_allocation.delta_foliage_mass * community.stem_traits.p_foliage_for_reproductive_tissue - ) + ).squeeze() n_for_roots_growth = ( (1 / self.model_constants.root_turnover_c_n_ratio) * stem_allocation.delta_foliage_mass * community.stem_traits.zeta - ) + ).squeeze() # Update the stochiometry values - stochiometry.n_surplus -= ( - n_for_foliage_growth - + n_for_stem_growth - + n_for_rt_growth - + n_for_roots_growth + stochiometry.n_surplus = ( + stochiometry.n_surplus + - ( + stochiometry.n_for_foliage_growth( + stem_allocation.delta_foliage_mass + ) + + n_for_stem_growth + + n_for_rt_growth + + n_for_roots_growth + ).squeeze() + ) + stochiometry.n_foliage = ( + stochiometry.n_foliage + + stochiometry.n_for_foliage_growth(stem_allocation.delta_foliage_mass) ) - stochiometry.n_foliage = stochiometry.n_foliage + n_for_foliage_growth stochiometry.n_wood = stochiometry.n_wood + n_for_stem_growth stochiometry.n_reproductive_tissue = ( stochiometry.n_reproductive_tissue + n_for_rt_growth @@ -844,51 +852,43 @@ def allocate_gpp(self) -> None: stochiometry.n_roots = stochiometry.n_roots + n_for_roots_growth # Balance the N surplus/deficit with the symbiote carbon supply - self.data["plant_n_uptake_arbuscular"] = xr.full_like( - self.data["elevation"], 0 - ) - self.data["plant_n_uptake_ecto"] = xr.full_like(self.data["elevation"], 0) - n_available_from_symbiotes = ( + n_weighted_avg = np.dot( self.data["ecto_supply_limit_n"][cell_id] - + self.data["arbuscular_supply_limit_n"][cell_id] + + self.data["arbuscular_supply_limit_n"][cell_id], + cohorts.n_individuals, ) - for cohort in range(0, len(stochiometry.n_surplus)): - if stochiometry.n_surplus[cohort] >= 0: - # Take no N from the symbiotes and keep the surplus - self.data["plant_n_uptake_arbuscular"][cell_id] += 0 - self.data["plant_n_uptake_ecto"][cell_id] += 0 - elif stochiometry.n_surplus[cohort] + n_available_from_symbiotes >= 0: - # Take N from arbuscular and ecto symbiotes in proportion to their - # supply, and reset the surplus to 0 - self.data["plant_n_uptake_arbuscular"] = stochiometry.n_surplus * ( - self.data["arbs_supply_limit_n"] - / ( - self.data["ecto_supply_limit_n"] - + self.data["arbs_supply_limit_n"] - ) - ) - self.data["plant_n_uptake_ecto"] = ( - stochiometry.n_surplus - self.data["plant_n_uptake_arbuscular"] - ) - stochiometry.n_surplus = 0 + n_avaiable_per_cohort = n_weighted_avg / sum(n_weighted_avg) + n_available_per_stem = np.divide( + n_avaiable_per_cohort, cohorts.n_individuals + ) + + stochiometry.n_surplus = stochiometry.n_surplus + n_available_per_stem + + # Cohort by cohort, distribute the surplus/deficit across the tissue types + for cohort in range(len(cohorts.n_individuals)): + if stochiometry.n_surplus[cohort] < 0: + # Distribute deficit across the tissue types + stochiometry.distribute_deficit(cohort) + + elif ( + stochiometry.n_surplus[cohort] > 0 + and stochiometry.n_tissue_deficit[cohort] > 0 + ): + # Distribute the surplus across the tissue types + stochiometry.distribute_surplus(cohort) + else: - # NITROGEN DEFICIT!! - # Take all the N from the symbiotes, and decrease the deficit - self.data["plant_n_uptake_arbuscular"] = self.data[ - "arbuscular_supply_limit_n" - ] - self.data["plant_n_uptake_ecto"] = self.data["ecto_supply_limit_n"] - stochiometry.n_surplus = ( - stochiometry.n_surplus + n_available_from_symbiotes - ) - # TODO: what to do about the remaining deficit???? + # NO ADJUSTMENT REQUIRED - there is a surplus in the store, but + # there is no deficit in the tissue types. + pass # Update community allometry with new dbh values community.stem_allometry = StemAllometry( stem_traits=community.stem_traits, at_dbh=cohorts.dbh_values ) - self.allocations[cell_id] = stem_allocation + # self.allocations[cell_id] = stem_allocation + self.update_cn_ratios() def apply_mortality(self) -> None: """Apply mortality to plant cohorts. @@ -922,6 +922,15 @@ def apply_mortality(self) -> None: ) def update_cn_ratios(self) -> None: + """Update the C:N and C:P ratios of plant tissues. + + This function updates the C:N and C:P ratios of various plant tissues, including + deadwood, leaf turnover, plant reproductive tissue turnover, and root turnover. + + Warning: + At present, this function just sets values to original constants. + """ + # C:N and C:P ratios self.data["deadwood_c_n_ratio"] = xr.full_like( self.data["elevation"], self.model_constants.deadwood_c_n_ratio @@ -1028,7 +1037,7 @@ def partition_reproductive_tissue( return n_propagules, non_propagule_mass def temp_populate_symbiote_vars(self) -> None: - """Jacob has these vars in an open PR. Will delete before merging with develop.""" + """Jacob has these vars in an open PR. Will delete before merging.""" self.data["ecto_supply_limit_n"] = xr.full_like(self.data["elevation"], 1) self.data["arbuscular_supply_limit_n"] = xr.full_like(self.data["elevation"], 1) diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index 2340b0c18..85f624c84 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -17,8 +17,7 @@ from numpy.typing import NDArray from pyrealm.demography.community import Community from pyrealm.demography.core import CohortMethods, PandasExporter -from pyrealm.demography.tmodel import StemAllometry, StemAllocation -from pyrealm.demography.community import Community +from pyrealm.demography.tmodel import StemAllocation from virtual_ecosystem.models.plants.constants import PlantsConsts @@ -35,12 +34,8 @@ class StemStochiometry(CohortMethods, PandasExporter): # Init vars plant_constants: PlantsConsts """The plant constants used in the model.""" - num_cohorts: int - """The number of cohorts in the community.""" community: Community """The Community object for the stochiometry.""" - stem_allometry: StemAllometry - """The StemAllometry object for the community.""" # Post-init vars n_reproductive_tissue: NDArray[np.float64] = field(init=False) @@ -81,57 +76,64 @@ def __post_init__( # Initialise the arrays self.n_reproductive_tissue = np.full( - self.num_cohorts, + self.community.number_of_cohorts, self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio - * self.stem_allometry.reproductive_tissue_mass, + * self.community.stem_allometry.reproductive_tissue_mass, ) self.p_reproductive_tissue = np.full( - self.num_cohorts, + self.community.number_of_cohorts, self.plant_constants.plant_reproductive_tissue_turnover_c_p_ratio - * self.stem_allometry.reproductive_tissue_mass, + * self.community.stem_allometry.reproductive_tissue_mass, ) self.n_foliage = np.full( - self.num_cohorts, - self.plant_constants.foliage_c_n_ratio * self.stem_allometry.foliage_mass, + self.community.number_of_cohorts, + self.plant_constants.foliage_c_n_ratio + * self.community.stem_allometry.foliage_mass, ) self.p_foliage = np.full( - self.num_cohorts, - self.plant_constants.foliage_c_p_ratio * self.stem_allometry.foliage_mass, + self.community.number_of_cohorts, + self.plant_constants.foliage_c_p_ratio + * self.community.stem_allometry.foliage_mass, ) self.n_wood = np.full( - self.num_cohorts, - self.plant_constants.deadwood_c_n_ratio * self.stem_allometry.stem_mass, + self.community.number_of_cohorts, + self.plant_constants.deadwood_c_n_ratio + * self.community.stem_allometry.stem_mass, ) self.p_wood = np.full( - self.num_cohorts, - self.plant_constants.deadwood_c_p_ratio * self.stem_allometry.stem_mass, + self.community.number_of_cohorts, + self.plant_constants.deadwood_c_p_ratio + * self.community.stem_allometry.stem_mass, ) self.n_roots = np.full( - self.num_cohorts, + self.community.number_of_cohorts, self.plant_constants.root_turnover_c_n_ratio * self.community.stem_traits.zeta - * self.stem_allometry.foliage_mass, + * self.community.stem_allometry.foliage_mass, ) self.p_roots = np.full( - self.num_cohorts, + self.community.number_of_cohorts, self.plant_constants.root_turnover_c_p_ratio * self.community.stem_traits.zeta - * self.stem_allometry.foliage_mass, + * self.community.stem_allometry.foliage_mass, ) - self.n_surplus = np.full(self.num_cohorts, 0.0) - self.p_surplus = np.full(self.num_cohorts, 0.0) + self.n_surplus = np.full(self.community.number_of_cohorts, 0.0) + self.p_surplus = np.full(self.community.number_of_cohorts, 0.0) + @property def total_n(self) -> NDArray[np.float64]: """Calculate the total nitrogen mass for each cohort. Returns: The total nitrogen mass for each cohort. """ - return self.n_foliage + self.n_wood + self.n_roots + self.n_reproductive_tissue + return ( + self.n_foliage + self.n_wood + self.n_roots + self.n_reproductive_tissue + ).squeeze() def total_p(self) -> NDArray[np.float64]: """Calculate the total phosphorous mass for each cohort. @@ -141,27 +143,74 @@ def total_p(self) -> NDArray[np.float64]: """ return self.p_foliage + self.p_wood + self.p_roots + self.p_reproductive_tissue - def n_deficit(self) -> NDArray[np.float64]: + @property + def n_foliage_deficit(self) -> NDArray[np.float64]: + """Calculate the nitrogen deficit for foliage. + + Returns: + The nitrogen deficit for foliage. + """ + return ( + self.plant_constants.foliage_c_n_ratio + * self.community.stem_allometry.foliage_mass + - self.n_foliage + ).squeeze() + + @property + def n_wood_deficit(self) -> NDArray[np.float64]: + """Calculate the nitrogen deficit for wood. + + Returns: + The nitrogen deficit for wood. + """ + return ( + self.plant_constants.deadwood_c_n_ratio + * self.community.stem_allometry.stem_mass + - self.n_wood + ).squeeze() + + @property + def n_roots_deficit(self) -> NDArray[np.float64]: + """Calculate the nitrogen deficit for roots. + + Returns: + The nitrogen deficit for roots. + """ + return ( + self.plant_constants.root_turnover_c_n_ratio + * self.community.stem_traits.zeta + * self.community.stem_allometry.foliage_mass + - self.n_roots + ).squeeze() + + @property + def n_reproductive_tissue_deficit(self) -> NDArray[np.float64]: + """Calculate the nitrogen deficit for reproductive tissue. + + Returns: + The nitrogen deficit for reproductive tissue. + """ + return ( + self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio + * self.community.stem_allometry.reproductive_tissue_mass + - self.n_reproductive_tissue + ).squeeze() + + @property + def n_tissue_deficit(self) -> NDArray[np.float64]: """Calculate the nitrogen deficit for each cohort. Returns: The nitrogen deficit for each cohort. """ - ideal_n = ( - self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio - * self.stem_allometry.reproductive_tissue_mass - + self.plant_constants.foliage_c_n_ratio * self.stem_allometry.foliage_mass - + self.plant_constants.deadwood_c_n_ratio * self.stem_allometry.stem_mass - + ( - self.plant_constants.root_turnover_c_n_ratio - * self.community.stem_traits.zeta - * self.stem_allometry.foliage_mass - ) + return ( + self.n_foliage_deficit + + self.n_wood_deficit + + self.n_roots_deficit + + self.n_reproductive_tissue_deficit ) - return ideal_n - self.total_n() - def n_for_growth(self, allocation: StemAllocation) -> NDArray[np.float64]: """Calculate the nitrogen required for growth for each cohort. @@ -180,7 +229,7 @@ def n_for_growth(self, allocation: StemAllocation) -> NDArray[np.float64]: n_needed_for_rt_growth = ( allocation.delta_foliage_mass * (1 / self.plant_consts["plant_reproductive_tissue_turnover_c_n_ratio"]) - * self.stem_allometry.p_foliage_for_reproductive_tissue + * self.community.stem_allometry.p_foliage_for_reproductive_tissue ) n_needed_for_growth = ( @@ -194,16 +243,105 @@ def n_for_growth(self, allocation: StemAllocation) -> NDArray[np.float64]: @property def cn_ratio_foliage(self) -> NDArray[np.float64]: """Get the carbon to nitrogen ratio for foliage.""" - return self.stem_allometry.foliage_mass / self.n_foliage + return self.community.stem_allometry.foliage_mass / self.n_foliage @property def cn_ratio_roots(self) -> NDArray[np.float64]: """Get the carbon to nitrogen ratio for roots.""" return ( - self.community.stem_traits.zeta * self.stem_allometry.foliage_mass + self.community.stem_traits.zeta * self.community.stem_allometry.foliage_mass ) / self.n_roots @property def cn_ratio_reproductive_tissue(self) -> NDArray[np.float64]: """Get the carbon to nitrogen ratio for reproductive tissue.""" - return self.stem_allometry.reproductive_tissue_mass / self.n_reproductive_tissue + return ( + self.community.stem_allometry.reproductive_tissue_mass + / self.n_reproductive_tissue + ) + + def n_for_foliage_growth(self, delta_foliage_mass) -> NDArray[np.float64]: + """Get the nitrogen needed for foliage growth.""" + return ( + (1 / self.plant_constants.foliage_c_n_ratio) * delta_foliage_mass + ).squeeze() + + def distribute_deficit(self, cohort: int) -> None: + """Distribute the nitrogen deficit across the tissue types. + + Args: + cohort: The cohort to reconcile deficit. + """ + deficit = self.n_surplus[cohort] * -1 + + self.n_foliage[cohort] = self.n_foliage[cohort] - ( + deficit * self.n_foliage[cohort] / self.total_n[cohort] + ) + + self.n_wood[cohort] = self.n_wood[cohort] - ( + deficit * self.n_wood[cohort] / self.total_n[cohort] + ) + self.n_roots[cohort] -= deficit * self.n_roots[cohort] / self.total_n[cohort] + self.n_reproductive_tissue[cohort] -= ( + deficit * self.n_reproductive_tissue[cohort] / self.total_n[cohort] + ) + self.n_surplus[cohort] += deficit + + def distribute_surplus(self, cohort: int) -> None: + """Distribute the nitrogen surplus across the tissue types for a single cohort. + + Args: + cohort: The cohort to reconcile surplus. + """ + if self.n_surplus[cohort] > self.n_tissue_deficit[cohort]: + # If there is sufficient surplus N to cover the existing deficit, the + # amount of the deficit is subtracted from the surplus which persists until + # the next update. All tissue types are updated to the ideal N ratios. + self.n_surplus[cohort] = ( + self.n_surplus[cohort] - self.n_tissue_deficit[cohort] + ) + self.n_foliage[cohort] = ( + self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio + * self.community.stem_allometry.foliage_mass.squeeze()[cohort] + ) + self.n_wood[cohort] = ( + # self.plant_constants.deadwood_c_n_ratio * + self.community.stem_allometry.stem_mass.squeeze()[cohort] + ) + self.n_roots[cohort] = ( + self.plant_constants.root_turnover_c_n_ratio + * self.community.stem_traits.zeta[cohort] + * self.community.stem_allometry.foliage_mass.squeeze()[cohort] + ) + self.n_reproductive_tissue[cohort] = ( + self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio + * self.community.stem_allometry.reproductive_tissue_mass.squeeze()[ + cohort + ] + ) + else: + # If there is not enough surplus N to cover the deficit, the surplus is + # distributed across the tissue types in proportion to the N deficit. + # The surplus is then set to zero. + + self.n_foliage[cohort] += ( + self.n_surplus[cohort] + * self.n_foliage_deficit[cohort] + / self.n_tissue_deficit[cohort] + ) + self.n_wood[cohort] += ( + self.n_surplus[cohort] + * self.n_wood_deficit[cohort] + / self.n_tissue_deficit[cohort] + ) + self.n_roots[cohort] += ( + self.n_surplus[cohort] + * self.n_roots_deficit[cohort] + / self.n_tissue_deficit[cohort] + ) + self.n_reproductive_tissue[cohort] += ( + self.n_surplus[cohort] + * self.n_reproductive_tissue_deficit[cohort] + / self.n_tissue_deficit[cohort] + ) + self.n_surplus[cohort] = 0.0 From 0c481ba5391f12e74fc89739ebfcc25507bf955e Mon Sep 17 00:00:00 2001 From: Sally Date: Thu, 8 May 2025 11:35:39 +0100 Subject: [PATCH 04/22] Deleting debug print statements. --- virtual_ecosystem/main.py | 10 ++-------- 1 file changed, 2 insertions(+), 8 deletions(-) diff --git a/virtual_ecosystem/main.py b/virtual_ecosystem/main.py index 3bd238831..f17adeef1 100644 --- a/virtual_ecosystem/main.py +++ b/virtual_ecosystem/main.py @@ -161,7 +161,7 @@ def ve_run( ) if progress: print("* Saved model inital state") - print("checkpoint3") + # If no path for saving continuous data is specified, fall back on using out_path if "out_folder_continuous" not in config["core"]["data_output_options"]: config["core"]["data_output_options"]["out_folder_continuous"] = str(out_path) @@ -187,38 +187,32 @@ def ve_run( pbar = tqdm(total=core_components.model_timing.n_updates) time_index = 0 current_time = core_components.model_timing.start_time + while current_time < core_components.model_timing.end_time: - print("starting update ", time_index) LOGGER.info(f"Starting update {time_index}: {current_time}") current_time += core_components.model_timing.update_interval # Run update() method for every model for model in models_update.values(): - print(model) LOGGER.info(f"Updating model {model.model_name}") model.update(time_index) - print("done with ", model.model_name) # With updates complete increment the time_index time_index += 1 - print("yo yo yo") # Append updated data to the continuous data file if config["core"]["data_output_options"]["save_continuous_data"]: - print("next checkpoint") outfile_path = data.output_current_state( variables_to_save, config["core"]["data_output_options"], time_index ) continuous_data_files.append(outfile_path) - print("almost to end") pbar.update(n=1) pbar.close() if progress: print("* Simulation completed") - print("checkpoint4") # Merge all files together based on a list if config["core"]["data_output_options"]["save_continuous_data"]: merge_continuous_data_files( From 53e1ec6ce5fee7cae801039082e5eecffe5b5bfd Mon Sep 17 00:00:00 2001 From: Sally Date: Thu, 8 May 2025 11:37:18 +0100 Subject: [PATCH 05/22] Returning to original main.py --- virtual_ecosystem/main.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/virtual_ecosystem/main.py b/virtual_ecosystem/main.py index f17adeef1..f11309b2d 100644 --- a/virtual_ecosystem/main.py +++ b/virtual_ecosystem/main.py @@ -187,7 +187,6 @@ def ve_run( pbar = tqdm(total=core_components.model_timing.n_updates) time_index = 0 current_time = core_components.model_timing.start_time - while current_time < core_components.model_timing.end_time: LOGGER.info(f"Starting update {time_index}: {current_time}") @@ -207,12 +206,14 @@ def ve_run( variables_to_save, config["core"]["data_output_options"], time_index ) continuous_data_files.append(outfile_path) + pbar.update(n=1) pbar.close() if progress: print("* Simulation completed") + # Merge all files together based on a list if config["core"]["data_output_options"]["save_continuous_data"]: merge_continuous_data_files( From b5d599b0232245187c65de01ed47409c96a46a2e Mon Sep 17 00:00:00 2001 From: Sally Date: Tue, 3 Jun 2025 15:15:04 +0100 Subject: [PATCH 06/22] Updated stochiometry with abstracted class structure. --- tests/models/plants/test_stochiometry.py | 151 +++++ .../models/plants/plants_model.py | 147 +++-- .../models/plants/stochiometry.py | 540 ++++++++++-------- 3 files changed, 530 insertions(+), 308 deletions(-) create mode 100644 tests/models/plants/test_stochiometry.py diff --git a/tests/models/plants/test_stochiometry.py b/tests/models/plants/test_stochiometry.py new file mode 100644 index 000000000..812d5c903 --- /dev/null +++ b/tests/models/plants/test_stochiometry.py @@ -0,0 +1,151 @@ +"""Tests for the Stochiometry class.""" + +import numpy as np + + +def test_Stochiometry__init__(fxt_plants_model): + """Fixture for the Stochiometry class.""" + + from virtual_ecosystem.models.plants.stochiometry import ( + FoliageTissue, + ReproductiveTissue, + RootTissue, + StemStochiometry, + WoodTissue, + ) + + plant_consts = fxt_plants_model.model_constants + for cell_id in fxt_plants_model.communities.keys(): + community = fxt_plants_model.communities[cell_id] + + n_stochiometry = StemStochiometry( + tissues=[ + FoliageTissue( + element_name="Nitrogen", + community=community, + ideal_ratio=np.full( + community.number_of_cohorts, plant_consts.foliage_c_n_ratio + ), + actual_element_mass=community.stem_allometry.foliage_mass + * plant_consts.foliage_c_n_ratio, + reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, + ), + RootTissue( + element_name="Nitrogen", + community=community, + ideal_ratio=np.full( + community.number_of_cohorts, + plant_consts.root_turnover_c_n_ratio, + ), + actual_element_mass=plant_consts.root_turnover_c_n_ratio + * community.stem_traits.zeta + * community.stem_allometry.foliage_mass, + ), + WoodTissue( + element_name="Nitrogen", + community=community, + ideal_ratio=np.full( + community.number_of_cohorts, plant_consts.deadwood_c_n_ratio + ), + actual_element_mass=plant_consts.deadwood_c_n_ratio + * community.stem_allometry.stem_mass, + ), + ReproductiveTissue( + element_name="Nitrogen", + community=community, + ideal_ratio=np.full( + community.number_of_cohorts, + plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, + ), + actual_element_mass=community.stem_allometry.reproductive_tissue_mass + * plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, + ), + ], + n_cohorts=community.number_of_cohorts, + community=community, + ) + assert n_stochiometry is not None + + +def test_Stochiometry_FoliageTissue(fxt_plants_model): + """Test the foliage stochiometry.""" + from virtual_ecosystem.models.plants.stochiometry import FoliageTissue + + plant_consts = fxt_plants_model.model_constants + for cell_id in fxt_plants_model.communities.keys(): + community = fxt_plants_model.communities[cell_id] + + foliage_tissue = FoliageTissue( + element_name="Nitrogen", + community=community, + ideal_ratio=np.full( + community.number_of_cohorts, plant_consts.foliage_c_n_ratio + ), + actual_element_mass=community.stem_allometry.foliage_mass + * plant_consts.foliage_c_n_ratio, + reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, + ) + assert foliage_tissue is not None + + +def test_Stochiometry_ReproductiveTissue(fxt_plants_model): + """Test the reproductive tissue stochiometry.""" + from virtual_ecosystem.models.plants.stochiometry import ReproductiveTissue + + plant_consts = fxt_plants_model.model_constants + for cell_id in fxt_plants_model.communities.keys(): + community = fxt_plants_model.communities[cell_id] + + reproductive_tissue = ReproductiveTissue( + element_name="Nitrogen", + community=community, + ideal_ratio=np.full( + community.number_of_cohorts, + plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, + ), + actual_element_mass=community.stem_allometry.reproductive_tissue_mass + * plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, + ) + assert reproductive_tissue is not None + + +def test_Stochiometry_RootTissue(fxt_plants_model): + """Test the reproductive tissue stochiometry.""" + from virtual_ecosystem.models.plants.stochiometry import RootTissue + + plant_consts = fxt_plants_model.model_constants + for cell_id in fxt_plants_model.communities.keys(): + community = fxt_plants_model.communities[cell_id] + + root_tissue = RootTissue( + element_name="Nitrogen", + community=community, + ideal_ratio=np.full( + community.number_of_cohorts, plant_consts.root_turnover_c_n_ratio + ), + actual_element_mass=plant_consts.root_turnover_c_n_ratio + * community.stem_traits.zeta + * community.stem_allometry.foliage_mass, + ) + + assert root_tissue is not None + + +def test_Stochiometry_WoodTissue(fxt_plants_model): + """Test the reproductive tissue stochiometry.""" + from virtual_ecosystem.models.plants.stochiometry import WoodTissue + + plant_consts = fxt_plants_model.model_constants + for cell_id in fxt_plants_model.communities.keys(): + community = fxt_plants_model.communities[cell_id] + + wood_tissue = WoodTissue( + element_name="Nitrogen", + community=community, + ideal_ratio=np.full( + community.number_of_cohorts, plant_consts.deadwood_c_n_ratio + ), + actual_element_mass=plant_consts.deadwood_c_n_ratio + * community.stem_allometry.stem_mass, + ) + assert wood_tissue is not None diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index f9c478443..3d446378f 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -29,7 +29,13 @@ from virtual_ecosystem.models.plants.communities import PlantCommunities from virtual_ecosystem.models.plants.constants import PlantsConsts from virtual_ecosystem.models.plants.functional_types import get_flora_from_config -from virtual_ecosystem.models.plants.stochiometry import StemStochiometry +from virtual_ecosystem.models.plants.stochiometry import ( + FoliageTissue, + ReproductiveTissue, + RootTissue, + StemStochiometry, + WoodTissue, +) class PlantsModel( @@ -302,7 +308,58 @@ def _setup( # Initialize the stochiometries of each cohort self.stochiometries = { cell_id: StemStochiometry( - plant_constants=self.model_constants, + tissues=[ + FoliageTissue( + element_name="Nitrogen", + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].number_of_cohorts, + model_constants.foliage_c_n_ratio, + ), + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.foliage_mass + * model_constants.foliage_c_n_ratio, + reclaim_ratio=np.full( + self.communities[cell_id].number_of_cohorts, + model_constants.leaf_turnover_c_n_ratio, + ), + ), + RootTissue( + element_name="Nitrogen", + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].number_of_cohorts, + model_constants.root_turnover_c_n_ratio, + ), + actual_element_mass=model_constants.root_turnover_c_n_ratio + * self.communities[cell_id].stem_traits.zeta + * self.communities[cell_id].stem_allometry.foliage_mass, + ), + WoodTissue( + element_name="Nitrogen", + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].number_of_cohorts, + model_constants.deadwood_c_n_ratio, + ), + actual_element_mass=model_constants.deadwood_c_n_ratio + * self.communities[cell_id].stem_allometry.stem_mass, + ), + ReproductiveTissue( + element_name="Nitrogen", + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].number_of_cohorts, + model_constants.plant_reproductive_tissue_turnover_c_n_ratio, + ), + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.reproductive_tissue_mass + * self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio, # noqa: E501 + ), + ], + n_cohorts=self.communities[cell_id].number_of_cohorts, community=self.communities[cell_id], ) for cell_id in self.communities.keys() @@ -622,7 +679,7 @@ def estimate_gpp(self, time_index: int) -> None: # - Conversion factor from µmol H2O to m^3 (1.08015*10^-11) # - Concentration of dissolved N (kg m^-3) # - Kg to g (1000) - self.stochiometries[cell_id].n_surplus += ( + self.stochiometries[cell_id].element_surplus += ( self.per_stem_transpiration[cell_id] * (1.8015 * pow(10.0, -11)) * ( @@ -757,39 +814,7 @@ def allocate_gpp(self) -> None: ) # SECOND: ALLOCATE N TO REGROW WHAT WAS LOST TO TURNOVER - # Nitrogen is lost from the tree in the form of turnover, and so an - # equivalent allocation of N is required to replace what was lost. To - # represent this process, N is allocated from the surplus store in the same - # quantities as turnover. This uses the current ratios so that the tissue - # C:N ratios are maintained. - stochiometry.n_surplus = ( - stochiometry.n_surplus - - ( - ( - stem_allocation.foliage_turnover - * (1 / stochiometry.cn_ratio_foliage) - ) - - ( - stem_allocation.fine_root_turnover - * (1 / stochiometry.cn_ratio_roots) - ) - - ( - stem_allocation.reproductive_tissue_turnover - * (1 / stochiometry.cn_ratio_reproductive_tissue) - ) - ).squeeze() - ) - # N is lost in the form of senseced leaves, but some of the N emobodied in - # living leaves is returned to the surplus store before the leaves are shed. - stochiometry.n_surplus = stochiometry.n_surplus + ( - ( - stem_allocation.foliage_turnover - * ( - 1 / stochiometry.cn_ratio_foliage - - 1 / self.model_constants.leaf_turnover_c_n_ratio - ) - ).squeeze() - ) + stochiometry.account_for_element_loss_turnover(stem_allocation) # THIRD, ALLOCATE GPP TO ACTIVE NUTRIENT PATHWAYS: # Allocate the topsliced GPP to root exudates with remainder as active @@ -813,43 +838,9 @@ def allocate_gpp(self) -> None: # TODO: dimension mismatch (1d vs 2d array) - check in pyrealm cohorts.dbh_values = cohorts.dbh_values + stem_allocation.delta_dbh - # Allocate N to growth and update stochiometry values - - n_for_stem_growth = ( - 1 / self.model_constants.deadwood_c_n_ratio - ) * stem_allocation.delta_stem_mass.squeeze() - n_for_rt_growth = ( - (1 / self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio) - * stem_allocation.delta_foliage_mass - * community.stem_traits.p_foliage_for_reproductive_tissue - ).squeeze() - n_for_roots_growth = ( - (1 / self.model_constants.root_turnover_c_n_ratio) - * stem_allocation.delta_foliage_mass - * community.stem_traits.zeta - ).squeeze() - - # Update the stochiometry values - stochiometry.n_surplus = ( - stochiometry.n_surplus - - ( - stochiometry.n_for_foliage_growth( - stem_allocation.delta_foliage_mass - ) - + n_for_stem_growth - + n_for_rt_growth - + n_for_roots_growth - ).squeeze() - ) - stochiometry.n_foliage = ( - stochiometry.n_foliage - + stochiometry.n_for_foliage_growth(stem_allocation.delta_foliage_mass) - ) - stochiometry.n_wood = stochiometry.n_wood + n_for_stem_growth - stochiometry.n_reproductive_tissue = ( - stochiometry.n_reproductive_tissue + n_for_rt_growth - ) - stochiometry.n_roots = stochiometry.n_roots + n_for_roots_growth + # Subtract the N/P required from growth from the element store, and + # redistribute it to the individual tisuses. + stochiometry.account_for_growth(stem_allocation) # Balance the N surplus/deficit with the symbiote carbon supply n_weighted_avg = np.dot( @@ -862,17 +853,19 @@ def allocate_gpp(self) -> None: n_avaiable_per_cohort, cohorts.n_individuals ) - stochiometry.n_surplus = stochiometry.n_surplus + n_available_per_stem + stochiometry.element_surplus = ( + stochiometry.element_surplus + n_available_per_stem + ) # Cohort by cohort, distribute the surplus/deficit across the tissue types for cohort in range(len(cohorts.n_individuals)): - if stochiometry.n_surplus[cohort] < 0: + if stochiometry.element_surplus[cohort] < 0: # Distribute deficit across the tissue types stochiometry.distribute_deficit(cohort) elif ( - stochiometry.n_surplus[cohort] > 0 - and stochiometry.n_tissue_deficit[cohort] > 0 + stochiometry.element_surplus[cohort] > 0 + and stochiometry.tissue_deficit[cohort] > 0 ): # Distribute the surplus across the tissue types stochiometry.distribute_surplus(cohort) diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index 85f624c84..46698c322 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -19,252 +19,363 @@ from pyrealm.demography.core import CohortMethods, PandasExporter from pyrealm.demography.tmodel import StemAllocation -from virtual_ecosystem.models.plants.constants import PlantsConsts - @dataclass -class StemStochiometry(CohortMethods, PandasExporter): - """A class holding the ratios of Carbon to Nitrogen and Phosphorous for stems. +class Tissue(PandasExporter, CohortMethods): + """A dataclass to hold tissue stochiometry data for a set of plant cohorts. - This class holds the current ratios across tissue type for a community object, which - in essence is a series of cohorts. It acts in parallel with SteAllometry, a class - attribute of Community. + This class holds the current quantitiy of a given element (generally N or P) for a + specific plant tissue type (generally foliage, wood, roots or reproductive tissue). + The class also holds the ideal ratio of the element for that tissue type. They hold + an entry for each cohort in the data class. """ - # Init vars - plant_constants: PlantsConsts - """The plant constants used in the model.""" + element_name: str + """The name of the element type.""" + tissue_name: str + """The name of the tissue type.""" community: Community - """The Community object for the stochiometry.""" - - # Post-init vars - n_reproductive_tissue: NDArray[np.float64] = field(init=False) - """Per stem Nitrogen mass of reproductive tissue for each cohort.""" - p_reproductive_tissue: NDArray[np.float64] = field(init=False) - """Per stem Phosphorous mass of reproductive tissue for each cohort.""" - - n_foliage: NDArray[np.float64] = field(init=False) - """Per stem Nitrogen mass of foliage for each cohort.""" - p_foliage: NDArray[np.float64] = field(init=False) - """Per stem Phosphorous mass of foliage for each cohort.""" - - n_wood: NDArray[np.float64] = field(init=False) - """Per stem Nitrogen mass of wood for each cohort.""" - p_wood: NDArray[np.float64] = field(init=False) - """Per stem Phosphorous mass of wood for each cohort.""" - - n_roots: NDArray[np.float64] = field(init=False) - """Per stem Nitrogen mass of roots for each cohort.""" - p_roots: NDArray[np.float64] = field(init=False) - """Per stem Phosphorous mass of roots for each cohort.""" - - n_surplus: NDArray[np.float64] = field(init=False) - """Per stem Nitrogen surplus (or deficit) for each cohort.""" - p_surplus: NDArray[np.float64] = field(init=False) - """Per stem Phosphorous surplus (or deficit) for each cohort.""" - - def __post_init__( - self, - ) -> None: - """Initialise the stochiometry class. + """The community object that the tissue is associated with.""" + # Should this be stored only in stochiometry and not in tissue? - TODO: Where do these nutrients come from? + ideal_ratio: NDArray[np.float64] + """The ideal ratio of the element for the tissue type.""" + actual_element_mass: NDArray[np.float64] + """The actual mass of the element for the tissue type.""" - Args: - community: The Community object that parallels the Stochiometry. + def __post_init__(self) -> None: + """Initialize the Tissue object.""" + self.actual_element_mass = self.actual_element_mass.squeeze() + + @property + def carbon_mass(self) -> NDArray[np.float64]: + """Get the carbon mass for the tissue type. + + This method should be implemented by subclasses to return the carbon mass for + the specific tissue type. + + Returns: + The carbon mass for the specified tissue. """ + # This method should be implemented by subclasses + raise NotImplementedError("Carbon mass must be defined in subclasses.") - # Initialise the arrays - self.n_reproductive_tissue = np.full( - self.community.number_of_cohorts, - self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio - * self.community.stem_allometry.reproductive_tissue_mass, - ) - self.p_reproductive_tissue = np.full( - self.community.number_of_cohorts, - self.plant_constants.plant_reproductive_tissue_turnover_c_p_ratio - * self.community.stem_allometry.reproductive_tissue_mass, - ) + @property + def deficit(self) -> NDArray[np.float64]: + """Calculate the element deficit for the tissue type. - self.n_foliage = np.full( - self.community.number_of_cohorts, - self.plant_constants.foliage_c_n_ratio - * self.community.stem_allometry.foliage_mass, - ) - self.p_foliage = np.full( - self.community.number_of_cohorts, - self.plant_constants.foliage_c_p_ratio - * self.community.stem_allometry.foliage_mass, - ) + Returns: + The element deficit for the specified tissue. + """ + return self.ideal_ratio * self.carbon_mass - self.actual_element_mass - self.n_wood = np.full( - self.community.number_of_cohorts, - self.plant_constants.deadwood_c_n_ratio - * self.community.stem_allometry.stem_mass, - ) - self.p_wood = np.full( - self.community.number_of_cohorts, - self.plant_constants.deadwood_c_p_ratio - * self.community.stem_allometry.stem_mass, - ) + def element_needed_for_growth( + self, allocation: StemAllocation + ) -> NDArray[np.float64]: + """Calculate the element needed for growth for the tissue type. - self.n_roots = np.full( - self.community.number_of_cohorts, - self.plant_constants.root_turnover_c_n_ratio - * self.community.stem_traits.zeta - * self.community.stem_allometry.foliage_mass, + Returns: + The element needed for growth for the specified tissue. + """ + raise NotImplementedError( + "Element needed for growth must be defined in subclasses." ) - self.p_roots = np.full( - self.community.number_of_cohorts, - self.plant_constants.root_turnover_c_p_ratio - * self.community.stem_traits.zeta - * self.community.stem_allometry.foliage_mass, + + def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: + """Calculate the element lost to turnover for the tissue type. + + Returns: + The element lost to turnover for the specified tissue. + """ + raise NotImplementedError( + "Element needed for growth must be defined in subclasses." ) - self.n_surplus = np.full(self.community.number_of_cohorts, 0.0) - self.p_surplus = np.full(self.community.number_of_cohorts, 0.0) + @property + def Cx_ratio(self) -> NDArray[np.float64]: + """Get the carbon to element ratio for the tissue type. + + Returns: + The carbon to element ratio for the specified tissue. + """ + return self.carbon_mass / self.actual_element_mass + + +class FoliageTissue(Tissue): + """A class to hold foliage stochiometry data for a set of plant cohorts.""" + + # reclaim_ratio: NDArray[np.float64] + """The ratio of the element that can be r eclaimed from the sensced tissue.""" + + def __init__( + self, + element_name: str, + community: Community, + ideal_ratio: NDArray[np.float64], + actual_element_mass: NDArray[np.float64], + reclaim_ratio: NDArray[np.float64], + ): + super().__init__( + element_name=element_name, + tissue_name="Foliage", + community=community, + ideal_ratio=ideal_ratio, + actual_element_mass=actual_element_mass, + ) + self.reclaim_ratio = reclaim_ratio + """The ratio of the element that can be reclaimed from the senesced tissue.""" @property - def total_n(self) -> NDArray[np.float64]: - """Calculate the total nitrogen mass for each cohort. + def carbon_mass(self) -> NDArray[np.float64]: + """Get the carbon mass for foliage tissue. Returns: - The total nitrogen mass for each cohort. + The carbon mass for foliage tissue. """ - return ( - self.n_foliage + self.n_wood + self.n_roots + self.n_reproductive_tissue - ).squeeze() + return self.community.stem_allometry.foliage_mass.squeeze() - def total_p(self) -> NDArray[np.float64]: - """Calculate the total phosphorous mass for each cohort. + def element_needed_for_growth( + self, allocation: StemAllocation + ) -> NDArray[np.float64]: + """Calculate the nitrogen needed for growth for foliage tissue. Returns: - The total phosphorous mass for each cohort. + The nitrogen needed for growth for foliage tissue. """ - return self.p_foliage + self.p_wood + self.p_roots + self.p_reproductive_tissue + return (allocation.delta_foliage_mass * (1 / self.ideal_ratio)).squeeze() - @property - def n_foliage_deficit(self) -> NDArray[np.float64]: - """Calculate the nitrogen deficit for foliage. + def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: + """Calculate the element mass lost to turnover for foliage tissue. Returns: - The nitrogen deficit for foliage. + The nitrogen lost to turnover for foliage tissue. """ return ( - self.plant_constants.foliage_c_n_ratio - * self.community.stem_allometry.foliage_mass - - self.n_foliage + allocation.foliage_turnover + * ((1 / self.reclaim_ratio) - (1 / self.Cx_ratio)) ).squeeze() + +class ReproductiveTissue(Tissue): + """Holds reproductive tissue stochiometry data for a set of plant cohorts.""" + + def __init__( + self, + element_name: str, + community: Community, + ideal_ratio: NDArray[np.float64], + actual_element_mass: NDArray[np.float64], + ): + super().__init__( + element_name=element_name, + tissue_name="Reproductive", + community=community, + ideal_ratio=ideal_ratio, + actual_element_mass=actual_element_mass, + ) + @property - def n_wood_deficit(self) -> NDArray[np.float64]: - """Calculate the nitrogen deficit for wood. + def carbon_mass(self) -> NDArray[np.float64]: + """Get the carbon mass for reproductive tissue. + + Returns: + The carbon mass for reproductive tissue. + """ + return self.community.stem_allometry.reproductive_tissue_mass.squeeze() + + def element_needed_for_growth( + self, allocation: StemAllocation + ) -> NDArray[np.float64]: + """Calculate the nitrogen needed for growth for reproductive tissue. Returns: - The nitrogen deficit for wood. + The nitrogen needed for growth for reproductive tissue. """ return ( - self.plant_constants.deadwood_c_n_ratio - * self.community.stem_allometry.stem_mass - - self.n_wood + allocation.delta_foliage_mass + * (1 / self.ideal_ratio) + * self.community.stem_traits.p_foliage_for_reproductive_tissue ).squeeze() + def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: + """Calculate the element lost to turnover for reproductive tissue. + + Returns: + The element lost to turnover for reproductive tissue. + """ + return (allocation.reproductive_tissue_turnover * (1 / self.Cx_ratio)).squeeze() + + +class WoodTissue(Tissue): + """A class to hold wood stochiometry data for a set of plant cohorts.""" + + def __init__( + self, + element_name: str, + community: Community, + ideal_ratio: NDArray[np.float64], + actual_element_mass: NDArray[np.float64], + ): + super().__init__( + element_name=element_name, + tissue_name="Wood", + community=community, + ideal_ratio=ideal_ratio, + actual_element_mass=actual_element_mass, + ) + @property - def n_roots_deficit(self) -> NDArray[np.float64]: - """Calculate the nitrogen deficit for roots. + def carbon_mass(self) -> NDArray[np.float64]: + """Get the carbon mass for wood tissue. Returns: - The nitrogen deficit for roots. + The carbon mass for wood tissue. """ - return ( - self.plant_constants.root_turnover_c_n_ratio - * self.community.stem_traits.zeta - * self.community.stem_allometry.foliage_mass - - self.n_roots - ).squeeze() + return self.community.stem_allometry.stem_mass.squeeze() + + def element_needed_for_growth( + self, allocation: StemAllocation + ) -> NDArray[np.float64]: + """Calculate the nitrogen needed for growth for wood tissue. + + Returns: + The nitrogen needed for growth for wood tissue. + """ + return (allocation.delta_stem_mass * (1 / self.ideal_ratio)).squeeze() + + def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: + """Assume no wood tissue is lost. + + Returns: + The element lost to turnover for wood tissue. + """ + return np.zeros(self.community.number_of_cohorts) + + +class RootTissue(Tissue): + """A class to hold root stochiometry data for a set of plant cohorts.""" + + def __init__( + self, + element_name: str, + community: Community, + ideal_ratio: NDArray[np.float64], + actual_element_mass: NDArray[np.float64], + ): + super().__init__( + element_name=element_name, + tissue_name="Roots", + community=community, + ideal_ratio=ideal_ratio, + actual_element_mass=actual_element_mass, + ) @property - def n_reproductive_tissue_deficit(self) -> NDArray[np.float64]: - """Calculate the nitrogen deficit for reproductive tissue. + def carbon_mass(self) -> NDArray[np.float64]: + """Get the carbon mass for root tissue. Returns: - The nitrogen deficit for reproductive tissue. + The carbon mass for root tissue. """ return ( - self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio - * self.community.stem_allometry.reproductive_tissue_mass - - self.n_reproductive_tissue + self.community.stem_allometry.foliage_mass * self.community.stem_traits.zeta ).squeeze() - @property - def n_tissue_deficit(self) -> NDArray[np.float64]: - """Calculate the nitrogen deficit for each cohort. + def element_needed_for_growth( + self, allocation: StemAllocation + ) -> NDArray[np.float64]: + """Calculate the nitrogen needed for growth for root tissue. Returns: - The nitrogen deficit for each cohort. + The nitrogen needed for growth for root tissue. """ - return ( - self.n_foliage_deficit - + self.n_wood_deficit - + self.n_roots_deficit - + self.n_reproductive_tissue_deficit - ) + allocation.delta_foliage_mass + * (1 / self.ideal_ratio) + * self.community.stem_traits.zeta + ).squeeze() - def n_for_growth(self, allocation: StemAllocation) -> NDArray[np.float64]: - """Calculate the nitrogen required for growth for each cohort. + def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: + """Calculate the element lost to turnover for root tissue. Returns: - The nitrogen available for growth for each cohort. + The element lost to turnover for root tissue. """ + return (allocation.fine_root_turnover * (1 / self.Cx_ratio)).squeeze() - n_needed_for_foliage_growth = allocation.delta_foliage_mass * ( - 1 / self.plant_consts["foliage_c_n_ratio"] - ) - n_needed_for_stem_growth = allocation.delta_stem_mass * ( - 1 / self.plant_consts["deadwood_c_n_ratio"] - ) +@dataclass +class StemStochiometry(CohortMethods, PandasExporter): + """A class holding the ratios of Carbon to Nitrogen and Phosphorous for stems. - n_needed_for_rt_growth = ( - allocation.delta_foliage_mass - * (1 / self.plant_consts["plant_reproductive_tissue_turnover_c_n_ratio"]) - * self.community.stem_allometry.p_foliage_for_reproductive_tissue - ) + This class holds the current ratios across tissue type for a community object, which + in essence is a series of cohorts. It acts in parallel with StemAllometry, a class + attribute of Community. + """ - n_needed_for_growth = ( - n_needed_for_foliage_growth - + n_needed_for_stem_growth - + n_needed_for_rt_growth - ) + tissues: list[Tissue] + """Tissues for the associated stems.""" + n_cohorts: np.int64 + """The number of cohorts represented by the Stochiometry.""" + community: Community + """The community object that the stochiometry is associated with.""" + element_surplus: NDArray[np.float64] = field(init=False) + """The surplus of the element per cohort.""" - return n_needed_for_growth + def __post_init__(self) -> None: + """Initialize the element surplus for each cohort.""" + self.element_surplus = np.zeros(self.n_cohorts, dtype=np.float64) @property - def cn_ratio_foliage(self) -> NDArray[np.float64]: - """Get the carbon to nitrogen ratio for foliage.""" - return self.community.stem_allometry.foliage_mass / self.n_foliage + def total_element_mass(self) -> NDArray[np.float64]: + """Calculate the total element mass for each cohort. - @property - def cn_ratio_roots(self) -> NDArray[np.float64]: - """Get the carbon to nitrogen ratio for roots.""" - return ( - self.community.stem_traits.zeta * self.community.stem_allometry.foliage_mass - ) / self.n_roots + Returns: + The total nitrogen mass for each cohort. + """ + mass = np.zeros(self.n_cohorts) + for tissue in self.tissues: + mass += tissue.actual_element_mass + return mass @property - def cn_ratio_reproductive_tissue(self) -> NDArray[np.float64]: - """Get the carbon to nitrogen ratio for reproductive tissue.""" - return ( - self.community.stem_allometry.reproductive_tissue_mass - / self.n_reproductive_tissue - ) + def tissue_deficit(self) -> NDArray[np.float64]: + """Calculate the element deficit for a tissue type. - def n_for_foliage_growth(self, delta_foliage_mass) -> NDArray[np.float64]: - """Get the nitrogen needed for foliage growth.""" - return ( - (1 / self.plant_constants.foliage_c_n_ratio) * delta_foliage_mass - ).squeeze() + Returns: + The element deficit for the specified tissue. + """ + element_deficit = np.zeros(self.n_cohorts) + for tissue in self.tissues: + element_deficit += tissue.deficit + return element_deficit + + def account_for_growth(self, allocation: StemAllocation) -> None: + """Distribute the element needed for growth to each tissue type. + + This method updates the actual element mass for each tissue type based on the + element needed for growth calculated from the allocation. + + Args: + allocation: The allocation object containing the growth allocation data. + """ + for tissue in self.tissues: + tissue.actual_element_mass += tissue.element_needed_for_growth(allocation) + self.element_surplus -= tissue.element_needed_for_growth(allocation) + + def account_for_element_loss_turnover(self, allocation: StemAllocation) -> None: + """Calculate the total element lost to turnover for each cohort. + + Elements are lost from the tree in the form of turnover, and so an equivalent + amount of that element is required to replace what was lost. To represent this + process, the element is allocated from the surplus store in the same quantity + as turnover. This uses current ratios so that the C:x ratios are maintained. + + Returns: + The total element lost to turnover for each cohort. + """ + for tissue in self.tissues: + self.element_surplus -= tissue.element_turnover(allocation) def distribute_deficit(self, cohort: int) -> None: """Distribute the nitrogen deficit across the tissue types. @@ -272,20 +383,16 @@ def distribute_deficit(self, cohort: int) -> None: Args: cohort: The cohort to reconcile deficit. """ - deficit = self.n_surplus[cohort] * -1 + deficit = self.element_surplus[cohort] * -1 - self.n_foliage[cohort] = self.n_foliage[cohort] - ( - deficit * self.n_foliage[cohort] / self.total_n[cohort] - ) + for tissue in self.tissues: + tissue.actual_element_mass[cohort] = tissue.actual_element_mass[cohort] - ( + deficit + * tissue.actual_element_mass[cohort] + / self.total_element_mass[cohort] + ) - self.n_wood[cohort] = self.n_wood[cohort] - ( - deficit * self.n_wood[cohort] / self.total_n[cohort] - ) - self.n_roots[cohort] -= deficit * self.n_roots[cohort] / self.total_n[cohort] - self.n_reproductive_tissue[cohort] -= ( - deficit * self.n_reproductive_tissue[cohort] / self.total_n[cohort] - ) - self.n_surplus[cohort] += deficit + self.element_surplus[cohort] += deficit def distribute_surplus(self, cohort: int) -> None: """Distribute the nitrogen surplus across the tissue types for a single cohort. @@ -293,55 +400,26 @@ def distribute_surplus(self, cohort: int) -> None: Args: cohort: The cohort to reconcile surplus. """ - if self.n_surplus[cohort] > self.n_tissue_deficit[cohort]: + if self.element_surplus[cohort] > self.tissue_deficit[cohort]: # If there is sufficient surplus N to cover the existing deficit, the # amount of the deficit is subtracted from the surplus which persists until # the next update. All tissue types are updated to the ideal N ratios. - self.n_surplus[cohort] = ( - self.n_surplus[cohort] - self.n_tissue_deficit[cohort] - ) - self.n_foliage[cohort] = ( - self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio - * self.community.stem_allometry.foliage_mass.squeeze()[cohort] - ) - self.n_wood[cohort] = ( - # self.plant_constants.deadwood_c_n_ratio * - self.community.stem_allometry.stem_mass.squeeze()[cohort] - ) - self.n_roots[cohort] = ( - self.plant_constants.root_turnover_c_n_ratio - * self.community.stem_traits.zeta[cohort] - * self.community.stem_allometry.foliage_mass.squeeze()[cohort] - ) - self.n_reproductive_tissue[cohort] = ( - self.plant_constants.plant_reproductive_tissue_turnover_c_n_ratio - * self.community.stem_allometry.reproductive_tissue_mass.squeeze()[ - cohort - ] + self.element_surplus[cohort] = ( + self.element_surplus[cohort] - self.tissue_deficit[cohort] ) + for tissue in self.tissues: + tissue.actual_element_mass[cohort] = ( + tissue.ideal_ratio[cohort] * tissue.carbon_mass[cohort] + ) else: # If there is not enough surplus N to cover the deficit, the surplus is # distributed across the tissue types in proportion to the N deficit. # The surplus is then set to zero. - self.n_foliage[cohort] += ( - self.n_surplus[cohort] - * self.n_foliage_deficit[cohort] - / self.n_tissue_deficit[cohort] - ) - self.n_wood[cohort] += ( - self.n_surplus[cohort] - * self.n_wood_deficit[cohort] - / self.n_tissue_deficit[cohort] - ) - self.n_roots[cohort] += ( - self.n_surplus[cohort] - * self.n_roots_deficit[cohort] - / self.n_tissue_deficit[cohort] - ) - self.n_reproductive_tissue[cohort] += ( - self.n_surplus[cohort] - * self.n_reproductive_tissue_deficit[cohort] - / self.n_tissue_deficit[cohort] - ) - self.n_surplus[cohort] = 0.0 + for tissue in self.tissues: + tissue.actual_element_mass[cohort] += ( + self.element_surplus[cohort] + * tissue.deficit[cohort] + / self.tissue_deficit[cohort] + ) + self.element_surplus[cohort] = 0.0 From 432d7569041a6c74718b878e24520327d52fc43f Mon Sep 17 00:00:00 2001 From: Sally Date: Tue, 3 Jun 2025 15:27:33 +0100 Subject: [PATCH 07/22] Corrected spelling mistakes. --- virtual_ecosystem/models/plants/stochiometry.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index 46698c322..4b58444b9 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -1,6 +1,6 @@ """The :mod:`~virtual_ecosystem.models.plants.stochiometry` module contains the class for managing plant cohort stochiometry ratios. The carbon mass is stored in plant -alloemetry or allocation, so this class uses thoses as the anchor weights and stores +alloemetry or allocation, so this class uses those as the anchor weights and stores CN and CP ratios. The class holds current CN and CP ratios for foliage, reproductive tissue, wood, and @@ -24,7 +24,7 @@ class Tissue(PandasExporter, CohortMethods): """A dataclass to hold tissue stochiometry data for a set of plant cohorts. - This class holds the current quantitiy of a given element (generally N or P) for a + This class holds the current quantity of a given element (generally N or P) for a specific plant tissue type (generally foliage, wood, roots or reproductive tissue). The class also holds the ideal ratio of the element for that tissue type. They hold an entry for each cohort in the data class. From 90ffd56ac5d91255fd4acb3dd59f29ee888d18fc Mon Sep 17 00:00:00 2001 From: Sally Date: Tue, 3 Jun 2025 16:15:13 +0100 Subject: [PATCH 08/22] Fixing tests and updating pyrealm references to rc3 convention. --- tests/models/plants/test_plants_model.py | 13 ++++++---- tests/models/plants/test_stochiometry.py | 24 +++++++++---------- .../models/plants/plants_model.py | 12 +++++----- .../models/plants/stochiometry.py | 2 +- 4 files changed, 27 insertions(+), 24 deletions(-) diff --git a/tests/models/plants/test_plants_model.py b/tests/models/plants/test_plants_model.py index e82402f75..fe20ed377 100644 --- a/tests/models/plants/test_plants_model.py +++ b/tests/models/plants/test_plants_model.py @@ -394,10 +394,15 @@ def test_PlantsModel_calculate_mycorrhizal_uptakes(fxt_plants_model): fxt_plants_model.calculate_mycorrhizal_uptakes() # Check that all expected variables are generated and have the correct value - assert np.allclose(fxt_plants_model.data["plant_n_uptake_arbuscular"], 0.00216) - assert np.allclose(fxt_plants_model.data["plant_n_uptake_ecto"], 0.000805) - assert np.allclose(fxt_plants_model.data["plant_p_uptake_arbuscular"], 0.000117) - assert np.allclose(fxt_plants_model.data["plant_p_uptake_ecto"], 6.6e-5) + assert np.allclose(fxt_plants_model.data["plant_n_uptake_arbuscular"], 0.5) + assert np.allclose(fxt_plants_model.data["plant_n_uptake_ecto"], 0.5) + assert np.allclose(fxt_plants_model.data["plant_p_uptake_arbuscular"], 0.5) + assert np.allclose(fxt_plants_model.data["plant_p_uptake_ecto"], 0.5) + + # assert np.allclose(fxt_plants_model.data["plant_n_uptake_arbuscular"], 0.00216) + # assert np.allclose(fxt_plants_model.data["plant_n_uptake_ecto"], 0.000805) + # assert np.allclose(fxt_plants_model.data["plant_p_uptake_arbuscular"], 0.000117) + # assert np.allclose(fxt_plants_model.data["plant_p_uptake_ecto"], 6.6e-5) def test_PlantsModel_apply_mortality(fxt_plants_model): diff --git a/tests/models/plants/test_stochiometry.py b/tests/models/plants/test_stochiometry.py index 812d5c903..a94f352cf 100644 --- a/tests/models/plants/test_stochiometry.py +++ b/tests/models/plants/test_stochiometry.py @@ -24,7 +24,8 @@ def test_Stochiometry__init__(fxt_plants_model): element_name="Nitrogen", community=community, ideal_ratio=np.full( - community.number_of_cohorts, plant_consts.foliage_c_n_ratio + community.n_cohorts, + plant_consts.foliage_c_n_ratio, ), actual_element_mass=community.stem_allometry.foliage_mass * plant_consts.foliage_c_n_ratio, @@ -34,7 +35,7 @@ def test_Stochiometry__init__(fxt_plants_model): element_name="Nitrogen", community=community, ideal_ratio=np.full( - community.number_of_cohorts, + community.n_cohorts, plant_consts.root_turnover_c_n_ratio, ), actual_element_mass=plant_consts.root_turnover_c_n_ratio @@ -45,7 +46,8 @@ def test_Stochiometry__init__(fxt_plants_model): element_name="Nitrogen", community=community, ideal_ratio=np.full( - community.number_of_cohorts, plant_consts.deadwood_c_n_ratio + community.n_cohorts, + plant_consts.deadwood_c_n_ratio, ), actual_element_mass=plant_consts.deadwood_c_n_ratio * community.stem_allometry.stem_mass, @@ -54,14 +56,14 @@ def test_Stochiometry__init__(fxt_plants_model): element_name="Nitrogen", community=community, ideal_ratio=np.full( - community.number_of_cohorts, + community.n_cohorts, plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, ), actual_element_mass=community.stem_allometry.reproductive_tissue_mass * plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, ), ], - n_cohorts=community.number_of_cohorts, + n_cohorts=community.n_cohorts, community=community, ) assert n_stochiometry is not None @@ -78,9 +80,7 @@ def test_Stochiometry_FoliageTissue(fxt_plants_model): foliage_tissue = FoliageTissue( element_name="Nitrogen", community=community, - ideal_ratio=np.full( - community.number_of_cohorts, plant_consts.foliage_c_n_ratio - ), + ideal_ratio=np.full(community.n_cohorts, plant_consts.foliage_c_n_ratio), actual_element_mass=community.stem_allometry.foliage_mass * plant_consts.foliage_c_n_ratio, reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, @@ -100,7 +100,7 @@ def test_Stochiometry_ReproductiveTissue(fxt_plants_model): element_name="Nitrogen", community=community, ideal_ratio=np.full( - community.number_of_cohorts, + community.n_cohorts, plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, ), actual_element_mass=community.stem_allometry.reproductive_tissue_mass @@ -121,7 +121,7 @@ def test_Stochiometry_RootTissue(fxt_plants_model): element_name="Nitrogen", community=community, ideal_ratio=np.full( - community.number_of_cohorts, plant_consts.root_turnover_c_n_ratio + community.n_cohorts, plant_consts.root_turnover_c_n_ratio ), actual_element_mass=plant_consts.root_turnover_c_n_ratio * community.stem_traits.zeta @@ -142,9 +142,7 @@ def test_Stochiometry_WoodTissue(fxt_plants_model): wood_tissue = WoodTissue( element_name="Nitrogen", community=community, - ideal_ratio=np.full( - community.number_of_cohorts, plant_consts.deadwood_c_n_ratio - ), + ideal_ratio=np.full(community.n_cohorts, plant_consts.deadwood_c_n_ratio), actual_element_mass=plant_consts.deadwood_c_n_ratio * community.stem_allometry.stem_mass, ) diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index 966f50ede..7c0192f27 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -334,7 +334,7 @@ def _setup( element_name="Nitrogen", community=self.communities[cell_id], ideal_ratio=np.full( - self.communities[cell_id].number_of_cohorts, + self.communities[cell_id].n_cohorts, model_constants.foliage_c_n_ratio, ), actual_element_mass=self.communities[ @@ -342,7 +342,7 @@ def _setup( ].stem_allometry.foliage_mass * model_constants.foliage_c_n_ratio, reclaim_ratio=np.full( - self.communities[cell_id].number_of_cohorts, + self.communities[cell_id].n_cohorts, model_constants.leaf_turnover_c_n_ratio, ), ), @@ -350,7 +350,7 @@ def _setup( element_name="Nitrogen", community=self.communities[cell_id], ideal_ratio=np.full( - self.communities[cell_id].number_of_cohorts, + self.communities[cell_id].n_cohorts, model_constants.root_turnover_c_n_ratio, ), actual_element_mass=model_constants.root_turnover_c_n_ratio @@ -361,7 +361,7 @@ def _setup( element_name="Nitrogen", community=self.communities[cell_id], ideal_ratio=np.full( - self.communities[cell_id].number_of_cohorts, + self.communities[cell_id].n_cohorts, model_constants.deadwood_c_n_ratio, ), actual_element_mass=model_constants.deadwood_c_n_ratio @@ -371,7 +371,7 @@ def _setup( element_name="Nitrogen", community=self.communities[cell_id], ideal_ratio=np.full( - self.communities[cell_id].number_of_cohorts, + self.communities[cell_id].n_cohorts, model_constants.plant_reproductive_tissue_turnover_c_n_ratio, ), actual_element_mass=self.communities[ @@ -380,7 +380,7 @@ def _setup( * self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio, # noqa: E501 ), ], - n_cohorts=self.communities[cell_id].number_of_cohorts, + n_cohorts=self.communities[cell_id].n_cohorts, community=self.communities[cell_id], ) for cell_id in self.communities.keys() diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index 4b58444b9..23a6e0502 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -249,7 +249,7 @@ def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: Returns: The element lost to turnover for wood tissue. """ - return np.zeros(self.community.number_of_cohorts) + return np.zeros(self.community.n_cohorts) class RootTissue(Tissue): From 757d99d30078c4172d9d9acf5b0e4d43c6725ebc Mon Sep 17 00:00:00 2001 From: Sally Date: Wed, 4 Jun 2025 11:22:21 +0100 Subject: [PATCH 09/22] Sphinx targets issue resolved. --- docs/source/conf.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/docs/source/conf.py b/docs/source/conf.py index f42dc9bda..576a25b6b 100644 --- a/docs/source/conf.py +++ b/docs/source/conf.py @@ -163,6 +163,8 @@ class MyReferenceStyle(AuthorYearReferenceStyle): ("py:class", "Flora"), ("py:class", "PModelConst"), ("py:class", "CoreConst"), + ("py:class", "StemAllocation"), + ("py:class", "StemStochiometry"), ] intersphinx_mapping = { "numpy": ("https://numpy.org/doc/stable/", None), From f38b1166a876d118153fd73d04613c69dfdcb1e5 Mon Sep 17 00:00:00 2001 From: Sally Date: Wed, 4 Jun 2025 12:00:18 +0100 Subject: [PATCH 10/22] Remove temp mycorrizha functions. --- tests/models/plants/test_plants_model.py | 13 ++++--------- virtual_ecosystem/models/plants/plants_model.py | 10 ---------- 2 files changed, 4 insertions(+), 19 deletions(-) diff --git a/tests/models/plants/test_plants_model.py b/tests/models/plants/test_plants_model.py index fe20ed377..e82402f75 100644 --- a/tests/models/plants/test_plants_model.py +++ b/tests/models/plants/test_plants_model.py @@ -394,15 +394,10 @@ def test_PlantsModel_calculate_mycorrhizal_uptakes(fxt_plants_model): fxt_plants_model.calculate_mycorrhizal_uptakes() # Check that all expected variables are generated and have the correct value - assert np.allclose(fxt_plants_model.data["plant_n_uptake_arbuscular"], 0.5) - assert np.allclose(fxt_plants_model.data["plant_n_uptake_ecto"], 0.5) - assert np.allclose(fxt_plants_model.data["plant_p_uptake_arbuscular"], 0.5) - assert np.allclose(fxt_plants_model.data["plant_p_uptake_ecto"], 0.5) - - # assert np.allclose(fxt_plants_model.data["plant_n_uptake_arbuscular"], 0.00216) - # assert np.allclose(fxt_plants_model.data["plant_n_uptake_ecto"], 0.000805) - # assert np.allclose(fxt_plants_model.data["plant_p_uptake_arbuscular"], 0.000117) - # assert np.allclose(fxt_plants_model.data["plant_p_uptake_ecto"], 6.6e-5) + assert np.allclose(fxt_plants_model.data["plant_n_uptake_arbuscular"], 0.00216) + assert np.allclose(fxt_plants_model.data["plant_n_uptake_ecto"], 0.000805) + assert np.allclose(fxt_plants_model.data["plant_p_uptake_arbuscular"], 0.000117) + assert np.allclose(fxt_plants_model.data["plant_p_uptake_ecto"], 6.6e-5) def test_PlantsModel_apply_mortality(fxt_plants_model): diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index 7c0192f27..c9530b9a0 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -298,8 +298,6 @@ def _setup( self.flora = flora self.model_constants = model_constants - self.temp_populate_symbiote_vars() - # Adjust flora turnover rates to timestep # TODO: Pyrealm provides annual turnover rates. Dividing by the number of # updates_per_year to get monthly turnover values is naive and will @@ -1204,14 +1202,6 @@ def partition_reproductive_tissue( return n_propagules, non_propagule_mass - def temp_populate_symbiote_vars(self) -> None: - """Jacob has these vars in an open PR. Will delete before merging.""" - - self.data["ecto_supply_limit_n"] = xr.full_like(self.data["elevation"], 1) - self.data["arbuscular_supply_limit_n"] = xr.full_like(self.data["elevation"], 1) - self.data["ecto_supply_limit_p"] = xr.full_like(self.data["elevation"], 1) - self.data["arbuscular_supply_limit_p"] = xr.full_like(self.data["elevation"], 1) - def convert_to_litter_units( self, input_mass: NDArray[np.float64] ) -> NDArray[np.float64]: From 3875ee89302515bcabc03b4c193c249c37bf9b44 Mon Sep 17 00:00:00 2001 From: Sally Date: Wed, 4 Jun 2025 13:09:20 +0100 Subject: [PATCH 11/22] Update timeout. --- docs/source/using_the_ve/virtual_ecosystem_in_use.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/source/using_the_ve/virtual_ecosystem_in_use.md b/docs/source/using_the_ve/virtual_ecosystem_in_use.md index 0d06a8c89..4c7818877 100644 --- a/docs/source/using_the_ve/virtual_ecosystem_in_use.md +++ b/docs/source/using_the_ve/virtual_ecosystem_in_use.md @@ -1,6 +1,6 @@ --- execution: - timeout: 60 + timeout: 90 jupytext: formats: md:myst text_representation: From 5d73cb06f1857c74fa7143017853dd9a9fc1a81e Mon Sep 17 00:00:00 2001 From: Sally Date: Thu, 5 Jun 2025 16:37:52 +0100 Subject: [PATCH 12/22] Add stochiometry handling for Phosphorus, clean up some code. --- tests/models/plants/test_plants_model.py | 6 + tests/models/plants/test_stochiometry.py | 10 +- .../models/plants/plants_model.py | 262 +++++++++++------- .../models/plants/stochiometry.py | 20 +- 4 files changed, 181 insertions(+), 117 deletions(-) diff --git a/tests/models/plants/test_plants_model.py b/tests/models/plants/test_plants_model.py index e82402f75..c92881115 100644 --- a/tests/models/plants/test_plants_model.py +++ b/tests/models/plants/test_plants_model.py @@ -378,6 +378,10 @@ def test_PlantsModel_calculate_turnover_constant_override( def test_PlantsModel_calculate_nutrient_uptake(fxt_plants_model): """Test the calculate_nutrient_uptake method of the plants model.""" + # Provide transpiration values + fxt_plants_model.per_stem_transpiration = { + cell_id: np.array([10]) for cell_id in fxt_plants_model.communities.keys() + } # Check reset fxt_plants_model.calculate_nutrient_uptake() @@ -386,6 +390,8 @@ def test_PlantsModel_calculate_nutrient_uptake(fxt_plants_model): assert np.allclose(fxt_plants_model.data["plant_nitrate_uptake"], 7.5e-3) assert np.allclose(fxt_plants_model.data["plant_phosphorus_uptake"], 3.0e-5) + # TODO: add test for element uptake + def test_PlantsModel_calculate_mycorrhizal_uptakes(fxt_plants_model): """Test the calculate_mycorrhizal_uptakes method of the plants model.""" diff --git a/tests/models/plants/test_stochiometry.py b/tests/models/plants/test_stochiometry.py index a94f352cf..0cd1c033b 100644 --- a/tests/models/plants/test_stochiometry.py +++ b/tests/models/plants/test_stochiometry.py @@ -19,9 +19,9 @@ def test_Stochiometry__init__(fxt_plants_model): community = fxt_plants_model.communities[cell_id] n_stochiometry = StemStochiometry( + element="N", tissues=[ FoliageTissue( - element_name="Nitrogen", community=community, ideal_ratio=np.full( community.n_cohorts, @@ -32,7 +32,6 @@ def test_Stochiometry__init__(fxt_plants_model): reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, ), RootTissue( - element_name="Nitrogen", community=community, ideal_ratio=np.full( community.n_cohorts, @@ -43,7 +42,6 @@ def test_Stochiometry__init__(fxt_plants_model): * community.stem_allometry.foliage_mass, ), WoodTissue( - element_name="Nitrogen", community=community, ideal_ratio=np.full( community.n_cohorts, @@ -53,7 +51,6 @@ def test_Stochiometry__init__(fxt_plants_model): * community.stem_allometry.stem_mass, ), ReproductiveTissue( - element_name="Nitrogen", community=community, ideal_ratio=np.full( community.n_cohorts, @@ -63,7 +60,6 @@ def test_Stochiometry__init__(fxt_plants_model): * plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, ), ], - n_cohorts=community.n_cohorts, community=community, ) assert n_stochiometry is not None @@ -78,7 +74,6 @@ def test_Stochiometry_FoliageTissue(fxt_plants_model): community = fxt_plants_model.communities[cell_id] foliage_tissue = FoliageTissue( - element_name="Nitrogen", community=community, ideal_ratio=np.full(community.n_cohorts, plant_consts.foliage_c_n_ratio), actual_element_mass=community.stem_allometry.foliage_mass @@ -97,7 +92,6 @@ def test_Stochiometry_ReproductiveTissue(fxt_plants_model): community = fxt_plants_model.communities[cell_id] reproductive_tissue = ReproductiveTissue( - element_name="Nitrogen", community=community, ideal_ratio=np.full( community.n_cohorts, @@ -118,7 +112,6 @@ def test_Stochiometry_RootTissue(fxt_plants_model): community = fxt_plants_model.communities[cell_id] root_tissue = RootTissue( - element_name="Nitrogen", community=community, ideal_ratio=np.full( community.n_cohorts, plant_consts.root_turnover_c_n_ratio @@ -140,7 +133,6 @@ def test_Stochiometry_WoodTissue(fxt_plants_model): community = fxt_plants_model.communities[cell_id] wood_tissue = WoodTissue( - element_name="Nitrogen", community=community, ideal_ratio=np.full(community.n_cohorts, plant_consts.deadwood_c_n_ratio), actual_element_mass=plant_consts.deadwood_c_n_ratio diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index c9530b9a0..7a961bff6 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -212,7 +212,7 @@ def __init__( self.communities: PlantCommunities """An instance of PlantCommunities providing dictionary access keyed by cell id to PlantCommunity instances for each cell.""" - self.stochiometries: dict[int, StemStochiometry] + self.stochiometries: dict[int, dict[str, StemStochiometry]] """A dictionary keyed by cell id giving the stochiometry of each community.""" self.allocations: dict[int, StemAllocation] """A dictionary keyed by cell id giving the allocation of each community.""" @@ -326,61 +326,110 @@ def _setup( # Initialize the stochiometries of each cohort self.stochiometries = { - cell_id: StemStochiometry( - tissues=[ - FoliageTissue( - element_name="Nitrogen", - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.foliage_c_n_ratio, + cell_id: { + "N": StemStochiometry( + element="N", + tissues=[ + FoliageTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.foliage_c_n_ratio, + ), + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.foliage_mass + * model_constants.foliage_c_n_ratio, + reclaim_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.leaf_turnover_c_n_ratio, + ), ), - actual_element_mass=self.communities[ - cell_id - ].stem_allometry.foliage_mass - * model_constants.foliage_c_n_ratio, - reclaim_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.leaf_turnover_c_n_ratio, + RootTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.root_turnover_c_n_ratio, + ), + actual_element_mass=model_constants.root_turnover_c_n_ratio + * self.communities[cell_id].stem_traits.zeta + * self.communities[cell_id].stem_allometry.foliage_mass, ), - ), - RootTissue( - element_name="Nitrogen", - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.root_turnover_c_n_ratio, + WoodTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.deadwood_c_n_ratio, + ), + actual_element_mass=model_constants.deadwood_c_n_ratio + * self.communities[cell_id].stem_allometry.stem_mass, ), - actual_element_mass=model_constants.root_turnover_c_n_ratio - * self.communities[cell_id].stem_traits.zeta - * self.communities[cell_id].stem_allometry.foliage_mass, - ), - WoodTissue( - element_name="Nitrogen", - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.deadwood_c_n_ratio, + ReproductiveTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.plant_reproductive_tissue_turnover_c_n_ratio, + ), + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.reproductive_tissue_mass + * self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio, # noqa: E501 ), - actual_element_mass=model_constants.deadwood_c_n_ratio - * self.communities[cell_id].stem_allometry.stem_mass, - ), - ReproductiveTissue( - element_name="Nitrogen", - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.plant_reproductive_tissue_turnover_c_n_ratio, + ], + community=self.communities[cell_id], + ), + "P": StemStochiometry( + element="P", + tissues=[ + FoliageTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.foliage_c_p_ratio, + ), + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.foliage_mass + * model_constants.foliage_c_p_ratio, + reclaim_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.leaf_turnover_c_p_ratio, + ), ), - actual_element_mass=self.communities[ - cell_id - ].stem_allometry.reproductive_tissue_mass - * self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio, # noqa: E501 - ), - ], - n_cohorts=self.communities[cell_id].n_cohorts, - community=self.communities[cell_id], - ) + RootTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.root_turnover_c_p_ratio, + ), + actual_element_mass=model_constants.root_turnover_c_p_ratio + * self.communities[cell_id].stem_traits.zeta + * self.communities[cell_id].stem_allometry.foliage_mass, + ), + WoodTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.deadwood_c_p_ratio, + ), + actual_element_mass=model_constants.deadwood_c_p_ratio + * self.communities[cell_id].stem_allometry.stem_mass, + ), + ReproductiveTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.plant_reproductive_tissue_turnover_c_p_ratio, + ), + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.reproductive_tissue_mass + * self.model_constants.plant_reproductive_tissue_turnover_c_p_ratio, # noqa: E501 + ), + ], + community=self.communities[cell_id], + ), + } for cell_id in self.communities.keys() } @@ -450,14 +499,15 @@ def _update(self, time_index: int, **kwargs: Any) -> None: # Estimate the canopy GPP and growth with the updated this update self.calculate_light_use_efficiency() self.estimate_gpp(time_index=time_index) + + # Calculate uptake from each inorganic soil nutrient pool + self.calculate_nutrient_uptake() + self.allocate_gpp() # Calculate the turnover of each plant biomass pool self.calculate_turnover() - # Calculate uptake from each inorganic soil nutrient pool - self.calculate_nutrient_uptake() - # Calculate the rate at which plants take nutrients from mycorrhizal fungi self.calculate_mycorrhizal_uptakes() @@ -704,21 +754,6 @@ def estimate_gpp(self, time_index: int) -> None: * self.model_timing.update_interval_seconds ).sum(axis=0) - # Calculate N uptake (g N per stem) due to transpiration. Conversion units: - # - Per stem transpiration (µmol H2O per stem) - # - Conversion factor from µmol H2O to m^3 (1.08015*10^-11) - # - Concentration of dissolved N (kg m^-3) - # - Kg to g (1000) - self.stochiometries[cell_id].element_surplus += ( - self.per_stem_transpiration[cell_id] - * (1.8015 * pow(10.0, -11)) - * ( - self.data["dissolved_ammonium"][cell_id] - + self.data["dissolved_nitrate"][cell_id] - ).item() - * 1000 - ) - # Calculate the total transpiration per layer in mm m2 in mm, converted from # an initial value is in µmol m2 s1 transpiration[active_layers, cell_id] = ( @@ -767,7 +802,7 @@ def allocate_gpp(self) -> None: for cell_id in self.communities.keys(): community = self.communities[cell_id] cohorts = community.cohorts - stochiometry = self.stochiometries[cell_id] + stochiometries = self.stochiometries[cell_id] # Calculate the allocation of GPP per stem stem_allocation = StemAllocation( @@ -847,7 +882,8 @@ def allocate_gpp(self) -> None: ) # SECOND: ALLOCATE N TO REGROW WHAT WAS LOST TO TURNOVER - stochiometry.account_for_element_loss_turnover(stem_allocation) + for stochiometry in stochiometries.values(): + stochiometry.account_for_element_loss_turnover(stem_allocation) # THIRD, ALLOCATE GPP TO ACTIVE NUTRIENT PATHWAYS: # Allocate the topsliced GPP to root exudates with remainder as active @@ -879,9 +915,10 @@ def allocate_gpp(self) -> None: # Subtract the N/P required from growth from the element store, and # redistribute it to the individual tisuses. - stochiometry.account_for_growth(stem_allocation) + for stochiometry in stochiometries.values(): + stochiometry.account_for_growth(stem_allocation) - # Balance the N surplus/deficit with the symbiote carbon supply + # Balance the N & P surplus/deficit with the symbiote carbon supply n_weighted_avg = np.dot( self.data["ecto_supply_limit_n"][cell_id] + self.data["arbuscular_supply_limit_n"][cell_id], @@ -891,35 +928,48 @@ def allocate_gpp(self) -> None: n_available_per_stem = np.divide( n_avaiable_per_cohort, cohorts.n_individuals ) + stochiometries["N"].element_surplus = ( + stochiometries["N"].element_surplus + n_available_per_stem + ) + + p_weighted_avg = np.dot( + self.data["ecto_supply_limit_p"][cell_id] + + self.data["arbuscular_supply_limit_p"][cell_id], + cohorts.n_individuals, + ) + p_available_per_cohort = p_weighted_avg / sum(p_weighted_avg) + p_available_per_stem = np.divide( + p_available_per_cohort, cohorts.n_individuals + ) - stochiometry.element_surplus = ( - stochiometry.element_surplus + n_available_per_stem + stochiometries["P"].element_surplus = ( + stochiometries["P"].element_surplus + p_available_per_stem ) # Cohort by cohort, distribute the surplus/deficit across the tissue types for cohort in range(len(cohorts.n_individuals)): - if stochiometry.element_surplus[cohort] < 0: - # Distribute deficit across the tissue types - stochiometry.distribute_deficit(cohort) - - elif ( - stochiometry.element_surplus[cohort] > 0 - and stochiometry.tissue_deficit[cohort] > 0 - ): - # Distribute the surplus across the tissue types - stochiometry.distribute_surplus(cohort) - - else: - # NO ADJUSTMENT REQUIRED - there is a surplus in the store, but - # there is no deficit in the tissue types. - pass + for stochiometry in stochiometries.values(): + if stochiometry.element_surplus[cohort] < 0: + # Distribute deficit across the tissue types + stochiometry.distribute_deficit(cohort) + + elif ( + stochiometry.element_surplus[cohort] > 0 + and stochiometry.tissue_deficit[cohort] > 0 + ): + # Distribute the surplus across the tissue types + stochiometry.distribute_surplus(cohort) + + else: + # NO ADJUSTMENT REQUIRED - there is a surplus in the store, but + # there is no deficit in the tissue types. + pass # Update community allometry with new dbh values community.stem_allometry = StemAllometry( stem_traits=community.stem_traits, at_dbh=cohorts.dbh_values ) - # self.allocations[cell_id] = stem_allocation self.update_cn_ratios() def apply_mortality(self) -> None: @@ -959,6 +1009,8 @@ def update_cn_ratios(self) -> None: This function updates the C:N and C:P ratios of various plant tissues, including deadwood, leaf turnover, plant reproductive tissue turnover, and root turnover. + # TODO: Update this to use the Stochiometry class values. + Warning: At present, this function just sets values to original constants. """ @@ -1032,7 +1084,9 @@ def calculate_nutrient_uptake(self) -> None: """Calculate uptake of soil nutrients by the plant community. This function calculates the rate a which plants take up inorganic nutrients - (ammonium, nitrate, and labile phosphorus) from the soil. + (ammonium, nitrate, and labile phosphorus) from the soil. The function then + assigns the N/P uptake values to the respective community through the + stochiometry class. Warning: At present, this function just calculates uptake based on an entirely made @@ -1040,11 +1094,33 @@ def calculate_nutrient_uptake(self) -> None: """ # Assume plants can take 0.1% of the available nutrient per day - # TODO: This is now set by transpiration? self.data["plant_ammonium_uptake"] = self.data["dissolved_ammonium"] * 0.01 self.data["plant_nitrate_uptake"] = self.data["dissolved_nitrate"] * 0.01 self.data["plant_phosphorus_uptake"] = self.data["dissolved_phosphorus"] * 0.01 + # Calculate N/P uptake (g N/P per stem) due to transpiration. Multiply: + # - Per stem transpiration (µmol H2O per stem) + # - Conversion factor from µmol H2O to m^3 (1.08015*10^-11) + # - Concentration of N/P uptake (kg m^-3) + # - Kg to g (1000) + + for cell_id in self.communities.keys(): + self.stochiometries[cell_id]["N"].element_surplus += ( + self.per_stem_transpiration[cell_id] + * (1.8015 * pow(10.0, -11)) + * ( + self.data["plant_ammonium_uptake"][cell_id] + + self.data["plant_nitrate_uptake"][cell_id] + ).item() + * 1000 + ) + self.stochiometries[cell_id]["P"].element_surplus += ( + self.per_stem_transpiration[cell_id] + * (1.8015 * pow(10.0, -11)) + * (self.data["plant_phosphorus_uptake"][cell_id]).item() + * 1000 + ) + def calculate_mycorrhizal_uptakes(self) -> None: """Calculate the rate at which plants take nutrients from mycorrhizal fungi. diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index 23a6e0502..390941150 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -30,8 +30,6 @@ class Tissue(PandasExporter, CohortMethods): an entry for each cohort in the data class. """ - element_name: str - """The name of the element type.""" tissue_name: str """The name of the tissue type.""" community: Community @@ -109,14 +107,12 @@ class FoliageTissue(Tissue): def __init__( self, - element_name: str, community: Community, ideal_ratio: NDArray[np.float64], actual_element_mass: NDArray[np.float64], reclaim_ratio: NDArray[np.float64], ): super().__init__( - element_name=element_name, tissue_name="Foliage", community=community, ideal_ratio=ideal_ratio, @@ -161,13 +157,11 @@ class ReproductiveTissue(Tissue): def __init__( self, - element_name: str, community: Community, ideal_ratio: NDArray[np.float64], actual_element_mass: NDArray[np.float64], ): super().__init__( - element_name=element_name, tissue_name="Reproductive", community=community, ideal_ratio=ideal_ratio, @@ -211,13 +205,11 @@ class WoodTissue(Tissue): def __init__( self, - element_name: str, community: Community, ideal_ratio: NDArray[np.float64], actual_element_mass: NDArray[np.float64], ): super().__init__( - element_name=element_name, tissue_name="Wood", community=community, ideal_ratio=ideal_ratio, @@ -257,13 +249,11 @@ class RootTissue(Tissue): def __init__( self, - element_name: str, community: Community, ideal_ratio: NDArray[np.float64], actual_element_mass: NDArray[np.float64], ): super().__init__( - element_name=element_name, tissue_name="Roots", community=community, ideal_ratio=ideal_ratio, @@ -313,10 +303,10 @@ class StemStochiometry(CohortMethods, PandasExporter): attribute of Community. """ + element: str + """The name of the element (e.g., N or P).""" tissues: list[Tissue] """Tissues for the associated stems.""" - n_cohorts: np.int64 - """The number of cohorts represented by the Stochiometry.""" community: Community """The community object that the stochiometry is associated with.""" element_surplus: NDArray[np.float64] = field(init=False) @@ -324,7 +314,7 @@ class StemStochiometry(CohortMethods, PandasExporter): def __post_init__(self) -> None: """Initialize the element surplus for each cohort.""" - self.element_surplus = np.zeros(self.n_cohorts, dtype=np.float64) + self.element_surplus = np.zeros(self.community.n_cohorts, dtype=np.float64) @property def total_element_mass(self) -> NDArray[np.float64]: @@ -333,7 +323,7 @@ def total_element_mass(self) -> NDArray[np.float64]: Returns: The total nitrogen mass for each cohort. """ - mass = np.zeros(self.n_cohorts) + mass = np.zeros(self.community.n_cohorts) for tissue in self.tissues: mass += tissue.actual_element_mass return mass @@ -345,7 +335,7 @@ def tissue_deficit(self) -> NDArray[np.float64]: Returns: The element deficit for the specified tissue. """ - element_deficit = np.zeros(self.n_cohorts) + element_deficit = np.zeros(self.community.n_cohorts) for tissue in self.tissues: element_deficit += tissue.deficit return element_deficit From 532ee2802dfd84d352a5aeca24bb2c00e4a4be75 Mon Sep 17 00:00:00 2001 From: Sally Date: Mon, 16 Jun 2025 13:46:24 +0100 Subject: [PATCH 13/22] Updated foliage helper functions test. --- tests/models/plants/test_stochiometry.py | 183 +++++++++++++---------- 1 file changed, 107 insertions(+), 76 deletions(-) diff --git a/tests/models/plants/test_stochiometry.py b/tests/models/plants/test_stochiometry.py index 0cd1c033b..c88be749a 100644 --- a/tests/models/plants/test_stochiometry.py +++ b/tests/models/plants/test_stochiometry.py @@ -1,6 +1,113 @@ """Tests for the Stochiometry class.""" +from contextlib import nullcontext as does_not_raise + import numpy as np +import pytest + + +@pytest.mark.parametrize( + "expected_exception", + [ + (does_not_raise()), + ], +) +def test_FoliageTissue__init__(fxt_plants_model, expected_exception): + """Test the foliage stochiometry.""" + from virtual_ecosystem.models.plants.stochiometry import FoliageTissue + + plant_consts = fxt_plants_model.model_constants + for cell_id in fxt_plants_model.communities.keys(): + community = fxt_plants_model.communities[cell_id] + + tissue_model = FoliageTissue( + community=community, + ideal_ratio=np.full(community.n_cohorts, plant_consts.foliage_c_n_ratio), + actual_element_mass=community.stem_allometry.foliage_mass + * plant_consts.foliage_c_n_ratio, + reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, + ) + + with expected_exception: + if expected_exception is does_not_raise(): + assert isinstance(tissue_model, FoliageTissue) + + +@pytest.mark.parametrize( + argnames=( + "classname, carbon_mass, deficit, element_needed_for_growth, " + "element_turnover, Cx_ratio" + ), + argvalues=[ + pytest.param( + "FoliageTissue", + [1.00834120e01, 3.27620602e-01, 3.85363977e-03], + [0, 0, 0], + [-0.35486191, 0.02330898, 0.87583064], + [-9.83256242e02, -3.19470237e01, -3.75777103e-01], + [0.06666667, 0.06666667, 0.06666667], + ), + ], +) +def test_Tissue_class_functions( + fxt_plants_model, + classname, + carbon_mass, + deficit, + element_needed_for_growth, + element_turnover, + Cx_ratio, +): + """Test the carbon mass calculation in FoliageTissue.""" + from pyrealm.demography.tmodel import StemAllocation + + from virtual_ecosystem.models.plants.stochiometry import FoliageTissue + + plant_consts = fxt_plants_model.model_constants + fxt_plants_model.per_stem_gpp = { + cell_id: np.array([55]) for cell_id in fxt_plants_model.communities.keys() + } + cell_id = 0 + community = fxt_plants_model.communities[cell_id] + tissue_model = FoliageTissue( + community=community, + ideal_ratio=np.full(community.n_cohorts, plant_consts.foliage_c_n_ratio), + actual_element_mass=community.stem_allometry.foliage_mass + * plant_consts.foliage_c_n_ratio, + reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, + ) + + stem_allocation = StemAllocation( + stem_traits=community.stem_traits, + stem_allometry=community.stem_allometry, + whole_crown_gpp=fxt_plants_model.per_stem_gpp[cell_id], + ) + + # Assert that the calculated carbon mass matches the expected value + assert np.allclose( + tissue_model.carbon_mass, + carbon_mass, + ) + + assert np.allclose( + tissue_model.deficit, + deficit, + ) + + assert np.allclose( + tissue_model.element_needed_for_growth(stem_allocation), + element_needed_for_growth, + ) + + assert np.allclose( + tissue_model.element_turnover(stem_allocation), + element_turnover, + ) + + assert np.allclose( + tissue_model.Cx_ratio, + Cx_ratio, + ) def test_Stochiometry__init__(fxt_plants_model): @@ -63,79 +170,3 @@ def test_Stochiometry__init__(fxt_plants_model): community=community, ) assert n_stochiometry is not None - - -def test_Stochiometry_FoliageTissue(fxt_plants_model): - """Test the foliage stochiometry.""" - from virtual_ecosystem.models.plants.stochiometry import FoliageTissue - - plant_consts = fxt_plants_model.model_constants - for cell_id in fxt_plants_model.communities.keys(): - community = fxt_plants_model.communities[cell_id] - - foliage_tissue = FoliageTissue( - community=community, - ideal_ratio=np.full(community.n_cohorts, plant_consts.foliage_c_n_ratio), - actual_element_mass=community.stem_allometry.foliage_mass - * plant_consts.foliage_c_n_ratio, - reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, - ) - assert foliage_tissue is not None - - -def test_Stochiometry_ReproductiveTissue(fxt_plants_model): - """Test the reproductive tissue stochiometry.""" - from virtual_ecosystem.models.plants.stochiometry import ReproductiveTissue - - plant_consts = fxt_plants_model.model_constants - for cell_id in fxt_plants_model.communities.keys(): - community = fxt_plants_model.communities[cell_id] - - reproductive_tissue = ReproductiveTissue( - community=community, - ideal_ratio=np.full( - community.n_cohorts, - plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, - ), - actual_element_mass=community.stem_allometry.reproductive_tissue_mass - * plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, - ) - assert reproductive_tissue is not None - - -def test_Stochiometry_RootTissue(fxt_plants_model): - """Test the reproductive tissue stochiometry.""" - from virtual_ecosystem.models.plants.stochiometry import RootTissue - - plant_consts = fxt_plants_model.model_constants - for cell_id in fxt_plants_model.communities.keys(): - community = fxt_plants_model.communities[cell_id] - - root_tissue = RootTissue( - community=community, - ideal_ratio=np.full( - community.n_cohorts, plant_consts.root_turnover_c_n_ratio - ), - actual_element_mass=plant_consts.root_turnover_c_n_ratio - * community.stem_traits.zeta - * community.stem_allometry.foliage_mass, - ) - - assert root_tissue is not None - - -def test_Stochiometry_WoodTissue(fxt_plants_model): - """Test the reproductive tissue stochiometry.""" - from virtual_ecosystem.models.plants.stochiometry import WoodTissue - - plant_consts = fxt_plants_model.model_constants - for cell_id in fxt_plants_model.communities.keys(): - community = fxt_plants_model.communities[cell_id] - - wood_tissue = WoodTissue( - community=community, - ideal_ratio=np.full(community.n_cohorts, plant_consts.deadwood_c_n_ratio), - actual_element_mass=plant_consts.deadwood_c_n_ratio - * community.stem_allometry.stem_mass, - ) - assert wood_tissue is not None From e84f911d46c8d2997117c14ec312d55ed0f79c9c Mon Sep 17 00:00:00 2001 From: Sally Date: Mon, 16 Jun 2025 17:24:32 +0100 Subject: [PATCH 14/22] Make Stochiometry test fixture and member tests. --- tests/models/plants/test_stochiometry.py | 209 +++++++++++++++++------ 1 file changed, 153 insertions(+), 56 deletions(-) diff --git a/tests/models/plants/test_stochiometry.py b/tests/models/plants/test_stochiometry.py index c88be749a..b36e661df 100644 --- a/tests/models/plants/test_stochiometry.py +++ b/tests/models/plants/test_stochiometry.py @@ -1,18 +1,10 @@ """Tests for the Stochiometry class.""" -from contextlib import nullcontext as does_not_raise - import numpy as np import pytest -@pytest.mark.parametrize( - "expected_exception", - [ - (does_not_raise()), - ], -) -def test_FoliageTissue__init__(fxt_plants_model, expected_exception): +def test_FoliageTissue__init__(fxt_plants_model): """Test the foliage stochiometry.""" from virtual_ecosystem.models.plants.stochiometry import FoliageTissue @@ -28,9 +20,7 @@ def test_FoliageTissue__init__(fxt_plants_model, expected_exception): reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, ) - with expected_exception: - if expected_exception is does_not_raise(): - assert isinstance(tissue_model, FoliageTissue) + assert isinstance(tissue_model, FoliageTissue) @pytest.mark.parametrize( @@ -110,7 +100,8 @@ def test_Tissue_class_functions( ) -def test_Stochiometry__init__(fxt_plants_model): +@pytest.fixture +def fxt_stochiometry_model(fxt_plants_model): """Fixture for the Stochiometry class.""" from virtual_ecosystem.models.plants.stochiometry import ( @@ -122,51 +113,157 @@ def test_Stochiometry__init__(fxt_plants_model): ) plant_consts = fxt_plants_model.model_constants - for cell_id in fxt_plants_model.communities.keys(): - community = fxt_plants_model.communities[cell_id] + cell_id = 0 # Assuming we are testing for the first cell + community = fxt_plants_model.communities[cell_id] - n_stochiometry = StemStochiometry( - element="N", - tissues=[ - FoliageTissue( - community=community, - ideal_ratio=np.full( - community.n_cohorts, - plant_consts.foliage_c_n_ratio, - ), - actual_element_mass=community.stem_allometry.foliage_mass - * plant_consts.foliage_c_n_ratio, - reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, + n_stochiometry = StemStochiometry( + element="N", + tissues=[ + FoliageTissue( + community=community, + ideal_ratio=np.full( + community.n_cohorts, + plant_consts.foliage_c_n_ratio, ), - RootTissue( - community=community, - ideal_ratio=np.full( - community.n_cohorts, - plant_consts.root_turnover_c_n_ratio, - ), - actual_element_mass=plant_consts.root_turnover_c_n_ratio - * community.stem_traits.zeta - * community.stem_allometry.foliage_mass, + actual_element_mass=community.stem_allometry.foliage_mass + * plant_consts.foliage_c_n_ratio, + reclaim_ratio=plant_consts.leaf_turnover_c_n_ratio, + ), + RootTissue( + community=community, + ideal_ratio=np.full( + community.n_cohorts, + plant_consts.root_turnover_c_n_ratio, ), - WoodTissue( - community=community, - ideal_ratio=np.full( - community.n_cohorts, - plant_consts.deadwood_c_n_ratio, - ), - actual_element_mass=plant_consts.deadwood_c_n_ratio - * community.stem_allometry.stem_mass, + actual_element_mass=plant_consts.root_turnover_c_n_ratio + * community.stem_traits.zeta + * community.stem_allometry.foliage_mass, + ), + WoodTissue( + community=community, + ideal_ratio=np.full( + community.n_cohorts, + plant_consts.deadwood_c_n_ratio, ), - ReproductiveTissue( - community=community, - ideal_ratio=np.full( - community.n_cohorts, - plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, - ), - actual_element_mass=community.stem_allometry.reproductive_tissue_mass - * plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, + actual_element_mass=plant_consts.deadwood_c_n_ratio + * community.stem_allometry.stem_mass, + ), + ReproductiveTissue( + community=community, + ideal_ratio=np.full( + community.n_cohorts, + plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, ), - ], - community=community, - ) - assert n_stochiometry is not None + actual_element_mass=community.stem_allometry.reproductive_tissue_mass + * plant_consts.plant_reproductive_tissue_turnover_c_n_ratio, + ), + ], + community=community, + ) + return n_stochiometry + + +def test_Stochiometry__init__(fxt_stochiometry_model): + """Test the Stochiometry class initialization.""" + + from virtual_ecosystem.models.plants.stochiometry import StemStochiometry + + assert isinstance(fxt_stochiometry_model, StemStochiometry) + + +def test_Stochiometry_total_element_mass(fxt_stochiometry_model): + """Test the total_element_mass method of the Stochiometry class.""" + + # Calculate the total element mass + total_mass = fxt_stochiometry_model.total_element_mass + + assert np.allclose(total_mass, [1.31887368e05, 4.35408760e02, 5.93227761e-01]) + + +def test_Stochiometry_tissue_deficit(fxt_stochiometry_model): + """Test the tissue_deficit method of the Stochiometry class.""" + + # Calculate the tissue deficit + tissue_deficit = fxt_stochiometry_model.tissue_deficit + + expected_deficit = np.array([0.0, 0.0, 0.0]) + assert np.allclose(tissue_deficit, expected_deficit) + + +@pytest.mark.skip( + reason="Negative growth problem. Return to this test when that problem is fixed." +) +def test_Stochiometry_account_for_growth(fxt_plants_model, fxt_stochiometry_model): + """Test the account_for_growth method of the Stochiometry class.""" + + from virtual_ecosystem.models.plants.stochiometry import StemAllocation + + # Create a StemAllocation object + cell_id = 0 # Assuming we are testing for the first cell + community = fxt_plants_model.communities[cell_id] + fxt_plants_model.per_stem_gpp = { + cell_id: np.array([55]) for cell_id in fxt_plants_model.communities.keys() + } + stem_allocation = StemAllocation( + stem_traits=community.stem_traits, + stem_allometry=community.stem_allometry, + whole_crown_gpp=fxt_plants_model.per_stem_gpp[cell_id], + ) + + assert np.allclose( + fxt_stochiometry_model.total_element_mass, + [1.31887368e05, 4.35408760e02, 5.93227761e-01], + ) + + # Account for growth + fxt_stochiometry_model.account_for_growth(stem_allocation) + + print("new mass") + + assert np.all( + fxt_stochiometry_model.total_element_mass + >= [1.31887368e05, 4.35408760e02, 5.93227761e-01] + ) + + +def test_Stochiometry_account_for_element_loss_turnover( + fxt_plants_model, fxt_stochiometry_model +): + """Test the account_for_element_loss_turnover method of the Stochiometry class.""" + + from virtual_ecosystem.models.plants.stochiometry import StemAllocation + + # Create a StemAllocation object + cell_id = 0 # Assuming we are testing for the first cell + community = fxt_plants_model.communities[cell_id] + fxt_plants_model.per_stem_gpp = { + cell_id: np.array([55]) for cell_id in fxt_plants_model.communities.keys() + } + stem_allocation = StemAllocation( + stem_traits=community.stem_traits, + stem_allometry=community.stem_allometry, + whole_crown_gpp=fxt_plants_model.per_stem_gpp[cell_id], + ) + + assert np.allclose(fxt_stochiometry_model.element_surplus, [0.0, 0.0, 0.0]) + fxt_stochiometry_model.account_for_element_loss_turnover(stem_allocation) + + # The surplus should always be negative after accounting for turnover + assert np.all(fxt_stochiometry_model.element_surplus <= [0.0, 0.0, 0.0]) + + +def test_Stochiometry_distribute_deficit(fxt_stochiometry_model): + """Test the distribute_deficit method of the Stochiometry class.""" + + # Distribute the deficit + fxt_stochiometry_model.distribute_deficit(0) + + # TODO: finish this test + + +def test_Stochiometry_distrubte_surplus(fxt_stochiometry_model): + """Test the distribute_surplus method of the Stochiometry class.""" + + fxt_stochiometry_model.distribute_surplus(0) + + # TODO: finish this test From 72d8d201981ab8af590cfa0f06267984adbaaac4 Mon Sep 17 00:00:00 2001 From: Sally Date: Tue, 17 Jun 2025 15:29:45 +0100 Subject: [PATCH 15/22] Add asserts for final stochiometry tests. --- tests/models/plants/test_stochiometry.py | 16 ++++++++++++---- virtual_ecosystem/models/plants/stochiometry.py | 10 ++++++++-- 2 files changed, 20 insertions(+), 6 deletions(-) diff --git a/tests/models/plants/test_stochiometry.py b/tests/models/plants/test_stochiometry.py index b36e661df..f7b723688 100644 --- a/tests/models/plants/test_stochiometry.py +++ b/tests/models/plants/test_stochiometry.py @@ -253,17 +253,25 @@ def test_Stochiometry_account_for_element_loss_turnover( def test_Stochiometry_distribute_deficit(fxt_stochiometry_model): - """Test the distribute_deficit method of the Stochiometry class.""" + """Test the distribute_deficit method of the Stochiometry class. + + NOTE: This method should have a more robust test with a more deliberate test value. + This will be easier to implement once growth is functioning. + """ # Distribute the deficit fxt_stochiometry_model.distribute_deficit(0) - # TODO: finish this test + assert np.allclose(fxt_stochiometry_model.element_surplus[0], [0]) def test_Stochiometry_distrubte_surplus(fxt_stochiometry_model): - """Test the distribute_surplus method of the Stochiometry class.""" + """Test the distribute_surplus method of the Stochiometry class. + + NOTE: This method should have a more robust test with a more deliberate test value. + This will be easier to implement once growth is functioning. + """ fxt_stochiometry_model.distribute_surplus(0) - # TODO: finish this test + assert fxt_stochiometry_model.element_surplus[0] >= 0 diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index 390941150..65a5aaa1d 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -368,7 +368,13 @@ def account_for_element_loss_turnover(self, allocation: StemAllocation) -> None: self.element_surplus -= tissue.element_turnover(allocation) def distribute_deficit(self, cohort: int) -> None: - """Distribute the nitrogen deficit across the tissue types. + """Distribute the element deficit across the tissue types. + + During the update, the information about a surplus/deficit of element are stored + in the element_surplus. If there is a deficit (represented by a negative element + surplus), this method distributes the deficit across the tissue types. Then, the + element surplus is reset to 0. The deficit is distributed in proportion to the + total element mass of each tissue type for that cohort. Args: cohort: The cohort to reconcile deficit. @@ -382,7 +388,7 @@ def distribute_deficit(self, cohort: int) -> None: / self.total_element_mass[cohort] ) - self.element_surplus[cohort] += deficit + self.element_surplus[cohort] = 0 def distribute_surplus(self, cohort: int) -> None: """Distribute the nitrogen surplus across the tissue types for a single cohort. From 2e1231baab2ad7a80f01076597db976ce14a09fb Mon Sep 17 00:00:00 2001 From: Sally Matson Date: Tue, 1 Jul 2025 11:32:47 +0100 Subject: [PATCH 16/22] Update virtual_ecosystem/models/plants/stochiometry.py Co-authored-by: arne-scheire <97458577+arne-exe@users.noreply.github.com> --- virtual_ecosystem/models/plants/stochiometry.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index 65a5aaa1d..e73983580 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -103,7 +103,7 @@ class FoliageTissue(Tissue): """A class to hold foliage stochiometry data for a set of plant cohorts.""" # reclaim_ratio: NDArray[np.float64] - """The ratio of the element that can be r eclaimed from the sensced tissue.""" + """The ratio of the element that can be reclaimed from the senesced tissue.""" def __init__( self, From 90fd88eb2a52ad42ce73fcea3b6bbc130ff4d1a4 Mon Sep 17 00:00:00 2001 From: Sally Matson Date: Tue, 1 Jul 2025 15:42:42 +0100 Subject: [PATCH 17/22] Apply suggestions from code review Co-authored-by: arne-scheire <97458577+arne-exe@users.noreply.github.com> Co-authored-by: David Orme --- .../models/plants/plants_model.py | 43 +++++++------------ .../models/plants/stochiometry.py | 6 +-- 2 files changed, 18 insertions(+), 31 deletions(-) diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index 7a961bff6..9e0f4ea9b 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -914,37 +914,24 @@ def allocate_gpp(self) -> None: cohorts.dbh_values = cohorts.dbh_values + stem_allocation.delta_dbh # Subtract the N/P required from growth from the element store, and - # redistribute it to the individual tisuses. + # redistribute it to the individual tissues. for stochiometry in stochiometries.values(): stochiometry.account_for_growth(stem_allocation) # Balance the N & P surplus/deficit with the symbiote carbon supply - n_weighted_avg = np.dot( - self.data["ecto_supply_limit_n"][cell_id] - + self.data["arbuscular_supply_limit_n"][cell_id], - cohorts.n_individuals, - ) - n_avaiable_per_cohort = n_weighted_avg / sum(n_weighted_avg) - n_available_per_stem = np.divide( - n_avaiable_per_cohort, cohorts.n_individuals - ) - stochiometries["N"].element_surplus = ( - stochiometries["N"].element_surplus + n_available_per_stem - ) - - p_weighted_avg = np.dot( - self.data["ecto_supply_limit_p"][cell_id] - + self.data["arbuscular_supply_limit_p"][cell_id], - cohorts.n_individuals, - ) - p_available_per_cohort = p_weighted_avg / sum(p_weighted_avg) - p_available_per_stem = np.divide( - p_available_per_cohort, cohorts.n_individuals - ) - - stochiometries["P"].element_surplus = ( - stochiometries["P"].element_surplus + p_available_per_stem - ) + for element in ["N", "P"]: + element_weighted_avg = np.dot( + self.data["ecto_supply_limit_" + element.lower()][cell_id] + + self.data["arbuscular_supply_limit_" + element.lower()][cell_id], + cohorts.n_individuals, + ) + element_available_per_cohort = ( + element_weighted_avg / sum(element_weighted_avg) + ) + element_available_per_stem = np.divide( + element_available_per_cohort, cohorts.n_individuals + ) + stochiometries[element].element_surplus += element_available_per_stem # Cohort by cohort, distribute the surplus/deficit across the tissue types for cohort in range(len(cohorts.n_individuals)): @@ -1107,7 +1094,7 @@ def calculate_nutrient_uptake(self) -> None: for cell_id in self.communities.keys(): self.stochiometries[cell_id]["N"].element_surplus += ( self.per_stem_transpiration[cell_id] - * (1.8015 * pow(10.0, -11)) + * 1.8015e-11 * ( self.data["plant_ammonium_uptake"][cell_id] + self.data["plant_nitrate_uptake"][cell_id] diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index e73983580..ccdf3e119 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -1,12 +1,12 @@ """The :mod:`~virtual_ecosystem.models.plants.stochiometry` module contains the class for managing plant cohort stochiometry ratios. The carbon mass is stored in plant -alloemetry or allocation, so this class uses those as the anchor weights and stores +allometry or allocation, so this class uses those as the anchor weights and stores CN and CP ratios. The class holds current CN and CP ratios for foliage, reproductive tissue, wood, and -roots on the cohort level. Each tissue also has an idea CN and CP ratio, which is used +roots on the cohort level. Each tissue also has an ideal CN and CP ratio, which is used as a comparison in the case of any nutrient deficit. Senesced leaves also have fixed CN -and CP ratios, which is used for leaf turnover. +and CP ratios, which are used for leaf turnover. In the future, the ideal CN and CP ratios will be PFT traits. """ # noqa: D205 From 048b02dc286e716c317e4b0339cecd722f2a4630 Mon Sep 17 00:00:00 2001 From: Sally Date: Wed, 2 Jul 2025 15:57:40 +0100 Subject: [PATCH 18/22] Comments from PR, including changing Tissue to ABC. --- tests/models/plants/test_stochiometry.py | 2 +- .../models/plants/plants_model.py | 209 +++++++++--------- .../models/plants/stochiometry.py | 167 +++++--------- 3 files changed, 162 insertions(+), 216 deletions(-) diff --git a/tests/models/plants/test_stochiometry.py b/tests/models/plants/test_stochiometry.py index f7b723688..bed19863f 100644 --- a/tests/models/plants/test_stochiometry.py +++ b/tests/models/plants/test_stochiometry.py @@ -186,7 +186,7 @@ def test_Stochiometry_tissue_deficit(fxt_stochiometry_model): # Calculate the tissue deficit tissue_deficit = fxt_stochiometry_model.tissue_deficit - expected_deficit = np.array([0.0, 0.0, 0.0]) + expected_deficit = np.array([1.01616593e03, 3.30162938e01, 3.88354401e-01]) assert np.allclose(tissue_deficit, expected_deficit) diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index 9e0f4ea9b..187158040 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -324,114 +324,118 @@ def _setup( data=self.data, flora=self.flora, grid=self.grid ) - # Initialize the stochiometries of each cohort - self.stochiometries = { - cell_id: { - "N": StemStochiometry( - element="N", - tissues=[ - FoliageTissue( - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.foliage_c_n_ratio, - ), - actual_element_mass=self.communities[ - cell_id - ].stem_allometry.foliage_mass - * model_constants.foliage_c_n_ratio, - reclaim_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.leaf_turnover_c_n_ratio, - ), + # Initialize the stochiometries of each cohort. Each StemStochiometry object + # contains a list of StemTissue objects, which are the tissues that make up the + # stochiometry of the stem. The initial values for N and P are based on the + # ideal stochiometric ratios defined in the PlantsConsts class. + # TODO: #697 - these need to be configurable + self.stochiometries = {} + for cell_id in self.communities.keys(): + self.stochiometries[cell_id] = {} + self.stochiometries[cell_id]["N"] = StemStochiometry( + element="N", + tissues=[ + FoliageTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.foliage_c_n_ratio, ), - RootTissue( - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.root_turnover_c_n_ratio, - ), - actual_element_mass=model_constants.root_turnover_c_n_ratio - * self.communities[cell_id].stem_traits.zeta - * self.communities[cell_id].stem_allometry.foliage_mass, + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.foliage_mass + * model_constants.foliage_c_n_ratio, + reclaim_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.leaf_turnover_c_n_ratio, ), - WoodTissue( - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.deadwood_c_n_ratio, - ), - actual_element_mass=model_constants.deadwood_c_n_ratio - * self.communities[cell_id].stem_allometry.stem_mass, + ), + RootTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.root_turnover_c_n_ratio, ), - ReproductiveTissue( - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.plant_reproductive_tissue_turnover_c_n_ratio, - ), - actual_element_mass=self.communities[ - cell_id - ].stem_allometry.reproductive_tissue_mass - * self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio, # noqa: E501 + actual_element_mass=model_constants.root_turnover_c_n_ratio + * self.communities[cell_id].stem_traits.zeta + * self.communities[cell_id].stem_allometry.foliage_mass + * self.communities[cell_id].stem_traits.sla, + ), + WoodTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.deadwood_c_n_ratio, ), - ], - community=self.communities[cell_id], - ), - "P": StemStochiometry( - element="P", - tissues=[ - FoliageTissue( - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.foliage_c_p_ratio, - ), - actual_element_mass=self.communities[ - cell_id - ].stem_allometry.foliage_mass - * model_constants.foliage_c_p_ratio, - reclaim_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.leaf_turnover_c_p_ratio, - ), + actual_element_mass=model_constants.deadwood_c_n_ratio + * self.communities[cell_id].stem_allometry.stem_mass, + ), + ReproductiveTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.plant_reproductive_tissue_turnover_c_n_ratio, ), - RootTissue( - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.root_turnover_c_p_ratio, - ), - actual_element_mass=model_constants.root_turnover_c_p_ratio - * self.communities[cell_id].stem_traits.zeta - * self.communities[cell_id].stem_allometry.foliage_mass, + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.reproductive_tissue_mass + * self.model_constants.plant_reproductive_tissue_turnover_c_n_ratio, # noqa: E501 + ), + ], + community=self.communities[cell_id], + ) + self.stochiometries[cell_id]["P"] = StemStochiometry( + element="P", + tissues=[ + FoliageTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.foliage_c_p_ratio, ), - WoodTissue( - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.deadwood_c_p_ratio, - ), - actual_element_mass=model_constants.deadwood_c_p_ratio - * self.communities[cell_id].stem_allometry.stem_mass, + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.foliage_mass + * model_constants.foliage_c_p_ratio, + reclaim_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.leaf_turnover_c_p_ratio, ), - ReproductiveTissue( - community=self.communities[cell_id], - ideal_ratio=np.full( - self.communities[cell_id].n_cohorts, - model_constants.plant_reproductive_tissue_turnover_c_p_ratio, - ), - actual_element_mass=self.communities[ - cell_id - ].stem_allometry.reproductive_tissue_mass - * self.model_constants.plant_reproductive_tissue_turnover_c_p_ratio, # noqa: E501 + ), + RootTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.root_turnover_c_p_ratio, ), - ], - community=self.communities[cell_id], - ), - } - for cell_id in self.communities.keys() - } + actual_element_mass=model_constants.root_turnover_c_p_ratio + * self.communities[cell_id].stem_traits.zeta + * self.communities[cell_id].stem_allometry.foliage_mass + * self.communities[cell_id].stem_traits.sla, + ), + WoodTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.deadwood_c_p_ratio, + ), + actual_element_mass=model_constants.deadwood_c_p_ratio + * self.communities[cell_id].stem_allometry.stem_mass, + ), + ReproductiveTissue( + community=self.communities[cell_id], + ideal_ratio=np.full( + self.communities[cell_id].n_cohorts, + model_constants.plant_reproductive_tissue_turnover_c_p_ratio, + ), + actual_element_mass=self.communities[ + cell_id + ].stem_allometry.reproductive_tissue_mass + * self.model_constants.plant_reproductive_tissue_turnover_c_p_ratio, # noqa: E501 + ), + ], + community=self.communities[cell_id], + ) # This is widely used internally so store it as an attribute. self._canopy_layer_indices = self.layer_structure.index_canopy @@ -925,8 +929,8 @@ def allocate_gpp(self) -> None: + self.data["arbuscular_supply_limit_" + element.lower()][cell_id], cohorts.n_individuals, ) - element_available_per_cohort = ( - element_weighted_avg / sum(element_weighted_avg) + element_available_per_cohort = element_weighted_avg / sum( + element_weighted_avg ) element_available_per_stem = np.divide( element_available_per_cohort, cohorts.n_individuals @@ -1090,6 +1094,7 @@ def calculate_nutrient_uptake(self) -> None: # - Conversion factor from µmol H2O to m^3 (1.08015*10^-11) # - Concentration of N/P uptake (kg m^-3) # - Kg to g (1000) + # TODO: scale by atmospheric pressure and temperature (#927) for cell_id in self.communities.keys(): self.stochiometries[cell_id]["N"].element_surplus += ( diff --git a/virtual_ecosystem/models/plants/stochiometry.py b/virtual_ecosystem/models/plants/stochiometry.py index ccdf3e119..7e63e04d6 100644 --- a/virtual_ecosystem/models/plants/stochiometry.py +++ b/virtual_ecosystem/models/plants/stochiometry.py @@ -11,6 +11,7 @@ In the future, the ideal CN and CP ratios will be PFT traits. """ # noqa: D205 +from abc import ABC, abstractmethod from dataclasses import dataclass, field import numpy as np @@ -21,7 +22,7 @@ @dataclass -class Tissue(PandasExporter, CohortMethods): +class Tissue(ABC): """A dataclass to hold tissue stochiometry data for a set of plant cohorts. This class holds the current quantity of a given element (generally N or P) for a @@ -30,11 +31,9 @@ class Tissue(PandasExporter, CohortMethods): an entry for each cohort in the data class. """ - tissue_name: str - """The name of the tissue type.""" community: Community """The community object that the tissue is associated with.""" - # Should this be stored only in stochiometry and not in tissue? + # TODO: consider where best to store shared attributes like community. ideal_ratio: NDArray[np.float64] """The ideal ratio of the element for the tissue type.""" @@ -42,22 +41,9 @@ class Tissue(PandasExporter, CohortMethods): """The actual mass of the element for the tissue type.""" def __post_init__(self) -> None: - """Initialize the Tissue object.""" + """Post-initialization to properly format the actual element mass.""" self.actual_element_mass = self.actual_element_mass.squeeze() - @property - def carbon_mass(self) -> NDArray[np.float64]: - """Get the carbon mass for the tissue type. - - This method should be implemented by subclasses to return the carbon mass for - the specific tissue type. - - Returns: - The carbon mass for the specified tissue. - """ - # This method should be implemented by subclasses - raise NotImplementedError("Carbon mass must be defined in subclasses.") - @property def deficit(self) -> NDArray[np.float64]: """Calculate the element deficit for the tissue type. @@ -67,28 +53,6 @@ def deficit(self) -> NDArray[np.float64]: """ return self.ideal_ratio * self.carbon_mass - self.actual_element_mass - def element_needed_for_growth( - self, allocation: StemAllocation - ) -> NDArray[np.float64]: - """Calculate the element needed for growth for the tissue type. - - Returns: - The element needed for growth for the specified tissue. - """ - raise NotImplementedError( - "Element needed for growth must be defined in subclasses." - ) - - def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: - """Calculate the element lost to turnover for the tissue type. - - Returns: - The element lost to turnover for the specified tissue. - """ - raise NotImplementedError( - "Element needed for growth must be defined in subclasses." - ) - @property def Cx_ratio(self) -> NDArray[np.float64]: """Get the carbon to element ratio for the tissue type. @@ -98,29 +62,29 @@ def Cx_ratio(self) -> NDArray[np.float64]: """ return self.carbon_mass / self.actual_element_mass + @property + @abstractmethod + def carbon_mass(self) -> NDArray[np.float64]: + """Calculate the carbon mass for the tissue type.""" + + @abstractmethod + def element_needed_for_growth( + self, allocation: StemAllocation + ) -> NDArray[np.float64]: + """Calculate the element needed for growth for the tissue type.""" + + @abstractmethod + def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: + """Calculate the element lost to turnover for the tissue type.""" + +@dataclass class FoliageTissue(Tissue): """A class to hold foliage stochiometry data for a set of plant cohorts.""" - # reclaim_ratio: NDArray[np.float64] + reclaim_ratio: NDArray[np.float64] """The ratio of the element that can be reclaimed from the senesced tissue.""" - def __init__( - self, - community: Community, - ideal_ratio: NDArray[np.float64], - actual_element_mass: NDArray[np.float64], - reclaim_ratio: NDArray[np.float64], - ): - super().__init__( - tissue_name="Foliage", - community=community, - ideal_ratio=ideal_ratio, - actual_element_mass=actual_element_mass, - ) - self.reclaim_ratio = reclaim_ratio - """The ratio of the element that can be reclaimed from the senesced tissue.""" - @property def carbon_mass(self) -> NDArray[np.float64]: """Get the carbon mass for foliage tissue. @@ -133,10 +97,10 @@ def carbon_mass(self) -> NDArray[np.float64]: def element_needed_for_growth( self, allocation: StemAllocation ) -> NDArray[np.float64]: - """Calculate the nitrogen needed for growth for foliage tissue. + """Calculate the element quantity needed for growth for foliage tissue. Returns: - The nitrogen needed for growth for foliage tissue. + The element quantity needed for growth for foliage tissue. """ return (allocation.delta_foliage_mass * (1 / self.ideal_ratio)).squeeze() @@ -144,7 +108,7 @@ def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: """Calculate the element mass lost to turnover for foliage tissue. Returns: - The nitrogen lost to turnover for foliage tissue. + The element quantity lost to turnover for foliage tissue. """ return ( allocation.foliage_turnover @@ -152,22 +116,10 @@ def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: ).squeeze() +@dataclass class ReproductiveTissue(Tissue): """Holds reproductive tissue stochiometry data for a set of plant cohorts.""" - def __init__( - self, - community: Community, - ideal_ratio: NDArray[np.float64], - actual_element_mass: NDArray[np.float64], - ): - super().__init__( - tissue_name="Reproductive", - community=community, - ideal_ratio=ideal_ratio, - actual_element_mass=actual_element_mass, - ) - @property def carbon_mass(self) -> NDArray[np.float64]: """Get the carbon mass for reproductive tissue. @@ -180,10 +132,10 @@ def carbon_mass(self) -> NDArray[np.float64]: def element_needed_for_growth( self, allocation: StemAllocation ) -> NDArray[np.float64]: - """Calculate the nitrogen needed for growth for reproductive tissue. + """Calculate the element needed for growth for reproductive tissue. Returns: - The nitrogen needed for growth for reproductive tissue. + The element quantity needed for growth for reproductive tissue. """ return ( allocation.delta_foliage_mass @@ -200,22 +152,10 @@ def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: return (allocation.reproductive_tissue_turnover * (1 / self.Cx_ratio)).squeeze() +@dataclass class WoodTissue(Tissue): """A class to hold wood stochiometry data for a set of plant cohorts.""" - def __init__( - self, - community: Community, - ideal_ratio: NDArray[np.float64], - actual_element_mass: NDArray[np.float64], - ): - super().__init__( - tissue_name="Wood", - community=community, - ideal_ratio=ideal_ratio, - actual_element_mass=actual_element_mass, - ) - @property def carbon_mass(self) -> NDArray[np.float64]: """Get the carbon mass for wood tissue. @@ -228,10 +168,10 @@ def carbon_mass(self) -> NDArray[np.float64]: def element_needed_for_growth( self, allocation: StemAllocation ) -> NDArray[np.float64]: - """Calculate the nitrogen needed for growth for wood tissue. + """Calculate the element needed for growth for wood tissue. Returns: - The nitrogen needed for growth for wood tissue. + The element needed for growth for wood tissue. """ return (allocation.delta_stem_mass * (1 / self.ideal_ratio)).squeeze() @@ -244,22 +184,10 @@ def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: return np.zeros(self.community.n_cohorts) +@dataclass class RootTissue(Tissue): """A class to hold root stochiometry data for a set of plant cohorts.""" - def __init__( - self, - community: Community, - ideal_ratio: NDArray[np.float64], - actual_element_mass: NDArray[np.float64], - ): - super().__init__( - tissue_name="Roots", - community=community, - ideal_ratio=ideal_ratio, - actual_element_mass=actual_element_mass, - ) - @property def carbon_mass(self) -> NDArray[np.float64]: """Get the carbon mass for root tissue. @@ -268,21 +196,31 @@ def carbon_mass(self) -> NDArray[np.float64]: The carbon mass for root tissue. """ return ( - self.community.stem_allometry.foliage_mass * self.community.stem_traits.zeta + self.community.stem_allometry.foliage_mass + * self.community.stem_traits.zeta + * self.community.stem_traits.sla ).squeeze() def element_needed_for_growth( self, allocation: StemAllocation ) -> NDArray[np.float64]: - """Calculate the nitrogen needed for growth for root tissue. + """Calculate the element needed for growth for root tissue. + + The calculation is the NC ratio (1 / CN ratio) multiplied by the change in root + mass (change in foliage mass * zeta * SLA). + + Delta foliage mass (g C) + Zeta: Ratio of fine-root mass to foliage area (kg C / m2) + SLA: Specific leaf area (m2 / kg C) Returns: - The nitrogen needed for growth for root tissue. + The element needed for growth for root tissue. """ return ( - allocation.delta_foliage_mass - * (1 / self.ideal_ratio) + (1 / self.ideal_ratio) + * allocation.delta_foliage_mass * self.community.stem_traits.zeta + * self.community.stem_traits.sla ).squeeze() def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: @@ -296,17 +234,20 @@ def element_turnover(self, allocation: StemAllocation) -> NDArray[np.float64]: @dataclass class StemStochiometry(CohortMethods, PandasExporter): - """A class holding the ratios of Carbon to Nitrogen and Phosphorous for stems. + """A class holding elemental weights for a set of plant cohorts and tissues. This class holds the current ratios across tissue type for a community object, which in essence is a series of cohorts. It acts in parallel with StemAllometry, a class attribute of Community. + + The class is designed to be element-agnostic, so it can be used for any element as + required. """ element: str - """The name of the element (e.g., N or P).""" + """The name of the element.""" tissues: list[Tissue] - """Tissues for the associated stems.""" + """Tissues for the associated cohorts.""" community: Community """The community object that the stochiometry is associated with.""" element_surplus: NDArray[np.float64] = field(init=False) @@ -321,7 +262,7 @@ def total_element_mass(self) -> NDArray[np.float64]: """Calculate the total element mass for each cohort. Returns: - The total nitrogen mass for each cohort. + The total element mass for each cohort. """ mass = np.zeros(self.community.n_cohorts) for tissue in self.tissues: @@ -391,7 +332,7 @@ def distribute_deficit(self, cohort: int) -> None: self.element_surplus[cohort] = 0 def distribute_surplus(self, cohort: int) -> None: - """Distribute the nitrogen surplus across the tissue types for a single cohort. + """Distribute the element surplus across the tissue types for a single cohort. Args: cohort: The cohort to reconcile surplus. From 80589d60de540e682a1c636292fc4ef02a401896 Mon Sep 17 00:00:00 2001 From: Sally Matson Date: Thu, 3 Jul 2025 10:58:35 +0100 Subject: [PATCH 19/22] Increasing timeout --- docs/source/using_the_ve/virtual_ecosystem_in_use.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/source/using_the_ve/virtual_ecosystem_in_use.md b/docs/source/using_the_ve/virtual_ecosystem_in_use.md index 05a1c3f0d..aee9e1a81 100644 --- a/docs/source/using_the_ve/virtual_ecosystem_in_use.md +++ b/docs/source/using_the_ve/virtual_ecosystem_in_use.md @@ -1,6 +1,6 @@ --- execution: - timeout: 90 + timeout: 120 jupytext: formats: md:myst text_representation: From 61b786ad926df81935fd2292166bc7b0f1d2e292 Mon Sep 17 00:00:00 2001 From: Sally Date: Thu, 3 Jul 2025 14:54:41 +0100 Subject: [PATCH 20/22] Fix dbh bug. --- virtual_ecosystem/models/plants/plants_model.py | 5 ----- 1 file changed, 5 deletions(-) diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index f3a03da92..6f822a6c3 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -965,11 +965,6 @@ def allocate_gpp(self) -> None: ) ) - # FOURTH, ALLOCATE TO GROWTH: - # Grow the plants by increasing the stem dbh - # TODO: dimension mismatch (1d vs 2d array) - check in pyrealm - cohorts.dbh_values = cohorts.dbh_values + stem_allocation.delta_dbh - # Subtract the N/P required from growth from the element store, and # redistribute it to the individual tissues. for stochiometry in stochiometries.values(): From 9c9566af85bc0313388a50b731a9fd9bb91abd84 Mon Sep 17 00:00:00 2001 From: Sally Date: Thu, 3 Jul 2025 14:57:27 +0100 Subject: [PATCH 21/22] Updating comments --- virtual_ecosystem/models/plants/plants_model.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/virtual_ecosystem/models/plants/plants_model.py b/virtual_ecosystem/models/plants/plants_model.py index 6f822a6c3..208354665 100644 --- a/virtual_ecosystem/models/plants/plants_model.py +++ b/virtual_ecosystem/models/plants/plants_model.py @@ -818,7 +818,7 @@ def allocate_gpp(self) -> None: turnover values. """ - # First, initialize all turnover variables to 0 with the proper dimensions. + # Initialize all turnover variables to 0 with the proper dimensions. # Most variables are merged across PFTs and cohorts - one pool per cell. self.data["leaf_turnover"] = xr.full_like(self.data["elevation"], 0) self.data["root_turnover"] = xr.full_like(self.data["elevation"], 0) @@ -854,7 +854,7 @@ def allocate_gpp(self) -> None: whole_crown_gpp=self.per_stem_gpp[cell_id], ) - # FIRST, ALLOCATE TO TURNOVER: + # ALLOCATE TO TURNOVER: # Grow the plants by increasing the stem dbh # TODO: dimension mismatch (1d vs 2d array) - check in pyrealm # HACK: The current code prevents stems shrinking to zero and below. This is @@ -938,11 +938,11 @@ def allocate_gpp(self) -> None: canopy_non_propagule_mass * cohort_n_stems ) - # SECOND: ALLOCATE N TO REGROW WHAT WAS LOST TO TURNOVER + # ALLOCATE N TO REGROW WHAT WAS LOST TO TURNOVER for stochiometry in stochiometries.values(): stochiometry.account_for_element_loss_turnover(stem_allocation) - # THIRD, ALLOCATE GPP TO ACTIVE NUTRIENT PATHWAYS: + # ALLOCATE GPP TO ACTIVE NUTRIENT PATHWAYS: # Allocate the topsliced GPP to root exudates with remainder as active # nutrient pathways self.data["root_carbohydrate_exudation"][cell_id] = ( From a847c3899d34841c4a2cdb9ab206fde64e061643 Mon Sep 17 00:00:00 2001 From: Sally Date: Mon, 7 Jul 2025 15:57:13 +0100 Subject: [PATCH 22/22] Static timeout error. --- docs/source/using_the_ve/virtual_ecosystem_in_static_mode.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/source/using_the_ve/virtual_ecosystem_in_static_mode.md b/docs/source/using_the_ve/virtual_ecosystem_in_static_mode.md index d2ca0e894..cafb967fc 100644 --- a/docs/source/using_the_ve/virtual_ecosystem_in_static_mode.md +++ b/docs/source/using_the_ve/virtual_ecosystem_in_static_mode.md @@ -1,6 +1,6 @@ --- execution: - timeout: 90 + timeout: 120 jupytext: formats: md:myst main_language: python