当输出在scipy.integrate.odeint中达到0时停止集成

我编写了一个代码,用于查看带有拖动的物体的抛射物运动 . 我正在使用scipy的odeint来执行前向Euler方法 . 集成将一直运行,直到达到时间限制 . 我想在达到这个极限时或者当ry = 0的输出时(即射弹已落地)停止积分 .

def f(y, t):
    # takes vector y(t) and time t and returns function y_dot = F(t,y) as a vector
    # y(t) = [vx(t), vy(t), rx(t), ry(t)]

    # magnitude of velocity calculated
    v = np.sqrt(y[0] ** 2 + y[1] **2)

    # define new drag constant k
    k = cd * rho * v * A / (2 * m)

    return [-k * y[0], -k * y[1] - g, y[0], y[1]] 


def solve_f(v0, ang, h0, tstep, tmax=1000):
    # uses the forward Euler time integration method to solve the ODE using f function
    # create vector for time steps based on args
    t = np.linspace(0, tmax, tmax / tstep)

    # convert angle from degrees to radians
    ang = ang * np.pi / 180

    return odeint(f, [v0 * np.cos(ang), v0 * np.sin(ang), 0, h0], t)

solution = solve_f(v0, ang, h0, tstep)

我已经尝试了几个循环和类似的尝试在ry = 0时停止集成 . 并在下面找到了这个问题,但是无法用odeint实现类似的东西 . 解决方案[:,3]是ry的输出列 . 使用odeint有一种简单的方法吗?

https://stackoverflow.com/questions/26738676/does-scipy-integrate-ode-set-solout-work

回答(1)

2 years ago

结帐 scipy.integrate.ode here . 它比 odeint 更灵活,有助于您想要做的事情 .

一个使用垂直镜头的简单示例,集成到接地:

from scipy.integrate import ode, odeint
import scipy.constants as SPC

def f(t, y):
    return [y[1], -SPC.g]

v0 = 10
y0 = 0

r = ode(f)
r.set_initial_value([y0, v0], 0)

dt = 0.1
while r.successful() and r.y[0] >= 0:
    print('time: {:.3f}, y: {:.3f}, vy: {:.3f}'.format(r.t + dt, *r.integrate(r.t + dt)))

每次调用 r.integrate 时, r 将存储当前时间和y值 . 如果要存储它们,可以将它们传递给列表 .