Simulación en tiempo real suave con scipy solve_ivp.

Tengo un sistema de EDOs que hasta ahora he resuelto usando solve_ivp.

scipy.integrate.solve_ivp(fun=modelo, t_span=(0.0, t_final), y0=[s0])

Mi problema es que quiero resolver las EDOs en una simulación en ejecución, donde constantemente se inyectan nuevos valores en la simulación y se muestran los resultados. La simulación puede durar varias horas. Mi enfoque fue llamar a solve_ivp repetidamente, más o menos como sigue (esto me permite mostrar resultados intermedios y obtener nuevos datos, que no se muestran aquí):

t = 0.0
s = s0
while t < t_final:
resultado = scipy.integrate.solve_ivp(fun=modelo, t_span=(t, t + t_paso), y0=[s])
s = resultado.y[0][-1]
t += t_paso

Escribí algunos casos de prueba, donde resolví las EDOs analíticamente y por razones que aún no entiendo completamente, llamar repetidamente a solveivp siempre se acerca más a la solución analítica (sin ajustar manualmente los diferentes parámetros de solveivp). Mi pregunta es más en el sentido de ¿hay algo mal en este enfoque, o hay alguna funcionalidad en scipy o en otro paquete que pueda ser más adecuada para mis necesidades?