如何为SciPy odeint求解器传入数组输入信号而非定义函数

方法1：用插值适配odeint #

由于odeint求解时会自适应使用用户指定时间数组之外的中间步长，无法直接传入离散输入数组，我们可以通过线性插值将输入数组转换为能响应任意时间点的函数，适配odeint的调用逻辑。

示例代码 #

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

# 先生成时间数组和处理后的输入信号数组（这里演示加窗处理）
t = np.linspace(0, 10, 101)
input_array = np.sin(2 * np.pi * 50 * t)
input_array *= np.hamming(len(input_array))  # 加汉明窗

# 定义插值函数，根据任意时间值返回对应输入
def get_input(t_val, t_grid, input_grid):
    return np.interp(t_val, t_grid, input_grid)

def pend(y, t, b, c, t_grid, input_grid):
    theta, omega = y
    u = get_input(t, t_grid, input_grid)
    dydt = [omega, -b*omega - c*np.sin(theta) - u]
    return dydt

b = 0.25
c = 5.0
y0 = [0, 0.0]
# 将时间网格和输入数组通过args传入odeint
sol = odeint(pend, y0, t, args=(b, c, t, input_array))

plt.plot(t, sol, label=['theta(t)','omega(t)'])
plt.legend(loc='best')
plt.xlabel('t')
plt.grid()
plt.show()

方法2：使用solve_ivp（更灵活的新一代求解器） #

scipy.integrate.solve_ivp是SciPy官方推荐的现代ODE求解器，接口设计更灵活，同样可以通过插值适配离散输入数组，且参数顺序更符合常规习惯（t在前，y在后）。

示例代码 #

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

t = np.linspace(0, 10, 101)
input_array = np.sin(2 * np.pi * 50 * t)
input_array *= np.hamming(len(input_array))  # 加窗处理

def pend(t, y, b, c, t_grid, input_grid):
    theta, omega = y
    u = np.interp(t, t_grid, input_grid)
    dydt = [omega, -b*omega - c*np.sin(theta) - u]
    return dydt

b = 0.25
c = 5.0
y0 = [0, 0.0]
# 调用solve_ivp，指定求解时间区间和输出时间点
sol = solve_ivp(pend, [t[0], t[-1]], y0, args=(b, c, t, input_array), t_eval=t)

plt.plot(sol.t, sol.y.T, label=['theta(t)','omega(t)'])
plt.legend(loc='best')
plt.xlabel('t')
plt.grid()
plt.show()

关键说明 #

• 直接传数组不可行的原因：ODE求解器会自适应调整步长，并非仅使用用户提供的离散时间点，必须通过插值获取任意时间点的输入值。
• 如果输入是脉冲类离散信号，可以替换插值逻辑为np.where或阶跃判断，但插值是通用适配方案。

内容的提问来源于stack exchange，提问作者Nuopel