2013-03-02 83 views
4

我正在做一個非常簡單的概率計算,從一組A-Z(具有相應的概率x,y,z)獲取X,Y,Z的子集。在sympy中分解多邊形

因爲非常沉重的公式而且,爲了處理他們,我想簡化(或收集因素 - 我不知道確切的定義)使用sympy這些多項式。

所以..有這個

import sympy as sp 

x, y, z = sp.symbols('x y z') 

expression = (
    x * (1 - x) * y * (1 - x - y) * z + 
    x * (1 - x) * z * (1 - x - z) * y + 

    y * (1 - y) * x * (1 - y - x) * z + 
    y * (1 - y) * z * (1 - y - z) * x + 

    z * (1 - z) * y * (1 - z - y) * x + 
    z * (1 - z) * x * (1 - z - x) * y 
) 

我想要得到的東西(一個非常簡單的概率計算得到X,Y,從集AZ的的Z子集與相應的概率X,Y,Z的表達)這樣

x * y * z * (6 * (1 - x - y - z) + (x + y) ** 2 + (y + z) ** 2 + (x + z) ** 2) 

聚,在改寫的方式有儘可能少的操作(+-***,...)儘可能


我試過factor()collect()simplify()。但結果與我的預期不同。主要是我得到

2*x*y*z*(x**2 + x*y + x*z - 3*x + y**2 + y*z - 3*y + z**2 - 3*z + 3) 

我知道sympy可以結合多項式成簡單的形式:

sp.factor(x**2 + 2*x*y + y**2) # gives (x + y)**2 

但如何讓sympy到從表達式組合多項式以上?


如果在sympy中這是不可能的任務,可能還有其他的選擇嗎?

回答

3

把這些方法放在一起恰好給出了一個很好的答案。看看這個策略是否比你產生的方程更經常地工作,或者如果顧名思義這是一次幸運的結果,這將是很有趣的。

def iflfactor(eq): 
    """Return the "I'm feeling lucky" factored form of eq.""" 
    e = Mul(*[horner(e) if e.is_Add else e for e in 
     Mul.make_args(factor_terms(expand(eq)))]) 
    r, e = cse(e) 
    s = [ri[0] for ri in r] 
    e = Mul(*[collect(ei.expand(), s) if ei.is_Add else ei for ei in 
     Mul.make_args(e[0])]).subs(r) 
    return e 

>>> iflfactor(eq) # using your equation as eq 
2*x*y*z*(x**2 + x*y + y**2 + (z - 3)*(x + y + z) + 3) 
>>> _.count_ops() 
15 

BTW,factor_terms和gcd_terms之間的不同之處在於factor_terms將更加努力,拉出來的常用術語,同時保留表達時,原來的結構非常像,你會做手工(即尋找常用術語添加可以拉出)。

>>> factor_terms(x/(z+z*y)+x/z) 
x*(1 + 1/(y + 1))/z 
>>> gcd_terms(x/(z+z*y)+x/z) 
x*(y*z + 2*z)/(z*(y*z + z)) 

對於它的價值,

克里斯

+0

尼斯組合,但它很難理解。你能解釋算法/想法嗎?順便說一句,我發現了一個明顯的錯誤:'x *(1-x)* y *(1-x-y)* z + ...' - >'x /(1 - x)* y /(1 - x - y )* z + ...',對於這樣的情商你的組合不起作用(我想這是因爲明顯的事情,但因爲我不知道algorythm ....) – akaRem 2013-03-04 21:54:09

1

據我所知,沒有什麼功能可以做到這一點。我相信這實際上是一個非常難的問題。有關它的一些討論見Reduce the number of operations on a simple expression

然而,SymPy中有很多簡化函數可供您嘗試。一個你沒有提到過的結果是gcd_terms,它沒有進行擴展就將符號化的gcd分解出來。它給出了

>>> gcd_terms(expression) 
x*y*z*((-x + 1)*(-x - y + 1) + (-x + 1)*(-x - z + 1) + (-y + 1)*(-x - y + 1) + (-y + 1)*(-y - z + 1) + (-z + 1)*(-x - z + 1) + (-z + 1)*(-y - z + 1)) 

另一個有用的功能是.count_ops,它計算表達式中的操作數。例如

>>> expression.count_ops() 
47 
>>> factor(expression).count_ops() 
22 
>>> e = x * y * z * (6 * (1 - x - y - z) + (x + y) ** 2 + (y + z) ** 2 + (x + z) ** 2) 
>>> e.count_ops() 
18 

(注意:e.count_ops()是不一樣的,你算你自己,因爲SymPy自動6*(1 - x - y - z)分發到6 - 6*x - 6*y - 6*z)。

其它有用的功能:

  • cse:執行上表達的公共子表達式消除。有時你可以簡化個別部分,然後再將它們放在一起。這也有助於避免重複的計算。

  • horner:將Horner scheme應用於多項式。如果多項式在一個變量中,這將最小化操作次數。

  • factor_terms:類似於gcd_terms。其實我並不完全清楚它們有什麼不同。

注意,默認情況下,simplify會嘗試一些簡化,並返回由count_ops最小的一個。