【BZOJ5119】【CTT2017】生成树计数 DP 分治FFT 斯特林数
2017-12-12 10:02
465 查看
CTT=清华集训
对于一种连边方案T,设第i个点的度数为di,那么这棵树的价值为:
val(T)=(∏i=1nadiidmi)(∑i=1ndmi)
求所有生成树的价值和mod998244353
n≤30000,m≤30
先把式子化简:
ans=(∏i=1nadiidmi)(∑i=1ndmi)=∑i=1ndmi(∏j=1nadjjdmj)=∑i=1nadiid2mi(∏j=1,j!=inadjjdmj)
就是得到排列后选一个数i,把贡献乘上dmi
考虑最简单的做法,记fi,j为前i个数,占了j个空,没有额外的选择一个数的贡献,gi,j就是选了一个数。
很容易得到DP式。这个DP是O(n3)的,可以用FFT优化到O(n2logn)
注意到m很小,还记得那个乘方转组合数和斯特林数的套路吗?
先在prufer数列后面补上1~n,这样写了数字i的格子个数就是di
dmi就是用m种颜色染di个格子,每种颜色只能染一个格子,每个格子可以染多种颜色的方案数。
因为染了颜色的格子很少,所以枚举写了数字i的格子中染了颜色的格子数量k,如果染色的格子全在前面,方案数就是(n−2−jk)×S(m,k)×k!,贡献就是ak+1i。后面有就是(n−2−jk)×S(m,k+1)×(k+1)!
额外选一个数就把m改成2m
fi,j,gi,j的意义和以前类似,只不过改成了染色的格子个数是j。
最后枚举染色的格子个数i,没染色的格子可以随便填数,要乘上(∑ni=1ai)n−2−i
时间复杂度:O(n2m)
把组合数什么的全部拆开后就能得到一个比较漂亮的式子(其中F,G,A,B都是多项式):
FiGi=Fi−1A=Gi−1A+Fi−1B
可以用分治FFT优化。
时间复杂度:O(nmlog2n)
题目大意
有n个点,点权为ai,你要连接一条边,使该图变成一颗树。对于一种连边方案T,设第i个点的度数为di,那么这棵树的价值为:
val(T)=(∏i=1nadiidmi)(∑i=1ndmi)
求所有生成树的价值和mod998244353
n≤30000,m≤30
题解
很容易想到prufer序列先把式子化简:
ans=(∏i=1nadiidmi)(∑i=1ndmi)=∑i=1ndmi(∏j=1nadjjdmj)=∑i=1nadiid2mi(∏j=1,j!=inadjjdmj)
就是得到排列后选一个数i,把贡献乘上dmi
考虑最简单的做法,记fi,j为前i个数,占了j个空,没有额外的选择一个数的贡献,gi,j就是选了一个数。
很容易得到DP式。这个DP是O(n3)的,可以用FFT优化到O(n2logn)
注意到m很小,还记得那个乘方转组合数和斯特林数的套路吗?
先在prufer数列后面补上1~n,这样写了数字i的格子个数就是di
dmi就是用m种颜色染di个格子,每种颜色只能染一个格子,每个格子可以染多种颜色的方案数。
因为染了颜色的格子很少,所以枚举写了数字i的格子中染了颜色的格子数量k,如果染色的格子全在前面,方案数就是(n−2−jk)×S(m,k)×k!,贡献就是ak+1i。后面有就是(n−2−jk)×S(m,k+1)×(k+1)!
额外选一个数就把m改成2m
fi,j,gi,j的意义和以前类似,只不过改成了染色的格子个数是j。
最后枚举染色的格子个数i,没染色的格子可以随便填数,要乘上(∑ni=1ai)n−2−i
时间复杂度:O(n2m)
把组合数什么的全部拆开后就能得到一个比较漂亮的式子(其中F,G,A,B都是多项式):
FiGi=Fi−1A=Gi−1A+Fi−1B
可以用分治FFT优化。
时间复杂度:O(nmlog2n)
代码
暴力
#include<cstdio> #include<cstring> #include<algorithm> #include<cstdlib> #include<ctime> #include<utility> #include<cmath> #include<functional> using namespace std; typedef long long ll; typedef unsigned long long ull; typedef pair<int,int> pii; typedef pair<ll,ll> pll; void sort(int &a,int &b) { if(a>b) swap(a,b); } void open(const char *s) { #ifndef ONLINE_JUDGE char str[100]; sprintf(str,"%s.in",s); freopen(str,"r",stdin); sprintf(str,"%s.out",s); freopen(str,"w",stdout); #endif } int rd() { int s=0,c; while((c=getchar())<'0'||c>'9'); do { s=s*10+c-'0'; } while((c=getchar())>='0'&&c<='9'); return s; } int upmin(int &a,int b) { if(b<a) { a=b; return 1; } return 0; } int upmax(int &a,int b) { if(b>a) { a=b; return 1; } return 0; } const ll p=998244353; ll h[62][62]; ll f[2][60010]; ll g[2][60010]; ll a[30010]; ll d[30010][62]; ll fac[30010]; ll ifac[30010]; ll inv[30010]; void add(ll &a,ll b) { a=(a+b)%p; } ll fp(ll a,ll b) { ll s=1; for(;b;b>>=1,a=a*a%p) if(b&1) s=s*a%p; return s; } ll getc(int x,int y) { if(x<y) return 0; return fac[x]*ifac[y]%p*ifac[x-y]%p; } ll s1[60010][62]; ll s2[60010][62]; ll s3[60010][62]; ll s4[60010][62]; int main() { freopen("a.in","r",stdin); freopen("a2.out","w",stdout); int n,m; scanf("%d%d",&n,&m); h[0][0]=1; int i,j,k; fac[0]=fac[1]=ifac[0]=ifac[1]=inv[0]=inv[1]=1; for(i=2;i<=2*m||i<=n;i++) { inv[i]=-p/i*inv[p%i]%p; fac[i]=fac[i-1]*i%p; ifac[i]=ifac[i-1]*inv[i]%p; } for(i=1;i<=2*m;i++) for(j=1;j<=2*m;j++) h[i][j]=(h[i-1][j-1]+h[i-1][j]*j)%p; for(i=0;i<=2*m;i++) for(j=0;j<=2*m;j++) h[i][j]=h[i][j]*fac[j]%p; for(i=1;i<=n;i++) { scanf("%lld",&a[i]); d[i][0]=1; for(j=1;j<=2*m;j++) d[i][j]=d[i][j-1]*a[i]%p; } for(i=1;i<=n;i++) for(j=0;j<=2*m;j++) { s1[i][j]=ifac[j]*d[i][j]%p*((h[m][j]+h[m][j+1])%p)%p; s2[i][j]=ifac[j]*d[i][j]%p*((h[2*m][j]+h[2*m][j+1])%p)%p; } f[0][0]=1; int t=0; for(i=0;i<n;i++) { t^=1; memset(f[t],0,sizeof f[t]); memset(g[t],0,sizeof g[t]); for(j=0;j<=n-2;j++) { for(k=0;k<=m&&k<=n-2-j;k++) { add(f[t][j+k],f[t^1][j]*s1[i+1][k]); add(g[t][j+k],g[t^1][j]*s1[i+1][k]); } for(k=0;k<=2*m&&k<=n-2-j;k++) add(g[t][j+k],f[t^1][j]*s2[i+1][k]); } } ll ans=0,s=0; for(i=1;i<=n;i++) s=(s+a[i])%p; for(i=0;i<=n-2;i++) ans=(ans+g[t][i]*ifac[n-2-i]%p*fp(s,n-2-i))%p; for(i=1;i<=n;i++) ans=ans*a[i]%p; ans=ans*fac[n-2]%p; ans=(ans+p)%p; printf("%lld\n",ans); return 0; }
正解
#include<cstdio> #include<cstring> #include<algorithm> #include<cstdlib> #include<ctime> #include<utility> #include<cmath> #include<functional> using namespace std; typedef long long ll; typedef unsigned long long ull; typedef pair<int,int> pii; typedef pair<ll,ll> pll; void sort(int &a,int &b) { if(a>b) swap(a,b); } void open(const char *s) { #ifndef ONLINE_JUDGE char str[100]; sprintf(str,"%s.in",s); freopen(str,"r",stdin); sprintf(str,"%s.out",s); freopen(str,"w",stdout); #endif } int rd() { int s=0,c; while((c=getchar())<'0'||c>'9'); do { s=s*10+c-'0'; } while((c=getchar())>='0'&&c<='9'); return s; } int upmin(int &a,int b) { if(b<a) { a=b; return 1; } return 0; } int upmax(int &a,int b) { if(b>a) { a=b; return 1; } return 0; } const int p=998244353; int fp(ll a,ll b) { ll s=1; for(;b;b>>=1,a=1ll*a*a%p) if(b&1) s=1ll*s*a%p; return s; } namespace ntt { const int g=3; int w1[150000]; int w2[150000]; int rev[150000]; int n; void init(int m) { n=m; int i; for(i=2;i<=n;i<<=1) { w1[i]=fp(g,(p-1)/i); w2[i]=fp(w1[i],p-2); } rev[0]=0; for(i=1;i<n;i++) rev[i]=(rev[i>>1]>>1)|(i&1?n>>1:0); } void ntt(int *a,int t) { int u,v,w,wn; int i,j,k; for(i=0;i<n;i++) if(rev[i]<i) swap(a[i],a[rev[i]]); for(i=2;i<=n;i<<=1) { wn=(t==1?w1[i]:w2[i]); for(j=0;j<n;j+=i) { w=1; for(k=j;k<j+i/2;k++) { u=a[k]; v=1ll*a[k+i/2]*w%p; a[k]=(u+v)%p; a[k+i/2]=(u-v)%p; w=1ll*w*wn%p; } } } if(t==-1) { int inv=fp(n,p-2); for(i=0;i<n;i++) a[i]=1ll*a[i]*inv%p; } } void copy(int *a,const vector<int> &b,int len) { int i; for(i=0;i<=len;i++) a[i]=b[i]; for(i=len+1;i<n;i++) a[i]=0; } vector<int> back(int *a,int len) { vector<int> s; int i; for(i=0;i<=len;i++) s.push_back(a[i]); return s; } } int n,m; int h[62][62]; int f[2][60010]; int g[2][60010]; int a[30010]; int d[30010][62]; int fac[30010]; int ifac[30010]; int inv[30010]; void add(int &a,int b) { a=(a+b)%p; } int s1[60010][62]; int s2[60010][62]; auto gao(vector<int> &a,vector<int> &b,vector<int> &c,vector<int> &d) { static int a1[150000],a2[150000],a3[150000],a4[150000],s1[150000],s2[150000]; int len=1; int n1=a.size()-1; int n2=b.size()-1; int n3=c.size()-1; int n4=d.size()-1; while(len<=2*n2||len<=2*n4) len<<=1; ntt::init(len); ntt::copy(a1,a,n1); ntt::copy(a2,b,n2); ntt::copy(a3,c,n3); ntt::copy(a4,d,n4); ntt::ntt(a1,1); ntt::ntt(a2,1); ntt::ntt(a3,1); ntt::ntt(a4,1); int i; for(i=0;i<len;i++) { s1[i]=1ll*a1[i]*a3[i]%p; s2[i]=(1ll*a1[i]*a4[i]+1ll*a2[i]*a3[i])%p; } ntt::ntt(s1,-1); ntt::ntt(s2,-1); return make_pair(ntt::back(s1,min(n-2,n1+n3)),ntt::back(s2,min(n-2,max(n1+n4,n2+n3)))); } auto solve(int l,int r) { if(l==r) { vector<int> a(m+1),b(2*m+1); int i; for(i=0;i<=2*m;i++) { if(i<=m) a[i]=s1[l][i]; b[i]=s2[l][i]; } return make_pair(a,b); } int mid=(l+r)>>1; auto L=solve(l,mid); auto R=solve(mid+1,r); return gao(L.first,L.second,R.first,R.second); } int main() { open("a"); scanf("%d%d",&n,&m); h[0][0]=1; int i,j,k; fac[0]=fac[1]=ifac[0]=ifac[1]=inv[0]=inv[1]=1; for(i=2;i<=2*m||i<=n;i++) { inv[i]=1ll*-p/i*inv[p%i]%p; fac[i]=1ll*fac[i-1]*i%p; ifac[i]=1ll*ifac[i-1]*inv[i]%p; } for(i=1;i<=2*m;i++) for(j=1;j<=2*m;j++) h[i][j]=(h[i-1][j-1]+1ll*h[i-1][j]*j)%p; for(i=0;i<=2*m;i++) for(j=0;j<=2*m;j++) h[i][j]=1ll*h[i][j]*fac[j]%p; for(i=1;i<=n;i++) { scanf("%d",&a[i]); d[i][0]=1; for(j=1;j<=2*m;j++) d[i][j]=1ll*d[i][j-1]*a[i]%p; } for(i=1;i<=n;i++) for(j=0;j<=2*m;j++) { s1[i][j]=1ll*ifac[j]*d[i][j]%p*((1ll*h[m][j]+h[m][j+1])%p)%p; s2[i][j]=1ll*ifac[j]*d[i][j]%p*((1ll*h[2*m][j]+h[2*m][j+1])%p)%p; } // f[0][0]=1; int t=0; // for(i=0;i<n;i++) // { // t^=1; // memset(f[t],0,sizeof f[t]); // memset(g[t],0,sizeof g[t]); // for(j=0;j<=n-2;j++) // { // for(k=0;k<=m;k++) // { // add(f[t][j+k],f[t^1][j]*s1[i+1][k]); // add(g[t][j+k],g[t^1][j]*s1[i+1][k]); // } // for(k=0;k<=2*m;k++) // add(g[t][j+k],f[t^1][j]*s2[i+1][k]); // } // } auto dp=solve(1,n); int ans=0,s=0; for(i=1;i<=n;i++) s=(s+a[i])%p; for(i=0;i<=n-2;i++) // ans=(ans+g[t][i]*ifac[n-2-i]%p*fp(s,n-2-i))%p; ans=(ans+1ll*dp.second[i]*ifac[n-2-i]%p*fp(s,n-2-i))%p; for(i=1;i<=n;i++) ans=1ll*ans*a[i]%p; ans=1ll*ans*fac[n-2]%p; ans=(ans+p)%p; printf("%d\n",ans); return 0; }
相关文章推荐
- BZOJ 4555: [Tjoi2016&Heoi2016]求和 [分治FFT 组合计数 | 多项式求逆]
- BZOJ 3193: [JLOI2013]地形生成【计数dp
- 【NOI2007/BZOJ1494】生成树计数 插头DP
- [DP 倍增] BZOJ 4818 [Sdoi2017]序列计数
- BZOJ 2152: 聪聪可可(点分治/树形DP)
- 【codevs1359】【BZOJ1833】数字计数,进击的学弟与数位DP
- bzoj1016 最小生成树计数 (性质+克鲁斯卡尔)[省选计划系列题]
- [bzoj1016]:[JSOI2008]最小生成树计数
- 【BZOJ 1016】 1016: [JSOI2008]最小生成树计数 (DFS|矩阵树定理)
- 【BZOJ3294】[Cqoi2011]放棋子【计数DP】
- 最小生成树计数 BZOJ 1016
- BZOJ2159 Crash 的文明世界 【第二类斯特林数 + 树形dp】
- BZOJ 1833 count 数字计数 (数位DP)
- bzoj 4503: 两个串 (FFT+DP)
- [BZOJ1016][JSOI2008]最小生成树计数
- BZOJ1016: [JSOI2008]最小生成树计数
- bzoj 1833 [ZJOI2010]count 数字计数(数位DP)
- 【BZOJ1030】[JSOI2007]文本生成器【AC自动机】【计数DP】
- 【bzoj1016】[JSOI2008]最小生成树计数
- BZOJ 1801: [Ahoi2009]chess 中国象棋 [DP 组合计数]