Python实现RK2（中点法）迭代异常求助：数组y赋值问题排查

嘿，我来帮你捋一捋RK2（中点法）实现里的常见坑！你提到数组y只有前两个值正常，后续全异常，大概率是迭代公式写错了，或者数组更新的逻辑出了问题。

先明确RK2中点法的正确步骤 #

中点法的核心是用区间中点的斜率来更新值，正确的迭代公式是：

1. 计算当前点的斜率：k1 = f(t_n, y_n)
2. 计算中点的斜率：k2 = f(t_n + h/2, y_n + (h*k1)/2)
3. 更新下一个点的值：y_{n+1} = y_n + h*k2

很多人容易搞混的点：

• 把中点的步长写成h而不是h/2
• 错误用Heun法的公式（y_{n+1} = y_n + h*(k1+k2)/2）来替代中点法
• 微分方程函数f的参数顺序搞反（比如把f(y, t)写成f(t, y)）

给你一个可运行的正确示例代码 #

假设我们求解的是微分方程dy/dt = -y（解析解是y=e^(-t)），下面是正确的RK2中点法实现：

import numpy as np
import matplotlib.pyplot as plt

# 定义微分方程 dy/dt = f(t, y)
def f(t, y):
    return -y

def rk2_midpoint(f, t_start, y_start, step_size, total_steps):
    # 初始化时间数组和结果数组
    t = np.linspace(t_start, t_start + total_steps*step_size, total_steps + 1)
    y = np.zeros_like(t)
    y[0] = y_start  # 初始值
    
    # 迭代计算每个点
    for n in range(total_steps):
        # 计算k1
        k1 = f(t[n], y[n])
        # 计算中点的t和y
        t_mid = t[n] + step_size / 2
        y_mid = y[n] + (step_size * k1) / 2
        # 计算k2
        k2 = f(t_mid, y_mid)
        # 更新下一个y值
        y[n+1] = y[n] + step_size * k2
    
    return t, y

# 测试用参数
t0 = 0
y0 = 1
h = 0.1
steps = 20

# 运行RK2
t, y_rk2 = rk2_midpoint(f, t0, y0, h, steps)
# 计算解析解
y_exact = np.exp(-t)

# 绘图对比
plt.plot(t, y_rk2, marker='o', label='RK2 Midpoint')
plt.plot(t, y_exact, linestyle='--', label='Exact Solution')
plt.xlabel('t')
plt.ylabel('y')
plt.legend()
plt.show()

你可以对照检查自己的代码 #

• 公式正确性：对比你的k1、k2计算和更新步骤，是不是和上面一致？
• 数组索引：你的循环是不是从0到total_steps-1，确保每个y[n+1]都被赋值？有没有出现索引越界或者漏更新的情况？
• 微分方程函数：确认f的参数顺序、逻辑是否正确，比如有没有把自变量和因变量搞反？
• 数组初始化：你的y数组是不是和t数组长度一致？有没有默认值没被覆盖的情况？

如果你把自己的代码贴出来，还能更精准地定位问题，但先从上面这几点排查应该能解决大部分问题~

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