Botan  2.1.0
Crypto and TLS for C++11
mp_karat.cpp
Go to the documentation of this file.
1 /*
2 * Multiplication and Squaring
3 * (C) 1999-2010 Jack Lloyd
4 * 2016 Matthias Gierlings
5 *
6 * Botan is released under the Simplified BSD License (see license.txt)
7 */
8 
9 #include <botan/internal/mp_core.h>
10 #include <botan/internal/mp_asmi.h>
11 #include <botan/mem_ops.h>
12 
13 namespace Botan {
14 
15 namespace {
16 
17 const size_t KARATSUBA_MULTIPLY_THRESHOLD = 32;
18 const size_t KARATSUBA_SQUARE_THRESHOLD = 32;
19 
20 /*
21 * Simple O(N^2) Multiplication
22 */
23 void basecase_mul(word z[],
24  const word x[], size_t x_size,
25  const word y[], size_t y_size)
26  {
27  const size_t x_size_8 = x_size - (x_size % 8);
28 
29  clear_mem(z, x_size + y_size);
30 
31  for(size_t i = 0; i != y_size; ++i)
32  {
33  const word y_i = y[i];
34 
35  word carry = 0;
36 
37  for(size_t j = 0; j != x_size_8; j += 8)
38  carry = word8_madd3(z + i + j, x + j, y_i, carry);
39 
40  for(size_t j = x_size_8; j != x_size; ++j)
41  z[i+j] = word_madd3(x[j], y_i, z[i+j], &carry);
42 
43  z[x_size+i] = carry;
44  }
45  }
46 
47 /*
48 * Karatsuba Multiplication Operation
49 */
50 void karatsuba_mul(word z[], const word x[], const word y[], size_t N,
51  word workspace[])
52  {
53  if(N < KARATSUBA_MULTIPLY_THRESHOLD || N % 2)
54  {
55  if(N == 6)
56  return bigint_comba_mul6(z, x, y);
57  else if(N == 8)
58  return bigint_comba_mul8(z, x, y);
59  else if(N == 16)
60  return bigint_comba_mul16(z, x, y);
61  else
62  return basecase_mul(z, x, N, y, N);
63  }
64 
65  const size_t N2 = N / 2;
66 
67  const word* x0 = x;
68  const word* x1 = x + N2;
69  const word* y0 = y;
70  const word* y1 = y + N2;
71  word* z0 = z;
72  word* z1 = z + N;
73 
74  const int32_t cmp0 = bigint_cmp(x0, N2, x1, N2);
75  const int32_t cmp1 = bigint_cmp(y1, N2, y0, N2);
76 
77  clear_mem(workspace, 2*N);
78 
79  /*
80  * If either of cmp0 or cmp1 is zero then z0 or z1 resp is zero here,
81  * resulting in a no-op - z0*z1 will be equal to zero so we don't need to do
82  * anything, clear_mem above already set the correct result.
83  *
84  * However we ignore the result of the comparisons and always perform the
85  * subtractions and recursively multiply to avoid the timing channel.
86  */
87 
88  //if(cmp0 && cmp1)
89  {
90  if(cmp0 > 0)
91  bigint_sub3(z0, x0, N2, x1, N2);
92  else
93  bigint_sub3(z0, x1, N2, x0, N2);
94 
95  if(cmp1 > 0)
96  bigint_sub3(z1, y1, N2, y0, N2);
97  else
98  bigint_sub3(z1, y0, N2, y1, N2);
99 
100  karatsuba_mul(workspace, z0, z1, N2, workspace+N);
101  }
102 
103  karatsuba_mul(z0, x0, y0, N2, workspace+N);
104  karatsuba_mul(z1, x1, y1, N2, workspace+N);
105 
106  const word ws_carry = bigint_add3_nc(workspace + N, z0, N, z1, N);
107  word z_carry = bigint_add2_nc(z + N2, N, workspace + N, N);
108 
109  z_carry += bigint_add2_nc(z + N + N2, N2, &ws_carry, 1);
110  bigint_add2_nc(z + N + N2, N2, &z_carry, 1);
111 
112  if((cmp0 == cmp1) || (cmp0 == 0) || (cmp1 == 0))
113  bigint_add2(z + N2, 2*N-N2, workspace, N);
114  else
115  bigint_sub2(z + N2, 2*N-N2, workspace, N);
116  }
117 
118 /*
119 * Karatsuba Squaring Operation
120 */
121 void karatsuba_sqr(word z[], const word x[], size_t N, word workspace[])
122  {
123  if(N < KARATSUBA_SQUARE_THRESHOLD || N % 2)
124  {
125  if(N == 6)
126  return bigint_comba_sqr6(z, x);
127  else if(N == 8)
128  return bigint_comba_sqr8(z, x);
129  else if(N == 16)
130  return bigint_comba_sqr16(z, x);
131  else
132  return basecase_mul(z, x, N, x, N);
133  }
134 
135  const size_t N2 = N / 2;
136 
137  const word* x0 = x;
138  const word* x1 = x + N2;
139  word* z0 = z;
140  word* z1 = z + N;
141 
142  const int32_t cmp = bigint_cmp(x0, N2, x1, N2);
143 
144  clear_mem(workspace, 2*N);
145 
146  // See comment in karatsuba_mul
147 
148  //if(cmp)
149  {
150  if(cmp > 0)
151  bigint_sub3(z0, x0, N2, x1, N2);
152  else
153  bigint_sub3(z0, x1, N2, x0, N2);
154 
155  karatsuba_sqr(workspace, z0, N2, workspace+N);
156  }
157 
158  karatsuba_sqr(z0, x0, N2, workspace+N);
159  karatsuba_sqr(z1, x1, N2, workspace+N);
160 
161  const word ws_carry = bigint_add3_nc(workspace + N, z0, N, z1, N);
162  word z_carry = bigint_add2_nc(z + N2, N, workspace + N, N);
163 
164  z_carry += bigint_add2_nc(z + N + N2, N2, &ws_carry, 1);
165  bigint_add2_nc(z + N + N2, N2, &z_carry, 1);
166 
167  /*
168  * This is only actually required if cmp is != 0, however
169  * if cmp==0 then workspace[0:N] == 0 and avoiding the jump
170  * hides a timing channel.
171  */
172  bigint_sub2(z + N2, 2*N-N2, workspace, N);
173  }
174 
175 /*
176 * Pick a good size for the Karatsuba multiply
177 */
178 size_t karatsuba_size(size_t z_size,
179  size_t x_size, size_t x_sw,
180  size_t y_size, size_t y_sw)
181  {
182  if(x_sw > x_size || x_sw > y_size || y_sw > x_size || y_sw > y_size)
183  return 0;
184 
185  if(((x_size == x_sw) && (x_size % 2)) ||
186  ((y_size == y_sw) && (y_size % 2)))
187  return 0;
188 
189  const size_t start = (x_sw > y_sw) ? x_sw : y_sw;
190  const size_t end = (x_size < y_size) ? x_size : y_size;
191 
192  if(start == end)
193  {
194  if(start % 2)
195  return 0;
196  return start;
197  }
198 
199  for(size_t j = start; j <= end; ++j)
200  {
201  if(j % 2)
202  continue;
203 
204  if(2*j > z_size)
205  return 0;
206 
207  if(x_sw <= j && j <= x_size && y_sw <= j && j <= y_size)
208  {
209  if(j % 4 == 2 &&
210  (j+2) <= x_size && (j+2) <= y_size && 2*(j+2) <= z_size)
211  return j+2;
212  return j;
213  }
214  }
215 
216  return 0;
217  }
218 
219 /*
220 * Pick a good size for the Karatsuba squaring
221 */
222 size_t karatsuba_size(size_t z_size, size_t x_size, size_t x_sw)
223  {
224  if(x_sw == x_size)
225  {
226  if(x_sw % 2)
227  return 0;
228  return x_sw;
229  }
230 
231  for(size_t j = x_sw; j <= x_size; ++j)
232  {
233  if(j % 2)
234  continue;
235 
236  if(2*j > z_size)
237  return 0;
238 
239  if(j % 4 == 2 && (j+2) <= x_size && 2*(j+2) <= z_size)
240  return j+2;
241  return j;
242  }
243 
244  return 0;
245  }
246 
247 }
248 
249 /*
250 * Multiplication Algorithm Dispatcher
251 */
252 void bigint_mul(BigInt& z, const BigInt& x, const BigInt& y, word workspace[])
253  {
254  const size_t x_sig_words = x.sig_words();
255  const size_t y_sig_words = y.sig_words();
256 
257  clear_mem(z.mutable_data(), z.size());
258 
259  if(x_sig_words == 1)
260  {
261  bigint_linmul3(z.mutable_data(), y.data(), y_sig_words, x.data()[0]);
262  }
263  else if(y_sig_words == 1)
264  {
265  bigint_linmul3(z.mutable_data(), x.data(), x_sig_words, y.data()[0]);
266  }
267  else if(x_sig_words <= 4 && x.size() >= 4 &&
268  y_sig_words <= 4 && y.size() >= 4 && z.size() >= 8)
269  {
270  bigint_comba_mul4(z.mutable_data(), x.data(), y.data());
271  }
272  else if(x_sig_words <= 6 && x.size() >= 6 &&
273  y_sig_words <= 6 && y.size() >= 6 && z.size() >= 12)
274  {
275  bigint_comba_mul6(z.mutable_data(), x.data(), y.data());
276  }
277  else if(x_sig_words <= 8 && x.size() >= 8 &&
278  y_sig_words <= 8 && y.size() >= 8 && z.size() >= 16)
279  {
280  bigint_comba_mul8(z.mutable_data(), x.data(), y.data());
281  }
282  else if(x_sig_words <= 9 && x.size() >= 9 &&
283  y_sig_words <= 9 && y.size() >= 9 && z.size() >= 18)
284  {
285  bigint_comba_mul9(z.mutable_data(), x.data(), y.data());
286  }
287  else if(x_sig_words <= 16 && x.size() >= 16 &&
288  y_sig_words <= 16 && y.size() >= 16 && z.size() >= 32)
289  {
290  bigint_comba_mul16(z.mutable_data(), x.data(), y.data());
291  }
292  else if(x_sig_words < KARATSUBA_MULTIPLY_THRESHOLD ||
293  y_sig_words < KARATSUBA_MULTIPLY_THRESHOLD ||
294  !workspace)
295  {
296  basecase_mul(z.mutable_data(), x.data(), x_sig_words, y.data(), y_sig_words);
297  }
298  else
299  {
300  const size_t N = karatsuba_size(z.size(), x.size(), x_sig_words, y.size(), y_sig_words);
301 
302  if(N)
303  karatsuba_mul(z.mutable_data(), x.data(), y.data(), N, workspace);
304  else
305  basecase_mul(z.mutable_data(), x.data(), x_sig_words, y.data(), y_sig_words);
306  }
307  }
308 
309 /*
310 * Squaring Algorithm Dispatcher
311 */
312 void bigint_sqr(word z[], size_t z_size, word workspace[],
313  const word x[], size_t x_size, size_t x_sw)
314  {
315  BOTAN_ASSERT(z_size/2 >= x_sw, "Output size is sufficient");
316 
317  if(x_sw == 1)
318  {
319  bigint_linmul3(z, x, x_sw, x[0]);
320  }
321  else if(x_sw <= 4 && x_size >= 4 && z_size >= 8)
322  {
323  bigint_comba_sqr4(z, x);
324  }
325  else if(x_sw <= 6 && x_size >= 6 && z_size >= 12)
326  {
327  bigint_comba_sqr6(z, x);
328  }
329  else if(x_sw <= 8 && x_size >= 8 && z_size >= 16)
330  {
331  bigint_comba_sqr8(z, x);
332  }
333  else if(x_sw == 9 && x_size >= 9 && z_size >= 18)
334  {
335  bigint_comba_sqr9(z, x);
336  }
337  else if(x_sw <= 16 && x_size >= 16 && z_size >= 32)
338  {
339  bigint_comba_sqr16(z, x);
340  }
341  else if(x_size < KARATSUBA_SQUARE_THRESHOLD || !workspace)
342  {
343  basecase_mul(z, x, x_sw, x, x_sw);
344  }
345  else
346  {
347  const size_t N = karatsuba_size(z_size, x_size, x_sw);
348 
349  if(N)
350  karatsuba_sqr(z, x, N, workspace);
351  else
352  basecase_mul(z, x, x_sw, x, x_sw);
353  }
354  }
355 
356 }
size_t sig_words() const
Definition: bigint.h:393
int32_t bigint_cmp(const word x[], size_t x_size, const word y[], size_t y_size)
Definition: mp_core.cpp:378
void clear_mem(T *ptr, size_t n)
Definition: mem_ops.h:57
word bigint_sub2(word x[], size_t x_size, const word y[], size_t y_size)
Definition: mp_core.cpp:157
word bigint_add3_nc(word z[], const word x[], size_t x_size, const word y[], size_t y_size)
Definition: mp_core.cpp:113
word word_madd3(word a, word b, word c, word *d)
Definition: mp_madd.h:105
void bigint_comba_mul4(word z[8], const word x[4], const word y[4])
Definition: mp_comba.cpp:50
word * mutable_data()
Definition: bigint.h:419
word bigint_sub3(word z[], const word x[], size_t x_size, const word y[], size_t y_size)
Definition: mp_core.cpp:198
size_t size() const
Definition: bigint.h:387
void bigint_comba_mul9(word z[18], const word x[9], const word y[9])
Definition: mp_comba.cpp:474
#define BOTAN_ASSERT(expr, assertion_made)
Definition: assert.h:27
void bigint_linmul3(word z[], const word x[], size_t x_size, word y)
Definition: mp_core.cpp:240
void bigint_comba_sqr16(word z[32], const word x[16])
Definition: mp_comba.cpp:598
word word8_madd3(word z[8], const word x[8], word y, word carry)
Definition: mp_asmi.h:681
void bigint_sqr(word z[], size_t z_size, word workspace[], const word x[], size_t x_size, size_t x_sw)
Definition: mp_karat.cpp:312
const word * data() const
Definition: bigint.h:425
word bigint_add2_nc(word x[], size_t x_size, const word y[], size_t y_size)
Definition: mp_core.cpp:90
Definition: alg_id.cpp:13
void bigint_comba_sqr9(word z[18], const word x[9])
Definition: mp_comba.cpp:386
void bigint_comba_mul8(word z[16], const word x[8], const word y[8])
Definition: mp_comba.cpp:283
void bigint_comba_sqr8(word z[16], const word x[8])
Definition: mp_comba.cpp:208
void bigint_add2(word x[], size_t x_size, const word y[], size_t y_size)
Definition: mp_core.cpp:138
void bigint_comba_mul16(word z[32], const word x[16], const word y[16])
Definition: mp_comba.cpp:805
void bigint_mul(BigInt &z, const BigInt &x, const BigInt &y, word workspace[])
Definition: mp_karat.cpp:252
void bigint_comba_mul6(word z[12], const word x[6], const word y[6])
Definition: mp_comba.cpp:141
void bigint_comba_sqr4(word z[8], const word x[4])
Definition: mp_comba.cpp:17
void bigint_comba_sqr6(word z[12], const word x[6])
Definition: mp_comba.cpp:89