Solver for Ordinary Differential Equations (ODE) for COMPLEX variables
Solves the initial value problem for stiff or nonstiff systems of ordinary differential equations (ODE) in the form:
dy/dt = f(t,y)
where dy and y are complex variables.
The R function zvode provides an interface to the FORTRAN ODE
solver of the same name, written by Peter N. Brown, Alan C. Hindmarsh
and George D. Byrne.
zvode(y, times, func, parms, rtol = 1e-6, atol = 1e-6, jacfunc = NULL, jactype = "fullint", mf = NULL, verbose = FALSE, tcrit = NULL, hmin = 0, hmax = NULL, hini = 0, ynames = TRUE, maxord = NULL, bandup = NULL, banddown = NULL, maxsteps = 5000, dllname = NULL, initfunc = dllname, initpar = parms, rpar = NULL, ipar = NULL, nout = 0, outnames = NULL, forcings = NULL, initforc = NULL, fcontrol = NULL, ...)
y |
the initial (state) values for the ODE system. If |
times |
time sequence for which output is wanted; the first
value of |
func |
either an R-function that computes the values of the
derivatives in the ODE system (the model definition) at time
If The return value of If |
parms |
vector or list of parameters used in |
rtol |
relative error tolerance, either a scalar or an array as
long as |
atol |
absolute error tolerance, either a scalar or an array as
long as |
jacfunc |
if not In some circumstances, supplying
If the Jacobian is a full matrix, If the Jacobian is banded, |
jactype |
the structure of the Jacobian, one of
|
mf |
the "method flag" passed to function |
verbose |
if TRUE: full output to the screen, e.g. will
print the |
tcrit |
if not |
hmin |
an optional minimum value of the integration stepsize. In special situations this parameter may speed up computations with the cost of precision. Don't use hmin if you don't know why! |
hmax |
an optional maximum value of the integration stepsize. If
not specified, hmax is set to the largest difference in
|
hini |
initial step size to be attempted; if 0, the initial step size is determined by the solver. |
ynames |
logical; if |
maxord |
the maximum order to be allowed. |
bandup |
number of non-zero bands above the diagonal, in case the Jacobian is banded. |
banddown |
number of non-zero bands below the diagonal, in case the Jacobian is banded. |
maxsteps |
maximal number of steps per output interval taken by the solver. |
dllname |
a string giving the name of the shared library
(without extension) that contains all the compiled function or
subroutine definitions refered to in |
initfunc |
if not |
initpar |
only when ‘dllname’ is specified and an
initialisation function |
rpar |
only when ‘dllname’ is specified: a vector with
double precision values passed to the DLL-functions whose names are
specified by |
ipar |
only when ‘dllname’ is specified: a vector with
integer values passed to the dll-functions whose names are specified
by |
nout |
only used if |
outnames |
only used if ‘dllname’ is specified and
|
forcings |
only used if ‘dllname’ is specified: a list with
the forcing function data sets, each present as a two-columned matrix,
with (time,value); interpolation outside the interval
[min( See forcings or package vignette |
initforc |
if not |
fcontrol |
A list of control parameters for the forcing functions.
forcings or package vignette |
... |
additional arguments passed to |
see vode, the double precision version, for details.
A matrix of class deSolve with up to as many rows as elements
in times and as many columns as elements in y plus the
number of "global" values returned in the next elements of the return
from func,
plus and additional column for the time value. There will be a row
for each element in times unless the FORTRAN routine ‘zvode’
returns with an unrecoverable error. If y has a names
attribute, it will be used to label the columns of the output value.
From version 1.10.4, the default of atol was changed from 1e-8 to 1e-6, to be consistent with the other solvers.
The following text is adapted from the zvode.f source code:
When using zvode for a stiff system, it should only be used for
the case in which the function f is analytic, that is, when each f(i)
is an analytic function of each y(j). Analyticity means that the
partial derivative df(i)/dy(j) is a unique complex number, and this
fact is critical in the way zvode solves the dense or banded linear
systems that arise in the stiff case. For a complex stiff ODE system
in which f is not analytic, zvode is likely to have convergence
failures, and for this problem one should instead use ode on the
equivalent real system (in the real and imaginary parts of y).
Karline Soetaert <karline.soetaert@nioz.nl>
P. N. Brown, G. D. Byrne, and A. C. Hindmarsh, 1989. VODE: A Variable
Coefficient ODE Solver, SIAM J. Sci. Stat. Comput., 10, pp. 1038-1051.
Also, LLNL Report UCRL-98412, June 1988.
G. D. Byrne and A. C. Hindmarsh, 1975. A Polyalgorithm for the Numerical Solution of Ordinary Differential Equations. ACM Trans. Math. Software, 1, pp. 71-96.
A. C. Hindmarsh and G. D. Byrne, 1977. EPISODE: An Effective Package for the Integration of Systems of Ordinary Differential Equations. LLNL Report UCID-30112, Rev. 1.
G. D. Byrne and A. C. Hindmarsh, 1976. EPISODEB: An Experimental Package for the Integration of Systems of Ordinary Differential Equations with Banded Jacobians. LLNL Report UCID-30132, April 1976.
A. C. Hindmarsh, 1983. ODEPACK, a Systematized Collection of ODE Solvers. in Scientific Computing, R. S. Stepleman et al., eds., North-Holland, Amsterdam, pp. 55-64.
K. R. Jackson and R. Sacks-Davis, 1980. An Alternative Implementation of Variable Step-Size Multistep Formulas for Stiff ODEs. ACM Trans. Math. Software, 6, pp. 295-318.
Netlib: http://www.netlib.org
vode for the double precision version
## =======================================================================
## Example 1 - very simple example
## df/dt = 1i*f, where 1i is the imaginary unit
## The initial value is f(0) = 1 = 1+0i
## =======================================================================
ZODE <- function(Time, f, Pars) {
df <- 1i*f
return(list(df))
}
pars <- NULL
yini <- c(f = 1+0i)
times <- seq(0, 2*pi, length = 100)
out <- zvode(func = ZODE, y = yini, parms = pars, times = times,
atol = 1e-10, rtol = 1e-10)
# The analytical solution to this ODE is the exp-function:
# f(t) = exp(1i*t)
# = cos(t)+1i*sin(t) (due to Euler's equation)
analytical.solution <- exp(1i * times)
## compare numerical and analytical solution
tail(cbind(out[,2], analytical.solution))
## =======================================================================
## Example 2 - example in "zvode.f",
## df/dt = 1i*f (same as above ODE)
## dg/dt = -1i*g*g*f (an additional ODE depending on f)
##
## Initial values are
## g(0) = 1/2.1 and
## z(0) = 1
## =======================================================================
ZODE2<-function(Time,State,Pars) {
with(as.list(State), {
df <- 1i * f
dg <- -1i * g*g * f
return(list(c(df, dg)))
})
}
yini <- c(f = 1 + 0i, g = 1/2.1 + 0i)
times <- seq(0, 2*pi, length = 100)
out <- zvode(func = ZODE2, y = yini, parms = NULL, times = times,
atol = 1e-10, rtol = 1e-10)
## The analytical solution is
## f(t) = exp(1i*t) (same as above)
## g(t) = 1/(f(t) + 1.1)
analytical <- cbind(f = exp(1i * times), g = 1/(exp(1i * times) + 1.1))
## compare numerical solution and the two analytical ones:
tail(cbind(out[,2], analytical[,1]))Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.