G14 高斯消元法 • EH's Blog

题目情境#

题目链接：https://www.luogu.com.cn/problem/P3389 ↗

题目描述

给定一个线性方程组，对其求解。

$\begin{cases} a_{1, 1} x_1 + a_{1, 2} x_2 + \cdots + a_{1, n} x_n = b_1 \\ a_{2, 1} x_1 + a_{2, 2} x_2 + \cdots + a_{2, n} x_n = b_2 \\ \cdots \\ a_{n,1} x_1 + a_{n, 2} x_2 + \cdots + a_{n, n} x_n = b_n \end{cases}$

输入格式

第一行，一个正整数 $n$ 。

第二至 $n+1$ 行，每行 $n+1$ 个整数，为 $a_1, a_2, \dots ,a_n$ 和 $b$ ，代表一组方程。

输出格式

共 $n$ 行，每行一个数，第 $i$ 行为 $x_i$ （四舍五入保留 $2$ 位小数）。

如果不存在唯一解或无解，在第一行输出 No Solution.

输入

plaintext

输出

-0.97
5.18
-2.39

plaintext

说明/提示

数据范围： $1 \leq n \leq 100, \left | a_i \right| \leq {10}^4 , \left |b \right| \leq {10}^4$ 。保证数据若有解则所有解均满足 $|x_i|\le 10^3$ ，且 $x_i\pm 10^{-6}$ 和 $x_i$ 四舍五入后的结果相同（即不会因为较小的精度误差导致四舍五入后的结果不同）。

算法解析：#

高斯消元：逐列选主元（当前列绝对值最大），与当前行交换，主元行归一化，再用行变换消去下方行的当前列；回代求各变量。若某列找不到非零主元且常数项非零则无解（No Solution），否则唯一解。复杂度 $O(n^3)$ 。

Python代码实现#


def gauss():
    c, r = 0, 0
    for c in range(n):
        t = r
        for i in range(r, n):
            if abs(a[i][c]) > abs(a[t][c]):
                t = i
        if abs(a[t][c]) < 1e-6:
            continue
        for i in range(c, n + 1):
            a[t][i], a[r][i] = a[r][i], a[t][i]
        for i in range(n, c - 1, -1):
            a[r][i] /= a[r][c]
        for i in range(r + 1, n):
            if abs(a[i][c]) > 1e-6:
                for j in range(n, c - 1, -1):
                    a[i][j] -= a[i][c] * a[r][j]
        r += 1
    if r < n:
        for i in range(r, n):
            if abs(a[i][n]) > 1e-6:
                return 2
        return 1
    for i in range(n - 1, -1, -1):
        for j in range(i + 1, n):
            a[i][n] -= a[i][j] * a[j][n]
    return 0

n = int(input())
a = []
for i in range(n):
    a.append(list(map(float, input().split())))
t = gauss()
if t:
    print("No Solution")
else:
    for i in range(n):
        print("{:.2f}".format(a[i][n]))

python

C++代码实现#

#include <iostream>
#include <algorithm>
#include <cmath>
using namespace std;

const int N=110;
const double eps=1e-6;
int n;
double a[N][N]; //增广矩阵

int gauss(){
  int c,r;//当前列，当前行
  for(c=r=0; c<n; c++){
    //1.找到c列的最大行t
    int t=r;
    for(int i=r; i<n; i++)
      if(fabs(a[i][c])>fabs(a[t][c])) t=i;
    if(fabs(a[t][c])<eps) continue; //c列已0化

    //2.把最大行换到上面
    for(int i=c; i<n+1; i++)swap(a[t][i],a[r][i]);

    //3.把当前行r的第一个数，变成1
    for(int i=n; i>=c; i--)a[r][i]/=a[r][c];

    //4.把当前列c下面的所有数，全部消成0
    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) //a[i][n]==bi
        return 2; //左边=0，右边≠0,无解
    return 1; //0==0，无穷多解
  }
  //5.唯一解，从下往上回代，得到方程的解
  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;
}
int main(){
  cin >> n;
  for(int i=0; i<n; i++)
    for(int j=0; j<=n; j++)
      cin >> a[i][j];
  int t=gauss();
  if(t) puts("No Solution");
  else for(int i=0; i<n; i++)
        printf("%.2lf\n",a[i][n]);
}