...
Storage Routing links model the storage and movement of water through a length of a river (drainage link) using a hydrologic routing method. They can represent the travel time of water through a reach, the attenuation of flow rates due to channel shape, and roughness and reach processes.
...
Where S(x) is the link storage at flow rate x.
There are two choices in Source for defining the relationship in Equation 9:
...
Relationship in Equation 9 in MUSICX is defined using a functional relationship (a power function)
...
.
...
Using a Functional Relationship
The storage function has the following form:
...
As noted in the Assumptions and Constraints section, above, MUSICX can divide the storage representing the flow routing in a routing link into a number of equal divisions. This means that the flow in a routing link passes through a cascade of storages. The user inputs a value of the storage-delay constant, k, which applies to each division on the link for use in the K(
) relationship (Equation 16) for each division.
Solution Methodology
The solution methodology for Storage Routing links is discussed below
Ordering phase
In the ordering phase, the wave travel time must be estimated before the flows for the current time step are known. The estimated travel time is called “Order Travel Time” and its calculation depends on the type of routing used:
- For Storage Routing links with linear Muskingum routing (), the wave travel time is constant and known. The Order Travel Time is hence equal to the wave travel time and is calculated automatically by Source.
- For Storage Routing links with non-linear and variable parameter Muskingum routing, the travel time depends on the flow. The user must supply an “Average Regulated Flow”, which Source uses to calculate the corresponding order travel time.
Further information on how ordering, resource assessment and ownership are handled in links is provided in the chapters on Rules-Based Ordering - SRG, Resource Assessment - SRG and Ownership - SRG.
Flow phase
It should be noted that if x=1 then a direct solution for O can be made (as =
) when the division is flowing. Otherwise an iterative solution is required, for which a modified Newton-Raphson technique is used with two convergence criteria. If either one is satisfied the procedure is considered to have converged. The criteria are a change in index flow of less than 1.0e-8 m3/s (equivalent to 0.864 litres per day) and mass balance error less than 0.001 m3 (equivalent to 1 litre).
- Assume
= 0 then min = x•. Determine the storage associated with the minimum index flow S(min) and the mass balance error based on Equation 31 MBmin=fMB (min,0). If MBmin >-0.001 then the stream has ceased to flow and the storage is based on Equation 32; - Determine the maximum possible index flow:
- If max ≤ min then the stream has ceased to flow as fluxes exceed inflow and storage.
- Determine the storage associated with the maximum index flow S(max) and the mass balance error based on Equation 31 MBmax=fMB (max,O). If MBmax < 0.001 then the outflow is the maximum possible outflow;
- Initial estimate
.
Determine the storage associated with the index flow S(), O and the mass balance error based on Equation 31 MB = fMB (, 
). If |MB| < 0.001 a solution has been found; - Use a modified Newton-Raphson solver limited to 20 iterations
- If MB > 0.0 then
max = and MBmax = MB
else
min = and MBmin = MB
- If |-halve| < 1.0e-8 a solution has been found =halve.
- If |-bisect| < 1.0e-8 a solution has been found =bisect.
- Determine the storage associated with the halve index flow S(halve), and the mass balance error MBhalve = fMB (halve,). If |MBhalve| < 0.001 a solution has been found =halve.
- Determine the storage associated with the bisect index flow S(bisect), and the mass balance error MBbisect = fMB (bisect,). If |MBbisect|<0.001 a solution has been found =bisect.
- Determine the slope of the mass balance index flow function
- Determine the Newton-Raphson correction =
- NR = - correction
- If min < NR < max then
- If correction < 1.0e-8 then a solution has been found = NR.
- Determine the storage associated with the index flow S(NT) and the mass balance error MBNR=fMB(NR, ). If |MBNR | < 0.001 a solution has been found = NR.
- If MBNR > 0.0 and min< NR < max then max = NR and MBmax = MBNR else min = NR and MBmin = MBNR
- If MBbisect > 0.0 and min< bisect < max then max = bisect and MBmax = MBbisect else min = bisect and MBmin = MBbisect.
- If MBhalve > 0.0 and min< halve < max then max = halve and MBmax = MBhalve else min = bisect and MBmin = MBhalve
- If min has not been changed
= max
else if max has not been changed
= min
else if |MB|min ≤ MBmax
= min
max and MBmax set to previous max
else
= max
min and MBmin set to previous min.
- Determine the storage associated with the index flow S(), and the mass balance error MB = fMB (, ).
- If S() ≤ Dmax a solution has been found with zero outflow
- If |MB| < 0.001 a solution has been found
- Calculate the outflow:
In reaches with multiple divisions the outflow from the reach becomes the inflow into the downstream division. The outflow from the final division is the outflow for the routing reach.
Data
Refer to the Source User Guide for detailed data requirements and formats.
Input data
...
Input data
The information that users may provide is summarised in the Table below. Note that there is also a global requirement to specify the model time-step, dt. In addition, of course, data on the inflow to the routing link is required for each model time-step.
Table 1. Summary of parameters requiring user input
...
Parameters
...
Routing- Power Function
(Linear and Non-linear)
...
Routing- Piecewise
required for each model time-step.
Parameters | ||
|---|---|---|
Length of routing reach | Yes | Yes |
Initialisation Type (governs whether an initial value of storage or an initial value of flow is entered) | Yes | Yes |
Initial value of storage | Yes | Yes |
Initial value of flow (note this is mandatory for Lagged Flow)Yes | Yes | |
Lag time | No | No |
Number of routing divisions in the linkYes | Yes | |
Representative flow rate (to calculate link delivery time) | Yes | Yes |
x , the weighting factor in routing | Yes | Yes |
The storage delay constant, k, if m≠1 or the Muskingum K, if m=1, in the routing equation | Yes | No |
m, the exponent in the routing equation | Yes | No |
Flow versus travel time relationship | No | Yes |
Governing data for fluxes (reach processes) - details below | Yes | Yes |
Data needs for lateral fluxes (reach processes)
...