【线性代数 C++】求逆矩阵

news2024/11/18 8:11:48

对于 $n$ 阶矩阵 $A$ ，如果有 $n$ 阶矩阵 $B$ ，使 $A B = B A = E$ ，则说 $A$ 是可逆的，并把 $B$ 称为 $A$ 的逆矩阵.
$A$ 的逆矩阵记作 $A^{-1}$ ，则 $B=A^{-1}$ .
若 $\begin{vmatrix}A\end{vmatrix} \neq 0$ ，则 $A$ 可逆，且 $A^{-1}= \frac{1}{\begin{vmatrix}A \end{vmatrix}}A^*$ .
上式中， $\begin{vmatrix}A\end{vmatrix}$ 为 $A$ 的行列式， $A^*$ 为 $A$ 的伴随矩阵.

根据上述内容，为求逆矩阵，需要分别求得矩阵行列式的值和矩阵的逆矩阵.

1 矩阵行列式求值

1.1 定义

由 $n$ 阶方阵 $A$ 的元素构成的行列式（各元素位置不变），称为 $A$ 的行列式，记作 $\begin{vmatrix}A \end{vmatrix}$ 或 $d e t A$ .

1.2 C++代码

为便于管理函数，建立行列式类CDeterminant.

求行列式值的两种途径：一是根据行列式定义、二是利用代数余子式.
在CDeterminant类中增加两个成员函数GetCetValByDef()和GetDetValByRem()，相关代码参考【线性代数|行列式定义及其值】和【线性代数 | C++】行列式按行（列）展开。

2 求伴随矩阵

2.1 定义

行列式 $\begin{vmatrix}A\end{vmatrix}$ 各元素的代数余子式 $A_{ij}$ 构成的如下矩阵 $A^*=\begin{bmatrix}A_{11} & A_{21} & \cdots & A_{n1} \\A_{12} & A_{22} & \cdots & A_{n2} \\ \vdots & \vdots && \vdots \\A_{1n} & A_{2n} &\cdots & A_{nn}\end{bmatrix},$ 称为矩阵 $A$ 的伴随矩阵.
注意：伴随矩阵 $A^*$ 中元素的位置与对应行列式 $\begin{vmatrix}A\end{vmatrix}$ 中元素的位置呈转置关系.

2.2 C++代码

2.2.1 求余子式

由定义可知，求伴随矩阵，需先求行列式各元素的代数余子式.
求代数余子式，需先求余子式.
故，在CDeterminant类中增加一个求行列式余子式的成员函数GetDetRem().

//在CDeterminant.h声明成员函数
static bool GetDetRem
(
    const vector<vector<double>> &vvDetInput, //原始行列式
    int i,  //待求余子式元素的行号
    int j,  //待求余子式元素的列号
    vector<vector<double>> &vvDetRet    //求得的余子式
);

//在CDeterminant.cpp中定义成员函数
bool CDeterminant::GetDetRem
(
    const vector<vector<double>> &vvDetInput, 
    int i, 
    int j,
    vector<vector<double>> &vvDetRet
)
{
    if (false == IsDet(vvDetInput))//形参是否符合行列式格式要求
        return false;

    vvDetRet.clear();
    vvDetRet = vvDetInput;

    //删除元素所在的行
    vvDetRet.erase(vvDetRet.cbegin() + i);
    for (int i = 0; i < vvDetRet.size(); i++)
    {
        //删除元素所在的列
        vvDetRet[i].erase(vvDetRet[i].cbegin() + j);
    }
    
    return true;
}

2.2.2 求伴随矩阵

在CMatrix类中添加GetAdjointMat函数，用于求伴随矩阵.其思路是：
1. 按列逐行取得矩阵的元素；
2. 求相应的余子式
3. 对余子式求值
4. 求代数余子式
5. 将代数余子式的值，按行填入新矩阵
6. 得到伴随矩阵

//在CMatrix.h中声明成员函数
static bool GetAdjointMat
(
const vector<vector<double>> &vvMatInput,
vector<vector<double>> &vvMatRet
);

//在CMatrix.cpp中定义成员函数
bool CMatrix::GetAdjointMat
(
  const vector<vector<double>> &vvMatInput,
  vector<vector<double>> &vvMatRet
)
{
  //判断vector变量是否符合行列式格式
  if (false == CDeterminant::IsDet(vvMatInput))
    return false;
  
  vvMatRet.clear();
  vector<double> vTemp;//临时存储伴随矩阵的行元素
  vector<vector<double>> vvDetTemp;//临时存储余子式
  double iDetValTemp;//临时存储余子式的值
  //按列循环
  for (int i = 0; i < vvMatInput[0].size(); i++)
  {
    vTemp.clear();
    //按行循环
    for (int j = 0; j < vvMatInput.size(); j++)
    {
      //求余子式
      vvDetTemp.clear();
      CDeterminant::GetDetRem(vvMatInput, j, i, vvDetTemp);
      //余子式求值
      iDetValTemp = 0;
      CDeterminant::GetDetValByDef(vvDetTemp, iDetValTemp);
      //求代数余子式
      iDetValTemp = (pow(-1, i + j) * iDetValTemp);
      //代数余子式的值填入新矩阵的行元素中
      vTemp.push_back(iDetValTemp);
    }  
    //行元素填入新矩阵中
    vvMatRet.push_back(vTemp);  
  }  
  
  return true;
}

3 求逆矩阵

利用前面求得的 $\begin{vmatrix}A \end{vmatrix}$ 及 $A^*$ ，按照公式 $A^{-1}= \frac{1}{\begin{vmatrix}A \end{vmatrix}}A^*$ ，可以求得 $A$ 的逆矩阵.
在CMatrix类中增加GetInverseMat函数，以实现上述功能.

//在CMatrix.h中声明函数
static bool GetInverseMat
(
const vector<vector<double>> &vvMatInput,
vector<vector<double>> &vvInverseMat
);

//在CMatrix.cpp中定义函数
bool CMatrix::GetInverseMat
(
  const vector<vector<double>> &vvMatInput,
  vector<vector<double>> &vvInverseMat
)
{
  if (false == CDeterminant::IsDet(vvMatInput))
    return false;

  //方阵行列式值不等于0时，方阵可逆
  double fDetVal;
  CDeterminant::GetDetValByDef(vvMatInput, fDetVal);
  if (0 == fDetVal)
    return false;  
  fDetVal = 1.0 / fDetVal;

  vector<vector<double>> vvMatTemp;
  GetAdjointMat(vvMatInput, vvMatTemp);//求伴随矩阵

  vvInverseMat.clear();
  MatMulti(fDetVal, vvMatTemp, vvInverseMat);//求数与矩阵相乘

  return true;
}

4 测试

//在test.cpp中测试
#include <iostream>
#include <iomanip>
#include <vector>

#include "CMatrix.h"

using namespace std;

bool PrintMat
(
  const vector<vector<double>> &vvMat
)
{
  for (int i = 0; i < vvMat.size(); i++)
  {
    for (int j = 0; j < vvMat[i].size(); j++)
    {
      cout << setw(5) << vvMat[i][j];
    }
    cout << endl;
  }
  
  return true;
}

int main()
{
  vector<vector<double>> vvMatA{{ 1, 2, 3},
                                { 2, 2, 1},
                                { 3, 4, 3}};
  
  vector<vector<double>> vvMatRet;

  if (false == CMatrix::GetInverseMat(vvMatA, vvMatRet))
  {
    cout << "计算失败" << endl;
  }
  else
  {
    PrintMat(vvMatRet);
  } 

  return 0;
}