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