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 }