目录
列主元高斯消去法原理
在基本高斯消去法的消元过程中并没有考虑任何数值方面的问题,事实上这方面的问题是常见的,也是不能忽略的,即当主元,且很小时,高斯消去法虽然能执行下去,但用作为主元计算行乘数时,会扩大误差,导致结果不可靠,甚至严重失真。
基本高斯消去法的求解过程如下(具体原理参考https://blog.csdn.net/weixin_41788456/article/details/102485139):
设Ax=b,,若A的所有顺序主子式均不为零,则基本高斯消元无需换行进行到底,得到唯一解,其消元和回代的计算公式为:
(1)消元计算 对于
(2)回代计算
在上面的求解过程中可以看出,适当的改变算法可以有效的提高精度,其基本思想是,在高斯消去法的每一步先选择一个绝对值最大的元素,在进行行行交换、消元,即变换到第k步时,从第k列的及以下的各元素中选出绝对值最大者,然后通过行变换将它交换到主元素的位置上,再用其消去主对线以下的其他元素,最后变为同解的上三角形方程组,这种方法称为列主元高斯消去法。
列主元高斯消去法流程图
C++程序源代码
//列主元高斯消去法实现
//开发人员:chenshuai 开发日期:2019.11.24 邮箱:chenshuai0614@hrbeu.edu.cn
#include "pch.h"
#include <iostream>//基本数据流输入/输出
#include <iomanip> //参数化输入/输出
#include <vector>//STL动态数组容器
using namespace std;
//************************
//列主元高斯消去法公式
//***********************
vector<double> cpe_gaussian_elimination(vector<vector<double>>a, vector<double>b); //高斯消去法求解线性方程组AX=B
vector<double> cpe_gaussian_elimination(vector<vector<double>>a, vector<double>b)
{
int n = size(b);
vector<double>x; //定义方程组解
x.resize(n);
vector<double>mi_k; //定义消去过程中的中间变量
mi_k.resize(n);
double sum,max=0,c;
//n-1步消元
for (int k = 0; k < n - 1; k++)
{
//列主元,找到绝对值最大的主元
if (a[k][k] == 0 || fabs(a[k][k]) < epslion)//列主元的条件就是主元为0或者主元小于某一个值
{
max = fabs(a[k][k]);
for (int i = k; i < n; i++)
{
if (max < fabs(a[i][k]))
{
max = fabs(a[i][k]);
}
}
//交换该行
for (int i = k; i < n; i++)
{
if (max == fabs(a[i][k]))
{
for (int j = k; j < n; j++)
{
c = a[k][j];
a[k][j] = a[i][j];
a[i][j] = c;
}
c = b[k];
b[k] = b[i];
b[i] = c;
}
}
//求出第i次初等行变换系数
for (int j = k + 1; j < n; j++)
{
mi_k[j] = a[j][k] / a[k][k];
}
for (int i = k + 1; i < n; i++)
{
for (int j = 0; j < n; j++)
{
a[i][j] = a[i][j] - mi_k[i] * a[k][j];
}
b[i] = b[i] - mi_k[i] * b[k];
}
} //回代过程
x[n - 1] = b[n - 1] / a[n - 1][n - 1];
for (int i = n - 2; i >= 0; i--)
{
sum = 0;
for (int j = i + 1; j < n; j++)
{
sum = sum + a[i][j] * x[j];
}
x[i] = (b[i] - sum) / a[i][i];
}
return x;
}
int main()
{
int n = 2;
vector<vector<double>>a;
a.resize(n, vector<double>(n));
vector<double>b(n);
vector<double>x;
a[0] = { 0.007,-0.8 }, a[1] = {-0.1 ,10};
b = { 0.7,10};
x = cpe_gaussian_elimination(a, b);
cout << "解为:" << endl;
for (int i = 0; i < n; i++)
cout <<"x["<<i<<"]="<< fixed << setprecision(2) << setw(5) << x[i] << endl;
}
实例
用列主元高斯消去法求解线性方程组: