Population models, chemical kinetics, epidemics and mechanics are all ODEs you integrate numerically.
An ordinary differential equation gives the rate of change of a state: dy/dt = f(t, y). integrate.solve_ivp(fun, t_span, y0) integrates it from t_span[0] to t_span[1] starting at y0.
fun(t, y)takes time first and the state vector, and returns the derivative vector of the same length. Higher-order equations are rewritten as first-order systems: forx'' = -xuse the state[x, v]with derivative[v, -x].t_eval=np.linspace(...)chooses output times; otherwise you get the solver's own steps.dense_output=Truegivessol.sol(t)for any t.args=(...)passes parameters tofun.
The result has t, y (shape (n_states, n_times), one row per variable), success and message.
Methods: the default "RK45" suits most non-stiff problems. If it crawls with tiny steps, the problem is likely stiff (very different time scales, as in chemical kinetics); switch to "LSODA" (switches automatically), "BDF" or "Radau".
Accuracy is controlled by rtol (default 1e-3) and atol (1e-6). The defaults are loose; tighten them for long simulations.
Events are functions of (t, y) whose zero crossings are recorded in sol.t_events; set event.terminal = True to stop there, for example when a ball hits the ground.
Example
import numpy as np
from scipy.integrate import solve_ivp
def falling(t, state, drag): # state = [height, velocity]
h, v = state
return [v, -9.81 - drag * v]
def hit_ground(t, state, drag):
return state[0] # event when height crosses zero
hit_ground.terminal = True
hit_ground.direction = -1
sol = solve_ivp(falling, (0, 10), [100.0, 0.0], args=(0.1,),
events=hit_ground, rtol=1e-8, atol=1e-8)
print(sol.success, sol.y.shape[0]) # True 2
print(round(sol.t_events[0][0], 3)) # 4.882 (seconds until impact)
print(round(sol.y_events[0][0][1], 2)) # -37.89 (impact velocity, m/s)Remember
fun(t, y)returns dy/dt; rewrite higher-order equations as systems- Switch to LSODA, BDF or Radau for stiff problems
- Tighten
rtol/atoland use events to stop at conditions
Watch out for
- Writing the derivative function as
fun(y, t), the oldodeintorder.