高斯消元法是用来解线性方程组的一种算法
板子
for(reg i = 1;i < n;i++)
{
int r = i;
while(!a[r][i]) r++;
swap(a[r],a[i]);
double pass = a[i][i];
for(reg k = 1;k <= n;k++) a[i][k] /= pass;
for(reg j = 1;j < n;j++)
{
if(j == i) continue;
double div = a[j][i];
for(reg k = 1;k <= n;k++)
a[j][k] -= a[i][k] * div;
}
}
用途
这里主要说下用途
做很多关于期望的题时
\[dp[i]=f(dp[j])\]
\[dp[j]=f(dp[i])\]
出现了以上形式的\(dp\)方程
明显直接\(dp\)是不行的
但如果 我们能求出
\(n\)个关于\(dp[i]=f(dp[j])(i,j\in[1,n],s.t.i,j\in Z)\)
就可以高斯消元来求这些未知量
例题
#include <stack>
#include <cstdio>
#include <vector>
#include <string>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;
#define reg register int
#define isdigit(x) ('0' <= x&&x <= '9')
template<typename T>
inline T Read(T Type)
{
T x = 0,f = 1;
char a = getchar();
while(!isdigit(a)) {if(a == '-') f = -1;a = getchar();}
while(isdigit(a)) {x = (x << 1) + (x << 3) + (a ^ '0');a = getchar();}
return x * f;
}
const int MAXN = 510,MAXM = MAXN * MAXN;
struct node
{
int v,next;
}edge[MAXM << 1];
double ans,a[MAXN][MAXN],p[MAXN],edge_p[MAXM];
int cnt,n,m,head[MAXN],deg[MAXN];
inline void addedge(int u,int v)
{
deg[u]++;
edge[++cnt].v = v;
edge[cnt].next = head[u];
head[u] = cnt;
}
inline void gauss()
{
a[1][n] = -1;
for(reg i = 1;i < n;i++)
{
for(reg e = head[i];e;e = edge[e].next)
{
int v = edge[e].v;
if(v == n) continue;
a[i][v] += 1.0 / deg[v];
}
a[i][i] -= 1;
}
for(reg i = 1;i < n;i++)
{
int r = i;
while(!a[r][i]) r++;
swap(a[r],a[i]);
double pass = a[i][i];
for(reg k = 1;k <= n;k++) a[i][k] /= pass;
for(reg j = 1;j < n;j++)
{
if(j == i) continue;
double div = a[j][i];
for(reg k = 1;k <= n;k++)
a[j][k] -= a[i][k] * div;
}
}
for(reg i = 1;i < n;i++) p[i] = a[i][n] / a[i][i];
cnt = 0;
for(reg i = 1;i <= n;i++)
for(reg e = head[i];e;e = edge[e].next)
{
int v = edge[e].v;
if(v < i) continue;
edge_p[++cnt] = (i == n?0:p[i] / deg[i]) + (v == n?0:p[v] / deg[v]);
}
}
int main()
{
n = Read(1),m = Read(1);
for(reg i = 1;i <= m;i++)
{
int u = Read(1),v = Read(1);
addedge(u,v),addedge(v,u);
}
gauss();
sort(edge_p + 1,edge_p + 1 + cnt);
for(reg i = 1;i <= m;i++) ans += edge_p[i] * (m - i + 1);
printf("%.3f",ans);
return 0;
}
#include <stack>
#include <cstdio>
#include <vector>
#include <string>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;
#define reg register int
#define isdigit(x) ('0' <= x&&x <= '9')
template<typename T>
inline T Read(T Type)
{
T x = 0,f = 1;
char a = getchar();
while(!isdigit(a)) {if(a == '-') f = -1;a = getchar();}
while(isdigit(a)) {x = (x << 1) + (x << 3) + (a ^ '0');a = getchar();}
return x * f;
}
const int MAXN = 110,MAXM = 10010;
struct node
{
int v,w,next;
}edge[MAXM << 1];
double ans,a[MAXN][MAXN];
int n,m,maxt,cnt,head[MAXN],deg[MAXN];
inline void addedge(int u,int v,int w)
{
edge[++cnt].v = v;
edge[cnt].w = w;
edge[cnt].next = head[u];
head[u] = cnt;
}
inline int getlen(int x) {int res = 0;while(x) x >>= 1,res++;return res;}
inline int getx(int x,int t) {return (x >> (t - 1) & 1);}
inline double fabs(double x) {if(x < 0) return -x;return x;}
inline void gauss(int t)
{
for(reg i = 1;i <= n;i++)
for(reg j = 1;j <= n + 1;j++)
a[i][j] = 0;
for(reg i = 1;i < n;i++)
{
for(reg j = head[i];j;j = edge[j].next)
{
int v = edge[j].v,w = getx(edge[j].w,t);
if(w) a[i][v] -= 1.0,a[i][n + 1] -= 1.0;
else a[i][v] += 1;
}
a[i][i] -= deg[i];
}
a[n][n] = 1;
for(reg i = 1;i <= n;i++)
{
int r = i;
while(!a[r][i]) r++;
swap(a[r],a[i]);
double pass = a[i][i];
for(reg j = 1;j <= n + 1;j++) a[i][j] /= pass;
for(reg j = 1;j <= n;j++)
{
if(j == i) continue;
double pass = a[j][i];
for(reg k = 1;k <= n + 1;k++)
a[j][k] -= a[i][k] * pass;
}
}
ans += a[1][n + 1] / a[1][1] * (1 << (t - 1));
}
int main()
{
n = Read(1),m = Read(1);
for(reg i = 1;i <= m;i++)
{
int u = Read(1),v = Read(1),w = Read(1);
addedge(u,v,w),deg[u]++;
maxt = max(maxt,getlen(w));
if(u != v) addedge(v,u,w),deg[v]++;
}
for(reg i = 1;i <= maxt;i++) gauss(i);
printf("%.3f",ans);
return 0;
}
#include <stack>
#include <cstdio>
#include <vector>
#include <string>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;
#define reg register int
#define isdigit(x) ('0' <= x&&x <= '9')
template<typename T>
inline T Read(T Type)
{
T x = 0,f = 1;
char a = getchar();
while(!isdigit(a)) {if(a == '-') f = -1;a = getchar();}
while(isdigit(a)) {x = (x << 1) + (x << 3) + (a ^ '0');a = getchar();}
return x * f;
}
const int MAXN = 15;
double pow[MAXN],a[MAXN][MAXN];
inline double fabs(double x) {if(x < 0) return -x;return x;}
struct Matrix
{
int n,m;
double c[MAXN][MAXN];
}T;
int main()
{
int n = Read(1);
for(reg i = 1;i <= n + 1;i++)
{
for(reg j = 1;j <= n;j++)
{
scanf("%lf",&a[i][j]);
pow[i] += a[i][j] * a[i][j];
if(i != 1) T.c[i - 1][j] = 2 * (a[i - 1][j] - a[i][j]);
}
if(i != 1) T.c[i - 1][n + 1] = pow[i - 1] - pow[i];
}
for(reg i = 1;i <= n;i++)
{
int r = i;
for(reg j = i + 1;j <= n;j++)
if(fabs(T.c[r][i]) > fabs(T.c[j][i])) r = j;
swap(T.c[r],T.c[i]);
double pass = T.c[i][i];
for(reg j = 1;j <= n + 1;j++) T.c[i][j] /= pass;
for(reg j = 1;j <= n;j++)
{
if(i == j) continue;
double pass = T.c[j][i];
for(reg k = 1;k <= n + 1;k++) T.c[j][k] -= T.c[i][k] * pass;
}
}
for(reg i = 1;i <= n;i++)
printf("%.3f ",T.c[i][n + 1]);
return 0;
}