发布时间：2025-05-10

多项式

cdq分治FFT

给定序列 $g_{1\dots n - 1}$，求序列 $f_{0\dots n - 1}$。

其中 $f_i=\sum_{j=1}^if_{i-j}g_j$，边界为 $f_0=1$。

答案对 $998244353$ 取模。

考虑递归计算计算$f_l$到$f_{mid}$，接着计算后半部分贡献

$$w_x =\sum_{i=l}^{mid}f_i \times g_{x-i}\\ w_x=\sum_{i=l}^{x-1}f_i \times g_{x-i}\\ w_x=\sum_{i=0}^{x-l-1}f_{i+l} \times g_{x-i-l} 设a_i=f_{i+l},b_i=g_{i+1} w_x=\sum_{i=0}^{x-l-1}a_i \times b_{x-l-1-i}$$

namespace NTT{
    static constexpr int P=998244353;
    int rev[N],A[N],B[N],bit=0,tot,g=3,gi;
    int qpow(int a,int b){
        int ans=1;
        for(;b;b>>=1){
            if(b&1)ans=ans*a%P;
            a=a*a%P;
        }
        return ans;
    }
    int inv(int n){
        return qpow(n,P-2);
    }
    void ntt(int *a,int inv){
        for(int i=0;i<tot;i++){
            if(i<rev[i])swap(a[i],a[rev[i]]);
        }
        for(int mid=1;mid<tot;mid<<=1){
            int g1=qpow(inv==1?g:gi,(P-1)/(mid*2));
            for(int i=0;i<tot;i+=mid*2){
                int gk=1;
                for(int j=0;j<mid;j++,gk=gk*g1%P){
                    auto x=a[i+j],y=gk*a[i+j+mid]%P;
                    a[i+j]=(x+y)%P,a[i+j+mid]=(x-y+P)%P;
                }
            }
        }
        if(inv==-1){
            int tmp=NTT::inv(tot);
            for(int i=0;i<tot;i++)a[i]=a[i]*tmp%P;
        }
    }
    void init(int n){
        bit=0;
        while((1<<bit)<n)bit++;
        tot=1<<bit;
        for(int i=0;i<tot;i++)rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1));
        gi=inv(g);
    }
    void mul(int *a,int *b){
        ntt(a,1),ntt(b,1);
        for(int i=0;i<tot;i++)a[i]=a[i]*b[i]%P;
        ntt(a,-1);
    }
    void cdq_ntt(int *a,int *b,int l,int r){
        if(l==r)return;
        int mid=l+r>>1,len=r-l+1;
        cdq_ntt(a,b,l,mid);
        init(len);
        for(int i=l;i<=mid;i++)A[i-l]=a[i];
        for(int i=1;i<=r-l;i++)B[i-1]=b[i];
        mul(A,B);
        for(int i=mid+1;i<=r;i++)(a[i]+=A[i-l-1])%=P;
        for(int i=0;i<tot;i++)A[i]=B[i]=0;
        cdq_ntt(a,b,mid+1,r);
    }
}
int f[N],g[N];
void solve(){
    int n;
    cin>>n;
    f[0]=1;
    for(int i=1;i<n;i++)cin>>g[i];
    NTT::cdq_ntt(f,g,0,n-1);
    for(int i=0;i<n;i++)cout<<f[i]<<" \n"[i==n-1];
}

板子集合

namespace Poly{
    using poly=vector<int>;
    vector<int> rev;
    const int P=998244353;
    int qpow(int a,int b){
        int ans=1;
        for(;b;b>>=1){
            if(b&1)ans=ans*a%P;
            a=a*a%P;
        }
        return ans;
    }
    int bit=0,tot,g=3,gi=qpow(g,P-2);
    void ntt(poly &a,int inv){
        for(int i=0;i<tot;i++){
            if(i<rev[i])swap(a[i],a[rev[i]]);
        }
        for(int mid=1;mid<tot;mid<<=1){
            int g1=qpow(inv==1?g:gi,(P-1)/(mid*2));
            for(int i=0;i<tot;i+=mid*2){
                int gk=1;
                for(int j=0;j<mid;j++,gk=gk*g1%P){
                    auto x=a[i+j],y=gk*a[i+j+mid]%P;
                    a[i+j]=(x+y)%P,a[i+j+mid]=(x-y+P)%P;
                }
            }
        }
        if(inv==-1){
            int invn=qpow(tot,P-2);
            for(int i=0;i<tot;i++)a[i]=a[i]*invn%P;
        }
    }
    void init(int n){
        bit=0;
        while((1<<bit)<n)bit++;
        tot=1<<bit;
        rev.resize(tot);
        for(int i=0;i<tot;i++)rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1));
    }
    // 多项式乘法
    poly mul(poly a,poly b){
        int n=a.size(),m=b.size();
        init(n+m+1);
        a.resize(tot),b.resize(tot);
        ntt(a,1),ntt(b,1);
        for(int i=0;i<tot;i++)a[i]=a[i]*b[i]%P;
        ntt(a,-1);
        for(int i=n+m;i<tot;i++)a[i]=0;
        return a;
    }
    // 多项式求逆
    poly inv(const poly &f,int n){
        if(n==1)return {qpow(f[0],P-2)};
        auto g=inv(f,n+1>>1);
        static vector<int> c;
        init(2*n+1);
        c.resize(tot);
        g.resize(tot);
        for(int i=0;i<n;i++)c[i]=f[i];
        for(int i=n;i<tot;i++)c[i]=0;
        ntt(c,1),ntt(g,1);
        for(int i=0;i<tot;i++)g[i]=g[i]*(2-c[i]*g[i]%P+P)%P;
        ntt(g,-1);
        for(int i=n;i<tot;i++)g[i]=0;
        return g;
    }
    // 多项式求导
    poly der(const poly &f){
        int n=f.size();
        poly g(n-1);
        for(int i=1;i<n;i++)g[i-1]=f[i]*i%P;
        return g;
    }
    // 多项式求积分
    poly inte(const poly &f){
        int n=f.size();
        poly g(n+1);
        g[0]=0;
        for(int i=1;i<=n;i++)g[i]=f[i-1]*qpow(i,P-2)%P;
        return g;
    }
    // 多项式求值
    int eval(const poly &f,int x){
        int n=f.size(),ans=0;
        for(int i=n-1;i>=0;i--)ans=(ans*x+f[i])%P;
        return ans;
    }
    // 多项式翻转Ar(x)=x^n*A(1/x)
    poly reverse(const poly &f,int n){
        poly g(n);
        for(int i=0;i<n;i++)g[i]=f[n-i-1];
        return g;
    }
    // 多项式加法
    poly add(const poly &f,const poly &g){
        poly c(max(f.size(),g.size()));
        for(int i=0;i<max(f.size(),g.size());i++){
            if(i<f.size())c[i]=f[i];
            if(i<g.size())c[i]=(c[i]+g[i])%P;
        }
        return c;
    }
    // 多项式减法
    poly sub(const poly &f,const poly &g){
        poly c(max(f.size(),g.size()));
        for(int i=0;i<max(f.size(),g.size());i++){
            if(i<f.size())c[i]=f[i];
            if(i<g.size())c[i]=(c[i]-g[i]+P)%P;
        }
        return c;
    }
    // 多项式除法
    pair<poly,poly> div(const poly &f,const poly &g){
        int n=f.size(),m=g.size(),len=n-m+1;
        if(n<m)return {{},f};
        poly fr=reverse(f,n),gr=reverse(g,m);
        gr.resize(len);
        for(int i=len;i<m;i++)gr[i]=0;
        gr=inv(gr,len);
        poly q=mul(fr,gr);
        q=reverse(q,len);
        poly r=sub(f,mul(q,g));
        return {q,r};
    }
    // cipolla
    int cipolla(int x){
        if(qpow(x,P-1>>1)==P-1)return -1;
        while(1){
            int a=rng()%P;
            int w=(a*a%P-x+P)%P;
            auto mul=[&](PII a,PII b){
                return PII{(a.F*b.F%P+w*a.S%P*b.S%P)%P,(a.F*b.S%P+a.S*b.F%P)%P};
            };
            auto pow=[&](PII a,int b){
                PII res{1,0};
                for(;b;b>>=1){
                    if(b&1)res=mul(res,a);
                    a=mul(a,a);
                }
                return res.F;
            };
            if(qpow(w,P-1>>1)==P-1){
                int res=pow(PII{a,1},P+1>>1);
                return min(res,P-res);
            }
        }
    }
    // 多项式开根(a0!=1)
    poly sqrt(const poly &f,int n){
        if(n==1){
            if(f[0]==1)return {1};
            return {cipolla(f[0])};
        }
        auto g=sqrt(f,n+1>>1);
        static vector<int> c;
        static const int inv2=qpow(2,P-2);
        init(2*n+1);
        c.resize(tot);
        g.resize(tot);
        for(int i=0;i<n;i++)c[i]=f[i];
        for(int i=n;i<tot;i++)c[i]=0;
        auto ig=inv(g,n);
        ntt(g,1),ntt(ig,1),ntt(c,1);
        for(int i=0;i<tot;i++)g[i]=(g[i]*g[i]%P+c[i])%P*ig[i]%P*inv2%P;
        ntt(g,-1);
        for(int i=n;i<tot;i++)g[i]=0;
        return g;
    }   
    // 多项式对数ln
    poly ln(const poly &f){
        return inte(mul(der(f),inv(f,f.size())));
    }
    // 多项式指数exp
    poly exp(const poly &f,int n){
        if(n==1)return {1};
        poly g=exp(f,n+1>>1);
        g.resize(n);
        poly lng=ln(g);
        init(2*n+1);
        g.resize(tot);
        lng.resize(tot);
        for(int i=0;i<n;i++)lng[i]=(f[i]-lng[i]+P)%P;
        for(int i=n;i<tot;i++)lng[i]=g[i]=0;
        lng[0]++;
        ntt(lng,1);
        ntt(g,1);
        for(int i=0;i<tot;i++)g[i]=g[i]*lng[i]%P;
        ntt(g,-1);
        for(int i=n;i<tot;i++)g[i]=0;
        return g;
    }
    // 多项式快速幂
    poly qpow(const poly &f,int b){
        int n=f.size();
        poly res{1},a=f;
        res[0]=1;
        init(n<<1);
        res.resize(tot);
        a.resize(tot);
        for(;b;b>>=1){
            ntt(a,1);
            if(b&1){
                ntt(res,1);
                for(int i=0;i<tot;i++)res[i]=res[i]*a[i]%P;
                ntt(res,-1);
                for(int i=n;i<tot;i++)res[i]=0;
            }
            for(int i=0;i<tot;i++)a[i]=a[i]*a[i]%P;
            ntt(a,-1);
            for(int i=n;i<tot;i++)a[i]=0;
        }
        return res;
    }
    // 多项式多点求值
    vector<int> evals(const poly &f,vector<int> &x){
        vector<poly> T;
        int n=f.size(),m=x.size();
        T.assign(m<<2,{});
        function<void(int,int,int)> build=[&](int p,int l,int r){
            if(l==r){
                T[p]={(P-x[l])%P,1};
                return;
            }
            int mid=l+r>>1;
            build(p<<1,l,mid),build(p<<1|1,mid+1,r);
            T[p]=mul(T[p<<1],T[p<<1|1]);
            while(T[p].size()>1&&T[p].back()==0)T[p].pop_back();
        };
        vector<int> res(m);
        function<void(int,int,int,poly)> solve=[&](int p,int l,int r,poly f){
            if(r-l+1<=64){
                for(int i=l;i<=r;i++){
                    res[i]=eval(f,x[i]);
                }
                return;
            }
            auto [q1,r1]=div(f,T[p<<1]);
            auto [q2,r2]=div(f,T[p<<1|1]);
            while(r1.size()>1&&r1.back()==0)r1.pop_back();
            while(r2.size()>1&&r2.back()==0)r2.pop_back();
            int mid=l+r>>1;
            solve(p<<1,l,mid,r1);
            solve(p<<1|1,mid+1,r,r2);
        };
        build(1,0,m-1);
        if(f.size()>=T[1].size()){
            auto [q,r]=div(f,T[1]);
            solve(1,0,m-1,r);
        }
        else solve(1,0,m-1,f);
        return res;
    }
}

modint板子

struct mint {
    constexpr static int mod = 998244353;
	int x;
	inline mint(int o = 0) { x = o%mod; }
	inline mint & operator = (int o) { return x = o%mod, *this; }
	inline mint & operator += (mint o) { return (x += o.x) >= mod && (x -= mod), *this; }
	inline mint & operator -= (mint o) { return (x -= o.x) < 0 && (x += mod), *this; }
	inline mint & operator *= (mint o) { return x = (long long) x * o.x % mod, *this; }
    inline mint& qpow(int b){
        mint w = *this;
        mint ret(1);
        for(; b; b >>= 1, w *= w) if(b & 1) ret *= w;
        return x = ret.x, *this;
    }
	inline mint & operator /= (mint o) { return *this *= (o .qpow( (mod - 2))); }
	friend inline mint operator + (mint a, mint b) { return a += b; }
	friend inline mint operator - (mint a, mint b) { return a -= b; }
	friend inline mint operator * (mint a, mint b) { return a *= b; }
	friend inline mint operator / (mint a, mint b) { return a /= b; }
    friend ostream& operator<<(ostream& os, const mint& m) {
        return os << m.x;
    }
};

拉格朗日插值

$n+1$个数确定$n$次多项式, $f(k)=\sum_{i=0}^n y_i \prod \frac{k-x_j}{x_i-x_j}$

O(n^2)版本

int lagrange(int k,vector<int> &X,vector<int> &Y){
    int ans=0;
    int n=X.size();
    for(int i=0;i<n;i++){
        int fz=1,fm=1;
        for(int j=0;j<n;j++){
            if(i==j)continue;
            fz=fz*(k-X[j])%mod;
            fm=fm*(X[i]-X[j])%mod;
        }
        ans=(ans+Y[i]*fz%mod*qpow(fm,mod-2)%mod)%mod;
    }
    return (ans+mod)%mod;
}

O(n)版本

• 如果$x$取值连续，那么把$x_i$换成$i$,$f(k)=\sum_{i=0}^ny_i \prod \frac{k-i}{i-j}$
• 记$pre_i=\prod_{j=0}^i k-j,suf_i=\prod_{j=i}^n k-j$
• $f(k)=\sum_{i=0}^ny_i\frac{pre_{i-1}suf_{i+1}}{fac_{i-1}fac_{n-i}}$,$n-i$为奇数时分母取负

int lagrange(int k,vector<int> &X,vector<int> &Y,int n){
    int ans=0;
    vector<int> pre(n+10),suf(n+10),ifac(n+10);
    pre[0]=(k-X[0])%mod,suf[n+1]=1;
    for(int i=1;i<=n;i++)pre[i]=pre[i-1]*(k-X[i])%mod;
    for(int i=n;i>=0;i--)suf[i]=suf[i+1]*(k-X[i])%mod;
    ifac[0]=ifac[1]=1;
    for(int i=2;i<=n;i++)ifac[i]=-1ll*mod/i*ifac[mod%i]%mod;
    for(int i=2;i<=n;i++)ifac[i]=ifac[i]*ifac[i-1]%mod;
    for(int i=0;i<=n;i++){
        int tmp=Y[i]*(i==0?1:pre[i-1])%mod*suf[i+1]%mod;
        tmp=tmp*ifac[i]%mod*ifac[n-i]%mod;
        if(n-i&1)tmp=mod-tmp;
        ans=(ans+tmp)%mod;
    }
    return (ans+mod)%mod;
}

求$\sum_{i=1}^n i^k$

复杂度$O(k)$

void solve(){
    int n,k;
    cin>>n>>k;
    vector<int> X(k+2),Y(k+2);
    iota(X.begin(),X.end(),1);
    Y[0]=1;
    for(int i=1;i<=k+1;i++)Y[i]=(Y[i-1]+qpow(i+1,k))%mod;
    cout<<lagrange(n,X,Y,k+1);
}

重心拉格朗日插值法

• $g=\prod _{i=0}^nk-x_i,t_i=\frac{y_i}{\prod x_i-x_j}$
• $f(k)=g \sum_{i=0}^n\frac{t_i}{k-x_i}$ 每新加入一个点只需要计算$t_i$即可，预处理仍然是$O(n^2)$,但是插入是$O(n)$的

多项式乘法逆

给定一个多项式 $F(x)$ ，请求出一个多项式 $G(x)$， 满足 $F(x) * G(x) \equiv 1 \pmod{x^n}$。系数对 $998244353$ 取模。

牛顿迭代法 已知一个复合多项式 $G(F(x)) \equiv 0(\mod x^n)$，假设 $F_0(x) \equiv F(x)(\mod x^{\frac{n}{2}}),F_1(x) \equiv F(x)(\mod x^n)$

泰勒公式在$F_0(x)$展开$G(F(x))$,$G(F(x))=\sum_{i=0}^{\infty}\frac{G^{(i)}F_0(x)}{i!}(F(x)-F_0(x))^i$

$G(F_1(x)) \equiv G(F_0(x))+G'(F_0(x))(F_1(x)-F_0(x))(\mod x^n)$

$F_1(x)=F_0(x)-\frac{G(F_0(x))}{G'(F_0(x))}$

将$H(t)=F(x)-\frac{1}{t}$套进去

$G_1(x)=G_0(x)-\frac{F(x)-\frac{1}{G_0(x)}}{\frac{1}{G_0(x)^2}}$

$G_1(x)=G_0(x)(2-G_0(x)F(x))$

poly inv(poly &f,int n){
    if(n==1)return {qpow(f[0],P-2)};
    auto g=inv(f,n+1>>1);
    static vector<int> c;
    init(2*n+1);
    c.resize(tot);
    g.resize(tot);
    for(int i=0;i<n;i++)c[i]=f[i];
    for(int i=n;i<tot;i++)c[i]=0;
    ntt(c,1),ntt(g,1);
    for(int i=0;i<tot;i++)g[i]=g[i]*(2-c[i]*g[i]%P+P)%P;
    ntt(g,-1);
    for(int i=n;i<tot;i++)g[i]=0;
    return g;
}
void solve(){
    int n;
    cin>>n;
    poly f(n+1);
    for(int i=0;i<n;i++)cin>>f[i];
    auto g=inv(f,n);
    for(int i=0;i<n;i++)cout<<g[i]<<" ";
}

多项式求导

poly der(poly &f){
    int n=f.size();
    poly g(n-1);
    for(int i=1;i<n;i++)g[i-1]=f[i]*i%P;
    return g;
}

多项式求积分

poly inte(poly &f){
    int n=f.size();
    poly g(n+1);
    g[0]=0;
    for(int i=1;i<=n;i++)g[i]=f[i-1]*qpow(i,P-2)%P;
    return g;
}

多项式除法

给定一个 $n$ 次多项式 $F(x)$ 和一个 $m$ 次多项式 $G(x)$ ，请求出多项式 $Q(x)$, $R(x)$，满足以下条件：

• $Q(x)$ 次数为 $n-m$，$R(x)$ 次数小于 $m$
• $F(x) = Q(x) * G(x) + R(x)$

考虑一种翻转操作，$A_R(x)=x^nA(\frac{1}{x})$，即$A_{n-i}=A_{i}$

$F(\frac{1}{x})=Q(\frac{1}{x})G(\frac{1}{x})+R(\frac{1}{x})$

$x^nF(\frac{1}{x})=x^{n-m}Q(\frac{1}{x}) * x^mG(\frac{1}{x})+x^{n-m+1}*x^{m-1}R(\frac{1}{x})$

$F_R(x)=Q_R(x)G_R(x)+x^{n-m+1}R_R(x)$

$F_R(x) \equiv Q_R(x)G_R(x)(\mod x^{n-m+1})$

$Q_R(x) \equiv F_R(x)G_R^{-1}(x)(\mod x^{n-m+1})$

最后$R(x)=F(x)-G(x)Q(x)$

pair<poly,poly> div(const poly &f,const poly &g){
    int n=f.size(),m=g.size(),len=n-m+1;
    poly fr=reverse(f,n),gr=reverse(g,m);
    gr.resize(len);
    for(int i=len;i<m;i++)gr[i]=0;
    gr=inv(gr,len);
    poly q=mul(fr,gr,n,len);
    q=reverse(q,len);
    poly r=sub(f,mul(q,g,len,m));
    return {q,r};
}
void solve(){
    int n,m;
    cin>>n>>m;
    poly f(n+1),g(m+1);
    for(int i=0;i<=n;i++)cin>>f[i];
    for(int i=0;i<=m;i++)cin>>g[i];
    auto [q,r]=div(f,g);
    for(int i=0;i<=n-m;i++)cout<<q[i]<<" \n"[i==n-m];
    for(int i=0;i<m;i++)cout<<r[i]<<" \n"[i==m-1];
}

多项式开根(a0=1)

给定一个 $n-1$ 次多项式 $A(x)$，求一个在 ${} \bmod x^n$ 意义下的多项式 $B(x)$，使得 $B^2(x) \equiv A(x) \pmod{x^n}$。若有多解，请取零次项系数较小的作为答案。

多项式的系数在 ${}\bmod 998244353$ 的意义下进行运算。

牛顿迭代法，构造式子$H(t)=t^2-F(x)$

带入式子得$G_1(x)=G_0(x)-\frac{G_0(x)^2-F(x)}{2G_0(x)}=\frac{G_0(x)^2+F(x)}{2G_0(x)}=\frac{1}{2}(G_0(x)+\frac{F(x)}{G_0(x)})$

poly sqrt(const poly &f,int n){
    if(n==1)return {1};
    auto g=sqrt(f,n+1>>1);
    static vector<int> c;
    static const int inv2=qpow(2,P-2);
    init(2*n+1);
    c.resize(tot);
    g.resize(tot);
    for(int i=0;i<n;i++)c[i]=f[i];
    for(int i=n;i<tot;i++)c[i]=0;
    auto ig=inv(g,n);
    ntt(g,1),ntt(ig,1),ntt(c,1);
    for(int i=0;i<tot;i++)g[i]=(g[i]*g[i]%P+c[i])%P*ig[i]%P*inv2%P;
    ntt(g,-1);
    for(int i=n;i<tot;i++)g[i]=0;
    return g;
}  
void solve(){
    int n;
    cin>>n;
    poly f(n+1);
    for(int i=0;i<n;i++)cin>>f[i];
    auto g=sqrt(f,n+1);
    for(int i=0;i<n;i++)cout<<g[i]<<" ";
}

多项式开根a0!=1

需要检验a0是否有二次剩余，即$x^2 \equiv a_0(\mod 998244353)$是否成立

int cipolla(int x){
    if(qpow(x,P-1>>1)==P-1)return -1;
    while(1){
        int a=rng()%P;
        int w=(a*a%P-x+P)%P;
        auto mul=[&](PII a,PII b){
            return PII{(a.F*b.F%P+w*a.S%P*b.S%P)%P,(a.F*b.S%P+a.S*b.F%P)%P};
        };
        auto pow=[&](PII a,int b){
            PII res{1,0};
            for(;b;b>>=1){
                if(b&1)res=mul(res,a);
                a=mul(a,a);
            }
            return res.F;
        };
        if(qpow(w,P-1>>1)==P-1){
            int res=pow(PII{a,1},P+1>>1);
            return min(res,P-res);
        }
    }
}
// 多项式开根(a0!=1)
poly sqrt(const poly &f,int n){
    if(n==1){
        if(f[0]==1)return {1};
        return {cipolla(f[0])};
    }
    auto g=sqrt(f,n+1>>1);
    static vector<int> c;
    static const int inv2=qpow(2,P-2);
    init(2*n+1);
    c.resize(tot);
    g.resize(tot);
    for(int i=0;i<n;i++)c[i]=f[i];
    for(int i=n;i<tot;i++)c[i]=0;
    auto ig=inv(g,n);
    ntt(g,1),ntt(ig,1),ntt(c,1);
    for(int i=0;i<tot;i++)g[i]=(g[i]*g[i]%P+c[i])%P*ig[i]%P*inv2%P;
    ntt(g,-1);
    for(int i=n;i<tot;i++)g[i]=0;
    return g;
}

多项式对数函数（多项式 ln）(a_0=1)

给出 $n-1$ 次多项式 $A(x)$，求一个 $\bmod{\:x^n}$ 下的多项式 $B(x)$，满足 $B(x) \equiv \ln A(x)$.

在 $\text{mod } 998244353$ 下进行，且 $a_i \in [0, 998244353) \cap \mathbb{Z}$

保证 $a_0 = 1$.

令 $f(x)=\ln x$,$G(x) \equiv f(F(x))(\mod x^n)$

两边求导 $G'(x)=f'(F(x))F'(x)=\frac{F'(x)}{F(x)}$

求出来后再积回去

poly ln(const poly &f,int n){
    return inte(mul(der(f),inv(f,n)));
}

多项式指数函数（多项式 exp）(a_0=0)

给出 $n-1$ 次多项式 $A(x)$，求一个 $\bmod{\:x^n}$ 下的多项式 $B(x)$，满足 $B(x) \equiv \text e^{A(x)}$。系数对 $998244353$ 取模。

保证 $a_0 = 0$.

牛顿迭代法，$H(t)=\ln t-F(x)$,$G_1(x)\equiv G_0(x)(1-\ln G_0(x)+F(x))$

poly exp(const poly &f,int n){
    if(n==1)return {1};
    poly g=exp(f,n+1>>1);
    g.resize(n);
    poly lng=ln(g);
    init(2*n+1);
    g.resize(tot);
    lng.resize(tot);
    for(int i=0;i<n;i++)lng[i]=(f[i]-lng[i]+P)%P;
    for(int i=n;i<tot;i++)lng[i]=g[i]=0;
    lng[0]++;
    ntt(lng,1);
    ntt(g,1);
    for(int i=0;i<tot;i++)g[i]=g[i]*lng[i]%P;
    ntt(g,-1);
    for(int i=n;i<tot;i++)g[i]=0;
    return g;
}

多项式多点求值

// 多项式多点求值
vector<int> evals(const poly &f,vector<int> &x){
    vector<poly> T;
    int n=f.size(),m=x.size();
    T.assign(m<<2,{});
    function<void(int,int,int)> build=[&](int p,int l,int r){
        if(l==r){
            T[p]={(P-x[l])%P,1};
            return;
        }
        int mid=l+r>>1;
        build(p<<1,l,mid),build(p<<1|1,mid+1,r);
        T[p]=mul(T[p<<1],T[p<<1|1]);
        while(T[p].size()>1&&T[p].back()==0)T[p].pop_back();
    };
    vector<int> res(m);
    function<void(int,int,int,poly)> solve=[&](int p,int l,int r,poly f){
        if(r-l+1<=64){
            for(int i=l;i<=r;i++){
                res[i]=eval(f,x[i]);
            }
            return;
        }
        auto [q1,r1]=div(f,T[p<<1]);
        auto [q2,r2]=div(f,T[p<<1|1]);
        while(r1.size()>1&&r1.back()==0)r1.pop_back();
        while(r2.size()>1&&r2.back()==0)r2.pop_back();
        int mid=l+r>>1;
        solve(p<<1,l,mid,r1);
        solve(p<<1|1,mid+1,r,r2);
    };
    build(1,0,m-1);
    if(f.size()>=T[1].size()){
        auto [q,r]=div(f,T[1]);
        solve(1,0,m-1,r);
    }
    else solve(1,0,m-1,f);
    return res;
}
void solve(){
    int n,m;
    cin>>n>>m;
    n++;
    vector<int> f(n),a(m);
    for(int i=0;i<n;i++)cin>>f[i];
    for(int i=0;i<m;i++)cin>>a[i];
    auto ans=evals(f,a);
    for(int i=0;i<m;i++)cout<<ans[i]<<"\n";
}

多项式快速幂(a_0=1)

• 首先对k直接取模，然后直接算，复杂度$O(n \log n \log (k \mod 998244353))$

poly qpow(const poly &f,int b){
    int n=f.size();
    poly res{1},a=f;
    res[0]=1;
    init(n<<1);
    res.resize(tot);
    a.resize(tot);
    ntt(res,1),ntt(a,1);
    for(;b;b>>=1){
        if(b&1){
            for(int i=0;i<tot;i++)res[i]=res[i]*a[i]%P;
            ntt(res,-1);
            for(int i=n;i<tot;i++)res[i]=0;
            ntt(res,1);
        }
        for(int i=0;i<tot;i++)a[i]=a[i]*a[i]%P;
        ntt(a,-1);
        for(int i=n;i<tot;i++)a[i]=0;
        ntt(a,1);
    }
    ntt(res,-1);
    return res;
}

• 都取对数$\ln B(x)=k \ln A(x)$，先算出右边，再exp回去
• 复杂度$O(n \log n)$

poly qpow(const poly &f,int b){
    poly lnf=ln(f);
    for(int i=0;i<lnf.size();i++)lnf[i]=lnf[i]*b%P;
    return exp(lnf,f.size());
}

多项式复合逆

令 $n-1$ 次多项式 $F(x)=\sum\limits _{i=0}^{n-1} a_ix^i$。

给定 $n$ 和 $F(x)$ 的各项系数，要求一个 $n-1$ 次多项式 $G(x)$ 满足：

$$G(F(x))\equiv x\pmod{x^n}$$

求 $G(x)$ 的各项系数对 $998244353$ 取模的结果。

保证 $a_0=0$，$a_1\neq 0$。

生成函数

概念

$g(x)=a_0+a_1x+a_2x^2+...+a_nx^n(-1 \leq x \leq 1)$

常用公式

$A(x)=\sum_{n=0}^{\infty}a_nx^n$

1. 求$\sum_{n=0}^K$的方法为，求$\frac{A(x)}{1-x}$的$x^K$的系数

​ 因为$\frac{1}{1-x}=1+x+x^2+...+x^n$很容易想到

2. 二项式展开：$(a+b)^n=\sum_{k=0}^n\binom{n}{k}a^{n-k}b^k$
3. 负二项式展开：$(1-x)^{-n}=\sum_{k=0}^{\infty}\binom{n+k-1}{k}x^k(|x| \lt 1)$
4. 广义形式：$\binom{-n}{k}=(-1)^k\binom{n+k-1}{k}$,因此$(a+b)^{-n}=\sum_{k=0}^{\infty}(-1)^k\binom{n+k-1}{k}a^{-n-k}b^k$

泰勒展开

普通生成函数（OCF)

1. $\frac{1}{1-x}=1+x+x^2+x^3+...=\sum_{n=0}^{\infty}x^n(a_n=1)$(几何级数)
2. $\frac{1}{1-x^2}=1+x^2+x^4+...$
3. $\frac{1}{1-x^3}=1+x^3+x^6+...$
4. $\frac{1}{(1-x)^2}=1+2x+3x^2+...$
5. $\frac{x}{(1-x)^2}=0+x+2x^2+...(a_n=n)$(自然数序列)
6. $\frac{x(1+x)}{(1-x)^3}=0+x+4x^2+9x^3+...(a_n=n^2)$(平方数序列)
7. $(1+x)^n=\sum_{k=0}^n\binom{n}{k}x^k(a_k=\binom{n}{k})$(二项式系数)
8. $\frac{1}{(1-k)^n}=\sum_{n=0}^{\infty}\binom{n+k-1}{k-1}x^n(a_n=允许重复的选法)$
9. $\frac{x}{1-x-x^2}=0+x+x^2+2x^3+3x^4+...(F_0=0,F_1=1)$(斐波那契数列)
10. $\frac{-\ln(1-x)}{1-x}=\sum_{n=1}^{\infty}H_nx^n(H_n=1+\frac{1}{2}+...+\frac{1}{n})$(调和数)

指数生成函数(EGF)

1. $e^x=\sum_{n=0}^{\infty}\frac{x^n}{n!}$
2. $e^{-x}=1-\frac{x^1}{1}+\frac{x^2}{2!}-\frac{x^3}{3!}+...$
3. $\frac{e^x+e^{-x}}{2}=1+\frac{x^2}{2!}+\frac{x^4}{4!}+...$
4. $\frac{e^x-e^{-x}}{2}=x+\frac{x^3}{3!}+\frac{x^5}{5!}+...$
5. $\frac{1}{1-x}=\sum_{n=0}^{\infty}n! \cdot \frac{x^n}{n!}(a_n=n!)$(排列数)
6. $-ln(1-x)=\sum_{n=1}^{\infty}\frac{(n-1)!x^n}{n!}(a_n=(n-1)!)$(循环排列数)
7. $e^{e^x-1}=\sum_{n=0}^{\infty}B_n \cdot \frac{x^n}{n!}(B_n=贝尔数)$(贝尔数集合划分)
8. $\frac{(e^x-1)^k}{k!}=\sum_{n=k}^{\infty}S(n,k)\cdot \frac{x^n}{n!}(S(n,k)=第二类斯特林数)$
9. $\sin(x)=\sum_{n=0}^{\infty}(-1)^n\frac{x^{2n+1}}{(2n+1)!}$
10. $\cos(x)=\sum_{n=0}^{\infty}(-1)^n\frac{x^{2n}}{(2n)!}$

有穷序列的生成函数

1. $1+x+x^2=\frac{1-x^3}{1-x}$
2. $1+x+x^2+x^3=\frac{1-x^4}{1-x}$
3. $1+x^2+x^4+x^{2m}=\frac{1-x^{2(m+1)}}{1-x^2}$

P10780 BZOJ3028 食物

题目描述

明明这次又要出去旅游了，和上次不同的是，他这次要去宇宙探险！我们暂且不讨论他有多么 NC，他又幻想了他应该带一些什么东西。理所当然的，你当然要帮他计算携带 $n$ 件物品的方案数。

他这次又准备带一些受欢迎的食物，如：蜜桃多啦，鸡块啦，承德汉堡等等。

当然，他又有一些稀奇古怪的限制：

每种食物的限制如下：

• 承德汉堡：偶数个；
• 可乐：$0$ 个或 $1$ 个；
• 鸡腿：$0$ 个，$1$ 个或 $2$ 个；
• 蜜桃多：奇数个；
• 鸡块：$4$ 的倍数个；
• 包子：$0$ 个，$1$ 个，$2$ 个或 $3$ 个；
• 土豆片炒肉：不超过一个；
• 面包：$3$ 的倍数个；

注意，这里我们懒得考虑明明对于带的食物该怎么搭配着吃，也认为每种食物都是以『个』为单位（反正是幻想嘛），只要总数加起来是 $n$ 就算一种方案。因此，对于给出的 $n$，你需要计算出方案数，并对 $10007$ 取模。

提示

• 对于 $40\%$ 的数据，但是 $1\leq n\leq 10^5$；
• 对于所有数据，$1\leq n\leq 10^{500}$；

思路

• 承德汉堡的生成函数：$1+x^2+x^4+...+x^n=\frac{1}{1-x^2}$
• 可乐为：$1+x=\frac{1-x^2}{1-x}$
• 鸡腿为：$1+x+x^2=\frac{1-x^3}{1-x}$
• 蜜桃为：$x+x^3+x^5+...+x^n=\frac{x}{1-x^2}$
• 鸡块为：$1+x^4+x^8+...x^n=\frac{1}{1-x^4}$
• 包子为：$1+x+x^2+x^3+x^4=\frac{1-x^4}{1-x}$
• 土豆片炒肉为：$1+x=\frac{1-x^2}{1-x}$
• 面包为：$1+x^3+x^6+...=\frac{1}{1-x^3}$
• 全部相乘得到 $f(x)=\frac{x}{(1-x)^4}$
• 将式子展开，得到的 $x^n$ 的系数即为方案数
• $f(x)=x(\sum _{n=0}^{∞}C_{n+4-1}^{4-1}x^n)= \sum _{n=0}^{∞}C_{n+3}^{3}x^{n+1}$
• 答案即为 $C_{n+2}^3=\frac{(n+2)(n+1)n}{6}$

多重集组合数问题

有$n$种物品，每种物品有$a_i$个，问取$m$个物品的组合数？

那么构造$n$个生成函数，第$i$个为$1+x^1+x^2+...+x^{a_i}$

所有函数相乘，指数为$m$的系数即为答案

abc248c dice_sum

给定三个整数N,M,K，求多少种序列：

1. A长度为N
2. $1 \leq A_i \leq M(1 \leq i \leq N)$
3. $\sum_{i=1}^NA_i \leq K$

$f(x)=x+x^2+...+x^M=x\frac{1-x^M}{1-x}$

那么要求的是$f(x)^N$的前k项和，转化为求$\frac{f(x)^N}{1-x}$的$x^K$项

$\frac{f(x)^N}{1-x}=x^N(1-x^M)^N(1-x)^{-(N+1)}$

实际上求$(1-x^M)^N(1-x)^{-(N+1)}$的$x^{K-N}$项

$$(1-x^M)^N=\sum_{k=0}^N\binom{N}{k}(-1)^kx^{Mk}\\ (1-x)^{-(N+1)}=\sum_{k=0}^{\infty}\binom{-(N+1)}{k}(-1)^{k}x^k=\sum\binom{N+k}{k}(-1)^{2k}x^k=\sum\binom{N+k}{k}x^k\\ ((1-x^M)^N(1-x)^{-(N+1)})[x^{K-N}]=\sum_{i=0}\binom{N}{i}(-1)^i\binom{N+(K-N-iM)}{K-N-iM}$$

void solve(){
    init();
    int n,m,k;
    cin>>n>>m>>k;
    int ans=0;
    for(int i=0;i*m<=k-n;i++){
        int tmp=C(n,i)*C(k-i*m,k-n-i*m)%mod;
        if(i&1)tmp*=-1;
        ans=(ans+tmp+mod)%mod;
    }
    cout<<ans<<"\n";
}

更无脑的方法是直接分治ntt，复杂度$O(m\log m \log n)$,如果是$x^N(1-x^M)^N(1-x)^{-(N+1)}$，那么可以快速幂+多项式求逆+ntt，复杂度$(O(m\log m))$(过于抽象没有实现)

指数生成函数

$F(x)=\sum_{n \geq 0}a_n\frac{x^n}{n!}$

序列$<1,1,1,...>$的指数生成函数为$1+\frac{x}{1!}+\frac{x^2}{2!}+...=\sum\frac{x^n}{n!}=e^x$

序列$<1,p,p^2,..>$的指数生成函数为$\sum p^n \frac{x^n}{n!}=e^{px}$

多重集排列数问题

有$n$个物品，每种物品$a_i$个，问$m$个物品的排列数

那么答案就是$(1+\frac{x^1}{1}+\frac{x^2}{2}+...\frac{x^{a_i}}{a_i!})$所有的和，$\frac{x^m}{m!}$的系数

POJ3734 Blocks

长度$n$的序列，四种颜色，红色和绿色数量必须为偶数，求方案数

红色和绿色：$1+\frac{x^2}{2!}+\frac{x^4}{4!}+...=\frac{e^x+e^{-x}}{2}$

蓝色黄色：$1+\frac{x}{1!}+\frac{x^2}{2!}+...=e^x$

$(\frac{e^x+e^{-x}}{2})^2e^{2x}=\frac{e^{4x}+2e^{2x}+1}{4}$

因为$\sum_{i=0}\frac{k^ix^i}{i!}=e^{kx}$

因此$x^n$系数为$\frac{4^n+2\cdot 2^n}{4n!}$,$\frac{x^n}{n!}$系数为$\frac{4^n+2\cdot 2^n}{4}$

牛客17316(挑选队友)(分治FFT)

共有$n$位选手，$m$个群，$k$名选手被选择，问方案数

$s_i$表示该群的人数，可以发现该群的$x^k$的多项式系数为$C_{s_i}^k$，$m$个数卷积即可

使用分治fft

$T(n)=2T(\frac{n}{2})+O(n\log n)$,复杂度为$O(n \log^2 n)$

void solve(){
    int n,m,k;
    cin>>n>>m>>k;
    comb.init(n);
    vector<int> s(n+1);
    for(int i=1;i<=m;i++)cin>>s[i];
    auto dfs=[&](auto &&dfs,int l,int r){
        if(l==r){
            poly f(s[l]+1);
            for(int i=1;i<=s[l];i++)f[i]=comb.binom(s[l],i);
            return f;
        }
        int mid=l+r>>1;
        return mul(dfs(dfs,l,mid),dfs(dfs,mid+1,r));
    };
    auto res=dfs(dfs,1,m);
    cout<<res[k]<<"\n";
}

注意mul函数处的init(n+m+1)改成init(n+m)，否则复杂度是错的

abc234f Reordering(指数型)

给定字符串$S$，求他的所有排列能生成的子序列个数

26个字母，每个字母出现次数为$a_i$个，该字母的生成函数为$1+x+\frac{x^2}{2}+...+\frac{x^{a_i}}{a_i!}$

所有26个字母相乘

答案为$\sum_{i=1}^{n}i!\cdot f_i$

void solve(){
    string s;
    cin>>s;
    int n=s.size();
    comb.init(n);
    poly res{1};
    vector<int> a(26);
    for(auto c:s)a[c-'a']++;
    for(int i=0;i<26;i++){
        if(!a[i])continue;
        poly tmp(a[i]+1,1);
        for(int j=1;j<=a[i];j++)tmp[j]=comb.inv[j];
        res.resize(n+1);
        res=mul(res,tmp);
    }
    int ans=0;
    for(int i=1;i<=n;i++)ans=(ans+comb.fac[i]*res[i])%mod;
    cout<<ans<<"\n";
}

CF632F

$n$种物品，价值为$a_i$，求拿走恰好$k$种物品的所有可能价值，$1 \leq n,a_i,k \leq 1000$

定义$f_{i,j}$表示选择了$i$个物品，价值为$j$的方案数，$f_{i,j}=\min f_{i-1,j-a_k}$转移$O(n^4)$

换个思路，定义$f_i$表示价值为$i$的最小物品数

假设价值最小值为$x$，那么所有物品$a_i=a_i-x$，那么此时即使$f_i \lt k$,可以添加价值为0的物品数，答案就是$i+xk$

另一种思路：

定义$g_i$为价值的出现次数，$g(x)=1+b_1x^1+b_2x^2+...$

那么$f_{i,j}=f_{i-1,j-k}*g_k$是一个卷积形式，NTT求解(不卡模数的情况下)

CF1096G

一个$n$位数，每位可以是给出的$k$个数码的一个数，可以有前导$0$，输出前$\frac{n}{2}$和后$\frac{n}{2}$相等的方案数

设出现的数位集合为$S$，那么生成函数$F(X)=\sum_{i \in S}x^i$，答案为$\sum [x^i]F^{\frac{n}{2}}(x)$

GYM102920H

三条水平线a,b,c，每个水平线有$n$个点，高度相差1，求三点共线数量。

$a_i+c_j=2b_k$，构造生成函数，$F(x)=\sum_{i \in S}x^i$，$S$是值的集合，答案就是$A(x)C(x)$中$x^{2k}$的系数，前提$b_k=1$

做一遍FFT即可

HDU7328 Snake

有$n$个长度为1的蛇，求集合个数(比如$n=2$集合有{1->2},{2->1},{1,2},{2,1}，最终有$m$条蛇，每条蛇长度最多为$k$

容斥做法：

考虑没有排列的因素，$f(i)$表示至少$i$条蛇大于$k$的方案

$ans=f_0-f_1+f_2....$

• 对于$i$，选出$i$个集合等于$k$的，$\binom{m}{i}$
• 剩余$n-ik$分成$m$份，$\binom{n-ik-1}{m-1}$

考虑排序，$\frac{ans \cdot n!}{m!}$全排列，但是集合不排

void solve(){
    int n,m,k;
    cin>>n>>m>>k;
    vector<int> f(n+1);
    for(int i=0;i<=m;i++){
        f[i]=comb.binom(m,i)*comb.binom(n-i*k-1,m-1)%mod;
    }
    int ans=0;
    for(int i=0;i<=m;i++){
        if(i&1)ans=(ans+mod-f[i])%mod;
        else ans=(ans+f[i])%mod;
    }
    cout<<ans*comb.fac[n]%mod*comb.inv[m]%mod<<"\n";
}

生成函数做法：

每个集合的生成函数$F(x)=x^1+x^2+...+x^k=x \frac{1-x^k}{1-x}$

答案就是$F^m(x)[x^n]$，化简和abc248c一样

最后答案为$\frac{ans \cdot n!}{m!}$

2017CCPC 杭州Marriage(容斥+生成函数)

有$n$个家庭，每个家庭有$a_i$个男孩和$b_i$个女孩，$\sum a_i=\sum b_i$，除了同家庭，其余男女配对，问方案数

考虑容斥，$f_i$表示至少$k$对同家庭的

设$p=\min(a_i,b_i)$为当前家庭的贡献，设贡献$i$对，方案为$c_i=\binom{a}{i} \binom{b}{i}i!$，生成函数为$\sum_{i=0}^p c_ix^i$

那么所有生成函数相乘，$f_i=[x^i]$

容斥一下计算答案即可

复杂度$O(n \log^2n)$

分治NTT需要用优先队列启发式合并实现

void solve(){
    int n;
    cin>>n;
    vector<int> a(n+1),b(n+1);
    for(int i=1;i<=n;i++)cin>>a[i]>>b[i];
    int s=accumulate(a.begin()+1,a.end(),0LL);
    int t=0;
    for(int i=1;i<=n;i++)t+=min(a[i],b[i]);
    vector<poly> f(n+1);
    for(int i=1;i<=n;i++){
        int p=min(a[i],b[i]);
        f[i].resize(p+1);
        for(int j=0;j<=p;j++)f[i][j]=comb.binom(a[i],j)*comb.binom(b[i],j)%mod*comb.fac[j]%mod;
    }
    auto cmp=[&](int x,int y){
        return f[x].size()>f[y].size();
    };
    priority_queue<int,vector<int>,decltype(cmp)> pq(cmp);
    for(int i=1;i<=n;i++)pq.push(i);
    while(pq.size()>1){
        auto x=pq.top();pq.pop();
        auto y=pq.top();pq.pop();
        int z=f.size();
        auto g=mul(f[x],f[y]);
        g.resize(min((int)g.size(),t+1));
        f.push_back(g);
        pq.push(z);
    }
    auto F=f.back();
    int ans=0;
    for(int i=0;i<=t;i++){
        int tmp=F[i]*comb.fac[s-i]%mod;
        if(i&1)tmp=mod-tmp;
        ans=(ans+tmp)%mod;
    }
    cout<<ans<<'\n';
}

CCPC桂林B

求$A_i=(\sum_{j=i}^{\min(i+L-1,n)}A_j) \mod P$，执行$m$次

将数组$A$倒序，构造数组$G(x)=x^0+x^1+x^2+...+x^L$

每执行一次操作即$A(x)*G(x)$

答案是$A(x)*G^m(x)$

void solve(){
    int _=0;
    cin>>_;
    for(int T=1;T<=_;T++){
        int n,L,m;
        cin>>n>>L>>m;
        poly f(n);
        for(int i=0;i<n;i++)cin>>f[i];
        ranges::reverse(f);
        poly g(L,1);
        g=qpow(g,m,n);
        
        f=mul(f,g);
        f.resize(n);
        cout<<"Case "<<T<<": ";
        for(int i=n-1;i>=0;i--)cout<<f[i]<<" ";
        cout<<"\n";
    }
}

2018沈阳M

多重背包问题，每个位置有数量为$a$价值为$b$的物品

多次在线询问，求$[L,R]$选够价值$c$个物品的方案

每个位置的生成函数为$1+x^b+x^{2b}+...+x^{ab}=\frac{1-x^{(a+1)b}}{1-x^b}$

答案为$\prod_{i=l}^r\frac{1-x^{(a+1)b}}{1-x^b}=\prod_{i=1}^r\frac{x^{(a+1)b}-1}{x^b-1}\prod_{i=1}^{l-1}\frac{x^b-1}{x^{(a+1)b}-1}$

void solve(){
    int _;
    cin>>_;
    for(int T=1;T<=_;T++){
        cout<<"Case #"<<T<<":\n";
        int n,m;
        cin>>n>>m;
        const int M=2001;
        vector<int> f(M),g(M);
        vector sf(n+1,vector<int>(M)),sg(n+1,vector<int>(M));
        f[0]=g[0]=sf[0][0]=sg[0][0]=1;
        for(int i=1;i<M;i++)sf[0][i]=1;
        auto _01=[&](vector<int> &dp,int x){
            for(int i=M-1;i>=x;i--)dp[i]=(dp[i]-dp[i-x]+mod)%mod;
        };
        auto _complete=[&](vector<int> &dp,int x){
            for(int i=x;i+x<M;i++)dp[i]=(dp[i]+dp[i-x])%mod;
        };
        for(int i=1;i<=n;i++){
            int a,b;
            cin>>a>>b;
            _01(f,(a+1)*b);
            _01(g,b);
            _complete(f,b);
            _complete(g,(a+1)*b);
            sf[i][0]=f[0];
            sg[i][0]=g[0];
            for(int j=1;j<M;j++){
                sf[i][j]=(sf[i][j-1]+f[j])%mod;
                sg[i][j]=g[j];
            }
        }
        int ans=0;
        while(m--){
            int l,r,c;
            cin>>l>>r>>c;
            l=(l+ans)%n+1,r=(r+ans)%n+1;
            if(l>r)swap(l,r);
            ans=0;
            for(int i=0;i<=c;i++)ans=(ans+sf[r][i]*sg[l-1][c-i]%mod)%mod;
            cout<<ans<<"\n";
        }
    }
}

2020ICPC澳门