source: trunk/src/bignum.c @ 452

Last change on this file since 452 was 435, checked in by katerina, 9 years ago

Minor bugfixes.

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