Objects/longobject.c

   1 /***********************************************************
   2 Copyright 1991-1995 by Stichting Mathematisch Centrum, Amsterdam,
   3 The Netherlands.
   4
   5                         All Rights Reserved
   6
   7 Permission to use, copy, modify, and distribute this software and its
   8 documentation for any purpose and without fee is hereby granted,
   9 provided that the above copyright notice appear in all copies and that
  10 both that copyright notice and this permission notice appear in
  11 supporting documentation, and that the names of Stichting Mathematisch
  12 Centrum or CWI not be used in advertising or publicity pertaining to
  13 distribution of the software without specific, written prior permission.
  14
  15 STICHTING MATHEMATISCH CENTRUM DISCLAIMS ALL WARRANTIES WITH REGARD TO
  16 THIS SOFTWARE, INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND
  17 FITNESS, IN NO EVENT SHALL STICHTING MATHEMATISCH CENTRUM BE LIABLE
  18 FOR ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
  19 WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
  20 ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT
  21 OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
  22
  23 ******************************************************************/
  24
  25 /* Long (arbitrary precision) integer object implementation */
  26
  27 /* XXX The functional organization of this file is terrible */
  28
  29 #include "allobjects.h"
  30 #include "longintrepr.h"
  31 #include "mymath.h"
  32 #include <assert.h>
  33 #include <ctype.h>
  34
  35 #define ABS(x) ((x) < 0 ? -(x) : (x))
  36
  37 /* Forward */
  38 static longobject *long_normalize PROTO((longobject *));
  39 static longobject *mul1 PROTO((longobject *, wdigit));
  40 static longobject *muladd1 PROTO((longobject *, wdigit, wdigit));
  41 static longobject *divrem1 PROTO((longobject *, wdigit, digit *));
  42 static object *long_format PROTO((object *aa, int base));
  43
  44 static int ticker;      /* XXX Could be shared with ceval? */
  45
  46 #define SIGCHECK(block) \
  47         if (--ticker < 0) { \
  48                 ticker = 100; \
  49                 if (sigcheck()) { block; } \
  50         }
  51
  52 /* Normalize (remove leading zeros from) a long int object.
  53    Doesn't attempt to free the storage--in most cases, due to the nature
  54    of the algorithms used, this could save at most be one word anyway. */
  55
  56 static longobject *
  57 long_normalize(v)
  58         register longobject *v;
  59 {
  60         int j = ABS(v->ob_size);
  61         register int i = j;
  62
  63         while (i > 0 && v->ob_digit[i-1] == 0)
  64                 --i;
  65         if (i != j)
  66                 v->ob_size = (v->ob_size < 0) ? -(i) : i;
  67         return v;
  68 }
  69
  70 /* Allocate a new long int object with size digits.
  71    Return NULL and set exception if we run out of memory. */
  72
  73 longobject *
  74 alloclongobject(size)
  75         int size;
  76 {
  77         return NEWVAROBJ(longobject, &Longtype, size);
  78 }
  79
  80 /* Create a new long int object from a C long int */
  81
  82 object *
  83 newlongobject(ival)
  84         long ival;
  85 {
  86         /* Assume a C long fits in at most 3 'digits' */
  87         /* XXX On 64 bit machines this isn't true!!! */
  88         longobject *v = alloclongobject(3);
  89         if (v != NULL) {
  90                 if (ival < 0) {
  91                         ival = -ival;
  92                         v->ob_size = -(v->ob_size);
  93                 }
  94                 v->ob_digit[0] = ival & MASK;
  95                 v->ob_digit[1] = (ival >> SHIFT) & MASK;
  96                 v->ob_digit[2] = ((unsigned long)ival >> (2*SHIFT)) & MASK;
  97                 v = long_normalize(v);
  98         }
  99         return (object *)v;
 100 }
 101
 102 /* Create a new long int object from a C double */
 103
 104 object *
 105 #ifdef MPW
 106 dnewlongobject(double dval)
 107 #else
 108 dnewlongobject(dval)
 109         double dval;
 110 #endif /* MPW */
 111 {
 112         longobject *v;
 113         double frac;
 114         int i, ndig, expo, neg;
 115         neg = 0;
 116         if (dval < 0.0) {
 117                 neg = 1;
 118                 dval = -dval;
 119         }
 120         frac = frexp(dval, &expo); /* dval = frac*2**expo; 0.0 <= frac < 1.0 */
 121         if (expo <= 0)
 122                 return newlongobject(0L);
 123         ndig = (expo-1) / SHIFT + 1; /* Number of 'digits' in result */
 124         v = alloclongobject(ndig);
 125         if (v == NULL)
 126                 return NULL;
 127         frac = ldexp(frac, (expo-1) % SHIFT + 1);
 128         for (i = ndig; --i >= 0; ) {
 129                 long bits = (long)frac;
 130                 v->ob_digit[i] = bits;
 131                 frac = frac - (double)bits;
 132                 frac = ldexp(frac, SHIFT);
 133         }
 134         if (neg)
 135                 v->ob_size = -(v->ob_size);
 136         return (object *)v;
 137 }
 138
 139 /* Get a C long int from a long int object.
 140    Returns -1 and sets an error condition if overflow occurs. */
 141
 142 long
 143 getlongvalue(vv)
 144         object *vv;
 145 {
 146         register longobject *v;
 147         long x, prev;
 148         int i, sign;
 149
 150         if (vv == NULL || !is_longobject(vv)) {
 151                 err_badcall();
 152                 return -1;
 153         }
 154         v = (longobject *)vv;
 155         i = v->ob_size;
 156         sign = 1;
 157         x = 0;
 158         if (i < 0) {
 159                 sign = -1;
 160                 i = -(i);
 161         }
 162         while (--i >= 0) {
 163                 prev = x;
 164                 x = (x << SHIFT) + v->ob_digit[i];
 165                 if ((x >> SHIFT) != prev) {
 166                         err_setstr(OverflowError,
 167                                 "long int too long to convert");
 168                         return -1;
 169                 }
 170         }
 171         return x * sign;
 172 }
 173
 174 /* Get a C double from a long int object.  No overflow check. */
 175
 176 double
 177 dgetlongvalue(vv)
 178         object *vv;
 179 {
 180         register longobject *v;
 181         double x;
 182         double multiplier = (double) (1L << SHIFT);
 183         int i, sign;
 184
 185         if (vv == NULL || !is_longobject(vv)) {
 186                 err_badcall();
 187                 return -1;
 188         }
 189         v = (longobject *)vv;
 190         i = v->ob_size;
 191         sign = 1;
 192         x = 0.0;
 193         if (i < 0) {
 194                 sign = -1;
 195                 i = -(i);
 196         }
 197         while (--i >= 0) {
 198                 x = x*multiplier + (double)v->ob_digit[i];
 199         }
 200         return x * sign;
 201 }
 202
 203 /* Multiply by a single digit, ignoring the sign. */
 204
 205 static longobject *
 206 mul1(a, n)
 207         longobject *a;
 208         wdigit n;
 209 {
 210         return muladd1(a, n, (digit)0);
 211 }
 212
 213 /* Multiply by a single digit and add a single digit, ignoring the sign. */
 214
 215 static longobject *
 216 muladd1(a, n, extra)
 217         longobject *a;
 218         wdigit n;
 219         wdigit extra;
 220 {
 221         int size_a = ABS(a->ob_size);
 222         longobject *z = alloclongobject(size_a+1);
 223         twodigits carry = extra;
 224         int i;
 225
 226         if (z == NULL)
 227                 return NULL;
 228         for (i = 0; i < size_a; ++i) {
 229                 carry += (twodigits)a->ob_digit[i] * n;
 230                 z->ob_digit[i] = carry & MASK;
 231                 carry >>= SHIFT;
 232         }
 233         z->ob_digit[i] = carry;
 234         return long_normalize(z);
 235 }
 236
 237 /* Divide a long integer by a digit, returning both the quotient
 238    (as function result) and the remainder (through *prem).
 239    The sign of a is ignored; n should not be zero. */
 240
 241 static longobject *
 242 divrem1(a, n, prem)
 243         longobject *a;
 244         wdigit n;
 245         digit *prem;
 246 {
 247         int size = ABS(a->ob_size);
 248         longobject *z;
 249         int i;
 250         twodigits rem = 0;
 251
 252         assert(n > 0 && n <= MASK);
 253         z = alloclongobject(size);
 254         if (z == NULL)
 255                 return NULL;
 256         for (i = size; --i >= 0; ) {
 257                 rem = (rem << SHIFT) + a->ob_digit[i];
 258                 z->ob_digit[i] = rem/n;
 259                 rem %= n;
 260         }
 261         *prem = rem;
 262         return long_normalize(z);
 263 }
 264
 265 /* Convert a long int object to a string, using a given conversion base.
 266    Return a string object.
 267    If base is 8 or 16, add the proper prefix '0' or '0x'.
 268    External linkage: used in bltinmodule.c by hex() and oct(). */
 269
 270 static object *
 271 long_format(aa, base)
 272         object *aa;
 273         int base;
 274 {
 275         register longobject *a = (longobject *)aa;
 276         stringobject *str;
 277         int i;
 278         int size_a = ABS(a->ob_size);
 279         char *p;
 280         int bits;
 281         char sign = '\0';
 282
 283         if (a == NULL || !is_longobject(a)) {
 284                 err_badcall();
 285                 return NULL;
 286         }
 287         assert(base >= 2 && base <= 36);
 288
 289         /* Compute a rough upper bound for the length of the string */
 290         i = base;
 291         bits = 0;
 292         while (i > 1) {
 293                 ++bits;
 294                 i >>= 1;
 295         }
 296         i = 6 + (size_a*SHIFT + bits-1) / bits;
 297         str = (stringobject *) newsizedstringobject((char *)0, i);
 298         if (str == NULL)
 299                 return NULL;
 300         p = GETSTRINGVALUE(str) + i;
 301         *p = '\0';
 302         *--p = 'L';
 303         if (a->ob_size < 0)
 304                 sign = '-';
 305
 306         INCREF(a);
 307         do {
 308                 digit rem;
 309                 longobject *temp = divrem1(a, (digit)base, &rem);
 310                 if (temp == NULL) {
 311                         DECREF(a);
 312                         DECREF(str);
 313                         return NULL;
 314                 }
 315                 if (rem < 10)
 316                         rem += '0';
 317                 else
 318                         rem += 'A'-10;
 319                 assert(p > GETSTRINGVALUE(str));
 320                 *--p = rem;
 321                 DECREF(a);
 322                 a = temp;
 323                 SIGCHECK({
 324                         DECREF(a);
 325                         DECREF(str);
 326                         return NULL;
 327                 })
 328         } while (ABS(a->ob_size) != 0);
 329         DECREF(a);
 330         if (base == 8) {
 331                 if (size_a != 0)
 332                         *--p = '0';
 333         }
 334         else if (base == 16) {
 335                 *--p = 'x';
 336                 *--p = '0';
 337         }
 338         else if (base != 10) {
 339                 *--p = '#';
 340                 *--p = '0' + base%10;
 341                 if (base > 10)
 342                         *--p = '0' + base/10;
 343         }
 344         if (sign)
 345                 *--p = sign;
 346         if (p != GETSTRINGVALUE(str)) {
 347                 char *q = GETSTRINGVALUE(str);
 348                 assert(p > q);
 349                 do {
 350                 } while ((*q++ = *p++) != '\0');
 351                 q--;
 352                 resizestring((object **)&str, (int) (q - GETSTRINGVALUE(str)));
 353         }
 354         return (object *)str;
 355 }
 356
 357 #if 0
 358 /* Convert a string to a long int object, in a given base.
 359    Base zero implies a default depending on the number.
 360    External linkage: used in compile.c and stropmodule.c. */
 361
 362 object *
 363 long_scan(str, base)
 364         char *str;
 365         int base;
 366 {
 367         return long_escan(str, (char **)NULL, base);
 368 }
 369 #endif
 370
 371 object *
 372 long_escan(str, pend, base)
 373         char *str;
 374         char **pend;
 375         int base;
 376 {
 377         int sign = 1;
 378         longobject *z;
 379
 380         if (base != 0 && base < 2 || base > 36) {
 381                 err_setstr(ValueError, "invalid base for long literal");
 382                 return NULL;
 383         }
 384         while (*str != '\0' && isspace(Py_CHARMASK(*str)))
 385                 str++;
 386         if (*str == '+')
 387                 ++str;
 388         else if (*str == '-') {
 389                 ++str;
 390                 sign = -1;
 391         }
 392         while (*str != '\0' && isspace(Py_CHARMASK(*str)))
 393                 str++;
 394         if (base == 0) {
 395                 if (str[0] != '0')
 396                         base = 10;
 397                 else if (str[1] == 'x' || str[1] == 'X')
 398                         base = 16;
 399                 else
 400                         base = 8;
 401         }
 402         if (base == 16 && str[0] == '0' && (str[1] == 'x' || str[1] == 'X'))
 403                 str += 2;
 404         z = alloclongobject(0);
 405         for ( ; z != NULL; ++str) {
 406                 int k = -1;
 407                 longobject *temp;
 408
 409                 if (*str <= '9')
 410                         k = *str - '0';
 411                 else if (*str >= 'a')
 412                         k = *str - 'a' + 10;
 413                 else if (*str >= 'A')
 414                         k = *str - 'A' + 10;
 415                 if (k < 0 || k >= base)
 416                         break;
 417                 temp = muladd1(z, (digit)base, (digit)k);
 418                 DECREF(z);
 419                 z = temp;
 420         }
 421         if (sign < 0 && z != NULL && z->ob_size != 0)
 422                 z->ob_size = -(z->ob_size);
 423         if (pend)
 424                 *pend = str;
 425         return (object *) z;
 426 }
 427
 428 static longobject *x_divrem PROTO((longobject *, longobject *, longobject **));
 429 static object *long_pos PROTO((longobject *));
 430 static long_divrem PROTO((longobject *, longobject *,
 431         longobject **, longobject **));
 432
 433 /* Long division with remainder, top-level routine */
 434
 435 static int
 436 long_divrem(a, b, pdiv, prem)
 437         longobject *a, *b;
 438         longobject **pdiv;
 439         longobject **prem;
 440 {
 441         int size_a = ABS(a->ob_size), size_b = ABS(b->ob_size);
 442         longobject *z;
 443
 444         if (size_b == 0) {
 445                 err_setstr(ZeroDivisionError, "long division or modulo");
 446                 return -1;
 447         }
 448         if (size_a < size_b ||
 449                         size_a == size_b &&
 450                         a->ob_digit[size_a-1] < b->ob_digit[size_b-1]) {
 451                 /* |a| < |b|. */
 452                 *pdiv = alloclongobject(0);
 453                 INCREF(a);
 454                 *prem = (longobject *) a;
 455                 return 0;
 456         }
 457         if (size_b == 1) {
 458                 digit rem = 0;
 459                 z = divrem1(a, b->ob_digit[0], &rem);
 460                 if (z == NULL)
 461                         return -1;
 462                 *prem = (longobject *) newlongobject((long)rem);
 463         }
 464         else {
 465                 z = x_divrem(a, b, prem);
 466                 if (z == NULL)
 467                         return -1;
 468         }
 469         /* Set the signs.
 470            The quotient z has the sign of a*b;
 471            the remainder r has the sign of a,
 472            so a = b*z + r. */
 473         if ((a->ob_size < 0) != (b->ob_size < 0))
 474                 z->ob_size = -(z->ob_size);
 475         if (a->ob_size < 0 && (*prem)->ob_size != 0)
 476                 (*prem)->ob_size = -((*prem)->ob_size);
 477         *pdiv = z;
 478         return 0;
 479 }
 480
 481 /* Unsigned long division with remainder -- the algorithm */
 482
 483 static longobject *
 484 x_divrem(v1, w1, prem)
 485         longobject *v1, *w1;
 486         longobject **prem;
 487 {
 488         int size_v = ABS(v1->ob_size), size_w = ABS(w1->ob_size);
 489         digit d = (twodigits)BASE / (w1->ob_digit[size_w-1] + 1);
 490         longobject *v = mul1(v1, d);
 491         longobject *w = mul1(w1, d);
 492         longobject *a;
 493         int j, k;
 494
 495         if (v == NULL || w == NULL) {
 496                 XDECREF(v);
 497                 XDECREF(w);
 498                 return NULL;
 499         }
 500
 501         assert(size_v >= size_w && size_w > 1); /* Assert checks by div() */
 502         assert(v->ob_refcnt == 1); /* Since v will be used as accumulator! */
 503         assert(size_w == ABS(w->ob_size)); /* That's how d was calculated */
 504
 505         size_v = ABS(v->ob_size);
 506         a = alloclongobject(size_v - size_w + 1);
 507
 508         for (j = size_v, k = a->ob_size-1; a != NULL && k >= 0; --j, --k) {
 509                 digit vj = (j >= size_v) ? 0 : v->ob_digit[j];
 510                 twodigits q;
 511                 stwodigits carry = 0;
 512                 int i;
 513
 514                 SIGCHECK({
 515                         DECREF(a);
 516                         a = NULL;
 517                         break;
 518                 })
 519                 if (vj == w->ob_digit[size_w-1])
 520                         q = MASK;
 521                 else
 522                         q = (((twodigits)vj << SHIFT) + v->ob_digit[j-1]) /
 523                                 w->ob_digit[size_w-1];
 524
 525                 while (w->ob_digit[size_w-2]*q >
 526                                 ((
 527                                         ((twodigits)vj << SHIFT)
 528                                         + v->ob_digit[j-1]
 529                                         - q*w->ob_digit[size_w-1]
 530                                                                 ) << SHIFT)
 531                                 + v->ob_digit[j-2])
 532                         --q;
 533
 534                 for (i = 0; i < size_w && i+k < size_v; ++i) {
 535                         twodigits z = w->ob_digit[i] * q;
 536                         digit zz = z >> SHIFT;
 537                         carry += v->ob_digit[i+k] - z + ((twodigits)zz << SHIFT);
 538                         v->ob_digit[i+k] = carry & MASK;
 539                         carry = (carry >> SHIFT) - zz;
 540                 }
 541
 542                 if (i+k < size_v) {
 543                         carry += v->ob_digit[i+k];
 544                         v->ob_digit[i+k] = 0;
 545                 }
 546
 547                 if (carry == 0)
 548                         a->ob_digit[k] = q;
 549                 else {
 550                         assert(carry == -1);
 551                         a->ob_digit[k] = q-1;
 552                         carry = 0;
 553                         for (i = 0; i < size_w && i+k < size_v; ++i) {
 554                                 carry += v->ob_digit[i+k] + w->ob_digit[i];
 555                                 v->ob_digit[i+k] = carry & MASK;
 556                                 carry >>= SHIFT;
 557                         }
 558                 }
 559         } /* for j, k */
 560
 561         if (a == NULL)
 562                 *prem = NULL;
 563         else {
 564                 a = long_normalize(a);
 565                 *prem = divrem1(v, d, &d);
 566                 /* d receives the (unused) remainder */
 567                 if (*prem == NULL) {
 568                         DECREF(a);
 569                         a = NULL;
 570                 }
 571         }
 572         DECREF(v);
 573         DECREF(w);
 574         return a;
 575 }
 576
 577 /* Methods */
 578
 579 /* Forward */
 580 static void long_dealloc PROTO((object *));
 581 static object *long_repr PROTO((object *));
 582 static int long_compare PROTO((longobject *, longobject *));
 583 static long long_hash PROTO((longobject *));
 584
 585 static object *long_add PROTO((longobject *, longobject *));
 586 static object *long_sub PROTO((longobject *, longobject *));
 587 static object *long_mul PROTO((longobject *, longobject *));
 588 static object *long_div PROTO((longobject *, longobject *));
 589 static object *long_mod PROTO((longobject *, longobject *));
 590 static object *long_divmod PROTO((longobject *, longobject *));
 591 static object *long_pow PROTO((longobject *, longobject *, longobject *));
 592 static object *long_neg PROTO((longobject *));
 593 static object *long_pos PROTO((longobject *));
 594 static object *long_abs PROTO((longobject *));
 595 static int long_nonzero PROTO((longobject *));
 596 static object *long_invert PROTO((longobject *));
 597 static object *long_lshift PROTO((longobject *, longobject *));
 598 static object *long_rshift PROTO((longobject *, longobject *));
 599 static object *long_and PROTO((longobject *, longobject *));
 600 static object *long_xor PROTO((longobject *, longobject *));
 601 static object *long_or PROTO((longobject *, longobject *));
 602
 603 static void
 604 long_dealloc(v)
 605         object *v;
 606 {
 607         DEL(v);
 608 }
 609
 610 static object *
 611 long_repr(v)
 612         object *v;
 613 {
 614         return long_format(v, 10);
 615 }
 616
 617 static int
 618 long_compare(a, b)
 619         longobject *a, *b;
 620 {
 621         int sign;
 622
 623         if (a->ob_size != b->ob_size) {
 624                 if (ABS(a->ob_size) == 0 && ABS(b->ob_size) == 0)
 625                         sign = 0;
 626                 else
 627                         sign = a->ob_size - b->ob_size;
 628         }
 629         else {
 630                 int i = ABS(a->ob_size);
 631                 while (--i >= 0 && a->ob_digit[i] == b->ob_digit[i])
 632                         ;
 633                 if (i < 0)
 634                         sign = 0;
 635                 else {
 636                         sign = (int)a->ob_digit[i] - (int)b->ob_digit[i];
 637                         if (a->ob_size < 0)
 638                                 sign = -sign;
 639                 }
 640         }
 641         return sign < 0 ? -1 : sign > 0 ? 1 : 0;
 642 }
 643
 644 static long
 645 long_hash(v)
 646         longobject *v;
 647 {
 648         long x;
 649         int i, sign;
 650
 651         /* This is designed so that Python ints and longs with the
 652            same value hash to the same value, otherwise comparisons
 653            of mapping keys will turn out weird */
 654         i = v->ob_size;
 655         sign = 1;
 656         x = 0;
 657         if (i < 0) {
 658                 sign = -1;
 659                 i = -(i);
 660         }
 661         while (--i >= 0) {
 662                 /* Force a 32-bit circular shift */
 663                 x = ((x << SHIFT) & ~MASK) | ((x >> (32-SHIFT)) & MASK);
 664                 x += v->ob_digit[i];
 665         }
 666         x = x * sign;
 667         if (x == -1)
 668                 x = -2;
 669         return x;
 670 }
 671
 672
 673 /* Add the absolute values of two long integers. */
 674
 675 static longobject *x_add PROTO((longobject *, longobject *));
 676 static longobject *
 677 x_add(a, b)
 678         longobject *a, *b;
 679 {
 680         int size_a = ABS(a->ob_size), size_b = ABS(b->ob_size);
 681         longobject *z;
 682         int i;
 683         digit carry = 0;
 684
 685         /* Ensure a is the larger of the two: */
 686         if (size_a < size_b) {
 687                 { longobject *temp = a; a = b; b = temp; }
 688                 { int size_temp = size_a; size_a = size_b; size_b = size_temp; }
 689         }
 690         z = alloclongobject(size_a+1);
 691         if (z == NULL)
 692                 return NULL;
 693         for (i = 0; i < size_b; ++i) {
 694                 carry += a->ob_digit[i] + b->ob_digit[i];
 695                 z->ob_digit[i] = carry & MASK;
 696                 /* The following assumes unsigned shifts don't
 697                    propagate the sign bit. */
 698                 carry >>= SHIFT;
 699         }
 700         for (; i < size_a; ++i) {
 701                 carry += a->ob_digit[i];
 702                 z->ob_digit[i] = carry & MASK;
 703                 carry >>= SHIFT;
 704         }
 705         z->ob_digit[i] = carry;
 706         return long_normalize(z);
 707 }
 708
 709 /* Subtract the absolute values of two integers. */
 710
 711 static longobject *x_sub PROTO((longobject *, longobject *));
 712 static longobject *
 713 x_sub(a, b)
 714         longobject *a, *b;
 715 {
 716         int size_a = ABS(a->ob_size), size_b = ABS(b->ob_size);
 717         longobject *z;
 718         int i;
 719         int sign = 1;
 720         digit borrow = 0;
 721
 722         /* Ensure a is the larger of the two: */
 723         if (size_a < size_b) {
 724                 sign = -1;
 725                 { longobject *temp = a; a = b; b = temp; }
 726                 { int size_temp = size_a; size_a = size_b; size_b = size_temp; }
 727         }
 728         else if (size_a == size_b) {
 729                 /* Find highest digit where a and b differ: */
 730                 i = size_a;
 731                 while (--i >= 0 && a->ob_digit[i] == b->ob_digit[i])
 732                         ;
 733                 if (i < 0)
 734                         return alloclongobject(0);
 735                 if (a->ob_digit[i] < b->ob_digit[i]) {
 736                         sign = -1;
 737                         { longobject *temp = a; a = b; b = temp; }
 738                 }
 739                 size_a = size_b = i+1;
 740         }
 741         z = alloclongobject(size_a);
 742         if (z == NULL)
 743                 return NULL;
 744         for (i = 0; i < size_b; ++i) {
 745                 /* The following assumes unsigned arithmetic
 746                    works module 2**N for some N>SHIFT. */
 747                 borrow = a->ob_digit[i] - b->ob_digit[i] - borrow;
 748                 z->ob_digit[i] = borrow & MASK;
 749                 borrow >>= SHIFT;
 750                 borrow &= 1; /* Keep only one sign bit */
 751         }
 752         for (; i < size_a; ++i) {
 753                 borrow = a->ob_digit[i] - borrow;
 754                 z->ob_digit[i] = borrow & MASK;
 755                 borrow >>= SHIFT;
 756         }
 757         assert(borrow == 0);
 758         if (sign < 0)
 759                 z->ob_size = -(z->ob_size);
 760         return long_normalize(z);
 761 }
 762
 763 static object *
 764 long_add(a, b)
 765         longobject *a;
 766         longobject *b;
 767 {
 768         longobject *z;
 769
 770         if (a->ob_size < 0) {
 771                 if (b->ob_size < 0) {
 772                         z = x_add(a, b);
 773                         if (z != NULL && z->ob_size != 0)
 774                                 z->ob_size = -(z->ob_size);
 775                 }
 776                 else
 777                         z = x_sub(b, a);
 778         }
 779         else {
 780                 if (b->ob_size < 0)
 781                         z = x_sub(a, b);
 782                 else
 783                         z = x_add(a, b);
 784         }
 785         return (object *)z;
 786 }
 787
 788 static object *
 789 long_sub(a, b)
 790         longobject *a;
 791         longobject *b;
 792 {
 793         longobject *z;
 794
 795         if (a->ob_size < 0) {
 796                 if (b->ob_size < 0)
 797                         z = x_sub(a, b);
 798                 else
 799                         z = x_add(a, b);
 800                 if (z != NULL && z->ob_size != 0)
 801                         z->ob_size = -(z->ob_size);
 802         }
 803         else {
 804                 if (b->ob_size < 0)
 805                         z = x_add(a, b);
 806                 else
 807                         z = x_sub(a, b);
 808         }
 809         return (object *)z;
 810 }
 811
 812 static object *
 813 long_mul(a, b)
 814         longobject *a;
 815         longobject *b;
 816 {
 817         int size_a;
 818         int size_b;
 819         longobject *z;
 820         int i;
 821
 822         size_a = ABS(a->ob_size);
 823         size_b = ABS(b->ob_size);
 824         z = alloclongobject(size_a + size_b);
 825         if (z == NULL)
 826                 return NULL;
 827         for (i = 0; i < z->ob_size; ++i)
 828                 z->ob_digit[i] = 0;
 829         for (i = 0; i < size_a; ++i) {
 830                 twodigits carry = 0;
 831                 twodigits f = a->ob_digit[i];
 832                 int j;
 833
 834                 SIGCHECK({
 835                         DECREF(z);
 836                         return NULL;
 837                 })
 838                 for (j = 0; j < size_b; ++j) {
 839                         carry += z->ob_digit[i+j] + b->ob_digit[j] * f;
 840                         z->ob_digit[i+j] = carry & MASK;
 841                         carry >>= SHIFT;
 842                 }
 843                 for (; carry != 0; ++j) {
 844                         assert(i+j < z->ob_size);
 845                         carry += z->ob_digit[i+j];
 846                         z->ob_digit[i+j] = carry & MASK;
 847                         carry >>= SHIFT;
 848                 }
 849         }
 850         if (a->ob_size < 0)
 851                 z->ob_size = -(z->ob_size);
 852         if (b->ob_size < 0)
 853                 z->ob_size = -(z->ob_size);
 854         return (object *) long_normalize(z);
 855 }
 856
 857 /* The / and % operators are now defined in terms of divmod().
 858    The expression a mod b has the value a - b*floor(a/b).
 859    The long_divrem function gives the remainder after division of
 860    |a| by |b|, with the sign of a.  This is also expressed
 861    as a - b*trunc(a/b), if trunc truncates towards zero.
 862    Some examples:
 863          a       b      a rem b         a mod b
 864          13      10      3               3
 865         -13      10     -3               7
 866          13     -10      3              -7
 867         -13     -10     -3              -3
 868    So, to get from rem to mod, we have to add b if a and b
 869    have different signs.  We then subtract one from the 'div'
 870    part of the outcome to keep the invariant intact. */
 871
 872 static int l_divmod PROTO((longobject *, longobject *,
 873         longobject **, longobject **));
 874 static int
 875 l_divmod(v, w, pdiv, pmod)
 876         longobject *v;
 877         longobject *w;
 878         longobject **pdiv;
 879         longobject **pmod;
 880 {
 881         longobject *div, *mod;
 882
 883         if (long_divrem(v, w, &div, &mod) < 0)
 884                 return -1;
 885         if (mod->ob_size < 0 && w->ob_size > 0 ||
 886                                 mod->ob_size > 0 && w->ob_size < 0) {
 887                 longobject *temp;
 888                 longobject *one;
 889                 temp = (longobject *) long_add(mod, w);
 890                 DECREF(mod);
 891                 mod = temp;
 892                 if (mod == NULL) {
 893                         DECREF(div);
 894                         return -1;
 895                 }
 896                 one = (longobject *) newlongobject(1L);
 897                 if (one == NULL ||
 898                     (temp = (longobject *) long_sub(div, one)) == NULL) {
 899                         DECREF(mod);
 900                         DECREF(div);
 901                         XDECREF(one);
 902                         return -1;
 903                 }
 904                 DECREF(one);
 905                 DECREF(div);
 906                 div = temp;
 907         }
 908         *pdiv = div;
 909         *pmod = mod;
 910         return 0;
 911 }
 912
 913 static object *
 914 long_div(v, w)
 915         longobject *v;
 916         longobject *w;
 917 {
 918         longobject *div, *mod;
 919         if (l_divmod(v, w, &div, &mod) < 0)
 920                 return NULL;
 921         DECREF(mod);
 922         return (object *)div;
 923 }
 924
 925 static object *
 926 long_mod(v, w)
 927         longobject *v;
 928         longobject *w;
 929 {
 930         longobject *div, *mod;
 931         if (l_divmod(v, w, &div, &mod) < 0)
 932                 return NULL;
 933         DECREF(div);
 934         return (object *)mod;
 935 }
 936
 937 static object *
 938 long_divmod(v, w)
 939         longobject *v;
 940         longobject *w;
 941 {
 942         object *z;
 943         longobject *div, *mod;
 944         if (l_divmod(v, w, &div, &mod) < 0)
 945                 return NULL;
 946         z = newtupleobject(2);
 947         if (z != NULL) {
 948                 settupleitem(z, 0, (object *) div);
 949                 settupleitem(z, 1, (object *) mod);
 950         }
 951         else {
 952                 DECREF(div);
 953                 DECREF(mod);
 954         }
 955         return z;
 956 }
 957
 958 static object *
 959 long_pow(a, b, c)
 960         longobject *a;
 961         longobject *b;
 962         longobject *c;
 963 {
 964         longobject *z, *div, *mod;
 965         int size_b, i;
 966
 967         size_b = b->ob_size;
 968         if (size_b < 0) {
 969                 err_setstr(ValueError, "long integer to the negative power");
 970                 return NULL;
 971         }
 972         z = (longobject *)newlongobject(1L);
 973         INCREF(a);
 974         for (i = 0; i < size_b; ++i) {
 975                 digit bi = b->ob_digit[i];
 976                 int j;
 977
 978                 for (j = 0; j < SHIFT; ++j) {
 979                         longobject *temp;
 980
 981                         if (bi & 1) {
 982                                 temp = (longobject *)long_mul(z, a);
 983                                 DECREF(z);
 984                                 if ((object*)c!=None && temp!=NULL) {
 985                                         l_divmod(temp, c, &div, &mod);
 986                                         XDECREF(div);
 987                                         DECREF(temp);
 988                                         temp = mod;
 989                                 }
 990                                 z = temp;
 991                                 if (z == NULL)
 992                                         break;
 993                         }
 994                         bi >>= 1;
 995                         if (bi == 0 && i+1 == size_b)
 996                                 break;
 997                         temp = (longobject *)long_mul(a, a);
 998                         DECREF(a);
 999                         if ((object*)c!=None && temp!=NULL) {
1000                                 l_divmod(temp, c, &div, &mod);
1001                                 XDECREF(div);
1002                                 DECREF(temp);
1003                                 temp = mod;
1004                         }
1005                         a = temp;
1006                         if (a == NULL) {
1007                                 DECREF(z);
1008                                 z = NULL;
1009                                 break;
1010                         }
1011                 }
1012                 if (a == NULL || z == NULL)
1013                         break;
1014         }
1015         XDECREF(a);
1016         if ((object*)c!=None && z!=NULL) {
1017                         l_divmod(z, c, &div, &mod);
1018                         XDECREF(div);
1019                         DECREF(z);
1020                         z=mod;
1021         }
1022         return (object *)z;
1023 }
1024
1025 static object *
1026 long_invert(v)
1027         longobject *v;
1028 {
1029         /* Implement ~x as -(x+1) */
1030         longobject *x;
1031         longobject *w;
1032         w = (longobject *)newlongobject(1L);
1033         if (w == NULL)
1034                 return NULL;
1035         x = (longobject *) long_add(v, w);
1036         DECREF(w);
1037         if (x == NULL)
1038                 return NULL;
1039         if (x->ob_size != 0)
1040                 x->ob_size = -(x->ob_size);
1041         return (object *)x;
1042 }
1043
1044 static object *
1045 long_pos(v)
1046         longobject *v;
1047 {
1048         INCREF(v);
1049         return (object *)v;
1050 }
1051
1052 static object *
1053 long_neg(v)
1054         longobject *v;
1055 {
1056         longobject *z;
1057         int i, n;
1058         n = ABS(v->ob_size);
1059         if (n == 0) {
1060                 /* -0 == 0 */
1061                 INCREF(v);
1062                 return (object *) v;
1063         }
1064         z = alloclongobject(ABS(n));
1065         if (z == NULL)
1066                 return NULL;
1067         for (i = 0; i < n; i++)
1068                 z->ob_digit[i] = v->ob_digit[i];
1069         z->ob_size = -(v->ob_size);
1070         return (object *)z;
1071 }
1072
1073 static object *
1074 long_abs(v)
1075         longobject *v;
1076 {
1077         if (v->ob_size < 0)
1078                 return long_neg(v);
1079         else {
1080                 INCREF(v);
1081                 return (object *)v;
1082         }
1083 }
1084
1085 static int
1086 long_nonzero(v)
1087         longobject *v;
1088 {
1089         return ABS(v->ob_size) != 0;
1090 }
1091
1092 static object *
1093 long_rshift(a, b)
1094         longobject *a;
1095         longobject *b;
1096 {
1097         longobject *z;
1098         long shiftby;
1099         int newsize, wordshift, loshift, hishift, i, j;
1100         digit lomask, himask;
1101
1102         if (a->ob_size < 0) {
1103                 /* Right shifting negative numbers is harder */
1104                 longobject *a1, *a2, *a3;
1105                 a1 = (longobject *) long_invert(a);
1106                 if (a1 == NULL) return NULL;
1107                 a2 = (longobject *) long_rshift(a1, b);
1108                 DECREF(a1);
1109                 if (a2 == NULL) return NULL;
1110                 a3 = (longobject *) long_invert(a2);
1111                 DECREF(a2);
1112                 return (object *) a3;
1113         }
1114
1115         shiftby = getlongvalue((object *)b);
1116         if (shiftby == -1L && err_occurred())
1117                 return NULL;
1118         if (shiftby < 0) {
1119                 err_setstr(ValueError, "negative shift count");
1120                 return NULL;
1121         }
1122         wordshift = shiftby / SHIFT;
1123         newsize = ABS(a->ob_size) - wordshift;
1124         if (newsize <= 0) {
1125                 z = alloclongobject(0);
1126                 return (object *)z;
1127         }
1128         loshift = shiftby % SHIFT;
1129         hishift = SHIFT - loshift;
1130         lomask = ((digit)1 << hishift) - 1;
1131         himask = MASK ^ lomask;
1132         z = alloclongobject(newsize);
1133         if (z == NULL)
1134                 return NULL;
1135         if (a->ob_size < 0)
1136                 z->ob_size = -(z->ob_size);
1137         for (i = 0, j = wordshift; i < newsize; i++, j++) {
1138                 z->ob_digit[i] = (a->ob_digit[j] >> loshift) & lomask;
1139                 if (i+1 < newsize)
1140                         z->ob_digit[i] |=
1141                           (a->ob_digit[j+1] << hishift) & himask;
1142         }
1143         return (object *) long_normalize(z);
1144 }
1145
1146 static object *
1147 long_lshift(a, b)
1148         longobject *a;
1149         longobject *b;
1150 {
1151         longobject *z;
1152         long shiftby;
1153         int newsize, wordshift, loshift, hishift, i, j;
1154         digit lomask, himask;
1155
1156         shiftby = getlongvalue((object *)b);
1157         if (shiftby == -1L && err_occurred())
1158                 return NULL;
1159         if (shiftby < 0) {
1160                 err_setstr(ValueError, "negative shift count");
1161                 return NULL;
1162         }
1163         if (shiftby > MASK) {
1164                 err_setstr(ValueError, "outrageous left shift count");
1165                 return NULL;
1166         }
1167         if (shiftby % SHIFT == 0) {
1168                 wordshift = shiftby / SHIFT;
1169                 loshift = 0;
1170                 hishift = SHIFT;
1171                 newsize = ABS(a->ob_size) + wordshift;
1172                 lomask = MASK;
1173                 himask = 0;
1174         }
1175         else {
1176                 wordshift = shiftby / SHIFT + 1;
1177                 loshift = SHIFT - shiftby%SHIFT;
1178                 hishift = shiftby % SHIFT;
1179                 newsize = ABS(a->ob_size) + wordshift;
1180                 lomask = ((digit)1 << hishift) - 1;
1181                 himask = MASK ^ lomask;
1182         }
1183         z = alloclongobject(newsize);
1184         if (z == NULL)
1185                 return NULL;
1186         if (a->ob_size < 0)
1187                 z->ob_size = -(z->ob_size);
1188         for (i = 0; i < wordshift; i++)
1189                 z->ob_digit[i] = 0;
1190         for (i = wordshift, j = 0; i < newsize; i++, j++) {
1191                 if (i > 0)
1192                         z->ob_digit[i-1] |=
1193                                 (a->ob_digit[j] << hishift) & himask;
1194                 z->ob_digit[i] =
1195                         (a->ob_digit[j] >> loshift) & lomask;
1196         }
1197         return (object *) long_normalize(z);
1198 }
1199
1200
1201 /* Bitwise and/xor/or operations */
1202
1203 #define MAX(x, y) ((x) < (y) ? (y) : (x))
1204 #define MIN(x, y) ((x) > (y) ? (y) : (x))
1205
1206 static object *long_bitwise PROTO((longobject *, int, longobject *));
1207 static object *
1208 long_bitwise(a, op, b)
1209         longobject *a;
1210         int op; /* '&', '|', '^' */
1211         longobject *b;
1212 {
1213         digit maska, maskb; /* 0 or MASK */
1214         int negz;
1215         int size_a, size_b, size_z;
1216         longobject *z;
1217         int i;
1218         digit diga, digb;
1219         object *v;
1220
1221         if (a->ob_size < 0) {
1222                 a = (longobject *) long_invert(a);
1223                 maska = MASK;
1224         }
1225         else {
1226                 INCREF(a);
1227                 maska = 0;
1228         }
1229         if (b->ob_size < 0) {
1230                 b = (longobject *) long_invert(b);
1231                 maskb = MASK;
1232         }
1233         else {
1234                 INCREF(b);
1235                 maskb = 0;
1236         }
1237
1238         size_a = a->ob_size;
1239         size_b = b->ob_size;
1240         size_z = MAX(size_a, size_b);
1241         z = alloclongobject(size_z);
1242         if (a == NULL || b == NULL || z == NULL) {
1243                 XDECREF(a);
1244                 XDECREF(b);
1245                 XDECREF(z);
1246                 return NULL;
1247         }
1248
1249         negz = 0;
1250         switch (op) {
1251         case '^':
1252                 if (maska != maskb) {
1253                         maska ^= MASK;
1254                         negz = -1;
1255                 }
1256                 break;
1257         case '&':
1258                 if (maska && maskb) {
1259                         op = '|';
1260                         maska ^= MASK;
1261                         maskb ^= MASK;
1262                         negz = -1;
1263                 }
1264                 break;
1265         case '|':
1266                 if (maska || maskb) {
1267                         op = '&';
1268                         maska ^= MASK;
1269                         maskb ^= MASK;
1270                         negz = -1;
1271                 }
1272                 break;
1273         }
1274
1275         for (i = 0; i < size_z; ++i) {
1276                 diga = (i < size_a ? a->ob_digit[i] : 0) ^ maska;
1277                 digb = (i < size_b ? b->ob_digit[i] : 0) ^ maskb;
1278                 switch (op) {
1279                 case '&': z->ob_digit[i] = diga & digb; break;
1280                 case '|': z->ob_digit[i] = diga | digb; break;
1281                 case '^': z->ob_digit[i] = diga ^ digb; break;
1282                 }
1283         }
1284
1285         DECREF(a);
1286         DECREF(b);
1287         z = long_normalize(z);
1288         if (negz == 0)
1289                 return (object *) z;
1290         v = long_invert(z);
1291         DECREF(z);
1292         return v;
1293 }
1294
1295 static object *
1296 long_and(a, b)
1297         longobject *a;
1298         longobject *b;
1299 {
1300         return long_bitwise(a, '&', b);
1301 }
1302
1303 static object *
1304 long_xor(a, b)
1305         longobject *a;
1306         longobject *b;
1307 {
1308         return long_bitwise(a, '^', b);
1309 }
1310
1311 static object *
1312 long_or(a, b)
1313         longobject *a;
1314         longobject *b;
1315 {
1316         return long_bitwise(a, '|', b);
1317 }
1318
1319 static int
1320 long_coerce(pv, pw)
1321         object **pv;
1322         object **pw;
1323 {
1324         if (is_intobject(*pw)) {
1325                 *pw = newlongobject(getintvalue(*pw));
1326                 INCREF(*pv);
1327                 return 0;
1328         }
1329         return 1; /* Can't do it */
1330 }
1331
1332 static object *
1333 long_int(v)
1334         object *v;
1335 {
1336         long x;
1337         x = getlongvalue(v);
1338         if (err_occurred())
1339                 return NULL;
1340         return newintobject(x);
1341 }
1342
1343 static object *
1344 long_long(v)
1345         object *v;
1346 {
1347         INCREF(v);
1348         return v;
1349 }
1350
1351 static object *
1352 long_float(v)
1353         object *v;
1354 {
1355         return newfloatobject(dgetlongvalue(v));
1356 }
1357
1358 static object *
1359 long_oct(v)
1360         object *v;
1361 {
1362         return long_format(v, 8);
1363 }
1364
1365 static object *
1366 long_hex(v)
1367         object *v;
1368 {
1369         return long_format(v, 16);
1370 }
1371
1372
1373 #define UF (unaryfunc)
1374 #define BF (binaryfunc)
1375 #define TF (ternaryfunc)
1376 #define IF (inquiry)
1377
1378 static number_methods long_as_number = {
1379         BF long_add,    /*nb_add*/
1380         BF long_sub,    /*nb_subtract*/
1381         BF long_mul,    /*nb_multiply*/
1382         BF long_div,    /*nb_divide*/
1383         BF long_mod,    /*nb_remainder*/
1384         BF long_divmod, /*nb_divmod*/
1385         TF long_pow,    /*nb_power*/
1386         UF long_neg,    /*nb_negative*/
1387         UF long_pos,    /*tp_positive*/
1388         UF long_abs,    /*tp_absolute*/
1389         IF long_nonzero,/*tp_nonzero*/
1390         UF long_invert, /*nb_invert*/
1391         BF long_lshift, /*nb_lshift*/
1392         BF long_rshift, /*nb_rshift*/
1393         BF long_and,    /*nb_and*/
1394         BF long_xor,    /*nb_xor*/
1395         BF long_or,     /*nb_or*/
1396         (int (*) FPROTO((object **, object **)))
1397         (coercion)long_coerce, /*nb_coerce*/
1398         UF long_int,    /*nb_int*/
1399         UF long_long,   /*nb_long*/
1400         UF long_float,  /*nb_float*/
1401         UF long_oct,    /*nb_oct*/
1402         UF long_hex,    /*nb_hex*/
1403 };
1404
1405 typeobject Longtype = {
1406         OB_HEAD_INIT(&Typetype)
1407         0,
1408         "long int",
1409         sizeof(longobject) - sizeof(digit),
1410         sizeof(digit),
1411         (destructor)long_dealloc, /*tp_dealloc*/
1412         0,              /*tp_print*/
1413         0,              /*tp_getattr*/
1414         0,              /*tp_setattr*/
1415         (int (*) FPROTO((object *, object *)))
1416         (cmpfunc)long_compare, /*tp_compare*/
1417         (reprfunc)long_repr, /*tp_repr*/
1418         &long_as_number,/*tp_as_number*/
1419         0,              /*tp_as_sequence*/
1420         0,              /*tp_as_mapping*/
1421         (long (*) FPROTO((object *)))
1422         (hashfunc)long_hash, /*tp_hash*/
1423 };