跳到主要内容

高斯消元

  • 可在O(n3)O(n^3)求解线性方程组
    • 找到绝对值最大的一行
    • 将该行换到最上面
    • 将该行第一个数化成11
    • 把该列下边的数化成00
  • 注意部分过程需要从后往前遍历,否则出现写后读问题

代码

// a[N][N+1]是增广矩阵
int gauss()
{
int c, r;
// 按列遍历,化成行阶梯矩阵
for (c = 0, r = 0; c < n; c ++)
{
// 找到当前列绝对值最大元素所在的行(搜索第r行~最后一行)
int t = r;
for (int i = r + 1; i < n; i ++ )
if (fabs(a[i][c]) > fabs(a[t][c]))
t = i;

if (fabs(a[t][c]) < eps) continue; // 当前列全是0

for (int j = c; j <= n; j ++ ) swap(a[t][j], a[r][j]); // 将绝对值最大的行换到最顶端
for (int j = n; j >= c; j -- ) a[r][j] /= a[r][c]; // 将当前上的首位变成1,注意从后往前遍历
for (int i = r + 1; i < n; i ++ ) // 用当前行将下面所有的列消成0
if (fabs(a[i][c]) > eps)
for (int j = n; j >= c; j -- )
a[i][j] -= a[r][j] * a[i][c]; // 注意从后往前删,否则出现写后读错误

r ++ ;
}

if (r < n)
{
for (int i = r; i < n; i ++ )
if (fabs(a[i][n]) > eps)
return 2; // 出现0!=0,无解
return 1; // 都是0=0,有无穷多组解
}

// 化成单位阵,增广矩阵的扩展部分为方程的解
for (int i = n - 1; i >= 0; i -- )
for (int j = i + 1; j < n; j ++ )
a[i][n] -= a[i][j] * a[j][n];

return 0; // 有唯一解
}