求矩阵的秩

crushed-sw

已于 2023-11-27 15:51:16 修改

阅读量3.5k

点赞数 6

文章标签：矩阵线性代数 c++

于 2022-05-30 23:04:14 首次发布

本文链接：https://blog.csdn.net/qq_50952332/article/details/125055896

版权

求矩阵的秩

文章目录

求矩阵的秩
- 矩阵秩的定义：
- 矩阵秩的求解：
- - 定义法：
  - 阶梯型矩阵非零行数：
  - 全选主元高斯消去求秩

矩阵秩的求法也有几种，比如说它的定义和简便求法，让我们来复习一下：

矩阵秩的定义：

在 $\times n$ 的矩阵 $A$ 中有一个不等于0的r阶子式 $D$ ，且所有的r+1阶子式全等于0，则 $D$ 是该矩阵的最高阶非零子式。非零子式的最高阶数即叫做矩阵的秩，记作 $R (A)$ ，r是rank的缩写。

矩阵秩的求解：

定义法：

该方法是根据矩阵的秩的定义来求，如果找到k阶子式为0，而k-1阶子式不为0，那么k-1即该矩阵的秩，来看下面例子：
$\begin{pmatrix} 1 & 2 & 3 & 0\\ 0 & 1 & 0 & 1\\ 0 & 0 & 1 & 0\\ \end{pmatrix}$
针对矩阵 $A$ ，我们首先找到它的一个三阶子式是否为0，如我们可以找到：
$\begin{pmatrix} 1 & 2 & 3 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \\ \end{pmatrix}$
很显然，该子式等于-1不等于0，所以该矩阵的秩为3

阶梯型矩阵非零行数：

该方法便是我们求解矩阵的秩的时候常用的方法，用该方法求解可以分为两步：

1、现将原先的矩阵简化为阶梯型矩阵

2、数出新矩阵的非零行行数，即为该原矩阵的秩

例如下面这个矩阵：
$\begin{pmatrix} 1 & 2 & 4 & 1\\ 2 & 4 & 8 & 2\\ 3 & 6 & 2 & 0\\ \end{pmatrix}$
我们可以对其进行行变换：
$\begin{aligned} A&= \begin{pmatrix} 1 & 2 & 4 & 1\\ 2 & 4 & 8 & 2\\ 3 & 6 & 2 & 0\\ \end{pmatrix}\\\\ &\Rightarrow \begin{pmatrix} 1 & 2 & 4 & 1\\ 0 & 0 & 0 & 0\\ 0 & 0 & -10 & -3\\ \end{pmatrix}\\\\ &\Rightarrow \begin{pmatrix} 1 & 2 & 4 & 1\\ 0 & 0 & -10 & -3\\ 0 & 0 & 0 & 0\\ \end{pmatrix}\\ \end{aligned}$
因此可以得知该矩阵的秩为2。

全选主元高斯消去求秩

用以上的方法可能会有个小问题，就是消元之后会碰到该行元素都是0的情况，因此为了消除该不确定因素，确保该算法的稳定，我们使用全选主元高斯消去法来进行消元求秩，全选主元的过程如下所示：

对于矩阵消元过程中的第 $k$ 步，首先，在 $a_{kk}$ 的右下方（包括 $a_{kk}$ 自己）的 $n - k + 1$ 阶子矩阵中选取绝对值最大的元素，并将该元素所在的行号记录在 $I S (k)$ 中，而列号记录在 $J S (k)$ 中，然后通过行交换与列交换将该绝对值最大的元素交换到主对角线 $a_{kk}$ 的位置上，即做以下两步：
$\Leftrightarrow A(IS(k),l),l=1,2,...,n\\ A(l,k) \Leftrightarrow A(l,JS(k)),l=1,2,...,n\\$
经过全选主元过后，矩阵消元的过程是稳定的。但最后需要对结果进行恢复。恢复过程如下：
$\Leftrightarrow A(JS(k),l),l=1,2,...,n\\ A(l,k) \Leftrightarrow A(l,IS(k)),l=1,2,...,n\\$

还是以上那个例子：

$\begin{pmatrix} 1 & 2 & 4 & 1\\ 2 & 4 & 8 & 2\\ 3 & 6 & 2 & 0\\ \end{pmatrix}$
对 $a_{00}$ 进行全选主元，发现右下角子矩阵中最大的是8，因此保存行号与列号 $I S (0) = 2$ ， $J S (0) = 3$ 以便恢复，此时交换两元素：
$\begin{aligned} A&\stackrel{行变换后}{\Longrightarrow} \begin{pmatrix} 2 & 4 & 8 & 2\\ 1 & 2 & 4 & 1\\ 3 & 6 & 2 & 0\\ \end{pmatrix}\\ &\stackrel{列变换后}{\Longrightarrow} \begin{pmatrix} 8 & 4 & 2 & 2\\ 4 & 2 & 1 & 1\\ 2 & 6 & 3 & 0\\ \end{pmatrix} \end{aligned}$
此时再进行消元操作：
$\begin{aligned} A&\stackrel{归一化}{\Longrightarrow} \begin{pmatrix} 1 & \frac{1}{2} & \frac{1}{4} & \frac{1}{4}\\ 4 & 2 & 1 & 1\\ 2 & 6 & 3 & 0\\ \end{pmatrix}\\ &\stackrel{消元}{\Longrightarrow} \begin{pmatrix} 1 & \frac{1}{2} & \frac{1}{4} & \frac{1}{4}\\ 0 & 0 & 0 & 0\\ 0 & 5 & \frac{5}{2} & -\frac{1}{2}\\ \end{pmatrix} \end{aligned}$
因为此时的消元只为了求秩，因此不复原也可，若是为了求逆或求解矩阵，则必须还原。

此时再对 $a_{11}$ 全选主元，最大值为5，同样进行交换：
$\begin{aligned} A&\stackrel{行变换后}{\Longrightarrow} \begin{pmatrix} 1 & \frac{1}{2} & \frac{1}{4} & \frac{1}{4}\\ 0 & 5 & \frac{5}{2} & -\frac{1}{2}\\ 0 & 0 & 0 & 0\\ \end{pmatrix}\\ \end{aligned}$
此时再进行归一和消元：
$\begin{aligned} A&\stackrel{归一化}{\Longrightarrow} \begin{pmatrix} 1 & \frac{1}{2} & \frac{1}{4} & \frac{1}{4}\\ 0 & 1 & \frac{1}{2} & -\frac{1}{10}\\ 0 & 0 & 0 & 0\\ \end{pmatrix}\\ \end{aligned}$
此时对 $a_{22}$ 全选主元，最大值为0，此时便可知道，该矩阵的秩为2。

总结

来总结一下关于全选主元的方法：

设有 $m\times n$ 阶矩阵：
$\begin{pmatrix} a_{00} & a_{01} & \dots & a_{0, n-1}\\ a_{10} & a_{11} & \dots & a_{1, n-1}\\ \vdots & \vdots & \ddots & \vdots\\ a_{m-1,0} & a_{m-1,1} & \dots & a_{m-1, n-1}\\ \end{pmatrix}$
取 $k=min\{m,n\}$ 。对于 $r=0,1,\dots,k-1$ ，用全选主元高斯消去法将A变为上三角矩阵，知道某次 $a_{rr}=0$ 为止，矩阵 $A$ 的秩为 $r$ 。

对该算法进行接口设计 int rank(double* a, int m, int n)

接口与参数说明：

形参与函数类型	参数意义
`double* a`	存放矩阵 $A$ 的元素，返回时存放 $A$ 的逆矩阵
`int m`	矩阵A的行数
`int n`	矩阵A的列数
`int rank()`	返回A的秩

此时就可以轻松写出代码：

#include <cmath>

template <int M, int N>
int rank(double a[M][N]) {
    int i, j, k, l;
    int is, js, nn;
    double q, d;

    // 找出M和N最小值
    nn = M > N ? N : M;
    k = 0;

    for(l = 0; l != nn; ++l) {
        // 循环找出主元
        q = 0;
        for(i = l; i != M; ++i) {
            for(j = l; j != N; ++j) {
                d = fabs(a[i][j]);
                if(d > q) {
                    q = d;
                    is = i;
                    js = j;
                }
            }
        }

        // 主元为0，返回答案
        if(q + 1.0 == 1.0) {
            return k;
        }

        k++;
        // 主元行交换
        if(is != i) {
            for(j = l; j != N; ++j) {
                d = a[l][j];
                a[l][j] = a[is][j];
                a[is][j] = d;
            }
        }

        // 主元列交换
        if(js != j) {
            for(i = l; i != M; ++i) {
                d = a[i][js];
                a[i][js] = a[i][l];
                a[i][l] = d;
            }
        }

        // 对其他行数进行消元
        for(i = l+1; i != ; ++i) {
            d = a[i][l] / a[l][l];
            for(j = l+1; j != N; ++j) {
                a[i][j] = a[i][j] - d * a[l][j];
            }
        }
    }

    return k;
}

我们运用以上例子来进行检验：

#include "rank.hpp"

int main() {
    double a[3][4] = {{1, 2, 4, 1},
                      {2, 4, 8, 2},
                      {3, 6, 2, 0}};

    int r = rank<3, 4>(a);
    std::cout << r << std::endl;

    return 0;
}

程序输出如下：