【BZOJ5119】【CTT2017】生成树计数 DP 分治FFT 斯特林数

时间：2022-12-17 13:36:24

　　CTT=清华集训

题目大意

　　有 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;
}
标签：

相关文章