我想用python(sympy(求解一个由4个耦合微分方程组成的系统：
eqs = [Eq(cP1(t).diff(t), k1*cE1(t)**3), Eq(cE1(t).diff(t), -k1 * cE1(t)**3 + k6 * cE3(t)**2), Eq(cE2(t).diff(t), -k8 * cE2(t)), Eq(cE3(t).diff(t), k8 * cE2(t) - k6 * cE3(t)**2)]
当我试图用"；dsolve_system"：
solution = dsolve_system(eqs, ics={cP1(0): 0, cE1(0): cE1_0, cE2(0):cE2_0, cE3(0):cE3_0})
我得到的答案是：；NotImplementedError:
传递的ODE系统无法由dsolve_system解决">
有人知道问题出在哪里吗？或者在Sympy中有更好的方法来求解这个微分方程组吗
非常感谢！
显示完整的代码通常是友好的：
In [18]: cP1, cE1, cE2, cE3 = symbols('cP1, cE1:4', cls=Function)
In [19]: t, k1, k6, k8 = symbols('t, k1, k6, k8')
In [20]: eqs = [Eq(cP1(t).diff(t), k1*cE1(t)**3), Eq(cE1(t).diff(t), -k1 * cE1(t)**3 + k6 * cE3(t)**2),
...: Eq(cE2(t).diff(t), -k8 * cE2(t)), Eq(cE3(t).diff(t), k8 * cE2(t) - k6 * cE3(t)**2)]
...:
In [21]: for eq in eqs:
...: pprint(eq)
...:
d 3
──(cP₁(t)) = k₁⋅cE₁ (t)
dt
d 3 2
──(cE₁(t)) = - k₁⋅cE₁ (t) + k₆⋅cE₃ (t)
dt
d
──(cE₂(t)) = -k₈⋅cE₂(t)
dt
d 2
──(cE₃(t)) = - k₆⋅cE₃ (t) + k₈⋅cE₂(t)
dt
In [22]: dsolve(eqs)
---------------------------------------------------------------------------
NotImplementedError Traceback (most recent call last)
<ipython-input-22-69ab769a7261> in <module>
----> 1 dsolve(eqs)
~/current/sympy/sympy/sympy/solvers/ode/ode.py in dsolve(eq, func, hint, simplify, ics, xi, eta, x0, n, **kwargs)
609 "number of functions being equal to number of equations")
610 if match['type_of_equation'] is None:
--> 611 raise NotImplementedError
612 else:
613 if match['is_linear'] == True:
NotImplementedError:
这意味着dsolve还不能处理这种特定类型的系统。注意，在一般情况下，常微分方程的非线性系统不太可能有解析解(dsolve用于寻找解析解，如果你想要数值解，可以使用类似scipy的odeint(。
就非线性系统而言，这一个相对友好，因此可能会解决它。让我们看看。。。
首先，有一个守恒量(所有4个变量的总和(，我们可以用来消除一个方程。事实上，这并没有多大帮助，因为第一个方程已经与其他方程分离：如果我们知道cE1，我们可以积分，但如果其他变量已知，守恒量会更容易地给出。
该系统的结构如下：
cE2 ---> cE3 ---> cE1 ---> cP1
这意味着它可以作为一个常微分方程序列来求解，其中我们先求解cE2的第三个方程，然后求解cE3的第四个方程，再使用它来求解cE1等等。因此，我们可以将其简化为一个涉及多个非线性单常微分方程的序列的问题。这也不太可能有一个分析解决方案，但让我们试一试：
In [24]: dsolve(eqs[2], cE2(t))
Out[24]:
-k₈⋅t
cE₂(t) = C₁⋅ℯ
In [25]: cE2sol = dsolve(eqs[2], cE2(t)).rhs
In [26]: eqs[3].subs(cE2(t), cE2sol)
Out[26]:
d -k₈⋅t 2
──(cE₃(t)) = C₁⋅k₈⋅ℯ - k₆⋅cE₃ (t)
dt
在这一点上，原则上我们可以在这里求解cE3，但sympy没有任何方法来求解这个特定的非线性ODE，所以dsolve给出了一个级数解(我不认为这是你想要的(，唯一一个可能处理这个问题的解算器是lie_group，但它实际上失败了。
由于我们不能得到cE3的解的表达式，我们也不能继续求解cE1和cP1。失败的ODE是Riccati方程，但这种特殊类型的Ricatti方程尚未在dsolve中实现。看起来WolframAlpha用贝塞尔函数给出了一个答案：https://www.wolframalpha.com/input/?i=solve+dx%2Fdt+%3D+e%5E-t+-+x%5E2
从这个角度来看，我想我们不太可能求解下一个方程。在这一点上Wolfram Alpha也放弃了并说"；分类：非线性微分方程组"：https://www.wolframalpha.com/input/?i=solve+dx%2Fdt+%3D+e%5E-t+-+x%5E2%2C+dy%2Fdt++%3D+-y%5E3+%2B+x%5E2
我怀疑你的系统没有解析解，但你可以尝试数值解或其他更定性的分析。