Module for wrapping rkc.f.
Bases: odespy.solvers.Adaptive
Wrapper for rkc.f, a well-known Fortran ODE solver.
Besides the standard attributes of class Solver, class RKC also stores a dictionary statistics, which contains data and explanations from the execution of the RKC subroutine.
The Fortran source code can be obtained from netlib and contains more details. For convenience we quote here from rkc.f the main description of the method:
“ABSTRACT: RKC integrates initial value problems for systems of first order ordinary differential equations. It is based on a family of explicit Runge-Kutta-Chebyshev formulas of order two. The stability of members of the family increases quadratically in the number of stages m. An estimate of the spectral radius is used at each step to select the smallest m resulting in a stable integration. RKC is appropriate for the solution to modest accuracy of mildly stiff problems with eigenvalues of Jacobians that are close to the negative real axis. For such problems it has the advantages of explicit one-step methods and very low storage. If it should turn out that RKC is using m far beyond 100, the problem is not mildly stiff and alternative methods should be considered. Answers can be obtained cheaply anywhere in the interval of integration by means of a continuous extension evaluated in the subroutine RKCINT.
The initial value problems arising from semi-discretization of diffusion-dominated parabolic partial differential equations and of reaction-diffusion equations, especially in two and three spatial variables, exemplify the problems for which RKC was designed.” (rkc.f)
This wrapper does not call RKCINT but runs RKC between each time interval specified by the time_points array sent to the solve method.
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) |
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
|
spcrad | Python function of (u, t) returning the spectral radius of the Jacobian in the rkc.f solver. |
spcrad_f77 | Fortran version of spcrad function. This subroutine should be defined in form:
|
jac_constant | Flag to show whether Jacobian is constant, 0 (false) or 1 (true) (default: 0) |
Methods
adjust_parameters() | |
advance() | Advance solution one time step. |
check_atol() | ATOL need to be supplied as scalar or vector of length NEQ. |
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() | Import extension module _rkc and check that it exists. |
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() |