Examples

Here is a quick description of the examples found in the demo directory. To understand the variables and classes in the code blocks, check out Quickstart Guide.

Linked Polders

This example demonstrates the aggregation of 10 identical polders connected in series:

Water flows in at 0.5 CMS into one end of the polders. The polders can only pump into their immediate neighbors. At the other end from the inflow, there is a polder that can pump water out of the network, but it can only do so at a rate of 0.45 CMS. All the polders are at the ideal max water volume, so none want to store the water. The polders must decide how to allocate the inflow. There are 10 of these polders. Their Modelica models look like this:

model Polder
  // Model Elements
  Deltares.ChannelFlow.SimpleRouting.BoundaryConditions.Inflow Inflow;
  Deltares.ChannelFlow.SimpleRouting.BoundaryConditions.Inflow Outflow;
  Deltares.ChannelFlow.SimpleRouting.Storage.Storage Storage(V(min = 8000.0, max=42000.0));
  // Inputs
  input Real Inflow_Q(fixed = false, min=-5.0, max=5.0) = Inflow.Q;
  input Real Outflow_Q(fixed = false, min=-5.0, max=5.0) = Outflow.Q;
  // Outputs
  output Real Storage_V = Storage.V;
  output Real incurred_damage = (30000.0 - Storage.V) / 12000.0;
equation
  connect(Inflow.QOut, Storage.QIn);
  connect(Storage.QOut, Outflow.QOut);
end Polder;

And their python class looks like this:

class Polder(
    MetaMixin,
    GoalProgrammingMixin,
    CSVMixin,
    ModelicaMixin,
    CollocatedIntegratedOptimizationProblem,
):
    goal_priority_dict = {"water_level_goal_0": 1}

    def path_goals(self):
        g = super().path_goals()

        # Add goal on water volume states
        g.append(
            StateGoalPlus(
                self,
                state="Storage.V",
                priority=self.goal_priority_dict["water_level_goal_0"],
                ts_min="Storage.V_min",
                ts_max="Storage.V_max",
            )
        )
        return g

Here is the aggregator class in a separate file:

from pathlib import Path

import casadi as ca

from rtc_tools_meta.meta_optimization_problem import MetaOptimizationProblem
from rtc_tools_meta.util import collect_model_instance

base_path = Path(__file__).parent.resolve()
all_polders = {
    f"polder{i}": collect_model_instance("Polder", base_path / f"polder{i}/src/opt.py")
    for i in range(1, 11)
}


class LinkedPolders(MetaOptimizationProblem):
    submodels = all_polders
    variable_mappings = [
        (f"polder{i}::Outflow_Q", f"polder{i+1}::Inflow_Q", "-") for i in range(1, 10)
    ]

    def meta_objective(self):
        return ca.sumsqr(
            ca.vertcat(
                *[
                    self.states_in(f"{submodel}::incurred_damage")
                    for submodel in self.submodels.keys()
                ]
            )
        )


if __name__ == "__main__":
    problem = LinkedPolders()
    problem.optimize()

Being egalitarian people, the operators decide that the high water levels must be shared as much as possible between the models. They give each model a damage function correlating the cost of high water with the water level, and penalize this expression with a second-order penalty. This means that when all the objectives from each polder are added together, IPOPT will try to reduce the largest damage function first. The net effect is that when IPOPT is done, the models will be sharing the inflows (and the damage) as evenly as possible, and as little as possible.

In the results plotted below, it is immediately apparent that the polders are taking equal damage. They do this by each storing an equal share of the inflow, while passing the rest on.

(Source code, svg, png)

Polders and Boezem

In this example we add a meta model on top of the previous example. The last polder drains into a boezem with no tolerance for flooding:

To make sure the flow into the boezem is sufficiently small, we aggregate the boezem model into the system. Here is the boezem as a Modelica model:

model Boezem
  // Model Elements
  Deltares.ChannelFlow.SimpleRouting.BoundaryConditions.Inflow Inflow;
  Deltares.ChannelFlow.SimpleRouting.BoundaryConditions.Inflow Polder_Inflow;
  Deltares.ChannelFlow.SimpleRouting.BoundaryConditions.Inflow Outflow;
  Deltares.ChannelFlow.SimpleRouting.Nodes.Node Node(nin=2, nout=1);
  Deltares.ChannelFlow.SimpleRouting.Branches.Integrator Channel(V(min = 400.0, max=600.0));
  // Inputs
  input Real Inflow_Q(fixed = true) = Inflow.Q;
  input Real Polder_Inflow_Q(fixed = false, min=-5.0, max=5.0) = Polder_Inflow.Q;
  input Real Outflow_Q(fixed = false, min=0.0, max=0.45) = Outflow.Q;
  // Outputs
  output Real Channel_V = Channel.V;
equation
  connect(Polder_Inflow.QOut, Node.QIn[1]);
  connect(Inflow.QOut, Node.QIn[2]);
  connect(Node.QOut[1], Channel.QIn);
  connect(Channel.QOut, Outflow.QOut);
end Boezem;

And its python class looks like this:

class Boezem(
    MetaMixin,
    GoalProgrammingMixin,
    CSVMixin,
    ModelicaMixin,
    CollocatedIntegratedOptimizationProblem,
):
    goal_priority_dict = {"water_level_goal": 1}

    def path_goals(self):
        g = super().path_goals()

        # Add RangeGoal on water volume states with a priority of 1
        g.append(
            StateGoalPlus(
                self,
                state="Channel.V",
                priority=self.goal_priority_dict["water_level_goal"],
                ts_min="Channel.V_min",
                ts_max="Channel.V_max",
            )
        )
        return g

The damage function tells IPOPT to find a solution that maximizes the amount of water coming from the LinkedPolders model, subject to the hard constraints volume constraints in the model.

Here is the meta model:

from pathlib import Path

import casadi as ca

from rtc_tools_meta.meta_optimization_problem import MetaOptimizationProblem
from rtc_tools_meta.util import collect_model_instance


base_path = Path(__file__).parent.resolve()

polders = collect_model_instance("LinkedPolders", base_path / "linked_polders/opt.py")
boezem = collect_model_instance("Boezem", base_path / "boezem/src/opt.py")


class PoldersBoezem(MetaOptimizationProblem):
    submodels = {"LinkedPolders": polders, "Boezem": boezem}
    variable_mappings = [
        ("LinkedPolders::polder10::Outflow_Q", "Boezem::Polder_Inflow_Q", "-")
    ]

    def meta_objective(self):
        return -ca.sum1(self.states_in("Boezem::Polder_Inflow_Q"))


if __name__ == "__main__":
    problem = PoldersBoezem()
    problem.optimize()

Note that to run, we simply call the optimize() method!

In the results plotted below, it is immediately apparent that the polders are taking equal damage. They do this by each storing an equal share of the inflow, while passing the rest on. Furthermore, The boezem does not take its share of damage- it only takes the flow that it can safely handle.

(Source code, svg, png)