多项式的各种操作

码了个100行结构体多项式代码，把除了复合和复合逆以外的多项式操作都写了。顺便把多项式各种操作的原理写一下。

多项式乘法

这个。。NTT

class Array{
    private:
        vector<int>a;
    public:
        Array(const int size,const int f):a(size,f){}
        void push(int n){a.push_back(n);}
        Array(int* l=NULL,int* r=NULL){while(l!=r)push(*l),++l;}
        inline int size(){return a.size();}
        inline int& operator [] (const int x){return a[x];}
        void resize(int n){a.resize(n);}
        void clear(){a.clear();}
        void swap(){reverse(a.begin(),a.end());}
};
int* NTT(const int len,Array& a,const bool Ty,int* r=NULL){
    if(!r){
        r=new int[len];
        r[0]=0;int L=log2(len);
        for(int i=0;i<len;i++)
            r[i]=(r[i>>1]>>1)|((i&1)<<L-1);
    }
    for(int i=0;i<len;i++)
        if(i<r[i])swap(a[i],a[r[i]]);
    for(int i=1;i<len;i<<=1){
        int T=poww(Ty?g:332748118,(mod-1)/(i<<1),mod);
        for(int W=i<<1,j=0;j<len;j+=W){
            ll omega=1;
            for(int k=0;k<i;++k,omega=omega*T%mod){
                ll x(a[j+k]),y(omega*a[i+j+k]%mod);
                a[j+k]=x+y;(a[j+k]>mod)&&(a[j+k]-=mod);
                a[i+j+k]=x-y+mod;(a[i+j+k]>mod)&&(a[i+j+k]-=mod);
            }
        }
    }
    return r;
}
Array Mul(Array x,Array y){
    int n=x.size()-1,m=y.size()-1;
    int limit=1;
    while(limit<=n+m)limit<<=1;
    Array ans;
    x.resize(limit+1);
    y.resize(limit+1);
    int *r;
    r=NTT(limit,x,1);
    NTT(limit,y,1,r);
    for(int i=0;i<=limit;i++) x[i]=1ll*x[i]*y[i]%mod;
    NTT(limit,x,0,r);
    int tem=poww(limit,mod-2,mod);
    for(int i=0;i<=n+m;i++) x[i]=(ll)x[i]*tem%mod;
    x.resize(n+m+1);
    return x;
}

牛顿迭代

就是一个在实数域和复数域上近似求解方程的方法。

公式为：
$$已知方程f(x)=0，x_{n+1}=x_n-\frac{f(x_n)}{f^{\prime }(x_n)}$$
这个$x_n$不断迭代就会以平方收敛，靠近方程的解，也就是说牛顿迭代每次迭代都会使有效数字翻一倍

那么，把这个方法拓展到多项式上，就可以做到乘法逆元、exp的运算，迭代logN次就能得到正确答案。

乘法逆元

设关于$f(x)$的方程$f(x)g(x)=1$，即$f(x)g(x)-1=0$

把这个式子套牛顿迭代，把$f_n(x)$替代$x_n$：

$f_{n+1}(x)=f_n(x)-\frac{f_n(x)g(x)-1}{g(x)}$

化简得：

$f_{n+1}(x)=2f_n(x)-g(x)f_n^2(x)$

然后递归就好了，递归logN次，同时注意递归边界，即$f(x)$ 为0次多项式时的解，就是 $g(x)$ 0次项系数的模意义下的乘法逆元。

设这个复杂度为$T(N)$，那么$T(N)=T(N/2)+O(NlogN)$，可以证明（我不会证），$T(N)$其实就是$O(NlogN)$的，但是常数更大。

void Rev(Array &x,Array y){
    int n=x.size()-1,m=y.size()-1;
    int limit=1;
    while(limit<=n+m)limit<<=1;
    Array ans;
    x.resize(limit+1);
    y.resize(limit+1);
    int *r;
    r=NTT(limit,x,1);
    NTT(limit,y,1,r);
    for(int i=0;i<=limit;i++) x[i]=(ll)(2ll-(ll)x[i]*y[i]%mod+mod)%mod*y[i]%mod;
    NTT(limit,x,0,r);
    int tem=poww(limit,mod-2,mod);
    for(int i=0;i<=n+m;i++) x[i]=(ll)x[i]*tem%mod;
    x.resize(n+m+1);
}
Array Inv(Array a){
    int N=a.size();
    if(N==1)return Array(1,poww(a[0],mod-2,mod));
    Array b=a;b.resize(N+1>>1);
    b=Inv(b);b.resize(N);
    Rev(a,b);
    a.resize(N);
    return a;
}

ln

其实这个。。就是一个式子：
$$ln(f(x))=\int \frac{f^{\prime }(x)}{f(x)}dx$$
但是它需要逆元运算，所以复杂度是O(NlogN)

Array Drv(Array a){ //求导
	for(int i=0;i<a.size()-1;i++) a[i]=1ll*a[i+1]*(i+1)%mod;
	a[a.size()-1]=0;
	return a;
}
Array Idf(Array a){ //不定积分
	for(int i=a.size()-1;i>0;i--) a[i]=1ll*a[i-1]*poww(i,mod-2,mod)%mod;
	a[0]=0;
	return a;
}
Array ln(Array a){
	return Idf(Mul(Drv(a),Inv(a)));
}

exp

同样是牛顿迭代法，推式子即可：
$f_{n+1}(x)=f_n(x)-\frac{\ln (f_n(x))-g(x)}{\frac{1}{f_n(x)}}$
$f_{n+1}(x)=f_n(x)(1-\ln (f_n(x)+g(x))$

然后做法就和之前一样了，这个复杂度也是O(NlogN)，但是是套上逆元的O(NlogN)，也就是常数是逆元的常数的平方，常数很大。

Array Exp(Array a){
    int N=a.size();
    if(N==1) return Array(1,1);
    Array b=a;b.resize(N+1>>1);
    b=Exp(b);b.resize(N);
    Array b1=ln(b);
    for(int i=0;i<a.size();i++) b1[i]=1ll*(a[i]-b1[i]+mod)%mod;
    b1[0]=1ll*(b1[0]+1)%mod;
    a=Mul(b,b1);
    return a;
}

平方根

我一开始想法是求ln，然后除2再exp回去，但是T了（略略略），所以还是用牛顿迭代。

需要注意的是，这里的边界不太好处理，要处理模意义的平方根，所以用了bsgs算法。

int bsgs(int y,int z,int p){
	int m=sqrt(p)+1;
	map<int,int> mp;
	for(int i=0,t=z;i<m;i++,t=1ll*y*t%p) mp[t]=i;
	for(int i=1,tt=poww(y,m,p),t=tt;i<=m+1;i++,t=1ll*t*tt%p){
		int k=mp[t];if(k==0) continue;
		return i*m-k;
	}
}
int sqrtt(int x){
	return min(poww(3,bsgs(3,x,mod)/2,mod),mod-poww(3,bsgs(3,x,mod)/2,mod));
}
Array Sqrt(Array a){
	int N=a.size();
	if(N==1) return Array(1,sqrtt(a[0]));
	Array b=a;b.resize(N+1>>1);
	b=Sqrt(b);b.resize(N);
	a=Mul(a,Inv(b));
	a.resize(N);
	int inv2=poww(2,mod-2,mod);
	for(int i=0;i<a.size();i++) a[i]=1ll*(a[i]+b[i])*inv2%mod;
	return a;
}