题解 P3803 【【模板】多项式乘法(FFT)】

小菜鸟

2019-11-05 20:26:48

Solution

这题...我可能是唯一一个用三种语言各写一遍的人 ~~然而Pascal卡不过去~~ 并且我应该是第一个用C写FFT通过的( 所以这篇题解的重点不在数学,而在如何用C愉快地食用FFT( ~~既然前面有P党卡过去发题解那这个应该也行吧~~ --- 先重复一下做法:把两个多项式变换成点值表示,直接$O(n)$乘起来,然后变回系数表示就是了。。。 --- 然后进正题 由于C++有面向对象、函数重载等一系列特性,复数运算很容易以优美的方式呈现 然而C既不能重载函数,也难以面向对象... 然而在研究`std::complex`的时候,我们发现它的头文件`include`了一个`complex.h`库... 于是上网学习,发现C也有复数库(其实甚至不需要库,在C99中复数已经是内建类型了) 于是就极度舒适( --- FFT需要的基本用法: ```c #include<complex.h>//为了方便而包含,实际可以省略,但声明变量需将complex改为_Complex double complex a;//声明基于double的复数变量 creal(a)//取a的实部 cimag(a)//取a的虚部 a=1+2*_Complex_I//_Complex_I即是单位虚数i //另外不担心重名可以用I代替_Complex_I //然后不考虑可移植性的话,在GCC下提供1.0fi,1.0i,1.0li等写法 //同时从C++14开始,类似的1.0if,1.0i,1.0il是标准的虚数字面量 ``` 然后大量的C数学函数均有针对复数特化的泛型函数,在名字前写前缀`c`就是了 于是有了神奇的复数,我们就可以愉快地搞FFT了( --- 然后是不开IO优化卡不过去的代码 ```cpp #include<stdio.h> #include<ctype.h> #include<math.h> #include<complex.h> char gc(void) { static char buf[1<<16],*p1=buf,*p2=buf; if(p1==p2) { p2=(p1=buf)+fread(buf,1,1<<16,stdin); if(p1==p2)return EOF; } return *p1++; } // #define gc getchar void read(int *x) { *x=0; char c=gc(); while(!isdigit(c))c=gc(); while(isdigit(c)) { *x=*x*10+c-48; c=gc(); } } void write(int x) { if(x>=10)write(x/10); putchar(x%10+48); } #define N (1<<21) #define PI 3.141592653589793238462643383 int n,m; double complex a[N],b[N]; int rev[N]; void getrev(int n) { for(int i=0;i<n;++i) { rev[i]=(rev[i>>1]>>1)|(i&1?n>>1:0); } } void fft(double complex *a,int n,int dft) { for(int i=0;i<n;++i) { if(i<rev[i]) { double complex temp=a[i]; a[i]=a[rev[i]]; a[rev[i]]=temp; } } for(int step=1;step<n;step<<=1) { const double complex wn=cos(PI/step)+sin(dft*PI/step)*_Complex_I; for(int j=0;j<n;j+=step<<1) { double complex wnk=1; for(int k=j;k<j+step;++k) { double complex x=a[k]; double complex y=a[k+step]*wnk; a[k]=x+y; a[k+step]=x-y; wnk*=wn; } } } if(dft==-1) { double inv=1./n; for(int i=0;i<n;++i)a[i]*=inv; } } int main(void) { read(&n),read(&m); for(int i=0;i<=n;++i) { int x; read(&x); a[i]=x; } for(int i=0;i<=m;++i) { int x; read(&x); b[i]=x; } int len=1; while(len<=n+m+1)len<<=1; getrev(len); fft(a,len,1); fft(b,len,1); for(int i=0;i<len;++i)a[i]*=b[i]; fft(a,len,-1); for(int i=0;i<n+m+1;++i) { write((int)(creal(a[i])+0.5)); putchar(' '); } } ``` ~~awsl为什么这题突然减小时限卡常数啊~~ 警告,若编译器有`__STDC_NO_COMPLEX__`宏,则不支持复数(gcc可以放心 最后,虽然我是个C++&面向对象的死忠粉,但还是要说一句: ## C牛逼!!! 另外为什么C的语法高亮萎掉了啊求修复(被迫写`cpp`QwQ --- ### 2020.10.14更新 初赛前颓废学了一下Rust 考虑到FFT代码相对简单于是就来练练手 --- 简要介绍Rust,一个安全性和效率都相当不错的新兴语言 由Mozilla的大佬创造维护 定位和C++类似,零开销,多范式 但具有更多安全的、流行的特性 基本上可以把UB扼杀在编译期 ~~而且拥有比GNU STL清晰一些的标准库~~ --- 这篇题解本来就不是讲思路而是讲代码的... 所以下面直接给出代码实现(PS下面有效率对比) ```rust const P:i64=998244353; const g:i64=3; const gi:i64=332748118; fn power(mut a:i64,mut b:i64)->i64 { let mut res:i64=1; while b!=0 { if (b&1)!=0 { res=res*a%P; } a=a*a%P; b>>=1; } return res; } fn getrev(rev:&mut Vec<usize>,len:usize) //Rust对于所有权的控制很严格,只能用可变引用来改变外部的东西 //Rust的引用与C++的指针类似 { rev.clear(); rev.push(0); for i in 1..len { rev.push((rev[i>>1]>>1)|(if (i&1)==0{0}else{len>>1})); } } fn ntt(arr:&mut Vec<i64>,n:usize,dft:bool,rev:&Vec<usize>) { for i in 0..n { if i<rev[i] { arr.swap(i,rev[i]);//在Rust中,一个对象不能同时存在两个可变引用 //所以直接std::mem::swap(i,rev[i])是unsafe的,需要用Vec自己的方法 } } let mut step=1; while step<n { let wn:i64=power(if dft{g}else{gi},(P-1)/(step as i64*2)); let mut j=0; while(j<n) { let mut wnk=1; for k in j..j+step { let x=arr[k]; let y=arr[k+step]*wnk%P; arr[k]=(x+y)%P; arr[k+step]=(x-y+P)%P; wnk=wnk*wn%P; } j+=step*2; } step*=2; } if !dft { let inv=power(n as i64,P-2); for i in 0..n { arr[i]=arr[i]*inv%P; } } } extern "C" { fn getchar()->u8; } fn read()->i64//rust没有原生的比较简洁的快读,不得不学习巨佬cosmicAC的 { let mut res=0; let mut ch; loop { unsafe { ch=getchar(); } if ch>=48 && ch<=57 { break; } } loop { res=res*10+(ch as i64)-48; unsafe { ch=getchar(); } if ch<48 || ch>57 { break; } } return res; } fn main() { let n:usize=read() as usize;//Rust在safe mode里不能定义全局变量 let m:usize=read() as usize; let n=n+1; let m=m+1; let mut a:Vec<i64>=Vec::new(); let mut b:Vec<i64>=Vec::new(); let mut len=1; while len<=n+m { len*=2; } for i in 0..n { let x=read(); a.push(x); } while a.len()<len { a.push(0); } for i in 0..m { let x=read(); b.push(x); } while b.len()<len { b.push(0); } let mut rev:Vec<usize>=Vec::new(); for i in 0..len { rev.push(0); } getrev(&mut rev,len); ntt(&mut a,len,true,&rev); ntt(&mut b,len,true,&rev); let mut res:Vec<i64>=Vec::new(); for i in 0..len { res.push(a[i]*b[i]%P); } ntt(&mut res,len,false,&rev); for i in 0..n+m-1 { print!("{} ",res[i]); } } ``` C++在使用普通快读的情况下是1.4s(有没有`fread`区别不大) Rust是1.8s... 看来GCC开发这么多年毕竟是比rustc成熟的233