1

I am running a numerical integration using scipy ode on several million cases and occasionally I get the error:

Excess work done on this call (perhaps wrong Dfun type).

I increased nsteps a fair bit (500,000), and I also tried specifying the BDF method for 'stiff' problems, as described in this question

r = ode(RHS).set_initial_value(state0, t_initial).set_f_params(Efield,qmp)
r.set_integrator('vode',nsteps=500000,method='bdf')

My question is: what happens when this occurs? Is that integration run thrown out? Or, can I force it to be thrown out? I don't care if I lose a few out of millions of runs, but I don't want the data to be contaminated.

Alex
  • 427
  • 6
  • 12

1 Answers1

1

The return value is suspect when this error occurs. As a test case, I tried integrating dy/dt = y from 0 to 1 with mere 10 steps:

from scipy.integrate import ode
r = ode(lambda t, y: y).set_initial_value(1, 0).set_integrator('lsoda', nsteps=10)
print(r.integrate(1))

This prints [1.40213975] with the warning "Excess work done on this call (perhaps wrong Dfun type)". The return value is clearly wrong: the exact solution is exp(1) = 2.718... It seems the process of integration indeed stops when nsteps is reached, and whatever y is found by that time is returned, even if it is not y(1).

You can detect this by calling r.successful() after r.integrate: if it returns False, the most recent step failed.

Increasing nsteps (50 is enough here) eliminates the warning and returns correct value.

Aside: consider using new SciPy API for ODE such as solve_ivp, which eliminate some of the micromanaging of the old ode method (there is no nsteps option in solve_ivp).

  • OK, so to clarify, the integration does not actually stop at nsteps, it keeps going until it's done and just throws a warning at nsteps? – Alex May 15 '18 at 18:48
  • Re-reading the code, it does terminate when reaching nsteps, so my answer was wrong. See the revision. –  May 15 '18 at 19:37
  • OK thanks. So then is there a way to throw that integration out? Or maybe I can catch that warning and throw it out manually? – Alex May 15 '18 at 20:27
  • On another note, I am still using the old ODE because I need to bail out based on position and not time (see my prev question: https://stackoverflow.com/questions/50047265/scipy-ode-integrate-to-unknown-limit-of-t-based-on-the-size-of-solution) – Alex May 15 '18 at 20:28
  • Just check `r.successful()` after each `r.integrate`. Edited in the answer. –  May 15 '18 at 21:54