-
Notifications
You must be signed in to change notification settings - Fork 23
Expand file tree
/
Copy pathfft.cpp
More file actions
51 lines (41 loc) · 1.3 KB
/
fft.cpp
File metadata and controls
51 lines (41 loc) · 1.3 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
typedef complex<double> base;
void fft (vector<base> & a, bool invert) {
int n = (int) a.size();
if (n == 1) return;
vector<base> a0 (n/2), a1 (n/2);
for (int i=0, j=0; i<n; i+=2, ++j) {
a0[j] = a[i];
a1[j] = a[i+1];
}
fft (a0, invert);
fft (a1, invert);
double ang = 2*PI/n * (invert ? -1 : 1);
base w (1), wn (cos(ang), sin(ang));
for (int i=0; i<n/2; ++i) {
a[i] = a0[i] + w * a1[i];
a[i+n/2] = a0[i] - w * a1[i];
if (invert)
a[i] /= 2, a[i+n/2] /= 2;
w *= wn;
}
}
/*
1. Build two different vectors with coefficient of the equation ( 3 + x + 2x^2 will be [3,1,2] )
2. Send them to multiply ( arr, brr, res )
3. res vector contains result of multiplication.
WARNING: FFT causes precision error. multiply() rounds them to integer.
*/
void multiply (const vector<int> & a, const vector<int> & b, vector<int> & res) {
vector<base> fa (a.begin(), a.end()), fb (b.begin(), b.end());
size_t n = 1;
while (n < max (a.size(), b.size())) n <<= 1;
n <<= 1;
fa.resize (n), fb.resize (n);
fft (fa, false), fft (fb, false);
for (size_t i=0; i<n; ++i)
fa[i] *= fb[i];
fft (fa, true);
res.resize (n);
for (size_t i=0; i<n; ++i)
res[i] = int (fa[i].real() + 0.5);
}