From 4812bcac34862f12bcabf3865b32acd99b323bcd Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Wed, 10 Apr 2024 13:42:18 +0200 Subject: [PATCH 01/22] Python side --- core/src/schema.jl | 10 ++++ python/ribasim/ribasim/config.py | 5 ++ python/ribasim/ribasim/model.py | 12 +++- .../ribasim/ribasim/nodes/discrete_control.py | 13 +++- python/ribasim/ribasim/schemas.py | 9 +++ .../ribasim_testmodels/__init__.py | 2 + .../ribasim_testmodels/discrete_control.py | 60 +++++++++++++++++++ 7 files changed, 106 insertions(+), 5 deletions(-) diff --git a/core/src/schema.jl b/core/src/schema.jl index d59190c3b..c8dd7094d 100644 --- a/core/src/schema.jl +++ b/core/src/schema.jl @@ -2,6 +2,7 @@ # The identifier is parsed as ribasim.nodetype.kind, no capitals or underscores are allowed. @schema "ribasim.discretecontrol.condition" DiscreteControlCondition @schema "ribasim.discretecontrol.logic" DiscreteControlLogic +@schema "ribasim.discretecontrol.compoundvariable" DiscreteControlCompoundVariable @schema "ribasim.basin.static" BasinStatic @schema "ribasim.basin.time" BasinTime @schema "ribasim.basin.profile" BasinProfile @@ -183,6 +184,15 @@ end node_id::Int32 end +@version DiscreteControlCompoundVariableV1 begin + node_id::Int32 + name::String + listen_node_type::String + listen_node_id::Int + variable::String + weight::Float64 +end + @version DiscreteControlConditionV1 begin node_id::Int32 listen_node_type::String diff --git a/python/ribasim/ribasim/config.py b/python/ribasim/ribasim/config.py index f9605111d..ecc49badc 100644 --- a/python/ribasim/ribasim/config.py +++ b/python/ribasim/ribasim/config.py @@ -21,6 +21,7 @@ BasinStaticSchema, BasinSubgridSchema, BasinTimeSchema, + DiscreteControlCompoundVariableSchema, DiscreteControlConditionSchema, DiscreteControlLogicSchema, FlowBoundaryStaticSchema, @@ -287,6 +288,10 @@ class ManningResistance(MultiNodeModel): class DiscreteControl(MultiNodeModel): + compound_variable: TableModel[DiscreteControlCompoundVariableSchema] = Field( + default_factory=TableModel[DiscreteControlCompoundVariableSchema], + json_schema_extra={"sort_keys": ["name"]}, + ) condition: TableModel[DiscreteControlConditionSchema] = Field( default_factory=TableModel[DiscreteControlConditionSchema], json_schema_extra={ diff --git a/python/ribasim/ribasim/model.py b/python/ribasim/ribasim/model.py index 34c808f73..472f7d4cc 100644 --- a/python/ribasim/ribasim/model.py +++ b/python/ribasim/ribasim/model.py @@ -321,11 +321,17 @@ def plot_control_listen(self, ax): df_listen_edge = pd.concat([df_listen_edge, to_add]) # Listen edges from DiscreteControl - condition = self.discrete_control.condition.df - if condition is not None: - to_add = condition[ + for table in ( + self.discrete_control.condition.df, + self.discrete_control.compound_variable.df, + ): + if table is None: + continue + + to_add = table[ ["node_id", "listen_node_id", "listen_node_type"] ].drop_duplicates() + to_add = to_add[to_add["listen_node_type"] != "compound"] to_add.columns = ["control_node_id", "listen_node_id", "listen_node_type"] to_add["control_node_type"] = "DiscreteControl" df_listen_edge = pd.concat([df_listen_edge, to_add]) diff --git a/python/ribasim/ribasim/nodes/discrete_control.py b/python/ribasim/ribasim/nodes/discrete_control.py index eeb0b28ca..02b22840c 100644 --- a/python/ribasim/ribasim/nodes/discrete_control.py +++ b/python/ribasim/ribasim/nodes/discrete_control.py @@ -1,9 +1,18 @@ from pandas import DataFrame from ribasim.input_base import TableModel -from ribasim.schemas import DiscreteControlConditionSchema, DiscreteControlLogicSchema +from ribasim.schemas import ( + DiscreteControlCompoundVariableSchema, + DiscreteControlConditionSchema, + DiscreteControlLogicSchema, +) -__all__ = ["Condition", "Logic"] +__all__ = ["Condition", "Logic", "CompoundVariable"] + + +class CompoundVariable(TableModel[DiscreteControlCompoundVariableSchema]): + def __init__(self, **kwargs): + super().__init__(df=DataFrame(dict(**kwargs))) class Condition(TableModel[DiscreteControlConditionSchema]): diff --git a/python/ribasim/ribasim/schemas.py b/python/ribasim/ribasim/schemas.py index 0cf033671..aa4c738c2 100644 --- a/python/ribasim/ribasim/schemas.py +++ b/python/ribasim/ribasim/schemas.py @@ -48,6 +48,15 @@ class BasinTimeSchema(_BaseSchema): urban_runoff: Series[float] = pa.Field(nullable=True) +class DiscreteControlCompoundVariableSchema(_BaseSchema): + node_id: Series[Int32] = pa.Field(nullable=False, default=0) + name: Series[str] = pa.Field(nullable=False) + listen_node_type: Series[str] = pa.Field(nullable=False) + listen_node_id: Series[Int32] = pa.Field(nullable=False, default=0) + variable: Series[str] = pa.Field(nullable=False) + weight: Series[float] = pa.Field(nullable=False) + + class DiscreteControlConditionSchema(_BaseSchema): node_id: Series[Int32] = pa.Field(nullable=False, default=0) listen_node_type: Series[str] = pa.Field(nullable=False) diff --git a/python/ribasim_testmodels/ribasim_testmodels/__init__.py b/python/ribasim_testmodels/ribasim_testmodels/__init__.py index ceeb1fad2..4f9ce86c3 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/__init__.py +++ b/python/ribasim_testmodels/ribasim_testmodels/__init__.py @@ -28,6 +28,7 @@ ) from ribasim_testmodels.bucket import bucket_model, leaky_bucket_model from ribasim_testmodels.discrete_control import ( + compound_variable_condition_model, flow_condition_model, level_boundary_condition_model, level_setpoint_with_minmax_model, @@ -64,6 +65,7 @@ "basic_model", "basic_transient_model", "bucket_model", + "compound_variable_condition_model", "discrete_control_of_pid_control_model", "dutch_waterways_model", "flow_boundary_time_model", diff --git a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py index aae7674c9..2f14e8e9a 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py +++ b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py @@ -390,3 +390,63 @@ def level_setpoint_with_minmax_model() -> Model: ) return model + + +def compound_variable_condition_model() -> Model: + """ + Set up a minimal model containing a condition on a compound variable + for discrete control. + """ + + model = Model( + starttime="2020-01-01", + endtime="2021-01-01", + crs="EPSG:28992", + ) + + model.basin.add( + Node(1, Point(1, 0)), + [ + basin.Profile(area=1000.0, level=[0.0, 1.0]), + basin.State(level=[1.0]), + ], + ) + model.flow_boundary.add( + Node(2, Point(0, 0)), [flow_boundary.Static(flow_rate=[0.0])] + ) + model.flow_boundary.add( + Node(3, Point(0, 1)), + [flow_boundary.Time(time=["2020-01-01", "2022-01-01"], flow_rate=[0.0, 2.0])], + ) + model.pump.add( + Node(4, Point(2, 0)), + [pump.Static(control_state=["Off", "On"], flow_rate=[0.0, 1.0])], + ) + model.terminal.add(Node(5, Point(3, 0))) + model.discrete_control.add( + Node(6, Point(1, 1)), + [ + discrete_control.CompoundVariable( + name="flow_mean", + listen_node_type="FlowBoundary", + listen_node_id=[2, 3], + variable="flow", + weight=0.5, + ), + discrete_control.Condition( + listen_node_type="compound", + listen_node_id=[0], # Irrelevant + variable="flow_mean", + greater_than=1.0, + ), + discrete_control.Logic(truth_state=["T", "F"], control_state=["On", "Off"]), + ], + ) + + model.edge.add(model.flow_boundary[2], model.basin[1]) + model.edge.add(model.flow_boundary[3], model.basin[1]) + model.edge.add(model.basin[1], model.pump[4]) + model.edge.add(model.pump[4], model.terminal[5]) + model.edge.add(model.discrete_control[6], model.pump[4]) + + return model From 44fcaaea381fbe62f95758447748094b2502085c Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Wed, 10 Apr 2024 16:05:31 +0200 Subject: [PATCH 02/22] Change name of new schema --- core/src/schema.jl | 4 ++-- python/ribasim/ribasim/config.py | 6 +++--- .../ribasim/ribasim/nodes/discrete_control.py | 6 +++--- python/ribasim/ribasim/schemas.py | 2 +- .../ribasim_testmodels/discrete_control.py | 2 +- ribasim_qgis/core/nodes.py | 21 +++++++++++++++++++ 6 files changed, 31 insertions(+), 10 deletions(-) diff --git a/core/src/schema.jl b/core/src/schema.jl index c8dd7094d..2d3896e82 100644 --- a/core/src/schema.jl +++ b/core/src/schema.jl @@ -2,7 +2,7 @@ # The identifier is parsed as ribasim.nodetype.kind, no capitals or underscores are allowed. @schema "ribasim.discretecontrol.condition" DiscreteControlCondition @schema "ribasim.discretecontrol.logic" DiscreteControlLogic -@schema "ribasim.discretecontrol.compoundvariable" DiscreteControlCompoundVariable +@schema "ribasim.discretecontrol.compoundvariable" DiscreteControlCompoundvariable @schema "ribasim.basin.static" BasinStatic @schema "ribasim.basin.time" BasinTime @schema "ribasim.basin.profile" BasinProfile @@ -184,7 +184,7 @@ end node_id::Int32 end -@version DiscreteControlCompoundVariableV1 begin +@version DiscreteControlCompoundvariableV1 begin node_id::Int32 name::String listen_node_type::String diff --git a/python/ribasim/ribasim/config.py b/python/ribasim/ribasim/config.py index ecc49badc..64a5019e9 100644 --- a/python/ribasim/ribasim/config.py +++ b/python/ribasim/ribasim/config.py @@ -21,7 +21,7 @@ BasinStaticSchema, BasinSubgridSchema, BasinTimeSchema, - DiscreteControlCompoundVariableSchema, + DiscreteControlCompoundvariableSchema, DiscreteControlConditionSchema, DiscreteControlLogicSchema, FlowBoundaryStaticSchema, @@ -288,8 +288,8 @@ class ManningResistance(MultiNodeModel): class DiscreteControl(MultiNodeModel): - compound_variable: TableModel[DiscreteControlCompoundVariableSchema] = Field( - default_factory=TableModel[DiscreteControlCompoundVariableSchema], + compoundvariable: TableModel[DiscreteControlCompoundvariableSchema] = Field( + default_factory=TableModel[DiscreteControlCompoundvariableSchema], json_schema_extra={"sort_keys": ["name"]}, ) condition: TableModel[DiscreteControlConditionSchema] = Field( diff --git a/python/ribasim/ribasim/nodes/discrete_control.py b/python/ribasim/ribasim/nodes/discrete_control.py index 02b22840c..cc3cfcaf9 100644 --- a/python/ribasim/ribasim/nodes/discrete_control.py +++ b/python/ribasim/ribasim/nodes/discrete_control.py @@ -2,15 +2,15 @@ from ribasim.input_base import TableModel from ribasim.schemas import ( - DiscreteControlCompoundVariableSchema, + DiscreteControlCompoundvariableSchema, DiscreteControlConditionSchema, DiscreteControlLogicSchema, ) -__all__ = ["Condition", "Logic", "CompoundVariable"] +__all__ = ["Condition", "Logic", "Compoundvariable"] -class CompoundVariable(TableModel[DiscreteControlCompoundVariableSchema]): +class Compoundvariable(TableModel[DiscreteControlCompoundvariableSchema]): def __init__(self, **kwargs): super().__init__(df=DataFrame(dict(**kwargs))) diff --git a/python/ribasim/ribasim/schemas.py b/python/ribasim/ribasim/schemas.py index aa4c738c2..213d97667 100644 --- a/python/ribasim/ribasim/schemas.py +++ b/python/ribasim/ribasim/schemas.py @@ -48,7 +48,7 @@ class BasinTimeSchema(_BaseSchema): urban_runoff: Series[float] = pa.Field(nullable=True) -class DiscreteControlCompoundVariableSchema(_BaseSchema): +class DiscreteControlCompoundvariableSchema(_BaseSchema): node_id: Series[Int32] = pa.Field(nullable=False, default=0) name: Series[str] = pa.Field(nullable=False) listen_node_type: Series[str] = pa.Field(nullable=False) diff --git a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py index 2f14e8e9a..ea9eb3899 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py +++ b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py @@ -426,7 +426,7 @@ def compound_variable_condition_model() -> Model: model.discrete_control.add( Node(6, Point(1, 1)), [ - discrete_control.CompoundVariable( + discrete_control.Compoundvariable( name="flow_mean", listen_node_type="FlowBoundary", listen_node_id=[2, 3], diff --git a/ribasim_qgis/core/nodes.py b/ribasim_qgis/core/nodes.py index 95edb13bd..70750e419 100644 --- a/ribasim_qgis/core/nodes.py +++ b/ribasim_qgis/core/nodes.py @@ -691,6 +691,27 @@ def attributes(cls) -> list[QgsField]: ] +class DiscreteControlCompoundvariable(Input): + @classmethod + def input_type(cls) -> str: + return "DiscreteControl / compoundvariable" + + @classmethod + def geometry_type(cls) -> str: + return "No Geometry" + + @classmethod + def attributes(cls) -> list[QgsField]: + return [ + QgsField("node_id", QVariant.Int), + QgsField("name", QVariant.String), + QgsField("listen_node_type", QVariant.String), + QgsField("listen_node_id", QVariant.Int), + QgsField("variable", QVariant.String), + QgsField("weight", QVariant.Double), + ] + + class DiscreteControlCondition(Input): @classmethod def input_type(cls) -> str: From 1bc18e7bc62ff2dffd48e30aa7137df30fd77c7a Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Wed, 10 Apr 2024 16:49:39 +0200 Subject: [PATCH 03/22] Read compound variable data in core --- core/src/parameter.jl | 7 +++-- core/src/read.jl | 47 ++++++++++++++++++++++++++++--- core/src/schema.jl | 1 + python/ribasim/ribasim/schemas.py | 1 + 4 files changed, 49 insertions(+), 7 deletions(-) diff --git a/core/src/parameter.jl b/core/src/parameter.jl index 16774df21..2d8c5bca0 100644 --- a/core/src/parameter.jl +++ b/core/src/parameter.jl @@ -455,9 +455,10 @@ record: Namedtuple with discrete control information for results """ struct DiscreteControl <: AbstractParameterNode node_id::Vector{NodeID} - listen_node_id::Vector{NodeID} - variable::Vector{String} - look_ahead::Vector{Float64} + listen_node_id::Vector{Vector{NodeID}} + variable::Vector{Vector{String}} + weight::Vector{Vector{Float64}} + look_ahead::Vector{Vector{Float64}} greater_than::Vector{Float64} condition_value::Vector{Bool} control_state::Dict{NodeID, Tuple{String, Float64}} diff --git a/core/src/read.jl b/core/src/read.jl index 20643c943..af4e9345b 100644 --- a/core/src/read.jl +++ b/core/src/read.jl @@ -542,9 +542,49 @@ function Basin(db::DB, config::Config, chunk_sizes::Vector{Int})::Basin ) end +function get_compound_variables(compound_variable, condition) + listen_node_id = Vector{NodeID}[] + variable = Vector{String}[] + weight = Vector{Float64}[] + look_ahead = Vector{Float64}[] + + for cond in condition + if cond.listen_node_type == "compound" + compound_variable_data = filter( + row -> (row.node_id, row.name) == (cond.node_id, cond.variable), + compound_variable, + ) + listen_node_id_data = + NodeID.( + compound_variable_data.listen_node_type, + compound_variable_data.listen_node_id, + ) + @assert !isempty(listen_node_id_data) "No compound variable data found for name $(cond.variable)." + variable_data = compound_variable_data.variable + weight_data = compound_variable_data.weight + look_ahead_data = coalesce.(compound_variable_data.look_ahead, 0.0) + else + listen_node_id_data = [NodeID(cond.listen_node_type, cond.listen_node_id)] + variable_data = [cond.variable] + weight_data = [1.0] + look_ahead_data = [coalesce(cond.look_ahead, 0.0)] + end + + push!(listen_node_id, listen_node_id_data) + push!(variable, variable_data) + push!(weight, weight_data) + push!(look_ahead, look_ahead_data) + end + return listen_node_id, variable, weight, look_ahead +end + function DiscreteControl(db::DB, config::Config)::DiscreteControl + compound_variable = load_structvector(db, config, DiscreteControlCompoundvariableV1) condition = load_structvector(db, config, DiscreteControlConditionV1) + listen_node_id, variable, weight, look_ahead = + get_compound_variables(compound_variable, condition) + condition_value = fill(false, length(condition.node_id)) control_state::Dict{NodeID, Tuple{String, Float64}} = Dict() @@ -557,7 +597,6 @@ function DiscreteControl(db::DB, config::Config)::DiscreteControl end logic = load_structvector(db, config, DiscreteControlLogicV1) - logic_mapping = Dict{Tuple{NodeID, String}, String}() for (node_id, truth_state, control_state_) in @@ -567,7 +606,6 @@ function DiscreteControl(db::DB, config::Config)::DiscreteControl end logic_mapping = expand_logic_mapping(logic_mapping) - look_ahead = coalesce.(condition.look_ahead, 0.0) record = ( time = Float64[], @@ -578,8 +616,9 @@ function DiscreteControl(db::DB, config::Config)::DiscreteControl return DiscreteControl( NodeID.(NodeType.DiscreteControl, condition.node_id), # Not unique - NodeID.(condition.listen_node_type, condition.listen_node_id), - condition.variable, + listen_node_id, + variable, + weight, look_ahead, condition.greater_than, condition_value, diff --git a/core/src/schema.jl b/core/src/schema.jl index 2d3896e82..95d92d38b 100644 --- a/core/src/schema.jl +++ b/core/src/schema.jl @@ -191,6 +191,7 @@ end listen_node_id::Int variable::String weight::Float64 + look_ahead::Union{Missing, Float64} end @version DiscreteControlConditionV1 begin diff --git a/python/ribasim/ribasim/schemas.py b/python/ribasim/ribasim/schemas.py index 213d97667..24dd99f76 100644 --- a/python/ribasim/ribasim/schemas.py +++ b/python/ribasim/ribasim/schemas.py @@ -55,6 +55,7 @@ class DiscreteControlCompoundvariableSchema(_BaseSchema): listen_node_id: Series[Int32] = pa.Field(nullable=False, default=0) variable: Series[str] = pa.Field(nullable=False) weight: Series[float] = pa.Field(nullable=False) + look_ahead: Series[float] = pa.Field(nullable=True) class DiscreteControlConditionSchema(_BaseSchema): From 0cec53786665261ab9771ebb93a1b3a87e07d3be Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Wed, 10 Apr 2024 17:04:53 +0200 Subject: [PATCH 04/22] Handle compound variables in discrete_control callback --- core/src/callback.jl | 9 +++++++-- .../ribasim_testmodels/discrete_control.py | 2 +- 2 files changed, 8 insertions(+), 3 deletions(-) diff --git a/core/src/callback.jl b/core/src/callback.jl index bd94043a5..ee49c6086 100644 --- a/core/src/callback.jl +++ b/core/src/callback.jl @@ -194,15 +194,20 @@ function discrete_control_condition(out, u, t, integrator) (; p) = integrator (; discrete_control) = p - for (i, (listen_node_id, variable, greater_than, look_ahead)) in enumerate( + for (i, (listen_node_ids, variables, weights, greater_than, look_aheads)) in enumerate( zip( discrete_control.listen_node_id, discrete_control.variable, + discrete_control.weight, discrete_control.greater_than, discrete_control.look_ahead, ), ) - value = get_value(p, listen_node_id, variable, look_ahead, u, t) + value = 0.0 + for (listen_node_id, variable, weight, look_ahead) in + zip(listen_node_ids, variables, weights, look_aheads) + value += weight * get_value(p, listen_node_id, variable, look_ahead, u, t) + end diff = value - greater_than out[i] = diff end diff --git a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py index ea9eb3899..55aa1849b 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py +++ b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py @@ -430,7 +430,7 @@ def compound_variable_condition_model() -> Model: name="flow_mean", listen_node_type="FlowBoundary", listen_node_id=[2, 3], - variable="flow", + variable="flow_rate", weight=0.5, ), discrete_control.Condition( From 78e8364036f1ac720cd03bb13ffe1b678e8ea368 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Wed, 10 Apr 2024 18:24:48 +0200 Subject: [PATCH 05/22] Fix part of tests --- core/src/callback.jl | 7 +++++-- core/src/validation.jl | 43 +++++++++++++++++++++------------------ core/test/control_test.jl | 4 ++-- 3 files changed, 30 insertions(+), 24 deletions(-) diff --git a/core/src/callback.jl b/core/src/callback.jl index ee49c6086..0ac0701d4 100644 --- a/core/src/callback.jl +++ b/core/src/callback.jl @@ -277,7 +277,10 @@ function discrete_control_affect_upcrossing!(integrator, condition_idx) # only possibly the du. Parameter changes can change the flow on an edge discontinuously, # giving the possibility of logical paradoxes where certain parameter changes immediately # undo the truth state that caused that parameter change. - is_basin = id_index(basin.node_id, discrete_control.listen_node_id[condition_idx])[1] + listen_node_ids = discrete_control.listen_node_id[condition_idx] + is_basin = + length(listen_node_ids) == 1 ? id_index(basin.node_id, only(listen_node_ids))[1] : + false # NOTE: The above no longer works when listen feature ids can be something other than node ids # I think the more durable option is to give all possible condition types a different variable string, # e.g. basin.level and level_boundary.level @@ -374,7 +377,7 @@ function discrete_control_affect!( discrete_control.logic_mapping[(discrete_control_node_id, truth_state)] else error( - "Control state specified for neither $truth_state_crossing_specific nor $truth_state for DiscreteControl node $discrete_control_node_id.", + "Control state specified for neither $truth_state_crossing_specific nor $truth_state for $discrete_control_node_id.", ) end diff --git a/core/src/validation.jl b/core/src/validation.jl index 202eec6a8..478fca586 100644 --- a/core/src/validation.jl +++ b/core/src/validation.jl @@ -552,29 +552,32 @@ function valid_discrete_control(p::Parameters, config::Config)::Bool end end end - for (Δt, var, node_id) in zip(look_ahead, variable, listen_node_id) - if !iszero(Δt) - node_type = node_id.type - # TODO: If more transient listen variables must be supported, this validation must be more specific - # (e.g. for some node some variables are transient, some not). - if node_type ∉ [NodeType.FlowBoundary, NodeType.LevelBoundary] - errors = true - @error "Look ahead supplied for non-timeseries listen variable '$var' from listen node $node_id." - else - if Δt < 0 + for (look_aheads, variables, listen_node_ids) in + zip(look_ahead, variable, listen_node_id) + for (Δt, var, node_id) in zip(look_aheads, variables, listen_node_ids) + if !iszero(Δt) + node_type = node_id.type + # TODO: If more transient listen variables must be supported, this validation must be more specific + # (e.g. for some node some variables are transient, some not). + if node_type ∉ [NodeType.FlowBoundary, NodeType.LevelBoundary] errors = true - @error "Negative look ahead supplied for listen variable '$var' from listen node $node_id." + @error "Look ahead supplied for non-timeseries listen variable '$var' from listen node $node_id." else - node = getfield(p, graph[node_id].type) - idx = if node_type == NodeType.Basin - id_index(node.node_id, node_id) - else - searchsortedfirst(node.node_id, node_id) - end - interpolation = getfield(node, Symbol(var))[idx] - if t_end + Δt > interpolation.t[end] + if Δt < 0 errors = true - @error "Look ahead for listen variable '$var' from listen node $node_id goes past timeseries end during simulation." + @error "Negative look ahead supplied for listen variable '$var' from listen node $node_id." + else + node = getfield(p, graph[node_id].type) + idx = if node_type == NodeType.Basin + id_index(node.node_id, node_id) + else + searchsortedfirst(node.node_id, node_id) + end + interpolation = getfield(node, Symbol(var))[idx] + if t_end + Δt > interpolation.t[end] + errors = true + @error "Look ahead for listen variable '$var' from listen node $node_id goes past timeseries end during simulation." + end end end end diff --git a/core/test/control_test.jl b/core/test/control_test.jl index fdb2e80f8..55d934a10 100644 --- a/core/test/control_test.jl +++ b/core/test/control_test.jl @@ -54,7 +54,7 @@ end p = model.integrator.p (; discrete_control, flow_boundary) = p - Δt = discrete_control.look_ahead[1] + Δt = discrete_control.look_ahead[1][1] t = Ribasim.tsaves(model) t_control = discrete_control.record.time[2] @@ -78,7 +78,7 @@ end p = model.integrator.p (; discrete_control, level_boundary) = p - Δt = discrete_control.look_ahead[1] + Δt = discrete_control.look_ahead[1][1] t = Ribasim.tsaves(model) t_control = discrete_control.record.time[2] From ff54bb34d22c2167661222a7ff8156a3642ea387 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Wed, 10 Apr 2024 22:28:46 +0200 Subject: [PATCH 06/22] Fix the last test in a hacky way --- core/src/callback.jl | 14 ++++++++++---- 1 file changed, 10 insertions(+), 4 deletions(-) diff --git a/core/src/callback.jl b/core/src/callback.jl index 0ac0701d4..7e6daec53 100644 --- a/core/src/callback.jl +++ b/core/src/callback.jl @@ -281,15 +281,16 @@ function discrete_control_affect_upcrossing!(integrator, condition_idx) is_basin = length(listen_node_ids) == 1 ? id_index(basin.node_id, only(listen_node_ids))[1] : false + # NOTE: The above no longer works when listen feature ids can be something other than node ids # I think the more durable option is to give all possible condition types a different variable string, # e.g. basin.level and level_boundary.level - if variable[condition_idx] == "level" && control_state_change && is_basin + if variable[condition_idx][1] == "level" && control_state_change && is_basin # Calling water_balance is expensive, but it is a sure way of getting # du for the basin of this level condition du = zero(u) water_balance!(du, u, p, t) - _, condition_basin_idx = id_index(basin.node_id, listen_node_id[condition_idx]) + _, condition_basin_idx = id_index(basin.node_id, listen_node_id[condition_idx][1]) if du[condition_basin_idx] < 0.0 condition_value[condition_idx] = false @@ -316,13 +317,18 @@ function discrete_control_affect_downcrossing!(integrator, condition_idx) # only possibly the du. Parameter changes can change the flow on an edge discontinuously, # giving the possibility of logical paradoxes where certain parameter changes immediately # undo the truth state that caused that parameter change. - if variable[condition_idx] == "level" && control_state_change + listen_node_ids = discrete_control.listen_node_id[condition_idx] + is_basin = + length(listen_node_ids) == 1 ? id_index(basin.node_id, only(listen_node_ids))[1] : + false + + if variable[condition_idx][1] == "level" && control_state_change && is_basin # Calling water_balance is expensive, but it is a sure way of getting # du for the basin of this level condition du = zero(u) water_balance!(du, u, p, t) has_index, condition_basin_idx = - id_index(basin.node_id, listen_node_id[condition_idx]) + id_index(basin.node_id, listen_node_id[condition_idx][1]) if has_index && du[condition_basin_idx] > 0.0 condition_value[condition_idx] = true From acca9496859abb49abc1b67ad178bad07350f0e6 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 11 Apr 2024 07:09:16 +0200 Subject: [PATCH 07/22] Add tests --- core/test/control_test.jl | 20 +++++++++++++++++++ .../ribasim_testmodels/discrete_control.py | 4 ++-- 2 files changed, 22 insertions(+), 2 deletions(-) diff --git a/core/test/control_test.jl b/core/test/control_test.jl index 55d934a10..e0913ba1d 100644 --- a/core/test/control_test.jl +++ b/core/test/control_test.jl @@ -206,3 +206,23 @@ end @test isapprox(level[t_idx_target_jump], target_high, atol = 1e-1) @test isapprox(level[end], target_low, atol = 1e-1) end + +@testitem "Compound condition" begin + using Ribasim: NodeID + + toml_path = normpath( + @__DIR__, + "../../generated_testmodels/compound_variable_condition/ribasim.toml", + ) + @test ispath(toml_path) + model = Ribasim.run(toml_path) + (; discrete_control) = model.integrator.p + (; listen_node_id, variable, weight, record) = discrete_control + + @test listen_node_id == [[NodeID(:FlowBoundary, 2), NodeID(:FlowBoundary, 3)]] + @test variable == [["flow_rate", "flow_rate"]] + @test weight == [[0.5, 0.5]] + @test record.time ≈ [0.0, model.integrator.sol.t[end] / 2] + @test record.truth_state == ["F", "T"] + @test record.control_state == ["Off", "On"] +end diff --git a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py index 55aa1849b..a81b3e64c 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py +++ b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py @@ -416,7 +416,7 @@ def compound_variable_condition_model() -> Model: ) model.flow_boundary.add( Node(3, Point(0, 1)), - [flow_boundary.Time(time=["2020-01-01", "2022-01-01"], flow_rate=[0.0, 2.0])], + [flow_boundary.Time(time=["2020-01-01", "2021-01-01"], flow_rate=[0.0, 2.0])], ) model.pump.add( Node(4, Point(2, 0)), @@ -437,7 +437,7 @@ def compound_variable_condition_model() -> Model: listen_node_type="compound", listen_node_id=[0], # Irrelevant variable="flow_mean", - greater_than=1.0, + greater_than=0.5, ), discrete_control.Logic(truth_state=["T", "F"], control_state=["On", "Off"]), ], From 72e65b57b07385670aeb086ace947f6722fc51d5 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 11 Apr 2024 07:15:24 +0200 Subject: [PATCH 08/22] small plotting fix --- python/ribasim/ribasim/model.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/python/ribasim/ribasim/model.py b/python/ribasim/ribasim/model.py index cbe9798e9..bb0b74f2c 100644 --- a/python/ribasim/ribasim/model.py +++ b/python/ribasim/ribasim/model.py @@ -323,7 +323,7 @@ def plot_control_listen(self, ax): # Listen edges from DiscreteControl for table in ( self.discrete_control.condition.df, - self.discrete_control.compound_variable.df, + self.discrete_control.compoundvariable.df, ): if table is None: continue From 44eab71344d1f39ac54689446009b02cca674f07 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 11 Apr 2024 07:21:08 +0200 Subject: [PATCH 09/22] small qgis fix --- ribasim_qgis/core/nodes.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/ribasim_qgis/core/nodes.py b/ribasim_qgis/core/nodes.py index 70750e419..7fe9390f3 100644 --- a/ribasim_qgis/core/nodes.py +++ b/ribasim_qgis/core/nodes.py @@ -709,6 +709,7 @@ def attributes(cls) -> list[QgsField]: QgsField("listen_node_id", QVariant.Int), QgsField("variable", QVariant.String), QgsField("weight", QVariant.Double), + QgsField("look_ahead", QVariant.Double), ] @@ -729,6 +730,7 @@ def attributes(cls) -> list[QgsField]: QgsField("listen_node_id", QVariant.Int), QgsField("variable", QVariant.String), QgsField("greater_than", QVariant.Double), + QgsField("look_ahead", QVariant.Double), ] From 25997fea2bb007cc77798e1e18617f3604aaa0a3 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 11 Apr 2024 07:35:17 +0200 Subject: [PATCH 10/22] Update schemas in usage.qmd --- docs/core/usage.qmd | 17 +++++++++++++++-- 1 file changed, 15 insertions(+), 2 deletions(-) diff --git a/docs/core/usage.qmd b/docs/core/usage.qmd index 3613fa409..4cd4b3788 100644 --- a/docs/core/usage.qmd +++ b/docs/core/usage.qmd @@ -204,6 +204,7 @@ name it must have in the database if it is stored there. - Terminal: Water sink without state or properties - `Terminal / static`: - (only node IDs) - DiscreteControl: Set parameters of other nodes based on model state conditions (e.g. basin level) + - `DiscreteControl / compoundvariable` : Linear combinations of variables, e.g. averages or differences of levels - `DisceteControl / condition`: Conditions of the form 'the level in the basin with node id `n` is bigger than 2.0 m' - `DisceteControl / logic`: Translates the truth value of a set of conditions to parameter values for a controlled node - PidControl: Controls the level in a basin by continuously controlling the flow rate of a connected pump or outlet. See also [Wikipedia](https://en.wikipedia.org/wiki/PID_controller) and [PID controller in equations](equations.qmd#sec-PID). @@ -625,6 +626,18 @@ node_id | Int32 | - | sorted DiscreteControl is implemented based on [VectorContinuousCallback](https://docs.sciml.ai/DiffEqDocs/stable/features/callback_functions/#VectorContinuousCallback). +## DiscreteControl / compoundvariable + +column | type | unit | restriction +----------------- | -------- | ------- | ----------- +node_id | Int32 | - | sorted +name | String | - | - +listen_node_type | String | - | known node type +listen_node_id | Int32 | - | sorted per node_id +variable | String | - | must be "level" or "flow_rate", sorted per listen_node_id +weight | Float64 | - | - +look_ahead | Float64 | $s$ | Only on transient boundary conditions, non-negative (optional, default 0). + ## DiscreteControl / condition The condition schema defines conditions of the form 'the discrete_control node with this node id listens to whether the given variable of the node with the given listen feature id is grater than the given value'. If the condition variable comes from a time-series, a look ahead $\Delta t$ can be supplied. @@ -632,9 +645,9 @@ The condition schema defines conditions of the form 'the discrete_control node w column | type | unit | restriction ----------------- | -------- | ------- | ----------- node_id | Int32 | - | sorted -listen_node_type | String | - | known node type +listen_node_type | String | - | known node type or "compound" listen_node_id | Int32 | - | sorted per node_id -variable | String | - | must be "level" or "flow_rate", sorted per listen_node_id +variable | String | - | must be "level", "flow_rate" or a compound variable name, sorted per listen_node_id greater_than | Float64 | various | sorted per variable look_ahead | Float64 | $s$ | Only on transient boundary conditions, non-negative (optional, default 0) From 3144ef1098b6528137492c6a6597354db88c5b74 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 11 Apr 2024 07:42:23 +0200 Subject: [PATCH 11/22] Update docstrings --- core/src/parameter.jl | 6 ++++-- core/src/read.jl | 4 ++++ 2 files changed, 8 insertions(+), 2 deletions(-) diff --git a/core/src/parameter.jl b/core/src/parameter.jl index 2d8c5bca0..231f0a9c4 100644 --- a/core/src/parameter.jl +++ b/core/src/parameter.jl @@ -445,8 +445,10 @@ end """ node_id: node ID of the DiscreteControl node; these are not unique but repeated by the amount of conditions of this DiscreteControl node -listen_node_id: the ID of the node being condition on -variable: the name of the variable in the condition +listen_node_id: the IDs of the nodes being condition on +variable: the names of the variables in the condition +weight: the weight of the variables in the condition +look_ahead: the look ahead of variables in the condition in seconds greater_than: The threshold value in the condition condition_value: The current value of each condition control_state: Dictionary: node ID => (control state, control state start) diff --git a/core/src/read.jl b/core/src/read.jl index af4e9345b..bed7abb0c 100644 --- a/core/src/read.jl +++ b/core/src/read.jl @@ -542,6 +542,10 @@ function Basin(db::DB, config::Config, chunk_sizes::Vector{Int})::Basin ) end +""" +Parse compound variables from the compound variable table, and turn +variables in the condition table into trivial compound variables +""" function get_compound_variables(compound_variable, condition) listen_node_id = Vector{NodeID}[] variable = Vector{String}[] From e842bdfd9ea586f7acc8f9858cc7a23bc572c401 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 11 Apr 2024 09:00:07 +0200 Subject: [PATCH 12/22] update sorting and schema explanations --- core/src/validation.jl | 2 ++ docs/core/usage.qmd | 13 +++++++++++-- python/ribasim/ribasim/config.py | 6 ++---- 3 files changed, 15 insertions(+), 6 deletions(-) diff --git a/core/src/validation.jl b/core/src/validation.jl index 478fca586..1280049dc 100644 --- a/core/src/validation.jl +++ b/core/src/validation.jl @@ -102,6 +102,7 @@ sort_by_id_state_level(row) = (row.node_id, row.control_state, row.level) sort_by_priority(row) = (row.node_id, row.priority) sort_by_priority_time(row) = (row.node_id, row.priority, row.time) sort_by_subgrid_level(row) = (row.subgrid_id, row.basin_level) +sort_by_id_name(row) = (row.node_id, row.name) # get the right sort by function given the Schema, with sort_by_id as the default sort_by_function(table::StructVector{<:Legolas.AbstractRecord}) = sort_by_id @@ -110,6 +111,7 @@ sort_by_function(table::StructVector{BasinProfileV1}) = sort_by_id_level sort_by_function(table::StructVector{UserDemandStaticV1}) = sort_by_priority sort_by_function(table::StructVector{UserDemandTimeV1}) = sort_by_priority_time sort_by_function(table::StructVector{BasinSubgridV1}) = sort_by_subgrid_level +sort_by_function(table::StructVector{DiscreteControlCompoundvariableV1}) = sort_by_id_name const TimeSchemas = Union{ BasinTimeV1, diff --git a/docs/core/usage.qmd b/docs/core/usage.qmd index 4cd4b3788..aa7c5f73a 100644 --- a/docs/core/usage.qmd +++ b/docs/core/usage.qmd @@ -628,10 +628,18 @@ DiscreteControl is implemented based on [VectorContinuousCallback](https://docs. ## DiscreteControl / compoundvariable +The compound variable schema defines linear combinations of variables which can be used in conditions. This means that +this schema defines new variables with the given name that look like +$$ +\text{weight}_1 * \text{variable}_1 + \text{weight}_2 * \text{variable}_2 + \ldots, +$$ + +which can be for instance an average or a difference of variables. If a variable comes from a time-series, a look ahead $\Delta t$ can be supplied. + column | type | unit | restriction ----------------- | -------- | ------- | ----------- node_id | Int32 | - | sorted -name | String | - | - +name | String | - | sorted per node_id listen_node_type | String | - | known node type listen_node_id | Int32 | - | sorted per node_id variable | String | - | must be "level" or "flow_rate", sorted per listen_node_id @@ -640,7 +648,8 @@ look_ahead | Float64 | $s$ | Only on transient boundary conditions, ## DiscreteControl / condition -The condition schema defines conditions of the form 'the discrete_control node with this node id listens to whether the given variable of the node with the given listen feature id is grater than the given value'. If the condition variable comes from a time-series, a look ahead $\Delta t$ can be supplied. +The condition schema defines conditions of the form 'the discrete_control node with this node id listens to whether the given variable of the node with the given listen node id is greater than the given value'. If the condition variable comes from a time-series, a look ahead $\Delta t$ can be supplied. +Compound variables from the schema above can also be used, in which case the listen node ID and look ahead in this schema are ignored. column | type | unit | restriction ----------------- | -------- | ------- | ----------- diff --git a/python/ribasim/ribasim/config.py b/python/ribasim/ribasim/config.py index 64a5019e9..3fed1f8a7 100644 --- a/python/ribasim/ribasim/config.py +++ b/python/ribasim/ribasim/config.py @@ -290,13 +290,11 @@ class ManningResistance(MultiNodeModel): class DiscreteControl(MultiNodeModel): compoundvariable: TableModel[DiscreteControlCompoundvariableSchema] = Field( default_factory=TableModel[DiscreteControlCompoundvariableSchema], - json_schema_extra={"sort_keys": ["name"]}, + json_schema_extra={"sort_keys": ["node_id", "name"]}, ) condition: TableModel[DiscreteControlConditionSchema] = Field( default_factory=TableModel[DiscreteControlConditionSchema], - json_schema_extra={ - "sort_keys": ["node_id", "listen_node_id", "variable", "greater_than"] - }, + json_schema_extra={"sort_keys": ["node_id"]}, ) logic: TableModel[DiscreteControlLogicSchema] = Field( default_factory=TableModel[DiscreteControlLogicSchema], From c975b0c0bac347fdf66398740fa7b5bc11950ddc Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 11 Apr 2024 11:19:54 +0200 Subject: [PATCH 13/22] Fix table sorting --- core/src/validation.jl | 5 +++++ python/ribasim/ribasim/config.py | 10 +++++++++- python/ribasim/tests/test_io.py | 1 + 3 files changed, 15 insertions(+), 1 deletion(-) diff --git a/core/src/validation.jl b/core/src/validation.jl index 1280049dc..d4b46b203 100644 --- a/core/src/validation.jl +++ b/core/src/validation.jl @@ -102,7 +102,10 @@ sort_by_id_state_level(row) = (row.node_id, row.control_state, row.level) sort_by_priority(row) = (row.node_id, row.priority) sort_by_priority_time(row) = (row.node_id, row.priority, row.time) sort_by_subgrid_level(row) = (row.subgrid_id, row.basin_level) +sort_by_condition(row) = + (row.node_id, row.listen_node_type, row.listen_node_id, row.variable, row.greater_than) sort_by_id_name(row) = (row.node_id, row.name) +sort_by_truth_state(row) = (row.node_id, row.truth_state) # get the right sort by function given the Schema, with sort_by_id as the default sort_by_function(table::StructVector{<:Legolas.AbstractRecord}) = sort_by_id @@ -112,6 +115,8 @@ sort_by_function(table::StructVector{UserDemandStaticV1}) = sort_by_priority sort_by_function(table::StructVector{UserDemandTimeV1}) = sort_by_priority_time sort_by_function(table::StructVector{BasinSubgridV1}) = sort_by_subgrid_level sort_by_function(table::StructVector{DiscreteControlCompoundvariableV1}) = sort_by_id_name +sort_by_function(table::StructVector{DiscreteControlConditionV1}) = sort_by_condition +sort_by_function(table::StructVector{DiscreteControlLogicV1}) = sort_by_truth_state const TimeSchemas = Union{ BasinTimeV1, diff --git a/python/ribasim/ribasim/config.py b/python/ribasim/ribasim/config.py index 3fed1f8a7..e00f20ed2 100644 --- a/python/ribasim/ribasim/config.py +++ b/python/ribasim/ribasim/config.py @@ -294,7 +294,15 @@ class DiscreteControl(MultiNodeModel): ) condition: TableModel[DiscreteControlConditionSchema] = Field( default_factory=TableModel[DiscreteControlConditionSchema], - json_schema_extra={"sort_keys": ["node_id"]}, + json_schema_extra={ + "sort_keys": [ + "node_id", + "listen_node_type", + "listen_node_id", + "variable", + "greater_than", + ] + }, ) logic: TableModel[DiscreteControlLogicSchema] = Field( default_factory=TableModel[DiscreteControlLogicSchema], diff --git a/python/ribasim/tests/test_io.py b/python/ribasim/tests/test_io.py index bd6fcc7a8..e9b724eaa 100644 --- a/python/ribasim/tests/test_io.py +++ b/python/ribasim/tests/test_io.py @@ -102,6 +102,7 @@ def test_sort(level_setpoint_with_minmax, tmp_path): assert table.df.iloc[0]["greater_than"] == 15.0 assert table._sort_keys == [ "node_id", + "listen_node_type", "listen_node_id", "variable", "greater_than", From 0ba1d6fdfb66884bb0b60afe95d7d85fc056b453 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Mon, 15 Apr 2024 15:34:40 +0200 Subject: [PATCH 14/22] Refactor variable input --- core/src/read.jl | 50 ++++--------------- core/src/schema.jl | 11 ++-- core/src/validation.jl | 8 +-- python/ribasim/ribasim/config.py | 18 ++++--- python/ribasim/ribasim/model.py | 2 +- .../ribasim/ribasim/nodes/discrete_control.py | 6 +-- python/ribasim/ribasim/schemas.py | 23 ++++----- .../ribasim_testmodels/allocation.py | 18 ++++--- .../ribasim_testmodels/discrete_control.py | 44 +++++++++------- .../ribasim_testmodels/dutch_waterways.py | 6 ++- .../ribasim_testmodels/invalid.py | 6 ++- .../ribasim_testmodels/pid_control.py | 6 ++- ribasim_qgis/core/nodes.py | 9 +--- 13 files changed, 94 insertions(+), 113 deletions(-) diff --git a/core/src/read.jl b/core/src/read.jl index bed7abb0c..b9bf0c8ec 100644 --- a/core/src/read.jl +++ b/core/src/read.jl @@ -542,52 +542,22 @@ function Basin(db::DB, config::Config, chunk_sizes::Vector{Int})::Basin ) end -""" -Parse compound variables from the compound variable table, and turn -variables in the condition table into trivial compound variables -""" -function get_compound_variables(compound_variable, condition) +function DiscreteControl(db::DB, config::Config)::DiscreteControl + compound_variable = load_structvector(db, config, DiscreteControlVariableV1) + condition = load_structvector(db, config, DiscreteControlConditionV1) + listen_node_id = Vector{NodeID}[] variable = Vector{String}[] weight = Vector{Float64}[] look_ahead = Vector{Float64}[] - for cond in condition - if cond.listen_node_type == "compound" - compound_variable_data = filter( - row -> (row.node_id, row.name) == (cond.node_id, cond.variable), - compound_variable, - ) - listen_node_id_data = - NodeID.( - compound_variable_data.listen_node_type, - compound_variable_data.listen_node_id, - ) - @assert !isempty(listen_node_id_data) "No compound variable data found for name $(cond.variable)." - variable_data = compound_variable_data.variable - weight_data = compound_variable_data.weight - look_ahead_data = coalesce.(compound_variable_data.look_ahead, 0.0) - else - listen_node_id_data = [NodeID(cond.listen_node_type, cond.listen_node_id)] - variable_data = [cond.variable] - weight_data = [1.0] - look_ahead_data = [coalesce(cond.look_ahead, 0.0)] - end - - push!(listen_node_id, listen_node_id_data) - push!(variable, variable_data) - push!(weight, weight_data) - push!(look_ahead, look_ahead_data) + for id in unique(condition.node_id) + group = filter(row -> row.node_id == id, compound_variable) + push!(listen_node_id, NodeID.(group.listen_node_type, group.listen_node_id)) + push!(variable, group.variable) + push!(weight, coalesce.(group.weight, 1.0)) + push!(look_ahead, coalesce.(group.look_ahead, 0.0)) end - return listen_node_id, variable, weight, look_ahead -end - -function DiscreteControl(db::DB, config::Config)::DiscreteControl - compound_variable = load_structvector(db, config, DiscreteControlCompoundvariableV1) - condition = load_structvector(db, config, DiscreteControlConditionV1) - - listen_node_id, variable, weight, look_ahead = - get_compound_variables(compound_variable, condition) condition_value = fill(false, length(condition.node_id)) control_state::Dict{NodeID, Tuple{String, Float64}} = Dict() diff --git a/core/src/schema.jl b/core/src/schema.jl index 95d92d38b..60147c38a 100644 --- a/core/src/schema.jl +++ b/core/src/schema.jl @@ -1,8 +1,8 @@ # These schemas define the name of database tables and the configuration file structure # The identifier is parsed as ribasim.nodetype.kind, no capitals or underscores are allowed. +@schema "ribasim.discretecontrol.variable" DiscreteControlVariable @schema "ribasim.discretecontrol.condition" DiscreteControlCondition @schema "ribasim.discretecontrol.logic" DiscreteControlLogic -@schema "ribasim.discretecontrol.compoundvariable" DiscreteControlCompoundvariable @schema "ribasim.basin.static" BasinStatic @schema "ribasim.basin.time" BasinTime @schema "ribasim.basin.profile" BasinProfile @@ -184,23 +184,18 @@ end node_id::Int32 end -@version DiscreteControlCompoundvariableV1 begin +@version DiscreteControlVariableV1 begin node_id::Int32 - name::String listen_node_type::String listen_node_id::Int variable::String - weight::Float64 + weight::Union{Missing, Float64} look_ahead::Union{Missing, Float64} end @version DiscreteControlConditionV1 begin node_id::Int32 - listen_node_type::String - listen_node_id::Int32 - variable::String greater_than::Float64 - look_ahead::Union{Missing, Float64} end @version DiscreteControlLogicV1 begin diff --git a/core/src/validation.jl b/core/src/validation.jl index d4b46b203..846cca5b4 100644 --- a/core/src/validation.jl +++ b/core/src/validation.jl @@ -104,7 +104,9 @@ sort_by_priority_time(row) = (row.node_id, row.priority, row.time) sort_by_subgrid_level(row) = (row.subgrid_id, row.basin_level) sort_by_condition(row) = (row.node_id, row.listen_node_type, row.listen_node_id, row.variable, row.greater_than) -sort_by_id_name(row) = (row.node_id, row.name) +sort_by_variable(row) = + (row.node_id, row.listen_node_type, row.listen_node_id, row.variable) +sort_by_id_greater_than(row) = (row.node_id, row.greater_than) sort_by_truth_state(row) = (row.node_id, row.truth_state) # get the right sort by function given the Schema, with sort_by_id as the default @@ -114,8 +116,8 @@ sort_by_function(table::StructVector{BasinProfileV1}) = sort_by_id_level sort_by_function(table::StructVector{UserDemandStaticV1}) = sort_by_priority sort_by_function(table::StructVector{UserDemandTimeV1}) = sort_by_priority_time sort_by_function(table::StructVector{BasinSubgridV1}) = sort_by_subgrid_level -sort_by_function(table::StructVector{DiscreteControlCompoundvariableV1}) = sort_by_id_name -sort_by_function(table::StructVector{DiscreteControlConditionV1}) = sort_by_condition +sort_by_function(table::StructVector{DiscreteControlVariableV1}) = sort_by_variable +sort_by_function(table::StructVector{DiscreteControlConditionV1}) = sort_by_id_greater_than sort_by_function(table::StructVector{DiscreteControlLogicV1}) = sort_by_truth_state const TimeSchemas = Union{ diff --git a/python/ribasim/ribasim/config.py b/python/ribasim/ribasim/config.py index e00f20ed2..dd7fbbf0d 100644 --- a/python/ribasim/ribasim/config.py +++ b/python/ribasim/ribasim/config.py @@ -21,9 +21,9 @@ BasinStaticSchema, BasinSubgridSchema, BasinTimeSchema, - DiscreteControlCompoundvariableSchema, DiscreteControlConditionSchema, DiscreteControlLogicSchema, + DiscreteControlVariableSchema, FlowBoundaryStaticSchema, FlowBoundaryTimeSchema, FlowDemandStaticSchema, @@ -288,18 +288,22 @@ class ManningResistance(MultiNodeModel): class DiscreteControl(MultiNodeModel): - compoundvariable: TableModel[DiscreteControlCompoundvariableSchema] = Field( - default_factory=TableModel[DiscreteControlCompoundvariableSchema], - json_schema_extra={"sort_keys": ["node_id", "name"]}, - ) - condition: TableModel[DiscreteControlConditionSchema] = Field( - default_factory=TableModel[DiscreteControlConditionSchema], + variable: TableModel[DiscreteControlVariableSchema] = Field( + default_factory=TableModel[DiscreteControlVariableSchema], json_schema_extra={ "sort_keys": [ "node_id", "listen_node_type", "listen_node_id", "variable", + ] + }, + ) + condition: TableModel[DiscreteControlConditionSchema] = Field( + default_factory=TableModel[DiscreteControlConditionSchema], + json_schema_extra={ + "sort_keys": [ + "node_id", "greater_than", ] }, diff --git a/python/ribasim/ribasim/model.py b/python/ribasim/ribasim/model.py index d8d6b08e0..8bbb1666a 100644 --- a/python/ribasim/ribasim/model.py +++ b/python/ribasim/ribasim/model.py @@ -323,7 +323,7 @@ def plot_control_listen(self, ax): # Listen edges from DiscreteControl for table in ( self.discrete_control.condition.df, - self.discrete_control.compoundvariable.df, + self.discrete_control.variable.df, ): if table is None: continue diff --git a/python/ribasim/ribasim/nodes/discrete_control.py b/python/ribasim/ribasim/nodes/discrete_control.py index cc3cfcaf9..df7b9f42b 100644 --- a/python/ribasim/ribasim/nodes/discrete_control.py +++ b/python/ribasim/ribasim/nodes/discrete_control.py @@ -2,15 +2,15 @@ from ribasim.input_base import TableModel from ribasim.schemas import ( - DiscreteControlCompoundvariableSchema, DiscreteControlConditionSchema, DiscreteControlLogicSchema, + DiscreteControlVariableSchema, ) -__all__ = ["Condition", "Logic", "Compoundvariable"] +__all__ = ["Condition", "Logic", "Variable"] -class Compoundvariable(TableModel[DiscreteControlCompoundvariableSchema]): +class Variable(TableModel[DiscreteControlVariableSchema]): def __init__(self, **kwargs): super().__init__(df=DataFrame(dict(**kwargs))) diff --git a/python/ribasim/ribasim/schemas.py b/python/ribasim/ribasim/schemas.py index 24dd99f76..260831282 100644 --- a/python/ribasim/ribasim/schemas.py +++ b/python/ribasim/ribasim/schemas.py @@ -48,23 +48,9 @@ class BasinTimeSchema(_BaseSchema): urban_runoff: Series[float] = pa.Field(nullable=True) -class DiscreteControlCompoundvariableSchema(_BaseSchema): - node_id: Series[Int32] = pa.Field(nullable=False, default=0) - name: Series[str] = pa.Field(nullable=False) - listen_node_type: Series[str] = pa.Field(nullable=False) - listen_node_id: Series[Int32] = pa.Field(nullable=False, default=0) - variable: Series[str] = pa.Field(nullable=False) - weight: Series[float] = pa.Field(nullable=False) - look_ahead: Series[float] = pa.Field(nullable=True) - - class DiscreteControlConditionSchema(_BaseSchema): node_id: Series[Int32] = pa.Field(nullable=False, default=0) - listen_node_type: Series[str] = pa.Field(nullable=False) - listen_node_id: Series[Int32] = pa.Field(nullable=False, default=0) - variable: Series[str] = pa.Field(nullable=False) greater_than: Series[float] = pa.Field(nullable=False) - look_ahead: Series[float] = pa.Field(nullable=True) class DiscreteControlLogicSchema(_BaseSchema): @@ -73,6 +59,15 @@ class DiscreteControlLogicSchema(_BaseSchema): control_state: Series[str] = pa.Field(nullable=False) +class DiscreteControlVariableSchema(_BaseSchema): + node_id: Series[Int32] = pa.Field(nullable=False, default=0) + listen_node_type: Series[str] = pa.Field(nullable=False) + listen_node_id: Series[Int32] = pa.Field(nullable=False, default=0) + variable: Series[str] = pa.Field(nullable=False) + weight: Series[float] = pa.Field(nullable=True) + look_ahead: Series[float] = pa.Field(nullable=True) + + class FlowBoundaryStaticSchema(_BaseSchema): node_id: Series[Int32] = pa.Field(nullable=False, default=0) active: Series[pa.BOOL] = pa.Field(nullable=True) diff --git a/python/ribasim_testmodels/ribasim_testmodels/allocation.py b/python/ribasim_testmodels/ribasim_testmodels/allocation.py index 1e102e8c0..71353cf93 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/allocation.py +++ b/python/ribasim_testmodels/ribasim_testmodels/allocation.py @@ -406,11 +406,13 @@ def fractional_flow_subnetwork_model() -> Model: model.discrete_control.add( Node(10, Point(-1, 2), subnetwork_id=2), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="FlowBoundary", listen_node_id=[1], variable="flow_rate", - greater_than=3e-3, + ), + discrete_control.Condition( + greater_than=[3e-3], ), discrete_control.Logic(truth_state=["F", "T"], control_state=["A", "B"]), ], @@ -500,11 +502,13 @@ def allocation_example_model() -> Model: model.discrete_control.add( Node(11, Point(4.5, 0.25), subnetwork_id=2), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="Basin", listen_node_id=[5], variable="level", - greater_than=0.52, + ), + discrete_control.Condition( + greater_than=[0.52], ), discrete_control.Logic( truth_state=["T", "F"], control_state=["divert", "close"] @@ -679,11 +683,13 @@ def main_network_with_subnetworks_model() -> Model: model.discrete_control.add( Node(33, Point(13, 5), subnetwork_id=5), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="Basin", listen_node_id=[25], variable="level", - greater_than=0.003, + ), + discrete_control.Condition( + greater_than=[0.003], ), discrete_control.Logic(truth_state=["F", "T"], control_state=["A", "B"]), ], diff --git a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py index a81b3e64c..1a6037d8a 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py +++ b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py @@ -54,10 +54,12 @@ def pump_discrete_control_model() -> Model: model.discrete_control.add( Node(5, Point(1, 1)), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="Basin", listen_node_id=[1, 3], variable="level", + ), + discrete_control.Condition( greater_than=[0.8, 0.4], ), discrete_control.Logic( @@ -69,10 +71,12 @@ def pump_discrete_control_model() -> Model: model.discrete_control.add( Node(6, Point(2, -1)), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="Basin", - listen_node_id=3, + listen_node_id=[3], variable="level", + ), + discrete_control.Condition( greater_than=[0.45], ), discrete_control.Logic( @@ -140,13 +144,15 @@ def flow_condition_model() -> Model: model.discrete_control.add( Node(5, Point(1, 1)), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="FlowBoundary", - listen_node_id=1, + listen_node_id=[1], variable="flow_rate", - greater_than=[20 / (86400)], look_ahead=60 * 86400, ), + discrete_control.Condition( + greater_than=[20 / (86400)], + ), discrete_control.Logic(truth_state=["T", "F"], control_state=["off", "on"]), ], ) @@ -203,13 +209,15 @@ def level_boundary_condition_model() -> Model: model.discrete_control.add( Node(6, Point(1.5, 1)), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="LevelBoundary", listen_node_id=[1], variable="level", - greater_than=6.0, look_ahead=60 * 86400, ), + discrete_control.Condition( + greater_than=[6.0], + ), discrete_control.Logic(truth_state=["T", "F"], control_state=["on", "off"]), ], ) @@ -275,11 +283,13 @@ def tabulated_rating_curve_control_model() -> Model: model.discrete_control.add( Node(4, Point(1, 1)), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="Basin", listen_node_id=[1], variable="level", - greater_than=0.5, + ), + discrete_control.Condition( + greater_than=[0.5], ), discrete_control.Logic( truth_state=["T", "F"], control_state=["low", "high"] @@ -342,10 +352,12 @@ def level_setpoint_with_minmax_model() -> Model: model.discrete_control.add( Node(7, Point(1, 0)), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="Basin", - listen_node_id=1, + listen_node_id=[1], variable="level", + ), + discrete_control.Condition( # min, setpoint, max greater_than=[5.0, 10.0, 15.0], ), @@ -426,18 +438,14 @@ def compound_variable_condition_model() -> Model: model.discrete_control.add( Node(6, Point(1, 1)), [ - discrete_control.Compoundvariable( - name="flow_mean", + discrete_control.Variable( listen_node_type="FlowBoundary", listen_node_id=[2, 3], variable="flow_rate", weight=0.5, ), discrete_control.Condition( - listen_node_type="compound", - listen_node_id=[0], # Irrelevant - variable="flow_mean", - greater_than=0.5, + greater_than=[0.5], ), discrete_control.Logic(truth_state=["T", "F"], control_state=["On", "Off"]), ], diff --git a/python/ribasim_testmodels/ribasim_testmodels/dutch_waterways.py b/python/ribasim_testmodels/ribasim_testmodels/dutch_waterways.py index fd990957b..b13bb4bcb 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/dutch_waterways.py +++ b/python/ribasim_testmodels/ribasim_testmodels/dutch_waterways.py @@ -121,10 +121,12 @@ def dutch_waterways_model() -> Model: model.discrete_control.add( Node(17, Point(183612, 441833), name="Controller Driel"), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="FlowBoundary", - listen_node_id=1, + listen_node_id=[1], variable="flow_rate", + ), + discrete_control.Condition( greater_than=[250, 275, 750, 800], ), discrete_control.Logic( diff --git a/python/ribasim_testmodels/ribasim_testmodels/invalid.py b/python/ribasim_testmodels/ribasim_testmodels/invalid.py index 906a7b485..29c9a8405 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/invalid.py +++ b/python/ribasim_testmodels/ribasim_testmodels/invalid.py @@ -148,16 +148,18 @@ def invalid_discrete_control_model() -> Model: model.discrete_control.add( Node(5, Point(1, 1)), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type=["Basin", "FlowBoundary", "FlowBoundary"], listen_node_id=[1, 4, 4], variable=["level", "flow_rate", "flow_rate"], - greater_than=[0.5, 1.5, 1.5], # Invalid: look_ahead can only be specified for timeseries variables. # Invalid: this look_ahead will go past the provided timeseries during simulation. # Invalid: look_ahead must be non-negative. look_ahead=[100.0, 40 * 24 * 60 * 60, -10.0], ), + discrete_control.Condition( + greater_than=[0.5, 1.5, 1.5], + ), # Invalid: DiscreteControl node #4 has 2 conditions so # truth states have to be of length 2 discrete_control.Logic(truth_state=["FFFF"], control_state=["foo"]), diff --git a/python/ribasim_testmodels/ribasim_testmodels/pid_control.py b/python/ribasim_testmodels/ribasim_testmodels/pid_control.py index d8c591c09..05ebc31b7 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/pid_control.py +++ b/python/ribasim_testmodels/ribasim_testmodels/pid_control.py @@ -129,11 +129,13 @@ def discrete_control_of_pid_control_model() -> Model: model.discrete_control.add( Node(7, Point(0, 1)), [ - discrete_control.Condition( + discrete_control.Variable( listen_node_type="LevelBoundary", listen_node_id=[1], variable="level", - greater_than=5.0, + ), + discrete_control.Condition( + greater_than=[5.0], ), discrete_control.Logic( truth_state=["T", "F"], control_state=["target_high", "target_low"] diff --git a/ribasim_qgis/core/nodes.py b/ribasim_qgis/core/nodes.py index 7fe9390f3..5c0393a81 100644 --- a/ribasim_qgis/core/nodes.py +++ b/ribasim_qgis/core/nodes.py @@ -691,10 +691,10 @@ def attributes(cls) -> list[QgsField]: ] -class DiscreteControlCompoundvariable(Input): +class DiscreteControlVariable(Input): @classmethod def input_type(cls) -> str: - return "DiscreteControl / compoundvariable" + return "DiscreteControl / variable" @classmethod def geometry_type(cls) -> str: @@ -704,7 +704,6 @@ def geometry_type(cls) -> str: def attributes(cls) -> list[QgsField]: return [ QgsField("node_id", QVariant.Int), - QgsField("name", QVariant.String), QgsField("listen_node_type", QVariant.String), QgsField("listen_node_id", QVariant.Int), QgsField("variable", QVariant.String), @@ -726,11 +725,7 @@ def geometry_type(cls) -> str: def attributes(cls) -> list[QgsField]: return [ QgsField("node_id", QVariant.Int), - QgsField("listen_node_type", QVariant.String), - QgsField("listen_node_id", QVariant.Int), - QgsField("variable", QVariant.String), QgsField("greater_than", QVariant.Double), - QgsField("look_ahead", QVariant.Double), ] From 17d2a1531e3ff970938e1d91e17b6e8f383d56d1 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Mon, 15 Apr 2024 17:00:33 +0200 Subject: [PATCH 15/22] Add compound_variable_id --- core/src/read.jl | 21 +++++++++++++------ core/src/schema.jl | 2 ++ core/src/validation.jl | 1 - docs/python/examples.ipynb | 8 +++---- python/ribasim/ribasim/config.py | 2 +- python/ribasim/ribasim/model.py | 11 +++------- python/ribasim/ribasim/schemas.py | 2 ++ .../ribasim_testmodels/allocation.py | 6 ++++++ .../ribasim_testmodels/discrete_control.py | 13 ++++++++++++ .../ribasim_testmodels/dutch_waterways.py | 2 ++ .../ribasim_testmodels/invalid.py | 2 ++ .../ribasim_testmodels/pid_control.py | 2 ++ 12 files changed, 52 insertions(+), 20 deletions(-) diff --git a/core/src/read.jl b/core/src/read.jl index b9bf0c8ec..a0aec5fda 100644 --- a/core/src/read.jl +++ b/core/src/read.jl @@ -546,17 +546,26 @@ function DiscreteControl(db::DB, config::Config)::DiscreteControl compound_variable = load_structvector(db, config, DiscreteControlVariableV1) condition = load_structvector(db, config, DiscreteControlConditionV1) + node_id = NodeID[] listen_node_id = Vector{NodeID}[] variable = Vector{String}[] weight = Vector{Float64}[] look_ahead = Vector{Float64}[] for id in unique(condition.node_id) - group = filter(row -> row.node_id == id, compound_variable) - push!(listen_node_id, NodeID.(group.listen_node_type, group.listen_node_id)) - push!(variable, group.variable) - push!(weight, coalesce.(group.weight, 1.0)) - push!(look_ahead, coalesce.(group.look_ahead, 0.0)) + group_id = filter(row -> row.node_id == id, compound_variable) + for group_variable in + StructVector.(IterTools.groupby(row -> row.compound_variable_id, group_id)) + first_row = first(group_variable) + push!(node_id, NodeID(NodeType.DiscreteControl, first_row.node_id)) + push!( + listen_node_id, + NodeID.(group_variable.listen_node_type, group_variable.listen_node_id), + ) + push!(variable, group_variable.variable) + push!(weight, coalesce.(group_variable.weight, 1.0)) + push!(look_ahead, coalesce.(group_variable.look_ahead, 0.0)) + end end condition_value = fill(false, length(condition.node_id)) @@ -589,7 +598,7 @@ function DiscreteControl(db::DB, config::Config)::DiscreteControl ) return DiscreteControl( - NodeID.(NodeType.DiscreteControl, condition.node_id), # Not unique + node_id, # Not unique listen_node_id, variable, weight, diff --git a/core/src/schema.jl b/core/src/schema.jl index 60147c38a..60eb87ac8 100644 --- a/core/src/schema.jl +++ b/core/src/schema.jl @@ -186,6 +186,7 @@ end @version DiscreteControlVariableV1 begin node_id::Int32 + compound_variable_id::Int listen_node_type::String listen_node_id::Int variable::String @@ -195,6 +196,7 @@ end @version DiscreteControlConditionV1 begin node_id::Int32 + compound_variable_id::Int greater_than::Float64 end diff --git a/core/src/validation.jl b/core/src/validation.jl index 846cca5b4..580b4d5ba 100644 --- a/core/src/validation.jl +++ b/core/src/validation.jl @@ -118,7 +118,6 @@ sort_by_function(table::StructVector{UserDemandTimeV1}) = sort_by_priority_time sort_by_function(table::StructVector{BasinSubgridV1}) = sort_by_subgrid_level sort_by_function(table::StructVector{DiscreteControlVariableV1}) = sort_by_variable sort_by_function(table::StructVector{DiscreteControlConditionV1}) = sort_by_id_greater_than -sort_by_function(table::StructVector{DiscreteControlLogicV1}) = sort_by_truth_state const TimeSchemas = Union{ BasinTimeV1, diff --git a/docs/python/examples.ipynb b/docs/python/examples.ipynb index 5a76e3ae7..f78357420 100644 --- a/docs/python/examples.ipynb +++ b/docs/python/examples.ipynb @@ -492,10 +492,10 @@ "model.discrete_control.add(\n", " Node(7, Point(1.0, 0.0)),\n", " [\n", + " discrete_control.Variable(\n", + " listen_node_id=[1], listen_node_type=\"Basin\", variable=\"level\"\n", + " ),\n", " discrete_control.Condition(\n", - " listen_node_id=[1, 1, 1],\n", - " listen_node_type=[\"Basin\", \"Basin\", \"Basin\"],\n", - " variable=[\"level\", \"level\", \"level\"],\n", " greater_than=[5.0, 10.0, 15.0],\n", " ),\n", " discrete_control.Logic(\n", @@ -691,7 +691,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Now run the model with `level_setpoint_with_minmax/ribasim.toml`.\n", + "Now run the model with `ribasim level_setpoint_with_minmax/ribasim.toml`.\n", "After running the model, read back the results:\n" ] }, diff --git a/python/ribasim/ribasim/config.py b/python/ribasim/ribasim/config.py index dd7fbbf0d..d837fc8fc 100644 --- a/python/ribasim/ribasim/config.py +++ b/python/ribasim/ribasim/config.py @@ -310,7 +310,7 @@ class DiscreteControl(MultiNodeModel): ) logic: TableModel[DiscreteControlLogicSchema] = Field( default_factory=TableModel[DiscreteControlLogicSchema], - json_schema_extra={"sort_keys": ["node_id", "truth_state"]}, + json_schema_extra={"sort_keys": ["node_id"]}, ) diff --git a/python/ribasim/ribasim/model.py b/python/ribasim/ribasim/model.py index 8bbb1666a..966a65888 100644 --- a/python/ribasim/ribasim/model.py +++ b/python/ribasim/ribasim/model.py @@ -321,14 +321,9 @@ def plot_control_listen(self, ax): df_listen_edge = pd.concat([df_listen_edge, to_add]) # Listen edges from DiscreteControl - for table in ( - self.discrete_control.condition.df, - self.discrete_control.variable.df, - ): - if table is None: - continue - - to_add = table[ + df_variable = self.discrete_control.variable.df + if df_variable is not None: + to_add = df_variable[ ["node_id", "listen_node_id", "listen_node_type"] ].drop_duplicates() to_add = to_add[to_add["listen_node_type"] != "compound"] diff --git a/python/ribasim/ribasim/schemas.py b/python/ribasim/ribasim/schemas.py index 260831282..2a8b5eab5 100644 --- a/python/ribasim/ribasim/schemas.py +++ b/python/ribasim/ribasim/schemas.py @@ -50,6 +50,7 @@ class BasinTimeSchema(_BaseSchema): class DiscreteControlConditionSchema(_BaseSchema): node_id: Series[Int32] = pa.Field(nullable=False, default=0) + compound_variable_id: Series[Int32] = pa.Field(nullable=False, default=0) greater_than: Series[float] = pa.Field(nullable=False) @@ -61,6 +62,7 @@ class DiscreteControlLogicSchema(_BaseSchema): class DiscreteControlVariableSchema(_BaseSchema): node_id: Series[Int32] = pa.Field(nullable=False, default=0) + compound_variable_id: Series[Int32] = pa.Field(nullable=False, default=0) listen_node_type: Series[str] = pa.Field(nullable=False) listen_node_id: Series[Int32] = pa.Field(nullable=False, default=0) variable: Series[str] = pa.Field(nullable=False) diff --git a/python/ribasim_testmodels/ribasim_testmodels/allocation.py b/python/ribasim_testmodels/ribasim_testmodels/allocation.py index 71353cf93..51980abc5 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/allocation.py +++ b/python/ribasim_testmodels/ribasim_testmodels/allocation.py @@ -410,9 +410,11 @@ def fractional_flow_subnetwork_model() -> Model: listen_node_type="FlowBoundary", listen_node_id=[1], variable="flow_rate", + compound_variable_id=1, ), discrete_control.Condition( greater_than=[3e-3], + compound_variable_id=1, ), discrete_control.Logic(truth_state=["F", "T"], control_state=["A", "B"]), ], @@ -506,9 +508,11 @@ def allocation_example_model() -> Model: listen_node_type="Basin", listen_node_id=[5], variable="level", + compound_variable_id=1, ), discrete_control.Condition( greater_than=[0.52], + compound_variable_id=1, ), discrete_control.Logic( truth_state=["T", "F"], control_state=["divert", "close"] @@ -687,9 +691,11 @@ def main_network_with_subnetworks_model() -> Model: listen_node_type="Basin", listen_node_id=[25], variable="level", + compound_variable_id=1, ), discrete_control.Condition( greater_than=[0.003], + compound_variable_id=1, ), discrete_control.Logic(truth_state=["F", "T"], control_state=["A", "B"]), ], diff --git a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py index 1a6037d8a..d8d911674 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py +++ b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py @@ -58,9 +58,11 @@ def pump_discrete_control_model() -> Model: listen_node_type="Basin", listen_node_id=[1, 3], variable="level", + compound_variable_id=[1, 2], ), discrete_control.Condition( greater_than=[0.8, 0.4], + compound_variable_id=[1, 2], ), discrete_control.Logic( truth_state=["FF", "TF", "FT", "TT"], @@ -75,9 +77,11 @@ def pump_discrete_control_model() -> Model: listen_node_type="Basin", listen_node_id=[3], variable="level", + compound_variable_id=1, ), discrete_control.Condition( greater_than=[0.45], + compound_variable_id=1, ), discrete_control.Logic( truth_state=["T", "F"], @@ -149,9 +153,11 @@ def flow_condition_model() -> Model: listen_node_id=[1], variable="flow_rate", look_ahead=60 * 86400, + compound_variable_id=1, ), discrete_control.Condition( greater_than=[20 / (86400)], + compound_variable_id=1, ), discrete_control.Logic(truth_state=["T", "F"], control_state=["off", "on"]), ], @@ -214,9 +220,11 @@ def level_boundary_condition_model() -> Model: listen_node_id=[1], variable="level", look_ahead=60 * 86400, + compound_variable_id=1, ), discrete_control.Condition( greater_than=[6.0], + compound_variable_id=1, ), discrete_control.Logic(truth_state=["T", "F"], control_state=["on", "off"]), ], @@ -287,6 +295,7 @@ def tabulated_rating_curve_control_model() -> Model: listen_node_type="Basin", listen_node_id=[1], variable="level", + compound_variable_id=1, ), discrete_control.Condition( greater_than=[0.5], @@ -356,10 +365,12 @@ def level_setpoint_with_minmax_model() -> Model: listen_node_type="Basin", listen_node_id=[1], variable="level", + compound_variable_id=1, ), discrete_control.Condition( # min, setpoint, max greater_than=[5.0, 10.0, 15.0], + compound_variable_id=1, ), discrete_control.Logic( truth_state=["FFF", "U**", "T*F", "**D", "TTT"], @@ -443,9 +454,11 @@ def compound_variable_condition_model() -> Model: listen_node_id=[2, 3], variable="flow_rate", weight=0.5, + compound_variable_id=1, ), discrete_control.Condition( greater_than=[0.5], + compound_variable_id=1, ), discrete_control.Logic(truth_state=["T", "F"], control_state=["On", "Off"]), ], diff --git a/python/ribasim_testmodels/ribasim_testmodels/dutch_waterways.py b/python/ribasim_testmodels/ribasim_testmodels/dutch_waterways.py index b13bb4bcb..9fb315a7f 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/dutch_waterways.py +++ b/python/ribasim_testmodels/ribasim_testmodels/dutch_waterways.py @@ -125,9 +125,11 @@ def dutch_waterways_model() -> Model: listen_node_type="FlowBoundary", listen_node_id=[1], variable="flow_rate", + compound_variable_id=1, ), discrete_control.Condition( greater_than=[250, 275, 750, 800], + compound_variable_id=1, ), discrete_control.Logic( truth_state=["FFFF", "U***", "T**F", "***D", "TTTT"], diff --git a/python/ribasim_testmodels/ribasim_testmodels/invalid.py b/python/ribasim_testmodels/ribasim_testmodels/invalid.py index 29c9a8405..463603bf6 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/invalid.py +++ b/python/ribasim_testmodels/ribasim_testmodels/invalid.py @@ -156,9 +156,11 @@ def invalid_discrete_control_model() -> Model: # Invalid: this look_ahead will go past the provided timeseries during simulation. # Invalid: look_ahead must be non-negative. look_ahead=[100.0, 40 * 24 * 60 * 60, -10.0], + compound_variable_id=[1, 2, 3], ), discrete_control.Condition( greater_than=[0.5, 1.5, 1.5], + compound_variable_id=[1, 2, 3], ), # Invalid: DiscreteControl node #4 has 2 conditions so # truth states have to be of length 2 diff --git a/python/ribasim_testmodels/ribasim_testmodels/pid_control.py b/python/ribasim_testmodels/ribasim_testmodels/pid_control.py index 05ebc31b7..4a543acc3 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/pid_control.py +++ b/python/ribasim_testmodels/ribasim_testmodels/pid_control.py @@ -133,9 +133,11 @@ def discrete_control_of_pid_control_model() -> Model: listen_node_type="LevelBoundary", listen_node_id=[1], variable="level", + compound_variable_id=1, ), discrete_control.Condition( greater_than=[5.0], + compound_variable_id=1, ), discrete_control.Logic( truth_state=["T", "F"], control_state=["target_high", "target_low"] From 75b56a077f379e1510fbf9223f1ffa80e6f717e0 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Tue, 16 Apr 2024 11:49:38 +0200 Subject: [PATCH 16/22] Fix bugs --- core/src/callback.jl | 106 +++++++++++++----- core/src/parameter.jl | 7 +- core/src/read.jl | 80 ++++++++++--- core/src/validation.jl | 12 +- core/test/control_test.jl | 12 +- python/ribasim/ribasim/config.py | 3 +- .../ribasim_testmodels/discrete_control.py | 1 + 7 files changed, 158 insertions(+), 63 deletions(-) diff --git a/core/src/callback.jl b/core/src/callback.jl index 7e6daec53..9dbb48acd 100644 --- a/core/src/callback.jl +++ b/core/src/callback.jl @@ -9,10 +9,17 @@ function set_initial_discrete_controlled_parameters!( (; p) = integrator (; discrete_control) = p - n_conditions = length(discrete_control.condition_value) + n_conditions = sum(length(vec) for vec in discrete_control.condition_value) condition_diffs = zeros(Float64, n_conditions) discrete_control_condition(condition_diffs, storage0, integrator.t, integrator) - discrete_control.condition_value .= (condition_diffs .> 0.0) + + idx_start = 1 + for (i, vec) in enumerate(discrete_control.condition_value) + l = length(vec) + idx_end = idx_start + l - 1 + discrete_control.condition_value[i] .= (condition_diffs[idx_start:idx_end] .> 0) + idx_start += l + end # For every discrete_control node find a condition_idx it listens to for discrete_control_node_id in unique(discrete_control.node_id) @@ -88,7 +95,7 @@ function create_callbacks( saved = SavedResults(saved_flow, saved_vertical_flux, saved_subgrid_level) - n_conditions = length(discrete_control.node_id) + n_conditions = sum(length(vec) for vec in discrete_control.greater_than; init = 0) if n_conditions > 0 discrete_control_cb = VectorContinuousCallback( discrete_control_condition, @@ -193,23 +200,25 @@ Listens for changes in condition truths. function discrete_control_condition(out, u, t, integrator) (; p) = integrator (; discrete_control) = p - - for (i, (listen_node_ids, variables, weights, greater_than, look_aheads)) in enumerate( - zip( - discrete_control.listen_node_id, - discrete_control.variable, - discrete_control.weight, - discrete_control.greater_than, - discrete_control.look_ahead, - ), + condition_idx = 0 + + for (listen_node_ids, variables, weights, greater_thans, look_aheads) in zip( + discrete_control.listen_node_id, + discrete_control.variable, + discrete_control.weight, + discrete_control.greater_than, + discrete_control.look_ahead, ) value = 0.0 for (listen_node_id, variable, weight, look_ahead) in zip(listen_node_ids, variables, weights, look_aheads) value += weight * get_value(p, listen_node_id, variable, look_ahead, u, t) end - diff = value - greater_than - out[i] = diff + for greater_than in greater_thans + condition_idx += 1 + diff = value - greater_than + out[condition_idx] = diff + end end end @@ -259,6 +268,22 @@ function get_value( return value end +function get_discrete_control_indices(discrete_control::DiscreteControl, condition_idx::Int) + (; greater_than) = discrete_control + condition_idx_now = 1 + + for (compound_variable_idx, vec) in enumerate(greater_than) + l = length(vec) + + if condition_idx_now + l > condition_idx + greater_than_idx = condition_idx - condition_idx_now + 1 + return compound_variable_idx, greater_than_idx + end + + condition_idx_now += l + end +end + """ An upcrossing means that a condition (always greater than) becomes true. """ @@ -267,7 +292,9 @@ function discrete_control_affect_upcrossing!(integrator, condition_idx) (; discrete_control, basin) = p (; variable, condition_value, listen_node_id) = discrete_control - condition_value[condition_idx] = true + compound_variable_idx, greater_than_idx = + get_discrete_control_indices(discrete_control, condition_idx) + condition_value[compound_variable_idx][greater_than_idx] = true control_state_change = discrete_control_affect!(integrator, condition_idx, true) @@ -277,7 +304,7 @@ function discrete_control_affect_upcrossing!(integrator, condition_idx) # only possibly the du. Parameter changes can change the flow on an edge discontinuously, # giving the possibility of logical paradoxes where certain parameter changes immediately # undo the truth state that caused that parameter change. - listen_node_ids = discrete_control.listen_node_id[condition_idx] + listen_node_ids = discrete_control.listen_node_id[compound_variable_idx] is_basin = length(listen_node_ids) == 1 ? id_index(basin.node_id, only(listen_node_ids))[1] : false @@ -285,15 +312,16 @@ function discrete_control_affect_upcrossing!(integrator, condition_idx) # NOTE: The above no longer works when listen feature ids can be something other than node ids # I think the more durable option is to give all possible condition types a different variable string, # e.g. basin.level and level_boundary.level - if variable[condition_idx][1] == "level" && control_state_change && is_basin + if variable[compound_variable_idx][1] == "level" && control_state_change && is_basin # Calling water_balance is expensive, but it is a sure way of getting # du for the basin of this level condition du = zero(u) water_balance!(du, u, p, t) - _, condition_basin_idx = id_index(basin.node_id, listen_node_id[condition_idx][1]) + _, condition_basin_idx = + id_index(basin.node_id, listen_node_id[compound_variable_idx][1]) if du[condition_basin_idx] < 0.0 - condition_value[condition_idx] = false + condition_value[compound_variable_idx][greater_than_idx] = false discrete_control_affect!(integrator, condition_idx, false) end end @@ -307,7 +335,9 @@ function discrete_control_affect_downcrossing!(integrator, condition_idx) (; discrete_control, basin) = p (; variable, condition_value, listen_node_id) = discrete_control - condition_value[condition_idx] = false + compound_variable_idx, greater_than_idx = + get_discrete_control_indices(discrete_control, condition_idx) + condition_value[compound_variable_idx][greater_than_idx] = false control_state_change = discrete_control_affect!(integrator, condition_idx, false) @@ -317,21 +347,23 @@ function discrete_control_affect_downcrossing!(integrator, condition_idx) # only possibly the du. Parameter changes can change the flow on an edge discontinuously, # giving the possibility of logical paradoxes where certain parameter changes immediately # undo the truth state that caused that parameter change. - listen_node_ids = discrete_control.listen_node_id[condition_idx] + compound_variable_idx, greater_than_idx = + get_discrete_control_indices(discrete_control, condition_idx) + listen_node_ids = discrete_control.listen_node_id[compound_variable_idx] is_basin = length(listen_node_ids) == 1 ? id_index(basin.node_id, only(listen_node_ids))[1] : false - if variable[condition_idx][1] == "level" && control_state_change && is_basin + if variable[compound_variable_idx][1] == "level" && control_state_change && is_basin # Calling water_balance is expensive, but it is a sure way of getting # du for the basin of this level condition du = zero(u) water_balance!(du, u, p, t) has_index, condition_basin_idx = - id_index(basin.node_id, listen_node_id[condition_idx][1]) + id_index(basin.node_id, listen_node_id[compound_variable_idx][1]) if has_index && du[condition_basin_idx] > 0.0 - condition_value[condition_idx] = true + condition_value[compound_variable_idx][greater_than_idx] = true discrete_control_affect!(integrator, condition_idx, true) end end @@ -349,20 +381,34 @@ function discrete_control_affect!( (; discrete_control, graph) = p # Get the discrete_control node that listens to this condition - discrete_control_node_id = discrete_control.node_id[condition_idx] + + compound_variable_idx, _ = get_discrete_control_indices(discrete_control, condition_idx) + discrete_control_node_id = discrete_control.node_id[compound_variable_idx] # Get the indices of all conditions that this control node listens to - condition_ids = discrete_control.node_id .== discrete_control_node_id + where_node_id = searchsorted(discrete_control.node_id, discrete_control_node_id) # Get the truth state for this discrete_control node - truth_values = [ifelse(b, "T", "F") for b in discrete_control.condition_value] - truth_state = join(truth_values[condition_ids], "") + truth_values = cat( + [ + [ifelse(b, "T", "F") for b in discrete_control.condition_value[i]] for + i in where_node_id + ]...; + dims = 1, + ) + truth_state = join(truth_values, "") # Get the truth specific about the latest crossing if !ismissing(upcrossing) - truth_values[condition_idx] = upcrossing ? "U" : "D" + truth_value_idx = + condition_idx - sum( + length(vec) for + vec in discrete_control.condition_value[1:(where_node_id.start - 1)]; + init = 0, + ) + truth_values[truth_value_idx] = upcrossing ? "U" : "D" end - truth_state_crossing_specific = join(truth_values[condition_ids], "") + truth_state_crossing_specific = join(truth_values, "") # What the local control state should be control_state_new = diff --git a/core/src/parameter.jl b/core/src/parameter.jl index 231f0a9c4..e898280ff 100644 --- a/core/src/parameter.jl +++ b/core/src/parameter.jl @@ -457,12 +457,15 @@ record: Namedtuple with discrete control information for results """ struct DiscreteControl <: AbstractParameterNode node_id::Vector{NodeID} + # Definition of compound variables listen_node_id::Vector{Vector{NodeID}} variable::Vector{Vector{String}} weight::Vector{Vector{Float64}} look_ahead::Vector{Vector{Float64}} - greater_than::Vector{Float64} - condition_value::Vector{Bool} + # Definition of conditions (one or more greater_than per compound variable) + greater_than::Vector{Vector{Float64}} + condition_value::Vector{BitVector} + # Definition of logic control_state::Dict{NodeID, Tuple{String, Float64}} logic_mapping::Dict{Tuple{NodeID, String}, String} record::@NamedTuple{ diff --git a/core/src/read.jl b/core/src/read.jl index a0aec5fda..62a7ddc89 100644 --- a/core/src/read.jl +++ b/core/src/read.jl @@ -542,33 +542,79 @@ function Basin(db::DB, config::Config, chunk_sizes::Vector{Int})::Basin ) end -function DiscreteControl(db::DB, config::Config)::DiscreteControl - compound_variable = load_structvector(db, config, DiscreteControlVariableV1) - condition = load_structvector(db, config, DiscreteControlConditionV1) - +function parse_variables_and_conditions(compound_variable, condition) node_id = NodeID[] listen_node_id = Vector{NodeID}[] variable = Vector{String}[] weight = Vector{Float64}[] look_ahead = Vector{Float64}[] + greater_than = Vector{Float64}[] + condition_value = BitVector[] + errors = false for id in unique(condition.node_id) - group_id = filter(row -> row.node_id == id, compound_variable) - for group_variable in - StructVector.(IterTools.groupby(row -> row.compound_variable_id, group_id)) - first_row = first(group_variable) - push!(node_id, NodeID(NodeType.DiscreteControl, first_row.node_id)) - push!( - listen_node_id, - NodeID.(group_variable.listen_node_type, group_variable.listen_node_id), + condition_group_id = filter(row -> row.node_id == id, condition) + variable_group_id = filter(row -> row.node_id == id, compound_variable) + for compound_variable_id in unique(condition_group_id.compound_variable_id) + condition_group_variable = filter( + row -> row.compound_variable_id == compound_variable_id, + condition_group_id, ) - push!(variable, group_variable.variable) - push!(weight, coalesce.(group_variable.weight, 1.0)) - push!(look_ahead, coalesce.(group_variable.look_ahead, 0.0)) + variable_group_variable = filter( + row -> row.compound_variable_id == compound_variable_id, + variable_group_id, + ) + discrete_control_id = NodeID(NodeType.DiscreteControl, id) + if isempty(variable_group_variable) + errors = true + @error "compound_variable_id $compound_variable_id for $discrete_control_id in condition table but not in variable table" + else + push!(node_id, discrete_control_id) + push!( + listen_node_id, + NodeID.( + variable_group_variable.listen_node_type, + variable_group_variable.listen_node_id, + ), + ) + push!(variable, variable_group_variable.variable) + push!(weight, coalesce.(variable_group_variable.weight, 1.0)) + push!(look_ahead, coalesce.(variable_group_variable.look_ahead, 0.0)) + push!(greater_than, condition_group_variable.greater_than) + push!( + condition_value, + BitVector(zeros(length(condition_group_variable.greater_than))), + ) + end end end + return node_id, + listen_node_id, + variable, + weight, + look_ahead, + greater_than, + condition_value, + !errors +end + +function DiscreteControl(db::DB, config::Config)::DiscreteControl + condition = load_structvector(db, config, DiscreteControlConditionV1) + compound_variable = load_structvector(db, config, DiscreteControlVariableV1) + + node_id, + listen_node_id, + variable, + weight, + look_ahead, + greater_than, + condition_value, + valid = parse_variables_and_conditions(compound_variable, condition) + + if !valid + error("Problems encountered when parsing DiscreteControl variables and conditions.") + end - condition_value = fill(false, length(condition.node_id)) control_state::Dict{NodeID, Tuple{String, Float64}} = Dict() rows = execute(db, "SELECT from_node_id, edge_type FROM Edge ORDER BY fid") @@ -603,7 +649,7 @@ function DiscreteControl(db::DB, config::Config)::DiscreteControl variable, weight, look_ahead, - condition.greater_than, + greater_than, condition_value, control_state, logic_mapping, diff --git a/core/src/validation.jl b/core/src/validation.jl index 580b4d5ba..9d6a621d2 100644 --- a/core/src/validation.jl +++ b/core/src/validation.jl @@ -102,12 +102,9 @@ sort_by_id_state_level(row) = (row.node_id, row.control_state, row.level) sort_by_priority(row) = (row.node_id, row.priority) sort_by_priority_time(row) = (row.node_id, row.priority, row.time) sort_by_subgrid_level(row) = (row.subgrid_id, row.basin_level) -sort_by_condition(row) = - (row.node_id, row.listen_node_type, row.listen_node_id, row.variable, row.greater_than) sort_by_variable(row) = (row.node_id, row.listen_node_type, row.listen_node_id, row.variable) -sort_by_id_greater_than(row) = (row.node_id, row.greater_than) -sort_by_truth_state(row) = (row.node_id, row.truth_state) +sort_by_condition(row) = (row.node_id, row.compound_variable_id, row.greater_than) # get the right sort by function given the Schema, with sort_by_id as the default sort_by_function(table::StructVector{<:Legolas.AbstractRecord}) = sort_by_id @@ -117,7 +114,7 @@ sort_by_function(table::StructVector{UserDemandStaticV1}) = sort_by_priority sort_by_function(table::StructVector{UserDemandTimeV1}) = sort_by_priority_time sort_by_function(table::StructVector{BasinSubgridV1}) = sort_by_subgrid_level sort_by_function(table::StructVector{DiscreteControlVariableV1}) = sort_by_variable -sort_by_function(table::StructVector{DiscreteControlConditionV1}) = sort_by_id_greater_than +sort_by_function(table::StructVector{DiscreteControlConditionV1}) = sort_by_condition const TimeSchemas = Union{ BasinTimeV1, @@ -499,7 +496,8 @@ Check: """ function valid_discrete_control(p::Parameters, config::Config)::Bool (; discrete_control, graph) = p - (; node_id, logic_mapping, look_ahead, variable, listen_node_id) = discrete_control + (; node_id, logic_mapping, look_ahead, variable, listen_node_id, greater_than) = + discrete_control t_end = seconds_since(config.endtime, config.starttime) errors = false @@ -512,7 +510,7 @@ function valid_discrete_control(p::Parameters, config::Config)::Bool truth_states_wrong_length = String[] # The number of conditions of this DiscreteControl node - n_conditions = length(searchsorted(node_id, id)) + n_conditions = sum(length(greater_than[i]) for i in searchsorted(node_id, id)) for (key, control_state) in logic_mapping id_, truth_state = key diff --git a/core/test/control_test.jl b/core/test/control_test.jl index e0913ba1d..21bcc653c 100644 --- a/core/test/control_test.jl +++ b/core/test/control_test.jl @@ -37,11 +37,11 @@ # Control times t_1 = discrete_control.record.time[3] t_1_index = findfirst(>=(t_1), t) - @test level[1, t_1_index] <= discrete_control.greater_than[1] + @test level[1, t_1_index] <= discrete_control.greater_than[1][1] t_2 = discrete_control.record.time[4] t_2_index = findfirst(>=(t_2), t) - @test level[2, t_2_index] >= discrete_control.greater_than[2] + @test level[2, t_2_index] >= discrete_control.greater_than[2][1] flow = get_tmp(graph[].flow, 0) @test all(iszero, flow) @@ -60,7 +60,7 @@ end t_control = discrete_control.record.time[2] t_control_index = searchsortedfirst(t, t_control) - greater_than = discrete_control.greater_than[1] + greater_than = discrete_control.greater_than[1][1] flow_t_control = flow_boundary.flow_rate[1](t_control) flow_t_control_ahead = flow_boundary.flow_rate[1](t_control + Δt) @@ -84,7 +84,7 @@ end t_control = discrete_control.record.time[2] t_control_index = searchsortedfirst(t, t_control) - greater_than = discrete_control.greater_than[1] + greater_than = discrete_control.greater_than[1][1] level_t_control = level_boundary.level[1](t_control) level_t_control_ahead = level_boundary.level[1](t_control + Δt) @@ -167,8 +167,8 @@ end t_in = discrete_control.record.time[3] t_none_2 = discrete_control.record.time[4] - level_min = greater_than[1] - setpoint = greater_than[2] + level_min = greater_than[1][1] + setpoint = greater_than[1][2] t_1_none_index = findfirst(>=(t_none_1), t) t_in_index = findfirst(>=(t_in), t) diff --git a/python/ribasim/ribasim/config.py b/python/ribasim/ribasim/config.py index d837fc8fc..30408fe05 100644 --- a/python/ribasim/ribasim/config.py +++ b/python/ribasim/ribasim/config.py @@ -304,13 +304,14 @@ class DiscreteControl(MultiNodeModel): json_schema_extra={ "sort_keys": [ "node_id", + "compound_variable_id", "greater_than", ] }, ) logic: TableModel[DiscreteControlLogicSchema] = Field( default_factory=TableModel[DiscreteControlLogicSchema], - json_schema_extra={"sort_keys": ["node_id"]}, + json_schema_extra={"sort_keys": ["node_id", "truth_state"]}, ) diff --git a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py index d8d911674..650e73876 100644 --- a/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py +++ b/python/ribasim_testmodels/ribasim_testmodels/discrete_control.py @@ -299,6 +299,7 @@ def tabulated_rating_curve_control_model() -> Model: ), discrete_control.Condition( greater_than=[0.5], + compound_variable_id=1, ), discrete_control.Logic( truth_state=["T", "F"], control_state=["low", "high"] From a65f6f69266caa7749715d5d495898b5bb701b57 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Tue, 16 Apr 2024 11:57:32 +0200 Subject: [PATCH 17/22] Small fixes --- core/src/callback.jl | 2 +- python/ribasim/tests/test_io.py | 4 +--- 2 files changed, 2 insertions(+), 4 deletions(-) diff --git a/core/src/callback.jl b/core/src/callback.jl index 9dbb48acd..05b49dda4 100644 --- a/core/src/callback.jl +++ b/core/src/callback.jl @@ -9,7 +9,7 @@ function set_initial_discrete_controlled_parameters!( (; p) = integrator (; discrete_control) = p - n_conditions = sum(length(vec) for vec in discrete_control.condition_value) + n_conditions = sum(length(vec) for vec in discrete_control.condition_value; init = 0) condition_diffs = zeros(Float64, n_conditions) discrete_control_condition(condition_diffs, storage0, integrator.t, integrator) diff --git a/python/ribasim/tests/test_io.py b/python/ribasim/tests/test_io.py index e9b724eaa..e6507d77e 100644 --- a/python/ribasim/tests/test_io.py +++ b/python/ribasim/tests/test_io.py @@ -102,9 +102,7 @@ def test_sort(level_setpoint_with_minmax, tmp_path): assert table.df.iloc[0]["greater_than"] == 15.0 assert table._sort_keys == [ "node_id", - "listen_node_type", - "listen_node_id", - "variable", + "compound_variable_id", "greater_than", ] table.sort() From 9765dc2d977947b0c3104eed8a1969a4ffe7ddbd Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Tue, 16 Apr 2024 12:47:46 +0200 Subject: [PATCH 18/22] Fix examples notebook --- docs/python/examples.ipynb | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/docs/python/examples.ipynb b/docs/python/examples.ipynb index f78357420..d225579d3 100644 --- a/docs/python/examples.ipynb +++ b/docs/python/examples.ipynb @@ -493,10 +493,13 @@ " Node(7, Point(1.0, 0.0)),\n", " [\n", " discrete_control.Variable(\n", - " listen_node_id=[1], listen_node_type=\"Basin\", variable=\"level\"\n", + " listen_node_id=[1],\n", + " listen_node_type=\"Basin\",\n", + " variable=\"level\",\n", + " compound_variable_id=1,\n", " ),\n", " discrete_control.Condition(\n", - " greater_than=[5.0, 10.0, 15.0],\n", + " greater_than=[5.0, 10.0, 15.0], compound_variable_id=1\n", " ),\n", " discrete_control.Logic(\n", " truth_state=[\"FFF\", \"U**\", \"T*F\", \"**D\", \"TTT\"],\n", @@ -1183,10 +1186,12 @@ "model.discrete_control.add(\n", " Node(11, Point(4.5, 0.25), subnetwork_id=1),\n", " [\n", - " discrete_control.Condition(\n", + " discrete_control.Variable(\n", " listen_node_id=[5],\n", " listen_node_type=[\"Basin\"],\n", " variable=[\"level\"],\n", + " ),\n", + " discrete_control.Condition(\n", " greater_than=[0.52],\n", " ),\n", " discrete_control.Logic(\n", From 267d412a48eb9aab3ee5887c6213be5858941ad8 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Tue, 16 Apr 2024 13:12:26 +0200 Subject: [PATCH 19/22] Update docstrings --- core/src/callback.jl | 24 ++++++------------------ core/src/parameter.jl | 15 +++++++-------- core/src/read.jl | 2 ++ core/src/util.jl | 16 ++++++++++++++++ 4 files changed, 31 insertions(+), 26 deletions(-) diff --git a/core/src/callback.jl b/core/src/callback.jl index 05b49dda4..fb1918831 100644 --- a/core/src/callback.jl +++ b/core/src/callback.jl @@ -13,11 +13,13 @@ function set_initial_discrete_controlled_parameters!( condition_diffs = zeros(Float64, n_conditions) discrete_control_condition(condition_diffs, storage0, integrator.t, integrator) + # Set the discrete control value (bool) per compound variable idx_start = 1 - for (i, vec) in enumerate(discrete_control.condition_value) + for (compound_variable_idx, vec) in enumerate(discrete_control.condition_value) l = length(vec) idx_end = idx_start + l - 1 - discrete_control.condition_value[i] .= (condition_diffs[idx_start:idx_end] .> 0) + discrete_control.condition_value[compound_variable_idx] .= + (condition_diffs[idx_start:idx_end] .> 0) idx_start += l end @@ -202,6 +204,7 @@ function discrete_control_condition(out, u, t, integrator) (; discrete_control) = p condition_idx = 0 + # Loop over compound variables for (listen_node_ids, variables, weights, greater_thans, look_aheads) in zip( discrete_control.listen_node_id, discrete_control.variable, @@ -214,6 +217,7 @@ function discrete_control_condition(out, u, t, integrator) zip(listen_node_ids, variables, weights, look_aheads) value += weight * get_value(p, listen_node_id, variable, look_ahead, u, t) end + # Loop over greater_than values for this compound_variable for greater_than in greater_thans condition_idx += 1 diff = value - greater_than @@ -268,22 +272,6 @@ function get_value( return value end -function get_discrete_control_indices(discrete_control::DiscreteControl, condition_idx::Int) - (; greater_than) = discrete_control - condition_idx_now = 1 - - for (compound_variable_idx, vec) in enumerate(greater_than) - l = length(vec) - - if condition_idx_now + l > condition_idx - greater_than_idx = condition_idx - condition_idx_now + 1 - return compound_variable_idx, greater_than_idx - end - - condition_idx_now += l - end -end - """ An upcrossing means that a condition (always greater than) becomes true. """ diff --git a/core/src/parameter.jl b/core/src/parameter.jl index e898280ff..9ea874d74 100644 --- a/core/src/parameter.jl +++ b/core/src/parameter.jl @@ -443,14 +443,13 @@ struct Terminal <: AbstractParameterNode end """ -node_id: node ID of the DiscreteControl node; these are not unique but repeated - by the amount of conditions of this DiscreteControl node -listen_node_id: the IDs of the nodes being condition on -variable: the names of the variables in the condition -weight: the weight of the variables in the condition -look_ahead: the look ahead of variables in the condition in seconds -greater_than: The threshold value in the condition -condition_value: The current value of each condition +node_id: node ID of the DiscreteControl node per compound variable (can contain repeats) +listen_node_id: the IDs of the nodes being condition on per compound variable +variable: the names of the variables in the condition per compound variable +weight: the weight of the variables in the condition per compound variable +look_ahead: the look ahead of variables in the condition in seconds per compound_variable +greater_than: The threshold values per compound variable +condition_value: The current truth value of each condition per compound_variable per greater_than control_state: Dictionary: node ID => (control state, control state start) logic_mapping: Dictionary: (control node ID, truth state) => control state record: Namedtuple with discrete control information for results diff --git a/core/src/read.jl b/core/src/read.jl index 62a7ddc89..bd07fd17d 100644 --- a/core/src/read.jl +++ b/core/src/read.jl @@ -552,9 +552,11 @@ function parse_variables_and_conditions(compound_variable, condition) condition_value = BitVector[] errors = false + # Loop over unique discrete_control node IDs (on which at least one condition is defined) for id in unique(condition.node_id) condition_group_id = filter(row -> row.node_id == id, condition) variable_group_id = filter(row -> row.node_id == id, compound_variable) + # Loop over compound variables for this node ID for compound_variable_id in unique(condition_group_id.compound_variable_id) condition_group_variable = filter( row -> row.compound_variable_id == compound_variable_id, diff --git a/core/src/util.jl b/core/src/util.jl index 9f28670cf..d84567632 100644 --- a/core/src/util.jl +++ b/core/src/util.jl @@ -709,3 +709,19 @@ function get_influx(basin::Basin, basin_idx::Int)::Float64 return precipitation[basin_idx] - evaporation[basin_idx] + drainage[basin_idx] - infiltration[basin_idx] end + +function get_discrete_control_indices(discrete_control::DiscreteControl, condition_idx::Int) + (; greater_than) = discrete_control + condition_idx_now = 1 + + for (compound_variable_idx, vec) in enumerate(greater_than) + l = length(vec) + + if condition_idx_now + l > condition_idx + greater_than_idx = condition_idx - condition_idx_now + 1 + return compound_variable_idx, greater_than_idx + end + + condition_idx_now += l + end +end From c662fe7f37ea30b32f17fee2f143250ca34c1819 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Tue, 16 Apr 2024 13:28:39 +0200 Subject: [PATCH 20/22] Update usage.qmd --- core/src/schema.jl | 6 +++--- docs/core/usage.qmd | 44 +++++++++++++++++--------------------- ribasim_qgis/core/nodes.py | 2 ++ 3 files changed, 25 insertions(+), 27 deletions(-) diff --git a/core/src/schema.jl b/core/src/schema.jl index 60eb87ac8..0875e48c5 100644 --- a/core/src/schema.jl +++ b/core/src/schema.jl @@ -186,9 +186,9 @@ end @version DiscreteControlVariableV1 begin node_id::Int32 - compound_variable_id::Int + compound_variable_id::Int32 listen_node_type::String - listen_node_id::Int + listen_node_id::Int32 variable::String weight::Union{Missing, Float64} look_ahead::Union{Missing, Float64} @@ -196,7 +196,7 @@ end @version DiscreteControlConditionV1 begin node_id::Int32 - compound_variable_id::Int + compound_variable_id::Int32 greater_than::Float64 end diff --git a/docs/core/usage.qmd b/docs/core/usage.qmd index aa7c5f73a..a736c41cd 100644 --- a/docs/core/usage.qmd +++ b/docs/core/usage.qmd @@ -204,8 +204,8 @@ name it must have in the database if it is stored there. - Terminal: Water sink without state or properties - `Terminal / static`: - (only node IDs) - DiscreteControl: Set parameters of other nodes based on model state conditions (e.g. basin level) - - `DiscreteControl / compoundvariable` : Linear combinations of variables, e.g. averages or differences of levels - - `DisceteControl / condition`: Conditions of the form 'the level in the basin with node id `n` is bigger than 2.0 m' + - `DiscreteControl / variable` : Linear combinations of variables, e.g. averages or differences of levels + - `DisceteControl / condition`: Conditions of the form 'the compound variable with ID `n` is bigger than 2.0' - `DisceteControl / logic`: Translates the truth value of a set of conditions to parameter values for a controlled node - PidControl: Controls the level in a basin by continuously controlling the flow rate of a connected pump or outlet. See also [Wikipedia](https://en.wikipedia.org/wiki/PID_controller) and [PID controller in equations](equations.qmd#sec-PID). - `PidControl / static`: The proportional, integral, and derivative parameters, the target value and which basin should be controlled @@ -626,39 +626,35 @@ node_id | Int32 | - | sorted DiscreteControl is implemented based on [VectorContinuousCallback](https://docs.sciml.ai/DiffEqDocs/stable/features/callback_functions/#VectorContinuousCallback). -## DiscreteControl / compoundvariable +## DiscreteControl / variable The compound variable schema defines linear combinations of variables which can be used in conditions. This means that -this schema defines new variables with the given name that look like +this schema defines new variables with the given `compound_variable_id` that look like $$ \text{weight}_1 * \text{variable}_1 + \text{weight}_2 * \text{variable}_2 + \ldots, $$ which can be for instance an average or a difference of variables. If a variable comes from a time-series, a look ahead $\Delta t$ can be supplied. -column | type | unit | restriction ------------------ | -------- | ------- | ----------- -node_id | Int32 | - | sorted -name | String | - | sorted per node_id -listen_node_type | String | - | known node type -listen_node_id | Int32 | - | sorted per node_id -variable | String | - | must be "level" or "flow_rate", sorted per listen_node_id -weight | Float64 | - | - -look_ahead | Float64 | $s$ | Only on transient boundary conditions, non-negative (optional, default 0). +column | type | unit | restriction +-------------------- | -------- | ------- | ----------- +node_id | Int32 | - | sorted +compound_variable_id | String | - | sorted per node_id +listen_node_type | String | - | known node type +listen_node_id | Int32 | - | sorted per node_id +variable | String | - | must be "level" or "flow_rate", sorted per listen_node_id +weight | Float64 | - | (optional, default 1.0) +look_ahead | Float64 | $s$ | Only on transient boundary conditions, non-negative (optional, default 0.0). ## DiscreteControl / condition -The condition schema defines conditions of the form 'the discrete_control node with this node id listens to whether the given variable of the node with the given listen node id is greater than the given value'. If the condition variable comes from a time-series, a look ahead $\Delta t$ can be supplied. -Compound variables from the schema above can also be used, in which case the listen node ID and look ahead in this schema are ignored. - -column | type | unit | restriction ------------------ | -------- | ------- | ----------- -node_id | Int32 | - | sorted -listen_node_type | String | - | known node type or "compound" -listen_node_id | Int32 | - | sorted per node_id -variable | String | - | must be "level", "flow_rate" or a compound variable name, sorted per listen_node_id -greater_than | Float64 | various | sorted per variable -look_ahead | Float64 | $s$ | Only on transient boundary conditions, non-negative (optional, default 0) +The condition schema defines conditions of the form 'the `discrete_control` node with this `node_id`listens to whether the variable given by the `node_id` and `compound_variable_id` is greater than `greater_than`'. If the condition variable comes from a time-series, a look ahead $\Delta t$ can be supplied. Multiple conditions with different `greater_than` values can be defined on the same compound_variable. + +column | type | unit | restriction +-------------------- | -------- | ------- | ----------- +node_id | Int32 | - | sorted +compound_variable_id | Int32 | - | - +greater_than | Float64 | various | sorted per variable ## DiscreteControl / logic diff --git a/ribasim_qgis/core/nodes.py b/ribasim_qgis/core/nodes.py index 5c0393a81..97682e780 100644 --- a/ribasim_qgis/core/nodes.py +++ b/ribasim_qgis/core/nodes.py @@ -704,6 +704,7 @@ def geometry_type(cls) -> str: def attributes(cls) -> list[QgsField]: return [ QgsField("node_id", QVariant.Int), + QgsField("compound_variable_id", QVariant.Int), QgsField("listen_node_type", QVariant.String), QgsField("listen_node_id", QVariant.Int), QgsField("variable", QVariant.String), @@ -725,6 +726,7 @@ def geometry_type(cls) -> str: def attributes(cls) -> list[QgsField]: return [ QgsField("node_id", QVariant.Int), + QgsField("compound_variable_id", QVariant.Int), QgsField("greater_than", QVariant.Double), ] From 3707f5083438ad477c7b3f7688da8bb0a6f7a973 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Tue, 16 Apr 2024 18:24:00 +0200 Subject: [PATCH 21/22] Fix examples notebook? --- docs/python/examples.ipynb | 2 ++ 1 file changed, 2 insertions(+) diff --git a/docs/python/examples.ipynb b/docs/python/examples.ipynb index d225579d3..32bbb556b 100644 --- a/docs/python/examples.ipynb +++ b/docs/python/examples.ipynb @@ -1187,11 +1187,13 @@ " Node(11, Point(4.5, 0.25), subnetwork_id=1),\n", " [\n", " discrete_control.Variable(\n", + " compound_variable_id=1,\n", " listen_node_id=[5],\n", " listen_node_type=[\"Basin\"],\n", " variable=[\"level\"],\n", " ),\n", " discrete_control.Condition(\n", + " compound_variable_id=1,\n", " greater_than=[0.52],\n", " ),\n", " discrete_control.Logic(\n", From a140316ea9c8d0cfe70d6d19004fa175fc513cf9 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Tue, 16 Apr 2024 19:24:01 +0200 Subject: [PATCH 22/22] Comments adressed --- core/src/validation.jl | 2 -- docs/core/usage.qmd | 4 ++-- python/ribasim/ribasim/model.py | 1 - 3 files changed, 2 insertions(+), 5 deletions(-) diff --git a/core/src/validation.jl b/core/src/validation.jl index cf42c6a4c..ae039a561 100644 --- a/core/src/validation.jl +++ b/core/src/validation.jl @@ -574,8 +574,6 @@ function valid_discrete_control(p::Parameters, config::Config)::Bool for (Δt, var, node_id) in zip(look_aheads, variables, listen_node_ids) if !iszero(Δt) node_type = node_id.type - # TODO: If more transient listen variables must be supported, this validation must be more specific - # (e.g. for some node some variables are transient, some not). if node_type ∉ [NodeType.FlowBoundary, NodeType.LevelBoundary] errors = true @error "Look ahead supplied for non-timeseries listen variable '$var' from listen node $node_id." diff --git a/docs/core/usage.qmd b/docs/core/usage.qmd index a736c41cd..88374f8ef 100644 --- a/docs/core/usage.qmd +++ b/docs/core/usage.qmd @@ -639,7 +639,7 @@ which can be for instance an average or a difference of variables. If a variable column | type | unit | restriction -------------------- | -------- | ------- | ----------- node_id | Int32 | - | sorted -compound_variable_id | String | - | sorted per node_id +compound_variable_id | Int32 | - | sorted per node_id listen_node_type | String | - | known node type listen_node_id | Int32 | - | sorted per node_id variable | String | - | must be "level" or "flow_rate", sorted per listen_node_id @@ -648,7 +648,7 @@ look_ahead | Float64 | $s$ | Only on transient boundary condition ## DiscreteControl / condition -The condition schema defines conditions of the form 'the `discrete_control` node with this `node_id`listens to whether the variable given by the `node_id` and `compound_variable_id` is greater than `greater_than`'. If the condition variable comes from a time-series, a look ahead $\Delta t$ can be supplied. Multiple conditions with different `greater_than` values can be defined on the same compound_variable. +The condition schema defines conditions of the form 'the `discrete_control` node with this `node_id`listens to whether the variable given by the `node_id` and `compound_variable_id` is greater than `greater_than`'. Multiple conditions with different `greater_than` values can be defined on the same compound_variable. column | type | unit | restriction -------------------- | -------- | ------- | ----------- diff --git a/python/ribasim/ribasim/model.py b/python/ribasim/ribasim/model.py index 966a65888..2b4f195ef 100644 --- a/python/ribasim/ribasim/model.py +++ b/python/ribasim/ribasim/model.py @@ -326,7 +326,6 @@ def plot_control_listen(self, ax): to_add = df_variable[ ["node_id", "listen_node_id", "listen_node_type"] ].drop_duplicates() - to_add = to_add[to_add["listen_node_type"] != "compound"] to_add.columns = ["control_node_id", "listen_node_id", "listen_node_type"] to_add["control_node_type"] = "DiscreteControl" df_listen_edge = pd.concat([df_listen_edge, to_add])