xref: /optee_os/lib/libmbedtls/mbedtls/library/bignum.c (revision 731185b11620a6de1f824278ab3c166c7853ef66)
1 /*
2  *  Multi-precision integer library
3  *
4  *  Copyright The Mbed TLS Contributors
5  *  SPDX-License-Identifier: Apache-2.0
6  *
7  *  Licensed under the Apache License, Version 2.0 (the "License"); you may
8  *  not use this file except in compliance with the License.
9  *  You may obtain a copy of the License at
10  *
11  *  http://www.apache.org/licenses/LICENSE-2.0
12  *
13  *  Unless required by applicable law or agreed to in writing, software
14  *  distributed under the License is distributed on an "AS IS" BASIS, WITHOUT
15  *  WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
16  *  See the License for the specific language governing permissions and
17  *  limitations under the License.
18  */
19 
20 /*
21  *  The following sources were referenced in the design of this Multi-precision
22  *  Integer library:
23  *
24  *  [1] Handbook of Applied Cryptography - 1997
25  *      Menezes, van Oorschot and Vanstone
26  *
27  *  [2] Multi-Precision Math
28  *      Tom St Denis
29  *      https://github.com/libtom/libtommath/blob/develop/tommath.pdf
30  *
31  *  [3] GNU Multi-Precision Arithmetic Library
32  *      https://gmplib.org/manual/index.html
33  *
34  */
35 
36 #include "common.h"
37 
38 #if defined(MBEDTLS_BIGNUM_C)
39 
40 #include "mbedtls/bignum.h"
41 #include "bignum_core.h"
42 #include "bn_mul.h"
43 #include "mbedtls/platform_util.h"
44 #include "mbedtls/error.h"
45 #include "constant_time_internal.h"
46 
47 #include <limits.h>
48 #include <string.h>
49 
50 #include "mbedtls/platform.h"
51 
52 #include <mempool.h>
53 #include <util.h>
54 
55 #define MPI_VALIDATE_RET(cond)                                       \
56     MBEDTLS_INTERNAL_VALIDATE_RET(cond, MBEDTLS_ERR_MPI_BAD_INPUT_DATA)
57 #define MPI_VALIDATE(cond)                                           \
58     MBEDTLS_INTERNAL_VALIDATE(cond)
59 
60 #define MPI_SIZE_T_MAX  ((size_t) -1)   /* SIZE_T_MAX is not standard */
61 
62 void *mbedtls_mpi_mempool;
63 
64 /* Implementation that should never be optimized out by the compiler */
65 static void mbedtls_mpi_zeroize(mbedtls_mpi_uint *v, size_t n)
66 {
67     mbedtls_platform_zeroize(v, ciL * n);
68 }
69 
70 /*
71  * Initialize one MPI
72  */
73 static void mpi_init(mbedtls_mpi *X, short use_mempool)
74 {
75     MPI_VALIDATE(X != NULL);
76 
77     X->s = 1;
78     X->use_mempool = use_mempool;
79     X->n = 0;
80     X->p = NULL;
81 }
82 
83 void mbedtls_mpi_init(mbedtls_mpi *X)
84 {
85     mpi_init(X, 0 /*use_mempool*/);
86 }
87 
88 void mbedtls_mpi_init_mempool(mbedtls_mpi *X)
89 {
90     mpi_init(X, !!mbedtls_mpi_mempool /*use_mempool*/);
91 }
92 
93 /*
94  * Unallocate one MPI
95  */
96 void mbedtls_mpi_free(mbedtls_mpi *X)
97 {
98     if (X == NULL) {
99         return;
100     }
101 
102     if (X->p != NULL) {
103         mbedtls_mpi_zeroize(X->p, X->n);
104         if(X->use_mempool)
105             mempool_free(mbedtls_mpi_mempool, X->p);
106         else
107         mbedtls_free(X->p);
108     }
109 
110     X->s = 1;
111     X->n = 0;
112     X->p = NULL;
113 }
114 
115 /*
116  * Enlarge to the specified number of limbs
117  */
118 int mbedtls_mpi_grow(mbedtls_mpi *X, size_t nblimbs)
119 {
120     mbedtls_mpi_uint *p;
121     MPI_VALIDATE_RET(X != NULL);
122 
123     if (nblimbs > MBEDTLS_MPI_MAX_LIMBS) {
124         return MBEDTLS_ERR_MPI_ALLOC_FAILED;
125     }
126 
127     if (X->n < nblimbs) {
128         if(X->use_mempool) {
129             p = mempool_alloc(mbedtls_mpi_mempool, nblimbs * ciL);
130             if(p == NULL)
131                 return MBEDTLS_ERR_MPI_ALLOC_FAILED;
132             memset(p, 0, nblimbs * ciL);
133         } else {
134                 p = (mbedtls_mpi_uint *) mbedtls_calloc(nblimbs, ciL);
135                 if (p == NULL)
136                     return MBEDTLS_ERR_MPI_ALLOC_FAILED;
137         }
138 
139         if (X->p != NULL) {
140             memcpy(p, X->p, X->n * ciL);
141             mbedtls_mpi_zeroize(X->p, X->n);
142             if (X->use_mempool)
143                 mempool_free(mbedtls_mpi_mempool, X->p);
144             else
145                 mbedtls_free(X->p);
146         }
147 
148         X->n = nblimbs;
149         X->p = p;
150     }
151 
152     return 0;
153 }
154 
155 /*
156  * Resize down as much as possible,
157  * while keeping at least the specified number of limbs
158  */
159 int mbedtls_mpi_shrink(mbedtls_mpi *X, size_t nblimbs)
160 {
161     mbedtls_mpi_uint *p;
162     size_t i;
163     MPI_VALIDATE_RET(X != NULL);
164 
165     if (nblimbs > MBEDTLS_MPI_MAX_LIMBS) {
166         return MBEDTLS_ERR_MPI_ALLOC_FAILED;
167     }
168 
169     /* Actually resize up if there are currently fewer than nblimbs limbs. */
170     if (X->n <= nblimbs) {
171         return mbedtls_mpi_grow(X, nblimbs);
172     }
173     /* After this point, then X->n > nblimbs and in particular X->n > 0. */
174 
175     for (i = X->n - 1; i > 0; i--) {
176         if (X->p[i] != 0) {
177             break;
178         }
179     }
180     i++;
181 
182     if (i < nblimbs) {
183         i = nblimbs;
184     }
185 
186     if (X->use_mempool) {
187         p = mempool_alloc(mbedtls_mpi_mempool, i * ciL);
188         if (p == NULL)
189             return MBEDTLS_ERR_MPI_ALLOC_FAILED;
190         memset(p, 0, i * ciL);
191     } else {
192         if ((p = (mbedtls_mpi_uint *) mbedtls_calloc(i, ciL)) == NULL)
193             return MBEDTLS_ERR_MPI_ALLOC_FAILED;
194     }
195 
196     if (X->p != NULL) {
197         memcpy(p, X->p, i * ciL);
198         mbedtls_mpi_zeroize(X->p, X->n);
199         if (X->use_mempool)
200             mempool_free(mbedtls_mpi_mempool, X->p);
201         else
202             mbedtls_free(X->p);
203     }
204 
205     X->n = i;
206     X->p = p;
207 
208     return 0;
209 }
210 
211 /* Resize X to have exactly n limbs and set it to 0. */
212 static int mbedtls_mpi_resize_clear(mbedtls_mpi *X, size_t limbs)
213 {
214     if (limbs == 0) {
215         mbedtls_mpi_free(X);
216         return 0;
217     } else if (X->n == limbs) {
218         memset(X->p, 0, limbs * ciL);
219         X->s = 1;
220         return 0;
221     } else {
222         mbedtls_mpi_free(X);
223         return mbedtls_mpi_grow(X, limbs);
224     }
225 }
226 
227 /*
228  * Copy the contents of Y into X.
229  *
230  * This function is not constant-time. Leading zeros in Y may be removed.
231  *
232  * Ensure that X does not shrink. This is not guaranteed by the public API,
233  * but some code in the bignum module relies on this property, for example
234  * in mbedtls_mpi_exp_mod().
235  */
236 int mbedtls_mpi_copy(mbedtls_mpi *X, const mbedtls_mpi *Y)
237 {
238     int ret = 0;
239     size_t i;
240     MPI_VALIDATE_RET(X != NULL);
241     MPI_VALIDATE_RET(Y != NULL);
242 
243     if (X == Y) {
244         return 0;
245     }
246 
247     if (Y->n == 0) {
248         if (X->n != 0) {
249             X->s = 1;
250             memset(X->p, 0, X->n * ciL);
251         }
252         return 0;
253     }
254 
255     for (i = Y->n - 1; i > 0; i--) {
256         if (Y->p[i] != 0) {
257             break;
258         }
259     }
260     i++;
261 
262     X->s = Y->s;
263 
264     if (X->n < i) {
265         MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, i));
266     } else {
267         memset(X->p + i, 0, (X->n - i) * ciL);
268     }
269 
270     memcpy(X->p, Y->p, i * ciL);
271 
272 cleanup:
273 
274     return ret;
275 }
276 
277 /*
278  * Swap the contents of X and Y
279  */
280 void mbedtls_mpi_swap(mbedtls_mpi *X, mbedtls_mpi *Y)
281 {
282     mbedtls_mpi T;
283     MPI_VALIDATE(X != NULL);
284     MPI_VALIDATE(Y != NULL);
285 
286     memcpy(&T,  X, sizeof(mbedtls_mpi));
287     memcpy(X,  Y, sizeof(mbedtls_mpi));
288     memcpy(Y, &T, sizeof(mbedtls_mpi));
289 }
290 
291 static inline mbedtls_mpi_uint mpi_sint_abs(mbedtls_mpi_sint z)
292 {
293     if (z >= 0) {
294         return z;
295     }
296     /* Take care to handle the most negative value (-2^(biL-1)) correctly.
297      * A naive -z would have undefined behavior.
298      * Write this in a way that makes popular compilers happy (GCC, Clang,
299      * MSVC). */
300     return (mbedtls_mpi_uint) 0 - (mbedtls_mpi_uint) z;
301 }
302 
303 /*
304  * Set value from integer
305  */
306 int mbedtls_mpi_lset(mbedtls_mpi *X, mbedtls_mpi_sint z)
307 {
308     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
309     MPI_VALIDATE_RET(X != NULL);
310 
311     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, 1));
312     memset(X->p, 0, X->n * ciL);
313 
314     X->p[0] = mpi_sint_abs(z);
315     X->s    = (z < 0) ? -1 : 1;
316 
317 cleanup:
318 
319     return ret;
320 }
321 
322 /*
323  * Get a specific bit
324  */
325 int mbedtls_mpi_get_bit(const mbedtls_mpi *X, size_t pos)
326 {
327     MPI_VALIDATE_RET(X != NULL);
328 
329     if (X->n * biL <= pos) {
330         return 0;
331     }
332 
333     return (X->p[pos / biL] >> (pos % biL)) & 0x01;
334 }
335 
336 /*
337  * Set a bit to a specific value of 0 or 1
338  */
339 int mbedtls_mpi_set_bit(mbedtls_mpi *X, size_t pos, unsigned char val)
340 {
341     int ret = 0;
342     size_t off = pos / biL;
343     size_t idx = pos % biL;
344     MPI_VALIDATE_RET(X != NULL);
345 
346     if (val != 0 && val != 1) {
347         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
348     }
349 
350     if (X->n * biL <= pos) {
351         if (val == 0) {
352             return 0;
353         }
354 
355         MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, off + 1));
356     }
357 
358     X->p[off] &= ~((mbedtls_mpi_uint) 0x01 << idx);
359     X->p[off] |= (mbedtls_mpi_uint) val << idx;
360 
361 cleanup:
362 
363     return ret;
364 }
365 
366 /*
367  * Return the number of less significant zero-bits
368  */
369 size_t mbedtls_mpi_lsb(const mbedtls_mpi *X)
370 {
371     size_t i, j, count = 0;
372     MBEDTLS_INTERNAL_VALIDATE_RET(X != NULL, 0);
373 
374     for (i = 0; i < X->n; i++) {
375         for (j = 0; j < biL; j++, count++) {
376             if (((X->p[i] >> j) & 1) != 0) {
377                 return count;
378             }
379         }
380     }
381 
382     return 0;
383 }
384 
385 /*
386  * Return the number of bits
387  */
388 size_t mbedtls_mpi_bitlen(const mbedtls_mpi *X)
389 {
390     return mbedtls_mpi_core_bitlen(X->p, X->n);
391 }
392 
393 /*
394  * Return the total size in bytes
395  */
396 size_t mbedtls_mpi_size(const mbedtls_mpi *X)
397 {
398     return (mbedtls_mpi_bitlen(X) + 7) >> 3;
399 }
400 
401 /*
402  * Convert an ASCII character to digit value
403  */
404 static int mpi_get_digit(mbedtls_mpi_uint *d, int radix, char c)
405 {
406     *d = 255;
407 
408     if (c >= 0x30 && c <= 0x39) {
409         *d = c - 0x30;
410     }
411     if (c >= 0x41 && c <= 0x46) {
412         *d = c - 0x37;
413     }
414     if (c >= 0x61 && c <= 0x66) {
415         *d = c - 0x57;
416     }
417 
418     if (*d >= (mbedtls_mpi_uint) radix) {
419         return MBEDTLS_ERR_MPI_INVALID_CHARACTER;
420     }
421 
422     return 0;
423 }
424 
425 /*
426  * Import from an ASCII string
427  */
428 int mbedtls_mpi_read_string(mbedtls_mpi *X, int radix, const char *s)
429 {
430     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
431     size_t i, j, slen, n;
432     int sign = 1;
433     mbedtls_mpi_uint d;
434     mbedtls_mpi T;
435     MPI_VALIDATE_RET(X != NULL);
436     MPI_VALIDATE_RET(s != NULL);
437 
438     if (radix < 2 || radix > 16) {
439         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
440     }
441 
442     mbedtls_mpi_init_mempool(&T);
443 
444     if (s[0] == 0) {
445         mbedtls_mpi_free(X);
446         return 0;
447     }
448 
449     if (s[0] == '-') {
450         ++s;
451         sign = -1;
452     }
453 
454     slen = strlen(s);
455 
456     if (radix == 16) {
457         if (slen > MPI_SIZE_T_MAX >> 2) {
458             return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
459         }
460 
461         n = BITS_TO_LIMBS(slen << 2);
462 
463         MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, n));
464         MBEDTLS_MPI_CHK(mbedtls_mpi_lset(X, 0));
465 
466         for (i = slen, j = 0; i > 0; i--, j++) {
467             MBEDTLS_MPI_CHK(mpi_get_digit(&d, radix, s[i - 1]));
468             X->p[j / (2 * ciL)] |= d << ((j % (2 * ciL)) << 2);
469         }
470     } else {
471         MBEDTLS_MPI_CHK(mbedtls_mpi_lset(X, 0));
472 
473         for (i = 0; i < slen; i++) {
474             MBEDTLS_MPI_CHK(mpi_get_digit(&d, radix, s[i]));
475             MBEDTLS_MPI_CHK(mbedtls_mpi_mul_int(&T, X, radix));
476             MBEDTLS_MPI_CHK(mbedtls_mpi_add_int(X, &T, d));
477         }
478     }
479 
480     if (sign < 0 && mbedtls_mpi_bitlen(X) != 0) {
481         X->s = -1;
482     }
483 
484 cleanup:
485 
486     mbedtls_mpi_free(&T);
487 
488     return ret;
489 }
490 
491 /*
492  * Helper to write the digits high-order first.
493  */
494 static int mpi_write_hlp(mbedtls_mpi *X, int radix,
495                          char **p, const size_t buflen)
496 {
497     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
498     mbedtls_mpi_uint r;
499     size_t length = 0;
500     char *p_end = *p + buflen;
501 
502     do {
503         if (length >= buflen) {
504             return MBEDTLS_ERR_MPI_BUFFER_TOO_SMALL;
505         }
506 
507         MBEDTLS_MPI_CHK(mbedtls_mpi_mod_int(&r, X, radix));
508         MBEDTLS_MPI_CHK(mbedtls_mpi_div_int(X, NULL, X, radix));
509         /*
510          * Write the residue in the current position, as an ASCII character.
511          */
512         if (r < 0xA) {
513             *(--p_end) = (char) ('0' + r);
514         } else {
515             *(--p_end) = (char) ('A' + (r - 0xA));
516         }
517 
518         length++;
519     } while (mbedtls_mpi_cmp_int(X, 0) != 0);
520 
521     memmove(*p, p_end, length);
522     *p += length;
523 
524 cleanup:
525 
526     return ret;
527 }
528 
529 /*
530  * Export into an ASCII string
531  */
532 int mbedtls_mpi_write_string(const mbedtls_mpi *X, int radix,
533                              char *buf, size_t buflen, size_t *olen)
534 {
535     int ret = 0;
536     size_t n;
537     char *p;
538     mbedtls_mpi T;
539     MPI_VALIDATE_RET(X    != NULL);
540     MPI_VALIDATE_RET(olen != NULL);
541     MPI_VALIDATE_RET(buflen == 0 || buf != NULL);
542 
543     if (radix < 2 || radix > 16) {
544         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
545     }
546 
547     n = mbedtls_mpi_bitlen(X);   /* Number of bits necessary to present `n`. */
548     if (radix >=  4) {
549         n >>= 1;                 /* Number of 4-adic digits necessary to present
550                                   * `n`. If radix > 4, this might be a strict
551                                   * overapproximation of the number of
552                                   * radix-adic digits needed to present `n`. */
553     }
554     if (radix >= 16) {
555         n >>= 1;                 /* Number of hexadecimal digits necessary to
556                                   * present `n`. */
557 
558     }
559     n += 1; /* Terminating null byte */
560     n += 1; /* Compensate for the divisions above, which round down `n`
561              * in case it's not even. */
562     n += 1; /* Potential '-'-sign. */
563     n += (n & 1);   /* Make n even to have enough space for hexadecimal writing,
564                      * which always uses an even number of hex-digits. */
565 
566     if (buflen < n) {
567         *olen = n;
568         return MBEDTLS_ERR_MPI_BUFFER_TOO_SMALL;
569     }
570 
571     p = buf;
572     mbedtls_mpi_init_mempool(&T);
573 
574     if (X->s == -1) {
575         *p++ = '-';
576         buflen--;
577     }
578 
579     if (radix == 16) {
580         int c;
581         size_t i, j, k;
582 
583         for (i = X->n, k = 0; i > 0; i--) {
584             for (j = ciL; j > 0; j--) {
585                 c = (X->p[i - 1] >> ((j - 1) << 3)) & 0xFF;
586 
587                 if (c == 0 && k == 0 && (i + j) != 2) {
588                     continue;
589                 }
590 
591                 *(p++) = "0123456789ABCDEF" [c / 16];
592                 *(p++) = "0123456789ABCDEF" [c % 16];
593                 k = 1;
594             }
595         }
596     } else {
597         MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&T, X));
598 
599         if (T.s == -1) {
600             T.s = 1;
601         }
602 
603         MBEDTLS_MPI_CHK(mpi_write_hlp(&T, radix, &p, buflen));
604     }
605 
606     *p++ = '\0';
607     *olen = p - buf;
608 
609 cleanup:
610 
611     mbedtls_mpi_free(&T);
612 
613     return ret;
614 }
615 
616 #if defined(MBEDTLS_FS_IO)
617 /*
618  * Read X from an opened file
619  */
620 int mbedtls_mpi_read_file(mbedtls_mpi *X, int radix, FILE *fin)
621 {
622     mbedtls_mpi_uint d;
623     size_t slen;
624     char *p;
625     /*
626      * Buffer should have space for (short) label and decimal formatted MPI,
627      * newline characters and '\0'
628      */
629     char s[MBEDTLS_MPI_RW_BUFFER_SIZE];
630 
631     MPI_VALIDATE_RET(X   != NULL);
632     MPI_VALIDATE_RET(fin != NULL);
633 
634     if (radix < 2 || radix > 16) {
635         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
636     }
637 
638     memset(s, 0, sizeof(s));
639     if (fgets(s, sizeof(s) - 1, fin) == NULL) {
640         return MBEDTLS_ERR_MPI_FILE_IO_ERROR;
641     }
642 
643     slen = strlen(s);
644     if (slen == sizeof(s) - 2) {
645         return MBEDTLS_ERR_MPI_BUFFER_TOO_SMALL;
646     }
647 
648     if (slen > 0 && s[slen - 1] == '\n') {
649         slen--; s[slen] = '\0';
650     }
651     if (slen > 0 && s[slen - 1] == '\r') {
652         slen--; s[slen] = '\0';
653     }
654 
655     p = s + slen;
656     while (p-- > s) {
657         if (mpi_get_digit(&d, radix, *p) != 0) {
658             break;
659         }
660     }
661 
662     return mbedtls_mpi_read_string(X, radix, p + 1);
663 }
664 
665 /*
666  * Write X into an opened file (or stdout if fout == NULL)
667  */
668 int mbedtls_mpi_write_file(const char *p, const mbedtls_mpi *X, int radix, FILE *fout)
669 {
670     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
671     size_t n, slen, plen;
672     /*
673      * Buffer should have space for (short) label and decimal formatted MPI,
674      * newline characters and '\0'
675      */
676     char s[MBEDTLS_MPI_RW_BUFFER_SIZE];
677     MPI_VALIDATE_RET(X != NULL);
678 
679     if (radix < 2 || radix > 16) {
680         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
681     }
682 
683     memset(s, 0, sizeof(s));
684 
685     MBEDTLS_MPI_CHK(mbedtls_mpi_write_string(X, radix, s, sizeof(s) - 2, &n));
686 
687     if (p == NULL) {
688         p = "";
689     }
690 
691     plen = strlen(p);
692     slen = strlen(s);
693     s[slen++] = '\r';
694     s[slen++] = '\n';
695 
696     if (fout != NULL) {
697         if (fwrite(p, 1, plen, fout) != plen ||
698             fwrite(s, 1, slen, fout) != slen) {
699             return MBEDTLS_ERR_MPI_FILE_IO_ERROR;
700         }
701     } else {
702         mbedtls_printf("%s%s", p, s);
703     }
704 
705 cleanup:
706 
707     return ret;
708 }
709 #endif /* MBEDTLS_FS_IO */
710 
711 /*
712  * Import X from unsigned binary data, little endian
713  *
714  * This function is guaranteed to return an MPI with exactly the necessary
715  * number of limbs (in particular, it does not skip 0s in the input).
716  */
717 int mbedtls_mpi_read_binary_le(mbedtls_mpi *X,
718                                const unsigned char *buf, size_t buflen)
719 {
720     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
721     const size_t limbs = CHARS_TO_LIMBS(buflen);
722 
723     /* Ensure that target MPI has exactly the necessary number of limbs */
724     MBEDTLS_MPI_CHK(mbedtls_mpi_resize_clear(X, limbs));
725 
726     MBEDTLS_MPI_CHK(mbedtls_mpi_core_read_le(X->p, X->n, buf, buflen));
727 
728 cleanup:
729 
730     /*
731      * This function is also used to import keys. However, wiping the buffers
732      * upon failure is not necessary because failure only can happen before any
733      * input is copied.
734      */
735     return ret;
736 }
737 
738 /*
739  * Import X from unsigned binary data, big endian
740  *
741  * This function is guaranteed to return an MPI with exactly the necessary
742  * number of limbs (in particular, it does not skip 0s in the input).
743  */
744 int mbedtls_mpi_read_binary(mbedtls_mpi *X, const unsigned char *buf, size_t buflen)
745 {
746     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
747     const size_t limbs = CHARS_TO_LIMBS(buflen);
748 
749     MPI_VALIDATE_RET(X != NULL);
750     MPI_VALIDATE_RET(buflen == 0 || buf != NULL);
751 
752     /* Ensure that target MPI has exactly the necessary number of limbs */
753     MBEDTLS_MPI_CHK(mbedtls_mpi_resize_clear(X, limbs));
754 
755     MBEDTLS_MPI_CHK(mbedtls_mpi_core_read_be(X->p, X->n, buf, buflen));
756 
757 cleanup:
758 
759     /*
760      * This function is also used to import keys. However, wiping the buffers
761      * upon failure is not necessary because failure only can happen before any
762      * input is copied.
763      */
764     return ret;
765 }
766 
767 /*
768  * Export X into unsigned binary data, little endian
769  */
770 int mbedtls_mpi_write_binary_le(const mbedtls_mpi *X,
771                                 unsigned char *buf, size_t buflen)
772 {
773     return mbedtls_mpi_core_write_le(X->p, X->n, buf, buflen);
774 }
775 
776 /*
777  * Export X into unsigned binary data, big endian
778  */
779 int mbedtls_mpi_write_binary(const mbedtls_mpi *X,
780                              unsigned char *buf, size_t buflen)
781 {
782     return mbedtls_mpi_core_write_be(X->p, X->n, buf, buflen);
783 }
784 
785 /*
786  * Left-shift: X <<= count
787  */
788 int mbedtls_mpi_shift_l(mbedtls_mpi *X, size_t count)
789 {
790     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
791     size_t i, v0, t1;
792     mbedtls_mpi_uint r0 = 0, r1;
793     MPI_VALIDATE_RET(X != NULL);
794 
795     v0 = count / (biL);
796     t1 = count & (biL - 1);
797 
798     i = mbedtls_mpi_bitlen(X) + count;
799 
800     if (X->n * biL < i) {
801         MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, BITS_TO_LIMBS(i)));
802     }
803 
804     ret = 0;
805 
806     /*
807      * shift by count / limb_size
808      */
809     if (v0 > 0) {
810         for (i = X->n; i > v0; i--) {
811             X->p[i - 1] = X->p[i - v0 - 1];
812         }
813 
814         for (; i > 0; i--) {
815             X->p[i - 1] = 0;
816         }
817     }
818 
819     /*
820      * shift by count % limb_size
821      */
822     if (t1 > 0) {
823         for (i = v0; i < X->n; i++) {
824             r1 = X->p[i] >> (biL - t1);
825             X->p[i] <<= t1;
826             X->p[i] |= r0;
827             r0 = r1;
828         }
829     }
830 
831 cleanup:
832 
833     return ret;
834 }
835 
836 /*
837  * Right-shift: X >>= count
838  */
839 int mbedtls_mpi_shift_r(mbedtls_mpi *X, size_t count)
840 {
841     MPI_VALIDATE_RET(X != NULL);
842     if (X->n != 0) {
843         mbedtls_mpi_core_shift_r(X->p, X->n, count);
844     }
845     return 0;
846 }
847 
848 /*
849  * Compare unsigned values
850  */
851 int mbedtls_mpi_cmp_abs(const mbedtls_mpi *X, const mbedtls_mpi *Y)
852 {
853     size_t i, j;
854     MPI_VALIDATE_RET(X != NULL);
855     MPI_VALIDATE_RET(Y != NULL);
856 
857     for (i = X->n; i > 0; i--) {
858         if (X->p[i - 1] != 0) {
859             break;
860         }
861     }
862 
863     for (j = Y->n; j > 0; j--) {
864         if (Y->p[j - 1] != 0) {
865             break;
866         }
867     }
868 
869     if (i == 0 && j == 0) {
870         return 0;
871     }
872 
873     if (i > j) {
874         return 1;
875     }
876     if (j > i) {
877         return -1;
878     }
879 
880     for (; i > 0; i--) {
881         if (X->p[i - 1] > Y->p[i - 1]) {
882             return 1;
883         }
884         if (X->p[i - 1] < Y->p[i - 1]) {
885             return -1;
886         }
887     }
888 
889     return 0;
890 }
891 
892 /*
893  * Compare signed values
894  */
895 int mbedtls_mpi_cmp_mpi(const mbedtls_mpi *X, const mbedtls_mpi *Y)
896 {
897     size_t i, j;
898     MPI_VALIDATE_RET(X != NULL);
899     MPI_VALIDATE_RET(Y != NULL);
900 
901     for (i = X->n; i > 0; i--) {
902         if (X->p[i - 1] != 0) {
903             break;
904         }
905     }
906 
907     for (j = Y->n; j > 0; j--) {
908         if (Y->p[j - 1] != 0) {
909             break;
910         }
911     }
912 
913     if (i == 0 && j == 0) {
914         return 0;
915     }
916 
917     if (i > j) {
918         return X->s;
919     }
920     if (j > i) {
921         return -Y->s;
922     }
923 
924     if (X->s > 0 && Y->s < 0) {
925         return 1;
926     }
927     if (Y->s > 0 && X->s < 0) {
928         return -1;
929     }
930 
931     for (; i > 0; i--) {
932         if (X->p[i - 1] > Y->p[i - 1]) {
933             return X->s;
934         }
935         if (X->p[i - 1] < Y->p[i - 1]) {
936             return -X->s;
937         }
938     }
939 
940     return 0;
941 }
942 
943 /*
944  * Compare signed values
945  */
946 int mbedtls_mpi_cmp_int(const mbedtls_mpi *X, mbedtls_mpi_sint z)
947 {
948     mbedtls_mpi Y;
949     mbedtls_mpi_uint p[1];
950     MPI_VALIDATE_RET(X != NULL);
951 
952     *p  = mpi_sint_abs(z);
953     Y.s = (z < 0) ? -1 : 1;
954     Y.n = 1;
955     Y.p = p;
956 
957     return mbedtls_mpi_cmp_mpi(X, &Y);
958 }
959 
960 /*
961  * Unsigned addition: X = |A| + |B|  (HAC 14.7)
962  */
963 int mbedtls_mpi_add_abs(mbedtls_mpi *X, const mbedtls_mpi *A, const mbedtls_mpi *B)
964 {
965     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
966     size_t j;
967     MPI_VALIDATE_RET(X != NULL);
968     MPI_VALIDATE_RET(A != NULL);
969     MPI_VALIDATE_RET(B != NULL);
970 
971     if (X == B) {
972         const mbedtls_mpi *T = A; A = X; B = T;
973     }
974 
975     if (X != A) {
976         MBEDTLS_MPI_CHK(mbedtls_mpi_copy(X, A));
977     }
978 
979     /*
980      * X must always be positive as a result of unsigned additions.
981      */
982     X->s = 1;
983 
984     for (j = B->n; j > 0; j--) {
985         if (B->p[j - 1] != 0) {
986             break;
987         }
988     }
989 
990     /* Exit early to avoid undefined behavior on NULL+0 when X->n == 0
991      * and B is 0 (of any size). */
992     if (j == 0) {
993         return 0;
994     }
995 
996     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, j));
997 
998     /* j is the number of non-zero limbs of B. Add those to X. */
999 
1000     mbedtls_mpi_uint *p = X->p;
1001 
1002     mbedtls_mpi_uint c = mbedtls_mpi_core_add(p, p, B->p, j);
1003 
1004     p += j;
1005 
1006     /* Now propagate any carry */
1007 
1008     while (c != 0) {
1009         if (j >= X->n) {
1010             MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, j + 1));
1011             p = X->p + j;
1012         }
1013 
1014         *p += c; c = (*p < c); j++; p++;
1015     }
1016 
1017 cleanup:
1018 
1019     return ret;
1020 }
1021 
1022 /*
1023  * Unsigned subtraction: X = |A| - |B|  (HAC 14.9, 14.10)
1024  */
1025 int mbedtls_mpi_sub_abs(mbedtls_mpi *X, const mbedtls_mpi *A, const mbedtls_mpi *B)
1026 {
1027     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
1028     size_t n;
1029     mbedtls_mpi_uint carry;
1030     MPI_VALIDATE_RET(X != NULL);
1031     MPI_VALIDATE_RET(A != NULL);
1032     MPI_VALIDATE_RET(B != NULL);
1033 
1034     for (n = B->n; n > 0; n--) {
1035         if (B->p[n - 1] != 0) {
1036             break;
1037         }
1038     }
1039     if (n > A->n) {
1040         /* B >= (2^ciL)^n > A */
1041         ret = MBEDTLS_ERR_MPI_NEGATIVE_VALUE;
1042         goto cleanup;
1043     }
1044 
1045     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, A->n));
1046 
1047     /* Set the high limbs of X to match A. Don't touch the lower limbs
1048      * because X might be aliased to B, and we must not overwrite the
1049      * significant digits of B. */
1050     if (A->n > n && A != X) {
1051         memcpy(X->p + n, A->p + n, (A->n - n) * ciL);
1052     }
1053     if (X->n > A->n) {
1054         memset(X->p + A->n, 0, (X->n - A->n) * ciL);
1055     }
1056 
1057     carry = mbedtls_mpi_core_sub(X->p, A->p, B->p, n);
1058     if (carry != 0) {
1059         /* Propagate the carry through the rest of X. */
1060         carry = mbedtls_mpi_core_sub_int(X->p + n, X->p + n, carry, X->n - n);
1061 
1062         /* If we have further carry/borrow, the result is negative. */
1063         if (carry != 0) {
1064             ret = MBEDTLS_ERR_MPI_NEGATIVE_VALUE;
1065             goto cleanup;
1066         }
1067     }
1068 
1069     /* X should always be positive as a result of unsigned subtractions. */
1070     X->s = 1;
1071 
1072 cleanup:
1073     return ret;
1074 }
1075 
1076 /* Common function for signed addition and subtraction.
1077  * Calculate A + B * flip_B where flip_B is 1 or -1.
1078  */
1079 static int add_sub_mpi(mbedtls_mpi *X,
1080                        const mbedtls_mpi *A, const mbedtls_mpi *B,
1081                        int flip_B)
1082 {
1083     int ret, s;
1084     MPI_VALIDATE_RET(X != NULL);
1085     MPI_VALIDATE_RET(A != NULL);
1086     MPI_VALIDATE_RET(B != NULL);
1087 
1088     s = A->s;
1089     if (A->s * B->s * flip_B < 0) {
1090         int cmp = mbedtls_mpi_cmp_abs(A, B);
1091         if (cmp >= 0) {
1092             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_abs(X, A, B));
1093             /* If |A| = |B|, the result is 0 and we must set the sign bit
1094              * to +1 regardless of which of A or B was negative. Otherwise,
1095              * since |A| > |B|, the sign is the sign of A. */
1096             X->s = cmp == 0 ? 1 : s;
1097         } else {
1098             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_abs(X, B, A));
1099             /* Since |A| < |B|, the sign is the opposite of A. */
1100             X->s = -s;
1101         }
1102     } else {
1103         MBEDTLS_MPI_CHK(mbedtls_mpi_add_abs(X, A, B));
1104         X->s = s;
1105     }
1106 
1107 cleanup:
1108 
1109     return ret;
1110 }
1111 
1112 /*
1113  * Signed addition: X = A + B
1114  */
1115 int mbedtls_mpi_add_mpi(mbedtls_mpi *X, const mbedtls_mpi *A, const mbedtls_mpi *B)
1116 {
1117     return add_sub_mpi(X, A, B, 1);
1118 }
1119 
1120 /*
1121  * Signed subtraction: X = A - B
1122  */
1123 int mbedtls_mpi_sub_mpi(mbedtls_mpi *X, const mbedtls_mpi *A, const mbedtls_mpi *B)
1124 {
1125     return add_sub_mpi(X, A, B, -1);
1126 }
1127 
1128 /*
1129  * Signed addition: X = A + b
1130  */
1131 int mbedtls_mpi_add_int(mbedtls_mpi *X, const mbedtls_mpi *A, mbedtls_mpi_sint b)
1132 {
1133     mbedtls_mpi B;
1134     mbedtls_mpi_uint p[1];
1135     MPI_VALIDATE_RET(X != NULL);
1136     MPI_VALIDATE_RET(A != NULL);
1137 
1138     p[0] = mpi_sint_abs(b);
1139     B.s = (b < 0) ? -1 : 1;
1140     B.n = 1;
1141     B.p = p;
1142 
1143     return mbedtls_mpi_add_mpi(X, A, &B);
1144 }
1145 
1146 /*
1147  * Signed subtraction: X = A - b
1148  */
1149 int mbedtls_mpi_sub_int(mbedtls_mpi *X, const mbedtls_mpi *A, mbedtls_mpi_sint b)
1150 {
1151     mbedtls_mpi B;
1152     mbedtls_mpi_uint p[1];
1153     MPI_VALIDATE_RET(X != NULL);
1154     MPI_VALIDATE_RET(A != NULL);
1155 
1156     p[0] = mpi_sint_abs(b);
1157     B.s = (b < 0) ? -1 : 1;
1158     B.n = 1;
1159     B.p = p;
1160 
1161     return mbedtls_mpi_sub_mpi(X, A, &B);
1162 }
1163 
1164 /*
1165  * Baseline multiplication: X = A * B  (HAC 14.12)
1166  */
1167 int mbedtls_mpi_mul_mpi(mbedtls_mpi *X, const mbedtls_mpi *A, const mbedtls_mpi *B)
1168 {
1169     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
1170     size_t i, j;
1171     mbedtls_mpi TA, TB;
1172     int result_is_zero = 0;
1173     MPI_VALIDATE_RET(X != NULL);
1174     MPI_VALIDATE_RET(A != NULL);
1175     MPI_VALIDATE_RET(B != NULL);
1176 
1177     mbedtls_mpi_init_mempool(&TA); mbedtls_mpi_init_mempool(&TB);
1178 
1179     if (X == A) {
1180         MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&TA, A)); A = &TA;
1181     }
1182     if (X == B) {
1183         MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&TB, B)); B = &TB;
1184     }
1185 
1186     for (i = A->n; i > 0; i--) {
1187         if (A->p[i - 1] != 0) {
1188             break;
1189         }
1190     }
1191     if (i == 0) {
1192         result_is_zero = 1;
1193     }
1194 
1195     for (j = B->n; j > 0; j--) {
1196         if (B->p[j - 1] != 0) {
1197             break;
1198         }
1199     }
1200     if (j == 0) {
1201         result_is_zero = 1;
1202     }
1203 
1204     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, i + j));
1205     MBEDTLS_MPI_CHK(mbedtls_mpi_lset(X, 0));
1206 
1207     for (size_t k = 0; k < j; k++) {
1208         /* We know that there cannot be any carry-out since we're
1209          * iterating from bottom to top. */
1210         (void) mbedtls_mpi_core_mla(X->p + k, i + 1,
1211                                     A->p, i,
1212                                     B->p[k]);
1213     }
1214 
1215     /* If the result is 0, we don't shortcut the operation, which reduces
1216      * but does not eliminate side channels leaking the zero-ness. We do
1217      * need to take care to set the sign bit properly since the library does
1218      * not fully support an MPI object with a value of 0 and s == -1. */
1219     if (result_is_zero) {
1220         X->s = 1;
1221     } else {
1222         X->s = A->s * B->s;
1223     }
1224 
1225 cleanup:
1226 
1227     mbedtls_mpi_free(&TB); mbedtls_mpi_free(&TA);
1228 
1229     return ret;
1230 }
1231 
1232 /*
1233  * Baseline multiplication: X = A * b
1234  */
1235 int mbedtls_mpi_mul_int(mbedtls_mpi *X, const mbedtls_mpi *A, mbedtls_mpi_uint b)
1236 {
1237     MPI_VALIDATE_RET(X != NULL);
1238     MPI_VALIDATE_RET(A != NULL);
1239 
1240     size_t n = A->n;
1241     while (n > 0 && A->p[n - 1] == 0) {
1242         --n;
1243     }
1244 
1245     /* The general method below doesn't work if b==0. */
1246     if (b == 0 || n == 0) {
1247         return mbedtls_mpi_lset(X, 0);
1248     }
1249 
1250     /* Calculate A*b as A + A*(b-1) to take advantage of mbedtls_mpi_core_mla */
1251     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
1252     /* In general, A * b requires 1 limb more than b. If
1253      * A->p[n - 1] * b / b == A->p[n - 1], then A * b fits in the same
1254      * number of limbs as A and the call to grow() is not required since
1255      * copy() will take care of the growth if needed. However, experimentally,
1256      * making the call to grow() unconditional causes slightly fewer
1257      * calls to calloc() in ECP code, presumably because it reuses the
1258      * same mpi for a while and this way the mpi is more likely to directly
1259      * grow to its final size.
1260      *
1261      * Note that calculating A*b as 0 + A*b doesn't work as-is because
1262      * A,X can be the same. */
1263     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(X, n + 1));
1264     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(X, A));
1265     mbedtls_mpi_core_mla(X->p, X->n, A->p, n, b - 1);
1266 
1267 cleanup:
1268     return ret;
1269 }
1270 
1271 /*
1272  * Unsigned integer divide - double mbedtls_mpi_uint dividend, u1/u0, and
1273  * mbedtls_mpi_uint divisor, d
1274  */
1275 static mbedtls_mpi_uint mbedtls_int_div_int(mbedtls_mpi_uint u1,
1276                                             mbedtls_mpi_uint u0,
1277                                             mbedtls_mpi_uint d,
1278                                             mbedtls_mpi_uint *r)
1279 {
1280 #if defined(MBEDTLS_HAVE_UDBL)
1281     mbedtls_t_udbl dividend, quotient;
1282 #else
1283     const mbedtls_mpi_uint radix = (mbedtls_mpi_uint) 1 << biH;
1284     const mbedtls_mpi_uint uint_halfword_mask = ((mbedtls_mpi_uint) 1 << biH) - 1;
1285     mbedtls_mpi_uint d0, d1, q0, q1, rAX, r0, quotient;
1286     mbedtls_mpi_uint u0_msw, u0_lsw;
1287     size_t s;
1288 #endif
1289 
1290     /*
1291      * Check for overflow
1292      */
1293     if (0 == d || u1 >= d) {
1294         if (r != NULL) {
1295             *r = ~(mbedtls_mpi_uint) 0u;
1296         }
1297 
1298         return ~(mbedtls_mpi_uint) 0u;
1299     }
1300 
1301 #if defined(MBEDTLS_HAVE_UDBL)
1302     dividend  = (mbedtls_t_udbl) u1 << biL;
1303     dividend |= (mbedtls_t_udbl) u0;
1304     quotient = dividend / d;
1305     if (quotient > ((mbedtls_t_udbl) 1 << biL) - 1) {
1306         quotient = ((mbedtls_t_udbl) 1 << biL) - 1;
1307     }
1308 
1309     if (r != NULL) {
1310         *r = (mbedtls_mpi_uint) (dividend - (quotient * d));
1311     }
1312 
1313     return (mbedtls_mpi_uint) quotient;
1314 #else
1315 
1316     /*
1317      * Algorithm D, Section 4.3.1 - The Art of Computer Programming
1318      *   Vol. 2 - Seminumerical Algorithms, Knuth
1319      */
1320 
1321     /*
1322      * Normalize the divisor, d, and dividend, u0, u1
1323      */
1324     s = mbedtls_mpi_core_clz(d);
1325     d = d << s;
1326 
1327     u1 = u1 << s;
1328     u1 |= (u0 >> (biL - s)) & (-(mbedtls_mpi_sint) s >> (biL - 1));
1329     u0 =  u0 << s;
1330 
1331     d1 = d >> biH;
1332     d0 = d & uint_halfword_mask;
1333 
1334     u0_msw = u0 >> biH;
1335     u0_lsw = u0 & uint_halfword_mask;
1336 
1337     /*
1338      * Find the first quotient and remainder
1339      */
1340     q1 = u1 / d1;
1341     r0 = u1 - d1 * q1;
1342 
1343     while (q1 >= radix || (q1 * d0 > radix * r0 + u0_msw)) {
1344         q1 -= 1;
1345         r0 += d1;
1346 
1347         if (r0 >= radix) {
1348             break;
1349         }
1350     }
1351 
1352     rAX = (u1 * radix) + (u0_msw - q1 * d);
1353     q0 = rAX / d1;
1354     r0 = rAX - q0 * d1;
1355 
1356     while (q0 >= radix || (q0 * d0 > radix * r0 + u0_lsw)) {
1357         q0 -= 1;
1358         r0 += d1;
1359 
1360         if (r0 >= radix) {
1361             break;
1362         }
1363     }
1364 
1365     if (r != NULL) {
1366         *r = (rAX * radix + u0_lsw - q0 * d) >> s;
1367     }
1368 
1369     quotient = q1 * radix + q0;
1370 
1371     return quotient;
1372 #endif
1373 }
1374 
1375 /*
1376  * Division by mbedtls_mpi: A = Q * B + R  (HAC 14.20)
1377  */
1378 int mbedtls_mpi_div_mpi(mbedtls_mpi *Q, mbedtls_mpi *R, const mbedtls_mpi *A,
1379                         const mbedtls_mpi *B)
1380 {
1381     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
1382     size_t i, n, t, k;
1383     mbedtls_mpi X, Y, Z, T1, T2;
1384     mbedtls_mpi_uint TP2[3];
1385     MPI_VALIDATE_RET(A != NULL);
1386     MPI_VALIDATE_RET(B != NULL);
1387 
1388     if (mbedtls_mpi_cmp_int(B, 0) == 0) {
1389         return MBEDTLS_ERR_MPI_DIVISION_BY_ZERO;
1390     }
1391 
1392     mbedtls_mpi_init_mempool(&X); mbedtls_mpi_init_mempool(&Y);
1393     mbedtls_mpi_init_mempool(&Z); mbedtls_mpi_init_mempool(&T1);
1394     /*
1395      * Avoid dynamic memory allocations for constant-size T2.
1396      *
1397      * T2 is used for comparison only and the 3 limbs are assigned explicitly,
1398      * so nobody increase the size of the MPI and we're safe to use an on-stack
1399      * buffer.
1400      */
1401     T2.s = 1;
1402     T2.n = sizeof(TP2) / sizeof(*TP2);
1403     T2.p = TP2;
1404 
1405     if (mbedtls_mpi_cmp_abs(A, B) < 0) {
1406         if (Q != NULL) {
1407             MBEDTLS_MPI_CHK(mbedtls_mpi_lset(Q, 0));
1408         }
1409         if (R != NULL) {
1410             MBEDTLS_MPI_CHK(mbedtls_mpi_copy(R, A));
1411         }
1412         return 0;
1413     }
1414 
1415     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&X, A));
1416     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&Y, B));
1417     X.s = Y.s = 1;
1418 
1419     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(&Z, A->n + 2));
1420     MBEDTLS_MPI_CHK(mbedtls_mpi_lset(&Z,  0));
1421     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(&T1, A->n + 2));
1422 
1423     k = mbedtls_mpi_bitlen(&Y) % biL;
1424     if (k < biL - 1) {
1425         k = biL - 1 - k;
1426         MBEDTLS_MPI_CHK(mbedtls_mpi_shift_l(&X, k));
1427         MBEDTLS_MPI_CHK(mbedtls_mpi_shift_l(&Y, k));
1428     } else {
1429         k = 0;
1430     }
1431 
1432     n = X.n - 1;
1433     t = Y.n - 1;
1434     MBEDTLS_MPI_CHK(mbedtls_mpi_shift_l(&Y, biL * (n - t)));
1435 
1436     while (mbedtls_mpi_cmp_mpi(&X, &Y) >= 0) {
1437         Z.p[n - t]++;
1438         MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&X, &X, &Y));
1439     }
1440     MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&Y, biL * (n - t)));
1441 
1442     for (i = n; i > t; i--) {
1443         if (X.p[i] >= Y.p[t]) {
1444             Z.p[i - t - 1] = ~(mbedtls_mpi_uint) 0u;
1445         } else {
1446             Z.p[i - t - 1] = mbedtls_int_div_int(X.p[i], X.p[i - 1],
1447                                                  Y.p[t], NULL);
1448         }
1449 
1450         T2.p[0] = (i < 2) ? 0 : X.p[i - 2];
1451         T2.p[1] = (i < 1) ? 0 : X.p[i - 1];
1452         T2.p[2] = X.p[i];
1453 
1454         Z.p[i - t - 1]++;
1455         do {
1456             Z.p[i - t - 1]--;
1457 
1458             MBEDTLS_MPI_CHK(mbedtls_mpi_lset(&T1, 0));
1459             T1.p[0] = (t < 1) ? 0 : Y.p[t - 1];
1460             T1.p[1] = Y.p[t];
1461             MBEDTLS_MPI_CHK(mbedtls_mpi_mul_int(&T1, &T1, Z.p[i - t - 1]));
1462         } while (mbedtls_mpi_cmp_mpi(&T1, &T2) > 0);
1463 
1464         MBEDTLS_MPI_CHK(mbedtls_mpi_mul_int(&T1, &Y, Z.p[i - t - 1]));
1465         MBEDTLS_MPI_CHK(mbedtls_mpi_shift_l(&T1,  biL * (i - t - 1)));
1466         MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&X, &X, &T1));
1467 
1468         if (mbedtls_mpi_cmp_int(&X, 0) < 0) {
1469             MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&T1, &Y));
1470             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_l(&T1, biL * (i - t - 1)));
1471             MBEDTLS_MPI_CHK(mbedtls_mpi_add_mpi(&X, &X, &T1));
1472             Z.p[i - t - 1]--;
1473         }
1474     }
1475 
1476     if (Q != NULL) {
1477         MBEDTLS_MPI_CHK(mbedtls_mpi_copy(Q, &Z));
1478         Q->s = A->s * B->s;
1479     }
1480 
1481     if (R != NULL) {
1482         MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&X, k));
1483         X.s = A->s;
1484         MBEDTLS_MPI_CHK(mbedtls_mpi_copy(R, &X));
1485 
1486         if (mbedtls_mpi_cmp_int(R, 0) == 0) {
1487             R->s = 1;
1488         }
1489     }
1490 
1491 cleanup:
1492 
1493     mbedtls_mpi_free(&X); mbedtls_mpi_free(&Y); mbedtls_mpi_free(&Z);
1494     mbedtls_mpi_free(&T1);
1495     mbedtls_platform_zeroize(TP2, sizeof(TP2));
1496 
1497     return ret;
1498 }
1499 
1500 /*
1501  * Division by int: A = Q * b + R
1502  */
1503 int mbedtls_mpi_div_int(mbedtls_mpi *Q, mbedtls_mpi *R,
1504                         const mbedtls_mpi *A,
1505                         mbedtls_mpi_sint b)
1506 {
1507     mbedtls_mpi B;
1508     mbedtls_mpi_uint p[1];
1509     MPI_VALIDATE_RET(A != NULL);
1510 
1511     p[0] = mpi_sint_abs(b);
1512     B.s = (b < 0) ? -1 : 1;
1513     B.n = 1;
1514     B.p = p;
1515 
1516     return mbedtls_mpi_div_mpi(Q, R, A, &B);
1517 }
1518 
1519 /*
1520  * Modulo: R = A mod B
1521  */
1522 int mbedtls_mpi_mod_mpi(mbedtls_mpi *R, const mbedtls_mpi *A, const mbedtls_mpi *B)
1523 {
1524     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
1525     MPI_VALIDATE_RET(R != NULL);
1526     MPI_VALIDATE_RET(A != NULL);
1527     MPI_VALIDATE_RET(B != NULL);
1528 
1529     if (mbedtls_mpi_cmp_int(B, 0) < 0) {
1530         return MBEDTLS_ERR_MPI_NEGATIVE_VALUE;
1531     }
1532 
1533     MBEDTLS_MPI_CHK(mbedtls_mpi_div_mpi(NULL, R, A, B));
1534 
1535     while (mbedtls_mpi_cmp_int(R, 0) < 0) {
1536         MBEDTLS_MPI_CHK(mbedtls_mpi_add_mpi(R, R, B));
1537     }
1538 
1539     while (mbedtls_mpi_cmp_mpi(R, B) >= 0) {
1540         MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(R, R, B));
1541     }
1542 
1543 cleanup:
1544 
1545     return ret;
1546 }
1547 
1548 /*
1549  * Modulo: r = A mod b
1550  */
1551 int mbedtls_mpi_mod_int(mbedtls_mpi_uint *r, const mbedtls_mpi *A, mbedtls_mpi_sint b)
1552 {
1553     size_t i;
1554     mbedtls_mpi_uint x, y, z;
1555     MPI_VALIDATE_RET(r != NULL);
1556     MPI_VALIDATE_RET(A != NULL);
1557 
1558     if (b == 0) {
1559         return MBEDTLS_ERR_MPI_DIVISION_BY_ZERO;
1560     }
1561 
1562     if (b < 0) {
1563         return MBEDTLS_ERR_MPI_NEGATIVE_VALUE;
1564     }
1565 
1566     /*
1567      * handle trivial cases
1568      */
1569     if (b == 1 || A->n == 0) {
1570         *r = 0;
1571         return 0;
1572     }
1573 
1574     if (b == 2) {
1575         *r = A->p[0] & 1;
1576         return 0;
1577     }
1578 
1579     /*
1580      * general case
1581      */
1582     for (i = A->n, y = 0; i > 0; i--) {
1583         x  = A->p[i - 1];
1584         y  = (y << biH) | (x >> biH);
1585         z  = y / b;
1586         y -= z * b;
1587 
1588         x <<= biH;
1589         y  = (y << biH) | (x >> biH);
1590         z  = y / b;
1591         y -= z * b;
1592     }
1593 
1594     /*
1595      * If A is negative, then the current y represents a negative value.
1596      * Flipping it to the positive side.
1597      */
1598     if (A->s < 0 && y != 0) {
1599         y = b - y;
1600     }
1601 
1602     *r = y;
1603 
1604     return 0;
1605 }
1606 
1607 static void mpi_montg_init(mbedtls_mpi_uint *mm, const mbedtls_mpi *N)
1608 {
1609     *mm = mbedtls_mpi_core_montmul_init(N->p);
1610 }
1611 
1612 void mbedtls_mpi_montg_init( mbedtls_mpi_uint *mm, const mbedtls_mpi *N )
1613 {
1614 	mpi_montg_init( mm, N );
1615 }
1616 
1617 /** Montgomery multiplication: A = A * B * R^-1 mod N  (HAC 14.36)
1618  *
1619  * \param[in,out]   A   One of the numbers to multiply.
1620  *                      It must have at least as many limbs as N
1621  *                      (A->n >= N->n), and any limbs beyond n are ignored.
1622  *                      On successful completion, A contains the result of
1623  *                      the multiplication A * B * R^-1 mod N where
1624  *                      R = (2^ciL)^n.
1625  * \param[in]       B   One of the numbers to multiply.
1626  *                      It must be nonzero and must not have more limbs than N
1627  *                      (B->n <= N->n).
1628  * \param[in]       N   The modulus. \p N must be odd.
1629  * \param           mm  The value calculated by `mpi_montg_init(&mm, N)`.
1630  *                      This is -N^-1 mod 2^ciL.
1631  * \param[in,out]   T   A bignum for temporary storage.
1632  *                      It must be at least twice the limb size of N plus 1
1633  *                      (T->n >= 2 * N->n + 1).
1634  *                      Its initial content is unused and
1635  *                      its final content is indeterminate.
1636  *                      It does not get reallocated.
1637  */
1638 static void mpi_montmul(mbedtls_mpi *A, const mbedtls_mpi *B,
1639                         const mbedtls_mpi *N, mbedtls_mpi_uint mm,
1640                         mbedtls_mpi *T)
1641 {
1642     mbedtls_mpi_core_montmul(A->p, A->p, B->p, B->n, N->p, N->n, mm, T->p);
1643 }
1644 
1645 void mbedtls_mpi_montmul(mbedtls_mpi *A, const mbedtls_mpi *B,
1646                          const mbedtls_mpi *N, mbedtls_mpi_uint mm,
1647                          mbedtls_mpi *T )
1648 {
1649     mpi_montmul( A, B, N, mm, T);
1650 }
1651 
1652 /*
1653  * Montgomery reduction: A = A * R^-1 mod N
1654  *
1655  * See mpi_montmul() regarding constraints and guarantees on the parameters.
1656  */
1657 static void mpi_montred(mbedtls_mpi *A, const mbedtls_mpi *N,
1658                         mbedtls_mpi_uint mm, mbedtls_mpi *T)
1659 {
1660     mbedtls_mpi_uint z = 1;
1661     mbedtls_mpi U;
1662 
1663     U.n = U.s = (int) z;
1664     U.p = &z;
1665 
1666     mpi_montmul(A, &U, N, mm, T);
1667 }
1668 
1669 void mbedtls_mpi_montred(mbedtls_mpi *A, const mbedtls_mpi *N,
1670                          mbedtls_mpi_uint mm, mbedtls_mpi *T)
1671 {
1672     mpi_montred(A, N, mm, T);
1673 }
1674 
1675 /**
1676  * Select an MPI from a table without leaking the index.
1677  *
1678  * This is functionally equivalent to mbedtls_mpi_copy(R, T[idx]) except it
1679  * reads the entire table in order to avoid leaking the value of idx to an
1680  * attacker able to observe memory access patterns.
1681  *
1682  * \param[out] R        Where to write the selected MPI.
1683  * \param[in] T         The table to read from.
1684  * \param[in] T_size    The number of elements in the table.
1685  * \param[in] idx       The index of the element to select;
1686  *                      this must satisfy 0 <= idx < T_size.
1687  *
1688  * \return \c 0 on success, or a negative error code.
1689  */
1690 static int mpi_select(mbedtls_mpi *R, const mbedtls_mpi *T, size_t T_size, size_t idx)
1691 {
1692     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
1693 
1694     for (size_t i = 0; i < T_size; i++) {
1695         MBEDTLS_MPI_CHK(mbedtls_mpi_safe_cond_assign(R, &T[i],
1696                                                      (unsigned char) mbedtls_ct_size_bool_eq(i,
1697                                                                                              idx)));
1698     }
1699 
1700 cleanup:
1701     return ret;
1702 }
1703 
1704 /*
1705  * Sliding-window exponentiation: X = A^E mod N  (HAC 14.85)
1706  */
1707 int mbedtls_mpi_exp_mod(mbedtls_mpi *X, const mbedtls_mpi *A,
1708                         const mbedtls_mpi *E, const mbedtls_mpi *N,
1709                         mbedtls_mpi *prec_RR)
1710 {
1711     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
1712     size_t window_bitsize;
1713     size_t i, j, nblimbs;
1714     size_t bufsize, nbits;
1715     size_t exponent_bits_in_window = 0;
1716     mbedtls_mpi_uint ei, mm, state;
1717     mbedtls_mpi RR, T, WW, Apos;
1718     mbedtls_mpi *W;
1719     const size_t array_size_W = 2 << MBEDTLS_MPI_WINDOW_SIZE;
1720     int neg;
1721 
1722     MPI_VALIDATE_RET(X != NULL);
1723     MPI_VALIDATE_RET(A != NULL);
1724     MPI_VALIDATE_RET(E != NULL);
1725     MPI_VALIDATE_RET(N != NULL);
1726 
1727     if (mbedtls_mpi_cmp_int(N, 0) <= 0 || (N->p[0] & 1) == 0) {
1728         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
1729     }
1730 
1731     if (mbedtls_mpi_cmp_int(E, 0) < 0) {
1732         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
1733     }
1734 
1735     if (mbedtls_mpi_bitlen(E) > MBEDTLS_MPI_MAX_BITS ||
1736         mbedtls_mpi_bitlen(N) > MBEDTLS_MPI_MAX_BITS) {
1737         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
1738     }
1739 
1740     /*
1741      * Init temps and window size
1742      */
1743     mpi_montg_init(&mm, N);
1744     mbedtls_mpi_init_mempool(&RR); mbedtls_mpi_init(&T);
1745     mbedtls_mpi_init_mempool(&Apos);
1746     mbedtls_mpi_init_mempool(&WW);
1747 
1748     i = mbedtls_mpi_bitlen(E);
1749 
1750     window_bitsize = (i > 671) ? 6 : (i > 239) ? 5 :
1751                      (i >  79) ? 4 : (i >  23) ? 3 : 1;
1752 
1753 #if (MBEDTLS_MPI_WINDOW_SIZE < 6)
1754     if (window_bitsize > MBEDTLS_MPI_WINDOW_SIZE) {
1755         window_bitsize = MBEDTLS_MPI_WINDOW_SIZE;
1756     }
1757 #endif
1758 
1759     const size_t w_table_used_size = (size_t) 1 << window_bitsize;
1760 
1761     W = mempool_alloc(mbedtls_mpi_mempool,
1762                       sizeof( mbedtls_mpi ) * array_size_W);
1763     if (W == NULL) {
1764         ret = MBEDTLS_ERR_MPI_ALLOC_FAILED;
1765         goto cleanup;
1766     }
1767     for (i = 0; i < array_size_W; i++)
1768         mbedtls_mpi_init_mempool(W + i);
1769     /*
1770      * This function is not constant-trace: its memory accesses depend on the
1771      * exponent value. To defend against timing attacks, callers (such as RSA
1772      * and DHM) should use exponent blinding. However this is not enough if the
1773      * adversary can find the exponent in a single trace, so this function
1774      * takes extra precautions against adversaries who can observe memory
1775      * access patterns.
1776      *
1777      * This function performs a series of multiplications by table elements and
1778      * squarings, and we want the prevent the adversary from finding out which
1779      * table element was used, and from distinguishing between multiplications
1780      * and squarings. Firstly, when multiplying by an element of the window
1781      * W[i], we do a constant-trace table lookup to obfuscate i. This leaves
1782      * squarings as having a different memory access patterns from other
1783      * multiplications. So secondly, we put the accumulator X in the table as
1784      * well, and also do a constant-trace table lookup to multiply by X.
1785      *
1786      * This way, all multiplications take the form of a lookup-and-multiply.
1787      * The number of lookup-and-multiply operations inside each iteration of
1788      * the main loop still depends on the bits of the exponent, but since the
1789      * other operations in the loop don't have an easily recognizable memory
1790      * trace, an adversary is unlikely to be able to observe the exact
1791      * patterns.
1792      *
1793      * An adversary may still be able to recover the exponent if they can
1794      * observe both memory accesses and branches. However, branch prediction
1795      * exploitation typically requires many traces of execution over the same
1796      * data, which is defeated by randomized blinding.
1797      *
1798      * To achieve this, we make a copy of X and we use the table entry in each
1799      * calculation from this point on.
1800      */
1801     const size_t x_index = 0;
1802     mbedtls_mpi_copy(&W[x_index], X);
1803 
1804     j = N->n + 1;
1805     /* All W[i] and X must have at least N->n limbs for the mpi_montmul()
1806      * and mpi_montred() calls later. Here we ensure that W[1] and X are
1807      * large enough, and later we'll grow other W[i] to the same length.
1808      * They must not be shrunk midway through this function!
1809      */
1810     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(&W[x_index], j));
1811 
1812     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(&T, j * 2));
1813 
1814     /*
1815      * Compensate for negative A (and correct at the end)
1816      */
1817     neg = (A->s == -1);
1818     if (neg) {
1819         MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&Apos, A));
1820         Apos.s = 1;
1821         A = &Apos;
1822     }
1823 
1824     /*
1825      * If 1st call, pre-compute R^2 mod N
1826      */
1827     if (prec_RR == NULL || prec_RR->p == NULL) {
1828         MBEDTLS_MPI_CHK(mbedtls_mpi_lset(&RR, 1));
1829         MBEDTLS_MPI_CHK(mbedtls_mpi_shift_l(&RR, N->n * 2 * biL));
1830         MBEDTLS_MPI_CHK(mbedtls_mpi_mod_mpi(&RR, &RR, N));
1831 
1832         if (prec_RR != NULL) {
1833             memcpy(prec_RR, &RR, sizeof(mbedtls_mpi));
1834         }
1835     } else {
1836         memcpy(&RR, prec_RR, sizeof(mbedtls_mpi));
1837     }
1838 
1839     MBEDTLS_MPI_CHK(mbedtls_mpi_grow(&W[1],  j));
1840 
1841     /*
1842      * W[1] = A * R^2 * R^-1 mod N = A * R mod N
1843      */
1844     if (mbedtls_mpi_cmp_mpi(A, N) >= 0) {
1845         MBEDTLS_MPI_CHK(mbedtls_mpi_mod_mpi(&W[1], A, N));
1846         /* This should be a no-op because W[1] is already that large before
1847          * mbedtls_mpi_mod_mpi(), but it's necessary to avoid an overflow
1848          * in mpi_montmul() below, so let's make sure. */
1849         MBEDTLS_MPI_CHK(mbedtls_mpi_grow(&W[1], N->n + 1));
1850     } else {
1851         MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&W[1], A));
1852     }
1853 
1854     /* Note that this is safe because W[1] always has at least N->n limbs
1855      * (it grew above and was preserved by mbedtls_mpi_copy()). */
1856     mpi_montmul(&W[1], &RR, N, mm, &T);
1857 
1858     /*
1859      * W[x_index] = R^2 * R^-1 mod N = R mod N
1860      */
1861     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&W[x_index], &RR));
1862     mpi_montred(&W[x_index], N, mm, &T);
1863 
1864 
1865     if (window_bitsize > 1) {
1866         /*
1867          * W[i] = W[1] ^ i
1868          *
1869          * The first bit of the sliding window is always 1 and therefore we
1870          * only need to store the second half of the table.
1871          *
1872          * (There are two special elements in the table: W[0] for the
1873          * accumulator/result and W[1] for A in Montgomery form. Both of these
1874          * are already set at this point.)
1875          */
1876         j = w_table_used_size / 2;
1877 
1878         MBEDTLS_MPI_CHK(mbedtls_mpi_grow(&W[j], N->n + 1));
1879         MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&W[j], &W[1]));
1880 
1881         for (i = 0; i < window_bitsize - 1; i++) {
1882             mpi_montmul(&W[j], &W[j], N, mm, &T);
1883         }
1884 
1885         /*
1886          * W[i] = W[i - 1] * W[1]
1887          */
1888         for (i = j + 1; i < w_table_used_size; i++) {
1889             MBEDTLS_MPI_CHK(mbedtls_mpi_grow(&W[i], N->n + 1));
1890             MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&W[i], &W[i - 1]));
1891 
1892             mpi_montmul(&W[i], &W[1], N, mm, &T);
1893         }
1894     }
1895 
1896     nblimbs = E->n;
1897     bufsize = 0;
1898     nbits   = 0;
1899     state   = 0;
1900 
1901     while (1) {
1902         if (bufsize == 0) {
1903             if (nblimbs == 0) {
1904                 break;
1905             }
1906 
1907             nblimbs--;
1908 
1909             bufsize = sizeof(mbedtls_mpi_uint) << 3;
1910         }
1911 
1912         bufsize--;
1913 
1914         ei = (E->p[nblimbs] >> bufsize) & 1;
1915 
1916         /*
1917          * skip leading 0s
1918          */
1919         if (ei == 0 && state == 0) {
1920             continue;
1921         }
1922 
1923         if (ei == 0 && state == 1) {
1924             /*
1925              * out of window, square W[x_index]
1926              */
1927             MBEDTLS_MPI_CHK(mpi_select(&WW, W, w_table_used_size, x_index));
1928             mpi_montmul(&W[x_index], &WW, N, mm, &T);
1929             continue;
1930         }
1931 
1932         /*
1933          * add ei to current window
1934          */
1935         state = 2;
1936 
1937         nbits++;
1938         exponent_bits_in_window |= (ei << (window_bitsize - nbits));
1939 
1940         if (nbits == window_bitsize) {
1941             /*
1942              * W[x_index] = W[x_index]^window_bitsize R^-1 mod N
1943              */
1944             for (i = 0; i < window_bitsize; i++) {
1945                 MBEDTLS_MPI_CHK(mpi_select(&WW, W, w_table_used_size,
1946                                            x_index));
1947                 mpi_montmul(&W[x_index], &WW, N, mm, &T);
1948             }
1949 
1950             /*
1951              * W[x_index] = W[x_index] * W[exponent_bits_in_window] R^-1 mod N
1952              */
1953             MBEDTLS_MPI_CHK(mpi_select(&WW, W, w_table_used_size,
1954                                        exponent_bits_in_window));
1955             mpi_montmul(&W[x_index], &WW, N, mm, &T);
1956 
1957             state--;
1958             nbits = 0;
1959             exponent_bits_in_window = 0;
1960         }
1961     }
1962 
1963     /*
1964      * process the remaining bits
1965      */
1966     for (i = 0; i < nbits; i++) {
1967         MBEDTLS_MPI_CHK(mpi_select(&WW, W, w_table_used_size, x_index));
1968         mpi_montmul(&W[x_index], &WW, N, mm, &T);
1969 
1970         exponent_bits_in_window <<= 1;
1971 
1972         if ((exponent_bits_in_window & ((size_t) 1 << window_bitsize)) != 0) {
1973             MBEDTLS_MPI_CHK(mpi_select(&WW, W, w_table_used_size, 1));
1974             mpi_montmul(&W[x_index], &WW, N, mm, &T);
1975         }
1976     }
1977 
1978     /*
1979      * W[x_index] = A^E * R * R^-1 mod N = A^E mod N
1980      */
1981     mpi_montred(&W[x_index], N, mm, &T);
1982 
1983     if (neg && E->n != 0 && (E->p[0] & 1) != 0) {
1984         W[x_index].s = -1;
1985         MBEDTLS_MPI_CHK(mbedtls_mpi_add_mpi(&W[x_index], N, &W[x_index]));
1986     }
1987 
1988     /*
1989      * Load the result in the output variable.
1990      */
1991     mbedtls_mpi_copy(X, &W[x_index]);
1992 
1993 cleanup:
1994 
1995     if (W)
1996         for (i = 0; i < array_size_W; i++)
1997             mbedtls_mpi_free(W + i);
1998     mempool_free(mbedtls_mpi_mempool , W);
1999 
2000     mbedtls_mpi_free(&T);
2001     mbedtls_mpi_free(&Apos);
2002     mbedtls_mpi_free(&WW);
2003 
2004     if (prec_RR == NULL || prec_RR->p == NULL) {
2005         mbedtls_mpi_free(&RR);
2006     }
2007 
2008     return ret;
2009 }
2010 
2011 /*
2012  * Greatest common divisor: G = gcd(A, B)  (HAC 14.54)
2013  */
2014 int mbedtls_mpi_gcd(mbedtls_mpi *G, const mbedtls_mpi *A, const mbedtls_mpi *B)
2015 {
2016     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
2017     size_t lz, lzt;
2018     mbedtls_mpi TA, TB;
2019 
2020     MPI_VALIDATE_RET(G != NULL);
2021     MPI_VALIDATE_RET(A != NULL);
2022     MPI_VALIDATE_RET(B != NULL);
2023 
2024     mbedtls_mpi_init_mempool(&TA); mbedtls_mpi_init_mempool(&TB);
2025 
2026     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&TA, A));
2027     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&TB, B));
2028 
2029     lz = mbedtls_mpi_lsb(&TA);
2030     lzt = mbedtls_mpi_lsb(&TB);
2031 
2032     /* The loop below gives the correct result when A==0 but not when B==0.
2033      * So have a special case for B==0. Leverage the fact that we just
2034      * calculated the lsb and lsb(B)==0 iff B is odd or 0 to make the test
2035      * slightly more efficient than cmp_int(). */
2036     if (lzt == 0 && mbedtls_mpi_get_bit(&TB, 0) == 0) {
2037         ret = mbedtls_mpi_copy(G, A);
2038         goto cleanup;
2039     }
2040 
2041     if (lzt < lz) {
2042         lz = lzt;
2043     }
2044 
2045     TA.s = TB.s = 1;
2046 
2047     /* We mostly follow the procedure described in HAC 14.54, but with some
2048      * minor differences:
2049      * - Sequences of multiplications or divisions by 2 are grouped into a
2050      *   single shift operation.
2051      * - The procedure in HAC assumes that 0 < TB <= TA.
2052      *     - The condition TB <= TA is not actually necessary for correctness.
2053      *       TA and TB have symmetric roles except for the loop termination
2054      *       condition, and the shifts at the beginning of the loop body
2055      *       remove any significance from the ordering of TA vs TB before
2056      *       the shifts.
2057      *     - If TA = 0, the loop goes through 0 iterations and the result is
2058      *       correctly TB.
2059      *     - The case TB = 0 was short-circuited above.
2060      *
2061      * For the correctness proof below, decompose the original values of
2062      * A and B as
2063      *   A = sa * 2^a * A' with A'=0 or A' odd, and sa = +-1
2064      *   B = sb * 2^b * B' with B'=0 or B' odd, and sb = +-1
2065      * Then gcd(A, B) = 2^{min(a,b)} * gcd(A',B'),
2066      * and gcd(A',B') is odd or 0.
2067      *
2068      * At the beginning, we have TA = |A| and TB = |B| so gcd(A,B) = gcd(TA,TB).
2069      * The code maintains the following invariant:
2070      *     gcd(A,B) = 2^k * gcd(TA,TB) for some k   (I)
2071      */
2072 
2073     /* Proof that the loop terminates:
2074      * At each iteration, either the right-shift by 1 is made on a nonzero
2075      * value and the nonnegative integer bitlen(TA) + bitlen(TB) decreases
2076      * by at least 1, or the right-shift by 1 is made on zero and then
2077      * TA becomes 0 which ends the loop (TB cannot be 0 if it is right-shifted
2078      * since in that case TB is calculated from TB-TA with the condition TB>TA).
2079      */
2080     while (mbedtls_mpi_cmp_int(&TA, 0) != 0) {
2081         /* Divisions by 2 preserve the invariant (I). */
2082         MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&TA, mbedtls_mpi_lsb(&TA)));
2083         MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&TB, mbedtls_mpi_lsb(&TB)));
2084 
2085         /* Set either TA or TB to |TA-TB|/2. Since TA and TB are both odd,
2086          * TA-TB is even so the division by 2 has an integer result.
2087          * Invariant (I) is preserved since any odd divisor of both TA and TB
2088          * also divides |TA-TB|/2, and any odd divisor of both TA and |TA-TB|/2
2089          * also divides TB, and any odd divisor of both TB and |TA-TB|/2 also
2090          * divides TA.
2091          */
2092         if (mbedtls_mpi_cmp_mpi(&TA, &TB) >= 0) {
2093             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_abs(&TA, &TA, &TB));
2094             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&TA, 1));
2095         } else {
2096             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_abs(&TB, &TB, &TA));
2097             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&TB, 1));
2098         }
2099         /* Note that one of TA or TB is still odd. */
2100     }
2101 
2102     /* By invariant (I), gcd(A,B) = 2^k * gcd(TA,TB) for some k.
2103      * At the loop exit, TA = 0, so gcd(TA,TB) = TB.
2104      * - If there was at least one loop iteration, then one of TA or TB is odd,
2105      *   and TA = 0, so TB is odd and gcd(TA,TB) = gcd(A',B'). In this case,
2106      *   lz = min(a,b) so gcd(A,B) = 2^lz * TB.
2107      * - If there was no loop iteration, then A was 0, and gcd(A,B) = B.
2108      *   In this case, lz = 0 and B = TB so gcd(A,B) = B = 2^lz * TB as well.
2109      */
2110 
2111     MBEDTLS_MPI_CHK(mbedtls_mpi_shift_l(&TB, lz));
2112     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(G, &TB));
2113 
2114 cleanup:
2115 
2116     mbedtls_mpi_free(&TA); mbedtls_mpi_free(&TB);
2117 
2118     return ret;
2119 }
2120 
2121 /*
2122  * Fill X with size bytes of random.
2123  * The bytes returned from the RNG are used in a specific order which
2124  * is suitable for deterministic ECDSA (see the specification of
2125  * mbedtls_mpi_random() and the implementation in mbedtls_mpi_fill_random()).
2126  */
2127 int mbedtls_mpi_fill_random(mbedtls_mpi *X, size_t size,
2128                             int (*f_rng)(void *, unsigned char *, size_t),
2129                             void *p_rng)
2130 {
2131     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
2132     const size_t limbs = CHARS_TO_LIMBS(size);
2133 
2134     MPI_VALIDATE_RET(X     != NULL);
2135     MPI_VALIDATE_RET(f_rng != NULL);
2136 
2137     /* Ensure that target MPI has exactly the necessary number of limbs */
2138     MBEDTLS_MPI_CHK(mbedtls_mpi_resize_clear(X, limbs));
2139     if (size == 0) {
2140         return 0;
2141     }
2142 
2143     ret = mbedtls_mpi_core_fill_random(X->p, X->n, size, f_rng, p_rng);
2144 
2145 cleanup:
2146     return ret;
2147 }
2148 
2149 int mbedtls_mpi_random(mbedtls_mpi *X,
2150                        mbedtls_mpi_sint min,
2151                        const mbedtls_mpi *N,
2152                        int (*f_rng)(void *, unsigned char *, size_t),
2153                        void *p_rng)
2154 {
2155     if (min < 0) {
2156         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
2157     }
2158     if (mbedtls_mpi_cmp_int(N, min) <= 0) {
2159         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
2160     }
2161 
2162     /* Ensure that target MPI has exactly the same number of limbs
2163      * as the upper bound, even if the upper bound has leading zeros.
2164      * This is necessary for mbedtls_mpi_core_random. */
2165     int ret = mbedtls_mpi_resize_clear(X, N->n);
2166     if (ret != 0) {
2167         return ret;
2168     }
2169 
2170     return mbedtls_mpi_core_random(X->p, min, N->p, X->n, f_rng, p_rng);
2171 }
2172 
2173 /*
2174  * Modular inverse: X = A^-1 mod N  (HAC 14.61 / 14.64)
2175  */
2176 int mbedtls_mpi_inv_mod(mbedtls_mpi *X, const mbedtls_mpi *A, const mbedtls_mpi *N)
2177 {
2178     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
2179     mbedtls_mpi G, TA, TU, U1, U2, TB, TV, V1, V2;
2180     MPI_VALIDATE_RET(X != NULL);
2181     MPI_VALIDATE_RET(A != NULL);
2182     MPI_VALIDATE_RET(N != NULL);
2183 
2184     if (mbedtls_mpi_cmp_int(N, 1) <= 0) {
2185         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
2186     }
2187 
2188     mbedtls_mpi_init_mempool(&TA); mbedtls_mpi_init_mempool(&TU);
2189     mbedtls_mpi_init_mempool(&U1); mbedtls_mpi_init_mempool(&U2);
2190     mbedtls_mpi_init_mempool(&G); mbedtls_mpi_init_mempool(&TB);
2191     mbedtls_mpi_init_mempool(&TV); mbedtls_mpi_init_mempool(&V1);
2192     mbedtls_mpi_init_mempool(&V2);
2193 
2194     MBEDTLS_MPI_CHK(mbedtls_mpi_gcd(&G, A, N));
2195 
2196     if (mbedtls_mpi_cmp_int(&G, 1) != 0) {
2197         ret = MBEDTLS_ERR_MPI_NOT_ACCEPTABLE;
2198         goto cleanup;
2199     }
2200 
2201     MBEDTLS_MPI_CHK(mbedtls_mpi_mod_mpi(&TA, A, N));
2202     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&TU, &TA));
2203     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&TB, N));
2204     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&TV, N));
2205 
2206     MBEDTLS_MPI_CHK(mbedtls_mpi_lset(&U1, 1));
2207     MBEDTLS_MPI_CHK(mbedtls_mpi_lset(&U2, 0));
2208     MBEDTLS_MPI_CHK(mbedtls_mpi_lset(&V1, 0));
2209     MBEDTLS_MPI_CHK(mbedtls_mpi_lset(&V2, 1));
2210 
2211     do {
2212         while ((TU.p[0] & 1) == 0) {
2213             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&TU, 1));
2214 
2215             if ((U1.p[0] & 1) != 0 || (U2.p[0] & 1) != 0) {
2216                 MBEDTLS_MPI_CHK(mbedtls_mpi_add_mpi(&U1, &U1, &TB));
2217                 MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&U2, &U2, &TA));
2218             }
2219 
2220             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&U1, 1));
2221             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&U2, 1));
2222         }
2223 
2224         while ((TV.p[0] & 1) == 0) {
2225             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&TV, 1));
2226 
2227             if ((V1.p[0] & 1) != 0 || (V2.p[0] & 1) != 0) {
2228                 MBEDTLS_MPI_CHK(mbedtls_mpi_add_mpi(&V1, &V1, &TB));
2229                 MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&V2, &V2, &TA));
2230             }
2231 
2232             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&V1, 1));
2233             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&V2, 1));
2234         }
2235 
2236         if (mbedtls_mpi_cmp_mpi(&TU, &TV) >= 0) {
2237             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&TU, &TU, &TV));
2238             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&U1, &U1, &V1));
2239             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&U2, &U2, &V2));
2240         } else {
2241             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&TV, &TV, &TU));
2242             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&V1, &V1, &U1));
2243             MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&V2, &V2, &U2));
2244         }
2245     } while (mbedtls_mpi_cmp_int(&TU, 0) != 0);
2246 
2247     while (mbedtls_mpi_cmp_int(&V1, 0) < 0) {
2248         MBEDTLS_MPI_CHK(mbedtls_mpi_add_mpi(&V1, &V1, N));
2249     }
2250 
2251     while (mbedtls_mpi_cmp_mpi(&V1, N) >= 0) {
2252         MBEDTLS_MPI_CHK(mbedtls_mpi_sub_mpi(&V1, &V1, N));
2253     }
2254 
2255     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(X, &V1));
2256 
2257 cleanup:
2258 
2259     mbedtls_mpi_free(&TA); mbedtls_mpi_free(&TU); mbedtls_mpi_free(&U1); mbedtls_mpi_free(&U2);
2260     mbedtls_mpi_free(&G); mbedtls_mpi_free(&TB); mbedtls_mpi_free(&TV);
2261     mbedtls_mpi_free(&V1); mbedtls_mpi_free(&V2);
2262 
2263     return ret;
2264 }
2265 
2266 #if defined(MBEDTLS_GENPRIME)
2267 
2268 static const int small_prime[] =
2269 {
2270     3,    5,    7,   11,   13,   17,   19,   23,
2271     29,   31,   37,   41,   43,   47,   53,   59,
2272     61,   67,   71,   73,   79,   83,   89,   97,
2273     101,  103,  107,  109,  113,  127,  131,  137,
2274     139,  149,  151,  157,  163,  167,  173,  179,
2275     181,  191,  193,  197,  199,  211,  223,  227,
2276     229,  233,  239,  241,  251,  257,  263,  269,
2277     271,  277,  281,  283,  293,  307,  311,  313,
2278     317,  331,  337,  347,  349,  353,  359,  367,
2279     373,  379,  383,  389,  397,  401,  409,  419,
2280     421,  431,  433,  439,  443,  449,  457,  461,
2281     463,  467,  479,  487,  491,  499,  503,  509,
2282     521,  523,  541,  547,  557,  563,  569,  571,
2283     577,  587,  593,  599,  601,  607,  613,  617,
2284     619,  631,  641,  643,  647,  653,  659,  661,
2285     673,  677,  683,  691,  701,  709,  719,  727,
2286     733,  739,  743,  751,  757,  761,  769,  773,
2287     787,  797,  809,  811,  821,  823,  827,  829,
2288     839,  853,  857,  859,  863,  877,  881,  883,
2289     887,  907,  911,  919,  929,  937,  941,  947,
2290     953,  967,  971,  977,  983,  991,  997, -103
2291 };
2292 
2293 /*
2294  * Small divisors test (X must be positive)
2295  *
2296  * Return values:
2297  * 0: no small factor (possible prime, more tests needed)
2298  * 1: certain prime
2299  * MBEDTLS_ERR_MPI_NOT_ACCEPTABLE: certain non-prime
2300  * other negative: error
2301  */
2302 static int mpi_check_small_factors(const mbedtls_mpi *X)
2303 {
2304     int ret = 0;
2305     size_t i;
2306     mbedtls_mpi_uint r;
2307 
2308     if ((X->p[0] & 1) == 0) {
2309         return MBEDTLS_ERR_MPI_NOT_ACCEPTABLE;
2310     }
2311 
2312     for (i = 0; small_prime[i] > 0; i++) {
2313         if (mbedtls_mpi_cmp_int(X, small_prime[i]) <= 0) {
2314             return 1;
2315         }
2316 
2317         MBEDTLS_MPI_CHK(mbedtls_mpi_mod_int(&r, X, small_prime[i]));
2318 
2319         if (r == 0) {
2320             return MBEDTLS_ERR_MPI_NOT_ACCEPTABLE;
2321         }
2322     }
2323 
2324 cleanup:
2325     return ret;
2326 }
2327 
2328 /*
2329  * Miller-Rabin pseudo-primality test  (HAC 4.24)
2330  */
2331 static int mpi_miller_rabin(const mbedtls_mpi *X, size_t rounds,
2332                             int (*f_rng)(void *, unsigned char *, size_t),
2333                             void *p_rng)
2334 {
2335     int ret, count;
2336     size_t i, j, k, s;
2337     mbedtls_mpi W, R, T, A, RR;
2338 
2339     MPI_VALIDATE_RET(X     != NULL);
2340     MPI_VALIDATE_RET(f_rng != NULL);
2341 
2342     mbedtls_mpi_init_mempool(&W); mbedtls_mpi_init_mempool(&R);
2343     mbedtls_mpi_init_mempool(&T); mbedtls_mpi_init_mempool(&A);
2344     mbedtls_mpi_init_mempool(&RR);
2345 
2346     /*
2347      * W = |X| - 1
2348      * R = W >> lsb( W )
2349      */
2350     MBEDTLS_MPI_CHK(mbedtls_mpi_sub_int(&W, X, 1));
2351     s = mbedtls_mpi_lsb(&W);
2352     MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&R, &W));
2353     MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&R, s));
2354 
2355     for (i = 0; i < rounds; i++) {
2356         /*
2357          * pick a random A, 1 < A < |X| - 1
2358          */
2359         count = 0;
2360         do {
2361             MBEDTLS_MPI_CHK(mbedtls_mpi_fill_random(&A, X->n * ciL, f_rng, p_rng));
2362 
2363             j = mbedtls_mpi_bitlen(&A);
2364             k = mbedtls_mpi_bitlen(&W);
2365             if (j > k) {
2366                 A.p[A.n - 1] &= ((mbedtls_mpi_uint) 1 << (k - (A.n - 1) * biL - 1)) - 1;
2367             }
2368 
2369             if (count++ > 300) {
2370                 ret = MBEDTLS_ERR_MPI_NOT_ACCEPTABLE;
2371                 goto cleanup;
2372             }
2373 
2374         } while (mbedtls_mpi_cmp_mpi(&A, &W) >= 0 ||
2375                  mbedtls_mpi_cmp_int(&A, 1)  <= 0);
2376 
2377         /*
2378          * A = A^R mod |X|
2379          */
2380         MBEDTLS_MPI_CHK(mbedtls_mpi_exp_mod(&A, &A, &R, X, &RR));
2381 
2382         if (mbedtls_mpi_cmp_mpi(&A, &W) == 0 ||
2383             mbedtls_mpi_cmp_int(&A,  1) == 0) {
2384             continue;
2385         }
2386 
2387         j = 1;
2388         while (j < s && mbedtls_mpi_cmp_mpi(&A, &W) != 0) {
2389             /*
2390              * A = A * A mod |X|
2391              */
2392             MBEDTLS_MPI_CHK(mbedtls_mpi_mul_mpi(&T, &A, &A));
2393             MBEDTLS_MPI_CHK(mbedtls_mpi_mod_mpi(&A, &T, X));
2394 
2395             if (mbedtls_mpi_cmp_int(&A, 1) == 0) {
2396                 break;
2397             }
2398 
2399             j++;
2400         }
2401 
2402         /*
2403          * not prime if A != |X| - 1 or A == 1
2404          */
2405         if (mbedtls_mpi_cmp_mpi(&A, &W) != 0 ||
2406             mbedtls_mpi_cmp_int(&A,  1) == 0) {
2407             ret = MBEDTLS_ERR_MPI_NOT_ACCEPTABLE;
2408             break;
2409         }
2410     }
2411 
2412 cleanup:
2413     mbedtls_mpi_free(&W); mbedtls_mpi_free(&R);
2414     mbedtls_mpi_free(&T); mbedtls_mpi_free(&A);
2415     mbedtls_mpi_free(&RR);
2416 
2417     return ret;
2418 }
2419 
2420 /*
2421  * Pseudo-primality test: small factors, then Miller-Rabin
2422  */
2423 int mbedtls_mpi_is_prime_ext(const mbedtls_mpi *X, int rounds,
2424                              int (*f_rng)(void *, unsigned char *, size_t),
2425                              void *p_rng)
2426 {
2427     int ret = MBEDTLS_ERR_ERROR_CORRUPTION_DETECTED;
2428     mbedtls_mpi XX;
2429     MPI_VALIDATE_RET(X     != NULL);
2430     MPI_VALIDATE_RET(f_rng != NULL);
2431 
2432     XX.s = 1;
2433     XX.n = X->n;
2434     XX.p = X->p;
2435 
2436     if (mbedtls_mpi_cmp_int(&XX, 0) == 0 ||
2437         mbedtls_mpi_cmp_int(&XX, 1) == 0) {
2438         return MBEDTLS_ERR_MPI_NOT_ACCEPTABLE;
2439     }
2440 
2441     if (mbedtls_mpi_cmp_int(&XX, 2) == 0) {
2442         return 0;
2443     }
2444 
2445     if ((ret = mpi_check_small_factors(&XX)) != 0) {
2446         if (ret == 1) {
2447             return 0;
2448         }
2449 
2450         return ret;
2451     }
2452 
2453     return mpi_miller_rabin(&XX, rounds, f_rng, p_rng);
2454 }
2455 
2456 /*
2457  * Prime number generation
2458  *
2459  * To generate an RSA key in a way recommended by FIPS 186-4, both primes must
2460  * be either 1024 bits or 1536 bits long, and flags must contain
2461  * MBEDTLS_MPI_GEN_PRIME_FLAG_LOW_ERR.
2462  */
2463 int mbedtls_mpi_gen_prime(mbedtls_mpi *X, size_t nbits, int flags,
2464                           int (*f_rng)(void *, unsigned char *, size_t),
2465                           void *p_rng)
2466 {
2467 #ifdef MBEDTLS_HAVE_INT64
2468 // ceil(2^63.5)
2469 #define CEIL_MAXUINT_DIV_SQRT2 0xb504f333f9de6485ULL
2470 #else
2471 // ceil(2^31.5)
2472 #define CEIL_MAXUINT_DIV_SQRT2 0xb504f334U
2473 #endif
2474     int ret = MBEDTLS_ERR_MPI_NOT_ACCEPTABLE;
2475     size_t k, n;
2476     int rounds;
2477     mbedtls_mpi_uint r;
2478     mbedtls_mpi Y;
2479 
2480     MPI_VALIDATE_RET(X     != NULL);
2481     MPI_VALIDATE_RET(f_rng != NULL);
2482 
2483     if (nbits < 3 || nbits > MBEDTLS_MPI_MAX_BITS) {
2484         return MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
2485     }
2486 
2487     mbedtls_mpi_init_mempool(&Y);
2488 
2489     n = BITS_TO_LIMBS(nbits);
2490 
2491     if ((flags & MBEDTLS_MPI_GEN_PRIME_FLAG_LOW_ERR) == 0) {
2492         /*
2493          * 2^-80 error probability, number of rounds chosen per HAC, table 4.4
2494          */
2495         rounds = ((nbits >= 1300) ?  2 : (nbits >=  850) ?  3 :
2496                   (nbits >=  650) ?  4 : (nbits >=  350) ?  8 :
2497                   (nbits >=  250) ? 12 : (nbits >=  150) ? 18 : 27);
2498     } else {
2499         /*
2500          * 2^-100 error probability, number of rounds computed based on HAC,
2501          * fact 4.48
2502          */
2503         rounds = ((nbits >= 1450) ?  4 : (nbits >=  1150) ?  5 :
2504                   (nbits >= 1000) ?  6 : (nbits >=   850) ?  7 :
2505                   (nbits >=  750) ?  8 : (nbits >=   500) ? 13 :
2506                   (nbits >=  250) ? 28 : (nbits >=   150) ? 40 : 51);
2507     }
2508 
2509     while (1) {
2510         MBEDTLS_MPI_CHK(mbedtls_mpi_fill_random(X, n * ciL, f_rng, p_rng));
2511         /* make sure generated number is at least (nbits-1)+0.5 bits (FIPS 186-4 §B.3.3 steps 4.4, 5.5) */
2512         if (X->p[n-1] < CEIL_MAXUINT_DIV_SQRT2) {
2513             continue;
2514         }
2515 
2516         k = n * biL;
2517         if (k > nbits) {
2518             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(X, k - nbits));
2519         }
2520         X->p[0] |= 1;
2521 
2522         if ((flags & MBEDTLS_MPI_GEN_PRIME_FLAG_DH) == 0) {
2523             ret = mbedtls_mpi_is_prime_ext(X, rounds, f_rng, p_rng);
2524 
2525             if (ret != MBEDTLS_ERR_MPI_NOT_ACCEPTABLE) {
2526                 goto cleanup;
2527             }
2528         } else {
2529             /*
2530              * A necessary condition for Y and X = 2Y + 1 to be prime
2531              * is X = 2 mod 3 (which is equivalent to Y = 2 mod 3).
2532              * Make sure it is satisfied, while keeping X = 3 mod 4
2533              */
2534 
2535             X->p[0] |= 2;
2536 
2537             MBEDTLS_MPI_CHK(mbedtls_mpi_mod_int(&r, X, 3));
2538             if (r == 0) {
2539                 MBEDTLS_MPI_CHK(mbedtls_mpi_add_int(X, X, 8));
2540             } else if (r == 1) {
2541                 MBEDTLS_MPI_CHK(mbedtls_mpi_add_int(X, X, 4));
2542             }
2543 
2544             /* Set Y = (X-1) / 2, which is X / 2 because X is odd */
2545             MBEDTLS_MPI_CHK(mbedtls_mpi_copy(&Y, X));
2546             MBEDTLS_MPI_CHK(mbedtls_mpi_shift_r(&Y, 1));
2547 
2548             while (1) {
2549                 /*
2550                  * First, check small factors for X and Y
2551                  * before doing Miller-Rabin on any of them
2552                  */
2553                 if ((ret = mpi_check_small_factors(X)) == 0 &&
2554                     (ret = mpi_check_small_factors(&Y)) == 0 &&
2555                     (ret = mpi_miller_rabin(X, rounds, f_rng, p_rng))
2556                     == 0 &&
2557                     (ret = mpi_miller_rabin(&Y, rounds, f_rng, p_rng))
2558                     == 0) {
2559                     goto cleanup;
2560                 }
2561 
2562                 if (ret != MBEDTLS_ERR_MPI_NOT_ACCEPTABLE) {
2563                     goto cleanup;
2564                 }
2565 
2566                 /*
2567                  * Next candidates. We want to preserve Y = (X-1) / 2 and
2568                  * Y = 1 mod 2 and Y = 2 mod 3 (eq X = 3 mod 4 and X = 2 mod 3)
2569                  * so up Y by 6 and X by 12.
2570                  */
2571                 MBEDTLS_MPI_CHK(mbedtls_mpi_add_int(X,  X, 12));
2572                 MBEDTLS_MPI_CHK(mbedtls_mpi_add_int(&Y, &Y, 6));
2573             }
2574         }
2575     }
2576 
2577 cleanup:
2578 
2579     mbedtls_mpi_free(&Y);
2580 
2581     return ret;
2582 }
2583 
2584 #endif /* MBEDTLS_GENPRIME */
2585 
2586 #if defined(MBEDTLS_SELF_TEST)
2587 
2588 #define GCD_PAIR_COUNT  3
2589 
2590 static const int gcd_pairs[GCD_PAIR_COUNT][3] =
2591 {
2592     { 693, 609, 21 },
2593     { 1764, 868, 28 },
2594     { 768454923, 542167814, 1 }
2595 };
2596 
2597 /*
2598  * Checkup routine
2599  */
2600 int mbedtls_mpi_self_test(int verbose)
2601 {
2602     int ret, i;
2603     mbedtls_mpi A, E, N, X, Y, U, V;
2604 
2605     mbedtls_mpi_init_mempool(&A); mbedtls_mpi_init_mempool(&E);
2606     mbedtls_mpi_init_mempool(&N); mbedtls_mpi_init_mempool(&X);
2607     mbedtls_mpi_init_mempool(&Y); mbedtls_mpi_init_mempool(&U);
2608     mbedtls_mpi_init_mempool(&V);
2609 
2610     MBEDTLS_MPI_CHK(mbedtls_mpi_read_string(&A, 16,
2611                                             "EFE021C2645FD1DC586E69184AF4A31E" \
2612                                             "D5F53E93B5F123FA41680867BA110131" \
2613                                             "944FE7952E2517337780CB0DB80E61AA" \
2614                                             "E7C8DDC6C5C6AADEB34EB38A2F40D5E6"));
2615 
2616     MBEDTLS_MPI_CHK(mbedtls_mpi_read_string(&E, 16,
2617                                             "B2E7EFD37075B9F03FF989C7C5051C20" \
2618                                             "34D2A323810251127E7BF8625A4F49A5" \
2619                                             "F3E27F4DA8BD59C47D6DAABA4C8127BD" \
2620                                             "5B5C25763222FEFCCFC38B832366C29E"));
2621 
2622     MBEDTLS_MPI_CHK(mbedtls_mpi_read_string(&N, 16,
2623                                             "0066A198186C18C10B2F5ED9B522752A" \
2624                                             "9830B69916E535C8F047518A889A43A5" \
2625                                             "94B6BED27A168D31D4A52F88925AA8F5"));
2626 
2627     MBEDTLS_MPI_CHK(mbedtls_mpi_mul_mpi(&X, &A, &N));
2628 
2629     MBEDTLS_MPI_CHK(mbedtls_mpi_read_string(&U, 16,
2630                                             "602AB7ECA597A3D6B56FF9829A5E8B85" \
2631                                             "9E857EA95A03512E2BAE7391688D264A" \
2632                                             "A5663B0341DB9CCFD2C4C5F421FEC814" \
2633                                             "8001B72E848A38CAE1C65F78E56ABDEF" \
2634                                             "E12D3C039B8A02D6BE593F0BBBDA56F1" \
2635                                             "ECF677152EF804370C1A305CAF3B5BF1" \
2636                                             "30879B56C61DE584A0F53A2447A51E"));
2637 
2638     if (verbose != 0) {
2639         mbedtls_printf("  MPI test #1 (mul_mpi): ");
2640     }
2641 
2642     if (mbedtls_mpi_cmp_mpi(&X, &U) != 0) {
2643         if (verbose != 0) {
2644             mbedtls_printf("failed\n");
2645         }
2646 
2647         ret = 1;
2648         goto cleanup;
2649     }
2650 
2651     if (verbose != 0) {
2652         mbedtls_printf("passed\n");
2653     }
2654 
2655     MBEDTLS_MPI_CHK(mbedtls_mpi_div_mpi(&X, &Y, &A, &N));
2656 
2657     MBEDTLS_MPI_CHK(mbedtls_mpi_read_string(&U, 16,
2658                                             "256567336059E52CAE22925474705F39A94"));
2659 
2660     MBEDTLS_MPI_CHK(mbedtls_mpi_read_string(&V, 16,
2661                                             "6613F26162223DF488E9CD48CC132C7A" \
2662                                             "0AC93C701B001B092E4E5B9F73BCD27B" \
2663                                             "9EE50D0657C77F374E903CDFA4C642"));
2664 
2665     if (verbose != 0) {
2666         mbedtls_printf("  MPI test #2 (div_mpi): ");
2667     }
2668 
2669     if (mbedtls_mpi_cmp_mpi(&X, &U) != 0 ||
2670         mbedtls_mpi_cmp_mpi(&Y, &V) != 0) {
2671         if (verbose != 0) {
2672             mbedtls_printf("failed\n");
2673         }
2674 
2675         ret = 1;
2676         goto cleanup;
2677     }
2678 
2679     if (verbose != 0) {
2680         mbedtls_printf("passed\n");
2681     }
2682 
2683     MBEDTLS_MPI_CHK(mbedtls_mpi_exp_mod(&X, &A, &E, &N, NULL));
2684 
2685     MBEDTLS_MPI_CHK(mbedtls_mpi_read_string(&U, 16,
2686                                             "36E139AEA55215609D2816998ED020BB" \
2687                                             "BD96C37890F65171D948E9BC7CBAA4D9" \
2688                                             "325D24D6A3C12710F10A09FA08AB87"));
2689 
2690     if (verbose != 0) {
2691         mbedtls_printf("  MPI test #3 (exp_mod): ");
2692     }
2693 
2694     if (mbedtls_mpi_cmp_mpi(&X, &U) != 0) {
2695         if (verbose != 0) {
2696             mbedtls_printf("failed\n");
2697         }
2698 
2699         ret = 1;
2700         goto cleanup;
2701     }
2702 
2703     if (verbose != 0) {
2704         mbedtls_printf("passed\n");
2705     }
2706 
2707     MBEDTLS_MPI_CHK(mbedtls_mpi_inv_mod(&X, &A, &N));
2708 
2709     MBEDTLS_MPI_CHK(mbedtls_mpi_read_string(&U, 16,
2710                                             "003A0AAEDD7E784FC07D8F9EC6E3BFD5" \
2711                                             "C3DBA76456363A10869622EAC2DD84EC" \
2712                                             "C5B8A74DAC4D09E03B5E0BE779F2DF61"));
2713 
2714     if (verbose != 0) {
2715         mbedtls_printf("  MPI test #4 (inv_mod): ");
2716     }
2717 
2718     if (mbedtls_mpi_cmp_mpi(&X, &U) != 0) {
2719         if (verbose != 0) {
2720             mbedtls_printf("failed\n");
2721         }
2722 
2723         ret = 1;
2724         goto cleanup;
2725     }
2726 
2727     if (verbose != 0) {
2728         mbedtls_printf("passed\n");
2729     }
2730 
2731     if (verbose != 0) {
2732         mbedtls_printf("  MPI test #5 (simple gcd): ");
2733     }
2734 
2735     for (i = 0; i < GCD_PAIR_COUNT; i++) {
2736         MBEDTLS_MPI_CHK(mbedtls_mpi_lset(&X, gcd_pairs[i][0]));
2737         MBEDTLS_MPI_CHK(mbedtls_mpi_lset(&Y, gcd_pairs[i][1]));
2738 
2739         MBEDTLS_MPI_CHK(mbedtls_mpi_gcd(&A, &X, &Y));
2740 
2741         if (mbedtls_mpi_cmp_int(&A, gcd_pairs[i][2]) != 0) {
2742             if (verbose != 0) {
2743                 mbedtls_printf("failed at %d\n", i);
2744             }
2745 
2746             ret = 1;
2747             goto cleanup;
2748         }
2749     }
2750 
2751     if (verbose != 0) {
2752         mbedtls_printf("passed\n");
2753     }
2754 
2755 cleanup:
2756 
2757     if (ret != 0 && verbose != 0) {
2758         mbedtls_printf("Unexpected error, return code = %08X\n", (unsigned int) ret);
2759     }
2760 
2761     mbedtls_mpi_free(&A); mbedtls_mpi_free(&E); mbedtls_mpi_free(&N); mbedtls_mpi_free(&X);
2762     mbedtls_mpi_free(&Y); mbedtls_mpi_free(&U); mbedtls_mpi_free(&V);
2763 
2764     if (verbose != 0) {
2765         mbedtls_printf("\n");
2766     }
2767 
2768     return ret;
2769 }
2770 
2771 #endif /* MBEDTLS_SELF_TEST */
2772 
2773 #endif /* MBEDTLS_BIGNUM_C */
2774