矩阵求逆
设逆矩阵为$P$,该矩阵为$A$,单位矩阵为$E$
则有$P*A=E$ $P*E=P$
因为做初等行变换等价于被对应的初等矩阵左乘。
在$A$化为$E$过程中,对$E$做相同操作,就可以得到$P$。
初始化另一个矩阵为单位矩阵,
将本矩阵用高斯消元尝试消为单位矩阵。
注意该过程不全与高斯消元相同,消为上三角后,应当有回带的过程。
在消元过程中同时对另一个矩阵做相同操作,最后得到的即是逆矩阵。
如果消元失败,即该矩阵不存在逆矩阵
复杂度$O(n^3)$。
1 #include<iostream> 2 #include<cstdio> 3 using namespace std; 4 const int mod=1e9+7; 5 const int N=405; 6 int n,a[N][2*N]; 7 int qpow(int base,int k) 8 { 9 int ans=1; 10 while(k) 11 { 12 if(k&1) ans=1ll*ans*base%mod; 13 base=1ll*base*base%mod; 14 k>>=1; 15 } 16 return ans; 17 } 18 int main() 19 { 20 scanf("%d",&n); 21 for(int i=1;i<=n;++i) 22 for(int j=1;j<=n;++j) scanf("%d",&a[i][j]); 23 for(int i=1;i<=n;++i) a[i][i+n]=1; 24 for(int i=1;i<=n;++i) 25 { 26 int id=-1; 27 for(int k=i;k<=n;++k) 28 if(a[k][i]) 29 { 30 id=k; 31 break; 32 } 33 if(id==-1) 34 { 35 puts("No Solution"); 36 return 0; 37 } 38 for(int k=i;k<=2*n;++k) swap(a[i][k],a[id][k]); 39 int inv=qpow(a[i][i],mod-2); 40 for(int k=i;k<=2*n;++k) a[i][k]=1ll*a[i][k]*inv%mod; 41 for(int j=i+1;j<=n;++j) 42 { 43 int r=a[j][i]; 44 for(int k=i;k<=2*n;++k) a[j][k]=(a[j][k]-1ll*r*a[i][k]%mod+mod)%mod; 45 } 46 } 47 for(int i=n;i;--i) 48 for(int j=i-1;j;--j) 49 { 50 int r=a[j][i]; 51 for(int k=i;k<=2*n;++k) a[j][k]=(a[j][k]-1ll*r*a[i][k]%mod+mod)%mod; 52 } 53 for(int i=1;i<=n;++i) 54 { 55 for(int j=1;j<=n;++j) printf("%d ",a[i][j+n]); 56 puts(""); 57 } 58 return 0; 59 }
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cmath> 5 using namespace std; 6 const int N=155; 7 int n,m,hp; 8 int w[N],du[N],road[N][N]; 9 double a[N][2*N]; 10 double b[N],dp[10010][N]; 11 int main() 12 { 13 scanf("%d%d%d",&n,&m,&hp); 14 for(int i=1;i<=n;++i) scanf("%d",&w[i]); 15 for(int i=1,x,y;i<=m;++i) 16 { 17 scanf("%d%d",&x,&y); 18 if(x==y) 19 { 20 if(x==n) continue; 21 if(!w[x]) a[x][y]-=1; 22 ++road[x][x]; 23 ++du[x]; 24 continue; 25 } 26 if(!w[x]&&y!=n) a[x][y]-=1; 27 if(!w[y]&&x!=n) a[y][x]-=1; 28 ++road[x][y]; ++road[y][x]; 29 ++du[x]; ++du[y]; 30 } 31 for(int i=1;i<=n;++i) 32 for(int j=1;j<=n;++j) 33 a[i][j]/=du[j]; 34 for(int i=1;i<=n;++i) a[i][i]+=1; 35 for(int i=1;i<=n;++i) a[i][i+n]=1; 36 for(int i=1;i<=n;++i) 37 { 38 int id=0; 39 for(int j=i;j<=n;++j) 40 if(!id||abs(a[j][i])>abs(a[id][i])) id=j; 41 for(int j=i;j<=2*n;++j) swap(a[id][j],a[i][j]); 42 double r=1/a[i][i]; 43 for(int j=i;j<=2*n;++j) a[i][j]*=r; 44 for(int j=i+1;j<=n;++j) 45 { 46 double r=a[j][i]; 47 for(int k=i;k<=2*n;++k) a[j][k]-=r*a[i][k]; 48 } 49 } 50 for(int i=n;i;--i) 51 { 52 for(int j=i-1;j;--j) 53 { 54 double r=a[j][i]; 55 for(int k=i;k<=2*n;++k) a[j][k]-=r*a[i][k]; 56 } 57 } 58 double ans=0; 59 for(int i=hp;i;--i) 60 { 61 memset(b,0,sizeof(b)); 62 if(i==hp) b[1]=1; 63 for(int j=1;j<=n;++j) 64 for(int k=1;k<n;++k) 65 if(w[j]+i<=hp&&w[j]) b[j]+=dp[w[j]+i][k]*road[k][j]/du[k]; 66 for(int j=1;j<=n;++j) 67 for(int k=1;k<=n;++k) 68 dp[i][j]+=a[j][k+n]*b[k]; 69 ans+=dp[i][n]; 70 } 71 printf("%.8lf\n",ans); 72 return 0; 73 }