Skip to content

高斯消元

用于求解线性方程组。时间复杂度为 ,其中 为方程组的未知数个数。

C++模板

cpp
const double eps = 1e-7;
double a[105][105], ans[105];

int main() {
    int n, m;
    cin >> n >> m;
    for (int i = 1; i <= n; i++)
        for (int j = 1; j <= m; j++)
            cin >> a[i][j]; // 输入n行m列的矩阵,最后一列为常数项
    int rank = 0; // 记录矩阵的秩
    for (int col = 1; col <= m - 1; col++) { // m-1列(变量数)
        // 寻找当前列主元
        int pivot = -1;
        for (int i = rank + 1; i <= n; i++)
            if (fabs(a[i][col]) > eps) {
                pivot = i;
                break;
            }
        if (pivot == -1)
            continue; // 自由变量
        swap(a[++rank], a[pivot]); // 移动到有效行
        // 归一化主元行
        double div = a[rank][col];
        for (int j = col; j <= m; j++)
            a[rank][j] /= div;
        // 消去其他行当前列
        for (int i = 1; i <= n; i++) {
            if (i != rank && fabs(a[i][col]) > eps) {
                double factor = a[i][col];
                for (int j = col; j <= m; j++)
                    a[i][j] -= factor * a[rank][j];
            }
        }
    }

    // 检查无解(矛盾方程判断)
    for (int i = rank + 1; i <= n; i++)
        if (fabs(a[i][m]) > eps) {
            cout << -1 << endl;
            return 0;
        }

    if (rank < m - 1) { // 无穷多解(新增特解构造)
        vector<int> main_vars; // 主变量列
        // 记录主变量位置
        for (int r = 1, c = 1; r <= rank; r++) {
            while (fabs(a[r][c]) < eps)
                c++;
            main_vars.push_back(c++);
        }
        // 自由变量设为0,构造特解
        fill(ans, ans + 105, 0.0);
        for (int r = 1; r <= rank; r++) {
            int pc = main_vars[r - 1];
            ans[pc] = a[r][m];
            for (int j = pc + 1; j <= m - 1; j++)
                ans[pc] -= a[r][j] * ans[j];
        }
    }
    else { // 唯一解(原回代逻辑)
        for (int i = rank; i >= 1; i--) {
            int pc = find_if(a[i] + 1, a[i] + m, [](double x) {
                return fabs(x) > eps;
            }) - a[i];
            ans[pc] = a[i][m];
            for (int j = pc + 1; j <= m - 1; j++)
                ans[pc] -= a[i][j] * ans[j];
        }
    }

    // 输出结果(变量数为m-1)
    for (int i = 1; i <= m - 1; i++)
        cout << ans[i] << endl;
    return 0;
}

若无解则输出-1,若有无穷多解则输出特解。