RK4adapt function

A function which uses the Fourth order Runge-Kutta method with adaptive step size to solve a system of ODE's.

This function simulates a discrete time Markov chain with transition matrix P, state space 0,1,..,n and and initial state i for nsteps transitions.

RK4adapt(dydt, t0, y0, t1, h0 = 1, tol = 1e-10, ...)

Arguments

• dydt: a function giving the gradient of y(t).
• t0: initial value of t.
• y0: initial value of y(t).
• t1: system solved up to time t1.
• h0: initial step size
• tol: tolerance for adapting step size.
• ...: pass arguments to function dydt.

Details

We assume that P is well defined transition matrix with rows summing to 1.

Returns

Returns a list with elements t, a vector giving times, and y, a matrix whose rows give the solution at successive times.

References

Jones, O.D., R. Maillardet, and A.P. Robinson. 2009. An Introduction to Scientific Programming and Simulation, Using R. Chapman And Hall/CRC.

Examples

LV <- function(t=NULL, y, a, b, g, e, K=Inf)
  c(a*y[1]*(1 - y[1]/K) - b*y[1]*y[2], g*b*y[1]*y[2] - e*y[2])

xy <- RK4adapt(LV, 0, c(100, 50), 200, 1, tol=1e-3, 
               a=0.05, K=Inf, b=0.0002, g=0.8, e=0.03)

par(mfrow = c(2,1))
plot(xy$y[,1], xy$y[,2], type='p', 
     xlab='prey', ylab='pred', main='RK4, adaptive h')
plot(xy$t, xy$y[,1], type='p', xlab='time', 
     ylab='prey circles pred triangles', main='RK4, adaptive h')
points(xy$t, xy$y[,2], pch=2)
par(mfrow=c(1,1))