高斯-约旦消元法详解

注：本文章假设读者已经学会基础的高斯消元法

引入

高斯约旦消元法是高斯消元法的一种，一般用于求解线性方程组。

对于一个线性方程组

\[ \begin{cases} x + 3y + 4z = 5\\ x + 4y + 7z = 3\\ 9x + 3y + 2z = 2 \end{cases}\]

我们可以将其写成一个增广矩阵的形式

\[\left[\begin{array}{ccc|c} 1 & 3 & 4 & 5\\ 1 & 4 & 7 & 3\\ 9 & 3 & 2 & 2\\ \end{array}\right]\]

我们知道，高斯消元是将增广矩阵等号左边部分化成一个阶梯化矩阵。这样可以轻易求出最后一个未知数的值，再回代入上一个方程，最终计算出每一个未知数的值。

而高斯约旦消元法没有回代的过程，而是旨在将非对角线上的系数都化为 \(0\)。

步骤

1. 寻找主元，并将主元所在行放在未操作行的第一个。

对于这个矩阵，第一步先找到第一列最大值作为主元，如图应该为 \(9\)。

\[\left[\begin{array}{ccc|c} 1 & 3 & 4 & 5\\ 1 & 4 & 7 & 3\\ 9 & 3 & 2 & 2\\ \end{array}\right]\]

换位之后该图变化为

\[\left[\begin{array}{ccc|c} 9 & 3 & 2 & 2\\ 1 & 4 & 7 & 3\\ 1 & 3 & 4 & 5\\ \end{array}\right]\]

2.进行消元

将第一列除当前行以外的行都进行消元（包括已操作过的行）。
消去第一列化为

\[\left[\begin{array}{ccc|c} 9 & 3 & 2 & 2\\ 0 & \frac{11}{3} & \frac{61}{9} & \frac{25}{9}\\ 0 & \frac{8}{3} & \frac{34}{9} & \frac{43}{9}\\ \end{array}\right]\]

第一二步轮流操作，直到形成对角线。

Code

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
#define MAXN 105

using namespace std;

int n;
double a[MAXN][MAXN];

void swap_line(int x,int y)
{ // 将两行交换 
	for(int i = 1;i <= n+1; i++)
		swap(a[x][i],a[y][i]);
	
	return ;
}

int find(int top/*top代表列数*/)// 寻找当前列的主元，返回主元所在行 
{
	// 寻找这一列中最大的系数为主元  
	int maxn = top;
	for(int i = top+1;i <= n; i++)
	{
		if( fabs(a[i][top]) > fabs(a[maxn][top]) )
			maxn = i;
	}
	return maxn;
}

int main()
{
	scanf("%d",&n);
	for(int i = 1;i <= n; i++)
	{
		for(int j = 1;j <= n+1; j++)
			scanf("%lf",&a[i][j]);
	}
	
	for(int i = 1;i <= n; i++)
	{
		int _main = find(i);// 这里 _main 表示的是主元所在行  
		swap_line(i,_main);
		if( a[i][i] == 0 )
		{// 一列全是 0 则表示方程组无解  
			printf("No Solution\n");
			return 0;
		}
		
		for(int j = 1;j <= n; j++)
		{// 消元的过程  
			if( i != j )
			{
				double tmp = a[j][i]/a[i][i];
				for(int k = i;k <= n+1;k++)
					a[j][k] -= a[i][k]*tmp;
			}
		}
	}
	
	for(int i = 1;i <= n; i++)
	{// 要除以该项系数  
		a[i][n+1] /= a[i][i];
	}
	
	for(int i = 1;i <= n; i++)
	{
		printf("%.2lf\n",a[i][n+1]);
	}
	
	return 0;
}