solve_ivp from scipy does not integrate the whole range of tspan

Question

I'm trying to use solve_ivp from scipy in Python to solve an IVP. I specified the tspan argument of solve_ivp to be (0,10), as shown below. However, for some reason, the solutions I get always stop around t=2.5.

from scipy.integrate import solve_ivp
import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize as optim

def dudt(t, u):
    return u*(1-u/12)-4*np.heaviside(-(t-5), 1)

ic = [2,4,6,8,10,12,14,16,18,20]

sol = solve_ivp(dudt, (0, 10), ic, t_eval=np.linspace(0, 10, 10000))

for solution in sol.y:
    y = [y for y in solution if y >= 0]
    t = sol.t[:len(y)]
    plt.plot(t, y)

Answer 1

You can examine the sol object to see why the integration failed. It provides the message 'Required step size is less than spacing between numbers.' This usually indicates an implementation error in the right-hand side function or a singularity in the ODE.

Your ODE is simple enough to find the exact solution. We can consider the scalar case because each component is completely independent. For the initial condition $u(0)=2$, the exact solution is $$ u(t) = 6-2 \sqrt{3} \tan \left(\frac{t}{2 \sqrt{3}}+\tan ^{-1}\left(\frac{2}{\sqrt{3}}\right)\right). $$ The tangent function has a singularity at $\frac{\pi}{2}$ which occurs when $t \approx 2.47$. solve_ivp cannot handle singularities (nor should it be expected to).

You will have to decide if your ODE is correct and whether it should have this singularity. Even without it, some care will be needed to handle the discontinuity at $t=5$.