一、开头
神犇MX:大家好,我是神犇MX,每次都AK,喜欢d人。
xyz32768:大家好,我是蒟蒻xyz32768,每次都爆0,总是被d。
神犇MX:考下你,给定一个
n
n
个节点条边的无向连通图,每到一个节点就会等概率随机地选择下一个与当前的点有边相连的节点走,求从
1
1
走到的期望步数。
xyz32768:呃,
3
3
个节点,有边,有
12
1
2
的概率走
2
2
步,的概率走
4
4
步,的概率走
6
6
步,的概率走
8
8
步,…,期望值为。
神犇MX:我说你xyz32768你真是菜啊,高斯消元一下不就行了!
xyz32768:什么?高斯消元?
二、概述
高斯消元是解 n n 元一次方程组的算法。高斯消元的基本思想就是加减消元法。
三、预备
先将个
n
n
元一次方程存进一个的矩阵
a
a
,对于任意的表示第
i
i
个方程内,第个未知数的系数,
a[i][n+1]
a
[
i
]
[
n
+
1
]
表示第
i
i
个方程等号右边的数值。举个例子:
其中 X1,X2,X3 X 1 , X 2 , X 3 是未知数,那么储存的矩阵 a a 就为:
四、求解
分
n
n
步。每一步的操作都是把变为
1
1
,矩阵第列的其他元素变为
0
0
(即将第个方程的
Xi
X
i
系数变为
1
1
,其他方程的系数变为
0
0
)。第步的做法为:
1、将第
i
i
个方程的等式两边都除以,使第
i
i
个方程的系数为
1
1
,也就是说,把矩阵第行的所有元素除以
a[i][i]
a
[
i
]
[
i
]
。
2、消元,也就是把其他方程的
Xi
X
i
系数变为
0
0
。这时候就可以加减消元了,也就是说:
对于任何一个将
a[j][k]
a
[
j
]
[
k
]
减去
a[i][k]∗a[j][i]
a
[
i
]
[
k
]
∗
a
[
j
]
[
i
]
。
实现时,经常要在第
i
i
个到第个方程中,使用
Xi
X
i
系数绝对值最大的方程作为第
i
i
个方程(将第个元素绝对值最大的行与第
i
i
行交换),以避免浮点数计算的精度问题。
最后步都结束后,所有的
a[i][i]
a
[
i
]
[
i
]
(对角线上的元素)都会变为
1
1
,除了对角线上的元素和第列的元素之外,其他元素都为
0
0
,这时候方程的解就是第列的元素。
特殊情况:如果矩阵的某一行出现
0
0
0
0
0
0
的情况,则方程有无穷多个解,如果出现
0
0
...
.
.
.
b
b
(),则方程无解。
五、实例
模拟下求解上面给出的方程:
初始状态:
将 X1 X 1 的绝对值最大的方程置为第 1 1 行:
第 1 1 个方程的系数化为 1 1 :
将第 2 2 行消元,即用减去 2X1−4013X2+3013X3=4413 2 X 1 − 40 13 X 2 + 30 13 X 3 = 44 13 ,矩阵变为:
将第 3 3 行消元,即
将 X2 X 2 的绝对值最大的方程置为第 2 2 行并将系数化为:
将第 1 1 行和第行消元:
第 3 3 个方程系数化为:
将 1 1 、行消元:
即 X1=−1,X2=2,X3=5 X 1 = − 1 , X 2 = 2 , X 3 = 5 。
六、应用
高斯消元在竞赛中的部分用途如:
1、求解期望方程
求解期望方程等。举最简单的一个例子:给定一个无向连通图,每到一个节点就会等概率随机地选择下一个与当前节点有边相连的节点走,求从
1
1
走到的期望步数。
这时,就可以设
Xi
X
i
为
i
i
到的期望步数,其中
Xn=0
X
n
=
0
,那么容易得到,
Xi=1di∑(1+Xj)
X
i
=
1
d
i
∑
(
1
+
X
j
)
(
di
d
i
为点
i
i
的度数,为与
i
i
相邻的节点)
移项后就可以得到个方程式。高斯消元后就能求出每个
Xi
X
i
了。
2、求解异或方程(补充中)
3、求矩阵的秩(补充中)
4、求解方阵的行列式
https://blog.csdn.net/xyz32768/article/details/81413569
5、求解可逆方阵的逆矩阵(补充中)
七、代码
高斯消元模板:
void Gauss() {
int i, j, k;
for (i = 1; i <= n; i++) {
int u = i;
for (j = i + 1; j <= n; j++)
if (fabs(a[j][i]) > fabs(a[u][i])) u = j;
if (u != i) for (j = i; j <= n + 1; j++)
swap(a[i][j], a[u][j]);
double tmp = a[i][i];
for (j = i; j <= n + 1; j++) a[i][j] /= tmp;
for (j = 1; j <= n; j++) if (j != i) {
tmp = a[j][i];
for (k = i; k <= n + 1; k++)
a[j][k] -= a[i][k] * tmp;
}
}
}
八、题目
[BZOJ1013][JSOI2008]球形空间产生器sphere:
http://www.lydsy.com/JudgeOnline/problem.php?id=1013
[BZOJ2337][HNOI2011]XOR和路径:
http://www.lydsy.com/JudgeOnline/problem.php?id=2337
[BZOJ2707][SDOI2012]走迷宫:
http://www.lydsy.com/JudgeOnline/problem.php?id=2707
[BZOJ3143][HNOI2013]游走:
http://www.lydsy.com/JudgeOnline/problem.php?id=3143
[BZOJ4820][SDOI2017]硬币游戏:
http://www.lydsy.com/JudgeOnline/problem.php?id=4820