多项式
多项式Ⅰ

update20181229

这篇屑文章咕咕了,我更改了一些代码写法(非递归->递归),从牛顿迭代的思想理解了一下,这里的一些思想可能比较繁琐,可以看看多项式Ⅱ

前言:本人只会瞎背板子,根本没怎么理解,而且只学了一两个最简单的(学一点更新一点

前置定义:多项式的度表示多项式的最高次数。

多项式求逆

  • 定义

    \[A(x)A^{-1}(x)\equiv 1\pmod {x^n} \]

    \(n\)\(A(x)\)的度,这里是两个多项式进行卷积之后对一个多项式(\(x^n\))取模

    • 为什么要取模?

      如果不取模,除了只有常数项的多项式以外,其他多项式的逆元都是无穷级数。

    • 一个取模的感性理解

      \(x^n\)取模相当于把次数大于\(n\)的项都扔掉

  • 倍增

    若求出了\(A(x)A'(x)\equiv 1\pmod {x^{\lceil\frac{n}{2}\rceil}}\)\(A'(x)\)

    \[A(x)A^{-1}(x)\equiv 1\pmod {x^n} \]

    \[A(x)A^{-1}(x)\equiv 1 \pmod {x^{\lceil\frac{n}{2}\rceil}} \]

    \[A(x)(A^{-1}(x)-A'(x))\equiv 0 \pmod {x^{\lceil\frac{n}{2}\rceil}} \]

    \[A^{-1}(x)-A'(x)\equiv 0 \pmod {x^{\lceil\frac{n}{2}\rceil}} \]

    \[(A^{-1}(x)-A'(x))^2\equiv 0 \pmod {x^n} \]

    \[(A^{-1}(x))^2-2A^{-1}(x)A'(x)+A'^2(x)\equiv 0 \pmod {x^n} \]

    \[A^{-1}(x)-2A'(x)+A(x)A'^2(x)\equiv 0 \pmod {x^n} \]

    \[A^{-1}(x)\equiv 2A'(x)-A(x)A'^2(x)\pmod {x^n} \]

    上面过程以已意会为主,谢谢。正经的觉得疑惑还是从定义去理解+讨论。

  • 细节

    注意每次乘法时取长度取多不取少,比如那个三个多项式乘开开四倍长度。

  • Code:

    #include <cstdio>
    #include <algorithm>
    #define ll long long
    const int N=(1<<18)+10;
    const int mod=998244353,G=3,Gi=332748118;
    #define add(x,y) ((x+y)%mod)
    #define Mul(x,y) (1ll*(x)*(y)%mod)
    int qp(int d,int k){int f=1;while(k){if(k&1)f=Mul(f,d);d=Mul(d,d);k>>=1;}return f;}
    int A[N],B[N],f[N],b[2][N],turn[N];
    void NTT(int *a,int len,int typ)
    {
        for(int i=0;i<len;i++)
            if(i<turn[i])
                std::swap(a[i],a[turn[i]]);
        for(int le=1;le<len;le<<=1)
        {
            int wn=qp(typ?G:Gi,(mod-1)/(le<<1));
            for(int p=0;p<len;p+=le<<1)
            {
                int w=1;
                for(int i=p;i<p+le;i++,w=Mul(w,wn))
                {
                    int tx=a[i],ty=Mul(w,a[i+le]);
                    a[i]=add(tx,ty);
                    a[i+le]=add(tx,mod-ty);
                }
            }
        }
        if(!typ)
        {
            int inv=qp(len,mod-2);
            for(int i=0;i<len;i++)
                a[i]=Mul(a[i],inv);
        }
    }
    void mul(int *a,int *b,int len)
    {
        int L=-1;for(int i=1;i<len;i<<=1,++L);
        for(int i=0;i<len;i++) turn[i]=turn[i>>1]>>1|(i&1)<<L,A[i]=B[i]=0;
        for(int i=0;i<len>>1;i++) A[i]=a[i],B[i]=b[i];
        NTT(A,len,1),NTT(B,len,1);
        for(int i=0;i<len;i++) A[i]=Mul(A[i],B[i]);
        NTT(A,len,0);
        for(int i=0;i<len;i++) a[i]=A[i];
    }
    void inv(int *a,int n)
    {
        int cur=0,len=2;
        b[cur][0]=qp(f[0],mod-2);
        while(len<=(n<<2))
        {
            cur^=1;
            for(int i=0;i<len>>1;i++) b[cur][i]=add(b[cur^1][i],b[cur^1][i]);
            mul(b[cur^1],b[cur^1],len);
            mul(b[cur^1],f,len);
            for(int i=0;i<len;i++) b[cur][i]=add(b[cur][i],mod-b[cur^1][i]);
            len<<=1;
        }
        for(int i=0;i<n;i++)
            a[i]=b[cur][i];
    }
    int main()
    {
        int n;scanf("%d",&n);
        for(int i=0;i<n;i++) scanf("%d",f+i);
        inv(f,n);
        for(int i=0;i<n;i++) printf("%d ",f[i]);
        return 0;
    }
    

多项式除法

\[A(x)=B(x)D(x)+R(x) \]

注:\(A(x)\)被除多项式,度为\(n\)\(B(x)\)除他的多项式,度为\(m\)\(D(x)\)商,度为\(n-m\),余数为\(R(x)\),度不大于\(m-1\)

  • 一种转换

    \(A^r(x)=x^nA(\frac{1}{x})\)

    感性理解:把它的系数颠倒了,可以\(O(n)\)处理

    • 有什么用?

      通过取模把余数扔掉。

  • \[A(x)=B(x)D(x)+R(x) \]

    \[x^nA(x)=x^nB(x)D(x)+x^nR(x) \]

    \[A^r(x)=B^r(x)D^r(x)+x^{n-m+1}R^r(x) \]

    \[A^r(x)\equiv B^r(x)D^r(x) \pmod {x^{n-m+1}} \]

    然后求个逆把\(D^r(x)\)弄出来,带回去就行了

  • Code:

    #include <cstdio>
    #include <algorithm>
    const int N=(1<<20)+10;
    const int mod=998244353,Ge=3,Gi=332748118;
    #define add(x,y) ((x+y)%mod)
    #define Mul(x,y) (1ll*(x)*(y)%mod)
    int qp(int d,int k){int f=1;while(k){if(k&1)f=Mul(f,d);d=Mul(d,d);k>>=1;}return f;}
    int F[N],G[N],b[2][N],X[N],Y[N],A[N],B[N],turn[N];
    void NTT(int *a,int len,int typ)
    {
        for(int i=0;i<len;i++)
            if(i<turn[i])
                std::swap(a[i],a[turn[i]]);
        for(int le=1;le<len;le<<=1)
        {
            int wn=qp(typ?Ge:Gi,(mod-1)/(le<<1));
            for(int p=0;p<len;p+=le<<1)
            {
                int w=1;
                for(int i=p;i<p+le;i++,w=Mul(w,wn))
                {
                    int tx=a[i],ty=Mul(w,a[i+le]);
                    a[i]=add(tx,ty);
                    a[i+le]=add(tx,mod-ty);
                }
            }
        }
        if(!typ)
        {
            int inv=qp(len,mod-2);
            for(int i=0;i<len;i++)
                a[i]=Mul(a[i],inv);
        }
    }
    void polymul(int *a,int *b,int len)
    {
        int L=-1;for(int i=1;i<len;i<<=1,++L);
        for(int i=0;i<len;i++) turn[i]=turn[i>>1]>>1|(i&1)<<L,A[i]=B[i]=0;
        for(int i=0;i<len>>1;i++) A[i]=a[i],B[i]=b[i];
        NTT(A,len,1),NTT(B,len,1);
        for(int i=0;i<len;i++) A[i]=Mul(A[i],B[i]);
        NTT(A,len,0);
        for(int i=0;i<len;i++) a[i]=A[i];
    }
    void polyinv(int *a,int n)
    {
        int cur=0,len=2;
        b[cur][0]=qp(a[0],mod-2);
        while(len<=(n<<2))
        {
            cur^=1;
            for(int i=0;i<len>>1;i++) b[cur][i]=add(b[cur^1][i],b[cur^1][i]);
            polymul(b[cur^1],b[cur^1],len);
            polymul(b[cur^1],a,len);
            for(int i=0;i<len;i++) b[cur][i]=add(b[cur][i],mod-b[cur^1][i]);
            len<<=1;
        }
        for(int i=0;i<n;i++)
            a[i]=b[cur][i];
    }
    void Reverse(int *a,int len)
    {
        for(int i=0;i<<1<len;i++) std::swap(a[i],a[len-i-1]);
    }
    void polydiv(int *a,int *b,int n,int m)
    {
        for(int i=0;i<n;i++) X[i]=a[i];
        for(int i=0;i<m;i++) Y[i]=b[i];
        Reverse(a,n),Reverse(b,m);
        polyinv(b,n-m+1);
        int len=1;while(len<=n<<2) len<<=1;
        polymul(a,b,len);
        Reverse(a,n-m+1);
        for(int i=n-m+1;i<len;i++) a[i]=0;
        polymul(Y,a,len);
        for(int i=0;i<m-1;i++) b[i]=add(X[i],mod-Y[i]);
    }
    int main()
    {
        int n,m;scanf("%d%d",&n,&m);++n,++m;
        for(int i=0;i<n;i++) scanf("%d",F+i);
        for(int i=0;i<m;i++) scanf("%d",G+i);
        polydiv(F,G,n,m);
        for(int i=0;i<n-m+1;i++) printf("%d ",F[i]);puts("");
        for(int i=0;i<m-1;i++) printf("%d ",G[i]);
        return 0;
    }
    

多项式开根

  • 像求逆那样进行倍增

    \[T^2\equiv F \pmod {2^n} \]

    \[T'^2\equiv F\pmod {2^{\lceil\frac{n}{2}\rceil}} \]

    \[(T+T')(T-T')\equiv 0 \pmod {2^{\lceil\frac{n}{2}\rceil}} \]

    \[(T-T')\equiv 0\pmod {2^{\lceil\frac{n}{2}\rceil}} \]

    \[T^2-2TT'+T'^2\equiv 0\pmod {2^n} \]

    \[F-2TT'+T'^2\equiv 0\pmod {2^n} \]

    \[T\equiv \frac{F+T'^2}{2T'}\pmod {2^n} \]

    或者再化简一步

    \[T\equiv \frac{F}{2T'}+\frac{T'}{2}\pmod {2^n} \]

  • Code:

    void polysqrt(int *a,int n)
    {
        int len=2,cur=0;
        sq[cur][0]=1;
        while(len<=(n<<2))
        {
            cur^=1;
            for(int i=0;i<len>>1;i++) sq[cur][i]=mul(sq[cur^1][i],i2);
            for(int i=0;i<len>>1;i++) sq[cur^1][i]=add(sq[cur^1][i],sq[cur^1][i]);
            polyinv(sq[cur^1],len);
            polymul(sq[cur^1],a,len);
            for(int i=0;i<len;i++) sq[cur][i]=add(sq[cur][i],sq[cur^1][i]);
            len<<=1;
        }
        for(int i=0;i<n;i++) a[i]=sq[cur][i];
    }
    

一定还会更新的...