【BZOJ1494】【NOI2007】生成树计数（动态规划，矩阵快速幂）
【BZOJ1494】【NOI2007】生成树计数（动态规划，矩阵快速幂）
题面
Description
最近，小栋在无向连通图的生成树个数计算方面有了惊人的进展，他发现：
·n个结点的环的生成树个数为n。
·n个结点的完全图的生成树个数为n^(n-2)。这两个发现让小栋欣喜若狂，由此更加坚定了他继续计算生成树个数的
想法，他要计算出各种各样图的生成树数目。一天，小栋和同学聚会，大家围坐在一张大圆桌周围。小栋看了看，
马上想到了生成树问题。如果把每个同学看成一个结点，邻座（结点间距离为1）的同学间连一条边，就变成了一
个环。可是，小栋对环的计数已经十分娴熟且不再感兴趣。于是，小栋又把图变了一下：不仅把邻座的同学之间连
一条边，还把相隔一个座位（结点间距离为2）的同学之间也连一条边，将结点间有边直接相连的这两种情况统称
为有边相连，如图1所示。
小栋以前没有计算过这类图的生成树个数，但是，他想起了老师讲过的计算任意图的生成树个数的一种通用方法：
构造一个n×n的矩阵A={aij},其中
其中di表示结点i的度数。与图1相应的A矩阵如下所示。为了计算图1所对应的生成数的个数，只要去掉矩阵A的最
后一行和最后一列，得到一个(n-1)×(n-1)的矩阵B，计算出矩阵B的行列式的值便可得到图1的生成树的个数所以
生成树的个数为|B|=3528。小栋发现利用通用方法，因计算过于复杂而很难算出来，而且用其他方法也难以找到更
简便的公式进行计算。于是，他将图做了简化，从一个地方将圆桌断开，这样所有的同学形成了一条链，连接距离
为1和距离为2的点。例如八个点的情形如下：
这样生成树的总数就减少了很多。小栋不停的思考，一直到聚会结束，终于找到了一种快捷的方法计算出这个图的
生成树个数。可是，如果把距离为3的点也连起来，小栋就不知道如何快捷计算了。现在，请你帮助小栋计算这类
图的生成树的数目。
Input
包含两个整数k,n，由一个空格分隔。k表示要将所有距离不超过k（含k）的结点连接起来，n表示有n个结点。
Output
输出一个整数，表示生成树的个数。由于答案可能比较大，所以你 只要输出答案除65521 的余数即可。
Sample Input
3 5
Sample Output
75
HINT
题解
这是一道神仙题啊。
鉴于我自己在网上研究了很久各种题解才知道怎么做。
我还是打算好好地把这道题目从头到尾写一写。
从部分分开始吧。
\(Task1\ 60pts:k\le 5,n\le 100\)
送分的一档
直接暴力构图，矩阵树定理直接算就好了
时间复杂度\(O(n^3)\)
#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<set>
#include<map>
#include<vector>
#include<queue>
using namespace std;
#define ll long long
#define RG register
#define MOD 65521
#define MAX 111
inline int read()
{
RG int x=0,t=1;RG char ch=getchar();
while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
if(ch=='-')t=-1,ch=getchar();
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return x*t;
}
int n,K,c[MAX][MAX],ans=1;
int main()
{
K=read();n=read();
for(int i=1;i<=n;++i)
for(int j=i+1;j<=min(n,i+K);++j)
c[i][i]++,c[j][j]++,c[i][j]--,c[j][i]--;
for(int i=2;i<=n;++i)
for(int j=i+1;j<=n;++j)
while(c[j][i])
{
int t=c[i][i]/c[j][i];
for(int k=i;k<=n;++k)c[i][k]=(c[i][k]-1ll*c[j][k]*t%MOD)%MOD,swap(c[i][k],c[j][k]);
ans=-ans;
}
for(int i=2;i<=n;++i)ans=1ll*ans*c[i][i]%MOD;
printf("%d\n",ans);
return 0;
}
\(Task2\ 80pts:k\le 5,n\le 10000\)
这档分其实写出来基本就会满分了
发现\(k\)的范围十分的小，也就是每个点的边并不会很多
换而言之，如果一个点要和前面的所有点处在一个生成树中
意味着它必须和前面\(k\)个点中的至少一个处在同一个联通块中
所以，我们只需要考虑当前点和它前面\(K\)个点的状态就行了
那么，我们的状态是什么？
我们思考一下，对于生成树而言，我们比较在意的只有两个方面：联通块和边
对于相邻的\(K\)个点，单独拿出来显然是一个完全图，因此不需要考虑边的问题
所以，我们的状态和联通块有关。
现在的问题又变成了如何记录前面\(k\)个点的联通块呢？
我们可以按照顺序编号，保证任意一个编号的联通块在出现之前，所有小于它的编号都已经出现过。（最小表示法吼啊）
（不一定要这样，只需要一种能够保证所有相同的联通块方案只会被算一次就行了）
编完号之后，我们可以把前面\(k\)个点所在的联通块给压成十进制
因为\(k\le 5\)，就只需要三个二进制位，所以可以用八进制来把联通块的情况给压成十进制位。
那么，不同的联通块方案数有多少呢？
相当于把\(n\)个球放进任意数量个集合，也就是\(Bell\)数，算出来是\(52\)种
我们不关心每个联通块之间的连接方案，我们只关心如何对于两种联通块之间进行转移
所以我们可以枚举当前点和前面\(k\)个点的连边方案（最多就\(k\)条边）
然后暴力（用并查集）判断是否成环，
同时，最前面的那个点也必须和当前这\(k\)个点中的一个在同一个联通块中
这样就可以转移到另外一个联通块的情况
再用上面的最小表示法把它的编号还原出来
这样证明可以从当前状态向后面的状态进行转移。
\(dp\)要的是初始状态和转移。
显然搞出来了转移，考虑初始状态。
显然不能一开始不满\(k\)个点
所以直接从\(k\)个点开始计算
因为每个联通块之间是完全图，所以可以暴力计算联通块大小为\(x\)时的方案数
那么\(k\)个点时，某个联通块情况的方案数就是\(\prod num_{size}\)
联通块大小为\(1,2\)时，有\(1\)种方法
联通块大小为\(3\)有\(3\)种
联通块大小为\(4\)有\(16\)种
联通块大小为\(5\)有\(125\)种。
这样子，我们就可以\(O(52*2^5*K^2+n*52^2)\)转移了
这样子可以过\(80pts\)
#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<set>
#include<map>
#include<vector>
#include<queue>
using namespace std;
#define ll long long
#define RG register
#define MOD 65521
#define MAX 55
inline ll read()
{
RG ll x=0,t=1;RG char ch=getchar();
while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
if(ch=='-')t=-1,ch=getchar();
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return x*t;
}
ll n;int K,cnt;
int p[1<<20],st[MAX];
bool check(int t)//检查一个状态是否合法
{
int tmp=1<<1;//因为第一个点一定属于一号联通快，所以先把一号联通快放进去检查
for(int i=3;i<K+K+K;i+=3)
{
for(int j=1;j<((t>>i)&7);++j)//检查比当前编号小的所有编号是否都已经出现过
if(!(tmp&(1<<j)))return false;
tmp|=1<<((t>>i)&7);//将当前编号也给放进来
}
return true;
}
void dfs(int x,int t)//暴力找出所有状态，每个编号利用3个二进制位存
{
if(x==K){if(check(t))p[t]=++cnt,st[cnt]=t;return;}
for(int i=1;i<=K;++i)dfs(x+1,t|(i<<(x+x+x)));
}
int fa[MAX],a[MAX];
int getf(int x){return x==fa[x]?x:fa[x]=getf(fa[x]);}
int f[11111][MAX],g[MAX][MAX];
int main()
{
K=read();n=read();dfs(1,1);
for(int i=1;i<=cnt;++i)
{
f[K][i]=1;memset(a,0,sizeof(a));
for(int j=0;j<K;++j)++a[(st[i]>>(j*3))&7];
for(int j=1;j<=K;++j)
if(a[j]==3)f[K][i]=3;
else if(a[j]==4)f[K][i]=16;
else if(a[j]==5)f[K][i]=125;
int t=st[i];
for(int s=0;s<(1<<K);++s)//暴力枚举当前点对于前面几个点的连边状态
{
for(int j=0;j<=K;++j)fa[j]=j;
for(int j=0;j<K;++j)//利用并查集维护联通性
for(int k=j+1;k<K;++k)
if(((t>>(3*j))&7)==((t>>(3*k))&7))
fa[getf(j)]=getf(k);
bool cir=false;
for(int j=0;j<K;++j)//检查当前点的连边
if(s&(1<<j))
{
if(getf(K)==getf(j)){cir=true;break;}//出现了环
fa[getf(K)]=getf(j);
}
if(cir)continue;//连边不合法
for(int j=1;j<=K;++j)//最前面的点必须和后面的一个点联通，否则就无法联通了
if(getf(0)==getf(j)){cir=true;break;}
if(!cir)continue;
int now=0,used=0;
for(int j=0;j<K;++j)//当前存在合法的联通方案，因此当前的状态可以转移到另外的一个状态上去
if(!(now&(7<<(j*3))))//当前点不在任意一个联通块中
{
now|=++used<<(j*3);//新的联通块
for(int k=j+1;k<K;++k)//把所有在一个联通块里的点丢到状态里去
if(getf(j+1)==getf(k+1))
now|=used<<(k*3);
}
g[i][p[now]]++;
}
}
for(int i=K+1;i<=n;++i)
for(int j=1;j<=cnt;++j)
for(int k=1;k<=cnt;++k)
if(g[j][k])f[i][k]=(f[i][k]+1ll*g[j][k]*f[i-1][j])%MOD;
printf("%d\n",f[n][1]);
return 0;
}
\(AC:n\le 10^{15}\)
明摆着\(log\)算法，
发现每次转移相同，直接矩阵快速幂就行了
时间复杂度\(O(logn·52^3)\)，预处理的时间写在上面了。。
#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<set>
#include<map>
#include<vector>
#include<queue>
using namespace std;
#define ll long long
#define RG register
#define MOD 65521
#define MAX 55
ll n;int K,cnt;
int p[1<<20],st[MAX];
struct Matrix
{
int s[MAX][MAX];
void clear(){memset(s,0,sizeof(s));}
void init(){clear();for(int i=1;i<=cnt;++i)s[i][i]=1;}
}G;
Matrix operator*(Matrix a,Matrix b)
{
Matrix ret;ret.clear();
for(int i=1;i<=cnt;++i)
for(int j=1;j<=cnt;++j)
for(int k=1;k<=cnt;++k)
ret.s[i][j]=(ret.s[i][j]+1ll*a.s[i][k]*b.s[k][j])%MOD;
return ret;
}
Matrix fpow(Matrix a,ll b)
{
Matrix s;s.init();
while(b){if(b&1)s=s*a;a=a*a;b>>=1;}
return s;
}
bool check(int t)//检查一个状态是否合法
{
int tmp=1<<1;//因为第一个点一定属于一号联通快，所以先把一号联通快放进去检查
for(int i=3;i<K+K+K;i+=3)
{
for(int j=1;j<((t>>i)&7);++j)//检查比当前编号小的所有编号是否都已经出现过
if(!(tmp&(1<<j)))return false;
tmp|=1<<((t>>i)&7);//将当前编号也给放进来
}
return true;
}
void dfs(int x,int t)//暴力找出所有状态，每个编号利用3个二进制位存
{
if(x==K){if(check(t))p[t]=++cnt,st[cnt]=t;return;}
for(int i=1;i<=K;++i)dfs(x+1,t|(i<<(x+x+x)));
}
int fa[MAX],a[MAX];
int getf(int x){return x==fa[x]?x:fa[x]=getf(fa[x]);}
int f[MAX],g[MAX][MAX];
int main()
{
scanf("%d%lld",&K,&n);dfs(1,1);
for(int i=1;i<=cnt;++i)
{
f[i]=1;memset(a,0,sizeof(a));
for(int j=0;j<K;++j)++a[(st[i]>>(j*3))&7];
for(int j=1;j<=K;++j)
if(a[j]==3)f[i]=3;
else if(a[j]==4)f[i]=16;
else if(a[j]==5)f[i]=125;
int t=st[i];
for(int s=0;s<(1<<K);++s)//暴力枚举当前点对于前面几个点的连边状态
{
for(int j=0;j<=K;++j)fa[j]=j;
for(int j=0;j<K;++j)//利用并查集维护联通性
for(int k=j+1;k<K;++k)
if(((t>>(3*j))&7)==((t>>(3*k))&7))
fa[getf(j)]=getf(k);
bool cir=false;
for(int j=0;j<K;++j)//检查当前点的连边
if(s&(1<<j))
{
if(getf(K)==getf(j)){cir=true;break;}//出现了环
fa[getf(K)]=getf(j);
}
if(cir)continue;//连边不合法
for(int j=1;j<=K;++j)//最前面的点必须和后面的一个点联通，否则就无法联通了
if(getf(0)==getf(j)){cir=true;break;}
if(!cir)continue;
int now=0,used=0;
for(int j=0;j<K;++j)//当前存在合法的联通方案，因此当前的状态可以转移到另外的一个状态上去
if(!(now&(7<<(j*3))))//当前点不在任意一个联通块中
{
now|=++used<<(j*3);//新的联通块
for(int k=j+1;k<K;++k)//把所有在一个联通块里的点丢到状态里去
if(getf(j+1)==getf(k+1))
now|=used<<(k*3);
}
g[i][p[now]]++;
}
}
for(int i=1;i<=cnt;++i)
for(int j=1;j<=cnt;++j)
G.s[i][j]=g[i][j];
G=fpow(G,n-K);
int ans=0;
for(int i=1;i<=cnt;++i)
ans=(ans+1ll*G.s[i][1]*f[i])%MOD;
printf("%d\n",ans);
return 0;
}
