矩阵初等行变换法则

任一行可以与另一行进行加减。
任一行可以乘或除以一个非零常数（除其实就是乘一个倒数）。
任两行可以交换位置。

线性方程组

形如

a_{1,1}x_1+a_{1,2}x_2+\dots+a_{1,n}x_n=b_1 \\ a_{2,1}x_1+a_{2,2}x_2+\dots+a_{2,n}x_n=b_n \\ \vdots \\ a_{n,1}x_1+a_{n,2}x_2+\dots+a_{n,n}x_n=b_n

其中，系数矩阵为

A=\left( \begin{matrix} &a_{1,1} &a_{1,2} &\dots &a_{1,n} \\ &a_{2,1} &a_{2,2} &\dots &a_{2,n} \\ &\vdots &\vdots &\vdots &\vdots \\ &a_{n,1} &a_{n,2} &\dots &a_{n,n} \\ \end{matrix} \right)

右值向量（矩阵）为

B=\left ( \begin{matrix} b_1\\ b_2\\ \vdots\\ b_n \end{matrix} \right)

解向量为

X=\left ( \begin{matrix} x_1\\ x_2\\ \vdots\\ x_n \end{matrix} \right)

因此，方程组可表示为

AX=B

矩阵的秩

矩阵可由初等行变换化为行最简形矩阵，所谓行最简型矩阵，即在阶梯形矩阵中，若非零行的第一个非零元素全是1，且非零行的第一个元素1所在列的其余元素全为零，就称该矩阵为行最简形矩阵。矩阵的秩就是行最简形矩阵非零行的个数，以 $r(M)$ 来表示矩阵 $M$ 的秩。如：

M=\left (\begin{matrix} &1 &0 &0 &-1 \\ &0 &1 &0 &-2 \\ &0 &0 &1 &2 \end{matrix} \right)

则 $r(M)=3$

高斯消元法(列主元法)

其实就是线性代数中的矩阵行化简算法。

思路

要解上述方程组，需要引入增广矩阵

(A\vdots b)= \left( \begin{matrix} &a_{1,1} &a_{1,2} &\dots &a_{1,n} &b_1 \\ &a_{2,1} &a_{2,2} &\dots &a_{2,n} &b_2 \\ &\vdots &\vdots &\vdots &\vdots &\vdots \\ &a_{n,1} &a_{n,2} &\dots &a_{n,n} &b_n \\ \end{matrix} \right)

其实就是在系数矩阵 $A$ 右侧添加右值向量 $b$

若 $r(A)=r(A\vdots b)$ 且 $r(A)=n$ ，则方程组有唯一解
若 $r(A)=r(A\vdots b)$ 且 $r(A)\ne n$ ，则方程组有无穷个解
若 $r(A)\ne r(A\vdots b)$ ，则方程组无解

算法思想

假设行数为 $1 \sim n$ ，列数为 $1\sim n+1$

化简矩阵

初始化当前行为 $i=1$
在 $i\sim n$ 行中寻找绝对值最大的 $a_{ii}$ 所在行 $j$ （最大系数可减小误差）
若 $a_{jj}=0$ 则说明 $r(A)\ne n$ ，无唯一解，返回 $false$
交换第 $i,j$ 两行，使得增广矩阵保持为上三角矩阵
第 $i$ 行所有元素除以系数 $a_{ii}$
若 $i=n$ 说明这是末尾行，结束矩阵化简，在求解向量后返回1
第 $i+1\sim n$ 行，减去 $a_{i+1,i}$ 倍第 $i$ 行，消除其余行的第 $i$ 列系数
$i=i+1$ ，跳回到第 $2$ 步，寻找下一行

求解向量

此时矩阵为

(A\vdots b)= \left( \begin{matrix} &1 &a_{1,2} &\dots &a_{1,n} &b_1 \\ & &1 &\dots &a_{2,n} &b_2 \\ & & &\ddots & &\vdots \\ & & & &1 &b_n \\ \end{matrix} \right)

即

x_1+a_{1,2}x_2+\dots +a_{1,n}x_n=b_1 \\ x_2+a_{2,3}x_3+\dots +a_{2,n}x_n=b_2 \\ x_n=b_n

此时有

x_n=b_n \\ x_{n-1}=b_{n-1}-a_{n-1,n}*b_n \\ \vdots \\ x_1=b_{1}-a_{1,n}*b_n-a_{1,n-1}*b_{n-1}-\dots a_{1,2}a_2

解向量为

X=\left( \begin{matrix} x_1\\ x_2\\ \vdots\\ x_n \end{matrix} \right)

算法模板

int gauss(double num[100][101],int n,double x[]){
	for(int i=0;i<n;i++){//循环n次,第i轮循环行为i~n-1，列为i~n
		int maxRow=i;//maxRow记录系数最大的行，作为被减行减小误差
		for(int j=i+1;j<n;j++){
			if(abs(num[j][i])>abs(num[maxRow][i])) maxRow=j;
		}
		if(abs(num[maxRow][i])<zero) return 0;//x系数为0则增广矩阵无唯一解,返回0
		if(maxRow!=i){//交换最大行到i行,使之保持为上三角矩阵
			for(int j=i;j<n+1;j++){
				swap(num[maxRow][j],num[i][j]);
			}
		}
		for(int j=n;j>=i;j--){//化最大行第一个系数为1
			num[i][j]/=num[i][i];//从后向前除以系数,否则需要临时变量记录[i][i]的系数
		}
		for(int j=i+1;j<n;j++){//被系数行减去
			for(int k=n;k>=i;k--){
				num[j][k]-=num[j][i]*num[i][k];//减去了系数行乘以对应系数
			}
		}
	}
	for(int i=n-1;i>=0;i--){//逆向求解向量
		x[i]=num[i][n];//赋初值使得ax=b
		for(int j=i+1;j<n;j++)
			x[i]-=num[i][j]*x[j];//减去其他解向量
	}
	return 1;
}

例题

题目链接

题目背景

Gauss消元

题目描述

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

输入格式

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

输出格式

共 $n$ 行，每行一个数，第 $i$ 行为 $x_i$ （保留2位小数）如果不存在唯一解，在第一行输出"No Solution".

输入输出样例

输入 #1

输出 #1

-0.97
5.18
-2.39

说明/提示

$1 \leq n \leq 100, \left | a_i \right| \leq {10}^4 , \left |b \right| \leq {10}^4$

AC代码

#include <iostream>
#include <algorithm>
#define zero 1e-10
using namespace std;
int gauss(double num[100][101],int n,double x[]){
	for(int i=0;i<n;i++){//循环n次,第i轮循环行为i~n-1，列为i~n
		int maxRow=i;//maxRow记录系数最大的行，作为被减行减小误差
		for(int j=i+1;j<n;j++){
			if(abs(num[j][i])>abs(num[maxRow][i])) maxRow=j;
		}
		if(abs(num[maxRow][i])<zero) return 0;//x系数为0则增广矩阵无唯一解,返回0
		if(maxRow!=i){//交换最大行到i行,使之保持为上三角矩阵
			for(int j=i;j<n+1;j++){
				swap(num[maxRow][j],num[i][j]);
			}
		}
		for(int j=n;j>=i;j--){//化最大行第一个系数为1
			num[i][j]/=num[i][i];//从后向前除以系数,否则需要临时变量记录[i][i]的系数
		}
		for(int j=i+1;j<n;j++){//被系数行减去
			for(int k=n;k>=i;k--){
				num[j][k]-=num[j][i]*num[i][k];//减去了系数行乘以对应系数
			}
		}
	}
	for(int i=n-1;i>=0;i--){//逆向求解向量
		x[i]=num[i][n];//赋初值使得ax=b
		for(int j=i+1;j<n;j++)
			x[i]-=num[i][j]*x[j];//减去其他解向量
	}
	return 1;
}
int main(){
	int n;
	double num[100][101];//矩阵大小是n*n+1
	double x[100];//存储解向量x
	scanf("%d",&n);
	for(int i=0;i<n;i++){
		for(int j=0;j<n+1;j++){
			scanf("%lf",&num[i][j]);
		}
	}
	if(gauss(num,n,x)){
		for(int i=0;i<n;i++){
			printf("%.2lf\n",x[i]);
		}
	}else{
		printf("No Solution");
	}
	return 0;
}
/*
3
1 3 4 5
1 4 7 3
9 3 2 2
*/

【C++】高斯消元算法