Currently, my problem requires me to solve a system a large system of non-linear ODEs (up to ~5000). So far, I have been using
scipy.integrate.odeint as my workhorse.
A simpler subset of my problem only involves ~50 non-linear ODEs. Solving the system over an interval $[0, 10]$ takes about 3 times longer than solving the system over $[0, 1]$. Not bad!
However, solving the system over $[0, 1000]$ takes about 1500 times longer than solving over $[0, 1]$. Not so nice!
An obvious reason why this might be the case is because at some points between $[0, 1000]$, the solution takes far longer to converge than usual -- so, the average time taken to solve the system at each time step is greater over $[0, 1000]$ then over $[0, 1]$. Is my intuition correct over here?
What other things might be issues?