具体原理就不讲了qwq,毕竟证明我也不太懂

FFT(快速傅立叶变换)&NTT(快速数论变换)

FFT

 //求多项式乘积
//要求多项式A和多项式B的积多项式C
//具体操作就是
//DFT(A),DFT(B)->暴力乘积->拉格朗日插值(即IDFT(C))->C
//其中DFT表示离散傅里叶变换
//通俗的来说就是用点值表示多项式
//使用神秘单位复数根将时间复杂度降至O(nlogn)
//ps:但是常数巨大
//pps:应用非常广泛,非常多题目都要fft or ntt优化,板子一定要背熟
#include<iostream>
#include<cstring>
#include<cstdio>
#include<cmath>
#define pw(n) (1<<n)
using namespace std;
const double pi=acos(-);
struct complex{
double a,b;
complex(double _a=,double _b=){
a=_a;
b=_b;
}
friend complex operator +(complex x,complex y){return complex(x.a+y.a,x.b+y.b);}
friend complex operator -(complex x,complex y){return complex(x.a-y.a,x.b-y.b);}
friend complex operator *(complex x,complex y){return complex(x.a*y.a-x.b*y.b,x.a*y.b+x.b*y.a);}
friend complex operator *(complex x,double y){return complex(x.a*y,x.b*y);}
friend complex operator /(complex x,double y){return complex(x.a/y,x.b/y);}
}a[],b[];
int n,m,bit,bitnum=,rev[pw()];
void getrev(int l){//Reverse
for(int i=;i<pw(l);i++){
rev[i]=(rev[i>>]>>)|((i&)<<(l-));
}
}
void FFT(complex *s,int op){
for(int i=;i<bit;i++)if(i<rev[i])swap(s[i],s[rev[i]]);
for(int i=;i<bit;i<<=){
complex w(cos(pi/i),op*sin(pi/i));
for(int p=i<<,j=;j<bit;j+=p){//Butterfly
complex wk(,);
for(int k=j;k<i+j;k++,wk=wk*w){
complex x=s[k],y=wk*s[k+i];
s[k]=x+y;
s[k+i]=x-y;
}
}
}
if(op==-){
for(int i=;i<=bit;i++){
s[i]=s[i]/(double)bit;
}
}
}
int main(){
scanf("%d%d",&n,&m);
for(int i=;i<=n;i++)scanf("%lf",&a[i].a);
for(int i=;i<=m;i++)scanf("%lf",&b[i].a);
m+=n;
for(bit=;bit<=m;bit<<=)bitnum++;
getrev(bitnum);
FFT(a,);
FFT(b,);
for(int i=;i<=bit;i++)a[i]=a[i]*b[i];
FFT(a,-);
for(int i=;i<=m;i++)printf("%d ",(int)(a[i].a+0.5));
return ;
}

NTT

 //费马数数论变换
//大家觉得998244353好还是1004535809好?^_^
#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdio>
#include<cmath>
#define pw(n) (1<<n)
using namespace std;
const int N=,P=,g=;//或P=1004535809
int n,m,bit,bitnum=,a[N+],b[N+],rev[N+];
void getrev(int l){
for(int i=;i<pw(l);i++){
rev[i]=(rev[i>>]>>)|((i&)<<(l-));
}
}
int fastpow(int a,int b){
int ans=;
for(;b;b>>=,a=1LL*a*a%P){
if(b&)ans=1LL*ans*a%P;
}
return ans;
}
void NTT(int *s,int op){
for(int i=;i<bit;i++)if(i<rev[i])swap(s[i],s[rev[i]]);
for(int i=;i<bit;i<<=){
int w=fastpow(g,(P-)/(i<<));
for(int p=i<<,j=;j<bit;j+=p){
int wk=;
for(int k=j;k<i+j;k++,wk=1LL*wk*w%P){
int x=s[k],y=1LL*s[k+i]*wk%P;
s[k]=(x+y)%P;
s[k+i]=(x-y+P)%P;
}
}
}
if(op==-){
reverse(s+,s+bit);
int inv=fastpow(bit,P-);
for(int i=;i<bit;i++)a[i]=1LL*a[i]*inv%P;
}
}
int main(){
scanf("%d%d",&n,&m);
for(int i=;i<=n;i++)scanf("%d",&a[i]);
for(int i=;i<=m;i++)scanf("%d",&b[i]);
m+=n;
for(bit=;bit<=m;bit<<=)bitnum++;
getrev(bitnum);
NTT(a,);
NTT(b,);
for(int i=;i<bit;i++)a[i]=1LL*a[i]*b[i]%P;
NTT(a,-);
for(int i=m;i>=;i--)printf("%d ",a[i]);
return ;
}

最新文章

  1. Mybatis的mapper文件中$和#的区别
  2. 【noip2014T3】
  3. No bootable device-insert boot disk and press any key
  4. UVa 10868 (物理) Bungee Jumping
  5. 小C的填数游戏
  6. 【2013长沙区域赛】部分题解 hdu4791—4801
  7. SQL 事务及实例演示
  8. [置顶] ios 一个不错的图片浏览分享框架demo
  9. Intent的属性及Intent-filter配置——Data、Type属性与intent-filter配置
  10. css 清除浮动的方法
  11. Windows搭建以太坊的私有链环境
  12. Kotlin 初窥门径[1]:基础概念
  13. ORACLE中seq$表更新频繁的分析
  14. js获取地址栏参数数据
  15. Java基础语法(上篇)
  16. HTTP认证方式与https简介
  17. centos 7.5+如何格式化硬盘
  18. javafx由浅到深的 认识(一)
  19. 1.5 A better alternative thing: React Native
  20. Android OpenGL ES 开发(N): OpenGL ES 2.0 机型兼容问题整理

热门文章

  1. POJ 1611 The Suspects【并查集】
  2. 使用layer.tips实现鼠标悬浮时触发事件提示消息实现
  3. java中的string trim具体有什么用处。。。
  4. bootstrap fileinput控件 + django后台上传、回显简单使用
  5. libvirtd.service
  6. [转载] Linux新手必看:浅谈如何学习linux
  7. myEclies项目导入Eclipse中常见问题
  8. class的基本操作方法
  9. MyBatis学习总结(18)——MyBatis与Hibernate区别
  10. UVALive 5545 Glass Beads