python - sympy 能解出这些多项式平方根方程吗?
问题描述
我尝试在 SymPy 中以 2 种不同形式求解一个方程,但存在错误和无限循环。这是在 Python 3.8.8 和 SymPy 1.7.1 上,在两个设备上的结果相同。相比之下,Mathematica 在一秒钟内解决了这两个问题。
方程
from sympy import *
x,xh,xhh,p,q,r = symbols('x xh xhh p q r')
eq1 = ((2*q**2-2*r**2)/(-q - r + xh) - 2*q - 2*r + x + xhh
+ sqrt(4*p*x + 4*p*xh + 4*q**2 + x**2 + 2*x*xh + xh**2)
+ sqrt(4*p*xh + 4*q**2 + xh**2 + xhh**2 + 2*xhh*(2*p + xh))
)
eq2 = (4*q**2 + 4*q*r - q*x - 2*q*xh - r*x - 2*r*xh + x*xh + xhh*(-q - r + xh)
+ (-q - r + xh)*sqrt(4*p*xh + 4*q**2 + xh**2 + xhh**2 + 2*xhh*(2*p + xh))
+ (-q - r + xh)*sqrt(4*p*x + 4*p*xh + 4*q**2 + x**2 + 2*x*xh + xh**2)
)
A,C,D,F,G,H = symbols('A C D F G H')
eq3 = A + xhh + C/D + sqrt(F) + sqrt(G + H*xhh + xhh**2)
解决尝试
解决第一种形式solve(eq1,xhh)
返回错误
Traceback (most recent call last):
File "<ipython-input-24-9923c6f02abd>", line 1, in <module>
solve(eq1,xhh)
File "/Users/jacobrichardson/anaconda/envs/env_sympy/lib/python3.8/site-packages/sympy/solvers/solvers.py", line 1097, in solve
solution = _solve(f[0], *symbols, **flags)
File "/Users/jacobrichardson/anaconda/envs/env_sympy/lib/python3.8/site-packages/sympy/solvers/solvers.py", line 1460, in _solve
gen = f_num.match(D.xreplace({d: w}))[w]
TypeError: 'NoneType' object is not subscriptable
在我的两台设备上都没有完成第二次solve(eq2,xhh)
运行一个多小时!
确实,在我的手动帮助下,第 3 种形式解决solve(eq3,xhh)
了给予
的问题[(-A**2*D**2 - 2*A*C*D - 2*A*D**2*sqrt(F) - C**2 - 2*C*D*sqrt(F) - D**2*F + D**2*G)/(D*(2*A*D + 2*C + 2*D*sqrt(F) - D*H))]
有什么方法可以让 SymPy 像Mathematicaeq1
和eq2
一样轻松方便地Solve[eq1,xhh]
求解方程Solve[eq2,xhh]
?
等价方程
print(simplify(eq1 * (-q - r + xh) - eq2)) # "0" equations 1 and 2 equivalent up to one singularity
print(simplify(eq1 -
eq3.subs([(A,-2*q - 2*r + x),(C,2*q**2 - 2*r**2),(D,-q - r + xh),
(F,4*p*x + 4*p*xh + 4*q**2 + x**2 + 2*x*xh + xh**2),
(G,4*p*xh + 4*q**2 + xh**2),(H,2*(2*p + xh))]))) # "0" equations 1 and 3 equivalent
解决方案
我试过solve(eq1, xhh)
了,我看到您在使用 sympy 1.7.1(最新版本)进行测试时显示的错误。我认为这是一个已经修复的错误,因为我没有看到最新的“主”版本的 sympy。事实上,这是修复它的拉取请求:
https://github.com/sympy/sympy/pull/20842
取而代之的是主人,我发现它在您报告时很慢eq2
。打断它,我可以看到可以禁用的检查代码很慢。它check=False
需要 4 秒eq1
和 5 秒eq2
:
In [2]: %time solve(eq1, xhh, check=False)
CPU times: user 4.21 s, sys: 45.3 ms, total: 4.25 s
Wall time: 4.33 s
Out[2]:
⎡ ______________________
⎢ 2 2 2 4 3 3 3 3 ╱ 2
⎢- 2⋅p⋅q ⋅x - 4⋅p⋅q⋅r⋅x + 4⋅p⋅q⋅x⋅xh - 2⋅p⋅r ⋅x + 4⋅p⋅r⋅x⋅xh - 2⋅p⋅x⋅xh - 8⋅q - 16⋅q ⋅r + 4⋅q ⋅x + 8⋅q ⋅xh + 4⋅q ⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q
⎢────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
⎢
⎢
⎣
_____________________ ___________________________________________
2 2 2 2 2 2 2 ╱ 2 2 2 2 2 2 2
+ x + 2⋅x⋅xh + xh - 8⋅q ⋅r + 8⋅q ⋅r⋅x + 16⋅q ⋅r⋅xh + 8⋅q ⋅r⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - q ⋅x - 7⋅q ⋅x⋅xh - q ⋅x⋅╲╱
─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
___________________________________________ ___________________________________________
╱ 2 2 2 2 2 2 ╱ 2 2 2 2 2
4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - 2⋅q ⋅xh - 6⋅q ⋅xh⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh + 4⋅q⋅r ⋅x + 8⋅q⋅r ⋅xh + 4⋅q
─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
____________
2 2 2 3 2 2 2 2 ╱
- 2⋅p⋅q - 4⋅p⋅q⋅r + 4⋅p⋅q⋅xh - 2⋅p⋅r + 4⋅p⋅r⋅xh - 2⋅p⋅xh - 4⋅q - 8⋅q ⋅r + q ⋅x + 5⋅q ⋅xh + q ⋅╲╱ 4⋅p⋅x + 4⋅p
___________________________________________ ___________________________________________
2 ╱ 2 2 2 2 ╱ 2 2 2
⋅r ⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - 2⋅q⋅r⋅x - 10⋅q⋅r⋅x⋅xh - 2⋅q⋅r⋅x⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - 4⋅q⋅r⋅
─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
_______________________________ ___________________________________________ _
2 2 2 2 ╱ 2 2 2 ╱
⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - 4⋅q⋅r + 2⋅q⋅r⋅x + 6⋅q⋅r⋅xh + 2⋅q⋅r⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - 2⋅q⋅x⋅xh - 2⋅q⋅xh⋅╲╱
___________________________________________ _________________________________________
2 ╱ 2 2 2 2 2 ╱ 2 2
xh - 8⋅q⋅r⋅xh⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh + 2⋅q⋅x ⋅xh + 4⋅q⋅x⋅xh + 2⋅q⋅x⋅xh⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh
─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
__________________________________________ ___________________________________________ ________
2 2 2 2 2 2 ╱ 2 2 2 ╱
4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh + r ⋅x + r ⋅xh + r ⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - 2⋅r⋅x⋅xh - 2⋅r⋅xh⋅╲╱ 4⋅p⋅x +
__ ___________________________________________ ___________________________________________
2 2 ╱ 2 2 2 2 2 2 2 ╱ 2 2 2 2
+ 2⋅q⋅xh ⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - r ⋅x - 3⋅r ⋅x⋅xh - r ⋅x⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - 2⋅r ⋅
─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
___________________________________ ___________________________________________
2 2 2 2 3 2 ╱ 2 2 2
4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh + x⋅xh - xh + xh ⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh
___________________________________________ __________________________________________
2 2 ╱ 2 2 2 2 2 ╱ 2 2 2
xh - 2⋅r ⋅xh⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh + 2⋅r⋅x ⋅xh + 4⋅r⋅x⋅xh + 2⋅r⋅x⋅xh⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh
─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
_ ___________________________________________ ___________________________________________⎤
2 ╱ 2 2 2 2 2 3 2 ╱ 2 2 2 ⎥
+ 2⋅r⋅xh ⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh - x ⋅xh - x⋅xh - x⋅xh ⋅╲╱ 4⋅p⋅x + 4⋅p⋅xh + 4⋅q + x + 2⋅x⋅xh + xh ⎥
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────⎥
⎥
⎥
⎦
推荐阅读
- python - Python:路径中的奇怪前缀
- c++ - MFC C++ 如何使编辑控件的背景透明?
- html - 如何在单个应用程序中拥有多个 Web 配置
- python - 用 Not In 列出理解
- php - Mysql选择结果作为文本变量
- java - Gif 循环运行而不是单次迭代
- python - 解析原始文本数据并在 Python 中提取特定值
- kubernetes - 如何使用 kubectl 通过 jsonpath 获取输出
- ios - 在 iOS 上为 React-Native 添加 Emoji 到 TextInput 的问题
- c# - DatabaseGenerated(DatabaseGeneratedOption.Identity) vs Key