给出两幅 (n(leq 400)) 个点的无向图 (G_1 ,G_2)，对于 (G_1) 的每一颗生成树，它的权值定义为有多少条边在 (G_2) 中出现。求 (G_1) 所有生成树的权值和。
Solution
很容易想到，设 (G_1) 中每条边的权值，这条边在 (G_2) 中出现则权值为 (1)，否则权值为 (0)。
现在就真的是求所有生成树的边权和的权值和了。
然而标准的 Matrix-Tree Theorem 求的是生成树的边权积的和。
现在我们定义每条只出现在 (G_1) 中的边边权为 (1)，同时出现在 (G_1,G_2) 中的边权为 (x)，则基尔霍夫矩阵的每个元素嗾使一个多项式，记为 (B(x))。
(det B(x)) 是一个 (n-1) 次多项式 (f(x) = sum a_i x^i)，那么其中 (a_i) 就是使用了 (i) 条公共变的生成树个数。
于是答案就是 (f'(1)=sum ia_i=det B(1) cdot sum_i sum_j (B^{-1}(1))_{i,j}cdot B'(1)_{i,j})
于是用 Gauss 消元法求行列式和逆矩阵即可
#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 405;
const int mod = 998244353;
namespace mat {
int f[N][N<<1],a[N][N],n;
inline void exgcd(int a,int b,int &x,int &y) {
if(!b) {
x=1,y=0;
return;
}
exgcd(b,a%b,x,y);
int t=x;
x=y,y=t-(a/b)*y;
}
inline int inv(int a,int b) {
int x,y;
return exgcd(a,b,x,y),(x%b+b)%b;
}
int getdet() {
int det=1;
int flag=0;
for(int i=1; i<=n; i++) {
for(int j=i+1; j<=n; j++) {
int x=i,y=j;
while(a[y][i]!=0) {
int t=a[x][i]*inv(a[y][i],mod)%mod;
for(int k=i; k<=n; k++) (a[x][k]-=t*a[y][k]%mod)%=mod;
swap(x,y);
}
if(x!=i) {
for(int k=1; k<=n; k++) {
swap(a[x][k],a[i][k]);
}
flag^=1;
}
}
if(a[i][i]==0) return 0;
det=det*a[i][i]%mod;
}
if(flag) det=-det;
det%=mod; det+=mod; det%=mod;
return det;
}
int solve() {
int m=n*2;
for(int i=1;i<=n;i++) {
for(int j=1;j<=n;j++) a[i][j]=f[i][j], f[i][j+n]=0;
}
int ret= getdet();
for(int i=1; i<=n; ++i) {
f[i][n+i]=1;
}
for(int i=1; i<=n; ++i) {
for(int j=i; j<=n; j++)
if(f[j][i]) {
for(int k=1; k<=m; k++)
swap(f[i][k],f[j][k]);
break;
}
if(!f[i][i]) {
return 0;
}
int r=inv(f[i][i],mod);
for(int j=i; j<=m; ++j)
f[i][j]=f[i][j]*r%mod;
for(int j=1; j<=n; ++j)
if(j!=i) {
r=f[j][i];
for(int k=i; k<=m; ++k)
f[j][k]=(f[j][k]-r*f[i][k]%mod+mod)%mod;
}
}
return ret;
}
}
int n,b[N][N],bd[N][N];
char g1[N][N],g2[N][N];
signed main() {
cin>>n;
for(int i=1;i<=n;i++) {
cin>>g1[i]+1;
}
for(int i=1;i<=n;i++) {
cin>>g2[i]+1;
}
for(int i=1;i<=n;i++) {
for(int j=1;j<i;j++) {
if(g1[i][j]=='1') b[i][j]=-1, b[j][i]=-1, b[i][i]++, b[j][j]++;
if(g1[i][j]=='1' && g2[i][j]=='1') bd[i][j]=-1, bd[j][i]=-1, bd[i][i]++, bd[j][j]++;
}
}
for(int i=1;i<n;i++) {
for(int j=1;j<n;j++) {
mat::f[i][j]=b[i][j];
}
}
--n;
mat::n=n;
int det=mat::solve();
int ans=0;
/*for(int i=1;i<=n;i++) {
for(int j=1;j<=n;j++) cout<<b[i][j]<<" ";
cout<<endl;
}
for(int i=1;i<=n;i++) {
for(int j=1;j<=n;j++) cout<<bd[i][j]<<" ";
cout<<endl;
}
for(int i=1;i<=n;i++) {
for(int j=1;j<=n;j++) cout<<mat::f[i][j+n]<<" ";
cout<<endl;
}*/
for(int i=1;i<=n;i++) {
for(int j=1;j<=n;j++) {
ans+=mat::f[i][j+n]*bd[i][j];
ans%=mod;
ans+=mod;
ans%=mod;
}
}
//cout<<ans<<" "<<det<<endl;
cout<<((ans*det)%mod+mod)%mod;
}