多项式八分之三家桶

因为目前半家桶里只会俩所以就成四分之一了。

UPD:又会了一个,现在是八分之三了。

UUPD:虽然四个齐了,但你不感觉八分之三很带劲吗(

UUUPD:就当多项式全家桶写就完了(虽然不可能写完

多项式乘法逆

设多项式 \(A(x)\)\(n\) 次多项式,求多项式 \(B\) ,使 \(A(x)*B(x)\equiv 1\pmod{x^n}\)

\(A(x)*C(x)\equiv 1\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\) ,那么

\[A(x)*B(x)\equiv 1\pmod{x^n}\to A(x)*B(x)\equiv 1\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\\ \; \\A(x)*(B(x)-C(x))\equiv 0\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\\ \; \\B(x)-C(x)\equiv 0\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\to (B(x)-C(x))^2\equiv 0\pmod{x^n}\\ \; \\B^2(x)-2B(x)C(x)+C^2(x)\equiv 0\pmod{x^n}\\ \; \\B(x)-2C(x)+A(x)C^2(x)\equiv 0\pmod{x^n}\to B(x)\equiv2C(x)-A(x)C^2(x)\pmod{x^n}
\]

递归求解。由主定理复杂度 \(O(n\log n)\)

洛谷P4238[模板]多项式乘法逆
##include<bits/stdc++.h>
using namespace std;

namespace IO{
    typedef long long LL;
    int read(){
        int x=0,f=0; char ch=getchar();
        while(ch<'0'||ch>'9'){ f|=ch=='-'; ch=getchar(); }
        while(ch>='0'&&ch<='9'){ x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
        return f?-x:x;
    } char outp[50];
    void write(int x,char sp,int len=0){
        if(x<0) putchar('-'), x=-x;
        do{ outp[len++]=x%10+'0'; x/=10; }while(x);
        for(int i=len-1;~i;i--) putchar(outp[i]); putchar(sp);
    }
    void ckmin(int& x,int y){ x=x<y?x:y; }
    void ckmax(int& x,int y){ x=x>y?x:y; }
} using namespace IO;

const int NN=400010,mod=998244353;
int n;

namespace Poly_Calculation{
    const LL rt=3,two=998244355;
    int l,len,r[NN];
    LL a[NN],b[NN],g[NN],tp[NN];
    void pls(LL& x,LL y){ (x+=y)>=mod?(x-=mod):(x<0?(x+=mod):0); }
    LL qpow(LL a,int b=mod-2,LL res=1){
        for(;b;b>>=1,a=a*a%mod)
            if(b&1) res=res*a%mod;
        return res;
    }
    void prework(int n){
        for(len=1,l=0;len<=n+n;len<<=1) ++l;
        for(int i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<l-1);
        g[0]=1; g[1]=qpow(rt,(mod-1)/len);
        for(int i=2;i<len;i++) g[i]=g[i-1]*g[1]%mod;
    }
    void ntt(LL* a,int n){
        for(int i=0;i<n;i++) if(i<r[i]) swap(a[i],a[r[i]]);
        for(int t=n>>1,d=1;d<n;d<<=1,t>>=1)
            for(int i=0;i<n;i+=(d<<1)) for(int j=0;j<d;j++){
                LL tmp=g[j*t]*a[i+j+d]%mod;
                pls(a[i+j+d]=a[i+j],-tmp); pls(a[i+j],tmp);
            }
    }
    void polyinv(LL* a,LL* res,int deg){
        if(deg==1) return void(res[0]=qpow(a[0]));
        polyinv(a,res,deg+1>>1); prework(deg);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        for(int i=deg;i<len;i++) tp[i]=0;
        ntt(res,len); ntt(tp,len);
        for(int i=0;i<len;i++) res[i]=(two-res[i]*tp[i]%mod)*res[i]%mod, g[i]=qpow(g[i]);
        ntt(res,len); LL inv=qpow(len);
        for(int i=0;i<deg;i++) res[i]=res[i]*inv%mod;
        for(int i=deg;i<len;i++) res[i]=0;
    }
} using namespace Poly_Calculation;

signed main(){
    n=read();
    for(int i=0;i<n;i++) a[i]=read();
    polyinv(a,b,n);
    for(int i=0;i<n;i++) write(b[i],' ');
    return puts(""),0;
}

多项式开根

已知 \(n\) 多项式 \(A(x)\) ,求多项式 \(B(x)\) ,满足 \(B^2(x)\equiv A(x)\pmod{x^n}\)

类似地,设 \(C^2(x)\equiv A(x)\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\) ,那么

\[B(x)\equiv C(x)\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\\ \; \\B(x)-C(x)\equiv 0\pmod{x^{\left\lceil\frac{n}{2}\right\rceil}}\to (B(x)-C(x))^2\equiv 0\pmod{x^n}\\ \; \\B^2(x)-2B(x)C(x)+C^2(x)\equiv 0\pmod{x^n}\\ \; \\A(x)-2B(x)C(x)+C^2(x)\equiv 0\pmod {x^n}\\ \; \\B(x)=\equiv\frac{C(x)+\frac{A(x)}{C(x)}}{2}\pmod{x^n}
\]

递归并求逆求解。

这份板子因为某种神秘力量常数过丑了。有时间的话应该会来改改吧。(

UPD:取模优化去掉开O2跑快了三秒。。不懂什么原理。虽然还是挺慢但至少能看了。先这样吧

洛谷P5205[模板]多项式开根
##include<bits/stdc++.h>
using namespace std;

namespace IO{
    typedef long long LL;
    int read(){
        int x=0,f=0; char ch=getchar();
        while(ch<'0'||ch>'9'){ f|=ch=='-'; ch=getchar(); }
        while(ch>='0'&&ch<='9'){ x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
        return f?-x:x;
    } char outp[50];
    void write(int x,char sp,int len=0){
        if(x<0) putchar('-'), x=-x;
        do{ outp[len++]=x%10+'0'; x/=10; }while(x);
        for(int i=len-1;~i;i--) putchar(outp[i]); putchar(sp);
    }
    void ckmin(int& x,int y){ x=x<y?x:y; }
    void ckmax(int& x,int y){ x=x>y?x:y; }
} using namespace IO;

const int NN=400010,mod=998244353;
int n;

namespace Poly_Calculation{
    const LL rt=3,two=998244355,inv2=499122177;
    int l,len,r[NN];
    LL a[NN],b[NN],g[NN],iv[NN],tp[NN];
    void pls(LL& x,LL y){ (x+=y)>=mod?(x-=mod):(x<0?(x+=mod):0); }
    LL qpow(LL a,int b=mod-2,LL res=1){
        for(;b;b>>=1,a=a*a%mod)
            if(b&1) res=res*a%mod;
        return res;
    }
    void prework(int n){
        for(len=1,l=0;len<=(n<<1);len<<=1) ++l;
        for(int i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<l-1);
        g[0]=1; g[1]=qpow(rt,(mod-1)/len);
        for(int i=2;i<len;i++) g[i]=g[i-1]*g[1]%mod;
    }
    void ntt(LL* a,int n){
        for(int i=0;i<n;i++) if(i<r[i]) swap(a[i],a[r[i]]);
        for(int t=n>>1,d=1;d<n;d<<=1,t>>=1)
            for(int i=0;i<n;i+=(d<<1)) for(int j=0;j<d;j++){
                LL tmp=g[j*t]*a[i+j+d]%mod;
                a[i+j+d]=(a[i+j]-tmp+mod)%mod;
                a[i+j]=(a[i+j]+tmp)%mod;
            }
    }
    void polyinv(LL* a,LL* res,int deg){
        if(deg==1) return void(res[0]=qpow(a[0]));
        polyinv(a,res,deg+1>>1); prework(deg);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        for(int i=deg;i<len;i++) tp[i]=0;
        ntt(res,len); ntt(tp,len);
        for(int i=0;i<len;i++) res[i]=(two-res[i]*tp[i]%mod)*res[i]%mod, g[i]=qpow(g[i]);
        ntt(res,len); LL inv=qpow(len);
        for(int i=0;i<deg;i++) res[i]=res[i]*inv%mod;
        for(int i=deg;i<len;i++) res[i]=0;
    }
    void polysqr(LL* a,LL *res,int deg){
        if(deg==1) return void(res[0]=1);
        polysqr(a,res,deg+1>>1);
        for(int i=0;i<len;i++) iv[i]=0;
        polyinv(res,iv,deg); prework(deg);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        for(int i=deg;i<len;i++) tp[i]=0;
        ntt(res,len); ntt(iv,len); ntt(tp,len);
        for(int i=0;i<len;i++) res[i]=(tp[i]*iv[i]%mod+res[i])*inv2%mod, g[i]=qpow(g[i]);
        ntt(res,len); LL inv=qpow(len);
        for(int i=0;i<deg;i++) res[i]=res[i]*inv%mod;
        for(int i=deg;i<len;i++) res[i]=0;
    }
} using namespace Poly_Calculation;

signed main(){
    n=read();
    for(int i=0;i<n;i++) a[i]=read();
    polysqr(a,b,n);
    for(int i=0;i<n;i++) write(b[i],' ');
    return puts(""),0;
}

多项式对数函数

已知 \(n\) 次多项式 \(A\) ,求多项式 \(B\) ,满足 \(B(x)\equiv \ln(A(x))\pmod{x^n}\)

\(\ln\) 很ex,求导化掉。

\[B'(x)=\ln'(A(x))A'(x)=\frac{A'(x)}{A(x)}
\]

求乘法逆之后求出导数再积分回去。

不定积分和求导是互逆的,

\[(x^a)'=ax^{x-1}\\ \; \\ \int x^adx=\frac{1}{a+1}x^{a+1}
\]

洛谷P4725[模板]多项式对数函数
##include<bits/stdc++.h>
using namespace std;

namespace IO{
    typedef long long LL;
    int read(){
        int x=0,f=0; char ch=getchar();
        while(ch<'0'||ch>'9'){ f|=ch=='-'; ch=getchar(); }
        while(ch>='0'&&ch<='9'){ x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
        return f?-x:x;
    } char outp[50];
    void write(int x,char sp,int len=0){
        if(x<0) putchar('-'), x=-x;
        do{ outp[len++]=x%10+'0'; x/=10; }while(x);
        for(int i=len-1;~i;i--) putchar(outp[i]); putchar(sp);
    }
    void ckmin(int& x,int y){ x=x<y?x:y; }
    void ckmax(int& x,int y){ x=x>y?x:y; }
} using namespace IO;

const int NN=600010,mod=998244353;
int n;

namespace Poly_Calculation{
    const LL rt=3,two=mod+2;
    int l,len,r[NN];
    LL a[NN],b[NN],g[NN],iv[NN],tp[NN];
    LL qpow(LL a,int b=mod-2,LL res=1){
        for(;b;b>>=1,a=a*a%mod)
            if(b&1) res=res*a%mod;
        return res;
    }
    void prework(int n){
        for(len=1,l=0;len<=(n<<1);len<<=1) ++l;
        for(int i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<l-1);
        g[0]=1; g[1]=qpow(rt,(mod-1)/len);
        for(int i=2;i<len;i++) g[i]=g[i-1]*g[1]%mod;
    }
    void ntt(LL* a,int n){
        for(int i=0;i<n;i++) if(r[i]>i) swap(a[i],a[r[i]]);
        for(int t=n>>1,d=1;d<n;d<<=1,t>>=1)
            for(int i=0;i<n;i+=(d<<1)) for(int j=0;j<d;j++){
                LL tmp=g[t*j]*a[i+j+d]%mod;
                a[i+j+d]=(a[i+j]-tmp+mod)%mod; a[i+j]=(a[i+j]+tmp)%mod;
            }
    }
    void deriv(LL* a,LL* res,int n){ for(int i=1;i<n;i++) res[i-1]=a[i]*i%mod; res[n-1]=0; }
    void integ(LL* a,LL* res,int n){ for(int i=1;i<n;i++) res[i]=a[i-1]*qpow(i)%mod; res[0]=0; }
    void polyinv(LL* a,LL* res,int deg){
        if(deg==1) return void(res[0]=qpow(a[0]));
        polyinv(a,res,deg+1>>1); prework(deg);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        for(int i=deg;i<len;i++) tp[i]=0;
        ntt(tp,len); ntt(res,len);
        for(int i=0;i<len;i++) res[i]=(two-res[i]*tp[i]%mod)*res[i]%mod, g[i]=qpow(g[i]);
        ntt(res,len); LL inv=qpow(len);
        for(int i=0;i<deg;i++) res[i]=res[i]*inv%mod;
        for(int i=deg;i<len;i++) res[i]=0;
    }
    void polyln(LL* a,LL* res,int deg){
        polyinv(a,iv,deg); prework(deg);
        for(int i=0;i<len;i++) tp[i]=0;
        deriv(a,tp,deg);
        ntt(iv,len); ntt(tp,len);
        for(int i=0;i<len;i++) tp[i]=tp[i]*iv[i]%mod, g[i]=qpow(g[i]);
        ntt(tp,len); LL inv=qpow(len);
        for(int i=0;i<len;i++) tp[i]=tp[i]*inv%mod;
        integ(tp,res,deg);
    }
} using namespace Poly_Calculation;

signed main(){
    n=read();
    for(int i=0;i<n;i++) a[i]=read();
    polyln(a,b,n);
    for(int i=0;i<n;i++) write(b[i],' ');
    return puts(""),0;
}

多项式指数函数

已知 \(n\) 次多项式 \(A(x)\) ,求 \(n\) 次多项式 \(B(x)\) ,满足 \(B(x)\equiv e^{A(x)}\pmod{x^n}\)

仍考虑递归求解。如果已知函数 \(G(x)\) ,求满足 \(G(F(z))\equiv 0(\mod z^n)\)\(n\) 次多项式 \(F(z)\) ,可以用多项式牛顿迭代。

洛谷题解第一篇的图

现在令 \(G(x)=\ln B(x)-A(x)\) ,令 \(B_0(x)\)\(\mod x^{\left\lceil\frac{n}{2}\right\rfloor}\) 意义下的解,那么就是要求 \(G\) 的一个零点。

\(G\) 求导,把 \(A(x)\) 看作常数,那么 \(G'(B(x))=\frac{1}{B(x)}\) ,代入多项式牛顿迭代,有

\[B(x)\equiv B_0(x)-\frac{G(F_0(x))}{G'(F_0(x))}\pmod{x^n}\\\;\\B(x)\equiv B_0(x)(1-\ln B_0(x)+A(x))\pmod{x^n}
\]

\(\ln\) 并递归求解。

洛谷P4726[模板]多项式指数函数
##include<bits/stdc++.h>
using namespace std;

namespace IO{
    typedef long long LL;
    int read(){
        int x=0,f=0; char ch=getchar();
        while(ch<'0'||ch>'9'){ f|=ch=='-'; ch=getchar(); }
        while(ch>='0'&&ch<='9'){ x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
        return f?-x:x;
    } char outp[50];
    void write(int x,char sp,int len=0){
        if(x<0) putchar('-'), x=-x;
        do{ outp[len++]=x%10+'0'; x/=10; }while(x);
        for(int i=len-1;~i;i--) putchar(outp[i]); putchar(sp);
    }
    void ckmin(int& x,int y){ x=x<y?x:y; }
    void ckmax(int& x,int y){ x=x>y?x:y; }
} using namespace IO;

const int NN=600010,mod=998244353;
int n;

namespace Poly_Calculation{
    const LL rt=3,one=mod+1,two=mod+2;
    int l,len,r[NN];
    LL a[NN],g[NN],b[NN],iv[NN],tp[NN],ln[NN];
    LL qpow(LL a,int b=mod-2,LL res=1){
        for(;b;b>>=1,a=a*a%mod)
            if(b&1) res=res*a%mod;
        return res;
    }
    void prework(int n){
        for(len=1,l=0;len<(n<<1);len<<=1) ++l;
        for(int i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<l-1);
        g[0]=1; g[1]=qpow(rt,(mod-1)/len);
        for(int i=2;i<len;i++) g[i]=g[i-1]*g[1]%mod;
    }
    void ntt(LL* a,int n){
        for(int i=0;i<n;i++) if(r[i]>i) swap(a[i],a[r[i]]);
        for(int t=n>>1,d=1;d<n;d<<=1,t>>=1)
            for(int i=0;i<n;i+=(d<<1)) for(int j=0;j<d;j++){
                LL tmp=g[t*j]*a[i+j+d]%mod;
                a[i+j+d]=(a[i+j]-tmp+mod)%mod; a[i+j]=(a[i+j]+tmp)%mod;
            }
    }
    void polyinv(LL* a,LL* res,int deg){
        if(deg==1) return void(res[0]=qpow(a[0]));
        polyinv(a,res,deg+1>>1); prework(deg);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        ntt(tp,len); ntt(res,len);
        for(int i=0;i<len;i++) res[i]=(two-res[i]*tp[i]%mod)*res[i]%mod, g[i]=qpow(g[i]);
        ntt(res,len); LL inv=qpow(len);
        for(int i=0;i<deg;i++) res[i]=res[i]*inv%mod, tp[i]=0;
        for(int i=deg;i<len;i++) res[i]=tp[i]=0;
    }
    void deriv(LL* a,LL* res,int n){ for(int i=1;i<n;i++) res[i-1]=a[i]*i%mod; res[n-1]=0; }
    void integ(LL* a,LL* res,int n){ for(int i=1;i<n;i++) res[i]=a[i-1]*qpow(i)%mod; res[0]=0; }
    void polyln(LL* a,LL* res,int deg){
        polyinv(a,iv,deg); deriv(a,tp,deg); prework(deg);
        ntt(iv,len); ntt(tp,len);
        for(int i=0;i<len;i++) tp[i]=tp[i]*iv[i]%mod, g[i]=qpow(g[i]);
        ntt(tp,len); LL inv=qpow(len);
        for(int i=0;i<deg;i++) tp[i]=tp[i]*inv%mod;
        integ(tp,res,deg);
        for(int i=0;i<len;i++) tp[i]=iv[i]=0;
    }
    void polyexp(LL* a,LL* res,int deg){
        if(deg==1) return void(res[0]=1);
        polyexp(a,res,deg+1>>1); polyln(res,ln,deg); prework(deg);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        ntt(tp,len); ntt(res,len); ntt(ln,len);
        for(int i=0;i<len;i++) res[i]=(one-ln[i]+tp[i])*res[i]%mod, g[i]=qpow(g[i]);
        ntt(res,len); LL inv=qpow(len);
        for(int i=0;i<deg;i++) res[i]=res[i]*inv%mod, tp[i]=ln[i]=0;
        for(int i=deg;i<len;i++) res[i]=tp[i]=ln[i]=0;
    }
} using namespace Poly_Calculation;

signed main(){
    n=read();
    for(int i=0;i<n;i++) a[i]=read();
    polyexp(a,b,n);
    for(int i=0;i<n;i++) write(b[i],' ');
    return puts(""),0;
}

多项式半家桶

namespace 套起来能算封装吗(

洛谷P4389付公主的背包
#include<bits/stdc++.h>
#define int long long
using namespace std;

namespace IO{
    typedef unsigned long long ULL;
    typedef double DB; typedef long long LL;
    int read(){
        int x=0,f=0; char ch=getchar();
        while(ch<'0'||ch>'9'){ f|=ch=='-'; ch=getchar(); }
        while(ch>='0'&&ch<='9'){ x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
        return f?-x:x;
    } char outp[50];
    void write(int x,char sp,int len=0){
        if(x<0) putchar('-'), x=-x;
        do{ outp[len++]=(x%10)^48; x/=10; }while(x);
        for(int i=len-1;~i;i--) putchar(outp[i]); putchar(sp);
    }
    void ckmin(int& x,int y){ x=x<y?x:y; }
    void ckmax(int& x,int y){ x=x>y?x:y; }
} using namespace IO;

const int NN=600010,mod=998244353;
int n,m,buc[NN],inv[NN];
int qpow(int a,int b=mod-2,int res=1){
    for(;b;b>>=1,a=a*a%mod)
        if(b&1) res=res*a%mod;
    return res;
}

namespace PolyCalculation{
    const int rt=3,one=mod+1,two=mod+2,inv2=mod+1>>1;
    int l,len,r[NN],a[NN],g[NN],tp[NN],iv[NN],ln[NN],res[NN];
    void revg(){ g[1]=qpow(g[1]); for(int i=2;i<len;i++) g[i]=g[i-1]*g[1]%mod; }
    void deriv(int* a,int* r,int n){ for(int i=1;i<n;i++) r[i-1]=i*a[i]%mod; r[n-1]=0; }
    void integ(int* a,int* r,int n){ for(int i=1;i<n;i++) r[i]=a[i-1]*qpow(i)%mod; r[0]=0; }
    void prework(int n){
        for(len=1,l=0;len<=n;len<<=1) ++l;
        for(int i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<l-1);
        g[0]=1; g[1]=qpow(rt,(mod-1)/len);
        for(int i=2;i<len;i++) g[i]=g[i-1]*g[1]%mod;
    }
    void ntt(int* a,int n){
        for(int i=0;i<n;i++) if(r[i]>i) swap(a[i],a[r[i]]);
        for(int t=n>>1,d=1;d<n;t>>=1,d<<=1)
            for(int i=0;i<n;i+=(d<<1)) for(int j=0;j<d;j++){
                int tmp=g[t*j]*a[i+j+d]%mod;
                a[i+j+d]=(a[i+j]-tmp+mod)%mod; a[i+j]=(a[i+j]+tmp)%mod;
            }
    }
    void polyinv(int* a,int* r,int deg){
        if(deg==1) return void(r[0]=qpow(a[0]));
        polyinv(a,r,deg+1>>1); prework(deg<<1);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        ntt(tp,len); ntt(r,len);
        for(int i=0;i<len;i++) r[i]=(two-r[i]*tp[i]%mod)*r[i]%mod;
        revg(); ntt(r,len); int inv=qpow(len);
        for(int i=0;i<deg;i++) r[i]=r[i]*inv%mod, tp[i]=0;
        for(int i=deg;i<len;i++) r[i]=tp[i]=0;
    }
    void polysqr(LL* a,LL *r,int deg){
        if(deg==1) return void(r[0]=1);
        polysqr(a,r,deg+1>>1); polyinv(r,iv,deg); prework(deg<<1);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        ntt(r,len); ntt(iv,len); ntt(tp,len);
        for(int i=0;i<len;i++) r[i]=(tp[i]*iv[i]%mod+r[i])*inv2%mod;
        revg(); ntt(r,len); LL inv=qpow(len);
        for(int i=0;i<deg;i++) r[i]=r[i]*inv%mod, iv[i]=tp[i]=0;
        for(int i=deg;i<len;i++) r[i]=0, iv[i]=tp[i]=0;
    }
    void polyln(int* a,int* r,int deg){
        polyinv(a,iv,deg); deriv(a,tp,deg);
        prework(deg<<1); ntt(tp,len); ntt(iv,len);
        for(int i=0;i<len;i++) tp[i]=tp[i]*iv[i]%mod;
        revg(); ntt(tp,len); int inv=qpow(len);
        for(int i=0;i<len;i++) tp[i]=tp[i]*inv%mod;
        integ(tp,r,deg);
        for(int i=0;i<len;i++) tp[i]=iv[i]=0;
    }
    void polyexp(int* a,int* r,int deg){
        if(deg==1) return void(r[0]=1);
        polyexp(a,r,deg+1>>1); polyln(r,ln,deg); prework(deg<<1);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        ntt(tp,len); ntt(ln,len); ntt(r,len);
        for(int i=0;i<len;i++) r[i]=(one-ln[i]+tp[i])*r[i]%mod;
        revg(); ntt(r,len); int inv=qpow(len);
        for(int i=0;i<deg;i++) r[i]=r[i]*inv%mod, ln[i]=tp[i]=0;
        for(int i=deg;i<len;i++) r[i]=ln[i]=tp[i]=0;
    }
} using namespace PolyCalculation;

signed main(){
    n=read(); m=read(); inv[0]=inv[1]=1;
    for(int i=1;i<=n;i++) ++buc[read()];
    for(int i=2;i<=m;i++) inv[i]=mod-mod/i*inv[mod%i]%mod;
    for(int i=1;i<=m;i++) if(buc[i])
        for(int j=i,k=1;j<=m;j+=i,k++)
            a[j]=(a[j]+buc[i]*inv[k]);
    polyexp(a,res,m+1);
    for(int i=1;i<=m;i++) write(res[i],'\n');
    return 0;
}

多项式另八分之五家桶(

任意模数多项式乘法

出题人毒瘤模数不符合 \(\tt{NTT}\) 要求时就要写这个。有两种写法,一是三模数 \(\tt{NTT}\) ,一是拆系数 \(\tt{FFT}\) 。前者要做变换的次数较多,而后者经过一系列NB优化后一共只用四次 \(\tt{FFT}\) ,因此学了后者。

三模数 \(\tt{NTT}\) 思路大体是用三个 \(\tt{NTT}\) 模数通过中国剩余定理合并出一个大于运算中最大值的超大模数,合并答案后取模。好像要 \(9\)\(\tt{NTT}\)

拆系数 \(\tt{FFT}\) (\(\tt{MTT}\))将多项式 \(F(x)=\sum_{i=0}^na_ix^i\) 拆为 \(F(x)=\sum_{i=0}^n(b_i+c_i)x^i\) ,其中 \(b_i=\left\lfloor\frac{a_i}{\sqrt{mod}}\right\rfloor\)\(c_i=a_i\mod\sqrt{mod}\) 。实际运算中 \(\sqrt{mod}\) 可以直接取 \(2^{15}\) 。化简式子后可以 \(7\)\(\tt{FFT}\) ,再发掘一些共轭与 \(\tt{DFT}\) 的性质后可以优化到 \(4\)\(\tt{FFT}\)具体可见洛谷题解

洛谷P4245[模板]任意模数多项式乘法
#include<bits/stdc++.h>
using namespace std;

namespace IO{
    typedef long double DB;
    typedef long long LL;
    int read(){
        int x=0,f=0; char ch=getchar();
        while(ch<'0'||ch>'9'){ f|=ch=='-'; ch=getchar(); }
        while(ch>='0'&&ch<='9'){ x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
        return f?-x:x;
    } char outp[50];
    void write(int x,char sp,int len=0){
        if(x<0) putchar('-'), x=-x;
        do{ outp[len++]=x%10+'0'; x/=10; }while(x);
        for(int i=len-1;~i;i--) putchar(outp[i]); putchar(sp);
    }
    void ckmin(int& x,int y){ x=x<y?x:y; }
    void ckmax(int& x,int y){ x=x>y?x:y; }
} using namespace IO;

const int NN=800010;
int n,m,mod,f[NN],g[NN],res[NN];

namespace Poly_Calculation{
    const DB PI=acos(-1.0);
    struct compx{
        DB r,i;
        compx(){} compx(DB x,DB y){ r=x; i=y; }
        compx conj(){ return compx(r,-i); }
        compx operator+(const compx& t)const{ return compx(r+t.r,i+t.i); }
        compx operator-(const compx& t)const{ return compx(r-t.r,i-t.i); }
        compx operator*(const compx& t)const{ return compx(r*t.r-i*t.i,r*t.i+i*t.r); }
    }a[NN],b[NN],c[NN],d[NN],w[NN];
    int l,len,r[NN];
    void prework(int n){
        for(len=1,l=0;len<n;len<<=1) ++l;
        for(int i=0;i<len;i++){
            r[i]=(r[i>>1]>>1)|((i&1)<<l-1);
            w[i]=compx(cos(2.0*i*PI/len),sin(2.0*i*PI/len));
        }
    }
    void fft(compx* a,int n){
        for(int i=0;i<n;i++) if(r[i]>i) swap(a[i],a[r[i]]);
        for(int t=n>>1,d=1;d<n;d<<=1,t>>=1)
            for(int i=0;i<n;i+=(d<<1)) for(int j=0;j<d;j++){
                compx tmp=w[j*t]*a[i+j+d];
                a[i+j+d]=a[i+j]-tmp; a[i+j]=a[i+j]+tmp;
            }
    }
    void mtt(int* f,int* g,int* res,int n){
        prework(n);
        for(int i=0;i<n;i++){
            a[i]=compx(f[i]>>15,f[i]&32767);
            c[i]=compx(g[i]>>15,g[i]&32767);
        }
        fft(a,len); fft(c,len);
        for(int i=1;i<len;i++) b[i]=a[len-i].conj(), d[i]=c[len-i].conj();
        b[0]=a[0].conj(); d[0]=c[0].conj();
        for(int i=0;i<len;i++){
            compx aa=(a[i]+b[i])*compx(0.5,0), bb=(a[i]-b[i])*compx(0,-0.5);
            compx cc=(c[i]+d[i])*compx(0.5,0), dd=(c[i]-d[i])*compx(0,-0.5);
            a[i]=aa*cc+compx(0,1)*(aa*dd+cc*bb); b[i]=bb*dd; w[i].i=-w[i].i;
        }
        fft(a,len); fft(b,len);
        for(int i=0;i<len;i++){
            LL x=(LL)(a[i].r/len+0.5)%mod,y=(LL)(a[i].i/len+0.5)%mod,z=(LL)(b[i].r/len+0.5)%mod;
            res[i]=((1ll<<30)*x%mod+(1ll<<15)*y%mod+z+mod)%mod;
        }
    }
} using namespace Poly_Calculation;

signed main(){
    n=read()+1; m=read()+1; mod=read();
    for(int i=0;i<n;i++) f[i]=read();
    for(int i=0;i<m;i++) g[i]=read();
    mtt(f,g,res,n+m);
    for(int i=0;i<n+m-1;i++) write(res[i],' ');
    return puts(""),0;
}

分治 \(\tt{FFT}\)

\(f_i=\sum f_jg_{i-j}\) 这样像卷积但又有些ex的式子。做CDQ分治就好。

洛谷P4721[模板]分治FFT
#include<bits/stdc++.h>
using namespace std;

namespace IO{
    #define int LL
    typedef long long LL;
    int read(){
        int x=0,f=0; char ch=getchar();
        while(ch<'0'||ch>'9'){ f|=ch=='-'; ch=getchar(); }
        while(ch>='0'&&ch<='9'){ x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
        return f?-x:x;
    } char outp[50];
    void write(int x,char sp,int len=0){
        if(x<0) putchar('-'), x=-x;
        do{ outp[len++]=x%10+'0'; x/=10; }while(x);
        for(int i=len-1;~i;i--) putchar(outp[i]); putchar(sp);
    }
    void ckmin(int& x,int y){ x=x<y?x:y; }
    void ckmax(int& x,int y){ x=x>y?x:y; }
} using namespace IO;

const int NN=600010,mod=998244353;
int n;
int qpow(int a,int b=mod-2,int res=1){
    for(;b;b>>=1,a=a*a%mod)
        if(b&1) res=res*a%mod;
    return res;
}

namespace Poly_Calculation{
    const int rt=3;
    int l,len,r[NN],a[NN],b[NN],g[NN],tp[NN],pt[NN];
    void prework(int n){
        for(len=1,l=0;len<=n;len<<=1) ++l;
        for(int i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<l-1);
        g[0]=1; g[1]=qpow(rt,(mod-1)/len);
        for(int i=2;i<len;i++) g[i]=g[i-1]*g[1]%mod;
    }
    void ntt(int* a,int n){
        for(int i=0;i<n;i++) if(r[i]>i) swap(a[i],a[r[i]]);
        for(int t=n>>1,d=1;d<n;d<<=1,t>>=1)
            for(int i=0;i<n;i+=(d<<1)) for(int j=0;j<d;j++){
                int tmp=g[j*t]*a[i+j+d]%mod;
                a[i+j+d]=(a[i+j]-tmp+mod)%mod; a[i+j]=(a[i+j]+tmp)%mod;
            }
    }
    void conquer(int* ff,int* gg,int l,int r){
        if(l==r){ if(!l) ff[0]=1; return; }
        int mid=l+r>>1; conquer(ff,gg,l,mid); prework(r-l+1);
        for(int i=l;i<=mid;i++) tp[i-l]=ff[i];
        for(int i=0;i<=r-l;i++) pt[i]=gg[i];
        ntt(tp,len); ntt(pt,len);
        for(int i=0;i<len;i++) tp[i]=tp[i]*pt[i]%mod, g[i]=qpow(g[i]);
        ntt(tp,len); int inv=qpow(len);
        for(int i=mid+1;i<=r;i++) ff[i]=(ff[i]+tp[i-l]*inv)%mod;
        for(int i=0;i<len;i++) tp[i]=pt[i]=0;
        conquer(ff,gg,mid+1,r);
    }
} using namespace Poly_Calculation;

signed main(){
    n=read();
    for(int i=1;i<n;i++) b[i]=read();
    conquer(a,b,0,n-1);
    for(int i=0;i<n;i++) write(a[i],' ');
    return puts(""),0;
}

多项式除法/多项式取模

已知 \(n\) 次多项式 \(F\)\(m\) 次多项式 \(G\) ,求 \(n-m\) 次多项式 \(Q\) 与次数不大于 \(m\) 的多项式 \(R\) ,满足 \(F(x)=G(x)Q(x)+R(x)\)

把多项式系数翻转,然后可以提出 \(x\) 的次幂,令 \(F\) 反转后为 \(F_r\) ,有

\[F_r(x)=Q_r(x)G_r(x)+x^{n-m+1}R_r(x)
\]

在模 \(x^{n-m+1}\) 意义下后面的东西为 \(0\) ,于是模 \(x^{n-m+1}\) 求逆后容易得出两个多项式。

洛谷P4512[模板]多项式除法
#include<bits/stdc++.h>
using namespace std;

namespace IO{
    #define int LL
    typedef long long LL;
    int read(){
        int x=0,f=0; char ch=getchar();
        while(ch<'0'||ch>'9'){ f|=ch=='-'; ch=getchar(); }
        while(ch>='0'&&ch<='9'){ x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
        return f?-x:x;
    } char outp[50];
    void write(int x,char sp,int len=0){
        if(x<0) putchar('-'), x=-x;
        do{ outp[len++]=x%10+'0'; x/=10; }while(x);
        for(int i=len-1;~i;i--) putchar(outp[i]); putchar(sp);
    }
    void ckmin(int& x,int y){ x=x<y?x:y; }
    void ckmax(int& x,int y){ x=x>y?x:y; }
} using namespace IO;

const int NN=400010,mod=998244353;
int n,m,f[NN],q[NN],res[NN];
int qpow(int a,int b=mod-2,int res=1){
    for(;b;b>>=1,a=a*a%mod)
        if(b&1) res=res*a%mod;
    return res;
}

namespace Poly_Calculation{
    const int rt=3,two=mod+2;
    int l,len,r[NN],g[NN],ff[NN],qq[NN],hh[NN],dd[NN],vv[NN],tp[NN],pt[NN];
    void revg(){ g[1]=qpow(g[1]); for(int i=2;i<len;i++) g[i]=g[i-1]*g[1]%mod; }
    void prework(int n){
        for(len=1,l=0;len<=n;len<<=1) ++l;
        for(int i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<l-1);
        g[0]=1; g[1]=qpow(rt,(mod-1)/len);
        for(int i=2;i<len;i++) g[i]=g[i-1]*g[1]%mod;
    }
    void ntt(int* a,int n){
        for(int i=0;i<n;i++) if(r[i]>i) swap(a[i],a[r[i]]);
        for(int t=n>>1,d=1;d<n;d<<=1,t>>=1)
            for(int i=0;i<n;i+=(d<<1)) for(int j=0;j<d;j++){
                int tmp=g[j*t]*a[i+j+d]%mod;
                a[i+j+d]=(a[i+j]-tmp+mod)%mod; a[i+j]=(a[i+j]+tmp)%mod;
            }
    }
    void mul(int* a,int* b,int* c,int n){
        for(int i=0;i<n;i++) tp[i]=a[i], pt[i]=b[i];
        prework(n); ntt(tp,len); ntt(pt,len);
        for(int i=0;i<len;i++) tp[i]=tp[i]*pt[i]%mod;
        revg(); ntt(tp,len); int inv=qpow(len);
        for(int i=0;i<len;i++) c[i]=tp[i]*inv%mod, tp[i]=pt[i]=0;
    }
    void polyinv(int* a,int* res,int deg){
        if(deg==1) return void(res[0]=qpow(a[0]));
        polyinv(a,res,deg+1>>1); prework(deg<<1);
        for(int i=0;i<deg;i++) tp[i]=a[i];
        ntt(tp,len); ntt(res,len);
        for(int i=0;i<len;i++) res[i]=(two-res[i]*tp[i]%mod)*res[i]%mod;
        revg(); ntt(res,len); int inv=qpow(len);
        for(int i=0;i<deg;i++) res[i]=res[i]*inv%mod, tp[i]=0;
        for(int i=deg;i<len;i++) res[i]=tp[i]=0;
    }
    void polymod(int* a,int* b,int* q,int* r,int n,int m){
        for(int i=0;i<=n;i++) q[n-i]=ff[i]=a[i];
        for(int i=0;i<=m;i++) dd[m-i]=hh[i]=b[i];
        for(int i=n-m+2;i<=m;i++) dd[i]=0;
        polyinv(dd,vv,n-m+1); mul(q,vv,q,n<<1); reverse(q,q+n-m+1);
        for(int i=n-m+1;i<=n;i++) q[i]=0; mul(hh,q,hh,n<<1);
        for(int i=0;i<m;i++) r[i]=(ff[i]+mod-hh[i])%mod;
        for(int i=0;i<len;i++) dd[i]=ff[i]=hh[i]=0;
    }
} using namespace Poly_Calculation;

signed main(){
    n=read(); m=read();
    for(int i=0;i<=n;i++) f[i]=read();
    for(int i=0;i<=m;i++) q[i]=read();
    polymod(f,q,qq,res,n,m);
    for(int i=0;i<=n-m;i++)write(qq[i],' ');puts("");
    for(int i=0;i<m;i++) write(res[i],' '); puts("");
    return 0;
}

原文链接: https://www.cnblogs.com/keeeeeeen/p/15834885.html

欢迎关注

微信关注下方公众号,第一时间获取干货硬货;公众号内回复【pdf】免费获取数百本计算机经典书籍

    多项式八分之三家桶

原创文章受到原创版权保护。转载请注明出处:https://www.ccppcoding.com/archives/187480

非原创文章文中已经注明原地址,如有侵权,联系删除

关注公众号【高性能架构探索】,第一时间获取最新文章

转载文章受原作者版权保护。转载请注明原作者出处!

(0)
上一篇 2023年2月12日 上午11:18
下一篇 2023年2月12日 上午11:18

相关推荐