For anyone interested I"ve posted this in Scientific Computing Stack Exchange. I mention there that for new projects, one should use the more modern solve_ivp, but stil I'm curious to understand the interaction between requested output step size and the internal integration step size. I thought they were decoupled but there seems to be only 9 nonzero values of dh when 11 output steps are requested, but almost 300 nonzero dh values for a million. 9 seems too low. https://scicomp.stackexchange.com/questions/42613/understanding-scipys-odein...