Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Make LevelDemand min_level and max_level optional #1430

Merged
merged 3 commits into from
Apr 29, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 3 additions & 2 deletions core/src/read.jl
Original file line number Diff line number Diff line change
Expand Up @@ -93,7 +93,7 @@ function parse_static_and_time(

errors = false
t_end = seconds_since(config.endtime, config.starttime)
trivial_timespan = [nextfloat(-Inf), prevfloat(Inf)]
trivial_timespan = [0.0, prevfloat(Inf)]
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Without this change we would get a NaN instead of Inf from DataInterpolations.


for (node_idx, node_id) in enumerate(node_ids)
if node_id in static_node_ids
Expand Down Expand Up @@ -812,7 +812,7 @@ function UserDemand(db::DB, config::Config)::UserDemand
realized_bmi = zeros(n_user)
demand = zeros(n_user, n_priority)
demand_reduced = zeros(n_user, n_priority)
trivial_timespan = [nextfloat(-Inf), prevfloat(Inf)]
trivial_timespan = [0.0, prevfloat(Inf)]
demand_itp = [
[LinearInterpolation(zeros(2), trivial_timespan) for i in eachindex(priorities)] for j in eachindex(node_ids)
]
Expand Down Expand Up @@ -877,6 +877,7 @@ function LevelDemand(db::DB, config::Config)::LevelDemand
static,
time,
time_interpolatables = [:min_level, :max_level],
defaults = (; min_level = -Inf, max_level = Inf),
)

if !valid
Expand Down
8 changes: 4 additions & 4 deletions core/src/schema.jl
Original file line number Diff line number Diff line change
Expand Up @@ -250,16 +250,16 @@ end

@version LevelDemandStaticV1 begin
node_id::Int32
min_level::Float64
max_level::Float64
min_level::Union{Missing, Float64}
max_level::Union{Missing, Float64}
priority::Int32
end

@version LevelDemandTimeV1 begin
node_id::Int32
time::DateTime
min_level::Float64
max_level::Float64
min_level::Union{Missing, Float64}
max_level::Union{Missing, Float64}
priority::Int32
end

Expand Down
14 changes: 5 additions & 9 deletions core/src/util.jl
Original file line number Diff line number Diff line change
Expand Up @@ -31,17 +31,10 @@ function get_storage_from_level(basin::Basin, state_idx::Int, level::Float64)::F
storage_discrete = basin.storage[state_idx]
area_discrete = basin.area[state_idx]
level_discrete = basin.level[state_idx]
bottom = first(level_discrete)

if level < bottom
node_id = basin.node_id.values[state_idx]
@error "The level $level of $node_id is lower than the bottom of this basin; $bottom."
return NaN
end

level_lower_index = searchsortedlast(level_discrete, level)

# If the level is equal to the bottom then the storage is 0
# If the level is at or below the bottom then the storage is 0
if level_lower_index == 0
return 0.0
end
Expand Down Expand Up @@ -77,7 +70,10 @@ function get_storages_from_levels(basin::Basin, levels::Vector)::Vector{Float64}

for (i, level) in enumerate(levels)
storage = get_storage_from_level(basin, i, level)
if isnan(storage)
bottom = first(basin.level[i])
node_id = basin.node_id.values[i]
if level < bottom
@error "The initial level ($level) of $node_id is below the bottom ($bottom)."
errors = true
end
storages[i] = storage
Expand Down
10 changes: 10 additions & 0 deletions core/test/allocation_test.jl
Original file line number Diff line number Diff line change
Expand Up @@ -313,6 +313,8 @@ end
end

@testitem "Allocation level control" begin
import JuMP

toml_path = normpath(@__DIR__, "../../generated_testmodels/level_demand/ribasim.toml")
@test ispath(toml_path)
model = Ribasim.run(toml_path)
Expand Down Expand Up @@ -377,6 +379,14 @@ end
stage_6_start_idx = findfirst(stage_6)
u_stage_6(τ) = storage[stage_6_start_idx]
@test storage[stage_6] ≈ u_stage_6.(t[stage_6]) rtol = 1e-4

# Isolated LevelDemand + Basin pair to test optional min_level
problem = allocation.allocation_models[2].problem
@test JuMP.value(only(problem[:F_basin_in])) == 0.0
@test JuMP.value(only(problem[:F_basin_out])) == 0.0
q = JuMP.normalized_rhs(only(problem[:basin_outflow]))
storage_surplus = 1000.0 # Basin #7 is 1000 m2 and 1 m above LevelDemand max_level
@test q ≈ storage_surplus / Δt_allocation
end

@testitem "Flow demand" begin
Expand Down
7 changes: 5 additions & 2 deletions core/test/utils_test.jl
Original file line number Diff line number Diff line change
Expand Up @@ -114,14 +114,17 @@ end
@test length(logger.logs) == 1
@test logger.logs[1].level == Error
@test logger.logs[1].message ==
"The level -1.0 of Basin #1 is lower than the bottom of this basin; 0.0."
"The initial level (-1.0) of Basin #1 is below the bottom (0.0)."

# Converting from storages to levels and back should return the same storages
storages = range(0.0, 2 * storage[end], 50)
levels = [Ribasim.get_area_and_level(basin, 1, s)[2] for s in storages]
storages_ = [Ribasim.get_storage_from_level(basin, 1, l) for l in levels]

@test storages ≈ storages_

# At or below bottom the storage is 0
@test Ribasim.get_storage_from_level(basin, 1, 0.0) == 0.0
@test Ribasim.get_storage_from_level(basin, 1, -1.0) == 0.0
end

@testitem "Expand logic_mapping" begin
Expand Down
11 changes: 7 additions & 4 deletions docs/core/usage.qmd
Original file line number Diff line number Diff line change
Expand Up @@ -479,14 +479,17 @@ min_level | Float64 | $m$ | -

A `LevelDemand` node associates a minimum and a maximum level with connected basins to be used by the allocation algorithm.
Below the minimum level the basin has a demand of the supplied priority,
above the maximum level the basin acts as a source, used by all nodes with demands in order of priority.
The same `LevelDemand` node can be used for basins in different subnetworks.
above the maximum level the basin has a surplus and acts as a source, used by all nodes with demands in order of priority.
The same `LevelDemand` node can be used for Basins in different subnetworks.

Both `min_level` and `max_level` are optional, to be able to handle only the demand or surplus side.
If both are missing, `LevelDemand` won't have any effects on allocation.

column | type | unit | restriction
------------- | ------- | ------------ | -----------
node_id | Int32 | - | sorted
min_level | Float64 | $m$ | -
max_level | Float64 | $m$ | -
min_level | Float64 | $m$ | (optional, default -Inf)
max_level | Float64 | $m$ | (optional, default Inf)
priority | Int32 | - | positive

## LevelDemand / time
Expand Down
8 changes: 4 additions & 4 deletions python/ribasim/ribasim/schemas.py
Original file line number Diff line number Diff line change
Expand Up @@ -115,16 +115,16 @@ class LevelBoundaryTimeSchema(_BaseSchema):

class LevelDemandStaticSchema(_BaseSchema):
node_id: Series[Int32] = pa.Field(nullable=False, default=0)
min_level: Series[float] = pa.Field(nullable=False)
max_level: Series[float] = pa.Field(nullable=False)
min_level: Series[float] = pa.Field(nullable=True)
max_level: Series[float] = pa.Field(nullable=True)
priority: Series[Int32] = pa.Field(nullable=False, default=0)


class LevelDemandTimeSchema(_BaseSchema):
node_id: Series[Int32] = pa.Field(nullable=False, default=0)
time: Series[Timestamp] = pa.Field(nullable=False)
min_level: Series[float] = pa.Field(nullable=False)
max_level: Series[float] = pa.Field(nullable=False)
min_level: Series[float] = pa.Field(nullable=True)
max_level: Series[float] = pa.Field(nullable=True)
priority: Series[Int32] = pa.Field(nullable=False, default=0)


Expand Down
14 changes: 13 additions & 1 deletion python/ribasim_testmodels/ribasim_testmodels/allocation.py
Original file line number Diff line number Diff line change
Expand Up @@ -883,7 +883,7 @@ def subnetworks_with_sources_model() -> Model:


def level_demand_model() -> Model:
"""Small model with a LevelDemand."""
"""Small model with LevelDemand nodes."""

model = Model(
starttime="2020-01-01",
Expand Down Expand Up @@ -922,12 +922,24 @@ def level_demand_model() -> Model:
[basin.Profile(area=1000.0, level=[0.0, 1.0]), basin.State(level=[0.5])],
)

# Isolated LevelDemand + Basin pair to test optional min_level
model.level_demand.add(
Node(6, Point(3, -1), subnetwork_id=3),
[level_demand.Static(max_level=[1.0], priority=1)],
)
model.basin.add(
Node(7, Point(3, 0), subnetwork_id=3),
[basin.Profile(area=1000.0, level=[0.0, 1.0]), basin.State(level=[2.0])],
)

model.edge.add(model.flow_boundary[1], model.basin[2], subnetwork_id=2)
model.edge.add(model.basin[2], model.user_demand[3])
model.edge.add(model.level_demand[4], model.basin[2])
model.edge.add(model.user_demand[3], model.basin[5])
model.edge.add(model.level_demand[4], model.basin[5])

model.edge.add(model.level_demand[6], model.basin[7])

return model


Expand Down
Loading