使用odeint求解ODE后调用quad积分时出现IndexError的排查与解决

我用odeint求解了v(t)和u(t)的耦合微分方程，代码如下：

f = odeint(ODEs, f0, t)
    
v = f[:,0]    
u = f[:,1]    

想要计算v(t)在[0,1]区间的积分，于是用scipy.integrate.quad写了这段代码：

x = lambda t: f[t,0]

delta_x=quad(x, 0,1)
print(delta_x)

运行后触发索引错误：

IndexError: only integers, slices (`:`), ellipsis (`...`), numpy.newaxis (`None`) and integer or boolean arrays are valid indices

Jupyter Notebook高亮了f[t,0]这一行。以下是最小可复现代码：

import numpy as np
from scipy.integrate import odeint, quad 

def ODEs(f,t):
    v = f[0]
    u = f[1]
    
    dvdt = u-v
    dudt = (v-u)/3
    
    return [dvdt, dudt]

f0 = [0.3, 0.2]

t = np.linspace(0,4,10000)
f = odeint(ODEs, f0, t)
v = f[:,0]    
u = f[:,1]  

x = lambda t: f[t,0]

delta_x=quad(x, 0,1)
print(delta_x)

错误原因 #

quad在数值积分过程中，会向传入的函数传递浮点数类型的采样点（比如0.1、0.5这类非整数），而f是numpy数组，只能接受整数、切片等合法索引类型。当lambda函数里用浮点数t去索引f[t,0]时，就会触发索引错误。

解决方法 #

方法一：用插值函数包装v(t) #

通过插值将离散的v(t)转换为连续可调用的函数，让quad可以正常调用：

from scipy.interpolate import interp1d

# 创建v(t)的线性插值函数，可选'cubic'等更高阶插值
v_continuous = interp1d(t, v, kind='linear')
delta_x = quad(v_continuous, 0, 1)
print(delta_x)

方法二：直接对离散点用数值积分 #

既然已经有了t和v的离散采样点，直接用梯形法或辛普森法积分更高效：

# 筛选出t在[0,1]区间内的子集
mask = (t >= 0) & (t <= 1)
t_segment = t[mask]
v_segment = v[mask]

# 梯形法积分
delta_x_trapz = np.trapz(v_segment, t_segment)
# 辛普森法积分（需scipy版本≥1.10）
from scipy.integrate import simpson
delta_x_simpson = simpson(v_segment, t_segment)

print("梯形法积分结果:", delta_x_trapz)
print("辛普森法积分结果:", delta_x_simpson)

方法三：通过查找最近索引匹配（不推荐） #

如果硬要保留lambda的写法，可以通过查找与传入浮点数最接近的t的索引，但这种方法精度依赖采样密度，不如插值法：

x = lambda t_val: v[np.argmin(np.abs(t - t_val))]
delta_x = quad(x, 0, 1)
print(delta_x)

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