您的位置:首页 > 其它

高斯消元

2017-01-19 19:00 260 查看
bool equal(double source, double target)
{
if (abs(source - target) <=  0.000001)
return true;
else return false;
}
vector<double> Gauss(vector<vector<double> > A, vector<double> b,int &state)
{
//state -1 无解 0 唯一解 1 无穷解
//N个方程,M个未知数
int N = A.size(), M = A[0].size();
state = 0;
//处理上三角矩阵
for (size_t i = 0; i < M; i++)
{
bool find = false;
//令对角线上的元素 !=0
for (size_t j = i; j < N; j++)
{
if (A[j][i] != 0)
{
swap(A[i], A[j]);
find = true;
break;
}
}
//多解的情况
if (!find)
{
state = 1;
continue;
}
//将当前列第i行后的所有行都化为0(非上三角部分的化为0)
for (size_t j = i + 1; j < N; j++)
{
//double aji = A[j][i],aii= A[i][i];
double tmp = A[j][i] / A[i][i];
for (size_t k = 0; k < M; k++)
A[j][k] -= A[i][k] * tmp;
//A[j][k] -= A[i][k] * aji / aii;
b[j]-=b[i] * tmp;
}
}
//判断无解
for (size_t i = 0; i < N; i++)
{
bool allZero=true;
for (size_t j = 0; j < M; j++)
{
if (!equal(A[i][j],0))
{
allZero = false;
break;
}
}
if (allZero&&!equal(b[i],0))
{
state = -1;
return vector<double>();
}
}
//判断无穷解
if (state == 1) return vector<double>();
//求值过程
vector<double> result(M, 0);
for (int i = M - 1; i >= 0; i--)
{
for (size_t j = i + 1; j < M; j++)
{
b[i] -= A[i][j] * result[j];
A[i][j] = 0;
}
result[i] = b[i] / A[i][i];
}
return result;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签:  gauss