monotone

monotone Mtn Source Tree

Root/botan/mp_mul.cpp

1/*************************************************
2* Karatsuba Multiplication Source File *
3* (C) 1999-2007 The Botan Project *
4*************************************************/
5
6#include <botan/mp_core.h>
7#include <botan/mem_ops.h>
8
9namespace Botan {
10
11namespace {
12
13/*************************************************
14* Simple O(N^2) Multiplication *
15*************************************************/
16void bigint_simple_mul(word z[], const word x[], u32bit x_size,
17 const word y[], u32bit y_size)
18 {
19 clear_mem(z, x_size + y_size);
20
21 for(u32bit j = 0; j != x_size; ++j)
22 z[j+y_size] = bigint_mul_add_words(z + j, y, y_size, x[j]);
23 }
24
25/*************************************************
26* Karatsuba Multiplication Operation *
27*************************************************/
28void karatsuba_mul(word z[], const word x[], const word y[], u32bit N,
29 word workspace[])
30 {
31 const u32bit KARATSUBA_MUL_LOWER_SIZE = BOTAN_KARAT_MUL_THRESHOLD;
32
33 if(N == 6)
34 bigint_comba_mul6(z, x, y);
35 else if(N == 8)
36 bigint_comba_mul8(z, x, y);
37 else if(N < KARATSUBA_MUL_LOWER_SIZE || N % 2)
38 bigint_simple_mul(z, x, N, y, N);
39 else
40 {
41 const u32bit N2 = N / 2;
42
43 const word* x0 = x;
44 const word* x1 = x + N2;
45 const word* y0 = y;
46 const word* y1 = y + N2;
47 word* z0 = z;
48 word* z1 = z + N;
49
50 const s32bit cmp0 = bigint_cmp(x0, N2, x1, N2);
51 const s32bit cmp1 = bigint_cmp(y1, N2, y0, N2);
52
53 clear_mem(workspace, 2*N);
54
55 if(cmp0 && cmp1)
56 {
57 if(cmp0 > 0)
58 bigint_sub3(z0, x0, N2, x1, N2);
59 else
60 bigint_sub3(z0, x1, N2, x0, N2);
61
62 if(cmp1 > 0)
63 bigint_sub3(z1, y1, N2, y0, N2);
64 else
65 bigint_sub3(z1, y0, N2, y1, N2);
66
67 karatsuba_mul(workspace, z0, z1, N2, workspace+N);
68 }
69
70 karatsuba_mul(z0, x0, y0, N2, workspace+N);
71 karatsuba_mul(z1, x1, y1, N2, workspace+N);
72
73 word carry = bigint_add3_nc(workspace+N, z0, N, z1, N);
74 carry += bigint_add2_nc(z + N2, N, workspace + N, N);
75 bigint_add2_nc(z + N + N2, N2, &carry, 1);
76
77 if((cmp0 == cmp1) || (cmp0 == 0) || (cmp1 == 0))
78 bigint_add2(z + N2, 2*N-N2, workspace, N);
79 else
80 bigint_sub2(z + N2, 2*N-N2, workspace, N);
81 }
82 }
83
84/*************************************************
85* Pick a good size for the Karatsuba multiply *
86*************************************************/
87u32bit karatsuba_size(u32bit z_size,
88 u32bit x_size, u32bit x_sw,
89 u32bit y_size, u32bit y_sw)
90 {
91 if(x_sw > x_size || x_sw > y_size || y_sw > x_size || y_sw > y_size)
92 return 0;
93
94 if(((x_size == x_sw) && (x_size % 2)) ||
95 ((y_size == y_sw) && (y_size % 2)))
96 return 0;
97
98 const u32bit start = (x_sw > y_sw) ? x_sw : y_sw;
99 const u32bit end = (x_size < y_size) ? x_size : y_size;
100
101 if(start == end)
102 {
103 if(start % 2)
104 return 0;
105 return start;
106 }
107
108 for(u32bit j = start; j <= end; ++j)
109 {
110 if(j % 2)
111 continue;
112
113 if(2*j > z_size)
114 return 0;
115
116 if(x_sw <= j && j <= x_size && y_sw <= j && j <= y_size)
117 {
118 if(j % 4 == 2 &&
119 (j+2) <= x_size && (j+2) <= y_size && 2*(j+2) <= z_size)
120 return j+2;
121 return j;
122 }
123 }
124
125 return 0;
126 }
127
128/*************************************************
129* Handle small operand multiplies *
130*************************************************/
131void handle_small_mul(word z[], u32bit z_size,
132 const word x[], u32bit x_size, u32bit x_sw,
133 const word y[], u32bit y_size, u32bit y_sw)
134 {
135 if(x_sw == 1) bigint_linmul3(z, y, y_sw, x[0]);
136 else if(y_sw == 1) bigint_linmul3(z, x, x_sw, y[0]);
137
138 else if(x_sw <= 4 && x_size >= 4 &&
139 y_sw <= 4 && y_size >= 4 && z_size >= 8)
140 bigint_comba_mul4(z, x, y);
141
142 else if(x_sw <= 6 && x_size >= 6 &&
143 y_sw <= 6 && y_size >= 6 && z_size >= 12)
144 bigint_comba_mul6(z, x, y);
145
146 else if(x_sw <= 8 && x_size >= 8 &&
147 y_sw <= 8 && y_size >= 8 && z_size >= 16)
148 bigint_comba_mul8(z, x, y);
149
150 else
151 bigint_simple_mul(z, x, x_sw, y, y_sw);
152 }
153
154}
155
156/*************************************************
157* Multiplication Algorithm Dispatcher *
158*************************************************/
159void bigint_mul(word z[], u32bit z_size, word workspace[],
160 const word x[], u32bit x_size, u32bit x_sw,
161 const word y[], u32bit y_size, u32bit y_sw)
162 {
163 if(x_size <= 8 || y_size <= 8)
164 {
165 handle_small_mul(z, z_size, x, x_size, x_sw, y, y_size, y_sw);
166 return;
167 }
168
169 const u32bit N = karatsuba_size(z_size, x_size, x_sw, y_size, y_sw);
170
171 if(N)
172 {
173 clear_mem(workspace, 2*N);
174 karatsuba_mul(z, x, y, N, workspace);
175 }
176 else
177 bigint_simple_mul(z, x, x_sw, y, y_sw);
178 }
179
180}

Archive Download this file

Branches

Tags

Quick Links:     www.monotone.ca    -     Downloads    -     Documentation    -     Wiki    -     Code Forge    -     Build Status