Skip to content

Commit 99cfccc

Browse files
author
Marco Dal Molin
committed
Improvements to the numerical routines. Version 1.3.0
1 parent 6f750e0 commit 99cfccc

20 files changed

Lines changed: 750 additions & 63 deletions

File tree

setup.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -3,7 +3,7 @@
33

44
setuptools.setup(
55
name='superflexpy',
6-
version='1.2.1',
6+
version='1.3.0',
77
author='Marco Dal Molin, Fabrizio Fenicia, Dmitri Kavetski',
88
author_email='marco.dalmolin.1991@gmail.com',
99
description='Framework for building hydrological models',

superflexpy/framework/element.py

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -504,6 +504,8 @@ class ODEsElement(StateParameterizedElement):
504504
float
505505
Maximum value of the state. Used, sometimes, by the numerical solver
506506
to search for the solution.
507+
list(floats)
508+
Values of the derivatives of the fluxes w.r.t. the states.
507509
"""
508510

509511
def __init__(self, parameters, states, approximation, id):

superflexpy/implementation/__init__.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -23,4 +23,4 @@
2323
DESIGNED BY: Marco Dal Molin, Fabrizio Fenicia, Dmitri Kavetski
2424
"""
2525

26-
from . import computation, elements, models
26+
from . import elements, models, numerical_approximators, root_finders

superflexpy/implementation/elements/gr4j.py

Lines changed: 26 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -237,11 +237,16 @@ def _flux_function_python(S, S0, ind, P, x1, alpha, beta, ni, PET, dt):
237237
- ((x1[ind]**(1 - beta[ind])) / ((beta[ind] - 1) * dt[ind])) * (ni[ind]**(beta[ind] - 1)) * (S**beta[ind]) # Perc
238238
],
239239
0.0,
240-
S0 + P[ind] * (1 - (S / x1[ind])**alpha[ind]) * dt[ind]
240+
S0 + P[ind] * (1 - (S / x1[ind])**alpha[ind]) * dt[ind],
241+
[
242+
- (P[ind] * alpha[ind] / x1[ind]) * ((S / x1[ind])**(alpha[ind] - 1)),
243+
- (PET[ind] / x1[ind]) * (2 - alpha[ind] * ((S / x1[ind])**(alpha[ind] - 1))),
244+
- beta[ind] * ((x1[ind]**(1 - beta[ind])) / ((beta[ind] - 1) * dt[ind])) * (ni[ind]**(beta[ind] - 1)) * (S**(beta[ind] - 1))
245+
]
241246
)
242247

243248
@staticmethod
244-
@nb.jit('Tuple((UniTuple(f8, 3), f8, f8))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
249+
@nb.jit('Tuple((UniTuple(f8, 3), f8, f8, UniTuple(f8, 3)))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
245250
nopython=True)
246251
def _flux_function_numba(S, S0, ind, P, x1, alpha, beta, ni, PET, dt):
247252

@@ -252,7 +257,12 @@ def _flux_function_numba(S, S0, ind, P, x1, alpha, beta, ni, PET, dt):
252257
- ((x1[ind]**(1 - beta[ind])) / ((beta[ind] - 1) * dt[ind])) * (ni[ind]**(beta[ind] - 1)) * (S**beta[ind]) # Perc
253258
),
254259
0.0,
255-
S0 + P[ind] * (1 - (S / x1[ind])**alpha[ind]) * dt[ind]
260+
S0 + P[ind] * (1 - (S / x1[ind])**alpha[ind]) * dt[ind],
261+
(
262+
- (P[ind] * alpha[ind] / x1[ind]) * ((S / x1[ind])**(alpha[ind] - 1)),
263+
- (PET[ind] / x1[ind]) * (2 - alpha[ind] * ((S / x1[ind])**(alpha[ind] - 1))),
264+
- beta[ind] * ((x1[ind]**(1 - beta[ind])) / ((beta[ind] - 1) * dt[ind])) * (ni[ind]**(beta[ind] - 1)) * (S**(beta[ind] - 1))
265+
)
256266
)
257267

258268

@@ -365,11 +375,16 @@ def _flux_function_python(S, S0, ind, P, x2, x3, gamma, omega, dt):
365375
- (x2[ind] * (S / x3[ind])**omega[ind]), # F
366376
],
367377
0.0,
368-
S0 + P[ind] * dt[ind]
378+
S0 + P[ind] * dt[ind],
379+
[
380+
0.0,
381+
- ((x3[ind]**(1 - gamma[ind])) / ((gamma[ind] - 1) * dt[ind])) * (S**(gamma[ind] - 1)) * gamma[ind],
382+
- (omega[ind] * x2[ind] * ((S / x3[ind])**(omega[ind] - 1))) / x3[ind]
383+
]
369384
)
370385

371386
@staticmethod
372-
@nb.jit('Tuple((UniTuple(f8, 3), f8, f8))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
387+
@nb.jit('Tuple((UniTuple(f8, 3), f8, f8, UniTuple(f8, 3)))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
373388
nopython=True)
374389
def _flux_function_numba(S, S0, ind, P, x2, x3, gamma, omega, dt):
375390

@@ -380,7 +395,12 @@ def _flux_function_numba(S, S0, ind, P, x2, x3, gamma, omega, dt):
380395
- (x2[ind] * (S / x3[ind])**omega[ind]), # F
381396
),
382397
0.0,
383-
S0 + P[ind] * dt[ind]
398+
S0 + P[ind] * dt[ind],
399+
(
400+
0.0,
401+
- ((x3[ind]**(1 - gamma[ind])) / ((gamma[ind] - 1) * dt[ind])) * (S**(gamma[ind] - 1)) * gamma[ind],
402+
- (omega[ind] * x2[ind] * ((S / x3[ind])**(omega[ind] - 1)))
403+
)
384404
)
385405

386406

superflexpy/implementation/elements/hbv.py

Lines changed: 25 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -29,6 +29,7 @@
2929

3030
from ...framework.element import ODEsElement
3131
import numba as nb
32+
import numpy as np
3233

3334

3435
class PowerReservoir(ODEsElement):
@@ -135,11 +136,15 @@ def _fluxes_function_python(S, S0, ind, P, k, alpha, dt):
135136
- k[ind] * S**alpha[ind],
136137
],
137138
0.0,
138-
S0 + P[ind] * dt[ind]
139+
S0 + P[ind] * dt[ind],
140+
[
141+
0.0,
142+
- k[ind] * alpha[ind] * S**(alpha[ind] - 1)
143+
]
139144
)
140145

141146
@staticmethod
142-
@nb.jit('Tuple((UniTuple(f8, 2), f8, f8))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:])',
147+
@nb.jit('Tuple((UniTuple(f8, 2), f8, f8, UniTuple(f8, 2)))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:])',
143148
nopython=True)
144149
def _fluxes_function_numba(S, S0, ind, P, k, alpha, dt):
145150
# This method is used only when solving the equation
@@ -150,7 +155,11 @@ def _fluxes_function_numba(S, S0, ind, P, k, alpha, dt):
150155
- k[ind] * S**alpha[ind],
151156
),
152157
0.0,
153-
S0 + P[ind] * dt[ind]
158+
S0 + P[ind] * dt[ind],
159+
(
160+
0.0,
161+
- k[ind] * alpha[ind] * S**(alpha[ind] - 1)
162+
)
154163
)
155164

156165

@@ -293,11 +302,16 @@ def _fluxes_function_python(S, S0, ind, P, Smax, Ce, m, beta, PET, dt):
293302
- P[ind] * (S / Smax[ind])**beta[ind],
294303
],
295304
0.0,
296-
S0 + P[ind] * dt[ind]
305+
S0 + P[ind] * dt[ind],
306+
[
307+
0.0,
308+
- (Ce[ind] * PET[ind] * m[ind] * (m[ind] + 1) * Smax[ind])/((S + m[ind] * Smax[ind])**2),
309+
- (P[ind] * beta[ind] / Smax[ind]) * (S / Smax[ind])**(beta[ind] - 1),
310+
]
297311
)
298312

299313
@staticmethod
300-
@nb.jit('Tuple((UniTuple(f8, 3), f8, f8))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
314+
@nb.jit('Tuple((UniTuple(f8, 3), f8, f8,UniTuple(f8, 3)))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
301315
nopython=True)
302316
def _fluxes_function_numba(S, S0, ind, P, Smax, Ce, m, beta, PET, dt):
303317
# TODO: handle time variable parameters (Smax) -> overflow
@@ -309,5 +323,10 @@ def _fluxes_function_numba(S, S0, ind, P, Smax, Ce, m, beta, PET, dt):
309323
- P[ind] * (S / Smax[ind])**beta[ind],
310324
),
311325
0.0,
312-
S0 + P[ind] * dt[ind]
326+
S0 + P[ind] * dt[ind],
327+
(
328+
0.0,
329+
- (Ce[ind] * PET[ind] * m[ind] * (m[ind] + 1) * Smax[ind])/((S + m[ind] * Smax[ind])**2),
330+
- (P[ind] * beta[ind] / Smax[ind]) * (S / Smax[ind])**(beta[ind] - 1),
331+
)
313332
)

superflexpy/implementation/elements/hymod.py

Lines changed: 24 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -183,11 +183,16 @@ def _fluxes_function_python(S, S0, ind, P, Smax, m, beta, PET, dt):
183183
- P[ind] * (1 - (1 - (S / Smax[ind]))**beta[ind]),
184184
],
185185
0.0,
186-
S0 + P[ind] * dt[ind]
186+
S0 + P[ind] * dt[ind],
187+
[
188+
0.0,
189+
- PET[ind] * m[ind] * Smax[ind] * (1 + m[ind]) / ((S + Smax[ind] * m[ind])**2),
190+
- P[ind] * beta[ind] * ((1 - (S / Smax[ind]))**(beta[ind] - 1)) / Smax[ind],
191+
]
187192
)
188193

189194
@staticmethod
190-
@nb.jit('Tuple((UniTuple(f8, 3), f8, f8))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
195+
@nb.jit('Tuple((UniTuple(f8, 3), f8, f8, UniTuple(f8, 3)))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
191196
nopython=True)
192197
def _fluxes_function_numba(S, S0, ind, P, Smax, m, beta, PET, dt):
193198
# TODO: handle time variable parameters (Smax) -> overflow
@@ -199,7 +204,12 @@ def _fluxes_function_numba(S, S0, ind, P, Smax, m, beta, PET, dt):
199204
- P[ind] * (1 - (1 - (S / Smax[ind]))**beta[ind]),
200205
),
201206
0.0,
202-
S0 + P[ind] * dt[ind]
207+
S0 + P[ind] * dt[ind],
208+
(
209+
0.0,
210+
- PET[ind] * m[ind] * Smax[ind] * (1 + m[ind]) / ((S + Smax[ind] * m[ind])**2),
211+
- P[ind] * beta[ind] * ((1 - (S / Smax[ind]))**(beta[ind] - 1)) / Smax[ind],
212+
)
203213
)
204214

205215

@@ -307,11 +317,15 @@ def _fluxes_function_python(S, S0, ind, P, k, dt):
307317
- k[ind] * S,
308318
],
309319
0.0,
310-
S0 + P[ind] * dt[ind]
320+
S0 + P[ind] * dt[ind],
321+
[
322+
0.0,
323+
- k[ind]
324+
]
311325
)
312326

313327
@staticmethod
314-
@nb.jit('Tuple((UniTuple(f8, 2), f8, f8))(optional(f8), f8, i4, f8[:], f8[:], f8[:])',
328+
@nb.jit('Tuple((UniTuple(f8, 2), f8, f8, UniTuple(f8, 2)))(optional(f8), f8, i4, f8[:], f8[:], f8[:])',
315329
nopython=True)
316330
def _fluxes_function_numba(S, S0, ind, P, k, dt):
317331
# This method is used only when solving the equation
@@ -322,5 +336,9 @@ def _fluxes_function_numba(S, S0, ind, P, k, dt):
322336
- k[ind] * S,
323337
),
324338
0.0,
325-
S0 + P[ind] * dt[ind]
339+
S0 + P[ind] * dt[ind],
340+
(
341+
0.0,
342+
- k[ind]
343+
)
326344
)

superflexpy/implementation/elements/thur_model_hess.py

Lines changed: 10 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -76,7 +76,7 @@ def __init__(self, parameters, states, approximation, id):
7676
if approximation.architecture == 'numba':
7777
self._fluxes = [self._flux_function_numba]
7878
elif approximation.architecture == 'python':
79-
self._fluxes = [self._flux_function_numba]
79+
self._fluxes = [self._flux_function_python]
8080

8181
def set_input(self, input):
8282
"""
@@ -165,10 +165,14 @@ def _flux_function_python(S, S0, ind, snow, T, t0, k, m, dt):
165165
],
166166
0.0,
167167
S0 + snow[ind] * dt[ind],
168+
[
169+
0.0,
170+
- melt_potential * np.exp(-(S / m[ind])) / m[ind]
171+
]
168172
)
169173

170174
@staticmethod
171-
@nb.jit('Tuple((UniTuple(f8, 2), f8, f8))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
175+
@nb.jit('Tuple((UniTuple(f8, 2), f8, f8, UniTuple(f8, 2)))(optional(f8), f8, i4, f8[:], f8[:], f8[:], f8[:], f8[:], f8[:])',
172176
nopython=True)
173177
def _flux_function_numba(S, S0, ind, snow, T, t0, k, m, dt):
174178

@@ -185,6 +189,10 @@ def _flux_function_numba(S, S0, ind, snow, T, t0, k, m, dt):
185189
),
186190
0.0,
187191
S0 + snow[ind] * dt[ind],
192+
(
193+
0.0,
194+
- melt_potential * np.exp(-(S / m[ind])) / m[ind]
195+
)
188196
)
189197

190198

superflexpy/implementation/models/gr4j.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -25,8 +25,8 @@
2525
This file implements a version of the model GR4J
2626
"""
2727

28-
from superflexpy.implementation.computation.pegasus_root_finding import PegasusPython
29-
from superflexpy.implementation.computation.implicit_euler import ImplicitEulerPython
28+
from superflexpy.implementation.root_finders.pegasus import PegasusPython
29+
from superflexpy.implementation.numerical_approximators.implicit_euler import ImplicitEulerPython
3030
from superflexpy.implementation.elements.gr4j import InterceptionFilter, ProductionStore, UnitHydrograph1, UnitHydrograph2, RoutingStore, FluxAggregator
3131
from superflexpy.implementation.elements.structure_elements import Transparent, Splitter, Junction
3232
from superflexpy.framework.unit import Unit

superflexpy/implementation/models/hymod.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -25,8 +25,8 @@
2525
This file implements a version of the model Hymod
2626
"""
2727

28-
from superflexpy.implementation.computation.pegasus_root_finding import PegasusPython
29-
from superflexpy.implementation.computation.implicit_euler import ImplicitEulerPython
28+
from superflexpy.implementation.root_finders.pegasus import PegasusPython
29+
from superflexpy.implementation.numerical_approximators.implicit_euler import ImplicitEulerPython
3030
from superflexpy.implementation.elements.hymod import UpperZone, LinearReservoir
3131
from superflexpy.implementation.elements.structure_elements import Junction, Splitter, Transparent
3232
from superflexpy.framework.unit import Unit

superflexpy/implementation/models/m4_sf_2011.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -32,8 +32,8 @@
3232
Water Resour. Res., 47, W11511, doi:10.1029/2011WR010748
3333
"""
3434

35-
from superflexpy.implementation.computation.pegasus_root_finding import PegasusPython
36-
from superflexpy.implementation.computation.implicit_euler import ImplicitEulerPython
35+
from superflexpy.implementation.root_finders.pegasus import PegasusPython
36+
from superflexpy.implementation.numerical_approximators.implicit_euler import ImplicitEulerPython
3737
from superflexpy.implementation.elements.hbv import UnsaturatedReservoir, PowerReservoir
3838
from superflexpy.framework.unit import Unit
3939

0 commit comments

Comments
 (0)