Mercurial > dropbear
comparison libtommath/bn_mp_div.c @ 1655:f52919ffd3b1
update ltm to 1.1.0 and enable FIPS 186.4 compliant key-generation (#79)
* make key-generation compliant to FIPS 186.4
* fix includes in tommath_class.h
* update fuzzcorpus instead of error-out
* fixup fuzzing make-targets
* update Makefile.in
* apply necessary patches to ltm sources
* clean-up not required ltm files
* update to vanilla ltm 1.1.0
this already only contains the required files
* remove set/get double
author | Steffen Jaeckel <s_jaeckel@gmx.de> |
---|---|
date | Mon, 16 Sep 2019 15:50:38 +0200 |
parents | 8bba51a55704 |
children | 1051e4eea25a |
comparison
equal
deleted
inserted
replaced
1654:cc0fc5131c5c | 1655:f52919ffd3b1 |
---|---|
1 #include <tommath_private.h> | 1 #include "tommath_private.h" |
2 #ifdef BN_MP_DIV_C | 2 #ifdef BN_MP_DIV_C |
3 /* LibTomMath, multiple-precision integer library -- Tom St Denis | 3 /* LibTomMath, multiple-precision integer library -- Tom St Denis |
4 * | 4 * |
5 * LibTomMath is a library that provides multiple-precision | 5 * LibTomMath is a library that provides multiple-precision |
6 * integer arithmetic as well as number theoretic functionality. | 6 * integer arithmetic as well as number theoretic functionality. |
7 * | 7 * |
8 * The library was designed directly after the MPI library by | 8 * The library was designed directly after the MPI library by |
9 * Michael Fromberger but has been written from scratch with | 9 * Michael Fromberger but has been written from scratch with |
10 * additional optimizations in place. | 10 * additional optimizations in place. |
11 * | 11 * |
12 * The library is free for all purposes without any express | 12 * SPDX-License-Identifier: Unlicense |
13 * guarantee it works. | |
14 * | |
15 * Tom St Denis, [email protected], http://libtom.org | |
16 */ | 13 */ |
17 | 14 |
18 #ifdef BN_MP_DIV_SMALL | 15 #ifdef BN_MP_DIV_SMALL |
19 | 16 |
20 /* slower bit-bang division... also smaller */ | 17 /* slower bit-bang division... also smaller */ |
21 int mp_div(mp_int * a, mp_int * b, mp_int * c, mp_int * d) | 18 int mp_div(const mp_int *a, const mp_int *b, mp_int *c, mp_int *d) |
22 { | 19 { |
23 mp_int ta, tb, tq, q; | 20 mp_int ta, tb, tq, q; |
24 int res, n, n2; | 21 int res, n, n2; |
25 | 22 |
26 /* is divisor zero ? */ | 23 /* is divisor zero ? */ |
27 if (mp_iszero (b) == MP_YES) { | 24 if (mp_iszero(b) == MP_YES) { |
28 return MP_VAL; | 25 return MP_VAL; |
29 } | 26 } |
30 | 27 |
31 /* if a < b then q=0, r = a */ | 28 /* if a < b then q=0, r = a */ |
32 if (mp_cmp_mag (a, b) == MP_LT) { | 29 if (mp_cmp_mag(a, b) == MP_LT) { |
33 if (d != NULL) { | 30 if (d != NULL) { |
34 res = mp_copy (a, d); | 31 res = mp_copy(a, d); |
35 } else { | 32 } else { |
36 res = MP_OKAY; | 33 res = MP_OKAY; |
37 } | 34 } |
38 if (c != NULL) { | 35 if (c != NULL) { |
39 mp_zero (c); | 36 mp_zero(c); |
40 } | 37 } |
41 return res; | 38 return res; |
42 } | 39 } |
43 | 40 |
44 /* init our temps */ | 41 /* init our temps */ |
45 if ((res = mp_init_multi(&ta, &tb, &tq, &q, NULL)) != MP_OKAY) { | 42 if ((res = mp_init_multi(&ta, &tb, &tq, &q, NULL)) != MP_OKAY) { |
46 return res; | 43 return res; |
47 } | 44 } |
48 | 45 |
49 | 46 |
50 mp_set(&tq, 1); | 47 mp_set(&tq, 1uL); |
51 n = mp_count_bits(a) - mp_count_bits(b); | 48 n = mp_count_bits(a) - mp_count_bits(b); |
52 if (((res = mp_abs(a, &ta)) != MP_OKAY) || | 49 if (((res = mp_abs(a, &ta)) != MP_OKAY) || |
53 ((res = mp_abs(b, &tb)) != MP_OKAY) || | 50 ((res = mp_abs(b, &tb)) != MP_OKAY) || |
54 ((res = mp_mul_2d(&tb, n, &tb)) != MP_OKAY) || | 51 ((res = mp_mul_2d(&tb, n, &tb)) != MP_OKAY) || |
55 ((res = mp_mul_2d(&tq, n, &tq)) != MP_OKAY)) { | 52 ((res = mp_mul_2d(&tq, n, &tq)) != MP_OKAY)) { |
56 goto LBL_ERR; | 53 goto LBL_ERR; |
57 } | 54 } |
58 | 55 |
59 while (n-- >= 0) { | 56 while (n-- >= 0) { |
60 if (mp_cmp(&tb, &ta) != MP_GT) { | 57 if (mp_cmp(&tb, &ta) != MP_GT) { |
61 if (((res = mp_sub(&ta, &tb, &ta)) != MP_OKAY) || | 58 if (((res = mp_sub(&ta, &tb, &ta)) != MP_OKAY) || |
62 ((res = mp_add(&q, &tq, &q)) != MP_OKAY)) { | 59 ((res = mp_add(&q, &tq, &q)) != MP_OKAY)) { |
63 goto LBL_ERR; | 60 goto LBL_ERR; |
64 } | 61 } |
65 } | 62 } |
66 if (((res = mp_div_2d(&tb, 1, &tb, NULL)) != MP_OKAY) || | 63 if (((res = mp_div_2d(&tb, 1, &tb, NULL)) != MP_OKAY) || |
67 ((res = mp_div_2d(&tq, 1, &tq, NULL)) != MP_OKAY)) { | 64 ((res = mp_div_2d(&tq, 1, &tq, NULL)) != MP_OKAY)) { |
68 goto LBL_ERR; | 65 goto LBL_ERR; |
69 } | 66 } |
70 } | 67 } |
71 | 68 |
72 /* now q == quotient and ta == remainder */ | 69 /* now q == quotient and ta == remainder */ |
73 n = a->sign; | 70 n = a->sign; |
74 n2 = (a->sign == b->sign) ? MP_ZPOS : MP_NEG; | 71 n2 = (a->sign == b->sign) ? MP_ZPOS : MP_NEG; |
75 if (c != NULL) { | 72 if (c != NULL) { |
76 mp_exch(c, &q); | 73 mp_exch(c, &q); |
77 c->sign = (mp_iszero(c) == MP_YES) ? MP_ZPOS : n2; | 74 c->sign = (mp_iszero(c) == MP_YES) ? MP_ZPOS : n2; |
78 } | 75 } |
79 if (d != NULL) { | 76 if (d != NULL) { |
80 mp_exch(d, &ta); | 77 mp_exch(d, &ta); |
81 d->sign = (mp_iszero(d) == MP_YES) ? MP_ZPOS : n; | 78 d->sign = (mp_iszero(d) == MP_YES) ? MP_ZPOS : n; |
82 } | 79 } |
83 LBL_ERR: | 80 LBL_ERR: |
84 mp_clear_multi(&ta, &tb, &tq, &q, NULL); | 81 mp_clear_multi(&ta, &tb, &tq, &q, NULL); |
85 return res; | 82 return res; |
86 } | 83 } |
87 | 84 |
98 * case that y has fewer than three digits, etc.. | 95 * case that y has fewer than three digits, etc.. |
99 * | 96 * |
100 * The overall algorithm is as described as | 97 * The overall algorithm is as described as |
101 * 14.20 from HAC but fixed to treat these cases. | 98 * 14.20 from HAC but fixed to treat these cases. |
102 */ | 99 */ |
103 int mp_div (mp_int * a, mp_int * b, mp_int * c, mp_int * d) | 100 int mp_div(const mp_int *a, const mp_int *b, mp_int *c, mp_int *d) |
104 { | 101 { |
105 mp_int q, x, y, t1, t2; | 102 mp_int q, x, y, t1, t2; |
106 int res, n, t, i, norm, neg; | 103 int res, n, t, i, norm, neg; |
107 | 104 |
108 /* is divisor zero ? */ | 105 /* is divisor zero ? */ |
109 if (mp_iszero (b) == MP_YES) { | 106 if (mp_iszero(b) == MP_YES) { |
110 return MP_VAL; | 107 return MP_VAL; |
111 } | 108 } |
112 | 109 |
113 /* if a < b then q=0, r = a */ | 110 /* if a < b then q=0, r = a */ |
114 if (mp_cmp_mag (a, b) == MP_LT) { | 111 if (mp_cmp_mag(a, b) == MP_LT) { |
115 if (d != NULL) { | 112 if (d != NULL) { |
116 res = mp_copy (a, d); | 113 res = mp_copy(a, d); |
117 } else { | 114 } else { |
118 res = MP_OKAY; | 115 res = MP_OKAY; |
119 } | 116 } |
120 if (c != NULL) { | 117 if (c != NULL) { |
121 mp_zero (c); | 118 mp_zero(c); |
122 } | 119 } |
123 return res; | 120 return res; |
124 } | 121 } |
125 | 122 |
126 if ((res = mp_init_size (&q, a->used + 2)) != MP_OKAY) { | 123 if ((res = mp_init_size(&q, a->used + 2)) != MP_OKAY) { |
127 return res; | 124 return res; |
128 } | 125 } |
129 q.used = a->used + 2; | 126 q.used = a->used + 2; |
130 | 127 |
131 if ((res = mp_init (&t1)) != MP_OKAY) { | 128 if ((res = mp_init(&t1)) != MP_OKAY) { |
132 goto LBL_Q; | 129 goto LBL_Q; |
133 } | 130 } |
134 | 131 |
135 if ((res = mp_init (&t2)) != MP_OKAY) { | 132 if ((res = mp_init(&t2)) != MP_OKAY) { |
136 goto LBL_T1; | 133 goto LBL_T1; |
137 } | 134 } |
138 | 135 |
139 if ((res = mp_init_copy (&x, a)) != MP_OKAY) { | 136 if ((res = mp_init_copy(&x, a)) != MP_OKAY) { |
140 goto LBL_T2; | 137 goto LBL_T2; |
141 } | 138 } |
142 | 139 |
143 if ((res = mp_init_copy (&y, b)) != MP_OKAY) { | 140 if ((res = mp_init_copy(&y, b)) != MP_OKAY) { |
144 goto LBL_X; | 141 goto LBL_X; |
145 } | 142 } |
146 | 143 |
147 /* fix the sign */ | 144 /* fix the sign */ |
148 neg = (a->sign == b->sign) ? MP_ZPOS : MP_NEG; | 145 neg = (a->sign == b->sign) ? MP_ZPOS : MP_NEG; |
149 x.sign = y.sign = MP_ZPOS; | 146 x.sign = y.sign = MP_ZPOS; |
150 | 147 |
151 /* normalize both x and y, ensure that y >= b/2, [b == 2**DIGIT_BIT] */ | 148 /* normalize both x and y, ensure that y >= b/2, [b == 2**DIGIT_BIT] */ |
152 norm = mp_count_bits(&y) % DIGIT_BIT; | 149 norm = mp_count_bits(&y) % DIGIT_BIT; |
153 if (norm < (int)(DIGIT_BIT-1)) { | 150 if (norm < (DIGIT_BIT - 1)) { |
154 norm = (DIGIT_BIT-1) - norm; | 151 norm = (DIGIT_BIT - 1) - norm; |
155 if ((res = mp_mul_2d (&x, norm, &x)) != MP_OKAY) { | 152 if ((res = mp_mul_2d(&x, norm, &x)) != MP_OKAY) { |
156 goto LBL_Y; | 153 goto LBL_Y; |
157 } | 154 } |
158 if ((res = mp_mul_2d (&y, norm, &y)) != MP_OKAY) { | 155 if ((res = mp_mul_2d(&y, norm, &y)) != MP_OKAY) { |
159 goto LBL_Y; | 156 goto LBL_Y; |
160 } | 157 } |
161 } else { | 158 } else { |
162 norm = 0; | 159 norm = 0; |
163 } | 160 } |
164 | 161 |
165 /* note hac does 0 based, so if used==5 then its 0,1,2,3,4, e.g. use 4 */ | 162 /* note hac does 0 based, so if used==5 then its 0,1,2,3,4, e.g. use 4 */ |
166 n = x.used - 1; | 163 n = x.used - 1; |
167 t = y.used - 1; | 164 t = y.used - 1; |
168 | 165 |
169 /* while (x >= y*b**n-t) do { q[n-t] += 1; x -= y*b**{n-t} } */ | 166 /* while (x >= y*b**n-t) do { q[n-t] += 1; x -= y*b**{n-t} } */ |
170 if ((res = mp_lshd (&y, n - t)) != MP_OKAY) { /* y = y*b**{n-t} */ | 167 if ((res = mp_lshd(&y, n - t)) != MP_OKAY) { /* y = y*b**{n-t} */ |
171 goto LBL_Y; | |
172 } | |
173 | |
174 while (mp_cmp (&x, &y) != MP_LT) { | |
175 ++(q.dp[n - t]); | |
176 if ((res = mp_sub (&x, &y, &x)) != MP_OKAY) { | |
177 goto LBL_Y; | 168 goto LBL_Y; |
178 } | 169 } |
179 } | 170 |
180 | 171 while (mp_cmp(&x, &y) != MP_LT) { |
181 /* reset y by shifting it back down */ | 172 ++(q.dp[n - t]); |
182 mp_rshd (&y, n - t); | 173 if ((res = mp_sub(&x, &y, &x)) != MP_OKAY) { |
183 | 174 goto LBL_Y; |
184 /* step 3. for i from n down to (t + 1) */ | 175 } |
185 for (i = n; i >= (t + 1); i--) { | 176 } |
186 if (i > x.used) { | 177 |
187 continue; | 178 /* reset y by shifting it back down */ |
188 } | 179 mp_rshd(&y, n - t); |
189 | 180 |
190 /* step 3.1 if xi == yt then set q{i-t-1} to b-1, | 181 /* step 3. for i from n down to (t + 1) */ |
191 * otherwise set q{i-t-1} to (xi*b + x{i-1})/yt */ | 182 for (i = n; i >= (t + 1); i--) { |
192 if (x.dp[i] == y.dp[t]) { | 183 if (i > x.used) { |
193 q.dp[(i - t) - 1] = ((((mp_digit)1) << DIGIT_BIT) - 1); | 184 continue; |
194 } else { | 185 } |
195 mp_word tmp; | 186 |
196 tmp = ((mp_word) x.dp[i]) << ((mp_word) DIGIT_BIT); | 187 /* step 3.1 if xi == yt then set q{i-t-1} to b-1, |
197 tmp |= ((mp_word) x.dp[i - 1]); | 188 * otherwise set q{i-t-1} to (xi*b + x{i-1})/yt */ |
198 tmp /= ((mp_word) y.dp[t]); | 189 if (x.dp[i] == y.dp[t]) { |
199 if (tmp > (mp_word) MP_MASK) { | 190 q.dp[(i - t) - 1] = ((mp_digit)1 << (mp_digit)DIGIT_BIT) - (mp_digit)1; |
200 tmp = MP_MASK; | 191 } else { |
201 } | 192 mp_word tmp; |
202 q.dp[(i - t) - 1] = (mp_digit) (tmp & (mp_word) (MP_MASK)); | 193 tmp = (mp_word)x.dp[i] << (mp_word)DIGIT_BIT; |
203 } | 194 tmp |= (mp_word)x.dp[i - 1]; |
204 | 195 tmp /= (mp_word)y.dp[t]; |
205 /* while (q{i-t-1} * (yt * b + y{t-1})) > | 196 if (tmp > (mp_word)MP_MASK) { |
206 xi * b**2 + xi-1 * b + xi-2 | 197 tmp = MP_MASK; |
207 | 198 } |
208 do q{i-t-1} -= 1; | 199 q.dp[(i - t) - 1] = (mp_digit)(tmp & (mp_word)MP_MASK); |
200 } | |
201 | |
202 /* while (q{i-t-1} * (yt * b + y{t-1})) > | |
203 xi * b**2 + xi-1 * b + xi-2 | |
204 | |
205 do q{i-t-1} -= 1; | |
206 */ | |
207 q.dp[(i - t) - 1] = (q.dp[(i - t) - 1] + 1uL) & (mp_digit)MP_MASK; | |
208 do { | |
209 q.dp[(i - t) - 1] = (q.dp[(i - t) - 1] - 1uL) & (mp_digit)MP_MASK; | |
210 | |
211 /* find left hand */ | |
212 mp_zero(&t1); | |
213 t1.dp[0] = ((t - 1) < 0) ? 0u : y.dp[t - 1]; | |
214 t1.dp[1] = y.dp[t]; | |
215 t1.used = 2; | |
216 if ((res = mp_mul_d(&t1, q.dp[(i - t) - 1], &t1)) != MP_OKAY) { | |
217 goto LBL_Y; | |
218 } | |
219 | |
220 /* find right hand */ | |
221 t2.dp[0] = ((i - 2) < 0) ? 0u : x.dp[i - 2]; | |
222 t2.dp[1] = ((i - 1) < 0) ? 0u : x.dp[i - 1]; | |
223 t2.dp[2] = x.dp[i]; | |
224 t2.used = 3; | |
225 } while (mp_cmp_mag(&t1, &t2) == MP_GT); | |
226 | |
227 /* step 3.3 x = x - q{i-t-1} * y * b**{i-t-1} */ | |
228 if ((res = mp_mul_d(&y, q.dp[(i - t) - 1], &t1)) != MP_OKAY) { | |
229 goto LBL_Y; | |
230 } | |
231 | |
232 if ((res = mp_lshd(&t1, (i - t) - 1)) != MP_OKAY) { | |
233 goto LBL_Y; | |
234 } | |
235 | |
236 if ((res = mp_sub(&x, &t1, &x)) != MP_OKAY) { | |
237 goto LBL_Y; | |
238 } | |
239 | |
240 /* if x < 0 then { x = x + y*b**{i-t-1}; q{i-t-1} -= 1; } */ | |
241 if (x.sign == MP_NEG) { | |
242 if ((res = mp_copy(&y, &t1)) != MP_OKAY) { | |
243 goto LBL_Y; | |
244 } | |
245 if ((res = mp_lshd(&t1, (i - t) - 1)) != MP_OKAY) { | |
246 goto LBL_Y; | |
247 } | |
248 if ((res = mp_add(&x, &t1, &x)) != MP_OKAY) { | |
249 goto LBL_Y; | |
250 } | |
251 | |
252 q.dp[(i - t) - 1] = (q.dp[(i - t) - 1] - 1uL) & MP_MASK; | |
253 } | |
254 } | |
255 | |
256 /* now q is the quotient and x is the remainder | |
257 * [which we have to normalize] | |
209 */ | 258 */ |
210 q.dp[(i - t) - 1] = (q.dp[(i - t) - 1] + 1) & MP_MASK; | 259 |
211 do { | 260 /* get sign before writing to c */ |
212 q.dp[(i - t) - 1] = (q.dp[(i - t) - 1] - 1) & MP_MASK; | 261 x.sign = (x.used == 0) ? MP_ZPOS : a->sign; |
213 | 262 |
214 /* find left hand */ | 263 if (c != NULL) { |
215 mp_zero (&t1); | 264 mp_clamp(&q); |
216 t1.dp[0] = ((t - 1) < 0) ? 0 : y.dp[t - 1]; | 265 mp_exch(&q, c); |
217 t1.dp[1] = y.dp[t]; | 266 c->sign = neg; |
218 t1.used = 2; | 267 } |
219 if ((res = mp_mul_d (&t1, q.dp[(i - t) - 1], &t1)) != MP_OKAY) { | 268 |
220 goto LBL_Y; | 269 if (d != NULL) { |
221 } | 270 if ((res = mp_div_2d(&x, norm, &x, NULL)) != MP_OKAY) { |
222 | 271 goto LBL_Y; |
223 /* find right hand */ | 272 } |
224 t2.dp[0] = ((i - 2) < 0) ? 0 : x.dp[i - 2]; | 273 mp_exch(&x, d); |
225 t2.dp[1] = ((i - 1) < 0) ? 0 : x.dp[i - 1]; | 274 } |
226 t2.dp[2] = x.dp[i]; | 275 |
227 t2.used = 3; | 276 res = MP_OKAY; |
228 } while (mp_cmp_mag(&t1, &t2) == MP_GT); | 277 |
229 | 278 LBL_Y: |
230 /* step 3.3 x = x - q{i-t-1} * y * b**{i-t-1} */ | 279 mp_clear(&y); |
231 if ((res = mp_mul_d (&y, q.dp[(i - t) - 1], &t1)) != MP_OKAY) { | 280 LBL_X: |
232 goto LBL_Y; | 281 mp_clear(&x); |
233 } | 282 LBL_T2: |
234 | 283 mp_clear(&t2); |
235 if ((res = mp_lshd (&t1, (i - t) - 1)) != MP_OKAY) { | 284 LBL_T1: |
236 goto LBL_Y; | 285 mp_clear(&t1); |
237 } | 286 LBL_Q: |
238 | 287 mp_clear(&q); |
239 if ((res = mp_sub (&x, &t1, &x)) != MP_OKAY) { | 288 return res; |
240 goto LBL_Y; | |
241 } | |
242 | |
243 /* if x < 0 then { x = x + y*b**{i-t-1}; q{i-t-1} -= 1; } */ | |
244 if (x.sign == MP_NEG) { | |
245 if ((res = mp_copy (&y, &t1)) != MP_OKAY) { | |
246 goto LBL_Y; | |
247 } | |
248 if ((res = mp_lshd (&t1, (i - t) - 1)) != MP_OKAY) { | |
249 goto LBL_Y; | |
250 } | |
251 if ((res = mp_add (&x, &t1, &x)) != MP_OKAY) { | |
252 goto LBL_Y; | |
253 } | |
254 | |
255 q.dp[(i - t) - 1] = (q.dp[(i - t) - 1] - 1UL) & MP_MASK; | |
256 } | |
257 } | |
258 | |
259 /* now q is the quotient and x is the remainder | |
260 * [which we have to normalize] | |
261 */ | |
262 | |
263 /* get sign before writing to c */ | |
264 x.sign = (x.used == 0) ? MP_ZPOS : a->sign; | |
265 | |
266 if (c != NULL) { | |
267 mp_clamp (&q); | |
268 mp_exch (&q, c); | |
269 c->sign = neg; | |
270 } | |
271 | |
272 if (d != NULL) { | |
273 if ((res = mp_div_2d (&x, norm, &x, NULL)) != MP_OKAY) { | |
274 goto LBL_Y; | |
275 } | |
276 mp_exch (&x, d); | |
277 } | |
278 | |
279 res = MP_OKAY; | |
280 | |
281 LBL_Y:mp_clear (&y); | |
282 LBL_X:mp_clear (&x); | |
283 LBL_T2:mp_clear (&t2); | |
284 LBL_T1:mp_clear (&t1); | |
285 LBL_Q:mp_clear (&q); | |
286 return res; | |
287 } | 289 } |
288 | 290 |
289 #endif | 291 #endif |
290 | 292 |
291 #endif | 293 #endif |
292 | 294 |
293 /* ref: $Format:%D$ */ | 295 /* ref: HEAD -> master, tag: v1.1.0 */ |
294 /* git commit: $Format:%H$ */ | 296 /* git commit: 08549ad6bc8b0cede0b357a9c341c5c6473a9c55 */ |
295 /* commit time: $Format:%ai$ */ | 297 /* commit time: 2019-01-28 20:32:32 +0100 */ |