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