快速傅里叶变换 (fast Fourier transform)
xn={x0,x1,…xn-1} (num:N)
旋转因子系数:
d=2pik/N
旋转因子
wk(n)=(cos(dn)+isin(dn)) n=[0,N-1]
y(k)= sum(x(n)wk(n),0,N-1)
y(k)={y(0),y(1),…y(N-1)} 傅里叶级数
x(n)wk(n)的级数是:
1.d=2pik/N 这个系数决定转动圈数,不懂看第二个再说
2.y(k)=x(0)wk(0)+x(1)wk(2)+…+x(n-1)wk(n-1)
旋转因子:cos(t)+isin(t)
t=2pikn/N
T=2pi就是一圈,那么可以得出步进量:2pi*k/N, 圈数:k
在一圈中:wk(n)的单位圆上t<=pi/2,是他的展开区 t<=pi是他的对称区
#include <bits/stdc++.h>
#define fu(x , y , z) for(int x = y ; x <= z ; x ++)
#define fd(x , y , z) for(int x = y ; x >= z ; x --)
#define LL long long
using namespace std;
const int N = 4e6 + 5;
const double pi = acos (-1.0);
struct node {
double x , y;
} a[N] , b[N];
int n , m , len = 1 , r[N] , l;
node operator + (node a, node b) { return (node){a.x + b.x , a.y + b.y};}
node operator - (node a, node b) { return (node){a.x - b.x , a.y - b.y};}
node operator * (node a, node b) { return (node){a.x * b.x - a.y * b.y , a.x * b.y + a.y * b.x};}
int read () {
int val = 0 , fu = 1;
char ch = getchar ();
while (ch < '0' || ch > '9') {
if (ch == '-') fu = -1;
ch = getchar ();
}
while (ch >= '0' && ch <= '9') {
val = val * 10 + (ch - '0');
ch = getchar ();
}
return val * fu;
}
void fft (node *A , int inv) {
for (int i = 0 ; i < len ; i ++)
if (i < r[i])
swap (A[i] , A[r[i]]);
for (int mid = 1 ; mid < len ; mid <<= 1) { //mid={1,2,4,8,16,32,64,...}
node wn = (node){cos (1.0 * pi / mid) , inv * sin (1.0 * pi / mid)};
for (int R = mid << 1 , j = 0 ; j < len ; j += R) {
node w = (node){1 , 0};
for (int k = 0 ; k < mid ; k ++ , w = w * wn) {
node x = A[j + k] , y = w * A[j + mid + k];
A[j + k] = x + y;
A[j + mid + k] = x - y;
}
}
}
}
int main () {
n = read () , m = read ();
fu (i , 0 , n) a[i].x = read ();
fu (i , 0 , m) b[i].x = read ();
while (len <= n + m) len <<= 1 , l ++;
for (int i = 0 ; i < len ; i ++)
r[i] = (r[i >> 1] >> 1) | ((i & 1) << (l - 1));
fft (a , 1);
fft (b , 1);
fu (i , 0 , len)
a[i] = a[i] * b[i];
fft (a , -1);
// for (int i = 0 ; i <= n + m ; i ++) cout << a[i].x << " " << a[i].y << "\n";
fu (i , 0 , n + m)
printf ("%d " , (int)(a[i].x / len + 0.5));
return 0;
}
计算2^n的逆排长度:2的6次方:64=0b100-0000;逆排长度0b11-1111
while (len <= n + m) len <<= 1 , l ++;
产生逆排:0,32,
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]]);
#include <bits/stdc++.h>
#define fu(x , y , z) for(int x = y ; x <= z ; x ++)
using namespace std;
const int N = 4e6 + 5;
const double pi = acos (-1.0);
int n , m1 , m2 , rev[N];
complex<double> a[N] , b[N];
void fft (complex<double> *a , int type) {
fu (i , 0 , n - 1)
if (i < rev[i]) swap (a[i] , a[rev[i]]);
for (int j = 1 ; j < n ; j <<= 1) {
complex<double> W(cos (pi / j) , sin (pi / j) * type);
for (int k = 0 ; k < n ; k += (j << 1)) {
complex<double> w(1.0 , 0.0);
fu (i , 0 , j - 1) {
complex<double> ye , yo;
ye = a[i + k] , yo = a[i + j + k] * w;
a[i + k] = ye + yo;
a[i + k] = ye + yo;
a[i + j + k] = ye - yo;
w *= W;
}
}
}
}
int main () {
scanf ("%d%d" , &m1 , &m2);
fu (i , 0 , m1) cin >> a[i];
fu (i , 0 , m2) cin >> b[i];
n = m1 + m2;
int t = 0;
while (n >= (1 << t))
t ++;
n = (1 << t);
fu (i , 0 , n - 1)
rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (n >> 1) : 0);
fft (a , 1) , fft (b , 1);
fu (i , 0 , n)
a[i] *= b[i];
fft (a , -1);
fu (i , 0 , m1 + m2)
printf ("%d " , (int)(a[i].real() / (double)n + 0.5));
return 0;
}
f0=(a0,a1,…an-2)
f1=(a1,a3,…an-1)
f(wnk)=f0(wnk/2)+wnkf1(wkn/2)
f(wn(k+n/2)=f0(wnk/2)-wnkf1(wkn/2)
#include<cstdio>
2 #include<iostream>
3 #include<cmath>
4 #include<cstring>
5 #include<algorithm>
6 #include<cstdlib>
7 using namespace std;
8 const int mod=1e9+7;
9 const double pi=acos(-1);
10 struct cn
11 {
12 double x,y;
13 cn (double x=0,double y=0):x(x),y(y) {}
14 }a[300005],b[300005],c[300005];
15 cn operator + (const cn &a,const cn &b) {return cn(a.x+b.x,a.y+b.y);}
16 cn operator - (const cn &a,const cn &b) {return cn(a.x-b.x,a.y-b.y);}
17 cn operator * (const cn &a,const cn &b) {return cn(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}
18 void fft(cn a[],int n,int l,int f)
19 {
20 int rev[n+5];
21 rev[0]=0;
22 for (int i=1; i<n; i++){
23 rev[i]=(rev[i>>1]>>1)|((i&1)<<l-1);
24 if (i<rev[i]) swap(a[i],a[rev[i]]);
25 }
26 for (int i=1; i<n; i<<=1){
27 cn wi(cos(pi/i),f*sin(pi/i));
28 for (int j=0; j<n; j+=i*2){
29 cn w(1,0);
30 for (int k=0; k<i; k++){
31 cn x=a[j+k],y=w*a[j+k+i];
32 a[j+k]=x+y;
33 a[j+k+i]=x-y;
34 w=w*wi;
35 }
36 }
37 }
38 if (f==-1)
39 for (int i=0; i<n; i++){
40 a[i].x/=n; a[i].y/=n;
41 }
42 }
43 int main()
44 {
45 int n,m;
46 scanf("%d%d",&n,&m); n++; m++;
47 for (int i=0; i<n; i++) scanf("%lf",&a[i].x);
48 for (int i=0; i<m; i++) scanf("%lf",&b[i].x);
49 int l=0,N=1;
50 while (N<n+m-1) N<<=1,l++;
51 fft(a,N,l,1);
52 fft(b,N,l,1);
53 for (int i=0; i<N; i++) c[i]=a[i]*b[i];
54 fft(c,N,l,-1);
55 for (int i=0; i<n+m-1; i++) printf("%d ",(int)(c[i].x+0.5));
56 return 0;
57 }