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