矩阵求逆

1. 矩阵求逆原理介绍

矩阵求逆的原理有很多,此处仅介绍两种:利用伴随矩阵 和 利用行变换
前者具有较好的精度,后者具有较好的计算速率

1.1 利用伴随矩阵求逆

(1) 代数余子式

一个 \(n \times n\) 矩阵 \(A\)\(A\)\((i,j)\) 处的代数余子式为 \(A\) 去掉第 \(i\) 行和第 \(j\) 列后剩下的矩阵的行列式的值再乘上 \((-1)^{i+j}\)

(2) 伴随矩阵

一个 \(n \times n\) 矩阵 \(A\)\(A\) 的伴随矩阵 \(C_A\) 为一个 \(n \times n\) 矩阵,\(C_A\)\((i,j)\) 处的值为 \(A\)\((i,j)\) 处的代数余子式

(3) 转置矩阵

一个 \(n \times n\) 矩阵 \(A\)\(A\) 的转置矩阵 \(A^T\) 为一个 \(n \times n\) 矩阵,\(A^T\)\((i,j)\) 处的值为 \(A\)\((j,i)\) 处的值

(4) 逆矩阵

一个 \(n \times n\) 矩阵 \(A\)\(A\) 的逆矩阵 \(A^{-1}\)\(\frac{1}{d} \cdot (C_A)^T\)。其中 d 是矩阵 \(A\) 的行列式,\(C_A\)\(A\) 的伴随矩阵,\((C_A)^T\)\(C_A\) 的转置矩阵

1.2 利用高斯行变换

一个 \(n \times n\) 的矩阵 \(A\),构造一个 \(n \times n\) 的单位阵作为其扩展矩阵,拼接后成为一个 \(n \times 2n\) 的矩阵 \([A_{n \times n}|I_{n \times n}]\)

再通过初等行变换,将上述矩阵转变为 \([I_{n \times n}|P_{n \times n}]\) 的形式,则有 \(A^{-1} = P\)

2. 矩阵行列式

2.1 使用余子式

对于 \(n \times n\) 的矩阵 \(A\),其行列式 \(det(A)\) 的值等于 任意一行的元素 依次乘以 对应位置的余子式 的结果之和

2.2 利用行变换

通过初等行变换将矩阵变成对角矩阵(除了对角线之外其它元素均为 0 的矩阵),对角线元素的乘积即为行列式的值

3. 矩阵求逆代码(python)

两种求逆方法各有优劣,利用伴随矩阵的方法在对整数矩阵求逆时精度更好,会损失精度的除法运算只在最后一步中使用;利用高斯行变换求解时可能会因为行变换的除法而损失一定的精度

但就时间复杂度来说,伴随矩阵所消耗的时间复杂度远大于使用高斯行变换的方式

3.1 利用伴随矩阵

def _determinant(matrix):
    """递归计算矩阵行列式"""
    if len(matrix) == 1:
        return matrix[0][0]
    else:
        d = 0
        matrix_row, matrix_col = len(matrix), len(matrix[0])
        for col in range(matrix_col):
            # 删去第 0 行 col 列
            m = [[matrix[i][j] for j in range(matrix_col) if j != col] for i in range(matrix_row) if i != 0]
            d = d + (-1) ** (col & 1) * matrix[0][col] * _determinant(m)
        return d


def _cofactor(matrix):
    """计算伴随矩阵"""
    matrix_row, matrix_col = len(matrix), len(matrix[0])
    res = [[0 for _ in range(matrix_col)] for _ in range(matrix_row)]

    for row in range(matrix_row):
        for col in range(matrix_col):
            # 删去第 i 行 j 列
            m = [[matrix[i][j] for j in range(matrix_col) if j != col] for i in range(matrix_row) if i != row]
            # 计算余子式
            res[row][col] = (-1) ** ((row + col) & 1) * _determinant(m)
    return res


def _traverse(matrix):
    """矩阵转置"""
    return [[matrix[col][row] for col in range(len(matrix[0]))] for row in range(len(matrix))]


def inverse(matrix):
    d = _determinant(matrix)  # 行列式
    if d == 0:  # 行列式等于0,不存在逆矩阵
        return None
    else:
        cofactor = _cofactor(matrix)  # 获取伴随矩阵
        res = _traverse(cofactor)  # 伴随矩阵的转置
        # 将每个元素都除以 d
        for row in range(len(res)):
            for col in range(len(res[0])):
                res[row][col] = res[row][col] / d
        return res

3.2 利用高斯行变换

def inverse(matrix):
    row, col = len(matrix), len(matrix[0])  # 矩阵的行列
    t_matrix = [[matrix[r][c] for c in range(col)] for r in range(row)]
    e_matrix = [[0 if c != r else 1 for c in range(col)] for r in range(row)]  # 扩展矩阵

    for i in range(row):
        # 寻找第i列不为0的行
        for r in range(i, row):
            if t_matrix[r][i] != 0:
                if i != r:
                    t_matrix[i], t_matrix[r] = t_matrix[r], t_matrix[i]
                    e_matrix[i], e_matrix[r] = e_matrix[r], e_matrix[i]
                break
        else:  # 找不到对应的行,没有逆矩阵
            return None

        # 对当前行的变换
        temp = t_matrix[i][i]
        for c in range(col):
            t_matrix[i][c] /= temp
            e_matrix[i][c] /= temp
        # 对其它行的变换
        for r in range(row):
            if r != i:
                temp = t_matrix[r][i]
                for c in range(col):
                    e_matrix[r][c] = e_matrix[r][c] - e_matrix[i][c] * temp
                    t_matrix[r][c] = t_matrix[r][c] - t_matrix[i][c] * temp
    return e_matrix

4. 最后

上文提过,利用伴随矩阵求逆的方式会消耗极大的时间,在求解 6 阶以上的矩阵的时候就已经开始明显变慢;而行变换的方式会产生精度问题,在求解高阶矩阵的逆矩阵是很可能会由于精度产生错误

例如求解如下矩阵的逆矩阵:

m = [[7, 18, 16, 9, 13],
     [21, 0, 21, 1, 22],
     [1, 20, 11, 2, 17],
     [1, 0, 14, 20, 5],
     [10, 13, 11, 22, 5]]

这是一个 5 阶矩阵,当使用行变换求解时就会由于除法精度问题而产生错误

为了避免精度的问题,可以采用 python 的分数库 Fraction 来计算除法,或者使用其它的求逆方法(例如 QR 分解)

posted @ 2021-07-15 11:48  kentle  阅读(1380)  评论(2编辑  收藏  举报