The optimal step size in the Lotka-Volterra equations

Question

Given a Lotka-Volterra predator-prey system,

$x'(t) = ax - bxy$ and $y'(t) = cxy - dy$,

and the initial values (for example, $x(0) = 600$ and $y(0) = 400$), I need to estimate the optimal step size for the Lotka-Volterra predator-prey system when solved using the Python function integrate.odeint() with t in the range 0, 50. Specifically, I want to determine the best number of values of t between 0 and 50 to obtain accurate results. I have noticed that typically about 1000 values of t are included in this case, but I'm not sure why. Is there a formula or a method to calculate the optimal number of values for t in this context?

Answer 1

All the solvers in scipy use variable step size for the integration. They also provide an interpolation formula, so that you can get values at arbitrary times.

So your question is not posed in the right direction. The output step-size does not control the step size of the internal integration steps. (There is a minor influence due to the implementation strategy, it is observable but does not change the accuracy.) The internal step size and thus the accuracy of the solution is controlled via the error threshold parameters atol and rtol. Their default values are good to produce accurate plots over a variety of standard problems like Volterra-Lotka. Try out how the values at the long end of your interval change if you vary these parameters.

The evaluation time points that you pass with the time array need to have enough sampling density that the bows in the plot look nice. Or so that they are sufficient for any evaluation of the returned function value table afterwards, statistical or otherwise.

With the extended return objects of odeint you can get information on the internal steps that were taken, but in the normal use of the function that is not possible. solve_ivp has standard modes that return the internal steps only, with or without interpolation function.