以下是你的方程式:
In [45]: from sympy import solve
...:
...: from sympy import symbols
...: x,y,a,b,c,d,g,h = symbols('x y a b c d g h')
...: a = 0.12
...: b = -2.031529100521498e-30
...: c = b
...: d = 1
...: g = 11
...: h = 0
...: F = (y + a/2)*x*(2*x**2 - 2*y*a + 2*d**2 + b + c) - (x*y*g**2)
...: G = (x**2 - y*a + b + d**2)*(x**2 - y*a + c+d**2) - 4*((y+a/2)**2)*x**2 - (h**2 + (x**2)*(g*
...: *2))
我要把浮点数转换成有理数，因为浮点数和多项式一般不能很好地混合.我只想在这里指出，您的值b非常小，以至于实际上是零(在标准的64位浮点中).如果希望b为非零，则需要从一开始就避免浮点数，因为在执行上面所示的代码后，浮点数已经消失了.
In [46]: F, G = nsimplify(F), nsimplify(G) # don't use floats with polynomials
In [47]: F
Out[47]:
⎛ 2 6⋅y ⎞
-121⋅x⋅y + x⋅(y + 3/50)⋅⎜2⋅x - ─── + 2⎟
⎝ 25 ⎠
In [48]: G
Out[48]:
2
2 2 2 ⎛ 2 3⋅y ⎞
- 4⋅x ⋅(y + 3/50) - 121⋅x + ⎜x - ─── + 1⎟
⎝ 25 ⎠
方程式F分解，所以让我们依次取这两个因子:
In [56]: F.factor()
Out[56]:
⎛ 2 2 2 ⎞
x⋅⎝1250⋅x ⋅y + 75⋅x - 150⋅y - 74384⋅y + 75⎠
─────────────────────────────────────────────
625
In [57]: _, F1, F2 = F.factor().args
In [58]: F1
Out[58]: x
In [59]: solve([F1, G], [x, y])
Out[59]: [(0, 25/3)]
我相信这就是你提到的解决方案.您认为这是不正确的，但对于代码中所示的方程式来说，这绝对是正确的.
这是最简单的一个.现在，对于F2，我们计算Groebner基数:
In [61]: gb = groebner([F2, G], [x, y], method='f5b')
In [62]: print(gb[0])
x**2 + 16*y**4/121 + 595216*y**3/9075 + 892608*y**2/75625 + 1844017554*y/1890625 - 75643/75625
In [63]: print(gb[1])
y**5 + 74411*y**4/150 + 148777*y**3/1250 + 1845583341*y**2/250000 + 691424307*y/781250 - 113451/125000
我们可以看到，gb[1]是y中的一个五次数.这个五次曲线是不可约的，并且可能没有根的公式(Abel-Ruffini)，但当系数是显式有理数时，SymPy可以将根表示为RootOf:
In [66]: r1, r2, r3, r4, r5 = Poly(gb[-1]).all_roots()
In [67]: r1
Out[67]:
⎛ 5 4 3 2
CRootOf⎝18750000⋅x + 9301375000⋅x + 2231655000⋅x + 138418750575⋅x + 16594183368⋅x - 17017650, 0
⎞
⎠
这些根中的每一个都可以在第一个方程中用来求解x，例如:
In [69]: print(solve([gb[0],y-r1], [x, y]))
[(-sqrt(-2250000*CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0)**4 - 200836800*CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0)**2 + 17019675 - 16596157986*CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0) - 1116030000*CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0)**3)/4125, CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0)), (sqrt(-2250000*CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0)**4 - 200836800*CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0)**2 + 17019675 - 16596157986*CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0) - 1116030000*CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0)**3)/4125, CRootOf(18750000*x**5 + 9301375000*x**4 + 2231655000*x**3 + 138418750575*x**2 + 16594183368*x - 17017650, 0))]
还有其他方法可以做到这一点，但这就给出了RootOf方面的所有解决方案.我不知道这是不是您想要的，但您可能希望SymPy返回一些不可能返回的东西.以下是全套解决方案:
In [71]: sx1, sx2 = solve(gb[0], x)
In [72]: sols = solve([F1, G], [x, y], dict=True)
In [73]: sols.extend([{x:sx1.subs(y, r), y:r} for r in [r1,r2,r3,r4,r5]])
In [74]: sols.extend([{x:sx2.subs(y, r), y:r} for r in [r1,r2,r3,r4,r5]])
这是它们在数字上的样子:
In [76]: for s in sols:
...: print(s[x].n(3), s[y].n(3))
...:
0 8.33
-0.0610 -496.
-10.9 -0.121
-0.0917 0.00102
-7.71 + 0.091*I -0.045 - 3.86*I
-7.71 - 0.091*I -0.045 + 3.86*I
0.0610 -496.
10.9 -0.121
0.0917 0.00102
7.71 - 0.091*I -0.045 - 3.86*I
7.71 + 0.091*I -0.045 + 3.86*I
这个问题很模糊，你是否有兴趣得到a，b等符号形式的解决方案，而不是明确的数字.gb[-1]是一个不可约的五次曲线，这一事实使得符号系数的情况不太可能存在任何有意义的显式解:
https://en.wikipedia.org/wiki/Abel%E2%80%93Ruffini_theorem个
一般来说，当面对这类问题时，最好考虑不同的方法.显式解可能是不可能的(就像这里)，或者太复杂而没有多大用处(例如，当使用Cardano的四次公式时).当我给出这个建议时，似乎很少人注意到它，但实际上，在这样的情况下，寻求explicit个解决方案要好得多.原始的多项式方程已经是解的更好的表示，而不是任何基本公式所能提供的:对于符号工作，使用这些多项式作为系统解的implicit表示通常更好.我不能建议如何做到这一点，因为我不知道你下一步打算做什么(如果你想要答案，请问一个新的问题).
