1 /********************************************************************/
2 /*                                                                  */
3 /*  big_gmp.c     Functions for bigInteger using the gmp library.   */
4 /*  Copyright (C) 1989 - 2019  Thomas Mertes                        */
5 /*                                                                  */
6 /*  This file is part of the Seed7 Runtime Library.                 */
7 /*                                                                  */
8 /*  The Seed7 Runtime Library is free software; you can             */
9 /*  redistribute it and/or modify it under the terms of the GNU     */
10 /*  Lesser General Public License as published by the Free Software */
11 /*  Foundation; either version 2.1 of the License, or (at your      */
12 /*  option) any later version.                                      */
13 /*                                                                  */
14 /*  The Seed7 Runtime Library is distributed in the hope that it    */
15 /*  will be useful, but WITHOUT ANY WARRANTY; without even the      */
16 /*  implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR */
17 /*  PURPOSE.  See the GNU Lesser General Public License for more    */
18 /*  details.                                                        */
19 /*                                                                  */
20 /*  You should have received a copy of the GNU Lesser General       */
21 /*  Public License along with this program; if not, write to the    */
22 /*  Free Software Foundation, Inc., 51 Franklin Street,             */
23 /*  Fifth Floor, Boston, MA  02110-1301, USA.                       */
24 /*                                                                  */
25 /*  Module: Seed7 Runtime Library                                   */
26 /*  File: seed7/src/big_gmp.c                                       */
27 /*  Changes: 2008, 2009, 2010, 2013 - 2019  Thomas Mertes           */
28 /*  Content: Functions for bigInteger using the gmp library.        */
29 /*                                                                  */
30 /********************************************************************/
31 
32 #define LOG_FUNCTIONS 0
33 #define VERBOSE_EXCEPTIONS 0
34 
35 #include "version.h"
36 
37 #if BIGINT_LIB == BIG_GMP_LIBRARY
38 #include "stdlib.h"
39 #include "stdio.h"
40 #include "string.h"
41 #include "limits.h"
42 #include "gmp.h"
43 
44 #include "common.h"
45 #include "data_rtl.h"
46 #include "heaputl.h"
47 #include "striutl.h"
48 #include "int_rtl.h"
49 #include "rtl_err.h"
50 
51 #undef EXTERN
52 #define EXTERN
53 #include "big_drv.h"
54 
55 
56 /* The integers in the gmp library are limited in size.    */
57 /* If the size of an integer gets too large the message:   */
58 /* gmp: overflow in mpz type                               */
59 /* is written and the program is terminated with abort().  */
60 /* The setting AVOID_OVERFLOW_IN_MPZ_TYPE activates code   */
61 /* to avoid the gmp overflow error for shift operations.   */
62 /* Instead the exception OVERFLOW_IN_MPZ_ERROR is raised.  */
63 #define AVOID_OVERFLOW_IN_MPZ_TYPE 1
64 #define GMP_MAX_LIMB_SIZE INT_MAX
65 #define OVERFLOW_IN_MPZ_ERROR MEMORY_ERROR
66 
67 #define USE_CUSTOM_ALLOCATION 0
68 
69 #define HEAP_ALLOC_BIG     (bigIntType) malloc(sizeof(__mpz_struct))
70 #define HEAP_FREE_BIG(var) free(var)
71 
72 #if WITH_BIGINT_FREELIST
73 
74 static freeListElemType flist = NULL;
75 static unsigned int flist_allowed = 100;
76 
77 #define POP_BIG_OK      flist != NULL
78 #define PUSH_BIG_OK     flist_allowed > 0
79 
80 #define POP_BIG(var)    {var = (bigIntType) flist; flist = flist->next; flist_allowed++; }
81 #define PUSH_BIG(var)   {((freeListElemType) var)->next = flist; flist = (freeListElemType) var; flist_allowed--; }
82 
83 #define ALLOC_BIG(var)  if (POP_BIG_OK) POP_BIG(var) else var = HEAP_ALLOC_BIG;
84 #define FREE_BIG(var)   if (PUSH_BIG_OK) PUSH_BIG(var) else HEAP_FREE_BIG(var);
85 
86 #else
87 
88 #define ALLOC_BIG(var)  var = HEAP_ALLOC_BIG
89 #define FREE_BIG(var)   HEAP_FREE_BIG(var)
90 
91 #endif
92 
93 
94 
95 #if USE_CUSTOM_ALLOCATION
alloc_func(size_t size)96 static void *alloc_func (size_t size)
97 
98   {
99     void *memory;
100 
101   /* alloc_func */
102     /* printf("alloc_func(" FMT_U_MEM ")\n", (memSizeType) size); */
103     memory = malloc(size);
104     if (unlikely(memory == NULL)) {
105       raise_error(MEMORY_ERROR);
106       printf(" ***** MEMORY_ERROR *****\n");
107       exit(0);
108     } /* if */
109     return memory;
110   } /* alloc_func */
111 
112 
113 
realloc_func(void * ptr,size_t old_size,size_t new_size)114 static void *realloc_func (void *ptr, size_t old_size, size_t new_size)
115 
116   {
117     void *memory;
118 
119   /* realloc_func */
120     /* printf("realloc_func(*, " FMT_U_MEM ", " FMT_U_MEM ")\n",
121        (memSizeType) old_size, (memSizeType) new_size); */
122     memory = realloc(ptr, new_size);
123     if (unlikely(memory == NULL)) {
124       raise_error(MEMORY_ERROR);
125       printf(" ***** MEMORY_ERROR *****\n");
126       exit(0);
127     } /* if */
128     return memory;
129   } /* realloc_func */
130 #endif
131 
132 
133 
134 /**
135  *  Setup bigInteger computations.
136  *  This function must be called before doing any bigInteger computations.
137  */
setupBig(void)138 void setupBig (void)
139 
140   { /* setupBig */
141 #if USE_CUSTOM_ALLOCATION
142     mp_set_memory_functions(alloc_func, realloc_func, NULL);
143 #endif
144   } /* setupBig */
145 
146 
147 
bigHexCStri(const const_bigIntType big1)148 cstriType bigHexCStri (const const_bigIntType big1)
149 
150   {
151     size_t sizeInBytes;
152     size_t bytesExported;
153     size_t pos;
154     int sign;
155     unsigned int carry;
156     const_cstriType stri_ptr;
157     ustriType buffer;
158     memSizeType charIndex;
159     size_t result_size;
160     cstriType result;
161 
162   /* bigHexCStri */
163     if (big1 != NULL) {
164       sizeInBytes = (mpz_sizeinbase(big1, 2) + CHAR_BIT - 1) / CHAR_BIT;
165       buffer = (ustriType) malloc(sizeInBytes);
166       if (unlikely(buffer == NULL)) {
167         raise_error(MEMORY_ERROR);
168         result = NULL;
169       } else {
170         mpz_export(buffer, &bytesExported, 1, 1, 0, 0, big1);
171         sign = mpz_sgn(big1);
172         if (sign < 0) {
173           carry = 1;
174           pos = bytesExported;
175           while (pos > 0) {
176             pos--;
177             carry += ~buffer[pos] & 0xFF;
178             buffer[pos] = (ucharType) (carry & 0xFF);
179             carry >>= CHAR_BIT;
180           } /* while */
181         } /* if */
182         result_size = 3 + sizeInBytes * 2;
183         if ((sign > 0 && buffer[0] > BYTE_MAX) ||
184             (sign < 0 && buffer[0] <= BYTE_MAX)) {
185           result_size += 2;
186         } /* if */
187         if (unlikely(!ALLOC_CSTRI(result, result_size))) {
188           raise_error(MEMORY_ERROR);
189         } else {
190           if (sign == 0) {
191             strcpy(result, "16#00");
192           } else {
193             strcpy(result, "16#");
194             charIndex = 3;
195             if (sign < 0) {
196               if (buffer[0] <= BYTE_MAX) {
197                 strcpy(&result[charIndex], "ff");
198                 charIndex += 2;
199               } /* if */
200             } else {
201               if (buffer[0] > BYTE_MAX) {
202                 strcpy(&result[charIndex], "00");
203                 charIndex += 2;
204               } /* if */
205             } /* for */
206             for (pos = 0; pos < bytesExported; pos++) {
207               sprintf(&result[charIndex], "%02x", buffer[pos]);
208               charIndex += 2;
209             } /* for */
210           } /* if */
211         } /* if */
212         free(buffer);
213       } /* if */
214     } else {
215       stri_ptr = " *NULL_BIGINT* ";
216       result_size = STRLEN(" *NULL_BIGINT* ");
217       if (unlikely(!ALLOC_CSTRI(result, result_size))) {
218         raise_error(MEMORY_ERROR);
219         return NULL;
220       } else {
221         strcpy(result, stri_ptr);
222       } /* if */
223     } /* if */
224     return result;
225   } /* bigHexCStri */
226 
227 
228 
229 /**
230  *  Compute the absolute value of a 'bigInteger' number.
231  *  @return the absolute value.
232  */
bigAbs(const const_bigIntType big1)233 bigIntType bigAbs (const const_bigIntType big1)
234 
235   {
236     bigIntType absoluteValue;
237 
238   /* bigAbs */
239     ALLOC_BIG(absoluteValue);
240     mpz_init(absoluteValue);
241     mpz_abs(absoluteValue, big1);
242     return absoluteValue;
243   } /* bigAbs */
244 
245 
246 
247 /**
248  *  Compute the absolute value of a 'bigInteger' number.
249  *  Big1 is assumed to be a temporary value which is reused.
250  *  @return the absolute value.
251  */
bigAbsTemp(bigIntType big1)252 bigIntType bigAbsTemp (bigIntType big1)
253 
254   { /* bigAbsTemp */
255     mpz_abs(big1, big1);
256     return big1;
257   } /* bigAbsTemp */
258 
259 
260 
261 /**
262  *  Add two 'bigInteger' numbers.
263  *  @return the sum of the two numbers.
264  */
bigAdd(const_bigIntType summand1,const_bigIntType summand2)265 bigIntType bigAdd (const_bigIntType summand1, const_bigIntType summand2)
266 
267   {
268     bigIntType sum;
269 
270   /* bigAdd */
271     ALLOC_BIG(sum);
272     mpz_init(sum);
273     mpz_add(sum, summand1, summand2);
274     return sum;
275   } /* bigAdd */
276 
277 
278 
279 /**
280  *  Increment a 'bigInteger' variable by a delta.
281  *  Adds delta to *big_variable.
282  *  @param delta The delta to be added to *big_variable.
283  */
bigAddAssign(bigIntType * const big_variable,const const_bigIntType delta)284 void bigAddAssign (bigIntType *const big_variable, const const_bigIntType delta)
285 
286   { /* bigAddAssign */
287     mpz_add(*big_variable, *big_variable, delta);
288   } /* bigAddAssign */
289 
290 
291 
292 /**
293  *  Increment a 'bigInteger' variable by a delta.
294  *  Adds delta to *big_variable.
295  *  @param delta The delta to be added to *big_variable.
296  *         Delta must be in the range of a long int.
297  */
bigAddAssignSignedDigit(bigIntType * const big_variable,const intType delta)298 void bigAddAssignSignedDigit (bigIntType *const big_variable, const intType delta)
299 
300   { /* bigAddAssignSignedDigit */
301     if (delta < 0) {
302       mpz_sub_ui(*big_variable, *big_variable, (unsigned long int) -delta);
303     } else {
304       mpz_add_ui(*big_variable, *big_variable, (unsigned long int) delta);
305     } /* if */
306   } /* bigAddAssignSignedDigit */
307 
308 
309 
310 /**
311  *  Add two 'bigInteger' numbers.
312  *  Summand1 is assumed to be a temporary value which is reused.
313  *  @return the sum of the two numbers in 'summand1'.
314  */
bigAddTemp(bigIntType summand1,const const_bigIntType summand2)315 bigIntType bigAddTemp (bigIntType summand1, const const_bigIntType summand2)
316 
317   { /* bigAddTemp */
318     mpz_add(summand1, summand1, summand2);
319     return summand1;
320   } /* bigAddTemp */
321 
322 
323 
bigAnd(const_bigIntType big1,const_bigIntType big2)324 bigIntType bigAnd (const_bigIntType big1, const_bigIntType big2)
325 
326   {
327     bigIntType result;
328 
329   /* bigAnd */
330     ALLOC_BIG(result);
331     mpz_init(result);
332     mpz_and(result, big1, big2);
333     return result;
334   } /* bigAnd */
335 
336 
337 
338 /**
339  *  Number of bits in the minimum two's-complement representation.
340  *  The high bits equivalent to the sign bit are not part of the
341  *  minimum two's-complement representation.
342  *  @return the number of bits.
343  *  @exception RANGE_ERROR The result does not fit into an integer.
344  */
bigBitLength(const const_bigIntType big1)345 intType bigBitLength (const const_bigIntType big1)
346 
347   {
348     int sign;
349     mpz_t help;
350     intType bitLength;
351 
352   /* bigBitLength */
353     sign = mpz_sgn(big1);
354     if (sign < 0) {
355       if (mpz_cmp_si(big1, -1) == 0) {
356         bitLength = 0;
357       } else {
358         mpz_init(help);
359         mpz_add_ui(help, big1, 1);
360         bitLength = (intType) mpz_sizeinbase(help, 2);
361         mpz_clear(help);
362       } /* if */
363     } else if (sign == 0) {
364       bitLength = 0;
365     } else {
366       bitLength = (intType) mpz_sizeinbase(big1, 2);
367     } /* if */
368     return bitLength;
369   } /* bigBitLength */
370 
371 
372 
373 /**
374  *  Compare two 'bigInteger' numbers.
375  *  @return -1, 0 or 1 if the first argument is considered to be
376  *          respectively less than, equal to, or greater than the
377  *          second.
378  */
bigCmp(const const_bigIntType big1,const const_bigIntType big2)379 intType bigCmp (const const_bigIntType big1, const const_bigIntType big2)
380 
381   {
382     intType signumValue;
383 
384   /* bigCmp */
385     signumValue = mpz_cmp(big1, big2);
386     if (signumValue < 0) {
387       signumValue = -1;
388     } else if (signumValue > 0) {
389       signumValue = 1;
390     } /* if */
391     return signumValue;
392   } /* bigCmp */
393 
394 
395 
396 /**
397  *  Reinterpret the generic parameters as bigIntType and call bigCmp.
398  *  Function pointers in C programs generated by the Seed7 compiler
399  *  may point to this function. This assures correct behaviour even
400  *  if sizeof(genericType) != sizeof(bigIntType).
401  *  @return -1, 0 or 1 if the first argument is considered to be
402  *          respectively less than, equal to, or greater than the
403  *          second.
404  */
bigCmpGeneric(const genericType value1,const genericType value2)405 intType bigCmpGeneric (const genericType value1, const genericType value2)
406 
407   { /* bigCmpGeneric */
408     return bigCmp(((const_rtlObjectType *) &value1)->value.bigIntValue,
409                   ((const_rtlObjectType *) &value2)->value.bigIntValue);
410   } /* bigCmpGeneric */
411 
412 
413 
414 /**
415  *  Compare 'big1' with the bigdigit 'number'.
416  *  The range of 'number' is restricted and it is the job of the
417  *  compiler to assure that 'number' is within the allowed range.
418  *  @param number Number that must be in the range of
419  *         a long int.
420  *  @return -1, 0 or 1 if the first argument is considered to be
421  *          respectively less than, equal to, or greater than the
422  *          second.
423  */
bigCmpSignedDigit(const const_bigIntType big1,intType number)424 intType bigCmpSignedDigit (const const_bigIntType big1, intType number)
425 
426   {
427     intType signumValue;
428 
429   /* bigCmpSignedDigit */
430     signumValue = mpz_cmp_si(big1, number);
431     if (signumValue < 0) {
432       signumValue = -1;
433     } else if (signumValue > 0) {
434       signumValue = 1;
435     } /* if */
436     return signumValue;
437   } /* bigCmpSignedDigit */
438 
439 
440 
441 /**
442  *  Assign source to *dest.
443  *  A copy function assumes that *dest contains a legal value.
444  */
bigCpy(bigIntType * const dest,const const_bigIntType source)445 void bigCpy (bigIntType *const dest, const const_bigIntType source)
446 
447   { /* bigCpy */
448     mpz_set(*dest, source);
449   } /* bigCpy */
450 
451 
452 
453 /**
454  *  Reinterpret the generic parameters as bigIntType and call bigCpy.
455  *  Function pointers in C programs generated by the Seed7 compiler
456  *  may point to this function. This assures correct behaviour even
457  *  if sizeof(genericType) != sizeof(bigIntType).
458  */
bigCpyGeneric(genericType * const dest,const genericType source)459 void bigCpyGeneric (genericType *const dest, const genericType source)
460 
461   { /* bigCpyGeneric */
462     bigCpy(&((rtlObjectType *) dest)->value.bigIntValue,
463            ((const_rtlObjectType *) &source)->value.bigIntValue);
464   } /* bigCpyGeneric */
465 
466 
467 
468 /**
469  *  Return a copy of source, that can be assigned to a new destination.
470  *  It is assumed that the destination of the assignment is undefined.
471  *  Create functions can be used to initialize Seed7 constants.
472  *  @return a copy of source.
473  */
bigCreate(const const_bigIntType source)474 bigIntType bigCreate (const const_bigIntType source)
475 
476   {
477     bigIntType result;
478 
479   /* bigCreate */
480     ALLOC_BIG(result);
481     mpz_init_set(result, source);
482     return result;
483   } /* bigCreate */
484 
485 
486 
487 /**
488  *  Generic Create function to be used via function pointers.
489  *  Function pointers in C programs generated by the Seed7 compiler
490  *  may point to this function. This assures correct behaviour even
491  *  if sizeof(genericType) != sizeof(bigIntType).
492  */
bigCreateGeneric(const genericType source)493 genericType bigCreateGeneric (const genericType source)
494 
495   {
496     rtlObjectType result;
497 
498   /* bigCreateGeneric */
499     INIT_GENERIC_PTR(result.value.genericValue);
500     result.value.bigIntValue =
501         bigCreate(((const_rtlObjectType *) &source)->value.bigIntValue);
502     return result.value.genericValue;
503   } /* bigCreateGeneric */
504 
505 
506 
507 /**
508  *  Decrement a 'bigInteger' variable.
509  *  Decrements *big_variable by 1.
510  */
bigDecr(bigIntType * const big_variable)511 void bigDecr (bigIntType *const big_variable)
512 
513   { /* bigDecr */
514     mpz_sub_ui(*big_variable, *big_variable, 1);
515   } /* bigDecr */
516 
517 
518 
519 /**
520  *  Free the memory referred by 'old_bigint'.
521  *  After bigDestr is left 'old_bigint' refers to not existing memory.
522  *  The memory where 'old_bigint' is stored can be freed afterwards.
523  */
bigDestr(const const_bigIntType old_bigint)524 void bigDestr (const const_bigIntType old_bigint)
525 
526   { /* bigDestr */
527     if (old_bigint != NULL) {
528       mpz_clear((bigIntType) old_bigint);
529       FREE_BIG((bigIntType) old_bigint);
530     } /* if */
531   } /* bigDestr */
532 
533 
534 
535 /**
536  *  Generic Destr function to be used via function pointers.
537  *  Function pointers in C programs generated by the Seed7 compiler
538  *  may point to this function. This assures correct behaviour even
539  *  if sizeof(genericType) != sizeof(bigIntType).
540  */
bigDestrGeneric(const genericType old_value)541 void bigDestrGeneric (const genericType old_value)
542 
543   { /* bigDestrGeneric */
544     bigDestr(((const_rtlObjectType *) &old_value)->value.bigIntValue);
545   } /* bigDestrGeneric */
546 
547 
548 
549 /**
550  *  Integer division truncated towards zero.
551  *  The remainder of this division is computed with bigRem.
552  *  @return the quotient of the integer division.
553  *  @exception NUMERIC_ERROR If a division by zero occurs.
554  */
bigDiv(const const_bigIntType dividend,const const_bigIntType divisor)555 bigIntType bigDiv (const const_bigIntType dividend, const const_bigIntType divisor)
556 
557   {
558     bigIntType quotient;
559 
560   /* bigDiv */
561     logFunction(printf("bigDiv(%s,", bigHexCStri(dividend));
562                 printf("%s)\n", bigHexCStri(divisor)););
563     if (unlikely(mpz_sgn(divisor) == 0)) {
564       logError(printf("bigDiv(%s, %s): Division by zero.\n",
565                       bigHexCStri(dividend), bigHexCStri(divisor)););
566       raise_error(NUMERIC_ERROR);
567       quotient = NULL;
568     } else {
569       ALLOC_BIG(quotient);
570       mpz_init(quotient);
571       mpz_tdiv_q(quotient, dividend, divisor);
572     } /* if */
573     logFunction(printf("bigDiv --> %s\n", bigHexCStri(quotient)););
574     return quotient;
575   } /* bigDiv */
576 
577 
578 
579 /**
580  *  Integer division truncated towards zero.
581  *  The memory for the quotient is requested and the normalized
582  *  quotient is returned. The memory for the remainder is
583  *  requested and the normalized remainder is assigned to
584  *  *remainderAddr.
585  *  @return the quotient of the integer division.
586  *  @exception NUMERIC_ERROR If a division by zero occurs.
587  */
bigDivRem(const const_bigIntType dividend,const const_bigIntType divisor,bigIntType * remainderAddr)588 bigIntType bigDivRem (const const_bigIntType dividend, const const_bigIntType divisor,
589     bigIntType *remainderAddr)
590 
591   {
592     bigIntType quotient;
593 
594   /* bigDivRem */
595     logFunction(printf("bigDivRem(%s,", bigHexCStri(dividend));
596                 printf("%s)\n", bigHexCStri(divisor)););
597     if (unlikely(mpz_sgn(divisor) == 0)) {
598       logError(printf("bigDiv(%s, %s): Division by zero.\n",
599                       bigHexCStri(dividend), bigHexCStri(divisor)););
600       raise_error(NUMERIC_ERROR);
601       quotient = NULL;
602     } else {
603       ALLOC_BIG(quotient);
604       mpz_init(quotient);
605       ALLOC_BIG(*remainderAddr);
606       mpz_init(*remainderAddr);
607       mpz_tdiv_qr(quotient, *remainderAddr, dividend, divisor);
608     } /* if */
609     logFunction(printf("bigDivRem --> %s", bigHexCStri(quotient));
610                 printf(" (%s)\n", bigHexCStri(*remainderAddr)););
611     return quotient;
612   } /* bigDivRem */
613 
614 
615 
616 /**
617  *  Check if two 'bigInteger' numbers are equal.
618  *  @return TRUE if both numbers are equal,
619  *          FALSE otherwise.
620  */
bigEq(const const_bigIntType big1,const const_bigIntType big2)621 boolType bigEq (const const_bigIntType big1, const const_bigIntType big2)
622 
623   { /* bigEq */
624     return mpz_cmp(big1, big2) == 0;
625   } /* bigEq */
626 
627 
628 
629 /**
630  *  Check if 'big1' is equal to the bigdigit 'number'.
631  *  The range of 'number' is restricted and it is the job of the
632  *  compiler to assure that 'number' is within the allowed range.
633  *  @param number Number that must be in the range of
634  *         a long int.
635  *  @return TRUE if 'big1' and 'number' are equal,
636  *          FALSE otherwise.
637  */
bigEqSignedDigit(const const_bigIntType big1,intType number)638 boolType bigEqSignedDigit (const const_bigIntType big1, intType number)
639 
640   { /* bigEqSignedDigit */
641     return mpz_cmp_si(big1, number) == 0;
642   } /* bigEqSignedDigit */
643 
644 
645 
646 /**
647  *  Convert a byte buffer (interpreted as big-endian) to a bigInteger.
648  *  @param size Size of the byte buffer to be converted (in bytes).
649  *  @param buffer Byte buffer to be converted. The bytes are interpreted
650  *         as binary big-endian representation with a base of 256.
651  *  @param isSigned Defines if 'buffer' is interpreted as signed value.
652  *         If 'isSigned' is TRUE the twos-complement representation
653  *         is used. In this case the result is negative if the most
654  *         significant byte (the first byte) has an ordinal > BYTE_MAX (=127).
655  *  @return a bigInteger created from the big-endian bytes.
656  */
bigFromByteBufferBe(const memSizeType size,const const_ustriType buffer,const boolType isSigned)657 bigIntType bigFromByteBufferBe (const memSizeType size,
658     const const_ustriType buffer, const boolType isSigned)
659 
660   {
661     memSizeType pos;
662     ustriType negated_buffer;
663     unsigned int carry;
664     bigIntType result;
665 
666   /* bigFromByteBufferBe */
667     logFunction(printf("bigFromByteBufferBe(" FMT_U_MEM ", 0x" FMT_X_MEM ", %d)\n",
668                        size, (memSizeType) buffer, isSigned););
669     ALLOC_BIG(result);
670     mpz_init(result);
671     if (isSigned && size != 0 && buffer[0] > BYTE_MAX) {
672       negated_buffer = (ustriType) malloc(size);
673       carry = 1;
674       pos = size;
675       while (pos > 0) {
676         pos--;
677         carry += ~buffer[pos] & 0xFF;
678         negated_buffer[pos] = (ucharType) (carry & 0xFF);
679         carry >>= CHAR_BIT;
680       } /* for */
681       mpz_import(result, (size_t) size, 1, 1, 0, 0, negated_buffer);
682       free(negated_buffer);
683       mpz_neg(result, result);
684     } else {
685       mpz_import(result, (size_t) size, 1, 1, 0, 0, buffer);
686     } /* if */
687     logFunction(printf("bigFromByteBufferBe --> %s\n", bigHexCStri(result)););
688     return result;
689   } /* bigFromByteBufferBe */
690 
691 
692 
693 /**
694  *  Convert a byte buffer (interpreted as little-endian) to a bigInteger.
695  *  @param size Size of the byte buffer to be converted (in bytes).
696  *  @param buffer Byte buffer to be converted. The bytes are interpreted
697  *         as binary little-endian representation with a base of 256.
698  *  @param isSigned Defines if 'buffer' is interpreted as signed value.
699  *         If 'isSigned' is TRUE the twos-complement representation
700  *         is used. In this case the result is negative if the most
701  *         significant byte (the last byte) has an ordinal > BYTE_MAX (=127).
702  *  @return a bigInteger created from the little-endian bytes.
703  */
bigFromByteBufferLe(const memSizeType size,const const_ustriType buffer,const boolType isSigned)704 bigIntType bigFromByteBufferLe (const memSizeType size,
705     const const_ustriType buffer, const boolType isSigned)
706 
707   {
708     memSizeType pos;
709     ustriType negated_buffer;
710     unsigned int carry;
711     bigIntType result;
712 
713   /* bigFromByteBufferLe */
714     logFunction(printf("bigFromByteBufferLe(" FMT_U_MEM ", 0x" FMT_X_MEM ", %d)\n",
715                        size, (memSizeType) buffer, isSigned););
716     ALLOC_BIG(result);
717     mpz_init(result);
718     if (isSigned && size != 0 && buffer[size - 1] > BYTE_MAX) {
719       negated_buffer = (ustriType) malloc(size);
720       carry = 1;
721       pos = 0;
722       while (pos < size) {
723         carry += ~buffer[pos] & 0xFF;
724         negated_buffer[pos] = (ucharType) (carry & 0xFF);
725         carry >>= CHAR_BIT;
726         pos++;
727       } /* for */
728       mpz_import(result, (size_t) size, -1, 1, 0, 0, negated_buffer);
729       free(negated_buffer);
730       mpz_neg(result, result);
731     } else {
732       mpz_import(result, (size_t) size, -1, 1, 0, 0, buffer);
733     } /* if */
734     logFunction(printf("bigFromByteBufferLe --> %s\n", bigHexCStri(result)););
735     return result;
736   } /* bigFromByteBufferLe */
737 
738 
739 
740 /**
741  *  Convert a bstring (interpreted as big-endian) to a bigInteger.
742  *  @param bstri Bstring to be converted. The bytes are interpreted
743  *         as binary big-endian representation with a base of 256.
744  *  @param isSigned Defines if 'bstri' is interpreted as signed value.
745  *         If 'isSigned' is TRUE the twos-complement representation
746  *         is used. In this case the result is negative if the most
747  *         significant byte (the first byte) has an ordinal > BYTE_MAX (=127).
748  *  @return a bigInteger created from the big-endian bytes.
749  */
bigFromBStriBe(const const_bstriType bstri,const boolType isSigned)750 bigIntType bigFromBStriBe (const const_bstriType bstri, const boolType isSigned)
751 
752   { /* bigFromBStriBe */
753     logFunction(printf("bigFromBStriBe(\"%s\", %d)\n",
754                        bstriAsUnquotedCStri(bstri), isSigned););
755     return bigFromByteBufferBe(bstri->size, bstri->mem, isSigned);
756   } /* bigFromBStriBe */
757 
758 
759 
760 /**
761  *  Convert a bstring (interpreted as little-endian) to a bigInteger.
762  *  @param bstri Bstring to be converted. The bytes are interpreted
763  *         as binary little-endian representation with a base of 256.
764  *  @param isSigned Defines if 'bstri' is interpreted as signed value.
765  *         If 'isSigned' is TRUE the twos-complement representation
766  *         is used. In this case the result is negative if the most
767  *         significant byte (the last byte) has an ordinal > BYTE_MAX (=127).
768  *  @return a bigInteger created from the little-endian bytes.
769  */
bigFromBStriLe(const const_bstriType bstri,const boolType isSigned)770 bigIntType bigFromBStriLe (const const_bstriType bstri, const boolType isSigned)
771 
772   { /* bigFromBStriLe */
773     logFunction(printf("bigFromBStriLe(\"%s\", %d)\n",
774                        bstriAsUnquotedCStri(bstri), isSigned););
775     return bigFromByteBufferLe(bstri->size, bstri->mem, isSigned);
776   } /* bigFromBStriLe */
777 
778 
779 
780 /**
781  *  Convert an int32Type number to 'bigInteger'.
782  *  @return the bigInteger result of the conversion.
783  */
bigFromInt32(int32Type number)784 bigIntType bigFromInt32 (int32Type number)
785 
786   {
787     bigIntType result;
788 
789   /* bigFromInt32 */
790     ALLOC_BIG(result);
791     mpz_init_set_si(result, number);
792     return result;
793   } /* bigFromInt32 */
794 
795 
796 
797 #ifdef INT64TYPE
798 /**
799  *  Convert an int64Type number to 'bigInteger'.
800  *  @return the bigInteger result of the conversion.
801  */
bigFromInt64(int64Type number)802 bigIntType bigFromInt64 (int64Type number)
803 
804   {
805     bigIntType result;
806 
807   /* bigFromInt64 */
808     ALLOC_BIG(result);
809 #if LONG_SIZE == 64
810     mpz_init_set_si(result, number);
811 #else
812     {
813       mpz_t help;
814 
815       mpz_init_set_si(result, (int32Type) ((number >> 32) & 0xFFFFFFFF));
816       mpz_mul_2exp(result, result, 32);
817       mpz_init_set_ui(help, (uint32Type) (number & 0xFFFFFFFF));
818       mpz_ior(result, result, help);
819       mpz_clear(help);
820     }
821 #endif
822     return result;
823   } /* bigFromInt64 */
824 #endif
825 
826 
827 
828 /**
829  *  Convert an uint32Type number to 'bigInteger'.
830  *  @return the bigInteger result of the conversion.
831  */
bigFromUInt32(uint32Type number)832 bigIntType bigFromUInt32 (uint32Type number)
833 
834   {
835     bigIntType result;
836 
837   /* bigFromUInt32 */
838     ALLOC_BIG(result);
839     mpz_init_set_ui(result, number);
840     return result;
841   } /* bigFromUInt32 */
842 
843 
844 
845 #ifdef INT64TYPE
846 /**
847  *  Convert an uint64Type number to 'bigInteger'.
848  *  @return the bigInteger result of the conversion.
849  */
bigFromUInt64(uint64Type number)850 bigIntType bigFromUInt64 (uint64Type number)
851 
852   {
853     mpz_t help;
854     bigIntType result;
855 
856   /* bigFromUInt64 */
857     logFunction(printf("bigFromUInt64(" FMT_U64 ")\n", number););
858     ALLOC_BIG(result);
859     mpz_init_set_ui(result, (uint32Type) ((number >> 32) & 0xFFFFFFFF));
860     mpz_mul_2exp(result, result, 32);
861     mpz_init_set_ui(help, (uint32Type) (number & 0xFFFFFFFF));
862     mpz_ior(result, result, help);
863     mpz_clear(help);
864     logFunction(printf("bigFromUInt64 --> %s\n", bigHexCStri(result)););
865     return result;
866   } /* bigFromUInt64 */
867 #endif
868 
869 
870 
871 /**
872  *  Compute the greatest common divisor of two 'bigInteger' numbers.
873  *  @return the greatest common divisor of the two numbers.
874  *          The greatest common divisor is positive or zero.
875  */
bigGcd(const const_bigIntType big1,const const_bigIntType big2)876 bigIntType bigGcd (const const_bigIntType big1,
877     const const_bigIntType big2)
878 
879   {
880     bigIntType result;
881 
882   /* bigGcd */
883     ALLOC_BIG(result);
884     mpz_init(result);
885     mpz_gcd(result, big1, big2);
886     return result;
887   } /* bigGcd */
888 
889 
890 
891 /**
892  *  Compute the hash value of a 'bigInteger' number.
893  *  @return the hash value.
894  */
bigHashCode(const const_bigIntType big1)895 intType bigHashCode (const const_bigIntType big1)
896 
897   {
898     size_t size;
899     intType result;
900 
901   /* bigHashCode */
902     size = mpz_size(big1);
903     if (size == 0) {
904       result = 0;
905     } else {
906       result = (intType) (mpz_getlimbn(big1, 0) << 5 ^ size << 3 ^ mpz_getlimbn(big1, size - 1));
907     } /* if */
908     return result;
909   } /* bigHashCode */
910 
911 
912 
913 /**
914  *  Increment a 'bigInteger' variable.
915  *  Increments *big_variable by 1.
916  */
bigIncr(bigIntType * const big_variable)917 void bigIncr (bigIntType *const big_variable)
918 
919   { /* bigIncr */
920     mpz_add_ui(*big_variable, *big_variable, 1);
921   } /* bigIncr */
922 
923 
924 
925 /**
926  *  Compute the exponentiation of a 'bigInteger' base with an integer exponent.
927  *  @return the result of the exponentiation.
928  *  @exception NUMERIC_ERROR If the exponent is negative.
929  */
bigIPow(const const_bigIntType base,intType exponent)930 bigIntType bigIPow (const const_bigIntType base, intType exponent)
931 
932   {
933     bigIntType power;
934 
935   /* bigIPow */
936     logFunction(printf("bigIPow(%s, " FMT_D ")\n",
937                        bigHexCStri(base), exponent););
938     if (unlikely(exponent < 0)) {
939       logError(printf("bigIPow(%s, " FMT_D "): "
940                       "Exponent is negative.\n",
941                       bigHexCStri(base), exponent););
942       raise_error(NUMERIC_ERROR);
943       power = NULL;
944     } else {
945       ALLOC_BIG(power);
946       mpz_init(power);
947       mpz_pow_ui(power, base, (uintType) exponent);
948     } /* if */
949     logFunction(printf("bigIPow --> %s\n", bigHexCStri(power)););
950     return power;
951   } /* bigIPow */
952 
953 
954 
955 /**
956  *  Compute the exponentiation of a bigdigit base with an integer exponent.
957  *  @param base Base that must be in the range of a long int.
958  *  @return the result of the exponentiation.
959  *  @exception NUMERIC_ERROR If the exponent is negative.
960  */
bigIPowSignedDigit(intType base,intType exponent)961 bigIntType bigIPowSignedDigit (intType base, intType exponent)
962 
963   {
964     bigIntType power;
965 
966   /* bigIPowSignedDigit */
967     if (unlikely(exponent < 0)) {
968       logError(printf("bigIPowSignedDigit(" FMT_D ", " FMT_D "): "
969                       "Exponent is negative.\n",
970                       base, exponent););
971       raise_error(NUMERIC_ERROR);
972       power = NULL;
973     } else {
974       ALLOC_BIG(power);
975       mpz_init_set_si(power, base);
976       mpz_pow_ui(power, power, (uintType) exponent);
977     } /* if */
978     return power;
979   } /* bigIPowSignedDigit */
980 
981 
982 
983 /**
984  *  Compute the truncated base 10 logarithm of a 'bigInteger' number.
985  *  The definition of 'log10' is extended by defining log10(0) = -1_.
986  *  @return the truncated base 10 logarithm.
987  *  @exception NUMERIC_ERROR The number is negative.
988  */
bigLog10(const const_bigIntType big1)989 bigIntType bigLog10 (const const_bigIntType big1)
990 
991   {
992     int sign;
993     char *cstri;
994     bigIntType logarithm;
995 
996   /* bigLog10 */
997     sign = mpz_sgn(big1);
998     if (unlikely(sign < 0)) {
999       logError(printf("bigLog10(%s): Number is negative.\n",
1000                       bigHexCStri(big1)););
1001       raise_error(NUMERIC_ERROR);
1002       logarithm = NULL;
1003     } else if (sign == 0) {
1004       ALLOC_BIG(logarithm);
1005       mpz_init_set_si(logarithm, -1);
1006     } else {
1007       cstri = mpz_get_str(NULL, 10, big1);
1008       ALLOC_BIG(logarithm);
1009       mpz_init_set_ui(logarithm, strlen(cstri) - 1);
1010       free(cstri);
1011     } /* if */
1012     return logarithm;
1013   } /* bigLog10 */
1014 
1015 
1016 
1017 /**
1018  *  Compute the truncated base 2 logarithm of a 'bigInteger' number.
1019  *  The definition of 'log2' is extended by defining log2(0) = -1_.
1020  *  @return the truncated base 2 logarithm.
1021  *  @exception NUMERIC_ERROR The number is negative.
1022  */
bigLog2(const const_bigIntType big1)1023 bigIntType bigLog2 (const const_bigIntType big1)
1024 
1025   {
1026     int sign;
1027     bigIntType logarithm;
1028 
1029   /* bigLog2 */
1030     sign = mpz_sgn(big1);
1031     if (unlikely(sign < 0)) {
1032       logError(printf("bigLog2(%s): Number is negative.\n",
1033                       bigHexCStri(big1)););
1034       raise_error(NUMERIC_ERROR);
1035       logarithm = NULL;
1036     } else if (sign == 0) {
1037       ALLOC_BIG(logarithm);
1038       mpz_init_set_si(logarithm, -1);
1039     } else {
1040       ALLOC_BIG(logarithm);
1041       mpz_init_set_ui(logarithm, mpz_sizeinbase(big1, 2) - 1);
1042     } /* if */
1043     return logarithm;
1044   } /* bigLog2 */
1045 
1046 
1047 
1048 /**
1049  *  Create a number from the lower bits of big1.
1050  *  This corresponds to the modulo if the dividend is a power of two:
1051  *   bigLowerBits(big1, bits)  corresponds to  big1 mod (2_ ** bits)
1052  *  @param bits Number of lower bits to select from big1.
1053  *  @return a number in the range 0 .. pred(2_ ** bits).
1054  *  @exception NUMERIC_ERROR The number of bits is negative.
1055  */
bigLowerBits(const const_bigIntType big1,const intType bits)1056 bigIntType bigLowerBits (const const_bigIntType big1, const intType bits)
1057 
1058   {
1059     bigIntType result;
1060 
1061   /* bigLowerBits */
1062     if (unlikely(bits < 0)) {
1063       logError(printf("bigLowerBits(%s, " FMT_D "): "
1064                       "Number of bits is negative.\n",
1065                       bigHexCStri(big1), bits););
1066       raise_error(NUMERIC_ERROR);
1067       result = NULL;
1068     } else {
1069       ALLOC_BIG(result);
1070       mpz_init(result);
1071       mpz_fdiv_r_2exp(result, big1, (uintType) bits);
1072     } /* if */
1073     return result;
1074   } /* bigLowerBits */
1075 
1076 
1077 
1078 /**
1079  *  Create a number from the lower bits of big1.
1080  *  Big1 is assumed to be a temporary value which is reused.
1081  *  This corresponds to the modulo if the dividend is a power of two:
1082  *   bigLowerBits(big1, bits)  corresponds to  big1 mod (2_ ** bits)
1083  *  @param bits Number of lower bits to select from big1.
1084  *  @return a number in the range 0 .. pred(2_ ** bits).
1085  *  @exception NUMERIC_ERROR The number of bits is negative.
1086  */
bigLowerBitsTemp(const bigIntType big1,const intType bits)1087 bigIntType bigLowerBitsTemp (const bigIntType big1, const intType bits)
1088 
1089   { /* bigLowerBitsTemp */
1090     if (unlikely(bits < 0)) {
1091       FREE_BIG(big1);
1092       logError(printf("bigLowerBitsTemp(%s, " FMT_D "): "
1093                       "Number of bits is negative.\n",
1094                       bigHexCStri(big1), bits););
1095       raise_error(NUMERIC_ERROR);
1096       return NULL;
1097     } else {
1098       mpz_fdiv_r_2exp(big1, big1, (uintType) bits);
1099     } /* if */
1100     return big1;
1101   } /* bigLowerBitsTemp */
1102 
1103 
1104 
bigLowerBits64(const const_bigIntType big1)1105 uint64Type bigLowerBits64 (const const_bigIntType big1)
1106 
1107   {
1108     mpz_t mod64;
1109     uint64Type result;
1110 
1111   /* bigLowerBits64 */
1112     logFunction(printf("bigLowerBits64(%s)\n", bigHexCStri(big1)););
1113     mpz_init(mod64);
1114     mpz_fdiv_r_2exp(mod64, big1, 64);
1115 #if LONG_SIZE == 64
1116     result = mpz_get_ui(mod64);
1117 #else
1118     {
1119       mpz_t help;
1120 
1121       mpz_init_set(help, mod64);
1122       mpz_fdiv_q_2exp(help, help, 32);
1123       result = (uint64Type) (mpz_get_ui(help) & 0xFFFFFFFF) << 32 |
1124                (uint64Type) (mpz_get_ui(mod64) & 0xFFFFFFFF);
1125       mpz_clear(help);
1126     }
1127 #endif
1128     logFunction(printf("bigLowerBits64(%s) --> " FMT_U64 "\n",
1129                        bigHexCStri(big1), result););
1130     return result;
1131   } /* bigLowerBits64 */
1132 
1133 
1134 
1135 /**
1136  *  Index of the lowest-order one bit.
1137  *  For A <> 0 this is equal to the number of lowest-order zero bits.
1138  *  @return the number of lowest-order zero bits or -1 for lowestSetBit(0).
1139  */
bigLowestSetBit(const const_bigIntType big1)1140 intType bigLowestSetBit (const const_bigIntType big1)
1141 
1142   { /* bigLowestSetBit */
1143     return (intType) mpz_scan1(big1, 0);
1144   } /* bigLowestSetBit */
1145 
1146 
1147 
1148 /**
1149  *  Shift a 'bigInteger' number left by lshift bits.
1150  *  If lshift is negative a right shift is done instead.
1151  *  A << B is equivalent to A * 2_ ** B if B >= 0 holds.
1152  *  A << B is equivalent to A mdiv 2_ ** -B if B < 0 holds.
1153  *  @return the left shifted number.
1154  */
bigLShift(const const_bigIntType big1,const intType lshift)1155 bigIntType bigLShift (const const_bigIntType big1, const intType lshift)
1156 
1157   {
1158     bigIntType result;
1159 
1160   /* bigLShift */
1161     logFunction(printf("bigLShift(%s, " FMT_D ")\n",
1162                        bigHexCStri(big1), lshift););
1163     ALLOC_BIG(result);
1164     mpz_init(result);
1165     if (lshift <= 0) {
1166       /* The unsigned value is negated to avoid a signed integer */
1167       /* overflow if the smallest signed integer is negated.     */
1168       mpz_fdiv_q_2exp(result, big1, -(uintType) lshift);
1169 #if AVOID_OVERFLOW_IN_MPZ_TYPE
1170     } else if (unlikely((uintType) lshift / GMP_LIMB_BITS + 1 >
1171                         GMP_MAX_LIMB_SIZE - mpz_size(big1) &&
1172                         mpz_sgn(big1) != 0)) {
1173       mpz_clear(result);
1174       FREE_BIG(result);
1175       raise_error(OVERFLOW_IN_MPZ_ERROR);
1176       result = NULL;
1177 #endif
1178     } else {
1179       mpz_mul_2exp(result, big1, (uintType) lshift);
1180     } /* if */
1181     logFunction(printf("bigLShift --> %s\n", bigHexCStri(result)););
1182     return result;
1183   } /* bigLShift */
1184 
1185 
1186 
1187 /**
1188  *  Shift a number left by lshift bits and assign the result back to number.
1189  *  If lshift is negative a right shift is done instead.
1190  */
bigLShiftAssign(bigIntType * const big_variable,intType lshift)1191 void bigLShiftAssign (bigIntType *const big_variable, intType lshift)
1192 
1193   {
1194     bigIntType big1;
1195 
1196   /* bigLShiftAssign */
1197     big1 = *big_variable;
1198     logFunction(printf("bigLShiftAssign(%s, " FMT_D ")\n",
1199                        bigHexCStri(big1), lshift););
1200     if (lshift < 0) {
1201       /* The unsigned value is negated to avoid a signed integer */
1202       /* overflow if the smallest signed integer is negated.     */
1203       mpz_fdiv_q_2exp(big1, big1, -(uintType) lshift);
1204     } else if (lshift == 0) {
1205       /* do nothing */
1206 #if AVOID_OVERFLOW_IN_MPZ_TYPE
1207     } else if (unlikely((uintType) lshift / GMP_LIMB_BITS + 1 >
1208                         GMP_MAX_LIMB_SIZE - mpz_size(big1) &&
1209                         mpz_sgn(big1) != 0)) {
1210       raise_error(OVERFLOW_IN_MPZ_ERROR);
1211 #endif
1212     } else {
1213       mpz_mul_2exp(big1, big1, (uintType) lshift);
1214     } /* if */
1215     logFunction(printf("bigLShiftAssign --> %s\n", bigHexCStri(big1)););
1216   } /* bigLShiftAssign */
1217 
1218 
1219 
1220 /**
1221  * Shift one left by 'lshift' bits.
1222  * If 'lshift' is positive or zero this corresponds to
1223  * the computation of a power of two:
1224  *  bigLShiftOne(lshift)  corresponds to  2_ ** lshift
1225  * If 'lshift' is negative the result is zero.
1226  * @return one shifted left by 'lshift'.
1227  */
bigLShiftOne(const intType lshift)1228 bigIntType bigLShiftOne (const intType lshift)
1229 
1230   {
1231     mpz_t one;
1232     bigIntType result;
1233 
1234   /* bigLShiftOne */
1235     ALLOC_BIG(result);
1236     if (lshift < 0) {
1237       mpz_init_set_ui(result, 0);
1238 #if AVOID_OVERFLOW_IN_MPZ_TYPE
1239     } else if (unlikely((uintType) lshift / GMP_LIMB_BITS + 1 >
1240                         GMP_MAX_LIMB_SIZE - 1)) {
1241       FREE_BIG(result);
1242       raise_error(OVERFLOW_IN_MPZ_ERROR);
1243       result = NULL;
1244 #endif
1245     } else {
1246       mpz_init(result);
1247       mpz_init_set_ui(one, 1);
1248       mpz_mul_2exp(result, one, (uintType) lshift);
1249       mpz_clear(one);
1250     } /* if */
1251     return result;
1252   } /* bigLShiftOne */
1253 
1254 
1255 
1256 /**
1257  *  Exponentiation if the base is a power of two.
1258  *  @param log2base Logarithm of the actual base ( =log2(base) )
1259  *  @return (2 ** log2base) ** exponent
1260  *  @exception NUMERIC_ERROR If log2base or exponent is negative.
1261  */
bigLog2BaseIPow(const intType log2base,const intType exponent)1262 bigIntType bigLog2BaseIPow (const intType log2base, const intType exponent)
1263 
1264   {
1265     uintType high_shift;
1266     uintType low_shift;
1267     bigIntType power;
1268 
1269   /* bigLog2BaseIPow */
1270     if (unlikely(log2base < 0 || exponent < 0)) {
1271       logError(printf("bigLog2BaseIPow(" FMT_D ", " FMT_D "): "
1272                       "Log2base or exponent is negative.\n",
1273                       log2base, exponent););
1274       raise_error(NUMERIC_ERROR);
1275       power = NULL;
1276     } else if (likely(log2base == 1)) {
1277       power = bigLShiftOne(exponent);
1278     } else if (log2base <= 10 && exponent <= MAX_DIV_10) {
1279       power = bigLShiftOne(log2base * exponent);
1280     } else {
1281       low_shift = uintMult((uintType) log2base, (uintType) exponent, &high_shift);
1282       if (unlikely(high_shift != 0 || (intType) low_shift < 0)) {
1283         raise_error(MEMORY_ERROR);
1284         power = NULL;
1285       } else {
1286         power = bigLShiftOne((intType) low_shift);
1287       } /* if */
1288     } /* if */
1289     return power;
1290   } /* bigLog2BaseIPow */
1291 
1292 
1293 
1294 /**
1295  *  Integer division truncated towards negative infinity.
1296  *  The modulo (remainder) of this division is computed with bigMod.
1297  *  Therefore this division is called modulo division (MDiv).
1298  *  @return the quotient of the integer division.
1299  *  @exception NUMERIC_ERROR If a division by zero occurs.
1300  */
bigMDiv(const const_bigIntType dividend,const const_bigIntType divisor)1301 bigIntType bigMDiv (const const_bigIntType dividend, const const_bigIntType divisor)
1302 
1303   {
1304     bigIntType quotient;
1305 
1306   /* bigMDiv */
1307     logFunction(printf("bigMDiv(%s,", bigHexCStri(dividend));
1308                 printf("%s)\n", bigHexCStri(divisor)););
1309     if (unlikely(mpz_sgn(divisor) == 0)) {
1310       logError(printf("bigMDiv(%s, %s): Division by zero.\n",
1311                       bigHexCStri(dividend), bigHexCStri(divisor)););
1312       raise_error(NUMERIC_ERROR);
1313       quotient = NULL;
1314     } else {
1315       ALLOC_BIG(quotient);
1316       mpz_init(quotient);
1317       mpz_fdiv_q(quotient, dividend, divisor);
1318     } /* if */
1319     logFunction(printf("bigMDiv --> %s\n", bigHexCStri(quotient)););
1320     return quotient;
1321   } /* bigMDiv */
1322 
1323 
1324 
1325 /**
1326  *  Compute the modulo (remainder) of the integer division bigMDiv.
1327  *  The modulo has the same sign as the divisor.
1328  *  @return the modulo of the integer division.
1329  *  @exception NUMERIC_ERROR If a division by zero occurs.
1330  */
bigMod(const const_bigIntType dividend,const const_bigIntType divisor)1331 bigIntType bigMod (const const_bigIntType dividend, const const_bigIntType divisor)
1332 
1333   {
1334     bigIntType modulo;
1335 
1336   /* bigMod */
1337     logFunction(printf("bigMod(%s,", bigHexCStri(dividend));
1338                 printf("%s)\n", bigHexCStri(divisor)););
1339     if (unlikely(mpz_sgn(divisor) == 0)) {
1340       logError(printf("bigMod(%s, %s): Division by zero.\n",
1341                       bigHexCStri(dividend), bigHexCStri(divisor)););
1342       raise_error(NUMERIC_ERROR);
1343       modulo = NULL;
1344     } else {
1345       ALLOC_BIG(modulo);
1346       mpz_init(modulo);
1347       mpz_fdiv_r(modulo, dividend, divisor);
1348     } /* if */
1349     logFunction(printf("bigMod --> %s\n", bigHexCStri(modulo)););
1350     return modulo;
1351   } /* bigMod */
1352 
1353 
1354 
1355 /**
1356  *  Multiply two 'bigInteger' numbers.
1357  *  @return the product of the two numbers.
1358  */
bigMult(const const_bigIntType factor1,const const_bigIntType factor2)1359 bigIntType bigMult (const const_bigIntType factor1, const const_bigIntType factor2)
1360 
1361   {
1362     bigIntType product;
1363 
1364   /* bigMult */
1365     ALLOC_BIG(product);
1366     mpz_init(product);
1367     mpz_mul(product, factor1, factor2);
1368     return product;
1369   } /* bigMult */
1370 
1371 
1372 
1373 /**
1374  *  Multiply a 'bigInteger' number by a factor and assign the result back to number.
1375  */
bigMultAssign(bigIntType * const big_variable,const_bigIntType factor)1376 void bigMultAssign (bigIntType *const big_variable, const_bigIntType factor)
1377 
1378   { /* bigMultAssign */
1379     mpz_mul(*big_variable, *big_variable, factor);
1380   } /* bigMultAssign */
1381 
1382 
1383 
1384 /**
1385  *  Multiply factor1 with the bigdigit factor2.
1386  *  The range of factor2 is restricted and it is the job of the
1387  *  compiler to assure that factor2 is within the allowed range.
1388  *  @param factor2 Multiplication factor that must be
1389  *         in the range of a long int.
1390  *  @return the product of factor1 * factor2.
1391  */
bigMultSignedDigit(const_bigIntType factor1,intType factor2)1392 bigIntType bigMultSignedDigit (const_bigIntType factor1, intType factor2)
1393 
1394   {
1395     bigIntType product;
1396 
1397   /* bigMultSignedDigit */
1398     ALLOC_BIG(product);
1399     mpz_init(product);
1400     mpz_mul_si(product, factor1, factor2);
1401     return product;
1402   } /* bigMultSignedDigit */
1403 
1404 
1405 
1406 /**
1407  *  Minus sign, negate a 'bigInteger' number.
1408  *  @return the negated value of the number.
1409  */
bigNegate(const const_bigIntType big1)1410 bigIntType bigNegate (const const_bigIntType big1)
1411 
1412   {
1413     bigIntType negatedValue;
1414 
1415   /* bigNegate */
1416     ALLOC_BIG(negatedValue);
1417     mpz_init(negatedValue);
1418     mpz_neg(negatedValue, big1);
1419     return negatedValue;
1420   } /* bigNegate */
1421 
1422 
1423 
1424 /**
1425  *  Minus sign, negate a 'bigInteger' number.
1426  *  Big1 is assumed to be a temporary value which is reused.
1427  *  @return the negated value of the number.
1428  */
bigNegateTemp(bigIntType big1)1429 bigIntType bigNegateTemp (bigIntType big1)
1430 
1431   { /* bigNegateTemp */
1432     mpz_neg(big1, big1);
1433     return big1;
1434   } /* bigNegateTemp */
1435 
1436 
1437 
1438 /**
1439  *  Determine if a 'bigInteger' number is odd.
1440  *  @return TRUE if the number is odd,
1441  *          FALSE otherwise.
1442  */
bigOdd(const const_bigIntType big1)1443 boolType bigOdd (const const_bigIntType big1)
1444 
1445   { /* bigOdd */
1446     return mpz_odd_p(big1);
1447   } /* bigOdd */
1448 
1449 
1450 
bigOr(const_bigIntType big1,const_bigIntType big2)1451 bigIntType bigOr (const_bigIntType big1, const_bigIntType big2)
1452 
1453   {
1454     bigIntType result;
1455 
1456   /* bigOr */
1457     ALLOC_BIG(result);
1458     mpz_init(result);
1459     mpz_ior(result, big1, big2);
1460     return result;
1461   } /* bigOr */
1462 
1463 
1464 
1465 /**
1466  *  Convert a string to a 'bigInteger' number.
1467  *  The string must contain an integer literal consisting of an
1468  *  optional + or - sign, followed by a sequence of digits. Other
1469  *  characters as well as leading or trailing whitespace characters are
1470  *  not allowed. The sequence of digits is taken to be decimal.
1471  *  @return the 'bigInteger' result of the conversion.
1472  *  @exception RANGE_ERROR If the string is empty or does not contain
1473  *             an integer literal.
1474  *  @exception MEMORY_ERROR  Not enough memory to represent the result.
1475  */
bigParse(const const_striType stri)1476 bigIntType bigParse (const const_striType stri)
1477 
1478   {
1479     cstriType cstri;
1480     int mpz_result;
1481     errInfoType err_info = OKAY_NO_ERROR;
1482     bigIntType result;
1483 
1484   /* bigParse */
1485     cstri = stri_to_cstri(stri, &err_info);
1486     if (unlikely(cstri == NULL)) {
1487       logError(printf("bigParse: stri_to_cstri(\"%s\", *) failed:\n"
1488                       "err_info=%d\n",
1489                       striAsUnquotedCStri(stri), err_info););
1490       raise_error(err_info);
1491       result = NULL;
1492     } else if (strpbrk(cstri, " \f\n\r\t\v") != NULL) {
1493       logError(printf("bigParse(\"%s\"): String contains whitespace characters.\n",
1494                       striAsUnquotedCStri(stri)););
1495       raise_error(RANGE_ERROR);
1496       result = NULL;
1497     } else {
1498       ALLOC_BIG(result);
1499       if (cstri[0] == '+' && cstri[1] != '-') {
1500         mpz_result = mpz_init_set_str(result, &cstri[1], 10);
1501       } else {
1502         mpz_result = mpz_init_set_str(result, cstri, 10);
1503       } /* if */
1504       free_cstri(cstri, stri);
1505       if (unlikely(mpz_result != 0)) {
1506         mpz_clear(result);
1507         FREE_BIG(result);
1508         logError(printf("bigParse(\"%s\"): "
1509                         "mpz_init_set_str(*, \"%s\", 10) failed.\n",
1510                         striAsUnquotedCStri(stri),
1511                         cstri[0] == '+' && cstri[1] != '-' ? &cstri[1] : cstri););
1512         raise_error(RANGE_ERROR);
1513         result = NULL;
1514       } /* if */
1515     } /* if */
1516     return result;
1517   } /* bigParse */
1518 
1519 
1520 
1521 /**
1522  *  Convert a numeric string, with a specified radix, to a 'bigInteger'.
1523  *  The numeric string must contain the representation of an integer
1524  *  in the specified radix. It consists of an optional + or - sign,
1525  *  followed by a sequence of digits in the specified radix. Digit values
1526  *  from 10 upward can be encoded with upper or lower case letters.
1527  *  E.g.: 10 can be encoded with A or a, 11 with B or b, etc. Other
1528  *  characters as well as leading or trailing whitespace characters
1529  *  are not allowed.
1530  *  @param base Radix of the integer in the 'stri' parameter.
1531  *  @return the 'bigInteger' result of the conversion.
1532  *  @exception RANGE_ERROR If base < 2 or base > 36 holds or
1533  *             the string is empty or it does not contain an integer
1534  *             literal with the specified base.
1535  *  @exception MEMORY_ERROR  Not enough memory to represent the result.
1536  */
bigParseBased(const const_striType stri,intType base)1537 bigIntType bigParseBased (const const_striType stri, intType base)
1538 
1539   {
1540     cstriType cstri;
1541     int mpz_result;
1542     errInfoType err_info = OKAY_NO_ERROR;
1543     bigIntType result;
1544 
1545   /* bigParseBased */
1546     if (unlikely(stri->size == 0 || base < 2 || base > 36)) {
1547       logError(printf("bigParseBased(\"%s\", " FMT_D "): "
1548                       "String empty or base not in allowed range.\n",
1549                       striAsUnquotedCStri(stri), base););
1550       raise_error(RANGE_ERROR);
1551       result = NULL;
1552     } else {
1553       cstri = stri_to_cstri(stri, &err_info);
1554       if (unlikely(cstri == NULL)) {
1555         logError(printf("bigParseBased: stri_to_cstri(\"%s\", *) failed:\n"
1556                         "err_info=%d\n",
1557                         striAsUnquotedCStri(stri), err_info););
1558         raise_error(err_info);
1559         result = NULL;
1560       } else if (strpbrk(cstri, " \f\n\r\t\v") != NULL) {
1561         logError(printf("bigParseBased(\"%s\", " FMT_D "): "
1562                         "String contains whitespace characters.\n",
1563                         striAsUnquotedCStri(stri), base););
1564         raise_error(RANGE_ERROR);
1565         result = NULL;
1566       } else {
1567         ALLOC_BIG(result);
1568         if (cstri[0] == '+' && cstri[1] != '-') {
1569           mpz_result = mpz_init_set_str(result, &cstri[1], (int) base);
1570         } else {
1571           mpz_result = mpz_init_set_str(result, cstri, (int) base);
1572         } /* if */
1573         free_cstri(cstri, stri);
1574         if (unlikely(mpz_result != 0)) {
1575           mpz_clear(result);
1576           FREE_BIG(result);
1577           logError(printf("bigParseBased(\"%s\", " FMT_D "): "
1578                           "mpz_init_set_str(*, \"%s\", " FMT_D ") failed.\n",
1579                           striAsUnquotedCStri(stri), base,
1580                           cstri[0] == '+' && cstri[1] != '-' ? &cstri[1] : cstri,
1581                           base););
1582           raise_error(RANGE_ERROR);
1583           result = NULL;
1584         } /* if */
1585       } /* if */
1586     } /* if */
1587     return result;
1588   } /* bigParseBased */
1589 
1590 
1591 
1592 /**
1593  *  Predecessor of a 'bigInteger' number.
1594  *  pred(A) is equivalent to A-1 .
1595  *  @return big1 - 1 .
1596  */
bigPred(const const_bigIntType big1)1597 bigIntType bigPred (const const_bigIntType big1)
1598 
1599   {
1600     bigIntType predecessor;
1601 
1602   /* bigPred */
1603     ALLOC_BIG(predecessor);
1604     mpz_init(predecessor);
1605     mpz_sub_ui(predecessor, big1, 1);
1606     return predecessor;
1607   } /* bigPred */
1608 
1609 
1610 
1611 /**
1612  *  Returns a signed big integer decremented by 1.
1613  *  Big1 is assumed to be a temporary value which is reused.
1614  */
bigPredTemp(bigIntType big1)1615 bigIntType bigPredTemp (bigIntType big1)
1616 
1617   { /* bigPredTemp */
1618     mpz_sub_ui(big1, big1, 1);
1619     return big1;
1620   } /* bigPredTemp */
1621 
1622 
1623 
1624 /**
1625  *  Convert a big integer number to a string using a radix.
1626  *  The conversion uses the numeral system with the given base.
1627  *  Digit values from 10 upward are encoded with letters.
1628  *  For negative numbers a minus sign is prepended.
1629  *  @param upperCase Decides about the letter case.
1630  *  @return the string result of the conversion.
1631  *  @exception RANGE_ERROR If base < 2 or base > 36 holds.
1632  *  @exception MEMORY_ERROR Not enough memory to represent the result.
1633  */
bigRadix(const const_bigIntType big1,intType base,boolType upperCase)1634 striType bigRadix (const const_bigIntType big1, intType base,
1635     boolType upperCase)
1636 
1637   {
1638     char *cstri;
1639     striType result;
1640 
1641   /* bigRadix */
1642     logFunction(printf("bigRadix(%s, " FMT_D ", %d)\n",
1643                        bigHexCStri(big1), base, upperCase););
1644     if (unlikely(base < 2 || base > 36)) {
1645       logError(printf("bigRadix(%s, " FMT_D ", %d): "
1646                       "Base not in allowed range.\n",
1647                       bigHexCStri(big1), base, upperCase););
1648       raise_error(RANGE_ERROR);
1649       result = NULL;
1650     } else {
1651       if (upperCase) {
1652         cstri = mpz_get_str(NULL, -(int) base, big1);
1653       } else {
1654         cstri = mpz_get_str(NULL, (int) base, big1);
1655       } /* if */
1656       result = cstri_to_stri(cstri);
1657       free(cstri);
1658       if (unlikely(result == NULL)) {
1659         raise_error(MEMORY_ERROR);
1660       } /* if */
1661     } /* if */
1662     logFunction(printf("bigRadix --> \"%s\"\n", striAsUnquotedCStri(result)););
1663     return result;
1664   } /* bigRadix */
1665 
1666 
1667 
1668 /**
1669  *  Compute pseudo-random number in the range [low, high].
1670  *  The random values are uniform distributed.
1671  *  @return a random number such that low <= rand(low, high) and
1672  *          rand(low, high) <= high holds.
1673  *  @exception RANGE_ERROR The range is empty (low > high holds).
1674  */
bigRand(const const_bigIntType low,const const_bigIntType high)1675 bigIntType bigRand (const const_bigIntType low,
1676     const const_bigIntType high)
1677 
1678   {
1679     static boolType seed_necessary = TRUE;
1680     static gmp_randstate_t state;
1681     mpz_t range_limit;
1682     bigIntType randomNumber;
1683 
1684   /* bigRand */
1685     mpz_init(range_limit);
1686     mpz_sub(range_limit, high, low);
1687     if (unlikely(mpz_sgn(range_limit) < 0)) {
1688       logError(printf("bigRand(%s, %s): "
1689                       "The range is empty (low > high holds).\n",
1690                       bigHexCStri(low), bigHexCStri(high)););
1691       raise_error(RANGE_ERROR);
1692       randomNumber = 0;
1693     } else {
1694       mpz_add_ui(range_limit, range_limit, 1);
1695       ALLOC_BIG(randomNumber);
1696       mpz_init(randomNumber);
1697       if (seed_necessary) {
1698         gmp_randinit_default(state);
1699         seed_necessary = FALSE;
1700       } /* if */
1701       mpz_urandomm(randomNumber, state, range_limit);
1702       mpz_add(randomNumber, randomNumber, low);
1703     } /* if */
1704     mpz_clear(range_limit);
1705     return randomNumber;
1706   } /* bigRand */
1707 
1708 
1709 
1710 /**
1711  *  Compute the remainder of the integer division bigDiv.
1712  *  The remainder has the same sign as the dividend.
1713  *  @return the remainder of the integer division.
1714  *  @exception NUMERIC_ERROR If a division by zero occurs.
1715  */
bigRem(const const_bigIntType dividend,const const_bigIntType divisor)1716 bigIntType bigRem (const const_bigIntType dividend, const const_bigIntType divisor)
1717 
1718   {
1719     bigIntType remainder;
1720 
1721   /* bigRem */
1722     logFunction(printf("bigRem(%s,", bigHexCStri(dividend));
1723                 printf("%s)\n", bigHexCStri(divisor)););
1724     if (unlikely(mpz_sgn(divisor) == 0)) {
1725       logError(printf("bigRem(%s, %s): Division by zero.\n",
1726                       bigHexCStri(dividend), bigHexCStri(divisor)););
1727       raise_error(NUMERIC_ERROR);
1728       remainder = NULL;
1729     } else {
1730       ALLOC_BIG(remainder);
1731       mpz_init(remainder);
1732       mpz_tdiv_r(remainder, dividend, divisor);
1733     } /* if */
1734     logFunction(printf("bigRem --> %s\n", bigHexCStri(remainder)););
1735     return remainder;
1736   } /* bigRem */
1737 
1738 
1739 
1740 /**
1741  *  Shift a 'bigInteger' number right by rshift bits.
1742  *  If rshift is negative a left shift is done instead.
1743  *  A >> B is equivalent to A mdiv 2_ ** B if B >= 0 holds.
1744  *  A >> B is equivalent to A * 2_ ** -B if B < 0 holds.
1745  *  @return the right shifted number.
1746  */
bigRShift(const const_bigIntType big1,const intType rshift)1747 bigIntType bigRShift (const const_bigIntType big1, const intType rshift)
1748 
1749   {
1750     bigIntType result;
1751 
1752   /* bigRShift */
1753     logFunction(printf("bigRShift(%s, " FMT_D ")\n",
1754                        bigHexCStri(big1), rshift););
1755     ALLOC_BIG(result);
1756     mpz_init(result);
1757     if (rshift >= 0) {
1758       mpz_fdiv_q_2exp(result, big1, (uintType) rshift);
1759 #if AVOID_OVERFLOW_IN_MPZ_TYPE
1760     } else if (unlikely((-(uintType) rshift) / GMP_LIMB_BITS + 1 >
1761                         GMP_MAX_LIMB_SIZE - mpz_size(big1) &&
1762                         mpz_sgn(big1) != 0)) {
1763       mpz_clear(result);
1764       FREE_BIG(result);
1765       raise_error(OVERFLOW_IN_MPZ_ERROR);
1766       result = NULL;
1767 #endif
1768     } else {
1769       /* The unsigned value is negated to avoid a signed integer */
1770       /* overflow if the smallest signed integer is negated.     */
1771       mpz_mul_2exp(result, big1, -(uintType) rshift);
1772     } /* if */
1773     logFunction(printf("bigRShift --> %s\n", bigHexCStri(result)););
1774     return result;
1775   } /* bigRShift */
1776 
1777 
1778 
1779 /**
1780  *  Shift a number right by rshift bits and assign the result back to number.
1781  *  If rshift is negative a left shift is done instead.
1782  */
bigRShiftAssign(bigIntType * const big_variable,intType rshift)1783 void bigRShiftAssign (bigIntType *const big_variable, intType rshift)
1784 
1785   {
1786     bigIntType big1;
1787 
1788   /* bigRShiftAssign */
1789     big1 = *big_variable;
1790     logFunction(printf("bigRShiftAssign(%s, " FMT_D ")\n",
1791                        bigHexCStri(big1), rshift););
1792     if (rshift > 0) {
1793       mpz_fdiv_q_2exp(big1, big1, (uintType) rshift);
1794     } else if (rshift == 0) {
1795       /* do nothing */
1796 #if AVOID_OVERFLOW_IN_MPZ_TYPE
1797     } else if (unlikely((-(uintType) rshift) / GMP_LIMB_BITS + 1 >
1798                         GMP_MAX_LIMB_SIZE - mpz_size(big1) &&
1799                         mpz_sgn(big1) != 0)) {
1800       raise_error(OVERFLOW_IN_MPZ_ERROR);
1801 #endif
1802     } else {
1803       /* The unsigned value is negated to avoid a signed integer */
1804       /* overflow if the smallest signed integer is negated.     */
1805       mpz_mul_2exp(big1, big1, -(uintType) rshift);
1806     } /* if */
1807     logFunction(printf("bigRShiftAssign --> %s\n", bigHexCStri(big1)););
1808   } /* bigRShiftAssign */
1809 
1810 
1811 
1812 /**
1813  *  Compute the subtraction of two 'bigInteger' numbers.
1814  *  @return the difference of the two numbers.
1815  */
bigSbtr(const const_bigIntType minuend,const const_bigIntType subtrahend)1816 bigIntType bigSbtr (const const_bigIntType minuend, const const_bigIntType subtrahend)
1817 
1818   {
1819     bigIntType difference;
1820 
1821   /* bigSbtr */
1822     ALLOC_BIG(difference);
1823     mpz_init(difference);
1824     mpz_sub(difference, minuend, subtrahend);
1825     return difference;
1826   } /* bigSbtr */
1827 
1828 
1829 
1830 /**
1831  *  Decrement a 'bigInteger' variable by a delta.
1832  *  Subtracts delta from *big_variable.
1833  */
bigSbtrAssign(bigIntType * const big_variable,const const_bigIntType delta)1834 void bigSbtrAssign (bigIntType *const big_variable, const const_bigIntType delta)
1835 
1836   { /* bigSbtrAssign */
1837     mpz_sub(*big_variable, *big_variable, delta);
1838   } /* bigSbtrAssign */
1839 
1840 
1841 
1842 /**
1843  *  Compute the subtraction of two 'bigInteger' numbers.
1844  *  Minuend is assumed to be a temporary value which is reused.
1845  *  @return the difference of the two numbers in 'minuend'.
1846  */
bigSbtrTemp(bigIntType minuend,const_bigIntType subtrahend)1847 bigIntType bigSbtrTemp (bigIntType minuend, const_bigIntType subtrahend)
1848 
1849   { /* bigSbtrTemp */
1850     mpz_sub(minuend, minuend, subtrahend);
1851     return minuend;
1852   } /* bigSbtrTemp */
1853 
1854 
1855 
1856 /**
1857  *  Compute the square of a 'bigInteger'.
1858  *  This function is used by the compiler to optimize
1859  *  multiplication and exponentiation operations.
1860  *  @return the square of big1.
1861  */
bigSquare(const_bigIntType big1)1862 bigIntType bigSquare (const_bigIntType big1)
1863 
1864   {
1865     bigIntType result;
1866 
1867   /* bigSquare */
1868     ALLOC_BIG(result);
1869     mpz_init(result);
1870     mpz_mul(result, big1, big1);
1871     return result;
1872   } /* bigSquare */
1873 
1874 
1875 
1876 /**
1877  *  Convert a 'bigInteger' number to a string.
1878  *  The number is converted to a string with decimal representation.
1879  *  For negative numbers a minus sign is prepended.
1880  *  @return the string result of the conversion.
1881  *  @exception MEMORY_ERROR  Not enough memory to represent the result.
1882  */
bigStr(const const_bigIntType big1)1883 striType bigStr (const const_bigIntType big1)
1884 
1885   {
1886     char *cstri;
1887     striType result;
1888 
1889   /* bigStr */
1890     logFunction(printf("bigStr(%s)\n", bigHexCStri(big1)););
1891     cstri = mpz_get_str(NULL, 10, big1);
1892     result = cstri_to_stri(cstri);
1893     free(cstri);
1894     if (unlikely(result == NULL)) {
1895       raise_error(MEMORY_ERROR);
1896     } /* if */
1897     logFunction(printf("bigStr --> \"%s\"\n", striAsUnquotedCStri(result)););
1898     return result;
1899   } /* bigStr */
1900 
1901 
1902 
1903 /**
1904  *  Successor of a 'bigInteger' number.
1905  *  succ(A) is equivalent to A+1 .
1906  *  @return big1 + 1 .
1907  */
bigSucc(const const_bigIntType big1)1908 bigIntType bigSucc (const const_bigIntType big1)
1909 
1910   {
1911     bigIntType successor;
1912 
1913   /* bigSucc */
1914     ALLOC_BIG(successor);
1915     mpz_init(successor);
1916     mpz_add_ui(successor, big1, 1);
1917     return successor;
1918   } /* bigSucc */
1919 
1920 
1921 
1922 /**
1923  *  Successor of a 'bigInteger' number.
1924  *  Big1 is assumed to be a temporary value which is reused.
1925  *  @return big1 + 1 .
1926  */
bigSuccTemp(bigIntType big1)1927 bigIntType bigSuccTemp (bigIntType big1)
1928 
1929   { /* bigSuccTemp */
1930     mpz_add_ui(big1, big1, 1);
1931     return big1;
1932   } /* bigSuccTemp */
1933 
1934 
1935 
1936 /**
1937  *  Convert a 'bigInteger' into a big-endian 'bstring'.
1938  *  The result uses binary representation with a base of 256.
1939  *  @param big1 BigInteger number to be converted.
1940  *  @param isSigned Determines the signedness of the result.
1941  *         If 'isSigned' is TRUE the result is encoded with the
1942  *         twos-complement representation. In this case a negative
1943  *         'big1' is converted to a result where the most significant
1944  *         byte (the first byte) has an ordinal > BYTE_MAX (=127).
1945  *  @return a bstring with the big-endian representation.
1946  *  @exception RANGE_ERROR If 'big1' is negative and 'isSigned' is FALSE.
1947  *  @exception MEMORY_ERROR Not enough memory to represent the result.
1948  */
bigToBStriBe(const const_bigIntType big1,const boolType isSigned)1949 bstriType bigToBStriBe (const const_bigIntType big1, const boolType isSigned)
1950 
1951   {
1952     size_t sizeInBytes;
1953     size_t bytesExported;
1954     size_t pos;
1955     int sign;
1956     unsigned int carry;
1957     ustriType buffer;
1958     memSizeType charIndex;
1959     memSizeType result_size;
1960     bstriType result;
1961 
1962   /* bigToBStriBe */
1963     sizeInBytes = (mpz_sizeinbase(big1, 2) + CHAR_BIT - 1) / CHAR_BIT;
1964     buffer = (ustriType) malloc(sizeInBytes);
1965     if (unlikely(buffer == NULL)) {
1966       raise_error(MEMORY_ERROR);
1967       result = NULL;
1968     } else {
1969       mpz_export(buffer, &bytesExported, 1, 1, 0, 0, big1);
1970       sign = mpz_sgn(big1);
1971       if (isSigned) {
1972         if (sign < 0) {
1973           carry = 1;
1974           pos = bytesExported;
1975           while (pos > 0) {
1976             pos--;
1977             carry += ~buffer[pos] & 0xFF;
1978             buffer[pos] = (ucharType) (carry & 0xFF);
1979             carry >>= CHAR_BIT;
1980           } /* while */
1981         } /* if */
1982         result_size = sizeInBytes;
1983         if ((sign > 0 && buffer[0] > BYTE_MAX) ||
1984             (sign < 0 && buffer[0] <= BYTE_MAX)) {
1985           result_size++;
1986         } /* if */
1987       } else {
1988         if (unlikely(sign < 0)) {
1989           logError(printf("bigToBStriBe(%s, %d): "
1990                           "Number is negative and 'isSigned' is FALSE.\n",
1991                           bigHexCStri(big1), isSigned););
1992           raise_error(RANGE_ERROR);
1993           free(buffer);
1994           return NULL;
1995         } else {
1996           result_size = sizeInBytes;
1997         } /* if */
1998       } /* if */
1999       if (unlikely(!ALLOC_BSTRI_CHECK_SIZE(result, result_size))) {
2000         raise_error(MEMORY_ERROR);
2001       } else {
2002         result->size = result_size;
2003         if (sign == 0) {
2004           result->mem[0] = 0;
2005         } else {
2006           charIndex = 0;
2007           if (isSigned) {
2008             if (sign < 0) {
2009               if (buffer[0] <= BYTE_MAX) {
2010                 result->mem[charIndex] = UBYTE_MAX;
2011                 charIndex++;
2012               } /* if */
2013             } else {
2014               if (buffer[0] > BYTE_MAX) {
2015                 result->mem[charIndex] = 0;
2016                 charIndex++;
2017               } /* if */
2018             } /* if */
2019           } /* if */
2020           memcpy(&result->mem[charIndex], buffer, bytesExported);
2021         } /* if */
2022       } /* if */
2023       free(buffer);
2024     } /* if */
2025     return result;
2026   } /* bigToBStriBe */
2027 
2028 
2029 
2030 /**
2031  *  Convert a 'bigInteger' into a little-endian 'bstring'.
2032  *  The result uses binary representation with a base of 256.
2033  *  @param big1 BigInteger number to be converted.
2034  *  @param isSigned Determines the signedness of the result.
2035  *         If 'isSigned' is TRUE the result is encoded with the
2036  *         twos-complement representation. In this case a negative
2037  *         'big1' is converted to a result where the most significant
2038  *         byte (the last byte) has an ordinal > BYTE_MAX (=127).
2039  *  @return a bstring with the little-endian representation.
2040  *  @exception RANGE_ERROR If 'big1' is negative and 'isSigned' is FALSE.
2041  *  @exception MEMORY_ERROR Not enough memory to represent the result.
2042  */
bigToBStriLe(const const_bigIntType big1,const boolType isSigned)2043 bstriType bigToBStriLe (const const_bigIntType big1, const boolType isSigned)
2044 
2045   {
2046     size_t sizeInBytes;
2047     size_t bytesExported;
2048     size_t pos;
2049     int sign;
2050     unsigned int carry;
2051     ustriType buffer;
2052     memSizeType result_size;
2053     bstriType result;
2054 
2055   /* bigToBStriLe */
2056     sizeInBytes = (mpz_sizeinbase(big1, 2) + CHAR_BIT - 1) / CHAR_BIT;
2057     buffer = (ustriType) malloc(sizeInBytes);
2058     if (unlikely(buffer == NULL)) {
2059       raise_error(MEMORY_ERROR);
2060       result = NULL;
2061     } else {
2062       mpz_export(buffer, &bytesExported, -1, 1, 0, 0, big1);
2063       sign = mpz_sgn(big1);
2064       if (isSigned) {
2065         if (sign < 0) {
2066           carry = 1;
2067           pos = 0;
2068           while (pos < bytesExported) {
2069             carry += ~buffer[pos] & 0xFF;
2070             buffer[pos] = (ucharType) (carry & 0xFF);
2071             carry >>= CHAR_BIT;
2072             pos++;
2073           } /* while */
2074         } /* if */
2075         result_size = sizeInBytes;
2076         if ((sign > 0 && buffer[bytesExported - 1] > BYTE_MAX) ||
2077             (sign < 0 && buffer[bytesExported - 1] <= BYTE_MAX)) {
2078           result_size++;
2079         } /* if */
2080       } else {
2081         if (unlikely(sign < 0)) {
2082           logError(printf("bigToBStriLe(%s, %d): "
2083                           "Number is negative and 'isSigned' is FALSE.\n",
2084                           bigHexCStri(big1), isSigned););
2085           raise_error(RANGE_ERROR);
2086           free(buffer);
2087           return NULL;
2088         } else {
2089           result_size = sizeInBytes;
2090         } /* if */
2091       } /* if */
2092       if (unlikely(!ALLOC_BSTRI_CHECK_SIZE(result, result_size))) {
2093         raise_error(MEMORY_ERROR);
2094       } else {
2095         result->size = result_size;
2096         if (sign == 0) {
2097           result->mem[0] = 0;
2098         } else {
2099           memcpy(result->mem, buffer, bytesExported);
2100           if (isSigned) {
2101             if (sign < 0) {
2102               if (buffer[bytesExported - 1] <= BYTE_MAX) {
2103                 result->mem[bytesExported] = UBYTE_MAX;
2104               } /* if */
2105             } else {
2106               if (buffer[bytesExported - 1] > BYTE_MAX) {
2107                 result->mem[bytesExported] = 0;
2108               } /* if */
2109             } /* if */
2110           } /* if */
2111         } /* if */
2112       } /* if */
2113       free(buffer);
2114     } /* if */
2115     return result;
2116   } /* bigToBStriLe */
2117 
2118 
2119 
2120 /**
2121  *  Convert a 'bigInteger' to an 'int16Type' number.
2122  *  @return the int16Type result of the conversion.
2123  *  @param big1 BigInteger to be converted.
2124  *  @param err_info Unchanged if the function succeeds or
2125  *                  RANGE_ERROR The number is too small or too big
2126  *                  to fit into a int16Type value.
2127  */
bigToInt16(const const_bigIntType big1,errInfoType * err_info)2128 int16Type bigToInt16 (const const_bigIntType big1, errInfoType *err_info)
2129 
2130   {
2131     long int result;
2132 
2133   /* bigToInt16 */
2134     if (unlikely(!mpz_fits_slong_p(big1))) {
2135       logError(printf("bigToInt16(%s): mpz_fits_slong_p failed.\n",
2136                       bigHexCStri(big1)););
2137       *err_info = RANGE_ERROR;
2138       return 0;
2139     } else {
2140       result = mpz_get_si(big1);
2141 #if LONG_SIZE > 16
2142       if (unlikely(result < INT16TYPE_MIN || result > INT16TYPE_MAX)) {
2143         *err_info = RANGE_ERROR;
2144         return 0;
2145       } /* if */
2146 #endif
2147       return (int16Type) result;
2148     } /* if */
2149   } /* bigToInt16 */
2150 
2151 
2152 
2153 /**
2154  *  Convert a 'bigInteger' to an 'int32Type' number.
2155  *  @return the int32Type result of the conversion.
2156  *  @param big1 BigInteger to be converted.
2157  *  @param err_info Only used if err_info is not NULL.
2158  *                  Unchanged if the function succeeds or
2159  *                  RANGE_ERROR The number is too small or too big
2160  *                  to fit into a int32Type value.
2161  *  @exception RANGE_ERROR If err_info is NULL and the number is
2162  *             too small or too big to fit into a int32Type value.
2163  */
bigToInt32(const const_bigIntType big1,errInfoType * err_info)2164 int32Type bigToInt32 (const const_bigIntType big1, errInfoType *err_info)
2165 
2166   {
2167     long int result;
2168 
2169   /* bigToInt32 */
2170     logFunction(printf("bigToInt32(%s)\n", bigHexCStri(big1)););
2171     if (unlikely(!mpz_fits_slong_p(big1))) {
2172       logError(printf("bigToInt32(%s): mpz_fits_slong_p failed.\n",
2173                       bigHexCStri(big1)););
2174       if (err_info == NULL) {
2175         raise_error(RANGE_ERROR);
2176       } else {
2177         *err_info = RANGE_ERROR;
2178       } /* if */
2179       result = 0;
2180     } else {
2181       result = mpz_get_si(big1);
2182 #if LONG_SIZE > 32
2183       if (unlikely(result < INT32TYPE_MIN || result > INT32TYPE_MAX)) {
2184         if (err_info == NULL) {
2185           raise_error(RANGE_ERROR);
2186         } else {
2187           *err_info = RANGE_ERROR;
2188         } /* if */
2189         result = 0;
2190       } /* if */
2191 #endif
2192     } /* if */
2193     logFunction(printf("bigToInt32 --> " FMT_D32 "\n", (int32Type) result););
2194     return (int32Type) result;
2195   } /* bigToInt32 */
2196 
2197 
2198 
2199 #ifdef INT64TYPE
2200 /**
2201  *  Convert a 'bigInteger' to an 'int64Type' number.
2202  *  @return the int64Type result of the conversion.
2203  *  @param big1 BigInteger to be converted.
2204  *  @param err_info Only used if err_info is not NULL.
2205  *                  Unchanged if the function succeeds or
2206  *                  RANGE_ERROR The number is too small or too big
2207  *                  to fit into a int64Type value.
2208  *  @exception RANGE_ERROR If err_info is NULL and the number is
2209  *             too small or too big to fit into a int64Type value.
2210  */
bigToInt64(const const_bigIntType big1,errInfoType * err_info)2211 int64Type bigToInt64 (const const_bigIntType big1, errInfoType *err_info)
2212 
2213   {
2214     int64Type result;
2215 
2216   /* bigToInt64 */
2217     logFunction(printf("bigToInt64(%s)\n", bigHexCStri(big1)););
2218 #if LONG_SIZE == 64
2219     if (unlikely(!mpz_fits_slong_p(big1))) {
2220       logError(printf("bigToInt64(%s): mpz_fits_slong_p failed.\n",
2221                       bigHexCStri(big1)););
2222       if (err_info == NULL) {
2223         raise_error(RANGE_ERROR);
2224       } else {
2225         *err_info = RANGE_ERROR;
2226       } /* if */
2227       result = 0;
2228     } else {
2229       result = mpz_get_si(big1);
2230     } /* if */
2231 #else
2232     {
2233       mpz_t help;
2234 
2235       mpz_init_set(help, big1);
2236       mpz_fdiv_q_2exp(help, help, 32);
2237       if (unlikely(!mpz_fits_slong_p(help))) {
2238         logError(printf("bigToInt64(%s): mpz_fits_slong_p failed.\n",
2239                         bigHexCStri(big1)););
2240         if (err_info == NULL) {
2241           raise_error(RANGE_ERROR);
2242         } else {
2243           *err_info = RANGE_ERROR;
2244         } /* if */
2245         result = 0;
2246       } else {
2247         result = (int64Type) (mpz_get_si(help) & 0xFFFFFFFF) << 32 |
2248                  (int64Type) (mpz_get_si(big1) & 0xFFFFFFFF);
2249         mpz_clear(help);
2250       } /* if */
2251     }
2252 #endif
2253     logFunction(printf("bigToInt64 --> " FMT_D64 "\n", result););
2254     return result;
2255   } /* bigToInt64 */
2256 
2257 
2258 
2259 /**
2260  *  Convert a 'bigInteger' to an 'uint64Type' number.
2261  *  @return the uint64Type result of the conversion.
2262  *  @exception RANGE_ERROR The number is negative or too big to fit
2263  *             into a uint64Type value.
2264  */
bigToUInt64(const const_bigIntType big1)2265 uint64Type bigToUInt64 (const const_bigIntType big1)
2266 
2267   {
2268     uint64Type result;
2269 
2270   /* bigToUInt64 */
2271     logFunction(printf("bigToUInt64(%s)\n", bigHexCStri(big1)););
2272 #if LONG_SIZE == 64
2273     if (unlikely(!mpz_fits_ulong_p(big1))) {
2274       logError(printf("bigToUInt64(%s): mpz_fits_ulong_p failed.\n",
2275                       bigHexCStri(big1)););
2276       raise_error(RANGE_ERROR);
2277       result = 0;
2278     } else {
2279       result = mpz_get_ui(big1);
2280     } /* if */
2281 #else
2282     {
2283       mpz_t help;
2284 
2285       mpz_init_set(help, big1);
2286       mpz_fdiv_q_2exp(help, help, 32);
2287       if (unlikely(!mpz_fits_ulong_p(help))) {
2288         logError(printf("bigToUInt64(%s): mpz_fits_ulong_p failed.\n",
2289                         bigHexCStri(big1)););
2290         raise_error(RANGE_ERROR);
2291         result = 0;
2292       } else {
2293         result = (uint64Type) (mpz_get_ui(help) & 0xFFFFFFFF) << 32 |
2294                  (uint64Type) (mpz_get_ui(big1) & 0xFFFFFFFF);
2295         mpz_clear(help);
2296       } /* if */
2297     }
2298 #endif
2299     logFunction(printf("bigToUInt64(%s) --> " FMT_U64 "\n",
2300                        bigHexCStri(big1), result););
2301     return result;
2302   } /* bigToUInt64 */
2303 #endif
2304 
2305 
2306 
bigXor(const_bigIntType big1,const_bigIntType big2)2307 bigIntType bigXor (const_bigIntType big1, const_bigIntType big2)
2308 
2309   {
2310     bigIntType result;
2311 
2312   /* bigXor */
2313     ALLOC_BIG(result);
2314     mpz_init(result);
2315     mpz_xor(result, big1, big2);
2316     return result;
2317   } /* bigXor */
2318 
2319 
2320 
bigZero(void)2321 bigIntType bigZero (void)
2322 
2323   {
2324     bigIntType result;
2325 
2326   /* bigZero */
2327     ALLOC_BIG(result);
2328     mpz_init_set_ui(result, 0);
2329     return result;
2330   } /* bigZero */
2331 
2332 #endif
2333