source: trunk/src/bignum.c @ 544

Last change on this file since 544 was 544, checked in by katerina, 3 years ago

Fix for ticket #436 (new gcc compiler options, including LTO).

File size: 35.5 KB
Line 
1/* Implementation of bignums by Henrik.Johansson@Nexus.Comm.SE in 1991
2 * version 1.2
3 */
4
5/* Modifications in this file (Nov 1999, Rainer Wichmann):
6 * - a few compiler warnings fixed:
7 *   gcc insists on 'if ( ( b_eq_qr = ((b == q) || (b == r)) ) )'
8 *   instead of     'if (   b_eq_qr = ((b == q) || (b == r))   )'
9 * - code that is not used in samhain is #defined out by:
10 *    #if (defined (SH_WITH_CLIENT) || defined (SH_WITH_SERVER))
11 *   or
12 *    #ifdef UNUSED_CODE
13 * - the error message (error_string[]) is enclosed in N_(),_() macros
14 * - the following four #includes have been added
15 */
16#include "config_xor.h"
17
18
19#include <string.h>
20
21#ifdef HAVE_MEMORY_H
22#include <memory.h>
23#endif
24
25#include "samhain.h"
26#include "bignum.h"
27
28
29#if (!defined(HAVE_LIBGMP) || !defined(HAVE_GMP_H)) && (defined (SH_WITH_CLIENT) || defined (SH_WITH_SERVER))
30
31#define DIGIT BIGNUM_DIGIT
32#define DIGIT2 BIGNUM_TWO_DIGITS
33
34#ifndef TRUE
35#define TRUE (1 == 1)
36#endif
37#ifndef FALSE
38#define FALSE (1 != 1)
39#endif
40
41
42#define MIN_ALLOC ((sizeof(long) / sizeof(DIGIT)) << 1)
43
44/* Don't expect the "sign" field to have the right value (BIG_SIGN_0) at
45 * all times, so compare real content of bignum with zerop.
46 * Maybe it would be better to be sure at all times that the sign is right?
47 */
48#define zerop(BIG) ((*(BIG)->dp == 0) && ((BIG)->dgs_used == 1))
49#define uonep(BIG) ((*(BIG)->dp == 1) && ((BIG)->dgs_used == 1))
50#define NEGATE_SIGN(SGN) (-(SGN))
51
52#define DIGIT_BITS BIGNUM_DIGIT_BITS
53#define DIGIT2_BITS BIGNUM_TWO_DIGITS_BITS
54
55#define DIGIT_PART(N) ((DIGIT)((N) & ((((DIGIT2)1) << DIGIT_BITS) - 1)))
56#define RESULT_MINUSP(RES) (((RES) & ((DIGIT2)1 << (DIGIT2_BITS - 1))) != 0)
57#define SHIFT_DIGIT_DOWN(N) (DIGIT_PART((N) >> DIGIT_BITS))
58#define SHIFT_DIGIT_UP(N) ((DIGIT2)((N) << DIGIT_BITS))
59
60#define DIGIT_BASE (((DIGIT2)1) << DIGIT_BITS)
61#define MAX_DIGIT ((((DIGIT2)1) << DIGIT_BITS) - 1)
62
63#define uint unsigned int
64#define ulong unsigned long
65
66char *last_string    = NULL;
67char *big_end_string = NULL;
68
69ulong length_last_string;
70static char big_base_digits1[] =
71N_("00112233445566778899AaBbCcDdEeFfGgHhIiJjKkLlMmNnOoPpQqRrSsTtUuVvWwXxYyZz");
72static char error_string[] = N_("(big_errno != 0, something has gone wrong)");
73static char big_base_digits[73] = "\0";
74
75bignum big_one;
76bignum tmp_string, tmp_add, tmp_mul, tmp_q, tmp_r, tmp_rand, tmp_round;
77DIGIT *tmp_digits;
78ulong  tmp_ulong;
79
80/* #define COUNT_ALLOC */
81#ifdef COUNT_ALLOC
82ulong dgs_alloc = 0, digits_freed = 0;
83#endif
84
85int big_errno = BIG_OK;
86
87double log2tbl[] =
88{
89    -1.0, -1.0,
90    1.000000000000000, 1.584961891174317, 2.000000000000002, 2.321928024291994,
91    2.584960937500003, 2.807353973388675, 3.000000000000004, 3.169923782348637,
92    3.321928024291997, 3.459430694580083, 3.584960937500005, 3.700439453125006,
93    3.807353973388678, 3.906888961791999, 4.000000000000014, 4.087459564208999,
94    4.169921875000016, 4.247924804687517, 4.321926116943377, 4.392314910888691,
95    4.459430694580098, 4.523559570312520, 4.584960937500021, 4.643856048584007,
96    4.700439453125023, 4.754886627197290, 4.807353973388697, 4.857978820800807,
97    4.906887054443386, 4.954193115234403, 5.000000000000028, 5.044391632080107,
98    5.087459564209015, 5.129280090332062, 5.169921875000032
99};
100
101struct digit_blocks
102{
103    int digCnt;
104    DIGIT dig_base;
105} big_block[37];                /* The code uses index 2 to 36 */
106
107
108/* ----------------------------------------------------------------------
109 * All static utility functions are placed first in this file.
110 */
111
112#ifndef HAVE_STRCHR
113static char * strchr (char * str_ptr, char ch)     
114{
115    do
116    {
117        if (*str_ptr == ch)
118        {
119            return str_ptr;
120        }
121    } while (*str_ptr++ != '\0');
122    return NULL;
123}
124#endif
125
126#if 0
127
128extern int printf();
129
130static void
131print_digits(prefix, a)
132char *prefix;
133bignum *a;
134{
135    unsigned long i;
136    ulong tmp;
137
138    printf("%s", prefix);
139    if (a->sign == BIG_SIGN_MINUS)
140    {
141        printf("- ");
142    }
143    else
144    {
145        if (a->sign == BIG_SIGN_PLUS)
146        {
147            printf("+ ");
148        }
149        else
150        {
151            printf("+/- ");
152        }
153    }
154    for (i = a->dgs_used - 1; i > 0; i--)
155    {
156        tmp = a->dp[i];
157        printf("%lu, ", tmp);
158    }
159    tmp = *a->dp;
160    printf("%lu\n", tmp);
161}
162
163#else
164
165#define print_digits(prefix, big) /* */
166
167#endif
168
169static void
170init_digit_blocks(void)
171{
172    uint digcnt, base;
173    DIGIT maxdigit, curdigit, tmp;
174
175    for (base = 2; base <= 36; base++)
176    {
177        tmp = ((1 << (DIGIT_BITS - 1)) / base);
178        maxdigit = tmp * 2;
179        curdigit = 1;
180        digcnt = 0;
181        while (curdigit < maxdigit)
182        {
183            digcnt++;
184            curdigit *= base;
185        }
186        big_block[base].digCnt = digcnt;
187        big_block[base].dig_base = curdigit;
188    }
189}
190
191#ifdef COUNT_ALLOC
192static void
193free_digits(DIGIT * dp, ulong count)
194{
195  /* Mar 4 2000 R. Wichmann: check for dp == NULL */
196  if (!dp) 
197    return;
198  digits_freed += count;
199  free((char *)dp);
200  dp = NULL;
201}
202#else
203/* June 5 2000 R. Wichmann: check for dp == NULL */
204/* #define free_digits(DP,CNT) free((char *)DP) */
205static void
206free_digits(DIGIT * dp, ulong  count)
207{
208  /* Mar 4 2000 R. Wichmann: check for dp == NULL */
209  if (!dp) 
210    return;
211  (void) count;
212  free((char *)dp);
213  dp = NULL;
214}
215#endif
216
217static DIGIT *
218allocate_digits(ulong alloclen)
219{
220    DIGIT *digit_ptr;
221
222    if (big_errno != BIG_OK)
223    {
224        return NULL;
225    }
226#ifdef COUNT_ALLOC
227    dgs_alloc += alloclen;
228#endif
229    digit_ptr = calloc(alloclen, sizeof(DIGIT));
230    if (digit_ptr == NULL)
231    {
232        big_errno = BIG_MEMERR;
233        return NULL;
234    }
235    return digit_ptr;
236}
237
238static bigerr_t
239newsize(DIGIT **dp_p, ulong *cursize_p, ulong minsz, ulong newsz)
240{
241    if (*cursize_p >= minsz)
242    {
243        return big_errno;
244    }
245    free_digits(*dp_p, *cursize_p);
246    if (newsz < MIN_ALLOC)
247    {
248        newsz = MIN_ALLOC;
249    }
250    if ((*dp_p = allocate_digits(newsz)) == NULL)
251    {
252        return big_errno;
253    }
254    *cursize_p = newsz;
255    return big_errno;
256}
257
258/* `memcpy' uses an `int' as counter.  If an int can't hold a big enough
259 * counter, then `digits_cpy' becomes a function (and a little slower,
260 * probably.  Unfortunately `memcpy' takes a signed counter, but _that_
261 * size of numbers are almost too big!  Or isn't it?
262 */
263#ifdef MEMCPY_LONG_COUNTER
264/* extern char *memcpy(); */
265#define digits_cpy(dst, src, count) memcpy((char *)(dst), (char *)(src), \
266                                           (count) * sizeof(DIGIT))
267#else
268static void
269digits_cpy(DIGIT *dst, DIGIT *src, ulong count)
270{
271    while (count-- > 0)
272    {
273        *dst++ = *src++;
274    }
275}
276#endif
277
278static DIGIT *
279copy_digits(DIGIT *src, ulong cp_count, ulong alloclen)
280{
281    DIGIT *dst;
282
283    if (big_errno != BIG_OK)
284    {
285        return NULL;
286    }
287    if ((dst = allocate_digits(alloclen)) == NULL)
288    {
289        return NULL;
290    }
291    digits_cpy(dst, src, cp_count);
292    return dst;
293}
294
295static void
296add_digit(bignum * a, DIGIT d)
297{
298    DIGIT2 res = d;
299    DIGIT *last_digit, *digit_ptr = a->dp, *digvect;
300   
301    last_digit = digit_ptr + a->dgs_used;
302    while (digit_ptr < last_digit)
303    {
304        res = *digit_ptr + res;
305        *digit_ptr = DIGIT_PART(res);
306        res = SHIFT_DIGIT_DOWN(res);
307        digit_ptr++;
308        if (res == 0)
309        {
310            break;
311        }
312    }
313    if (res != 0)
314    {
315        if (a->dgs_used < a->dgs_alloc)
316        {
317            *digit_ptr = DIGIT_PART(res);
318        }
319        else
320        {
321            if ((digvect = copy_digits(a->dp, a->dgs_used,
322                                       a->dgs_used + 4)) == NULL)
323            {
324                return;         /* big_errno will be set to last error */
325            }
326            digvect[a->dgs_used] = DIGIT_PART(res);
327            free_digits(a->dp, a->dgs_alloc);
328            a->dgs_alloc = a->dgs_used + 4;
329            a->dp = digvect;
330        }
331        a->dgs_used++;
332    }
333}
334
335static DIGIT
336vect_div_digit(DIGIT *digvect, ulong *len, DIGIT d)
337{
338    DIGIT *digit_ptr = digvect + *len - 1, newval;
339    DIGIT2 rest = 0;
340
341    if (d == 0)
342    {
343        big_errno = BIG_DIV_ZERO;
344        return -1;
345    }
346    if (d == 1)
347    {
348        return 0;
349    }
350    while (digit_ptr >= digvect)
351    {
352        rest = (DIGIT2)SHIFT_DIGIT_UP(rest);
353        newval = DIGIT_PART((*digit_ptr + rest) / d);
354        rest = (*digit_ptr + rest) % d;
355        *digit_ptr = newval;
356        digit_ptr--;
357    }
358    if (*len > 1)
359    {
360        if (digvect[*len - 1] == 0)
361        {
362            *len -= 1;
363        }
364    }
365    return DIGIT_PART(rest);
366}
367
368static DIGIT
369udiv_digit(bignum *a, DIGIT d)
370{
371    DIGIT res;
372
373    res = vect_div_digit(a->dp, &a->dgs_used, d);
374    if (zerop(a))
375    {
376        a->sign = BIG_SIGN_0;
377    }
378    return res;
379}
380
381static DIGIT
382vect_mul_digit(DIGIT *digvect, ulong len, DIGIT x)
383{
384    DIGIT *digit_ptr = digvect + len;
385    DIGIT2 res = 0;
386
387    while (digvect < digit_ptr)
388    {
389        res += *digvect * (DIGIT2)x;
390        *digvect = DIGIT_PART(res);
391        res = SHIFT_DIGIT_DOWN(res);
392        digvect++;
393    }
394    return DIGIT_PART(res);
395}
396
397static void
398umul_digit(bignum *a, DIGIT x)
399{
400    DIGIT ovfl, *digvect;
401
402    ovfl = vect_mul_digit(a->dp, a->dgs_used, x);
403    if (ovfl != 0)
404    {
405        if (a->dgs_used < a->dgs_alloc)
406        {
407            a->dp[a->dgs_used] = ovfl;
408        }
409        else
410        {
411            digvect = copy_digits(a->dp,
412                                  a->dgs_used,
413                                  a->dgs_used + 4);
414            digvect[a->dgs_used] = ovfl;
415            free_digits(a->dp, a->dgs_alloc);
416            a->dgs_alloc = a->dgs_used + 4;
417            a->dp = digvect;
418        }
419        a->dgs_used++;
420    }
421}
422
423static int
424ucompare_digits(bignum *a, bignum *b)
425{
426    DIGIT *a_ptr, *b_ptr;
427    int retval = 0;
428
429    if (a->dgs_used  == b->dgs_used)
430    {
431        a_ptr = a->dp + a->dgs_used - 1;
432        b_ptr = b->dp + b->dgs_used - 1;
433        while ((*a_ptr == *b_ptr) && (a_ptr >= a->dp))
434        {
435            a_ptr--;
436            b_ptr--;
437        }
438        if (a_ptr < a->dp)
439        {
440            return retval;
441        }
442        else
443        {
444            if (retval == 0)
445                retval = (*a_ptr > *b_ptr) ? 1 : -1;
446        }
447        return retval;
448    }
449    return (a->dgs_used > b->dgs_used) ? 1 : -1;
450}
451
452static void
453uadd_digits(bignum *a, bignum *b, bignum *c)
454{
455    DIGIT *dp_x, *dp_y, *dp_z, *res_dp, *end_x, *end_y;
456    ulong len_x, len_y;
457    DIGIT2 res = 0;
458
459    if (a->dgs_used > b->dgs_used)
460    {
461        dp_x = a->dp;
462        len_x = a->dgs_used;
463        dp_y = b->dp;
464        len_y = b->dgs_used;
465    }
466    else
467    {
468        dp_x = b->dp;
469        len_x = b->dgs_used;
470        dp_y = a->dp;
471        len_y = a->dgs_used;
472    }
473    end_x = dp_x + len_x;
474    end_y = dp_y + len_y;
475    if (c->dgs_alloc >= len_x)
476    {
477        dp_z = c->dp;
478    }
479    else
480    {
481        if (newsize(&tmp_add.dp, &tmp_add.dgs_alloc,
482                    len_x, len_x + 4) != BIG_OK)
483        {
484            return;
485        }
486        dp_z = tmp_add.dp;
487    }
488    res_dp = dp_z;
489    while (dp_y < end_y)
490    {   
491        res += ((DIGIT2)*dp_x++) + *dp_y++;
492        *res_dp++ = DIGIT_PART(res);
493        res = SHIFT_DIGIT_DOWN(res);
494    }
495    while (dp_x < end_x)
496    {
497        res += *dp_x++;
498        *res_dp++ = DIGIT_PART(res);
499        res = SHIFT_DIGIT_DOWN(res);
500    }
501    if (res != 0)
502    {
503        *res_dp++ = DIGIT_PART(res);
504    }
505
506    if (dp_z != c->dp)
507    {
508        tmp_digits = c->dp;
509        c->dp = tmp_add.dp;
510        tmp_add.dp = tmp_digits;
511
512        tmp_ulong = c->dgs_alloc;
513        c->dgs_alloc = tmp_add.dgs_alloc;
514        tmp_add.dgs_alloc = tmp_ulong;
515    }
516    c->dgs_used = res_dp - c->dp;
517}
518
519static void
520usub_digits(bignum *a, bignum *b, bignum *c)
521{
522    DIGIT *dp_x, *dp_y, *dp_z, *res_dp, *end_x, *end_y;
523    ulong len_x, len_y;
524    DIGIT2 res = 0;
525
526    dp_x = a->dp;
527    len_x = a->dgs_used;
528    dp_y = b->dp;
529    len_y = b->dgs_used;
530
531    end_x = dp_x + len_x;
532    end_y = dp_y + len_y;
533    if (c->dgs_alloc >= len_x)
534    {
535        dp_z = c->dp;
536    }
537    else
538    {
539        if (newsize(&tmp_add.dp, &tmp_add.dgs_alloc,
540                    len_x, len_x + 2) != BIG_OK)
541        {
542            return;
543        }
544        dp_z = tmp_add.dp;
545    }
546    res_dp = dp_z;
547    while (dp_y < end_y)
548    {   
549        res = ((DIGIT2)*dp_x++) - *dp_y++ - RESULT_MINUSP(res);
550        *res_dp++ = DIGIT_PART(res);
551    }
552    while (dp_x < end_x)
553    {
554        res = *dp_x++ - RESULT_MINUSP(res);
555        *res_dp++ = DIGIT_PART(res);
556    }
557#ifdef BIG_CHECK_LIMITS
558    if (RESULT_MINUSP(res) != 0)
559    {
560        big_errno = BIG_ALGERR;
561        return;
562    }
563#endif
564    while ((*--res_dp == 0) && (res_dp > dp_z))
565    {
566        /* Move pointer down until we reach a non-zero */
567    }
568    if (dp_z != c->dp)
569    {
570        tmp_digits = c->dp;
571        c->dp = tmp_add.dp;
572        tmp_add.dp = tmp_digits;
573
574        tmp_ulong = tmp_add.dgs_alloc;
575        tmp_add.dgs_alloc = c->dgs_alloc;
576        c->dgs_alloc = tmp_ulong;
577    }
578    c->dgs_used = res_dp - dp_z + 1;
579}
580
581/*
582 *   This (pseudo) random number generator is not very good.  It has a long
583 * period [ 2 ^ (DIGIT_BITS * 2) (?)] before it starts with the same sequence
584 * again, but the lower bits are "less random" than they should be.  I've
585 * solved it by using word of double length, and returning the upper bits.
586 * Please mail me if you know how to make it into a "more random" generator.
587 *   One important thing though: it will have to reasonably fast.
588 *   The good thing with this one is that it is very portable, but that doesn't
589 * help much when you want _good_ random numbers.
590 *                                              Henrik.Johansson@Nexus.Comm.SE
591 */
592#ifdef UNUSED_CODE
593static DIGIT
594rand(void)
595{
596    static DIGIT2 x1 = 33, x2 = 45; /* Just give them two starting values */
597
598    if (x2 = BIG_RAND_A2 * x2 + BIG_RAND_C2, x2 == 45)
599    {
600        x1 = BIG_RAND_A1 * x1 + BIG_RAND_C1;
601    }
602    return SHIFT_DIGIT_DOWN(x1 + x2); /* Skip least significant bits */
603}
604/* UNUSED_CODE */
605#endif
606
607/* ----------------------------------------------------------------------
608 * All external functions comes here.
609 */
610
611bigerr_t
612big_init_pkg()
613{
614    strcpy (big_base_digits, _(big_base_digits1));      /* known to fit  */
615    init_digit_blocks();
616    big_create(&tmp_string);
617    big_create(&tmp_add);
618    big_create(&tmp_mul);
619    big_create(&tmp_q);
620    big_create(&tmp_r);
621    big_create(&tmp_rand);
622    big_create(&big_one);
623    big_set_long((long)1, &big_one);
624    length_last_string = 10;
625    if ((last_string = calloc(1,length_last_string)) == NULL)
626    {
627        big_errno = BIG_MEMERR;
628    }
629    return big_errno;
630}
631
632void
633big_release_pkg()
634{
635    big_destroy(&tmp_string);
636    big_destroy(&tmp_add);
637    big_destroy(&tmp_mul);
638    big_destroy(&tmp_q);
639    big_destroy(&tmp_r);
640    big_destroy(&tmp_round);
641    big_destroy(&big_one);
642    if (last_string != NULL)
643      free(last_string);
644#ifdef COUNT_ALLOC
645    printf("Allocated digits: %lu\n", dgs_alloc);
646    printf("Freed digits:     %lu\n", digits_freed);
647#endif
648}
649
650bigerr_t
651big_create(a)
652bignum *a;
653{
654  /* Make sure that big_create alway returns a NULL bignum
655   * that can safely be destroyed.
656   */
657    if (a != NULL)
658      {
659        a->dp = NULL;
660      }
661    if (big_errno != BIG_OK)
662    {
663        return big_errno;
664    }
665    if (a == NULL)
666    {
667      big_errno = BIG_ARGERR;
668      return big_errno;
669    }
670
671    a->sign = BIG_SIGN_0;
672    a->dgs_used = 1;
673    if ((a->dp = allocate_digits((long)sizeof(long))) == NULL)
674    {
675        return big_errno;
676    }
677    a->dgs_alloc = sizeof(long);
678    *a->dp = 0;
679    return big_errno;
680}
681
682void
683big_destroy(a)
684bignum *a;
685{
686  /* Mar 4, 2000 R. Wichmann: check for a == NULL */
687  /* free_digits() checks for a->dp == NULL       */
688  if (a != NULL)
689    free_digits(a->dp, a->dgs_alloc);
690}
691
692ulong
693big_bitcount(a)
694bignum *a;
695{
696    int bits = 0;
697    DIGIT high_digit;
698
699    if (big_errno != BIG_OK)
700    {
701        return 0;
702    }
703    high_digit = a->dp[a->dgs_used - 1];
704    while (high_digit != 0)
705    {
706        bits++;
707        high_digit >>= 1;
708    }
709    return DIGIT_BITS * (a->dgs_used - 1) + bits;
710}
711
712bigerr_t
713big_set_big(a, b)
714bignum *a;
715bignum *b;
716{
717    if ((big_errno != BIG_OK) || (a == b))
718    {
719        return big_errno;
720    }
721    if (newsize(&b->dp, &b->dgs_alloc, a->dgs_used, a->dgs_used) != BIG_OK)
722    {
723        return big_errno;
724    }
725    b->dgs_used = a->dgs_used;
726    b->sign = a->sign;
727    digits_cpy(b->dp, a->dp, a->dgs_used);
728    return big_errno;
729}
730
731
732void
733big_set_ulong(n, a)
734ulong n;
735bignum *a;
736{
737    unsigned int i;
738
739    if (big_errno != BIG_OK)
740    {
741        return;
742    }
743    if (n == 0)
744    {
745        *a->dp = 0;
746        a->dgs_used = 1;
747        a->sign = BIG_SIGN_0;
748    }
749    else
750    {
751        a->dgs_used = 0;
752        for (i = 0; i < sizeof(long) / sizeof(DIGIT); i++)
753        {
754            if (n == 0)
755            {
756                break;
757            }
758            a->dgs_used++;
759            a->dp[i] = DIGIT_PART(n);
760            n = SHIFT_DIGIT_DOWN(n);
761        }
762        a->sign = BIG_SIGN_PLUS;
763    }
764    print_digits("set_(u)long: a = ", a);
765}
766
767void
768big_set_long(n, a)
769long n;
770bignum *a;
771{
772    ulong m;
773
774    m = (n < 0) ? - n : n;
775    big_set_ulong(m, a);
776    if (!(n >= 0))
777    {
778        a->sign = BIG_SIGN_MINUS;
779    }
780}
781
782
783bigerr_t
784big_set_string(numstr, base, a)
785char *numstr;
786int base;
787bignum *a;
788{
789    char *src, *maxdigit, *chrptr;
790    DIGIT dig_base, dig_sum, last_base;
791    int cnt, maxcnt;
792
793    if (big_errno != BIG_OK)
794    {
795        return big_errno;
796    }
797
798    big_end_string = numstr;
799    if ((base < 2) || (base > 36) || (numstr == NULL) || (a == NULL))
800    {
801        big_errno = BIG_ARGERR;
802        return big_errno;
803    }
804
805    src            = numstr;
806
807    maxdigit = big_base_digits + (base << 1);
808
809    maxcnt = big_block[base].digCnt;
810    dig_base = big_block[base].dig_base;
811    a->dgs_used = 1;
812    *a->dp = 0;
813    a->sign = BIG_SIGN_PLUS;    /* Assume it will be positive */
814    while (strchr(" \t\n\r", *src) != NULL) /* Skip whitespaces */
815    {
816        src++;
817    }
818    if ((*src == '-') || (*src == '+')) /* First non-white is a sign? */
819    {
820        a->sign = (*src == '-') ? BIG_SIGN_MINUS : BIG_SIGN_PLUS;
821        src++;
822    }
823    chrptr = strchr(big_base_digits, *src);
824    if ((chrptr == NULL) || (chrptr >= maxdigit)) /* Next chr not a digit? */
825    {
826        big_end_string = src;
827        big_errno = BIG_ARGERR;
828        return big_errno;
829    }
830    while (*src == '0')         /* Next chr a '0'? */
831    {
832        src++;                  /* Read past all '0'es */
833    }
834    chrptr = strchr(big_base_digits, *src);
835    if ((chrptr == NULL) || (chrptr >= maxdigit)) /* Next char not a digit */
836    {
837        big_end_string = src;
838        a->sign = BIG_SIGN_0;   /* It was just a plain 0 */
839        return big_errno;
840    }
841    dig_sum = 0;
842    cnt = 0;
843    while ((chrptr = strchr(big_base_digits, *src)),
844           (chrptr != NULL) && (chrptr < maxdigit))
845    {
846        dig_sum = dig_sum * base + ((chrptr - big_base_digits) >> 1);
847        if (++cnt >= maxcnt)
848        {
849            umul_digit(a, dig_base);
850            add_digit(a, dig_sum);
851            dig_sum = 0;
852            cnt = 0;
853        }
854        src++;
855    }
856    if (cnt > 0)
857    {
858        last_base = base;
859        while (--cnt > 0)
860        {
861            last_base *= base;
862        }
863        umul_digit(a, last_base);
864        add_digit(a, dig_sum);
865    }
866    big_end_string = src;
867    return big_errno;
868}
869
870
871int
872big_ulong(a, n)
873bignum *a;
874ulong *n;
875{
876    ulong old_n;
877    DIGIT *dp;
878
879    if (big_errno != BIG_OK)
880    {
881        return FALSE;
882    }
883    if (a->dgs_used > sizeof(ulong) / sizeof(DIGIT))
884    {
885        return FALSE;
886    }
887    dp = a->dp + a->dgs_used - 1;
888    old_n = *n = *dp--;
889    while ((dp >= a->dp) && (old_n < *n))
890    {
891        old_n = *n;
892        *n = SHIFT_DIGIT_UP(*n) + *dp--;
893    } 
894    if (old_n >= *n)
895    {
896        return FALSE;
897    }
898    return FALSE;
899}
900
901int
902big_long(a, n)
903bignum *a;
904long *n;
905{
906    long old_n;
907    DIGIT *dp;
908
909    if (a->dgs_used > sizeof(ulong) / sizeof(DIGIT))
910    {
911        return FALSE;
912    }
913    dp = a->dp + a->dgs_used - 1;
914    old_n = *n = *dp--;
915    while ((dp >= a->dp) && (old_n <= *n))
916    {
917        old_n = *n;
918        *n = SHIFT_DIGIT_UP(*n) + *dp--;
919    } 
920    if (old_n >= *n)
921    {
922        return FALSE;
923    }
924    if (a->sign == BIG_SIGN_MINUS)
925    {
926      *n = -*n;
927    }
928    return FALSE;
929}
930
931
932char *
933big_string(a, base)
934bignum *a;
935int base;
936{
937    ulong bit_length, str_length;
938    char *dst;
939    DIGIT dig_base, rem;
940    int cnt, maxcnt;
941
942    if (big_errno != BIG_OK)
943    {
944        return _(error_string);
945    }
946    big_set_big(a, &tmp_string);
947                                /* Need more room than minimum here... */
948    bit_length = tmp_string.dgs_used * DIGIT_BITS;
949    /*    bit_length = big_bitcount(&tmp_string); */
950    str_length = (ulong)(bit_length / log2tbl[base] + 4);
951    if (str_length > length_last_string)
952    {
953        if (last_string != NULL)
954          free(last_string);
955        if ((last_string = calloc(1,str_length)) == NULL)
956        {
957            big_errno = BIG_MEMERR;
958            return NULL;
959        }
960        length_last_string = str_length;
961    }
962           
963    dst = last_string + length_last_string - 1;
964    *dst-- = '\0';
965    maxcnt = big_block[base].digCnt;
966    dig_base = big_block[base].dig_base;
967    while (tmp_string.dgs_used > 1)
968    {
969        rem = udiv_digit(&tmp_string, dig_base);
970        for (cnt = 0; cnt < maxcnt; cnt++)
971        {
972            *dst-- = big_base_digits[(rem % base) << 1];
973            rem /= base;
974        }
975    }
976    rem = *tmp_string.dp;
977    do
978    {
979        *dst-- = big_base_digits[(rem % base) << 1];
980        rem /= base;
981    } while (rem != 0);
982
983    if (a->sign == BIG_SIGN_MINUS)
984    {
985        *dst = '-';
986    }
987    else
988    {
989        dst++;
990    }
991    return dst;
992}
993
994bigerr_t
995big_negate(a, b)
996bignum *a;
997bignum *b;
998{
999    big_set_big(a, b);
1000    b->sign = NEGATE_SIGN(a->sign);
1001    return big_errno;
1002}
1003
1004int
1005big_sign(a)
1006bignum *a;
1007{
1008    return a->sign;
1009}
1010
1011bigerr_t
1012big_abs(a, b)
1013bignum *a;
1014bignum *b;
1015{
1016    big_set_big(a, b);
1017    if (a->sign == BIG_SIGN_MINUS)
1018    {
1019        b->sign = NEGATE_SIGN(a->sign);
1020    }
1021    return big_errno;
1022}
1023
1024int
1025big_compare(a, b)
1026bignum *a;
1027bignum *b;
1028{
1029    if (a->sign == b->sign)
1030    {
1031        if (a->sign == 0)
1032        {
1033            return 0;
1034        }
1035        return
1036            (a->sign == BIG_SIGN_MINUS)
1037            ? -ucompare_digits(a, b)
1038            : ucompare_digits(a, b);
1039    }
1040    return b->sign - a->sign;
1041}
1042
1043int
1044big_lessp(a, b)
1045bignum *a;
1046bignum *b;
1047{
1048    return big_compare(a, b) < 0;
1049}
1050
1051int
1052big_leqp(a, b)
1053bignum *a;
1054bignum *b;
1055{
1056    return !(big_compare(a, b) > 0);
1057}
1058
1059int
1060big_equalp(a, b)
1061bignum *a;
1062bignum *b;
1063{
1064    return big_compare(a, b) == 0;
1065}
1066
1067int
1068big_geqp(a, b)
1069bignum *a;
1070bignum *b;
1071{
1072    return !(big_compare(a, b) < 0);
1073}
1074
1075int
1076big_greaterp(a, b)
1077bignum *a;
1078bignum *b;
1079{
1080    return big_compare(a, b) > 0;
1081}
1082
1083int
1084big_zerop(a)
1085bignum *a;
1086{
1087    return a->sign == BIG_SIGN_0;
1088}
1089
1090int
1091big_evenp(a)
1092bignum *a;
1093{
1094    return ((*a->dp & 0x01) == 0);
1095}
1096
1097int
1098big_oddp(a)
1099bignum *a;
1100{
1101    return ((*a->dp & 0x01) == 1);
1102}
1103
1104bigerr_t
1105big_add(a, b, c)
1106bignum *a;
1107bignum *b;
1108bignum *c;
1109{
1110    int cmp;
1111
1112    if (big_errno != BIG_OK)
1113    {
1114        return big_errno;
1115    }
1116    print_digits("add:\ta = ", a);
1117    print_digits("\tb = ", b);
1118    if (a->sign == b->sign)
1119    {
1120        uadd_digits(a, b, c);
1121        c->sign = a->sign;
1122    }
1123    else
1124    {
1125        cmp = ucompare_digits(a, b);
1126        if (cmp < 0)
1127        {
1128            usub_digits(b, a, c);
1129            if (zerop(c))
1130            {
1131                c->sign = BIG_SIGN_0;
1132            }
1133            else
1134            {
1135                c->sign = b->sign;
1136            }
1137        }
1138        else if (cmp > 0)
1139        {
1140            usub_digits(a, b, c);
1141            if (zerop(c))
1142            {
1143                c->sign = BIG_SIGN_0;
1144            }
1145            else
1146            {
1147                c->sign = a->sign;
1148            }
1149        }
1150        else
1151        {
1152            c->dgs_used = 1;
1153            *c->dp = 0;
1154            c->sign = BIG_SIGN_0;
1155        }
1156    }
1157    print_digits("\tc = ", c);
1158    return big_errno;
1159}
1160
1161bigerr_t
1162big_sub(a, b, c)
1163bignum *a;
1164bignum *b;
1165bignum *c;
1166{
1167    int cmp;
1168
1169    if (big_errno != BIG_OK)
1170    {
1171        return big_errno;
1172    }
1173    print_digits("sub:\ta = ", a);
1174    print_digits("\tb = ", b);
1175    if (a->sign == BIG_SIGN_0)
1176    {
1177        big_set_big(b, c);
1178        big_negate(c, c);
1179        print_digits("\tc = ", c);
1180        return big_errno;
1181    }
1182    if (b->sign == BIG_SIGN_0)
1183    {
1184        big_set_big(a, c);
1185        print_digits("\tc = ", c);
1186        return big_errno;
1187    }
1188
1189    cmp = ucompare_digits(a, b);
1190    if (cmp <= 0)
1191    {
1192        if (a->sign != b->sign)
1193        {
1194            uadd_digits(a, b, c);
1195            c->sign = a->sign;
1196        }
1197        else
1198        {
1199            usub_digits(b, a, c);
1200            c->sign = (zerop(c) ? BIG_SIGN_0 : NEGATE_SIGN(a->sign));
1201        }
1202    }
1203    else if (cmp > 0)
1204    {
1205        if (a->sign != b->sign)
1206        {
1207            uadd_digits(a, b, c);
1208            c->sign = a->sign;
1209        }
1210        else
1211        {
1212            usub_digits(a, b, c);
1213            c->sign = (zerop(c) ? BIG_SIGN_0 : b->sign);
1214        }
1215    }
1216    else
1217    {
1218        c->dgs_used = 1;
1219        *c->dp = 0;
1220        c->sign = BIG_SIGN_0;
1221    }
1222    print_digits("\tc = ", c);
1223    return big_errno;
1224}
1225
1226bigerr_t
1227big_mul(a, b, c)
1228bignum *a;
1229bignum *b;
1230bignum *c;
1231{
1232    DIGIT *dp_x, *dp_xstart, *dp_xend;
1233    DIGIT *dp_y, *dp_ystart, *dp_yend;
1234    DIGIT *dp_z, *dp_zstart, *dp_zend, *dp_zsumstart;
1235    ulong len_x, len_y /* , len_z */;
1236    DIGIT2 res;
1237    DIGIT tmp_res;           /* Without use of this, gcc (v 1.39) generates */
1238                             /* erroneous code on a sun386 machine */
1239                             /* Should be removed with #ifdef's, when */
1240                             /* not on a sun386, but... */
1241
1242    if (big_errno != BIG_OK)
1243    {
1244        return big_errno;
1245    }
1246    print_digits("mul:\ta = ", a);
1247    print_digits("\tb = ", b);
1248
1249    if (zerop(a) || zerop(b))
1250    {
1251        c->sign = BIG_SIGN_0;
1252        c->dgs_used = 1;
1253        *c->dp = 0;
1254        print_digits("(a=0 || b=0)c = ", c);
1255        return big_errno;
1256    }
1257    if (uonep(a))
1258    {
1259        big_set_big(b, c);
1260        c->sign = (a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS;
1261        print_digits("(abs(a)=1)c = ", c);
1262        return big_errno;
1263    }
1264    if (uonep(b))
1265    {
1266        big_set_big(a, c);
1267        c->sign = (a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS;
1268        print_digits("(abs(b)=1)c = ", c);
1269        return big_errno;
1270    }
1271
1272    if (a->dgs_used < b->dgs_used)
1273    {
1274        dp_xstart = a->dp;
1275        len_x = a->dgs_used;
1276        dp_ystart = b->dp;
1277        len_y = b->dgs_used;
1278    }
1279    else
1280    {
1281        dp_xstart = b->dp;
1282        len_x = b->dgs_used;
1283        dp_ystart = a->dp;
1284        len_y = a->dgs_used;
1285    }
1286    if ((c == a) || (c == b))
1287    {
1288        if (newsize(&tmp_mul.dp, &tmp_mul.dgs_alloc,
1289                    len_x + len_y, len_x + len_y + 2) != BIG_OK)
1290        {
1291            return big_errno;
1292        }
1293        dp_zsumstart = tmp_mul.dp;
1294        /* len_z = tmp_mul.dgs_alloc; */
1295    }   
1296    else
1297    {
1298        if (newsize(&c->dp, &c->dgs_alloc,
1299                    len_x + len_y, len_x + len_y + 2) != BIG_OK)
1300        {
1301            return big_errno;
1302        }
1303        dp_zsumstart = c->dp;
1304        /* len_z = c->dgs_alloc; */
1305    }
1306
1307    dp_xend = dp_xstart + len_x;
1308    dp_yend = dp_ystart + len_y;
1309    dp_zend = dp_zsumstart + len_y;
1310
1311    for (dp_z = dp_zsumstart; dp_z < dp_zend; dp_z++)
1312    {
1313        *dp_z = 0;              /* Zero out rightmost digits */
1314    }
1315
1316    dp_zstart = dp_zsumstart;
1317
1318    for (dp_x = dp_xstart; dp_x < dp_xend; dp_x++)
1319    {
1320        dp_z = dp_zstart;
1321        tmp_res = 0;
1322        for (dp_y = dp_ystart; dp_y < dp_yend; dp_y++)
1323        {
1324            res = (DIGIT2)(*dp_x) * (*dp_y) + (*dp_z) + tmp_res;
1325            *dp_z++ = DIGIT_PART(res);
1326            tmp_res = SHIFT_DIGIT_DOWN(res);
1327        }
1328        *dp_z = tmp_res;
1329        dp_zstart++;
1330    }
1331    if (dp_zsumstart != c->dp)
1332    {
1333        tmp_digits = c->dp;
1334        c->dp = tmp_mul.dp;
1335        tmp_mul.dp = tmp_digits;
1336
1337        tmp_ulong = c->dgs_alloc;
1338        c->dgs_alloc = tmp_mul.dgs_alloc;
1339        tmp_mul.dgs_alloc = tmp_ulong;
1340    }
1341    if (*dp_z == 0)
1342    {
1343        dp_z--;
1344    }
1345    c->dgs_used = dp_z - dp_zsumstart + 1;
1346    c->sign = a->sign * b->sign;
1347    print_digits("\tc = ", c);
1348    return big_errno;
1349}
1350
1351bigerr_t
1352big_trunc(a, b, q, r)
1353bignum *a;
1354bignum *b;
1355bignum *q;
1356bignum *r;
1357{
1358    DIGIT *v_end, *q_end, *r_end, *src, *dst;
1359    DIGIT norm, qhat, t1, t2, t3, v1, v2, u1, u2, u3;
1360    DIGIT2 temp, res, carry;
1361    long a_l, b_l, q_l, r_l;
1362#if 0
1363    ulong i;
1364#endif
1365    ulong j, m, n;
1366    int cmp, q_eq_a, q_eq_b, r_eq_a, r_eq_b;
1367
1368    if (big_errno != BIG_OK)
1369    {
1370        return big_errno;
1371    }
1372    print_digits("div:\ta = ", a);
1373    print_digits("\tb = ", b);
1374    if (zerop(b))
1375    {
1376        big_errno = BIG_DIV_ZERO;
1377        return big_errno;
1378    }
1379   
1380    if (q == r)
1381    {
1382        big_errno = BIG_ARGERR;
1383        return big_errno;
1384    }
1385
1386    if (b->dgs_used == 1)
1387    {
1388        big_set_big(a, q);
1389        q->sign = ((a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS);
1390        *r->dp = udiv_digit(q, *b->dp);
1391        r->dgs_used = 1;
1392        r->sign = (zerop(r) ? BIG_SIGN_0 : a->sign);
1393        print_digits("\t3:q = ", q);
1394        print_digits("\t  r = ", r);
1395        return big_errno;
1396    }
1397   
1398    if (big_long(a, &a_l))      /* Pretend it is a signed value */
1399    {
1400        big_long(b, &b_l);      /* |a| < |b| so this will succeed */
1401        q_l = a_l / b_l;        /* Compute with unsigned operators */
1402        r_l = a_l % b_l;
1403        big_set_long((long)q_l, q);
1404        big_set_long((long)r_l, r);
1405        print_digits("\t4:q = ", q);
1406        print_digits("\t  r = ", r);
1407        return big_errno;
1408    }
1409
1410    cmp = ucompare_digits(a, b); /* Unsigned compare, that is... */
1411    if (cmp < 0)
1412    {
1413        big_set_big(a, r);      /* r = a (don't care about big_errno here) */
1414        q->sign = BIG_SIGN_0;   /* q = 0 */
1415        *q->dp = 0;
1416        q->dgs_used = 1;
1417        print_digits("\t1:q = ", q);
1418        print_digits("\t  r = ", r);
1419        return big_errno;
1420    }
1421    else
1422    if (cmp == 0)
1423    {
1424        q->sign = ((a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS);
1425        *q->dp = 1;     /* q = 1 */
1426        q->dgs_used = 1;
1427        r->sign = BIG_SIGN_0;   /* r = 0 */
1428        *r->dp = 0;
1429        r->dgs_used = 1;
1430        print_digits("\t2:q = ", q);
1431        print_digits("\t  r = ", r);
1432        return big_errno;
1433    }
1434
1435    q_eq_a = (q == a);
1436    q_eq_b = (q == b);
1437    if (q_eq_a || q_eq_b)
1438    {
1439        q = &tmp_q;
1440    }
1441
1442    r_eq_a = (r == a);
1443    r_eq_b = (r == b);
1444    if (r_eq_a || r_eq_b)
1445    {
1446        r = &tmp_r;
1447    }
1448
1449    if (newsize(&r->dp, &r->dgs_alloc, /* At least one more dig in r */
1450                a->dgs_used + 1, a->dgs_used + 2) != BIG_OK)
1451    {
1452        return big_errno;
1453    }
1454    big_set_big(a, r);
1455    r->dp[a->dgs_used] = 0; /* In case no overflow in mult. */
1456   
1457    n = b->dgs_used;
1458    v_end = b->dp + n - 1;
1459    norm = DIGIT_PART((DIGIT_BASE / ((DIGIT2)*v_end + 1)));
1460    if (norm != 1)
1461    {
1462        umul_digit(r, norm);
1463        umul_digit(b, norm);
1464        print_digits("r = ", r);
1465        print_digits("b = ", b);
1466    }
1467    m = a->dgs_used + 1 - b->dgs_used;
1468    r_end = r->dp + a->dgs_used;
1469   
1470    if (newsize(&q->dp, &q->dgs_alloc, m, m + 2) != BIG_OK)
1471    {
1472        return big_errno;
1473    }
1474    q_end = q->dp + m - 1;
1475   
1476    v1 = *v_end;
1477    v2 = *(v_end - 1);
1478   
1479    for (j = 0; j < m; j++)                   /* m steps through division */
1480    {
1481        u1 = *r_end;                            /*  routine */
1482        u2 = *(r_end - 1);
1483        u3 = *(r_end - 2);
1484        qhat = ((u1 == v1) ?
1485                MAX_DIGIT :
1486                DIGIT_PART(((DIGIT2)u1 * DIGIT_BASE + u2) / v1));
1487        while (1)
1488        {
1489            t3 = DIGIT_PART(temp = (DIGIT2)qhat * v2);
1490            t2 = DIGIT_PART(temp = SHIFT_DIGIT_DOWN(temp) + v1 * (DIGIT2)qhat);
1491            t1 = DIGIT_PART(SHIFT_DIGIT_DOWN(temp));
1492#if 0
1493            printf("t1 = %lu, ", (ulong)t1);
1494            printf("t2 = %lu, ", (ulong)t2);
1495            printf("t3 = %lu\n", (ulong)t3);
1496#endif
1497            if (t1 < u1) break;
1498            if (t1 > u1) {--qhat; continue; }
1499            if (t2 < u2) break;
1500            if (t2 > u2) { --qhat; continue; }
1501            if (t3 <= u3) break;
1502            qhat--;
1503        }
1504       
1505        /* This is a tricky one - multiply and subtract simultaneously */
1506        carry = 1;
1507        res = 0;
1508        src = b->dp;
1509        dst = r->dp + m - j - 1;
1510        while (src <= v_end)
1511        {
1512            res = (DIGIT2)qhat * *(src++) + SHIFT_DIGIT_DOWN(res);
1513            carry += (DIGIT2)(*dst) + MAX_DIGIT - DIGIT_PART(res);
1514            *(dst++) = DIGIT_PART(carry);
1515            carry = DIGIT_PART(SHIFT_DIGIT_DOWN(carry));
1516        }
1517        carry += (DIGIT2)(*dst) + MAX_DIGIT - SHIFT_DIGIT_DOWN(res);
1518        *dst = DIGIT_PART(carry);
1519        carry = DIGIT_PART(SHIFT_DIGIT_DOWN(carry));
1520       
1521        if (carry == 0)
1522        {
1523            qhat--;
1524            src = b->dp;
1525            dst = r->dp + m - j - 1;
1526            while (dst <= r_end)
1527            {
1528                carry = (DIGIT2)(*dst) + *src++ + SHIFT_DIGIT_DOWN(carry);
1529                *dst++ = DIGIT_PART(carry);
1530            }
1531            *dst = 0;
1532        }
1533        *(q_end - j) = DIGIT_PART(qhat);
1534        r_end--;
1535    }
1536    r->sign = a->sign;
1537#if 0
1538    i = r->dgs_used;
1539#endif
1540    while ((*r_end == 0) && (r_end > r->dp))
1541    {
1542        r_end--;
1543    }
1544    if (r_end == r->dp)
1545    {
1546        r->dgs_used = 1;
1547        r->sign = BIG_SIGN_0;
1548    }
1549    else
1550    {
1551        r->dgs_used = r_end - r->dp + 1;
1552        r->sign = a->sign;
1553    }
1554    if (norm != 1)
1555    {
1556        udiv_digit(b, norm);
1557        udiv_digit(r, norm);
1558    }
1559    while ((*q_end == 0) && (q_end > q->dp)) /* This is not needed!(?) */
1560    {
1561        q_end--;
1562    }
1563    /* was ifdef'd out already in original */
1564#if 0
1565    i = m - 1;
1566    while ((i > 0) && (q->dp[i--] == 0))
1567    {
1568        /* Loop through all zeroes */
1569    }
1570#endif
1571    q->dgs_used = q_end - q->dp + 1;
1572    q->sign =  ((a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS);
1573   
1574    if (q_eq_a)
1575    {
1576        big_set_big(q, a);
1577    }
1578    else
1579    if (q_eq_b)
1580    {
1581        big_set_big(q, b);
1582    }
1583
1584    if (r_eq_b)
1585    {
1586        big_set_big(r, b);
1587    }
1588    else
1589    if (r_eq_a)
1590    {
1591        big_set_big(r, a);
1592    }
1593
1594    print_digits("\t5:q = ", q);
1595    print_digits("\t  r = ", r);
1596    return big_errno;
1597}
1598
1599bigerr_t
1600big_floor(a, b, q, r)
1601bignum *a;
1602bignum *b;
1603bignum *q;
1604bignum *r;
1605{
1606    int b_eq_qr, sign_eq;
1607
1608    if ( ( b_eq_qr = ((b == q) || (b == r)) ) )
1609    {
1610        big_set_big(b, &tmp_mul);
1611    }
1612    sign_eq = a->sign == b->sign;
1613    big_trunc(a, b, q, r);
1614    if (sign_eq)
1615    {
1616        return big_errno;
1617    }
1618    if (!zerop(r))
1619    {
1620        if (b_eq_qr)
1621        {
1622            big_add(r, &tmp_mul, r);
1623        }
1624        else
1625        {
1626            big_add(r, b, r);
1627        }
1628        print_digits("big_one = ", &big_one);
1629        big_sub(q, &big_one, q);
1630    }
1631    return big_errno;
1632}
1633
1634bigerr_t
1635big_ceil(a, b, q, r)
1636bignum *a;
1637bignum *b;
1638bignum *q;
1639bignum *r;
1640{
1641    int b_eq_qr, sign_diff;
1642
1643    if ( ( b_eq_qr = ((b == q) || (b == r)) ) )
1644    {
1645        big_set_big(b, &tmp_mul);
1646    }
1647    sign_diff = a->sign != b->sign;
1648    big_trunc(a, b, q, r);
1649    if (sign_diff)
1650    {
1651        return big_errno;
1652    }
1653    if (!zerop(r))
1654    {
1655        if (b_eq_qr)
1656        {
1657            big_sub(r, &tmp_mul, r);
1658        }
1659        else
1660        {
1661            big_sub(r, b, r);
1662        }
1663        big_add(q, &big_one, q);
1664    }
1665    return big_errno;
1666}
1667
1668#ifdef UNUSED_CODE
1669/* This one doesn't work to 100%.  I was a little braindamaged when I wrote
1670 * this, but I'll eventually fix it.   Zzzzzzz.
1671 */
1672bigerr_t
1673big_round(a, b, q, r)
1674bignum *a;
1675bignum *b;
1676bignum *q;
1677bignum *r;
1678{
1679    int b_eq_qr, b_neg_p, a_sgn_neq_b_sgn;
1680
1681    if ( ( b_eq_qr = ((b == q) || (b == r)) ) )
1682    {
1683        big_set_big(b, &tmp_round);
1684    }
1685    b_neg_p = b->sign == BIG_SIGN_MINUS;
1686    a_sgn_neq_b_sgn = a->sign != b->sign;
1687    big_trunc(a, b, q, r);
1688
1689    big_set_big(r, &tmp_add);
1690    umul_digit(&tmp_add, 2);
1691    if (ucompare_digits(&tmp_add, b) > 0) /* |2 * r| > |b| */
1692    {
1693        if (q->sign == BIG_SIGN_0)
1694        {
1695            if (a_sgn_neq_b_sgn)
1696            {
1697                big_sub(q, &big_one, q);
1698            }
1699            else
1700            {
1701                big_add(q, &big_one, q);
1702            }
1703        }
1704        else
1705        {
1706            if (q->sign == BIG_SIGN_MINUS)
1707            {
1708                big_sub(q, &big_one, q);
1709            }
1710            else
1711            {
1712                big_add(q, &big_one, q);
1713            }
1714        }
1715        if (b_eq_qr)
1716        {
1717            if (q->sign == BIG_SIGN_PLUS)
1718            {
1719                big_sub(r, &tmp_round, r);
1720            }
1721            else
1722            {
1723                big_add(r, &tmp_round, r);
1724            }
1725        }
1726        else
1727        {
1728            if (q->sign == BIG_SIGN_PLUS)
1729            {
1730                big_sub(r, b, r);
1731            }
1732            else
1733            {
1734                big_add(r, b, r);
1735            }
1736        }
1737    }
1738    return big_errno;
1739}
1740
1741bigerr_t
1742big_random(a, b)
1743bignum *a;
1744bignum *b;
1745{
1746    unsigned long i;
1747    int a_sgn = a->sign;
1748
1749    if (big_errno != BIG_OK)
1750    {
1751        return big_errno;
1752    }
1753    if (zerop(a))               /* a = 0 -> big_random => 0 (special case) */
1754    {
1755        *b->dp = 0;
1756        b->dgs_used = 1;
1757        b->sign = a_sgn;
1758        return big_errno;
1759    }
1760    if (newsize(&tmp_rand.dp, &tmp_rand.dgs_alloc,
1761                a->dgs_used + 1, a->dgs_used + 1) != BIG_OK)
1762    {
1763        return big_errno;
1764    }
1765    for (i = 0; i <= a->dgs_used; i++)
1766    {
1767        tmp_rand.dp[i] = rand();
1768    }
1769    while (tmp_rand.dp[a->dgs_used] == 0) /* Make sure high digit is non-0 */
1770    {
1771        tmp_rand.dp[a->dgs_used] = rand();
1772    }
1773    tmp_rand.dgs_used = a->dgs_used + 1;
1774    tmp_rand.sign = BIG_SIGN_PLUS;
1775    a->sign = BIG_SIGN_PLUS;
1776    big_trunc(&tmp_rand, a, &tmp_q, b); /* Dangerous to use tmp_q here... */
1777    a->sign = a_sgn;
1778    b->sign = zerop(&tmp_rand) ? BIG_SIGN_0 : a_sgn;
1779    return big_errno;
1780}
1781
1782/* UNUSED_CODE */
1783#endif
1784
1785/* ----------------------------------------------------------------------
1786 * External functions that do not need to know anything about the internal
1787 * representation of a bignum.
1788 */
1789
1790#ifdef UNUSED_CODE
1791
1792int
1793big_expt(a, z, x)
1794bignum *a;
1795unsigned long z;
1796bignum *x;
1797{
1798    bignum b;
1799
1800    big_create(&b);
1801    big_set_big(a, &b);
1802    big_set_long((long)1, x);
1803   
1804    while ((z != 0) && (big_errno == BIG_OK))
1805    {
1806        while ((z & 0x01) == 0)
1807        {
1808            z >>= 1;
1809            big_mul(&b, &b, &b);
1810        }
1811        z -= 1;
1812        big_mul(x, &b, x);
1813    }
1814
1815    big_destroy(&b);
1816    return big_errno;
1817}
1818
1819/* UNUSED_CODE */
1820#endif
1821
1822int
1823big_exptmod(a_in, z_in, n, x)
1824bignum *a_in;
1825bignum *z_in;
1826bignum *n;
1827bignum *x;
1828{
1829    bignum a, z, b0, b1, b2, dummy;
1830
1831    big_create(&a);
1832    big_create(&z);
1833    big_create(&b0);
1834    big_create(&b1);
1835    big_create(&b2);
1836    big_create(&dummy);
1837
1838    big_set_big(a_in, &a);
1839    big_set_big(z_in, &z);
1840    big_set_long((long)1, x);
1841    big_set_long((long)0, &b0);
1842    big_set_long((long)1, &b1);
1843    big_set_long((long)2, &b2);
1844
1845    /* No foolproof testing on big_errno - it really ought to be done */
1846    while ((big_compare(&z, &b0) != 0) && (big_errno == BIG_OK))
1847    {
1848        while (big_evenp(&z) && (big_errno == BIG_OK))
1849        {
1850            big_trunc(&z, &b2, &z, &dummy);
1851            big_mul(&a, &a, &a);
1852            big_trunc(&a, n, &dummy, &a);
1853        }
1854        big_sub(&z, &b1, &z);
1855        big_mul(x, &a, x);
1856        big_trunc(x, n, &dummy, x);
1857    }
1858
1859    big_destroy(&dummy);
1860    big_destroy(&b2);
1861    big_destroy(&b1);
1862    big_destroy(&b0);
1863    big_destroy(&z);
1864    big_destroy(&a);
1865
1866    return big_errno;
1867}
1868
1869#ifdef UNUSED_CODE
1870
1871int
1872big_gcd(a, b, g)
1873bignum *a;
1874bignum *b;
1875bignum *g;
1876{
1877    bignum a1, b1, tmp;
1878
1879    if (big_zerop(b))
1880    {
1881        big_abs(a, g);
1882        return big_errno;
1883    }
1884
1885    big_create(&a1);
1886    big_create(&b1);
1887    big_create(&tmp);
1888
1889    big_abs(a, &a1);
1890    big_abs(b, &b1);
1891
1892    while (!big_zerop(&b1) && (big_errno == BIG_OK))
1893    {
1894        big_floor(&a1, &b1, &tmp, &a1);
1895        if (big_zerop(&a1))
1896        {
1897            break;
1898        }
1899        big_floor(&b1, &a1, &tmp, &b1);
1900    }
1901
1902    if (big_zerop(&a1))
1903    {
1904        big_set_big(&b1, g);
1905    }
1906    else
1907    {
1908        big_set_big(&a1, g);
1909    }
1910
1911    big_destroy(&tmp);
1912    big_destroy(&b1);
1913    big_destroy(&a1);
1914
1915    return big_errno;
1916}
1917/* UNUSED_CODE */
1918#endif
1919
1920/* #if (defined (SH_WITH_CLIENT) || defined (SH_WITH_SERVER)) */
1921#endif
1922
Note: See TracBrowser for help on using the repository browser.