odespy.radau5

class odespy.radau5.Radau5(f, **kwargs)[source]

Bases: odespy.solvers.Solver

This is a wrapper for radau5.f, a file containing the the radau5 FORTRAN subtroutine for the numerical solution of a stiff ODEs or differential algebraic (DAE) systems, involving first order equations on the form

\[Mu' = f(u,t).\]

The ODE system can be linearly implicit (mass-matrix \(M\) is not identity matrix) or explicit.

The method used is an implicit Runge-Kutta method (Radau IIA) of order 5 with step size control. The file radau5.f can be downloaded from http://www.unige.ch/~hairer/software.html. Details about the implementation can be found in the book “Solving Ordinary Differential Equations II, Stiff and Differential-Algebraic Problems”, by E. Hairer and G. Wanner, 2nd ed., 1996, published by Springer-Verlag (ISBN: 3-540-60452-9).

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)
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_f77_radau5

Fortran subroutine for jac to be provided to Radau5. This subroutine should be defined as (note that rpar,ipar are not used):

      subroutine jac_f77_radau5
     &           (neq,t,u,dfu,ldfu,rpar,ipar)
Cf2py intent(hide) neq,rpar,ipar
Cf2py intent(in)   t,u,ldfu
Cf2py intent(out) dfu
      integer neq,ipar,ldfu
      double precision t,u,dfu,rpar
      dimension u(neq),dfu(ldfu,neq)
      dfu(1,1) = ...
      dfu(1,2) = ...
      ...
      return
      end
atol Absolute tolerance for solution. (default: 1e-08)
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.
rtol Relative tolerance for solution. (default: 1e-06)
nsteps Max no of internal solver steps per time step. (default: 1000)
ml Number of lower non-zero diagonals in a banded Jacobian.
mu Number of upper non-zero diagonals in a banded Jacobian.
first_step Suggested first time step size for an adaptive algorithm.
jac Jacobian of right-hand side function f (df/du). (default: None)
safety Safety factor on new step selection (default: 0.9)
max_step Maximum step size for an adaptive algorithm.
nsteps Max no of internal solver steps per time step. (default: 1000)
mas User-supplied function to calculate the mass matrix M in the linearly implicit ODE system M*u’ = f(u, t). This function is used only when M is not I (the identity). M can be supplied as full matrix or banded matrix. The mas function returns M and takes any number of user-defined parameters as arguments, given through the parameter mas_args. If M is sparse, a banded matrix with dimensions (mlmas+mumas+1,neq) is returned, where mlmas and mumas are the lower and upper half-bandwidth of the banded matrix, and neq is the number of equations.
mas_f77

Fortran subroutine for mas. This subroutine has the signature (note that rpar,ipar are not used):

      subroutine mas_f77(neq,mas,lmas,rpar,ipar)
Cf2py intent(hide)   rpar,ipar
Cf2py intent(in)     neq,lmas
Cf2py intent(out)    mas
      integer neq,lmas,ipar
      double precision mas(lmas,neq),rpar
      mas(1,1) = ...
      ...
      return
      end
mlmas Lower bandwith of mass matrix M
mumas Upper bandwith of mass matrix M
mas_args Extra positional arguments to the mas function (it is called as mas(*mas_args). (default: ())

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)
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_f77_radau5

Fortran subroutine for jac to be provided to Radau5. This subroutine should be defined as (note that rpar,ipar are not used):

      subroutine jac_f77_radau5
     &           (neq,t,u,dfu,ldfu,rpar,ipar)
Cf2py intent(hide) neq,rpar,ipar
Cf2py intent(in)   t,u,ldfu
Cf2py intent(out) dfu
      integer neq,ipar,ldfu
      double precision t,u,dfu,rpar
      dimension u(neq),dfu(ldfu,neq)
      dfu(1,1) = ...
      dfu(1,2) = ...
      ...
      return
      end
atol Absolute tolerance for solution. (default: 1e-08)
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.
rtol Relative tolerance for solution. (default: 1e-06)
nsteps Max no of internal solver steps per time step. (default: 1000)
ml Number of lower non-zero diagonals in a banded Jacobian.
mu Number of upper non-zero diagonals in a banded Jacobian.
first_step Suggested first time step size for an adaptive algorithm.
jac Jacobian of right-hand side function f (df/du). (default: None)
safety Safety factor on new step selection (default: 0.9)
max_step Maximum step size for an adaptive algorithm.
nsteps Max no of internal solver steps per time step. (default: 1000)
mas User-supplied function to calculate the mass matrix M in the linearly implicit ODE system M*u’ = f(u, t). This function is used only when M is not I (the identity). M can be supplied as full matrix or banded matrix. The mas function returns M and takes any number of user-defined parameters as arguments, given through the parameter mas_args. If M is sparse, a banded matrix with dimensions (mlmas+mumas+1,neq) is returned, where mlmas and mumas are the lower and upper half-bandwidth of the banded matrix, and neq is the number of equations.
mas_f77

Fortran subroutine for mas. This subroutine has the signature (note that rpar,ipar are not used):

      subroutine mas_f77(neq,mas,lmas,rpar,ipar)
Cf2py intent(hide)   rpar,ipar
Cf2py intent(in)     neq,lmas
Cf2py intent(out)    mas
      integer neq,lmas,ipar
      double precision mas(lmas,neq),rpar
      mas(1,1) = ...
      ...
      return
      end
mlmas Lower bandwith of mass matrix M
mumas Upper bandwith of mass matrix M
mas_args Extra positional arguments to the mas function (it is called as mas(*mas_args). (default: ())

Methods

adjust_parameters()
advance() This function intends to one step forward with Radau5.
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.
func_wrappers()
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() Import extension module _radau5 and check that it exists.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_dummy_functions()
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
set_tol() Tolerance parameters, rtol & atol, can be provided as either two scalars
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.radau5'
adjust_parameters()[source]
advance()[source]

This function intends to one step forward with Radau5.

func_wrappers()[source]
initialize()[source]

Import extension module _radau5 and check that it exists.

initialize_for_solve()[source]
set_dummy_functions()[source]
set_tol()[source]

Tolerance parameters, rtol & atol, can be provided as either two scalars or two vectors.

validate_data()[source]
class odespy.radau5.Radau5Explicit(f, **kwargs)[source]

Bases: odespy.radau5.Radau5

Radau5 solver for explcit ODE problem.

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)
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_f77_radau5

Fortran subroutine for jac to be provided to Radau5. This subroutine should be defined as (note that rpar,ipar are not used):

      subroutine jac_f77_radau5
     &           (neq,t,u,dfu,ldfu,rpar,ipar)
Cf2py intent(hide) neq,rpar,ipar
Cf2py intent(in)   t,u,ldfu
Cf2py intent(out) dfu
      integer neq,ipar,ldfu
      double precision t,u,dfu,rpar
      dimension u(neq),dfu(ldfu,neq)
      dfu(1,1) = ...
      dfu(1,2) = ...
      ...
      return
      end
atol Absolute tolerance for solution. (default: 1e-08)
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.
rtol Relative tolerance for solution. (default: 1e-06)
nsteps Max no of internal solver steps per time step. (default: 1000)
ml Number of lower non-zero diagonals in a banded Jacobian.
mu Number of upper non-zero diagonals in a banded Jacobian.
first_step Suggested first time step size for an adaptive algorithm.
jac Jacobian of right-hand side function f (df/du). (default: None)
safety Safety factor on new step selection (default: 0.9)
max_step Maximum step size for an adaptive algorithm.
nsteps Max no of internal solver steps per time step. (default: 1000)

Methods

adjust_parameters()
advance() This function intends to one step forward with Radau5.
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.
func_wrappers()
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() Import extension module _radau5 and check that it exists.
initialize_for_solve()
set([strict]) Assign values to one or more parameters, specified as keyword arguments.
set_dummy_functions()
set_initial_condition(U0) Function set_initial_condition() is used to set initial value of
set_tol() Tolerance parameters, rtol & atol, can be provided as either two scalars
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.radau5'
odespy.radau5.Radau5Implicit

alias of Radau5