線性方程組求解（高斯、高斯若爾當）+應用（求逆，求行列式）

高斯消去法

#include <iostream>
#include <vector>
#include <math.h>
using namespace std;
typedef vector<vector<double>> Matrix;
void dispMatrix(Matrix m);      //輸出矩陣值
void dispRes(vector<double> r);        //輸出解向量
int MajorRow(Matrix &m, int i); //選取列主元
double det(Matrix m);
void normlize(Matrix &m, int i); //行列式第i行歸一
vector<double> Gauss(Matrix &m);
int main()
{
    ios::sync_with_stdio("false");
    int n;
    cin >> n;
    Matrix M1(n, vector<double>(n + 1, 0));
    for (int i = 0; i < n; i++)
        for (int j = 0; j < n + 1; j++) //這裡需要注意是輸入n*n 還是n*n+1 根據需要改寫
            cin >> M1[i][j];

    cout << det(M1) << endl;
    vector<double> res = Gauss(M1);
    dispRes(res);
    system("pause");
    return 0;
}
/* 
3
2 -1 3 1
4 2 5 4
1 2 0 7*/
void dispMatrix(Matrix m)
{
    cout << "Matrix:\n";
    for (auto i : m)
    {
        for (auto ele : i)
            cout << ele << " ";
        cout << endl;
    }
}
void dispRes(vector<double> r)
{
    cout << "res:\n";
    for (auto ele : r)
        cout << ele << " ";
    cout << endl;
}
int MajorRow(Matrix &m, int i)
{
    int n = (int)m.size();
    double max = 0;
    int j, row;                   //第i列的最大值在第j行
    for (row = i; row < n; row++) //從第i行往後找列主元！
        if (fabs(m[row][i]) > max)
        {
            max = fabs(m[row][i]);
            j = row;
        }
    if (j != i)
        return j;
    else
        return -1; //是最大的值
}
double det(Matrix m)
{
    int n = (int)m.size();
    double res = 1.0;
    for (int i = 0; i < n - 1; i++)
    {
        int p = MajorRow(m, i);
        if (p != -1)
            res *= -1; //交換行 deta需要乘-1
        for (int j = i + 1; j < n; j++)
        {
            double mi = m[j][i] / m[i][i];
            for (int k = i; k < n; k++) //注意 這裡是n 因為是行列式
                m[j][k] -= mi * m[i][k];
        }
        res *= m[i][i];
    }
    res *= m[n - 1][n - 1];
    return res;
}
void normlize(Matrix &m, int i)
{
    for (auto ele : m[i])
        ele /= m[i][i];
}
vector<double> Gauss(Matrix &m)
{
    int n = (int)m.size();
    vector<double> res(n);
    for (int i = 0; i < n - 1; i++)
    {
        int p = MajorRow(m, i);
        if (p != -1) //i不是列主元
            swap(m[p], m[i]);
        for (int j = i + 1; j < n; j++)
        {
            double mi = m[j][i] / m[i][i];
            for (int k = i; k < n + 1; k++)
                m[j][k] -= mi * m[i][k];
        }
    }
    // dispMatrix(m, n);
    // 回代
    for (int i = n - 1; i >= 0; i--)
    {
        res[i] = m[i][n];
        for (int j = n - 1; j > i; j--)
            res[i] -= m[i][j] * res[j];
        res[i] /= m[i][i];
    }
    return res;
}
           

高斯若爾當消去法

#include <iostream>
#include <vector>
#include <math.h>
using namespace std;
typedef vector<vector<double>> Matrix;

void dispMatrix(Matrix m);                         //輸出矩陣值
void dispRes(vector<double> r);                    //輸出解向量
int MajorRow(Matrix &m, int i);                    //選取列主元
void normlize(Matrix &m, int i);                   //行列式第i行歸一
vector<double> GaussJordan_elimination(Matrix &m); //高斯若爾當消去主函數 傳回解向量
vector<vector<double>> inverse(Matrix &m);         //求逆
Matrix eye(int n);
/* 
todo:
3
2 2 3 3
4 7 7 1
-2 4 5 -7
3
1 2 3
2 1 2
1 3 4
3
2 -1 3 1
4 2 5 4
1 2 0 7
*/
int main()
{
    ios::sync_with_stdio("false");

    int n;
    cin >> n;
    Matrix M1(n, vector<double>(n + 1, 0));
    for (int i = 0; i < n; i++)
        for (int j = 0; j < n + 1; j++) //這裡需要注意是輸入n*n 還是n*n+1 根據需要改寫
            cin >> M1[i][j];

    system("pause");
    return 0;
}
void dispMatrix(Matrix m)
{
    cout << "Matrix:\n";
    for (auto i : m)
    {
        for (auto ele : i)
            cout << ele << " ";
        cout << endl;
    }
}
void dispRes(vector<double> r)
{
    cout << "res:\n";
    for (auto ele : r)
        cout << ele << " ";
    cout << endl;
}
int MajorRow(Matrix &m, int i)
{
    int n = (int)m.size();
    double max = 0;
    int j, row;                   //第i列的最大值在第j行
    for (row = i; row < n; row++) //從第i行往後找列主元！
        if (fabs(m[row][i]) > max)
        {
            max = fabs(m[row][i]);
            j = row;
        }
    if (j != i)
        return j;
    else
        return -1; //是最大的值
}
void normlize(Matrix &m, int k)
{
    int n = (int)m.size();
    if (m[k][k] == 1)
        return; //不需要歸一
    for (int i = k + 1; i < n + 1; i++)
        m[k][i] /= m[k][k];
    m[k][k] = 1;
}
// 高斯若爾當消去
vector<double> GaussJordan_elimination(Matrix &m)
{
    int n = (int)m.size();
    for (int i = 0; i < n; i++)
    {
        int p = MajorRow(m, i);
        if (p != -1) //i不是列主元
            swap(m[p], m[i]);
        normlize(m, i);
        for (int j = 0; j < n; j++)
        {
            if (i == j)
                continue;
            double mi = m[j][i];
            for (int k = i; k < n + 1; k++)
                m[j][k] -= mi * m[i][k];
        }
    }
    dispMatrix(m);
    vector<double> res;
    for (auto ele : m)
        res.push_back(ele[n]);
    return res;
}
Matrix eye(int n)
{
    Matrix m(n, vector<double>(n, 0));
    for (int i = 0; i < n; i++)
        m[i][i] = 1;
    return m;
}
//矩陣求逆
Matrix inverse(Matrix &m)
{
    int n = (int)m.size();
    Matrix res = eye(n);
    for (int i = 0; i < n; i++)
    {
        // 列選主元
        int p = MajorRow(m, i);
        if (p != -1) //i不是列主元
        {
            swap(m[p], m[i]);
            swap(res[p], res[i]);
        }
        // 歸一
        if (m[i][i] != 1)
        {
            for (int p = 0; p < n; p++)
                res[i][p] /= m[i][i];
            normlize(m, i);
        }
        //
        for (int j = 0; j < n; j++)
        {
            if (i == j)
                continue;
            double mi = m[j][i];
            for (int k = i; k < n; k++)
                m[j][k] -= mi * m[i][k];
            for (int k = 0; k < n; k++)
                res[j][k] -= mi * res[i][k];
        }
    }
    return res;
}