求支持自适应步长的Python单步ODE求解器

核心解决办法：用迭代式自适应求解器 #

Python的scipy.integrate模块里有迭代型的自适应步长ODE求解器（比如RK45、DOP853），这类求解器支持手动调用step()方法执行单步自适应积分——完全符合你的需求：只计算求解器自身选择的下一个时间点，不会浪费算力去算后续多余的步骤。

具体代码示例 #

下面是封装成你理想调用形式的代码，专门针对step_num=1的场景（适配你每步更新ODE函数的需求）：

from scipy.integrate import RK45
from math import sin, sqrt
import matplotlib.pyplot as plt
import numpy as np

def solve(ode_func, t_span, init_cond, step_num):
    t0, t_final = t_span
    y0 = init_cond
    # 初始化自适应求解器
    solver = RK45(ode_func, t0, y0, t_final)
    
    sol_t = [t0]
    sol_y = [y0]
    
    for _ in range(step_num):
        # 执行一步自适应积分
        solver.step()
        # 记录当前步的结果
        sol_t.append(solver.t)
        sol_y.append(solver.y.copy())
        
        # 如果已经到达时间终点，提前终止
        if solver.status == 'finished':
            break
    
    # 返回初始点之外的结果（对应t₁到tₙ）
    return sol_t[1:], sol_y[1:]

# 测试用ODE函数，和你Matlab示例中的一致
def ode_func(t, y):
    return 0.25*(1+sin(t)) - 0.1*sqrt(y)

# 初始化参数
current_t = 0
current_y = [1.0]
t_max = 6
sol_t = [current_t]
sol_y = [current_y.copy()]

# 模拟你的循环场景：每步可更新ODE，每次只计算1步
while current_t < t_max:
    # 这里可根据需求动态更新ode_func，比如基于当前解调整方程参数
    # 示例：ode_func = new_equation(current_t, current_y)
    
    # 仅计算1步自适应积分
    next_t, next_y = solve(ode_func, (current_t, t_max), current_y, step_num=1)
    current_t = next_t[0]
    current_y = next_y[0]
    
    sol_t.append(current_t)
    sol_y.append(current_y)
    
    # 此处放入你的自定义操作
    # do something here

# 生成和Matlab示例一致的茎状图
plt.stem(sol_t, np.array(sol_y).flatten())
plt.xlabel('时间')
plt.ylabel('解')
plt.show()

重点说明 #

• 效率优势：RK45这类求解器的step()方法会自动计算当前最优自适应步长，直接推进到下一个时间点，不像你Matlab的实现那样要先计算整个时间窗再丢弃多余结果，完全避免算力浪费。
• 动态更新支持：因为每次只执行一步，你可以在每步之后任意修改ode_func（比如调整方程参数），再传入下一次求解，完美适配你的场景。
• 可选求解器：除了RK45，还可以使用DOP853（更高精度）、Radau（适配刚性ODE）等，用法完全一致，只需替换求解器类即可。

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