scipy - 使用 scipy.integrate.solve_ivp 时如何处理微分函数中由 t 索引的参数

Question

我想使用solve_ivp来模拟一个演化出以下微分函数的系统：

dy(t) = A D(t) + B y(t)

至于我，t 是从 0 开始计数的整数，我做到了

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp

def fun(t, y, D, A, B):
    dD1 = A*D[int(t)] - B*y
    return dD1#[dD1, dD2]

D = np.zeros((1024))
D[10:34] = 12
tspan = np.linspace(0, np.size(D)-1, np.size(D))

sol = solve_ivp(lambda t, y: fun(t, y, D, A=0.8, B=0.025),
                [tspan[0], tspan[-1]], 
                [0], method='RK45', t_eval=tspan)

plt.figure()
plt.plot(sol.t, sol.y[0,])

它给了我预期的结果，在 t=34 处有一个 y 的峰值。 但是，如果我通过简单地替换来移动 D 具有非零值的间隔

D[10:34] = 12

和

D[420:444] = 12

我希望我可以得到相同的形状，但只能移动，而结果给了我全部 0...此外，有时我会遇到警告消息：

RuntimeWarning：在 double_scalars max(1, SAFETY * error_norm ** (-1 / (order + 1))) 中遇到除以零

Answer 1

估计是找到原因了，单纯是步长造成的，所以传参数

first_step=np.size(D)/100, max_step=np.size(D)/100

solve_ivp 函数解决了这个问题。