一、矩阵相乘

1.1numpy实现代码

python中矩阵相乘，即两个二维数组做矩阵乘法，可以使用numpy实现，也可以使用sympy实现，以numpy实现为例，numpy使用np.dot做矩阵乘法运算，代码如下：

import numpy as np

a01=np.random.randint(10,size=(4,3))
print(a01)
array([[6, 4, 9],
       [1, 4, 8],
       [8, 0, 0],
       [5, 6, 8]])
a02=np.random.randint(10,size=(3,4))
print(a02)
array([[5, 7, 7, 9],
       [6, 4, 6, 2],
       [5, 2, 8, 2]])

a03=np.dot(a01,a02)
print(a03)
array([[ 99,  76, 138,  80],
       [ 69,  39,  95,  33],
       [ 40,  56,  56,  72],
       [101,  75, 135,  73]])

1.2结果以分数表示

numpy默认的计算结果是以小数形式呈现，在一些特殊应用背景中需要将结果转换为分数形式，示例代码如下：

import numpy as np
from fractions import Fraction

a01=np.random.randint(10,size=(4,3)).astype(float)
a01[0,0]=0.5
print(a01)
array([[0.5, 4. , 3. ],
       [0. , 4. , 4. ],
       [4. , 3. , 6. ],
       [0. , 2. , 6. ]])
a02=np.random.randint(10,size=(3,4))
print(a02)
array([[5, 7, 7, 9],
       [6, 4, 6, 2],
       [5, 2, 8, 2]])

a03=np.dot(a01,a02)
print(a03)#a03为小数形式
array([[41.5, 25.5, 51.5, 18.5],
       [44. , 24. , 56. , 16. ],
       [68. , 52. , 94. , 54. ],
       [42. , 20. , 60. , 16. ]])

a04=np.vectorize(lambda x: Fraction.from_float(x).limit_denominator())(a03)
print(a04)#a04为分数形式
array([[Fraction(83, 2), Fraction(51, 2), Fraction(103, 2), Fraction(37, 2)],
       [Fraction(44, 1), Fraction(24, 1), Fraction(56, 1), Fraction(16, 1)],
       [Fraction(68, 1), Fraction(52, 1), Fraction(94, 1), Fraction(54, 1)],
       [Fraction(42, 1), Fraction(20, 1), Fraction(60, 1), Fraction(16, 1)]], dtype=object)

二、矩阵求逆

矩阵求逆操作即求矩阵的逆矩阵，逆矩阵与原矩阵相乘的结果为单位矩阵。

2.1numpy求解

numpy中矩阵求逆通过np.linalg.inv函数实现，下面给出求逆的小数实现与分数实现，代码如下：

import numpy as np
from fractions import Fraction

#numpy实现
a01=np.random.randint(5,size=(2,2))
print(a01)
array([[4, 0],
       [3, 4]])
a02 = np.linalg.inv(a01.astype(float))#求逆矩阵
print(a02)
array([[ 0.25  ,  0.    ],
       [-0.1875,  0.25  ]])

a03=np.vectorize(lambda x: Fraction.from_float(x).limit_denominator())(a02)
print(a03)#分数形式
array([[Fraction(1, 4), Fraction(0, 1)],
       [Fraction(-3, 16), Fraction(1, 4)]], dtype=object)

2.2sympy求解

sympy中矩阵求逆通过inv函数实现，从求解结果分析，sympy会自动给出结果的简洁形式，如果可表示为分数则用分数表示，在展示效果上优于numpy。代码如下：

from sympy import Matrix

#sympy实现
a04=Matrix([[4,0],[3,4]])
print(a04)
Matrix([
[4, 0],
[3, 4]])
print(a04.inv())
Matrix([
[  1/4,   0],
[-3/16, 1/4]])

三、矩阵特征值与特征向量求解

矩阵特征值求解，可用于PCA主成分分析等数据降维算法中，也可以用于图像处理中的图像压缩算法等。

3.1numpy求解

numpy求解矩阵特征值与特征向量，通过np.linalg.eig函数实现，代码如下：

import numpy as np

#numpy实现
a01=np.array([[-1,2,0],[0,3,0],[2,1,-1]])
print(a01)
array([[-1,  2,  0],
       [ 0,  3,  0],
       [ 2,  1, -1]])
eigenvalue, featurevector = np.linalg.eig(a01)
print(eigenvalue)
array([-1., -1.,  3.])
print(featurevector)
array([[ 0.00000000e+00,  1.11022302e-16,  4.08248290e-01],
       [ 0.00000000e+00,  0.00000000e+00,  8.16496581e-01],
       [ 1.00000000e+00, -1.00000000e+00,  4.08248290e-01]])

3.2sympy求解

sympy求解矩阵特征值与特征向量，通过eigenvals、eigenvects两个函数实现，从求解结果看，sympy的求解结果比numpy求解结果更为直观。代码如下：

from sympy import Matrix

#sympy实现
a02=Matrix([[-1,2,0],[0,3,0],[2,1,-1]])
print(a02)
Matrix([
[-1, 2,  0],
[ 0, 3,  0],
[ 2, 1, -1]])
print(a02.eigenvals())
{3: 1, -1: 2}#3: 1意为特征值为3，个数为1个，-1: 2意为特征值为-1，个数为2个
print(a02.eigenvects())
[(-1, 2, [Matrix([
[0],
[0],
[1]])]), (3, 1, [Matrix([
[1],
[2],
[1]])])]

四、矩阵约当标准型与转换矩阵求解

矩阵的约当标准型与转换矩阵求解，可用于矩阵幂的运算简化，大幅度提高计算速度，也可以用于微分方程的求解等等。

矩阵约当标准型与转换矩阵使用numpy求解较为繁琐，因此使用sympy求解，代码如下：

from sympy import Matrix

a01=Matrix([[-1,2,0],[0,3,0],[2,1,-1]])
print(a01)
Matrix([
[-1, 2,  0],
[ 0, 3,  0],
[ 2, 1, -1]])
p_mat,j_mat=a01.jordan_form()#p_mat为转换矩阵，j_mat为约当标准型
print(p_mat)
Matrix([
[0, 1, 1],
[0, 0, 2],
[2, 0, 1]])
print(j_mat)
Matrix([
[-1,  1, 0],
[ 0, -1, 0],
[ 0,  0, 3]])

注意sympy求解的约当标准型的1放在上三角，有些数学教材中将1放在下三角。

五、矩阵奇异值分解

矩阵奇异值分解的用途与矩阵特征值求解类似，可用于PCA主成分分析，图像压缩，也可用于语音识别与信号处理等领域。sympy不直接支持矩阵奇异值的求解，因此使用numpy求解，函数为np.linalg.svd，代码如下：

import numpy as np

a02=np.array([[4,0],[3,0],[0,0]])
print(a02)
array([[4, 0],
       [3, 0],
       [0, 0]])
U, S, Vt = np.linalg.svd(a02)
print(U)#左奇异向量组成的酉矩阵
array([[-0.8, -0.6,  0. ],
       [-0.6,  0.8,  0. ],
       [ 0. ,  0. ,  1. ]])
print(S)#奇异值
array([5., 0.])
print(Vt)#右奇异向量组成的酉矩阵
array([[-1., -0.],
       [ 0.,  1.]])

奇异向量组成的酉矩阵并不唯一。

六、矩阵方程组求解

矩阵方程组求解比较常用，使用sympy求解，求解过程为先定义系数矩阵、结果矩阵和变量symbols，然后使用linsolve函数求解，代码如下：

from sympy import Matrix, symbols, linsolve

A = Matrix([[2,3,1],[4,2,3], [7,1,-1]])
print(A)
Matrix([
[2, 3,  1],
[4, 2,  3],
[7, 1, -1]])
B = Matrix([[4],[17],[1]])
print(B)
Matrix([
[ 4],
[17],
[ 1]])
x, y, z = symbols('x y z')
result=linsolve((A,B),x,y,z)
print(result)
{(1, -1, 5)}

A = Matrix([[1,1,1,1],[4,3,5,-1], [2,1,3,-3]])
print(A)
Matrix([
[1, 1, 1,  1],
[4, 3, 5, -1],
[2, 1, 3, -3]])
B = Matrix([[-1],[-1],[1]])
print(B)
Matrix([
[-1],
[-1],
[ 1]])
x, y, z = symbols('x y z')
result=linsolve((A,B),x,y,z)
print(result)
{(2 - 2*z, z - 3, z)}

从示例代码可看出，sympy可以解参数矩阵非奇异时的方程组，也可以解参数矩阵奇异时的方程组，并且给出通解。