feat!: Add a variable API and optimisation problems

pywr-book/py-listings/optimisation/model.json

@@ -0,0 +1,112 @@
+{
+  "metadata": {
+    "title": "Simple problem",
+    "description": "A simple model with a variable demand, used to test optimisation problems."
+  },
+  "time": {
+    "start": "2021-01-01",
+    "end": "2021-01-31",
+    "timestep": {
+      "type": "Days",
+      "days": 1
+    }
+  },
+  "network": {
+    "nodes": [
+      {
+        "meta": {
+          "name": "input1"
+        },
+        "type": "Input",
+        "max_flow": {
+          "type": "Literal",
+          "value": 8.0
+        }
+      },
+      {
+        "meta": {
+          "name": "link1"
+        },
+        "type": "Link"
+      },
+      {
+        "meta": {
+          "name": "output1"
+        },
+        "type": "Output",
+        "cost": {
+          "type": "Literal",
+          "value": -10.0
+        },
+        "max_flow": {
+          "type": "Parameter",
+          "name": "demand"
+        }
+      }
+    ],
+    "edges": [
+      {
+        "from_node": "input1",
+        "to_node": "link1"
+      },
+      {
+        "from_node": "link1",
+        "to_node": "output1"
+      }
+    ],
+    "parameters": [
+      {
+        "meta": {
+          "name": "demand"
+        },
+        "type": "Constant",
+        "value": {
+          "type": "Literal",
+          "value": 10.0
+        }
+      }
+    ],
+    "metric_sets": [
+      {
+        "meta": {"name": "demand"},
+        "metrics": [
+          {
+            "type": "Parameter",
+            "name": "demand"
+          }
+        ]
+      },
+      {
+        "meta": {"name": "delivered"},
+        "metrics": [
+          {
+            "type": "Node",
+            "name": "output1"
+          }
+        ]
+      }
+    ],
+    "outputs": [
+      {
+        "meta": {"name": "mean-demand"},
+        "type": "Memory",
+        "metric_set": "demand",
+        "aggregation": {
+          "time": {
+            "type": "Mean"
+          }
+        }
+      },
+      {
+        "meta": {"name": "mean-delivered"},
+        "type": "Memory",
+        "metric_set": "delivered",
+        "aggregation": {
+          "time": {
+            "type": "Mean"
+          }
+        }
+      }
+    ]
+  }
+}

pywr-book/py-listings/optimisation/problem.json

@@ -0,0 +1,31 @@
+{
+  "variables": [
+    {
+      "name": "demand",
+      "config": {
+        "type": "ActivationFunction",
+        "activation": {
+          "type": "Unit",
+          "min": 0.0,
+          "max": 20.0
+        }
+      }
+    }
+  ],
+  "objectives": [
+    {
+      "output": "mean-demand",
+      "sense": "Minimise"
+    },
+    {
+      "output": "mean-delivered",
+      "sense": "Maximise"
+    }
+  ],
+  "constraints": [
+    {
+      "output": "mean-delivered",
+      "lower_bound": 5.0
+    }
+  ]
+}

pywr-book/py-listings/optimisation/random_search.py

@@ -0,0 +1,87 @@
+"""Search for good values of a model's variable parameters using random sampling.
+
+This example uses a simple random search to show how a `ModelProblem` can be used by an
+optimisation algorithm. In practice, a library of optimisation algorithms (e.g. multi-objective
+evolutionary algorithms) would be used to choose the values to evaluate.
+"""
+
+import random
+from pathlib import Path
+
+from pywr import Evaluation, ModelProblem, ModelSchema, ProblemConfig
+
+
+# ANCHOR: build
+def build_problem() -> ModelProblem:
+    here = Path(__file__).parent
+
+    schema = ModelSchema.from_path(here / "model.json")
+    config = ProblemConfig.from_path(here / "problem.json")
+
+    return schema.build_problem(config, data_path=here)
+
+
+# ANCHOR_END: build
+
+
+def dominates(a: list[float], b: list[float]) -> bool:
+    """Whether objective values `a` dominate `b`, when all objectives are minimised."""
+    return all(x <= y for x, y in zip(a, b, strict=True)) and any(
+        x < y for x, y in zip(a, b, strict=True)
+    )
+
+
+# ANCHOR: search
+def random_search(
+    problem: ModelProblem, num_samples: int, seed: int = 42
+) -> list[tuple[list[float], Evaluation]]:
+    """Evaluate random values, and return the feasible values that are not dominated."""
+    rng = random.Random(seed)
+    lower, _ = problem.lower_bounds()
+    upper, _ = problem.upper_bounds()
+
+    samples = [
+        [rng.uniform(lb, ub) for lb, ub in zip(lower, upper, strict=True)]
+        for _ in range(num_samples)
+    ]
+
+    # Each sample is evaluated by running the model; the samples are run in parallel.
+    evaluations = problem.evaluate_batch(samples)
+
+    # Negate the objectives that are maximised, so that all objectives are minimised.
+    signs = [1.0 if sense == "Minimise" else -1.0 for sense in problem.objective_senses]
+
+    def minimised(evaluation: Evaluation) -> list[float]:
+        return [s * o for s, o in zip(signs, evaluation.objectives, strict=True)]
+
+    feasible = [
+        (sample, evaluation)
+        for sample, evaluation in zip(samples, evaluations, strict=True)
+        if evaluation.is_feasible
+    ]
+
+    return [
+        (sample, evaluation)
+        for sample, evaluation in feasible
+        if not any(
+            dominates(minimised(other), minimised(evaluation)) for _, other in feasible
+        )
+    ]
+
+
+# ANCHOR_END: search
+
+
+if __name__ == "__main__":
+    problem = build_problem()
+    print(problem)
+
+    results = random_search(problem, num_samples=50)
+    for sample, evaluation in sorted(results, key=lambda result: result[0]):
+        objectives = ", ".join(
+            f"{name}={value:.2f}"
+            for name, value in zip(
+                problem.objective_names, evaluation.objectives, strict=True
+            )
+        )
+        print(f"demand={sample[0]:.2f}: {objectives}")

pywr-book/src/SUMMARY.md

@@ -18,6 +18,7 @@
 - [Custom parameters](./custom_parameters.md)
 - [River routing and attenuation](./river-routing.md)
 - [Multi-models](./multi_models.md)
+- [Optimisation](./optimisation.md)
 - [Migration guide](./migration_guide.md)
 - [Model management](./model-management/README.md)
     - [Version control](./model-management/version-control.md)

pywr-book/src/optimisation.md

@@ -0,0 +1,89 @@
+# Optimisation
+
+Pywr models can be used with optimisation algorithms, such as multi-objective evolutionary
+algorithms, to search for good values of some of a model's parameters. The algorithm repeatedly
+chooses values for these _variable_ parameters, runs the model with them, and uses the results to
+calculate the _objectives_ (values to minimise or maximise) and _constraints_ (values that must be
+within bounds).
+
+Pywr does not include any optimisation algorithms. Instead, it provides a way to define a problem
+for a model, and to evaluate the problem for given variable values, which can be used with an
+optimisation library of your choice.
+
+## Defining a problem
+
+A problem is defined in a JSON file that is separate from the model, so that a model can be used
+with different problems. The file has three sections:
+
+- `variables`: the parameters whose values can be varied, and how. The values of the variables
+  are given as a list, and are split between the parameters in the order they are listed here.
+- `objectives`: the outputs whose aggregated values are minimised or maximised.
+- `constraints` (optional): the outputs whose aggregated values must be within a lower and/or
+  upper bound.
+
+Objectives and constraints refer to outputs by name. The outputs must be able to aggregate their
+results to a single value, such as a `Memory` output with a time and scenario aggregation.
+
+The example below defines a problem for a model with a single demand. The demand is variable
+between 0.0 and 20.0. The objectives are to minimise the mean demand and to maximise the mean
+delivered flow, and a constraint requires the mean delivered flow to be at least 5.0.
+
+```json
+{{#include ../py-listings/optimisation/problem.json}}
+```
+
+### Variable parameters
+
+Each variable parameter has a `config` whose `type` depends on the type of the parameter:
+
+| Config type          | Parameter types        | Variable values                                                           |
+|----------------------|------------------------|---------------------------------------------------------------------------|
+| `ActivationFunction` | `Constant`, `Offset`   | One value, transformed by the activation function.                        |
+| `MonthlyProfile`     | `MonthlyProfile`       | Twelve values, one for each month, within the given bounds.               |
+| `RbfProfile`         | `RbfProfile`           | One value per point, and optionally one day of the year per point.        |
+
+Activation functions transform the variable value before it is used by the parameter. For
+example, the `Unit` function uses the value directly within its `min` and `max`, while the
+`BinaryStep` function returns one of two values depending on whether the variable value is
+positive. This allows, for example, a binary-like decision to be represented by a continuous
+variable.
+
+Most variable values are floats. The day of the year of an `RbfProfile`'s points are integers,
+and are given separately from the float values.
+
+## Evaluating a problem in Python
+
+A problem is built from a model schema and a `ProblemConfig`:
+
+[//]: # (@formatter:off)
+
+```python,ignore
+{{ #include ../py-listings/optimisation/random_search.py:build}}
+```
+
+[//]: # (@formatter:on)
+
+The `ModelProblem` gives the number of variable values (`size`), their bounds (`lower_bounds()`
+and `upper_bounds()`), and the names and senses of the objectives and constraints. Its
+`evaluate` method runs the model with the given variable values, and returns an `Evaluation` of
+the objective and constraint values. `evaluate_batch` evaluates several sets of values in
+parallel, which is useful for population-based algorithms.
+
+The objective values are returned as calculated, so objectives that are maximised are not
+negated. `Evaluation.constraint_violations` gives the amount by which each constraint is outside
+its bounds, which is zero when the constraint is satisfied.
+
+The listing below uses a simple random search to find the values that give the best trade-off
+between the two objectives.
+
+[//]: # (@formatter:off)
+
+```python,ignore
+{{ #include ../py-listings/optimisation/random_search.py:search}}
+```
+
+[//]: # (@formatter:on)
+
+> **Note:** Each evaluation sets up and runs the model from scratch, including all of its outputs.
+> Outputs that write files, such as CSV or HDF5 outputs, will rewrite their files on every
+> evaluation. Models used for optimisation should normally only use `Memory` outputs.

pywr-core/src/agg_funcs/mod.rs

@@ -185,12 +185,13 @@ impl AggFuncU64 {
 
 #[cfg(test)]
 mod tests {
-    use crate::agg_funcs::{AggFuncF64, AggFuncU64, py};
+    #[cfg(feature = "pyo3")]
+    use crate::agg_funcs::py;
+    use crate::agg_funcs::{AggFuncF64, AggFuncU64};
     use float_cmp::assert_approx_eq;
-    use pyo3::{Py, PyAny};
     #[cfg(feature = "pyo3")]
     use pyo3::{
-        Python,
+        Py, PyAny, Python,
         ffi::c_str,
         prelude::{PyAnyMethods, PyModule},
         types::{PyDict, PyTuple},

pywr-core/src/lib.rs

@@ -7,8 +7,10 @@ pub mod edge;
 pub mod metric;
 pub mod models;
 pub mod network;
+pub mod network_variable_config;
 pub mod node;
 pub mod parameters;
+pub mod problem;
 pub mod recorders;
 pub mod scenario;
 pub mod solvers;