參考（大部分證明摘自）：https://oi.men.ci/fft-notes/
【簡介】
快速傅里葉變換（FFT）是一種可以在$O(nlogn)$時間內完成的離散傅里葉變換（DFT）算法，在OI中主要用於加速向量卷積/多項式乘法運算。
【前置技能】
【引入】
有兩個多項式$A(x)$和$B(x)$，求$C(x)=A(x)*B(x)$。
$A(x)=\sum_{i=0}^{n-1}a_ix^i$
$B(x)=\sum_{i=0}^{n-1}b_ix^i$
$C(x)=A(x)*B(x)=\sum_{i=0}^{2n-2}(\sum_{j=0}^i a_j*b_{i-j})x^i$
用朴素的做法復雜度顯然是$O(n^2)$的。
如果再推推？
$C(x)=A(x)*B(x)$
$=\sum_{i=0}^{2n-2}(\sum_{j=0}^i a_j*b_{i-j})x^i$
$=\sum_{i=0}^{2n-2}\sum_{j=0}^i(a_j*x^j)(b_{i-j}*x^{i-j})$
$=\sum_{j=0}^{2n-2}(a_j*x^j)\sum_{i=0}^{2n-2}(b_i*x^i)$
式子左邊和右邊分別是多項式$A(x),B(x)$的點值表示，求一個多項式的點值表示和通過點值表示求多項式都可以通過FFT和IDFT做到$O(nlogn)$，所以我們就可以在$O(nlogn)$的時間進行多項式乘法運算了。
【系數表示法】
設$A(x)$表示一個$n-1$次多項式，所有項的系數組成的$n$維向量$(a_0, a_1, a_2, ..., a_{n-1})$唯一確定了一個多項式，這種表示法叫系數表示法。
$A(x)=\sum_{i=0}^{n-1}a_ix^i$
【點值表示法】
如果將一組互不相同的$n$個$x$代入多項式，得到$n$個值，組成的$n$維向量為$(x_0, x_1, x_2, ..., x_{n-1}),(y_0, y_1, y_2, ..., y_{n-1})$，唯一確定了一個多項式（證明在下方），這種表示法叫點值表示法。
$y_i=A(x_i)=\sum_{j=0}^{n-1} a_jx_i^j$
一個$n-1$次多項式在$n$個不同點的取值唯一確定一個多項式的證明：
若命題不成立，則存在兩個不同的$n-1$次多項式$A(x),B(x)$滿足$\forall i\subseteq [0,n-1]，A(x_i)=B(x_i)$。
令$C(x)=A(x)-B(x)$，則$C(x)$也是一個$n-1$次多項式，且滿足$\forall i\subseteq [0,n-1]，C(x_i)=0$。
即$C(x)$有$n$個根，與代數基本定理（一個$n-1$次多項式在復數域上至多只有$n-1$個根）相矛盾，故原命題成立。
【卷積】
$(f*g)(n)$稱為$f,g$的卷積，其離散的定義為：
$(f*g)(n)=\sum_{-\infty }^{\infty}f(\tau)g(n-\tau)$
特征為$\tau+(n-\tau)=n$
【復數、復平面】相關知識見：http://www.ruanyifeng.com/blog/2012/09/imaginary_number.html
【單位根】
在復平面上以原點為起點，單位圓的$n$等分點為終點，作$n$個向量，所得的幅角為正且最小的向量對應的復數$\omega _n$，稱為$n$次單位根，則其余的$n-1$個向量對應的復數分別為$\omega _n^2,\omega _n^3...,\omega _n^n$，顯然$\omega _n^0=\omega _n^n=1$。
因為單位根的幅角為$\frac{1}{n}$，所以有$\omega _n^k=cos k\frac{2\pi}{n}+i sin k\frac{2\pi}{n}$。
【單位根的性質】
①$\omega _{2n}^{2k}=\omega _n^k$
②$\omega _n^{k+\frac{n}{2}}=-\omega _n^k$
結合單位圓較易理解，證明略
【快速傅里葉變換(FFT)】
點值向量$(A(\omega_n^0),A(\omega_n^1),A(\omega_n^2),...,A(\omega_n^{n-1}))$被稱為其系數向量$(a_0, a_1, a_2, ..., a_{n-1})$的離散傅里葉變換。
把$A(x)$按奇偶次項分組。
令
$A_1(x)=(a_0+a_2x+a_4x^2+...+a_{n-2}x^{\frac{n}{2}-1})$
$A_2(x)=(a_1+a_3x+a_5x^2+...+a_{n-1}x^{\frac{n}{2}-1})$
$A(x)=A_1(x^2)+xA_2(x^2)$
若$k<\frac{n}{2}$，求$A(\omega _n^k)$：
$A(\omega _n^k)=A_1(\omega _n^{2k})+\omega _n^kA_2(\omega _n^{2k})$
$=A_1(\omega _{\frac{n}{2}}^k)+\omega _n^kA_2(\omega _{\frac{n}{2}}^k)$
求$A(\omega _n^{k+\frac{n}{2}})$：
$A(\omega_n^{k+\frac{n}{2}})=A_1(\omega _n^{2k+n})+\omega _n^{k+\frac{n}{2}}A_2(\omega _n^{2k+n})$
$=A_1(\omega _n^{2k}\times \omega_n^n)-\omega_n^kA_2(\omega_n^{2k}\times \omega_n^n)$
$=A_1(\omega_n^{2k})-\omega_n^kA_2(\omega_n^{2k})$
$=A_1(\omega_{\frac{n}{2}}^k)-\omega_n^kA_2(\omega_{\frac{n}{2}}^k)$
根據上方的推導，我們顯然可以使用一個分治的算法來求得所有的點值，效率$O(nlogn)$。
這就是Cooley-Tukey算法。
【傅里葉逆變換(IDFT)】
將點值表示的多項式轉換為系數表示，可以使用傅里葉逆變換。
設$(y_0, y_1, y_2, ..., y_{n-1})$為$(a_0, a_1, a_2, ..., a_{n-1})$的離散傅里葉變換，設另一個向量$(c_0, c_1, c_2, ..., c_{n-1})$滿足
$c_k=\sum_{i=0}^{n-1}y_i(\omega _n^{-k})^i$
$=\sum_{i=0}^{n-1}\sum_{j=0}^{n-1}(a_j(\omega _n^i)^j)(\omega_n^{-k})^i$
$=\sum_{i=0}^{n-1}\sum_{j=0}^{n-1}(a_j(\omega _n^j)^i(\omega_n^{-k})^i)$
$=\sum_{i=0}^{n-1}\sum_{j=0}^{n-1}a_j(\omega_n^{j-k})^i$
$=\sum_{j=0}^{n-1}a_j(\sum_{i=0}^{n-1}(\omega_n^{j-k})^i)$
設
$S(\omega_n^k)=1+\omega_n^k+(\omega_n^k)^2+...+(\omega_n^k)^{n-1}$
$k=0$時，顯然$S(\omega_n^k)=n$
$k\neq 0$時，則等式兩邊乘上$\omega_n^k$
$\omega_n^kS(\omega_n^k)=\omega_n^k+(\omega_n^k)^2+(\omega_n^k)^3+...+(\omega_n^k)^{n}$
兩式相減，得
$\omega_n^kS(\omega_n^k)-S(\omega_n^k)=(\omega_n^k)^n-1$
$S(\omega_n^k)=\frac{(\omega_n^k)^n-1}{\omega_n^k-1}$
分子為0，分母不為0，故
$S(\omega_n^k)=0$
考慮上式
$c_k=\sum_{j=0}^{n-1}a_j(\sum_{i=0}^{n-1}(\omega_n^{j-k})^i)$
$=\sum_{j=0}^{n-1}a_jS(\omega_n^{j-k})$
當$j=k$時，$S(\omega_n^{j-k})=n$，否則為0，故
$c_i=na_i$
$a_i=\frac{1}{n}c_i$
所以使用單位根的倒數作為$x$，做一次FFT，對得到的每個數除以$n$，即可得到點值表示對應的系數表示。
【實現】
【迭代FFT】
遞歸版FFT速度十分慢，所以一般使用迭代版的FFT。
觀察規律可知FFT分治到終止時，某個編號的初始下標是終止下標二進制位翻轉后的值，然后就可以模擬分治進行迭代FFT了。
【蝴蝶操作】
運用一個臨時變量省去一個中轉數組。
若此時$a(k)$為$A_1(\omega_{\frac{n}{2}}^k)$，$a(k+\frac{n}{2})$為$A_2(\omega_{\frac{n}{2}}^k)$，用$a(k)+\omega_n^k\times a(\frac{n}{2}+k)$更新$a(k)$，用$a(k)-\omega_n^k\times a(\frac{n}{2}+k)$更新$a(\frac{n}{2}+k)$時，用一個臨時變量$t$進行以下操作。
$t=\omega_n^k\times a(\frac{n}{2}+k)$
$a(\frac{n}{2}+k)=a(k)-t$
$a(k)=a(k)+t$
【注意】最后求系數取整時為減少誤差一般+0.5以四舍五入
【流程】
將數從初始位置換到終止位置。
模擬分治進行迭代FFT。
【例題】
裸的大數乘法，顯然兩個大數可以看成兩個多項式
#include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> #include<algorithm> using namespace std; const int maxn=500010, inf=1e9; const double pi=acos(-1); struct cpx{double r, i; cpx(double _r=0, double _i=0):r(_r),i(_i){};}a[maxn], b[maxn], c[maxn]; cpx operator + (cpx a, cpx b){return cpx(a.r+b.r, a.i+b.i);} cpx operator - (cpx a, cpx b){return cpx(a.r-b.r, a.i-b.i);} cpx operator * (cpx a, cpx b){return cpx(a.r*b.r-a.i*b.i, a.r*b.i+a.i*b.r);} int n, mx; int xs[maxn]; char s[maxn]; inline void fft(int N, cpx a[], int f) { for(int i=0, j=0;i<N;i++) { if(i<j) swap(a[i], a[j]); for(int k=N>>1;(j^=k)<k;k>>=1); } for(int i=2;i<=N;i<<=1) for(int j=0, m=i>>1;j<N;j+=i) { cpx nw=cpx(cos(2.*pi/i), sin(2.*pi/i)*f), w=cpx(1, 0); for(int k=0;k<m;k++) { cpx t=a[j+k+m]*w; a[j+k+m]=a[j+k]-t; a[j+k]=a[j+k]+t; w=w*nw; } } } int main() { scanf("%d", &n); for(mx=1;mx<(n<<1);mx<<=1); scanf("%s", s+1); for(int i=0;i<n;i++) a[i]=s[n-i]-'0'; scanf("%s", s+1); for(int i=0;i<n;i++) b[i]=s[n-i]-'0'; fft(mx, a, 1); fft(mx, b, 1); for(int i=0;i<mx;i++) c[i]=a[i]*b[i]; fft(mx, c, -1); for(int i=0;i<mx;i++) xs[i]=(int)(c[i].r/mx+0.5); for(int i=0;i<mx;i++) if(xs[i]>=10) { xs[i+1]+=xs[i]/10, xs[i]%=10; if(i==mx-1) mx++; } mx--; while(!xs[mx] && mx>0) mx--; for(int i=mx;~i;i--) printf("%d", xs[i]); }
$c[k]=\sum_{i=k}^{n-1} a[i]*b[i-k]$不是卷積的形式，但化一化把$a$倒過來就是卷積的形式了。
$c[k]=ans[n-k-1]=\sum_{i=0}^{n-k-1}a[n-k-i-1]*b[i]$，所以FFT加速后輸出$ans[n-1]$~$ans[0]$即可。
#include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> #include<algorithm> using namespace std; const int maxn=500010, inf=1e9; const double pi=acos(-1); struct cpx{double r, i; cpx(double _r=0, double _i=0):r(_r),i(_i){};}a[maxn], b[maxn], c[maxn]; cpx operator + (cpx a, cpx b){return cpx(a.r+b.r, a.i+b.i);} cpx operator - (cpx a, cpx b){return cpx(a.r-b.r, a.i-b.i);} cpx operator * (cpx a, cpx b){return cpx(a.r*b.r-a.i*b.i, a.r*b.i+a.i*b.r);} int n, x, mx; int xs[maxn]; inline void read(int &k) { int f=1; k=0; char c=getchar(); while(c<'0' || c>'9') c=='-'&&(f=-1), c=getchar(); while(c<='9' && c>='0') k=k*10+c-'0', c=getchar(); k*=f; } inline void fft(int N, cpx a[], int f) { for(int i=0, j=0;i<N;i++) { if(i<j) swap(a[i], a[j]); for(int k=N>>1;(j^=k)<k;k>>=1); } for(int i=2;i<=N;i<<=1) for(int j=0, m=i>>1;j<N;j+=i) { cpx nw=cpx(cos(2.*pi/i), sin(2.*pi/i)*f), w=cpx(1, 0); for(int k=0;k<m;k++) { cpx t=a[j+k+m]*w; a[j+k+m]=a[j+k]-t; a[j+k]=a[j+k]+t; w=w*nw; } } } int main() { read(n); for(mx=1;mx<(n<<1);mx<<=1); for(int i=0;i<n;i++) read(x), a[n-i-1]=cpx(x), read(x), b[i]=cpx(x); fft(mx, a, 1); fft(mx, b, 1); for(int i=0;i<mx;i++) c[i]=a[i]*b[i]; fft(mx, c, -1); for(int i=0;i<n;i++) xs[i]=(int)(c[i].r/mx+0.5); for(int i=n-1;~i;i--) printf("%d\n", xs[i]); }
$E_i=\sum_{j=0}^{i-1}\frac{q_j}{(i-j)^2}-\sum_{j=i+1}^{n-1}\frac{q_j}{(j-i)^2}$
$A_i=q_i$
$B_i=\frac{1}{i^2}$
$E_i=\sum_{j=0}^{i-1}A_j*B_{i-j}-\sum_{j=i+1}^{n-1}A_j*B_{j-i}$
注意$B_0$取不到，需要設為0，然后$j=i$就可以取了，左邊直接FFT。
右半部分不是卷積的形式，需要轉化：
$\sum_{j=i+1}^{n-1}A_j*B_{j-i}$
$\sum_{j=1}^{n-i-1}A_{i+j}*B_{j}$
$\sum_{j=1}^{n-i-1}A’_{n-j-i-1}*B_{j}=Ans2_{n-i-1}$
$E_i=Ans1_i-Ans2_{n-i-1}$
#include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> #include<algorithm> #define ll long long using namespace std; const int maxn=500010, inf=1e9; const double pi=acos(-1); struct cpx{double r, i; cpx(double _r=0, double _i=0):r(_r),i(_i){};} a1[maxn], a2[maxn], b[maxn], ans1[maxn], ans2[maxn]; cpx operator + (cpx a, cpx b){return cpx(a.r+b.r, a.i+b.i);} cpx operator - (cpx a, cpx b){return cpx(a.r-b.r, a.i-b.i);} cpx operator * (cpx a, cpx b){return cpx(a.r*b.r-a.i*b.i, a.r*b.i+a.i*b.r);} int n, mx; double x; inline void fft(int N, cpx a[], int f) { for(int i=0, j=0;i<N;i++) { if(i<j) swap(a[i], a[j]); for(int k=N>>1;(j^=k)<k;k>>=1); } for(int i=2;i<=N;i<<=1) for(int j=0, m=i>>1;j<N;j+=i) { cpx nw=cpx(cos(2.*pi/i), sin(2.*pi/i)*f), w=cpx(1, 0); for(int k=0;k<m;k++) { cpx t=a[j+k+m]*w; a[j+k+m]=a[j+k]-t; a[j+k]=a[j+k]+t; w=w*nw; } } } int main() { scanf("%d", &n); for(mx=1;mx<(n<<1);mx<<=1); for(int i=0;i<n;i++) scanf("%lf", &x), a1[i]=cpx(x), a2[n-i-1]=cpx(x); for(int i=1;i<n;i++) b[i]=cpx(1.0/i/i); fft(mx, a1, 1); fft(mx, a2, 1); fft(mx, b, 1); for(int i=0;i<mx;i++) ans1[i]=a1[i]*b[i]; for(int i=0;i<mx;i++) ans2[i]=a2[i]*b[i]; fft(mx, ans1, -1); fft(mx, ans2, -1); for(int i=0;i<n;i++) printf("%.3lf\n", ans1[i].r/mx-ans2[n-i-1].r/mx); }
好厲害的題，完全沒看出是FFT
不連續回文子序列=回文子序列-連續回文子序列
后者顯然可以用manacher求得
對於回文子序列，只要知道了一個位置左右兩邊有多少對關於此位置對稱的相同的字符$f_i$，用$2^{f_i}-1$就能得到答案了。
求$f_i$我們同樣可以使用manacher的思想，先在原串$S$中插入一下分隔符，如$aba$改為$\$\# a\# b\# a\# $得到新串$S'$。
考慮兩個相等的字符對哪個位置有貢獻，顯然這三個位置組成了一個等差數列，故$x+y=2mid$，因為$x$和$y$是字符，且字符在$S'$中下標為偶數，所以有$\frac{x}{2}+\frac{y}{2}=mid$，而$\frac{x}{2}$在原串剛好是字符位置$+1$的位置。那么就能知道，對於$S_x=S_y$，他們對$x+y+2$有貢獻。
對$a$和$b$的貢獻分別算，算$a$時把字符是$a$的位置貢獻$gx_i$設為1，$b$的位置的貢獻$gx_i$設為0，算$b$貢獻時則反之。
那么我們就能得到一個求$f_i$的公式：
可以發現分子是一個卷積的形式，直接上FFT就好了。
#include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> #include<algorithm> #define MOD(x) ((x)>=mod?(x)-mod:(x)) #define ll long long using namespace std; const int maxn=500010, inf=1e9, mod=1e9+7; const double pi=acos(-1); struct cpx{double r, i;cpx(double _r=0, double _i=0):r(_r), i(_i){};} a[maxn], b[maxn]; cpx operator + (cpx a, cpx b){return cpx(a.r+b.r, a.i+b.i);} cpx operator - (cpx a, cpx b){return cpx(a.r-b.r, a.i-b.i);} cpx operator * (cpx a, cpx b){return cpx(a.r*b.r-a.i*b.i, a.r*b.i+a.i*b.r);} int n, N, ans, mid, mx, mxx; int p[maxn], mi[maxn], f[maxn]; char ch[maxn], s[maxn]; inline void fft(int N, cpx a[], int f) { for(int i=0, j=0;i<N;i++) { if(i<j) swap(a[i], a[j]); for(int k=N>>1;(j^=k)<k;k>>=1); } for(int i=2;i<=N;i<<=1) for(int j=0, m=i>>1;j<N;j+=i) { cpx nw=cpx(cos(2.*pi/i), sin(2.*pi/i)*f), w=cpx(1, 0); for(int k=0;k<m;k++) { cpx t=a[j+k+m]*w; a[j+k+m]=a[j+k]-t; a[j+k]=a[j+k]+t; w=w*nw; } } if(f==-1) for(int i=0;i<N;i++) a[i].r/=N; } inline void getans(int ty) { for(int i=0;i<n;i++) a[i]=cpx(ch[i]-'a'==ty, 0); for(int i=n;i<mx;i++) a[i]=cpx(0, 0); fft(mx, a, 1); for(int i=0;i<mx;i++) b[i]=a[i]*a[i]; fft(mx, b, -1); for(int i=2;i<N;i++) f[i]+=(int)(b[i-2].r+0.5); } int main() { scanf("%s", ch); n=strlen(ch); for(mx=1;mx<(n<<1);mx<<=1); for(int i=1;i<=n;i++) s[i<<1]=ch[i-1], s[(i<<1)|1]='#'; N=(n<<1)+2; s[1]=s[N]='#'; ans=mid=mxx=0; for(int i=1;i<=N;i++) { if(mxx>i) p[i]=min(p[(mid<<1)-i], mxx-i); else p[i]=1; while(s[i-p[i]]==s[i+p[i]]) p[i]++; if(p[i]+i>mxx) mxx=p[i]+i, mid=i; ans=MOD(ans+mod-(p[i]>>1)); } mi[0]=1; for(int i=1;i<N;i++) mi[i]=1ll*(mi[i-1]<<1)%mod; getans(0); getans(1); for(int i=2;i<N;i++) ans=MOD(ans+mi[(f[i]+1)>>1]-1); printf("%d\n", ans); }
母函數+容斥用FFT優化...
對於每一個物品$i$，使$A[a_i]$++$,B[2*a_i]$++$,C[3*a_i]$++，即$A,B,C$三個多項式代表選$1$個，$2$個，$3$個同樣物品的方案。
對於選$1$個物品的方案，顯然就是$A$了。
對於選$2$個物品的方案，$A^2$會多算$(x_i, x_i)$的情況，於是需要減去$B$，故方案為$\frac{A^2-B}{2}$。
對於選$3$個物品的方案，$A^3$會多算$(x_i, x_i, x_i),(x_i, x_i, y_i),(x_i, y_i, x_i), (y_i, x_i, x_i)$的情況，后三者的方案數剛好是$A*B*3$，但是$A*B*3$還包括$3$個$(x_i, x_i, x_i)$，而我們只需要減去一個，所以要多加上$2*C$，故方案為$\frac{A^3-A*B*3+2*C}{6}$。
#include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> #include<algorithm> #define ll long long using namespace std; const int maxn=500010, inf=1e9; const double pi=acos(-1); struct cpx{double r, i; cpx(double _r=0, double _i=0):r(_r), i(_i){};} a[maxn], b[maxn], c[maxn], a2[maxn], a3[maxn], ab[maxn]; cpx operator + (cpx a, cpx b){return cpx(a.r+b.r, a.i+b.i);} cpx operator - (cpx a, cpx b){return cpx(a.r-b.r, a.i-b.i);} cpx operator * (cpx a, cpx b){return cpx(a.r*b.r-a.i*b.i, a.r*b.i+a.i*b.r);} int n, mx, mxx, x; int ans[maxn]; inline void read(int &k) { int f=1; k=0; char c=getchar(); while(c<'0' || c>'9') c=='-'&&(f=-1), c=getchar(); while(c<='9' && c>='0') k=k*10+c-'0', c=getchar(); k*=f; } inline void fft(int N, cpx a[], int f) { for(int i=0, j=0;i<N;i++) { if(i<j) swap(a[i], a[j]); for(int k=N>>1;(j^=k)<k;k>>=1); } for(int i=2;i<=N;i<<=1) for(int j=0, m=i>>1;j<N;j+=i) { cpx nw=cpx(cos(2.*pi/i), sin(2.*pi/i)*f), w=cpx(1, 0); for(int k=0;k<m;k++) { cpx t=a[j+k+m]*w; a[j+k+m]=a[j+k]-t; a[j+k]=a[j+k]+t; w=w*nw; } } if(f==-1) for(int i=0;i<N;i++) a[i].r/=N; } int main() { read(n); for(int i=0;i<n;i++) read(x), ans[x]++, a[x].r++, b[x<<1].r++, c[x*3].r++, mxx=max(mxx, 3*x); for(mx=1;mx<mxx;mx<<=1); fft(mx, a, 1); fft(mx, b, 1); fft(mx, c, 1); for(int i=0;i<mx;i++) a[i]=(a[i]*a[i]*a[i]-cpx(3)*a[i]*b[i]+cpx(2)*c[i])*cpx(1.0/6)+(a[i]*a[i]-b[i])*cpx(1.0/2)+a[i]; fft(mx, a, -1); for(int i=1;i<mx;i++) if((int)(a[i].r+0.5))printf("%d %d\n", i, (int)(a[i].r+0.5)); }
設旋轉了$j$位，亮度增加了$c$，則差異值為$\sum_{i=0}^{n-1}(x_{i+j}-y_i+c)^2$
把式子拆開之后可以發現除了$-x_{i+j}*y_i$之外都是常數項和一次項，都可以直接預處理。
對$x$倍增一次，然后把$x$倒過來就可以對每種顏色求$-x_{i+j}*y_i$了。
$\sum_{i=0}^{n-j-1}x_{i+j}*y_i$
$\sum_{i=0}^{n-j-1}x_{n-i-j-1}*y_i$
然后枚舉旋轉次數和顏色即可。
233
#include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> #include<algorithm> #define ll long long using namespace std; const int maxn=500010; const ll inf=1e18; const double pi=acos(-1); struct cpx{double r, i;cpx(double _r=0, double _i=0):r(_r), i(_i){};} a[maxn], b[maxn], c[maxn]; cpx operator + (cpx a, cpx b){return cpx(a.r+b.r, a.i+b.i);} cpx operator - (cpx a, cpx b){return cpx(a.r-b.r, a.i-b.i);} cpx operator * (cpx a, cpx b){return cpx(a.r*b.r-a.i*b.i, a.r*b.i+a.i*b.r);} int n, m, x, mx; ll sqrx, sumx, sqry, sumy; ll xy[maxn], ans; inline void read(int &k) { int f=1; k=0; char c=getchar(); while(c<'0' || c>'9') c=='-'&&(f=-1), c=getchar(); while(c<='9' && c>='0') k=k*10+c-'0', c=getchar(); k*=f; } inline void fft(int N, cpx a[], int f) { for(int i=0, j=0;i<N;i++) { if(i<j) swap(a[i], a[j]); for(int k=N>>1;(j^=k)<k;k>>=1); } for(int i=2;i<=N;i<<=1) for(int j=0, m=i>>1;j<N;j+=i) { cpx nw=cpx(cos(2.*pi/i), sin(2.*pi/i)*f), w=cpx(1, 0); for(int k=0;k<m;k++) { cpx t=a[j+k+m]*w; a[j+k+m]=a[j+k]-t; a[j+k]=a[j+k]+t; w=w*nw; } } if(f==-1) for(int i=0;i<N;i++) a[i].r/=N; } int main() { read(n); read(m); for(int i=0;i<n;i++) read(x), sqrx+=1ll*x*x, sumx+=x, a[i]=a[i+n]=cpx(x); for(int i=0;i<n;i++) read(x), sqry+=1ll*x*x, sumy+=x, b[i]=cpx(x); n<<=1; for(mx=1;mx<(n<<1);mx<<=1); reverse(a, a+n); fft(mx, a, 1); fft(mx, b, 1); for(int i=0;i<mx;i++) c[i]=a[i]*b[i]; fft(mx, c, -1); for(int i=0;i<n;i++) xy[i]=(ll)(c[n-i-1].r+0.5); n>>=1; ans=inf; for(int i=0;i<n;i++) for(int j=-m;j<=m;j++) ans=min(ans, sqrx+sqry+1ll*j*j*n+2ll*(-xy[i]+sumx*j-sumy*j)); printf("%lld\n", ans); }
例7：codeforces 528D. Fuzzy Search
字符串匹配問題使用FFT的話，應該對於每一種字符分別統計。把母串對於這種字符能夠匹配的位置$i$，令$A_i=1$，否則為$0$，對模式串該字符位置$j$，令$B_j=1$，否則為$0$。
那么求每個位置能否匹配，等價於求：
$\sum_{l=0}^{n-1}\sum_{i=0}^{l+m-1}A_{l+i}\cdot B_{i}$
$\sum_{l=0}^{n-1}\sum_{i=0}^{l+m-1}A_{l+i}\cdot B_{m-i-1}$
卷積形式直接上FFT即可。
雖然$i$超過$m-1$時候$B$會越界，但是FFT處理的時候實際上不會越界，完全可以放心搞...
對於這題預處理出母串每個位置能匹配哪些字符即可然后直接上這種做法最后判斷每個位置能匹配的次數是否為模式串長度即可。
#include<iostream> #include<cstring> #include<cstdlib> #include<cstdio> #include<cmath> #include<algorithm> #define ll long long using namespace std; const int maxn=1000010, inf=1e9; const double pi=acos(-1); struct cpx{double r, i;cpx(double _r=0, double _i=0):r(_r), i(_i){};} a[maxn], b[maxn], c[maxn]; cpx operator + (cpx a, cpx b){return cpx(a.r+b.r, a.i+b.i);} cpx operator - (cpx a, cpx b){return cpx(a.r-b.r, a.i-b.i);} cpx operator * (cpx a, cpx b){return cpx(a.r*b.r-a.i*b.i, a.i*b.r+a.r*b.i);} int n, m, k, N, ANS; int ok[maxn][4], cnt[4], ans[maxn]; char s[maxn], t[maxn]; inline void read(int &k) { int f=1; k=0; char c=getchar(); while(c<'0' || c>'9') c=='-'&&(f=-1), c=getchar(); while(c<='9' && c>='0') k=k*10+c-'0', c=getchar(); k*=f; } inline int col(char ch) { if(ch=='A') return 0; if(ch=='T') return 1; if(ch=='G') return 2; return 3; } inline void fft(int N, cpx a[], int f) { for(int i=0, j=0;i<N;i++) { if(i<j) swap(a[i], a[j]); for(int k=N>>1;(j^=k)<k;k>>=1); } for(int i=2;i<=N;i<<=1) for(int j=0, m=i>>1;j<N;j+=i) { cpx nw=cpx(cos(2.*pi/i), sin(2.*pi/i)*f), w=cpx(1, 0); for(int k=0;k<m;k++) { cpx t=a[j+k+m]*w; a[j+k+m]=a[j+k]-t; a[j+k]=a[j+k]+t; w=w*nw; } } if(f==-1) for(int i=0;i<N;i++) a[i].r/=N; } int main() { read(n); read(m); read(k); scanf("%s", s); scanf("%s", t); for(int i=0;i<n;i++) { cnt[col(s[i])]++; if(i-k-1>=0) cnt[col(s[i-k-1])]--; for(int j=0;j<4;j++) if(cnt[j]) ok[i][j]=1; } cnt[0]=cnt[1]=cnt[2]=cnt[3]=0; for(int i=n-1;~i;i--) { cnt[col(s[i])]++; if(i+k+1<n) cnt[col(s[i+k+1])]--; for(int j=0;j<4;j++) if(cnt[j]) ok[i][j]=1; } for(N=1;N<(n+m);N<<=1); for(int i=0;i<4;i++) { memset(a, 0, sizeof(a)); memset(b, 0, sizeof(b)); for(int j=0;j<n;j++) a[j]=cpx(ok[j][i], 0); for(int j=0;j<m;j++) b[m-j-1]=cpx(col(t[j])==i, 0); fft(N, a, 1); fft(N, b, 1); for(int j=0;j<N;j++) c[j]=a[j]*b[j]; fft(N, c, -1); for(int j=0;j<n;j++) ans[j]+=(int)(c[j+m-1].r+0.5); } for(int i=0;i<n;i++) if(ans[i]==m) ANS++; printf("%d\n", ANS); }
例8：