需要背的板子
这是一个大坑。
先把写到的搬上来：
fread/fwrite 模板（注：需要在最后加上 flush，并且不支持负数）：
namespace IO{
const int SIZE=(1<<21)+1;
char ibuf[SIZE],*iS,*iT,obuf[SIZE],*oS=obuf,*oT=oS+SIZE-1;
#define gc() (iS==iT?(iT=(iS=ibuf)+fread(ibuf,1,SIZE,stdin),(iS==iT?EOF:*(iS++))):*(iS++))
void flush(){
fwrite(obuf,1,oS-obuf,stdout);
oS=obuf;
}
void putc(char x){
*oS++=x;
if(oS==oT){
flush();
}
}
template<typename Elem>
void read(Elem &a){
a=0;
char c=gc();
while(c<'0'||c>'9'){
c=gc();
}
while(c>='0'&&c<='9'){
a=(a<<1)+(a<<3)+(c^48);
c=gc();
}
}
template<typename Elem>
void write(Elem a){
if(a<10){
putc(a+'0');
return;
}
write(a/10);
putc(a%10+'0');
}
#undef gc
}
using IO::read;using IO::write;using IO::putc;
快速乘：
ll quick_mul(ll a,ll b,ll Mod){
ll z=(long double)a/Mod*b;
ll ans=(ull)a*b-(ull)z*Mod;
return (ans+Mod)%Mod;
}
FWT_or
void FWT_or(int *a,int flag,int len){
for(int i=1;i<len;i<<=1){
for(int j=0;j<len;j+=(i<<1)){
for(int k=0;k<i;k++){
if(flag==1){
a[i+j+k]=(a[i+j+k]+a[j+k])%Mod;
}
else{
a[i+j+k]=(a[i+j+k]-a[j+k]+Mod)%Mod;
}
}
}
}
}
FWT_and
void FWT_and(int *a,int flag,int len){
for(int i=1;i<len;i<<=1){
for(int j=0;j<len;j+=(i<<1)){
for(int k=0;k<i;k++){
if(flag==1){
a[j+k]=(a[i+j+k]+a[j+k])%Mod;
}
else{
a[j+k]=(a[j+k]-a[i+j+k]+Mod)%Mod;
}
}
}
}
}
FWT_xor
void FWT_xor(int *a,int flag,int len){
for(int i=1;i<len;i<<=1){
for(int j=0;j<len;j+=(i<<1)){
for(int k=0;k<i;k++){
int Nx=a[j+k],Ny=a[i+j+k];
a[j+k]=(Nx+Ny)%Mod;
a[i+j+k]=(Nx-Ny+Mod)%Mod;
if(flag==-1){
a[j+k]=1ll*a[j+k]*inv_2%Mod;
a[i+j+k]=1ll*a[i+j+k]*inv_2%Mod;
}
}
}
}
}
FFT
struct Complex{
double a,r;
friend Complex operator +(Complex a,Complex b){
Complex ans;
ans.a=a.a+b.a;
ans.r=a.r+b.r;
return ans;
}
friend Complex operator -(Complex a,Complex b){
Complex ans;
ans.a=a.a-b.a;
ans.r=a.r-b.r;
return ans;
}
friend Complex operator *(Complex a,Complex b){
Complex ans;
ans.a=a.a*b.a-a.r*b.r;
ans.r=a.a*b.r+a.r*b.a;
return ans;
}
friend Complex operator /(Complex a,int b){
a.a/=b;
a.r/=b;
return a;
}
Complex(double _a=0.0,double _r=0.0){
a=_a;
r=_r;
}
};
void FFT(Complex *a,int flag,int n){
static int R[Maxn+5],last_len=0;
int len=1,L=0;
while(len<n){
len<<=1;
L++;
}
if(len!=last_len){
last_len=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 j=1;j<len;j<<=1){
Complex T(cos(M_PI/j),flag*sin(M_PI/j));
for(int k=0;k<len;k+=(j<<1)){
Complex t(1,0);
for(int l=0;l<j;l++,t=t*T){
Complex Nx=a[k+l],Ny=t*a[j+k+l];
a[k+l]=Nx+Ny;
a[j+k+l]=Nx-Ny;
}
}
}
if(flag==-1){
for(int i=0;i<len;i++){
a[i]=a[i]/len;
}
}
}
MTT（假的）
namespace Poly{
typedef long double ld;
const ld PI=acos(-1.0);
struct Complex{
ld r,i;
Complex(ld _r=0,ld _i=0){
r=_r;
i=_i;
}
Complex conj(){
return Complex(r,-i);
}
friend Complex operator +(Complex a,Complex b){
return Complex(a.r+b.r,a.i+b.i);
}
friend Complex operator -(Complex a,Complex b){
return Complex(a.r-b.r,a.i-b.i);
}
friend Complex operator *(Complex a,Complex b){
return Complex(a.r*b.r-a.i*b.i,a.r*b.i+a.i*b.r);
}
Complex operator /=(int n){
r/=n;
i/=n;
return (*this);
}
};
typedef long long ll;
Complex a[Maxn+5],b[Maxn+5],c[Maxn+5],d[Maxn+5];
void FFT(Complex *a,int flag,int n){
static int R[Maxn+5],last_len=0;
int len=1,L=0;
while(len<n){
len<<=1;
L++;
}
if(len!=last_len){
last_len=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 j=1;j<len;j<<=1){
Complex T=Complex(cos(PI/j),flag*sin(PI/j));
for(int k=0;k<len;k+=(j<<1)){
Complex t=Complex(1,0);
for(int l=0;l<j;l++,t=t*T){
Complex Nx=a[k+l],Ny=t*a[j+k+l];
a[k+l]=Nx+Ny;
a[j+k+l]=Nx-Ny;
}
}
}
if(flag==-1){
for(int i=0;i<len;i++){
a[i]/=len;
}
}
}
void MTT(int *f,int *g,int *h,int n,int p){
for(int i=0;i<n;i++){
a[i].r=(f[i]>>15),a[i].i=(f[i]&32767);
c[i].r=(g[i]>>15),c[i].i=(g[i]&32767);
}
FFT(a,1,n);
FFT(c,1,n);
b[0]=a[0].conj();
for(int i=1;i<n;i++){
b[i]=a[n-i].conj();
}
d[0]=c[0].conj();
for(int i=1;i<n;i++){
d[i]=c[n-i].conj();
}
for(int i=0;i<n;i++){
Complex tmp_a=(a[i]+b[i])*Complex(0.5,0),
tmp_b=(a[i]-b[i])*Complex(0,-0.5),
tmp_c=(c[i]+d[i])*Complex(0.5,0),
tmp_d=(c[i]-d[i])*Complex(0,-0.5);
a[i]=tmp_a*tmp_c+Complex(0,1)*(tmp_a*tmp_d+tmp_b*tmp_c);
b[i]=tmp_b*tmp_d;
}
FFT(a,-1,n);
FFT(b,-1,n);
for(int i=0;i<n;i++){
int tmp_a=((ll)(a[i].r+0.5))%p,tmp_b=((ll)(a[i].i+0.5))%p,tmp_c=((ll)(b[i].r+0.5))%p;
h[i]=((1ll*tmp_a*(1<<30)+1ll*tmp_b*(1<<15)+tmp_c)%p+p)%p;
}
}
}
多项式（存在大量数组未清空的情况，保证思路正确，不保证实现正确）：
namespace Poly{
int add(int x,int y){
return x+y>=Mod?x+y-Mod:x+y;
}
int del(int x,int y){
return x-y<0?x-y+Mod:x-y;
}
int inv[Maxn+5],pow_G[Maxn+5];
void init_G(int len){
int t=quick_power(G,(Mod-1)/len,Mod);
pow_G[len>>1]=1;
for(int i=(len>>1)+1;i<len;i++){
pow_G[i]=1ll*pow_G[i-1]*t%Mod;
}
for(int i=(len>>1)-1;i>0;i--){
pow_G[i]=pow_G[i<<1];
}
}
void init(){
int len=1;
while(len<=Maxn){
len<<=1;
}
len>>=1;
init_G(len);
inv[1]=1;
for(int i=2;i<=Maxn;i++){
inv[i]=1ll*(Mod-Mod/i)*inv[Mod%i]%Mod;
}
}
void NTT(int *a,int flag,int n){
static int R[Maxn+5],last_len=0;
int len=1,L=0;
while(len<n){
len<<=1;
L++;
}
if(last_len!=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]){
std::swap(a[i],a[R[i]]);
}
}
for(int j=1;j<len;j<<=1){
int *t=pow_G+j;
for(int k=0;k<len;k+=(j<<1)){
for(int l=0;l<j;l++){
int Nx=a[k+l],Ny=1ll*t[l]*a[j+k+l]%Mod;
a[k+l]=add(Nx,Ny);
a[j+k+l]=del(Nx,Ny);
}
}
}
if(flag==-1){
std::reverse(a+1,a+len);
for(int i=0,t=inv[len];i<len;i++){
a[i]=1ll*a[i]*t%Mod;
}
}
}
void poly_inv(int *a,int *b,int len){
static int A[Maxn+5],B[Maxn+5];
if(len==1){
b[0]=find_inv(a[0]);
return;
}
poly_inv(a,b,len>>1);
for(int i=(len>>1);i<len;i++){
b[i]=0;
}
for(int i=0;i<len;i++){
A[i]=a[i];
B[i]=b[i];
}
for(int i=len;i<(len<<1);i++){
A[i]=B[i]=0;
}
NTT(A,1,len<<1);
NTT(B,1,len<<1);
for(int i=0;i<(len<<1);i++){
A[i]=1ll*A[i]*B[i]%Mod*B[i]%Mod;
}
NTT(A,-1,len<<1);
for(int i=0;i<len;i++){
b[i]=del(add(b[i],b[i]),A[i]);
}
}
void poly_dif(int *a,int len){
for(int i=1;i<len;i++){
a[i-1]=1ll*i*a[i]%Mod;
}
a[len-1]=0;
}
void poly_int(int *a,int len){
for(int i=len-2;i>=0;i--){
a[i+1]=1ll*inv[i+1]*a[i]%Mod;
}
a[0]=0;
}
void poly_ln(int *a,int *b,int len){
static int A[Maxn+5],B[Maxn+5];
for(int i=0;i<len;i++){
A[i]=a[i];
}
poly_dif(A,len),poly_inv(a,B,len);
NTT(A,1,len<<1);
NTT(B,1,len<<1);
for(int i=0;i<(len<<1);i++){
A[i]=1ll*A[i]*B[i]%Mod;
}
NTT(A,-1,len<<1);
poly_int(A,len);
for(int i=0;i<len;i++){
b[i]=A[i];
}
}
void poly_exp(int *a,int *b,int len){
static int A[Maxn+5],B[Maxn+5];
if(len==1){
b[0]=1;
return;
}
poly_exp(a,b,len>>1);
for(int i=(len>>1);i<len;i++){
b[i]=0;
}
poly_ln(b,A,len);
for(int i=0;i<len;i++){
B[i]=del(a[i],A[i]);
A[i]=b[i];
}
B[0]=add(B[0],1);
NTT(A,1,len<<1);
NTT(B,1,len<<1);
for(int i=0;i<(len<<1);i++){
A[i]=1ll*A[i]*B[i]%Mod;
}
NTT(A,-1,len<<1);
for(int i=0;i<len;i++){
b[i]=A[i];
}
}
namespace Bin{
int bin[Maxn<<5|5];
int *pos=bin;
int* new_array(int len){
int *ans=pos;
pos+=len;
return ans;
}
}
void poly_mul(int *a,int len_a,int *b,int len_b,int *c){
static int A[Maxn+5],B[Maxn+5];
int len=1;
while(len<=(len_a+len_b)){
len<<=1;
}
for(int i=0;i<=len_a;i++){
A[i]=a[i];
}
for(int i=0;i<=len_b;i++){
B[i]=b[i];
}
for(int i=len_a+1;i<len;i++){
A[i]=0;
}
for(int i=len_b+1;i<len;i++){
B[i]=0;
}
NTT(A,1,len);
NTT(B,1,len);
for(int i=0;i<len;i++){
A[i]=1ll*A[i]*B[i]%Mod;
}
NTT(A,-1,len);
for(int i=0;i<=len_a+len_b;i++){
c[i]=A[i];
}
}
int *g[Maxn+5];
void poly_calc(int p,int l,int r,int *a){
g[p]=Bin::new_array(r-l+2);
if(l==r){
g[p][1]=1;
g[p][0]=(a[l]==0?0:Mod-a[l]);
return;
}
int mid=(l+r)>>1;
poly_calc(p<<1,l,mid,a);
poly_calc(p<<1|1,mid+1,r,a);
poly_mul(g[p<<1],mid-l+1,g[p<<1|1],r-mid,g[p]);
}
int *t[Maxn+5];
void poly_get_mul(int *a,int len_a,int *b,int len_b,int *c){
static int A[Maxn+5],B[Maxn+5];
for(int i=0;i<=len_a;i++){
A[i]=a[i];
}
for(int i=0;i<=len_b;i++){
B[i]=b[len_b-i];
}
int len=1;
while(len<=len_a){
len<<=1;
}
for(int i=len_a+1;i<len;i++){
A[i]=0;
}
for(int i=len_b+1;i<len;i++){
B[i]=0;
}
NTT(A,1,len);
NTT(B,1,len);
for(int i=0;i<len;i++){
A[i]=1ll*A[i]*B[i]%Mod;
}
NTT(A,-1,len);
for(int i=len_b;i<=len_a;i++){
c[i-len_b]=A[i];
}
}
void poly_get(int p,int l,int r,int *a){
if(l==r){
a[l]=t[p][0];
return;
}
int mid=(l+r)>>1;
t[p<<1]=Bin::new_array(mid-l+1);
t[p<<1|1]=Bin::new_array(r-mid);
poly_get_mul(t[p],r-l,g[p<<1|1],r-mid,t[p<<1]);
poly_get_mul(t[p],r-l,g[p<<1],mid-l+1,t[p<<1|1]);
poly_get(p<<1,l,mid,a);
poly_get(p<<1|1,mid+1,r,a);
}
void poly_work(int *f,int n,int *a,int m,int *ans){
static int A[Maxn+5],B[Maxn+5];
poly_calc(1,1,m,a);
int len=1;
while(len<=m){
len<<=1;
}
for(int i=0;i<=m;i++){
A[i]=g[1][m-i];
}
poly_inv(A,B,len);
for(int i=m+1;i<len;i++){
B[i]=0;
}
for(int i=0;i<=n;i++){
A[i]=f[n-i];
}
poly_mul(A,n,B,m,A);
t[1]=Bin::new_array(m);
for(int i=0;i<m;i++){
t[1][i]=A[i];
}
for(int i=0;i<=m;i++){
B[i]=0;
}
poly_get(1,1,m,B);
for(int i=1;i<=m;i++){
ans[i]=(f[0]+1ll*a[i]*B[i])%Mod;
}
}
}
using Poly::NTT;
using Poly::poly_inv;
using Poly::poly_ln;
using Poly::poly_exp;
\(O(\frac{n\log^2 n}{\log\log n})\) 的多项式 exp：
namespace Poly{
const int Maxf=(1<<18);
int pow_G[Maxf+5],inv[Maxf+5];
void init(){
int len=1;
while(len<=Maxf){
len<<=1;
}
len>>=1;
pow_G[len>>1]=1;
int t=quick_power(G,(Mod-1)/len,Mod);
for(int i=(len>>1)+1;i<len;i++){
pow_G[i]=1ll*pow_G[i-1]*t%Mod;
}
for(int i=(len>>1)-1;i>=0;i--){
pow_G[i]=pow_G[i<<1];
}
inv[1]=1;
for(int i=2;i<=Maxf;i++){
inv[i]=1ll*(Mod-Mod/i)*inv[Mod%i]%Mod;
}
}
void NTT(int *a,int flag,int n){
static int R[Maxf+5],last_len=0;
static ull A[Maxf+5];
int len=1,L=0;
while(len<n){
len<<=1;
L++;
}
if(len!=last_len){
last_len=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++){
A[i]=a[R[i]];
}
for(int j=1;j<len;j<<=1){
int *t=pow_G+j;
for(int k=0;k<len;k+=(j<<1)){
for(int l=0;l<j;l++){
ull Nx=A[k+l],Ny=A[j+k+l]*t[l]%Mod;
A[k+l]=Nx+Ny;
A[j+k+l]=Nx+Mod-Ny;
}
}
}
for(int i=0;i<len;i++){
a[i]=A[i]%Mod;
}
if(flag==-1){
std::reverse(a+1,a+len);
for(int i=0,t=find_inv(len);i<len;i++){
a[i]=1ll*a[i]*t%Mod;
}
}
}
void poly_exp(int *g,int *f,int l,int r){
if(r-l+1<=64){
for(int i=l;i<=r&&i<=n;i++){
if(i==0){
f[i]=1;
continue;
}
ull sum=0;
for(int j=l;j<i;j++){
sum+=1ll*f[j]*g[i-j];
if(((j-l)&15)==15){
sum%=Mod;
}
}
f[i]=(f[i]+sum)%Mod;
f[i]=1ll*f[i]*inv[i]%Mod;
}
return;
}
int B=16;
int len=(r-l)/B+1;
int pos[B];
int L=1;
while(L<=(len<<1)){
L<<=1;
}
pos[0]=0;
for(int i=1;i<B;i++){
pos[i]=pos[i-1]+L;
}
int *_tmp_f=new int[L*B],*_tmp_g=new int[L*B];
int *tmp_f[B],*tmp_g[B];
for(int i=0;i<B;i++){
tmp_f[i]=_tmp_f+i*L;
tmp_g[i]=_tmp_g+i*L;
for(int j=0;j<L;j++){
tmp_g[i][j]=g[j+i*len];
tmp_f[i][j]=0;
}
NTT(tmp_g[i],1,L);
}
for(int i=0,pos_l=l,pos_r=l+len-1;i<B;i++,pos_l+=len,pos_r=std::min(r,pos_r+len)){
NTT(tmp_f[i],-1,L);
for(int j=pos_l;j<=pos_r;j++){
f[j]=(f[j]+tmp_f[i][j-pos_l+len])%Mod;
}
poly_exp(g,f,pos_l,pos_r);
for(int j=0;j<L;j++){
tmp_f[i][j]=0;
}
for(int j=pos_l;j<=pos_r;j++){
tmp_f[i][j-pos_l]=f[j];
}
NTT(tmp_f[i],1,L);
for(int j=i+1;j<B;j++){
for(int k=0;k<L;k++){
tmp_f[j][k]=(tmp_f[j][k]+1ll*tmp_f[i][k]*tmp_g[j-i-1][k])%Mod;
}
}
}
delete _tmp_f;
delete _tmp_g;
}
void poly_exp(int *a,int *b,int n){
static int A[Maxf+5];
int len=1;
while(len<=n){
len<<=1;
}
for(int i=0;i<=n;i++){
A[i]=1ll*a[i]*i%Mod;
}
for(int i=n+1;i<len;i++){
A[i]=0;
}
poly_exp(A,b,0,n);
}
}
Berlekamp-Massey 算法（欸不对呀这玩意儿貌似不用背）：
namespace Subtask_BM{
std::vector<int> r[Maxn+5];
int r_tot;
int fail_pos[Maxn+5],delta_val[Maxn+5];
std::vector<int> solve(std::vector<int> a){
r_tot=0;
for(int i=0;i<(int)a.size();i++){
int sum=0;
for(int j=0;j<(int)r[r_tot].size();j++){
sum=(sum+1ll*a[i-j-1]*r[r_tot][j])%Mod;
}
if(a[i]==sum){
continue;
}
delta_val[r_tot]=(a[i]-sum+Mod)%Mod;
fail_pos[r_tot]=i;
if(r_tot==0){
r_tot++;
r[r_tot].clear();
r[r_tot].resize(i+1,0);
continue;
}
int min_id=0,min_val=Inf;
for(int j=0;j<r_tot;j++){
if(i-fail_pos[j]+(int)r[j].size()<min_val){
min_val=i-fail_pos[j]+r[j].size();
min_id=j;
}
}
r_tot++;
r[r_tot].clear();
r[r_tot].resize(std::max(min_val,(int)r[r_tot-1].size()),0);
for(int j=0;j<(int)r[r_tot-1].size();j++){
r[r_tot][j]=r[r_tot-1][j];
}
int tmp=1ll*delta_val[r_tot-1]*find_inv(delta_val[min_id])%Mod;
r[r_tot][i-fail_pos[min_id]-1]=(r[r_tot][i-fail_pos[min_id]-1]+tmp)%Mod;
for(int j=0;j<(int)r[min_id].size();j++){
r[r_tot][i-fail_pos[min_id]+j]=(r[r_tot][i-fail_pos[min_id]+j]-1ll*tmp*r[min_id][j]%Mod+Mod)%Mod;
}
}
return r[r_tot];
}
}
本作品采用知识共享署名-非商业性使用-相同方式共享 4.0 国际许可协议进行许可。
