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