当我不知道解析导数时求解 Python 中的微分方程

Question

我正在尝试求解 Python 中的 first-order ODE：

其中 Gamma 和 u 是方阵。我并不是一直都明确知道 u(t)，但我确实在进行早期计算时在离散时间步知道它。

我在网上找到的每个 Python 求解器的例子（例如 this one 对 scipy.integrate.odeint 和 scipy.integrate.ode）都知道导数的表达式是解析的函数时间.

有没有办法在不知道导数的解析表达式的情况下调用这些（或其他微分方程求解器）？

现在，我已经编写了自己的 Runge-Kutta 求解器并用 numba 对它进行了 jit 处理。

Answer 1

您可以使用 SciPy interpolation methods, such as interp1d 中的任何一个，根据您的离散数据创建一个可调用函数，并将其传递给 odeint。三次样条插值，

f = interp1d(x, y, kind='cubic')

应该够好了。

Answer 2

Is there a way to call these (or other differential equation solvers) without knowing an analytic expression for the derivative?

是的，none 您提到的求解器（大多数其他求解器）都需要导数的解析表达式。相反，他们调用您提供的函数，该函数必须评估给定时间和状态的导数。因此，您的代码大致如下所示：

def my_derivative(time,flat_Gamma):
    Gamma = flat_Gamma.reshape(dim_1,dim_2)
    u = get_u_from_time(time)
    dGamma_dt = u.dot(Gamma)
    return dGamma_dt.flatten()

from scipy.integrate import ode
my_integrator = ode(my_derivative)
…

你的情况的困难在于你必须确保 get_u_from_time 为每次调用它提供适当的结果。可能最可靠和最简单的解决方案是使用插值法（参见）。

您也可以尝试将您的集成步骤与您拥有的数据相匹配，但至少对于 scipy.integrate.odeint 和 scipy.integrate.ode 这将非常乏味，因为所有集成商都使用不方便的内部步骤以此目的。例如，the fifth-order Dormand–Prince method (DoPri5) uses internal steps of 1/5, 3/10, 4/5, 8/9, and 1. This means that if you have temporally equidistant data for u, you would need 90 data points for each integration step (as 1/90 is the greatest common divisor of the internal steps). The only integrator that could make this remotely feasible is the Bogacki–Shampine integrator (RK23) from cipy.integrate.solve_ivp 内部步长为 1/2、3/4 和 1。

当我不知道解析导数时求解 Python 中的微分方程

Solve differential equation in Python when I don't know the derivative analytically

python

numpy

scipy

differential-equations

numba