Skip to content

Commit

Permalink
Precalculate resistance neighbors
Browse files Browse the repository at this point in the history
  • Loading branch information
@visr
visr committed Apr 29, 2024
1 parent 00b9896 commit 225d912
Show file tree
Hide file tree
Showing 3 changed files with 42 additions and 30 deletions.
20 changes: 12 additions & 8 deletions core/src/parameter.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -241,19 +241,18 @@ struct TabulatedRatingCurve{C} <: AbstractParameterNode
end

"""
Requirements:
* from: must be (Basin,) node
* to: must be (Basin,) node
node_id: node ID of the LinearResistance node
inflow_id: node ID across the incoming flow edge
outflow_id: node ID across the outgoing flow edge
active: whether this node is active and thus contributes flows
resistance: the resistance to flow; `Q_unlimited = Δh/resistance`
max_flow_rate: the maximum flow rate allowed through the node; `Q = clamp(Q_unlimited, -max_flow_rate, max_flow_rate)`
control_mapping: dictionary from (node_id, control_state) to resistance and/or active state
"""
struct LinearResistance <: AbstractParameterNode
node_id::Vector{NodeID}
inflow_id::Vector{NodeID}
outflow_id::Vector{NodeID}
active::BitVector
resistance::Vector{Float64}
max_flow_rate::Vector{Float64}
Expand All @@ -263,8 +262,11 @@ end
"""
This is a simple Manning-Gauckler reach connection.
* Length describes the reach length.
* roughness describes Manning's n in (SI units).
node_id: node ID of the ManningResistance node
inflow_id: node ID across the incoming flow edge
outflow_id: node ID across the outgoing flow edge
length: reach length
manning_n: roughness; Manning's n in (SI units).
The profile is described by a trapezoid:
Expand All @@ -288,13 +290,15 @@ Requirements:
* from: must be (Basin,) node
* to: must be (Basin,) node
* length > 0
* roughess > 0
* manning_n > 0
* profile_width >= 0
* profile_slope >= 0
* (profile_width == 0) xor (profile_slope == 0)
"""
struct ManningResistance <: AbstractParameterNode
node_id::Vector{NodeID}
inflow_id::Vector{NodeID}
outflow_id::Vector{NodeID}
active::BitVector
length::Vector{Float64}
manning_n::Vector{Float64}
Expand Down
20 changes: 14 additions & 6 deletions core/src/read.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -230,7 +230,7 @@ function initialize_allocation!(p::Parameters, config::Config)::Nothing
return nothing
end

function LinearResistance(db::DB, config::Config)::LinearResistance
function LinearResistance(db::DB, config::Config, graph::MetaGraph)::LinearResistance
static = load_structvector(db, config, LinearResistanceStaticV1)
defaults = (; max_flow_rate = Inf, active = true)
parsed_parameters, valid =
Expand All @@ -242,8 +242,12 @@ function LinearResistance(db::DB, config::Config)::LinearResistance
)
end

node_id = NodeID.(NodeType.LinearResistance, parsed_parameters.node_id)

return LinearResistance(
NodeID.(NodeType.LinearResistance, parsed_parameters.node_id),
node_id,
inflow_id.(Ref(graph), node_id),
outflow_id.(Ref(graph), node_id),
BitVector(parsed_parameters.active),
parsed_parameters.resistance,
parsed_parameters.max_flow_rate,
Expand Down Expand Up @@ -321,7 +325,7 @@ function TabulatedRatingCurve(db::DB, config::Config)::TabulatedRatingCurve
return TabulatedRatingCurve(node_ids, active, interpolations, time, control_mapping)
end

function ManningResistance(db::DB, config::Config)::ManningResistance
function ManningResistance(db::DB, config::Config, graph::MetaGraph)::ManningResistance
static = load_structvector(db, config, ManningResistanceStaticV1)
parsed_parameters, valid =
parse_static_and_time(db, config, "ManningResistance"; static)
Expand All @@ -330,8 +334,12 @@ function ManningResistance(db::DB, config::Config)::ManningResistance
error("Errors occurred when parsing ManningResistance data.")
end

node_id = NodeID.(NodeType.ManningResistance, parsed_parameters.node_id)

return ManningResistance(
NodeID.(NodeType.ManningResistance, parsed_parameters.node_id),
node_id,
inflow_id.(Ref(graph), node_id),
outflow_id.(Ref(graph), node_id),
BitVector(parsed_parameters.active),
parsed_parameters.length,
parsed_parameters.manning_n,
Expand Down Expand Up @@ -1029,8 +1037,8 @@ function Parameters(db::DB, config::Config)::Parameters
error("Invalid edge(s) found.")
end

linear_resistance = LinearResistance(db, config)
manning_resistance = ManningResistance(db, config)
linear_resistance = LinearResistance(db, config, graph)
manning_resistance = ManningResistance(db, config, graph)
tabulated_rating_curve = TabulatedRatingCurve(db, config)
fractional_flow = FractionalFlow(db, config)
level_boundary = LevelBoundary(db, config)
Expand Down
32 changes: 16 additions & 16 deletions core/src/solve.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -336,24 +336,24 @@ function formulate_flow!(
(; graph) = p
(; node_id, active, resistance, max_flow_rate) = linear_resistance
for (i, id) in enumerate(node_id)
basin_a_id = inflow_id(graph, id)
basin_b_id = outflow_id(graph, id)
inflow_id = linear_resistance.inflow_id[i]
outflow_id = linear_resistance.outflow_id[i]

if active[i]
h_a = get_level(p, basin_a_id, t; storage)
h_b = get_level(p, basin_b_id, t; storage)
h_a = get_level(p, inflow_id, t; storage)
h_b = get_level(p, outflow_id, t; storage)
q_unlimited = (h_a - h_b) / resistance[i]
q = clamp(q_unlimited, -max_flow_rate[i], max_flow_rate[i])

# add reduction_factor on highest level
if q > 0
q *= low_storage_factor(storage, p.basin.node_id, basin_a_id, 10.0)
q *= low_storage_factor(storage, p.basin.node_id, inflow_id, 10.0)
else
q *= low_storage_factor(storage, p.basin.node_id, basin_b_id, 10.0)
q *= low_storage_factor(storage, p.basin.node_id, outflow_id, 10.0)
end

set_flow!(graph, basin_a_id, id, q)
set_flow!(graph, id, basin_b_id, q)
set_flow!(graph, inflow_id, id, q)
set_flow!(graph, id, outflow_id, q)
end
end
return nothing
Expand Down Expand Up @@ -438,17 +438,17 @@ function formulate_flow!(
(; node_id, active, length, manning_n, profile_width, profile_slope) =
manning_resistance
for (i, id) in enumerate(node_id)
basin_a_id = inflow_id(graph, id)
basin_b_id = outflow_id(graph, id)
inflow_id = manning_resistance.inflow_id[i]
outflow_id = manning_resistance.outflow_id[i]

if !active[i]
continue
end

h_a = get_level(p, basin_a_id, t; storage)
h_b = get_level(p, basin_b_id, t; storage)
bottom_a = basin_bottom(basin, basin_a_id)
bottom_b = basin_bottom(basin, basin_b_id)
h_a = get_level(p, inflow_id, t; storage)
h_b = get_level(p, outflow_id, t; storage)
bottom_a = basin_bottom(basin, inflow_id)
bottom_b = basin_bottom(basin, outflow_id)
slope = profile_slope[i]
width = profile_width[i]
n = manning_n[i]
Expand Down Expand Up @@ -478,8 +478,8 @@ function formulate_flow!(

q = q_sign * A / n * R_h^(2 / 3) * sqrt(Δh / L * 2 / π * atan(k * Δh) + eps)

set_flow!(graph, basin_a_id, id, q)
set_flow!(graph, id, basin_b_id, q)
set_flow!(graph, inflow_id, id, q)
set_flow!(graph, id, outflow_id, q)
end
return nothing
end
Expand Down

0 comments on commit 225d912

Please sign in to comment.