多项式Ⅱ

前言

因为之前的 Ⅰ 写的太屑了,所以又写了篇更屑的Ⅱ,其实是不同时间的个人垃圾理解了..

牛顿迭代

  • 要求解\(f(x)=0\)

    \[\frac{f(x_0)}{x_0-x}=f'(x_0)\\ x=x_0-\frac{f(x_0)}{f'(x_0)} \]

    不断往下迭代可以找到近似解。

多项式牛顿迭代

和普通的牛顿迭代其实一样

\[F=F_0-\frac{G(F_0)}{G'(F_0)} \]

如何理解近似解?

\(\bmod x^i\)开始向上\(\bmod x^{i+1}...\)

  • 多项式求逆

    \[A*B\equiv 1\pmod {x^n} \]

    构造

    \[G(B)=A*B-1 \]

    牛顿迭代找根

\[ \begin{aligned} B&=B_0-\frac{AB_0-1}{A}\\ &=B_0-B_0(AB_0-1)\\ &=2B_0-AB_0^2 \end{aligned} \]

向上倍增即可。

  • 多项式开根

    \[B^2\equiv A\pmod {x^n} \]

    构造

    \[G(B)=A-B^2 \]

    牛顿迭代找根

    \[\begin{aligned} B&=B_0-\frac{A-B_0^2}{-2B_0}\\ &=\frac{A}{2B_0}+\frac{B_0}{2} \end{aligned} \]

    向上倍增即可

  • 多项式exp

    \[B\equiv e^A \pmod {x^n} \]

    并没有弄清楚这个咋定义的,因此还不会直接用。

    \[\ln B=A \]

    构造

    \[G(B)=\ln B-A \]

    \[\begin{aligned} B&=B_0-\frac{\ln B_0-A}{\frac{1}{B_0}}\\ &=B_0(1-\ln B_0+A) \end{aligned} \]

多项式取ln

  • 多项式积分和求导

    \[\begin{aligned} F(x)&=\sum_{i=0}^na_ix^i\\ \frac{\mathrm dF}{\mathrm dx}&=\sum_{i=0}^{n-1}a_{i+1}(i+1)x^i\\ \int F(x)&=\sum_{i=1}^n\frac{a_{i-1}}{i}x^i \end{aligned} \]

  • 多项式取ln

    \[\begin{aligned} B&=\ln A\pmod {x^n}\\ B'&=\ln'(A)A'\pmod {x^n}\\ &=\frac{A'}{A}\pmod {x^n} \end{aligned} \]

  • 多项式快速幂

    \[F^k=e^{k\ln n}\pmod {x^n} \]

    常数无敌的\(n\log n\)

拉格朗日反演和快速插值什么的,我还没学...


放一个exp的代码

提几个细节

  1. 不知道为什么递归版效率优于非递归版,而递归版又好写一点,因此以后就写递归版了
  2. 每个需要倍增的,每次要外面开两个自己的数组,不许改传进来的值,注意这两个数组一定需要清空
  3. 要注意一下每次NTT的长度什么的
  4. 不要封装卷积乘法,因为慢,可以把多个多项式乘起来全部先转换成点表达,然后一起乘,再转回系数,这样需要注意长度。
#include <cstdio>
#include <algorithm>
const int N=(1<<18)+10;
const int mod=998244353,Gi=332748118;
#define mul(a,b) (1ll*(a)*(b)%mod)
#define add(a,b) ((a+b)%mod)
int invA[N],invB[N],lnA[N],lnB[N],exA[N],exB[N],ans[N],f[N],turn[N],n;
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;}
void NTT(int *a,int typ,int len)
{
    int L=-1;for(int i=1;i<len;i<<=1) ++L;
    for(int i=1;i<len;i++)
    {
        turn[i]=turn[i>>1]>>1|(i&1)<<L;
        if(i<turn[i]) std::swap(a[i],a[turn[i]]);
    }
    for(int le=1;le<len;le<<=1)
    {
        int wn=qp(typ?3: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 polyinv(int *a,int *b,int len)
{
    if(len==1) {b[0]=qp(a[0],mod-2);return;}
    polyinv(a,b,len>>1);
    for(int i=0;i<len<<1;i++) invA[i]=invB[i]=0;//
    for(int i=0;i<len;i++) invA[i]=a[i],invB[i]=b[i];
    NTT(invA,1,len<<1),NTT(invB,1,len<<1);
    for(int i=0;i<len<<1;i++) invA[i]=mul(invB[i],add(2,mod-mul(invA[i],invB[i])));
    NTT(invA,0,len<<1);
    for(int i=0;i<len;i++) b[i]=invA[i];
}
void direv(int *a,int n){for(int i=0;i<n;i++)a[i]=mul(a[i+1],i+1);a[n]=0;}
void inter(int *a,int n){for(int i=n;i;i--)a[i]=mul(a[i-1],qp(i,mod-2));a[0]=0;}
void polyln(int *a,int *b,int len)
{
    for(int i=0;i<len;i++) lnA[i]=a[i],lnB[i]=b[i];
    polyinv(lnA,lnB,len);
    direv(lnA,len-1);
    NTT(lnA,1,len<<1),NTT(lnB,1,len<<1);
    for(int i=0;i<len<<1;i++) lnA[i]=mul(lnA[i],lnB[i]);
    NTT(lnA,0,len<<1);
    inter(lnA,len-1);
    for(int i=0;i<len;i++) b[i]=lnA[i];
}
void polyexp(int *a,int *b,int len)
{
    if(len==1){b[0]=1;return;}//默认a[0]=0
    polyexp(a,b,len>>1);
    for(int i=0;i<len<<1;i++) exA[i]=exB[i]=0;//
    polyln(b,exA,len);
    for(int i=0;i<len;i++) exA[i]=add(a[i]+(i==0),mod-exA[i]);
    for(int i=0;i<len;i++) exB[i]=b[i];
    NTT(exA,1,len<<1),NTT(exB,1,len<<1);
    for(int i=0;i<len<<1;i++) exA[i]=mul(exA[i],exB[i]);
    NTT(exA,0,len<<1);
    for(int i=0;i<len;i++) b[i]=exA[i];
}
int main()
{
    scanf("%d",&n);
    for(int i=0;i<n;i++) scanf("%d",f+i);
    int len=1;while(len<n) len<<=1;
    polyexp(f,ans,len);
    for(int i=0;i<n;i++) printf("%d ",ans[i]);
    return 0;
}

2018.12.29

posted @ 2018-12-29 15:09  露迭月  阅读(207)  评论(0编辑  收藏  举报