多项式入门题

· · 题解

题目简述

无标号无根树计数。n\leq2\times10^5

解题思路

charp1. 约定与预备知识

charp 2. Euler 变换

上面我们的式子是有标号计数,对于无标号计数可以考虑先欧拉变换。我们只关注每一个大小相同的组内部的总体分配情况,而不是单个组。我们可以尝试枚举每种方案选择的次数,不难得到这样的生成函数:

\xi(F)=\prod_{i=1}^{+\infty} (\dfrac{1}{1-x^i})^{f_i}

发现这时候原来在系数的 f 变换到了指数上去,感觉不太好用的样子。我们可以考虑化简一下,对 \xi(F) 使用 \lnf_i 变回系数:

\begin{aligned} \ln \xi(F)&=\sum_{i=1}^{+\infty}(-f_i\times \ln(1-x^i))\\&=\sum_{i=1}^{+\infty}-f_i\times∫\dfrac{-ix^{i-1}}{1-x^i}\\&=\sum_{i=1}^{+\infty}-f_i\times∫-ix^{i-1}(1+x^i+x^{2i}+\dots)\\&=\sum_{i=1}^{+\infty} f_i\times \sum_{j=1}^{+\infty}\dfrac1jx^{ij}\\&=\sum_{j=1}^{+\infty}\dfrac{F(x^j)}{j} \end{aligned}

好看了不少。

charp 3. 无标号有根树计数

假设无标号有根树的生成函数为 F。则 F=x·\xi(F)。这里的含义是当我们确定一个点作为根的时候

,那么相当于给其他元素分组,大小为 i 的组内部有 f_i 的情况。我们考虑解这个方程,首先对两边求导:

\begin{aligned}F'&=x\xi'(F)+\xi(F)\\&=x·\exp'(\sum_{j=1}^{+\infty} \dfrac{F(x^j)}{j})+\dfrac{F}{x}\\&=x·(\sum_{j=1}^{+\infty}\dfrac{F(x^j)}{j})'·\exp(\ln \xi(F))+\dfrac{F}{x}\\&=F\times(\sum_{j=1}^{+\infty}F'(x^j)x^{j-1})+\dfrac Fx\end{aligned}

两边同时乘上 x 会好看一些,变成:

ϑF=F\times(\sum_{j=1}^{+\infty}F'(x^j)x^j)+F

假设 G=\sum_{j=1}^{+\infty}F'(x^k)·x^k,不难得到 g_n=\sum_{d|n} d\times f_n,由上式

n\times f_n=\sum_{i=1}^{n-1} f_ig_{n-i}+f_n

即:

f_n=\dfrac{1}{n-1}\sum_{i=1}^{n-1} f_ig_{n-i}

使用半在线卷积求解即可。

charp 4. 实现细节

有些同学可能不会半在线卷积,这里简要提一嘴。首先根据 P4721 【模板】分治 FFT 的基本思路,我们可以 cdq 分治,这里 g 是固定的,所以对于区间 f[l,mid]f[mid+1,r] 的贡献只要卷上 g[0,r-l] 即可。但是本题当中有个问题,就是可能 f[mid+1,r] 会影响到 g[0,r-l],这时候我们只要考虑 mid+1\leq r-l 当且仅当 l=1,所以我们可以把 l=1 我们只能卷 f[l,mid],g[l,mid] 其他我们后面补回来,即当 l\not=1 时,卷上 f[l,mid],g[0,r-l]g[l,mid],f[0,r-l] 即可。

#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define MP make_pair
#define vll vector<long long>
const int MAXN=4e5+5;
const int MOD=998244353;
namespace polynomial{// yjl poly plank 2.0 ver?
int bfly[MAXN];ll inver[MAXN];
int clogg(int x){return (int)ceil(log2(x));} 
ll ksm(ll a,int b){ll res=1;while(b){if(b&1)res=res*a%MOD;a=a*a%MOD,b>>=1;}return res;}
void butterfly(int l){
    static int las;
    if(las!=l){
        las=l; 
        for(int i=1;i<(1<<l);i++)
            bfly[i]=(bfly[i>>1]>>1)|((i&1)<<l-1);
    } 
}
void NTT(vll &f,int l,int typ){
    butterfly(l);f.resize(1<<l);
    for(int i=0;i<(1<<l);i++)
        if(bfly[i]<i) swap(f[i],f[bfly[i]]);
    for(int i=0;i<l;i++){
        ll step=ksm(3,MOD-1+(MOD-1>>i+1)*typ);
        for(int j=0;j<(1<<l);j+=(1<<i+1)){
            ll cur=1;
            for(int k=j;k<j+(1<<i);k++){
                ll u=f[k],v=f[k+(1<<i)]*cur%MOD;
                f[k]=(u+v)%MOD;f[k+(1<<i)]=(u-v+MOD)%MOD;
                cur=cur*step%MOD;
            }
        }
    }
    if(typ==-1){
        ll val=ksm(1<<l,MOD-2);
        for(int i=0;i<(1<<l);i++)
            f[i]=val*f[i]%MOD;
    }
    return;
}
}using namespace polynomial;
int n;
vll f,g;
void solve(int l,int r){
    if(l==r){
        if(l!=1) f[l]=f[l]*ksm(l-1,MOD-2)%MOD;
        for(int i=l;i<=n;i+=l)
            g[i]+=f[l]*l%MOD,g[i]%=MOD;
        return;
    }
    int mid=l+r>>1;
    solve(l,mid);
    if(l==1){
        vll _f=vll(begin(f)+l,begin(f)+mid+1),_g=vll(begin(g),begin(g)+mid+1);
        int k=clogg(_f.size()+_g.size());
        NTT(_f,k,1);NTT(_g,k,1);
        for(int i=0;i<(1<<k);i++) _f[i]=_f[i]*_g[i]%MOD;
        NTT(_f,k,-1);
        for(int i=mid+1;i<=r;i++)
            f[i]+=_f[i-l],f[i]%=MOD;
        solve(mid+1,r);
        return;
    }
    vll _f=vll(begin(f)+l,begin(f)+mid+1),_g=vll(begin(g),begin(g)+r-l+1);
    int k=clogg(_f.size()+_g.size());
    NTT(_f,k,1);NTT(_g,k,1);
    for(int i=0;i<(1<<k);i++) _f[i]=_f[i]*_g[i]%MOD;
    NTT(_f,k,-1);
    for(int i=mid+1;i<=r;i++)
        f[i]+=_f[i-l],f[i]%=MOD;
    _f=vll(begin(f),begin(f)+r-l+1),_g=vll(begin(g)+l,begin(g)+mid+1);
    k=clogg(_f.size()+_g.size());
    NTT(_f,k,1);NTT(_g,k,1);
    for(int i=0;i<(1<<k);i++) _f[i]=_f[i]*_g[i]%MOD;
    NTT(_f,k,-1);
    for(int i=mid+1;i<=r;i++)
        f[i]+=_f[i-l],f[i]%=MOD;
    solve(mid+1,r);
    return;
}
int main(){
    ios::sync_with_stdio(false);
    cin>>n;f.resize(n+1),g.resize(n+1);
    f[1]=1;
    solve(1,n);
    ll ans=f[n];
    for(int i=1;i<=(n-1)/2;i++) ans-=f[i]*f[n-i]%MOD;
    if(n%2==0) ans-=f[n/2]*(f[n/2]-1)/2%MOD;
    cout<<(ans%MOD+MOD)%MOD;
    return 0;
}