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