如何用scipy.integrate.solve_ivp将微分方程组解代入另一方程组？

• 利用稠密输出函数获取连续的x(t)
  你在求解原方程组时已设置dense_output=True，这会让solve_ivp返回的sol对象包含一个可调用的sol.sol方法——它能对任意时刻t进行插值，返回该时刻的状态向量[x(t), y(t)]。这是获取连续x(t)的核心方式，避免了仅使用离散的sol.y[0]带来的精度局限。

• 修改伴随方程组函数，引入原解的插值方法
  调整adjunta函数，让它接收sol、omega1等必要参数，并通过sol.sol(t)获取当前时刻的x(t)值：

def adjunta(t, v, mu, omega1, sol):
    zx, zy = v
    # 获取当前t对应的x(t)
    x_t = sol.sol(t)[0]
    dzxdt = -zx*(mu/omega1)*(1 - x_t**2) - zy/(mu*omega1)
    dzydt = (zx*mu)/omega1
    return (dzxdt, dzydt)

• 求解伴随方程组
  调用solve_ivp求解伴随方程时，将mu、omega1和sol作为参数传入args，同时设置合适的初始条件和时间区间：

# 根据研究需求设置伴随方程的初始条件
z0 = [1.0, 0.0]
# 求解伴随方程组，时间区间可与原方程组保持一致
adj_sol = solve_ivp(adjunta, [-tend, tend], z0, args=[mu, omega1, sol], max_step=1e-2)

• 关键注意点
  • sol.sol(t)支持标量和数组形式的t，完全适配solve_ivp内部对导数的计算需求；
  • 若伴随方程的求解区间超出原方程组的求解范围，插值会进入外推模式，精度会下降，需确保原求解区间足够覆盖需求。

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