高斯消元解线性方程组

最新推荐文章于 2025-01-06 06:30:00 发布

W⁡angduoyu

最新推荐文章于 2025-01-06 06:30:00 发布

阅读量181

点赞数

分类专栏：数学学习总结文章标签：线性代数矩阵算法

本文链接：https://blog.youkuaiyun.com/ecjtu2020/article/details/121051797

版权

学习总结同时被 2 个专栏收录

39 篇文章

订阅专栏

数学

7 篇文章

订阅专栏

本文详细介绍了高斯消元法，用于解决包含多个方程和未知数的线性方程组。通过构建增广矩阵并进行初等行变换，将方程组转化为阶梯型矩阵，进而求解。此外，还展示了如何将非线性方程组通过巧妙转化，转化为线性方程组，应用高斯消元法求解。最后，给出了具体的浮点型和布尔型模板代码实现，并提供了一个实例来说明如何运用这种方法解决球面上点的坐标问题。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

高斯消元： $O(n^3)$

写出增广矩阵，将增广矩阵通过初等行变换变成阶梯型，然后往回带
适用于求解包含 $n$ 个方程， $n$ 个未知数的多元线性方程组
例如下式多元方程组：
$\begin{cases} a_{1,1}*x_1 + a_{1,2}*x_i +\quad \cdots \quad+ a_{1, n}*x_n = b_1\\ a_{2,1}*x_1 + a_{2,2}*x_i +\quad \cdots \quad+ a_{2, n}*x_n = b_2 \\ \qquad \qquad \vdots \quad \\ a_{n,1}*x_1 + a_{n,2}*x_i +\quad \cdots \quad+ a_{n, n}*x_n = b_n\end{cases}$
写出其增广矩阵：
$\begin{bmatrix} a_{1,1}&a_{1,2} \quad \cdots \quad a_{1, n} \quad b_1\\a_{2,1}&a_{2,2} \quad \cdots \quad a_{2, n} \quad b_2 \\ \qquad \qquad \vdots \quad \\ a_{n,1}&a_{n,2} \quad \cdots \quad a_{n, n} \quad b_n \end{bmatrix}\quad$

通过初等行变换，变成阶梯型矩阵（取 $n = 4$ ）为例：
$\begin{bmatrix} a_{1,1}&a_{1,2}& a_{1,3} &a_{1, 4} & b_1\\ 0&a_{2,2}& a_{2,3} &a_{2, 4} & b_2 \\ 0&0& a_{3,3} &a_{3, 4} & b_3 \\ 0&0&0 &a_{4, 4} & b_4 \end{bmatrix}\quad$

高斯消元模板：

int gauss() // 浮点值
{
    int c, r;
    for(c = 0, r = 0; c < n; ++ c)
    {
        int t = r; // 找到最大行
        for(int i = r + 1; i < n; ++i)
        {
            if(fabs(a[i][c]) > fabs(a[t][c]))
                t = i;
        }
        if(fabs(a[t][c]) < eps) continue;
        for(int i = c; i <= n; ++i) swap(a[t][i], a[r][i]); // 把最大行换上去
        for(int i = n; i >= c; -- i)  a[r][i] /= a[r][c]; // 把这一行第一个数变为1
        for(int i = r + 1; i < n; ++i)
        {
            if(fabs(a[i][c]) > eps) // 不是0
            {
                for(int j = n; j >= c; --j)
                {
                    a[i][j] -= a[r][j] * a[i][c];
                }
            }
        }
        r ++;
    }
    if(r < n) {
    // 当系数全为0，常数不为0时，无解
        for(int i = r; i < n; ++i)
            if(fabs(a[i][n]) > eps) return 2;
        return 1;
    }
    
    for(int i = n - 1; i >= 0; -- i)
    {
        for(int j = i + 1; j < n; ++ j)
        {
            a[i][n] -= a[j][n] * a[i][j];
        }
    }
    return 0;
}
// 布尔值
// 异或性质：
// a ^ b ^ c = x, e ^ f = y.
// a ^ b ^ c ^ e ^ f = x ^ y
int gauss() 
{
    int c, r;
    for (c = 0, r = 0; c < n; c ++ )
    {
        int t = r;
        for (int i = r; i < n; i ++ )  // 找非零行
            if (a[i][c])
                t = i;

        if (!a[t][c]) continue;

        for (int i = c; i <= n; i ++ ) swap(a[r][i], a[t][i]);  // 将非零行换到最顶端
        for (int i = r + 1; i < n; i ++ )  // 用当前行将下面所有的列消成0
            if (a[i][c])
                for (int j = n; j >= c; j -- )
                    a[i][j] ^= a[r][j];

        r ++ ;
    }

    if (r < n)
    {
        for (int i = r; i < n; i ++ )
            if (a[i][n])
                return 2;  // 无解
        return 1;  // 有多组解
    }

    for (int i = n - 1; i >= 0; i -- )
        for (int j = i + 1; j < n; j ++ )
            a[i][n] ^= a[i][j] * a[j][n];

    return 0;  // 有唯一解
}

题意：

在这里插入图片描述

Sol：

球体上所有点到球心的距离是相等的，因次只需求出一个点 $(x_1, x_2, x_3 \cdots, x_n)$ ，使得：
$\sum_{j=1}^n(a_{i,j} - x_j)^2 = C$
其中： $C$ 是常数， $\in [1, n+1]$ ，球面上第 $i$ 个点的坐标为 $(a_{i,1}, a_{i, 2} \cdots, a_{i, n})$ 。
该方程组由 $n + 1$ 个 $n$ 元二次方程构成，不是线性方程组。但是，我们可以两两相减，变成 $n$ 个 $n$ 元一次方程组，同时消去了常数 $C$ 。
$\sum_{j=1}^n (a_{i,j}^2 - a_{i+1, j}^2 - 2x_j(a_{i,j} - a_{i+1,j})) = 0 \qquad i \in[1, n]$
将变量放一边，常量放一边有：
$\sum_{j=1}^n 2(a_{i,j} - a_{i+1,j})x_j = \sum_{j=1}^n (a_{i,j}^2 - a_{i+1,j}^2) \qquad i \in [1, n]$
写成增广矩阵：
$\begin{bmatrix} 2(a_{1,1} - a_{2,1}) & 2(a_{1,2} - a_{2,2}) \quad \cdots \quad 2(a_{1,n} - a_{2,n}) \quad \sum_{j=1}^n (a_{1,j}^2 - a_{2,j}^2)\\ 2(a_{2,1} - a_{3,1}) & 2(a_{2,2} - a_{3,2}) \quad \cdots \quad 2(a_{2,n} - a_{3,n}) \quad \sum_{j=1}^n (a_{2,j}^2 - a_{3,j}^2)\\ \qquad \qquad \vdots \quad \\ 2(a_{n,1} - a_{n+1,1}) & 2(a_{n,2} - a_{n+1,2}) \quad \cdots \quad 2(a_{n,n} - a_{n+1,n}) \quad \sum_{j=1}^n (a_{n,j}^2 - a_{n+1,j}^2) \end{bmatrix}\quad$
利用高斯消元求解

Code：

#include <bits/stdc++.h>
using namespace std;
const int N = 20;
const double eps = 1e-8;
int n;
double a[N][N], d[N][N];


int gauss()
{
    int c, r;
    for(r = 0, c = 0; c < n; ++ c)
    {
        int t = r;
        for(int i = r + 1; i < n; ++i) 
            if(fabs(a[i][c]) > fabs(a[t][c]))
                t = i;
        if(fabs(a[t][c]) < eps) continue;
        for(int i = c; i <= n; ++i) swap(a[t][i], a[r][i]);
        for(int i = n; i >= c; -- i) a[r][i] /= a[r][c];
        for(int i = r + 1; i < n; ++i)
        {
            if(fabs(a[i][c]) > eps)
            {
                for(int j = n; j >= c; -- j)
                {
                    a[i][j] -= a[i][c] * a[r][j];
                }
            }
        }
        r ++;
    }
    if(r < n)
    {
        for(int i = r; i < n; ++i)
            if(fabs(a[i][n]) > eps) return 2;
        return 1;
    }
    for(int i = n - 1; i >= 0; -- i)
        for(int j = i + 1; j < n; ++j)
            a[i][n] -= a[j][n] * a[i][j];
    return 0;
}

double pow2(double a, double b)
{
    return a * a - b * b;
}

int main()
{
    cin >> n;
    for(int i = 0; i <= n; ++i) 
        for(int j = 0; j < n; ++j)
            cin >> d[i][j];
    // 增广矩阵
    for(int i = 0; i < n; ++i)
        for(int j = 0; j <= n; ++j)
        {
            if(j == n)
            {
                for(int k = 0; k < n; ++k)
                {
                    a[i][j] += pow2(d[i][k], d[i + 1][k]);
                }
            }
            else
            {
                a[i][j] = 2 * (d[i][j] - d[i + 1][j]);
            }
        }
    int t = gauss();
    for(int i = 0; i < n; ++i) printf("%.3f ", a[i][n]);
}