如何解決多元方程組編程？

Question

我試圖解決多元方程，這是一些Java代碼的結果的系統。在運行前，變量的形式和數量都是已知的。一個例子是

(I) (e-a*d*e-b*d*e+2*b*d*f+2*b*d*e*g)/(-1+a*d+b*d)+f == 0 

(II) e*g+((f+e*g)*a*d)/(-1+a*d+b*d)==0 

(III) -e*h+((-f-e*g)*d)/(-1+a*d+b*d)==0 

(IV) -e*j+((-f-e*g)*c)/(-1+a*d+b*d)==0 

我嘗試使用Symja，它簡單地返回輸入，和SymPy，這會引發一個錯誤

ZeroDivisionError: polynomial division 

的變量是所有從區間[0,1]，和我需要所有解決方案。數學是能夠解決這個問題，但因爲它是商業軟件，我很遺憾不能在這個項目中使用它。

我將要在其上使用的軟件的任何建議表示感謝。我真的會喜歡SymPy上班，我不明白爲什麼它會拋出這個錯誤，想法表示讚賞。下面一MWE到SymPy錯誤：

from sympy.solvers import solve 
from sympy.abc import a,b,c,d,e,f,g,h,j 

lst = a,b,c,d,e,f,g,h,j 
sys = [(e-a*d*e-b*d*e+2*b*d*f+2*b*d*e*g)/(-1+a*d+b*d)+f,e*g+((f+e*g)*a*d)/(-1+a*d+b*d),-e*h+((-f-e*g)*d)/(-1+a*d+b*d),-e*j+((-f-e*g)*c)/(-1+a*d+b*d)] 

solution = solve(sys, lst) 
print solution 

Mathematica的版本是：

eqn = {(e - a*d*e - b*d*e + 2*b*d*f + 2*b*d*e*g)/(-1 + a*d + b*d) + f == 0, e*g + ((f + e*g)*a*d)/(-1 + a*d + b*d) == 0, -e*h + ((-f - e*g)*d)/(-1 + a*d + b*d) == 0, -e*j + ((-f - e*g)*c)/(-1 + a*d + b*d) == 0}; 
Simplify[Solve[eqn, {a, b, c, d, e, f, g, h, j}]] 

輸出：

{{e -> 0, f -> 0}, 
{c -> (1 - 2 a d - 3 b d) j, f -> ((-1 + 2 a d + b d) e)/(-1 + 2 a d + 3 b d), g -> (a d)/(1 - 2 a d - 3 b d), h -> d/(1 - 2 a d - 3 b d)}, 
{a -> 0, c -> j - 3 b d j, f -> ((-1 + b d) e)/(-1 + 3 b d), g -> 0, h -> d/(1 - 3 b d)}, 
{a -> (1 - b d)/(2 d), c -> -2 b d j, f -> 0, g -> 1/4 - 1/(4 b d), h -> -(1/(2 b))}} 

Answer 1

如果這個問題被關閉，這裏是在數學的嘗試版。但是它目前缺少限制變量到區間[0,1]的條件。

Solve[{ 
    (e - a d e - b d e + 2 b d f + 2 b d e g)/(-1 + a d + b d) + f == 0, 
    e g + ((f + e g) a d)/(-1 + a d + b d) == 0, 
    -e h + ((-f - e g) d)/(-1 + a d + b d) == 0, 
    -e j + ((-f - e g) c)/(-1 + a d + b d) == 0 
    }, 
{a, b, c, d, e, f, g, h, j}] 

Answer 2

好吧，我能解決自己使用SAGE的問題，與給出的輸出作爲數學一樣。

sage: a,b,c,d,e,f,g,h,j = var('a,b,c,d,e,f,g,h,j') 
sage: qe = [(e-a*d*e-b*d*e+2*b*d*f+2*b*d*e*g)+f*(-1+a*d+b*d),e*g*(-1+a*d+b*d)+((f+e*g)*a*d),-e*h*(-1+a*d+b*d)+((-f-e*g)*d),-e*j*(-1+a*d+b*d)+((-f-e*g)*c)] 
sage: print(solve(qe,a,b,c,d,e,f,g,h,j, solution_dict=True)) 

給人的輸出

[{g: r5, j: r7, b: r2, d: r4, e: 0, h: r6, c: r3, f: 0, a: r1}, 
{g: r12, j: r14, b: r8, d: r10, e: r11, h: r13, c: r9, f: -r11*r12, a: -(r10*r8 - 1)/r10}, 
{g: r15*r18, j: r19, b: -1/6*(4*r15*r16*r18 - r16 + 3*r17)*(4*r15*r16*r18^2/(4*r15*r16*r18 - r16 + 3*r17) + 2*r16*r18/(4*r15*r16*r18 - r16 + 3*r17) - r18)/(r16*r18^2), d: 2*r16*r18/(4*r15*r16*r18 - r16 + 3*r17), e: r16, h: r18, c: 2*r16*r19/(4*r15*r16*r18 - r16 + 3*r17), f: r17, a: r15}, 
{g: r20*r23, j: 0, b: -1/6*(4*r20*r21*r23 - r21 + 3*r22)*(4*r20*r21*r23^2/(4*r20*r21*r23 - r21 + 3*r22) + 2*r21*r23/(4*r20*r21*r23 - r21 + 3*r22) - r23)/(r21*r23^2), d: 2*r21*r23/(4*r20*r21*r23 - r21 + 3*r22), e: r21, h: r23, c: 0, f: r22, a: r20}, 
{g: -2*r24*r27 - 1, j: 0, b: r24, d: r25, e: r26, h: r27, c: 0, f: 2*r24*r26*r27 + r26, a: -(r24*r25 - 1)/r25}, 
{g: 0, j: r30, b: r29, d: 0, e: r31, h: 0, c: r30, f: r31, a: r28}, 
{g: 0, j: 0, b: r33, d: 0, e: r34, h: 0, c: 0, f: r34, a: r32}, 
{g: -1, j: 0, b: 0, d: r35, e: r36, h: r37, c: 0, f: r36, a: 1/r35}, 
{g: r40, j: 0, b: 0, d: r38, e: r39, h: 2*r38*r40 + r38, c: 0, f: r39, a: r40/(2*r38*r40 + r38)}, 
{g: -1/2, j: 0, b: -1/2*(r41 - r42)/(r41*r43), d: -2/3*r41*r43/(r41 - r42), e: r41, h: r43, c: 0, f: r42, a: -1/2/r43}] 

Answer 3

注意，你有9個變量和4個方程。因此，可以消除四個變量 - 你需要告訴Sympy其中。

下面的代碼首先乘以分母出和消除a,b,c,d：

import sympy as sy 
from IPython.display import display # for pretty printing 

# sy.init_printing() # LaTeX-like pretty printing for IPython 

a, b, c, d, e, f, g, h, j = sy.symbols("a, b, c, d, e, f, g, h, j", real=True) 
lst = a, b, c, d, e, f, g, h, j 
sys0 = sy.Matrix([(e-a*d*e-b*d*e+2*b*d*f+2*b*d*e*g)/(-1+a*d+b*d)+f, 
        e*g+((f+e*g)*a*d)/(-1+a*d+b*d), 
        -e*h+((-f-e*g)*d)/(-1+a*d+b*d), 
        -e*j+((-f-e*g)*c)/(-1+a*d+b*d)]) 

# Denominator can be factored out: 
den = a*d + b*d - 1 
sys1 = sy.simplify(sys0*den).expand() 
print("Factored out denominator:") 
display(sys1) 

# Elimnate four variables: 
sol1 = sy.solve(sys1, a, b, c, d, dict=True) 
print("Solutions:") 
display(sol1) 
print("Substituting back into the equation gives obviously 0, i.e.:") 
display(sy.simplify(sys1.subs(sol1[0])).T) 

print("The denominator != 0 results in:") 
den1 = sy.simplify(den.subs(sol1[0])) 
display(sy.solve(den1)) 

產生以下輸出：

Factored out denominator: 
Matrix([ 
[-a*d*e + a*d*f + 2*b*d*e*g - b*d*e + 3*b*d*f + e - f], 
[     2*a*d*e*g + a*d*f + b*d*e*g - e*g], 
[    -a*d*e*h - b*d*e*h - d*e*g - d*f + e*h], 
[    -a*d*e*j - b*d*e*j - c*e*g - c*f + e*j]]) 
Solutions: 
[{d: 2*e*h/(4*e*g - e + 3*f), 
    c: 2*e*j/(4*e*g - e + 3*f), 
    a: g/h, 
    b: (-e + f)/(2*e*h)}] 
Substituting back into the equation gives obviously 0, i.e.: 
Matrix([[0, 0, 0, 0]]) 
The denominator != 0 results in: 
[{e: -f/g}] 

所以所得到的標識是：

-f/g != e # denominator - only false, if f,e=0 since f,g,e>=0 
a = g/h 
b = (-e + f)/(2*e*h) 
c = 2*e*j/(4*e*g - e + 3*f) 
d = 2*e*h/(4*e*g - e + 3*f) 

據我所知，Sympy無法處理multiv還存在着不平等（儘管我想證明是錯誤的）。但結果是很簡單做手工：

0 <= g <= h 
0 <= f-e <= 2*e*h 
0 <= 2*e*j <= 4*e*g - e + 3*f 
0 <= 2*e*h <= 4*e*g - e + 3*f 
0 <= e,f,g,h,j <= 1 

的情況下f,e=0是有效的爲好。可以通過檢查sys0.subs(e,0).diff(f)不依賴於f來驗證。