高斯消元小结

Contents
  1. 1. Ⅰ. 高斯消元认知
  2. 2. Ⅱ. 高斯消元模版
    1. 2.1. 浮点数高斯消元
    2. 2.2. 模意义下高斯消元
    3. 2.3. XOR高斯消元
  3. 3. Ⅲ. 题目讲解
    1. 3.1. HDU 4200 Bad Wiring
    2. 3.2. SGU 275 To xor or not to xor
    3. 3.3. HDU 2262 Where is the canteen

Ⅰ. 高斯消元认知

  • $这东西英文名叫Gaussian Elimination$
  • $高斯消元对矩阵进行操作,对满足一定运算规律的方程利用矩阵初等变换进行消元$
  • $主要用来求解线性方程组、矩阵的秩、以及可逆方阵的逆矩阵$
  • $还可以用来线性空间的基向量,比如线性基(线性无关的极大子集)$
  • $PS:关于线性相关,具体可以读2014集训队论文《浅谈线性相关》$

Ⅱ. 高斯消元模版

浮点数高斯消元

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
double a[N][N], ans[N];
bool isFreeX[N];
void getAns(int n, int m, int r) {
// memset(ans, 0, sizeof ans); //not necessary
for(int i = r - 1; ~i; --i) {
for(int j = 0; j < m; ++j) {
if(!sgn(a[i][j])) continue;
ans[j] = a[i][m];
for(int k = j + 1; k < m; ++k)
ans[j] -= a[i][k] * ans[k];
ans[j] = ans[j] / a[i][j];
break;
}
}
}
int gauss(int n, int m) {
for(int i = 0; i < m; ++i) isFreeX[i] = false;
int r = 0, c = 0;
for(; r < n && c < m; ++r, ++c) {
int maxR = r; //row transform
for(int i = r + 1; i < n; ++i)
if(abs(a[i][c]) > abs(a[maxR][c])) maxR = i;
if(maxR != r) swap(a[maxR], a[r]);
if(!sgn(a[r][c])) { --r; isFreeX[c] = true; continue;}
//eliminate coefficient
for(int i = r + 1; i < n; ++i) {
if(sgn(a[i][c])) {
double delta = a[i][c] / a[r][c];
for(int j = c; j <= m; ++j)
a[i][j] -= delta * a[r][j];
}
}
}
for(int i = r; i < n; i++) if(sgn(a[i][m])) return -1;
getAns(n, m, r);
//at last, r is rank, m - r is the number of freeX
return r;
}

模意义下高斯消元

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
int a[N][N], ans[N];
bool isFreeX[N];
inline int inv(int x) {return 1;}
int getAns(int n, int m, int r) {
int ret = 0;
for(int i = r - 1; ~i; --i) {
for(int j = 0; j < m; ++j) {
if(!a[i][j]) continue;
ans[j] = a[i][m];
for(int k = j + 1; k < m; ++k) {
ans[j] -= a[i][k] * ans[k];
ans[j] %= MOD;
if(ans[j] < 0) ans[j] += MOD;
}
ans[j] = ans[j] * inv(a[i][j]) % MOD;
break;
}
}
for(int i = 0; i < m; ++i) ret += ans[i];
return ret;
}
int gauss(int n, int m) {
for(int i = 0; i < m; ++i) isFreeX[i] = false;
int r = 0, c = 0;
for(; r < n && c < m; ++r, ++c) {
int maxR = r; //row transform
for(int i = r + 1; i < n; ++i)
if(abs(a[i][c]) > abs(a[maxR][c])) maxR = i;
if(maxR != r) swap(a[maxR], a[r]);
if(!a[r][c]) { --r; isFreeX[c] = true; continue;}
//eliminate coefficient
for(int i = r + 1; i < n; ++i) {
if(a[i][c]) {
int delta = a[i][c] * inv(a[r][c]);
for(int j = c; j <= m; ++j) {
a[i][j] -= delta * a[r][j];
a[i][j] %= MOD;
if(a[i][j] < 0) a[i][j] += MOD;
}
}
}
}
for(int i = r; i < n; i++) if(a[i][m]) return -1;
//at last, r is rank, m - r is the number of freeX
return r;
}
int getMinAns(int n, int m, int r) {
int ret = INF, freeX = m - r;
for(int s = 0; s < 1 << freeX; ++s) {
if(__builtin_popcount(s) >= ret) continue;
int cnt = 0;
for(int j = 0; j < m; ++j) {
if(isFreeX[j]) {
ans[j] = s >> cnt & 1;
++cnt;
}
}
ret = min(ret, getAns(n, m, r));
}
return ret;
}

XOR高斯消元

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
int a[N][N], ans[N];
int xorGauss() {
int r = 0, c = 0;
for(; r < n && c < n; ++r, ++c) {
int p = r;
for(; p < n; ++p) if(a[p][c]) break;
if(p == n) {--r; continue;}
swap(a[p], a[r]);
for(int i = 0; i < n; ++i) {
if(i != r && a[i][c]) {
for(int j = c; j <= n; ++j)
a[i][j] ^= a[r][j];
}
}
}
for(int i = r; i < n; ++i) if(a[i][n]) return -1;
//memset(ans, 0, sizeof ans); //not necessary
for(int i = r - 1; ~i; --i)
for(int j = 0; j < n; ++j)
if(a[i][j]) {ans[j] = a[i][n]; break;}
return r;
}

Ⅲ. 题目讲解

HDU 4200 Bad Wiring

分析:
$开关问题,POJ的几个数据都比较弱,所以就选了这个$
$建立模2意义下的方程组,然后高斯消元,要求最小解,二进制枚举自由变元回代即可$
代码:

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
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
//
// Created by TaoSama on 2016-07-28
// Copyright (c) 2016 TaoSama. All rights reserved.
//
#pragma comment(linker, "/STACK:102400000,102400000")
#include <algorithm>
#include <cctype>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <iomanip>
#include <iostream>
#include <map>
#include <queue>
#include <string>
#include <set>
#include <vector>
using namespace std;
#define pr(x) cout << #x << " = " << x << " "
#define prln(x) cout << #x << " = " << x << endl
const int N = 100 + 10, INF = 0x3f3f3f3f, MOD = 2;
int n, d;
int a[N][N], ans[N];
bool isFreeX[N];
inline int inv(int x) {return 1;}
int getAns(int n, int m, int r) {
int ret = 0;
for(int i = r - 1; ~i; --i) {
for(int j = 0; j < m; ++j) {
if(!a[i][j]) continue;
ans[j] = a[i][m];
for(int k = j + 1; k < m; ++k) {
ans[j] -= a[i][k] * ans[k];
ans[j] %= MOD;
if(ans[j] < 0) ans[j] += MOD;
}
ans[j] = ans[j] * inv(a[i][j]) % MOD;
break;
}
}
for(int i = 0; i < m; ++i) ret += ans[i];
return ret;
}
int gauss(int n, int m) {
for(int i = 0; i < m; ++i) isFreeX[i] = false;
int r = 0, c = 0;
for(; r < n && c < m; ++r, ++c) {
int maxR = r; //row transform
for(int i = r + 1; i < n; ++i)
if(abs(a[i][c]) > abs(a[maxR][c])) maxR = i;
if(maxR != r) swap(a[maxR], a[r]);
if(!a[r][c]) { --r; isFreeX[c] = true; continue;}
//eliminate coefficient
for(int i = r + 1; i < n; ++i) {
if(a[i][c]) {
int delta = a[i][c] * inv(a[r][c]);
for(int j = c; j <= m; ++j) {
a[i][j] -= delta * a[r][j];
a[i][j] %= MOD;
if(a[i][j] < 0) a[i][j] += MOD;
}
}
}
}
for(int i = r; i < n; i++) if(a[i][m]) return -1;
//at last, r is rank, m - r is the number of freeX
return r;
}
int getMinAns(int n, int m, int r) {
int ret = INF, freeX = m - r;
for(int s = 0; s < 1 << freeX; ++s) {
if(__builtin_popcount(s) >= ret) continue;
int cnt = 0;
for(int j = 0; j < m; ++j) {
if(isFreeX[j]) {
ans[j] = s >> cnt & 1;
++cnt;
}
}
ret = min(ret, getAns(n, m, r));
}
return ret;
}
int main() {
#ifdef LOCAL
freopen("C:\\Users\\TaoSama\\Desktop\\in.txt", "r", stdin);
// freopen("C:\\Users\\TaoSama\\Desktop\\out.txt","w",stdout);
#endif
ios_base::sync_with_stdio(0);
clock_t _ = clock();
int t; scanf("%d", &t);
while(t--) {
scanf("%d%d", &n, &d);
memset(a, 0, sizeof a);
for(int i = 0; i < n; ++i) {
int v; scanf("%d", &v);
a[i][n] = -v + MOD;
for(int j = max(0, i - d); j <= min(n - 1, i + d); ++j)
a[i][j] = 1;
}
int r = gauss(n, n);
if(~r) printf("%d\n", getMinAns(n, n, r));
else puts("impossible");
}
#ifdef LOCAL
printf("\nTime cost: %.2fs\n", 1.0 * (clock() - _) / CLOCKS_PER_SEC);
#endif
return 0;
}

SGU 275 To xor or not to xor

分析:
$求N个数的子集的最大异或和,高斯消元后求出的是极大线性无关组即线性基$
$线性基异或起来显然就是最大异或和了$
代码:

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
//
// Created by TaoSama on 2016-07-28
// Copyright (c) 2016 TaoSama. All rights reserved.
//
#pragma comment(linker, "/STACK:102400000,102400000")
#include <algorithm>
#include <cctype>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <iomanip>
#include <iostream>
#include <map>
#include <queue>
#include <string>
#include <set>
#include <vector>
using namespace std;
#define pr(x) cout << #x << " = " << x << " "
#define prln(x) cout << #x << " = " << x << endl
const int N = 1e5 + 10, INF = 0x3f3f3f3f, MOD = 1e9 + 7;
typedef long long LL;
int n;
LL a[N];
int xorGauss() {
int r = 0, c = 63;
for(; r < n && ~c; ++r, --c) {
int p = r;
for(; p < n; ++p) if(a[p] >> c & 1) break;
if(p == n) {--r; continue;}
swap(a[p], a[r]);
for(int i = 0; i < n; ++i)
if(i != r && a[i] >> c & 1) a[i] ^= a[r];
}
return r;
}
int main() {
#ifdef LOCAL
freopen("C:\\Users\\TaoSama\\Desktop\\in.txt", "r", stdin);
// freopen("C:\\Users\\TaoSama\\Desktop\\out.txt","w",stdout);
#endif
ios_base::sync_with_stdio(0);
clock_t _ = clock();
scanf("%d", &n);
for(int i = 0; i < n; ++i) scanf("%lld", a + i);
int r = xorGauss();
LL ans = 0;
for(int i = 0; i < r; ++i) ans ^= a[i];
printf("%lld\n", ans);
#ifdef LOCAL
printf("\nTime cost: %.2fs\n", 1.0 * (clock() - _) / CLOCKS_PER_SEC);
#endif
return 0;
}

HDU 2262 Where is the canteen

分析:
$高斯消元求解期望dp$
$f[x][y]:=(x, y)到终点的期望次数,显然终点的期望是0,即f[ex][ey]=0$
$期望倒推,f[x][y]=(\sum p_i\times f[nx][ny]) + 1$
$bfs求出每个点的可达性,之后对于每个点建立期望方程,然后高斯消元求解即可$
代码:

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
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
//
// Created by TaoSama on 2016-07-28
// Copyright (c) 2016 TaoSama. All rights reserved.
//
#pragma comment(linker, "/STACK:102400000,102400000")
#include <algorithm>
#include <cctype>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <iomanip>
#include <iostream>
#include <map>
#include <queue>
#include <string>
#include <set>
#include <vector>
using namespace std;
#define pr(x) cout << #x << " = " << x << " "
#define prln(x) cout << #x << " = " << x << endl
const int N = 15 * 15 + 10, INF = 0x3f3f3f3f, MOD = 1e9 + 7;
const double EPS = 1e-8;
int sgn(double x) {
return x < -EPS ? -1 : x > EPS;
}
int n, m;
char s[N][N];
double a[N][N], ans[N];
bool isFreeX[N];
void getAns(int n, int m, int r) {
// memset(ans, 0, sizeof ans); //not necessary
for(int i = r - 1; ~i; --i) {
for(int j = 0; j < m; ++j) {
if(!sgn(a[i][j])) continue;
ans[j] = a[i][m];
for(int k = j + 1; k < m; ++k)
ans[j] -= a[i][k] * ans[k];
ans[j] = ans[j] / a[i][j];
break;
}
}
}
int gauss(int n, int m) {
for(int i = 0; i < m; ++i) isFreeX[i] = false;
int r = 0, c = 0;
for(; r < n && c < m; ++r, ++c) {
int maxR = r; //row transform
for(int i = r + 1; i < n; ++i)
if(abs(a[i][c]) > abs(a[maxR][c])) maxR = i;
if(maxR != r) swap(a[maxR], a[r]);
if(!sgn(a[r][c])) { --r; isFreeX[c] = true; continue;}
//eliminate coefficient
for(int i = r + 1; i < n; ++i) {
if(sgn(a[i][c])) {
double delta = a[i][c] / a[r][c];
for(int j = c; j <= m; ++j)
a[i][j] -= delta * a[r][j];
}
}
}
for(int i = r; i < n; i++) if(sgn(a[i][m])) return -1;
getAns(n, m, r);
//at last, r is rank, m - r is the number of freeX
return r;
}
inline int ID(int x, int y) {
return x * m + y;
}
bool can[N][N];
const int d[][2] = { -1, 0, 0, -1, 1, 0, 0, 1};
bool bfs(int sx, int sy) {
memset(can, 0, sizeof can);
queue<pair<int, int> > q;
q.push({sx, sy});
can[sx][sy] = true;
bool ok = false;
while(q.size()) {
auto u = q.front(); q.pop();
int x, y; tie(x, y) = u;
if(s[x][y] == '$') ok = true;
for(int i = 0; i < 4; ++i) {
int nx = x + d[i][0], ny = y + d[i][1];
if(nx < 0 || nx >= n || ny < 0 || ny >= m) continue;
if(s[nx][ny] == '#' || can[nx][ny]) continue;
can[nx][ny] = true;
q.push({nx, ny});
}
}
return ok;
}
int main() {
#ifdef LOCAL
freopen("C:\\Users\\TaoSama\\Desktop\\in.txt", "r", stdin);
// freopen("C:\\Users\\TaoSama\\Desktop\\out.txt","w",stdout);
#endif
ios_base::sync_with_stdio(0);
clock_t _ = clock();
while(scanf("%d%d", &n, &m) == 2) {
int sx, sy; sx = sy = -1;
for(int i = 0; i < n; ++i) {
scanf("%s", s[i]);
if(~sx) continue;
for(int j = 0; j < m; ++j)
if(s[i][j] == '@') sx = i, sy = j;
}
if(!bfs(sx, sy)) {puts("-1"); continue;}
//E_i = sum(E_j)/cnt + 1 => sum(E_j) - cnt*E_i = -cnt
memset(a, 0, sizeof a);
for(int i = 0; i < n; ++i) {
for(int j = 0; j < m; ++j) {
if(s[i][j] == '#') continue; //obstacle
if(s[i][j] == '$') { //destination
a[ID(i, j)][ID(i, j)] = 1;
a[ID(i, j)][n * m] = 0; //E_exit = 0;
continue;
}
int cnt = 0;
for(int k = 0; k < 4; ++k) {
int x = i + d[k][0], y = j + d[k][1];
if(x < 0 || x >= n || y < 0 || y >= m || !can[x][y]) continue;
++cnt;
a[ID(i, j)][ID(x, y)] = 1;
}
a[ID(i, j)][ID(i, j)] = -cnt;
a[ID(i, j)][n * m] = -cnt;
}
}
// for(int i = 0; i < n * m; ++i)
// for(int j = 0; j <= n * m; ++j)
// printf("%.3f%c", a[i][j], " \n"[j == n * m]);
gauss(n * m, n * m);
double result = ans[ID(sx, sy)];
printf("%.6f\n", result);
}
#ifdef LOCAL
printf("\nTime cost: %.2fs\n", 1.0 * (clock() - _) / CLOCKS_PER_SEC);
#endif
return 0;
}


1. 除非注明,本博文即为原创,转载请注明链接地址
2. 本博文只代表博主当时的观点或结论,请不要恶意攻击
3. 如果本文帮到了您,不妨点一下 下面分享到 按钮,让更多的人看到