网站建设一般收多少定金,wordpress的用户,软件跟网站开发的区别,沈阳京科医院是正规医院吗文章目录前言系数表示法和点值表示法单位根离散傅立叶变换#xff08;DFT#xff09;位逆序置换#xff08;蝴蝶变换#xff09;离散傅立叶逆变换#xff08;IDFT#xff09;代码NTT代码所谓快速傅立叶变换#xff0c;就是傅立叶发明的一种快速的变换 #xff08;逃DFT位逆序置换蝴蝶变换离散傅立叶逆变换IDFT代码NTT代码所谓快速傅立叶变换就是傅立叶发明的一种快速的变换 逃
前言
FFT用于在nlognnlognnlogn的复杂度内求两个多项式相乘。 算是多项式的入门吧。 十一被rabbit的数学打蒙了所以开始碰这些东西。 不知道以后会不会直接鸽掉。 update on 2022.1.4:三个月后终于回来填坑了。 很妙的算法。 代码写成迭代精炼之后很优美。
系数表示法和点值表示法
系数表示法就是我们常用的函数的表示形式。 对于一个n次的函数点值表示法则是给出了n1个互异的点通过这n1个互异的点就可以唯一确定的描述出这个函数可以想想高斯消元。
这东西有啥用 考虑如果两个多项式都表示成了点值形式且选取的点的横坐标相同我把它们的纵坐标对应乘起来就能得到它们乘积的多项式的点值表示 乘起来的复杂度降到了O(n)O(n)O(n)
但是这个点值表示法求起来暴力似乎还是n方啊… 所以快速傅立叶变换其实解决的就是系数表示法和点值表示法快速互相转换的问题
单位根
规定 若复数ω\omegaω满足ωn1\omega^n1ωn1就称ω\omegaω是n次单位根 对于一个n项的多项式我们让它的点值表示取的横坐标为ωnk0≤kn)\omega_n^k0\leq k n)ωnk0≤kn) 然后就会有很多很妙的性质
ωnmkmωnk\omega_{nm}^{km}\omega_n^kωnmkmωnkωnknωnk\omega_{n}^{kn}\omega_n^kωnknωnkωnkn2−wnk\omega_n^{k\frac{n}{2}}-w_n^kωnk2n−wnk
都可以结合几何意义较为显然的证明。
离散傅立叶变换DFT
那么回到问题。 首先我们要把多项式填零直到 n2kn2^kn2k 的长度。 接下来我们要求 f(x)a0a1xa2x2a3x3...an−1xn−1f(x)a_0a_1xa_2x^2a_3x^3...a_{n-1}x^{n-1}f(x)a0a1xa2x2a3x3...an−1xn−1 让我们把这个函数分奇偶拆成两个函数 A(x)a0a2xa4x2a6x3...an−1xn2A(x)a_0a_2xa_4x^2a_6x^3...a_{n-1}x^{\frac{n}{2}}A(x)a0a2xa4x2a6x3...an−1x2n B(x)a1a3xa5x2a7x3...an−2xn2B(x)a_1a_3xa_5x^2a_7x^3...a_{n-2}x^{\frac{n}{2}}B(x)a1a3xa5x2a7x3...an−2x2n 那么就有 f(x)A(x2)xB(x2)f(x)A(x^2)xB(x^2)f(x)A(x2)xB(x2) 我们把 xωnkx\omega_n^kxωnk 带入。 分情况讨论 设0≤kn20\le k\dfrac{n}{2}0≤k2n 当指数 n2 \dfrac{n}{2}2n 时 f(ωnk)A((ωnk)2)ωnkB((ωnk)2)f(\omega_n^k)A((\omega_n^k)^2)\omega_n^kB((\omega_n^k)^2)f(ωnk)A((ωnk)2)ωnkB((ωnk)2) A(ωn2k)ωnkB(ωn2k)A(\omega_n^{2k})\omega_n^kB(\omega_n^{2k})A(ωn2k)ωnkB(ωn2k) A(ωn2k)ωnkB(ωn2k)A(\omega_{\frac{n}{2}}^{k})\omega_n^kB(\omega_{\frac{n}{2}}^{k})A(ω2nk)ωnkB(ω2nk) 类似的当指数 ≥n2\ge\dfrac{n}{2}≥2n 时 f(ωnkn2)A(ωn2kn)ωnkn2B(ωn2kn)f(\omega_n^{k\frac{n}{2}})A(\omega_n^{2kn})\omega_n^{k\frac{n}{2}}B(\omega_n^{2kn})f(ωnk2n)A(ωn2kn)ωnk2nB(ωn2kn) A(ωn2k)−ωnkB(ωn2k)A(\omega_n^{2k})-\omega_n^{k}B(\omega_n^{2k})A(ωn2k)−ωnkB(ωn2k) A(ωn2k)−ωnkB(ωn2k)A(\omega_{\frac{n}{2}}^{k})-\omega_n^{k}B(\omega_{\frac{n}{2}}^{k})A(ω2nk)−ωnkB(ω2nk) 这样我们就可以在 O(len)O(len)O(len) 的复杂度内合并两个长度为 lenlenlen 的多项式就可以利用分治把复杂度降到 O(nlogn)O(n\log n)O(nlogn)。
位逆序置换蝴蝶变换
这个也很妙 有点线性求逆元那味了 通过观察发现系数 i 调换位置之后会去往它的二进制表示反过来的位置记为rir_iri 而上面那个东西可以通过 r[i]r[i1]1[i1]∗(n1)r[i]r[i1]1[i\1]*(n1) r[i]r[i1]1[i1]∗(n1) 来线性的求 这个还是挺好理解的画一画就行了
离散傅立叶逆变换IDFT
我们变完之后当然还要变回来。 我会 n3n^3n3 高斯消元 设之前DFT得到的点值表示 bkf(ωnk)∑i0n−1ai(ωnk)ib_kf(\omega_n^k)\sum_{i0}^{n-1}a_i(\omega_{n}^k)^ibkf(ωnk)∑i0n−1ai(ωnk)i 构造一个函数 g(x)b0b1xb2x2b3x3...bn−1xn−1g(x)b_0b_1xb_2x^2b_3x^3...b_{n-1}x_{n-1}g(x)b0b1xb2x2b3x3...bn−1xn−1 然后我们用 ωn0,ωn−1,ωn−2,...,ωn−n\omega_n^{0},\omega_n^{-1},\omega_n^{-2},...,\omega_n^{-n}ωn0,ωn−1,ωn−2,...,ωn−n 作为横坐标带入。不难发现刚才 DFT 需要的性质依然成立所以可以用同样的方法求解 那么我们在第 kkk 位得到的新的点值就是 g(ωn−k)∑i0n−1bi(ωn−k)ig(\omega_n^{-k})\sum_{i0}^{n-1}b_i(\omega_n^{-k})^ig(ωn−k)i0∑n−1bi(ωn−k)i ∑i0n−1∑j0n−1aj(ωni)j(ωn−k)i\sum_{i0}^{n-1}\sum_{j0}^{n-1}a_j(\omega_{n}^i)^j(\omega_n^{-k})^ii0∑n−1j0∑n−1aj(ωni)j(ωn−k)i ∑j0j−1aj∑i0n−1(ωnj−k)i\sum_{j0}^{j-1}a_j\sum_{i0}^{n-1}(\omega_{n}^{j-k})^ij0∑j−1aji0∑n−1(ωnj−k)i 考虑 ∑i0n−1(ωnj−k)i\sum_{i0}^{n-1}(\omega_{n}^{j-k})^i∑i0n−1(ωnj−k)i这一项当 jkjkjk 时显然等于0当 j≠kj\ne kjk 时根据等比公式有 ∑i0n−1(ωnj−k)i(wnj−k)n−1wnj−k−1\sum_{i0}^{n-1}(\omega_{n}^{j-k})^i\frac{(w_n^{j-k})^n-1}{w_n^{j-k}-1}i0∑n−1(ωnj−k)iwnj−k−1(wnj−k)n−1 (wnn)j−k−1wnj−k−11j−k−1wnj−k−10\frac{(w_n^{n})^{j-k}-1}{w_n^{j-k}-1}\frac{1^{j-k}-1}{w_n^{j-k}-1}0wnj−k−1(wnn)j−k−1wnj−k−11j−k−10 所以原来的式子只有当 jkjkjk 的时候有值等于 nakna_knak。 所以我们直接 NTT 完除 nnn 即可。
代码
#includebits/stdc.h
using namespace std;
#define ll long long
#define il inline
const int N5e6100;
const int M150;
const int mod998244353;
const double piacos(-1.0);
inline ll read(){ll x0,f1;char cgetchar();while(!isdigit(c)){if(c-) f-1;cgetchar();}while(isdigit(c)){xx*10c-0;cgetchar();}return x*f;
}
int n,m,len,k;
struct node{double x,y;node(double a0,double b0){xa;yb;}
}A[N],B[N];
il node operator * (node a,node b){return (node){a.x*b.x-a.y*b.y,a.x*b.ya.y*b.x};
}
il node operator (node a,node b){return (node){a.xb.x,a.yb.y};
}
il node operator - (node a,node b){return (node){a.x-b.x,a.y-b.y};
}
int r[N];
il void fft(node *x,int flag){for(int i0;ilen;i){if(ir[i]) swap(x[i],x[r[i]]);}for(int l1;llen;l1){node o(cos(pi/l),flag*sin(pi/l));for(int st0;stlen;stl1){node t(1,0);for(int j0;jl;j,tt*o){node ux[stj],vt*x[stjl];x[stj]uv;x[stjl]u-v;}}}return;
}
int main(){nread();mread();for(int i0;in;i) A[i].xread();for(int i0;im;i) B[i].xread();lennm;while((1k)len) k;len1k;for(int i1;ilen;i) r[i](r[i1]1)|((i1)(k-1));fft(A,1);fft(B,1);for(int i0;ilen;i) A[i]A[i]*B[i];fft(A,-1);for(int i0;inm;i) printf(%d ,(int)(A[i].x/(len)0.5));return 0;
}
/**/NTT
NTT说简单也简单说难也难。 简单是因为板子几乎和FFT一模一样。 难是因为那个群的概念完全看不懂qwq。 然后我就选择调过了证明部分。 背下来背下来——宋队 总的来说写一写对背板子有帮助的我可以理解的。 对于一个质数P,若其存在原根基本全是3且可以表示为 Pk∗2n1Pk*2^n1Pk∗2n1 它就是一个ntt模数。 ggg为什么叫原根 因为它有一些关键性质 1.gk∗n1(modP)1. g^{k*n}1(mod P)1.gk∗n1(modP) 2.gi(modP)g^i(mod P)gi(modP) 在 0≤i≤p0\leq i\leq p0≤i≤p 的范围内取值两两不同。 设gkg^kgk 为单位根记为gng_ngn。 那么它就有很多熟悉的性质
gnn1(modP)g_n^n1(mod P)gnn1(modP)g2k2kgkkg_{2k}^{2k}g_k^kg2k2kgkkgnkn/2−gnkg_n^{kn/2}-g_n^kgnkn/2−gnk
没错推FFT时 ω\omegaω 需要有的性质它全有 所以就可以直接上FFT的板子啦 然后烦人的复数、浮点运算和精度问题全都拜拜了。 针不戳。 宋队诚不欺我。 逆变换还是倒数从FFT的求共轭改成求逆元。 就像喝水一样。
代码
il void ntt(ll *x,int flag){for(int i0;ilim;i){if(ir[i]) swap(x[i],x[r[i]]);}for(register int l1;llim;l1){register ll gksm(3,(mod-1)/(l1));if(flag-1) gksm(g,mod-2);for(register int st0;stlim;stl1){register ll now1;for(register int i0;il;i,nownow*g%mod){ll ux[sti],vx[stil]*now%mod;x[sti](uv)%mod;x[stil](u-vmod)%mod;}}}if(flag-1){register ll niksm(lim,mod-2);for(register int i0;ilim;i) x[i]x[i]*ni%mod;}
}
