66
77 .. _numerical_solver :
88
9+ Numerical implementation
10+ ========================
11+
912Numerical routines for solving ODEs
10- ===================================
13+ -----------------------------------
1114
1215:ref: `reservoirs ` are the most common elements in conceptual hydrological
1316models. Reservoirs are controlled by one (or more) ordinary differential
@@ -35,10 +38,13 @@ These steps can be performed using two SuperflexPy components:
3538:code: `NumericalApproximator ` and :code: `RootFinder `.
3639
3740SuperflexPy provides two built-in numerical approximators (implicit and explicit
38- Euler) and a root finder (Pegasus method).
41+ Euler) and a root finder (Pegasus method). These methods are best suited when
42+ dealing with smooth flux functions. If a user wants to experiment with
43+ discontinuous flux functions, other ODE solution algorithms should be
44+ considered.
3945
40- Other ODE solution algorithms, e.g. Runge-Kutta, can be accommodated
41- by extending the classes :code: `NumericalApproximator ` and/or
46+ The implementation of other ODE solution algorithms, e.g. Runge-Kutta, can be
47+ done by extending the classes :code: `NumericalApproximator ` and/or
4248:code: `RootFinder `. Even more generally, ODE solvers from external libraries
4349could be used within SuperflexPy by incorporating them in a class that respects
4450the interface expected by the :code: `NumericalApproximator `.
@@ -47,7 +53,7 @@ The following sections describe the standard approach to create customized
4753numerical approximators and root finders.
4854
4955Numerical approximator
50- ----------------------
56+ ......................
5157
5258A customized numerical approximator can be implemented by extending the class
5359:code: `NumericalApproximator ` and implementing two methods: :code: `_get_fluxes `
@@ -77,7 +83,7 @@ For further details, please see the implementation of :code:`ImplicitEuler` and
7783:code: `ExplicitEuler `.
7884
7985Root finder
80- -----------
86+ ...........
8187
8288A customized root finder can implemented by extending the class
8389:code: `RootFinder ` implementing the method :code: `solve `.
@@ -94,11 +100,35 @@ the fluxes, :code:`S0` is the initial state, :code:`dt` is the time step,
94100:code: `fluxes `, and :code: `ind ` is the index of the input arrays to use.
95101
96102The method :code: `solve ` is responsible for finding the numerical solution of
97- the approximated ODE.
103+ the approximated ODE. In case of failure, the method should either raise a
104+ :code: `RuntimeError ` (Python implementation) or return :code: `numpy.nan ` (this
105+ is not ideal but, since Numba does not support raising exceptions, is the
106+ solution that should be adopted with the Numba implementation).
98107
99108To understand better how this method works, please see the implementation of
100109:code: `Pegasus `.
101110
111+ Sequential solution of the elements
112+ -----------------------------------
113+
114+ The SuperflexPy framework is built on a model representation that maps to a
115+ directional acyclic graph. Model elements are solved sequentially from upstream
116+ to downstream, with the output from each element being used as input to
117+ its downstream elements.
118+
119+ When fixed-step solvers are used (e.g. implicit Euler), this
120+ "one-element-at-a-time" strategy is equivalent to applying the same (fixed-step)
121+ solver to the entire ODE system simultaneously.
122+
123+ When more advanced solvers use internal substepping, then the
124+ "one-element-at-a-time" strategy does introduces additional approximation error.
125+ This additional approximation error is due to treating the fluxes as constant
126+ over the time step, whereas the exact solution would have varying fluxes within
127+ the time step. However, in most practical applications, this "uniform flux"
128+ approximation is already applied to the meteorological inputs (precipitation and
129+ PET), hence applying it to internal fluxes does not represent a large additional
130+ approximation.
131+
102132Computational efficiency with Numpy and Numba
103133---------------------------------------------
104134
0 commit comments