Usage Examples

Julia

Use of this package can be found in the Streamfall.jl package here.

Examples of direct function calls through Julia's ccall method is shown below.

Julia pre-v1.5.3

ccall((:calc_outflow, "./lib/ihacres.so"), Cdouble,(Cdouble, Cdouble), 1.0, 1.0)

Julia v1.5.3 onwards

These versions include a @ccall macro which simplifies usage. Note: if on Windows, replace .so with .dll.

const IHACRES = "./lib/ihacres.so"
@ccall IHACRES.calc_outflow(10.0::Cdouble, 8.0::Cdouble)::Cdouble

Python

Use of the compiled DLL/so/dynlib is possible via the ctypes package, but declaring input/return types can get tedious (see the tutorial here).

Instead, the most straightforward approach is to use the nimporter package.

In future, a dedicated Python package may be made available.

1. Clone ihacres_nim

$ git clone https://github.com/ConnectedSystems/ihacres_nim.git

2. Install nimporter and ihacres-nim as a Python package

$ pip install nimporter
$ pip install -e .

Example usage is provided under tests/py_example.py

$ python tests/py_example.py

The contents of py_example.py are shown below.

import nimporter  # Required prior to any Nim module import

from ihacres import climate
from ihacres import flow
from ihacres import cmd as ihacres_cmd

nimporter.build_nim_extensions()  # Build Nim extension


def run_ihacres(cmd, rainfall, evaporation, inflow, quickflow, slowflow,
                gw_state, gw_exchange, extraction):
    '''Example run function.

    Parameters
    ----------
    cmd : catchment moisture deficit, $M_{k}$ in the literature.
    rainfall : precipitation ($P$)
    evaporation : ($E$)
    inflow : flow from previous node
    quickflow : quickflow at previous time step, $(V_{q-1})$
    slowflow : slowflow at previous time step, $(V_{s-1})$
    gw_state : groundwater storage index at $t-1$
    gw_exchange : volume flux, interaction with groundwater
    extraction : volume extracted from stream
    '''
    d, d2, e, f = run_ihacres.d, run_ihacres.d2, run_ihacres.e, run_ihacres.f
    a, b, alpha, s = run_ihacres.a, run_ihacres.b, run_ihacres.alpha, run_ihacres.s
    catchment_area = run_ihacres.catchment_area

    mf, U, r = ihacres_cmd.calc_ft_interim_cmd(cmd, rainfall, d, d2, alpha)
    ET = climate.calc_ET(e, evaporation, mf, f, d)
    cmd = ihacres_cmd.calc_cmd(cmd, rainfall, ET, U, r)

    Vq, Vs, outflow = flow.calc_ft_flows(quickflow, slowflow, U, r,
                                         catchment_area, a, b)

    # if node routes to another node
    gw_state, outflow = flow.routing(gw_state, s, inflow, outflow, extraction, gw_exchange)

    return {
        "flow": (gw_state, outflow),
        "state": (Vq, Vs, cmd)
    }

# Here we leverage the fact that everything in Python is an object
# including functions. We assign the model parameters to the function
# such that the function itself represents a node in a stream network.
# Of course, you could (should?) define a Class instead.
run_ihacres.d = 200.0  # flow threshold
run_ihacres.d2 = 2.0   # flow threshold, multiplier applied to `d`
run_ihacres.e = 0.1    # temperature to PET conversion factor
run_ihacres.f = 0.1    # plant stress threshold factor (applied to `d`). Determines effective rainfall.
run_ihacres.a = 54.35  # quickflow scaling factor
run_ihacres.b = 0.012  # slowflow scaling factor
run_ihacres.alpha = 0.727   #  effective rainfall scaling factor
run_ihacres.s = 2.5  # groundwater store factor
run_ihacres.catchment_area = 100.0  # area in km^2


# Initial conditions
cmd, quickflow, slowflow = 214.65, 0, 0

# Assume these are all 0.0
inflow = 0.0
gw_exchange = 0.0
extraction = 0.0

# Ideally, these would be read in via Numpy or Pandas
rainfall_ts = [70.0, 10.0, 0.0, 0.0, 200.0]
evaporation_ts = [2.0, 6.5, 7.0, 5.0, 1.0]

# Set up arrays to record state
outflow = [None] * len(rainfall_ts)
gw_state = [None] * (len(rainfall_ts)+1)
gw_state[0] = 0.0

# Run model (this would be in its own function)
for i in range(len(outflow)):
    progress = run_ihacres(cmd, rainfall_ts[i], evaporation_ts[i], inflow, quickflow, slowflow,
                           gw_state[i], gw_exchange, extraction)
    quickflow, slowflow, cmd = progress["state"]
    gw_state[i+1], outflow[i] = progress["flow"]

# Remove initial gw state to line up records
gw_state = gw_state[1:]

print(outflow)

The streamflow output should be:

    [2219.61106301161, 300.57003807505123, 36.753372670041266, 14.953875444140436, 8905.129670965001]