The WaterNetwork component uses the America's Water Network to describe how downstream flows depend on upstream flows and withdrawals.
The America's Water Network is described in detail in its working paper. The data contained in the network and its loading are described in the sister waternet notebook.
The component definition is very simple:
In [1]:
using Mimi
@defcomp WaterNetwork begin
gauges = Index()
# External
added = Parameter(index=[gauges, time], unit="1000 m^3") # Water added at node from runoff
removed = Parameter(index=[gauges, time], unit="1000 m^3") # Water removed from node
returned = Parameter(index=[gauges, time], unit="1000 m^3") # Water returns to a node from canals
inflows = Variable(index=[gauges, time], unit="1000 m^3") # Sum of upstream outflows
outflows = Variable(index=[gauges, time], unit="1000 m^3") # inflow + added - removed + returned
end
Out[1]:
All of these values are defined for each gauge, or node on the river network. They present a water balance:
$$outflows_i(t) = \sum_{j \in N(i)} inflows_i(t) + added_i(t) - removed_i(t) + returned_i(t)$$
Where $N(i)$ are the neighbors of gauge $i$, where neighbors are defined as the upstream sources of flow.
The added parameter is set to the runoff calculated by VIC and loaded in weather.jl.
In optimization mode, both the returned and removed parameters are set by the optimization to satisfy allocation demands.
In simulation mode,
removed parameter is set to the contents of data/cache/counties/extraction/withdrawals.jld, if available.returned parameter is set to the contents of data/cache/counties/extraction/returns.jld, if available.The removed and returned variables are generated by optimize-surface.jl, which finds plausible extractions to satisfy USGS demands.
The logic for optimization is complicated by the process of matching canal-based withdrawals to gauge-based flows. However, aside from this, the entire network follows two simple principes.
The following shows these gauge-specific allotments, and some interesting caveats.
First, the allotments take a while to calculate (which they are done by constraintoffset_waternetwork_outflows(m); in optimize.jl if need be), so they are cached in the file data/partialhouse2.jld.
In [3]:
using OptiMimi
cwro = deserialize(open("../data/partialhouse2.jld", "r"))
Out[3]:
This is a Linear Programming Hall file, used for constraint vectors.
Additional informatino on the gauges is stored in the waternet.RData file, generated outside of Julia.
In [4]:
using DataFrames
stations = read_rda("../data/waternet.RData", convertdataframes=true)["network"]
Out[4]:
Unfortunately, the allotment entries in the optimization problem are neither in the same order (they're in vertex order), nor have all of the entries that are in the water network.
In [5]:
println(nrow(stations))
println(length(cwro.f))
We need to use the wateridverts dictionary, which maps from gauge IDs to vertices in the network.
In [7]:
using Graphs
wateridverts = deserialize(open("../data/cache/wateridverts.jld", "r"))
wateridverts
#waternet = deserialize(open("../data/cache/waternet.jld", "r"))
Out[7]:
And now we can collect all of the information to produce a map.
In [8]:
lats = []
lons = []
runoff = []
for ii in 1:nrow(stations)
gaugeid = "$(stations[ii, :collection]).$(stations[ii, :colid])"
if gaugeid in keys(wateridverts)
runoff = [runoff; cwro.f[vertex_index(wateridverts[gaugeid])]]
lats = [lats; stations[ii, :lat]]
lons = [lons; stations[ii, :lon]]
end
end
runoff[isna(runoff)] = 0
Out[8]:
In [44]:
using RCall
cex = convert(Vector{Float64}, runoff) / 1e9
R"library(maps)"
R"map('usa', xlim=c(-140, -60), ylim=c(20, 50))"
R"points($lons, $lats, cex=$cex)"
Out[44]: