基础线性代数大记（二）三道高消题
笔者的睿智前言
\[我 \text{ TMD } 不想搬砖
\]
所以……
我要搬砖
这傻吊游戏， 不玩也罢
矩阵逆
将矩阵 \(A_{nn}\) 消元成上三角矩阵后， 如果主元的个数为 n， 那么该矩阵可逆。主元不为 n 的 n 阶向量组是线性相关的。
矩阵的初等变换可以通过矩阵乘法实现， 所以通过初等变换将一个矩阵变成单位矩阵， 等于将这个矩阵乘上了这个矩阵的逆矩阵， 记录下初等变换的所有细节， 对一个单位矩阵做同一操作， 就可以得到原来矩阵的逆矩阵了。
至于矩阵是否可逆， 可以用高斯消元判。
写的模板是洛谷的， 我天生自带大常数死得很惨， 开 O2 优化才能过。
#include<bits/stdc++.h>
typedef long long lo;
using namespace std;
const int N = 405;
const int mod = 1000000007;
int n;
lo poW(lo a, int b) {
lo res = 1ll;
for(;b;b>>=1,a=(a*a)%mod)if(b&1)res=(res*a)%mod;
return res%mod;
}
struct mat{
lo a[N][N];
mat(){memset(a,0,sizeof a);}
void Swap(int i,int j) {
for(int k=1;k<=n;++k) swap(a[i][k], a[j][k]);
}
void Mul(int x, lo k) {
for(int i=1;i<=n;++i) {
a[x][i] *= k;
a[x][i] %= mod;
}
}
void Ins(int x,int y,lo k) {
for(int i=1;i<=n;++i) {
a[x][i] += a[y][i]*k%mod;
a[x][i] %= mod;
}
}
void print() {
for(int i=1;i<=n;++i,putchar('\n'))
for(int j=1;j<=n;++j) cout << (a[i][j]+mod)%mod<< ' ';
putchar('\n');
}
} A,B;
signed main() {
scanf("%d", &n);
for(int i=1;i<=n;++i) {
for(int j=1;j<=n;++j) {
scanf("%lld", &A.a[i][j]);
A.a[i][j] %= mod;
}
}
for(int i=1;i<=n;++i) B.a[i][i]=1;
for(int i=1;i<=n;++i) {
if(A.a[i][i]%mod==0)
for(int j=i;j<=n;++j) if(A.a[j][i]%mod) {
B.Swap(i,j);
A.Swap(i,j);
break;
}
if(A.a[i][i]%mod==0) {
puts("No Solution");return 0;
}
A.a[i][i] = (A.a[i][i]%mod+mod)%mod;
B.Mul(i,poW(A.a[i][i],mod-2));
A.Mul(i,poW(A.a[i][i],mod-2));
for(int j=i+1;j<=n;++j) {
B.Ins(j, i, -A.a[j][i]);
A.Ins(j, i, -A.a[j][i]);
}
}
for(int i=n;i>=1;--i) {
for(int j=i-1;j>=1;--j) {
B.Ins(j,i,-A.a[j][i]);
A.Ins(j,i,-A.a[j][i]);
}
}
B.print();
return 0;
}
Broken Robot
列方程， 高消即可。
注意特判 m=1 的情况， 这种数据不能被列方程解决。
#include<bits/stdc++.h>
using namespace std;
const int N = 1003;
int n,m,x,y;
double a[N][N];
void Gauss()
{
for(int i=1;i<=m;++i) {
//i,i i,i+1
double tmp = a[i][i];
a[i][i+1]/=tmp, a[i][i]/=tmp;
if(i==m) break;
a[i][m+1]/=tmp;
tmp = a[i+1][i];
a[i+1][i]-=tmp*a[i][i], a[i+1][i+1]-=tmp*a[i][i+1], a[i+1][m+1]-=tmp*a[i][m+1];
}
for(int i=m;i>1;--i) {
a[i-1][m+1] -= a[i-1][i]*a[i][m+1];
}
}
int main()
{
scanf("%d%d%d%d", &n,&m,&x,&y);
if(m==1) {
printf("%.4lf", (double)2 * (n - 1));
return 0;
}
n -= (x-1);
// 2/3 f[i,1] - 1/3 f[i,2] = 1/3 f[i+1,1] + 1
//-1/4 f[i,j-1] + 3/4 f[i,j] - 1/4 f[i,j+1] = 1/4 f[i+1,j] + 1
// f[i,m] = 1/3 f[i,m] + 1/3 f[i,m-1] + 1/3 f[i+1,m] + 1
for(int T=1;T<n;++T) {
a[1][m+1]=a[1][m+1]/3.0+1, a[m][m+1]=a[m][m+1]/3.0+1;
for(int j=2;j<m;++j) a[j][m+1] = a[j][m+1]/4.0 + 1;
a[1][1]=2.0/3.0, a[1][2]=-1.0/3.0;
a[m][m]=2.0/3.0, a[m][m-1]=-1.0/3.0;
for(int j=2;j<m;++j) {
a[j][j] = 3.0/4.0;
a[j][j-1] = a[j][j+1] = -1.0/4.0;
}
Gauss();
}
printf("%.4lf", a[y][m+1]);
return 0;
}
//Author:XuHt
#include <cstdio>
#include <iostream>
using namespace std;
const int N = 1006;
int n, m, x, y;
double f[N][N], d[N][N];
void work() {
for (int i = 1; i <= m; i++) {
double w = 1.0 / d[i][i];
d[i][i] *= w;
d[i][m+1] *= w;
if (i == m) break;
d[i][i+1] *= w;
w = d[i+1][i] / d[i][i];
d[i+1][i] -= w * d[i][i];
d[i+1][i+1] -= w * d[i][i+1];
d[i+1][m+1] -= w * d[i][m+1];
}
for (int i = m - 1; i; i--)
d[i][m+1] -= d[i+1][m+1] * d[i][i+1];
}
int main() {
cin >> n >> m >> x >> y;
for (int i = n - 1; i >= x; i--) {
d[1][1] = d[m][m] = -2 / 3.0;
d[1][2] = d[m][m-1] = 1 / 3.0;
for (int j = 2; j < m; j++) {
d[j][m+1] = -f[i+1][j] / 4.0 - 1;
d[j][j] = -3 / 4.0;
d[j][j-1] = d[j][j+1] = 1 / 4.0;
}
if (m == 1) d[1][1] = -1 / 2.0;
d[1][m+1] = -f[i+1][1] / 3.0 - 1;
d[m][m+1] = -f[i+1][m] / 3.0 - 1;
if (m == 1) d[m][m+1] = -f[i+1][m] / 2.0 - 1;
work();
for (int j = 1; j <= m; j++) f[i][j] = d[j][m+1];
}
printf("%.10f\n", f[x][y]);
return 0;
}
HNOI2011 XOR和路径
根据期望的线性性， 求出每一位的期望然后相加。
\(f(i,x)\) 表示第 \(i\) 位， x 到 n 路径 Xor 值不为 0 的概率。
\[f(i,x) = \frac{1}{k(x)}\sum\limits_{s\mid w(x,s){\rm xor}2^i=1}(1-f(i,s)) + \frac{1}{k(x)}\sum\limits_{s\mid w(x,s){\rm xor}2^i=0}f(i,s)
\]
#include<bits/stdc++.h>
using namespace std;
const int N = 103;
const int M = 10003;
int n,m;
int ct, hd[N], nt[M<<1], vr[M<<1], w[M<<1];
void ad(int x,int y,int z) {
vr[++ct]=y, w[ct]=z, nt[ct]=hd[x], hd[x]=ct;
}
struct mat{
double a[N][N];
mat(){memset(a,0,sizeof a);}
void Swap(int x,int y) {
for(int i=1;i<=n+1;++i) swap(a[x][i],a[y][i]);
}
void Mul(int x,double k) {
for(int i=1;i<=n+1;++i) a[x][i] *= k;
}
void Ins(int x,int y,double k) {
for(int i=1;i<=n+1;++i) a[x][i] += a[y][i]*k;
}
};
double sol(int t) {
mat A;
for(int x=1;x<=n;++x) {
A.a[x][x] = 1;
if(x==n) continue;
int sum=0;
for(int i=hd[x];i;i=nt[i]) {
++sum;
}
for(int i=hd[x];i;i=nt[i]) {
int y = vr[i];
if((w[i]>>t)&1) {
A.a[x][n+1] += 1.0/sum;
A.a[x][y] += 1.0/sum;
} else {
A.a[x][y] += -1.0/sum;
}
}
}
for(int i=1;i<=n;++i) {
int mx=i; for(int j=i+1;j<=n;++j)if(fabs(A.a[j][i])>fabs(A.a[mx][i]))mx=j;
if(mx^i) A.Swap(i,mx);
A.Mul(i,1/A.a[i][i]);
for(int j=i+1;j<=n;++j) {
A.Ins(j,i,-A.a[j][i]);
}
}
for(int i=n;i>=1;--i) {
double tmp = A.a[i][i];
A.a[i][i]/=tmp, A.a[i][n+1]/=tmp;
for(int j=i-1;j>=1;--j) A.Ins(j,i,-A.a[j][i]);
}
return A.a[1][n+1];
}
int main()
{
cin>>n>>m;
for(int i=0;i<m;++i) {
int u,v,w; scanf("%d%d%d",&u,&v,&w);
ad(u,v,w); if(u^v)ad(v,u,w);
}
double Ans = 0;
for(int T=0;T<=30;++T) Ans += sol(T)*(1<<T);
printf("%.3lf", Ans);
return 0;
}
