Ask Your Question
0

integration - StopIntegration support?

asked 2018-11-17 20:26:25 +0100

PatB gravatar image

updated 2018-11-18 19:18:43 +0100

slelievre gravatar image

Is there a way to specify a stopping rule or condition like r[t] < rMin when integrating?

Here's an example in Mathematica:

s = NDSolve[{r''[t] == -GM / r[t]^2, r[0] = 1000 * s0, r'[0] = 0,
             WhenEvent[r[t] < rMin, "**StopIntegration**"]}, r, {t, 0, 86400 * 100}]
edit retag flag offensive close merge delete

1 Answer

Sort by ยป oldest newest most voted
2

answered 2018-11-19 13:54:18 +0100

eric_g gravatar image

updated 2018-11-19 14:03:03 +0100

You can do it by setting the attribute terminal to True in the event function of solve_ivp:

from scipy.integrate import solve_ivp
def rhs(t, (r, v)):     # r'' = -1/r^2 rewritten as the system
    return (v, -1/r^2)  # r' = v, v' = -1/r^2
r_min = 200.
def small_radius(t, (r, v)): 
    return 0 if r < r_min else 1
small_radius.terminal = True
t0, r0, v0 = 0, 1000., 0  # initial conditions
t_max = 50000
sol = solve_ivp(rhs, (t0, t_max), (r0, v0), events=small_radius, rtol=1e-8)
sol.message

gives

'A termination event occurred.'

You can plot the obtained solution by

line(zip(sol.t, sol.y[0]), axes_labels=[r'$t$', r'$r(t)$'])

and check that it indeed stops at r=r_min: image description

For more details, see the documentation of solve_ivp.

edit flag offensive delete link more

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

Stats

Asked: 2018-11-17 20:26:25 +0100

Seen: 970 times

Last updated: Nov 19 '18