Hydrological operators#

This section details fluxes and states computation for a given cell \(x\in\Omega\) also denoted \(x \in \mathcal{T}_\Omega\) and for a given time step \(t\in[1..N_t]\) considering a regular temporal grid of time step \(\Delta t\). For this cell at time step \(t\), we denote by \(P(t)\) and \(E(t)\) the local total rainfall and evapotranspiration.

GR#

The operators used here come from GR (Génie Rural) models of the literature:

Interception#

Given an interception reservoir \(\mathcal{I}\) of maximum capacity \(c_i\). If potential evapotranspiration \(E\) is greater than the sum of liquid precipitation \(P\) and initial water level of the interception reservoir \(h_i\), then the interception reservoir is emptied. Conversely, if the sum of liquid precipitation and initial level of the interception reservoir is greater than potential evapotranspiration, the interception reservoir is filled in depending on it’s available storage:

\begin{eqnarray} &E_i(t)& &=& &\min \left[ E(t), \; P(t) + h_i(t-1) \right]& \\ &P_{n}(t)& &=& &\max \left[ 0, \; P(t) + h_i(t-1) - c_i - E_i(t) \right]& \\ &h_i(t)& &=& &h_i(t-1) + P(t) - E_i(t) - P_{n}(t)& \end{eqnarray}

where \(P_{n}\) corresponds to the remaining rainfall amount (throughfall) inflowing next flow operators.

Production#

The production component ensures the role of non linear runoff production function.

Initially proposed for a minimal complexity description of catchment water balance functioning, based on empirical modeling, the GR model [Edijatno and Michel, 1989] considers a production reservoir \(\mathcal{P}\) of maximum depth \(c_p\) and water level \(h_p\). The neutralized rainfall and evaporation are respectively denoted \(P_{n}\) and \(E_n\), with \(E_n = E - E_i\):

\begin{eqnarray} dh_p = \begin{cases} &\left( 1 - \left( \frac{h_p}{c_p} \right) ^ 2 \right) dP_{n} &\text{if} \; P_{n} > 0 \\ &-\frac{h_p}{c_p} \left(2 - \frac{h_p}{c_p} \right) dE_n &\text{otherwise} \end{cases} \end{eqnarray}

Assuming a stepwise approximation of the inputs \(P_n\) and \(E_n\) the temporal integration of these ordinary differential equations, enabling analytical solutions, gives the infiltrating rainfall \(P_s\) and the actual evapotranspiration \(E_s\) from the reservoir store :

\begin{eqnarray} &P_s(t)& &=& &c_p \left( 1 - \left( \frac{h_p(t - 1)}{c_p} \right) ^ 2 \right) \; \frac{\tanh \left( \frac{P_{n}(t)}{c_p} \right) }{1 + \left( \frac{h_p(t - 1)}{c_p} \right) \; \tanh \left( \frac{P_{n}(t)}{c_p} \right)}& \\ &E_s(t)& &=& &h_p(t - 1) \left( 2 - \frac{h_p(t - 1)}{c_p} \right) \; \frac{\tanh \left( \frac{E_{n}(t)}{c_p} \right) }{1 + \left( 1 - \frac{h_p(t - 1)}{c_p} \right) \; \tanh \left( \frac{E_{n}(t)}{c_p} \right)}& \end{eqnarray}

\(h_p\) is the water level of the production reservoir and \(P_s\) and \(E_s\) are the amount of water gained or lost over \(\Delta t\) and used to update \(h_p\) at time step t:

\begin{eqnarray} h_p \left( t \right) = h_p(t - 1) + P_s(t) - E_s(t) \end{eqnarray}

Non conservative exchange#

A non-conservative exchange function, representing deep percolation or inter-catchment groudwater flow for instance, is expressed following [Edijatno, 1991]. Given a power law transfer storage \(\mathcal{T} _{ft}\) of capacity \(c_{ft}\) and water level \(h_{ft}\) the exchange term writes:

\begin{eqnarray} F(t) = exc \left( \frac{h_{ft}(t - 1)}{c_{ft}} \right) ^ {7 / 2} \end{eqnarray}

Transfer#

Lateral flows within pixels is represented with the following transfer formulations.

Single power law transfer storage and direct branch#

Transfer within pixels can be first represented by spliting the runoff \(P_r\) into \(Q9 = 0.9 P_r\) inflowing one branch containing a transfer reservoir \(\mathcal{T} _{ft}\) of capacity \(c_{ft}\) and water level \(h_{ft}\), and the remaining \(Q1 = 0.1 P_r\) inflowing a direct branch - i.e. without reservoir. Each transfer branch is also inflowed by the exchange term \(F\).

At the begining of the time step, the level of the power law transfer reservoir \(\mathcal{T} _{ft}\) is updated as:

\begin{eqnarray} h_{ft} \left(t ^ * \right) = \max \left( \epsilon, \; h_{ft}(t - 1) + 0.9 P_r(t) + F(t) \right) \end{eqnarray}

With \(\epsilon>0\), a fixed small constant.

Next, the outflow discharge from the transfer reservoir writes:

\begin{eqnarray} Q_{ft}(t) = h_{ft} \left( t ^ * \right) - \left[ h_{ft} \left( t ^ * \right) ^ {-4} + c_{ft} ^ {-4} \right] ^ {-1 / 4} \end{eqnarray}

The level of the transfer reservoir is updated as:

\begin{eqnarray} h_{ft}(t)=h_{ft} \left( t ^ * \right) - Q_{ft}(t) \end{eqnarray}

In the branch without reservoir the outflow discharge writes:

\begin{eqnarray} Q_d(t) = \max \left[0, \; 0.1 P_r(t) + F(t) \right] \end{eqnarray}

With flux \(Q_t\) inflowing the routing part equal to:

\begin{eqnarray} Q_t(t) = Q_{ft}(t) + Q_d(t) \end{eqnarray}

Double power law transfer storages and direct branch#

Transfer within pixels is represented similarly as in Single power law transfer storage and direct branch but with a second power law transfer reservoir. Again, the runoff \(P_r\) is splitted into \(Q9 = 0.9 P_r\) and \(Q1 = 0.1 P_r\), the latter \(Q1\) inflowing a direct branch - i.e. without reservoir. In the reservoirs branch, the inflow \(Q9\) is separated a second time into \(40 \%\) and \(60 \%\) respectively inflowing two transfer reservoirs \(\mathcal{T} _{ft}\) and \(\mathcal{T} _{st}\).

Again, the exchange term \(F\) is applied to the direct branch and to the reservoir \(\mathcal{T} _{ft}\)

\begin{eqnarray} &h_{ft} \left(t ^ * \right)& &=& &\max \left( \epsilon, \; h_{ft}(t - 1) + 0.9 \times 0.6 \; P_r(t) + F(t) \right)& \\ &h_{st} \left(t ^ * \right)& &=& &\max \left( \epsilon, \; h_{st}(t - 1) + 0.9 \times 0.4 \; P_r(t) \right)& \end{eqnarray}

Next, the outflow discharges from the reservoirs writes:

\begin{eqnarray} &Q_{ft}\left(t \right)& &=& & h_{ft} \left( t ^ * \right) - \left[ h_{ft} \left( t ^ * \right) ^ {-4} + c_{ft} ^ {-4} \right] ^ {-1 / 4} & \\ &Q_{st}\left(t \right)& &=& & h_{st} \left( t ^ * \right) - \left[ h_{st} \left( t ^ * \right) ^ {-4} + c_{st} ^ {-4} \right] ^ {-1 / 4}& \end{eqnarray}

In the branch without reservoir:

\begin{eqnarray} Q_d(t) = \max \left[0, \; 0.1 P_r(t) + F(t) \right] \end{eqnarray}

The level of the transfer reservoirs is updated as:

\begin{eqnarray} &h_{ft}(t)& &=& &h_{ft} \left( t ^ * \right) - Q_{ft}(t)& \\ &h_{st}(t)& &=& &h_{st} \left( t ^ * \right) - Q_{st}(t)& \end{eqnarray}

With flux \(Q_t\) inflowing the routing part equal to:

\begin{eqnarray} Q_t(t) = Q_{ft}(t) + Q_{st}(t) + Q_d(t) \end{eqnarray}

Generic#

Surface routing#

Surface runoff is conveyed from pixels to pixel to the outlet of the basin, following the drainage plan \(\mathcal{D}_{\Omega}\). Several routing models of different complexity are available.

Linear Reservoir#

Given \(N_{\text{xup}}\) upstream cells within \(\Omega\) flowing into cell \(x\) as imposed by the flow direction map \(\mathcal{D}_{\Omega}\), the upstream runoff is:

\begin{eqnarray} Q_{up}(x,t) = \sum _{k} ^ {N_{\text{xup}}} Q(k, t) \end{eqnarray}

With for a given upstream cell \(k\):

\begin{eqnarray} Q(k,t) = h_{lr}(k,t) \left[ 1 - \exp \left( - \frac{\Delta t}{60 \; lr} \right) \right] + Q_{up}(k,t) + Q_t(k,t) \end{eqnarray}

Updating the level in the routing storage:

\begin{eqnarray} h_{lr}(k,t) = h_{lr}(k,t-1) - Q(k,t) \end{eqnarray}