Hi,
to diagnose your problem, I tried the model with your original
parameters, using a Runge-Kutta fixed step integrator and smaller time
steps. You may try to make them even smaller.
rk4 has no warranted accuracy, so it is less reliable than "lsoda" etc.
However, it can be useful for debugging, because it runs constantly
through, even if NA or NaN's occur. We see that the states grow
drastically to extremely high (or even negative) values.
This points indeed to unrealistic parameter values, or to a
model-misspecification, which means that necessary feedback mechanisms
are missing, so that i12 can grow, even if nothing is left (i1, i2) from
which it can grow from.
Regards, Thomas
init <- c(i1=10, i2=10, i12=0)
parameters <- c(alpha1=0.7, alpha2=0.5, beta1=0.5, beta2=0.3,
gamma1=0.5, gamma2=0.5, delta=0.5, N=100)
times <- seq(0, 10, by=.01)
simul <- as.data.frame(ode(y = init, times = times, func = system,
parms = parameters, method="rk4"))
head(simul, 200)
time i1 i2 i12
1 0.00 1.000000e+01 1.000000e+01 0.000000e+00
2 0.01 1.685572e+01 1.433028e+01 6.451953e-01
3 0.02 2.567008e+01 1.875164e+01 4.259949e+00
4 0.03 3.251911e+01 2.091655e+01 4.215409e+01
5 0.04 -1.270280e+00 -1.988876e+00 1.581789e+02
6 0.05 -8.108876e+02 -4.992572e+02 2.631133e+04
7 0.06 -4.166359e+73 -2.975971e+73 6.427467e+98
8 0.07 NaN NaN NaN
9 0.08 NaN NaN NaN
...
--
Dr. Thomas Petzoldt
Technische Universitaet Dresden
Faculty of Environmental Sciences
Institute of Hydrobiology
01062 Dresden, Germany
E-Mail: thomas.petzoldt at tu-dresden.de
http://tu-dresden.de/Members/thomas.petzoldt