odespy.solvers

This module contains the base class Solver package and the implementations of many subclasses.

A new non-adaptive time-stepping scheme can often be implemented by just subclassing Solver and writing the advance method to define the time-stepping scheme. Here is how the 4th-order Runge-Kutta is implemented:

class RK4(Solver):
    quick_description = "Explicit 4th-order Runge-Kutta method"

    def advance(self):
        u, f, n, t = self.u, self.f, self.n, self.t
        dt = t[n+1] - t[n]
        dt2 = dt/2.0
        K1 = dt*f(u[n], t[n])
        K2 = dt*f(u[n] + 0.5*K1, t[n] + dt2)
        K3 = dt*f(u[n] + 0.5*K2, t[n] + dt2)
        K4 = dt*f(u[n] + K3, t[n] + dt)
        u_new = u[n] + (1/6.0)*(K1 + 2*K2 + 2*K3 + K4)
        return u_new

The quick_description string is needed for the class to appear in the package’s table of contents.

Implicit solver: See examples in the solver.py file, class BackwardEuler, for instance.

Adaptive solver: See examples in the solver.py file, class RKFehlberg, for instance.

Wrapping other packages. This section is for developers who intend to wrap an existing package and integrate with the present one..

If the original package is not written in Python language, developers need to apply some specific tools (like F2PY, Cython, or SWIG) to create an extension module to make the package accessible from Python code.

On the other hand, if the original package is also a Python software, 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).

By an attempt to import these necessary modules (often set in method initialize()), we can check whether the necessary dependencies are installed properly.

Definition of parameters and their properties. Each solver has a set of specific parameters depending on its underlying method. For example, adaptive solvers will be more likely to apply (many) attributes for step-size control. When integrating a new method, first search in the _parameters dictionary in solver.py if some parameters can be reused for the new solver. New parameters are added to the _parameters dictionary, either in the solver.py or by importing solver and updating solver._parameters, see rkc.py for an example.

Each solver class lists the required and optional parameters it needs in the class variables _optional_parameters and _required_parameters. Very often these lists are inherited, or just a few new parameters are added to the list already defined in the superclass.

Sometimes values of parameters (or other properties) need to be changed in a solver, e.g., because there are certain relations between various parameters. Appropriate adjustments and checks are done in the method initialize_for_solve, which is called in the beginning of the “solve” process (before any computations take place). Many classes provide examples on this, e.g., class RKC in rkc.py. This particular class shows how to generate input parameters to the Fortran code that are not given by the user, but automatically derived in the Python code from other data.

Another method that is called from solve is validate_data. The solver class can use this method to check for consistency of data structures before starting the numerical computations. As an example, the class Lsodes in the odepack.py file imposes relations between the input data: two input integer arrays ia and ja must be input simultaneously, len(ia) == neq + 1, and ia[neq] = len(ja). These checks are done in Python before calling the Fortran solver.

The solve and advance methods. Simple methods can just implement advance to bring the solution one step forward in time. The existing solve method derived from the superclass is general enough to administer the whole solution process.

Adaptive solvers will typically think of advance as bringing the solution to the next user-desired time level, using an unknown set of smaller time steps whose sizes must be determined. Then there will hence be a time loop within advance, while the outer time loop in solve takes care of the stepping between the user-desired time levels. (The simplest methods just takes one step between the user-desired time levels.)

When wrapping solvers in external software, it is occasionally not feasible to implement advance because one wants to rely on the 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 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.

class odespy.solvers.AdamsBashMoulton2(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Two-step (3rd-order) Adams-Bashforth method:

predictor = u[n] + dt/12.*(23.*f(u[n], t[n]) - 16*f(u[n-1], t[n-1]) +
                    5*f(u[n-2], t[n-2]))
corrector = u[n] + dt/12.*(8.*f(u[n], t[n]) - f(u[n-1], t[n-1]) +
                    5*f(predictor, t[n+1]))

for constant time step dt.

RK2 is used as default solver for first steps.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
start_method Method for the first steps in multi-step solvers. (default: RK2)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data()
__module__ = 'odespy.solvers'
advance()[source]
initialize_for_solve()[source]
quick_description = 'Explicit 2nd-order Adams-Bashforth-Moulton method'
validate_data()[source]
class odespy.solvers.AdamsBashMoulton3(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Three-step (4th-order) Adams-Bashforth method:

predictor = u[n] + dt/24.*(55.*f(u[n], t[n]) - 59*f(u[n-1], t[n-1]) +
                           37*f(u[n-2], t[n-2]) - 9*f(u[n-3], t[n-3]))
corrector = u[n] + dt/24.*(19.*f(u[n], t[n]) - 5*f(u[n-1], t[n-1]) +
                           f(u[n-2], t[n-2]) + 9*f(predictor, t[n+1]))

for constant time step dt.

RK2 is used as default solver for first steps.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
start_method Method for the first steps in multi-step solvers. (default: RK2)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data()
__module__ = 'odespy.solvers'
advance()[source]
initialize_for_solve()[source]
quick_description = 'Explicit 3rd-order Adams-Bashforth-Moulton method'
validate_data()[source]
class odespy.solvers.AdamsBashforth2(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Second-order Adams-Bashforth method:

u[n+1] = u[n] + dt/2.*(3*f(u[n], t[n]) - f(u[n-1], t[n-1]))

for constant time step dt.

RK2 is used as default solver in first step.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
start_method Method for the first steps in multi-step solvers. (default: RK2)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() Check that the time steps are constant.
__module__ = 'odespy.solvers'
advance()[source]
initialize_for_solve()[source]
quick_description = 'Explicit 2nd-order Adams-Bashforth method'
validate_data()[source]

Check that the time steps are constant.

class odespy.solvers.AdamsBashforth3(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Third-order Adams-Bashforth method:

u[n+1] = u[n] + dt/12.*(23*f(u[n], t[n]) - 16*f(u[n-1], t[n-1])
                        + 5*f(u[n-2], t[n-2]))

for constant time step dt.

RK2 is used as default solver for first steps.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
start_method Method for the first steps in multi-step solvers. (default: RK2)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data()
__module__ = 'odespy.solvers'
advance()[source]
initialize_for_solve()[source]
quick_description = 'Explicit 3rd-order Adams-Bashforth method'
validate_data()[source]
class odespy.solvers.AdamsBashforth4(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Fourth-order Adams-Bashforth method:

u[n+1] = u[n] + dt/24.*(55.*f(u[n], t[n]) - 59*f(u[n-1], t[n-1]) +
                        37*f(u[n-2], t[n-2]) - 9*f(u[n-3], t[n-3]))

for constant time step dt.

RK2 is used as default solver for first steps.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
start_method Method for the first steps in multi-step solvers. (default: RK2)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data()
__module__ = 'odespy.solvers'
advance()[source]
initialize_for_solve()[source]
quick_description = 'Explicit 4th-order Adams-Bashforth method'
validate_data()[source]
class odespy.solvers.Adaptive(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Superclass for adaptive solvers. Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance() Advance solution one time step.
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
initialize_for_solve()[source]
class odespy.solvers.AdaptiveResidual(f, **kwargs)[source]

Bases: odespy.solvers.Adaptive

Designed for educational purposes to demonstrate a possible adaptive strategy.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.
solver Name of solver class in solvers that need an extra solver (e.g., AdaptiveResidual). (default: RK4)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance() Advance solution one time step.
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
residual(u_n, u_next, t_n, t_next)
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate])
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__init__(f, **kwargs)[source]
__module__ = 'odespy.solvers'
quick_description = 'Very simple adaptive strategy based on the residual'
residual(u_n, u_next, t_n, t_next)[source]
solve(time_points, terminate=None)[source]
class odespy.solvers.Backward2Step(f, **kwargs)[source]

Bases: odespy.solvers.SolverImplicit

Implicit Backward Euler method with 2 steps:

u[n+1] = u[n]*4/3 - u[n-1]/3 + (t[n+1-t[n-1]])*f(t[n+1], u[n+1])/3

The 1st-order Backward Euler method is used for the first step.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
h_in_fd_jac h in finite difference approximation of the Jacobian. (default: 0.0001)
nonlinear_solver Newton or Picard nonlinear solver. (default: Picard)
max_iter Max no of iterations in nonlinear solver (default: 25)
eps_iter Max error measure in nonlinear solver (default: 0.0001)
relaxation relaxation argument (r): new_solution = r*solution + (1-r)*old_solution (default: 1.0)

Methods

Newton_system(ukp1)
Picard_update(ukp1)
adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
Newton_system(ukp1)[source]
Picard_update(ukp1)[source]
__module__ = 'odespy.solvers'
quick_description = 'Implicit 2nd-order Backward Euler method'
class odespy.solvers.BackwardEuler(f, **kwargs)[source]

Bases: odespy.solvers.SolverImplicit

Implicit Backward Euler method:

u[n+1] = u[n] + dt*f(t[n+1], u[n+1])

The nonlinear system is solved by Newton or Picard iteration.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
h_in_fd_jac h in finite difference approximation of the Jacobian. (default: 0.0001)
nonlinear_solver Newton or Picard nonlinear solver. (default: Picard)
max_iter Max no of iterations in nonlinear solver (default: 25)
eps_iter Max error measure in nonlinear solver (default: 0.0001)
relaxation relaxation argument (r): new_solution = r*solution + (1-r)*old_solution (default: 1.0)

Methods

Newton_system(ukp1)
Picard_update(ukp1)
adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
Newton_system(ukp1)[source]
Picard_update(ukp1)[source]
__module__ = 'odespy.solvers'
quick_description = 'Implicit 1st-order Backward Euler method'
class odespy.solvers.Dop853(f, **kwargs)[source]

Bases: odespy.solvers.ode_scipy

Wrapper for scipy.integrate.ode.dop853, which applies the Dormand&Prince method of order 8(5,3), based on the Fortran implementation by Hairer and Wanner. See http://www.unige.ch/~hairer/software.html.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
nsteps Max no of internal solver steps per time step. (default: 1000)
ifactor Maximum factor for increasing the step size (default: 2)
dfactor Maximum factor for decreasing the step size (default: 0.5)
beta Beta argument for stabilized step size control in Dormand&Prince methods from scipy
safety Safety factor on new step selection (default: 0.9)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize()
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
quick_description = 'Adaptive Dormand & Prince method of order 8(5,3) (scipy)'
class odespy.solvers.Dopri5(f, **kwargs)[source]

Bases: odespy.solvers.ode_scipy

Wrapper for scipy.integrate.ode.dopri5, which applies the Dormand&Prince method of order 5(4), based on the Fortran implementation by Hairer and Wanner. See http://www.unige.ch/~hairer/software.html.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
nsteps Max no of internal solver steps per time step. (default: 1000)
ifactor Maximum factor for increasing the step size (default: 2)
dfactor Maximum factor for decreasing the step size (default: 0.5)
beta Beta argument for stabilized step size control in Dormand&Prince methods from scipy
safety Safety factor on new step selection (default: 0.9)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize()
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
quick_description = 'Dormand & Prince method of order 5(4) (scipy)'
odespy.solvers.Euler

alias of ForwardEuler

class odespy.solvers.ForwardEuler(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Forward Euler scheme:

u[n+1] = u[n] + dt*f(u[n], t[n])

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
quick_description = 'The simple explicit (forward) Euler scheme'
class odespy.solvers.Heun(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Heun’s method, also known as an RungeKutta2 or Trapezoidal method. Basically, it is a central difference method, with one iteration and the Forward Euler scheme as start value. In this sense, it is a predictor-corrector method.

Scheme:

u[n+1] = u[n] + 0.5*dt*(f(u[n],t[n]) + f(u[n]+dt*f(u[n],t[n]),t[n+1]))

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
quick_description = "Heun's explicit method (similar to RK2)"
class odespy.solvers.Leapfrog(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Leapfrog scheme:

u[n+1] = u[n-1] + dt2*f(u[n], t[n])

with:

dt2 = t[n+1] - t[n-1]

Forward Euler is used for the first step.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
quick_description = 'Standard explicit Leapfrog scheme'
class odespy.solvers.LeapfrogFiltered(f, **kwargs)[source]

Bases: odespy.solvers.Solver

The standard Leapfrog scheme reads:

u[n+1] = u[k-1] + dt2*f(u[n], t[n])

with:

dt2 = t[n+1] - t[k-1]

Forward Euler is used for the first step. Since Leapfrog gives oscillatory solutions, this class applies a common filtering technique:

u[n] = u[n] + gamma*(u[n-1] - 2*u[n] + u[n+1])

with gamma=0.6 as in the NCAR Climate Model.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
quick_description = 'Filtered Leapfrog scheme'
class odespy.solvers.MidpointImplicit(f, **kwargs)[source]

Bases: odespy.solvers.SolverImplicit

Midpoint Implicit method:

u[n+1] = u[n] + dt*f((u[n+1] + u[n])/2., t[n] + dt/2.)

The nonlinear system is solved by Picard or Newton iteration.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
h_in_fd_jac h in finite difference approximation of the Jacobian. (default: 0.0001)
nonlinear_solver Newton or Picard nonlinear solver. (default: Picard)
max_iter Max no of iterations in nonlinear solver (default: 25)
eps_iter Max error measure in nonlinear solver (default: 0.0001)
relaxation relaxation argument (r): new_solution = r*solution + (1-r)*old_solution (default: 1.0)

Methods

Newton_system(ukp1)
Picard_update(ukp1)
adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
Newton_system(ukp1)[source]
Picard_update(ukp1)[source]
__module__ = 'odespy.solvers'
quick_description = 'Implicit 2nd-order Midpoint method'
class odespy.solvers.MidpointIter(f, **kwargs)[source]

Bases: odespy.solvers.Solver

A midpoint/central difference method with max_iter fixed-point iterations to solve the nonlinear system. The Forward Euler scheme is recovered if max_iter=1 and f(u,t) is independent of t. For max_iter=2 we have the Heun/RK2 scheme.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
max_iter Max no of iterations in nonlinear solver (default: 25)
eps_iter Max error measure in nonlinear solver (default: 0.0001)

Methods

adjust_parameters()
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
adjust_parameters()[source]
advance()[source]
quick_description = 'Explicit 2nd-order iterated Midpoint method'
class odespy.solvers.MySolver(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Users can define a solver with supplying a function myadvance(), and make use of all possible parameters in this module:

myadvance(MySolver_instance)  -->  return u_new

Example:

def myadvance_(ms):
    f, u, t, n, atol = ms.f, ms.u, ms.t, ms.n, ms.atol
    # All class attributes can be obtained
    u_new = ...
    return u_new

def f(u,t):
    udot = ...
    return udot

method = MySolver(f, myadvance=myadvance_)
method.set_initial_condition(u0)
u,t = method.solve(time_points)

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.
myadvance User supplied function to advance current solution one step forward. See documents of class MySolver.

Optional input arguments:

Name Description
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
nsteps Max no of internal solver steps per time step. (default: 1000)
verbose Integer reflecting output of intermediate quantities. (default: 0)
method_order Method order for user-defined method if known.A integer for 1-level methods, or a pair of integer for 2-levels methods.
g_f77

Fortran subroutine for g. This subroutine has the signature:

      subroutine g_f77(neq, t, u, ng, groot)
Cf2py intent(hide) neq
Cf2py optional, intent(hide) ng
Cf2py intent(in) t, u
Cf2py intent(out) groot
      integer neq, ng
      double precision t, u, groot
      dimension u(neq), groot(ng)
      groot = ...
      return
      end
use_special Switch for using special times
ia Integer array with length neq+1 which contains starting locations in ja of the descriptions for columns 1...neq. ia(1) == 1. The last element ia[neq+1] should equal to the total number of nonzero locations assumed. For each column j = 1...neq, the values of the row index i in column j, where a nonzero element may occur, are given by i == ja(k) where ia(j) <=’ k < ia(j+1).
dfactor Maximum factor for decreasing the step size (default: 0.5)
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
u_exact Function of t returning exact solution. (default: None)
mb Block size. Describe the block-tridiagonal form of matrix A together with nb.
myadvance User supplied function to advance current solution one step forward. See documents of class MySolver.
res User-supplied function to calculate the residual vector, defined by r = g(u,t) - A(u,t) * s. Used in the linearly implicit solvers: Lsodi, Lsodis, Lsoibt. The res function has the signature res(u,t,s,ires) and returns the tuple (r,ires), where ires is an int. On input, ires indicates how ODEPACK would use the returned array “r”: ires=1 means the full residual exactly, ires=2 means that r is used only to compute the Jacobian dr/du by finite differences. res should set the flag ires if it encounters a halt condition or illegal input. Otherwise, it should not be reset. On output, the value 1 or -1 represents a normal return.
nb Number of blocks in the main diagonal. In each of the nb block-rows of the matrix P (each consisting of mb consecutive rows), the nonzero elements are to lie in three consecutive mb by mb blocks. In block-rows 2 through nb-1, these are centered about the main diagonal. In block rows 1 and nb, they are the diagonal blocks and the two blocks adjacent to the diagonal block. (Thus block positions (1,3) and (nb,nb-2) can be nonzero.) Require: mb>=1, nb>=4, mb*nb==neq.
h_in_fd_jac h in finite difference approximation of the Jacobian. (default: 0.0001)
ng No of components in constraint function g.
first_step Suggested first time step size for an adaptive algorithm.
rtol Relative tolerance for solution. (default: 1e-06)
safety Safety factor on new step selection (default: 0.9)
adams_or_bdf Method in vode or solvers in odepack: “adams” or “bdf”. (default: adams)
stiff Boolean flag to indicate stiffness.
butcher_tableau 2d-array which contains the butcher table for user- supplied Runge-Kutta method. (n,n) array for 1-level Runge-Kutta methods.(n+1,n) array for 2-level Runge- Kutta methods.
ode_method solver type: “adams” or “bdf” (default: adams)
theta Weight in [0,1] used for “theta-rule” finite difference approx. (default: 0.5)
order Maximum order used by the integrator (<= 12 for “adams”, <= 5 for “bdf”). (default: 4)
complex_valued True if f is complex valued. (default: False)
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_constant Flag to show whether Jacobian is constant, 0 (false) or 1 (true) (default: 0)
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
min_step Minimum step size for an adaptive algorithm.
max_iter Max no of iterations in nonlinear solver (default: 25)
f_f77

Fortran subroutine for f. This subroutine has the signature:

      subroutine f_f77(neq,t,u,udot)
Cf2py intent(hide)   neq
Cf2py intent(out)    udot
      integer neq
      double precision t,u(neq),udot(neq)
      udot = ...
      return
      end
jac_banded Function for Jacobian in banded matrix form. Used in Lsode, Lsoda, Lsodar. jac_banded(u,t,ml,mu) returns df/du as an array of size neq times ml+mu+1.
beta Beta argument for stabilized step size control in Dormand&Prince methods from scipy
atol Absolute tolerance for solution. (default: 1e-08)
relaxation relaxation argument (r): new_solution = r*solution + (1-r)*old_solution (default: 1.0)
odelab_solver Name of Solver class in odelab (default: RungeKutta4)
solver Name of solver class in solvers that need an extra solver (e.g., AdaptiveResidual). (default: RK4)
ja Integer array containing the row indices where the nonzero elements occur, in columnwise order. Describes the sparsity matrix structure together with ia. In Lsodes, ia and ja describe the structure of Jacobian matrix; while in Lsodis, ia and ja are used to describe the structure of matrix A.
init_step Fixed step size for time mesh.
max_step Maximum step size for an adaptive algorithm.
g Callable object to define constraint functions. g(u, t) returns a vector of the values of the constraints (left-hand sides in the constraint equations).
f Right-hand side f(u,t) defining the ODE.
start_method Method for the first steps in multi-step solvers. (default: RK2)
ml Number of lower non-zero diagonals in a banded Jacobian.
ifactor Maximum factor for increasing the step size (default: 2)
eps_iter Max error measure in nonlinear solver (default: 0.0001)
mu Number of upper non-zero diagonals in a banded Jacobian.
jac_f77

Fortran subroutine for jac. This subroutine has the signature:

      subroutine jac_f77
     1 (neq, t, u, ml, mu, pd, nrowpd)
Cf2py intent(hide) neq,ml,mu,nrowpd
Cf2py intent(out) pd
      integer neq,ml,mu,nrowpd
      double precision t,u,pd
      dimension u(neq),pd(neq,neq)
      pd = ...
      return
      end
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
ydoti Real array for the initial value of dy/dt. (default: [])
strictdt Uniform time mesh vs exact dt spacings (default: True)
nonlinear_solver Newton or Picard nonlinear solver. (default: Picard)
specialtimes List of special times to use during iteration
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
class odespy.solvers.RK2(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Standard Runge-Kutta 2nd method:

u[n+1] = u[n] + dt*f(u[n] + 0.5*(dt*f(u[n],t[n])),t[n] + 0.5*dt)

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
quick_description = 'Explicit 2nd-order Runge-Kutta method'
class odespy.solvers.RK3(f, **kwargs)[source]

Bases: odespy.solvers.Solver

RungeKutta3 method:

u[n+1] = u[n] + (1/6.0)*(K1 + 4*K2 + K3)

where:

K1 = dt*f(u[n], t[n])
K2 = dt*f(u[n] + 0.5*K1, t[n] + 0.5*dt)
K3 = dt*f(u[n] - K1 + 2*K2, t[n] + dt)

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
quick_description = 'Explicit 3rd-order Runge-Kutta method'
class odespy.solvers.RK34(*args, **kwargs)[source]

Bases: odespy.solvers.Adaptive

NOTE: This class does not work!

Adaptive 4th-order Runge-Kutta method. For each time level t[n], the method takes many adaptively chosen (small) time steps to hit the next target level t[n+1]. All computed u and t values are available as self.u_adaptive and self.t_adaptive, if desired.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance() Advance from t[n] to t[n+1] in (small) adaptive steps.
advance_intermediate(u, t, dt) Advance the solution an intermediate addaptive step.
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__init__(*args, **kwargs)[source]
__module__ = 'odespy.solvers'
advance()[source]

Advance from t[n] to t[n+1] in (small) adaptive steps.

Adaptivity logic: Define h as the local adaptive time step and dt = t[1]-t[0]. Start with h as first_step if first_step <= dt, otherwise h = dt. At each user-specified time level: step h forward, estimate error, accept solution if error less than tolerance, adjust h, retry with new h until solution is less than error. Then proceed with a new step. Continue until next user-specified time level.

advance_intermediate(u, t, dt)[source]

Advance the solution an intermediate addaptive step. Parmeters: u at time t, dt is current step size.

initialize_for_solve()[source]
quick_description = 'Adaptive 4th-order Runge-Kutta method'
class odespy.solvers.RK4(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Standard RK4 method:

u[n+1] = u[n] + (1/6.0)*(K1 + 2*K2 + 2*K3 + K4)

where:

K1 = dt*f(u[n], t[n])
K2 = dt*f(u[n] + 0.5*K1, t[n] + 0.5*dt)
K3 = dt*f(u[n] + 0.5*K2, t[n] + 0.5*dt)
K4 = dt*f(u[n] + K3, t[n] + dt)

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
quick_description = 'Explicit 4th-order Runge-Kutta method'
class odespy.solvers.RKFehlberg(f, **kwargs)[source]

Bases: odespy.solvers.Adaptive

The classical adaptive Runge-Kutta-Fehlberg method of order 4-5. The method is also available in more sophisticated form as class Fehlberg (in the RungeKutta module).

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
quick_description = 'Adaptive Runge-Kutta-Fehlberg (4,5) method'
class odespy.solvers.Solver(f, **kwargs)[source]

Superclass for numerical methods solving ODE problem

u’(t) = f(u, t), u(0) = U0

where u and U0 are scalars (for scalar ODEs) or vectors (for systems of ODEs).

Attributes stored in this class:

Name Description
u array of point values of the solution function
t array of time values: u[i] corresponds to t[i]
n the most recently computed solution is u[n+1]
f function wrapping the user’s right-hand side f(u, t), used in all algorithms
users_f the user’s original function implementing f(u, t)
PRM an attribute for each optional and required parameter

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance() Advance solution one time step.
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve() Setting values of internal attributes to be used in iteration.
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__init__(f, **kwargs)[source]

f is the right-hand side function of the ODE u’ = f(u,t). The legal keyword arguments (in kwargs) are documented in the tables in the doc string of this class. The f function must return a float or complex object in case of a scalar ODE and a list or array of float or complex objects in case of a system of ODEs.

This constructor makes a dictionary self._parameters holding all the required and optional parameters for this solver (fetched from the global _parameters dictionary in this module). The method adjust_parameters (implemented in subclasses) is called to adjust default parameter settings if needed. Then all keys in self._parameters become class attributes, filled with default values. Thereafter, all keyword arguments (in kwargs) with None as value are removed as keyword arguments. The next step is to call set(**kwargs), i.e., use the keyword arguments to modify the values of the attributes that represent the parameters in this solver. Finally, the constructor calls the method initialize (to be implemeneted in subclasses, e.g., for importing necessary modules for the solver).

Instead of supplying keyword arguments to this constructor, the user can at any time call the set method with keyword arguments in order to specify parameters.

__module__ = 'odespy.solvers'
__repr__()[source]

Return solvername(f=..., param1=..., etc.).

__str__()[source]

Return solver name, plus parameters that are different from their default values.

adjust_parameters()[source]

This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary. The method is called from the constructor.

Further adjustments of self._parameters can be done in initialize_for_solve when all data for the solver are available.

advance()[source]

Advance solution one time step.

check_conditional_parameters()[source]

This function is used to check whether conditional parameters are provided when specified condition fulfilled.

This function is not intended for simple solvers. So it is not called automatically in current solvers.py. But for some complicated solvers as ones in ODEPACK, they are very useful and convenient.

Future developers can apply these functions at appropriate locations with corresponding property-setting in adjust_parameters().

For example, in Lsode_ODEPACK, when iter_method is set to 4, it indicates that ODEPACK would apply user-supplied banded Jacoabian function in corrector iteration. Then we need to confirm either ‘jac_banded’ or ‘jac_fortran’ is supplied. Besides, ‘ml’ & ‘mu’ are also necessary for iteration with banded Jacobian matrix. Thus in order to confirm sufficient conditional inputs, we set parameters[‘iter_method’][‘condition_list’] =

{...,‘4’: ((‘jac_banded’,’jac_fortran’),ml,mu),...}

In this function, we would check all the legal parameters with specified condition-list, and make sure all the conditional parameters with current value is supplied.

check_extra(**kwargs)[source]

A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter. This method runs the user-given function(s) on the relevant set of parameters.

check_input_range(**kwargs)[source]

Check whether all existing inputs are in right specified range.

check_input_types(**kwargs)[source]

Check whether all existing inputs are of right specified type.

compile_string_functions(f, **kwargs)[source]

Compile functions which are supplied as Fortran strings.

constant_time_step()[source]

Check if self.t has a uniform partition.

get(parameter_name=None, print_info=False)[source]

Return value of specified input parameters. If parameter_name is None, return dict of all inputs.

get_parameter_info(print_info=False)[source]

Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).

If print_info is True, the self._parameters dict is pretty printed, otherwise it is returned.

has_u_t_all()[source]

Return True if self.u_all and self.t_all, defined in subclasses derived from Adaptive, are computed. This function for testing is placed in class Solver so that any class can check of intermediate time steps have been computed and are available.

initialize()[source]

Subclass-specific initialization. Called from constructor. Typical use: import modules needed in methods in the class and provide error messages if modules are not installed.

initialize_for_solve()[source]

Setting values of internal attributes to be used in iteration.

These internal attributes are ususally dependent on the values of other attributes. For example, for Rkc, self.itol should be initialized here as a flag to indicate whether self.atol is supplied as scalar or sequence.

In subclasses, this function can be extended when required.

set(strict=False, **kwargs)[source]

Assign values to one or more parameters, specified as keyword arguments.

The legal parameters that can be set are contained in the dict self._parameters.

If strict is true, only registered parameter names are accepted, otherwise unregistered parameters are ignored.

The provided parameters (keyword arguments in kwargs) are first checked for legal type and legal range.

Types and ranges of attributes are defined in self._parameters, which is initialized with default settings and optionally modified in the adjust_parameters method.

set_initial_condition(U0)[source]

Function set_initial_condition() is used to set initial value of independent variables.

solve(time_points, terminate=None)[source]

Compute discrete solution u of the ODE problem at time points specified in the array time_points. An optional user-supplied function terminate(u, t, step_no) can be supplied to terminate the solution process (terminate returns True or False) at some time earlier than time_points[-1].

Most classes in this solver hierarchy inherit this solve method and implement their special advance method to advance the solution one step. Some solver classes will implement their own solve method, for instance if they wrap some underlying software that has a suitable solve functionality.

The algorithm steps in this solve method goes as follows. The initialize_for_solve method is called to initialize various data needed in the solution process (self. u, for instance). Thereafter, validate_data is called to perform a consistency check on data. We are then ready for the core of the method: the time loop.

Output:
u : array to hold solution values corresponding to points t : array to hold time values.Usually same as time_points
switch_to(solver_target, print_info=False, **kwargs)[source]

Create a new solver instance which switch to another subclass with same values of common attributes.

solver_target is either as a string (class name) or a class, i.e., RK4 or 'RK4'. The kwargs arguments are optional parameters to reset/supplement values of argmenters in the solver we switch to. The instance of the target solver is returned.

Example:

>>> import odespy
>>> f = lambda u,t: -u
>>> time_points = np.linspace(0.,2.,11)
>>> exact_u = np.exp(-time_points)
>>> m1 = odespy.RK2(f)
>>> m1.set_initial_condition(1.)
>>> u1, t = m1.solve(time_points)
>>> print 'Normarized error with RK2 is %g' % np.linalg.norm(u1 - exact_u)
Normarized error with RK2 is 0.0077317
>>> m2 = m1.switch_to(odespy.RKFehlberg, rtol=1e-18)
>>> u2, t = m2.solve(time_points)
>>> print 'Normarized error with RKFehlberg is %g' % np.linalg.norm(u2 - exact_u)
Normarized error with RKFehlberg is 8.55517e-08
validate_data()[source]

This function is used for extra checking and validating of attributes before the computations start.

This version checks that the time_points is correctly set up. The function also check that all required parameters are initialized. Subclass versions may introduce additional tests and help.

class odespy.solvers.SolverImplicit(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Super class for implicit methods for ODEs. Existing solvers are: BackwardEuler, Backward2Step, ThetaRule

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
h_in_fd_jac h in finite difference approximation of the Jacobian. (default: 0.0001)
nonlinear_solver Newton or Picard nonlinear solver. (default: Picard)
max_iter Max no of iterations in nonlinear solver (default: 25)
eps_iter Max error measure in nonlinear solver (default: 0.0001)
relaxation relaxation argument (r): new_solution = r*solution + (1-r)*old_solution (default: 1.0)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
initialize_for_solve()[source]
class odespy.solvers.ThetaRule(f, **kwargs)[source]

Bases: odespy.solvers.SolverImplicit

This class implements the implicit theta-rule method for solving ODEs. The method includes a parameter theta used to weight the previous and new time step when evaluating the right-hand side:

u[n+1] = u[n] + dt*(theta*f(u[n+1],t[n+1]) + (1 - theta)*f(u[n],t[n]))

Here, theta is a float in [0,1].

The nonlinear system is solved by Picard or Newton iteration.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
h_in_fd_jac h in finite difference approximation of the Jacobian. (default: 0.0001)
nonlinear_solver Newton or Picard nonlinear solver. (default: Picard)
max_iter Max no of iterations in nonlinear solver (default: 25)
eps_iter Max error measure in nonlinear solver (default: 0.0001)
relaxation relaxation argument (r): new_solution = r*solution + (1-r)*old_solution (default: 1.0)
theta Weight in [0,1] used for “theta-rule” finite difference approx. (default: 0.5)

Methods

Newton_system(ukp1)
Picard_update(ukp1)
adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize() Subclass-specific initialization.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
Newton_system(ukp1)[source]
Picard_update(ukp1)[source]
__module__ = 'odespy.solvers'
quick_description = 'Unified Forward/Backward Euler and Midpoint methods'
odespy.solvers.Trapezoidal

alias of Heun

class odespy.solvers.Vode(f, **kwargs)[source]

Bases: odespy.solvers.ode_scipy

Wrapper for scipy.integrate.ode.vode, which is a wrapper for vode.f, which intends to solve initial value problems of stiff or nonstiff type. The well-known vode.f solver applies backward differential formulae for iteration.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
nsteps Max no of internal solver steps per time step. (default: 1000)
adams_or_bdf Method in vode or solvers in odepack: “adams” or “bdf”. (default: adams)
min_step Minimum step size for an adaptive algorithm.
order Maximum order used by the integrator (<= 12 for “adams”, <= 5 for “bdf”). (default: 4)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize()
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
initialize_for_solve()[source]
quick_description = 'Adams/BDF Vode adaptive method (vode.f wrapper)'
odespy.solvers.approximate_Jacobian(f, u0, t0, h)[source]

Compute approximate Jacobian of fucntion f at current point (u0,t0). Method: forward finite difference approximation with step size h. Output: a two-dimensional array holding the Jacobian matrix.

odespy.solvers.compile_f77(*function_strings, **kwargs)[source]

Compile a list of strings, each string containing a Fortran 77 subroutine (typically for f, jac, etc.). Return the extension module if keyword argument return_module=True, otherwise return the subroutines/functions corresponding to the given strings (default). The name of the extension module is tmp_callback and the corresponding file (to be cleaned up) is tmp_callback.so.

odespy.solvers.list_all_solvers()[source]

Return all solver classes in this package, excluding superclasses.

odespy.solvers.list_available_solvers()[source]

Return all available solver classes in this package.

odespy.solvers.list_not_suitable_complex_solvers()[source]

List all solvers that do not work for complex-valued ODEs.

class odespy.solvers.lsoda_scipy(f, **kwargs)[source]

Bases: odespy.solvers.Adaptive

Wrapper of the scipy.integrate.odeint solver so that it can be called from the odespy user interface.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
ml Number of lower non-zero diagonals in a banded Jacobian.
mu Number of upper non-zero diagonals in a banded Jacobian.

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance() Advance solution one time step.
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize()
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate])
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
initialize()[source]
quick_description = 'Wrapper of lsoda (scipy.integrate.odeint)'
solve(time_points, terminate=None)[source]
class odespy.solvers.ode_scipy(f, **kwargs)[source]

Bases: odespy.solvers.Adaptive

Superclass wrapper for scipy.integrate.ode classes. Existing solvers in subclasses are: Vode, Dopri5, Dop853.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
nsteps Max no of internal solver steps per time step. (default: 1000)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance()
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize()
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) Compute discrete solution u of the ODE problem at time points
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
advance()[source]
initialize()[source]
initialize_for_solve()[source]
class odespy.solvers.odefun_sympy(f, **kwargs)[source]

Bases: odespy.solvers.Solver

Wrapper for the sympy.mpmath.odefun method, which applies a high-order Taylor series method to solve ODEs.

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance() Advance solution one time step.
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize()
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) The complete solve method must be overridded in this class since sympy.mpmath.odefun is such a solve method.
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
__module__ = 'odespy.solvers'
initialize()[source]
initialize_for_solve()[source]
quick_description = 'Very accurate high order Taylor method (from SymPy)'
solve(time_points, terminate=None)[source]

The complete solve method must be overridded in this class since sympy.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)).

class odespy.solvers.odelab(f, **kwargs)[source]

Bases: odespy.solvers.Adaptive

Odespy wrapper for the odelab package. Typical use:

f = lambda u: -u   # u' = -u
import odespy, numpy
solver = odespy.odelab(f, odelab_solver='Butcher')
solver.set_initial_condition(1.0)
u, t = solver.solve(numpy.linspace(0, 3, 21))

Odelab offers the following solvers (for parameter odelab_solver): ‘ExplicitEuler’, ‘ExplicitTrapezoidal’, ‘ImplicitEuler’, ‘RungeKutta34’, ‘RungeKutta4’, ‘ImplicitEuler’, ‘AdamsBashforth1’, ‘Heun’, ‘Kutta38’, ‘Kutta4’

Note: ode15s is actually scipy’s vode code called with BDF order=5. Odelab’s DAE methods are not yet supported (could be by proper wrapping of f and g in an odelab.system.System object).

Required input arguments:

Name Description
f Right-hand side f(u,t) defining the ODE.
odelab_solver Name of Solver class in odelab (default: RungeKutta4)

Optional input arguments:

Name Description
f_args Extra positional arguments to f: f(u, t, *f_args, **f_kwargs). (default: ())
f_kwargs Extra keyword arguments to f: f(u, t, *f_args, **f_kwargs). (default: {})
complex_valued True if f is complex valued. (default: False)
disk_storage Indicates whether u is stored in memory or in file. If string, it is the filename; if False or “”, u is kept in memory; if True, a default filename tmp_odspy.dat is used. (default: False)
verbose Integer reflecting output of intermediate quantities. (default: 0)
u_exact Function of t returning exact solution. (default: None)
rtol Relative tolerance for solution. (default: 1e-06)
atol Absolute tolerance for solution. (default: 1e-08)
first_step Suggested first time step size for an adaptive algorithm.
min_step Minimum step size for an adaptive algorithm.
max_step Maximum step size for an adaptive algorithm.
jac Jacobian of right-hand side function f (df/du). (default: None)
jac_args Extra positional arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: ())
jac_kwargs Extra keyword arguments to jac: jac(u, t, *jac_args,**jac_kwargs). (default: {})

Methods

adjust_parameters() This method allows subclasses to adjust (modify or add) entries in the self._parameters dictionary.
advance() Advance solution one time step.
check_conditional_parameters() This function is used to check whether conditional parameters are provided when specified condition fulfilled.
check_extra(**kwargs) A parameter may have a keyword extra_check for user-given functions that performs consistency checks on the parameter.
check_input_range(**kwargs) Check whether all existing inputs are in right specified range.
check_input_types(**kwargs) Check whether all existing inputs are of right specified type.
compile_string_functions(f, **kwargs) Compile functions which are supplied as Fortran strings.
constant_time_step() Check if self.t has a uniform partition.
get([parameter_name, print_info]) Return value of specified input parameters.
get_parameter_info([print_info]) Return a dictionary containing all properties of all legal parameters in current subclass (i.e., the parameters in self._parameters).
has_u_t_all() Return True if self.u_all and self.t_all, defined in
initialize()
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
solve(time_points[, terminate]) The complete solve method must be overridded in this class since odelab.
switch_to(solver_target[, print_info]) Create a new solver instance which switch to another subclass with same values of common attributes.
validate_data() This function is used for extra checking and validating of attributes before the computations start.
DAE_solvers = ['SymplicticEuler', 'MultiRKDAE', 'RKDAE', 'Spark', 'Spark2', 'Spark3', 'Spark4']
__module__ = 'odespy.solvers'
initialize()[source]
initialize_for_solve()[source]
not_valid = ['SymplecticEuler', 'LDIRK343', 'LobattoIIIA', 'LobattoIIIB', 'LobattoIIIC', 'LobattoIIICs', 'LobattoIIID', 'MultiRKDAE', 'RKDAE', 'RadauIIA', 'RungeKutta', 'Spark', 'Spark2', 'Spark3', 'Spark4', 'AdamsBashforth', 'AdamsBashforth2', 'AdamsBashforth2e', 'Butcher', 'Butcher1', 'Butcher3', 'ExplicitGeneralLinear', 'GeneralLinear', 'Kutta', 'McLachlan', 'NonHolonomicEnergy', 'NonHolonomicEnergy0', 'NonHolonomicEnergyEMP', 'NonHolonomicLeapFrog', 'SymplecticEuler', 'ABLawson', 'ABLawson2', 'ABLawson3', 'ABLawson4', 'ABNorset4', 'Exponential', 'GenLawson45', 'HochOst4', 'Lawson4', 'LawsonEuler', 'Phi', 'Polynomial', 'RKMK4T']
quick_description = 'interface to all solvers in odelab'
solve(time_points, terminate=None)[source]

The complete solve method must be overridded in this class since odelab. is such a solve method.

solver = 'RKMK4T'
solvers = ['ExplicitEuler', 'ExplicitTrapezoidal', 'ImplicitEuler', 'RungeKutta34', 'RungeKutta4', 'ImplicitEuler', 'AdamsBashforth1', 'Heun', 'Kutta38', 'Kutta4']
odespy.solvers.table_of_parameters(classname)[source]

Return a table (in reST format) of the required parameters in a class and a table of the optional parameters. The returned string is typially appended to the doc string of a solver class so that the user can easily see which parameters that must and can be provided.

odespy.solvers.typeset_toc(toc)[source]