首页 > 代码库 > 【高斯消元】兼 【期望dp】例题
【高斯消元】兼 【期望dp】例题
【总览】
高斯消元基本思想是将方程式的系数和常数化为矩阵,通过将矩阵通过行变换成为阶梯状(三角形),然后从小往上逐一求解。
如:$3X_1 + 2X_2 + 1X_3 = 3$
$ X_2 + 2X_3 = 1$
$2X_1 + X_3 = 0$
化为矩阵为:--->----->----->
然后就可以通过最后一行直接求出$X_3 = ...$,将其带回第二行,算出$X_2$,同理算出$X_1$。
代码很好理解:
inline void gauss(){ int i, j, k, l; for(i = 1; i <= n; i++){ l = i; for(j = i + 1; j <= n; j++) if(fabs(matrix[j][i]) > fabs(matrix[l][i])) l = j; if(l != i) for(j = i; j <= n + 1; j++) swap(matrix[i][j], matrix[l][j]); for(j = i + 1; j <= n; j++){ double tmp = matrix[j][i] / matrix[i][i]; for(k = i; k <= n + 1; k++) matrix[j][k] -= matrix[i][k] * tmp; } } for(i = n; i >= 1; i--){ double t = matrix[i][n + 1]; for(j = n; j > i; j--) t -= ans[j] * matrix[i][j]; ans[i] = t / matrix[i][i]; }}
高斯消元最常应用在 期望DP 中。下面是几道例题。
“期望dp讲解及例题”
【BZOJ1013】球形空间产生器sphere
由给出的$n + 1$个坐标,可以列出 $n$个方程,剩下的模板。
【CODE】
#include<iostream>#include<cstring>#include<string>#include<cstdio>#include<cstdlib>#include<algorithm>#include<vector>#include<cmath>using namespace std;const int N = 20;double matrix[N][N], last[N], t, ans[N];int n;inline int read(){ int i = 0, f = 1; char ch = getchar(); for(; (ch < ‘0‘ || ch > ‘9‘) && ch != ‘-‘; ch = getchar()); if(ch == ‘-‘) f = -1, ch = getchar(); for(; ch >= ‘0‘ && ch <= ‘9‘; ch = getchar()) i = (i << 3) + (i << 1) + (ch -‘0‘); return i * f;}inline void wr(int x){ if(x < 0) putchar(‘-‘), x = -x; if(x > 9) wr(x / 10); putchar(x % 10 + ‘0‘);}inline void gauss(){ int i, j, l, k; for(i = 1; i <= n; i++){ l = i; for(j = i + 1; j <= n; j++) if(fabs(matrix[j][i]) > fabs(matrix[l][i])) l = j; if(l != i) for(j = i; j <= n + 1; j++) swap(matrix[i][j], matrix[l][j]); for(j = i + 1; j <= n; j++){ double tmp = matrix[j][i] / matrix[i][i]; for(k = i; k <= n + 1; k++) matrix[j][k] -= matrix[i][k] * tmp; } } for(i = n; i >= 1; i--){ double tmp = matrix[i][n + 1]; for(j = n; j > i; j--) tmp -= ans[j] * matrix[i][j]; ans[i] = tmp / matrix[i][i]; }}int main(){ n = read(); for(int i = 1; i <= n; i++) scanf("%lf", &last[i]); for(int i = 1; i <= n; i++){ t = 0; for(int j = 1; j <= n; j++){ double tmp; scanf("%lf", &tmp); matrix[i][j] = 2 * (tmp - last[j]); t += tmp * tmp - last[j] * last[j]; last[j] = tmp; } matrix[i][n + 1] = t; } gauss(); for(int i = 1; i <= n; i++){ if(i < n) printf("%.3lf ", ans[i]); else printf("%.3lf\n", ans[i]); } return 0;}
【BZOJ3143】游走
因为要求期望的最小值,那么走的次数多的边肯定要让花费(编号)尽可能小,所以可以先求出从每个点出发次数的期望值$E_i$,那么对于一条边而言,走这条边的期望次数就是$E_i / degree[i] + E_j / degree[j]$,只要排一遍序就好。
求点的期望:$E_i = \sum (E_{son[i]} / degree[i])$
【code】
#include<iostream>#include<cstdio>#include<cstdlib>#include<cstring>#include<string>#include<cmath>#include<algorithm>#include<vector>#define eps 1e-10using namespace std;const int N = 600;double matrix[N][N], ans[N], ret;int n, m, num, degree[N];int st[500000], ed[500000];double gg[500000];inline void addEdge(const int &u, const int &v){ degree[u]++; degree[v]++; st[++num] = u, ed[num] = v;}inline int read(){ int i = 0, f = 1; char ch = getchar(); for(; (ch < ‘0‘ || ch > ‘9‘) && ch != ‘-‘; ch = getchar()); if(ch == ‘-‘) f = -1, ch = getchar(); for(; ch >= ‘0‘ && ch <= ‘9‘; ch = getchar()) i = (i << 3) + (i << 1) + (ch -‘0‘); return i * f;}inline void gauss(){ int i, j, l, k; for(i = 1; i <= n; i++){ l = i; for(j = i + 1; j <= n; j++) if(fabs(matrix[j][i]) > fabs(matrix[l][i])) l = j; if(l != i) for(j = i; j <= n + 1; j++) swap(matrix[i][j], matrix[l][j]); for(j = i + 1; j <= n; j++){ double tmp = matrix[j][i] / matrix[i][i]; for(k = i; k <= n + 1; k++) matrix[j][k] -= matrix[i][k] * tmp; } } for(i = n; i >= 1; i--){ double tmp = matrix[i][n + 1]; for(j = n; j > i; j--) tmp -= ans[j] * matrix[i][j]; ans[i] = tmp / matrix[i][i]; }}inline bool cmp (double a, double b){ return a > b;}int main(){ n = read(), m = read(); for(int i = 1; i <= m; i++){ int u = read(), v = read(); addEdge(u, v); } int i, j; for(i = 1; i <= m; i++){ matrix[st[i]][ed[i]] += 1.0 /degree[ed[i]]; matrix[ed[i]][st[i]] += 1.0 /degree[st[i]]; } for(i = 1; i <= n; i++) matrix[n][i] = 0; for(i = 1; i <= n; i++) matrix[i][i] = -1.0; matrix[1][n + 1] = -1.0; gauss(); for(i = 1; i <= m; i++) gg[i] = ans[st[i]] / degree[st[i]] + ans[ed[i]] / degree[ed[i]]; sort(gg + 1, gg + m + 1, cmp); for(i = 1; i <= m; i++) ret += gg[i] * i; printf("%.3f\n", ret); return 0;}
【bzoj2337】XOR和路径
学到了!看见求异或和$----->$按位计算:即一位一位的计算答案每一位上为1的期望值,这样就可以轻松统计出答案。
每一位都要重新构造矩阵求期望,设$a[i]$表示从$i$到$n$的路径异或和(这一位)为$1$的期望概率(总是≤$1$)
对于当前第$i + 1$位,若$(dis >> i) & 1$(这一位为1),那么要异或和为$1$,要求他从关联点异或和为$0$转移来,
同理,若这一位为$0$,要求从$1$转移来。即:$$a[i] = \sum a[son[i]](dis这一位为0) / degree[i] + \sum a[son[i]](dis这一位为1) / degree[i]$$。
【code】
#include<iostream>#include<cstring>#include<string>#include<algorithm>#include<cstdio>#include<cstdlib>#include<cmath>using namespace std;const int N = 105, M = 10005;int n, m;int ecnt, st[M << 1], ed[M << 1], len[M << 1], degree[N];double matrix[N][N], ans[N], ret;inline void addEdge(const int &u, const int &v, const int &l){ st[++ecnt] = u, ed[ecnt] = v, len[ecnt] = l; degree[u]++; if(u != v) st[++ecnt] = v, ed[ecnt] = u, len[ecnt] = l, degree[v]++;}inline void gauss(){ int i, j, k, l; for(i = 1; i <= n; i++){ l = i; for(j = i + 1; j <= n; j++) if(fabs(matrix[j][i]) > fabs(matrix[l][i])) l = j; if(l != i) for(j = i; j <= n + 1; j++) swap(matrix[i][j], matrix[l][j]); for(j = i + 1; j <= n; j++){ double tmp = matrix[j][i] / matrix[i][i]; for(k = i; k <= n + 1; k++) matrix[j][k] -= matrix[i][k] * tmp; } } for(i = n; i >= 1; i--){ double t = matrix[i][n + 1]; for(j = n; j > i; j--) t -= ans[j] * matrix[i][j]; ans[i] = t / matrix[i][i]; }}int main(){ scanf("%d%d", &n, &m); int i, j, k; for(i = 1; i <= m; i++){ int u, v, w; scanf("%d%d%d", &u, &v, &w); addEdge(u, v, w); } for(i = 0; i <= 30; i++){ memset(matrix, 0, sizeof matrix); memset(ans, 0, sizeof ans); for(j = 1; j <= n; j++) matrix[j][j] = 1; for(j = 1; j <= ecnt; j++){ int l = len[j], u = st[j], v = ed[j]; if(u == n) continue; if((l >> i) & 1){ matrix[u][v] += 1.0 / degree[u]; matrix[u][n + 1] += 1.0 / degree[u]; } else matrix[u][v] -= 1.0 / degree[u]; } gauss(); ret += ans[1] * (1 << i); } printf("%.3f\n", ret); return 0;}
【高斯消元】兼 【期望dp】例题