001    /* java.math.BigInteger -- Arbitary precision integers
002       Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2005, 2006  Free Software Foundation, Inc.
003    
004    This file is part of GNU Classpath.
005    
006    GNU Classpath is free software; you can redistribute it and/or modify
007    it under the terms of the GNU General Public License as published by
008    the Free Software Foundation; either version 2, or (at your option)
009    any later version.
010     
011    GNU Classpath is distributed in the hope that it will be useful, but
012    WITHOUT ANY WARRANTY; without even the implied warranty of
013    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
014    General Public License for more details.
015    
016    You should have received a copy of the GNU General Public License
017    along with GNU Classpath; see the file COPYING.  If not, write to the
018    Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
019    02110-1301 USA.
020    
021    Linking this library statically or dynamically with other modules is
022    making a combined work based on this library.  Thus, the terms and
023    conditions of the GNU General Public License cover the whole
024    combination.
025    
026    As a special exception, the copyright holders of this library give you
027    permission to link this library with independent modules to produce an
028    executable, regardless of the license terms of these independent
029    modules, and to copy and distribute the resulting executable under
030    terms of your choice, provided that you also meet, for each linked
031    independent module, the terms and conditions of the license of that
032    module.  An independent module is a module which is not derived from
033    or based on this library.  If you modify this library, you may extend
034    this exception to your version of the library, but you are not
035    obligated to do so.  If you do not wish to do so, delete this
036    exception statement from your version. */
037    
038    
039    package java.math;
040    
041    import gnu.java.math.MPN;
042    
043    import java.io.IOException;
044    import java.io.ObjectInputStream;
045    import java.io.ObjectOutputStream;
046    import java.util.Random;
047    
048    /**
049     * Written using on-line Java Platform 1.2 API Specification, as well
050     * as "The Java Class Libraries", 2nd edition (Addison-Wesley, 1998) and
051     * "Applied Cryptography, Second Edition" by Bruce Schneier (Wiley, 1996).
052     * 
053     * Based primarily on IntNum.java BitOps.java by Per Bothner (per@bothner.com)
054     * (found in Kawa 1.6.62).
055     *
056     * @author Warren Levy (warrenl@cygnus.com)
057     * @date December 20, 1999.
058     * @status believed complete and correct.
059     */
060    public class BigInteger extends Number implements Comparable<BigInteger>
061    {
062      /** All integers are stored in 2's-complement form.
063       * If words == null, the ival is the value of this BigInteger.
064       * Otherwise, the first ival elements of words make the value
065       * of this BigInteger, stored in little-endian order, 2's-complement form. */
066      private transient int ival;
067      private transient int[] words;
068    
069      // Serialization fields.
070      private int bitCount = -1;
071      private int bitLength = -1;
072      private int firstNonzeroByteNum = -2;
073      private int lowestSetBit = -2;
074      private byte[] magnitude;
075      private int signum;
076      private static final long serialVersionUID = -8287574255936472291L;
077    
078    
079      /** We pre-allocate integers in the range minFixNum..maxFixNum. 
080       * Note that we must at least preallocate 0, 1, and 10.  */
081      private static final int minFixNum = -100;
082      private static final int maxFixNum = 1024;
083      private static final int numFixNum = maxFixNum-minFixNum+1;
084      private static final BigInteger[] smallFixNums = new BigInteger[numFixNum];
085    
086      static
087      {
088        for (int i = numFixNum;  --i >= 0; )
089          smallFixNums[i] = new BigInteger(i + minFixNum);
090      }
091    
092      /**
093       * The constant zero as a BigInteger.
094       * @since 1.2
095       */
096      public static final BigInteger ZERO = smallFixNums[0 - minFixNum];
097    
098      /**
099       * The constant one as a BigInteger.
100       * @since 1.2
101       */
102      public static final BigInteger ONE = smallFixNums[1 - minFixNum];
103    
104      /**
105       * The constant ten as a BigInteger.
106       * @since 1.5
107       */
108      public static final BigInteger TEN = smallFixNums[10 - minFixNum];
109    
110      /* Rounding modes: */
111      private static final int FLOOR = 1;
112      private static final int CEILING = 2;
113      private static final int TRUNCATE = 3;
114      private static final int ROUND = 4;
115    
116      /** When checking the probability of primes, it is most efficient to
117       * first check the factoring of small primes, so we'll use this array.
118       */
119      private static final int[] primes =
120        {   2,   3,   5,   7,  11,  13,  17,  19,  23,  29,  31,  37,  41,  43,
121           47,  53,  59,  61,  67,  71,  73,  79,  83,  89,  97, 101, 103, 107,
122          109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181,
123          191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251 };
124    
125      /** HAC (Handbook of Applied Cryptography), Alfred Menezes & al. Table 4.4. */
126      private static final int[] k =
127          {100,150,200,250,300,350,400,500,600,800,1250, Integer.MAX_VALUE};
128      private static final int[] t =
129          { 27, 18, 15, 12,  9,  8,  7,  6,  5,  4,   3, 2};
130    
131      private BigInteger()
132      {
133      }
134    
135      /* Create a new (non-shared) BigInteger, and initialize to an int. */
136      private BigInteger(int value)
137      {
138        ival = value;
139      }
140    
141      public BigInteger(String val, int radix)
142      {
143        BigInteger result = valueOf(val, radix);
144        this.ival = result.ival;
145        this.words = result.words;
146      }
147    
148      public BigInteger(String val)
149      {
150        this(val, 10);
151      }
152    
153      /* Create a new (non-shared) BigInteger, and initialize from a byte array. */
154      public BigInteger(byte[] val)
155      {
156        if (val == null || val.length < 1)
157          throw new NumberFormatException();
158    
159        words = byteArrayToIntArray(val, val[0] < 0 ? -1 : 0);
160        BigInteger result = make(words, words.length);
161        this.ival = result.ival;
162        this.words = result.words;
163      }
164    
165      public BigInteger(int signum, byte[] magnitude)
166      {
167        if (magnitude == null || signum > 1 || signum < -1)
168          throw new NumberFormatException();
169    
170        if (signum == 0)
171          {
172            int i;
173            for (i = magnitude.length - 1; i >= 0 && magnitude[i] == 0; --i)
174              ;
175            if (i >= 0)
176              throw new NumberFormatException();
177            return;
178          }
179    
180        // Magnitude is always positive, so don't ever pass a sign of -1.
181        words = byteArrayToIntArray(magnitude, 0);
182        BigInteger result = make(words, words.length);
183        this.ival = result.ival;
184        this.words = result.words;
185    
186        if (signum < 0)
187          setNegative();
188      }
189    
190      public BigInteger(int numBits, Random rnd)
191      {
192        if (numBits < 0)
193          throw new IllegalArgumentException();
194    
195        init(numBits, rnd);
196      }
197    
198      private void init(int numBits, Random rnd)
199      {
200        int highbits = numBits & 31;
201        // minimum number of bytes to store the above number of bits
202        int highBitByteCount = (highbits + 7) / 8;
203        // number of bits to discard from the last byte
204        int discardedBitCount = highbits % 8;
205        if (discardedBitCount != 0)
206          discardedBitCount = 8 - discardedBitCount;
207        byte[] highBitBytes = new byte[highBitByteCount];
208        if (highbits > 0)
209          {
210            rnd.nextBytes(highBitBytes);
211            highbits = (highBitBytes[highBitByteCount - 1] & 0xFF) >>> discardedBitCount;
212            for (int i = highBitByteCount - 2; i >= 0; i--)
213              highbits = (highbits << 8) | (highBitBytes[i] & 0xFF);
214          }
215        int nwords = numBits / 32;
216    
217        while (highbits == 0 && nwords > 0)
218          {
219            highbits = rnd.nextInt();
220            --nwords;
221          }
222        if (nwords == 0 && highbits >= 0)
223          {
224            ival = highbits;
225          }
226        else
227          {
228            ival = highbits < 0 ? nwords + 2 : nwords + 1;
229            words = new int[ival];
230            words[nwords] = highbits;
231            while (--nwords >= 0)
232              words[nwords] = rnd.nextInt();
233          }
234      }
235    
236      public BigInteger(int bitLength, int certainty, Random rnd)
237      {
238        this(bitLength, rnd);
239    
240        // Keep going until we find a probable prime.
241        BigInteger result;
242        while (true)
243          {
244            // ...but first ensure that BI has bitLength bits
245            result = setBit(bitLength - 1);
246            this.ival = result.ival;
247            this.words = result.words;
248            if (isProbablePrime(certainty))
249              return;
250    
251            init(bitLength, rnd);
252          }
253      }
254    
255      /** 
256       *  Return a BigInteger that is bitLength bits long with a
257       *  probability < 2^-100 of being composite.
258       *
259       *  @param bitLength length in bits of resulting number
260       *  @param rnd random number generator to use
261       *  @throws ArithmeticException if bitLength < 2
262       *  @since 1.4
263       */
264      public static BigInteger probablePrime(int bitLength, Random rnd)
265      {
266        if (bitLength < 2)
267          throw new ArithmeticException();
268    
269        return new BigInteger(bitLength, 100, rnd);
270      }
271    
272      /** Return a (possibly-shared) BigInteger with a given long value. */
273      public static BigInteger valueOf(long val)
274      {
275        if (val >= minFixNum && val <= maxFixNum)
276          return smallFixNums[(int) val - minFixNum];
277        int i = (int) val;
278        if ((long) i == val)
279          return new BigInteger(i);
280        BigInteger result = alloc(2);
281        result.ival = 2;
282        result.words[0] = i;
283        result.words[1] = (int)(val >> 32);
284        return result;
285      }
286    
287      /** Make a canonicalized BigInteger from an array of words.
288       * The array may be reused (without copying). */
289      private static BigInteger make(int[] words, int len)
290      {
291        if (words == null)
292          return valueOf(len);
293        len = BigInteger.wordsNeeded(words, len);
294        if (len <= 1)
295          return len == 0 ? ZERO : valueOf(words[0]);
296        BigInteger num = new BigInteger();
297        num.words = words;
298        num.ival = len;
299        return num;
300      }
301    
302      /** Convert a big-endian byte array to a little-endian array of words. */
303      private static int[] byteArrayToIntArray(byte[] bytes, int sign)
304      {
305        // Determine number of words needed.
306        int[] words = new int[bytes.length/4 + 1];
307        int nwords = words.length;
308    
309        // Create a int out of modulo 4 high order bytes.
310        int bptr = 0;
311        int word = sign;
312        for (int i = bytes.length % 4; i > 0; --i, bptr++)
313          word = (word << 8) | (bytes[bptr] & 0xff);
314        words[--nwords] = word;
315    
316        // Elements remaining in byte[] are a multiple of 4.
317        while (nwords > 0)
318          words[--nwords] = bytes[bptr++] << 24 |
319                            (bytes[bptr++] & 0xff) << 16 |
320                            (bytes[bptr++] & 0xff) << 8 |
321                            (bytes[bptr++] & 0xff);
322        return words;
323      }
324    
325      /** Allocate a new non-shared BigInteger.
326       * @param nwords number of words to allocate
327       */
328      private static BigInteger alloc(int nwords)
329      {
330        BigInteger result = new BigInteger();
331        if (nwords > 1)
332        result.words = new int[nwords];
333        return result;
334      }
335    
336      /** Change words.length to nwords.
337       * We allow words.length to be upto nwords+2 without reallocating.
338       */
339      private void realloc(int nwords)
340      {
341        if (nwords == 0)
342          {
343            if (words != null)
344              {
345                if (ival > 0)
346                  ival = words[0];
347                words = null;
348              }
349          }
350        else if (words == null
351                 || words.length < nwords
352                 || words.length > nwords + 2)
353          {
354            int[] new_words = new int [nwords];
355            if (words == null)
356              {
357                new_words[0] = ival;
358                ival = 1;
359              }
360            else
361              {
362                if (nwords < ival)
363                  ival = nwords;
364                System.arraycopy(words, 0, new_words, 0, ival);
365              }
366            words = new_words;
367          }
368      }
369    
370      private boolean isNegative()
371      {
372        return (words == null ? ival : words[ival - 1]) < 0;
373      }
374    
375      public int signum()
376      {
377        if (ival == 0 && words == null)
378          return 0;
379        int top = words == null ? ival : words[ival-1];
380        return top < 0 ? -1 : 1;
381      }
382    
383      private static int compareTo(BigInteger x, BigInteger y)
384      {
385        if (x.words == null && y.words == null)
386          return x.ival < y.ival ? -1 : x.ival > y.ival ? 1 : 0;
387        boolean x_negative = x.isNegative();
388        boolean y_negative = y.isNegative();
389        if (x_negative != y_negative)
390          return x_negative ? -1 : 1;
391        int x_len = x.words == null ? 1 : x.ival;
392        int y_len = y.words == null ? 1 : y.ival;
393        if (x_len != y_len)
394          return (x_len > y_len) != x_negative ? 1 : -1;
395        return MPN.cmp(x.words, y.words, x_len);
396      }
397    
398      /** @since 1.2 */
399      public int compareTo(BigInteger val)
400      {
401        return compareTo(this, val);
402      }
403    
404      public BigInteger min(BigInteger val)
405      {
406        return compareTo(this, val) < 0 ? this : val;
407      }
408    
409      public BigInteger max(BigInteger val)
410      {
411        return compareTo(this, val) > 0 ? this : val;
412      }
413    
414      private boolean isZero()
415      {
416        return words == null && ival == 0;
417      }
418    
419      private boolean isOne()
420      {
421        return words == null && ival == 1;
422      }
423    
424      /** Calculate how many words are significant in words[0:len-1].
425       * Returns the least value x such that x>0 && words[0:x-1]==words[0:len-1],
426       * when words is viewed as a 2's complement integer.
427       */
428      private static int wordsNeeded(int[] words, int len)
429      {
430        int i = len;
431        if (i > 0)
432          {
433            int word = words[--i];
434            if (word == -1)
435              {
436                while (i > 0 && (word = words[i - 1]) < 0)
437                  {
438                    i--;
439                    if (word != -1) break;
440                  }
441              }
442            else
443              {
444                while (word == 0 && i > 0 && (word = words[i - 1]) >= 0)  i--;
445              }
446          }
447        return i + 1;
448      }
449    
450      private BigInteger canonicalize()
451      {
452        if (words != null
453            && (ival = BigInteger.wordsNeeded(words, ival)) <= 1)
454          {
455            if (ival == 1)
456              ival = words[0];
457            words = null;
458          }
459        if (words == null && ival >= minFixNum && ival <= maxFixNum)
460          return smallFixNums[ival - minFixNum];
461        return this;
462      }
463    
464      /** Add two ints, yielding a BigInteger. */
465      private static BigInteger add(int x, int y)
466      {
467        return valueOf((long) x + (long) y);
468      }
469    
470      /** Add a BigInteger and an int, yielding a new BigInteger. */
471      private static BigInteger add(BigInteger x, int y)
472      {
473        if (x.words == null)
474          return BigInteger.add(x.ival, y);
475        BigInteger result = new BigInteger(0);
476        result.setAdd(x, y);
477        return result.canonicalize();
478      }
479    
480      /** Set this to the sum of x and y.
481       * OK if x==this. */
482      private void setAdd(BigInteger x, int y)
483      {
484        if (x.words == null)
485          {
486            set((long) x.ival + (long) y);
487            return;
488          }
489        int len = x.ival;
490        realloc(len + 1);
491        long carry = y;
492        for (int i = 0;  i < len;  i++)
493          {
494            carry += ((long) x.words[i] & 0xffffffffL);
495            words[i] = (int) carry;
496            carry >>= 32;
497          }
498        if (x.words[len - 1] < 0)
499          carry--;
500        words[len] = (int) carry;
501        ival = wordsNeeded(words, len + 1);
502      }
503    
504      /** Destructively add an int to this. */
505      private void setAdd(int y)
506      {
507        setAdd(this, y);
508      }
509    
510      /** Destructively set the value of this to a long. */
511      private void set(long y)
512      {
513        int i = (int) y;
514        if ((long) i == y)
515          {
516            ival = i;
517            words = null;
518          }
519        else
520          {
521            realloc(2);
522            words[0] = i;
523            words[1] = (int) (y >> 32);
524            ival = 2;
525          }
526      }
527    
528      /** Destructively set the value of this to the given words.
529      * The words array is reused, not copied. */
530      private void set(int[] words, int length)
531      {
532        this.ival = length;
533        this.words = words;
534      }
535    
536      /** Destructively set the value of this to that of y. */
537      private void set(BigInteger y)
538      {
539        if (y.words == null)
540          set(y.ival);
541        else if (this != y)
542          {
543            realloc(y.ival);
544            System.arraycopy(y.words, 0, words, 0, y.ival);
545            ival = y.ival;
546          }
547      }
548    
549      /** Add two BigIntegers, yielding their sum as another BigInteger. */
550      private static BigInteger add(BigInteger x, BigInteger y, int k)
551      {
552        if (x.words == null && y.words == null)
553          return valueOf((long) k * (long) y.ival + (long) x.ival);
554        if (k != 1)
555          {
556            if (k == -1)
557              y = BigInteger.neg(y);
558            else
559              y = BigInteger.times(y, valueOf(k));
560          }
561        if (x.words == null)
562          return BigInteger.add(y, x.ival);
563        if (y.words == null)
564          return BigInteger.add(x, y.ival);
565        // Both are big
566        if (y.ival > x.ival)
567          { // Swap so x is longer then y.
568            BigInteger tmp = x;  x = y;  y = tmp;
569          }
570        BigInteger result = alloc(x.ival + 1);
571        int i = y.ival;
572        long carry = MPN.add_n(result.words, x.words, y.words, i);
573        long y_ext = y.words[i - 1] < 0 ? 0xffffffffL : 0;
574        for (; i < x.ival;  i++)
575          {
576            carry += ((long) x.words[i] & 0xffffffffL) + y_ext;
577            result.words[i] = (int) carry;
578            carry >>>= 32;
579          }
580        if (x.words[i - 1] < 0)
581          y_ext--;
582        result.words[i] = (int) (carry + y_ext);
583        result.ival = i+1;
584        return result.canonicalize();
585      }
586    
587      public BigInteger add(BigInteger val)
588      {
589        return add(this, val, 1);
590      }
591    
592      public BigInteger subtract(BigInteger val)
593      {
594        return add(this, val, -1);
595      }
596    
597      private static BigInteger times(BigInteger x, int y)
598      {
599        if (y == 0)
600          return ZERO;
601        if (y == 1)
602          return x;
603        int[] xwords = x.words;
604        int xlen = x.ival;
605        if (xwords == null)
606          return valueOf((long) xlen * (long) y);
607        boolean negative;
608        BigInteger result = BigInteger.alloc(xlen + 1);
609        if (xwords[xlen - 1] < 0)
610          {
611            negative = true;
612            negate(result.words, xwords, xlen);
613            xwords = result.words;
614          }
615        else
616          negative = false;
617        if (y < 0)
618          {
619            negative = !negative;
620            y = -y;
621          }
622        result.words[xlen] = MPN.mul_1(result.words, xwords, xlen, y);
623        result.ival = xlen + 1;
624        if (negative)
625          result.setNegative();
626        return result.canonicalize();
627      }
628    
629      private static BigInteger times(BigInteger x, BigInteger y)
630      {
631        if (y.words == null)
632          return times(x, y.ival);
633        if (x.words == null)
634          return times(y, x.ival);
635        boolean negative = false;
636        int[] xwords;
637        int[] ywords;
638        int xlen = x.ival;
639        int ylen = y.ival;
640        if (x.isNegative())
641          {
642            negative = true;
643            xwords = new int[xlen];
644            negate(xwords, x.words, xlen);
645          }
646        else
647          {
648            negative = false;
649            xwords = x.words;
650          }
651        if (y.isNegative())
652          {
653            negative = !negative;
654            ywords = new int[ylen];
655            negate(ywords, y.words, ylen);
656          }
657        else
658          ywords = y.words;
659        // Swap if x is shorter then y.
660        if (xlen < ylen)
661          {
662            int[] twords = xwords;  xwords = ywords;  ywords = twords;
663            int tlen = xlen;  xlen = ylen;  ylen = tlen;
664          }
665        BigInteger result = BigInteger.alloc(xlen+ylen);
666        MPN.mul(result.words, xwords, xlen, ywords, ylen);
667        result.ival = xlen+ylen;
668        if (negative)
669          result.setNegative();
670        return result.canonicalize();
671      }
672    
673      public BigInteger multiply(BigInteger y)
674      {
675        return times(this, y);
676      }
677    
678      private static void divide(long x, long y,
679                                 BigInteger quotient, BigInteger remainder,
680                                 int rounding_mode)
681      {
682        boolean xNegative, yNegative;
683        if (x < 0)
684          {
685            xNegative = true;
686            if (x == Long.MIN_VALUE)
687              {
688                divide(valueOf(x), valueOf(y),
689                       quotient, remainder, rounding_mode);
690                return;
691              }
692            x = -x;
693          }
694        else
695          xNegative = false;
696    
697        if (y < 0)
698          {
699            yNegative = true;
700            if (y == Long.MIN_VALUE)
701              {
702                if (rounding_mode == TRUNCATE)
703                  { // x != Long.Min_VALUE implies abs(x) < abs(y)
704                    if (quotient != null)
705                      quotient.set(0);
706                    if (remainder != null)
707                      remainder.set(x);
708                  }
709                else
710                  divide(valueOf(x), valueOf(y),
711                          quotient, remainder, rounding_mode);
712                return;
713              }
714            y = -y;
715          }
716        else
717          yNegative = false;
718    
719        long q = x / y;
720        long r = x % y;
721        boolean qNegative = xNegative ^ yNegative;
722    
723        boolean add_one = false;
724        if (r != 0)
725          {
726            switch (rounding_mode)
727              {
728              case TRUNCATE:
729                break;
730              case CEILING:
731              case FLOOR:
732                if (qNegative == (rounding_mode == FLOOR))
733                  add_one = true;
734                break;
735              case ROUND:
736                add_one = r > ((y - (q & 1)) >> 1);
737                break;
738              }
739          }
740        if (quotient != null)
741          {
742            if (add_one)
743              q++;
744            if (qNegative)
745              q = -q;
746            quotient.set(q);
747          }
748        if (remainder != null)
749          {
750            // The remainder is by definition: X-Q*Y
751            if (add_one)
752              {
753                // Subtract the remainder from Y.
754                r = y - r;
755                // In this case, abs(Q*Y) > abs(X).
756                // So sign(remainder) = -sign(X).
757                xNegative = ! xNegative;
758              }
759            else
760              {
761                // If !add_one, then: abs(Q*Y) <= abs(X).
762                // So sign(remainder) = sign(X).
763              }
764            if (xNegative)
765              r = -r;
766            remainder.set(r);
767          }
768      }
769    
770      /** Divide two integers, yielding quotient and remainder.
771       * @param x the numerator in the division
772       * @param y the denominator in the division
773       * @param quotient is set to the quotient of the result (iff quotient!=null)
774       * @param remainder is set to the remainder of the result
775       *  (iff remainder!=null)
776       * @param rounding_mode one of FLOOR, CEILING, TRUNCATE, or ROUND.
777       */
778      private static void divide(BigInteger x, BigInteger y,
779                                 BigInteger quotient, BigInteger remainder,
780                                 int rounding_mode)
781      {
782        if ((x.words == null || x.ival <= 2)
783            && (y.words == null || y.ival <= 2))
784          {
785            long x_l = x.longValue();
786            long y_l = y.longValue();
787            if (x_l != Long.MIN_VALUE && y_l != Long.MIN_VALUE)
788              {
789                divide(x_l, y_l, quotient, remainder, rounding_mode);
790                return;
791              }
792          }
793    
794        boolean xNegative = x.isNegative();
795        boolean yNegative = y.isNegative();
796        boolean qNegative = xNegative ^ yNegative;
797    
798        int ylen = y.words == null ? 1 : y.ival;
799        int[] ywords = new int[ylen];
800        y.getAbsolute(ywords);
801        while (ylen > 1 && ywords[ylen - 1] == 0)  ylen--;
802    
803        int xlen = x.words == null ? 1 : x.ival;
804        int[] xwords = new int[xlen+2];
805        x.getAbsolute(xwords);
806        while (xlen > 1 && xwords[xlen-1] == 0)  xlen--;
807    
808        int qlen, rlen;
809    
810        int cmpval = MPN.cmp(xwords, xlen, ywords, ylen);
811        if (cmpval < 0)  // abs(x) < abs(y)
812          { // quotient = 0;  remainder = num.
813            int[] rwords = xwords;  xwords = ywords;  ywords = rwords;
814            rlen = xlen;  qlen = 1;  xwords[0] = 0;
815          }
816        else if (cmpval == 0)  // abs(x) == abs(y)
817          {
818            xwords[0] = 1;  qlen = 1;  // quotient = 1
819            ywords[0] = 0;  rlen = 1;  // remainder = 0;
820          }
821        else if (ylen == 1)
822          {
823            qlen = xlen;
824            // Need to leave room for a word of leading zeros if dividing by 1
825            // and the dividend has the high bit set.  It might be safe to
826            // increment qlen in all cases, but it certainly is only necessary
827            // in the following case.
828            if (ywords[0] == 1 && xwords[xlen-1] < 0)
829              qlen++;
830            rlen = 1;
831            ywords[0] = MPN.divmod_1(xwords, xwords, xlen, ywords[0]);
832          }
833        else  // abs(x) > abs(y)
834          {
835            // Normalize the denominator, i.e. make its most significant bit set by
836            // shifting it normalization_steps bits to the left.  Also shift the
837            // numerator the same number of steps (to keep the quotient the same!).
838    
839            int nshift = MPN.count_leading_zeros(ywords[ylen - 1]);
840            if (nshift != 0)
841              {
842                // Shift up the denominator setting the most significant bit of
843                // the most significant word.
844                MPN.lshift(ywords, 0, ywords, ylen, nshift);
845    
846                // Shift up the numerator, possibly introducing a new most
847                // significant word.
848                int x_high = MPN.lshift(xwords, 0, xwords, xlen, nshift);
849                xwords[xlen++] = x_high;
850              }
851    
852            if (xlen == ylen)
853              xwords[xlen++] = 0;
854            MPN.divide(xwords, xlen, ywords, ylen);
855            rlen = ylen;
856            MPN.rshift0 (ywords, xwords, 0, rlen, nshift);
857    
858            qlen = xlen + 1 - ylen;
859            if (quotient != null)
860              {
861                for (int i = 0;  i < qlen;  i++)
862                  xwords[i] = xwords[i+ylen];
863              }
864          }
865    
866        if (ywords[rlen-1] < 0)
867          {
868            ywords[rlen] = 0;
869            rlen++;
870          }
871    
872        // Now the quotient is in xwords, and the remainder is in ywords.
873    
874        boolean add_one = false;
875        if (rlen > 1 || ywords[0] != 0)
876          { // Non-zero remainder i.e. in-exact quotient.
877            switch (rounding_mode)
878              {
879              case TRUNCATE:
880                break;
881              case CEILING:
882              case FLOOR:
883                if (qNegative == (rounding_mode == FLOOR))
884                  add_one = true;
885                break;
886              case ROUND:
887                // int cmp = compareTo(remainder<<1, abs(y));
888                BigInteger tmp = remainder == null ? new BigInteger() : remainder;
889                tmp.set(ywords, rlen);
890                tmp = shift(tmp, 1);
891                if (yNegative)
892                  tmp.setNegative();
893                int cmp = compareTo(tmp, y);
894                // Now cmp == compareTo(sign(y)*(remainder<<1), y)
895                if (yNegative)
896                  cmp = -cmp;
897                add_one = (cmp == 1) || (cmp == 0 && (xwords[0]&1) != 0);
898              }
899          }
900        if (quotient != null)
901          {
902            quotient.set(xwords, qlen);
903            if (qNegative)
904              {
905                if (add_one)  // -(quotient + 1) == ~(quotient)
906                  quotient.setInvert();
907                else
908                  quotient.setNegative();
909              }
910            else if (add_one)
911              quotient.setAdd(1);
912          }
913        if (remainder != null)
914          {
915            // The remainder is by definition: X-Q*Y
916            remainder.set(ywords, rlen);
917            if (add_one)
918              {
919                // Subtract the remainder from Y:
920                // abs(R) = abs(Y) - abs(orig_rem) = -(abs(orig_rem) - abs(Y)).
921                BigInteger tmp;
922                if (y.words == null)
923                  {
924                    tmp = remainder;
925                    tmp.set(yNegative ? ywords[0] + y.ival : ywords[0] - y.ival);
926                  }
927                else
928                  tmp = BigInteger.add(remainder, y, yNegative ? 1 : -1);
929                // Now tmp <= 0.
930                // In this case, abs(Q) = 1 + floor(abs(X)/abs(Y)).
931                // Hence, abs(Q*Y) > abs(X).
932                // So sign(remainder) = -sign(X).
933                if (xNegative)
934                  remainder.setNegative(tmp);
935                else
936                  remainder.set(tmp);
937              }
938            else
939              {
940                // If !add_one, then: abs(Q*Y) <= abs(X).
941                // So sign(remainder) = sign(X).
942                if (xNegative)
943                  remainder.setNegative();
944              }
945          }
946      }
947    
948      public BigInteger divide(BigInteger val)
949      {
950        if (val.isZero())
951          throw new ArithmeticException("divisor is zero");
952    
953        BigInteger quot = new BigInteger();
954        divide(this, val, quot, null, TRUNCATE);
955        return quot.canonicalize();
956      }
957    
958      public BigInteger remainder(BigInteger val)
959      {
960        if (val.isZero())
961          throw new ArithmeticException("divisor is zero");
962    
963        BigInteger rem = new BigInteger();
964        divide(this, val, null, rem, TRUNCATE);
965        return rem.canonicalize();
966      }
967    
968      public BigInteger[] divideAndRemainder(BigInteger val)
969      {
970        if (val.isZero())
971          throw new ArithmeticException("divisor is zero");
972    
973        BigInteger[] result = new BigInteger[2];
974        result[0] = new BigInteger();
975        result[1] = new BigInteger();
976        divide(this, val, result[0], result[1], TRUNCATE);
977        result[0].canonicalize();
978        result[1].canonicalize();
979        return result;
980      }
981    
982      public BigInteger mod(BigInteger m)
983      {
984        if (m.isNegative() || m.isZero())
985          throw new ArithmeticException("non-positive modulus");
986    
987        BigInteger rem = new BigInteger();
988        divide(this, m, null, rem, FLOOR);
989        return rem.canonicalize();
990      }
991    
992      /** Calculate the integral power of a BigInteger.
993       * @param exponent the exponent (must be non-negative)
994       */
995      public BigInteger pow(int exponent)
996      {
997        if (exponent <= 0)
998          {
999            if (exponent == 0)
1000              return ONE;
1001              throw new ArithmeticException("negative exponent");
1002          }
1003        if (isZero())
1004          return this;
1005        int plen = words == null ? 1 : ival;  // Length of pow2.
1006        int blen = ((bitLength() * exponent) >> 5) + 2 * plen;
1007        boolean negative = isNegative() && (exponent & 1) != 0;
1008        int[] pow2 = new int [blen];
1009        int[] rwords = new int [blen];
1010        int[] work = new int [blen];
1011        getAbsolute(pow2);  // pow2 = abs(this);
1012        int rlen = 1;
1013        rwords[0] = 1; // rwords = 1;
1014        for (;;)  // for (i = 0;  ; i++)
1015          {
1016            // pow2 == this**(2**i)
1017            // prod = this**(sum(j=0..i-1, (exponent>>j)&1))
1018            if ((exponent & 1) != 0)
1019              { // r *= pow2
1020                MPN.mul(work, pow2, plen, rwords, rlen);
1021                int[] temp = work;  work = rwords;  rwords = temp;
1022                rlen += plen;
1023                while (rwords[rlen - 1] == 0)  rlen--;
1024              }
1025            exponent >>= 1;
1026            if (exponent == 0)
1027              break;
1028            // pow2 *= pow2;
1029            MPN.mul(work, pow2, plen, pow2, plen);
1030            int[] temp = work;  work = pow2;  pow2 = temp;  // swap to avoid a copy
1031            plen *= 2;
1032            while (pow2[plen - 1] == 0)  plen--;
1033          }
1034        if (rwords[rlen - 1] < 0)
1035          rlen++;
1036        if (negative)
1037          negate(rwords, rwords, rlen);
1038        return BigInteger.make(rwords, rlen);
1039      }
1040    
1041      private static int[] euclidInv(int a, int b, int prevDiv)
1042      {
1043        if (b == 0)
1044          throw new ArithmeticException("not invertible");
1045    
1046        if (b == 1)
1047            // Success:  values are indeed invertible!
1048            // Bottom of the recursion reached; start unwinding.
1049            return new int[] { -prevDiv, 1 };
1050    
1051        int[] xy = euclidInv(b, a % b, a / b);      // Recursion happens here.
1052        a = xy[0]; // use our local copy of 'a' as a work var
1053        xy[0] = a * -prevDiv + xy[1];
1054        xy[1] = a;
1055        return xy;
1056      }
1057    
1058      private static void euclidInv(BigInteger a, BigInteger b,
1059                                    BigInteger prevDiv, BigInteger[] xy)
1060      {
1061        if (b.isZero())
1062          throw new ArithmeticException("not invertible");
1063    
1064        if (b.isOne())
1065          {
1066            // Success:  values are indeed invertible!
1067            // Bottom of the recursion reached; start unwinding.
1068            xy[0] = neg(prevDiv);
1069            xy[1] = ONE;
1070            return;
1071          }
1072    
1073        // Recursion happens in the following conditional!
1074    
1075        // If a just contains an int, then use integer math for the rest.
1076        if (a.words == null)
1077          {
1078            int[] xyInt = euclidInv(b.ival, a.ival % b.ival, a.ival / b.ival);
1079            xy[0] = new BigInteger(xyInt[0]);
1080            xy[1] = new BigInteger(xyInt[1]);
1081          }
1082        else
1083          {
1084            BigInteger rem = new BigInteger();
1085            BigInteger quot = new BigInteger();
1086            divide(a, b, quot, rem, FLOOR);
1087            // quot and rem may not be in canonical form. ensure
1088            rem.canonicalize();
1089            quot.canonicalize();
1090            euclidInv(b, rem, quot, xy);
1091          }
1092    
1093        BigInteger t = xy[0];
1094        xy[0] = add(xy[1], times(t, prevDiv), -1);
1095        xy[1] = t;
1096      }
1097    
1098      public BigInteger modInverse(BigInteger y)
1099      {
1100        if (y.isNegative() || y.isZero())
1101          throw new ArithmeticException("non-positive modulo");
1102    
1103        // Degenerate cases.
1104        if (y.isOne())
1105          return ZERO;
1106        if (isOne())
1107          return ONE;
1108    
1109        // Use Euclid's algorithm as in gcd() but do this recursively
1110        // rather than in a loop so we can use the intermediate results as we
1111        // unwind from the recursion.
1112        // Used http://www.math.nmsu.edu/~crypto/EuclideanAlgo.html as reference.
1113        BigInteger result = new BigInteger();
1114        boolean swapped = false;
1115    
1116        if (y.words == null)
1117          {
1118            // The result is guaranteed to be less than the modulus, y (which is
1119            // an int), so simplify this by working with the int result of this
1120            // modulo y.  Also, if this is negative, make it positive via modulo
1121            // math.  Note that BigInteger.mod() must be used even if this is
1122            // already an int as the % operator would provide a negative result if
1123            // this is negative, BigInteger.mod() never returns negative values.
1124            int xval = (words != null || isNegative()) ? mod(y).ival : ival;
1125            int yval = y.ival;
1126    
1127            // Swap values so x > y.
1128            if (yval > xval)
1129              {
1130                int tmp = xval; xval = yval; yval = tmp;
1131                swapped = true;
1132              }
1133            // Normally, the result is in the 2nd element of the array, but
1134            // if originally x < y, then x and y were swapped and the result
1135            // is in the 1st element of the array.
1136            result.ival =
1137              euclidInv(yval, xval % yval, xval / yval)[swapped ? 0 : 1];
1138    
1139            // Result can't be negative, so make it positive by adding the
1140            // original modulus, y.ival (not the possibly "swapped" yval).
1141            if (result.ival < 0)
1142              result.ival += y.ival;
1143          }
1144        else
1145          {
1146            // As above, force this to be a positive value via modulo math.
1147            BigInteger x = isNegative() ? this.mod(y) : this;
1148    
1149            // Swap values so x > y.
1150            if (x.compareTo(y) < 0)
1151              {
1152                result = x; x = y; y = result; // use 'result' as a work var
1153                swapped = true;
1154              }
1155            // As above (for ints), result will be in the 2nd element unless
1156            // the original x and y were swapped.
1157            BigInteger rem = new BigInteger();
1158            BigInteger quot = new BigInteger();
1159            divide(x, y, quot, rem, FLOOR);
1160            // quot and rem may not be in canonical form. ensure
1161            rem.canonicalize();
1162            quot.canonicalize();
1163            BigInteger[] xy = new BigInteger[2];
1164            euclidInv(y, rem, quot, xy);
1165            result = swapped ? xy[0] : xy[1];
1166    
1167            // Result can't be negative, so make it positive by adding the
1168            // original modulus, y (which is now x if they were swapped).
1169            if (result.isNegative())
1170              result = add(result, swapped ? x : y, 1);
1171          }
1172        
1173        return result;
1174      }
1175    
1176      public BigInteger modPow(BigInteger exponent, BigInteger m)
1177      {
1178        if (m.isNegative() || m.isZero())
1179          throw new ArithmeticException("non-positive modulo");
1180    
1181        if (exponent.isNegative())
1182          return modInverse(m).modPow(exponent.negate(), m);
1183        if (exponent.isOne())
1184          return mod(m);
1185    
1186        // To do this naively by first raising this to the power of exponent
1187        // and then performing modulo m would be extremely expensive, especially
1188        // for very large numbers.  The solution is found in Number Theory
1189        // where a combination of partial powers and moduli can be done easily.
1190        //
1191        // We'll use the algorithm for Additive Chaining which can be found on
1192        // p. 244 of "Applied Cryptography, Second Edition" by Bruce Schneier.
1193        BigInteger s = ONE;
1194        BigInteger t = this;
1195        BigInteger u = exponent;
1196    
1197        while (!u.isZero())
1198          {
1199            if (u.and(ONE).isOne())
1200              s = times(s, t).mod(m);
1201            u = u.shiftRight(1);
1202            t = times(t, t).mod(m);
1203          }
1204    
1205        return s;
1206      }
1207    
1208      /** Calculate Greatest Common Divisor for non-negative ints. */
1209      private static int gcd(int a, int b)
1210      {
1211        // Euclid's algorithm, copied from libg++.
1212        int tmp;
1213        if (b > a)
1214          {
1215            tmp = a; a = b; b = tmp;
1216          }
1217        for(;;)
1218          {
1219            if (b == 0)
1220              return a;
1221            if (b == 1)
1222              return b;
1223            tmp = b;
1224                b = a % b;
1225                a = tmp;
1226              }
1227          }
1228    
1229      public BigInteger gcd(BigInteger y)
1230      {
1231        int xval = ival;
1232        int yval = y.ival;
1233        if (words == null)
1234          {
1235            if (xval == 0)
1236              return abs(y);
1237            if (y.words == null
1238                && xval != Integer.MIN_VALUE && yval != Integer.MIN_VALUE)
1239              {
1240                if (xval < 0)
1241                  xval = -xval;
1242                if (yval < 0)
1243                  yval = -yval;
1244                return valueOf(gcd(xval, yval));
1245              }
1246            xval = 1;
1247          }
1248        if (y.words == null)
1249          {
1250            if (yval == 0)
1251              return abs(this);
1252            yval = 1;
1253          }
1254        int len = (xval > yval ? xval : yval) + 1;
1255        int[] xwords = new int[len];
1256        int[] ywords = new int[len];
1257        getAbsolute(xwords);
1258        y.getAbsolute(ywords);
1259        len = MPN.gcd(xwords, ywords, len);
1260        BigInteger result = new BigInteger(0);
1261        result.ival = len;
1262        result.words = xwords;
1263        return result.canonicalize();
1264      }
1265    
1266      /**
1267       * <p>Returns <code>true</code> if this BigInteger is probably prime,
1268       * <code>false</code> if it's definitely composite. If <code>certainty</code>
1269       * is <code><= 0</code>, <code>true</code> is returned.</p>
1270       *
1271       * @param certainty a measure of the uncertainty that the caller is willing
1272       * to tolerate: if the call returns <code>true</code> the probability that
1273       * this BigInteger is prime exceeds <code>(1 - 1/2<sup>certainty</sup>)</code>.
1274       * The execution time of this method is proportional to the value of this
1275       * parameter.
1276       * @return <code>true</code> if this BigInteger is probably prime,
1277       * <code>false</code> if it's definitely composite.
1278       */
1279      public boolean isProbablePrime(int certainty)
1280      {
1281        if (certainty < 1)
1282          return true;
1283    
1284        /** We'll use the Rabin-Miller algorithm for doing a probabilistic
1285         * primality test.  It is fast, easy and has faster decreasing odds of a
1286         * composite passing than with other tests.  This means that this
1287         * method will actually have a probability much greater than the
1288         * 1 - .5^certainty specified in the JCL (p. 117), but I don't think
1289         * anyone will complain about better performance with greater certainty.
1290         *
1291         * The Rabin-Miller algorithm can be found on pp. 259-261 of "Applied
1292         * Cryptography, Second Edition" by Bruce Schneier.
1293         */
1294    
1295        // First rule out small prime factors
1296        BigInteger rem = new BigInteger();
1297        int i;
1298        for (i = 0; i < primes.length; i++)
1299          {
1300            if (words == null && ival == primes[i])
1301              return true;
1302    
1303            divide(this, smallFixNums[primes[i] - minFixNum], null, rem, TRUNCATE);
1304            if (rem.canonicalize().isZero())
1305              return false;
1306          }
1307    
1308        // Now perform the Rabin-Miller test.
1309    
1310        // Set b to the number of times 2 evenly divides (this - 1).
1311        // I.e. 2^b is the largest power of 2 that divides (this - 1).
1312        BigInteger pMinus1 = add(this, -1);
1313        int b = pMinus1.getLowestSetBit();
1314    
1315        // Set m such that this = 1 + 2^b * m.
1316        BigInteger m = pMinus1.divide(valueOf(2L << b - 1));
1317    
1318        // The HAC (Handbook of Applied Cryptography), Alfred Menezes & al. Note
1319        // 4.49 (controlling the error probability) gives the number of trials
1320        // for an error probability of 1/2**80, given the number of bits in the
1321        // number to test.  we shall use these numbers as is if/when 'certainty'
1322        // is less or equal to 80, and twice as much if it's greater.
1323        int bits = this.bitLength();
1324        for (i = 0; i < k.length; i++)
1325          if (bits <= k[i])
1326            break;
1327        int trials = t[i];
1328        if (certainty > 80)
1329          trials *= 2;
1330        BigInteger z;
1331        for (int t = 0; t < trials; t++)
1332          {
1333            // The HAC (Handbook of Applied Cryptography), Alfred Menezes & al.
1334            // Remark 4.28 states: "...A strategy that is sometimes employed
1335            // is to fix the bases a to be the first few primes instead of
1336            // choosing them at random.
1337            z = smallFixNums[primes[t] - minFixNum].modPow(m, this);
1338            if (z.isOne() || z.equals(pMinus1))
1339              continue;                     // Passes the test; may be prime.
1340    
1341            for (i = 0; i < b; )
1342              {
1343                if (z.isOne())
1344                  return false;
1345                i++;
1346                if (z.equals(pMinus1))
1347                  break;                    // Passes the test; may be prime.
1348    
1349                z = z.modPow(valueOf(2), this);
1350              }
1351    
1352            if (i == b && !z.equals(pMinus1))
1353              return false;
1354          }
1355        return true;
1356      }
1357    
1358      private void setInvert()
1359      {
1360        if (words == null)
1361          ival = ~ival;
1362        else
1363          {
1364            for (int i = ival;  --i >= 0; )
1365              words[i] = ~words[i];
1366          }
1367      }
1368    
1369      private void setShiftLeft(BigInteger x, int count)
1370      {
1371        int[] xwords;
1372        int xlen;
1373        if (x.words == null)
1374          {
1375            if (count < 32)
1376              {
1377                set((long) x.ival << count);
1378                return;
1379              }
1380            xwords = new int[1];
1381            xwords[0] = x.ival;
1382            xlen = 1;
1383          }
1384        else
1385          {
1386            xwords = x.words;
1387            xlen = x.ival;
1388          }
1389        int word_count = count >> 5;
1390        count &= 31;
1391        int new_len = xlen + word_count;
1392        if (count == 0)
1393          {
1394            realloc(new_len);
1395            for (int i = xlen;  --i >= 0; )
1396              words[i+word_count] = xwords[i];
1397          }
1398        else
1399          {
1400            new_len++;
1401            realloc(new_len);
1402            int shift_out = MPN.lshift(words, word_count, xwords, xlen, count);
1403            count = 32 - count;
1404            words[new_len-1] = (shift_out << count) >> count;  // sign-extend.
1405          }
1406        ival = new_len;
1407        for (int i = word_count;  --i >= 0; )
1408          words[i] = 0;
1409      }
1410    
1411      private void setShiftRight(BigInteger x, int count)
1412      {
1413        if (x.words == null)
1414          set(count < 32 ? x.ival >> count : x.ival < 0 ? -1 : 0);
1415        else if (count == 0)
1416          set(x);
1417        else
1418          {
1419            boolean neg = x.isNegative();
1420            int word_count = count >> 5;
1421            count &= 31;
1422            int d_len = x.ival - word_count;
1423            if (d_len <= 0)
1424              set(neg ? -1 : 0);
1425            else
1426              {
1427                if (words == null || words.length < d_len)
1428                  realloc(d_len);
1429                MPN.rshift0 (words, x.words, word_count, d_len, count);
1430                ival = d_len;
1431                if (neg)
1432                  words[d_len-1] |= -2 << (31 - count);
1433              }
1434          }
1435      }
1436    
1437      private void setShift(BigInteger x, int count)
1438      {
1439        if (count > 0)
1440          setShiftLeft(x, count);
1441        else
1442          setShiftRight(x, -count);
1443      }
1444    
1445      private static BigInteger shift(BigInteger x, int count)
1446      {
1447        if (x.words == null)
1448          {
1449            if (count <= 0)
1450              return valueOf(count > -32 ? x.ival >> (-count) : x.ival < 0 ? -1 : 0);
1451            if (count < 32)
1452              return valueOf((long) x.ival << count);
1453          }
1454        if (count == 0)
1455          return x;
1456        BigInteger result = new BigInteger(0);
1457        result.setShift(x, count);
1458        return result.canonicalize();
1459      }
1460    
1461      public BigInteger shiftLeft(int n)
1462      {
1463        return shift(this, n);
1464      }
1465    
1466      public BigInteger shiftRight(int n)
1467      {
1468        return shift(this, -n);
1469      }
1470    
1471      private void format(int radix, StringBuffer buffer)
1472      {
1473        if (words == null)
1474          buffer.append(Integer.toString(ival, radix));
1475        else if (ival <= 2)
1476          buffer.append(Long.toString(longValue(), radix));
1477        else
1478          {
1479            boolean neg = isNegative();
1480            int[] work;
1481            if (neg || radix != 16)
1482              {
1483                work = new int[ival];
1484                getAbsolute(work);
1485              }
1486            else
1487              work = words;
1488            int len = ival;
1489    
1490            if (radix == 16)
1491              {
1492                if (neg)
1493                  buffer.append('-');
1494                int buf_start = buffer.length();
1495                for (int i = len;  --i >= 0; )
1496                  {
1497                    int word = work[i];
1498                    for (int j = 8;  --j >= 0; )
1499                      {
1500                        int hex_digit = (word >> (4 * j)) & 0xF;
1501                        // Suppress leading zeros:
1502                        if (hex_digit > 0 || buffer.length() > buf_start)
1503                          buffer.append(Character.forDigit(hex_digit, 16));
1504                      }
1505                  }
1506              }
1507            else
1508              {
1509                int i = buffer.length();
1510                for (;;)
1511                  {
1512                    int digit = MPN.divmod_1(work, work, len, radix);
1513                    buffer.append(Character.forDigit(digit, radix));
1514                    while (len > 0 && work[len-1] == 0) len--;
1515                    if (len == 0)
1516                      break;
1517                  }
1518                if (neg)
1519                  buffer.append('-');
1520                /* Reverse buffer. */
1521                int j = buffer.length() - 1;
1522                while (i < j)
1523                  {
1524                    char tmp = buffer.charAt(i);
1525                    buffer.setCharAt(i, buffer.charAt(j));
1526                    buffer.setCharAt(j, tmp);
1527                    i++;  j--;
1528                  }
1529              }
1530          }
1531      }
1532    
1533      public String toString()
1534      {
1535        return toString(10);
1536      }
1537    
1538      public String toString(int radix)
1539      {
1540        if (words == null)
1541          return Integer.toString(ival, radix);
1542        if (ival <= 2)
1543          return Long.toString(longValue(), radix);
1544        int buf_size = ival * (MPN.chars_per_word(radix) + 1);
1545        StringBuffer buffer = new StringBuffer(buf_size);
1546        format(radix, buffer);
1547        return buffer.toString();
1548      }
1549    
1550      public int intValue()
1551      {
1552        if (words == null)
1553          return ival;
1554        return words[0];
1555      }
1556    
1557      public long longValue()
1558      {
1559        if (words == null)
1560          return ival;
1561        if (ival == 1)
1562          return words[0];
1563        return ((long)words[1] << 32) + ((long)words[0] & 0xffffffffL);
1564      }
1565    
1566      public int hashCode()
1567      {
1568        // FIXME: May not match hashcode of JDK.
1569        return words == null ? ival : (words[0] + words[ival - 1]);
1570      }
1571    
1572      /* Assumes x and y are both canonicalized. */
1573      private static boolean equals(BigInteger x, BigInteger y)
1574      {
1575        if (x.words == null && y.words == null)
1576          return x.ival == y.ival;
1577        if (x.words == null || y.words == null || x.ival != y.ival)
1578          return false;
1579        for (int i = x.ival; --i >= 0; )
1580          {
1581            if (x.words[i] != y.words[i])
1582              return false;
1583          }
1584        return true;
1585      }
1586    
1587      /* Assumes this and obj are both canonicalized. */
1588      public boolean equals(Object obj)
1589      {
1590        if (! (obj instanceof BigInteger))
1591          return false;
1592        return equals(this, (BigInteger) obj);
1593      }
1594    
1595      private static BigInteger valueOf(String s, int radix)
1596           throws NumberFormatException
1597      {
1598        int len = s.length();
1599        // Testing (len < MPN.chars_per_word(radix)) would be more accurate,
1600        // but slightly more expensive, for little practical gain.
1601        if (len <= 15 && radix <= 16)
1602          return valueOf(Long.parseLong(s, radix));
1603    
1604        int i, digit;
1605        boolean negative;
1606        byte[] bytes;
1607        char ch = s.charAt(0);
1608        if (ch == '-')
1609          {
1610            negative = true;
1611            i = 1;
1612            bytes = new byte[len - 1];
1613          }
1614        else
1615          {
1616            negative = false;
1617            i = 0;
1618            bytes = new byte[len];
1619          }
1620        int byte_len = 0;
1621        for ( ; i < len;  i++)
1622          {
1623            ch = s.charAt(i);
1624            digit = Character.digit(ch, radix);
1625            if (digit < 0)
1626              throw new NumberFormatException();
1627            bytes[byte_len++] = (byte) digit;
1628          }
1629        return valueOf(bytes, byte_len, negative, radix);
1630      }
1631    
1632      private static BigInteger valueOf(byte[] digits, int byte_len,
1633                                        boolean negative, int radix)
1634      {
1635        int chars_per_word = MPN.chars_per_word(radix);
1636        int[] words = new int[byte_len / chars_per_word + 1];
1637        int size = MPN.set_str(words, digits, byte_len, radix);
1638        if (size == 0)
1639          return ZERO;
1640        if (words[size-1] < 0)
1641          words[size++] = 0;
1642        if (negative)
1643          negate(words, words, size);
1644        return make(words, size);
1645      }
1646    
1647      public double doubleValue()
1648      {
1649        if (words == null)
1650          return (double) ival;
1651        if (ival <= 2)
1652          return (double) longValue();
1653        if (isNegative())
1654          return neg(this).roundToDouble(0, true, false);
1655          return roundToDouble(0, false, false);
1656      }
1657    
1658      public float floatValue()
1659      {
1660        return (float) doubleValue();
1661      }
1662    
1663      /** Return true if any of the lowest n bits are one.
1664       * (false if n is negative).  */
1665      private boolean checkBits(int n)
1666      {
1667        if (n <= 0)
1668          return false;
1669        if (words == null)
1670          return n > 31 || ((ival & ((1 << n) - 1)) != 0);
1671        int i;
1672        for (i = 0; i < (n >> 5) ; i++)
1673          if (words[i] != 0)
1674            return true;
1675        return (n & 31) != 0 && (words[i] & ((1 << (n & 31)) - 1)) != 0;
1676      }
1677    
1678      /** Convert a semi-processed BigInteger to double.
1679       * Number must be non-negative.  Multiplies by a power of two, applies sign,
1680       * and converts to double, with the usual java rounding.
1681       * @param exp power of two, positive or negative, by which to multiply
1682       * @param neg true if negative
1683       * @param remainder true if the BigInteger is the result of a truncating
1684       * division that had non-zero remainder.  To ensure proper rounding in
1685       * this case, the BigInteger must have at least 54 bits.  */
1686      private double roundToDouble(int exp, boolean neg, boolean remainder)
1687      {
1688        // Compute length.
1689        int il = bitLength();
1690    
1691        // Exponent when normalized to have decimal point directly after
1692        // leading one.  This is stored excess 1023 in the exponent bit field.
1693        exp += il - 1;
1694    
1695        // Gross underflow.  If exp == -1075, we let the rounding
1696        // computation determine whether it is minval or 0 (which are just
1697        // 0x0000 0000 0000 0001 and 0x0000 0000 0000 0000 as bit
1698        // patterns).
1699        if (exp < -1075)
1700          return neg ? -0.0 : 0.0;
1701    
1702        // gross overflow
1703        if (exp > 1023)
1704          return neg ? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
1705    
1706        // number of bits in mantissa, including the leading one.
1707        // 53 unless it's denormalized
1708        int ml = (exp >= -1022 ? 53 : 53 + exp + 1022);
1709    
1710        // Get top ml + 1 bits.  The extra one is for rounding.
1711        long m;
1712        int excess_bits = il - (ml + 1);
1713        if (excess_bits > 0)
1714          m = ((words == null) ? ival >> excess_bits
1715               : MPN.rshift_long(words, ival, excess_bits));
1716        else
1717          m = longValue() << (- excess_bits);
1718    
1719        // Special rounding for maxval.  If the number exceeds maxval by
1720        // any amount, even if it's less than half a step, it overflows.
1721        if (exp == 1023 && ((m >> 1) == (1L << 53) - 1))
1722          {
1723            if (remainder || checkBits(il - ml))
1724              return neg ? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
1725            else
1726              return neg ? - Double.MAX_VALUE : Double.MAX_VALUE;
1727          }
1728    
1729        // Normal round-to-even rule: round up if the bit dropped is a one, and
1730        // the bit above it or any of the bits below it is a one.
1731        if ((m & 1) == 1
1732            && ((m & 2) == 2 || remainder || checkBits(excess_bits)))
1733          {
1734            m += 2;
1735            // Check if we overflowed the mantissa
1736            if ((m & (1L << 54)) != 0)
1737              {
1738                exp++;
1739                // renormalize
1740                m >>= 1;
1741              }
1742            // Check if a denormalized mantissa was just rounded up to a
1743            // normalized one.
1744            else if (ml == 52 && (m & (1L << 53)) != 0)
1745              exp++;
1746          }
1747            
1748        // Discard the rounding bit
1749        m >>= 1;
1750    
1751        long bits_sign = neg ? (1L << 63) : 0;
1752        exp += 1023;
1753        long bits_exp = (exp <= 0) ? 0 : ((long)exp) << 52;
1754        long bits_mant = m & ~(1L << 52);
1755        return Double.longBitsToDouble(bits_sign | bits_exp | bits_mant);
1756      }
1757    
1758      /** Copy the abolute value of this into an array of words.
1759       * Assumes words.length >= (this.words == null ? 1 : this.ival).
1760       * Result is zero-extended, but need not be a valid 2's complement number.
1761       */
1762      private void getAbsolute(int[] words)
1763      {
1764        int len;
1765        if (this.words == null)
1766          {
1767            len = 1;
1768            words[0] = this.ival;
1769          }
1770        else
1771          {
1772            len = this.ival;
1773            for (int i = len;  --i >= 0; )
1774              words[i] = this.words[i];
1775          }
1776        if (words[len - 1] < 0)
1777          negate(words, words, len);
1778        for (int i = words.length;  --i > len; )
1779          words[i] = 0;
1780      }
1781    
1782      /** Set dest[0:len-1] to the negation of src[0:len-1].
1783       * Return true if overflow (i.e. if src is -2**(32*len-1)).
1784       * Ok for src==dest. */
1785      private static boolean negate(int[] dest, int[] src, int len)
1786      {
1787        long carry = 1;
1788        boolean negative = src[len-1] < 0;
1789        for (int i = 0;  i < len;  i++)
1790          {
1791            carry += ((long) (~src[i]) & 0xffffffffL);
1792            dest[i] = (int) carry;
1793            carry >>= 32;
1794          }
1795        return (negative && dest[len-1] < 0);
1796      }
1797    
1798      /** Destructively set this to the negative of x.
1799       * It is OK if x==this.*/
1800      private void setNegative(BigInteger x)
1801      {
1802        int len = x.ival;
1803        if (x.words == null)
1804          {
1805            if (len == Integer.MIN_VALUE)
1806              set(- (long) len);
1807            else
1808              set(-len);
1809            return;
1810          }
1811        realloc(len + 1);
1812        if (negate(words, x.words, len))
1813          words[len++] = 0;
1814        ival = len;
1815      }
1816    
1817      /** Destructively negate this. */
1818      private void setNegative()
1819      {
1820        setNegative(this);
1821      }
1822    
1823      private static BigInteger abs(BigInteger x)
1824      {
1825        return x.isNegative() ? neg(x) : x;
1826      }
1827    
1828      public BigInteger abs()
1829      {
1830        return abs(this);
1831      }
1832    
1833      private static BigInteger neg(BigInteger x)
1834      {
1835        if (x.words == null && x.ival != Integer.MIN_VALUE)
1836          return valueOf(- x.ival);
1837        BigInteger result = new BigInteger(0);
1838        result.setNegative(x);
1839        return result.canonicalize();
1840      }
1841    
1842      public BigInteger negate()
1843      {
1844        return neg(this);
1845      }
1846    
1847      /** Calculates ceiling(log2(this < 0 ? -this : this+1))
1848       * See Common Lisp: the Language, 2nd ed, p. 361.
1849       */
1850      public int bitLength()
1851      {
1852        if (words == null)
1853          return MPN.intLength(ival);
1854          return MPN.intLength(words, ival);
1855      }
1856    
1857      public byte[] toByteArray()
1858      {
1859        // Determine number of bytes needed.  The method bitlength returns
1860        // the size without the sign bit, so add one bit for that and then
1861        // add 7 more to emulate the ceil function using integer math.
1862        byte[] bytes = new byte[(bitLength() + 1 + 7) / 8];
1863        int nbytes = bytes.length;
1864    
1865        int wptr = 0;
1866        int word;
1867    
1868        // Deal with words array until one word or less is left to process.
1869        // If BigInteger is an int, then it is in ival and nbytes will be <= 4.
1870        while (nbytes > 4)
1871          {
1872            word = words[wptr++];
1873            for (int i = 4; i > 0; --i, word >>= 8)
1874              bytes[--nbytes] = (byte) word;
1875          }
1876    
1877        // Deal with the last few bytes.  If BigInteger is an int, use ival.
1878        word = (words == null) ? ival : words[wptr];
1879        for ( ; nbytes > 0; word >>= 8)
1880          bytes[--nbytes] = (byte) word;
1881    
1882        return bytes;
1883      }
1884    
1885      /** Return the boolean opcode (for bitOp) for swapped operands.
1886       * I.e. bitOp(swappedOp(op), x, y) == bitOp(op, y, x).
1887       */
1888      private static int swappedOp(int op)
1889      {
1890        return
1891        "\000\001\004\005\002\003\006\007\010\011\014\015\012\013\016\017"
1892        .charAt(op);
1893      }
1894    
1895      /** Do one the the 16 possible bit-wise operations of two BigIntegers. */
1896      private static BigInteger bitOp(int op, BigInteger x, BigInteger y)
1897      {
1898        switch (op)
1899          {
1900            case 0:  return ZERO;
1901            case 1:  return x.and(y);
1902            case 3:  return x;
1903            case 5:  return y;
1904            case 15: return valueOf(-1);
1905          }
1906        BigInteger result = new BigInteger();
1907        setBitOp(result, op, x, y);
1908        return result.canonicalize();
1909      }
1910    
1911      /** Do one the the 16 possible bit-wise operations of two BigIntegers. */
1912      private static void setBitOp(BigInteger result, int op,
1913                                   BigInteger x, BigInteger y)
1914      {
1915        if ((y.words != null) && (x.words == null || x.ival < y.ival))
1916          {
1917            BigInteger temp = x;  x = y;  y = temp;
1918            op = swappedOp(op);
1919          }
1920        int xi;
1921        int yi;
1922        int xlen, ylen;
1923        if (y.words == null)
1924          {
1925            yi = y.ival;
1926            ylen = 1;
1927          }
1928        else
1929          {
1930            yi = y.words[0];
1931            ylen = y.ival;
1932          }
1933        if (x.words == null)
1934          {
1935            xi = x.ival;
1936            xlen = 1;
1937          }
1938        else
1939          {
1940            xi = x.words[0];
1941            xlen = x.ival;
1942          }
1943        if (xlen > 1)
1944          result.realloc(xlen);
1945        int[] w = result.words;
1946        int i = 0;
1947        // Code for how to handle the remainder of x.
1948        // 0:  Truncate to length of y.
1949        // 1:  Copy rest of x.
1950        // 2:  Invert rest of x.
1951        int finish = 0;
1952        int ni;
1953        switch (op)
1954          {
1955          case 0:  // clr
1956            ni = 0;
1957            break;
1958          case 1: // and
1959            for (;;)
1960              {
1961                ni = xi & yi;
1962                if (i+1 >= ylen) break;
1963                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
1964              }
1965            if (yi < 0) finish = 1;
1966            break;
1967          case 2: // andc2
1968            for (;;)
1969              {
1970                ni = xi & ~yi;
1971                if (i+1 >= ylen) break;
1972                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
1973              }
1974            if (yi >= 0) finish = 1;
1975            break;
1976          case 3:  // copy x
1977            ni = xi;
1978            finish = 1;  // Copy rest
1979            break;
1980          case 4: // andc1
1981            for (;;)
1982              {
1983                ni = ~xi & yi;
1984                if (i+1 >= ylen) break;
1985                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
1986              }
1987            if (yi < 0) finish = 2;
1988            break;
1989          case 5: // copy y
1990            for (;;)
1991              {
1992                ni = yi;
1993                if (i+1 >= ylen) break;
1994                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
1995              }
1996            break;
1997          case 6:  // xor
1998            for (;;)
1999              {
2000                ni = xi ^ yi;
2001                if (i+1 >= ylen) break;
2002                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
2003              }
2004            finish = yi < 0 ? 2 : 1;
2005            break;
2006          case 7:  // ior
2007            for (;;)
2008              {
2009                ni = xi | yi;
2010                if (i+1 >= ylen) break;
2011                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
2012              }
2013            if (yi >= 0) finish = 1;
2014            break;
2015          case 8:  // nor
2016            for (;;)
2017              {
2018                ni = ~(xi | yi);
2019                if (i+1 >= ylen) break;
2020                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
2021              }
2022            if (yi >= 0)  finish = 2;
2023            break;
2024          case 9:  // eqv [exclusive nor]
2025            for (;;)
2026              {
2027                ni = ~(xi ^ yi);
2028                if (i+1 >= ylen) break;
2029                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
2030              }
2031            finish = yi >= 0 ? 2 : 1;
2032            break;
2033          case 10:  // c2
2034            for (;;)
2035              {
2036                ni = ~yi;
2037                if (i+1 >= ylen) break;
2038                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
2039              }
2040            break;
2041          case 11:  // orc2
2042            for (;;)
2043              {
2044                ni = xi | ~yi;
2045                if (i+1 >= ylen) break;
2046                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
2047              }
2048            if (yi < 0)  finish = 1;
2049            break;
2050          case 12:  // c1
2051            ni = ~xi;
2052            finish = 2;
2053            break;
2054          case 13:  // orc1
2055            for (;;)
2056              {
2057                ni = ~xi | yi;
2058                if (i+1 >= ylen) break;
2059                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
2060              }
2061            if (yi >= 0) finish = 2;
2062            break;
2063          case 14:  // nand
2064            for (;;)
2065              {
2066                ni = ~(xi & yi);
2067                if (i+1 >= ylen) break;
2068                w[i++] = ni;  xi = x.words[i];  yi = y.words[i];
2069              }
2070            if (yi < 0) finish = 2;
2071            break;
2072          default:
2073          case 15:  // set
2074            ni = -1;
2075            break;
2076          }
2077        // Here i==ylen-1; w[0]..w[i-1] have the correct result;
2078        // and ni contains the correct result for w[i+1].
2079        if (i+1 == xlen)
2080          finish = 0;
2081        switch (finish)
2082          {
2083          case 0:
2084            if (i == 0 && w == null)
2085              {
2086                result.ival = ni;
2087                return;
2088              }
2089            w[i++] = ni;
2090            break;
2091          case 1:  w[i] = ni;  while (++i < xlen)  w[i] = x.words[i];  break;
2092          case 2:  w[i] = ni;  while (++i < xlen)  w[i] = ~x.words[i];  break;
2093          }
2094        result.ival = i;
2095      }
2096    
2097      /** Return the logical (bit-wise) "and" of a BigInteger and an int. */
2098      private static BigInteger and(BigInteger x, int y)
2099      {
2100        if (x.words == null)
2101          return valueOf(x.ival & y);
2102        if (y >= 0)
2103          return valueOf(x.words[0] & y);
2104        int len = x.ival;
2105        int[] words = new int[len];
2106        words[0] = x.words[0] & y;
2107        while (--len > 0)
2108          words[len] = x.words[len];
2109        return make(words, x.ival);
2110      }
2111    
2112      /** Return the logical (bit-wise) "and" of two BigIntegers. */
2113      public BigInteger and(BigInteger y)
2114      {
2115        if (y.words == null)
2116          return and(this, y.ival);
2117        else if (words == null)
2118          return and(y, ival);
2119    
2120        BigInteger x = this;
2121        if (ival < y.ival)
2122          {
2123            BigInteger temp = this;  x = y;  y = temp;
2124          }
2125        int i;
2126        int len = y.isNegative() ? x.ival : y.ival;
2127        int[] words = new int[len];
2128        for (i = 0;  i < y.ival;  i++)
2129          words[i] = x.words[i] & y.words[i];
2130        for ( ; i < len;  i++)
2131          words[i] = x.words[i];
2132        return make(words, len);
2133      }
2134    
2135      /** Return the logical (bit-wise) "(inclusive) or" of two BigIntegers. */
2136      public BigInteger or(BigInteger y)
2137      {
2138        return bitOp(7, this, y);
2139      }
2140    
2141      /** Return the logical (bit-wise) "exclusive or" of two BigIntegers. */
2142      public BigInteger xor(BigInteger y)
2143      {
2144        return bitOp(6, this, y);
2145      }
2146    
2147      /** Return the logical (bit-wise) negation of a BigInteger. */
2148      public BigInteger not()
2149      {
2150        return bitOp(12, this, ZERO);
2151      }
2152    
2153      public BigInteger andNot(BigInteger val)
2154      {
2155        return and(val.not());
2156      }
2157    
2158      public BigInteger clearBit(int n)
2159      {
2160        if (n < 0)
2161          throw new ArithmeticException();
2162    
2163        return and(ONE.shiftLeft(n).not());
2164      }
2165    
2166      public BigInteger setBit(int n)
2167      {
2168        if (n < 0)
2169          throw new ArithmeticException();
2170    
2171        return or(ONE.shiftLeft(n));
2172      }
2173    
2174      public boolean testBit(int n)
2175      {
2176        if (n < 0)
2177          throw new ArithmeticException();
2178    
2179        return !and(ONE.shiftLeft(n)).isZero();
2180      }
2181    
2182      public BigInteger flipBit(int n)
2183      {
2184        if (n < 0)
2185          throw new ArithmeticException();
2186    
2187        return xor(ONE.shiftLeft(n));
2188      }
2189    
2190      public int getLowestSetBit()
2191      {
2192        if (isZero())
2193          return -1;
2194    
2195        if (words == null)
2196          return MPN.findLowestBit(ival);
2197        else
2198          return MPN.findLowestBit(words);
2199      }
2200    
2201      // bit4count[I] is number of '1' bits in I.
2202      private static final byte[] bit4_count = { 0, 1, 1, 2,  1, 2, 2, 3,
2203                                                 1, 2, 2, 3,  2, 3, 3, 4};
2204    
2205      private static int bitCount(int i)
2206      {
2207        int count = 0;
2208        while (i != 0)
2209          {
2210            count += bit4_count[i & 15];
2211            i >>>= 4;
2212          }
2213        return count;
2214      }
2215    
2216      private static int bitCount(int[] x, int len)
2217      {
2218        int count = 0;
2219        while (--len >= 0)
2220          count += bitCount(x[len]);
2221        return count;
2222      }
2223    
2224      /** Count one bits in a BigInteger.
2225       * If argument is negative, count zero bits instead. */
2226      public int bitCount()
2227      {
2228        int i, x_len;
2229        int[] x_words = words;
2230        if (x_words == null)
2231          {
2232            x_len = 1;
2233            i = bitCount(ival);
2234          }
2235        else
2236          {
2237            x_len = ival;
2238            i = bitCount(x_words, x_len);
2239          }
2240        return isNegative() ? x_len * 32 - i : i;
2241      }
2242    
2243      private void readObject(ObjectInputStream s)
2244        throws IOException, ClassNotFoundException
2245      {
2246        s.defaultReadObject();
2247        if (magnitude.length == 0 || signum == 0)
2248          {
2249            this.ival = 0;
2250            this.words = null;
2251          }
2252        else
2253          {
2254            words = byteArrayToIntArray(magnitude, signum < 0 ? -1 : 0);
2255            BigInteger result = make(words, words.length);
2256            this.ival = result.ival;
2257            this.words = result.words;        
2258          }    
2259      }
2260    
2261      private void writeObject(ObjectOutputStream s)
2262        throws IOException, ClassNotFoundException
2263      {
2264        signum = signum();
2265        magnitude = signum == 0 ? new byte[0] : toByteArray();
2266        s.defaultWriteObject();
2267      }
2268    }