Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Change sympy solver to use mpmath directly #15

Open
wants to merge 2 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion odespy/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -730,7 +730,7 @@ def f(self, u, t):
...
Vode
lsoda_scipy
odefun_sympy
odefun_mpmath
odelab

A similar function, ``list_available_solvers``, returns a list of the
Expand Down
4 changes: 2 additions & 2 deletions odespy/problems.py
Original file line number Diff line number Diff line change
Expand Up @@ -122,7 +122,7 @@ def __init__(self, a=1, b=0, c=2):
self.not_suitable_solvers = [
'AdaptiveResidual', 'Lsodi', 'Lsodis', 'Lsoibt']
if self.a < 0:
self.not_suitable_solvers += ['lsoda_scipy', 'odefun_sympy']
self.not_suitable_solvers += ['lsoda_scipy', 'odefun_mpmath']

def f(self, u, t):
return self.a + (u - self.a*t - self.b)**self.c
Expand Down Expand Up @@ -194,7 +194,7 @@ def __init__(self, a=1, b=0, A=1):
raise ValueError('a=%g is too small' % a)
self.not_suitable_solvers = [
'Lsodi', 'Lsodis', 'Lsoibt', 'MyRungeKutta',
'MySolver', 'Lsodes', 'SymPy_odefun']
'MySolver', 'Lsodes', 'mpmath_odefun']

def f(self, u, t):
return self.a*u + self.b
Expand Down
28 changes: 14 additions & 14 deletions odespy/solvers.py
Original file line number Diff line number Diff line change
Expand Up @@ -45,7 +45,7 @@ class ``BackwardEuler``, for instance.
developers need to install and import (in ``initialize``) the desired
package as a Python module and just write a class in the ``Solver``
hierarchy that calls the proper functions in the package. See the
classes ``odefun_sympy`` and ``ode_scipy`` (and its subclasses).
classes ``odefun_mpmath`` and ``ode_scipy`` (and its subclasses).

By an attempt to import these necessary modules (often set in method
initialize()), we can check whether the necessary dependencies are
Expand Down Expand Up @@ -106,7 +106,7 @@ class ``BackwardEuler``, for instance.
software's ability to solve the whole ODE problem. Then it is more
natural to write a new ``solve`` method (using ``Solver.solve`` as
model) and call up the solve functionality in the external
software. Class ``odefun_sympy`` provides an example. On the other
software. Class ``odefun_mpmath`` provides an example. On the other
hand, the wrappers for the ``scipy`` solvers (``vode`` for instance)
applies ``solve`` from the present package and the ``scipy`` functions
for doing one (adaptive) time step, called in the ``advance`` method.
Expand Down Expand Up @@ -993,7 +993,7 @@ def set_initial_condition(self, U0):
# Test first if U0 is sequence (len(U0) possible),
# and use that as indicator for system of ODEs.
# The below code should work for U0 having
# float,int,sympy.mpmath.mpi and other objects as elements.
# float,int,mpmath.mpi and other objects as elements.
try:
self.neq = len(U0)
U0 = np.asarray(U0) # (assume U0 is sequence)
Expand Down Expand Up @@ -2001,30 +2001,30 @@ def approx_Jacobian(f, u0, t0, h):
return J.transpose()


class odefun_sympy(Solver):
class odefun_mpmath(Solver):
"""
Wrapper for the sympy.mpmath.odefun method, which applies a high-order
Wrapper for the mpmath.odefun method, which applies a high-order
Taylor series method to solve ODEs.
"""
quick_description = "Very accurate high order Taylor method (from SymPy)"
quick_description = "Very accurate high order Taylor method (from mpmath)"

def initialize(self):
try:
import sympy
self.sympy = sympy
import mpmath
self.mpmath = mpmath
except ImportError:
raise ImportError,'sympy is not installed - needed for sympy_odefun'
raise ImportError('mpmath is not installed - needed for mpmath_odefun')

def initialize_for_solve(self):
# sympy.odefun requires f(t, u), not f(u, t, *args, **kwargs)
# mpmath.odefun requires f(t, u), not f(u, t, *args, **kwargs)
# (self.f already handles f_args, f_kwargs)
self.f4odefun = lambda t, u: self.f(u, t)
Solver.initialize_for_solve(self)

def solve(self, time_points, terminate=None):
"""
The complete solve method must be overridded in this class
since sympy.mpmath.odefun is such a solve method.
since mpmath.odefun is such a solve method.

The class stores an attribute ufunc (return from odefun)
which can be used to evaluate u at any time point (ufunc(t)).
Expand All @@ -2034,8 +2034,8 @@ def solve(self, time_points, terminate=None):
self.t = np.asarray(time_points)
self.initialize_for_solve()

self.sympy.mpmath.mp.dps = 15 # accuracy
self.ufunc = self.sympy.mpmath.odefun(
self.mpmath.mp.dps = 15 # accuracy
self.ufunc = self.mpmath.odefun(
self.f4odefun, time_points[0], self.U0)

# u and t to be returned are now computed by sampling self.ufunc
Expand Down Expand Up @@ -3062,7 +3062,7 @@ def list_not_suitable_complex_solvers():
'Lsoda', 'Lsodar', 'Lsode', 'Lsodes',
'Lsodi', 'Lsodis', 'Lsoibt',
'RKC', 'RKF45',
'odefun_sympy', 'lsoda_scipy',
'odefun_mpmath', 'lsoda_scipy',
'Radau5', 'Radau5Explicit', 'Radau5Implicit',
'EulerCromer',
]
8 changes: 4 additions & 4 deletions odespy/tests/test_basics.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,7 @@
f_kwargs=dict(a=-1., b=0.,),
# These solvers are not suitable for this ODE problem
exceptions=['Lsodi', 'Lsodis', 'Lsoibt', 'MyRungeKutta',
'MySolver', 'Lsodes', 'SymPy_odefun', 'AdaptiveResidual',
'MySolver', 'Lsodes', 'mpmath_odefun', 'AdaptiveResidual',
'Radau5', 'Radau5Explicit', 'Radau5Implicit', 'EulerCromer'],
time_points=np.linspace(0., 1., 10),
terminate=lambda u,t,step_number: u[step_number] <= 0.4,
Expand Down Expand Up @@ -78,7 +78,7 @@
Trapezoidal=0.0008237438335759,
Vode=0.0000226336696082,
lsoda_scipy=0.0000865549875511,
odefun_sympy=0.0000000000000000,
odefun_mpmath=0.0000000000000000,
),
)

Expand Down Expand Up @@ -136,7 +136,7 @@
Trapezoidal=0.0012657016106976,
Vode=0.0000019483075576,
lsoda_scipy=0.0000003263486278,
odefun_sympy=0.0000000000000000,
odefun_mpmath=0.0000000000000000,
),
)

Expand Down Expand Up @@ -191,7 +191,7 @@
Trapezoidal=0.0375041929371673,
Vode=0.0375037017349515,
lsoda_scipy=0.0375036977160605,
odefun_sympy=0.0375036952176233,
odefun_mpmath=0.0375036952176233,
),
)

Expand Down