Model verification¶
In this notebook the model is verified by comparing outcomes of two different examples.
[1]:
import pandas as pd
import networkx as nx
import matplotlib.pyplot as plt
import numpy as np
[2]:
from context import RiverNetwork
from RiverNetwork import RiverNetwork
[3]:
%load_ext autoreload
%autoreload 2
University of Colorado Boulder example¶
This is the most straightforward example that I could find. It is a simple one segment model. \(k\) and \(x\) are estimated to 2.3h and 0.151559154 respectively.
Loading data¶
The structure is loaded as well as the inflow. On the graph we can see the correct \(k\) and \(x\)
[4]:
structure1 = RiverNetwork('../data/single-segment-boulder.xlsx')
structure1.draw(figsize=(3,3))
/home/docs/checkouts/readthedocs.org/user_builds/rna/envs/latest/lib/python3.7/site-packages/networkx/drawing/nx_pylab.py:579: MatplotlibDeprecationWarning:
The iterable function was deprecated in Matplotlib 3.1 and will be removed in 3.3. Use np.iterable instead.
if not cb.iterable(width):
/home/docs/checkouts/readthedocs.org/user_builds/rna/envs/latest/lib/python3.7/site-packages/networkx/drawing/nx_pylab.py:676: MatplotlibDeprecationWarning:
The iterable function was deprecated in Matplotlib 3.1 and will be removed in 3.3. Use np.iterable instead.
if cb.iterable(node_size): # many node sizes
[5]:
inflow = np.array(pd.read_excel('../data/example-inflow-boulder.xlsx').Inflow)
structure1.set_shape('S.1',21,inflow-50)
structure1.draw_Qin(only_sources=True,figsize=(7,4))
[5]:
(<Figure size 672x384 with 1 Axes>,
<matplotlib.axes._subplots.AxesSubplot at 0x7fac3591ffd0>)
Results of flow propagation¶
The flow is calculated for the sink node.
[6]:
structure1.calc_flow_propagation(20)
structure1.draw_Qin(figsize=(7,5))
[6]:
(<Figure size 672x480 with 1 Axes>,
<matplotlib.axes._subplots.AxesSubplot at 0x7fac3584fb38>)
Inflow for node S.1
[7]:
print(structure1.get_Graph().nodes['S.1']['Qin'])
[ 50 100 200 325 450 600 700 780 790 775 750 680 590 500 420 350 300 250
225 200]
Inflow/outflow for sink E.1
[8]:
print(structure1.get_Graph().nodes['E.1']['Qin'])
[ 50. 53.09090167 78.40717029 135.73277234 220.66311342
323.48352768 442.45783876 552.45606895 645.89060264 703.74630509
731.26560713 734.58020111 706.75303009 653.56523866 585.97973616
513.94847587 443.98213789 382.16017306 326.70590002 283.67408341]
These results are then compared to the results of the webpage. The figure shows that the results are almost similar. The small differences in outflow can be explained by rounding errors.
Karahan example¶
In this example the Wilson dataset from the Karahan paper is used. Karahan compared different estimation techniques. It is a simple one segment model. The most interesting difference here is that \(\Delta t\) is not 1 but karahan uses a value of 6. For this dataset the \(x\) is estimated on 0.021 and \(k\) to 29.165. This can also be seen in the figure. The base load is set to 22.
https://onlinelibrary.wiley.com/doi/full/10.1002/cae.20394
Loading data¶
[9]:
structure2 = RiverNetwork('../data/single-segment-karahan.xlsx',dt=6)
structure2.draw(figsize=(3,3))
/home/docs/checkouts/readthedocs.org/user_builds/rna/envs/latest/lib/python3.7/site-packages/networkx/drawing/nx_pylab.py:579: MatplotlibDeprecationWarning:
The iterable function was deprecated in Matplotlib 3.1 and will be removed in 3.3. Use np.iterable instead.
if not cb.iterable(width):
/home/docs/checkouts/readthedocs.org/user_builds/rna/envs/latest/lib/python3.7/site-packages/networkx/drawing/nx_pylab.py:676: MatplotlibDeprecationWarning:
The iterable function was deprecated in Matplotlib 3.1 and will be removed in 3.3. Use np.iterable instead.
if cb.iterable(node_size): # many node sizes
[10]:
inflow = np.array(pd.read_excel('../data/example-inflow-karahan.xlsx').Inflow)
structure2.set_shape('S.1',21,inflow-22)
structure2.draw_Qin(only_sources=True,figsize=(7,4))
[10]:
(<Figure size 672x384 with 1 Axes>,
<matplotlib.axes._subplots.AxesSubplot at 0x7fac357d0c88>)
Results of flow propagation¶
The flow is calculated for the sink node.
[11]:
structure2.calc_flow_propagation(22)
structure2.draw_Qin(figsize=(7,5))
[11]:
(<Figure size 672x480 with 1 Axes>,
<matplotlib.axes._subplots.AxesSubplot at 0x7fac3588c470>)
The input and output data is shown:
Inflow for node S.1
[12]:
structure2.get_Graph().nodes['S.1']['Qin']
[12]:
array([ 22, 23, 35, 71, 103, 111, 109, 100, 86, 71, 59, 47, 39,
32, 28, 24, 22, 21, 20, 19, 19, 18])
Inflow/outflow for node E.1
[13]:
structure2.get_Graph().nodes['E.1']['Qin']
[13]:
array([22. , 21.86603704, 20.52301889, 19.07762521, 26.90355864,
43.58406737, 59.57916704, 72.31400781, 80.64823559, 83.90617013,
82.50290056, 78.6275652 , 72.32100985, 65.4854351 , 58.20961341,
51.69799411, 45.50437379, 40.1551022 , 35.82045353, 32.26373075,
29.16949189, 26.93105771])
These results are compared to the results of Karahan. And as can be seen in the following figure, the output of the model is the same for Procedure II.
Karahan, H. (2012). Predicting Muskingum flood routing parameters using spreadsheets. Computer Applications in Engineering Education, 20(2), 280-286.
Warning
It should be noted that Karahan seems to use an invalid value for \(\Delta t = 6\). According to theory, the minimum value should be \(2kx = 13\).
Verification of multiple river segments - IJssel¶
The only model known that incorporates multiple river segements is the model of Ciullo. In this paper the IJssel is modelled as multiple muskingum segments in sequence. His code is published on Github: https://github.com/quaquel/epa1361_open.
Ciullo, A., de Bruijn, K. M., Kwakkel, J. H., & Klijn, F. (2019). Accounting for the uncertain effects of hydraulic interactions in optimising embankments heights: Proof of principle for the IJssel River. Journal of Flood Risk Management, e12532. https://onlinelibrary.wiley.com/doi/pdf/10.1111/jfr3.12532
Parts of this code were extracted to understand the techniques used. Then my own code was run on the same data and yielded exactly the same results. These two notebooks are not commented but show the same results. The model is thus verified for multiple segments.
IJssel model verification script