Skip to content

Commit 787d407

Browse files
committed
Reacting to new requirements on parameter setting from Solve IVP and fixing some other issues with order of set_integrator, set_initial_value and set_parameters. Fixes Issue #29.
1 parent 94ed85c commit 787d407

File tree

3 files changed

+85
-7
lines changed

3 files changed

+85
-7
lines changed

jitcode/_jitcode.py

Lines changed: 18 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -571,6 +571,7 @@ def set_integrator(self,name,nsteps=10**6,interpolate=True,**integrator_params):
571571
* `"dopri5"` – Dormand’s and Prince’s explicit fifth-order method via `ode`
572572
* `"RK45"` – Dormand’s and Prince’s explicit fifth-order method via `solve_ivp`
573573
* `"dop853"` – DoP853 (explicit) via `ode`
574+
* `"DOP853"` – DoP853 (explicit) via `solve_ivp`
574575
* `"RK23"` – Bogacki’s and Shampine’s explicit third-order method via `solve_ivp`
575576
* `"BDF"` – Implicit backward-differentiation formula via `solve_ivp`
576577
* `"lsoda"` – LSODA (implicit) via `ode`
@@ -598,6 +599,7 @@ def set_integrator(self,name,nsteps=10**6,interpolate=True,**integrator_params):
598599
self.generate_functions()
599600

600601
if info["backend"] == "ode":
602+
self.initialise()
601603
self.integrator = ODE_wrapper(self.f,self.jac)
602604
self.integrator.set_integrator(
603605
name,
@@ -611,26 +613,33 @@ def set_integrator(self,name,nsteps=10**6,interpolate=True,**integrator_params):
611613
self.integrator = IVP(
612614
name,
613615
self.f,
614-
self.jac,
616+
initialiser = lambda: self.initialise(force=True),
617+
jac = self.jac,
615618
**integrator_params
616619
)
617620

618-
# Restore state and params, if applicable:
621+
# Restore state and time, if applicable:
619622
try:
620-
self.set_initial_value(old_integrator._y,old_integrator.t)
623+
old_y = old_integrator._y
624+
old_t = old_integrator.t
621625
except (AttributeError,RuntimeError):
622626
pass
627+
else:
628+
self.set_initial_value(old_y,old_t)
623629

624-
def initialise(self):
630+
def initialise(self,force=False):
625631
if self._initialise is not None:
626-
self._initialise(
632+
if len(self.control_par_values)==len(self.control_pars):
633+
self._initialise(
627634
*self.control_par_values,
628635
*[callback for _,callback,_ in self.callback_functions]
629636
)
637+
elif force:
638+
raise RuntimeError("Something needs parameters to be set. Try calling `set_parameters` earlier.")
630639

631640
def set_parameters(self,*args):
632641
"""
633-
Same as `set_f_params` and `set_jac_params` for SciPy’s ODE (both sets of parameters are set simultaneuosly, because they should be the same anyway).
642+
Same as `set_f_params` and `set_jac_params` for SciPy’s ODE (both sets of parameters are set simultaneuosly, because they should be the same anyway).
634643
635644
The parameters can be passed as different arguments or as a list or other sequence.
636645
"""
@@ -642,6 +651,9 @@ def set_parameters(self,*args):
642651
if len(args)>1:
643652
raise TypeError("Argument must either be a single sequence or multiple numbers.")
644653

654+
if len(self.control_par_values)!=len(self.control_pars):
655+
raise ValueError("Number of values does not match number of control parameters.")
656+
645657
self.initialise()
646658

647659
def set_f_params(self, *args):

jitcode/integrator_tools.py

Lines changed: 9 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -40,12 +40,19 @@ class IVP_wrapper(object):
4040
This is a wrapper around the integrators from scipy.integrate.solve_ivp making them work like scipy.integrate.ode (mostly).
4141
"""
4242

43-
def __init__(self,name,f,jac=None,**kwargs):
43+
def __init__(
44+
self,
45+
name,
46+
f, jac=None,
47+
initialiser = lambda:None,
48+
**kwargs
49+
):
4450
info = integrator_info(name)
4551
self.ivp_class = info["integrator"]
4652
self.f = f
4753
self.jac = jac
4854
self.wants_jac = info["wants_jac"]
55+
self.initialiser = initialiser
4956

5057
# Dictionary to be passed as arguments to the integrator and store stuff
5158
self.kwargs = {
@@ -76,6 +83,7 @@ def try_to_initiate(self):
7683
"t0" in self.kwargs and
7784
"y0" in self.kwargs
7885
):
86+
self.initialiser()
7987
self.backend = self.ivp_class(**self.kwargs)
8088

8189
def set_initial_value(self, initial_value, time=0.0):

tests/test_order.py

Lines changed: 58 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,58 @@
1+
from jitcode import jitcode
2+
from symengine import Symbol
3+
from itertools import permutations
4+
5+
p = Symbol('p')
6+
7+
f = [1/p]
8+
9+
for integrator in ["dopri5","RK45"]:
10+
def set_integrator(ODE):
11+
ODE.set_integrator(integrator)
12+
13+
def set_parameters(ODE):
14+
ODE.set_parameters(1)
15+
16+
def set_initial_value(ODE):
17+
ODE.set_initial_value([0], 0)
18+
19+
at_east_one_working_order = False
20+
for functions in permutations([set_integrator,set_parameters,set_initial_value]):
21+
ODE = jitcode(f, control_pars=[p], verbose=False)
22+
try:
23+
for function in functions:
24+
function(ODE)
25+
except RuntimeError as exception:
26+
assert str(exception).startswith("Something needs parameters to be set")
27+
else:
28+
assert abs(ODE.integrate(1.)-1) < 1e-8
29+
at_east_one_working_order = True
30+
31+
print( ".", end="", flush=True )
32+
assert at_east_one_working_order
33+
34+
f_2 = [1]
35+
36+
for integrator in ["dopri5","RK45"]:
37+
def set_integrator(ODE):
38+
ODE.set_integrator(integrator)
39+
40+
def set_initial_value(ODE):
41+
ODE.set_initial_value([0], 0)
42+
43+
at_east_one_working_order = False
44+
for functions in permutations([set_integrator,set_initial_value]):
45+
ODE = jitcode(f_2, verbose=False)
46+
try:
47+
for function in functions:
48+
function(ODE)
49+
except RuntimeError as exception:
50+
assert str(exception).startswith("Something needs parameters to be set")
51+
else:
52+
assert abs(ODE.integrate(1.)-1) < 1e-8
53+
at_east_one_working_order = True
54+
55+
print( ".", end="", flush=True )
56+
assert at_east_one_working_order
57+
58+
print("")

0 commit comments

Comments
 (0)