学了三天，看了各种资料，总算差不多学会了、、、数学太差学这种东西太难了、、感觉学这东西还是背板 模板基本都是贴黈力哒，大家不用看啦

多项式乘法就是把多项式用点值表示，能在O(n)时间内完成乘法运算

然后通过快速傅里叶变换O(nlgn)把多项式转换为点值表示，再通过逆离散傅里叶变换O(nlgn)把点值转换为多项式，就完成了多项式乘法

单位根是满足x^n=1的复数wn,它可以看做是复平面上x轴正方向绕逆时针方向旋转2π/n的复数

然后w(2m,2n)=w(m,n),w(m,2n)=-(m+n,2n)，可以利用这两个性质加快运算

然后假设我们已经得到A的偶数项之和A0和奇数项之和A1的点值表示

然后A(w(m,n))=A0(w(m,n/2))+w(m,n)*A1(w(m,n/2)) A(w(m+n/2,n))=A0((w,n/2))-w(m,n)*A1(w(m,n/2))

然后发现只要O(n)的时间就可以得到A的点值表示，总共有lgn层，快速傅里叶变换的时间复杂度就是O(nlgn)

然后可以通过神奇的蝴蝶变换写成非递归的FFT，不仅代码短，还非常快

然后逆离散傅里叶变换就是把w(m,n)翻转，然后对得出的答案全除以N就行啦。。反正就多一行代码，我也不知道什么原理、

FFT模板

#include<cstdio>
#include<cmath>
#include<iostream>
#define t 17
#define N 262144
using namespace std;
struct Complex{double x,y;}a1[N],a2[N],b[N],c1[N],c2[N],w,w0,x;
int n,m,i,g[N];double pie=acos(-1);
Complex operator +(Complex x,Complex y){Complex z;z.x=x.x+y.x;z.y=x.y+y.y;return z;}//复数加法 
Complex operator -(Complex x,Complex y){Complex z;z.x=x.x-y.x;z.y=x.y-y.y;return z;}//复数减法 
Complex operator *(Complex x,Complex y){Complex z;z.x=x.x*y.x-x.y*y.y;z.y=x.x*y.y+x.y*y.x;return z;}//复数乘法 
void FFT(Complex a[N],int f){//f=1为傅里叶变换，f=-1为逆傅里叶变换
	int i,j,k;
	for(i=0;i<N;i++)if(g[i]>i)swap(a[i],a[g[i]]);//蝴蝶变化 
	for(i=1;i<N;i<<=1){
		w.x=cos(pie/i);w.y=sin(pie/i)*f;//w(1,n)
		for(j=0;j<N;j+=i<<1){
			w0.x=1;w0.y=0;
			for(k=0;k<i;k++){
				x=a[j+k];
				a[j+k]=x+(w0*a[j+k+i]);//a[k]=a[k](0)+w(k,n)a[k](1) 
				a[j+k+i]=x-(w0*a[j+k+i]);//a[k+n/2]=a[k](0)-w(k,n)a[k](1)+a[k](0)+w(k+n/2,n)a[k](1)
				w0=w0*w;//w(now+1,n)
			}
		}
	}
	if(f<0)for(i=0;i<N;i++)a[i].x/=N;
}
int main(){
	scanf("%d",&n);
	for(i=0;i<N;i++)g[i]=g[i>>1]>>1|((i&1)<<t);//蝴蝶变化预处理 
	for(i=0;i<n;i++)scanf("%lf",&a1[i].x);
	for(i=0;i<n;i++)a2[i].x=a1[n-i-1].x;
	for(i=1;i<n;i++)b[i].x=1.0/i/i;
	FFT(a1,1);FFT(a2,1);FFT(b,1);
	for(i=0;i<N;i++)c1[i]=a1[i]*b[i];
	for(i=0;i<N;i++)c2[i]=a2[i]*b[i];//点值乘法 
	FFT(c1,-1);FFT(c2,-1);
	for(i=0;i<n;i++)printf("%.6lf\n",c1[i].x-c2[n-i-1].x);
}

BZOJ2179 高精度乘法

#include<cstdio>
#include<cmath>
#include<algorithm>
#define t 16
#define N 131072
#define pie acos(-1)
using namespace std;
struct Cp{double x,y;}A[N],B[N],C[N],x,w,w0;
int n,i,g[N],ans[N];char s[N];
Cp operator +(Cp x,Cp y){Cp z;z.x=x.x+y.x;z.y=x.y+y.y;return z;}
Cp operator -(Cp x,Cp y){Cp z;z.x=x.x-y.x;z.y=x.y-y.y;return z;}
Cp operator *(Cp x,Cp y){Cp z;z.x=x.x*y.x-x.y*y.y;z.y=x.x*y.y+x.y*y.x;return z;}
void FFT(Cp A[N],int f){
	int i,j,k;for(i=0;i<N;i++)if(g[i]>i)swap(A[i],A[g[i]]);
	for(i=1;i<N;i<<=1){
		w.x=cos(pie/i);w.y=sin(pie/i)*f;
		for(j=0;j<N;j+=i<<1){
			w0.x=1;w0.y=0;
			for(k=0;k<i;k++){
				x=A[j+k];
				A[j+k]=x+(w0*A[j+k+i]);
				A[j+k+i]=x-(w0*A[j+k+i]);
				w0=w0*w;
			}
		}
	}
	if(f<0)for(i=0;i<N;i++)A[i].x/=N;
}
int main(){
	scanf("%d",&n);
	scanf("%s",s+1);for(i=1;i<=n;i++)A[n-i].x=(double)s[i]-'0';
	scanf("%s",s+1);for(i=1;i<=n;i++)B[n-i].x=(double)s[i]-'0';
	for(i=1;i<N;i++)g[i]=g[i>>1]>>1|((i&1)<<t);
	FFT(A,1);FFT(B,1);
	for(i=0;i<N;i++)C[i]=A[i]*B[i];
	FFT(C,-1);
	for(i=0;i<=2*n-2;i++){
		ans[i]+=int(C[i].x+0.5);
		if(ans[i]>=10)ans[i+1]+=ans[i]/10,ans[i]%=10;
	}
	if(ans[2*n-1])printf("%d",ans[2*n-1]);
	for(i=2*n-2;i>=0;i--)printf("%d",ans[i]);
}

BZOJ2194

#include<cstdio>
#include<cmath>
#include<algorithm>
#define t 17
#define N 262144
#define pie acos(-1)
using namespace std;
struct Cp{double x,y;}A[N],B[N],C[N],x,w,w0;
int n,i,g[N];
Cp operator +(Cp x,Cp y){Cp z;z.x=x.x+y.x;z.y=x.y+y.y;return z;}
Cp operator -(Cp x,Cp y){Cp z;z.x=x.x-y.x;z.y=x.y-y.y;return z;}
Cp operator *(Cp x,Cp y){Cp z;z.x=x.x*y.x-x.y*y.y;z.y=x.x*y.y+x.y*y.x;return z;}
void FFT(Cp A[N],int f){
	int i,j,k;for(i=0;i<N;i++)if(g[i]>i)swap(A[i],A[g[i]]);
	for(i=1;i<N;i<<=1){
		w.x=cos(pie/i);w.y=sin(pie/i)*f;
		for(j=0;j<N;j+=i<<1){
			w0.x=1;w0.y=0;
			for(k=0;k<i;k++){
				x=A[j+k];
				A[j+k]=x+(w0*A[j+k+i]);
				A[j+k+i]=x-(w0*A[j+k+i]);
				w0=w0*w;
			}
		}
	}
	if(f<0)for(i=0;i<N;i++)A[i].x/=N;
}
int main(){
	scanf("%d",&n);
	for(i=0;i<n;i++)scanf("%lf%lf",&A[i].x,&B[n-i-1].x);
	for(i=1;i<N;i++)g[i]=g[i>>1]>>1|((i&1)<<t);
	FFT(A,1);FFT(B,1);
	for(i=0;i<N;i++)C[i]=A[i]*B[i];
	FFT(C,-1);
	for(i=n-1;i<=2*n-2;i++)printf("%d\n",int(C[i].x+0.5));
}

BZOJ3160 先不管连续的，枚举每个中间点，答案要加上2^f[i]-1，f[i]表示以i为中心选择相同字符的方案数，显然要使得a[j]==a[k]且j+k=i，那么先把所有a的值赋为1，弄一遍FFT，再把所有值赋为2，弄一遍FFT，就可以算出方案了，去掉连续的直接Manacher即可

#include<bits/stdc++.h>
#define N 262144
#define M 1000000007
#define pi acos(-1)
using namespace std;typedef complex<double>C;
int n,m,i,t,l,ans,mx,p,r[N],g[N],pw[N];C A[N],B[N];char s[N],w[N];
void FFT(C *a,int f){
	int i,j,k,p;for(i=0;i<N;i++)if(g[i]>i)swap(a[i],a[g[i]]);
	for(i=1;i<N;i<<=1){
		C e(cos(pi/i),sin(pi/i)*f);
		for(j=0;j<N;j+=i*2){
			C w(1,0);for(k=0;k<i;k++,w*=e){
				C x=a[j+k],y=w*a[j+k+i];
				a[j+k]=x+y;a[j+k+i]=x-y;
			}
		}
	}
	if(f<0)for(i=0;i<N;i++)a[i]/=N;
}
int main(){
	for(scanf("%s",s),n=strlen(s),i=0;i<n;i++)if(s[i]=='a')A[i]=1;else B[i]=1;
	for(pw[0]=i=1;i<N;i++)g[i]=g[i>>1]>>1|((i&1)<<17),pw[i]=pw[i-1]*2%M;
	for(FFT(A,1),FFT(B,1),i=0;i<N;i++)A[i]*=A[i],B[i]*=B[i];
	for(FFT(A,-1),FFT(B,-1),i=0;i<N;i++)ans=(ans+pw[(int)(A[i].real()+B[i].real()+1.5)/2]-1)%M;
	for(w[l++]='d',w[l++]='c',i=0;i<n;i++)w[l++]=s[i],w[l++]='c';
	for(i=0;i<l;i++){
		r[i]=mx>i?min(r[2*p-i],mx-i):1;
		for(;w[i+r[i]]==w[i-r[i]];r[i]++);ans=(ans-r[i]/2+M)%M;
		if(i+r[i]>mx)mx=i+r[i],p=i;
	}
	printf("%d",ans);
}

快速数论变换——NTT

然后FFT的精度不够高，如果遇到毛主力这样的出题人多半是没希望的，于是NTT出现啦

NTT要求模数M必须满足M=w·2^k+1，如果是3^k的话据黈力说要分治NTT，那多半就GG啦

然后弄出个原根，一般题目都是3吧。。设原根为g,然后单位根就是(g^((M-1)/m))mod M，但程序里就每次b++;w=powmod(s,(1<<23-b)*119);这么搞就行啦，我也不知道为什么，反正蛮神哒。。

#include<cstdio>
#include<algorithm>
#define N 262144
#define t 17
#define s 3//原根
#define mod 998244353 
using namespace std;
typedef long long ll;
int n,m,i,g[N];
ll A[N],B[N],C[N],x,w,w0;
ll powmod(ll x,int y){//快速幂 
    ll ans=1;
    for(;y;y>>=1){
        if(y&1)ans=ans*x%mod;
        x=x*x%mod;
    }
    return ans;
}
void NTT(ll A[N],int f){//快速数论变换 
    int i,j,k,b=0;
    for(i=0;i<N;i++)if(g[i]>i)swap(A[i],A[g[i]]);//蝴蝶变换 
    for(i=1;i<N;i<<=1){
        b++;
        w=powmod(s,(1<<23-b)*119);//998244352=119*2^23
        if(f<0)w=powmod(w,mod-2);
        for(j=0;j<N;j+=i<<1){
            w0=1;
            for(k=0;k<i;k++){
                x=A[j+k];
                A[j+k]=(x+A[j+k+i]*w0)%mod;
                A[j+k+i]=(x-A[j+k+i]*w0)%mod;
                if(A[j+k+i]<0)A[j+k+i]+=mod;
                w0=(w0*w)%mod;
            }
        }
    }
    if(f<0){//逆变换
        ll v=powmod(N,mod-2);
        for(i=0;i<N;i++)A[i]=A[i]*v%mod;
    }
}
int main(){
    scanf("%d%d",&n,&m);
    for(i=0;i<=n;i++)scanf("%lld",&A[i]);
    for(i=0;i<=m;i++)scanf("%lld",&B[i]);
    for(i=1;i<N;i++)g[i]=g[i>>1]>>1|((i&1)<<t);//预处理蝴蝶变换
    NTT(A,1);NTT(B,1);
    for(i=0;i<N;i++)C[i]=A[i]*B[i]%mod;//点值乘法
    NTT(C,-1);
    for(i=0;i<=n+m;i++)printf("%lld ",C[i]);
}

BZOJ3992：首先很容易得到dp方程f[i][j]表示前i个数乘积为j的方案数

i这一层显然可以快速面摞掉，可是这样还是m^2lgn的

由于M是质数，那么把1~M-1中的每个数用原根表示出来，这样转移的时候就变成了一个多项式乘法

最后的时间复杂度是O(mlgnlgm)

#include<cstdio>
#include<algorithm>
#define N 16384
#define M 1004535809
using namespace std;typedef long long LL;
int k,n,i,p,S,x,w,g,h[N],to[N];LL iv,sum,A[N],B[N];
LL pow(LL a,LL b,LL p){for(sum=1,a%=p;b;a=a*a%p,b>>=1)if(b&1)sum=sum*a%p;return sum;}
bool ok(int x){for(int i=2;i*i<=n;i++)if((n-1)%i==0&&pow(x,(n-1)/i,n)==1)return 0;return 1;}
int get(){for(int i=2;i<=n;i++)if(ok(i))return i;}
void NTT(LL A[N],int f){
    int i,j,k,b=0;LL x,y,w,w0;
    for(i=0;i<N;i++)if(h[i]>i)swap(A[i],A[h[i]]);
    for(i=1;i<N;i<<=1){
        b++;w=pow(3,(1<<21-b)*479,M);if(f<0)w=pow(w,M-2,M);
        for(j=0;j<N;j+=i<<1)
            for(w0=1,k=0;k<i;w0=(w0*w)%M,k++)x=A[j+k],y=(A[j+k+i]*w0)%M,A[j+k]=(x+y)%M,A[j+k+i]=(x-y+M)%M;
    }
    if(f<0)for(i=0;i<N;i++)A[i]=A[i]*iv%M;
}
int main(){
	for(iv=pow(N,M-2,M),i=1;i<N;i++)h[i]=h[i>>1]>>1|((i&1)<<13);
	for(scanf("%d%d%d%d",&k,&n,&p,&S),g=get(),w=1,i=0;i<n-1;w=w*g%n,i++)to[w]=i;
	for(;S--;){scanf("%d",&x);if(x)A[to[x]]++;}
	for(B[0]=1;k;k>>=1){
		if(NTT(A,1),k&1){
			for(NTT(B,1),i=0;i<N;i++)B[i]=B[i]*A[i]%M;
			for(NTT(B,-1),i=N-1;i>=n-1;i--)B[i-n+1]=(B[i-n+1]+B[i])%M,B[i]=0;
		}
		for(i=0;i<N;i++)A[i]=A[i]*A[i]%M;
		for(NTT(A,-1),i=N-1;i>=n-1;i--)A[i-n+1]=(A[i-n+1]+A[i])%M,A[i]=0;
	}
	printf("%lld",B[to[p]]);
}

多项式求逆 遇到黈力这样的多项式大师出的题，就会碰到多项式求逆

假设我们要求多项式A(x)在模xⁿ意义下的逆B(x)，我们可以用倍增的思想，假设我们已经求出了模x^(t/2)意义下的答案C(x)，我们要求出模x^t意义下的答案B(x)，则A(x)·C(x)%x^(t/2)=1,A(x)·B(x)%x^(t/2)=1

则根据膜法的一些性质，(B(x)-C(x))%x^(t/2)=0

再根据膜法的一些性质把上式展开，得到(B(x)^2-2 B(x)·C(x)+C(x)^2)%x^(t/2)0

然后两边同乘A(x)在根据膜法的性质得到B(x)%(x^t)=C(x)·(2-A(x)·C(x))

然后每层做一次时间复杂度为tlgt,有lgn层，总时间复杂度为nlgn

寿司晚宴 n^2很好推，然后发现可以多项式求逆做

#include<bits/stdc++.h>
#define N 262144
#define M 998244353
using namespace std;
typedef long long LL;int n,i,T,g[N];LL A[N],B[N],C[N],x,w,w0,v;
LL pow(LL a,LL b,LL p){LL v=1;for(;b;a=a*a%p,b>>=1)if(b&1)v=v*a%p;return v;}
void Init(int n){
    int t,i;for(t=0;n>>t+2;t++);
    for(i=1;i<n;i++)g[i]=g[i/2]/2|((i&1)<<t);
}
void NTT(LL A[N],int f,int n){ 
    int i,j,k,b=0;Init(n);
    for(i=0;i<n;i++)if(g[i]>i)swap(A[i],A[g[i]]);
    for(i=1;i<n;i<<=1){
        b++;w=pow(3,(1<<23-b)*119,M);
        if(f<0)w=pow(w,M-2,M);
        for(j=0;j<n;j+=i*2)
            for(w0=1,k=0;k<i;w0=w0*w%M,k++)
                x=A[j+k],A[j+k]=(x+A[j+k+i]*w0)%M,A[j+k+i]=((x-A[j+k+i]*w0)%M+M)%M;
    }
    if(f<0)for(v=pow(n,M-2,M),i=0;i<n;i++)A[i]=A[i]*v%M;
}
void Inv(LL A[N],LL B[N],int n){
    int i,j;for(i=0;i<n;i++)C[i]=A[i];
    for(i=n;i<n*2;i++)C[i]=0;
    NTT(C,1,n*2);NTT(B,1,n*2);
    for(i=0;i<n*2;i++)C[i]=((2-C[i]*B[i])%M+M)%M;
    for(i=0;i<n*2;i++)B[i]=B[i]*C[i]%M;
    for(NTT(B,-1,n*2),i=n;i<n*2;i++)B[i]=0;
}
int main(){
    n=100000;for(A[0]=i=1;i<=n;i++)A[i]=A[i-1]*i%M;
    B[0]=pow(A[0],M-2,M);
    for(i=2;i/2<=n;i<<=1)Inv(A,B,i);
    A[0]--;NTT(A,1,N);NTT(B,1,N);
    for(i=0;i<N;i++)C[i]=A[i]*B[i]%M;NTT(C,-1,N);
    for(scanf("%d",&T);T--;printf("%lld\n",C[n]))scanf("%d",&n);
}

• lbnDP优化算法
• HNOI2016
• CQOI2016
• 模板
• 可持久化数据结构

• BZOJ 1588: [HNOI2002]营业额统计 && 红黑树板子
• [洛谷 P4389]付公主的背包
• [LibreOJ 2058][TJOI2016]求和
• [BZOJ 3512]DZY Loves Math IV
• [BZOJ 4916]神犇与蒟蒻