Unpin the update interval in hydrology

.pre-commit-config.yaml

@@ -5,16 +5,16 @@ repos:
       - id: check-merge-conflict
       - id: debug-statements
   - repo: https://github.com/astral-sh/ruff-pre-commit
-    rev: v0.12.5
+    rev: v0.12.7
     hooks:
       - id: ruff-check   #  Run the linter.
         args: [--fix, --exit-non-zero-on-fix]
       - id: ruff-format  #  Run the formatter.
   - repo: https://github.com/pre-commit/mirrors-mypy
-    rev: "v1.17.0"
+    rev: "v1.17.1"
     hooks:
       - id: mypy
-        additional_dependencies: [types-jsonschema, xarray, types-tabulate, numpy==2.1.3]
+        additional_dependencies: [types-jsonschema, xarray, types-tabulate, numpy]
   - repo: https://github.com/igorshubovych/markdownlint-cli
     rev: v0.45.0
     hooks:

docs/source/virtual_ecosystem/implementation/hydrology_implementation.md

@@ -47,10 +47,15 @@ water enters and leaves the plant model.
 :::
 
 ```{note}
-Many of the underlying processes are problematic at a monthly timestep, which is
-currently the only supported update interval. As a short-term work around, the input
-precipitation is randomly distributed over 30 days and input canopy transpiration is
-divided by 30, and the return variables are monthly means or monthly accumulated values.
+Calculating hydrological processes at coarse time scales is problematic and so the
+hydrology model uses an internal daily time step to model hydrology. When the model is
+run with a coarser update interval - and hence the precipitation and transpiration
+inputs are totals over more than one day - the hydrology model partitions the input data
+ into daily values. Precipitation is randomly partitioned between days and the total
+ transpiration is evenly divided across days.
+
+The values returned by the hydrology model are then monthly means or monthly accumulated
+values.
 ```
 
 ## Required variables

pyproject.toml

@@ -45,6 +45,7 @@ tomli = {version = "^2.0.1", python = "<3.11"}
 tomli-w = "^1.0.0"
 tqdm = "^4.66.2"
 xarray = ">=2024.6,<2026.0"
+types-tabulate = "^0.9.0.20241207"
 
 [tool.poetry.group.types.dependencies]
 types-dataclasses = "^0.6.6"
@@ -134,7 +135,6 @@ convention = "google"
 
 [tool.mypy]
 ignore_missing_imports = true
-plugins = "numpy.typing.mypy_plugin"
 
 [tool.jupytext]
 # Stop jupytext from removing mystnb and other settings in MyST Notebook YAML headers

tests/models/hydrology/test_hydrology_model.py

@@ -223,26 +223,68 @@ def test_generate_hydrology_model(
 
 
 @pytest.mark.parametrize(
-    "update_interval, raises",
+    "update_interval, raises, expected_2d, expected_1d",
     [
         pytest.param(
             pint.Quantity(1, "month"),
             does_not_raise(),
-            id="updates correctly",
+            {
+                "soil_moisture": [
+                    [248.938056, 248.937037, 248.935933, 248.936385],
+                    [218.994795, 218.994795, 218.994795, 218.994795],
+                ],
+                "matric_potential": [
+                    [-56.432398, -56.438614, -56.44538, -56.442609],
+                    [-217.596626, -217.596626, -217.596626, -217.596626],
+                ],
+                "vertical_flow": [
+                    [0.00017, 0.00017, 0.00017, 0.00017],
+                    [0.000526, 0.000526, 0.000526, 0.000526],
+                ],
+            },
+            {
+                "total_river_discharge": [0, 0, 67002, 22095],
+                "surface_runoff": [20.343781, 20.66599, 20.896484, 20.443394],
+                "surface_runoff_accumulated": [0, 0, 1470, 330],
+                "soil_evaporation": [5.870856, 5.870856, 5.870856, 5.870856],
+            },
+            id="1 month",
         ),
         pytest.param(
             pint.Quantity(1, "week"),
-            pytest.raises(NotImplementedError),
-            id="incorrect update frequency",
+            does_not_raise(),
+            {
+                "soil_moisture": [
+                    [249.628102, 249.628102, 249.628102, 249.628102],
+                    [215.713193, 215.713193, 215.713193, 215.713193],
+                ],
+                "matric_potential": [
+                    [-52.085807, -52.085807, -52.085807, -52.085807],
+                    [-196.720556, -196.720556, -196.720556, -196.720556],
+                ],
+                "vertical_flow": [
+                    [0.000295, 0.000295, 0.000295, 0.000295],
+                    [0.000611, 0.000611, 0.000611, 0.000611],
+                ],
+            },
+            {
+                "total_river_discharge": [0, 0, 5767, 1910],
+                "surface_runoff": [163.019971, 163.019971, 163.019971, 163.019971],
+                "surface_runoff_accumulated": [0, 0, 3395, 1127],
+                "soil_evaporation": [1.388223, 1.388223, 1.388223, 1.388223],
+            },
+            id="1 week",
         ),
     ],
 )
 def test_setup(
+    fixture_core_components,
     dummy_climate_data,
     fixture_config,
     update_interval,
     raises,
-    fixture_core_components,
+    expected_2d,
+    expected_1d,
 ):
     """Test set up and update."""
     from virtual_ecosystem.core.core_components import CoreComponents
@@ -302,21 +344,6 @@ def test_setup(
             model.update(time_index=1, seed=42)
 
             # Test 2d variables
-            expected_2d = {
-                "soil_moisture": [
-                    [248.938056, 248.937037, 248.935933, 248.936385],
-                    [218.994795, 218.994795, 218.994795, 218.994795],
-                ],
-                "matric_potential": [
-                    [-56.432398, -56.438614, -56.44538, -56.442609],
-                    [-217.596626, -217.596626, -217.596626, -217.596626],
-                ],
-                "vertical_flow": [
-                    [0.00017, 0.00017, 0.00017, 0.00017],
-                    [0.000526, 0.000526, 0.000526, 0.000526],
-                ],
-            }
-
             for var_name, expected_vals in expected_2d.items():
                 exp_var = lyr_strct.from_template()
                 exp_var[soil_indices] = expected_vals
@@ -329,12 +356,6 @@ def test_setup(
                 )
 
             # Test one dimensional variables
-            expected_1d = {
-                "total_river_discharge": [0, 0, 67002, 22095],
-                "surface_runoff": [20.343781, 20.66599, 20.896484, 20.443394],
-                "surface_runoff_accumulated": [0, 0, 1470, 330],
-                "soil_evaporation": [5.870856, 5.870856, 5.870856, 5.870856],
-            }
 
             for var_name, expected_vals in expected_1d.items():
                 np.testing.assert_allclose(

virtual_ecosystem/core/grid.py

@@ -433,14 +433,14 @@ def get_distances(
         """
 
         if cell_from is None:
-            _cell_from = np.arange(self.n_cells)
+            _cell_from: NDArray[np.int_] = np.arange(self.n_cells)
         else:
             _cell_from = np.array(
                 [cell_from] if isinstance(cell_from, int) else cell_from
             )
 
         if cell_to is None:
-            _cell_to = np.arange(self.n_cells)
+            _cell_to: NDArray[np.int_] = np.arange(self.n_cells)
         else:
             _cell_to = np.array([cell_to] if isinstance(cell_to, int) else cell_to)
 

virtual_ecosystem/models/hydrology/above_ground.py

@@ -290,7 +290,7 @@ def calculate_soil_evaporation(
         soil evaporation, [mm per time interval], aerodynamic resistance surface [s m-1]
     """
 
-    output = {}
+    output: dict[str, NDArray[np.floating]] = {}
 
     # Available soil moisture
     soil_moisture_free = np.clip(
@@ -307,7 +307,7 @@ def calculate_soil_evaporation(
     saturation_vapour_pressure = calc_vp_sat(
         ta=temperature,
         core_const=pyrealm_const(),
-    )
+    ).astype(np.floating)  # pyrealm returning np.float64
 
     saturated_specific_humidity = (
         gas_constant_water_vapour * saturation_vapour_pressure

virtual_ecosystem/models/hydrology/hydrology_model.py

@@ -15,7 +15,6 @@ class as a child of the :class:`~virtual_ecosystem.core.base_model.BaseModel` cl
 .. TODO:: time step and model structure
 
     * find a way to load daily (precipitation) data and loop over daily time_index
-    * allow for different time steps (currently only 30 days)
     * potentially move `calculate_drainage_map` to core
     * add abiotic constants from config
 
@@ -32,7 +31,6 @@ class as a child of the :class:`~virtual_ecosystem.core.base_model.BaseModel` cl
 
 import numpy as np
 from numpy.typing import NDArray
-from pint import Quantity
 from pyrealm.constants import CoreConst as PyrealmConst
 from xarray import DataArray
 
@@ -419,13 +417,16 @@ def _update(self, time_index: int, **kwargs: Any) -> None:
         and a number of parameters that as described in detail in
         :class:`~virtual_ecosystem.models.hydrology.constants.HydroConsts`.
         """
-        # Determine number of days, currently only 30 days (=1 month)
-        if self.model_timing.update_interval_quantity != Quantity("1 month"):
-            to_raise = NotImplementedError("This time step is currently not supported.")
-            LOGGER.error(to_raise)
-            raise to_raise
-
-        days: int = 30
+        # Determine number of days
+        days_float: float = self.model_timing.update_interval_seconds / 86400
+        days: int = int(days_float // 1)
+
+        # Check if the number of days is exact and warn if not
+        if not np.allclose(days_float % 1, 0):
+            LOGGER.warning(
+                f"Update interval is not a whole number of days ({days_float}),"
+                f" partitioning inputs among {days} days."
+            )
 
         # Set seed for random rainfall generator
         seed: None | int = kwargs.pop("seed", None)

virtual_ecosystem/models/litter/chemistry.py

@@ -306,7 +306,7 @@ def calculate_N_mineralisation(
         self,
         decay_rates: dict[str, NDArray[np.floating]],
         active_microbe_depth: float,
-    ) -> dict[str, NDArray[np.floating]]:
+    ) -> NDArray[np.floating]:
         """Method to calculate the amount of nitrogen mineralised by litter decay.
 
         This function finds the nitrogen mineralisation rate of each litter pool, by
@@ -358,7 +358,7 @@ def calculate_P_mineralisation(
         self,
         decay_rates: dict[str, NDArray[np.floating]],
         active_microbe_depth: float,
-    ) -> dict[str, NDArray[np.floating]]:
+    ) -> NDArray[np.floating]:
         """Method to calculate the amount of phosphorus mineralised by litter decay.
 
         This function finds the phosphorus mineralisation rate of each litter pool, by