Loping151 Blog
线性方程组 P2455

线性方程组 P2455

Sun Sep 22 2024
11 分钟
算法 C++ 线性代数
1570 字

今天的LeetCode题目太弱智了，补一道别的。

题目#

已知 nn 元线性一次方程组。

{a1,1x1+a1,2x2+⋯+a1,nxn=b1a2,1x1+a2,2x2+⋯+a2,nxn=b2⋯an,1x1+an,2x2+⋯+an,nxn=bn\begin{cases} a_{1, 1} x_1 + a_{1, 2} x_2 + \cdots + a_{1, n} x_n = b_1 \\ a_{2, 1} x_1 + a_{2, 2} x_2 + \cdots + a_{2, n} x_n = b_2 \\ \cdots \\ a_{n,1} x_1 + a_{n, 2} x_2 + \cdots + a_{n, n} x_n = b_n \end{cases}

请根据输入的数据，编程输出方程组的解的情况。

输入格式#

第一行输入未知数的个数 nn。
接下来 nn 行，每行 n+1n + 1 个整数，表示每一个方程的系数及方程右边的值。

输出格式#

如果有唯一解，则输出解。你的结果被认为正确，当且仅当对于每一个 xix_i 而言结果值与标准答案值的绝对误差或者相对误差不超过 0.010.01。

如果方程组无解输出 −1-1； 如果有无穷多实数解，输出 00；

样例#

样例输入#

PLAINTEXT
1
2
3
4
3
2 -1 1 1
4 1 -1 5
1 1 1 0

样例输出#

PLAINTEXT
1
2
3
x1=1.00
x2=0.00
x3=-1.00

提示#

【数据范围】
对于 100%100\% 的数据，1≤n≤501 \le n \le 50。对于 ∀1≤i,j≤n\forall 1\le i, j \le n，有 ∣ai,j∣≤100|a_{i, j}| \le 100，∣bi∣≤300|b_i| \le 300。

分析与解答#

我去，线性代数。这就不得不掏出我当年的笔记了。

这下实名上网了，大家不要学我。

不过这个题目没有这么复杂，指定了是方阵，不需要考虑解系啥的。我隐约记得当时隔壁非荣誉好像当时就要求实现了来着。。

题目要求解，首先想到的是Crammer法则，利于递归计算。但是这个方法的时间复杂度是 O(n!)O(n!)，对于 n=50n=50 的数据显然不可行。别的方法，我想到的是SVD，毕竟最近在写一些物理仿真，经常SVD。但要精确解的话不行。说到这个，SVD以后有机会也得重温下。

仔细想想，好像只能暴力的消元法了？试试。

CPP
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
#include <iostream>
#include <cmath>
#include <iomanip>
using namespace std;

const double EPS = 1e-9;
const int MAXN = 50;

double a[MAXN][MAXN + 1];
double x[MAXN];
int n;

int gauss() {
    for (int i = 0; i < n; i++) {
        int max_row = i;
        for (int j = i + 1; j < n; j++) {
            if (fabs(a[j][i]) > fabs(a[max_row][i])) {
                max_row = j;
            }
        }

        if (fabs(a[max_row][i]) < EPS) {
            continue;
        }

        for (int j = i; j <= n; j++) {
            swap(a[i][j], a[max_row][j]);
        }

        for (int j = 0; j < n; j++) {
            if (j != i && fabs(a[j][i]) > EPS) {
                double factor = a[j][i] / a[i][i];
                for (int k = i; k <= n; k++) {
                    a[j][k] -= factor * a[i][k];
                }
            }
        }
    }

    for (int i = 0; i < n; i++) {
        bool all_zero = true;
        for (int j = 0; j < n; j++) {
            if (fabs(a[i][j]) > EPS) {
                all_zero = false;
                break;
            }
        }
        if (all_zero && fabs(a[i][n]) > EPS) {
            return -1;
        }
        if (all_zero && fabs(a[i][n]) < EPS) {
            return 0;
        }
    }

    for (int i = 0; i < n; i++) {
        x[i] = a[i][n] / a[i][i];
    }

    return 1;
}

int main() {
    cin >> n;
    for (int i = 0; i < n; i++) {
        for (int j = 0; j <= n; j++) {
            cin >> a[i][j];
        }
    }

    int result = gauss();

    if (result == -1) {
        cout << -1 << endl; 
    } else if (result == 0) {
        cout << 0 << endl;
    } else {
        for (int i = 0; i < n; i++) {
            cout << fixed << setprecision(2) << "x" << i + 1 << "=" << x[i] << endl;
        }
    }

    return 0;
}

属于是在GPT的帮助下进行复健了。完整的写cpp还是没能一下子找回感觉。

为了避免可能的精度问题，总是选择每一个元最大的系数去消元。使用一个小量控制为0的精度。

然后上面的代码拿了90分。。错了的3个，我非常赖皮的发现直接输出-1能过，一共有4个，说明无解不对。无解是通过有没有一行系数全0但b不为0判断的。

上面的代码，如果遇到某一列全0，会交给最后的判断去判断无解还是无穷解。仔细想的话，这样其实是变成了n−1n-1元方程组。最后这一行一定全是0，但是当b为0的时候不一定就是无穷解，还需要继续检查完。因为去掉这一行一列以后方程组解的情况和这一行无关。

CPP
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
    int ans = 1;

    for (int i = 0; i < n; i++) {
        bool all_zero = true;
        for (int j = 0; j < n; j++) {
            if (fabs(a[i][j]) > EPS) {
                all_zero = false;
                break;
            }
        }
        if (all_zero && fabs(a[i][n]) > EPS) {
            return -1;
        }
        if (all_zero && fabs(a[i][n]) < EPS) {
            asn = 0;
        }
    }

    if (ans == 1)
        for (int i = 0; i < n; i++) {
            x[i] = a[i][n] / a[i][i];
        }

    return ans;

这样一来又过了一个。剩下两个还是没有过。吐了。于是我仔细想了反例：

TXT
1
2
3
4
3 0 0 0 0
0 0 0 1 1
0 0 0 1 2
0 0 0 0 0

写下来以后很快就意识到问题了：如果有全零列，不应该直接跳过，而是原地消下一个元。

CPP
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
#include <iostream>
#include <cmath>
#include <iomanip>
using namespace std;

const double EPS = 1e-9;
const int MAXN = 50;

double a[MAXN][MAXN + 1];
double x[MAXN];
int n;
int skip = 0;

int gauss() {
    for (int i = 0; i < n; i++) {
        int max_row = i - skip;
        for (int j = i - skip + 1; j < n; j++) {
            if (fabs(a[j][i]) > fabs(a[max_row][i])) {
                max_row = j;
            }
        }

        if (fabs(a[max_row][i]) < EPS) {
            skip += 1;
            continue;
        }

        for (int j = i; j <= n; j++) {
            swap(a[i - skip][j], a[max_row][j]);
        }

        for (int j = 0; j < n; j++) {
            if (j != i - skip && fabs(a[j][i]) > EPS) {
                double factor = a[j][i] / a[i - skip][i];
                for (int k = i; k <= n; k++) {
                    a[j][k] -= factor * a[i - skip][k];
                }
            }
        }
    }

    int ans = 1;

    for (int i = 0; i < n; i++) {
        bool all_zero = true;
        for (int j = 0; j < n; j++) {
            if (fabs(a[i][j]) > EPS) {
                all_zero = false;
                break;
            }
        }
        if (all_zero && fabs(a[i][n]) > EPS) {
            return -1;
        }
        if (all_zero && fabs(a[i][n]) < EPS) {
            ans = 0;
        }
    }

    if (ans == 1)
        for (int i = 0; i < n; i++) {
            x[i] = a[i][n] / a[i][i];
        }

    return ans;
}

int main() {
    cin >> n;
    for (int i = 0; i < n; i++) {
        for (int j = 0; j <= n; j++) {
            cin >> a[i][j];
        }
    }

    int result = gauss();

    if (result == -1) {
        cout << -1 << endl; 
    } else if (result == 0) {
        cout << 0 << endl;
    } else {
        for (int i = 0; i < n; i++) {
            cout << fixed << setprecision(2) << "x" << i + 1 << "=" << x[i] << endl;
        }
    }

    return 0;
}

其实i-skip直接写成新变量rank更好，不然容易漏掉。。提交了无数次才发现，感觉可以自裁了。