|
|
1.1 root 1:
2: //
3: // C++ multiple precision integer arithmetic library.
4: //
5: // Integers are encoded as INT_REPs. An INT_REP contains a
6: // reference count, to be used by the user data structure
7: // INT. The beg pointer points to the beginning of the
8: // heap allocated array of longs; last points to beg plus
9: // the length of the array. Integers are represented in
10: // base 2^16, with the least significant digit written at
11: // beg, the next at beg+1, etc. The wt pointer points
12: // to the long after the most significant digit. Negatives
13: // are stored in twos complement notation, with the most
14: // significant digit a long -1.
15: // Examples:
16: // 12345678901(10) <=> 2dfdc1c35 <=> beg-> 1c35,dfdc,2
17: // 3950617284 <=> eb79a2c4 <=> beg-> a2c4,eb79
18: // -3950617284 <=> <=> beg-> 5d3c,1486,-1
19: // The representations are always normalized to the shortest
20: // equivalent representation. Thus, zero has length 0;
21: // beg-> 5d3c,0,0 becomes beg->5d3c, and beg-> ffff,-1 becomes
22: // beg-> -1.
23: //
24: // 2^16 was chosen as the base since it is the largest base
25: // in which single precision operations can be carried out
26: // using machine longs (assuming longs have >= 32 bits).
27: // (This is not entirely true, as an operation in division
28: // potentially produces a value of more than 32 bits; a trick
29: // is used to avoid this problem.) Most of the algorithms
30: // are inherently independent of BASE, as long as BASE^2
31: // fits in a long. However, there are some short cuts used
32: // that make use of the fact that the digits contain 16 bits
33: // and that a long corresponds to 2 digits.
34: //
35:
36: #include "bignum.h"
37: #include <ctype.h>
38:
39: extern "C" {
40: extern int strlen(const char*);
41: extern char *realloc (char *, int);
42: extern void exit(int);
43: }
44:
45: #define NUMHDS 64 /* No. of INT_REPs allocated in a block. */
46: #define BASE 65536 /* Base BASE arithmetic : 2^16 */
47: #define BlockLen 9 /* Maximum no. of decimal digits fitting in a long. */
48:
49: inline int Min(int a,int b) { return(a < b ? a : b); }
50: inline int Max(int a,int b) { return(a > b ? a : b); }
51:
52: #define equalsOne(r) ((r->wt == r->beg+1) && (*(r->beg) == 1))
53: #define MakeZero(t) t->beg=new long[0];t->wt=t->last=t->beg
54:
55: enum Mode {
56: decimal, hexadecimal, octal
57: };
58:
59: static INT_REP *hfree; // Pointer to free list of headers.
60: static INT_REP OneRep("1"); // Representation of 1.
61: static INT_REP TwoRep("2"); // Representation of 2.
62: static INT_REP *Remainder; // Remainder from division, if requested.
63: static INT_REP BigTen(1000000000L); // Largest power of 10 fitting in long.
64: INT_REP ZeroRep(0, 0); // Representation of 0.
65: static INT Zero; // Guarantees refcnt of ZeroRep always > 0.
66:
67:
68: /* IROutBuf:
69: * Construct for handling reverse strings of
70: * indefinite length.
71: */
72:
73: #define IRBufSize 256
74:
75: struct IROutBuf {
76: char *beg;
77: char *last;
78: char *wt;
79:
80: IROutBuf() { beg = new char[IRBufSize]; wt = beg; last = beg + IRBufSize; }
81: ~IROutBuf() { delete beg; }
82:
83: int Length () { return ((int)(wt - beg)); }
84:
85: void Flush (register iostream &i)
86: {
87: register char *p = wt;
88:
89: while (p > beg) {
90: i.put(*(--p));
91: }
92: wt = beg;
93: }
94:
95: void Flush (register OutFn outf)
96: {
97: register char *p = wt;
98:
99: while (p > beg) {
100: (*outf)(*(--p));
101: }
102: (*outf)('\0');
103: wt = beg;
104: }
105:
106: void Flush (register char *buf, int len)
107: {
108: register char *p = wt;
109: register char *endp = wt - len;
110:
111: while (p > endp) {
112: *buf++ = *(--p);
113: }
114: *buf = '\0';
115: wt = beg;
116: }
117:
118: void Put (char c)
119: {
120: int newsize, oldsize;
121:
122: if (wt == last) { // Need more storage.
123: oldsize = (int)(last - beg);
124: newsize = oldsize + IRBufSize;
125: beg = realloc (beg, newsize);
126: last = beg + newsize;
127: wt = beg + oldsize;
128: }
129: *wt++ = c;
130: }
131:
132: };
133: static IROutBuf obuf;
134:
135: /* Fatal:
136: * Fatal error; print message and exit.
137: */
138: static void
139: Fatal (char *msg)
140: {
141:
142: fprintf (stderr, "bignum - %s\n", msg);
143: exit (1);
144: }
145:
146: /* MakeLong:
147: * Convert decimal string to unsigned long.
148: * We assume that the decimal will fit.
149: */
150: unsigned long MakeLong (char *sp, char *ep)
151: {
152: unsigned long ret = 0;
153: int c;
154:
155: while (sp < ep) {
156: c = *sp++;
157: if(!isdigit(c))
158: Fatal ("MakeLong - illegal decimal digit.");
159:
160: ret = 10*ret + (c - '0');
161: }
162: return(ret);
163:
164: }
165:
166: /* HighBit:
167: * Given positive long v < 2^16,
168: * return number of bits val must be left
169: * shifted to be >= 2^15.
170: */
171: int HighBit (register long v)
172: {
173: register int ret = 0;
174:
175: while ((v & 0x8000) == 0) {
176: v <<= 1;
177: ret++;
178: }
179:
180: return (ret);
181: }
182:
183: /* Shift:
184: * Given an array of longs, all < 2^16, left shift them shift
185: * bits, putting numfill of the resultant words in the storage
186: * supplied by fill, using 16 bits per long.
187: */
188: void Shift (long *sp, long *ep, int shift, long *fill, int numfill)
189: {
190: long val;
191:
192: // Easy case; no shifting, just copy.
193: if (shift == 0) {
194: while (numfill--)
195: *fill++ = (sp >= ep ? *sp-- : 0);
196: return;
197: }
198:
199: *fill = ((*sp--) << shift) & 0xFFFF;
200: while (numfill--) {
201: val = (sp >= ep ? *sp-- : 0);
202: *fill++ |= (val >> (16 - shift)) & 0xFFFF;
203: if (numfill == 0)
204: break;
205: *fill = (val << shift) & 0xFFFF;
206: }
207:
208: }
209:
210: /* IRMorehd:
211: * Allocate more headers.
212: */
213: void INT_REP::IRMorehd ()
214: {
215: register INT_REP *h, *kk;
216:
217: hfree = h = (INT_REP *)new char [sizeof(INT_REP) * NUMHDS];
218: if(hfree == 0)
219: Fatal ("IRMorehd : no more memory");
220:
221: // Link the headers together using the wt pointer.
222: kk = h;
223: while (h < hfree + NUMHDS) {
224: h->wt = ((long *)++kk);
225: h++;
226: }
227: h--;
228: h->wt = ((long *)0);
229: }
230:
231: /* IRAllocRep:
232: * Return next available header;
233: * allocate more if necessary.
234: */
235: static inline INT_REP *IRAllocRep () {
236: register INT_REP *next;
237:
238: if (hfree == 0) {
239: INT_REP::IRMorehd ();
240: }
241: next = hfree;
242: hfree = (INT_REP *)(next->wt);
243: return (next);
244: }
245:
246: /* StripZero:
247: * Strip trailing zeros and reset wt pointer.
248: */
249: void INT_REP::StripZero()
250: {
251: register long *pp;
252:
253: if(wt == beg)
254: return;
255: pp = wt - 1;
256: while ((*pp == 0) && (pp > beg))
257: pp--;
258: if (*pp != 0)
259: pp++;
260: wt = pp;
261: }
262:
263: /* StripMinus:
264: * If this is negative,
265: * strip trailing -1's and reset wt pointer.
266: */
267: void INT_REP::StripMinus()
268: {
269: register long *pp;
270:
271: if(isNegative()) {
272: pp = wt - 1; // Points at -1.
273: if (pp > beg) {
274: pp--;
275: while ((*pp == (BASE-1)) && (pp > beg))
276: pp--;
277: if (*pp != (BASE-1))
278: pp++;
279: *pp++ = -1;
280: wt = pp;
281: }
282: }
283: }
284:
285: /* Twos:
286: * Return largest power of 2 dividing this.
287: * Assume this is nonzero.
288: */
289: int INT_REP::Twos()
290: {
291: register int ret = 0;
292: register long *p, val;
293:
294: p = beg;
295: while(*p == 0) {
296: ret += 16;
297: p++;
298: }
299:
300: val = *p;
301: while ((val & 0x1) == 0) {
302: val >>= 1;
303: ret++;
304: }
305:
306: return (ret);
307: }
308:
309: /* FillLong:
310: * Use long to fill this.
311: * Assume storage has already been allocated.
312: */
313: inline void INT_REP::FillLong(unsigned long val)
314: {
315: beg[0] = val & 0xFFFF;
316: beg[1] = (val >> 16) & 0xFFFF;
317: if (val & 0xFFFF0000)
318: wt = beg + 2;
319: else if (val & 0xFFFF)
320: wt = beg + 1;
321: else
322: wt = beg;
323: }
324:
325: /* Trunc:
326: * Use low order digits to fill a long.
327: * If digits are lost, set error = 1;
328: * else error = 0.
329: */
330: long INT_REP::Trunc(int &error)
331: {
332: long ret = 0;
333: int len;
334: register long *ptr;
335:
336: error = 0;
337: if(isZero())
338: return(ret);
339:
340: len = Length();
341: ptr = beg;
342:
343: if (isNegative()) {
344: if (len == 1)
345: ret = -1;
346: else {
347: ret = *ptr++;
348: ret |= (*ptr << 16) & 0x7FFFFFFF;
349: if ((len > 3) || ((*ptr & 0x8000) == 0))
350: error = 1;
351: }
352: }
353: else {
354: ret = *ptr++;
355: if (len > 1) {
356: ret |= *ptr << 16;
357: if ((len > 2) || (*ptr & 0x8000))
358: error = 1;
359: }
360: }
361:
362: return (ret);
363: }
364:
365: /* Comp:
366: * this < b => -1
367: * this == b => 0
368: * this > b => 1
369: */
370: int INT_REP::Comp(INT_REP *b)
371: {
372: INT_REP *diff;
373: int df;
374:
375: diff = Sub(b);
376: if (diff->isZero())
377: df = 0;
378: else if (diff->isNegative())
379: df = -1;
380: else
381: df = 1;
382:
383: delete diff;
384: return (df);
385: }
386:
387:
388: /* Add:
389: * Return new rep, which is the sum of
390: * this and b.
391: */
392: INT_REP *INT_REP::Add(INT_REP *b)
393: {
394: register INT_REP *p;
395: register long n;
396: register int carry;
397: int size;
398: long *ap, *bp, *pp;
399:
400: if ((this == 0) || (b == 0))
401: Fatal ("Add - Null pointers.");
402:
403: if (isZero()) {
404: p = new INT_REP(b);
405: return (p);
406: }
407: if (b->isZero()) {
408: p = new INT_REP(this);
409: return (p);
410: }
411:
412: size = Max (Length(), b->Length());
413: p = new INT_REP(size);
414:
415: ap = this->beg;
416: bp = b->beg;
417: pp = p->beg;
418: carry = 0;
419: while(--size >= 0){
420: n = (ap == wt ? 0 : *ap++) + (bp == b->wt ? 0 : *bp++) + carry;
421: if (n >= BASE) {
422: carry = 1;
423: n -= BASE;
424: }
425: else if (n < 0) {
426: carry = -1;
427: n += BASE;
428: }
429: else carry = 0;
430: *pp++ = n;
431: }
432: p->wt = pp;
433: if(carry != 0)
434: p->Putc(carry);
435:
436: // Remove trailing zeros.
437: p->StripZero();
438:
439: // Remove redundant -1's in front of terminating -1.
440: p->StripMinus();
441:
442: return(p);
443: }
444:
445: /* Sub:
446: * Return new rep, which is the difference of
447: * this and b.
448: */
449: INT_REP *INT_REP::Sub(INT_REP *b)
450: {
451: INT_REP *nb, *diff;
452:
453: if ((this == 0) || (b == 0))
454: Fatal ("Sub - Null pointers.");
455:
456: nb = new INT_REP (b);
457: nb->Chsign();
458:
459: diff = Add(nb);
460: delete nb;
461: return (diff);
462: }
463:
464: /* DivMod:
465: * Return new rep, which is:
466: * mode == 0 => quotient of this and ddivr;
467: * mode == 1 => remainder of quotient of this and ddivr;
468: * mode == 2 => quotient of this and ddivr, with the
469: * remainder stored in Remainder.
470: * Based on the multiple precision division algorithm
471: * in Knuth, v.2.
472: */
473: INT_REP *INT_REP::DivMod(INT_REP *ddivr, int mode)
474: {
475: register INT_REP *quot, *divd, *divr;
476: long *quotp, *divrp, *divdp;
477: int remsign, divsign, offset;
478: long carry, borrow;
479: long *ap, q;
480: unsigned long uval, delta;
481: long val;
482: int divlen, n, exp;
483: long varray[2], uarray[3];
484:
485: if ((this == 0) || (ddivr == 0))
486: Fatal ("DivMod - Null arguments.");
487:
488: if (ddivr->isZero())
489: Fatal ("DivMod : divide by zero.");
490:
491: // Make divisor and dividend positive.
492: Remainder = 0;
493: divsign = remsign = 0;
494: if (ddivr->isNegative ()) {
495: divr = new INT_REP(ddivr);
496: divr->Chsign();
497: divsign = ~divsign;
498: }
499: else
500: divr = ddivr;
501:
502: divd = new INT_REP(this);
503: if (divd->isNegative()) {
504: divd->Chsign();
505: divsign = ~divsign;
506: remsign = ~remsign;
507: }
508:
509: offset = divd->Length() - divr->Length();
510: if(offset < 0) {
511: quot = new INT_REP(0);
512: goto ddone;
513: }
514:
515: // Allocate quotient.
516: quot = new INT_REP (offset + 1);
517: quot->wt = quot->last;
518:
519: val = *(divr->wt - 1); // Most significant "digit" in divisor.
520: divlen = divr->Length();
521:
522: // If divisor is 1 digit, handle specially.
523: if (divlen == 1) {
524: quotp = quot->beg + offset;
525: divdp = divd->beg + offset;
526: carry = 0;
527:
528: while (offset-- >= 0) {
529: uval = carry*BASE + *divdp--;
530: *quotp-- = uval/val;
531: carry = uval%val;
532: }
533:
534: // Set remainder.
535: divd->beg[0] = carry;
536: divd->wt = divd->beg + 1;
537:
538: goto ddone;
539: }
540:
541: // General case.
542: // Initialize. We are essentially finding the least power of 2
543: // that we can multiply times the most significant digit of the divisor
544: // and make it >= BASE/2. The algorithm then guarantees that our
545: // "guess" for the quotient will be at most two greater than the
546: // correct answer.
547:
548: divd->Putc (0);
549: exp = HighBit (val); // 0 <= exp < 16.
550: Shift (divr->wt - 1, divr->beg, exp, varray, 2);
551: divdp = divd->wt - 1;
552: divd->wt -= offset;
553: quotp = quot->wt - 1;
554:
555: // fprintf(stderr, "exp = %d , v1 = %x, v2 = %x\n", exp, varray[0], varray[1]);
556: while (offset-- >= 0) {
557: // divd->Prt();
558: Shift (divdp, divd->beg, exp, uarray, 3);
559: // fprintf(stderr, "u1 = %x, u2 = %x, u3 = %x\n", uarray[0], uarray[1], uarray[2]);
560:
561: // Guess a quotient.
562: if (uarray[0] == varray[0])
563: q = BASE-1;
564: else
565: q = ((unsigned long)(BASE*uarray[0] + uarray[1]))/varray[0];
566: // fprintf(stderr, "initial q = %x\n",q);
567:
568: // Adjust.
569: // while q*v1 > (BASE*u0 + u1 - q*v0)BASE + u2, decrement q.
570: // To avoid an overflow problem, first check if
571: // (BASE*u0 + u1 - q*v0) >= BASE; if it is, the right side is
572: // >= BASE^2, while the left is < BASE^2 and we are done.
573: // Otherwise, the right side is < BASE^2, and we can do an
574: // ordinary arithmetic check.
575: while (1) {
576: delta = BASE*uarray[0] + uarray[1] - q*varray[0];
577: if (delta >= BASE)
578: break;
579: if ((unsigned long)(varray[1]*q) <= (unsigned long)(delta*BASE+uarray[2]))
580: break;
581: q--;
582: }
583: // fprintf(stderr, "final q = %x\n",q);
584:
585: // Multiply divisor by q, and subtract from dividend.
586: n = divlen;
587: carry = 0;
588: borrow = 0;
589: ap = divdp - n;
590: divrp = divr->beg;
591:
592: while (n >= 0) {
593: if (n > 0)
594: uval = q*(*divrp++) + carry;
595: else
596: uval = carry;
597:
598: if (uval >= BASE) {
599: carry = uval/BASE;
600: uval %= BASE;
601: }
602: else
603: carry = 0;
604:
605: val = (*ap) - uval - borrow;
606: // fprintf(stderr, "n = %d, borrow = %x, carry = %x\n",
607: // n, borrow, carry);
608: // fprintf(stderr, "*ap = %x, val = %x, uval = %x\n",
609: // *ap, val, uval);
610: if (val < 0) {
611: val += BASE;
612: borrow = 1;
613: }
614: else
615: borrow = 0;
616:
617: *ap++ = val;
618: n--;
619: }
620:
621: // If guess is still too high (remainder negative => had to borrow),
622: // subtract one from quotient and add one copy of divisor
623: // back into remainder.
624: if (borrow) {
625: *quotp-- = q - 1;
626: n = divlen;
627: carry = 0;
628: ap = divdp - n;
629: divrp = divr->beg;
630:
631: while (n >= 0) {
632: if (n > 0)
633: uval = (*ap) + (*divrp++) + carry;
634: else
635: uval = (*ap) + carry;
636: // fprintf(stderr, "borrow n = %d, carry = %x, *ap = %x, uval = %x\n",
637: // n, carry, *ap, uval);
638:
639: if (uval >= BASE) {
640: uval -= BASE;
641: carry = 1;
642: }
643: else
644: carry = 0;
645:
646: *ap++ = uval;
647: n--;
648: }
649: }
650: else
651: *quotp-- = q;
652:
653: divdp--;
654: }
655:
656: ddone:
657: if (divr != ddivr) delete divr;
658: switch (mode) {
659: case 0 : // Return quotient; delete remainder.
660: delete divd;
661: break;
662: case 1 : // Return remainder; delete quotient and pretend
663: // that the remainder is the quotient.
664: delete quot;
665: quot = divd;
666: divsign = remsign;
667: break;
668: case 2 : // Return quotient; store remainder in Remainder.
669: divd->StripZero();
670: if (remsign < 0) divd->Chsign(); // Use correct sign.
671: Remainder = divd;
672: break;
673: }
674:
675: // Remove prefix zeros.
676: quot->StripZero();
677: if (divsign < 0) quot->Chsign(); // Use correct sign.
678:
679: return (quot);
680: }
681:
682: /* Exp:
683: * Return new rep, which is the first argument
684: * raised to the second argument power.
685: */
686: INT_REP *INT_REP::Exp(INT_REP *exp)
687: {
688: INT_REP *r;
689: INT_REP *e;
690: INT_REP *p;
691: INT_REP *t, *cp, *e1;
692: int neg = 0;
693:
694: r = new INT_REP(&OneRep);
695:
696: // Check for special cases.
697: if((exp->isZero()) || equalsOne(this))
698: return (r);
699: if(isZero()) {
700: r->SetZero();
701: return (r);
702: }
703:
704: // Create copies of arguments.
705: p = new INT_REP(this);
706: e = new INT_REP(exp);
707:
708: // Check for negative exponent.
709: if (e->isNegative()) {
710: neg++;
711: e->Chsign();
712: }
713:
714: // Compute exponentiation.
715: while (e->isNonzero()) {
716: e1 = e->DivMod(&TwoRep,2);
717: delete e;
718: e = e1;
719:
720: if (Remainder->isNonzero()) {
721: e1 = p->Mult(r);
722: delete r;
723: r = e1;
724: }
725: delete Remainder;
726:
727: t = new INT_REP(p);
728: cp = p->Mult(t);
729: delete p;
730: delete t;
731: p = cp;
732: }
733:
734: delete e;
735: delete p;
736:
737: // Adjust for negative exponent.
738: if (neg) {
739: t = OneRep.Div(r);
740: delete r;
741: return (t);
742: }
743: else
744: return (r);
745:
746: }
747:
748: /* Mult:
749: * Return new rep, which is the product of
750: * this and ddivr.
751: */
752: INT_REP *INT_REP::Mult(INT_REP *q)
753: {
754: register INT_REP *mp,*mq,*mr;
755: int sign,offset,carry;
756: unsigned long mt;
757: long cq, cp, mcr;
758: long *pptr, *qptr, *rptr;
759:
760: if ((this == 0) || (q == 0))
761: Fatal ("Mult - Null pointers.");
762:
763: offset = sign = 0;
764: if (isNegative()) {
765: mp = new INT_REP (this);
766: mp->Chsign();
767: sign = ~sign;
768: }
769: else
770: mp = this;
771:
772: if (q->isNegative()) {
773: mq = new INT_REP (q);
774: mq->Chsign();
775: sign = ~sign;
776: }
777: else
778: mq = q;
779:
780: mr = new INT_REP (mp->Length() + mq->Length());
781: mr->Zero();
782: qptr = mq->beg;
783: while (qptr < mq->wt) {
784: cq = *qptr++;
785: pptr = mp->beg;
786: rptr = mr->beg + offset;
787: carry = 0;
788: while (pptr < mp->wt) {
789: cp = *pptr++;
790: mcr = (rptr == mr->wt ? 0 : *rptr);
791: mt = cp*cq + carry + mcr;
792: carry = mt>>16; // mt/BASE.
793: mr->Alterc(mt&0xFFFF, rptr++); // mt%BASE.
794: }
795: offset++;
796: if(carry != 0){
797: mcr = (rptr == mr->wt ? 0 : *rptr);
798: mr->Alterc(mcr+carry, rptr); // Why is mcr+carry < BASE?
799: }
800: }
801:
802: if(sign < 0){
803: mr->Chsign();
804: }
805: if(mp != this) delete mp;
806: if(mq != q) delete mq;
807:
808: return (mr);
809: }
810:
811: /* INT_REP:
812: * Create new rep, with
813: * size bytes of storage.
814: * If setRef == 0, do not initialize refcnt. This is used by
815: * ZeroRep in case it gets initialized after some INT with
816: * no arguments.
817: */
818: INT_REP::INT_REP (int size, int setRef)
819: {
820: long *ptr;
821:
822: if (this == 0)
823: this = IRAllocRep ();
824:
825: ptr = new long[(unsigned)size];
826: if (ptr == 0)
827: Fatal ("INT_REP - Out of memory.");
828:
829: wt = beg = ptr;
830: last = ptr + size;
831: if (setRef) refcnt = 0;
832: }
833:
834: /* INT_REP:
835: * Create new rep,
836: * using supplied long.
837: */
838: INT_REP::INT_REP (long val)
839: {
840: int neg = 0;
841: int size = 2; // Number of longs needed to store positive 32 bits.
842:
843: if (this == 0)
844: this = IRAllocRep ();
845:
846: if (val == 0) { // If zero, ...
847: MakeZero (this);
848: return;
849: }
850:
851: if (val < 0) {
852: neg++;
853: size++; // May need extra word for sign.
854: }
855: beg = new long[size];
856: if (beg == 0)
857: Fatal ("INT_REP - Out of memory.");
858: last = beg + size;
859: refcnt = 0;
860:
861: // Store bits.
862: FillLong((unsigned long)val);
863:
864: // If negative, add sign, and strip redundant signs.
865: if (neg) {
866: Putc(-1);
867: StripMinus();
868: }
869:
870: }
871:
872: /* INT_REP:
873: * Create new rep,
874: * duplicating old rep.
875: */
876: INT_REP::INT_REP (INT_REP *old)
877: {
878: long *ptr, *oldptr;
879: int size = old->Length();
880:
881: if (this == 0)
882: this = IRAllocRep ();
883:
884: ptr = new long[(unsigned)size];
885: if (ptr == 0)
886: Fatal ("INT_REP - Out of memory.");
887: beg = ptr;
888: wt = last = ptr + size;
889: refcnt = 0;
890:
891: oldptr = old->beg;
892: while (ptr < last) {
893: *ptr++ = *oldptr++;
894: }
895: }
896:
897: /* INT_REP:
898: * Create new rep,
899: * from ascii representation s.
900: */
901: INT_REP::INT_REP (char *s)
902: {
903: register char *optr;
904: short minus = 0;
905: Mode mode;
906: int size;
907:
908: if (this == 0)
909: this = IRAllocRep ();
910: refcnt = 0;
911:
912: // Handle NULL pointer or NULL string.
913: if ((s == (char *)0) || (*s == '\0')) {
914: MakeZero(this);
915: return;
916: }
917:
918: // Check for sign.
919: if (*s == '-') {
920: s++;
921: minus = 1;
922: }
923:
924: // Determine mode.
925: if (*s == '0') {
926: s++;
927: if (*s == '\0') {
928: MakeZero(this);
929: return;
930: }
931: else if ((*s == 'x') || (*s == 'X')) {
932: s++;
933: mode = hexadecimal;
934: }
935: else
936: mode = octal;
937: }
938: else
939: mode = decimal;
940:
941: // Find end of string.
942: optr = s;
943: while(*optr) optr++;
944: if (s == optr)
945: Fatal ("INT_REP : Improper string.");
946:
947: // Fill representation, using correct mode.
948: if (mode == decimal)
949: FillDec (s, optr);
950: else {
951: // Allocate enough space.
952: size = (strlen(s) + 3)/4;
953: beg = new long[size];
954: if (beg == 0)
955: Fatal ("INT_REP - Out of memory.");
956: last = beg + size;
957:
958: // Fill.
959: if (mode == hexadecimal)
960: FillHex(s, optr);
961: else
962: FillOct(s, optr);
963:
964: // Remove extraneous zeros.
965: StripZero();
966: }
967:
968: // Correct for sign.
969: if (minus)
970: Chsign ();
971: }
972:
973: /* ~INT_REP:
974: * destructor.
975: */
976: INT_REP::~INT_REP ()
977: {
978: if(beg) delete beg;
979:
980: // Put back on free list.
981: wt = (long *)hfree;
982: hfree = this;
983: this = 0;
984: }
985:
986: /* More:
987: * Obtain more storage for this,
988: * using realloc, and resetting
989: * pointers.
990: */
991: void INT_REP::More()
992: {
993: register unsigned size;
994: register long *p;
995:
996: if ((size = (last - beg) * 2) == 0)
997: size = 1;
998: p = (long *)realloc((char *)beg, (unsigned)(size*sizeof(long)));
999: if(p == 0)
1000: Fatal ("More: realloc failed");
1001:
1002: wt = p + (wt - beg);
1003: beg = p;
1004: last = p + size;
1005: }
1006:
1007: /* Zero:
1008: * Zero out all the storage.
1009: */
1010: void INT_REP::Zero()
1011: {
1012: register long *pp;
1013:
1014: for (pp = beg; pp < last;)
1015: *pp++ = 0;
1016: }
1017:
1018:
1019: /* Prt:
1020: * Print internal information;
1021: * for debugging.
1022: */
1023: void INT_REP::Prt()
1024: {
1025: long *ptr = beg;
1026:
1027: while (ptr != wt)
1028: fprintf(stderr, "%x,", *ptr++);
1029: fprintf (stderr, " : this = %x, beg = %x, last = %x, wt = %x\n",
1030: this, beg, last, wt);
1031: }
1032:
1033: /* PutDec:
1034: * Store decimal representation of this in obuf.
1035: * NOTE : The routine destroys this.
1036: */
1037: INT_REP *INT_REP::PutDec ()
1038: {
1039: register INT_REP *divd = this, *div, *mod;
1040: register long val;
1041: int count;
1042:
1043: while (1) {
1044: div = divd->DivMod (&BigTen, 2);
1045: mod = Remainder;
1046: // div->Prt();
1047: // mod->Prt();
1048:
1049: // Store remainder.
1050: val = 0;
1051: switch (mod->Length()) {
1052: case 2 :
1053: val = (mod->beg[1]) << 16;
1054: /* Fall through. */
1055: case 1 :
1056: val |= (mod->beg[0]);
1057: }
1058: // fprintf (stderr, "val = %d (0x%x)\n",val, val);
1059: delete mod;
1060: delete divd;
1061:
1062: count = BlockLen;
1063: while (val) {
1064: obuf.Put((char)(val%10 + '0'));
1065: val /= 10;
1066: count--;
1067: }
1068:
1069: if (div->isZero())
1070: return (div);
1071:
1072: while (count--)
1073: obuf.Put('0');
1074: divd = div;
1075: }
1076:
1077: }
1078:
1079: /* PutHex:
1080: * Store hexadecimal representation of this in obuf.
1081: */
1082: INT_REP *INT_REP::PutHex ()
1083: {
1084: register unsigned long *ptr = (unsigned long *)beg;
1085: register unsigned long val, digit;
1086: register int count;
1087:
1088: while (1) {
1089: val = *ptr++;
1090: count = 4; // 4 hex digits per long
1091:
1092: while (val) {
1093: if ((digit = (val & 0xF)) < 10)
1094: obuf.Put((char)(digit + '0'));
1095: else
1096: obuf.Put((char)(digit + ('A' - 10)));
1097: val >>= 4;
1098: count--;
1099: }
1100:
1101: if (ptr == (unsigned long *) wt)
1102: return (this);
1103:
1104: while (count--)
1105: obuf.Put ('0');
1106: }
1107: }
1108:
1109: /* PutOct:
1110: * Store octal representation of this in obuf.
1111: */
1112: INT_REP *INT_REP::PutOct ()
1113: {
1114: register unsigned long *ptr = (unsigned long *)beg;
1115: register unsigned long val = 0;
1116: register int count;
1117: register int numbits = 0;
1118:
1119: while (1) {
1120: val |= (*ptr++) << numbits;
1121: if (numbits == 2)
1122: count = 6;
1123: else
1124: count = 5;
1125: numbits += 16;
1126:
1127: while (val && count) {
1128: obuf.Put((char)((val & 0x7) + '0'));
1129: val >>= 3;
1130: numbits -= 3;
1131: count--;
1132: }
1133:
1134: if (ptr == (unsigned long *) wt)
1135: return (this);
1136:
1137: if (count) {
1138: numbits -= (count * 3);
1139: while (count--)
1140: obuf.Put ('0');
1141: }
1142: }
1143: }
1144:
1145: /* Print:
1146: * Print INT_REP.
1147: * We use the output mode supplied by the iostream.
1148: * After taking care of 0 and the minus sign,
1149: * we use the appropriate function to store the
1150: * representation in obuf.
1151: * Finally, we flush obuf onto the iostream.
1152: */
1153: void INT_REP::Print(iostream &i)
1154: {
1155: register INT_REP *divd;
1156: #ifdef HEX
1157: register long *rp;
1158: #endif HEX
1159:
1160: if (Length() == 0) {
1161: i.put('0');
1162: return;
1163: }
1164:
1165: divd = new INT_REP (this);
1166: if (isNegative()) {
1167: divd->Chsign();
1168: i.put('-');
1169: }
1170:
1171:
1172: // Print.
1173: #ifdef PRE6
1174: switch (i.convbase()) {
1175: case 0 :
1176: case 10 :
1177: divd = divd->PutDec ();
1178: break;
1179: case 8 :
1180: divd = divd->PutOct ();
1181: break;
1182: case 16 :
1183: divd = divd->PutHex ();
1184: break;
1185: }
1186: #else
1187: switch (i.flags()) {
1188: default :
1189: case ios::dec :
1190: divd = divd->PutDec ();
1191: break;
1192: case ios::oct :
1193: divd = divd->PutOct ();
1194: break;
1195: case ios::hex :
1196: divd = divd->PutHex ();
1197: break;
1198: }
1199: #endif
1200:
1201: delete divd;
1202: obuf.Flush(i);
1203: #ifdef HEX
1204: rp = wt - 1;
1205: while (rp >= beg)
1206: printf("%04x", *rp--);
1207: printf("\n");
1208: #endif HEX
1209:
1210: }
1211:
1212: /* Print:
1213: * Print INT_REP.
1214: * We use the output function supplied by the user.
1215: * The mode parameter determines the output mode.
1216: * After taking care of 0 and the minus sign,
1217: * we use the appropriate function to store the
1218: * representation in obuf.
1219: * Finally, we flush obuf using the user's function.
1220: */
1221: void INT_REP::Print(OutFn outf, int mode)
1222: {
1223: register INT_REP *divd;
1224:
1225: if (Length() == 0) {
1226: (*outf)('0');
1227: (*outf)('\0');
1228: return;
1229: }
1230:
1231: divd = new INT_REP (this);
1232: if (isNegative()) {
1233: divd->Chsign();
1234: (*outf)('-');
1235: }
1236:
1237: // Print.
1238: switch (mode) {
1239: case 0 :
1240: divd = divd->PutDec ();
1241: break;
1242: case 1 :
1243: divd = divd->PutOct ();
1244: break;
1245: case 2 :
1246: divd = divd->PutHex ();
1247: break;
1248: }
1249:
1250: delete divd;
1251: obuf.Flush(outf);
1252:
1253: }
1254:
1255: /* Store:
1256: * Store representation of this in given buffer,
1257: * null terminated.
1258: * Buflen tells how many characters are in the user's
1259: * buffer; -1 => unspecified.
1260: * Representation is base 10 <=> mode = 0
1261: * base 8 <=> mode = 1
1262: * base 16 <=> mode = 2
1263: * Return number of (non-null) characters.
1264: */
1265: int INT_REP::Store (char *buf, int buflen, int mode)
1266: {
1267: register INT_REP *divd;
1268: int length;
1269:
1270: if (buflen == 0)
1271: return (0);
1272: else if (buflen == 1) {
1273: *buf = '\0';
1274: return (0);
1275: }
1276:
1277: if (isZero()) {
1278: *buf++ = '0';
1279: *buf = '\0';
1280: return (1);
1281: }
1282:
1283: divd = new INT_REP (this);
1284: if (isNegative())
1285: divd->Chsign();
1286:
1287: // Print.
1288: switch (mode) {
1289: case 0 :
1290: divd = divd->PutDec ();
1291: break;
1292: case 1 :
1293: divd = divd->PutOct ();
1294: break;
1295: case 2 :
1296: divd = divd->PutHex ();
1297: break;
1298: }
1299: delete divd;
1300:
1301: // If necessary, add minus sign.
1302: if (isNegative())
1303: obuf.Put('-');
1304:
1305: // Compute number of characters to store.
1306: length = obuf.Length ();
1307: if ((buflen > 0) && (length + 1 > buflen))
1308: length = buflen - 1;
1309:
1310: obuf.Flush(buf, length);
1311:
1312: return (length);
1313:
1314: }
1315:
1316: /* Chsign:
1317: * Change sign of the rep.
1318: */
1319: void INT_REP::Chsign ()
1320: {
1321: register long carry;
1322: register long ct;
1323: register long *readp;
1324:
1325: if (Length () == 0) return;
1326:
1327: carry = 0;
1328: readp = beg;
1329: while (readp != wt) {
1330: ct = BASE - *readp - carry;
1331:
1332: if(ct >= BASE) {
1333: ct -= BASE;
1334: if (ct == 1) {
1335: *readp = 1;
1336: return;
1337: }
1338: else if (carry == 1) {
1339: wt = readp;
1340: return;
1341: }
1342: else
1343: carry = 0;
1344: }
1345: else
1346: carry = 1;
1347:
1348: *readp++ = ct;
1349: }
1350: readp--;
1351: if (*readp == BASE-1)
1352: *readp = -1;
1353: else
1354: Putc(-1);
1355:
1356: }
1357:
1358: /* FillDec:
1359: * Use decimal representation to fill this.
1360: * sp points to first character, ep points
1361: * to terminating NULL. Use standard decimal
1362: * to binary conversion. Assume no
1363: * storage has already been allocated.
1364: */
1365: void INT_REP::FillDec(char *sp, char *ep)
1366: {
1367: char *pp;
1368: INT_REP *u, *newi, *product;
1369: unsigned long val;
1370:
1371: // Create utility INT_REP.
1372: u = new INT_REP(2);
1373:
1374: // Determine leftmost block of characters.
1375: pp = sp + ((int)(ep - sp))%BlockLen;
1376:
1377: // Convert to an INT_REP.
1378: val = MakeLong(sp, pp);
1379: u->FillLong (val);
1380: sp = pp;
1381: newi = new INT_REP (u);
1382:
1383: while (sp < ep) {
1384: product = newi->Mult(&BigTen);
1385: delete newi;
1386: val = MakeLong(sp, sp+BlockLen);
1387: u->FillLong (val);
1388: sp += BlockLen;
1389: newi = product->Add(u);
1390: delete product;
1391: }
1392:
1393: // Copy newi into this.
1394: beg = newi->beg;
1395: last = newi->last;
1396: wt = newi->wt;
1397: newi->beg = 0;
1398:
1399: // Delete temporaries.
1400: delete newi;
1401: delete u;
1402:
1403: }
1404:
1405: /* FillHex:
1406: * Use hex representation to fill this.
1407: * sp points to first character, ep points
1408: * to terminating NULL. Assume adequate
1409: * storage has already been allocated.
1410: */
1411: void INT_REP::FillHex(char *sp, char *ep)
1412: {
1413: long *lptr = beg;
1414: long word = 0;
1415: long c;
1416: short shift = 0;
1417:
1418: while (ep > sp) {
1419: c = *(--ep);
1420: if(!isxdigit(c))
1421: Fatal ("FillHex - Illegal hex character.");
1422: if (isdigit(c))
1423: c -= '0';
1424: else if (islower(c))
1425: c -= 'a' - 10;
1426: else
1427: c -= 'A' - 10;
1428:
1429: word |= c << shift;
1430: shift += 4;
1431: if (shift == 16) {
1432: shift = 0;
1433: *lptr++ = word;
1434: word = 0;
1435: }
1436: }
1437:
1438: // Add incomplete word.
1439: if (shift)
1440: *lptr++ = word;
1441:
1442: // Set wt.
1443: wt = lptr;
1444: }
1445:
1446: /* FillOct:
1447: * Use octal representation to fill this.
1448: * sp points to first character, ep points
1449: * to terminating NULL. Assume adequate
1450: * storage has already been allocated.
1451: */
1452: void INT_REP::FillOct(char *sp, char *ep)
1453: {
1454: long *lptr = beg;
1455: long word = 0;
1456: long c;
1457: short shift = 0;
1458:
1459: while (ep > sp) {
1460: c = *(--ep);
1461: if(!isdigit(c) || (c > '7'))
1462: Fatal ("FillOct - Illegal octal character.");
1463:
1464: c -= '0';
1465: word |= c << shift;
1466: shift += 3;
1467: if (shift >= 16) {
1468: *lptr++ = word & 0xFFFF;
1469: shift -= 16;
1470: word = c >> (3 - shift);
1471: }
1472: }
1473:
1474: // Add incomplete word.
1475: if (shift)
1476: *lptr++ = word;
1477:
1478: // Set wt.
1479: wt = lptr;
1480:
1481: }
1482:
1483: /* And:
1484: * Bit-wise and.
1485: */
1486: INT_REP *INT_REP::And(INT_REP *b)
1487: {
1488: INT_REP *newi;
1489: int length;
1490: long *ap, *bp, *pp;
1491:
1492: if (isZero() || b->isZero())
1493: return (&ZeroRep);
1494:
1495: length = Min(Length(), b->Length());
1496: newi = new INT_REP(length);
1497:
1498: ap = beg;
1499: bp = b->beg;
1500: pp = newi->beg;
1501:
1502: while (length-- > 0)
1503: *pp++ = (*ap++) & (*bp++);
1504: newi->wt = pp;
1505:
1506: // Strip trailing zeros.
1507: newi->StripZero();
1508:
1509: return (newi);
1510: }
1511:
1512: /* Or:
1513: * Bit-wise or.
1514: */
1515: INT_REP *INT_REP::Or(INT_REP *b)
1516: {
1517: INT_REP *newi;
1518: int length;
1519: long *ap, *bp, *pp;
1520:
1521: if (isZero()) {
1522: newi = new INT_REP(b);
1523: return newi;
1524: }
1525: if (b->isZero()) {
1526: newi = new INT_REP(this);
1527: return newi;
1528: }
1529:
1530: if (isNegative()) {
1531: if(b->isNegative())
1532: length = Min(Length(), b->Length());
1533: else
1534: length = Length();
1535: }
1536: else {
1537: if(b->isNegative())
1538: length = b->Length();
1539: else
1540: length = Max(Length(), b->Length());
1541: }
1542: newi = new INT_REP(length);
1543:
1544: ap = beg;
1545: bp = b->beg;
1546: pp = newi->beg;
1547:
1548: while (length-- > 0)
1549: *pp++ = (ap < wt ? *ap++ : 0) & (bp < b->wt ? *bp++ : 0);
1550: newi->wt = pp;
1551:
1552: // Strip possible trailing -1's.
1553: newi->StripMinus();
1554:
1555: return (newi);
1556:
1557: }
1558:
1559: /* Xor:
1560: * Bit-wise xor.
1561: */
1562: INT_REP *INT_REP::Xor(INT_REP *b)
1563: {
1564: INT_REP *newi;
1565: INT_REP *pos, *neg;
1566: long *ap, *bp, *pp;
1567: long extval;
1568: int length;
1569:
1570: if (isZero()) {
1571: newi = new INT_REP(b);
1572: return newi;
1573: }
1574: if (b->isZero()) {
1575: newi = new INT_REP(this);
1576: return newi;
1577: }
1578:
1579: if(isNegative() == b->isNegative()) { // Same sign.
1580: length = Max(Length(), b->Length());
1581: newi = new INT_REP(length);
1582:
1583: ap = beg;
1584: bp = b->beg;
1585: pp = newi->beg;
1586: if(isNegative())
1587: extval = -1;
1588: else
1589: extval = 0;
1590:
1591: while (length-- > 0)
1592: *pp++ = ((ap < wt ? *ap++ : extval) ^ (bp < b->wt ? *bp++ : extval)) & 0xFFFF;
1593: newi->wt = pp;
1594:
1595: // Remove trailing zeros.
1596: newi->StripZero();
1597: }
1598: else { // Different sign.
1599: // Set neg pointing to negative value, pos pointing to positive.
1600: if (isNegative()) {
1601: neg = this;
1602: pos = b;
1603: }
1604: else {
1605: neg = b;
1606: pos = this;
1607: }
1608:
1609: length = Max(pos->Length()+1, neg->Length());
1610: newi = new INT_REP(length);
1611:
1612: ap = pos->beg;
1613: bp = neg->beg;
1614: pp = newi->beg;
1615:
1616: while (--length > 0)
1617: *pp++ = ((ap < pos->wt ? *ap++ : 0) ^ (bp < neg->wt ? *bp++ : -1)) & 0xFFFF;
1618: *pp++ = -1;
1619: newi->wt = pp;
1620:
1621: // Remove trailing -1's.
1622: newi->StripMinus();
1623: }
1624:
1625: return (newi);
1626: }
1627:
1628: /* Not:
1629: * Bit-wise not.
1630: */
1631: INT_REP *INT_REP::Not()
1632: {
1633: int length;
1634: long *ap, *pp;
1635: INT_REP *newi;
1636:
1637: if (isZero()) {
1638: newi = new INT_REP(1);
1639: newi->Putc(-1);
1640: return (newi);
1641: }
1642:
1643: if (isNegative()) {
1644: length = Length() - 1;
1645: newi = new INT_REP(length);
1646:
1647: ap = beg;
1648: pp = newi->beg;
1649: while(length--)
1650: *pp++ = (~(*ap++)) & 0xFFFF;
1651: newi->wt = pp;
1652: }
1653: else {
1654: length = Length();
1655: newi = new INT_REP(length+1);
1656:
1657: ap = beg;
1658: pp = newi->beg;
1659: while(length--)
1660: *pp++ = (~(*ap++)) & 0xFFFF;
1661: *pp++ = -1;
1662: newi->wt = pp;
1663: }
1664:
1665: return (newi);
1666: }
1667:
1668: /* LeftShift:
1669: * Bit-wise left shift.
1670: */
1671: INT_REP *INT_REP::LeftShift(INT_REP *b)
1672: {
1673: long shift;
1674: INT_REP *bb, *newi;
1675: int error;
1676:
1677: shift = b->Trunc(error);
1678: if (!error) // Shift can be represented as a long.
1679: return(LeftShift(shift));
1680:
1681: if (isZero()) {
1682: newi = new INT_REP(this);
1683: return newi;
1684: }
1685: if (b->isNegative()) {
1686: bb = new INT_REP(b);
1687: bb->Chsign();
1688: newi = RightShift(bb);
1689: delete bb;
1690: return (newi);
1691: }
1692:
1693: bb = TwoRep.Exp(b);
1694: newi = Mult(bb);
1695: delete bb;
1696: return (newi);
1697:
1698: }
1699:
1700: INT_REP *INT_REP::LeftShift(long shift)
1701: {
1702: INT_REP *newi;
1703: int length;
1704: int bitshift;
1705: long val;
1706: long *ap, *pp;
1707:
1708: // Special cases.
1709: if ((shift == 0) || isZero()) {
1710: newi = new INT_REP(this);
1711: return newi;
1712: }
1713: if (shift < 0) {
1714: return (RightShift(-1*shift));
1715: }
1716:
1717: length = Length() + ((shift + 15)>>4); // (shift + 15)/16.
1718: bitshift = shift & 0xF; // shift%16.
1719: newi = new INT_REP(length);
1720: newi->wt = newi->last;
1721:
1722: // Initialize.
1723: ap = wt-1;
1724: pp = newi->last - 1;
1725:
1726: // Shift.
1727: if(bitshift) {
1728: if (isNegative())
1729: *pp = -1;
1730: else
1731: *pp = 0;
1732:
1733: while (ap >= beg) {
1734: val = *(ap--) & 0xFFFF;
1735: *pp-- |= val >> (16 - bitshift);
1736: *pp = (val << bitshift) & 0xFFFF;
1737: }
1738: }
1739: else {
1740: while (ap >= beg) {
1741: *pp-- = *ap--;
1742: }
1743: pp++;
1744: }
1745:
1746: // Fill with 0's.
1747: while (pp > newi->beg)
1748: *(--pp) = 0;
1749:
1750:
1751: if (isNegative())
1752: newi->StripMinus();
1753: else
1754: newi->StripZero();
1755:
1756: return (newi);
1757: }
1758:
1759: /* RightShift:
1760: * Bit-wise right shift.
1761: */
1762: INT_REP *INT_REP::RightShift(INT_REP *b)
1763: {
1764: long shift;
1765: INT_REP *bb, *newi;
1766: int error;
1767:
1768: shift = b->Trunc(error);
1769: if (!error) // Shift can be represented as a long.
1770: return(RightShift(shift));
1771:
1772: if (isZero()) {
1773: newi = new INT_REP(this);
1774: return newi;
1775: }
1776: if (b->isNegative()) {
1777: bb = new INT_REP(b);
1778: bb->Chsign();
1779: newi = LeftShift(bb);
1780: delete bb;
1781: return (newi);
1782: }
1783:
1784: bb = TwoRep.Exp(b);
1785: newi = Div(bb);
1786: delete bb;
1787: return (newi);
1788:
1789: }
1790:
1791: INT_REP *INT_REP::RightShift(long shift)
1792: {
1793: INT_REP *newi;
1794: int length;
1795: int offset, bitshift;
1796: long val;
1797: long *ap, *pp;
1798:
1799: // Special cases.
1800: if ((shift == 0) || isZero()) {
1801: newi = new INT_REP(this);
1802: return newi;
1803: }
1804: if (shift < 0) {
1805: return (LeftShift(-1*shift));
1806: }
1807:
1808: offset = shift>>4; // shift/16.
1809: // If offset >= Length(), shifting produces either 0 or -1.
1810: if (offset >= Length()) {
1811: if (isNegative()) {
1812: newi = new INT_REP(1);
1813: newi->Putc(-1);
1814: return newi;
1815: }
1816: else {
1817: newi = new INT_REP(0);
1818: return newi;
1819: }
1820: }
1821:
1822: length = Length() - offset;
1823: bitshift = shift & 0xF; // shift%16.
1824: newi = new INT_REP(length);
1825: newi->wt = newi->last;
1826:
1827: // Initialize.
1828: ap = beg + offset;
1829: pp = newi->beg;
1830:
1831: // Shift.
1832: if (bitshift) {
1833: *pp = (*ap++ >> bitshift) & 0xFFFF;
1834:
1835: while (ap < wt) {
1836: val = *ap++ & 0xFFFF;
1837: *pp++ |= (val << (16 - bitshift)) & 0xFFFF;
1838: *pp = (val >> bitshift);
1839: }
1840:
1841: if (isNegative()) {
1842: *pp = -1;
1843: newi->StripMinus();
1844: }
1845: else
1846: newi->StripZero();
1847: }
1848: else {
1849: while (ap < wt)
1850: *pp++ = *ap++;
1851: }
1852:
1853: return (newi);
1854: }
1855:
1856: /* GCD:
1857: * Returns greatest common divisor.
1858: */
1859: INT_REP *INT_REP::GCD(INT_REP *b)
1860: {
1861: INT_REP *newi;
1862: INT_REP *t, *u, *v, *tmp;
1863: int k, l;
1864:
1865: // If this or b is 0.
1866: if (isZero()) {
1867: newi = new INT_REP (b);
1868: if (b->isNegative())
1869: newi->Chsign();
1870: return (newi);
1871: }
1872: if (b->isZero()) {
1873: newi = new INT_REP (this);
1874: if (this->isNegative())
1875: newi->Chsign();
1876: return (newi);
1877: }
1878:
1879: // Make this and b positive.
1880: if (isNegative()) {
1881: u = new INT_REP (this);
1882: u->Chsign();
1883: }
1884: else
1885: u = this;
1886: if (b->isNegative()) {
1887: v = new INT_REP (b);
1888: v->Chsign();
1889: }
1890: else
1891: v = b;
1892:
1893: // Compute highest power of 2 dividing both this and b;
1894: // divide out this power.
1895: k = Min (Twos(), b->Twos());
1896: tmp = u->RightShift(k);
1897: if (u != this) delete u;
1898: u = tmp;
1899: tmp = v->RightShift(k);
1900: if (v != b) delete v;
1901: v = tmp;
1902:
1903: // Initialize.
1904: if (u->isOdd()) {
1905: t = new INT_REP (v);
1906: t->Chsign();
1907: }
1908: else
1909: t = new INT_REP (u);
1910:
1911: // Loop.
1912: while (t->isNonzero()) {
1913:
1914: // Remove powers of 2 from t.
1915: l = t->Twos();
1916: tmp = t->RightShift (l);
1917: delete t;
1918: t = tmp;
1919:
1920: if (t->isNegative()) {
1921: delete v;
1922: v = new INT_REP (t);
1923: v->Chsign();
1924: }
1925: else {
1926: delete u;
1927: u = new INT_REP (t);
1928: }
1929:
1930: delete t;
1931: t = u->Sub(v);
1932: }
1933:
1934: // Add in power of 2 computed above.
1935: newi = u->LeftShift (k);
1936:
1937: delete u;
1938: delete v;
1939:
1940: return (newi);
1941: }
1942:
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.