首页 > 解决方案 > 用 Python 求解积分-微分耦合方程组

问题描述

我正在尝试使用 Python 对本文中描述的方程组进行数值求解,Eqs。30 和 31,简化形式如下:

在此处输入图像描述

其中G(k)D(k)是一些已知函数,与 无关 Y。当然,所有的量也是的函数t。作者评论说,由于各种函数表现出依赖性,因此需要数值解。

例如,我通常按照此处此处所示的方式实现此类耦合方程的解,但现在额外的 k 依赖性让我有点困惑。

有什么建议么?非常感谢。

标签: pythonscipy

解决方案


IDESolver是由Josh Karpel创建的通用数值积分微分方程求解器。它的最新版本允许用户解决多维、耦合的 IDE。从提供的示例中,IDE 像

在此处输入图像描述

使用解析解(sin x, cos x),可以使用以下代码解决:

import numpy as np
import matplotlib.pyplot as plt
from idesolver import IDESolver

solver = IDESolver(
            x=np.linspace(0, 7, 100),
            y_0=[0, 1],
            c=lambda x, y: [0.5 * (y[1] + 1), -0.5 * y[0]],
            d=lambda x: -0.5,
            f=lambda y: y,
            lower_bound=lambda x: 0,
            upper_bound=lambda x: x,
)

solver.solve()

fig = plt.figure(dpi = 600)
ax = fig.add_subplot(111)

exact = [np.sin(solver.x), np.cos(solver.x)]

ax.plot(solver.x, solver.y[1], label = 'IDESolver Solution', linestyle = '-', linewidth = 3)
ax.plot(solver.x, exact[1], label = 'Analytic Solution', linestyle = ':', linewidth = 3)

ax.legend(loc = 'best')
ax.grid(True)

ax.set_title(f'Solution for Global Error Tolerance = {solver.global_error_tolerance}')
ax.set_xlabel(r'$x$')
ax.set_ylabel(r'$y(x)$')

plt.show()

推荐阅读