Python中如何解方程[HTML代码]

　Numpy 求解线性方程组

　　例如我们要解一个这样的二元一次方程组：

　　x + 2y = 3

　　4x ＋ 5y = 6

　　当然我们可以手动写出解析解，然后写一个函数来求解，这实际上只是用 Python 来单纯做“数值计算”. 但实际上，numpy.linalg.solve 可以直接求解线性方程组.

　　一般地，我们设解线性方程组形如 Ax=b，其中 A 是系数矩阵，b 是一维(n 维也可以，这个下面会提到)，x 是未知变量. 再拿上面地最简单的二元一次方程组为例，我们用 numpy.linalg.solve 可以这样写：

　　In [1]: import numpy as np

　　...: A = np.mat('1,2; 4,5') # 构造系数矩阵 A

　　...: b = np.mat('3,6').T # 构造转置矩阵 b （这里必须为列向量）

　　...: r = np.linalg.solve(A,b) # 调用 solve 函数求解

　　...: print r

　　...:

　　Out[1]: [[-1.]

　　[ 2.]]

　　那么前面提到的“ n 维”情形是什么呢？实际上就是同时求解多组形式相同的二元一次方程组，例如我们想同时求解这样两组：

　　x + 2y = 3

　　4x ＋ 5y = 6

　　和

　　x + 2y = 7

　　4x ＋ 5y = 8

　　就可以这样写：

　　In [2]: import numpy as np

　　...: A = np.mat('1,2; 4,5') # 构造系数矩阵 A

　　...: b = np.array([[3,6], [7,8]]).T # 构造转置矩阵 b （这里必须为列向量），

　　...: 注意这里用的是 array

　　...: r = np.linalg.solve(A,b) # 调用 solve 函数求解

　　...: print r

　　...:

　　Out[2]: [[-1. -6.33333333]

　　[ 2. 6.66666667]]

　　SciPy 求解非线性方程组

　　先看官方文档的介绍：

　　scipy.optimize.fsolve(func, x0, args=(), fprime=None, full_output=0, col_deriv=0, xtol=1.49012e-08, maxfev=0, band=None, epsfcn=None, factor=100, diag=None)[source]

　　一般来说，我们只需要用到 func 和 x0 就够了. func 是自己构造的函数，也就是需要求解的方程组的左端（右端为 0），而 x0 则是给定的初值.

　　我们来看一个具体的例子，求解：

　　x + 2y + 3z - 6 = 0

　　5 * (x ** 2) + 6 * (y ** 2) + 7 * (z ** 2) - 18 = 0

　　9 * (x ** 3) + 10 * (y ** 3) + 11 * (z ** 3) - 30 = 0

　　就可以这么写：

　　In [3]: from scipy.optimize import fsolve

　　...:

　　...: def func(i):

　　...: x, y, z = i[0], i[1], i[2]

　　...: return [

　　...: x + 2 * y + 3 * z - 6,

　　...: 5 * (x ** 2) + 6 * (y ** 2) + 7 * (z ** 2) - 18,

　　...: 9 * (x ** 3) + 10 * (y ** 3) + 11 * (z ** 3) - 30

　　...:

　　...: r = fsolve(func,[0, 0, 0])

　　...: print r

　　...:

　　Out[3]: [ 1.00000001 0.99999998 1.00000001]

　　当然，SciPy 也可以用来求解线性方程组，这是因为 scipy.optimize.fsolve 本质上是最小二乘法来逼近真实结果.

　　SymPy 通吃一切

　　例如求解一个：

　　x + 2 * (x ** 2) + 3 * (x ** 3) - 6 = 0

　　直接就是：

　　In [4]: from sympy import *

　　...: x = symbols('x')

　　...: solve(x + 2 * (x ** 2) + 3 * (x ** 3) - 6, x)

　　Out[4]: [1, -5/6 - sqrt(47)*I/6, -5/6 + sqrt(47)*I/6]

　　另外，Wayne Shi 的这篇 使用 Python 解数学方程 ，就重点讲述了 SymPy 解线性方程组的方法，所以我也就不再赘述了。