【C/C++】高斯列主元消元法求解线性方程组

初始增广矩阵：

矩阵当前状态：

3.0000 -1.0000 1.0000 4.0000
1.0000 1.0000 1.0000 6.0000
2.0000 3.0000 -1.0000 12.0000

方程组的解：
x1 = 1.0000
x2 = 3.0000
x3 = 2.0000

验证结果：
3x1 - x2 + x3 = 4.0000 (预期4)
x1 + x2 + x3 = 6.0000 (预期6)
2x1 + 3x2 - x3 = 12.0000 (预期12)

#include<stdio.h>
#include<math.h>

#define MAX_MATRIX 10

/**
 * @brief SwapRow 进行行交换
 * @param m 待计算的矩阵
 *        row 待交行的行
 *        max_row 待交换的另一行
 *        n 矩阵行数
 */
static void SwapRow(double m[][MAX_MATRIX], int row, int max_row, int n) {
	double swap;
	for (int k = row; k <= n; k++) {
		swap = m[row][k];
		m[row][k] = m[max_row][k];
		m[max_row][k] = swap;
	}
}

/**
 * @brief 组上三角矩阵
 * @param m 待计算的矩阵
 *        n 矩阵行数
 */
static void SelectColE(double m[][MAX_MATRIX], int n) {
	int max_row_e = 0;  //主元所在行
	double ratio = 0;   //消元因数
	for (int j = 0; j < n; j++) {
		max_row_e = j;
		for (int i = j; i < n; i++) {
			if (fabs(m[i][j]) > fabs(m[max_row_e][j])) {
				max_row_e = i;
			}
		}
		if (max_row_e != j) {
			SwapRow(m, j, max_row_e, n);   //与最大主元所在行交换
		}
		//消元
		for (int i = j + 1; i < n; i++) {
			ratio = m[i][j] / m[j][j];
			for (int k = j; k < n + 1; k++) {
				m[i][k] -= m[j][k] * ratio;
			}
		}
	}
}

/**
 * @brief: Gauss 高斯列主元消元法求解线性方程(A*X = B)
 * @param: m 由于A|B组成的增广矩阵,X为待求的解
 *         n 求解的元数,n要小于MAX_MATRIX
 * @result:所求结果存放在m[][n]中
 */
void Gauss(double m[][MAX_MATRIX], int n) {
	SelectColE(m, n);   // 列选主元并消元成上三角
	// 回代求解,结果存在m[][n]中
	for(int i = n - 1; i >= 0; i--) {
		for(int j = i + 1; j < n; j++) {
			m[i][n] -= m[i][j] * m[j][n];
		}
		m[i][n] /= m[i][i];
	}
}

double a[3][MAX_MATRIX] = {
	{3,-1, 1, 4},  //A|B
	{1, 1, 1, 6},
	{2, 3,-1, 12}
};

int main(int argc ,char **argv) {
	Gauss(a, 3);
	printf("%f,%f,%f\r\n",a[0][3], a[1][3], a[2][3]);
	return 0;
}