0

I'm integrating my ode with scipy.integrate.solve_ivp. I have a vector of six elements, three positions and three velocities, which is the **x** of my system. I'd like the integration to stop when the module of the position (being math.sqrt(x[0]**2+x[1]**2+x[2]**2) ) is higher than I certain value. I know this can be done by passing an event function to solve_ivp. However, it's not clear to me how to pass the function to solve_ivp. I'd like also to include the option terminal=True, to stop the integration when the condition is met. Can you please show me how to use it?

I defined the event function as:

def my_fun(t,x):     
    r=math.sqrt(x[0]**2+x[1]**2+x[2]**2)     
        if 10000-r<=0:
             x=0
     return x

I tried the solution suggested here How to pass parameters to event functions in`scipy.integrate.solve_ivp`? But I get an error from scipy,

  File "C:\Users\vsartoretto\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\scipy\integrate\_ivp\base.py", line 249, in __call__
    t = np.asarray(t)
TypeError: only size-1 arrays can be converted to Python scalars

this is how I'm calling it

sol = solve_ivp(fun=lambda t,x:my_eq(t,x,param), t_span=integr_time, y0=np.concatenate((r,v)), method='LSODA',
                    events=my_fun, atol=1e-12, rtol=1e-13)

UPDATE: if I remove the events function everything works nominally

vittosa95
  • 1
  • 1

1 Answers1

0

the issue was in the way I defined the event function. This one works fine: def alt_reached(t,y): return (mh.sqrt(y[0]**2+y[1]**2+y[2]**2) - (6378.1363+531))

vittosa95
  • 1
  • 1
  • Yes, the event function is a scalar valued function that has the event location as root. You could alternatively take `return sum(y**2)-r**2`, as `y` will be a numpy array when called from the solver. – Lutz Lehmann Jun 24 '23 at 17:33