zvode {deSolve}  R Documentation 
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 = 1e6, atol = 1e6, 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 Rfunction 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 nonzero bands above the diagonal, in case the Jacobian is banded. 
banddown 
number of nonzero 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 DLLfunctions whose names are
specified by 
ipar 
only when ‘dllname’ is specified: a vector with
integer values passed to the dllfunctions 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 twocolumned 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 1e8 to 1e6, 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. 10381051.
Also, LLNL Report UCRL98412, 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. 7196.
A. C. Hindmarsh and G. D. Byrne, 1977. EPISODE: An Effective Package for the Integration of Systems of Ordinary Differential Equations. LLNL Report UCID30112, 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 UCID30132, April 1976.
A. C. Hindmarsh, 1983. ODEPACK, a Systematized Collection of ODE Solvers. in Scientific Computing, R. S. Stepleman et al., eds., NorthHolland, Amsterdam, pp. 5564.
K. R. Jackson and R. SacksDavis, 1980. An Alternative Implementation of Variable StepSize Multistep Formulas for Stiff ODEs. ACM Trans. Math. Software, 6, pp. 295318.
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 = 1e10, rtol = 1e10) # The analytical solution to this ODE is the expfunction: # 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 = 1e10, rtol = 1e10) ## 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]))