本文已参与「新人创作礼」活动，一起开启掘金创作之路。

高斯消元

前置知识：增广矩阵是什么，矩阵行列变换，对角阵等

高斯消元主要用于多元线性方程组的求解，只对系数进行操作，不对未知数进行操作。

高斯消元知识点可以参考以下链接：

www.tonyyin.top/2021/09/gau…

oi-wiki.org/math/gauss/

blog.nowcoder.net/n/5705bab81…

我习惯把高斯消元分为几步：

化为上三角矩阵
- 找主元（主元即对角线上位置的元素）
  即按列来找，找每一列中绝对值最大的一行，将该行和当前访问到的行进行交换
- 交换
- 主元归一
  即将当前的对角元素变为1，对应该行也要做相同的变换
- 消元
化为对角阵（主对角线都为1，其余位置都为0）

例1 球形空间产生器

题目链接：www.acwing.com/problem/con…

有一个球形空间产生器能够在 n 维空间中产生一个坚硬的球体。
现在，你被困在了这个 n 维球体中，你只知道球面上 n+1 个点的坐标，你需要以最快的速度确定这个 n 维球体的球心坐标，以便于摧毁这个球形空间产生器。

思路

因为是n维空间，一个点的坐标有n个值。可以设球心坐标为 $(x_1,x_2,...,x_n)$ 表示n维情况下的坐标,设半径为 $R$ 给出了n+1个点的坐标，设第i个点的坐标为 $(a_{i,1},a_{i,2},...,a_{i,n})$ 球形可以想到球上的点到球心的距离相等，进而可以列出方程组

$\begin{cases} (a_{1,1}-x_1)^2+...+(a_{1,n}-x_n)^2=R^2【式1】\\ (a_{2,1}-x_1)^2+...+(a_{2,n}-x_n)^2=R^2【式2】\\ ......\\ (a_{n+1,1}-x_1)^2 +...+(a_{n+1,n}-x_n)^2=R^2【式n+1】 \end{cases}$

可以发现并不是满足线性的条件，但是我们却可以把它化成线性的情况

从第二个式子开始，依次减去第一个式子，然后我们就可以得到

$\begin{cases} 2(a_{2,1}-a_{1,1})x_1+...+2(a_{2,n}-a_{1,n})x_n=(a_{2,1}^2+...+a_{2,n}^2)-(a_{1,1}^2+...+a_{n,n}^2)\\ 2(a_{3,1}-a_{1,1})x_1+...+2(a_{3,n}-a_{1,n})x_n=(a_{3,1}^2+...+a_{3,n}^2)-(a_{1,1}^2+...+a_{n,n}^2)\\ ......\\ 2(a_{n+1,1}-a_{1,1})x_1+...+2(a_{n+1,n}-a_{1,n})x_n=(a_{n+1,1}^2+...+a_{n+1,n}^2)-(a_{1,1}^2+...+a_{n,n}^2) \end{cases}$

这就是线性的了，而且共有n个式子，有n个未知数，刚好可以解出来，式子太长我们做一下代换

$b_{1,1}=2(a_{2,1}-a_{1,1})$

$b_{n,1}=2(a_{n+1,1}-a_{1,1})$

$b_{1,n+1}=(a_{2,1}^2+...+a_{2,n}^2)-(a_{1,1}^2+...+a_{n,n}^2)$

$b_{n,n+1}=(a_{n+1,1}^2+...+a_{n+1,n}^2)-(a_{1,1}^2+...+a_{n,n}^2)$

其它的大家可以类推

然后就变成了这个矩阵，对这个矩阵进行高斯消元

$\begin{bmatrix} b_{1,1}&...&b_{1,n+1}\\ .&...&.\\ .&...&.\\ b_{n,1}&...&b_{n,n+1} \end{bmatrix}$

消元步骤：

化为上三角矩阵
- 找主元（主元即对角线上位置的元素）找该列绝对值最大的一个数，将那行和本行交换
- 交换
- 主元归一

即将当前的对角元素变为1，对应该行也要做相同的变换 - 消元即对角元素下面的值都要变为0，那么需要下面的每一行的所有值减去对角元素那一行元素的一个倍数

化为对角阵（主对角线都为1，其余位置都为0）

代码：

#include<bits/stdc++.h>
using namespace std;
const int N = 15;
int n;
double a[N][N],b[N][N];

void gauss()
{
	//化为上三角矩阵
	for(int r=1,c=1;c<=n;r++,c++)//主循环每次查一个对角元素。按列查，主要要看对角元素，r=c
	{
		//找主元，找该列绝对值最大的一个数，将那行和本行交换
		int t = r;
		for(int i=r+1;i<=n;i++)
			if(fabs(b[i][c]) > fabs(b[t][c]))
				t = i;
		//交换，行的交换
		for(int i=c;i<=n+1;i++)	swap(b[r][i],b[t][i]);
		//主元归一，将对角元素化为1，即该行所有值除以该行对应的对角元素
		for(int i=n+1;i>=c;--i)	b[r][i] /= b[r][c];
		//消元，即对角元素下面的值都要变为0，那么需要下面的每一行的所有值减去对角元素那一行元素的一个倍数（b[i][c]）
		for(int i=r+1;i<=n;i++)
			for(int j=n+1;j>=c;j--)
				b[i][j] -= b[r][j]*b[i][c];
	}
	//化为对角阵（只有对角线有元素且为1）
	for(int i=n;i>=1;i--)//按列消，对角元素行=列
		for(int j=i-1;j>=1;j--)//消去第i行上面所有的元素
		{
			b[j][n+1] -= b[i][n+1]*b[j][i];//只对右边的值做变化
			b[j][i] = 0;//最后一定消完为0
		} 
}
int main()
{
	scanf("%d",&n);
	for(int i=0;i<n+1;i++)
		for(int j=1;j<=n;j++)
			scanf("%lf",&a[i][j]);
	for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++)
		{
			b[i][j] = 2 * (a[i][j]-a[0][j]);
			b[i][n+1] += a[i][j]*a[i][j] - a[0][j]*a[0][j];
		}
	gauss();
	for(int i=1;i<=n;i++)
		printf("%.3lf ",b[i][n+1]);
	return 0;
}

高斯消元及实现

高斯消元

例1 球形空间产生器

思路

代码：