Python 入门基础知识 - 科学计算 - 矩阵运算和解线性方程

矩阵运算和解线性方程

使用SciPy可以完成矩阵分解、线性方程组求解和多项式求根等数学运算。SciPy中提供的函数名与MATLAB中的函数名大部分相同，用法也差不多。熟悉SciPy的学者可以很快熟悉SciPy。

矩阵运算

除了基本的矩阵乘除运算以外，还可以使用SciPy做矩阵分解运算、求逆运算等。SciPy中的linalg模块提供了和线性代数相关的函数，可以用于对矩阵进行运算。常用的linalg模块中常用的函数如下

inv 求解矩阵的逆
det 求解方阵的行列式
norm 求解向量的模
lstsq 求最小二乘法解
eig 求解特征值和特征向量
orth 求解矩阵的标准正交集
lu 矩阵的LU分解
qr 矩阵的QR分解
cholesky 矩阵的Cholesky分解

如下实例演示，在Python的交互式命令行中使用SciPy中的linalg模块对矩阵进行运算。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
>>> import scipy # 导入SciPy模块
>>> from scipy import linalg 导入linalg模块
>>> a = scipy.mat('[1 2 3; 2 2 1; 3 4 3]') 生成矩阵
>>> a
matrix([[1, 2, 3],
        [2, 2, 1],
        [3, 4, 3]])
>>> b = linalg.inv(a) # 求a的逆矩阵
>>> b
array([[ 1. ,  3. , -2. ],
       [-1.5, -3. ,  2.5],
       [ 1. ,  1. , -1. ]])
>>> a * b # 求a * b ，有误差，应为对角矩阵
matrix([[  1.00000000e+00,   0.00000000e+00,   0.00000000e+00],
        [ -2.22044605e-16,   1.00000000e+00,   4.44089210e-16],
        [  0.00000000e+00,   0.00000000e+00,   1.00000000e+00]])
>>> linalg.det(a) # 求a的行列式值，有误差，应为2
1.9999999999999996
>>> a = scipy.mat('[1 2 3]')
>>> linalg.norm(a) # 求a的模
3.7416573867739413
>>> a = scipy.mat('[1 -1; -1 1; 2 -2; -3 1]') # 建立系数矩阵
>>> a # 检验系数矩阵
matrix([[ 1, -1],
        [-1,  1],
        [ 2, -2],
        [-3,  1]])
>>> b = scipy.mat('[1;2;3;4]') # 建立常数项向量
>>> b
matrix([[1],
        [2],
        [3],
        [4]])
>>> x,y,z,w = linalg.lstsq(a, b) # 求最小二乘法解
>>> x # x为解
array([[-2.41666667],
       [-3.25      ]])
>>> y
array([ 9.83333333])
>>> z
2
>>> w
array([ 4.56605495,  1.07291295])
>>> a = scipy.mat('[-1 1 0; -4 3 0; 1 0 2]')
>>> a
matrix([[-1,  1,  0],
        [-4,  3,  0],
        [ 1,  0,  2]])
>>> x,y = linalg.eig(a) # 求a的特征值
>>> x
array([ 2.+0.j,  1.+0.j,  1.+0.j]) # a的特征值
>>> y # a的特征向量
array([[ 0.        ,  0.40824829,  0.40824829],
       [ 0.        ,  0.81649658,  0.81649658],
       [ 1.        , -0.40824829, -0.40824829]])
>>> a = scipy.mat('[]1 2 3; 0 1 2; 2 4 1')
>>> a
matrix([[1, 2, 3],
        [0, 1, 2],
        [2, 4, 1]])
>>> x,y,z = linalg.lu(a) # 求a的LU分解
>>> x
array([[ 0.,  0.,  1.],
       [ 0.,  1.,  0.],
       [ 1.,  0.,  0.]])
>>> y # L矩阵
array([[ 1. ,  0. ,  0. ],
       [ 0. ,  1. ,  0. ],
       [ 0.5,  0. ,  1. ]])
>>> z # U矩阵
array([[ 2. ,  4. ,  1. ],
       [ 0. ,  1. ,  2. ],
       [ 0. ,  0. ,  2.5]])
>>> a = scipy.mat('[16 4 8; 4 5 -4; 8 -4 22]')
>>> a
matrix([[16,  4,  8],
        [ 4,  5, -4],
        [ 8, -4, 22]])
>>> linalg.cholesky(a) # 求a的Cholesky分解
array([[ 4.,  1.,  2.],
       [ 0.,  2., -3.],
       [ 0.,  0.,  3.]])  

SciPy中的linalg提供了基本的解线性方程组的函数，只要给出方程的系数矩阵和常数项向量，就可以获得方程组的解。对于非线性方程组可以采用数值解法。linalg模块中的solve函数可以用于解线性方程组，有如下所示的方程组

2X1 + X2 - 5X3 + X4 = 8
X1 - 3X2 - 6X4 = 9
2X2 - X3 + 2X4 = -5
X1 + 4X2 - 7X3 + 6X4 = 0

代码实现演示如下，在Python命令行中，使用linalg模块中的solve函数解上述线性方程组

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
>>> import scipy
>>> from scipy import linalg
>>> a = scipy.mat('[2 1 -5 1; 1 3 0 -6; 0 2 -1 2;1 4 -7 6]') # 建立系数矩阵
>>> a
matrix([[ 2,  1, -5,  1],
        [ 1,  3,  0, -6],
        [ 0,  2, -1,  2],
        [ 1,  4, -7,  6]])
>>> b = scipy.mat('[8;9;-5;0]') # 建立常数项向量
>>> b
matrix([[ 8],
        [ 9],
        [-5],
        [ 0]])
>>> linalg.solve(a, b)
array([[-3.26086957],
       [-1.56521739],
       [-3.7826087 ],
       [-2.82608696]])  

运行实例环境

1
2
# _*_ coding: utf-8 -*-
# version 2.7.13  