矩阵树杂题
哥哥要是这般态度，倒不如直接不理我的好，倒显得我无理取闹了些
\(~~~~\) 好矩阵树，赏我吧。
「2017 山东一轮集训 Day5」苹果树
题意
\(~~~~\) 把 \(n\) 个带权的点连成一棵树，认为 \(val_i=-1\) 的点是坏的，并认为一棵树的权值为所有大小 \(\geq2\) 的非坏点连通块的点权和。求有多少树的点权 $\leq $ 定值 \(lim\).
\(~~~~\) \(1\leq n\leq 40,1\leq lim\leq 10^9\)。
题解
\(~~~~\) 把这题分成两个部分看：第一部分很容易想到应求好点中选若干个使权值 \(\leq lim\) 的方案数。
\(~~~~\) 具体解决的话可以用 \(\text{meet in middle}\).也就是分成两部分搜索得到各自的解，然后卷起来得到最终解，但注意到出现的权值可能很大，所以我们在用到选取 \(k\) 个时再枚举两边各自选的数量然后双指针求解。
\(~~~~\) 第二部分就是用矩阵树求解答案，考虑使已经选出的点之间连边成树。假设我们选 \(k\) 个点计入权值贡献，那所有点可以分为三部分：\(\text{1.}\)好点且选入了 \(k\) 个点 \(\text{2.}\)好点但没选入 \(k\) 个点 \(\text{3.}\)坏点。然后连上 \(1-1~~~1-3~~~2-3~~~3-3\) 的边，但注意到这样只会保证是至多选入 \(k\) 个点，因为 \(1\) 类点可能也会成为大小为 \(1\) 的连通块，所以容斥掉更小的答案，系数是 \(k\) 个点选出 \(i\) 个，即 \(C_{k}^i\)。
代码
查看代码
#include <cstdio>
#include <vector>
#include <cstring>
#include <iostream>
#include <algorithm>
#define ll long long
using namespace std;
ll n,lim,Sta=1;
ll arr[105],f[105];
ll G[105][105];
ll C[105][105],Ans[105];
const ll MOD=1e9+7;
vector<ll>V[2][105];
void dfs(ll s,ll t,ll cnt,ll val,ll op)
{
if(s==t+1)
{
if(val<=lim) V[op][cnt].push_back(val);
return;
}
dfs(s+1,t,cnt,val,op); dfs(s+1,t,cnt+1,val+arr[s],op);
}
ll Merge(ll x,ll y)
{
ll ret=0;
for(ll i=0,j=(ll)V[1][y].size()-1;i<(ll)V[0][x].size();i++)
{
while(j>=0&&V[0][x][i]+V[1][y][j]>lim) j--;
ret+=j+1;
}
return ret;
}
ll qpow(ll a,ll b)
{
ll ret=1;
while(b)
{
if(b&1) ret=ret*a%MOD;
b>>=1;a=a*a%MOD;
}
return ret;
}
ll Abs(ll x){return x>0?x:-x;}
ll Gauss(ll N)
{
ll ret=1,sig=1;
for(int i=1;i<=N;i++)
{
for(int j=i+1;j<=N;j++)
{
while(G[j][i])
{
ll div=G[i][i]*qpow(G[j][i],MOD-2)%MOD;
for(int k=i;k<=N;k++) G[i][k]=(G[i][k]-G[j][k]*div%MOD+MOD)%MOD;
swap(G[i],G[j]);
sig*=-1;
}
}
ret=ret*G[i][i]%MOD;
}
return (ret*sig+MOD)%MOD;
}
bool cmp(ll x,ll y){return x>y;}
int main() {
scanf("%lld %lld",&n,&lim); Sta=n;
for(ll i=1;i<=n;i++) scanf("%lld",&arr[i]);
sort(arr+1,arr+1+n,cmp);
for(ll i=Sta;i>=1&&arr[i]==-1;i--,Sta=i);
ll mid=Sta>>1;
dfs(1,mid,0,0,0); dfs(mid+1,Sta,0,0,1);
for(ll i=1;i<=mid;i++) sort(V[0][i].begin(),V[0][i].end());
for(ll i=1;i<=Sta-mid;i++) sort(V[1][i].begin(),V[1][i].end());
for(ll i=0;i<=n;i++)
{
C[i][0]=C[i][i]=1;
for(ll j=1;j<i;j++) C[i][j]=(C[i-1][j-1]+C[i-1][j])%MOD;
}
ll FAns=0;
for(ll k=0;k<=n;k++)
{
ll Sum=0;
for(ll i=0;i<=k;i++) Sum=(Sum+Merge(i,k-i))%MOD;
memset(G,0,sizeof(G));
for(ll i=1;i<=k;i++)
{
for(ll j=1;j<=k;j++) G[i][j]--,G[i][i]++;
for(ll j=Sta+1;j<=n;j++) G[i][j]--,G[i][i]++;
}
for(ll i=k+1;i<=Sta;i++)
for(ll j=Sta+1;j<=n;j++) G[i][j]--,G[i][i]++;
for(ll i=Sta+1;i<=n;i++)
for(ll j=1;j<=n;j++) G[i][j]--,G[i][i]++;
Ans[k]=Gauss(n-1);
for(ll i=0;i<k;i++) Ans[k]=(Ans[k]-(C[k][i]*Ans[i]%MOD)+MOD)%MOD;
FAns=(FAns+(Ans[k]*Sum%MOD))%MOD;
}
printf("%lld",FAns);
return 0;
}
「省选联考 2020 A 卷」 作业题
题意
\(~~~~\) 定义一颗树的权值为 边权和 $\times $ 所有边权的 \(\gcd\)，求一张图所有生成树的权值和。
\(~~~~\) \(1\leq n\leq 30,1\leq m\leq \frac{n(n-1)}{2},1\leq w_i\leq 152501\) 。
题解
\(~~~~\) 首先显然要化式子，\(\gcd\) 太不优美了。而由欧拉反演：\(n=\sum_{d|n} \varphi(d)\) ，所以题目求的变成：
\(~~~~\) 套路地枚举因数 \(d\) ，变为：
\(~~~~\) 其中 \(T_d\) 表示生成树内所有边的权值均能被 \(d\) 整除的所有生成树。
\(~~~~\) 所以现在求的是对于某个因数所有生成树的边权和。
\(~~~~\) 然后又是一个套路（鉴定为：博主矩阵树不合格）对于最常见的矩阵树我们维护的实际是边权积，那我们改造成维护对 \(x^2\) 取模的多项式，考虑 \((ax+b)\times (cx+d) \equiv (ad+bc)x+bd \pmod {x^2}\) 那令 \(b=d=1\) 我们就可以维护和了。
\(~~~~\) 然后注意如果一个因数都没有出现 \(n-1\) 次那就不用跑一遍矩阵树（否则你会T成暴力分）
代码
查看代码
#include <cstdio>
#include <vector>
#include <cstring>
#include <algorithm>
#define ll long long
using namespace std;
const int MOD=998244353;
ll qpow(ll a,ll b)
{
ll ret=1;
while(b)
{
if(b&1) ret=ret*a%MOD;
b>>=1;a=a*a%MOD;
}
return ret;
}
struct Poly{
ll a,b;
Poly(){}
Poly(ll A,ll B){a=A,b=B;}
Poly operator +(const Poly ano){return Poly((a+ano.a)%MOD,(b+ano.b)%MOD);}
Poly operator -(const Poly ano){return Poly((a-ano.a+MOD)%MOD,(b-ano.b+MOD)%MOD);}
Poly operator *(const Poly ano){return Poly((1ll*a*ano.b%MOD+1ll*b*ano.a%MOD)%MOD,1ll*b*ano.b%MOD);}
Poly operator /(const Poly ano){return Poly(1ll*((1ll*a*ano.b%MOD-1ll*b*ano.a%MOD+MOD)%MOD)*qpow(1ll*ano.b*ano.b%MOD,MOD-2)%MOD,1ll*b*qpow(ano.b,MOD-2)%MOD);}
}G[35][35];
int n,m;
ll Gauss(ll N)
{
Poly ret=Poly(0,1);int sig=1;
for(int i=1;i<=N;i++)
{
if(!G[i][i].b)
{
int pos=0;
for(int j=i+1;j<=N;j++)
if(G[j][i].b){pos=j;sig*=-1;break;}
swap(G[i],G[pos]);
}
Poly Inv=Poly(0,1)/G[i][i];
for(int j=i+1;j<=N;j++)
{
Poly Div=G[j][i]*Inv;
for(int k=i;k<=N;k++) G[j][k]=G[j][k]-Div*G[i][k];
}
ret=ret*G[i][i];
}
return (ret.a*sig+MOD)%MOD;
}
int u[505],v[505],w[505];
int Solve(int t)
{
memset(G,0,sizeof(G));
for(int i=1;i<=m;i++)
{
if(w[i]%t) continue;
G[u[i]][u[i]]=G[u[i]][u[i]]+Poly(w[i],1);
G[v[i]][v[i]]=G[v[i]][v[i]]+Poly(w[i],1);
G[u[i]][v[i]]=G[u[i]][v[i]]-Poly(w[i],1);
G[v[i]][u[i]]=G[v[i]][u[i]]-Poly(w[i],1);
}
return Gauss(n-1);
};
int phi[200005],Prime[200005],cnt;
void makePrime(int N)
{
phi[1]=1;
for(int i=2;i<=N;i++)
{
if(!phi[i]) Prime[++cnt]=i,phi[i]=i-1;
for(int j=1;j<=cnt&&i*Prime[j]<=N;j++)
{
if(!(i%Prime[j])){phi[i*Prime[j]]=phi[i]*Prime[j];break;}
phi[i*Prime[j]]=phi[i]*(Prime[j]-1);
}
}
}
int buc[200005];
void Divide(int x)
{
for(int i=1;i*i<=x;i++)
{
if(x%i) continue;
buc[i]++;
if(x/i!=i) buc[x/i]++;
}
}
int main() {
makePrime(152501);
scanf("%d %d",&n,&m);
for(int i=1;i<=m;i++) scanf("%d %d %d",&u[i],&v[i],&w[i]),Divide(w[i]);
ll Ans=0;
for(int Factor=1;Factor<=152501;Factor++)
{if(buc[Factor]<n-1) continue;Ans=(Ans+1ll*phi[Factor]*Solve(Factor)%MOD)%MOD;}
printf("%lld",Ans);
return 0;
}
浙公网安备 33010602011771号