Quickstart (.py)¶

Note - this script can also be opened in interactive Python if you wanted to play around. On the GitHub it is in docs/demo/scripts

1. Introduction

2. Data

3. Create nodes

4. Create arcs

5. Create model

6. Run model

7. What next?

Introduction¶

WSIMOD is a Python package that lets you create nodes that represent physical things in the water cycle and enables them talk to each other. In this demo, we will create some nodes, some arcs (the things that link nodes), and make a model that lets us simulate the flow of water through the water cycle.

We will create a simple catchment model that contains a hydrological node, a sewer node, a groundwater node, and a river node. Using these nodes, we will simulate a mixed urban-rural catchment.

Imports and forcing data¶

Import packages

In [1]:

import os

import pandas as pd
from matplotlib import pyplot as plt

from wsimod.core import constants
from wsimod.orchestration.model import Model

import os import pandas as pd from matplotlib import pyplot as plt from wsimod.core import constants from wsimod.orchestration.model import Model

Load input data, in this example we use precipitation, temperature, and evapotranspiration (et0).

In [2]:

# Select the root path for the data folder. Use the appropriate value for your case.
data_folder = os.path.join(os.path.abspath(""), "docs", "demo", "data")


input_fid = os.path.join(data_folder, "processed", "timeseries_data.csv")
input_data = pd.read_csv(input_fid)
input_data.loc[input_data.variable == "precipitation", "value"] *= constants.MM_TO_M
input_data.date = pd.to_datetime(input_data.date)
input_data = input_data.loc[input_data.site == "oxford_land"]
dates = input_data.date.drop_duplicates()
print(input_data.sample(10))

# Select the root path for the data folder. Use the appropriate value for your case. data_folder = os.path.join(os.path.abspath(""), "docs", "demo", "data") input_fid = os.path.join(data_folder, "processed", "timeseries_data.csv") input_data = pd.read_csv(input_fid) input_data.loc[input_data.variable == "precipitation", "value"] *= constants.MM_TO_M input_data.date = pd.to_datetime(input_data.date) input_data = input_data.loc[input_data.site == "oxford_land"] dates = input_data.date.drop_duplicates() print(input_data.sample(10))

              site       date       variable      value
79995  oxford_land 2012-12-03    temperature   5.509375
76644  oxford_land 2011-09-21  precipitation   0.000800
79521  oxford_land 2011-08-17    temperature  16.578571
76888  oxford_land 2012-05-22  precipitation   0.000000
76661  oxford_land 2011-10-08  precipitation   0.000100
76663  oxford_land 2011-10-10  precipitation   0.000000
78889  oxford_land 2009-11-23    temperature   9.789286
77346  oxford_land 2009-08-28            et0   0.002000
77971  oxford_land 2011-05-15            et0   0.002000
77885  oxford_land 2011-02-18            et0   0.002000

Input data is stored in dicts, where a key is a variable on a given day.

In [3]:

land_inputs = input_data.set_index(["variable", "date"]).value.to_dict()

example_date = pd.to_datetime("2009-03-03")

print(land_inputs[("precipitation", example_date)])
print(land_inputs[("et0", example_date)])
print(land_inputs[("temperature", example_date)])

land_inputs = input_data.set_index(["variable", "date"]).value.to_dict() example_date = pd.to_datetime("2009-03-03") print(land_inputs[("precipitation", example_date)]) print(land_inputs[("et0", example_date)]) print(land_inputs[("temperature", example_date)])

0.015
0.002
7.625

Create nodes¶

Nodes can be defined as dictionaries of parameters. Different nodes require different parameters, you can see the documentation in the API to understand what parameters can be set. Although every node should have a type (type_) and a name.

In [4]:

sewer = {"type_": "Sewer", "capacity": 0.04, "name": "my_sewer"}

surface1 = {
    "type_": "ImperviousSurface",
    "surface": "urban",
    "area": 10,
    "pollutant_load": {"phosphate": 1e-7},
}

surface2 = {
    "type_": "PerviousSurface",
    "surface": "rural",
    "area": 100,
    "depth": 0.5,
    "pollutant_load": {"phosphate": 1e-7},
}

land = {
    "type_": "Land",
    "data_input_dict": land_inputs,
    "surfaces": [surface1, surface2],
    "name": "my_land",
}

gw = {"type_": "Groundwater", "area": 100, "capacity": 100, "name": "my_groundwater"}

node = {"type_": "Node", "name": "my_river"}

waste = {"type_": "Waste", "name": "my_outlet"}

sewer = {"type_": "Sewer", "capacity": 0.04, "name": "my_sewer"} surface1 = { "type_": "ImperviousSurface", "surface": "urban", "area": 10, "pollutant_load": {"phosphate": 1e-7}, } surface2 = { "type_": "PerviousSurface", "surface": "rural", "area": 100, "depth": 0.5, "pollutant_load": {"phosphate": 1e-7}, } land = { "type_": "Land", "data_input_dict": land_inputs, "surfaces": [surface1, surface2], "name": "my_land", } gw = {"type_": "Groundwater", "area": 100, "capacity": 100, "name": "my_groundwater"} node = {"type_": "Node", "name": "my_river"} waste = {"type_": "Waste", "name": "my_outlet"}

Create arcs¶

Arcs can also be created as dictionaries, they don't typically need any numerical parameters (although there are some exceptions in the case study demo). Though they do need to specify the in_port (where the arc starts), and out_port (where it finishes). It's also handy to give each arc a name.

In [5]:

urban_drainage = {
    "type_": "Arc",
    "in_port": "my_land",
    "out_port": "my_sewer",
    "name": "urban_drainage",
}

percolation = {
    "type_": "Arc",
    "in_port": "my_land",
    "out_port": "my_groundwater",
    "name": "percolation",
}

runoff = {
    "type_": "Arc",
    "in_port": "my_land",
    "out_port": "my_river",
    "name": "runoff",
}

storm_outflow = {
    "type_": "Arc",
    "in_port": "my_sewer",
    "out_port": "my_river",
    "name": "storm_outflow",
}

baseflow = {
    "type_": "Arc",
    "in_port": "my_groundwater",
    "out_port": "my_river",
    "name": "baseflow",
}

catchment_outflow = {
    "type_": "Arc",
    "in_port": "my_river",
    "out_port": "my_outlet",
    "name": "catchment_outflow",
}

urban_drainage = { "type_": "Arc", "in_port": "my_land", "out_port": "my_sewer", "name": "urban_drainage", } percolation = { "type_": "Arc", "in_port": "my_land", "out_port": "my_groundwater", "name": "percolation", } runoff = { "type_": "Arc", "in_port": "my_land", "out_port": "my_river", "name": "runoff", } storm_outflow = { "type_": "Arc", "in_port": "my_sewer", "out_port": "my_river", "name": "storm_outflow", } baseflow = { "type_": "Arc", "in_port": "my_groundwater", "out_port": "my_river", "name": "baseflow", } catchment_outflow = { "type_": "Arc", "in_port": "my_river", "out_port": "my_outlet", "name": "catchment_outflow", }

Create model¶

We can create a model object and add dates

In [6]:

my_model = Model()
my_model.dates = dates

my_model = Model() my_model.dates = dates

Add nodes in a list. The Model object will create the nodes from the dict entries.

In [7]:

my_model.add_nodes([sewer, land, gw, node, waste])

my_model.add_nodes([sewer, land, gw, node, waste])

Add arcs in a list.

In [8]:

my_model.add_arcs(
    [urban_drainage, percolation, runoff, storm_outflow, baseflow, catchment_outflow]
)

my_model.add_arcs( [urban_drainage, percolation, runoff, storm_outflow, baseflow, catchment_outflow] )

Run model¶

Run the model

In [9]:

flows, _, _, _ = my_model.run()
flows = pd.DataFrame(flows)

flows, _, _, _ = my_model.run() flows = pd.DataFrame(flows)

  0%|                                                                                                                                                                                                                      | 0/1456 [00:00<?, ?it/s]

 13%|██████████████████████████▍                                                                                                                                                                               | 191/1456 [00:00<00:00, 1895.99it/s]

 26%|████████████████████████████████████████████████████▊                                                                                                                                                     | 381/1456 [00:00<00:00, 1729.80it/s]

 38%|████████████████████████████████████████████████████████████████████████████▉                                                                                                                             | 555/1456 [00:00<00:00, 1657.95it/s]

 50%|████████████████████████████████████████████████████████████████████████████████████████████████████▍                                                                                                     | 724/1456 [00:00<00:00, 1666.16it/s]

 62%|████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████▎                                                                             | 896/1456 [00:00<00:00, 1680.54it/s]

 73%|███████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████                                                      | 1065/1456 [00:00<00:00, 1664.01it/s]

 85%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████                               | 1232/1456 [00:00<00:00, 1656.74it/s]

 96%|████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████▉        | 1398/1456 [00:00<00:00, 1646.68it/s]

100%|█████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 1456/1456 [00:00<00:00, 1668.23it/s]

Plot results

In [10]:

f, axs_ = plt.subplots(3, 2, figsize=(10, 10))

for axs, variable in zip(axs_.T, ["flow", "phosphate"]):
    flows_plot = flows.pivot(index="time", columns="arc", values=variable)
    input_data.pivot(index="date", columns="variable", values="value")[
        ["precipitation"]
    ].plot(ax=axs[0], xlabel="", xticks=[], sharex=True)

    flows_plot[["runoff", "urban_drainage", "percolation"]].plot(
        ax=axs[1], xlabel="", xticks=[], sharex=True
    )
    flows_plot[["catchment_outflow", "storm_outflow", "baseflow"]].plot(ax=axs[2])
axs_[1, 0].set_title("Flows (m3/d)")
axs_[1, 1].set_title("Phosphate (kg/d)")
f.tight_layout()

f, axs_ = plt.subplots(3, 2, figsize=(10, 10)) for axs, variable in zip(axs_.T, ["flow", "phosphate"]): flows_plot = flows.pivot(index="time", columns="arc", values=variable) input_data.pivot(index="date", columns="variable", values="value")[ ["precipitation"] ].plot(ax=axs[0], xlabel="", xticks=[], sharex=True) flows_plot[["runoff", "urban_drainage", "percolation"]].plot( ax=axs[1], xlabel="", xticks=[], sharex=True ) flows_plot[["catchment_outflow", "storm_outflow", "baseflow"]].plot(ax=axs[2]) axs_[1, 0].set_title("Flows (m3/d)") axs_[1, 1].set_title("Phosphate (kg/d)") f.tight_layout()

What next?¶

If you are hydrologically minded then you might be interested in a more detailed tutorial for our Land node. If you want an overview to see how many different nodes in WSIMOD work, then the Oxford case study might be more interesting.