用 sympy 区分方程组

Question

我希望对常微分方程组的时间求导数。例如，起点如下：

from sympy import *

t = Symbol('t')
s = Function('s')(t)
i = Function('i')(t)
beta = Symbol('beta')
gamma = Symbol('gamma')

eqs = [Eq(s.diff(t), -beta*s*i), Eq(i.diff(t), beta*s*i-gamma*i)]

我想获得 s 和 i 的二阶导数以及更高的时间导数（作为 s、i、beta 和 gamma 的函数）。

在我看来，自然的符号解决方案似乎是将 diff 应用于列表或列表中的方程式。这些似乎都不起作用，在 Whosebug 和 google 上搜索后，我没有任何迹象表明我还能尝试什么。

版本信息：

Sympy 1.10.1

Python 3.8.13

Answer 1

如您所见：

无法将导数应用于 Python 列表：diff 不知道如何处理 list 类型的对象，因为它们不是符号对象。
无法对 Equality 对象 (Eq) 应用导数，因为这些对象表示等式关系而不是数学方程，因此它们不支持数学运算。 :|

最好的方法是将等式写成 (Left Hand Side) - (Right Hand Side):

eq1 = s.diff(t) + beta*s*i
eq2 = i.diff(t) - beta*s*i-gamma*i

然后，你就可以区分它们了。例如：

eq1.diff(t)
# output: beta*i(t)*Derivative(s(t), t) + beta*s(t)*Derivative(i(t), t) + Derivative(s(t), (t, 2))

Answer 2

你已经有了方程的一阶导数。你可以再次区分它以获得二阶导数：

In [10]: for eq in eqs:
    ...:     pprint(eq.rhs)
    ...: 
-β⋅i(t)⋅s(t)
β⋅i(t)⋅s(t) - γ⋅i(t)

In [11]: for eq in eqs:
    ...:     pprint(eq.rhs.diff(t))
    ...: 
         d                 d       
- β⋅i(t)⋅──(s(t)) - β⋅s(t)⋅──(i(t))
         dt                dt      
       d                 d            d       
β⋅i(t)⋅──(s(t)) + β⋅s(t)⋅──(i(t)) - γ⋅──(i(t))
       dt                dt           dt

如果你想使用微分方程从那里消除 s'(t) 和 i'(t) 那么你可以这样做：

In [12]: rep = {eq.lhs: eq.rhs for eq in eqs}

In [13]: for eq in eqs:
    ...:     pprint(eq.rhs.diff(t).subs(rep))
    ...: 
 2  2                                        
β ⋅i (t)⋅s(t) - β⋅(β⋅i(t)⋅s(t) - γ⋅i(t))⋅s(t)
   2  2                                                                   
- β ⋅i (t)⋅s(t) + β⋅(β⋅i(t)⋅s(t) - γ⋅i(t))⋅s(t) - γ⋅(β⋅i(t)⋅s(t) - γ⋅i(t))

用 sympy 区分方程组

Differentiating a system of equations with sympy

python

symbolic-math

sympy