/
Karatsuba.cpp
75 lines (54 loc) · 1.81 KB
/
Karatsuba.cpp
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
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
// functon Karatsuba (and stupid as well) computes c += a * b, not c = a * b
using hvect = vector<modulo<>>::iterator;
using hcvect = vector<modulo<>>::const_iterator;
void add(hcvect abegin, hcvect aend, hvect ans)
{
for (auto it = abegin; it != aend; ++it, ++ans)
*ans += *it;
}
void sub(hcvect abegin, hcvect aend, hvect ans)
{
for (auto it = abegin; it != aend; ++it, ++ans)
*ans -= *it;
}
void stupid(int siz, hcvect abegin, hcvect bbegin, hvect ans)
{
for (int i = 0; i < siz; i++)
for (int j = 0; j < siz; j++)
*(ans + i + j) += *(abegin + i) * *(bbegin + j);
}
void Karatsuba(size_t siz, hcvect abegin, hcvect bbegin, hvect ans, hvect small, hvect big, hvect sum)
{
assert((siz & (siz - 1)) == 0);
if (siz <= 32)
{
stupid(siz, abegin, bbegin, ans);
return;
}
auto amid = abegin + siz / 2, aend = abegin + siz;
auto bmid = bbegin + siz / 2, bend = bbegin + siz;
auto smid = sum + siz / 2, send = sum + siz;
fill(small, small + siz, 0);
Karatsuba(siz / 2, abegin, bbegin, small, small + siz, big + siz, sum);
fill(big, big + siz, 0);
Karatsuba(siz / 2, amid, bmid, big, small + siz, big + siz, sum);
copy(abegin, amid, sum);
add(amid, aend, sum);
copy(bbegin, bmid, sum + siz / 2);
add(bmid, bend, sum + siz / 2);
Karatsuba(siz / 2, sum, smid, ans + siz / 2, small + siz, big + siz, send);
add(small, small + siz, ans);
sub(small, small + siz, ans + siz / 2);
add(big, big + siz, ans + siz);
sub(big, big + siz, ans + siz / 2);
}
void mult(vector<modulo<>> a, vector<modulo<>> b, vector<modulo<>> &c)
{
a.resize(up(max(a.size(), b.size())), 0);
b.resize(a.size(), 0);
c.resize(max(c.size(), a.size() * 2), 0);
vector<modulo<>> small(2 * a.size());
auto big = small;
auto sum = small;
Karatsuba(a.size(), a.begin(), b.begin(), c.begin(), small.begin(), big.begin(), sum.begin());
}