Read some interpolation tables into the Model struct

core/src/Ribasim.jl

@@ -50,6 +50,7 @@ include("validation.jl")
 include("solve.jl")
 include("config.jl")
 using .config
+include("export.jl")
 include("utils.jl")
 include("lib.jl")
 include("io.jl")

core/src/bmi.jl

@@ -28,7 +28,7 @@ function BMI.initialize(T::Type{Model}, config::Config)::Model
     # All data from the GeoPackage that we need during runtime is copied into memory,
     # so we can directly close it again.
     db = SQLite.DB(gpkg_path)
-    local parameters, state, n, tstops
+    local parameters, state, n, tstops, level_exporters
     try
         parameters = Parameters(db, config)
 
@@ -82,6 +82,8 @@ function BMI.initialize(T::Type{Model}, config::Config)::Model
         # use state
         state = load_structvector(db, config, BasinStateV1)
         n = length(get_ids(db, "Basin"))
+
+        level_exporters = create_level_exporters(db, config, basin)
     finally
         # always close the GeoPackage, also in case of an error
         close(db)
@@ -150,7 +152,7 @@ function BMI.initialize(T::Type{Model}, config::Config)::Model
 
     set_initial_discrete_controlled_parameters!(integrator, storage)
 
-    return Model(integrator, config, saved_flow)
+    return Model(integrator, config, saved_flow, level_exporters)
 end
 
 """

core/src/export.jl

@@ -0,0 +1,58 @@
+# This module exports a water level:
+#
+# * the water level of the original hydrodynamic model before lumping.
+# * a differently aggregated water level, used for e.g. coupling to MODFLOW.
+#
+# The second is arguably easier to interpret.
+
+"""
+basin_level: a view on Ribasim's basin level.
+level: the interpolated water level
+tables: the interpolator callables
+
+All members of this struct have length n_elem.
+"""
+struct LevelExporter
+    basin_index::Vector{Int}
+    interpolations::Vector{ScalarInterpolation}
+    level::Vector{Float64}
+end
+
+function LevelExporter(tables, node_to_basin::Dict{Int, Int})::LevelExporter
+    basin_ids = Int[]
+    interpolations = ScalarInterpolation[]
+
+    for group in IterTools.groupby(row -> row.element_id, tables)
+        node_id = first(getproperty.(group, :node_id))
+        basin_level = getproperty.(group, :basin_level)
+        element_level = getproperty.(group, :level)
+        # Ensure it doesn't extrapolate before the first value.
+        new_interp = LinearInterpolation([element_level[1], element_level...], [prevfloat(basin_level[1]), basin_level...])
+        push!(basin_ids, node_to_basin[node_id])
+        push!(interpolations, new_interp)
+    end
+
+    return LevelExporter(basin_ids, interpolations, fill(NaN, length(basin_ids)))
+end
+
+function create_level_exporters(db::DB, config::Config, basin::Basin)::Dict{String, LevelExporter}
+    node_to_basin = Dict(node_id => index for (index, node_id) in enumerate(basin.node_id))
+    tables = load_structvector(db, config, LevelExporterStaticV1)
+    level_exporters = Dict{String, LevelExporter}()
+    if length(tables) > 0
+        for group in IterTools.groupby(row -> row.name, tables)
+            name = first(getproperty.(group, :name))
+            level_exporters[name] = LevelExporter(group, node_to_basin)
+        end
+    end
+    return level_exporters
+end
+
+"""
+Compute a new water level for each external element.
+"""
+function update!(exporter::LevelExporter, basin_level)
+    for (i, (index, interp)) in enumerate(zip(exporter.basin_index, exporter.interpolations))
+        exporter.level[i] = interp(basin_level[index])
+    end
+end

core/src/lib.jl

@@ -12,12 +12,14 @@ struct Model{T}
     integrator::T
     config::Config
     saved_flow::SavedValues{Float64, Vector{Float64}}
+    level_exporters::Dict{String, LevelExporter}
     function Model(
         integrator::T,
         config,
         saved_flow,
+        level_exporters,
     ) where {T <: SciMLBase.AbstractODEIntegrator}
-        new{T}(integrator, config, saved_flow)
+        new{T}(integrator, config, saved_flow, level_exporters)
     end
 end
 

core/src/validation.jl

@@ -24,6 +24,7 @@
 @schema "ribasim.outlet.static" OutletStatic
 @schema "ribasim.user.static" UserStatic
 @schema "ribasim.user.time" UserTime
+@schema "ribasim.levelexporter.static" LevelExporterStatic
 
 const delimiter = " / "
 tablename(sv::Type{SchemaVersion{T, N}}) where {T, N} = join(nodetype(sv), delimiter)
@@ -308,6 +309,14 @@ end
     priority::Int
 end
 
+@version LevelExporterStaticV1 begin
+    name::String
+    element_id::Int
+    node_id::Int
+    basin_level::Float64
+    level::Float64
+end
+
 function variable_names(s::Any)
     filter(x -> !(x in (:node_id, :control_state)), fieldnames(s))
 end

docs/schema/Config.schema.json

@@ -75,6 +75,12 @@
                 "timing": false
             }
         },
+        "fractional_flow": {
+            "$ref": "https://deltares.github.io/Ribasim/schema/fractional_flow.schema.json",
+            "default": {
+                "static": null
+            }
+        },
         "terminal": {
             "$ref": "https://deltares.github.io/Ribasim/schema/terminal.schema.json",
             "default": {
@@ -156,8 +162,8 @@
                 "static": null
             }
         },
-        "fractional_flow": {
-            "$ref": "https://deltares.github.io/Ribasim/schema/fractional_flow.schema.json",
+        "level_exporter": {
+            "$ref": "https://deltares.github.io/Ribasim/schema/level_exporter.schema.json",
             "default": {
                 "static": null
             }
@@ -174,6 +180,7 @@
         "output",
         "solver",
         "logging",
+        "fractional_flow",
         "terminal",
         "pid_control",
         "level_boundary",
@@ -186,6 +193,6 @@
         "discrete_control",
         "outlet",
         "linear_resistance",
-        "fractional_flow"
+        "level_exporter"
     ]
 }

docs/schema/LevelExporterStatic.schema.json

@@ -0,0 +1,42 @@
+{
+    "$schema": "https://json-schema.org/draft/2020-12/schema",
+    "$id": "https://deltares.github.io/Ribasim/schema/LevelExporterStatic.schema.json",
+    "title": "LevelExporterStatic",
+    "description": "A LevelExporterStatic object based on Ribasim.LevelExporterStaticV1",
+    "type": "object",
+    "properties": {
+        "name": {
+            "format": "default",
+            "type": "string"
+        },
+        "element_id": {
+            "format": "default",
+            "type": "integer"
+        },
+        "node_id": {
+            "format": "default",
+            "type": "integer"
+        },
+        "basin_level": {
+            "format": "double",
+            "type": "number"
+        },
+        "level": {
+            "format": "double",
+            "type": "number"
+        },
+        "remarks": {
+            "description": "a hack for pandera",
+            "type": "string",
+            "format": "default",
+            "default": ""
+        }
+    },
+    "required": [
+        "name",
+        "element_id",
+        "node_id",
+        "basin_level",
+        "level"
+    ]
+}

docs/schema/level_exporter.schema.json

@@ -0,0 +1,23 @@
+{
+    "$schema": "https://json-schema.org/draft/2020-12/schema",
+    "$id": "https://deltares.github.io/Ribasim/schema/level_exporter.schema.json",
+    "title": "level_exporter",
+    "description": "A level_exporter object based on Ribasim.config.level_exporter",
+    "type": "object",
+    "properties": {
+        "static": {
+            "format": "default",
+            "anyOf": [
+                {
+                    "type": "null"
+                },
+                {
+                    "type": "string"
+                }
+            ],
+            "default": null
+        }
+    },
+    "required": [
+    ]
+}

docs/schema/root.schema.json

@@ -37,6 +37,9 @@
         "LevelBoundaryTime": {
             "$ref": "LevelBoundaryTime.schema.json"
         },
+        "LevelExporterStatic": {
+            "$ref": "LevelExporterStatic.schema.json"
+        },
         "LinearResistanceStatic": {
             "$ref": "LinearResistanceStatic.schema.json"
         },

python/ribasim/ribasim/__init__.py

@@ -5,6 +5,7 @@
 from ribasim.config import Config, Logging, Solver
 from ribasim.geometry.edge import Edge
 from ribasim.geometry.node import Node
+from ribasim.level_exporter import LevelExporter
 from ribasim.model import Model
 from ribasim.node_types.basin import Basin
 from ribasim.node_types.discrete_control import DiscreteControl
@@ -42,4 +43,5 @@
     "DiscreteControl",
     "PidControl",
     "User",
+    "LevelExporter",
 ]

python/ribasim/ribasim/config.py

@@ -38,6 +38,10 @@ class Logging(BaseModel):
     timing: bool = False
 
 
+class FractionalFlow(BaseModel):
+    static: Optional[str] = None
+
+
 class Terminal(BaseModel):
     static: Optional[str] = None
 
@@ -95,7 +99,7 @@ class LinearResistance(BaseModel):
     static: Optional[str] = None
 
 
-class FractionalFlow(BaseModel):
+class LevelExporter(BaseModel):
     static: Optional[str] = None
 
 
@@ -142,6 +146,9 @@ class Config(BaseModel):
             {"verbosity": {"level": 0}, "timing": False}
         )
     )
+    fractional_flow: FractionalFlow = Field(
+        default_factory=lambda: FractionalFlow.parse_obj({"static": None})
+    )
     terminal: Terminal = Field(
         default_factory=lambda: Terminal.parse_obj({"static": None})
     )
@@ -180,6 +187,6 @@ class Config(BaseModel):
     linear_resistance: LinearResistance = Field(
         default_factory=lambda: LinearResistance.parse_obj({"static": None})
     )
-    fractional_flow: FractionalFlow = Field(
-        default_factory=lambda: FractionalFlow.parse_obj({"static": None})
+    level_exporter: LevelExporter = Field(
+        default_factory=lambda: LevelExporter.parse_obj({"static": None})
     )

python/ribasim/ribasim/level_exporter.py

@@ -0,0 +1,17 @@
+from pandera.typing import DataFrame
+
+from ribasim.input_base import TableModel
+from ribasim.schemas import LevelExporterStaticSchema  # type: ignore
+
+
+class LevelExporter(TableModel):
+    """The level exporter export Ribasim water levels."""
+
+    static: DataFrame[LevelExporterStaticSchema]
+
+    def sort(self):
+        self.static.sort_values(
+            ["name", "element_id", "node_id", "basin_level"],
+            ignore_index=True,
+            inplace=True,
+        )

python/ribasim/ribasim/model.py

@@ -21,6 +21,7 @@
 # Do not import from ribasim namespace: will create import errors.
 # E.g. not: from ribasim import Basin
 from ribasim.input_base import TableModel
+from ribasim.level_exporter import LevelExporter
 from ribasim.node_types.basin import Basin
 from ribasim.node_types.discrete_control import DiscreteControl
 from ribasim.node_types.flow_boundary import FlowBoundary
@@ -104,6 +105,7 @@ class Model(BaseModel):
     discrete_control: Optional[DiscreteControl]
     pid_control: Optional[PidControl]
     user: Optional[User]
+    level_exporter: Optional[LevelExporter]
     starttime: datetime.datetime
     endtime: datetime.datetime
     solver: Optional[Solver]

python/ribasim/ribasim/models.py

@@ -103,6 +103,15 @@ class LevelBoundaryTime(BaseModel):
     remarks: str = Field("", description="a hack for pandera")
 
 
+class LevelExporterStatic(BaseModel):
+    name: str
+    element_id: int
+    node_id: int
+    basin_level: float
+    level: float
+    remarks: str = Field("", description="a hack for pandera")
+
+
 class LinearResistanceStatic(BaseModel):
     node_id: int
     active: Optional[bool] = None
@@ -229,6 +238,7 @@ class Root(BaseModel):
     FractionalFlowStatic: Optional[FractionalFlowStatic] = None
     LevelBoundaryStatic: Optional[LevelBoundaryStatic] = None
     LevelBoundaryTime: Optional[LevelBoundaryTime] = None
+    LevelExporterStatic: Optional[LevelExporterStatic] = None
     LinearResistanceStatic: Optional[LinearResistanceStatic] = None
     ManningResistanceStatic: Optional[ManningResistanceStatic] = None
     Node: Optional[Node] = None

python/ribasim_testmodels/ribasim_testmodels/trivial.py

@@ -96,13 +96,32 @@ def trivial_model() -> ribasim.Model:
         )
     )
 
+    # Create a level exporter from one basin to three elements. Scale one to one, but:
+    #
+    # 1. start at -1.0
+    # 2. start at 0.0
+    # 3. start at 1.0
+    #
+    level_exporter = ribasim.LevelExporter(
+        static=pd.DataFrame(
+            data={
+                "name": "primary-system",
+                "element_id": [1, 1, 2, 2, 3, 3],
+                "node_id": [1, 1, 1, 1, 1, 1],
+                "basin_level": [0.0, 1.0, 0.0, 1.0, 0.0, 1.0],
+                "level": [-1.0, 0.0, 0.0, 1.0, 1.0, 2.0],
+            }
+        )
+    )
+
     model = ribasim.Model(
         modelname="trivial",
         node=node,
         edge=edge,
         basin=basin,
         terminal=terminal,
         tabulated_rating_curve=rating_curve,
+        level_exporter=level_exporter,
         starttime="2020-01-01 00:00:00",
         endtime="2021-01-01 00:00:00",
     )