sympy —— Python 符号运算

0. 素因子分解

>> import sympy
>> from sympy.ntheory import factorint
>> factorint(100)
{2: 1, 3: 1, 5: 2}

1. sympy：Python 自带符号运算库

1. from sympy import *

2. 自带常量E和复数单位I

3. 简单测试

>>>x = Symbol('x', real=True)
>>>y = Symbol('y', real=True) 
>>>(x+y)**2
(x + y)**2
>>>expand((x+y)**2, real=True)
x**2 + 2*x*y + y**2

比如我们使用 sympy 来验证如下的等式：

$$\begin{split} x^2+y^2=&(x+(2xy)^{1/2}+y)(x-(2xy)^{1/2}+y)\\ \end{split}$$

>> from sympy import Symbol, expand
>> from fractions import Fraction

>> x = Symbol('x')
>> y = Symbol('y')
>> expand((x+(2*x*y)**Fraction(1, 2)+y)*(x-(2*x*y)**Fraction(1, 2)+y))
x**2 + y**2

该等式的简单证明如下：

$$\begin{split} x^2+y^2=&\left(x+y\right)^2-2xy\\ =&\left(x+y+\sqrt{2xy}\right)\left(x+y-\sqrt{2xy}\right) \end{split}$$

2. sympy：多项式运算

$$x\cdot \left(x^2+x-1\right)^2$$

>> from sympy import poly
>> from sympy.abc import x
>> poly(x*(x**2+x-1)**2)
Poly(x**5 + 2*x**4 - x**3 - 2*x**2 + x, x, domain='ZZ')

Polynomials Manipulation Module Reference

3. 解方程

使用 Python 解数学方程

• 解一元一次方程：$3x-6=0$

  >> from sympy import Symbol, solve
  >> x = Symbol('x')
  >> solve(x*3-6, x)
  [2]

• 解 二元一次方程组：

  $$\left\{ \begin{split} y=1-x\\ 3x+2y=5 \end{split} \right.$$

  >> x = Symbol('x')
  >> y = Symbol('y')
  >> f1 = y+x-1
  >> f2 = 3*x+2*y-5
  >> solve([f1, f2], [x, y])
  {x: 3, y: -2}

4. python 仿真

class array(object):
    def __init__(self, value, name=None):
        self.value = value
        if name:
            self.grad = lambda g: {name : g}
            # 定义梯度函数，返回一个字典

    # '+' 加号运算符重载
    def __add__(self, other):
        assert isinstance(other, int)
        ret = array(self.value+other)
        ret.grad = self.grad
        return ret

    # '*' 乘法运算符重载
    def __mul__(self, other):
        assert isinstance(other, array)
        ret = array(self.value*other.value)
        def grad(g):
            x = self.grad(g*other.value)
            x.update(other.grad(self.value*g))
                #  也即x*y，对x（self）的导数等于y（other），反之亦然
            return x
        ret.grad = grad
        return ret
x = array(1, 'x')
y = array(2, 'y')
w = x*y
z = w + 1
print(z.value)
print(z.grad(1))
            # z.grad是成员函数，而不是一个变量，或者一个字典