1 /*
2 * Copyright (c) 1996, 2014, Oracle and/or its affiliates. All rights reserved.
3 * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
4 *
5 * This code is free software; you can redistribute it and/or modify it
6 * under the terms of the GNU General Public License version 2 only, as
7 * published by the Free Software Foundation. Oracle designates this
8 * particular file as subject to the "Classpath" exception as provided
9 * by Oracle in the LICENSE file that accompanied this code.
10 *
11 * This code is distributed in the hope that it will be useful, but WITHOUT
12 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14 * version 2 for more details (a copy is included in the LICENSE file that
15 * accompanied this code).
16 *
17 * You should have received a copy of the GNU General Public License version
18 * 2 along with this work; if not, write to the Free Software Foundation,
19 * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
20 *
21 * Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA
22 * or visit www.oracle.com if you need additional information or have any
23 * questions.
24 */
25
26 /*
27 * Portions Copyright (c) 1995 Colin Plumb. All rights reserved.
28 */
29
30 package java.math;
31
32 import java.io.IOException;
33 import java.io.ObjectInputStream;
34 import java.io.ObjectOutputStream;
35 import java.io.ObjectStreamField;
36 import java.util.Arrays;
37 import java.util.Random;
38 import java.util.concurrent.ThreadLocalRandom;
39 import sun.misc.DoubleConsts;
40 import sun.misc.FloatConsts;
41
42 /**
43 * Immutable arbitrary-precision integers. All operations behave as if
44 * BigIntegers were represented in two's-complement notation (like Java's
45 * primitive integer types). BigInteger provides analogues to all of Java's
46 * primitive integer operators, and all relevant methods from java.lang.Math.
47 * Additionally, BigInteger provides operations for modular arithmetic, GCD
48 * calculation, primality testing, prime generation, bit manipulation,
49 * and a few other miscellaneous operations.
50 *
51 * <p>Semantics of arithmetic operations exactly mimic those of Java's integer
52 * arithmetic operators, as defined in <i>The Java Language Specification</i>.
53 * For example, division by zero throws an {@code ArithmeticException}, and
54 * division of a negative by a positive yields a negative (or zero) remainder.
55 * All of the details in the Spec concerning overflow are ignored, as
56 * BigIntegers are made as large as necessary to accommodate the results of an
57 * operation.
58 *
59 * <p>Semantics of shift operations extend those of Java's shift operators
60 * to allow for negative shift distances. A right-shift with a negative
61 * shift distance results in a left shift, and vice-versa. The unsigned
62 * right shift operator ({@code >>>}) is omitted, as this operation makes
63 * little sense in combination with the "infinite word size" abstraction
64 * provided by this class.
65 *
66 * <p>Semantics of bitwise logical operations exactly mimic those of Java's
67 * bitwise integer operators. The binary operators ({@code and},
68 * {@code or}, {@code xor}) implicitly perform sign extension on the shorter
69 * of the two operands prior to performing the operation.
70 *
71 * <p>Comparison operations perform signed integer comparisons, analogous to
72 * those performed by Java's relational and equality operators.
73 *
74 * <p>Modular arithmetic operations are provided to compute residues, perform
75 * exponentiation, and compute multiplicative inverses. These methods always
76 * return a non-negative result, between {@code 0} and {@code (modulus - 1)},
77 * inclusive.
78 *
79 * <p>Bit operations operate on a single bit of the two's-complement
80 * representation of their operand. If necessary, the operand is sign-
81 * extended so that it contains the designated bit. None of the single-bit
82 * operations can produce a BigInteger with a different sign from the
83 * BigInteger being operated on, as they affect only a single bit, and the
84 * "infinite word size" abstraction provided by this class ensures that there
85 * are infinitely many "virtual sign bits" preceding each BigInteger.
86 *
87 * <p>For the sake of brevity and clarity, pseudo-code is used throughout the
88 * descriptions of BigInteger methods. The pseudo-code expression
89 * {@code (i + j)} is shorthand for "a BigInteger whose value is
90 * that of the BigInteger {@code i} plus that of the BigInteger {@code j}."
91 * The pseudo-code expression {@code (i == j)} is shorthand for
92 * "{@code true} if and only if the BigInteger {@code i} represents the same
93 * value as the BigInteger {@code j}." Other pseudo-code expressions are
94 * interpreted similarly.
95 *
96 * <p>All methods and constructors in this class throw
97 * {@code NullPointerException} when passed
98 * a null object reference for any input parameter.
99 *
100 * BigInteger must support values in the range
101 * -2<sup>{@code Integer.MAX_VALUE}</sup> (exclusive) to
102 * +2<sup>{@code Integer.MAX_VALUE}</sup> (exclusive)
103 * and may support values outside of that range.
104 *
105 * The range of probable prime values is limited and may be less than
106 * the full supported positive range of {@code BigInteger}.
107 * The range must be at least 1 to 2<sup>500000000</sup>.
108 *
109 * @implNote
110 * BigInteger constructors and operations throw {@code ArithmeticException} when
111 * the result is out of the supported range of
112 * -2<sup>{@code Integer.MAX_VALUE}</sup> (exclusive) to
113 * +2<sup>{@code Integer.MAX_VALUE}</sup> (exclusive).
114 *
115 * @see BigDecimal
116 * @author Josh Bloch
117 * @author Michael McCloskey
118 * @author Alan Eliasen
119 * @author Timothy Buktu
120 * @since 1.1
121 */
122
123 public class BigInteger extends Number implements Comparable<BigInteger> {
124 /**
125 * The signum of this BigInteger: -1 for negative, 0 for zero, or
126 * 1 for positive. Note that the BigInteger zero <i>must</i> have
127 * a signum of 0. This is necessary to ensures that there is exactly one
128 * representation for each BigInteger value.
129 */
130 final int signum;
131
132 /**
133 * The magnitude of this BigInteger, in <i>big-endian</i> order: the
134 * zeroth element of this array is the most-significant int of the
135 * magnitude. The magnitude must be "minimal" in that the most-significant
136 * int ({@code mag[0]}) must be non-zero. This is necessary to
137 * ensure that there is exactly one representation for each BigInteger
138 * value. Note that this implies that the BigInteger zero has a
139 * zero-length mag array.
140 */
141 final int[] mag;
142
143 // The following fields are stable variables. A stable variable's value
144 // changes at most once from the default zero value to a non-zero stable
145 // value. A stable value is calculated lazily on demand.
146
147 /**
148 * One plus the bitCount of this BigInteger. This is a stable variable.
149 *
150 * @see #bitCount
151 */
152 private int bitCountPlusOne;
153
154 /**
155 * One plus the bitLength of this BigInteger. This is a stable variable.
156 * (either value is acceptable).
157 *
158 * @see #bitLength()
159 */
160 private int bitLengthPlusOne;
161
162 /**
163 * Two plus the lowest set bit of this BigInteger. This is a stable variable.
164 *
165 * @see #getLowestSetBit
166 */
167 private int lowestSetBitPlusTwo;
168
169 /**
170 * Two plus the index of the lowest-order int in the magnitude of this
171 * BigInteger that contains a nonzero int. This is a stable variable. The
172 * least significant int has int-number 0, the next int in order of
173 * increasing significance has int-number 1, and so forth.
174 *
175 * <p>Note: never used for a BigInteger with a magnitude of zero.
176 *
177 * @see #firstNonzeroIntNum()
178 */
179 private int firstNonzeroIntNumPlusTwo;
180
181 /**
182 * This mask is used to obtain the value of an int as if it were unsigned.
183 */
184 final static long LONG_MASK = 0xffffffffL;
185
186 /**
187 * This constant limits {@code mag.length} of BigIntegers to the supported
188 * range.
189 */
190 private static final int MAX_MAG_LENGTH = Integer.MAX_VALUE / Integer.SIZE + 1; // (1 << 26)
191
192 /**
193 * Bit lengths larger than this constant can cause overflow in searchLen
194 * calculation and in BitSieve.singleSearch method.
195 */
196 private static final int PRIME_SEARCH_BIT_LENGTH_LIMIT = 500000000;
197
198 /**
199 * The threshold value for using Karatsuba multiplication. If the number
200 * of ints in both mag arrays are greater than this number, then
201 * Karatsuba multiplication will be used. This value is found
202 * experimentally to work well.
203 */
204 private static final int KARATSUBA_THRESHOLD = 80;
205
206 /**
207 * The threshold value for using 3-way Toom-Cook multiplication.
208 * If the number of ints in each mag array is greater than the
209 * Karatsuba threshold, and the number of ints in at least one of
210 * the mag arrays is greater than this threshold, then Toom-Cook
211 * multiplication will be used.
212 */
213 private static final int TOOM_COOK_THRESHOLD = 240;
214
215 /**
216 * The threshold value for using Karatsuba squaring. If the number
217 * of ints in the number are larger than this value,
218 * Karatsuba squaring will be used. This value is found
219 * experimentally to work well.
220 */
221 private static final int KARATSUBA_SQUARE_THRESHOLD = 128;
222
223 /**
224 * The threshold value for using Toom-Cook squaring. If the number
225 * of ints in the number are larger than this value,
226 * Toom-Cook squaring will be used. This value is found
227 * experimentally to work well.
228 */
229 private static final int TOOM_COOK_SQUARE_THRESHOLD = 216;
230
231 /**
232 * The threshold value for using Burnikel-Ziegler division. If the number
233 * of ints in the divisor are larger than this value, Burnikel-Ziegler
234 * division may be used. This value is found experimentally to work well.
235 */
236 static final int BURNIKEL_ZIEGLER_THRESHOLD = 80;
237
238 /**
239 * The offset value for using Burnikel-Ziegler division. If the number
240 * of ints in the divisor exceeds the Burnikel-Ziegler threshold, and the
241 * number of ints in the dividend is greater than the number of ints in the
242 * divisor plus this value, Burnikel-Ziegler division will be used. This
243 * value is found experimentally to work well.
244 */
245 static final int BURNIKEL_ZIEGLER_OFFSET = 40;
246
247 /**
248 * The threshold value for using Schoenhage recursive base conversion. If
249 * the number of ints in the number are larger than this value,
250 * the Schoenhage algorithm will be used. In practice, it appears that the
251 * Schoenhage routine is faster for any threshold down to 2, and is
252 * relatively flat for thresholds between 2-25, so this choice may be
253 * varied within this range for very small effect.
254 */
255 private static final int SCHOENHAGE_BASE_CONVERSION_THRESHOLD = 20;
256
257 /**
258 * The threshold value for using squaring code to perform multiplication
259 * of a {@code BigInteger} instance by itself. If the number of ints in
260 * the number are larger than this value, {@code multiply(this)} will
261 * return {@code square()}.
262 */
263 private static final int MULTIPLY_SQUARE_THRESHOLD = 20;
264
265 // Constructors
266
267 /**
268 * Translates a byte sub-array containing the two's-complement binary
269 * representation of a BigInteger into a BigInteger. The sub-array is
270 * specified via an offset into the array and a length. The sub-array is
271 * assumed to be in <i>big-endian</i> byte-order: the most significant
272 * byte is the element at index {@code off}.
273 *
274 * An {@code IndexOutOfBoundsException} is thrown if the length of the array
275 * {@code val} is non-zero and either {@code off} is negative, {@code len}
276 * is negative, or {@code off+len} is greater than the length of
277 * {@code val}.
278 *
279 * @param val byte array containing a sub-array which is the big-endian
280 * two's-complement binary representation of a BigInteger.
281 * @param off the start offset of the binary representation.
282 * @param len the number of bytes to use.
283 * @throws NumberFormatException {@code val} is zero bytes long.
284 */
285 public BigInteger(byte[] val, int off, int len) {
286 if (val.length == 0) {
287 throw new NumberFormatException("Zero length BigInteger");
288 } else if ((off < 0) || (off > val.length) || (len < 0) ||
289 ((off + len) > val.length) || ((off + len) < 0)) {
290 throw new IndexOutOfBoundsException();
291 }
292
293 if (val[off] < 0) {
294 mag = makePositive(val, off, len);
295 signum = -1;
296 } else {
297 mag = stripLeadingZeroBytes(val, off, len);
298 signum = (mag.length == 0 ? 0 : 1);
299 }
300 if (mag.length >= MAX_MAG_LENGTH) {
301 checkRange();
302 }
303 }
304
305 /**
306 * Translates a byte array containing the two's-complement binary
307 * representation of a BigInteger into a BigInteger. The input array is
308 * assumed to be in <i>big-endian</i> byte-order: the most significant
309 * byte is in the zeroth element.
310 *
311 * @param val big-endian two's-complement binary representation of a
312 * BigInteger.
313 * @throws NumberFormatException {@code val} is zero bytes long.
314 */
315 public BigInteger(byte[] val) {
316 this(val, 0, val.length);
317 }
318
319 /**
320 * This private constructor translates an int array containing the
321 * two's-complement binary representation of a BigInteger into a
322 * BigInteger. The input array is assumed to be in <i>big-endian</i>
323 * int-order: the most significant int is in the zeroth element.
324 */
325 private BigInteger(int[] val) {
326 if (val.length == 0)
327 throw new NumberFormatException("Zero length BigInteger");
328
329 if (val[0] < 0) {
330 mag = makePositive(val);
331 signum = -1;
332 } else {
333 mag = trustedStripLeadingZeroInts(val);
334 signum = (mag.length == 0 ? 0 : 1);
335 }
336 if (mag.length >= MAX_MAG_LENGTH) {
337 checkRange();
338 }
339 }
340
341 /**
342 * Translates the sign-magnitude representation of a BigInteger into a
343 * BigInteger. The sign is represented as an integer signum value: -1 for
344 * negative, 0 for zero, or 1 for positive. The magnitude is a sub-array of
345 * a byte array in <i>big-endian</i> byte-order: the most significant byte
346 * is the element at index {@code off}. A zero value of the length
347 * {@code len} is permissible, and will result in a BigInteger value of 0,
348 * whether signum is -1, 0 or 1.
349 *
350 * An {@code IndexOutOfBoundsException} is thrown if the length of the array
351 * {@code magnitude} is non-zero and either {@code off} is negative,
352 * {@code len} is negative, or {@code off+len} is greater than the length of
353 * {@code magnitude}.
354 *
355 * @param signum signum of the number (-1 for negative, 0 for zero, 1
356 * for positive).
357 * @param magnitude big-endian binary representation of the magnitude of
358 * the number.
359 * @throws NumberFormatException {@code signum} is not one of the three
360 * legal values (-1, 0, and 1), or {@code signum} is 0 and
361 * {@code magnitude} contains one or more non-zero bytes.
362 */
363 public BigInteger(int signum, byte[] magnitude, int off, int len) {
364 if ((off < 0) || (off > magnitude.length) || (len < 0)
365 || ((off + len) > magnitude.length) || ((off + len) < 0)) {
366 throw new IndexOutOfBoundsException();
367 }
368 this.mag = stripLeadingZeroBytes(magnitude, off, len);
369
370 if (signum < -1 || signum > 1)
371 throw(new NumberFormatException("Invalid signum value"));
372
373 if (this.mag.length == 0) {
374 this.signum = 0;
375 } else {
376 if (signum == 0)
377 throw(new NumberFormatException("signum-magnitude mismatch"));
378 this.signum = signum;
379 }
380 if (mag.length >= MAX_MAG_LENGTH) {
381 checkRange();
382 }
383 }
384
385 /**
386 * Translates the sign-magnitude representation of a BigInteger into a
387 * BigInteger. The sign is represented as an integer signum value: -1 for
388 * negative, 0 for zero, or 1 for positive. The magnitude is a byte array
389 * in <i>big-endian</i> byte-order: the most significant byte is the
390 * zeroth element. A zero-length magnitude array is permissible, and will
391 * result in a BigInteger value of 0, whether signum is -1, 0 or 1.
392 *
393 * @param signum signum of the number (-1 for negative, 0 for zero, 1
394 * for positive).
395 * @param magnitude big-endian binary representation of the magnitude of
396 * the number.
397 * @throws NumberFormatException {@code signum} is not one of the three
398 * legal values (-1, 0, and 1), or {@code signum} is 0 and
399 * {@code magnitude} contains one or more non-zero bytes.
400 */
401 public BigInteger(int signum, byte[] magnitude) {
402 this(signum, magnitude, 0, magnitude.length);
403 }
404
405 /**
406 * A constructor for internal use that translates the sign-magnitude
407 * representation of a BigInteger into a BigInteger. It checks the
408 * arguments and copies the magnitude so this constructor would be
409 * safe for external use.
410 */
411 private BigInteger(int signum, int[] magnitude) {
412 this.mag = stripLeadingZeroInts(magnitude);
413
414 if (signum < -1 || signum > 1)
415 throw(new NumberFormatException("Invalid signum value"));
416
417 if (this.mag.length == 0) {
418 this.signum = 0;
419 } else {
420 if (signum == 0)
421 throw(new NumberFormatException("signum-magnitude mismatch"));
422 this.signum = signum;
423 }
424 if (mag.length >= MAX_MAG_LENGTH) {
425 checkRange();
426 }
427 }
428
429 /**
430 * Translates the String representation of a BigInteger in the
431 * specified radix into a BigInteger. The String representation
432 * consists of an optional minus or plus sign followed by a
433 * sequence of one or more digits in the specified radix. The
434 * character-to-digit mapping is provided by {@code
435 * Character.digit}. The String may not contain any extraneous
436 * characters (whitespace, for example).
437 *
438 * @param val String representation of BigInteger.
439 * @param radix radix to be used in interpreting {@code val}.
440 * @throws NumberFormatException {@code val} is not a valid representation
441 * of a BigInteger in the specified radix, or {@code radix} is
442 * outside the range from {@link Character#MIN_RADIX} to
443 * {@link Character#MAX_RADIX}, inclusive.
444 * @see Character#digit
445 */
446 public BigInteger(String val, int radix) {
447 int cursor = 0, numDigits;
448 final int len = val.length();
449
450 if (radix < Character.MIN_RADIX || radix > Character.MAX_RADIX)
451 throw new NumberFormatException("Radix out of range");
452 if (len == 0)
453 throw new NumberFormatException("Zero length BigInteger");
454
455 // Check for at most one leading sign
456 int sign = 1;
457 int index1 = val.lastIndexOf('-');
458 int index2 = val.lastIndexOf('+');
459 if (index1 >= 0) {
460 if (index1 != 0 || index2 >= 0) {
461 throw new NumberFormatException("Illegal embedded sign character");
462 }
463 sign = -1;
464 cursor = 1;
465 } else if (index2 >= 0) {
466 if (index2 != 0) {
467 throw new NumberFormatException("Illegal embedded sign character");
468 }
469 cursor = 1;
470 }
471 if (cursor == len)
472 throw new NumberFormatException("Zero length BigInteger");
473
474 // Skip leading zeros and compute number of digits in magnitude
475 while (cursor < len &&
476 Character.digit(val.charAt(cursor), radix) == 0) {
477 cursor++;
478 }
479
480 if (cursor == len) {
481 signum = 0;
482 mag = ZERO.mag;
483 return;
484 }
485
486 numDigits = len - cursor;
487 signum = sign;
488
489 // Pre-allocate array of expected size. May be too large but can
490 // never be too small. Typically exact.
491 long numBits = ((numDigits * bitsPerDigit[radix]) >>> 10) + 1;
492 if (numBits + 31 >= (1L << 32)) {
493 reportOverflow();
494 }
495 int numWords = (int) (numBits + 31) >>> 5;
496 int[] magnitude = new int[numWords];
497
498 // Process first (potentially short) digit group
499 int firstGroupLen = numDigits % digitsPerInt[radix];
500 if (firstGroupLen == 0)
501 firstGroupLen = digitsPerInt[radix];
502 String group = val.substring(cursor, cursor += firstGroupLen);
503 magnitude[numWords - 1] = Integer.parseInt(group, radix);
504 if (magnitude[numWords - 1] < 0)
505 throw new NumberFormatException("Illegal digit");
506
507 // Process remaining digit groups
508 int superRadix = intRadix[radix];
509 int groupVal = 0;
510 while (cursor < len) {
511 group = val.substring(cursor, cursor += digitsPerInt[radix]);
512 groupVal = Integer.parseInt(group, radix);
513 if (groupVal < 0)
514 throw new NumberFormatException("Illegal digit");
515 destructiveMulAdd(magnitude, superRadix, groupVal);
516 }
517 // Required for cases where the array was overallocated.
518 mag = trustedStripLeadingZeroInts(magnitude);
519 if (mag.length >= MAX_MAG_LENGTH) {
520 checkRange();
521 }
522 }
523
524 /*
525 * Constructs a new BigInteger using a char array with radix=10.
526 * Sign is precalculated outside and not allowed in the val.
527 */
528 BigInteger(char[] val, int sign, int len) {
529 int cursor = 0, numDigits;
530
531 // Skip leading zeros and compute number of digits in magnitude
532 while (cursor < len && Character.digit(val[cursor], 10) == 0) {
533 cursor++;
534 }
535 if (cursor == len) {
536 signum = 0;
537 mag = ZERO.mag;
538 return;
539 }
540
541 numDigits = len - cursor;
542 signum = sign;
543 // Pre-allocate array of expected size
544 int numWords;
545 if (len < 10) {
546 numWords = 1;
547 } else {
548 long numBits = ((numDigits * bitsPerDigit[10]) >>> 10) + 1;
549 if (numBits + 31 >= (1L << 32)) {
550 reportOverflow();
551 }
552 numWords = (int) (numBits + 31) >>> 5;
553 }
554 int[] magnitude = new int[numWords];
555
556 // Process first (potentially short) digit group
557 int firstGroupLen = numDigits % digitsPerInt[10];
558 if (firstGroupLen == 0)
559 firstGroupLen = digitsPerInt[10];
560 magnitude[numWords - 1] = parseInt(val, cursor, cursor += firstGroupLen);
561
562 // Process remaining digit groups
563 while (cursor < len) {
564 int groupVal = parseInt(val, cursor, cursor += digitsPerInt[10]);
565 destructiveMulAdd(magnitude, intRadix[10], groupVal);
566 }
567 mag = trustedStripLeadingZeroInts(magnitude);
568 if (mag.length >= MAX_MAG_LENGTH) {
569 checkRange();
570 }
571 }
572
573 // Create an integer with the digits between the two indexes
574 // Assumes start < end. The result may be negative, but it
575 // is to be treated as an unsigned value.
576 private int parseInt(char[] source, int start, int end) {
577 int result = Character.digit(source[start++], 10);
578 if (result == -1)
579 throw new NumberFormatException(new String(source));
580
581 for (int index = start; index < end; index++) {
582 int nextVal = Character.digit(source[index], 10);
583 if (nextVal == -1)
584 throw new NumberFormatException(new String(source));
585 result = 10*result + nextVal;
586 }
587
588 return result;
589 }
590
591 // bitsPerDigit in the given radix times 1024
592 // Rounded up to avoid underallocation.
593 private static long bitsPerDigit[] = { 0, 0,
594 1024, 1624, 2048, 2378, 2648, 2875, 3072, 3247, 3402, 3543, 3672,
595 3790, 3899, 4001, 4096, 4186, 4271, 4350, 4426, 4498, 4567, 4633,
596 4696, 4756, 4814, 4870, 4923, 4975, 5025, 5074, 5120, 5166, 5210,
597 5253, 5295};
598
599 // Multiply x array times word y in place, and add word z
600 private static void destructiveMulAdd(int[] x, int y, int z) {
601 // Perform the multiplication word by word
602 long ylong = y & LONG_MASK;
603 long zlong = z & LONG_MASK;
604 int len = x.length;
605
606 long product = 0;
607 long carry = 0;
608 for (int i = len-1; i >= 0; i--) {
609 product = ylong * (x[i] & LONG_MASK) + carry;
610 x[i] = (int)product;
611 carry = product >>> 32;
612 }
613
614 // Perform the addition
615 long sum = (x[len-1] & LONG_MASK) + zlong;
616 x[len-1] = (int)sum;
617 carry = sum >>> 32;
618 for (int i = len-2; i >= 0; i--) {
619 sum = (x[i] & LONG_MASK) + carry;
620 x[i] = (int)sum;
621 carry = sum >>> 32;
622 }
623 }
624
625 /**
626 * Translates the decimal String representation of a BigInteger into a
627 * BigInteger. The String representation consists of an optional minus
628 * sign followed by a sequence of one or more decimal digits. The
629 * character-to-digit mapping is provided by {@code Character.digit}.
630 * The String may not contain any extraneous characters (whitespace, for
631 * example).
632 *
633 * @param val decimal String representation of BigInteger.
634 * @throws NumberFormatException {@code val} is not a valid representation
635 * of a BigInteger.
636 * @see Character#digit
637 */
638 public BigInteger(String val) {
639 this(val, 10);
640 }
641
642 /**
643 * Constructs a randomly generated BigInteger, uniformly distributed over
644 * the range 0 to (2<sup>{@code numBits}</sup> - 1), inclusive.
645 * The uniformity of the distribution assumes that a fair source of random
646 * bits is provided in {@code rnd}. Note that this constructor always
647 * constructs a non-negative BigInteger.
648 *
649 * @param numBits maximum bitLength of the new BigInteger.
650 * @param rnd source of randomness to be used in computing the new
651 * BigInteger.
652 * @throws IllegalArgumentException {@code numBits} is negative.
653 * @see #bitLength()
654 */
655 public BigInteger(int numBits, Random rnd) {
656 this(1, randomBits(numBits, rnd));
657 }
658
659 private static byte[] randomBits(int numBits, Random rnd) {
660 if (numBits < 0)
661 throw new IllegalArgumentException("numBits must be non-negative");
662 int numBytes = (int)(((long)numBits+7)/8); // avoid overflow
663 byte[] randomBits = new byte[numBytes];
664
665 // Generate random bytes and mask out any excess bits
666 if (numBytes > 0) {
667 rnd.nextBytes(randomBits);
668 int excessBits = 8*numBytes - numBits;
669 randomBits[0] &= (1 << (8-excessBits)) - 1;
670 }
671 return randomBits;
672 }
673
674 /**
675 * Constructs a randomly generated positive BigInteger that is probably
676 * prime, with the specified bitLength.
677 *
678 * <p>It is recommended that the {@link #probablePrime probablePrime}
679 * method be used in preference to this constructor unless there
680 * is a compelling need to specify a certainty.
681 *
682 * @param bitLength bitLength of the returned BigInteger.
683 * @param certainty a measure of the uncertainty that the caller is
684 * willing to tolerate. The probability that the new BigInteger
685 * represents a prime number will exceed
686 * (1 - 1/2<sup>{@code certainty}</sup>). The execution time of
687 * this constructor is proportional to the value of this parameter.
688 * @param rnd source of random bits used to select candidates to be
689 * tested for primality.
690 * @throws ArithmeticException {@code bitLength < 2} or {@code bitLength} is too large.
691 * @see #bitLength()
692 */
693 public BigInteger(int bitLength, int certainty, Random rnd) {
694 BigInteger prime;
695
696 if (bitLength < 2)
697 throw new ArithmeticException("bitLength < 2");
698 prime = (bitLength < SMALL_PRIME_THRESHOLD
699 ? smallPrime(bitLength, certainty, rnd)
700 : largePrime(bitLength, certainty, rnd));
701 signum = 1;
702 mag = prime.mag;
703 }
704
705 // Minimum size in bits that the requested prime number has
706 // before we use the large prime number generating algorithms.
707 // The cutoff of 95 was chosen empirically for best performance.
708 private static final int SMALL_PRIME_THRESHOLD = 95;
709
710 // Certainty required to meet the spec of probablePrime
711 private static final int DEFAULT_PRIME_CERTAINTY = 100;
712
713 /**
714 * Returns a positive BigInteger that is probably prime, with the
715 * specified bitLength. The probability that a BigInteger returned
716 * by this method is composite does not exceed 2<sup>-100</sup>.
717 *
718 * @param bitLength bitLength of the returned BigInteger.
719 * @param rnd source of random bits used to select candidates to be
720 * tested for primality.
721 * @return a BigInteger of {@code bitLength} bits that is probably prime
722 * @throws ArithmeticException {@code bitLength < 2} or {@code bitLength} is too large.
723 * @see #bitLength()
724 * @since 1.4
725 */
726 public static BigInteger probablePrime(int bitLength, Random rnd) {
727 if (bitLength < 2)
728 throw new ArithmeticException("bitLength < 2");
729
730 return (bitLength < SMALL_PRIME_THRESHOLD ?
731 smallPrime(bitLength, DEFAULT_PRIME_CERTAINTY, rnd) :
732 largePrime(bitLength, DEFAULT_PRIME_CERTAINTY, rnd));
733 }
734
735 /**
736 * Find a random number of the specified bitLength that is probably prime.
737 * This method is used for smaller primes, its performance degrades on
738 * larger bitlengths.
739 *
740 * This method assumes bitLength > 1.
741 */
742 private static BigInteger smallPrime(int bitLength, int certainty, Random rnd) {
743 int magLen = (bitLength + 31) >>> 5;
744 int temp[] = new int[magLen];
745 int highBit = 1 << ((bitLength+31) & 0x1f); // High bit of high int
746 int highMask = (highBit << 1) - 1; // Bits to keep in high int
747
748 while (true) {
749 // Construct a candidate
750 for (int i=0; i < magLen; i++)
751 temp[i] = rnd.nextInt();
752 temp[0] = (temp[0] & highMask) | highBit; // Ensure exact length
753 if (bitLength > 2)
754 temp[magLen-1] |= 1; // Make odd if bitlen > 2
755
756 BigInteger p = new BigInteger(temp, 1);
757
758 // Do cheap "pre-test" if applicable
759 if (bitLength > 6) {
760 long r = p.remainder(SMALL_PRIME_PRODUCT).longValue();
761 if ((r%3==0) || (r%5==0) || (r%7==0) || (r%11==0) ||
762 (r%13==0) || (r%17==0) || (r%19==0) || (r%23==0) ||
763 (r%29==0) || (r%31==0) || (r%37==0) || (r%41==0))
764 continue; // Candidate is composite; try another
765 }
766
767 // All candidates of bitLength 2 and 3 are prime by this point
768 if (bitLength < 4)
769 return p;
770
771 // Do expensive test if we survive pre-test (or it's inapplicable)
772 if (p.primeToCertainty(certainty, rnd))
773 return p;
774 }
775 }
776
777 private static final BigInteger SMALL_PRIME_PRODUCT
778 = valueOf(3L*5*7*11*13*17*19*23*29*31*37*41);
779
780 /**
781 * Find a random number of the specified bitLength that is probably prime.
782 * This method is more appropriate for larger bitlengths since it uses
783 * a sieve to eliminate most composites before using a more expensive
784 * test.
785 */
786 private static BigInteger largePrime(int bitLength, int certainty, Random rnd) {
787 BigInteger p;
788 p = new BigInteger(bitLength, rnd).setBit(bitLength-1);
789 p.mag[p.mag.length-1] &= 0xfffffffe;
790
791 // Use a sieve length likely to contain the next prime number
792 int searchLen = getPrimeSearchLen(bitLength);
793 BitSieve searchSieve = new BitSieve(p, searchLen);
794 BigInteger candidate = searchSieve.retrieve(p, certainty, rnd);
795
796 while ((candidate == null) || (candidate.bitLength() != bitLength)) {
797 p = p.add(BigInteger.valueOf(2*searchLen));
798 if (p.bitLength() != bitLength)
799 p = new BigInteger(bitLength, rnd).setBit(bitLength-1);
800 p.mag[p.mag.length-1] &= 0xfffffffe;
801 searchSieve = new BitSieve(p, searchLen);
802 candidate = searchSieve.retrieve(p, certainty, rnd);
803 }
804 return candidate;
805 }
806
807 /**
808 * Returns the first integer greater than this {@code BigInteger} that
809 * is probably prime. The probability that the number returned by this
810 * method is composite does not exceed 2<sup>-100</sup>. This method will
811 * never skip over a prime when searching: if it returns {@code p}, there
812 * is no prime {@code q} such that {@code this < q < p}.
813 *
814 * @return the first integer greater than this {@code BigInteger} that
815 * is probably prime.
816 * @throws ArithmeticException {@code this < 0} or {@code this} is too large.
817 * @since 1.5
818 */
819 public BigInteger nextProbablePrime() {
820 if (this.signum < 0)
821 throw new ArithmeticException("start < 0: " + this);
822
823 // Handle trivial cases
824 if ((this.signum == 0) || this.equals(ONE))
825 return TWO;
826
827 BigInteger result = this.add(ONE);
828
829 // Fastpath for small numbers
830 if (result.bitLength() < SMALL_PRIME_THRESHOLD) {
831
832 // Ensure an odd number
833 if (!result.testBit(0))
834 result = result.add(ONE);
835
836 while (true) {
837 // Do cheap "pre-test" if applicable
838 if (result.bitLength() > 6) {
839 long r = result.remainder(SMALL_PRIME_PRODUCT).longValue();
840 if ((r%3==0) || (r%5==0) || (r%7==0) || (r%11==0) ||
841 (r%13==0) || (r%17==0) || (r%19==0) || (r%23==0) ||
842 (r%29==0) || (r%31==0) || (r%37==0) || (r%41==0)) {
843 result = result.add(TWO);
844 continue; // Candidate is composite; try another
845 }
846 }
847
848 // All candidates of bitLength 2 and 3 are prime by this point
849 if (result.bitLength() < 4)
850 return result;
851
852 // The expensive test
853 if (result.primeToCertainty(DEFAULT_PRIME_CERTAINTY, null))
854 return result;
855
856 result = result.add(TWO);
857 }
858 }
859
860 // Start at previous even number
861 if (result.testBit(0))
862 result = result.subtract(ONE);
863
864 // Looking for the next large prime
865 int searchLen = getPrimeSearchLen(result.bitLength());
866
867 while (true) {
868 BitSieve searchSieve = new BitSieve(result, searchLen);
869 BigInteger candidate = searchSieve.retrieve(result,
870 DEFAULT_PRIME_CERTAINTY, null);
871 if (candidate != null)
872 return candidate;
873 result = result.add(BigInteger.valueOf(2 * searchLen));
874 }
875 }
876
877 private static int getPrimeSearchLen(int bitLength) {
878 if (bitLength > PRIME_SEARCH_BIT_LENGTH_LIMIT + 1) {
879 throw new ArithmeticException("Prime search implementation restriction on bitLength");
880 }
881 return bitLength / 20 * 64;
882 }
883
884 /**
885 * Returns {@code true} if this BigInteger is probably prime,
886 * {@code false} if it's definitely composite.
887 *
888 * This method assumes bitLength > 2.
889 *
890 * @param certainty a measure of the uncertainty that the caller is
891 * willing to tolerate: if the call returns {@code true}
892 * the probability that this BigInteger is prime exceeds
893 * {@code (1 - 1/2<sup>certainty</sup>)}. The execution time of
894 * this method is proportional to the value of this parameter.
895 * @return {@code true} if this BigInteger is probably prime,
896 * {@code false} if it's definitely composite.
897 */
898 boolean primeToCertainty(int certainty, Random random) {
899 int rounds = 0;
900 int n = (Math.min(certainty, Integer.MAX_VALUE-1)+1)/2;
901
902 // The relationship between the certainty and the number of rounds
903 // we perform is given in the draft standard ANSI X9.80, "PRIME
904 // NUMBER GENERATION, PRIMALITY TESTING, AND PRIMALITY CERTIFICATES".
905 int sizeInBits = this.bitLength();
906 if (sizeInBits < 100) {
907 rounds = 50;
908 rounds = n < rounds ? n : rounds;
909 return passesMillerRabin(rounds, random);
910 }
911
912 if (sizeInBits < 256) {
913 rounds = 27;
914 } else if (sizeInBits < 512) {
915 rounds = 15;
916 } else if (sizeInBits < 768) {
917 rounds = 8;
918 } else if (sizeInBits < 1024) {
919 rounds = 4;
920 } else {
921 rounds = 2;
922 }
923 rounds = n < rounds ? n : rounds;
924
925 return passesMillerRabin(rounds, random) && passesLucasLehmer();
926 }
927
928 /**
929 * Returns true iff this BigInteger is a Lucas-Lehmer probable prime.
930 *
931 * The following assumptions are made:
932 * This BigInteger is a positive, odd number.
933 */
934 private boolean passesLucasLehmer() {
935 BigInteger thisPlusOne = this.add(ONE);
936
937 // Step 1
938 int d = 5;
939 while (jacobiSymbol(d, this) != -1) {
940 // 5, -7, 9, -11, ...
941 d = (d < 0) ? Math.abs(d)+2 : -(d+2);
942 }
943
944 // Step 2
945 BigInteger u = lucasLehmerSequence(d, thisPlusOne, this);
946
947 // Step 3
948 return u.mod(this).equals(ZERO);
949 }
950
951 /**
952 * Computes Jacobi(p,n).
953 * Assumes n positive, odd, n>=3.
954 */
955 private static int jacobiSymbol(int p, BigInteger n) {
956 if (p == 0)
957 return 0;
958
959 // Algorithm and comments adapted from Colin Plumb's C library.
960 int j = 1;
961 int u = n.mag[n.mag.length-1];
962
963 // Make p positive
964 if (p < 0) {
965 p = -p;
966 int n8 = u & 7;
967 if ((n8 == 3) || (n8 == 7))
968 j = -j; // 3 (011) or 7 (111) mod 8
969 }
970
971 // Get rid of factors of 2 in p
972 while ((p & 3) == 0)
973 p >>= 2;
974 if ((p & 1) == 0) {
975 p >>= 1;
976 if (((u ^ (u>>1)) & 2) != 0)
977 j = -j; // 3 (011) or 5 (101) mod 8
978 }
979 if (p == 1)
980 return j;
981 // Then, apply quadratic reciprocity
982 if ((p & u & 2) != 0) // p = u = 3 (mod 4)?
983 j = -j;
984 // And reduce u mod p
985 u = n.mod(BigInteger.valueOf(p)).intValue();
986
987 // Now compute Jacobi(u,p), u < p
988 while (u != 0) {
989 while ((u & 3) == 0)
990 u >>= 2;
991 if ((u & 1) == 0) {
992 u >>= 1;
993 if (((p ^ (p>>1)) & 2) != 0)
994 j = -j; // 3 (011) or 5 (101) mod 8
995 }
996 if (u == 1)
997 return j;
998 // Now both u and p are odd, so use quadratic reciprocity
999 assert (u < p);
1000 int t = u; u = p; p = t;
1001 if ((u & p & 2) != 0) // u = p = 3 (mod 4)?
1002 j = -j;
1003 // Now u >= p, so it can be reduced
1004 u %= p;
1005 }
1006 return 0;
1007 }
1008
1009 private static BigInteger lucasLehmerSequence(int z, BigInteger k, BigInteger n) {
1010 BigInteger d = BigInteger.valueOf(z);
1011 BigInteger u = ONE; BigInteger u2;
1012 BigInteger v = ONE; BigInteger v2;
1013
1014 for (int i=k.bitLength()-2; i >= 0; i--) {
1015 u2 = u.multiply(v).mod(n);
1016
1017 v2 = v.square().add(d.multiply(u.square())).mod(n);
1018 if (v2.testBit(0))
1019 v2 = v2.subtract(n);
1020
1021 v2 = v2.shiftRight(1);
1022
1023 u = u2; v = v2;
1024 if (k.testBit(i)) {
1025 u2 = u.add(v).mod(n);
1026 if (u2.testBit(0))
1027 u2 = u2.subtract(n);
1028
1029 u2 = u2.shiftRight(1);
1030 v2 = v.add(d.multiply(u)).mod(n);
1031 if (v2.testBit(0))
1032 v2 = v2.subtract(n);
1033 v2 = v2.shiftRight(1);
1034
1035 u = u2; v = v2;
1036 }
1037 }
1038 return u;
1039 }
1040
1041 /**
1042 * Returns true iff this BigInteger passes the specified number of
1043 * Miller-Rabin tests. This test is taken from the DSA spec (NIST FIPS
1044 * 186-2).
1045 *
1046 * The following assumptions are made:
1047 * This BigInteger is a positive, odd number greater than 2.
1048 * iterations<=50.
1049 */
1050 private boolean passesMillerRabin(int iterations, Random rnd) {
1051 // Find a and m such that m is odd and this == 1 + 2**a * m
1052 BigInteger thisMinusOne = this.subtract(ONE);
1053 BigInteger m = thisMinusOne;
1054 int a = m.getLowestSetBit();
1055 m = m.shiftRight(a);
1056
1057 // Do the tests
1058 if (rnd == null) {
1059 rnd = ThreadLocalRandom.current();
1060 }
1061 for (int i=0; i < iterations; i++) {
1062 // Generate a uniform random on (1, this)
1063 BigInteger b;
1064 do {
1065 b = new BigInteger(this.bitLength(), rnd);
1066 } while (b.compareTo(ONE) <= 0 || b.compareTo(this) >= 0);
1067
1068 int j = 0;
1069 BigInteger z = b.modPow(m, this);
1070 while (!((j == 0 && z.equals(ONE)) || z.equals(thisMinusOne))) {
1071 if (j > 0 && z.equals(ONE) || ++j == a)
1072 return false;
1073 z = z.modPow(TWO, this);
1074 }
1075 }
1076 return true;
1077 }
1078
1079 /**
1080 * This internal constructor differs from its public cousin
1081 * with the arguments reversed in two ways: it assumes that its
1082 * arguments are correct, and it doesn't copy the magnitude array.
1083 */
1084 BigInteger(int[] magnitude, int signum) {
1085 this.signum = (magnitude.length == 0 ? 0 : signum);
1086 this.mag = magnitude;
1087 if (mag.length >= MAX_MAG_LENGTH) {
1088 checkRange();
1089 }
1090 }
1091
1092 /**
1093 * This private constructor is for internal use and assumes that its
1094 * arguments are correct.
1095 */
1096 private BigInteger(byte[] magnitude, int signum) {
1097 this.signum = (magnitude.length == 0 ? 0 : signum);
1098 this.mag = stripLeadingZeroBytes(magnitude, 0, magnitude.length);
1099 if (mag.length >= MAX_MAG_LENGTH) {
1100 checkRange();
1101 }
1102 }
1103
1104 /**
1105 * Throws an {@code ArithmeticException} if the {@code BigInteger} would be
1106 * out of the supported range.
1107 *
1108 * @throws ArithmeticException if {@code this} exceeds the supported range.
1109 */
1110 private void checkRange() {
1111 if (mag.length > MAX_MAG_LENGTH || mag.length == MAX_MAG_LENGTH && mag[0] < 0) {
1112 reportOverflow();
1113 }
1114 }
1115
1116 private static void reportOverflow() {
1117 throw new ArithmeticException("BigInteger would overflow supported range");
1118 }
1119
1120 //Static Factory Methods
1121
1122 /**
1123 * Returns a BigInteger whose value is equal to that of the
1124 * specified {@code long}. This "static factory method" is
1125 * provided in preference to a ({@code long}) constructor
1126 * because it allows for reuse of frequently used BigIntegers.
1127 *
1128 * @param val value of the BigInteger to return.
1129 * @return a BigInteger with the specified value.
1130 */
1131 public static BigInteger valueOf(long val) {
1132 // If -MAX_CONSTANT < val < MAX_CONSTANT, return stashed constant
1133 if (val == 0)
1134 return ZERO;
1135 if (val > 0 && val <= MAX_CONSTANT)
1136 return posConst[(int) val];
1137 else if (val < 0 && val >= -MAX_CONSTANT)
1138 return negConst[(int) -val];
1139
1140 return new BigInteger(val);
1141 }
1142
1143 /**
1144 * Constructs a BigInteger with the specified value, which may not be zero.
1145 */
1146 private BigInteger(long val) {
1147 if (val < 0) {
1148 val = -val;
1149 signum = -1;
1150 } else {
1151 signum = 1;
1152 }
1153
1154 int highWord = (int)(val >>> 32);
1155 if (highWord == 0) {
1156 mag = new int[1];
1157 mag[0] = (int)val;
1158 } else {
1159 mag = new int[2];
1160 mag[0] = highWord;
1161 mag[1] = (int)val;
1162 }
1163 }
1164
1165 /**
1166 * Returns a BigInteger with the given two's complement representation.
1167 * Assumes that the input array will not be modified (the returned
1168 * BigInteger will reference the input array if feasible).
1169 */
1170 private static BigInteger valueOf(int val[]) {
1171 return (val[0] > 0 ? new BigInteger(val, 1) : new BigInteger(val));
1172 }
1173
1174 // Constants
1175
1176 /**
1177 * Initialize static constant array when class is loaded.
1178 */
1179 private final static int MAX_CONSTANT = 16;
1180 private static BigInteger posConst[] = new BigInteger[MAX_CONSTANT+1];
1181 private static BigInteger negConst[] = new BigInteger[MAX_CONSTANT+1];
1182
1183 /**
1184 * The cache of powers of each radix. This allows us to not have to
1185 * recalculate powers of radix^(2^n) more than once. This speeds
1186 * Schoenhage recursive base conversion significantly.
1187 */
1188 private static volatile BigInteger[][] powerCache;
1189
1190 /** The cache of logarithms of radices for base conversion. */
1191 private static final double[] logCache;
1192
1193 /** The natural log of 2. This is used in computing cache indices. */
1194 private static final double LOG_TWO = Math.log(2.0);
1195
1196 static {
1197 for (int i = 1; i <= MAX_CONSTANT; i++) {
1198 int[] magnitude = new int[1];
1199 magnitude[0] = i;
1200 posConst[i] = new BigInteger(magnitude, 1);
1201 negConst[i] = new BigInteger(magnitude, -1);
1202 }
1203
1204 /*
1205 * Initialize the cache of radix^(2^x) values used for base conversion
1206 * with just the very first value. Additional values will be created
1207 * on demand.
1208 */
1209 powerCache = new BigInteger[Character.MAX_RADIX+1][];
1210 logCache = new double[Character.MAX_RADIX+1];
1211
1212 for (int i=Character.MIN_RADIX; i <= Character.MAX_RADIX; i++) {
1213 powerCache[i] = new BigInteger[] { BigInteger.valueOf(i) };
1214 logCache[i] = Math.log(i);
1215 }
1216 }
1217
1218 /**
1219 * The BigInteger constant zero.
1220 *
1221 * @since 1.2
1222 */
1223 public static final BigInteger ZERO = new BigInteger(new int[0], 0);
1224
1225 /**
1226 * The BigInteger constant one.
1227 *
1228 * @since 1.2
1229 */
1230 public static final BigInteger ONE = valueOf(1);
1231
1232 /**
1233 * The BigInteger constant two. (Not exported.)
1234 */
1235 private static final BigInteger TWO = valueOf(2);
1236
1237 /**
1238 * The BigInteger constant -1. (Not exported.)
1239 */
1240 private static final BigInteger NEGATIVE_ONE = valueOf(-1);
1241
1242 /**
1243 * The BigInteger constant ten.
1244 *
1245 * @since 1.5
1246 */
1247 public static final BigInteger TEN = valueOf(10);
1248
1249 // Arithmetic Operations
1250
1251 /**
1252 * Returns a BigInteger whose value is {@code (this + val)}.
1253 *
1254 * @param val value to be added to this BigInteger.
1255 * @return {@code this + val}
1256 */
1257 public BigInteger add(BigInteger val) {
1258 if (val.signum == 0)
1259 return this;
1260 if (signum == 0)
1261 return val;
1262 if (val.signum == signum)
1263 return new BigInteger(add(mag, val.mag), signum);
1264
1265 int cmp = compareMagnitude(val);
1266 if (cmp == 0)
1267 return ZERO;
1268 int[] resultMag = (cmp > 0 ? subtract(mag, val.mag)
1269 : subtract(val.mag, mag));
1270 resultMag = trustedStripLeadingZeroInts(resultMag);
1271
1272 return new BigInteger(resultMag, cmp == signum ? 1 : -1);
1273 }
1274
1275 /**
1276 * Package private methods used by BigDecimal code to add a BigInteger
1277 * with a long. Assumes val is not equal to INFLATED.
1278 */
1279 BigInteger add(long val) {
1280 if (val == 0)
1281 return this;
1282 if (signum == 0)
1283 return valueOf(val);
1284 if (Long.signum(val) == signum)
1285 return new BigInteger(add(mag, Math.abs(val)), signum);
1286 int cmp = compareMagnitude(val);
1287 if (cmp == 0)
1288 return ZERO;
1289 int[] resultMag = (cmp > 0 ? subtract(mag, Math.abs(val)) : subtract(Math.abs(val), mag));
1290 resultMag = trustedStripLeadingZeroInts(resultMag);
1291 return new BigInteger(resultMag, cmp == signum ? 1 : -1);
1292 }
1293
1294 /**
1295 * Adds the contents of the int array x and long value val. This
1296 * method allocates a new int array to hold the answer and returns
1297 * a reference to that array. Assumes x.length > 0 and val is
1298 * non-negative
1299 */
1300 private static int[] add(int[] x, long val) {
1301 int[] y;
1302 long sum = 0;
1303 int xIndex = x.length;
1304 int[] result;
1305 int highWord = (int)(val >>> 32);
1306 if (highWord == 0) {
1307 result = new int[xIndex];
1308 sum = (x[--xIndex] & LONG_MASK) + val;
1309 result[xIndex] = (int)sum;
1310 } else {
1311 if (xIndex == 1) {
1312 result = new int[2];
1313 sum = val + (x[0] & LONG_MASK);
1314 result[1] = (int)sum;
1315 result[0] = (int)(sum >>> 32);
1316 return result;
1317 } else {
1318 result = new int[xIndex];
1319 sum = (x[--xIndex] & LONG_MASK) + (val & LONG_MASK);
1320 result[xIndex] = (int)sum;
1321 sum = (x[--xIndex] & LONG_MASK) + (highWord & LONG_MASK) + (sum >>> 32);
1322 result[xIndex] = (int)sum;
1323 }
1324 }
1325 // Copy remainder of longer number while carry propagation is required
1326 boolean carry = (sum >>> 32 != 0);
1327 while (xIndex > 0 && carry)
1328 carry = ((result[--xIndex] = x[xIndex] + 1) == 0);
1329 // Copy remainder of longer number
1330 while (xIndex > 0)
1331 result[--xIndex] = x[xIndex];
1332 // Grow result if necessary
1333 if (carry) {
1334 int bigger[] = new int[result.length + 1];
1335 System.arraycopy(result, 0, bigger, 1, result.length);
1336 bigger[0] = 0x01;
1337 return bigger;
1338 }
1339 return result;
1340 }
1341
1342 /**
1343 * Adds the contents of the int arrays x and y. This method allocates
1344 * a new int array to hold the answer and returns a reference to that
1345 * array.
1346 */
1347 private static int[] add(int[] x, int[] y) {
1348 // If x is shorter, swap the two arrays
1349 if (x.length < y.length) {
1350 int[] tmp = x;
1351 x = y;
1352 y = tmp;
1353 }
1354
1355 int xIndex = x.length;
1356 int yIndex = y.length;
1357 int result[] = new int[xIndex];
1358 long sum = 0;
1359 if (yIndex == 1) {
1360 sum = (x[--xIndex] & LONG_MASK) + (y[0] & LONG_MASK) ;
1361 result[xIndex] = (int)sum;
1362 } else {
1363 // Add common parts of both numbers
1364 while (yIndex > 0) {
1365 sum = (x[--xIndex] & LONG_MASK) +
1366 (y[--yIndex] & LONG_MASK) + (sum >>> 32);
1367 result[xIndex] = (int)sum;
1368 }
1369 }
1370 // Copy remainder of longer number while carry propagation is required
1371 boolean carry = (sum >>> 32 != 0);
1372 while (xIndex > 0 && carry)
1373 carry = ((result[--xIndex] = x[xIndex] + 1) == 0);
1374
1375 // Copy remainder of longer number
1376 while (xIndex > 0)
1377 result[--xIndex] = x[xIndex];
1378
1379 // Grow result if necessary
1380 if (carry) {
1381 int bigger[] = new int[result.length + 1];
1382 System.arraycopy(result, 0, bigger, 1, result.length);
1383 bigger[0] = 0x01;
1384 return bigger;
1385 }
1386 return result;
1387 }
1388
1389 private static int[] subtract(long val, int[] little) {
1390 int highWord = (int)(val >>> 32);
1391 if (highWord == 0) {
1392 int result[] = new int[1];
1393 result[0] = (int)(val - (little[0] & LONG_MASK));
1394 return result;
1395 } else {
1396 int result[] = new int[2];
1397 if (little.length == 1) {
1398 long difference = ((int)val & LONG_MASK) - (little[0] & LONG_MASK);
1399 result[1] = (int)difference;
1400 // Subtract remainder of longer number while borrow propagates
1401 boolean borrow = (difference >> 32 != 0);
1402 if (borrow) {
1403 result[0] = highWord - 1;
1404 } else { // Copy remainder of longer number
1405 result[0] = highWord;
1406 }
1407 return result;
1408 } else { // little.length == 2
1409 long difference = ((int)val & LONG_MASK) - (little[1] & LONG_MASK);
1410 result[1] = (int)difference;
1411 difference = (highWord & LONG_MASK) - (little[0] & LONG_MASK) + (difference >> 32);
1412 result[0] = (int)difference;
1413 return result;
1414 }
1415 }
1416 }
1417
1418 /**
1419 * Subtracts the contents of the second argument (val) from the
1420 * first (big). The first int array (big) must represent a larger number
1421 * than the second. This method allocates the space necessary to hold the
1422 * answer.
1423 * assumes val >= 0
1424 */
1425 private static int[] subtract(int[] big, long val) {
1426 int highWord = (int)(val >>> 32);
1427 int bigIndex = big.length;
1428 int result[] = new int[bigIndex];
1429 long difference = 0;
1430
1431 if (highWord == 0) {
1432 difference = (big[--bigIndex] & LONG_MASK) - val;
1433 result[bigIndex] = (int)difference;
1434 } else {
1435 difference = (big[--bigIndex] & LONG_MASK) - (val & LONG_MASK);
1436 result[bigIndex] = (int)difference;
1437 difference = (big[--bigIndex] & LONG_MASK) - (highWord & LONG_MASK) + (difference >> 32);
1438 result[bigIndex] = (int)difference;
1439 }
1440
1441 // Subtract remainder of longer number while borrow propagates
1442 boolean borrow = (difference >> 32 != 0);
1443 while (bigIndex > 0 && borrow)
1444 borrow = ((result[--bigIndex] = big[bigIndex] - 1) == -1);
1445
1446 // Copy remainder of longer number
1447 while (bigIndex > 0)
1448 result[--bigIndex] = big[bigIndex];
1449
1450 return result;
1451 }
1452
1453 /**
1454 * Returns a BigInteger whose value is {@code (this - val)}.
1455 *
1456 * @param val value to be subtracted from this BigInteger.
1457 * @return {@code this - val}
1458 */
1459 public BigInteger subtract(BigInteger val) {
1460 if (val.signum == 0)
1461 return this;
1462 if (signum == 0)
1463 return val.negate();
1464 if (val.signum != signum)
1465 return new BigInteger(add(mag, val.mag), signum);
1466
1467 int cmp = compareMagnitude(val);
1468 if (cmp == 0)
1469 return ZERO;
1470 int[] resultMag = (cmp > 0 ? subtract(mag, val.mag)
1471 : subtract(val.mag, mag));
1472 resultMag = trustedStripLeadingZeroInts(resultMag);
1473 return new BigInteger(resultMag, cmp == signum ? 1 : -1);
1474 }
1475
1476 /**
1477 * Subtracts the contents of the second int arrays (little) from the
1478 * first (big). The first int array (big) must represent a larger number
1479 * than the second. This method allocates the space necessary to hold the
1480 * answer.
1481 */
1482 private static int[] subtract(int[] big, int[] little) {
1483 int bigIndex = big.length;
1484 int result[] = new int[bigIndex];
1485 int littleIndex = little.length;
1486 long difference = 0;
1487
1488 // Subtract common parts of both numbers
1489 while (littleIndex > 0) {
1490 difference = (big[--bigIndex] & LONG_MASK) -
1491 (little[--littleIndex] & LONG_MASK) +
1492 (difference >> 32);
1493 result[bigIndex] = (int)difference;
1494 }
1495
1496 // Subtract remainder of longer number while borrow propagates
1497 boolean borrow = (difference >> 32 != 0);
1498 while (bigIndex > 0 && borrow)
1499 borrow = ((result[--bigIndex] = big[bigIndex] - 1) == -1);
1500
1501 // Copy remainder of longer number
1502 while (bigIndex > 0)
1503 result[--bigIndex] = big[bigIndex];
1504
1505 return result;
1506 }
1507
1508 /**
1509 * Returns a BigInteger whose value is {@code (this * val)}.
1510 *
1511 * @implNote An implementation may offer better algorithmic
1512 * performance when {@code val == this}.
1513 *
1514 * @param val value to be multiplied by this BigInteger.
1515 * @return {@code this * val}
1516 */
1517 public BigInteger multiply(BigInteger val) {
1518 if (val.signum == 0 || signum == 0)
1519 return ZERO;
1520
1521 int xlen = mag.length;
1522
1523 if (val == this && xlen > MULTIPLY_SQUARE_THRESHOLD) {
1524 return square();
1525 }
1526
1527 int ylen = val.mag.length;
1528
1529 if ((xlen < KARATSUBA_THRESHOLD) || (ylen < KARATSUBA_THRESHOLD)) {
1530 int resultSign = signum == val.signum ? 1 : -1;
1531 if (val.mag.length == 1) {
1532 return multiplyByInt(mag,val.mag[0], resultSign);
1533 }
1534 if (mag.length == 1) {
1535 return multiplyByInt(val.mag,mag[0], resultSign);
1536 }
1537 int[] result = multiplyToLen(mag, xlen,
1538 val.mag, ylen, null);
1539 result = trustedStripLeadingZeroInts(result);
1540 return new BigInteger(result, resultSign);
1541 } else {
1542 if ((xlen < TOOM_COOK_THRESHOLD) && (ylen < TOOM_COOK_THRESHOLD)) {
1543 return multiplyKaratsuba(this, val);
1544 } else {
1545 return multiplyToomCook3(this, val);
1546 }
1547 }
1548 }
1549
1550 private static BigInteger multiplyByInt(int[] x, int y, int sign) {
1551 if (Integer.bitCount(y) == 1) {
1552 return new BigInteger(shiftLeft(x,Integer.numberOfTrailingZeros(y)), sign);
1553 }
1554 int xlen = x.length;
1555 int[] rmag = new int[xlen + 1];
1556 long carry = 0;
1557 long yl = y & LONG_MASK;
1558 int rstart = rmag.length - 1;
1559 for (int i = xlen - 1; i >= 0; i--) {
1560 long product = (x[i] & LONG_MASK) * yl + carry;
1561 rmag[rstart--] = (int)product;
1562 carry = product >>> 32;
1563 }
1564 if (carry == 0L) {
1565 rmag = java.util.Arrays.copyOfRange(rmag, 1, rmag.length);
1566 } else {
1567 rmag[rstart] = (int)carry;
1568 }
1569 return new BigInteger(rmag, sign);
1570 }
1571
1572 /**
1573 * Package private methods used by BigDecimal code to multiply a BigInteger
1574 * with a long. Assumes v is not equal to INFLATED.
1575 */
1576 BigInteger multiply(long v) {
1577 if (v == 0 || signum == 0)
1578 return ZERO;
1579 if (v == BigDecimal.INFLATED)
1580 return multiply(BigInteger.valueOf(v));
1581 int rsign = (v > 0 ? signum : -signum);
1582 if (v < 0)
1583 v = -v;
1584 long dh = v >>> 32; // higher order bits
1585 long dl = v & LONG_MASK; // lower order bits
1586
1587 int xlen = mag.length;
1588 int[] value = mag;
1589 int[] rmag = (dh == 0L) ? (new int[xlen + 1]) : (new int[xlen + 2]);
1590 long carry = 0;
1591 int rstart = rmag.length - 1;
1592 for (int i = xlen - 1; i >= 0; i--) {
1593 long product = (value[i] & LONG_MASK) * dl + carry;
1594 rmag[rstart--] = (int)product;
1595 carry = product >>> 32;
1596 }
1597 rmag[rstart] = (int)carry;
1598 if (dh != 0L) {
1599 carry = 0;
1600 rstart = rmag.length - 2;
1601 for (int i = xlen - 1; i >= 0; i--) {
1602 long product = (value[i] & LONG_MASK) * dh +
1603 (rmag[rstart] & LONG_MASK) + carry;
1604 rmag[rstart--] = (int)product;
1605 carry = product >>> 32;
1606 }
1607 rmag[0] = (int)carry;
1608 }
1609 if (carry == 0L)
1610 rmag = java.util.Arrays.copyOfRange(rmag, 1, rmag.length);
1611 return new BigInteger(rmag, rsign);
1612 }
1613
1614 /**
1615 * Multiplies int arrays x and y to the specified lengths and places
1616 * the result into z. There will be no leading zeros in the resultant array.
1617 */
1618 private int[] multiplyToLen(int[] x, int xlen, int[] y, int ylen, int[] z) {
1619 int xstart = xlen - 1;
1620 int ystart = ylen - 1;
1621
1622 if (z == null || z.length < (xlen+ ylen))
1623 z = new int[xlen+ylen];
1624
1625 long carry = 0;
1626 for (int j=ystart, k=ystart+1+xstart; j >= 0; j--, k--) {
1627 long product = (y[j] & LONG_MASK) *
1628 (x[xstart] & LONG_MASK) + carry;
1629 z[k] = (int)product;
1630 carry = product >>> 32;
1631 }
1632 z[xstart] = (int)carry;
1633
1634 for (int i = xstart-1; i >= 0; i--) {
1635 carry = 0;
1636 for (int j=ystart, k=ystart+1+i; j >= 0; j--, k--) {
1637 long product = (y[j] & LONG_MASK) *
1638 (x[i] & LONG_MASK) +
1639 (z[k] & LONG_MASK) + carry;
1640 z[k] = (int)product;
1641 carry = product >>> 32;
1642 }
1643 z[i] = (int)carry;
1644 }
1645 return z;
1646 }
1647
1648 /**
1649 * Multiplies two BigIntegers using the Karatsuba multiplication
1650 * algorithm. This is a recursive divide-and-conquer algorithm which is
1651 * more efficient for large numbers than what is commonly called the
1652 * "grade-school" algorithm used in multiplyToLen. If the numbers to be
1653 * multiplied have length n, the "grade-school" algorithm has an
1654 * asymptotic complexity of O(n^2). In contrast, the Karatsuba algorithm
1655 * has complexity of O(n^(log2(3))), or O(n^1.585). It achieves this
1656 * increased performance by doing 3 multiplies instead of 4 when
1657 * evaluating the product. As it has some overhead, should be used when
1658 * both numbers are larger than a certain threshold (found
1659 * experimentally).
1660 *
1661 * See: http://en.wikipedia.org/wiki/Karatsuba_algorithm
1662 */
1663 private static BigInteger multiplyKaratsuba(BigInteger x, BigInteger y) {
1664 int xlen = x.mag.length;
1665 int ylen = y.mag.length;
1666
1667 // The number of ints in each half of the number.
1668 int half = (Math.max(xlen, ylen)+1) / 2;
1669
1670 // xl and yl are the lower halves of x and y respectively,
1671 // xh and yh are the upper halves.
1672 BigInteger xl = x.getLower(half);
1673 BigInteger xh = x.getUpper(half);
1674 BigInteger yl = y.getLower(half);
1675 BigInteger yh = y.getUpper(half);
1676
1677 BigInteger p1 = xh.multiply(yh); // p1 = xh*yh
1678 BigInteger p2 = xl.multiply(yl); // p2 = xl*yl
1679
1680 // p3=(xh+xl)*(yh+yl)
1681 BigInteger p3 = xh.add(xl).multiply(yh.add(yl));
1682
1683 // result = p1 * 2^(32*2*half) + (p3 - p1 - p2) * 2^(32*half) + p2
1684 BigInteger result = p1.shiftLeft(32*half).add(p3.subtract(p1).subtract(p2)).shiftLeft(32*half).add(p2);
1685
1686 if (x.signum != y.signum) {
1687 return result.negate();
1688 } else {
1689 return result;
1690 }
1691 }
1692
1693 /**
1694 * Multiplies two BigIntegers using a 3-way Toom-Cook multiplication
1695 * algorithm. This is a recursive divide-and-conquer algorithm which is
1696 * more efficient for large numbers than what is commonly called the
1697 * "grade-school" algorithm used in multiplyToLen. If the numbers to be
1698 * multiplied have length n, the "grade-school" algorithm has an
1699 * asymptotic complexity of O(n^2). In contrast, 3-way Toom-Cook has a
1700 * complexity of about O(n^1.465). It achieves this increased asymptotic
1701 * performance by breaking each number into three parts and by doing 5
1702 * multiplies instead of 9 when evaluating the product. Due to overhead
1703 * (additions, shifts, and one division) in the Toom-Cook algorithm, it
1704 * should only be used when both numbers are larger than a certain
1705 * threshold (found experimentally). This threshold is generally larger
1706 * than that for Karatsuba multiplication, so this algorithm is generally
1707 * only used when numbers become significantly larger.
1708 *
1709 * The algorithm used is the "optimal" 3-way Toom-Cook algorithm outlined
1710 * by Marco Bodrato.
1711 *
1712 * See: http://bodrato.it/toom-cook/
1713 * http://bodrato.it/papers/#WAIFI2007
1714 *
1715 * "Towards Optimal Toom-Cook Multiplication for Univariate and
1716 * Multivariate Polynomials in Characteristic 2 and 0." by Marco BODRATO;
1717 * In C.Carlet and B.Sunar, Eds., "WAIFI'07 proceedings", p. 116-133,
1718 * LNCS #4547. Springer, Madrid, Spain, June 21-22, 2007.
1719 *
1720 */
1721 private static BigInteger multiplyToomCook3(BigInteger a, BigInteger b) {
1722 int alen = a.mag.length;
1723 int blen = b.mag.length;
1724
1725 int largest = Math.max(alen, blen);
1726
1727 // k is the size (in ints) of the lower-order slices.
1728 int k = (largest+2)/3; // Equal to ceil(largest/3)
1729
1730 // r is the size (in ints) of the highest-order slice.
1731 int r = largest - 2*k;
1732
1733 // Obtain slices of the numbers. a2 and b2 are the most significant
1734 // bits of the numbers a and b, and a0 and b0 the least significant.
1735 BigInteger a0, a1, a2, b0, b1, b2;
1736 a2 = a.getToomSlice(k, r, 0, largest);
1737 a1 = a.getToomSlice(k, r, 1, largest);
1738 a0 = a.getToomSlice(k, r, 2, largest);
1739 b2 = b.getToomSlice(k, r, 0, largest);
1740 b1 = b.getToomSlice(k, r, 1, largest);
1741 b0 = b.getToomSlice(k, r, 2, largest);
1742
1743 BigInteger v0, v1, v2, vm1, vinf, t1, t2, tm1, da1, db1;
1744
1745 v0 = a0.multiply(b0);
1746 da1 = a2.add(a0);
1747 db1 = b2.add(b0);
1748 vm1 = da1.subtract(a1).multiply(db1.subtract(b1));
1749 da1 = da1.add(a1);
1750 db1 = db1.add(b1);
1751 v1 = da1.multiply(db1);
1752 v2 = da1.add(a2).shiftLeft(1).subtract(a0).multiply(
1753 db1.add(b2).shiftLeft(1).subtract(b0));
1754 vinf = a2.multiply(b2);
1755
1756 // The algorithm requires two divisions by 2 and one by 3.
1757 // All divisions are known to be exact, that is, they do not produce
1758 // remainders, and all results are positive. The divisions by 2 are
1759 // implemented as right shifts which are relatively efficient, leaving
1760 // only an exact division by 3, which is done by a specialized
1761 // linear-time algorithm.
1762 t2 = v2.subtract(vm1).exactDivideBy3();
1763 tm1 = v1.subtract(vm1).shiftRight(1);
1764 t1 = v1.subtract(v0);
1765 t2 = t2.subtract(t1).shiftRight(1);
1766 t1 = t1.subtract(tm1).subtract(vinf);
1767 t2 = t2.subtract(vinf.shiftLeft(1));
1768 tm1 = tm1.subtract(t2);
1769
1770 // Number of bits to shift left.
1771 int ss = k*32;
1772
1773 BigInteger result = vinf.shiftLeft(ss).add(t2).shiftLeft(ss).add(t1).shiftLeft(ss).add(tm1).shiftLeft(ss).add(v0);
1774
1775 if (a.signum != b.signum) {
1776 return result.negate();
1777 } else {
1778 return result;
1779 }
1780 }
1781
1782
1783 /**
1784 * Returns a slice of a BigInteger for use in Toom-Cook multiplication.
1785 *
1786 * @param lowerSize The size of the lower-order bit slices.
1787 * @param upperSize The size of the higher-order bit slices.
1788 * @param slice The index of which slice is requested, which must be a
1789 * number from 0 to size-1. Slice 0 is the highest-order bits, and slice
1790 * size-1 are the lowest-order bits. Slice 0 may be of different size than
1791 * the other slices.
1792 * @param fullsize The size of the larger integer array, used to align
1793 * slices to the appropriate position when multiplying different-sized
1794 * numbers.
1795 */
1796 private BigInteger getToomSlice(int lowerSize, int upperSize, int slice,
1797 int fullsize) {
1798 int start, end, sliceSize, len, offset;
1799
1800 len = mag.length;
1801 offset = fullsize - len;
1802
1803 if (slice == 0) {
1804 start = 0 - offset;
1805 end = upperSize - 1 - offset;
1806 } else {
1807 start = upperSize + (slice-1)*lowerSize - offset;
1808 end = start + lowerSize - 1;
1809 }
1810
1811 if (start < 0) {
1812 start = 0;
1813 }
1814 if (end < 0) {
1815 return ZERO;
1816 }
1817
1818 sliceSize = (end-start) + 1;
1819
1820 if (sliceSize <= 0) {
1821 return ZERO;
1822 }
1823
1824 // While performing Toom-Cook, all slices are positive and
1825 // the sign is adjusted when the final number is composed.
1826 if (start == 0 && sliceSize >= len) {
1827 return this.abs();
1828 }
1829
1830 int intSlice[] = new int[sliceSize];
1831 System.arraycopy(mag, start, intSlice, 0, sliceSize);
1832
1833 return new BigInteger(trustedStripLeadingZeroInts(intSlice), 1);
1834 }
1835
1836 /**
1837 * Does an exact division (that is, the remainder is known to be zero)
1838 * of the specified number by 3. This is used in Toom-Cook
1839 * multiplication. This is an efficient algorithm that runs in linear
1840 * time. If the argument is not exactly divisible by 3, results are
1841 * undefined. Note that this is expected to be called with positive
1842 * arguments only.
1843 */
1844 private BigInteger exactDivideBy3() {
1845 int len = mag.length;
1846 int[] result = new int[len];
1847 long x, w, q, borrow;
1848 borrow = 0L;
1849 for (int i=len-1; i >= 0; i--) {
1850 x = (mag[i] & LONG_MASK);
1851 w = x - borrow;
1852 if (borrow > x) { // Did we make the number go negative?
1853 borrow = 1L;
1854 } else {
1855 borrow = 0L;
1856 }
1857
1858 // 0xAAAAAAAB is the modular inverse of 3 (mod 2^32). Thus,
1859 // the effect of this is to divide by 3 (mod 2^32).
1860 // This is much faster than division on most architectures.
1861 q = (w * 0xAAAAAAABL) & LONG_MASK;
1862 result[i] = (int) q;
1863
1864 // Now check the borrow. The second check can of course be
1865 // eliminated if the first fails.
1866 if (q >= 0x55555556L) {
1867 borrow++;
1868 if (q >= 0xAAAAAAABL)
1869 borrow++;
1870 }
1871 }
1872 result = trustedStripLeadingZeroInts(result);
1873 return new BigInteger(result, signum);
1874 }
1875
1876 /**
1877 * Returns a new BigInteger representing n lower ints of the number.
1878 * This is used by Karatsuba multiplication and Karatsuba squaring.
1879 */
1880 private BigInteger getLower(int n) {
1881 int len = mag.length;
1882
1883 if (len <= n) {
1884 return abs();
1885 }
1886
1887 int lowerInts[] = new int[n];
1888 System.arraycopy(mag, len-n, lowerInts, 0, n);
1889
1890 return new BigInteger(trustedStripLeadingZeroInts(lowerInts), 1);
1891 }
1892
1893 /**
1894 * Returns a new BigInteger representing mag.length-n upper
1895 * ints of the number. This is used by Karatsuba multiplication and
1896 * Karatsuba squaring.
1897 */
1898 private BigInteger getUpper(int n) {
1899 int len = mag.length;
1900
1901 if (len <= n) {
1902 return ZERO;
1903 }
1904
1905 int upperLen = len - n;
1906 int upperInts[] = new int[upperLen];
1907 System.arraycopy(mag, 0, upperInts, 0, upperLen);
1908
1909 return new BigInteger(trustedStripLeadingZeroInts(upperInts), 1);
1910 }
1911
1912 // Squaring
1913
1914 /**
1915 * Returns a BigInteger whose value is {@code (this<sup>2</sup>)}.
1916 *
1917 * @return {@code this<sup>2</sup>}
1918 */
1919 private BigInteger square() {
1920 if (signum == 0) {
1921 return ZERO;
1922 }
1923 int len = mag.length;
1924
1925 if (len < KARATSUBA_SQUARE_THRESHOLD) {
1926 int[] z = squareToLen(mag, len, null);
1927 return new BigInteger(trustedStripLeadingZeroInts(z), 1);
1928 } else {
1929 if (len < TOOM_COOK_SQUARE_THRESHOLD) {
1930 return squareKaratsuba();
1931 } else {
1932 return squareToomCook3();
1933 }
1934 }
1935 }
1936
1937 /**
1938 * Squares the contents of the int array x. The result is placed into the
1939 * int array z. The contents of x are not changed.
1940 */
1941 private static final int[] squareToLen(int[] x, int len, int[] z) {
1942 /*
1943 * The algorithm used here is adapted from Colin Plumb's C library.
1944 * Technique: Consider the partial products in the multiplication
1945 * of "abcde" by itself:
1946 *
1947 * a b c d e
1948 * * a b c d e
1949 * ==================
1950 * ae be ce de ee
1951 * ad bd cd dd de
1952 * ac bc cc cd ce
1953 * ab bb bc bd be
1954 * aa ab ac ad ae
1955 *
1956 * Note that everything above the main diagonal:
1957 * ae be ce de = (abcd) * e
1958 * ad bd cd = (abc) * d
1959 * ac bc = (ab) * c
1960 * ab = (a) * b
1961 *
1962 * is a copy of everything below the main diagonal:
1963 * de
1964 * cd ce
1965 * bc bd be
1966 * ab ac ad ae
1967 *
1968 * Thus, the sum is 2 * (off the diagonal) + diagonal.
1969 *
1970 * This is accumulated beginning with the diagonal (which
1971 * consist of the squares of the digits of the input), which is then
1972 * divided by two, the off-diagonal added, and multiplied by two
1973 * again. The low bit is simply a copy of the low bit of the
1974 * input, so it doesn't need special care.
1975 */
1976 int zlen = len << 1;
1977 if (z == null || z.length < zlen)
1978 z = new int[zlen];
1979
1980 // Store the squares, right shifted one bit (i.e., divided by 2)
1981 int lastProductLowWord = 0;
1982 for (int j=0, i=0; j < len; j++) {
1983 long piece = (x[j] & LONG_MASK);
1984 long product = piece * piece;
1985 z[i++] = (lastProductLowWord << 31) | (int)(product >>> 33);
1986 z[i++] = (int)(product >>> 1);
1987 lastProductLowWord = (int)product;
1988 }
1989
1990 // Add in off-diagonal sums
1991 for (int i=len, offset=1; i > 0; i--, offset+=2) {
1992 int t = x[i-1];
1993 t = mulAdd(z, x, offset, i-1, t);
1994 addOne(z, offset-1, i, t);
1995 }
1996
1997 // Shift back up and set low bit
1998 primitiveLeftShift(z, zlen, 1);
1999 z[zlen-1] |= x[len-1] & 1;
2000
2001 return z;
2002 }
2003
2004 /**
2005 * Squares a BigInteger using the Karatsuba squaring algorithm. It should
2006 * be used when both numbers are larger than a certain threshold (found
2007 * experimentally). It is a recursive divide-and-conquer algorithm that
2008 * has better asymptotic performance than the algorithm used in
2009 * squareToLen.
2010 */
2011 private BigInteger squareKaratsuba() {
2012 int half = (mag.length+1) / 2;
2013
2014 BigInteger xl = getLower(half);
2015 BigInteger xh = getUpper(half);
2016
2017 BigInteger xhs = xh.square(); // xhs = xh^2
2018 BigInteger xls = xl.square(); // xls = xl^2
2019
2020 // xh^2 << 64 + (((xl+xh)^2 - (xh^2 + xl^2)) << 32) + xl^2
2021 return xhs.shiftLeft(half*32).add(xl.add(xh).square().subtract(xhs.add(xls))).shiftLeft(half*32).add(xls);
2022 }
2023
2024 /**
2025 * Squares a BigInteger using the 3-way Toom-Cook squaring algorithm. It
2026 * should be used when both numbers are larger than a certain threshold
2027 * (found experimentally). It is a recursive divide-and-conquer algorithm
2028 * that has better asymptotic performance than the algorithm used in
2029 * squareToLen or squareKaratsuba.
2030 */
2031 private BigInteger squareToomCook3() {
2032 int len = mag.length;
2033
2034 // k is the size (in ints) of the lower-order slices.
2035 int k = (len+2)/3; // Equal to ceil(largest/3)
2036
2037 // r is the size (in ints) of the highest-order slice.
2038 int r = len - 2*k;
2039
2040 // Obtain slices of the numbers. a2 is the most significant
2041 // bits of the number, and a0 the least significant.
2042 BigInteger a0, a1, a2;
2043 a2 = getToomSlice(k, r, 0, len);
2044 a1 = getToomSlice(k, r, 1, len);
2045 a0 = getToomSlice(k, r, 2, len);
2046 BigInteger v0, v1, v2, vm1, vinf, t1, t2, tm1, da1;
2047
2048 v0 = a0.square();
2049 da1 = a2.add(a0);
2050 vm1 = da1.subtract(a1).square();
2051 da1 = da1.add(a1);
2052 v1 = da1.square();
2053 vinf = a2.square();
2054 v2 = da1.add(a2).shiftLeft(1).subtract(a0).square();
2055
2056 // The algorithm requires two divisions by 2 and one by 3.
2057 // All divisions are known to be exact, that is, they do not produce
2058 // remainders, and all results are positive. The divisions by 2 are
2059 // implemented as right shifts which are relatively efficient, leaving
2060 // only a division by 3.
2061 // The division by 3 is done by an optimized algorithm for this case.
2062 t2 = v2.subtract(vm1).exactDivideBy3();
2063 tm1 = v1.subtract(vm1).shiftRight(1);
2064 t1 = v1.subtract(v0);
2065 t2 = t2.subtract(t1).shiftRight(1);
2066 t1 = t1.subtract(tm1).subtract(vinf);
2067 t2 = t2.subtract(vinf.shiftLeft(1));
2068 tm1 = tm1.subtract(t2);
2069
2070 // Number of bits to shift left.
2071 int ss = k*32;
2072
2073 return vinf.shiftLeft(ss).add(t2).shiftLeft(ss).add(t1).shiftLeft(ss).add(tm1).shiftLeft(ss).add(v0);
2074 }
2075
2076 // Division
2077
2078 /**
2079 * Returns a BigInteger whose value is {@code (this / val)}.
2080 *
2081 * @param val value by which this BigInteger is to be divided.
2082 * @return {@code this / val}
2083 * @throws ArithmeticException if {@code val} is zero.
2084 */
2085 public BigInteger divide(BigInteger val) {
2086 if (val.mag.length < BURNIKEL_ZIEGLER_THRESHOLD ||
2087 mag.length - val.mag.length < BURNIKEL_ZIEGLER_OFFSET) {
2088 return divideKnuth(val);
2089 } else {
2090 return divideBurnikelZiegler(val);
2091 }
2092 }
2093
2094 /**
2095 * Returns a BigInteger whose value is {@code (this / val)} using an O(n^2) algorithm from Knuth.
2096 *
2097 * @param val value by which this BigInteger is to be divided.
2098 * @return {@code this / val}
2099 * @throws ArithmeticException if {@code val} is zero.
2100 * @see MutableBigInteger#divideKnuth(MutableBigInteger, MutableBigInteger, boolean)
2101 */
2102 private BigInteger divideKnuth(BigInteger val) {
2103 MutableBigInteger q = new MutableBigInteger(),
2104 a = new MutableBigInteger(this.mag),
2105 b = new MutableBigInteger(val.mag);
2106
2107 a.divideKnuth(b, q, false);
2108 return q.toBigInteger(this.signum * val.signum);
2109 }
2110
2111 /**
2112 * Returns an array of two BigIntegers containing {@code (this / val)}
2113 * followed by {@code (this % val)}.
2114 *
2115 * @param val value by which this BigInteger is to be divided, and the
2116 * remainder computed.
2117 * @return an array of two BigIntegers: the quotient {@code (this / val)}
2118 * is the initial element, and the remainder {@code (this % val)}
2119 * is the final element.
2120 * @throws ArithmeticException if {@code val} is zero.
2121 */
2122 public BigInteger[] divideAndRemainder(BigInteger val) {
2123 if (val.mag.length < BURNIKEL_ZIEGLER_THRESHOLD ||
2124 mag.length - val.mag.length < BURNIKEL_ZIEGLER_OFFSET) {
2125 return divideAndRemainderKnuth(val);
2126 } else {
2127 return divideAndRemainderBurnikelZiegler(val);
2128 }
2129 }
2130
2131 /** Long division */
2132 private BigInteger[] divideAndRemainderKnuth(BigInteger val) {
2133 BigInteger[] result = new BigInteger[2];
2134 MutableBigInteger q = new MutableBigInteger(),
2135 a = new MutableBigInteger(this.mag),
2136 b = new MutableBigInteger(val.mag);
2137 MutableBigInteger r = a.divideKnuth(b, q);
2138 result[0] = q.toBigInteger(this.signum == val.signum ? 1 : -1);
2139 result[1] = r.toBigInteger(this.signum);
2140 return result;
2141 }
2142
2143 /**
2144 * Returns a BigInteger whose value is {@code (this % val)}.
2145 *
2146 * @param val value by which this BigInteger is to be divided, and the
2147 * remainder computed.
2148 * @return {@code this % val}
2149 * @throws ArithmeticException if {@code val} is zero.
2150 */
2151 public BigInteger remainder(BigInteger val) {
2152 if (val.mag.length < BURNIKEL_ZIEGLER_THRESHOLD ||
2153 mag.length - val.mag.length < BURNIKEL_ZIEGLER_OFFSET) {
2154 return remainderKnuth(val);
2155 } else {
2156 return remainderBurnikelZiegler(val);
2157 }
2158 }
2159
2160 /** Long division */
2161 private BigInteger remainderKnuth(BigInteger val) {
2162 MutableBigInteger q = new MutableBigInteger(),
2163 a = new MutableBigInteger(this.mag),
2164 b = new MutableBigInteger(val.mag);
2165
2166 return a.divideKnuth(b, q).toBigInteger(this.signum);
2167 }
2168
2169 /**
2170 * Calculates {@code this / val} using the Burnikel-Ziegler algorithm.
2171 * @param val the divisor
2172 * @return {@code this / val}
2173 */
2174 private BigInteger divideBurnikelZiegler(BigInteger val) {
2175 return divideAndRemainderBurnikelZiegler(val)[0];
2176 }
2177
2178 /**
2179 * Calculates {@code this % val} using the Burnikel-Ziegler algorithm.
2180 * @param val the divisor
2181 * @return {@code this % val}
2182 */
2183 private BigInteger remainderBurnikelZiegler(BigInteger val) {
2184 return divideAndRemainderBurnikelZiegler(val)[1];
2185 }
2186
2187 /**
2188 * Computes {@code this / val} and {@code this % val} using the
2189 * Burnikel-Ziegler algorithm.
2190 * @param val the divisor
2191 * @return an array containing the quotient and remainder
2192 */
2193 private BigInteger[] divideAndRemainderBurnikelZiegler(BigInteger val) {
2194 MutableBigInteger q = new MutableBigInteger();
2195 MutableBigInteger r = new MutableBigInteger(this).divideAndRemainderBurnikelZiegler(new MutableBigInteger(val), q);
2196 BigInteger qBigInt = q.isZero() ? ZERO : q.toBigInteger(signum*val.signum);
2197 BigInteger rBigInt = r.isZero() ? ZERO : r.toBigInteger(signum);
2198 return new BigInteger[] {qBigInt, rBigInt};
2199 }
2200
2201 /**
2202 * Returns a BigInteger whose value is <tt>(this<sup>exponent</sup>)</tt>.
2203 * Note that {@code exponent} is an integer rather than a BigInteger.
2204 *
2205 * @param exponent exponent to which this BigInteger is to be raised.
2206 * @return <tt>this<sup>exponent</sup></tt>
2207 * @throws ArithmeticException {@code exponent} is negative. (This would
2208 * cause the operation to yield a non-integer value.)
2209 */
2210 public BigInteger pow(int exponent) {
2211 if (exponent < 0) {
2212 throw new ArithmeticException("Negative exponent");
2213 }
2214 if (signum == 0) {
2215 return (exponent == 0 ? ONE : this);
2216 }
2217
2218 BigInteger partToSquare = this.abs();
2219
2220 // Factor out powers of two from the base, as the exponentiation of
2221 // these can be done by left shifts only.
2222 // The remaining part can then be exponentiated faster. The
2223 // powers of two will be multiplied back at the end.
2224 int powersOfTwo = partToSquare.getLowestSetBit();
2225 long bitsToShift = (long)powersOfTwo * exponent;
2226 if (bitsToShift > Integer.MAX_VALUE) {
2227 reportOverflow();
2228 }
2229
2230 int remainingBits;
2231
2232 // Factor the powers of two out quickly by shifting right, if needed.
2233 if (powersOfTwo > 0) {
2234 partToSquare = partToSquare.shiftRight(powersOfTwo);
2235 remainingBits = partToSquare.bitLength();
2236 if (remainingBits == 1) { // Nothing left but +/- 1?
2237 if (signum < 0 && (exponent&1) == 1) {
2238 return NEGATIVE_ONE.shiftLeft(powersOfTwo*exponent);
2239 } else {
2240 return ONE.shiftLeft(powersOfTwo*exponent);
2241 }
2242 }
2243 } else {
2244 remainingBits = partToSquare.bitLength();
2245 if (remainingBits == 1) { // Nothing left but +/- 1?
2246 if (signum < 0 && (exponent&1) == 1) {
2247 return NEGATIVE_ONE;
2248 } else {
2249 return ONE;
2250 }
2251 }
2252 }
2253
2254 // This is a quick way to approximate the size of the result,
2255 // similar to doing log2[n] * exponent. This will give an upper bound
2256 // of how big the result can be, and which algorithm to use.
2257 long scaleFactor = (long)remainingBits * exponent;
2258
2259 // Use slightly different algorithms for small and large operands.
2260 // See if the result will safely fit into a long. (Largest 2^63-1)
2261 if (partToSquare.mag.length == 1 && scaleFactor <= 62) {
2262 // Small number algorithm. Everything fits into a long.
2263 int newSign = (signum <0 && (exponent&1) == 1 ? -1 : 1);
2264 long result = 1;
2265 long baseToPow2 = partToSquare.mag[0] & LONG_MASK;
2266
2267 int workingExponent = exponent;
2268
2269 // Perform exponentiation using repeated squaring trick
2270 while (workingExponent != 0) {
2271 if ((workingExponent & 1) == 1) {
2272 result = result * baseToPow2;
2273 }
2274
2275 if ((workingExponent >>>= 1) != 0) {
2276 baseToPow2 = baseToPow2 * baseToPow2;
2277 }
2278 }
2279
2280 // Multiply back the powers of two (quickly, by shifting left)
2281 if (powersOfTwo > 0) {
2282 if (bitsToShift + scaleFactor <= 62) { // Fits in long?
2283 return valueOf((result << bitsToShift) * newSign);
2284 } else {
2285 return valueOf(result*newSign).shiftLeft((int) bitsToShift);
2286 }
2287 }
2288 else {
2289 return valueOf(result*newSign);
2290 }
2291 } else {
2292 // Large number algorithm. This is basically identical to
2293 // the algorithm above, but calls multiply() and square()
2294 // which may use more efficient algorithms for large numbers.
2295 BigInteger answer = ONE;
2296
2297 int workingExponent = exponent;
2298 // Perform exponentiation using repeated squaring trick
2299 while (workingExponent != 0) {
2300 if ((workingExponent & 1) == 1) {
2301 answer = answer.multiply(partToSquare);
2302 }
2303
2304 if ((workingExponent >>>= 1) != 0) {
2305 partToSquare = partToSquare.square();
2306 }
2307 }
2308 // Multiply back the (exponentiated) powers of two (quickly,
2309 // by shifting left)
2310 if (powersOfTwo > 0) {
2311 answer = answer.shiftLeft(powersOfTwo*exponent);
2312 }
2313
2314 if (signum < 0 && (exponent&1) == 1) {
2315 return answer.negate();
2316 } else {
2317 return answer;
2318 }
2319 }
2320 }
2321
2322 /**
2323 * Returns a BigInteger whose value is the greatest common divisor of
2324 * {@code abs(this)} and {@code abs(val)}. Returns 0 if
2325 * {@code this == 0 && val == 0}.
2326 *
2327 * @param val value with which the GCD is to be computed.
2328 * @return {@code GCD(abs(this), abs(val))}
2329 */
2330 public BigInteger gcd(BigInteger val) {
2331 if (val.signum == 0)
2332 return this.abs();
2333 else if (this.signum == 0)
2334 return val.abs();
2335
2336 MutableBigInteger a = new MutableBigInteger(this);
2337 MutableBigInteger b = new MutableBigInteger(val);
2338
2339 MutableBigInteger result = a.hybridGCD(b);
2340
2341 return result.toBigInteger(1);
2342 }
2343
2344 /**
2345 * Package private method to return bit length for an integer.
2346 */
2347 static int bitLengthForInt(int n) {
2348 return 32 - Integer.numberOfLeadingZeros(n);
2349 }
2350
2351 /**
2352 * Left shift int array a up to len by n bits. Returns the array that
2353 * results from the shift since space may have to be reallocated.
2354 */
2355 private static int[] leftShift(int[] a, int len, int n) {
2356 int nInts = n >>> 5;
2357 int nBits = n&0x1F;
2358 int bitsInHighWord = bitLengthForInt(a[0]);
2359
2360 // If shift can be done without recopy, do so
2361 if (n <= (32-bitsInHighWord)) {
2362 primitiveLeftShift(a, len, nBits);
2363 return a;
2364 } else { // Array must be resized
2365 if (nBits <= (32-bitsInHighWord)) {
2366 int result[] = new int[nInts+len];
2367 System.arraycopy(a, 0, result, 0, len);
2368 primitiveLeftShift(result, result.length, nBits);
2369 return result;
2370 } else {
2371 int result[] = new int[nInts+len+1];
2372 System.arraycopy(a, 0, result, 0, len);
2373 primitiveRightShift(result, result.length, 32 - nBits);
2374 return result;
2375 }
2376 }
2377 }
2378
2379 // shifts a up to len right n bits assumes no leading zeros, 0<n<32
2380 static void primitiveRightShift(int[] a, int len, int n) {
2381 int n2 = 32 - n;
2382 for (int i=len-1, c=a[i]; i > 0; i--) {
2383 int b = c;
2384 c = a[i-1];
2385 a[i] = (c << n2) | (b >>> n);
2386 }
2387 a[0] >>>= n;
2388 }
2389
2390 // shifts a up to len left n bits assumes no leading zeros, 0<=n<32
2391 static void primitiveLeftShift(int[] a, int len, int n) {
2392 if (len == 0 || n == 0)
2393 return;
2394
2395 int n2 = 32 - n;
2396 for (int i=0, c=a[i], m=i+len-1; i < m; i++) {
2397 int b = c;
2398 c = a[i+1];
2399 a[i] = (b << n) | (c >>> n2);
2400 }
2401 a[len-1] <<= n;
2402 }
2403
2404 /**
2405 * Calculate bitlength of contents of the first len elements an int array,
2406 * assuming there are no leading zero ints.
2407 */
2408 private static int bitLength(int[] val, int len) {
2409 if (len == 0)
2410 return 0;
2411 return ((len - 1) << 5) + bitLengthForInt(val[0]);
2412 }
2413
2414 /**
2415 * Returns a BigInteger whose value is the absolute value of this
2416 * BigInteger.
2417 *
2418 * @return {@code abs(this)}
2419 */
2420 public BigInteger abs() {
2421 return (signum >= 0 ? this : this.negate());
2422 }
2423
2424 /**
2425 * Returns a BigInteger whose value is {@code (-this)}.
2426 *
2427 * @return {@code -this}
2428 */
2429 public BigInteger negate() {
2430 return new BigInteger(this.mag, -this.signum);
2431 }
2432
2433 /**
2434 * Returns the signum function of this BigInteger.
2435 *
2436 * @return -1, 0 or 1 as the value of this BigInteger is negative, zero or
2437 * positive.
2438 */
2439 public int signum() {
2440 return this.signum;
2441 }
2442
2443 // Modular Arithmetic Operations
2444
2445 /**
2446 * Returns a BigInteger whose value is {@code (this mod m}). This method
2447 * differs from {@code remainder} in that it always returns a
2448 * <i>non-negative</i> BigInteger.
2449 *
2450 * @param m the modulus.
2451 * @return {@code this mod m}
2452 * @throws ArithmeticException {@code m} ≤ 0
2453 * @see #remainder
2454 */
2455 public BigInteger mod(BigInteger m) {
2456 if (m.signum <= 0)
2457 throw new ArithmeticException("BigInteger: modulus not positive");
2458
2459 BigInteger result = this.remainder(m);
2460 return (result.signum >= 0 ? result : result.add(m));
2461 }
2462
2463 /**
2464 * Returns a BigInteger whose value is
2465 * <tt>(this<sup>exponent</sup> mod m)</tt>. (Unlike {@code pow}, this
2466 * method permits negative exponents.)
2467 *
2468 * @param exponent the exponent.
2469 * @param m the modulus.
2470 * @return <tt>this<sup>exponent</sup> mod m</tt>
2471 * @throws ArithmeticException {@code m} ≤ 0 or the exponent is
2472 * negative and this BigInteger is not <i>relatively
2473 * prime</i> to {@code m}.
2474 * @see #modInverse
2475 */
2476 public BigInteger modPow(BigInteger exponent, BigInteger m) {
2477 if (m.signum <= 0)
2478 throw new ArithmeticException("BigInteger: modulus not positive");
2479
2480 // Trivial cases
2481 if (exponent.signum == 0)
2482 return (m.equals(ONE) ? ZERO : ONE);
2483
2484 if (this.equals(ONE))
2485 return (m.equals(ONE) ? ZERO : ONE);
2486
2487 if (this.equals(ZERO) && exponent.signum >= 0)
2488 return ZERO;
2489
2490 if (this.equals(negConst[1]) && (!exponent.testBit(0)))
2491 return (m.equals(ONE) ? ZERO : ONE);
2492
2493 boolean invertResult;
2494 if ((invertResult = (exponent.signum < 0)))
2495 exponent = exponent.negate();
2496
2497 BigInteger base = (this.signum < 0 || this.compareTo(m) >= 0
2498 ? this.mod(m) : this);
2499 BigInteger result;
2500 if (m.testBit(0)) { // odd modulus
2501 result = base.oddModPow(exponent, m);
2502 } else {
2503 /*
2504 * Even modulus. Tear it into an "odd part" (m1) and power of two
2505 * (m2), exponentiate mod m1, manually exponentiate mod m2, and
2506 * use Chinese Remainder Theorem to combine results.
2507 */
2508
2509 // Tear m apart into odd part (m1) and power of 2 (m2)
2510 int p = m.getLowestSetBit(); // Max pow of 2 that divides m
2511
2512 BigInteger m1 = m.shiftRight(p); // m/2**p
2513 BigInteger m2 = ONE.shiftLeft(p); // 2**p
2514
2515 // Calculate new base from m1
2516 BigInteger base2 = (this.signum < 0 || this.compareTo(m1) >= 0
2517 ? this.mod(m1) : this);
2518
2519 // Caculate (base ** exponent) mod m1.
2520 BigInteger a1 = (m1.equals(ONE) ? ZERO :
2521 base2.oddModPow(exponent, m1));
2522
2523 // Calculate (this ** exponent) mod m2
2524 BigInteger a2 = base.modPow2(exponent, p);
2525
2526 // Combine results using Chinese Remainder Theorem
2527 BigInteger y1 = m2.modInverse(m1);
2528 BigInteger y2 = m1.modInverse(m2);
2529
2530 if (m.mag.length < MAX_MAG_LENGTH / 2) {
2531 result = a1.multiply(m2).multiply(y1).add(a2.multiply(m1).multiply(y2)).mod(m);
2532 } else {
2533 MutableBigInteger t1 = new MutableBigInteger();
2534 new MutableBigInteger(a1.multiply(m2)).multiply(new MutableBigInteger(y1), t1);
2535 MutableBigInteger t2 = new MutableBigInteger();
2536 new MutableBigInteger(a2.multiply(m1)).multiply(new MutableBigInteger(y2), t2);
2537 t1.add(t2);
2538 MutableBigInteger q = new MutableBigInteger();
2539 result = t1.divide(new MutableBigInteger(m), q).toBigInteger();
2540 }
2541 }
2542
2543 return (invertResult ? result.modInverse(m) : result);
2544 }
2545
2546 static int[] bnExpModThreshTable = {7, 25, 81, 241, 673, 1793,
2547 Integer.MAX_VALUE}; // Sentinel
2548
2549 /**
2550 * Returns a BigInteger whose value is x to the power of y mod z.
2551 * Assumes: z is odd && x < z.
2552 */
2553 private BigInteger oddModPow(BigInteger y, BigInteger z) {
2554 /*
2555 * The algorithm is adapted from Colin Plumb's C library.
2556 *
2557 * The window algorithm:
2558 * The idea is to keep a running product of b1 = n^(high-order bits of exp)
2559 * and then keep appending exponent bits to it. The following patterns
2560 * apply to a 3-bit window (k = 3):
2561 * To append 0: square
2562 * To append 1: square, multiply by n^1
2563 * To append 10: square, multiply by n^1, square
2564 * To append 11: square, square, multiply by n^3
2565 * To append 100: square, multiply by n^1, square, square
2566 * To append 101: square, square, square, multiply by n^5
2567 * To append 110: square, square, multiply by n^3, square
2568 * To append 111: square, square, square, multiply by n^7
2569 *
2570 * Since each pattern involves only one multiply, the longer the pattern
2571 * the better, except that a 0 (no multiplies) can be appended directly.
2572 * We precompute a table of odd powers of n, up to 2^k, and can then
2573 * multiply k bits of exponent at a time. Actually, assuming random
2574 * exponents, there is on average one zero bit between needs to
2575 * multiply (1/2 of the time there's none, 1/4 of the time there's 1,
2576 * 1/8 of the time, there's 2, 1/32 of the time, there's 3, etc.), so
2577 * you have to do one multiply per k+1 bits of exponent.
2578 *
2579 * The loop walks down the exponent, squaring the result buffer as
2580 * it goes. There is a wbits+1 bit lookahead buffer, buf, that is
2581 * filled with the upcoming exponent bits. (What is read after the
2582 * end of the exponent is unimportant, but it is filled with zero here.)
2583 * When the most-significant bit of this buffer becomes set, i.e.
2584 * (buf & tblmask) != 0, we have to decide what pattern to multiply
2585 * by, and when to do it. We decide, remember to do it in future
2586 * after a suitable number of squarings have passed (e.g. a pattern
2587 * of "100" in the buffer requires that we multiply by n^1 immediately;
2588 * a pattern of "110" calls for multiplying by n^3 after one more
2589 * squaring), clear the buffer, and continue.
2590 *
2591 * When we start, there is one more optimization: the result buffer
2592 * is implcitly one, so squaring it or multiplying by it can be
2593 * optimized away. Further, if we start with a pattern like "100"
2594 * in the lookahead window, rather than placing n into the buffer
2595 * and then starting to square it, we have already computed n^2
2596 * to compute the odd-powers table, so we can place that into
2597 * the buffer and save a squaring.
2598 *
2599 * This means that if you have a k-bit window, to compute n^z,
2600 * where z is the high k bits of the exponent, 1/2 of the time
2601 * it requires no squarings. 1/4 of the time, it requires 1
2602 * squaring, ... 1/2^(k-1) of the time, it reqires k-2 squarings.
2603 * And the remaining 1/2^(k-1) of the time, the top k bits are a
2604 * 1 followed by k-1 0 bits, so it again only requires k-2
2605 * squarings, not k-1. The average of these is 1. Add that
2606 * to the one squaring we have to do to compute the table,
2607 * and you'll see that a k-bit window saves k-2 squarings
2608 * as well as reducing the multiplies. (It actually doesn't
2609 * hurt in the case k = 1, either.)
2610 */
2611 // Special case for exponent of one
2612 if (y.equals(ONE))
2613 return this;
2614
2615 // Special case for base of zero
2616 if (signum == 0)
2617 return ZERO;
2618
2619 int[] base = mag.clone();
2620 int[] exp = y.mag;
2621 int[] mod = z.mag;
2622 int modLen = mod.length;
2623
2624 // Select an appropriate window size
2625 int wbits = 0;
2626 int ebits = bitLength(exp, exp.length);
2627 // if exponent is 65537 (0x10001), use minimum window size
2628 if ((ebits != 17) || (exp[0] != 65537)) {
2629 while (ebits > bnExpModThreshTable[wbits]) {
2630 wbits++;
2631 }
2632 }
2633
2634 // Calculate appropriate table size
2635 int tblmask = 1 << wbits;
2636
2637 // Allocate table for precomputed odd powers of base in Montgomery form
2638 int[][] table = new int[tblmask][];
2639 for (int i=0; i < tblmask; i++)
2640 table[i] = new int[modLen];
2641
2642 // Compute the modular inverse
2643 int inv = -MutableBigInteger.inverseMod32(mod[modLen-1]);
2644
2645 // Convert base to Montgomery form
2646 int[] a = leftShift(base, base.length, modLen << 5);
2647
2648 MutableBigInteger q = new MutableBigInteger(),
2649 a2 = new MutableBigInteger(a),
2650 b2 = new MutableBigInteger(mod);
2651
2652 MutableBigInteger r= a2.divide(b2, q);
2653 table[0] = r.toIntArray();
2654
2655 // Pad table[0] with leading zeros so its length is at least modLen
2656 if (table[0].length < modLen) {
2657 int offset = modLen - table[0].length;
2658 int[] t2 = new int[modLen];
2659 for (int i=0; i < table[0].length; i++)
2660 t2[i+offset] = table[0][i];
2661 table[0] = t2;
2662 }
2663
2664 // Set b to the square of the base
2665 int[] b = squareToLen(table[0], modLen, null);
2666 b = montReduce(b, mod, modLen, inv);
2667
2668 // Set t to high half of b
2669 int[] t = Arrays.copyOf(b, modLen);
2670
2671 // Fill in the table with odd powers of the base
2672 for (int i=1; i < tblmask; i++) {
2673 int[] prod = multiplyToLen(t, modLen, table[i-1], modLen, null);
2674 table[i] = montReduce(prod, mod, modLen, inv);
2675 }
2676
2677 // Pre load the window that slides over the exponent
2678 int bitpos = 1 << ((ebits-1) & (32-1));
2679
2680 int buf = 0;
2681 int elen = exp.length;
2682 int eIndex = 0;
2683 for (int i = 0; i <= wbits; i++) {
2684 buf = (buf << 1) | (((exp[eIndex] & bitpos) != 0)?1:0);
2685 bitpos >>>= 1;
2686 if (bitpos == 0) {
2687 eIndex++;
2688 bitpos = 1 << (32-1);
2689 elen--;
2690 }
2691 }
2692
2693 int multpos = ebits;
2694
2695 // The first iteration, which is hoisted out of the main loop
2696 ebits--;
2697 boolean isone = true;
2698
2699 multpos = ebits - wbits;
2700 while ((buf & 1) == 0) {
2701 buf >>>= 1;
2702 multpos++;
2703 }
2704
2705 int[] mult = table[buf >>> 1];
2706
2707 buf = 0;
2708 if (multpos == ebits)
2709 isone = false;
2710
2711 // The main loop
2712 while (true) {
2713 ebits--;
2714 // Advance the window
2715 buf <<= 1;
2716
2717 if (elen != 0) {
2718 buf |= ((exp[eIndex] & bitpos) != 0) ? 1 : 0;
2719 bitpos >>>= 1;
2720 if (bitpos == 0) {
2721 eIndex++;
2722 bitpos = 1 << (32-1);
2723 elen--;
2724 }
2725 }
2726
2727 // Examine the window for pending multiplies
2728 if ((buf & tblmask) != 0) {
2729 multpos = ebits - wbits;
2730 while ((buf & 1) == 0) {
2731 buf >>>= 1;
2732 multpos++;
2733 }
2734 mult = table[buf >>> 1];
2735 buf = 0;
2736 }
2737
2738 // Perform multiply
2739 if (ebits == multpos) {
2740 if (isone) {
2741 b = mult.clone();
2742 isone = false;
2743 } else {
2744 t = b;
2745 a = multiplyToLen(t, modLen, mult, modLen, a);
2746 a = montReduce(a, mod, modLen, inv);
2747 t = a; a = b; b = t;
2748 }
2749 }
2750
2751 // Check if done
2752 if (ebits == 0)
2753 break;
2754
2755 // Square the input
2756 if (!isone) {
2757 t = b;
2758 a = squareToLen(t, modLen, a);
2759 a = montReduce(a, mod, modLen, inv);
2760 t = a; a = b; b = t;
2761 }
2762 }
2763
2764 // Convert result out of Montgomery form and return
2765 int[] t2 = new int[2*modLen];
2766 System.arraycopy(b, 0, t2, modLen, modLen);
2767
2768 b = montReduce(t2, mod, modLen, inv);
2769
2770 t2 = Arrays.copyOf(b, modLen);
2771
2772 return new BigInteger(1, t2);
2773 }
2774
2775 /**
2776 * Montgomery reduce n, modulo mod. This reduces modulo mod and divides
2777 * by 2^(32*mlen). Adapted from Colin Plumb's C library.
2778 */
2779 private static int[] montReduce(int[] n, int[] mod, int mlen, int inv) {
2780 int c=0;
2781 int len = mlen;
2782 int offset=0;
2783
2784 do {
2785 int nEnd = n[n.length-1-offset];
2786 int carry = mulAdd(n, mod, offset, mlen, inv * nEnd);
2787 c += addOne(n, offset, mlen, carry);
2788 offset++;
2789 } while (--len > 0);
2790
2791 while (c > 0)
2792 c += subN(n, mod, mlen);
2793
2794 while (intArrayCmpToLen(n, mod, mlen) >= 0)
2795 subN(n, mod, mlen);
2796
2797 return n;
2798 }
2799
2800
2801 /*
2802 * Returns -1, 0 or +1 as big-endian unsigned int array arg1 is less than,
2803 * equal to, or greater than arg2 up to length len.
2804 */
2805 private static int intArrayCmpToLen(int[] arg1, int[] arg2, int len) {
2806 for (int i=0; i < len; i++) {
2807 long b1 = arg1[i] & LONG_MASK;
2808 long b2 = arg2[i] & LONG_MASK;
2809 if (b1 < b2)
2810 return -1;
2811 if (b1 > b2)
2812 return 1;
2813 }
2814 return 0;
2815 }
2816
2817 /**
2818 * Subtracts two numbers of same length, returning borrow.
2819 */
2820 private static int subN(int[] a, int[] b, int len) {
2821 long sum = 0;
2822
2823 while (--len >= 0) {
2824 sum = (a[len] & LONG_MASK) -
2825 (b[len] & LONG_MASK) + (sum >> 32);
2826 a[len] = (int)sum;
2827 }
2828
2829 return (int)(sum >> 32);
2830 }
2831
2832 /**
2833 * Multiply an array by one word k and add to result, return the carry
2834 */
2835 static int mulAdd(int[] out, int[] in, int offset, int len, int k) {
2836 long kLong = k & LONG_MASK;
2837 long carry = 0;
2838
2839 offset = out.length-offset - 1;
2840 for (int j=len-1; j >= 0; j--) {
2841 long product = (in[j] & LONG_MASK) * kLong +
2842 (out[offset] & LONG_MASK) + carry;
2843 out[offset--] = (int)product;
2844 carry = product >>> 32;
2845 }
2846 return (int)carry;
2847 }
2848
2849 /**
2850 * Add one word to the number a mlen words into a. Return the resulting
2851 * carry.
2852 */
2853 static int addOne(int[] a, int offset, int mlen, int carry) {
2854 offset = a.length-1-mlen-offset;
2855 long t = (a[offset] & LONG_MASK) + (carry & LONG_MASK);
2856
2857 a[offset] = (int)t;
2858 if ((t >>> 32) == 0)
2859 return 0;
2860 while (--mlen >= 0) {
2861 if (--offset < 0) { // Carry out of number
2862 return 1;
2863 } else {
2864 a[offset]++;
2865 if (a[offset] != 0)
2866 return 0;
2867 }
2868 }
2869 return 1;
2870 }
2871
2872 /**
2873 * Returns a BigInteger whose value is (this ** exponent) mod (2**p)
2874 */
2875 private BigInteger modPow2(BigInteger exponent, int p) {
2876 /*
2877 * Perform exponentiation using repeated squaring trick, chopping off
2878 * high order bits as indicated by modulus.
2879 */
2880 BigInteger result = ONE;
2881 BigInteger baseToPow2 = this.mod2(p);
2882 int expOffset = 0;
2883
2884 int limit = exponent.bitLength();
2885
2886 if (this.testBit(0))
2887 limit = (p-1) < limit ? (p-1) : limit;
2888
2889 while (expOffset < limit) {
2890 if (exponent.testBit(expOffset))
2891 result = result.multiply(baseToPow2).mod2(p);
2892 expOffset++;
2893 if (expOffset < limit)
2894 baseToPow2 = baseToPow2.square().mod2(p);
2895 }
2896
2897 return result;
2898 }
2899
2900 /**
2901 * Returns a BigInteger whose value is this mod(2**p).
2902 * Assumes that this {@code BigInteger >= 0} and {@code p > 0}.
2903 */
2904 private BigInteger mod2(int p) {
2905 if (bitLength() <= p)
2906 return this;
2907
2908 // Copy remaining ints of mag
2909 int numInts = (p + 31) >>> 5;
2910 int[] mag = new int[numInts];
2911 System.arraycopy(this.mag, (this.mag.length - numInts), mag, 0, numInts);
2912
2913 // Mask out any excess bits
2914 int excessBits = (numInts << 5) - p;
2915 mag[0] &= (1L << (32-excessBits)) - 1;
2916
2917 return (mag[0] == 0 ? new BigInteger(1, mag) : new BigInteger(mag, 1));
2918 }
2919
2920 /**
2921 * Returns a BigInteger whose value is {@code (this}<sup>-1</sup> {@code mod m)}.
2922 *
2923 * @param m the modulus.
2924 * @return {@code this}<sup>-1</sup> {@code mod m}.
2925 * @throws ArithmeticException {@code m} ≤ 0, or this BigInteger
2926 * has no multiplicative inverse mod m (that is, this BigInteger
2927 * is not <i>relatively prime</i> to m).
2928 */
2929 public BigInteger modInverse(BigInteger m) {
2930 if (m.signum != 1)
2931 throw new ArithmeticException("BigInteger: modulus not positive");
2932
2933 if (m.equals(ONE))
2934 return ZERO;
2935
2936 // Calculate (this mod m)
2937 BigInteger modVal = this;
2938 if (signum < 0 || (this.compareMagnitude(m) >= 0))
2939 modVal = this.mod(m);
2940
2941 if (modVal.equals(ONE))
2942 return ONE;
2943
2944 MutableBigInteger a = new MutableBigInteger(modVal);
2945 MutableBigInteger b = new MutableBigInteger(m);
2946
2947 MutableBigInteger result = a.mutableModInverse(b);
2948 return result.toBigInteger(1);
2949 }
2950
2951 // Shift Operations
2952
2953 /**
2954 * Returns a BigInteger whose value is {@code (this << n)}.
2955 * The shift distance, {@code n}, may be negative, in which case
2956 * this method performs a right shift.
2957 * (Computes <tt>floor(this * 2<sup>n</sup>)</tt>.)
2958 *
2959 * @param n shift distance, in bits.
2960 * @return {@code this << n}
2961 * @see #shiftRight
2962 */
2963 public BigInteger shiftLeft(int n) {
2964 if (signum == 0)
2965 return ZERO;
2966 if (n > 0) {
2967 return new BigInteger(shiftLeft(mag, n), signum);
2968 } else if (n == 0) {
2969 return this;
2970 } else {
2971 // Possible int overflow in (-n) is not a trouble,
2972 // because shiftRightImpl considers its argument unsigned
2973 return shiftRightImpl(-n);
2974 }
2975 }
2976
2977 /**
2978 * Returns a magnitude array whose value is {@code (mag << n)}.
2979 * The shift distance, {@code n}, is considered unnsigned.
2980 * (Computes <tt>this * 2<sup>n</sup></tt>.)
2981 *
2982 * @param mag magnitude, the most-significant int ({@code mag[0]}) must be non-zero.
2983 * @param n unsigned shift distance, in bits.
2984 * @return {@code mag << n}
2985 */
2986 private static int[] shiftLeft(int[] mag, int n) {
2987 int nInts = n >>> 5;
2988 int nBits = n & 0x1f;
2989 int magLen = mag.length;
2990 int newMag[] = null;
2991
2992 if (nBits == 0) {
2993 newMag = new int[magLen + nInts];
2994 System.arraycopy(mag, 0, newMag, 0, magLen);
2995 } else {
2996 int i = 0;
2997 int nBits2 = 32 - nBits;
2998 int highBits = mag[0] >>> nBits2;
2999 if (highBits != 0) {
3000 newMag = new int[magLen + nInts + 1];
3001 newMag[i++] = highBits;
3002 } else {
3003 newMag = new int[magLen + nInts];
3004 }
3005 int j=0;
3006 while (j < magLen-1)
3007 newMag[i++] = mag[j++] << nBits | mag[j] >>> nBits2;
3008 newMag[i] = mag[j] << nBits;
3009 }
3010 return newMag;
3011 }
3012
3013 /**
3014 * Returns a BigInteger whose value is {@code (this >> n)}. Sign
3015 * extension is performed. The shift distance, {@code n}, may be
3016 * negative, in which case this method performs a left shift.
3017 * (Computes <tt>floor(this / 2<sup>n</sup>)</tt>.)
3018 *
3019 * @param n shift distance, in bits.
3020 * @return {@code this >> n}
3021 * @see #shiftLeft
3022 */
3023 public BigInteger shiftRight(int n) {
3024 if (signum == 0)
3025 return ZERO;
3026 if (n > 0) {
3027 return shiftRightImpl(n);
3028 } else if (n == 0) {
3029 return this;
3030 } else {
3031 // Possible int overflow in {@code -n} is not a trouble,
3032 // because shiftLeft considers its argument unsigned
3033 return new BigInteger(shiftLeft(mag, -n), signum);
3034 }
3035 }
3036
3037 /**
3038 * Returns a BigInteger whose value is {@code (this >> n)}. The shift
3039 * distance, {@code n}, is considered unsigned.
3040 * (Computes <tt>floor(this * 2<sup>-n</sup>)</tt>.)
3041 *
3042 * @param n unsigned shift distance, in bits.
3043 * @return {@code this >> n}
3044 */
3045 private BigInteger shiftRightImpl(int n) {
3046 int nInts = n >>> 5;
3047 int nBits = n & 0x1f;
3048 int magLen = mag.length;
3049 int newMag[] = null;
3050
3051 // Special case: entire contents shifted off the end
3052 if (nInts >= magLen)
3053 return (signum >= 0 ? ZERO : negConst[1]);
3054
3055 if (nBits == 0) {
3056 int newMagLen = magLen - nInts;
3057 newMag = Arrays.copyOf(mag, newMagLen);
3058 } else {
3059 int i = 0;
3060 int highBits = mag[0] >>> nBits;
3061 if (highBits != 0) {
3062 newMag = new int[magLen - nInts];
3063 newMag[i++] = highBits;
3064 } else {
3065 newMag = new int[magLen - nInts -1];
3066 }
3067
3068 int nBits2 = 32 - nBits;
3069 int j=0;
3070 while (j < magLen - nInts - 1)
3071 newMag[i++] = (mag[j++] << nBits2) | (mag[j] >>> nBits);
3072 }
3073
3074 if (signum < 0) {
3075 // Find out whether any one-bits were shifted off the end.
3076 boolean onesLost = false;
3077 for (int i=magLen-1, j=magLen-nInts; i >= j && !onesLost; i--)
3078 onesLost = (mag[i] != 0);
3079 if (!onesLost && nBits != 0)
3080 onesLost = (mag[magLen - nInts - 1] << (32 - nBits) != 0);
3081
3082 if (onesLost)
3083 newMag = javaIncrement(newMag);
3084 }
3085
3086 return new BigInteger(newMag, signum);
3087 }
3088
3089 int[] javaIncrement(int[] val) {
3090 int lastSum = 0;
3091 for (int i=val.length-1; i >= 0 && lastSum == 0; i--)
3092 lastSum = (val[i] += 1);
3093 if (lastSum == 0) {
3094 val = new int[val.length+1];
3095 val[0] = 1;
3096 }
3097 return val;
3098 }
3099
3100 // Bitwise Operations
3101
3102 /**
3103 * Returns a BigInteger whose value is {@code (this & val)}. (This
3104 * method returns a negative BigInteger if and only if this and val are
3105 * both negative.)
3106 *
3107 * @param val value to be AND'ed with this BigInteger.
3108 * @return {@code this & val}
3109 */
3110 public BigInteger and(BigInteger val) {
3111 int[] result = new int[Math.max(intLength(), val.intLength())];
3112 for (int i=0; i < result.length; i++)
3113 result[i] = (getInt(result.length-i-1)
3114 & val.getInt(result.length-i-1));
3115
3116 return valueOf(result);
3117 }
3118
3119 /**
3120 * Returns a BigInteger whose value is {@code (this | val)}. (This method
3121 * returns a negative BigInteger if and only if either this or val is
3122 * negative.)
3123 *
3124 * @param val value to be OR'ed with this BigInteger.
3125 * @return {@code this | val}
3126 */
3127 public BigInteger or(BigInteger val) {
3128 int[] result = new int[Math.max(intLength(), val.intLength())];
3129 for (int i=0; i < result.length; i++)
3130 result[i] = (getInt(result.length-i-1)
3131 | val.getInt(result.length-i-1));
3132
3133 return valueOf(result);
3134 }
3135
3136 /**
3137 * Returns a BigInteger whose value is {@code (this ^ val)}. (This method
3138 * returns a negative BigInteger if and only if exactly one of this and
3139 * val are negative.)
3140 *
3141 * @param val value to be XOR'ed with this BigInteger.
3142 * @return {@code this ^ val}
3143 */
3144 public BigInteger xor(BigInteger val) {
3145 int[] result = new int[Math.max(intLength(), val.intLength())];
3146 for (int i=0; i < result.length; i++)
3147 result[i] = (getInt(result.length-i-1)
3148 ^ val.getInt(result.length-i-1));
3149
3150 return valueOf(result);
3151 }
3152
3153 /**
3154 * Returns a BigInteger whose value is {@code (~this)}. (This method
3155 * returns a negative value if and only if this BigInteger is
3156 * non-negative.)
3157 *
3158 * @return {@code ~this}
3159 */
3160 public BigInteger not() {
3161 int[] result = new int[intLength()];
3162 for (int i=0; i < result.length; i++)
3163 result[i] = ~getInt(result.length-i-1);
3164
3165 return valueOf(result);
3166 }
3167
3168 /**
3169 * Returns a BigInteger whose value is {@code (this & ~val)}. This
3170 * method, which is equivalent to {@code and(val.not())}, is provided as
3171 * a convenience for masking operations. (This method returns a negative
3172 * BigInteger if and only if {@code this} is negative and {@code val} is
3173 * positive.)
3174 *
3175 * @param val value to be complemented and AND'ed with this BigInteger.
3176 * @return {@code this & ~val}
3177 */
3178 public BigInteger andNot(BigInteger val) {
3179 int[] result = new int[Math.max(intLength(), val.intLength())];
3180 for (int i=0; i < result.length; i++)
3181 result[i] = (getInt(result.length-i-1)
3182 & ~val.getInt(result.length-i-1));
3183
3184 return valueOf(result);
3185 }
3186
3187
3188 // Single Bit Operations
3189
3190 /**
3191 * Returns {@code true} if and only if the designated bit is set.
3192 * (Computes {@code ((this & (1<<n)) != 0)}.)
3193 *
3194 * @param n index of bit to test.
3195 * @return {@code true} if and only if the designated bit is set.
3196 * @throws ArithmeticException {@code n} is negative.
3197 */
3198 public boolean testBit(int n) {
3199 if (n < 0)
3200 throw new ArithmeticException("Negative bit address");
3201
3202 return (getInt(n >>> 5) & (1 << (n & 31))) != 0;
3203 }
3204
3205 /**
3206 * Returns a BigInteger whose value is equivalent to this BigInteger
3207 * with the designated bit set. (Computes {@code (this | (1<<n))}.)
3208 *
3209 * @param n index of bit to set.
3210 * @return {@code this | (1<<n)}
3211 * @throws ArithmeticException {@code n} is negative.
3212 */
3213 public BigInteger setBit(int n) {
3214 if (n < 0)
3215 throw new ArithmeticException("Negative bit address");
3216
3217 int intNum = n >>> 5;
3218 int[] result = new int[Math.max(intLength(), intNum+2)];
3219
3220 for (int i=0; i < result.length; i++)
3221 result[result.length-i-1] = getInt(i);
3222
3223 result[result.length-intNum-1] |= (1 << (n & 31));
3224
3225 return valueOf(result);
3226 }
3227
3228 /**
3229 * Returns a BigInteger whose value is equivalent to this BigInteger
3230 * with the designated bit cleared.
3231 * (Computes {@code (this & ~(1<<n))}.)
3232 *
3233 * @param n index of bit to clear.
3234 * @return {@code this & ~(1<<n)}
3235 * @throws ArithmeticException {@code n} is negative.
3236 */
3237 public BigInteger clearBit(int n) {
3238 if (n < 0)
3239 throw new ArithmeticException("Negative bit address");
3240
3241 int intNum = n >>> 5;
3242 int[] result = new int[Math.max(intLength(), ((n + 1) >>> 5) + 1)];
3243
3244 for (int i=0; i < result.length; i++)
3245 result[result.length-i-1] = getInt(i);
3246
3247 result[result.length-intNum-1] &= ~(1 << (n & 31));
3248
3249 return valueOf(result);
3250 }
3251
3252 /**
3253 * Returns a BigInteger whose value is equivalent to this BigInteger
3254 * with the designated bit flipped.
3255 * (Computes {@code (this ^ (1<<n))}.)
3256 *
3257 * @param n index of bit to flip.
3258 * @return {@code this ^ (1<<n)}
3259 * @throws ArithmeticException {@code n} is negative.
3260 */
3261 public BigInteger flipBit(int n) {
3262 if (n < 0)
3263 throw new ArithmeticException("Negative bit address");
3264
3265 int intNum = n >>> 5;
3266 int[] result = new int[Math.max(intLength(), intNum+2)];
3267
3268 for (int i=0; i < result.length; i++)
3269 result[result.length-i-1] = getInt(i);
3270
3271 result[result.length-intNum-1] ^= (1 << (n & 31));
3272
3273 return valueOf(result);
3274 }
3275
3276 /**
3277 * Returns the index of the rightmost (lowest-order) one bit in this
3278 * BigInteger (the number of zero bits to the right of the rightmost
3279 * one bit). Returns -1 if this BigInteger contains no one bits.
3280 * (Computes {@code (this == 0? -1 : log2(this & -this))}.)
3281 *
3282 * @return index of the rightmost one bit in this BigInteger.
3283 */
3284 public int getLowestSetBit() {
3285 int lsb = lowestSetBitPlusTwo - 2;
3286 if (lsb == -2) { // lowestSetBit not initialized yet
3287 lsb = 0;
3288 if (signum == 0) {
3289 lsb -= 1;
3290 } else {
3291 // Search for lowest order nonzero int
3292 int i,b;
3293 for (i=0; (b = getInt(i)) == 0; i++)
3294 ;
3295 lsb += (i << 5) + Integer.numberOfTrailingZeros(b);
3296 }
3297 lowestSetBitPlusTwo = lsb + 2;
3298 }
3299 return lsb;
3300 }
3301
3302
3303 // Miscellaneous Bit Operations
3304
3305 /**
3306 * Returns the number of bits in the minimal two's-complement
3307 * representation of this BigInteger, <i>excluding</i> a sign bit.
3308 * For positive BigIntegers, this is equivalent to the number of bits in
3309 * the ordinary binary representation. (Computes
3310 * {@code (ceil(log2(this < 0 ? -this : this+1)))}.)
3311 *
3312 * @return number of bits in the minimal two's-complement
3313 * representation of this BigInteger, <i>excluding</i> a sign bit.
3314 */
3315 public int bitLength() {
3316 int n = bitLengthPlusOne - 1;
3317 if (n == -1) { // bitLength not initialized yet
3318 int[] m = mag;
3319 int len = m.length;
3320 if (len == 0) {
3321 n = 0; // offset by one to initialize
3322 } else {
3323 // Calculate the bit length of the magnitude
3324 int magBitLength = ((len - 1) << 5) + bitLengthForInt(mag[0]);
3325 if (signum < 0) {
3326 // Check if magnitude is a power of two
3327 boolean pow2 = (Integer.bitCount(mag[0]) == 1);
3328 for (int i=1; i< len && pow2; i++)
3329 pow2 = (mag[i] == 0);
3330
3331 n = (pow2 ? magBitLength -1 : magBitLength);
3332 } else {
3333 n = magBitLength;
3334 }
3335 }
3336 bitLengthPlusOne = n + 1;
3337 }
3338 return n;
3339 }
3340
3341 /**
3342 * Returns the number of bits in the two's complement representation
3343 * of this BigInteger that differ from its sign bit. This method is
3344 * useful when implementing bit-vector style sets atop BigIntegers.
3345 *
3346 * @return number of bits in the two's complement representation
3347 * of this BigInteger that differ from its sign bit.
3348 */
3349 public int bitCount() {
3350 int bc = bitCountPlusOne - 1;
3351 if (bc == -1) { // bitCount not initialized yet
3352 bc = 0; // offset by one to initialize
3353 // Count the bits in the magnitude
3354 for (int i=0; i < mag.length; i++)
3355 bc += Integer.bitCount(mag[i]);
3356 if (signum < 0) {
3357 // Count the trailing zeros in the magnitude
3358 int magTrailingZeroCount = 0, j;
3359 for (j=mag.length-1; mag[j] == 0; j--)
3360 magTrailingZeroCount += 32;
3361 magTrailingZeroCount += Integer.numberOfTrailingZeros(mag[j]);
3362 bc += magTrailingZeroCount - 1;
3363 }
3364 bitCountPlusOne = bc + 1;
3365 }
3366 return bc;
3367 }
3368
3369 // Primality Testing
3370
3371 /**
3372 * Returns {@code true} if this BigInteger is probably prime,
3373 * {@code false} if it's definitely composite. If
3374 * {@code certainty} is ≤ 0, {@code true} is
3375 * returned.
3376 *
3377 * @param certainty a measure of the uncertainty that the caller is
3378 * willing to tolerate: if the call returns {@code true}
3379 * the probability that this BigInteger is prime exceeds
3380 * (1 - 1/2<sup>{@code certainty}</sup>). The execution time of
3381 * this method is proportional to the value of this parameter.
3382 * @return {@code true} if this BigInteger is probably prime,
3383 * {@code false} if it's definitely composite.
3384 */
3385 public boolean isProbablePrime(int certainty) {
3386 if (certainty <= 0)
3387 return true;
3388 BigInteger w = this.abs();
3389 if (w.equals(TWO))
3390 return true;
3391 if (!w.testBit(0) || w.equals(ONE))
3392 return false;
3393
3394 return w.primeToCertainty(certainty, null);
3395 }
3396
3397 // Comparison Operations
3398
3399 /**
3400 * Compares this BigInteger with the specified BigInteger. This
3401 * method is provided in preference to individual methods for each
3402 * of the six boolean comparison operators ({@literal <}, ==,
3403 * {@literal >}, {@literal >=}, !=, {@literal <=}). The suggested
3404 * idiom for performing these comparisons is: {@code
3405 * (x.compareTo(y)} <<i>op</i>> {@code 0)}, where
3406 * <<i>op</i>> is one of the six comparison operators.
3407 *
3408 * @param val BigInteger to which this BigInteger is to be compared.
3409 * @return -1, 0 or 1 as this BigInteger is numerically less than, equal
3410 * to, or greater than {@code val}.
3411 */
3412 public int compareTo(BigInteger val) {
3413 if (signum == val.signum) {
3414 switch (signum) {
3415 case 1:
3416 return compareMagnitude(val);
3417 case -1:
3418 return val.compareMagnitude(this);
3419 default:
3420 return 0;
3421 }
3422 }
3423 return signum > val.signum ? 1 : -1;
3424 }
3425
3426 /**
3427 * Compares the magnitude array of this BigInteger with the specified
3428 * BigInteger's. This is the version of compareTo ignoring sign.
3429 *
3430 * @param val BigInteger whose magnitude array to be compared.
3431 * @return -1, 0 or 1 as this magnitude array is less than, equal to or
3432 * greater than the magnitude aray for the specified BigInteger's.
3433 */
3434 final int compareMagnitude(BigInteger val) {
3435 int[] m1 = mag;
3436 int len1 = m1.length;
3437 int[] m2 = val.mag;
3438 int len2 = m2.length;
3439 if (len1 < len2)
3440 return -1;
3441 if (len1 > len2)
3442 return 1;
3443 for (int i = 0; i < len1; i++) {
3444 int a = m1[i];
3445 int b = m2[i];
3446 if (a != b)
3447 return ((a & LONG_MASK) < (b & LONG_MASK)) ? -1 : 1;
3448 }
3449 return 0;
3450 }
3451
3452 /**
3453 * Version of compareMagnitude that compares magnitude with long value.
3454 * val can't be Long.MIN_VALUE.
3455 */
3456 final int compareMagnitude(long val) {
3457 assert val != Long.MIN_VALUE;
3458 int[] m1 = mag;
3459 int len = m1.length;
3460 if (len > 2) {
3461 return 1;
3462 }
3463 if (val < 0) {
3464 val = -val;
3465 }
3466 int highWord = (int)(val >>> 32);
3467 if (highWord == 0) {
3468 if (len < 1)
3469 return -1;
3470 if (len > 1)
3471 return 1;
3472 int a = m1[0];
3473 int b = (int)val;
3474 if (a != b) {
3475 return ((a & LONG_MASK) < (b & LONG_MASK))? -1 : 1;
3476 }
3477 return 0;
3478 } else {
3479 if (len < 2)
3480 return -1;
3481 int a = m1[0];
3482 int b = highWord;
3483 if (a != b) {
3484 return ((a & LONG_MASK) < (b & LONG_MASK))? -1 : 1;
3485 }
3486 a = m1[1];
3487 b = (int)val;
3488 if (a != b) {
3489 return ((a & LONG_MASK) < (b & LONG_MASK))? -1 : 1;
3490 }
3491 return 0;
3492 }
3493 }
3494
3495 /**
3496 * Compares this BigInteger with the specified Object for equality.
3497 *
3498 * @param x Object to which this BigInteger is to be compared.
3499 * @return {@code true} if and only if the specified Object is a
3500 * BigInteger whose value is numerically equal to this BigInteger.
3501 */
3502 public boolean equals(Object x) {
3503 // This test is just an optimization, which may or may not help
3504 if (x == this)
3505 return true;
3506
3507 if (!(x instanceof BigInteger))
3508 return false;
3509
3510 BigInteger xInt = (BigInteger) x;
3511 if (xInt.signum != signum)
3512 return false;
3513
3514 int[] m = mag;
3515 int len = m.length;
3516 int[] xm = xInt.mag;
3517 if (len != xm.length)
3518 return false;
3519
3520 for (int i = 0; i < len; i++)
3521 if (xm[i] != m[i])
3522 return false;
3523
3524 return true;
3525 }
3526
3527 /**
3528 * Returns the minimum of this BigInteger and {@code val}.
3529 *
3530 * @param val value with which the minimum is to be computed.
3531 * @return the BigInteger whose value is the lesser of this BigInteger and
3532 * {@code val}. If they are equal, either may be returned.
3533 */
3534 public BigInteger min(BigInteger val) {
3535 return (compareTo(val) < 0 ? this : val);
3536 }
3537
3538 /**
3539 * Returns the maximum of this BigInteger and {@code val}.
3540 *
3541 * @param val value with which the maximum is to be computed.
3542 * @return the BigInteger whose value is the greater of this and
3543 * {@code val}. If they are equal, either may be returned.
3544 */
3545 public BigInteger max(BigInteger val) {
3546 return (compareTo(val) > 0 ? this : val);
3547 }
3548
3549
3550 // Hash Function
3551
3552 /**
3553 * Returns the hash code for this BigInteger.
3554 *
3555 * @return hash code for this BigInteger.
3556 */
3557 public int hashCode() {
3558 int hashCode = 0;
3559
3560 for (int i=0; i < mag.length; i++)
3561 hashCode = (int)(31*hashCode + (mag[i] & LONG_MASK));
3562
3563 return hashCode * signum;
3564 }
3565
3566 /**
3567 * Returns the String representation of this BigInteger in the
3568 * given radix. If the radix is outside the range from {@link
3569 * Character#MIN_RADIX} to {@link Character#MAX_RADIX} inclusive,
3570 * it will default to 10 (as is the case for
3571 * {@code Integer.toString}). The digit-to-character mapping
3572 * provided by {@code Character.forDigit} is used, and a minus
3573 * sign is prepended if appropriate. (This representation is
3574 * compatible with the {@link #BigInteger(String, int) (String,
3575 * int)} constructor.)
3576 *
3577 * @param radix radix of the String representation.
3578 * @return String representation of this BigInteger in the given radix.
3579 * @see Integer#toString
3580 * @see Character#forDigit
3581 * @see #BigInteger(java.lang.String, int)
3582 */
3583 public String toString(int radix) {
3584 if (signum == 0)
3585 return "0";
3586 if (radix < Character.MIN_RADIX || radix > Character.MAX_RADIX)
3587 radix = 10;
3588
3589 // If it's small enough, use smallToString.
3590 if (mag.length <= SCHOENHAGE_BASE_CONVERSION_THRESHOLD)
3591 return smallToString(radix);
3592
3593 // Otherwise use recursive toString, which requires positive arguments.
3594 // The results will be concatenated into this StringBuilder
3595 StringBuilder sb = new StringBuilder();
3596 if (signum < 0) {
3597 toString(this.negate(), sb, radix, 0);
3598 sb.insert(0, '-');
3599 }
3600 else
3601 toString(this, sb, radix, 0);
3602
3603 return sb.toString();
3604 }
3605
3606 /** This method is used to perform toString when arguments are small. */
3607 private String smallToString(int radix) {
3608 if (signum == 0) {
3609 return "0";
3610 }
3611
3612 // Compute upper bound on number of digit groups and allocate space
3613 int maxNumDigitGroups = (4*mag.length + 6)/7;
3614 String digitGroup[] = new String[maxNumDigitGroups];
3615
3616 // Translate number to string, a digit group at a time
3617 BigInteger tmp = this.abs();
3618 int numGroups = 0;
3619 while (tmp.signum != 0) {
3620 BigInteger d = longRadix[radix];
3621
3622 MutableBigInteger q = new MutableBigInteger(),
3623 a = new MutableBigInteger(tmp.mag),
3624 b = new MutableBigInteger(d.mag);
3625 MutableBigInteger r = a.divide(b, q);
3626 BigInteger q2 = q.toBigInteger(tmp.signum * d.signum);
3627 BigInteger r2 = r.toBigInteger(tmp.signum * d.signum);
3628
3629 digitGroup[numGroups++] = Long.toString(r2.longValue(), radix);
3630 tmp = q2;
3631 }
3632
3633 // Put sign (if any) and first digit group into result buffer
3634 StringBuilder buf = new StringBuilder(numGroups*digitsPerLong[radix]+1);
3635 if (signum < 0) {
3636 buf.append('-');
3637 }
3638 buf.append(digitGroup[numGroups-1]);
3639
3640 // Append remaining digit groups padded with leading zeros
3641 for (int i=numGroups-2; i >= 0; i--) {
3642 // Prepend (any) leading zeros for this digit group
3643 int numLeadingZeros = digitsPerLong[radix]-digitGroup[i].length();
3644 if (numLeadingZeros != 0) {
3645 buf.append(zeros[numLeadingZeros]);
3646 }
3647 buf.append(digitGroup[i]);
3648 }
3649 return buf.toString();
3650 }
3651
3652 /**
3653 * Converts the specified BigInteger to a string and appends to
3654 * {@code sb}. This implements the recursive Schoenhage algorithm
3655 * for base conversions.
3656 * <p>
3657 * See Knuth, Donald, _The Art of Computer Programming_, Vol. 2,
3658 * Answers to Exercises (4.4) Question 14.
3659 *
3660 * @param u The number to convert to a string.
3661 * @param sb The StringBuilder that will be appended to in place.
3662 * @param radix The base to convert to.
3663 * @param digits The minimum number of digits to pad to.
3664 */
3665 private static void toString(BigInteger u, StringBuilder sb, int radix,
3666 int digits) {
3667 // If we're smaller than a certain threshold, use the smallToString
3668 // method, padding with leading zeroes when necessary.
3669 if (u.mag.length <= SCHOENHAGE_BASE_CONVERSION_THRESHOLD) {
3670 String s = u.smallToString(radix);
3671
3672 // Pad with internal zeros if necessary.
3673 // Don't pad if we're at the beginning of the string.
3674 if ((s.length() < digits) && (sb.length() > 0)) {
3675 for (int i=s.length(); i < digits; i++) {
3676 sb.append('0');
3677 }
3678 }
3679
3680 sb.append(s);
3681 return;
3682 }
3683
3684 int b, n;
3685 b = u.bitLength();
3686
3687 // Calculate a value for n in the equation radix^(2^n) = u
3688 // and subtract 1 from that value. This is used to find the
3689 // cache index that contains the best value to divide u.
3690 n = (int) Math.round(Math.log(b * LOG_TWO / logCache[radix]) / LOG_TWO - 1.0);
3691 BigInteger v = getRadixConversionCache(radix, n);
3692 BigInteger[] results;
3693 results = u.divideAndRemainder(v);
3694
3695 int expectedDigits = 1 << n;
3696
3697 // Now recursively build the two halves of each number.
3698 toString(results[0], sb, radix, digits-expectedDigits);
3699 toString(results[1], sb, radix, expectedDigits);
3700 }
3701
3702 /**
3703 * Returns the value radix^(2^exponent) from the cache.
3704 * If this value doesn't already exist in the cache, it is added.
3705 * <p>
3706 * This could be changed to a more complicated caching method using
3707 * {@code Future}.
3708 */
3709 private static BigInteger getRadixConversionCache(int radix, int exponent) {
3710 BigInteger[] cacheLine = powerCache[radix]; // volatile read
3711 if (exponent < cacheLine.length) {
3712 return cacheLine[exponent];
3713 }
3714
3715 int oldLength = cacheLine.length;
3716 cacheLine = Arrays.copyOf(cacheLine, exponent + 1);
3717 for (int i = oldLength; i <= exponent; i++) {
3718 cacheLine[i] = cacheLine[i - 1].pow(2);
3719 }
3720
3721 BigInteger[][] pc = powerCache; // volatile read again
3722 if (exponent >= pc[radix].length) {
3723 pc = pc.clone();
3724 pc[radix] = cacheLine;
3725 powerCache = pc; // volatile write, publish
3726 }
3727 return cacheLine[exponent];
3728 }
3729
3730 /* zero[i] is a string of i consecutive zeros. */
3731 private static String zeros[] = new String[64];
3732 static {
3733 zeros[63] =
3734 "000000000000000000000000000000000000000000000000000000000000000";
3735 for (int i=0; i < 63; i++)
3736 zeros[i] = zeros[63].substring(0, i);
3737 }
3738
3739 /**
3740 * Returns the decimal String representation of this BigInteger.
3741 * The digit-to-character mapping provided by
3742 * {@code Character.forDigit} is used, and a minus sign is
3743 * prepended if appropriate. (This representation is compatible
3744 * with the {@link #BigInteger(String) (String)} constructor, and
3745 * allows for String concatenation with Java's + operator.)
3746 *
3747 * @return decimal String representation of this BigInteger.
3748 * @see Character#forDigit
3749 * @see #BigInteger(java.lang.String)
3750 */
3751 public String toString() {
3752 return toString(10);
3753 }
3754
3755 /**
3756 * Returns a byte array containing the two's-complement
3757 * representation of this BigInteger. The byte array will be in
3758 * <i>big-endian</i> byte-order: the most significant byte is in
3759 * the zeroth element. The array will contain the minimum number
3760 * of bytes required to represent this BigInteger, including at
3761 * least one sign bit, which is {@code (ceil((this.bitLength() +
3762 * 1)/8))}. (This representation is compatible with the
3763 * {@link #BigInteger(byte[]) (byte[])} constructor.)
3764 *
3765 * @return a byte array containing the two's-complement representation of
3766 * this BigInteger.
3767 * @see #BigInteger(byte[])
3768 */
3769 public byte[] toByteArray() {
3770 int byteLen = bitLength()/8 + 1;
3771 byte[] byteArray = new byte[byteLen];
3772
3773 for (int i=byteLen-1, bytesCopied=4, nextInt=0, intIndex=0; i >= 0; i--) {
3774 if (bytesCopied == 4) {
3775 nextInt = getInt(intIndex++);
3776 bytesCopied = 1;
3777 } else {
3778 nextInt >>>= 8;
3779 bytesCopied++;
3780 }
3781 byteArray[i] = (byte)nextInt;
3782 }
3783 return byteArray;
3784 }
3785
3786 /**
3787 * Converts this BigInteger to an {@code int}. This
3788 * conversion is analogous to a
3789 * <i>narrowing primitive conversion</i> from {@code long} to
3790 * {@code int} as defined in section 5.1.3 of
3791 * <cite>The Java™ Language Specification</cite>:
3792 * if this BigInteger is too big to fit in an
3793 * {@code int}, only the low-order 32 bits are returned.
3794 * Note that this conversion can lose information about the
3795 * overall magnitude of the BigInteger value as well as return a
3796 * result with the opposite sign.
3797 *
3798 * @return this BigInteger converted to an {@code int}.
3799 * @see #intValueExact()
3800 */
3801 public int intValue() {
3802 int result = 0;
3803 result = getInt(0);
3804 return result;
3805 }
3806
3807 /**
3808 * Converts this BigInteger to a {@code long}. This
3809 * conversion is analogous to a
3810 * <i>narrowing primitive conversion</i> from {@code long} to
3811 * {@code int} as defined in section 5.1.3 of
3812 * <cite>The Java™ Language Specification</cite>:
3813 * if this BigInteger is too big to fit in a
3814 * {@code long}, only the low-order 64 bits are returned.
3815 * Note that this conversion can lose information about the
3816 * overall magnitude of the BigInteger value as well as return a
3817 * result with the opposite sign.
3818 *
3819 * @return this BigInteger converted to a {@code long}.
3820 * @see #longValueExact()
3821 */
3822 public long longValue() {
3823 long result = 0;
3824
3825 for (int i=1; i >= 0; i--)
3826 result = (result << 32) + (getInt(i) & LONG_MASK);
3827 return result;
3828 }
3829
3830 /**
3831 * Converts this BigInteger to a {@code float}. This
3832 * conversion is similar to the
3833 * <i>narrowing primitive conversion</i> from {@code double} to
3834 * {@code float} as defined in section 5.1.3 of
3835 * <cite>The Java™ Language Specification</cite>:
3836 * if this BigInteger has too great a magnitude
3837 * to represent as a {@code float}, it will be converted to
3838 * {@link Float#NEGATIVE_INFINITY} or {@link
3839 * Float#POSITIVE_INFINITY} as appropriate. Note that even when
3840 * the return value is finite, this conversion can lose
3841 * information about the precision of the BigInteger value.
3842 *
3843 * @return this BigInteger converted to a {@code float}.
3844 */
3845 public float floatValue() {
3846 if (signum == 0) {
3847 return 0.0f;
3848 }
3849
3850 int exponent = ((mag.length - 1) << 5) + bitLengthForInt(mag[0]) - 1;
3851
3852 // exponent == floor(log2(abs(this)))
3853 if (exponent < Long.SIZE - 1) {
3854 return longValue();
3855 } else if (exponent > Float.MAX_EXPONENT) {
3856 return signum > 0 ? Float.POSITIVE_INFINITY : Float.NEGATIVE_INFINITY;
3857 }
3858
3859 /*
3860 * We need the top SIGNIFICAND_WIDTH bits, including the "implicit"
3861 * one bit. To make rounding easier, we pick out the top
3862 * SIGNIFICAND_WIDTH + 1 bits, so we have one to help us round up or
3863 * down. twiceSignifFloor will contain the top SIGNIFICAND_WIDTH + 1
3864 * bits, and signifFloor the top SIGNIFICAND_WIDTH.
3865 *
3866 * It helps to consider the real number signif = abs(this) *
3867 * 2^(SIGNIFICAND_WIDTH - 1 - exponent).
3868 */
3869 int shift = exponent - FloatConsts.SIGNIFICAND_WIDTH;
3870
3871 int twiceSignifFloor;
3872 // twiceSignifFloor will be == abs().shiftRight(shift).intValue()
3873 // We do the shift into an int directly to improve performance.
3874
3875 int nBits = shift & 0x1f;
3876 int nBits2 = 32 - nBits;
3877
3878 if (nBits == 0) {
3879 twiceSignifFloor = mag[0];
3880 } else {
3881 twiceSignifFloor = mag[0] >>> nBits;
3882 if (twiceSignifFloor == 0) {
3883 twiceSignifFloor = (mag[0] << nBits2) | (mag[1] >>> nBits);
3884 }
3885 }
3886
3887 int signifFloor = twiceSignifFloor >> 1;
3888 signifFloor &= FloatConsts.SIGNIF_BIT_MASK; // remove the implied bit
3889
3890 /*
3891 * We round up if either the fractional part of signif is strictly
3892 * greater than 0.5 (which is true if the 0.5 bit is set and any lower
3893 * bit is set), or if the fractional part of signif is >= 0.5 and
3894 * signifFloor is odd (which is true if both the 0.5 bit and the 1 bit
3895 * are set). This is equivalent to the desired HALF_EVEN rounding.
3896 */
3897 boolean increment = (twiceSignifFloor & 1) != 0
3898 && ((signifFloor & 1) != 0 || abs().getLowestSetBit() < shift);
3899 int signifRounded = increment ? signifFloor + 1 : signifFloor;
3900 int bits = ((exponent + FloatConsts.EXP_BIAS))
3901 << (FloatConsts.SIGNIFICAND_WIDTH - 1);
3902 bits += signifRounded;
3903 /*
3904 * If signifRounded == 2^24, we'd need to set all of the significand
3905 * bits to zero and add 1 to the exponent. This is exactly the behavior
3906 * we get from just adding signifRounded to bits directly. If the
3907 * exponent is Float.MAX_EXPONENT, we round up (correctly) to
3908 * Float.POSITIVE_INFINITY.
3909 */
3910 bits |= signum & FloatConsts.SIGN_BIT_MASK;
3911 return Float.intBitsToFloat(bits);
3912 }
3913
3914 /**
3915 * Converts this BigInteger to a {@code double}. This
3916 * conversion is similar to the
3917 * <i>narrowing primitive conversion</i> from {@code double} to
3918 * {@code float} as defined in section 5.1.3 of
3919 * <cite>The Java™ Language Specification</cite>:
3920 * if this BigInteger has too great a magnitude
3921 * to represent as a {@code double}, it will be converted to
3922 * {@link Double#NEGATIVE_INFINITY} or {@link
3923 * Double#POSITIVE_INFINITY} as appropriate. Note that even when
3924 * the return value is finite, this conversion can lose
3925 * information about the precision of the BigInteger value.
3926 *
3927 * @return this BigInteger converted to a {@code double}.
3928 */
3929 public double doubleValue() {
3930 if (signum == 0) {
3931 return 0.0;
3932 }
3933
3934 int exponent = ((mag.length - 1) << 5) + bitLengthForInt(mag[0]) - 1;
3935
3936 // exponent == floor(log2(abs(this))Double)
3937 if (exponent < Long.SIZE - 1) {
3938 return longValue();
3939 } else if (exponent > Double.MAX_EXPONENT) {
3940 return signum > 0 ? Double.POSITIVE_INFINITY : Double.NEGATIVE_INFINITY;
3941 }
3942
3943 /*
3944 * We need the top SIGNIFICAND_WIDTH bits, including the "implicit"
3945 * one bit. To make rounding easier, we pick out the top
3946 * SIGNIFICAND_WIDTH + 1 bits, so we have one to help us round up or
3947 * down. twiceSignifFloor will contain the top SIGNIFICAND_WIDTH + 1
3948 * bits, and signifFloor the top SIGNIFICAND_WIDTH.
3949 *
3950 * It helps to consider the real number signif = abs(this) *
3951 * 2^(SIGNIFICAND_WIDTH - 1 - exponent).
3952 */
3953 int shift = exponent - DoubleConsts.SIGNIFICAND_WIDTH;
3954
3955 long twiceSignifFloor;
3956 // twiceSignifFloor will be == abs().shiftRight(shift).longValue()
3957 // We do the shift into a long directly to improve performance.
3958
3959 int nBits = shift & 0x1f;
3960 int nBits2 = 32 - nBits;
3961
3962 int highBits;
3963 int lowBits;
3964 if (nBits == 0) {
3965 highBits = mag[0];
3966 lowBits = mag[1];
3967 } else {
3968 highBits = mag[0] >>> nBits;
3969 lowBits = (mag[0] << nBits2) | (mag[1] >>> nBits);
3970 if (highBits == 0) {
3971 highBits = lowBits;
3972 lowBits = (mag[1] << nBits2) | (mag[2] >>> nBits);
3973 }
3974 }
3975
3976 twiceSignifFloor = ((highBits & LONG_MASK) << 32)
3977 | (lowBits & LONG_MASK);
3978
3979 long signifFloor = twiceSignifFloor >> 1;
3980 signifFloor &= DoubleConsts.SIGNIF_BIT_MASK; // remove the implied bit
3981
3982 /*
3983 * We round up if either the fractional part of signif is strictly
3984 * greater than 0.5 (which is true if the 0.5 bit is set and any lower
3985 * bit is set), or if the fractional part of signif is >= 0.5 and
3986 * signifFloor is odd (which is true if both the 0.5 bit and the 1 bit
3987 * are set). This is equivalent to the desired HALF_EVEN rounding.
3988 */
3989 boolean increment = (twiceSignifFloor & 1) != 0
3990 && ((signifFloor & 1) != 0 || abs().getLowestSetBit() < shift);
3991 long signifRounded = increment ? signifFloor + 1 : signifFloor;
3992 long bits = (long) ((exponent + DoubleConsts.EXP_BIAS))
3993 << (DoubleConsts.SIGNIFICAND_WIDTH - 1);
3994 bits += signifRounded;
3995 /*
3996 * If signifRounded == 2^53, we'd need to set all of the significand
3997 * bits to zero and add 1 to the exponent. This is exactly the behavior
3998 * we get from just adding signifRounded to bits directly. If the
3999 * exponent is Double.MAX_EXPONENT, we round up (correctly) to
4000 * Double.POSITIVE_INFINITY.
4001 */
4002 bits |= signum & DoubleConsts.SIGN_BIT_MASK;
4003 return Double.longBitsToDouble(bits);
4004 }
4005
4006 /**
4007 * Returns a copy of the input array stripped of any leading zero bytes.
4008 */
4009 private static int[] stripLeadingZeroInts(int val[]) {
4010 int vlen = val.length;
4011 int keep;
4012
4013 // Find first nonzero byte
4014 for (keep = 0; keep < vlen && val[keep] == 0; keep++)
4015 ;
4016 return java.util.Arrays.copyOfRange(val, keep, vlen);
4017 }
4018
4019 /**
4020 * Returns the input array stripped of any leading zero bytes.
4021 * Since the source is trusted the copying may be skipped.
4022 */
4023 private static int[] trustedStripLeadingZeroInts(int val[]) {
4024 int vlen = val.length;
4025 int keep;
4026
4027 // Find first nonzero byte
4028 for (keep = 0; keep < vlen && val[keep] == 0; keep++)
4029 ;
4030 return keep == 0 ? val : java.util.Arrays.copyOfRange(val, keep, vlen);
4031 }
4032
4033 /**
4034 * Returns a copy of the input array stripped of any leading zero bytes.
4035 */
4036 private static int[] stripLeadingZeroBytes(byte a[], int off, int len) {
4037 int indexBound = off + len;
4038 int keep;
4039
4040 // Find first nonzero byte
4041 for (keep = off; keep < indexBound && a[keep] == 0; keep++)
4042 ;
4043
4044 // Allocate new array and copy relevant part of input array
4045 int intLength = ((indexBound - keep) + 3) >>> 2;
4046 int[] result = new int[intLength];
4047 int b = indexBound - 1;
4048 for (int i = intLength-1; i >= 0; i--) {
4049 result[i] = a[b--] & 0xff;
4050 int bytesRemaining = b - keep + 1;
4051 int bytesToTransfer = Math.min(3, bytesRemaining);
4052 for (int j=8; j <= (bytesToTransfer << 3); j += 8)
4053 result[i] |= ((a[b--] & 0xff) << j);
4054 }
4055 return result;
4056 }
4057
4058 /**
4059 * Takes an array a representing a negative 2's-complement number and
4060 * returns the minimal (no leading zero bytes) unsigned whose value is -a.
4061 */
4062 private static int[] makePositive(byte a[], int off, int len) {
4063 int keep, k;
4064 int indexBound = off + len;
4065
4066 // Find first non-sign (0xff) byte of input
4067 for (keep=off; keep < indexBound && a[keep] == -1; keep++)
4068 ;
4069
4070
4071 /* Allocate output array. If all non-sign bytes are 0x00, we must
4072 * allocate space for one extra output byte. */
4073 for (k=keep; k < indexBound && a[k] == 0; k++)
4074 ;
4075
4076 int extraByte = (k == indexBound) ? 1 : 0;
4077 int intLength = ((indexBound - keep + extraByte) + 3) >>> 2;
4078 int result[] = new int[intLength];
4079
4080 /* Copy one's complement of input into output, leaving extra
4081 * byte (if it exists) == 0x00 */
4082 int b = indexBound - 1;
4083 for (int i = intLength-1; i >= 0; i--) {
4084 result[i] = a[b--] & 0xff;
4085 int numBytesToTransfer = Math.min(3, b-keep+1);
4086 if (numBytesToTransfer < 0)
4087 numBytesToTransfer = 0;
4088 for (int j=8; j <= 8*numBytesToTransfer; j += 8)
4089 result[i] |= ((a[b--] & 0xff) << j);
4090
4091 // Mask indicates which bits must be complemented
4092 int mask = -1 >>> (8*(3-numBytesToTransfer));
4093 result[i] = ~result[i] & mask;
4094 }
4095
4096 // Add one to one's complement to generate two's complement
4097 for (int i=result.length-1; i >= 0; i--) {
4098 result[i] = (int)((result[i] & LONG_MASK) + 1);
4099 if (result[i] != 0)
4100 break;
4101 }
4102
4103 return result;
4104 }
4105
4106 /**
4107 * Takes an array a representing a negative 2's-complement number and
4108 * returns the minimal (no leading zero ints) unsigned whose value is -a.
4109 */
4110 private static int[] makePositive(int a[]) {
4111 int keep, j;
4112
4113 // Find first non-sign (0xffffffff) int of input
4114 for (keep=0; keep < a.length && a[keep] == -1; keep++)
4115 ;
4116
4117 /* Allocate output array. If all non-sign ints are 0x00, we must
4118 * allocate space for one extra output int. */
4119 for (j=keep; j < a.length && a[j] == 0; j++)
4120 ;
4121 int extraInt = (j == a.length ? 1 : 0);
4122 int result[] = new int[a.length - keep + extraInt];
4123
4124 /* Copy one's complement of input into output, leaving extra
4125 * int (if it exists) == 0x00 */
4126 for (int i = keep; i < a.length; i++)
4127 result[i - keep + extraInt] = ~a[i];
4128
4129 // Add one to one's complement to generate two's complement
4130 for (int i=result.length-1; ++result[i] == 0; i--)
4131 ;
4132
4133 return result;
4134 }
4135
4136 /*
4137 * The following two arrays are used for fast String conversions. Both
4138 * are indexed by radix. The first is the number of digits of the given
4139 * radix that can fit in a Java long without "going negative", i.e., the
4140 * highest integer n such that radix**n < 2**63. The second is the
4141 * "long radix" that tears each number into "long digits", each of which
4142 * consists of the number of digits in the corresponding element in
4143 * digitsPerLong (longRadix[i] = i**digitPerLong[i]). Both arrays have
4144 * nonsense values in their 0 and 1 elements, as radixes 0 and 1 are not
4145 * used.
4146 */
4147 private static int digitsPerLong[] = {0, 0,
4148 62, 39, 31, 27, 24, 22, 20, 19, 18, 18, 17, 17, 16, 16, 15, 15, 15, 14,
4149 14, 14, 14, 13, 13, 13, 13, 13, 13, 12, 12, 12, 12, 12, 12, 12, 12};
4150
4151 private static BigInteger longRadix[] = {null, null,
4152 valueOf(0x4000000000000000L), valueOf(0x383d9170b85ff80bL),
4153 valueOf(0x4000000000000000L), valueOf(0x6765c793fa10079dL),
4154 valueOf(0x41c21cb8e1000000L), valueOf(0x3642798750226111L),
4155 valueOf(0x1000000000000000L), valueOf(0x12bf307ae81ffd59L),
4156 valueOf( 0xde0b6b3a7640000L), valueOf(0x4d28cb56c33fa539L),
4157 valueOf(0x1eca170c00000000L), valueOf(0x780c7372621bd74dL),
4158 valueOf(0x1e39a5057d810000L), valueOf(0x5b27ac993df97701L),
4159 valueOf(0x1000000000000000L), valueOf(0x27b95e997e21d9f1L),
4160 valueOf(0x5da0e1e53c5c8000L), valueOf( 0xb16a458ef403f19L),
4161 valueOf(0x16bcc41e90000000L), valueOf(0x2d04b7fdd9c0ef49L),
4162 valueOf(0x5658597bcaa24000L), valueOf( 0x6feb266931a75b7L),
4163 valueOf( 0xc29e98000000000L), valueOf(0x14adf4b7320334b9L),
4164 valueOf(0x226ed36478bfa000L), valueOf(0x383d9170b85ff80bL),
4165 valueOf(0x5a3c23e39c000000L), valueOf( 0x4e900abb53e6b71L),
4166 valueOf( 0x7600ec618141000L), valueOf( 0xaee5720ee830681L),
4167 valueOf(0x1000000000000000L), valueOf(0x172588ad4f5f0981L),
4168 valueOf(0x211e44f7d02c1000L), valueOf(0x2ee56725f06e5c71L),
4169 valueOf(0x41c21cb8e1000000L)};
4170
4171 /*
4172 * These two arrays are the integer analogue of above.
4173 */
4174 private static int digitsPerInt[] = {0, 0, 30, 19, 15, 13, 11,
4175 11, 10, 9, 9, 8, 8, 8, 8, 7, 7, 7, 7, 7, 7, 7, 6, 6, 6, 6,
4176 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 5};
4177
4178 private static int intRadix[] = {0, 0,
4179 0x40000000, 0x4546b3db, 0x40000000, 0x48c27395, 0x159fd800,
4180 0x75db9c97, 0x40000000, 0x17179149, 0x3b9aca00, 0xcc6db61,
4181 0x19a10000, 0x309f1021, 0x57f6c100, 0xa2f1b6f, 0x10000000,
4182 0x18754571, 0x247dbc80, 0x3547667b, 0x4c4b4000, 0x6b5a6e1d,
4183 0x6c20a40, 0x8d2d931, 0xb640000, 0xe8d4a51, 0x1269ae40,
4184 0x17179149, 0x1cb91000, 0x23744899, 0x2b73a840, 0x34e63b41,
4185 0x40000000, 0x4cfa3cc1, 0x5c13d840, 0x6d91b519, 0x39aa400
4186 };
4187
4188 /**
4189 * These routines provide access to the two's complement representation
4190 * of BigIntegers.
4191 */
4192
4193 /**
4194 * Returns the length of the two's complement representation in ints,
4195 * including space for at least one sign bit.
4196 */
4197 private int intLength() {
4198 return (bitLength() >>> 5) + 1;
4199 }
4200
4201 /* Returns sign bit */
4202 private int signBit() {
4203 return signum < 0 ? 1 : 0;
4204 }
4205
4206 /* Returns an int of sign bits */
4207 private int signInt() {
4208 return signum < 0 ? -1 : 0;
4209 }
4210
4211 /**
4212 * Returns the specified int of the little-endian two's complement
4213 * representation (int 0 is the least significant). The int number can
4214 * be arbitrarily high (values are logically preceded by infinitely many
4215 * sign ints).
4216 */
4217 private int getInt(int n) {
4218 if (n < 0)
4219 return 0;
4220 if (n >= mag.length)
4221 return signInt();
4222
4223 int magInt = mag[mag.length-n-1];
4224
4225 return (signum >= 0 ? magInt :
4226 (n <= firstNonzeroIntNum() ? -magInt : ~magInt));
4227 }
4228
4229 /**
4230 * Returns the index of the int that contains the first nonzero int in the
4231 * little-endian binary representation of the magnitude (int 0 is the
4232 * least significant). If the magnitude is zero, return value is undefined.
4233 *
4234 * <p>Note: never used for a BigInteger with a magnitude of zero.
4235 * @see #getInt.
4236 */
4237 private int firstNonzeroIntNum() {
4238 int fn = firstNonzeroIntNumPlusTwo - 2;
4239 if (fn == -2) { // firstNonzeroIntNum not initialized yet
4240 // Search for the first nonzero int
4241 int i;
4242 int mlen = mag.length;
4243 for (i = mlen - 1; i >= 0 && mag[i] == 0; i--)
4244 ;
4245 fn = mlen - i - 1;
4246 firstNonzeroIntNumPlusTwo = fn + 2; // offset by two to initialize
4247 }
4248 return fn;
4249 }
4250
4251 /** use serialVersionUID from JDK 1.1. for interoperability */
4252 private static final long serialVersionUID = -8287574255936472291L;
4253
4254 /**
4255 * Serializable fields for BigInteger.
4256 *
4257 * @serialField signum int
4258 * signum of this BigInteger
4259 * @serialField magnitude byte[]
4260 * magnitude array of this BigInteger
4261 * @serialField bitCount int
4262 * appears in the serialized form for backward compatibility
4263 * @serialField bitLength int
4264 * appears in the serialized form for backward compatibility
4265 * @serialField firstNonzeroByteNum int
4266 * appears in the serialized form for backward compatibility
4267 * @serialField lowestSetBit int
4268 * appears in the serialized form for backward compatibility
4269 */
4270 private static final ObjectStreamField[] serialPersistentFields = {
4271 new ObjectStreamField("signum", Integer.TYPE),
4272 new ObjectStreamField("magnitude", byte[].class),
4273 new ObjectStreamField("bitCount", Integer.TYPE),
4274 new ObjectStreamField("bitLength", Integer.TYPE),
4275 new ObjectStreamField("firstNonzeroByteNum", Integer.TYPE),
4276 new ObjectStreamField("lowestSetBit", Integer.TYPE)
4277 };
4278
4279 /**
4280 * Reconstitute the {@code BigInteger} instance from a stream (that is,
4281 * deserialize it). The magnitude is read in as an array of bytes
4282 * for historical reasons, but it is converted to an array of ints
4283 * and the byte array is discarded.
4284 * Note:
4285 * The current convention is to initialize the cache fields, bitCountPlusOne,
4286 * bitLengthPlusOne and lowestSetBitPlusTwo, to 0 rather than some other
4287 * marker value. Therefore, no explicit action to set these fields needs to
4288 * be taken in readObject because those fields already have a 0 value by
4289 * default since defaultReadObject is not being used.
4290 */
4291 private void readObject(java.io.ObjectInputStream s)
4292 throws java.io.IOException, ClassNotFoundException {
4293 // prepare to read the alternate persistent fields
4294 ObjectInputStream.GetField fields = s.readFields();
4295
4296 // Read the alternate persistent fields that we care about
4297 int sign = fields.get("signum", -2);
4298 byte[] magnitude = (byte[])fields.get("magnitude", null);
4299
4300 // Validate signum
4301 if (sign < -1 || sign > 1) {
4302 String message = "BigInteger: Invalid signum value";
4303 if (fields.defaulted("signum"))
4304 message = "BigInteger: Signum not present in stream";
4305 throw new java.io.StreamCorruptedException(message);
4306 }
4307 int[] mag = stripLeadingZeroBytes(magnitude, 0, magnitude.length);
4308 if ((mag.length == 0) != (sign == 0)) {
4309 String message = "BigInteger: signum-magnitude mismatch";
4310 if (fields.defaulted("magnitude"))
4311 message = "BigInteger: Magnitude not present in stream";
4312 throw new java.io.StreamCorruptedException(message);
4313 }
4314
4315 // Commit final fields via Unsafe
4316 UnsafeHolder.putSign(this, sign);
4317
4318 // Calculate mag field from magnitude and discard magnitude
4319 UnsafeHolder.putMag(this, mag);
4320 if (mag.length >= MAX_MAG_LENGTH) {
4321 try {
4322 checkRange();
4323 } catch (ArithmeticException e) {
4324 throw new java.io.StreamCorruptedException("BigInteger: Out of the supported range");
4325 }
4326 }
4327 }
4328
4329 // Support for resetting final fields while deserializing
4330 private static class UnsafeHolder {
4331 private static final sun.misc.Unsafe unsafe;
4332 private static final long signumOffset;
4333 private static final long magOffset;
4334 static {
4335 try {
4336 unsafe = sun.misc.Unsafe.getUnsafe();
4337 signumOffset = unsafe.objectFieldOffset
4338 (BigInteger.class.getDeclaredField("signum"));
4339 magOffset = unsafe.objectFieldOffset
4340 (BigInteger.class.getDeclaredField("mag"));
4341 } catch (Exception ex) {
4342 throw new ExceptionInInitializerError(ex);
4343 }
4344 }
4345
4346 static void putSign(BigInteger bi, int sign) {
4347 unsafe.putIntVolatile(bi, signumOffset, sign);
4348 }
4349
4350 static void putMag(BigInteger bi, int[] magnitude) {
4351 unsafe.putObjectVolatile(bi, magOffset, magnitude);
4352 }
4353 }
4354
4355 /**
4356 * Save the {@code BigInteger} instance to a stream. The magnitude of a
4357 * {@code BigInteger} is serialized as a byte array for historical reasons.
4358 * To maintain compatibility with older implementations, the integers
4359 * -1, -1, -2, and -2 are written as the values of the obsolete fields
4360 * {@code bitCount}, {@code bitLength}, {@code lowestSetBit}, and
4361 * {@code firstNonzeroByteNum}, respectively. These values are compatible
4362 * with older implementations, but will be ignored by current
4363 * implementations.
4364 */
4365 private void writeObject(ObjectOutputStream s) throws IOException {
4366 // set the values of the Serializable fields
4367 ObjectOutputStream.PutField fields = s.putFields();
4368 fields.put("signum", signum);
4369 fields.put("magnitude", magSerializedForm());
4370 // The values written for cached fields are compatible with older
4371 // versions, but are ignored in readObject so don't otherwise matter.
4372 fields.put("bitCount", -1);
4373 fields.put("bitLength", -1);
4374 fields.put("lowestSetBit", -2);
4375 fields.put("firstNonzeroByteNum", -2);
4376
4377 // save them
4378 s.writeFields();
4379 }
4380
4381 /**
4382 * Returns the mag array as an array of bytes.
4383 */
4384 private byte[] magSerializedForm() {
4385 int len = mag.length;
4386
4387 int bitLen = (len == 0 ? 0 : ((len - 1) << 5) + bitLengthForInt(mag[0]));
4388 int byteLen = (bitLen + 7) >>> 3;
4389 byte[] result = new byte[byteLen];
4390
4391 for (int i = byteLen - 1, bytesCopied = 4, intIndex = len - 1, nextInt = 0;
4392 i >= 0; i--) {
4393 if (bytesCopied == 4) {
4394 nextInt = mag[intIndex--];
4395 bytesCopied = 1;
4396 } else {
4397 nextInt >>>= 8;
4398 bytesCopied++;
4399 }
4400 result[i] = (byte)nextInt;
4401 }
4402 return result;
4403 }
4404
4405 /**
4406 * Converts this {@code BigInteger} to a {@code long}, checking
4407 * for lost information. If the value of this {@code BigInteger}
4408 * is out of the range of the {@code long} type, then an
4409 * {@code ArithmeticException} is thrown.
4410 *
4411 * @return this {@code BigInteger} converted to a {@code long}.
4412 * @throws ArithmeticException if the value of {@code this} will
4413 * not exactly fit in a {@code long}.
4414 * @see BigInteger#longValue
4415 * @since 1.8
4416 */
4417 public long longValueExact() {
4418 if (mag.length <= 2 && bitLength() <= 63)
4419 return longValue();
4420 else
4421 throw new ArithmeticException("BigInteger out of long range");
4422 }
4423
4424 /**
4425 * Converts this {@code BigInteger} to an {@code int}, checking
4426 * for lost information. If the value of this {@code BigInteger}
4427 * is out of the range of the {@code int} type, then an
4428 * {@code ArithmeticException} is thrown.
4429 *
4430 * @return this {@code BigInteger} converted to an {@code int}.
4431 * @throws ArithmeticException if the value of {@code this} will
4432 * not exactly fit in a {@code int}.
4433 * @see BigInteger#intValue
4434 * @since 1.8
4435 */
4436 public int intValueExact() {
4437 if (mag.length <= 1 && bitLength() <= 31)
4438 return intValue();
4439 else
4440 throw new ArithmeticException("BigInteger out of int range");
4441 }
4442
4443 /**
4444 * Converts this {@code BigInteger} to a {@code short}, checking
4445 * for lost information. If the value of this {@code BigInteger}
4446 * is out of the range of the {@code short} type, then an
4447 * {@code ArithmeticException} is thrown.
4448 *
4449 * @return this {@code BigInteger} converted to a {@code short}.
4450 * @throws ArithmeticException if the value of {@code this} will
4451 * not exactly fit in a {@code short}.
4452 * @see BigInteger#shortValue
4453 * @since 1.8
4454 */
4455 public short shortValueExact() {
4456 if (mag.length <= 1 && bitLength() <= 31) {
4457 int value = intValue();
4458 if (value >= Short.MIN_VALUE && value <= Short.MAX_VALUE)
4459 return shortValue();
4460 }
4461 throw new ArithmeticException("BigInteger out of short range");
4462 }
4463
4464 /**
4465 * Converts this {@code BigInteger} to a {@code byte}, checking
4466 * for lost information. If the value of this {@code BigInteger}
4467 * is out of the range of the {@code byte} type, then an
4468 * {@code ArithmeticException} is thrown.
4469 *
4470 * @return this {@code BigInteger} converted to a {@code byte}.
4471 * @throws ArithmeticException if the value of {@code this} will
4472 * not exactly fit in a {@code byte}.
4473 * @see BigInteger#byteValue
4474 * @since 1.8
4475 */
4476 public byte byteValueExact() {
4477 if (mag.length <= 1 && bitLength() <= 31) {
4478 int value = intValue();
4479 if (value >= Byte.MIN_VALUE && value <= Byte.MAX_VALUE)
4480 return byteValue();
4481 }
4482 throw new ArithmeticException("BigInteger out of byte range");
4483 }
4484 }