BZOJ 1494 [NOI2007]生成树计数
题目链接:http://www.lydsy.com/JudgeOnline/problem.php?id=1494
题目大意: 给定
n
个点的无向图 。节点编号1......n
第
i
个点与第j个点有边.当且仅当:∣∣i−j∣∣≤k
计算
n
个点时图的生成树数量。2≤k≤5 , 2≤n≤1015
答案 mod 65521
所以我们基于连通性进行递推。并通过构造矩阵进行快速幂。
为什么要基于连通性。如何表示状态。
首先:
当我们考虑加入一个新的点时。如果它不是我们加入的最后一个点。那么我们没有必要一定要使得加入这个点后 。建立一些边后。这个森林会成为一个树。
但我们需要使得更远点。与最后
k
个点有链接。
如果更远的点与最后k个点没有链接。新加入的点又不可能与之建立联系。
导致更远的点断开了链接。就不可能生成一个树。
其次。连通性其实就是特征。根据不同的连通性决定的一类森林进行转移。相同连通性的森林转移到另一连通性时。方案是一样的。
所以可以基于连通性进行递推。
那么如何表示后k 点的连通性呢。
我们用
k
进制数字来表示。[0,k−1]
后
k
个节点的第t 个节点所在集合编号就是这个数字的第
t
位数码。
可能会出现重复状态:
例如:0 1 11 2 2
这两个状态是一样的。
那么我们调整对状态表示的定义。对于确定的连通性。
某一点所在集合的编号如何确定是我们关心的。
那么不如所有集合中最早出现的点的名次作为集合编号。从
0
开始编号。
例如:2 1 2 3 0
他的状态就是:
0 1 0 2 3
我们计算出所有可能的状态。并按顺序对这些状态编号。从0开始。
那么剩下的是构造转移矩阵:
对于状态
a
向b 转移的方案数量我们记为:
M[b][a]
显然之前在一个联通块的点不可能分开。
所以如果在
a
中是连通的在b 中变为不连通是不可能的.这时候转移方案是0。
其次,如果在
a
中是不连通的。在b 中连通的。这时候。如何两个点与最新的点不在一个集合。那也是不可能的。之所以连通只可能是新的点建立 了链接。所以此时方案也是0
除了这两个情况之外。我们开始计算转移方案:
对于与新的点在一个联通块的点。之前状态所在的集合在
k
个点中占有的点数就是新点与之链接的方案数。根据乘法原理:M[b][a]=c1∗c2..∗cr
这表示让不相关的
r
个集合连通的方案只需要分开计算,相乘即可。
其次。构造出来的矩阵无法从0开始递推。
因为在n<k 时情况的特殊性。
我们需要计算。
k
个点时。在不同连通状态下的方案数。作为右侧向量。
从k 递推到
n
。
对于这一步。与之前的想法差不多。对于第i 个连通块的方案。
第一步:构造关于第i个连通块的子图。
第二步:构造关于这个子图的基尔霍夫矩阵。
第三步:对这个矩阵的 n−1 阶进行对角化或变为上下三角。
第四步:对角化后计算对角线的乘机即位方案数。(即行列式的值)
单独计算每个连通块的搭建方案后乘即可。
这样就可以愉快的矩阵快速幂了。
#include <stdio.h>
#include <string.h>
#include <algorithm>
#include <cmath>
#define MAXN 1000
using namespace std;
typedef long long LL;
const int mod=65521;
int T[100][6],dep;
int D[100];
bool vis[6];
int C[6];
int A[5];
int M[60][60];
int E[10][10];
int tmE[10][10];
int tma[6],tmb[6];
int cmp(int k,int a,int b) //a-->b
{
for(int i=0;i<k;i++)
{
tma[i]=T[a][i];
tmb[i]=T[b][i];
}
for(int i=0;i<k-1;i++)
for(int j=i+1;j<k-1;j++)
{
if(tma[i+1]==tma[j+1]&&tmb[i]!=tmb[j])return 0;
if(tma[i+1]!=tma[j+1]&&tmb[i]==tmb[j]&&tmb[i]!=tmb[k-1])return 0;
}
memset(vis,0,sizeof vis);
int ans=1;
for(int i=0;i<k-1;i++)
{
if(tmb[i]!=tmb[k-1])continue;
if(vis[tma[i+1]]) continue;
vis[tma[i+1]]=true;
int c=0;
for(int j=0;j<k;j++)
if(tma[j]==tma[i+1])c++;
ans=(LL)ans*c%mod;
}
return ans;
}
struct mat
{
int m[60][60];
mat()
{
memset(m,0,sizeof m);
}
void e()
{
for(int i=0;i<60;i++)m[i][i]=1;
}
mat operator *(const mat & a)const
{
mat b;
for(int i=0;i<dep;i++)
for(int j=0;j<dep;j++)
for(int k=0;k<dep;k++)
b.m[i][j]=((LL)b.m[i][j]+(LL)m[i][k]*a.m[k][j])%mod;
return b;
}
};
int Iv[mod+10];
int Kirchhof(int n)
{
n--;
int r=0;
for(int i=0; i<n ;i++)
{
for(int j=r; j<n ;j++)
if(tmE[j][i]>0)
{
for(int k=i;k<n;k++)
swap(tmE[j][k],tmE[r][k]);
break;
}
if(tmE[r][i]==0)continue;
for(int j=0;j<n;j++)
if(j!=r&&tmE[j][i]>0)
{
int tm=(LL)tmE[j][i]*Iv[tmE[r][i]]%mod;
for(int k=i;k<n;k++)
tmE[j][k]=(tmE[j][k]-(LL)tm*tmE[r][k]%mod+mod)%mod;
}
r++;
}
int ans=1;
for(int i=0;i<n;i++)
ans=(LL)ans*tmE[i][i]%mod;
return ans;
}
int B[100];
int clat(int k,int d)
{
int ans=1;
memset(vis,0,sizeof vis);
for(int i=0;i<k;i++)
{
if(vis[T[d][i]]) continue;
vis[T[d][i]]=true;
memset(tmE,0,sizeof tmE);
int u=0;
for(int j=0;j<k;j++)
if(T[d][j]==T[d][i])C[u++]=j;
for(int j=0;j<k;j++)
for(int h=0;h<k;h++)
{
if(T[d][j]!=T[d][i])continue;
if(T[d][h]!=T[d][i])continue;
int a=(int)(lower_bound(C,C+u,j)-C);
int b=(int)(lower_bound(C,C+u,h)-C);
tmE[a][b]=E[j][h];
}
for(int j=0;j<u;j++)
for(int h=0;h<u;h++)
if(j!=h)
tmE[j][j]=(tmE[j][j]-tmE[j][h]+mod)%mod;
ans=(LL)ans*Kirchhof(u)%mod;
}
return ans;
}
int slove(LL n,int k)
{
if(n<=(LL)k)
{
for(int i=0;i<n;i++)
for(int j=0;j<n;j++)
tmE[i][j]=E[i][j];
for(int i=0;i<n;i++)
for(int j=0;j<n;j++)
if(i!=j)
tmE[i][i]=(tmE[i][i]-tmE[i][j]+mod)%mod;
return Kirchhof((int)n);
}
mat a,tmp;
tmp.e();
for(int i=0;i<dep;i++)
for(int j=0;j<dep;j++)
a.m[i][j]=M[i][j];
n-=k;
while(n)
{
if(n&1)
tmp=tmp*a;
a=a*a;
n>>=1;
}
int ans=0;
for(int i=0;i<dep;i++)
ans=(ans+(LL)tmp.m[0][i]*clat(k,i)%mod)%mod;
return ans;
}
int main ()
{
Iv[1]=1;
for(int i=2;i<mod;i++) Iv[i]=mod-(LL)(mod/i)*Iv[mod%i]%mod;
LL n;
int k;
scanf("%d %lld",&k,&n);
int size=0,ten=1,bi=0;
while(bi<k)
{
ten*=10;
size=size*10+(k-1);
bi++;
}
ten/=10;
size+=1;
for(int i=0;i<size;i++)
{
for(int K=0,u=ten,d=i;K<k;K++)
{
A[K]=d/u;
d=d%u;
u/=10;
}
int cnt=0,flag=false;
for(int K=0;K<k;K++)
{
if(A[K]>cnt)
{
flag=true;
break;
}
if(A[K]==cnt) cnt++;
}
if(flag)continue;
for(int t=0;t<k;t++) T[dep][t]=A[t];
dep++;
}
for(int i=0;i<dep;i++)
for(int j=0;j<dep;j++)
M[i][j]=cmp(k,j,i);
for(int i=0;i<10;i++)
for(int j=0;j<10;j++)
{
if(i==j)continue;
if(abs(j-i)<=k)E[i][j]=mod-1;
}
printf("%d\n",slove(n,k));
return 0;
}