xref: /dragonfly/crypto/libressl/crypto/bn/bn_gcd.c (revision 72c33676)
1 /* $OpenBSD: bn_gcd.c,v 1.15 2017/01/29 17:49:22 beck Exp $ */
2 /* Copyright (C) 1995-1998 Eric Young (eay@cryptsoft.com)
3  * All rights reserved.
4  *
5  * This package is an SSL implementation written
6  * by Eric Young (eay@cryptsoft.com).
7  * The implementation was written so as to conform with Netscapes SSL.
8  *
9  * This library is free for commercial and non-commercial use as long as
10  * the following conditions are aheared to.  The following conditions
11  * apply to all code found in this distribution, be it the RC4, RSA,
12  * lhash, DES, etc., code; not just the SSL code.  The SSL documentation
13  * included with this distribution is covered by the same copyright terms
14  * except that the holder is Tim Hudson (tjh@cryptsoft.com).
15  *
16  * Copyright remains Eric Young's, and as such any Copyright notices in
17  * the code are not to be removed.
18  * If this package is used in a product, Eric Young should be given attribution
19  * as the author of the parts of the library used.
20  * This can be in the form of a textual message at program startup or
21  * in documentation (online or textual) provided with the package.
22  *
23  * Redistribution and use in source and binary forms, with or without
24  * modification, are permitted provided that the following conditions
25  * are met:
26  * 1. Redistributions of source code must retain the copyright
27  *    notice, this list of conditions and the following disclaimer.
28  * 2. Redistributions in binary form must reproduce the above copyright
29  *    notice, this list of conditions and the following disclaimer in the
30  *    documentation and/or other materials provided with the distribution.
31  * 3. All advertising materials mentioning features or use of this software
32  *    must display the following acknowledgement:
33  *    "This product includes cryptographic software written by
34  *     Eric Young (eay@cryptsoft.com)"
35  *    The word 'cryptographic' can be left out if the rouines from the library
36  *    being used are not cryptographic related :-).
37  * 4. If you include any Windows specific code (or a derivative thereof) from
38  *    the apps directory (application code) you must include an acknowledgement:
39  *    "This product includes software written by Tim Hudson (tjh@cryptsoft.com)"
40  *
41  * THIS SOFTWARE IS PROVIDED BY ERIC YOUNG ``AS IS'' AND
42  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
43  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
44  * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
45  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
46  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
47  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
48  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
49  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
50  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
51  * SUCH DAMAGE.
52  *
53  * The licence and distribution terms for any publically available version or
54  * derivative of this code cannot be changed.  i.e. this code cannot simply be
55  * copied and put under another distribution licence
56  * [including the GNU Public Licence.]
57  */
58 /* ====================================================================
59  * Copyright (c) 1998-2001 The OpenSSL Project.  All rights reserved.
60  *
61  * Redistribution and use in source and binary forms, with or without
62  * modification, are permitted provided that the following conditions
63  * are met:
64  *
65  * 1. Redistributions of source code must retain the above copyright
66  *    notice, this list of conditions and the following disclaimer.
67  *
68  * 2. Redistributions in binary form must reproduce the above copyright
69  *    notice, this list of conditions and the following disclaimer in
70  *    the documentation and/or other materials provided with the
71  *    distribution.
72  *
73  * 3. All advertising materials mentioning features or use of this
74  *    software must display the following acknowledgment:
75  *    "This product includes software developed by the OpenSSL Project
76  *    for use in the OpenSSL Toolkit. (http://www.openssl.org/)"
77  *
78  * 4. The names "OpenSSL Toolkit" and "OpenSSL Project" must not be used to
79  *    endorse or promote products derived from this software without
80  *    prior written permission. For written permission, please contact
81  *    openssl-core@openssl.org.
82  *
83  * 5. Products derived from this software may not be called "OpenSSL"
84  *    nor may "OpenSSL" appear in their names without prior written
85  *    permission of the OpenSSL Project.
86  *
87  * 6. Redistributions of any form whatsoever must retain the following
88  *    acknowledgment:
89  *    "This product includes software developed by the OpenSSL Project
90  *    for use in the OpenSSL Toolkit (http://www.openssl.org/)"
91  *
92  * THIS SOFTWARE IS PROVIDED BY THE OpenSSL PROJECT ``AS IS'' AND ANY
93  * EXPRESSED OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
94  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
95  * PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE OpenSSL PROJECT OR
96  * ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
97  * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
98  * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
99  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
100  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
101  * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
102  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
103  * OF THE POSSIBILITY OF SUCH DAMAGE.
104  * ====================================================================
105  *
106  * This product includes cryptographic software written by Eric Young
107  * (eay@cryptsoft.com).  This product includes software written by Tim
108  * Hudson (tjh@cryptsoft.com).
109  *
110  */
111 
112 #include <openssl/err.h>
113 
114 #include "bn_lcl.h"
115 
116 static BIGNUM *euclid(BIGNUM *a, BIGNUM *b);
117 static BIGNUM *BN_gcd_no_branch(BIGNUM *in, const BIGNUM *a, const BIGNUM *n,
118     BN_CTX *ctx);
119 
120 int
121 BN_gcd(BIGNUM *r, const BIGNUM *in_a, const BIGNUM *in_b, BN_CTX *ctx)
122 {
123 	BIGNUM *a, *b, *t;
124 	int ret = 0;
125 
126 	bn_check_top(in_a);
127 	bn_check_top(in_b);
128 
129 	BN_CTX_start(ctx);
130 	if ((a = BN_CTX_get(ctx)) == NULL)
131 		goto err;
132 	if ((b = BN_CTX_get(ctx)) == NULL)
133 		goto err;
134 
135 	if (BN_copy(a, in_a) == NULL)
136 		goto err;
137 	if (BN_copy(b, in_b) == NULL)
138 		goto err;
139 	a->neg = 0;
140 	b->neg = 0;
141 
142 	if (BN_cmp(a, b) < 0) {
143 		t = a;
144 		a = b;
145 		b = t;
146 	}
147 	t = euclid(a, b);
148 	if (t == NULL)
149 		goto err;
150 
151 	if (BN_copy(r, t) == NULL)
152 		goto err;
153 	ret = 1;
154 
155 err:
156 	BN_CTX_end(ctx);
157 	bn_check_top(r);
158 	return (ret);
159 }
160 
161 int
162 BN_gcd_ct(BIGNUM *r, const BIGNUM *in_a, const BIGNUM *in_b, BN_CTX *ctx)
163 {
164 	if (BN_gcd_no_branch(r, in_a, in_b, ctx) == NULL)
165 		return 0;
166 	return 1;
167 }
168 
169 int
170 BN_gcd_nonct(BIGNUM *r, const BIGNUM *in_a, const BIGNUM *in_b, BN_CTX *ctx)
171 {
172 	return BN_gcd(r, in_a, in_b, ctx);
173 }
174 
175 
176 static BIGNUM *
177 euclid(BIGNUM *a, BIGNUM *b)
178 {
179 	BIGNUM *t;
180 	int shifts = 0;
181 
182 	bn_check_top(a);
183 	bn_check_top(b);
184 
185 	/* 0 <= b <= a */
186 	while (!BN_is_zero(b)) {
187 		/* 0 < b <= a */
188 
189 		if (BN_is_odd(a)) {
190 			if (BN_is_odd(b)) {
191 				if (!BN_sub(a, a, b))
192 					goto err;
193 				if (!BN_rshift1(a, a))
194 					goto err;
195 				if (BN_cmp(a, b) < 0) {
196 					t = a;
197 					a = b;
198 					b = t;
199 				}
200 			}
201 			else		/* a odd - b even */
202 			{
203 				if (!BN_rshift1(b, b))
204 					goto err;
205 				if (BN_cmp(a, b) < 0) {
206 					t = a;
207 					a = b;
208 					b = t;
209 				}
210 			}
211 		}
212 		else			/* a is even */
213 		{
214 			if (BN_is_odd(b)) {
215 				if (!BN_rshift1(a, a))
216 					goto err;
217 				if (BN_cmp(a, b) < 0) {
218 					t = a;
219 					a = b;
220 					b = t;
221 				}
222 			}
223 			else		/* a even - b even */
224 			{
225 				if (!BN_rshift1(a, a))
226 					goto err;
227 				if (!BN_rshift1(b, b))
228 					goto err;
229 				shifts++;
230 			}
231 		}
232 		/* 0 <= b <= a */
233 	}
234 
235 	if (shifts) {
236 		if (!BN_lshift(a, a, shifts))
237 			goto err;
238 	}
239 	bn_check_top(a);
240 	return (a);
241 
242 err:
243 	return (NULL);
244 }
245 
246 
247 /* solves ax == 1 (mod n) */
248 static BIGNUM *BN_mod_inverse_no_branch(BIGNUM *in, const BIGNUM *a,
249     const BIGNUM *n, BN_CTX *ctx);
250 
251 static BIGNUM *
252 BN_mod_inverse_internal(BIGNUM *in, const BIGNUM *a, const BIGNUM *n, BN_CTX *ctx,
253     int ct)
254 {
255 	BIGNUM *A, *B, *X, *Y, *M, *D, *T, *R = NULL;
256 	BIGNUM *ret = NULL;
257 	int sign;
258 
259 	if (ct)
260 		return BN_mod_inverse_no_branch(in, a, n, ctx);
261 
262 	bn_check_top(a);
263 	bn_check_top(n);
264 
265 	BN_CTX_start(ctx);
266 	if ((A = BN_CTX_get(ctx)) == NULL)
267 		goto err;
268 	if ((B = BN_CTX_get(ctx)) == NULL)
269 		goto err;
270 	if ((X = BN_CTX_get(ctx)) == NULL)
271 		goto err;
272 	if ((D = BN_CTX_get(ctx)) == NULL)
273 		goto err;
274 	if ((M = BN_CTX_get(ctx)) == NULL)
275 		goto err;
276 	if ((Y = BN_CTX_get(ctx)) == NULL)
277 		goto err;
278 	if ((T = BN_CTX_get(ctx)) == NULL)
279 		goto err;
280 
281 	if (in == NULL)
282 		R = BN_new();
283 	else
284 		R = in;
285 	if (R == NULL)
286 		goto err;
287 
288 	BN_one(X);
289 	BN_zero(Y);
290 	if (BN_copy(B, a) == NULL)
291 		goto err;
292 	if (BN_copy(A, n) == NULL)
293 		goto err;
294 	A->neg = 0;
295 	if (B->neg || (BN_ucmp(B, A) >= 0)) {
296 		if (!BN_nnmod(B, B, A, ctx))
297 			goto err;
298 	}
299 	sign = -1;
300 	/* From  B = a mod |n|,  A = |n|  it follows that
301 	 *
302 	 *      0 <= B < A,
303 	 *     -sign*X*a  ==  B   (mod |n|),
304 	 *      sign*Y*a  ==  A   (mod |n|).
305 	 */
306 
307 	if (BN_is_odd(n) && (BN_num_bits(n) <= (BN_BITS <= 32 ? 450 : 2048))) {
308 		/* Binary inversion algorithm; requires odd modulus.
309 		 * This is faster than the general algorithm if the modulus
310 		 * is sufficiently small (about 400 .. 500 bits on 32-bit
311 		 * sytems, but much more on 64-bit systems) */
312 		int shift;
313 
314 		while (!BN_is_zero(B)) {
315 			/*
316 			 *      0 < B < |n|,
317 			 *      0 < A <= |n|,
318 			 * (1) -sign*X*a  ==  B   (mod |n|),
319 			 * (2)  sign*Y*a  ==  A   (mod |n|)
320 			 */
321 
322 			/* Now divide  B  by the maximum possible power of two in the integers,
323 			 * and divide  X  by the same value mod |n|.
324 			 * When we're done, (1) still holds. */
325 			shift = 0;
326 			while (!BN_is_bit_set(B, shift)) /* note that 0 < B */
327 			{
328 				shift++;
329 
330 				if (BN_is_odd(X)) {
331 					if (!BN_uadd(X, X, n))
332 						goto err;
333 				}
334 				/* now X is even, so we can easily divide it by two */
335 				if (!BN_rshift1(X, X))
336 					goto err;
337 			}
338 			if (shift > 0) {
339 				if (!BN_rshift(B, B, shift))
340 					goto err;
341 			}
342 
343 
344 			/* Same for  A  and  Y.  Afterwards, (2) still holds. */
345 			shift = 0;
346 			while (!BN_is_bit_set(A, shift)) /* note that 0 < A */
347 			{
348 				shift++;
349 
350 				if (BN_is_odd(Y)) {
351 					if (!BN_uadd(Y, Y, n))
352 						goto err;
353 				}
354 				/* now Y is even */
355 				if (!BN_rshift1(Y, Y))
356 					goto err;
357 			}
358 			if (shift > 0) {
359 				if (!BN_rshift(A, A, shift))
360 					goto err;
361 			}
362 
363 
364 			/* We still have (1) and (2).
365 			 * Both  A  and  B  are odd.
366 			 * The following computations ensure that
367 			 *
368 			 *     0 <= B < |n|,
369 			 *      0 < A < |n|,
370 			 * (1) -sign*X*a  ==  B   (mod |n|),
371 			 * (2)  sign*Y*a  ==  A   (mod |n|),
372 			 *
373 			 * and that either  A  or  B  is even in the next iteration.
374 			 */
375 			if (BN_ucmp(B, A) >= 0) {
376 				/* -sign*(X + Y)*a == B - A  (mod |n|) */
377 				if (!BN_uadd(X, X, Y))
378 					goto err;
379 				/* NB: we could use BN_mod_add_quick(X, X, Y, n), but that
380 				 * actually makes the algorithm slower */
381 				if (!BN_usub(B, B, A))
382 					goto err;
383 			} else {
384 				/*  sign*(X + Y)*a == A - B  (mod |n|) */
385 				if (!BN_uadd(Y, Y, X))
386 					goto err;
387 				/* as above, BN_mod_add_quick(Y, Y, X, n) would slow things down */
388 				if (!BN_usub(A, A, B))
389 					goto err;
390 			}
391 		}
392 	} else {
393 		/* general inversion algorithm */
394 
395 		while (!BN_is_zero(B)) {
396 			BIGNUM *tmp;
397 
398 			/*
399 			 *      0 < B < A,
400 			 * (*) -sign*X*a  ==  B   (mod |n|),
401 			 *      sign*Y*a  ==  A   (mod |n|)
402 			 */
403 
404 			/* (D, M) := (A/B, A%B) ... */
405 			if (BN_num_bits(A) == BN_num_bits(B)) {
406 				if (!BN_one(D))
407 					goto err;
408 				if (!BN_sub(M, A, B))
409 					goto err;
410 			} else if (BN_num_bits(A) == BN_num_bits(B) + 1) {
411 				/* A/B is 1, 2, or 3 */
412 				if (!BN_lshift1(T, B))
413 					goto err;
414 				if (BN_ucmp(A, T) < 0) {
415 					/* A < 2*B, so D=1 */
416 					if (!BN_one(D))
417 						goto err;
418 					if (!BN_sub(M, A, B))
419 						goto err;
420 				} else {
421 					/* A >= 2*B, so D=2 or D=3 */
422 					if (!BN_sub(M, A, T))
423 						goto err;
424 					if (!BN_add(D,T,B)) goto err; /* use D (:= 3*B) as temp */
425 						if (BN_ucmp(A, D) < 0) {
426 						/* A < 3*B, so D=2 */
427 						if (!BN_set_word(D, 2))
428 							goto err;
429 						/* M (= A - 2*B) already has the correct value */
430 					} else {
431 						/* only D=3 remains */
432 						if (!BN_set_word(D, 3))
433 							goto err;
434 						/* currently  M = A - 2*B,  but we need  M = A - 3*B */
435 						if (!BN_sub(M, M, B))
436 							goto err;
437 					}
438 				}
439 			} else {
440 				if (!BN_div_nonct(D, M, A, B, ctx))
441 					goto err;
442 			}
443 
444 			/* Now
445 			 *      A = D*B + M;
446 			 * thus we have
447 			 * (**)  sign*Y*a  ==  D*B + M   (mod |n|).
448 			 */
449 			tmp = A; /* keep the BIGNUM object, the value does not matter */
450 
451 			/* (A, B) := (B, A mod B) ... */
452 			A = B;
453 			B = M;
454 			/* ... so we have  0 <= B < A  again */
455 
456 			/* Since the former  M  is now  B  and the former  B  is now  A,
457 			 * (**) translates into
458 			 *       sign*Y*a  ==  D*A + B    (mod |n|),
459 			 * i.e.
460 			 *       sign*Y*a - D*A  ==  B    (mod |n|).
461 			 * Similarly, (*) translates into
462 			 *      -sign*X*a  ==  A          (mod |n|).
463 			 *
464 			 * Thus,
465 			 *   sign*Y*a + D*sign*X*a  ==  B  (mod |n|),
466 			 * i.e.
467 			 *        sign*(Y + D*X)*a  ==  B  (mod |n|).
468 			 *
469 			 * So if we set  (X, Y, sign) := (Y + D*X, X, -sign),  we arrive back at
470 			 *      -sign*X*a  ==  B   (mod |n|),
471 			 *       sign*Y*a  ==  A   (mod |n|).
472 			 * Note that  X  and  Y  stay non-negative all the time.
473 			 */
474 
475 			/* most of the time D is very small, so we can optimize tmp := D*X+Y */
476 			if (BN_is_one(D)) {
477 				if (!BN_add(tmp, X, Y))
478 					goto err;
479 			} else {
480 				if (BN_is_word(D, 2)) {
481 					if (!BN_lshift1(tmp, X))
482 						goto err;
483 				} else if (BN_is_word(D, 4)) {
484 					if (!BN_lshift(tmp, X, 2))
485 						goto err;
486 				} else if (D->top == 1) {
487 					if (!BN_copy(tmp, X))
488 						goto err;
489 					if (!BN_mul_word(tmp, D->d[0]))
490 						goto err;
491 				} else {
492 					if (!BN_mul(tmp, D,X, ctx))
493 						goto err;
494 				}
495 				if (!BN_add(tmp, tmp, Y))
496 					goto err;
497 			}
498 
499 			M = Y; /* keep the BIGNUM object, the value does not matter */
500 			Y = X;
501 			X = tmp;
502 			sign = -sign;
503 		}
504 	}
505 
506 	/*
507 	 * The while loop (Euclid's algorithm) ends when
508 	 *      A == gcd(a,n);
509 	 * we have
510 	 *       sign*Y*a  ==  A  (mod |n|),
511 	 * where  Y  is non-negative.
512 	 */
513 
514 	if (sign < 0) {
515 		if (!BN_sub(Y, n, Y))
516 			goto err;
517 	}
518 	/* Now  Y*a  ==  A  (mod |n|).  */
519 
520 	if (BN_is_one(A)) {
521 		/* Y*a == 1  (mod |n|) */
522 		if (!Y->neg && BN_ucmp(Y, n) < 0) {
523 			if (!BN_copy(R, Y))
524 				goto err;
525 		} else {
526 			if (!BN_nnmod(R, Y,n, ctx))
527 				goto err;
528 		}
529 	} else {
530 		BNerror(BN_R_NO_INVERSE);
531 		goto err;
532 	}
533 	ret = R;
534 
535 err:
536 	if ((ret == NULL) && (in == NULL))
537 		BN_free(R);
538 	BN_CTX_end(ctx);
539 	bn_check_top(ret);
540 	return (ret);
541 }
542 
543 BIGNUM *
544 BN_mod_inverse(BIGNUM *in, const BIGNUM *a, const BIGNUM *n, BN_CTX *ctx)
545 {
546 	int ct = ((BN_get_flags(a, BN_FLG_CONSTTIME) != 0) ||
547 	    (BN_get_flags(n, BN_FLG_CONSTTIME) != 0));
548 	return BN_mod_inverse_internal(in, a, n, ctx, ct);
549 }
550 
551 BIGNUM *
552 BN_mod_inverse_nonct(BIGNUM *in, const BIGNUM *a, const BIGNUM *n, BN_CTX *ctx)
553 {
554 	return BN_mod_inverse_internal(in, a, n, ctx, 0);
555 }
556 
557 BIGNUM *
558 BN_mod_inverse_ct(BIGNUM *in, const BIGNUM *a, const BIGNUM *n, BN_CTX *ctx)
559 {
560 	return BN_mod_inverse_internal(in, a, n, ctx, 1);
561 }
562 
563 /* BN_mod_inverse_no_branch is a special version of BN_mod_inverse.
564  * It does not contain branches that may leak sensitive information.
565  */
566 static BIGNUM *
567 BN_mod_inverse_no_branch(BIGNUM *in, const BIGNUM *a, const BIGNUM *n,
568     BN_CTX *ctx)
569 {
570 	BIGNUM *A, *B, *X, *Y, *M, *D, *T, *R = NULL;
571 	BIGNUM local_A, local_B;
572 	BIGNUM *pA, *pB;
573 	BIGNUM *ret = NULL;
574 	int sign;
575 
576 	bn_check_top(a);
577 	bn_check_top(n);
578 
579 	BN_CTX_start(ctx);
580 	if ((A = BN_CTX_get(ctx)) == NULL)
581 		goto err;
582 	if ((B = BN_CTX_get(ctx)) == NULL)
583 		goto err;
584 	if ((X = BN_CTX_get(ctx)) == NULL)
585 		goto err;
586 	if ((D = BN_CTX_get(ctx)) == NULL)
587 		goto err;
588 	if ((M = BN_CTX_get(ctx)) == NULL)
589 		goto err;
590 	if ((Y = BN_CTX_get(ctx)) == NULL)
591 		goto err;
592 	if ((T = BN_CTX_get(ctx)) == NULL)
593 		goto err;
594 
595 	if (in == NULL)
596 		R = BN_new();
597 	else
598 		R = in;
599 	if (R == NULL)
600 		goto err;
601 
602 	BN_one(X);
603 	BN_zero(Y);
604 	if (BN_copy(B, a) == NULL)
605 		goto err;
606 	if (BN_copy(A, n) == NULL)
607 		goto err;
608 	A->neg = 0;
609 
610 	if (B->neg || (BN_ucmp(B, A) >= 0)) {
611 		/* Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
612 	 	 * BN_div_no_branch will be called eventually.
613 	 	 */
614 		pB = &local_B;
615 		BN_with_flags(pB, B, BN_FLG_CONSTTIME);
616 		if (!BN_nnmod(B, pB, A, ctx))
617 			goto err;
618 	}
619 	sign = -1;
620 	/* From  B = a mod |n|,  A = |n|  it follows that
621 	 *
622 	 *      0 <= B < A,
623 	 *     -sign*X*a  ==  B   (mod |n|),
624 	 *      sign*Y*a  ==  A   (mod |n|).
625 	 */
626 
627 	while (!BN_is_zero(B)) {
628 		BIGNUM *tmp;
629 
630 		/*
631 		 *      0 < B < A,
632 		 * (*) -sign*X*a  ==  B   (mod |n|),
633 		 *      sign*Y*a  ==  A   (mod |n|)
634 		 */
635 
636 		/* Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
637 	 	 * BN_div_no_branch will be called eventually.
638 	 	 */
639 		pA = &local_A;
640 		BN_with_flags(pA, A, BN_FLG_CONSTTIME);
641 
642 		/* (D, M) := (A/B, A%B) ... */
643 		if (!BN_div_ct(D, M, pA, B, ctx))
644 			goto err;
645 
646 		/* Now
647 		 *      A = D*B + M;
648 		 * thus we have
649 		 * (**)  sign*Y*a  ==  D*B + M   (mod |n|).
650 		 */
651 		tmp = A; /* keep the BIGNUM object, the value does not matter */
652 
653 		/* (A, B) := (B, A mod B) ... */
654 		A = B;
655 		B = M;
656 		/* ... so we have  0 <= B < A  again */
657 
658 		/* Since the former  M  is now  B  and the former  B  is now  A,
659 		 * (**) translates into
660 		 *       sign*Y*a  ==  D*A + B    (mod |n|),
661 		 * i.e.
662 		 *       sign*Y*a - D*A  ==  B    (mod |n|).
663 		 * Similarly, (*) translates into
664 		 *      -sign*X*a  ==  A          (mod |n|).
665 		 *
666 		 * Thus,
667 		 *   sign*Y*a + D*sign*X*a  ==  B  (mod |n|),
668 		 * i.e.
669 		 *        sign*(Y + D*X)*a  ==  B  (mod |n|).
670 		 *
671 		 * So if we set  (X, Y, sign) := (Y + D*X, X, -sign),  we arrive back at
672 		 *      -sign*X*a  ==  B   (mod |n|),
673 		 *       sign*Y*a  ==  A   (mod |n|).
674 		 * Note that  X  and  Y  stay non-negative all the time.
675 		 */
676 
677 		if (!BN_mul(tmp, D, X, ctx))
678 			goto err;
679 		if (!BN_add(tmp, tmp, Y))
680 			goto err;
681 
682 		M = Y; /* keep the BIGNUM object, the value does not matter */
683 		Y = X;
684 		X = tmp;
685 		sign = -sign;
686 	}
687 
688 	/*
689 	 * The while loop (Euclid's algorithm) ends when
690 	 *      A == gcd(a,n);
691 	 * we have
692 	 *       sign*Y*a  ==  A  (mod |n|),
693 	 * where  Y  is non-negative.
694 	 */
695 
696 	if (sign < 0) {
697 		if (!BN_sub(Y, n, Y))
698 			goto err;
699 	}
700 	/* Now  Y*a  ==  A  (mod |n|).  */
701 
702 	if (BN_is_one(A)) {
703 		/* Y*a == 1  (mod |n|) */
704 		if (!Y->neg && BN_ucmp(Y, n) < 0) {
705 			if (!BN_copy(R, Y))
706 				goto err;
707 		} else {
708 			if (!BN_nnmod(R, Y, n, ctx))
709 				goto err;
710 		}
711 	} else {
712 		BNerror(BN_R_NO_INVERSE);
713 		goto err;
714 	}
715 	ret = R;
716 
717 err:
718 	if ((ret == NULL) && (in == NULL))
719 		BN_free(R);
720 	BN_CTX_end(ctx);
721 	bn_check_top(ret);
722 	return (ret);
723 }
724 
725 /*
726  * BN_gcd_no_branch is a special version of BN_mod_inverse_no_branch.
727  * that returns the GCD.
728  */
729 static BIGNUM *
730 BN_gcd_no_branch(BIGNUM *in, const BIGNUM *a, const BIGNUM *n,
731     BN_CTX *ctx)
732 {
733 	BIGNUM *A, *B, *X, *Y, *M, *D, *T, *R = NULL;
734 	BIGNUM local_A, local_B;
735 	BIGNUM *pA, *pB;
736 	BIGNUM *ret = NULL;
737 	int sign;
738 
739 	if (in == NULL)
740 		goto err;
741 	R = in;
742 
743 	bn_check_top(a);
744 	bn_check_top(n);
745 
746 	BN_CTX_start(ctx);
747 	if ((A = BN_CTX_get(ctx)) == NULL)
748 		goto err;
749 	if ((B = BN_CTX_get(ctx)) == NULL)
750 		goto err;
751 	if ((X = BN_CTX_get(ctx)) == NULL)
752 		goto err;
753 	if ((D = BN_CTX_get(ctx)) == NULL)
754 		goto err;
755 	if ((M = BN_CTX_get(ctx)) == NULL)
756 		goto err;
757 	if ((Y = BN_CTX_get(ctx)) == NULL)
758 		goto err;
759 	if ((T = BN_CTX_get(ctx)) == NULL)
760 		goto err;
761 
762 	BN_one(X);
763 	BN_zero(Y);
764 	if (BN_copy(B, a) == NULL)
765 		goto err;
766 	if (BN_copy(A, n) == NULL)
767 		goto err;
768 	A->neg = 0;
769 
770 	if (B->neg || (BN_ucmp(B, A) >= 0)) {
771 		/* Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
772 	 	 * BN_div_no_branch will be called eventually.
773 	 	 */
774 		pB = &local_B;
775 		BN_with_flags(pB, B, BN_FLG_CONSTTIME);
776 		if (!BN_nnmod(B, pB, A, ctx))
777 			goto err;
778 	}
779 	sign = -1;
780 	/* From  B = a mod |n|,  A = |n|  it follows that
781 	 *
782 	 *      0 <= B < A,
783 	 *     -sign*X*a  ==  B   (mod |n|),
784 	 *      sign*Y*a  ==  A   (mod |n|).
785 	 */
786 
787 	while (!BN_is_zero(B)) {
788 		BIGNUM *tmp;
789 
790 		/*
791 		 *      0 < B < A,
792 		 * (*) -sign*X*a  ==  B   (mod |n|),
793 		 *      sign*Y*a  ==  A   (mod |n|)
794 		 */
795 
796 		/* Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
797 	 	 * BN_div_no_branch will be called eventually.
798 	 	 */
799 		pA = &local_A;
800 		BN_with_flags(pA, A, BN_FLG_CONSTTIME);
801 
802 		/* (D, M) := (A/B, A%B) ... */
803 		if (!BN_div_ct(D, M, pA, B, ctx))
804 			goto err;
805 
806 		/* Now
807 		 *      A = D*B + M;
808 		 * thus we have
809 		 * (**)  sign*Y*a  ==  D*B + M   (mod |n|).
810 		 */
811 		tmp = A; /* keep the BIGNUM object, the value does not matter */
812 
813 		/* (A, B) := (B, A mod B) ... */
814 		A = B;
815 		B = M;
816 		/* ... so we have  0 <= B < A  again */
817 
818 		/* Since the former  M  is now  B  and the former  B  is now  A,
819 		 * (**) translates into
820 		 *       sign*Y*a  ==  D*A + B    (mod |n|),
821 		 * i.e.
822 		 *       sign*Y*a - D*A  ==  B    (mod |n|).
823 		 * Similarly, (*) translates into
824 		 *      -sign*X*a  ==  A          (mod |n|).
825 		 *
826 		 * Thus,
827 		 *   sign*Y*a + D*sign*X*a  ==  B  (mod |n|),
828 		 * i.e.
829 		 *        sign*(Y + D*X)*a  ==  B  (mod |n|).
830 		 *
831 		 * So if we set  (X, Y, sign) := (Y + D*X, X, -sign),  we arrive back at
832 		 *      -sign*X*a  ==  B   (mod |n|),
833 		 *       sign*Y*a  ==  A   (mod |n|).
834 		 * Note that  X  and  Y  stay non-negative all the time.
835 		 */
836 
837 		if (!BN_mul(tmp, D, X, ctx))
838 			goto err;
839 		if (!BN_add(tmp, tmp, Y))
840 			goto err;
841 
842 		M = Y; /* keep the BIGNUM object, the value does not matter */
843 		Y = X;
844 		X = tmp;
845 		sign = -sign;
846 	}
847 
848 	/*
849 	 * The while loop (Euclid's algorithm) ends when
850 	 *      A == gcd(a,n);
851 	 */
852 
853 	if (!BN_copy(R, A))
854 		goto err;
855 	ret = R;
856 err:
857 	if ((ret == NULL) && (in == NULL))
858 		BN_free(R);
859 	BN_CTX_end(ctx);
860 	bn_check_top(ret);
861 	return (ret);
862 }
863