odespy.problems

class odespy.problems.ComplexOscillator(w=1, U0=[1, 0])[source]

Bases: odespy.problems.Problem

Harmonic oscillator (u’’ + w*u = 0) expressed as a complex ODE: u’ = i*w*u.

Methods

default_parameters([resolution_per_period])
f(u, t)
get_initial_condition() Return vector of initial conditions.
jac(u, t)
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(w=1, U0=[1, 0])[source]
__module__ = 'odespy.problems'
default_parameters(resolution_per_period=20)[source]
f(u, t)[source]
jac(u, t)[source]
short_description = "Complex oscillator: u' = i*w*u"
u_exact(t)[source]
class odespy.problems.Exponential(a=1, b=0, A=1)[source]

Bases: odespy.problems.Problem

Methods

default_parameters()
f(u, t)
f_with_args(u, t, a, b)
f_with_kwargs(u, t[, a, b])
get_initial_condition() Return vector of initial conditions.
jac(u, t)
jac_with_args(u, t, a, b)
jac_with_kwargs(u, t[, a, b])
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(a=1, b=0, A=1)[source]
__module__ = 'odespy.problems'
default_parameters()[source]
f(u, t)[source]
f_with_args(u, t, a, b)[source]
f_with_kwargs(u, t, a=1, b=0)[source]
jac(u, t)[source]
jac_with_args(u, t, a, b)[source]
jac_with_kwargs(u, t, a=1, b=0)[source]
short_description = "Exponential solution: u' = a*u + b, u(0)=A"
u_exact(t)[source]
class odespy.problems.Gaussian0(c=2.0, s=1.0)[source]

Bases: odespy.problems.Problem

Methods

default_parameters() Compute suitable time_points, atol/rtol, etc.
f(u, t)
get_initial_condition() Return vector of initial conditions.
jac(u, t)
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(c=2.0, s=1.0)[source]
__module__ = 'odespy.problems'
f(u, t)[source]
jac(u, t)[source]
short_description = 'Gaussian function as solution'
u_exact(t)[source]
class odespy.problems.Gaussian1(c=2.0, s=1.0)[source]

Bases: odespy.problems.Problem

Methods

default_parameters() Compute suitable time_points, atol/rtol, etc.
f(u, t)
get_initial_condition() Return vector of initial conditions.
jac(u, t)
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(c=2.0, s=1.0)[source]
__module__ = 'odespy.problems'
f(u, t)[source]
jac(u, t)[source]
short_description = '1 + Gaussian function as solution'
u_exact(t)[source]
class odespy.problems.Linear1(a=1, b=0, c=2)[source]

Bases: odespy.problems.Problem

Linear solution of a nonlinear ODE, which is normally exactly reproduced by a numerical method.

Methods

default_parameters()
f(u, t)
get_initial_condition() Return vector of initial conditions.
jac(u, t)
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(a=1, b=0, c=2)[source]
__module__ = 'odespy.problems'
default_parameters()[source]
f(u, t)[source]
jac(u, t)[source]
short_description = "u' = a + (u - a*t - b)**c, u(0)=b"
u_exact(t)[source]
class odespy.problems.Linear2(a=1, b=0, c=2)[source]

Bases: odespy.problems.Linear1

Linear solution of nonlinear 2x2 system, which is normally exactly reproduced by a numerical method.

Methods

default_parameters()
f(u, t)
get_initial_condition() Return vector of initial conditions.
jac(u, t)
spcrad(u, t)
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(a=1, b=0, c=2)[source]
__module__ = 'odespy.problems'
f(u, t)[source]
jac(u, t)[source]
short_description = '2x2 nonlinear system with linear solution'
spcrad(u, t)[source]
u_exact(t)[source]
class odespy.problems.Linear2t(a=1, b=0, c=2)[source]

Bases: odespy.problems.Linear2

Linear solution of trivial 2x2 system, which is normally exactly reproduced by a numerical method.

Methods

default_parameters()
f(u, t)
get_initial_condition() Return vector of initial conditions.
jac(u, t)
spcrad(u, t)
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__module__ = 'odespy.problems'
f(u, t)[source]
jac(u, t)[source]
short_description = '2x2 trivial system with linear solution'
spcrad(u, t)[source]
class odespy.problems.LinearOscillator(m=1, k=1, x0=1, v0=0, F=None)[source]

Bases: odespy.problems.Problem

Methods

default_parameters([resolution_per_period])
f(u, t)
get_initial_condition() Return vector of initial conditions.
jac(u, t)
kinetic_energy(u)
potential_energy(u)
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(m=1, k=1, x0=1, v0=0, F=None)[source]
__module__ = 'odespy.problems'
default_parameters(resolution_per_period=20)[source]
f(u, t)[source]
jac(u, t)[source]
kinetic_energy(u)[source]
potential_energy(u)[source]
short_description = "Linear oscillator: m*u'' + k*u = F(t)"
u_exact(t)[source]
class odespy.problems.Logistic(a, R, A)[source]

Bases: odespy.problems.Problem

Methods

default_parameters() Compute suitable time_points, atol/rtol, etc.
f(u, t)
get_initial_condition() Return vector of initial conditions.
jac(u, t)
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(a, R, A)[source]
__module__ = 'odespy.problems'
f(u, t)[source]
jac(u, t)[source]
short_description = 'Logistic equation'
u_exact(t)[source]
class odespy.problems.Problem[source]

The user derives a subclass and implements the right-hand side function:

def f(self, u, t):
    """Python function for the right-hand side of the ODE u'=f(u,t)."""

Optionally, the Jacobian can be computed:

def jac(self, u, t):
    """Python function for the Jacobian of the right-hand side: df/du."""

Some solvers also allow constraings (DAE problems):

def constraints(self, u, t):
    """Python function for additional constraints: g(u,t)=0."""

Some problem classes will also define command-line arguments:

def define_command_line_arguments(self, parser):
    """
    Initialize an argparse object for reading command-line
    option-value pairs. `parser` is an ``argparse`` object.
    """

Other functions, for which there are default implementations in the superclass, are u_exact for returning the exact solution (if known), verify for checking that a numerical solution is correct (within a tolerance of the analytical solution),

See the tutorial for examples on subclasses of Problem.

Methods

default_parameters() Compute suitable time_points, atol/rtol, etc.
get_initial_condition() Return vector of initial conditions.
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t) Implementation of the exact solution.
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__contains__(attr)[source]

Return True if attr is a method in instance self.

__init__()[source]
__module__ = 'odespy.problems'
complex_ = False
default_parameters()[source]

Compute suitable time_points, atol/rtol, etc. for the particular problem. Useful for quick generation of test cases, demos, unit tests, etc.

get_initial_condition()[source]

Return vector of initial conditions. Not necessary to implement in subclass if the initial condition is stored in self.U0 (this method then returns that condition).

not_suitable_solvers = []
short_description = ''
stiff = False
terminate(u, t, step_number)[source]

Default terminate function, always returning False.

u_exact(t)[source]

Implementation of the exact solution. Return None if no exact solution is available.

verify(u, t, atol=None, rtol=None)[source]

Return True if u at time points t coincides with an exact solution within the prescribed tolerances. If one of the tolerances is None, return max computed error (infinity norm). Return None if the solution cannot be verified.

class odespy.problems.StiffSystem2x2(eps=0.01)[source]

Bases: odespy.problems.Problem

Methods

default_parameters() Compute suitable time_points, atol/rtol, etc.
f(u, t)
get_initial_condition() Return vector of initial conditions.
jac(u, t)
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(eps=0.01)[source]
__module__ = 'odespy.problems'
f(u, t)[source]
jac(u, t)[source]
short_description = "Potentially stiff 2x2 system u' = 1/eps*u"
stiff = True
u_exact(t)[source]
class odespy.problems.VanDerPolOscillator(U0=[2, 1], mu=3.0, f77=False)[source]

Bases: odespy.problems.Problem

Classical van der Pool oscillator:

\[y'' = \mu (1 - y^2) y' - y\]

with initial conditions \(y(0)=2, y'(0)=1\). The equation is rewritten as a system

\[\begin{split}u_0' &= u_1, \ u_1' &= \mu (1-u_0^2)u_1 - u_0\end{split}\]

with a Jacobian

\[\begin{split}\left(egin{array}{cc} 0 & 1\ -2\mu u_0 - 1 & \mu (1-u_0^2) \end{array}\end{split}\]

ight)

Methods

default_parameters([resolution_per_period])
f(u, t)
f_with_args(u, t, mu)
f_with_kwargs(u, t[, mu])
get_initial_condition() Return vector of initial conditions.
jac(u, t)
jac_with_args(u, t, mu)
jac_with_kwargs(u, t[, mu])
str_f_f77() Return f(u,t) as Fortran source code string.
str_jac_f77_fadau5() Return Jacobian Fortran source code string for Radau5.
str_jac_f77_lsode_dense() Return Fortran source for dense Jacobian matrix in LSODE format.
terminate(u, t, step_number) Default terminate function, always returning False.
u_exact(t)
verify(u, t[, atol, rtol]) Return True if u at time points t coincides with an exact
__init__(U0=[2, 1], mu=3.0, f77=False)[source]
__module__ = 'odespy.problems'
default_parameters(resolution_per_period=20)[source]
f(u, t)[source]
f_with_args(u, t, mu)[source]
f_with_kwargs(u, t, mu=3)[source]
jac(u, t)[source]
jac_with_args(u, t, mu)[source]
jac_with_kwargs(u, t, mu=3)[source]
short_description = 'Van der Pol oscillator'
str_f_f77()[source]

Return f(u,t) as Fortran source code string.

str_jac_f77_fadau5()[source]

Return Jacobian Fortran source code string for Radau5.

str_jac_f77_lsode_dense()[source]

Return Fortran source for dense Jacobian matrix in LSODE format.

u_exact(t)[source]
odespy.problems.convergence(u, t, u_ref=None, t_ref=None)[source]

Given a series of solutions and corresponding t arrays in u and t, use the analytical solution or a reference solution in u_ref and t_ref to compute errors. Compute pairwise convergence rates for errors on consecutive time meshes.

odespy.problems.default_oscillator(P, resolution_per_period=20)[source]
odespy.problems.tester(problems, methods, time_points=None, compare_tol=0.0001, solver_prm={})[source]

problems is a list of Problem subclass instances made ready. methods is a list of strings holding the method names.