1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/emul/compact/src/main/java/java/math/BigInteger.java Sat Sep 07 13:51:24 2013 +0200
1.3 @@ -0,0 +1,3186 @@
1.4 +/*
1.5 + * Copyright (c) 1996, 2007, Oracle and/or its affiliates. All rights reserved.
1.6 + * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
1.7 + *
1.8 + * This code is free software; you can redistribute it and/or modify it
1.9 + * under the terms of the GNU General Public License version 2 only, as
1.10 + * published by the Free Software Foundation. Oracle designates this
1.11 + * particular file as subject to the "Classpath" exception as provided
1.12 + * by Oracle in the LICENSE file that accompanied this code.
1.13 + *
1.14 + * This code is distributed in the hope that it will be useful, but WITHOUT
1.15 + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
1.16 + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
1.17 + * version 2 for more details (a copy is included in the LICENSE file that
1.18 + * accompanied this code).
1.19 + *
1.20 + * You should have received a copy of the GNU General Public License version
1.21 + * 2 along with this work; if not, write to the Free Software Foundation,
1.22 + * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
1.23 + *
1.24 + * Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA
1.25 + * or visit www.oracle.com if you need additional information or have any
1.26 + * questions.
1.27 + */
1.28 +
1.29 +/*
1.30 + * Portions Copyright (c) 1995 Colin Plumb. All rights reserved.
1.31 + */
1.32 +
1.33 +package java.math;
1.34 +
1.35 +import java.util.Random;
1.36 +import java.io.*;
1.37 +
1.38 +/**
1.39 + * Immutable arbitrary-precision integers. All operations behave as if
1.40 + * BigIntegers were represented in two's-complement notation (like Java's
1.41 + * primitive integer types). BigInteger provides analogues to all of Java's
1.42 + * primitive integer operators, and all relevant methods from java.lang.Math.
1.43 + * Additionally, BigInteger provides operations for modular arithmetic, GCD
1.44 + * calculation, primality testing, prime generation, bit manipulation,
1.45 + * and a few other miscellaneous operations.
1.46 + *
1.47 + * <p>Semantics of arithmetic operations exactly mimic those of Java's integer
1.48 + * arithmetic operators, as defined in <i>The Java Language Specification</i>.
1.49 + * For example, division by zero throws an {@code ArithmeticException}, and
1.50 + * division of a negative by a positive yields a negative (or zero) remainder.
1.51 + * All of the details in the Spec concerning overflow are ignored, as
1.52 + * BigIntegers are made as large as necessary to accommodate the results of an
1.53 + * operation.
1.54 + *
1.55 + * <p>Semantics of shift operations extend those of Java's shift operators
1.56 + * to allow for negative shift distances. A right-shift with a negative
1.57 + * shift distance results in a left shift, and vice-versa. The unsigned
1.58 + * right shift operator ({@code >>>}) is omitted, as this operation makes
1.59 + * little sense in combination with the "infinite word size" abstraction
1.60 + * provided by this class.
1.61 + *
1.62 + * <p>Semantics of bitwise logical operations exactly mimic those of Java's
1.63 + * bitwise integer operators. The binary operators ({@code and},
1.64 + * {@code or}, {@code xor}) implicitly perform sign extension on the shorter
1.65 + * of the two operands prior to performing the operation.
1.66 + *
1.67 + * <p>Comparison operations perform signed integer comparisons, analogous to
1.68 + * those performed by Java's relational and equality operators.
1.69 + *
1.70 + * <p>Modular arithmetic operations are provided to compute residues, perform
1.71 + * exponentiation, and compute multiplicative inverses. These methods always
1.72 + * return a non-negative result, between {@code 0} and {@code (modulus - 1)},
1.73 + * inclusive.
1.74 + *
1.75 + * <p>Bit operations operate on a single bit of the two's-complement
1.76 + * representation of their operand. If necessary, the operand is sign-
1.77 + * extended so that it contains the designated bit. None of the single-bit
1.78 + * operations can produce a BigInteger with a different sign from the
1.79 + * BigInteger being operated on, as they affect only a single bit, and the
1.80 + * "infinite word size" abstraction provided by this class ensures that there
1.81 + * are infinitely many "virtual sign bits" preceding each BigInteger.
1.82 + *
1.83 + * <p>For the sake of brevity and clarity, pseudo-code is used throughout the
1.84 + * descriptions of BigInteger methods. The pseudo-code expression
1.85 + * {@code (i + j)} is shorthand for "a BigInteger whose value is
1.86 + * that of the BigInteger {@code i} plus that of the BigInteger {@code j}."
1.87 + * The pseudo-code expression {@code (i == j)} is shorthand for
1.88 + * "{@code true} if and only if the BigInteger {@code i} represents the same
1.89 + * value as the BigInteger {@code j}." Other pseudo-code expressions are
1.90 + * interpreted similarly.
1.91 + *
1.92 + * <p>All methods and constructors in this class throw
1.93 + * {@code NullPointerException} when passed
1.94 + * a null object reference for any input parameter.
1.95 + *
1.96 + * @see BigDecimal
1.97 + * @author Josh Bloch
1.98 + * @author Michael McCloskey
1.99 + * @since JDK1.1
1.100 + */
1.101 +
1.102 +public class BigInteger extends Number implements Comparable<BigInteger> {
1.103 + /**
1.104 + * The signum of this BigInteger: -1 for negative, 0 for zero, or
1.105 + * 1 for positive. Note that the BigInteger zero <i>must</i> have
1.106 + * a signum of 0. This is necessary to ensures that there is exactly one
1.107 + * representation for each BigInteger value.
1.108 + *
1.109 + * @serial
1.110 + */
1.111 + final int signum;
1.112 +
1.113 + /**
1.114 + * The magnitude of this BigInteger, in <i>big-endian</i> order: the
1.115 + * zeroth element of this array is the most-significant int of the
1.116 + * magnitude. The magnitude must be "minimal" in that the most-significant
1.117 + * int ({@code mag[0]}) must be non-zero. This is necessary to
1.118 + * ensure that there is exactly one representation for each BigInteger
1.119 + * value. Note that this implies that the BigInteger zero has a
1.120 + * zero-length mag array.
1.121 + */
1.122 + final int[] mag;
1.123 +
1.124 + // These "redundant fields" are initialized with recognizable nonsense
1.125 + // values, and cached the first time they are needed (or never, if they
1.126 + // aren't needed).
1.127 +
1.128 + /**
1.129 + * One plus the bitCount of this BigInteger. Zeros means unitialized.
1.130 + *
1.131 + * @serial
1.132 + * @see #bitCount
1.133 + * @deprecated Deprecated since logical value is offset from stored
1.134 + * value and correction factor is applied in accessor method.
1.135 + */
1.136 + @Deprecated
1.137 + private int bitCount;
1.138 +
1.139 + /**
1.140 + * One plus the bitLength of this BigInteger. Zeros means unitialized.
1.141 + * (either value is acceptable).
1.142 + *
1.143 + * @serial
1.144 + * @see #bitLength()
1.145 + * @deprecated Deprecated since logical value is offset from stored
1.146 + * value and correction factor is applied in accessor method.
1.147 + */
1.148 + @Deprecated
1.149 + private int bitLength;
1.150 +
1.151 + /**
1.152 + * Two plus the lowest set bit of this BigInteger, as returned by
1.153 + * getLowestSetBit().
1.154 + *
1.155 + * @serial
1.156 + * @see #getLowestSetBit
1.157 + * @deprecated Deprecated since logical value is offset from stored
1.158 + * value and correction factor is applied in accessor method.
1.159 + */
1.160 + @Deprecated
1.161 + private int lowestSetBit;
1.162 +
1.163 + /**
1.164 + * Two plus the index of the lowest-order int in the magnitude of this
1.165 + * BigInteger that contains a nonzero int, or -2 (either value is acceptable).
1.166 + * The least significant int has int-number 0, the next int in order of
1.167 + * increasing significance has int-number 1, and so forth.
1.168 + * @deprecated Deprecated since logical value is offset from stored
1.169 + * value and correction factor is applied in accessor method.
1.170 + */
1.171 + @Deprecated
1.172 + private int firstNonzeroIntNum;
1.173 +
1.174 + /**
1.175 + * This mask is used to obtain the value of an int as if it were unsigned.
1.176 + */
1.177 + final static long LONG_MASK = 0xffffffffL;
1.178 +
1.179 + //Constructors
1.180 +
1.181 + /**
1.182 + * Translates a byte array containing the two's-complement binary
1.183 + * representation of a BigInteger into a BigInteger. The input array is
1.184 + * assumed to be in <i>big-endian</i> byte-order: the most significant
1.185 + * byte is in the zeroth element.
1.186 + *
1.187 + * @param val big-endian two's-complement binary representation of
1.188 + * BigInteger.
1.189 + * @throws NumberFormatException {@code val} is zero bytes long.
1.190 + */
1.191 + public BigInteger(byte[] val) {
1.192 + if (val.length == 0)
1.193 + throw new NumberFormatException("Zero length BigInteger");
1.194 +
1.195 + if (val[0] < 0) {
1.196 + mag = makePositive(val);
1.197 + signum = -1;
1.198 + } else {
1.199 + mag = stripLeadingZeroBytes(val);
1.200 + signum = (mag.length == 0 ? 0 : 1);
1.201 + }
1.202 + }
1.203 +
1.204 + /**
1.205 + * This private constructor translates an int array containing the
1.206 + * two's-complement binary representation of a BigInteger into a
1.207 + * BigInteger. The input array is assumed to be in <i>big-endian</i>
1.208 + * int-order: the most significant int is in the zeroth element.
1.209 + */
1.210 + private BigInteger(int[] val) {
1.211 + if (val.length == 0)
1.212 + throw new NumberFormatException("Zero length BigInteger");
1.213 +
1.214 + if (val[0] < 0) {
1.215 + mag = makePositive(val);
1.216 + signum = -1;
1.217 + } else {
1.218 + mag = trustedStripLeadingZeroInts(val);
1.219 + signum = (mag.length == 0 ? 0 : 1);
1.220 + }
1.221 + }
1.222 +
1.223 + /**
1.224 + * Translates the sign-magnitude representation of a BigInteger into a
1.225 + * BigInteger. The sign is represented as an integer signum value: -1 for
1.226 + * negative, 0 for zero, or 1 for positive. The magnitude is a byte array
1.227 + * in <i>big-endian</i> byte-order: the most significant byte is in the
1.228 + * zeroth element. A zero-length magnitude array is permissible, and will
1.229 + * result in a BigInteger value of 0, whether signum is -1, 0 or 1.
1.230 + *
1.231 + * @param signum signum of the number (-1 for negative, 0 for zero, 1
1.232 + * for positive).
1.233 + * @param magnitude big-endian binary representation of the magnitude of
1.234 + * the number.
1.235 + * @throws NumberFormatException {@code signum} is not one of the three
1.236 + * legal values (-1, 0, and 1), or {@code signum} is 0 and
1.237 + * {@code magnitude} contains one or more non-zero bytes.
1.238 + */
1.239 + public BigInteger(int signum, byte[] magnitude) {
1.240 + this.mag = stripLeadingZeroBytes(magnitude);
1.241 +
1.242 + if (signum < -1 || signum > 1)
1.243 + throw(new NumberFormatException("Invalid signum value"));
1.244 +
1.245 + if (this.mag.length==0) {
1.246 + this.signum = 0;
1.247 + } else {
1.248 + if (signum == 0)
1.249 + throw(new NumberFormatException("signum-magnitude mismatch"));
1.250 + this.signum = signum;
1.251 + }
1.252 + }
1.253 +
1.254 + /**
1.255 + * A constructor for internal use that translates the sign-magnitude
1.256 + * representation of a BigInteger into a BigInteger. It checks the
1.257 + * arguments and copies the magnitude so this constructor would be
1.258 + * safe for external use.
1.259 + */
1.260 + private BigInteger(int signum, int[] magnitude) {
1.261 + this.mag = stripLeadingZeroInts(magnitude);
1.262 +
1.263 + if (signum < -1 || signum > 1)
1.264 + throw(new NumberFormatException("Invalid signum value"));
1.265 +
1.266 + if (this.mag.length==0) {
1.267 + this.signum = 0;
1.268 + } else {
1.269 + if (signum == 0)
1.270 + throw(new NumberFormatException("signum-magnitude mismatch"));
1.271 + this.signum = signum;
1.272 + }
1.273 + }
1.274 +
1.275 + /**
1.276 + * Translates the String representation of a BigInteger in the
1.277 + * specified radix into a BigInteger. The String representation
1.278 + * consists of an optional minus or plus sign followed by a
1.279 + * sequence of one or more digits in the specified radix. The
1.280 + * character-to-digit mapping is provided by {@code
1.281 + * Character.digit}. The String may not contain any extraneous
1.282 + * characters (whitespace, for example).
1.283 + *
1.284 + * @param val String representation of BigInteger.
1.285 + * @param radix radix to be used in interpreting {@code val}.
1.286 + * @throws NumberFormatException {@code val} is not a valid representation
1.287 + * of a BigInteger in the specified radix, or {@code radix} is
1.288 + * outside the range from {@link Character#MIN_RADIX} to
1.289 + * {@link Character#MAX_RADIX}, inclusive.
1.290 + * @see Character#digit
1.291 + */
1.292 + public BigInteger(String val, int radix) {
1.293 + int cursor = 0, numDigits;
1.294 + final int len = val.length();
1.295 +
1.296 + if (radix < Character.MIN_RADIX || radix > Character.MAX_RADIX)
1.297 + throw new NumberFormatException("Radix out of range");
1.298 + if (len == 0)
1.299 + throw new NumberFormatException("Zero length BigInteger");
1.300 +
1.301 + // Check for at most one leading sign
1.302 + int sign = 1;
1.303 + int index1 = val.lastIndexOf('-');
1.304 + int index2 = val.lastIndexOf('+');
1.305 + if ((index1 + index2) <= -1) {
1.306 + // No leading sign character or at most one leading sign character
1.307 + if (index1 == 0 || index2 == 0) {
1.308 + cursor = 1;
1.309 + if (len == 1)
1.310 + throw new NumberFormatException("Zero length BigInteger");
1.311 + }
1.312 + if (index1 == 0)
1.313 + sign = -1;
1.314 + } else
1.315 + throw new NumberFormatException("Illegal embedded sign character");
1.316 +
1.317 + // Skip leading zeros and compute number of digits in magnitude
1.318 + while (cursor < len &&
1.319 + Character.digit(val.charAt(cursor), radix) == 0)
1.320 + cursor++;
1.321 + if (cursor == len) {
1.322 + signum = 0;
1.323 + mag = ZERO.mag;
1.324 + return;
1.325 + }
1.326 +
1.327 + numDigits = len - cursor;
1.328 + signum = sign;
1.329 +
1.330 + // Pre-allocate array of expected size. May be too large but can
1.331 + // never be too small. Typically exact.
1.332 + int numBits = (int)(((numDigits * bitsPerDigit[radix]) >>> 10) + 1);
1.333 + int numWords = (numBits + 31) >>> 5;
1.334 + int[] magnitude = new int[numWords];
1.335 +
1.336 + // Process first (potentially short) digit group
1.337 + int firstGroupLen = numDigits % digitsPerInt[radix];
1.338 + if (firstGroupLen == 0)
1.339 + firstGroupLen = digitsPerInt[radix];
1.340 + String group = val.substring(cursor, cursor += firstGroupLen);
1.341 + magnitude[numWords - 1] = Integer.parseInt(group, radix);
1.342 + if (magnitude[numWords - 1] < 0)
1.343 + throw new NumberFormatException("Illegal digit");
1.344 +
1.345 + // Process remaining digit groups
1.346 + int superRadix = intRadix[radix];
1.347 + int groupVal = 0;
1.348 + while (cursor < len) {
1.349 + group = val.substring(cursor, cursor += digitsPerInt[radix]);
1.350 + groupVal = Integer.parseInt(group, radix);
1.351 + if (groupVal < 0)
1.352 + throw new NumberFormatException("Illegal digit");
1.353 + destructiveMulAdd(magnitude, superRadix, groupVal);
1.354 + }
1.355 + // Required for cases where the array was overallocated.
1.356 + mag = trustedStripLeadingZeroInts(magnitude);
1.357 + }
1.358 +
1.359 + // Constructs a new BigInteger using a char array with radix=10
1.360 + BigInteger(char[] val) {
1.361 + int cursor = 0, numDigits;
1.362 + int len = val.length;
1.363 +
1.364 + // Check for leading minus sign
1.365 + int sign = 1;
1.366 + if (val[0] == '-') {
1.367 + if (len == 1)
1.368 + throw new NumberFormatException("Zero length BigInteger");
1.369 + sign = -1;
1.370 + cursor = 1;
1.371 + } else if (val[0] == '+') {
1.372 + if (len == 1)
1.373 + throw new NumberFormatException("Zero length BigInteger");
1.374 + cursor = 1;
1.375 + }
1.376 +
1.377 + // Skip leading zeros and compute number of digits in magnitude
1.378 + while (cursor < len && Character.digit(val[cursor], 10) == 0)
1.379 + cursor++;
1.380 + if (cursor == len) {
1.381 + signum = 0;
1.382 + mag = ZERO.mag;
1.383 + return;
1.384 + }
1.385 +
1.386 + numDigits = len - cursor;
1.387 + signum = sign;
1.388 +
1.389 + // Pre-allocate array of expected size
1.390 + int numWords;
1.391 + if (len < 10) {
1.392 + numWords = 1;
1.393 + } else {
1.394 + int numBits = (int)(((numDigits * bitsPerDigit[10]) >>> 10) + 1);
1.395 + numWords = (numBits + 31) >>> 5;
1.396 + }
1.397 + int[] magnitude = new int[numWords];
1.398 +
1.399 + // Process first (potentially short) digit group
1.400 + int firstGroupLen = numDigits % digitsPerInt[10];
1.401 + if (firstGroupLen == 0)
1.402 + firstGroupLen = digitsPerInt[10];
1.403 + magnitude[numWords - 1] = parseInt(val, cursor, cursor += firstGroupLen);
1.404 +
1.405 + // Process remaining digit groups
1.406 + while (cursor < len) {
1.407 + int groupVal = parseInt(val, cursor, cursor += digitsPerInt[10]);
1.408 + destructiveMulAdd(magnitude, intRadix[10], groupVal);
1.409 + }
1.410 + mag = trustedStripLeadingZeroInts(magnitude);
1.411 + }
1.412 +
1.413 + // Create an integer with the digits between the two indexes
1.414 + // Assumes start < end. The result may be negative, but it
1.415 + // is to be treated as an unsigned value.
1.416 + private int parseInt(char[] source, int start, int end) {
1.417 + int result = Character.digit(source[start++], 10);
1.418 + if (result == -1)
1.419 + throw new NumberFormatException(new String(source));
1.420 +
1.421 + for (int index = start; index<end; index++) {
1.422 + int nextVal = Character.digit(source[index], 10);
1.423 + if (nextVal == -1)
1.424 + throw new NumberFormatException(new String(source));
1.425 + result = 10*result + nextVal;
1.426 + }
1.427 +
1.428 + return result;
1.429 + }
1.430 +
1.431 + // bitsPerDigit in the given radix times 1024
1.432 + // Rounded up to avoid underallocation.
1.433 + private static long bitsPerDigit[] = { 0, 0,
1.434 + 1024, 1624, 2048, 2378, 2648, 2875, 3072, 3247, 3402, 3543, 3672,
1.435 + 3790, 3899, 4001, 4096, 4186, 4271, 4350, 4426, 4498, 4567, 4633,
1.436 + 4696, 4756, 4814, 4870, 4923, 4975, 5025, 5074, 5120, 5166, 5210,
1.437 + 5253, 5295};
1.438 +
1.439 + // Multiply x array times word y in place, and add word z
1.440 + private static void destructiveMulAdd(int[] x, int y, int z) {
1.441 + // Perform the multiplication word by word
1.442 + long ylong = y & LONG_MASK;
1.443 + long zlong = z & LONG_MASK;
1.444 + int len = x.length;
1.445 +
1.446 + long product = 0;
1.447 + long carry = 0;
1.448 + for (int i = len-1; i >= 0; i--) {
1.449 + product = ylong * (x[i] & LONG_MASK) + carry;
1.450 + x[i] = (int)product;
1.451 + carry = product >>> 32;
1.452 + }
1.453 +
1.454 + // Perform the addition
1.455 + long sum = (x[len-1] & LONG_MASK) + zlong;
1.456 + x[len-1] = (int)sum;
1.457 + carry = sum >>> 32;
1.458 + for (int i = len-2; i >= 0; i--) {
1.459 + sum = (x[i] & LONG_MASK) + carry;
1.460 + x[i] = (int)sum;
1.461 + carry = sum >>> 32;
1.462 + }
1.463 + }
1.464 +
1.465 + /**
1.466 + * Translates the decimal String representation of a BigInteger into a
1.467 + * BigInteger. The String representation consists of an optional minus
1.468 + * sign followed by a sequence of one or more decimal digits. The
1.469 + * character-to-digit mapping is provided by {@code Character.digit}.
1.470 + * The String may not contain any extraneous characters (whitespace, for
1.471 + * example).
1.472 + *
1.473 + * @param val decimal String representation of BigInteger.
1.474 + * @throws NumberFormatException {@code val} is not a valid representation
1.475 + * of a BigInteger.
1.476 + * @see Character#digit
1.477 + */
1.478 + public BigInteger(String val) {
1.479 + this(val, 10);
1.480 + }
1.481 +
1.482 + /**
1.483 + * Constructs a randomly generated BigInteger, uniformly distributed over
1.484 + * the range 0 to (2<sup>{@code numBits}</sup> - 1), inclusive.
1.485 + * The uniformity of the distribution assumes that a fair source of random
1.486 + * bits is provided in {@code rnd}. Note that this constructor always
1.487 + * constructs a non-negative BigInteger.
1.488 + *
1.489 + * @param numBits maximum bitLength of the new BigInteger.
1.490 + * @param rnd source of randomness to be used in computing the new
1.491 + * BigInteger.
1.492 + * @throws IllegalArgumentException {@code numBits} is negative.
1.493 + * @see #bitLength()
1.494 + */
1.495 + public BigInteger(int numBits, Random rnd) {
1.496 + this(1, randomBits(numBits, rnd));
1.497 + }
1.498 +
1.499 + private static byte[] randomBits(int numBits, Random rnd) {
1.500 + if (numBits < 0)
1.501 + throw new IllegalArgumentException("numBits must be non-negative");
1.502 + int numBytes = (int)(((long)numBits+7)/8); // avoid overflow
1.503 + byte[] randomBits = new byte[numBytes];
1.504 +
1.505 + // Generate random bytes and mask out any excess bits
1.506 + if (numBytes > 0) {
1.507 + rnd.nextBytes(randomBits);
1.508 + int excessBits = 8*numBytes - numBits;
1.509 + randomBits[0] &= (1 << (8-excessBits)) - 1;
1.510 + }
1.511 + return randomBits;
1.512 + }
1.513 +
1.514 + /**
1.515 + * Constructs a randomly generated positive BigInteger that is probably
1.516 + * prime, with the specified bitLength.
1.517 + *
1.518 + * <p>It is recommended that the {@link #probablePrime probablePrime}
1.519 + * method be used in preference to this constructor unless there
1.520 + * is a compelling need to specify a certainty.
1.521 + *
1.522 + * @param bitLength bitLength of the returned BigInteger.
1.523 + * @param certainty a measure of the uncertainty that the caller is
1.524 + * willing to tolerate. The probability that the new BigInteger
1.525 + * represents a prime number will exceed
1.526 + * (1 - 1/2<sup>{@code certainty}</sup>). The execution time of
1.527 + * this constructor is proportional to the value of this parameter.
1.528 + * @param rnd source of random bits used to select candidates to be
1.529 + * tested for primality.
1.530 + * @throws ArithmeticException {@code bitLength < 2}.
1.531 + * @see #bitLength()
1.532 + */
1.533 + public BigInteger(int bitLength, int certainty, Random rnd) {
1.534 + BigInteger prime;
1.535 +
1.536 + if (bitLength < 2)
1.537 + throw new ArithmeticException("bitLength < 2");
1.538 + // The cutoff of 95 was chosen empirically for best performance
1.539 + prime = (bitLength < 95 ? smallPrime(bitLength, certainty, rnd)
1.540 + : largePrime(bitLength, certainty, rnd));
1.541 + signum = 1;
1.542 + mag = prime.mag;
1.543 + }
1.544 +
1.545 + // Minimum size in bits that the requested prime number has
1.546 + // before we use the large prime number generating algorithms
1.547 + private static final int SMALL_PRIME_THRESHOLD = 95;
1.548 +
1.549 + // Certainty required to meet the spec of probablePrime
1.550 + private static final int DEFAULT_PRIME_CERTAINTY = 100;
1.551 +
1.552 + /**
1.553 + * Returns a positive BigInteger that is probably prime, with the
1.554 + * specified bitLength. The probability that a BigInteger returned
1.555 + * by this method is composite does not exceed 2<sup>-100</sup>.
1.556 + *
1.557 + * @param bitLength bitLength of the returned BigInteger.
1.558 + * @param rnd source of random bits used to select candidates to be
1.559 + * tested for primality.
1.560 + * @return a BigInteger of {@code bitLength} bits that is probably prime
1.561 + * @throws ArithmeticException {@code bitLength < 2}.
1.562 + * @see #bitLength()
1.563 + * @since 1.4
1.564 + */
1.565 + public static BigInteger probablePrime(int bitLength, Random rnd) {
1.566 + if (bitLength < 2)
1.567 + throw new ArithmeticException("bitLength < 2");
1.568 +
1.569 + // The cutoff of 95 was chosen empirically for best performance
1.570 + return (bitLength < SMALL_PRIME_THRESHOLD ?
1.571 + smallPrime(bitLength, DEFAULT_PRIME_CERTAINTY, rnd) :
1.572 + largePrime(bitLength, DEFAULT_PRIME_CERTAINTY, rnd));
1.573 + }
1.574 +
1.575 + /**
1.576 + * Find a random number of the specified bitLength that is probably prime.
1.577 + * This method is used for smaller primes, its performance degrades on
1.578 + * larger bitlengths.
1.579 + *
1.580 + * This method assumes bitLength > 1.
1.581 + */
1.582 + private static BigInteger smallPrime(int bitLength, int certainty, Random rnd) {
1.583 + int magLen = (bitLength + 31) >>> 5;
1.584 + int temp[] = new int[magLen];
1.585 + int highBit = 1 << ((bitLength+31) & 0x1f); // High bit of high int
1.586 + int highMask = (highBit << 1) - 1; // Bits to keep in high int
1.587 +
1.588 + while(true) {
1.589 + // Construct a candidate
1.590 + for (int i=0; i<magLen; i++)
1.591 + temp[i] = rnd.nextInt();
1.592 + temp[0] = (temp[0] & highMask) | highBit; // Ensure exact length
1.593 + if (bitLength > 2)
1.594 + temp[magLen-1] |= 1; // Make odd if bitlen > 2
1.595 +
1.596 + BigInteger p = new BigInteger(temp, 1);
1.597 +
1.598 + // Do cheap "pre-test" if applicable
1.599 + if (bitLength > 6) {
1.600 + long r = p.remainder(SMALL_PRIME_PRODUCT).longValue();
1.601 + if ((r%3==0) || (r%5==0) || (r%7==0) || (r%11==0) ||
1.602 + (r%13==0) || (r%17==0) || (r%19==0) || (r%23==0) ||
1.603 + (r%29==0) || (r%31==0) || (r%37==0) || (r%41==0))
1.604 + continue; // Candidate is composite; try another
1.605 + }
1.606 +
1.607 + // All candidates of bitLength 2 and 3 are prime by this point
1.608 + if (bitLength < 4)
1.609 + return p;
1.610 +
1.611 + // Do expensive test if we survive pre-test (or it's inapplicable)
1.612 + if (p.primeToCertainty(certainty, rnd))
1.613 + return p;
1.614 + }
1.615 + }
1.616 +
1.617 + private static final BigInteger SMALL_PRIME_PRODUCT
1.618 + = valueOf(3L*5*7*11*13*17*19*23*29*31*37*41);
1.619 +
1.620 + /**
1.621 + * Find a random number of the specified bitLength that is probably prime.
1.622 + * This method is more appropriate for larger bitlengths since it uses
1.623 + * a sieve to eliminate most composites before using a more expensive
1.624 + * test.
1.625 + */
1.626 + private static BigInteger largePrime(int bitLength, int certainty, Random rnd) {
1.627 + BigInteger p;
1.628 + p = new BigInteger(bitLength, rnd).setBit(bitLength-1);
1.629 + p.mag[p.mag.length-1] &= 0xfffffffe;
1.630 +
1.631 + // Use a sieve length likely to contain the next prime number
1.632 + int searchLen = (bitLength / 20) * 64;
1.633 + BitSieve searchSieve = new BitSieve(p, searchLen);
1.634 + BigInteger candidate = searchSieve.retrieve(p, certainty, rnd);
1.635 +
1.636 + while ((candidate == null) || (candidate.bitLength() != bitLength)) {
1.637 + p = p.add(BigInteger.valueOf(2*searchLen));
1.638 + if (p.bitLength() != bitLength)
1.639 + p = new BigInteger(bitLength, rnd).setBit(bitLength-1);
1.640 + p.mag[p.mag.length-1] &= 0xfffffffe;
1.641 + searchSieve = new BitSieve(p, searchLen);
1.642 + candidate = searchSieve.retrieve(p, certainty, rnd);
1.643 + }
1.644 + return candidate;
1.645 + }
1.646 +
1.647 + /**
1.648 + * Returns the first integer greater than this {@code BigInteger} that
1.649 + * is probably prime. The probability that the number returned by this
1.650 + * method is composite does not exceed 2<sup>-100</sup>. This method will
1.651 + * never skip over a prime when searching: if it returns {@code p}, there
1.652 + * is no prime {@code q} such that {@code this < q < p}.
1.653 + *
1.654 + * @return the first integer greater than this {@code BigInteger} that
1.655 + * is probably prime.
1.656 + * @throws ArithmeticException {@code this < 0}.
1.657 + * @since 1.5
1.658 + */
1.659 + public BigInteger nextProbablePrime() {
1.660 + if (this.signum < 0)
1.661 + throw new ArithmeticException("start < 0: " + this);
1.662 +
1.663 + // Handle trivial cases
1.664 + if ((this.signum == 0) || this.equals(ONE))
1.665 + return TWO;
1.666 +
1.667 + BigInteger result = this.add(ONE);
1.668 +
1.669 + // Fastpath for small numbers
1.670 + if (result.bitLength() < SMALL_PRIME_THRESHOLD) {
1.671 +
1.672 + // Ensure an odd number
1.673 + if (!result.testBit(0))
1.674 + result = result.add(ONE);
1.675 +
1.676 + while(true) {
1.677 + // Do cheap "pre-test" if applicable
1.678 + if (result.bitLength() > 6) {
1.679 + long r = result.remainder(SMALL_PRIME_PRODUCT).longValue();
1.680 + if ((r%3==0) || (r%5==0) || (r%7==0) || (r%11==0) ||
1.681 + (r%13==0) || (r%17==0) || (r%19==0) || (r%23==0) ||
1.682 + (r%29==0) || (r%31==0) || (r%37==0) || (r%41==0)) {
1.683 + result = result.add(TWO);
1.684 + continue; // Candidate is composite; try another
1.685 + }
1.686 + }
1.687 +
1.688 + // All candidates of bitLength 2 and 3 are prime by this point
1.689 + if (result.bitLength() < 4)
1.690 + return result;
1.691 +
1.692 + // The expensive test
1.693 + if (result.primeToCertainty(DEFAULT_PRIME_CERTAINTY, null))
1.694 + return result;
1.695 +
1.696 + result = result.add(TWO);
1.697 + }
1.698 + }
1.699 +
1.700 + // Start at previous even number
1.701 + if (result.testBit(0))
1.702 + result = result.subtract(ONE);
1.703 +
1.704 + // Looking for the next large prime
1.705 + int searchLen = (result.bitLength() / 20) * 64;
1.706 +
1.707 + while(true) {
1.708 + BitSieve searchSieve = new BitSieve(result, searchLen);
1.709 + BigInteger candidate = searchSieve.retrieve(result,
1.710 + DEFAULT_PRIME_CERTAINTY, null);
1.711 + if (candidate != null)
1.712 + return candidate;
1.713 + result = result.add(BigInteger.valueOf(2 * searchLen));
1.714 + }
1.715 + }
1.716 +
1.717 + /**
1.718 + * Returns {@code true} if this BigInteger is probably prime,
1.719 + * {@code false} if it's definitely composite.
1.720 + *
1.721 + * This method assumes bitLength > 2.
1.722 + *
1.723 + * @param certainty a measure of the uncertainty that the caller is
1.724 + * willing to tolerate: if the call returns {@code true}
1.725 + * the probability that this BigInteger is prime exceeds
1.726 + * {@code (1 - 1/2<sup>certainty</sup>)}. The execution time of
1.727 + * this method is proportional to the value of this parameter.
1.728 + * @return {@code true} if this BigInteger is probably prime,
1.729 + * {@code false} if it's definitely composite.
1.730 + */
1.731 + boolean primeToCertainty(int certainty, Random random) {
1.732 + int rounds = 0;
1.733 + int n = (Math.min(certainty, Integer.MAX_VALUE-1)+1)/2;
1.734 +
1.735 + // The relationship between the certainty and the number of rounds
1.736 + // we perform is given in the draft standard ANSI X9.80, "PRIME
1.737 + // NUMBER GENERATION, PRIMALITY TESTING, AND PRIMALITY CERTIFICATES".
1.738 + int sizeInBits = this.bitLength();
1.739 + if (sizeInBits < 100) {
1.740 + rounds = 50;
1.741 + rounds = n < rounds ? n : rounds;
1.742 + return passesMillerRabin(rounds, random);
1.743 + }
1.744 +
1.745 + if (sizeInBits < 256) {
1.746 + rounds = 27;
1.747 + } else if (sizeInBits < 512) {
1.748 + rounds = 15;
1.749 + } else if (sizeInBits < 768) {
1.750 + rounds = 8;
1.751 + } else if (sizeInBits < 1024) {
1.752 + rounds = 4;
1.753 + } else {
1.754 + rounds = 2;
1.755 + }
1.756 + rounds = n < rounds ? n : rounds;
1.757 +
1.758 + return passesMillerRabin(rounds, random) && passesLucasLehmer();
1.759 + }
1.760 +
1.761 + /**
1.762 + * Returns true iff this BigInteger is a Lucas-Lehmer probable prime.
1.763 + *
1.764 + * The following assumptions are made:
1.765 + * This BigInteger is a positive, odd number.
1.766 + */
1.767 + private boolean passesLucasLehmer() {
1.768 + BigInteger thisPlusOne = this.add(ONE);
1.769 +
1.770 + // Step 1
1.771 + int d = 5;
1.772 + while (jacobiSymbol(d, this) != -1) {
1.773 + // 5, -7, 9, -11, ...
1.774 + d = (d<0) ? Math.abs(d)+2 : -(d+2);
1.775 + }
1.776 +
1.777 + // Step 2
1.778 + BigInteger u = lucasLehmerSequence(d, thisPlusOne, this);
1.779 +
1.780 + // Step 3
1.781 + return u.mod(this).equals(ZERO);
1.782 + }
1.783 +
1.784 + /**
1.785 + * Computes Jacobi(p,n).
1.786 + * Assumes n positive, odd, n>=3.
1.787 + */
1.788 + private static int jacobiSymbol(int p, BigInteger n) {
1.789 + if (p == 0)
1.790 + return 0;
1.791 +
1.792 + // Algorithm and comments adapted from Colin Plumb's C library.
1.793 + int j = 1;
1.794 + int u = n.mag[n.mag.length-1];
1.795 +
1.796 + // Make p positive
1.797 + if (p < 0) {
1.798 + p = -p;
1.799 + int n8 = u & 7;
1.800 + if ((n8 == 3) || (n8 == 7))
1.801 + j = -j; // 3 (011) or 7 (111) mod 8
1.802 + }
1.803 +
1.804 + // Get rid of factors of 2 in p
1.805 + while ((p & 3) == 0)
1.806 + p >>= 2;
1.807 + if ((p & 1) == 0) {
1.808 + p >>= 1;
1.809 + if (((u ^ (u>>1)) & 2) != 0)
1.810 + j = -j; // 3 (011) or 5 (101) mod 8
1.811 + }
1.812 + if (p == 1)
1.813 + return j;
1.814 + // Then, apply quadratic reciprocity
1.815 + if ((p & u & 2) != 0) // p = u = 3 (mod 4)?
1.816 + j = -j;
1.817 + // And reduce u mod p
1.818 + u = n.mod(BigInteger.valueOf(p)).intValue();
1.819 +
1.820 + // Now compute Jacobi(u,p), u < p
1.821 + while (u != 0) {
1.822 + while ((u & 3) == 0)
1.823 + u >>= 2;
1.824 + if ((u & 1) == 0) {
1.825 + u >>= 1;
1.826 + if (((p ^ (p>>1)) & 2) != 0)
1.827 + j = -j; // 3 (011) or 5 (101) mod 8
1.828 + }
1.829 + if (u == 1)
1.830 + return j;
1.831 + // Now both u and p are odd, so use quadratic reciprocity
1.832 + assert (u < p);
1.833 + int t = u; u = p; p = t;
1.834 + if ((u & p & 2) != 0) // u = p = 3 (mod 4)?
1.835 + j = -j;
1.836 + // Now u >= p, so it can be reduced
1.837 + u %= p;
1.838 + }
1.839 + return 0;
1.840 + }
1.841 +
1.842 + private static BigInteger lucasLehmerSequence(int z, BigInteger k, BigInteger n) {
1.843 + BigInteger d = BigInteger.valueOf(z);
1.844 + BigInteger u = ONE; BigInteger u2;
1.845 + BigInteger v = ONE; BigInteger v2;
1.846 +
1.847 + for (int i=k.bitLength()-2; i>=0; i--) {
1.848 + u2 = u.multiply(v).mod(n);
1.849 +
1.850 + v2 = v.square().add(d.multiply(u.square())).mod(n);
1.851 + if (v2.testBit(0))
1.852 + v2 = v2.subtract(n);
1.853 +
1.854 + v2 = v2.shiftRight(1);
1.855 +
1.856 + u = u2; v = v2;
1.857 + if (k.testBit(i)) {
1.858 + u2 = u.add(v).mod(n);
1.859 + if (u2.testBit(0))
1.860 + u2 = u2.subtract(n);
1.861 +
1.862 + u2 = u2.shiftRight(1);
1.863 + v2 = v.add(d.multiply(u)).mod(n);
1.864 + if (v2.testBit(0))
1.865 + v2 = v2.subtract(n);
1.866 + v2 = v2.shiftRight(1);
1.867 +
1.868 + u = u2; v = v2;
1.869 + }
1.870 + }
1.871 + return u;
1.872 + }
1.873 +
1.874 + private static volatile Random staticRandom;
1.875 +
1.876 + private static Random getSecureRandom() {
1.877 + if (staticRandom == null) {
1.878 + staticRandom = new java.security.SecureRandom();
1.879 + }
1.880 + return staticRandom;
1.881 + }
1.882 +
1.883 + /**
1.884 + * Returns true iff this BigInteger passes the specified number of
1.885 + * Miller-Rabin tests. This test is taken from the DSA spec (NIST FIPS
1.886 + * 186-2).
1.887 + *
1.888 + * The following assumptions are made:
1.889 + * This BigInteger is a positive, odd number greater than 2.
1.890 + * iterations<=50.
1.891 + */
1.892 + private boolean passesMillerRabin(int iterations, Random rnd) {
1.893 + // Find a and m such that m is odd and this == 1 + 2**a * m
1.894 + BigInteger thisMinusOne = this.subtract(ONE);
1.895 + BigInteger m = thisMinusOne;
1.896 + int a = m.getLowestSetBit();
1.897 + m = m.shiftRight(a);
1.898 +
1.899 + // Do the tests
1.900 + if (rnd == null) {
1.901 + rnd = getSecureRandom();
1.902 + }
1.903 + for (int i=0; i<iterations; i++) {
1.904 + // Generate a uniform random on (1, this)
1.905 + BigInteger b;
1.906 + do {
1.907 + b = new BigInteger(this.bitLength(), rnd);
1.908 + } while (b.compareTo(ONE) <= 0 || b.compareTo(this) >= 0);
1.909 +
1.910 + int j = 0;
1.911 + BigInteger z = b.modPow(m, this);
1.912 + while(!((j==0 && z.equals(ONE)) || z.equals(thisMinusOne))) {
1.913 + if (j>0 && z.equals(ONE) || ++j==a)
1.914 + return false;
1.915 + z = z.modPow(TWO, this);
1.916 + }
1.917 + }
1.918 + return true;
1.919 + }
1.920 +
1.921 + /**
1.922 + * This internal constructor differs from its public cousin
1.923 + * with the arguments reversed in two ways: it assumes that its
1.924 + * arguments are correct, and it doesn't copy the magnitude array.
1.925 + */
1.926 + BigInteger(int[] magnitude, int signum) {
1.927 + this.signum = (magnitude.length==0 ? 0 : signum);
1.928 + this.mag = magnitude;
1.929 + }
1.930 +
1.931 + /**
1.932 + * This private constructor is for internal use and assumes that its
1.933 + * arguments are correct.
1.934 + */
1.935 + private BigInteger(byte[] magnitude, int signum) {
1.936 + this.signum = (magnitude.length==0 ? 0 : signum);
1.937 + this.mag = stripLeadingZeroBytes(magnitude);
1.938 + }
1.939 +
1.940 + //Static Factory Methods
1.941 +
1.942 + /**
1.943 + * Returns a BigInteger whose value is equal to that of the
1.944 + * specified {@code long}. This "static factory method" is
1.945 + * provided in preference to a ({@code long}) constructor
1.946 + * because it allows for reuse of frequently used BigIntegers.
1.947 + *
1.948 + * @param val value of the BigInteger to return.
1.949 + * @return a BigInteger with the specified value.
1.950 + */
1.951 + public static BigInteger valueOf(long val) {
1.952 + // If -MAX_CONSTANT < val < MAX_CONSTANT, return stashed constant
1.953 + if (val == 0)
1.954 + return ZERO;
1.955 + if (val > 0 && val <= MAX_CONSTANT)
1.956 + return posConst[(int) val];
1.957 + else if (val < 0 && val >= -MAX_CONSTANT)
1.958 + return negConst[(int) -val];
1.959 +
1.960 + return new BigInteger(val);
1.961 + }
1.962 +
1.963 + /**
1.964 + * Constructs a BigInteger with the specified value, which may not be zero.
1.965 + */
1.966 + private BigInteger(long val) {
1.967 + if (val < 0) {
1.968 + val = -val;
1.969 + signum = -1;
1.970 + } else {
1.971 + signum = 1;
1.972 + }
1.973 +
1.974 + int highWord = (int)(val >>> 32);
1.975 + if (highWord==0) {
1.976 + mag = new int[1];
1.977 + mag[0] = (int)val;
1.978 + } else {
1.979 + mag = new int[2];
1.980 + mag[0] = highWord;
1.981 + mag[1] = (int)val;
1.982 + }
1.983 + }
1.984 +
1.985 + /**
1.986 + * Returns a BigInteger with the given two's complement representation.
1.987 + * Assumes that the input array will not be modified (the returned
1.988 + * BigInteger will reference the input array if feasible).
1.989 + */
1.990 + private static BigInteger valueOf(int val[]) {
1.991 + return (val[0]>0 ? new BigInteger(val, 1) : new BigInteger(val));
1.992 + }
1.993 +
1.994 + // Constants
1.995 +
1.996 + /**
1.997 + * Initialize static constant array when class is loaded.
1.998 + */
1.999 + private final static int MAX_CONSTANT = 16;
1.1000 + private static BigInteger posConst[] = new BigInteger[MAX_CONSTANT+1];
1.1001 + private static BigInteger negConst[] = new BigInteger[MAX_CONSTANT+1];
1.1002 + static {
1.1003 + for (int i = 1; i <= MAX_CONSTANT; i++) {
1.1004 + int[] magnitude = new int[1];
1.1005 + magnitude[0] = i;
1.1006 + posConst[i] = new BigInteger(magnitude, 1);
1.1007 + negConst[i] = new BigInteger(magnitude, -1);
1.1008 + }
1.1009 + }
1.1010 +
1.1011 + /**
1.1012 + * The BigInteger constant zero.
1.1013 + *
1.1014 + * @since 1.2
1.1015 + */
1.1016 + public static final BigInteger ZERO = new BigInteger(new int[0], 0);
1.1017 +
1.1018 + /**
1.1019 + * The BigInteger constant one.
1.1020 + *
1.1021 + * @since 1.2
1.1022 + */
1.1023 + public static final BigInteger ONE = valueOf(1);
1.1024 +
1.1025 + /**
1.1026 + * The BigInteger constant two. (Not exported.)
1.1027 + */
1.1028 + private static final BigInteger TWO = valueOf(2);
1.1029 +
1.1030 + /**
1.1031 + * The BigInteger constant ten.
1.1032 + *
1.1033 + * @since 1.5
1.1034 + */
1.1035 + public static final BigInteger TEN = valueOf(10);
1.1036 +
1.1037 + // Arithmetic Operations
1.1038 +
1.1039 + /**
1.1040 + * Returns a BigInteger whose value is {@code (this + val)}.
1.1041 + *
1.1042 + * @param val value to be added to this BigInteger.
1.1043 + * @return {@code this + val}
1.1044 + */
1.1045 + public BigInteger add(BigInteger val) {
1.1046 + if (val.signum == 0)
1.1047 + return this;
1.1048 + if (signum == 0)
1.1049 + return val;
1.1050 + if (val.signum == signum)
1.1051 + return new BigInteger(add(mag, val.mag), signum);
1.1052 +
1.1053 + int cmp = compareMagnitude(val);
1.1054 + if (cmp == 0)
1.1055 + return ZERO;
1.1056 + int[] resultMag = (cmp > 0 ? subtract(mag, val.mag)
1.1057 + : subtract(val.mag, mag));
1.1058 + resultMag = trustedStripLeadingZeroInts(resultMag);
1.1059 +
1.1060 + return new BigInteger(resultMag, cmp == signum ? 1 : -1);
1.1061 + }
1.1062 +
1.1063 + /**
1.1064 + * Adds the contents of the int arrays x and y. This method allocates
1.1065 + * a new int array to hold the answer and returns a reference to that
1.1066 + * array.
1.1067 + */
1.1068 + private static int[] add(int[] x, int[] y) {
1.1069 + // If x is shorter, swap the two arrays
1.1070 + if (x.length < y.length) {
1.1071 + int[] tmp = x;
1.1072 + x = y;
1.1073 + y = tmp;
1.1074 + }
1.1075 +
1.1076 + int xIndex = x.length;
1.1077 + int yIndex = y.length;
1.1078 + int result[] = new int[xIndex];
1.1079 + long sum = 0;
1.1080 +
1.1081 + // Add common parts of both numbers
1.1082 + while(yIndex > 0) {
1.1083 + sum = (x[--xIndex] & LONG_MASK) +
1.1084 + (y[--yIndex] & LONG_MASK) + (sum >>> 32);
1.1085 + result[xIndex] = (int)sum;
1.1086 + }
1.1087 +
1.1088 + // Copy remainder of longer number while carry propagation is required
1.1089 + boolean carry = (sum >>> 32 != 0);
1.1090 + while (xIndex > 0 && carry)
1.1091 + carry = ((result[--xIndex] = x[xIndex] + 1) == 0);
1.1092 +
1.1093 + // Copy remainder of longer number
1.1094 + while (xIndex > 0)
1.1095 + result[--xIndex] = x[xIndex];
1.1096 +
1.1097 + // Grow result if necessary
1.1098 + if (carry) {
1.1099 + int bigger[] = new int[result.length + 1];
1.1100 + System.arraycopy(result, 0, bigger, 1, result.length);
1.1101 + bigger[0] = 0x01;
1.1102 + return bigger;
1.1103 + }
1.1104 + return result;
1.1105 + }
1.1106 +
1.1107 + /**
1.1108 + * Returns a BigInteger whose value is {@code (this - val)}.
1.1109 + *
1.1110 + * @param val value to be subtracted from this BigInteger.
1.1111 + * @return {@code this - val}
1.1112 + */
1.1113 + public BigInteger subtract(BigInteger val) {
1.1114 + if (val.signum == 0)
1.1115 + return this;
1.1116 + if (signum == 0)
1.1117 + return val.negate();
1.1118 + if (val.signum != signum)
1.1119 + return new BigInteger(add(mag, val.mag), signum);
1.1120 +
1.1121 + int cmp = compareMagnitude(val);
1.1122 + if (cmp == 0)
1.1123 + return ZERO;
1.1124 + int[] resultMag = (cmp > 0 ? subtract(mag, val.mag)
1.1125 + : subtract(val.mag, mag));
1.1126 + resultMag = trustedStripLeadingZeroInts(resultMag);
1.1127 + return new BigInteger(resultMag, cmp == signum ? 1 : -1);
1.1128 + }
1.1129 +
1.1130 + /**
1.1131 + * Subtracts the contents of the second int arrays (little) from the
1.1132 + * first (big). The first int array (big) must represent a larger number
1.1133 + * than the second. This method allocates the space necessary to hold the
1.1134 + * answer.
1.1135 + */
1.1136 + private static int[] subtract(int[] big, int[] little) {
1.1137 + int bigIndex = big.length;
1.1138 + int result[] = new int[bigIndex];
1.1139 + int littleIndex = little.length;
1.1140 + long difference = 0;
1.1141 +
1.1142 + // Subtract common parts of both numbers
1.1143 + while(littleIndex > 0) {
1.1144 + difference = (big[--bigIndex] & LONG_MASK) -
1.1145 + (little[--littleIndex] & LONG_MASK) +
1.1146 + (difference >> 32);
1.1147 + result[bigIndex] = (int)difference;
1.1148 + }
1.1149 +
1.1150 + // Subtract remainder of longer number while borrow propagates
1.1151 + boolean borrow = (difference >> 32 != 0);
1.1152 + while (bigIndex > 0 && borrow)
1.1153 + borrow = ((result[--bigIndex] = big[bigIndex] - 1) == -1);
1.1154 +
1.1155 + // Copy remainder of longer number
1.1156 + while (bigIndex > 0)
1.1157 + result[--bigIndex] = big[bigIndex];
1.1158 +
1.1159 + return result;
1.1160 + }
1.1161 +
1.1162 + /**
1.1163 + * Returns a BigInteger whose value is {@code (this * val)}.
1.1164 + *
1.1165 + * @param val value to be multiplied by this BigInteger.
1.1166 + * @return {@code this * val}
1.1167 + */
1.1168 + public BigInteger multiply(BigInteger val) {
1.1169 + if (val.signum == 0 || signum == 0)
1.1170 + return ZERO;
1.1171 +
1.1172 + int[] result = multiplyToLen(mag, mag.length,
1.1173 + val.mag, val.mag.length, null);
1.1174 + result = trustedStripLeadingZeroInts(result);
1.1175 + return new BigInteger(result, signum == val.signum ? 1 : -1);
1.1176 + }
1.1177 +
1.1178 + /**
1.1179 + * Package private methods used by BigDecimal code to multiply a BigInteger
1.1180 + * with a long. Assumes v is not equal to INFLATED.
1.1181 + */
1.1182 + BigInteger multiply(long v) {
1.1183 + if (v == 0 || signum == 0)
1.1184 + return ZERO;
1.1185 + if (v == BigDecimal.INFLATED)
1.1186 + return multiply(BigInteger.valueOf(v));
1.1187 + int rsign = (v > 0 ? signum : -signum);
1.1188 + if (v < 0)
1.1189 + v = -v;
1.1190 + long dh = v >>> 32; // higher order bits
1.1191 + long dl = v & LONG_MASK; // lower order bits
1.1192 +
1.1193 + int xlen = mag.length;
1.1194 + int[] value = mag;
1.1195 + int[] rmag = (dh == 0L) ? (new int[xlen + 1]) : (new int[xlen + 2]);
1.1196 + long carry = 0;
1.1197 + int rstart = rmag.length - 1;
1.1198 + for (int i = xlen - 1; i >= 0; i--) {
1.1199 + long product = (value[i] & LONG_MASK) * dl + carry;
1.1200 + rmag[rstart--] = (int)product;
1.1201 + carry = product >>> 32;
1.1202 + }
1.1203 + rmag[rstart] = (int)carry;
1.1204 + if (dh != 0L) {
1.1205 + carry = 0;
1.1206 + rstart = rmag.length - 2;
1.1207 + for (int i = xlen - 1; i >= 0; i--) {
1.1208 + long product = (value[i] & LONG_MASK) * dh +
1.1209 + (rmag[rstart] & LONG_MASK) + carry;
1.1210 + rmag[rstart--] = (int)product;
1.1211 + carry = product >>> 32;
1.1212 + }
1.1213 + rmag[0] = (int)carry;
1.1214 + }
1.1215 + if (carry == 0L)
1.1216 + rmag = java.util.Arrays.copyOfRange(rmag, 1, rmag.length);
1.1217 + return new BigInteger(rmag, rsign);
1.1218 + }
1.1219 +
1.1220 + /**
1.1221 + * Multiplies int arrays x and y to the specified lengths and places
1.1222 + * the result into z. There will be no leading zeros in the resultant array.
1.1223 + */
1.1224 + private int[] multiplyToLen(int[] x, int xlen, int[] y, int ylen, int[] z) {
1.1225 + int xstart = xlen - 1;
1.1226 + int ystart = ylen - 1;
1.1227 +
1.1228 + if (z == null || z.length < (xlen+ ylen))
1.1229 + z = new int[xlen+ylen];
1.1230 +
1.1231 + long carry = 0;
1.1232 + for (int j=ystart, k=ystart+1+xstart; j>=0; j--, k--) {
1.1233 + long product = (y[j] & LONG_MASK) *
1.1234 + (x[xstart] & LONG_MASK) + carry;
1.1235 + z[k] = (int)product;
1.1236 + carry = product >>> 32;
1.1237 + }
1.1238 + z[xstart] = (int)carry;
1.1239 +
1.1240 + for (int i = xstart-1; i >= 0; i--) {
1.1241 + carry = 0;
1.1242 + for (int j=ystart, k=ystart+1+i; j>=0; j--, k--) {
1.1243 + long product = (y[j] & LONG_MASK) *
1.1244 + (x[i] & LONG_MASK) +
1.1245 + (z[k] & LONG_MASK) + carry;
1.1246 + z[k] = (int)product;
1.1247 + carry = product >>> 32;
1.1248 + }
1.1249 + z[i] = (int)carry;
1.1250 + }
1.1251 + return z;
1.1252 + }
1.1253 +
1.1254 + /**
1.1255 + * Returns a BigInteger whose value is {@code (this<sup>2</sup>)}.
1.1256 + *
1.1257 + * @return {@code this<sup>2</sup>}
1.1258 + */
1.1259 + private BigInteger square() {
1.1260 + if (signum == 0)
1.1261 + return ZERO;
1.1262 + int[] z = squareToLen(mag, mag.length, null);
1.1263 + return new BigInteger(trustedStripLeadingZeroInts(z), 1);
1.1264 + }
1.1265 +
1.1266 + /**
1.1267 + * Squares the contents of the int array x. The result is placed into the
1.1268 + * int array z. The contents of x are not changed.
1.1269 + */
1.1270 + private static final int[] squareToLen(int[] x, int len, int[] z) {
1.1271 + /*
1.1272 + * The algorithm used here is adapted from Colin Plumb's C library.
1.1273 + * Technique: Consider the partial products in the multiplication
1.1274 + * of "abcde" by itself:
1.1275 + *
1.1276 + * a b c d e
1.1277 + * * a b c d e
1.1278 + * ==================
1.1279 + * ae be ce de ee
1.1280 + * ad bd cd dd de
1.1281 + * ac bc cc cd ce
1.1282 + * ab bb bc bd be
1.1283 + * aa ab ac ad ae
1.1284 + *
1.1285 + * Note that everything above the main diagonal:
1.1286 + * ae be ce de = (abcd) * e
1.1287 + * ad bd cd = (abc) * d
1.1288 + * ac bc = (ab) * c
1.1289 + * ab = (a) * b
1.1290 + *
1.1291 + * is a copy of everything below the main diagonal:
1.1292 + * de
1.1293 + * cd ce
1.1294 + * bc bd be
1.1295 + * ab ac ad ae
1.1296 + *
1.1297 + * Thus, the sum is 2 * (off the diagonal) + diagonal.
1.1298 + *
1.1299 + * This is accumulated beginning with the diagonal (which
1.1300 + * consist of the squares of the digits of the input), which is then
1.1301 + * divided by two, the off-diagonal added, and multiplied by two
1.1302 + * again. The low bit is simply a copy of the low bit of the
1.1303 + * input, so it doesn't need special care.
1.1304 + */
1.1305 + int zlen = len << 1;
1.1306 + if (z == null || z.length < zlen)
1.1307 + z = new int[zlen];
1.1308 +
1.1309 + // Store the squares, right shifted one bit (i.e., divided by 2)
1.1310 + int lastProductLowWord = 0;
1.1311 + for (int j=0, i=0; j<len; j++) {
1.1312 + long piece = (x[j] & LONG_MASK);
1.1313 + long product = piece * piece;
1.1314 + z[i++] = (lastProductLowWord << 31) | (int)(product >>> 33);
1.1315 + z[i++] = (int)(product >>> 1);
1.1316 + lastProductLowWord = (int)product;
1.1317 + }
1.1318 +
1.1319 + // Add in off-diagonal sums
1.1320 + for (int i=len, offset=1; i>0; i--, offset+=2) {
1.1321 + int t = x[i-1];
1.1322 + t = mulAdd(z, x, offset, i-1, t);
1.1323 + addOne(z, offset-1, i, t);
1.1324 + }
1.1325 +
1.1326 + // Shift back up and set low bit
1.1327 + primitiveLeftShift(z, zlen, 1);
1.1328 + z[zlen-1] |= x[len-1] & 1;
1.1329 +
1.1330 + return z;
1.1331 + }
1.1332 +
1.1333 + /**
1.1334 + * Returns a BigInteger whose value is {@code (this / val)}.
1.1335 + *
1.1336 + * @param val value by which this BigInteger is to be divided.
1.1337 + * @return {@code this / val}
1.1338 + * @throws ArithmeticException if {@code val} is zero.
1.1339 + */
1.1340 + public BigInteger divide(BigInteger val) {
1.1341 + MutableBigInteger q = new MutableBigInteger(),
1.1342 + a = new MutableBigInteger(this.mag),
1.1343 + b = new MutableBigInteger(val.mag);
1.1344 +
1.1345 + a.divide(b, q);
1.1346 + return q.toBigInteger(this.signum == val.signum ? 1 : -1);
1.1347 + }
1.1348 +
1.1349 + /**
1.1350 + * Returns an array of two BigIntegers containing {@code (this / val)}
1.1351 + * followed by {@code (this % val)}.
1.1352 + *
1.1353 + * @param val value by which this BigInteger is to be divided, and the
1.1354 + * remainder computed.
1.1355 + * @return an array of two BigIntegers: the quotient {@code (this / val)}
1.1356 + * is the initial element, and the remainder {@code (this % val)}
1.1357 + * is the final element.
1.1358 + * @throws ArithmeticException if {@code val} is zero.
1.1359 + */
1.1360 + public BigInteger[] divideAndRemainder(BigInteger val) {
1.1361 + BigInteger[] result = new BigInteger[2];
1.1362 + MutableBigInteger q = new MutableBigInteger(),
1.1363 + a = new MutableBigInteger(this.mag),
1.1364 + b = new MutableBigInteger(val.mag);
1.1365 + MutableBigInteger r = a.divide(b, q);
1.1366 + result[0] = q.toBigInteger(this.signum == val.signum ? 1 : -1);
1.1367 + result[1] = r.toBigInteger(this.signum);
1.1368 + return result;
1.1369 + }
1.1370 +
1.1371 + /**
1.1372 + * Returns a BigInteger whose value is {@code (this % val)}.
1.1373 + *
1.1374 + * @param val value by which this BigInteger is to be divided, and the
1.1375 + * remainder computed.
1.1376 + * @return {@code this % val}
1.1377 + * @throws ArithmeticException if {@code val} is zero.
1.1378 + */
1.1379 + public BigInteger remainder(BigInteger val) {
1.1380 + MutableBigInteger q = new MutableBigInteger(),
1.1381 + a = new MutableBigInteger(this.mag),
1.1382 + b = new MutableBigInteger(val.mag);
1.1383 +
1.1384 + return a.divide(b, q).toBigInteger(this.signum);
1.1385 + }
1.1386 +
1.1387 + /**
1.1388 + * Returns a BigInteger whose value is <tt>(this<sup>exponent</sup>)</tt>.
1.1389 + * Note that {@code exponent} is an integer rather than a BigInteger.
1.1390 + *
1.1391 + * @param exponent exponent to which this BigInteger is to be raised.
1.1392 + * @return <tt>this<sup>exponent</sup></tt>
1.1393 + * @throws ArithmeticException {@code exponent} is negative. (This would
1.1394 + * cause the operation to yield a non-integer value.)
1.1395 + */
1.1396 + public BigInteger pow(int exponent) {
1.1397 + if (exponent < 0)
1.1398 + throw new ArithmeticException("Negative exponent");
1.1399 + if (signum==0)
1.1400 + return (exponent==0 ? ONE : this);
1.1401 +
1.1402 + // Perform exponentiation using repeated squaring trick
1.1403 + int newSign = (signum<0 && (exponent&1)==1 ? -1 : 1);
1.1404 + int[] baseToPow2 = this.mag;
1.1405 + int[] result = {1};
1.1406 +
1.1407 + while (exponent != 0) {
1.1408 + if ((exponent & 1)==1) {
1.1409 + result = multiplyToLen(result, result.length,
1.1410 + baseToPow2, baseToPow2.length, null);
1.1411 + result = trustedStripLeadingZeroInts(result);
1.1412 + }
1.1413 + if ((exponent >>>= 1) != 0) {
1.1414 + baseToPow2 = squareToLen(baseToPow2, baseToPow2.length, null);
1.1415 + baseToPow2 = trustedStripLeadingZeroInts(baseToPow2);
1.1416 + }
1.1417 + }
1.1418 + return new BigInteger(result, newSign);
1.1419 + }
1.1420 +
1.1421 + /**
1.1422 + * Returns a BigInteger whose value is the greatest common divisor of
1.1423 + * {@code abs(this)} and {@code abs(val)}. Returns 0 if
1.1424 + * {@code this==0 && val==0}.
1.1425 + *
1.1426 + * @param val value with which the GCD is to be computed.
1.1427 + * @return {@code GCD(abs(this), abs(val))}
1.1428 + */
1.1429 + public BigInteger gcd(BigInteger val) {
1.1430 + if (val.signum == 0)
1.1431 + return this.abs();
1.1432 + else if (this.signum == 0)
1.1433 + return val.abs();
1.1434 +
1.1435 + MutableBigInteger a = new MutableBigInteger(this);
1.1436 + MutableBigInteger b = new MutableBigInteger(val);
1.1437 +
1.1438 + MutableBigInteger result = a.hybridGCD(b);
1.1439 +
1.1440 + return result.toBigInteger(1);
1.1441 + }
1.1442 +
1.1443 + /**
1.1444 + * Package private method to return bit length for an integer.
1.1445 + */
1.1446 + static int bitLengthForInt(int n) {
1.1447 + return 32 - Integer.numberOfLeadingZeros(n);
1.1448 + }
1.1449 +
1.1450 + /**
1.1451 + * Left shift int array a up to len by n bits. Returns the array that
1.1452 + * results from the shift since space may have to be reallocated.
1.1453 + */
1.1454 + private static int[] leftShift(int[] a, int len, int n) {
1.1455 + int nInts = n >>> 5;
1.1456 + int nBits = n&0x1F;
1.1457 + int bitsInHighWord = bitLengthForInt(a[0]);
1.1458 +
1.1459 + // If shift can be done without recopy, do so
1.1460 + if (n <= (32-bitsInHighWord)) {
1.1461 + primitiveLeftShift(a, len, nBits);
1.1462 + return a;
1.1463 + } else { // Array must be resized
1.1464 + if (nBits <= (32-bitsInHighWord)) {
1.1465 + int result[] = new int[nInts+len];
1.1466 + for (int i=0; i<len; i++)
1.1467 + result[i] = a[i];
1.1468 + primitiveLeftShift(result, result.length, nBits);
1.1469 + return result;
1.1470 + } else {
1.1471 + int result[] = new int[nInts+len+1];
1.1472 + for (int i=0; i<len; i++)
1.1473 + result[i] = a[i];
1.1474 + primitiveRightShift(result, result.length, 32 - nBits);
1.1475 + return result;
1.1476 + }
1.1477 + }
1.1478 + }
1.1479 +
1.1480 + // shifts a up to len right n bits assumes no leading zeros, 0<n<32
1.1481 + static void primitiveRightShift(int[] a, int len, int n) {
1.1482 + int n2 = 32 - n;
1.1483 + for (int i=len-1, c=a[i]; i>0; i--) {
1.1484 + int b = c;
1.1485 + c = a[i-1];
1.1486 + a[i] = (c << n2) | (b >>> n);
1.1487 + }
1.1488 + a[0] >>>= n;
1.1489 + }
1.1490 +
1.1491 + // shifts a up to len left n bits assumes no leading zeros, 0<=n<32
1.1492 + static void primitiveLeftShift(int[] a, int len, int n) {
1.1493 + if (len == 0 || n == 0)
1.1494 + return;
1.1495 +
1.1496 + int n2 = 32 - n;
1.1497 + for (int i=0, c=a[i], m=i+len-1; i<m; i++) {
1.1498 + int b = c;
1.1499 + c = a[i+1];
1.1500 + a[i] = (b << n) | (c >>> n2);
1.1501 + }
1.1502 + a[len-1] <<= n;
1.1503 + }
1.1504 +
1.1505 + /**
1.1506 + * Calculate bitlength of contents of the first len elements an int array,
1.1507 + * assuming there are no leading zero ints.
1.1508 + */
1.1509 + private static int bitLength(int[] val, int len) {
1.1510 + if (len == 0)
1.1511 + return 0;
1.1512 + return ((len - 1) << 5) + bitLengthForInt(val[0]);
1.1513 + }
1.1514 +
1.1515 + /**
1.1516 + * Returns a BigInteger whose value is the absolute value of this
1.1517 + * BigInteger.
1.1518 + *
1.1519 + * @return {@code abs(this)}
1.1520 + */
1.1521 + public BigInteger abs() {
1.1522 + return (signum >= 0 ? this : this.negate());
1.1523 + }
1.1524 +
1.1525 + /**
1.1526 + * Returns a BigInteger whose value is {@code (-this)}.
1.1527 + *
1.1528 + * @return {@code -this}
1.1529 + */
1.1530 + public BigInteger negate() {
1.1531 + return new BigInteger(this.mag, -this.signum);
1.1532 + }
1.1533 +
1.1534 + /**
1.1535 + * Returns the signum function of this BigInteger.
1.1536 + *
1.1537 + * @return -1, 0 or 1 as the value of this BigInteger is negative, zero or
1.1538 + * positive.
1.1539 + */
1.1540 + public int signum() {
1.1541 + return this.signum;
1.1542 + }
1.1543 +
1.1544 + // Modular Arithmetic Operations
1.1545 +
1.1546 + /**
1.1547 + * Returns a BigInteger whose value is {@code (this mod m}). This method
1.1548 + * differs from {@code remainder} in that it always returns a
1.1549 + * <i>non-negative</i> BigInteger.
1.1550 + *
1.1551 + * @param m the modulus.
1.1552 + * @return {@code this mod m}
1.1553 + * @throws ArithmeticException {@code m} ≤ 0
1.1554 + * @see #remainder
1.1555 + */
1.1556 + public BigInteger mod(BigInteger m) {
1.1557 + if (m.signum <= 0)
1.1558 + throw new ArithmeticException("BigInteger: modulus not positive");
1.1559 +
1.1560 + BigInteger result = this.remainder(m);
1.1561 + return (result.signum >= 0 ? result : result.add(m));
1.1562 + }
1.1563 +
1.1564 + /**
1.1565 + * Returns a BigInteger whose value is
1.1566 + * <tt>(this<sup>exponent</sup> mod m)</tt>. (Unlike {@code pow}, this
1.1567 + * method permits negative exponents.)
1.1568 + *
1.1569 + * @param exponent the exponent.
1.1570 + * @param m the modulus.
1.1571 + * @return <tt>this<sup>exponent</sup> mod m</tt>
1.1572 + * @throws ArithmeticException {@code m} ≤ 0 or the exponent is
1.1573 + * negative and this BigInteger is not <i>relatively
1.1574 + * prime</i> to {@code m}.
1.1575 + * @see #modInverse
1.1576 + */
1.1577 + public BigInteger modPow(BigInteger exponent, BigInteger m) {
1.1578 + if (m.signum <= 0)
1.1579 + throw new ArithmeticException("BigInteger: modulus not positive");
1.1580 +
1.1581 + // Trivial cases
1.1582 + if (exponent.signum == 0)
1.1583 + return (m.equals(ONE) ? ZERO : ONE);
1.1584 +
1.1585 + if (this.equals(ONE))
1.1586 + return (m.equals(ONE) ? ZERO : ONE);
1.1587 +
1.1588 + if (this.equals(ZERO) && exponent.signum >= 0)
1.1589 + return ZERO;
1.1590 +
1.1591 + if (this.equals(negConst[1]) && (!exponent.testBit(0)))
1.1592 + return (m.equals(ONE) ? ZERO : ONE);
1.1593 +
1.1594 + boolean invertResult;
1.1595 + if ((invertResult = (exponent.signum < 0)))
1.1596 + exponent = exponent.negate();
1.1597 +
1.1598 + BigInteger base = (this.signum < 0 || this.compareTo(m) >= 0
1.1599 + ? this.mod(m) : this);
1.1600 + BigInteger result;
1.1601 + if (m.testBit(0)) { // odd modulus
1.1602 + result = base.oddModPow(exponent, m);
1.1603 + } else {
1.1604 + /*
1.1605 + * Even modulus. Tear it into an "odd part" (m1) and power of two
1.1606 + * (m2), exponentiate mod m1, manually exponentiate mod m2, and
1.1607 + * use Chinese Remainder Theorem to combine results.
1.1608 + */
1.1609 +
1.1610 + // Tear m apart into odd part (m1) and power of 2 (m2)
1.1611 + int p = m.getLowestSetBit(); // Max pow of 2 that divides m
1.1612 +
1.1613 + BigInteger m1 = m.shiftRight(p); // m/2**p
1.1614 + BigInteger m2 = ONE.shiftLeft(p); // 2**p
1.1615 +
1.1616 + // Calculate new base from m1
1.1617 + BigInteger base2 = (this.signum < 0 || this.compareTo(m1) >= 0
1.1618 + ? this.mod(m1) : this);
1.1619 +
1.1620 + // Caculate (base ** exponent) mod m1.
1.1621 + BigInteger a1 = (m1.equals(ONE) ? ZERO :
1.1622 + base2.oddModPow(exponent, m1));
1.1623 +
1.1624 + // Calculate (this ** exponent) mod m2
1.1625 + BigInteger a2 = base.modPow2(exponent, p);
1.1626 +
1.1627 + // Combine results using Chinese Remainder Theorem
1.1628 + BigInteger y1 = m2.modInverse(m1);
1.1629 + BigInteger y2 = m1.modInverse(m2);
1.1630 +
1.1631 + result = a1.multiply(m2).multiply(y1).add
1.1632 + (a2.multiply(m1).multiply(y2)).mod(m);
1.1633 + }
1.1634 +
1.1635 + return (invertResult ? result.modInverse(m) : result);
1.1636 + }
1.1637 +
1.1638 + static int[] bnExpModThreshTable = {7, 25, 81, 241, 673, 1793,
1.1639 + Integer.MAX_VALUE}; // Sentinel
1.1640 +
1.1641 + /**
1.1642 + * Returns a BigInteger whose value is x to the power of y mod z.
1.1643 + * Assumes: z is odd && x < z.
1.1644 + */
1.1645 + private BigInteger oddModPow(BigInteger y, BigInteger z) {
1.1646 + /*
1.1647 + * The algorithm is adapted from Colin Plumb's C library.
1.1648 + *
1.1649 + * The window algorithm:
1.1650 + * The idea is to keep a running product of b1 = n^(high-order bits of exp)
1.1651 + * and then keep appending exponent bits to it. The following patterns
1.1652 + * apply to a 3-bit window (k = 3):
1.1653 + * To append 0: square
1.1654 + * To append 1: square, multiply by n^1
1.1655 + * To append 10: square, multiply by n^1, square
1.1656 + * To append 11: square, square, multiply by n^3
1.1657 + * To append 100: square, multiply by n^1, square, square
1.1658 + * To append 101: square, square, square, multiply by n^5
1.1659 + * To append 110: square, square, multiply by n^3, square
1.1660 + * To append 111: square, square, square, multiply by n^7
1.1661 + *
1.1662 + * Since each pattern involves only one multiply, the longer the pattern
1.1663 + * the better, except that a 0 (no multiplies) can be appended directly.
1.1664 + * We precompute a table of odd powers of n, up to 2^k, and can then
1.1665 + * multiply k bits of exponent at a time. Actually, assuming random
1.1666 + * exponents, there is on average one zero bit between needs to
1.1667 + * multiply (1/2 of the time there's none, 1/4 of the time there's 1,
1.1668 + * 1/8 of the time, there's 2, 1/32 of the time, there's 3, etc.), so
1.1669 + * you have to do one multiply per k+1 bits of exponent.
1.1670 + *
1.1671 + * The loop walks down the exponent, squaring the result buffer as
1.1672 + * it goes. There is a wbits+1 bit lookahead buffer, buf, that is
1.1673 + * filled with the upcoming exponent bits. (What is read after the
1.1674 + * end of the exponent is unimportant, but it is filled with zero here.)
1.1675 + * When the most-significant bit of this buffer becomes set, i.e.
1.1676 + * (buf & tblmask) != 0, we have to decide what pattern to multiply
1.1677 + * by, and when to do it. We decide, remember to do it in future
1.1678 + * after a suitable number of squarings have passed (e.g. a pattern
1.1679 + * of "100" in the buffer requires that we multiply by n^1 immediately;
1.1680 + * a pattern of "110" calls for multiplying by n^3 after one more
1.1681 + * squaring), clear the buffer, and continue.
1.1682 + *
1.1683 + * When we start, there is one more optimization: the result buffer
1.1684 + * is implcitly one, so squaring it or multiplying by it can be
1.1685 + * optimized away. Further, if we start with a pattern like "100"
1.1686 + * in the lookahead window, rather than placing n into the buffer
1.1687 + * and then starting to square it, we have already computed n^2
1.1688 + * to compute the odd-powers table, so we can place that into
1.1689 + * the buffer and save a squaring.
1.1690 + *
1.1691 + * This means that if you have a k-bit window, to compute n^z,
1.1692 + * where z is the high k bits of the exponent, 1/2 of the time
1.1693 + * it requires no squarings. 1/4 of the time, it requires 1
1.1694 + * squaring, ... 1/2^(k-1) of the time, it reqires k-2 squarings.
1.1695 + * And the remaining 1/2^(k-1) of the time, the top k bits are a
1.1696 + * 1 followed by k-1 0 bits, so it again only requires k-2
1.1697 + * squarings, not k-1. The average of these is 1. Add that
1.1698 + * to the one squaring we have to do to compute the table,
1.1699 + * and you'll see that a k-bit window saves k-2 squarings
1.1700 + * as well as reducing the multiplies. (It actually doesn't
1.1701 + * hurt in the case k = 1, either.)
1.1702 + */
1.1703 + // Special case for exponent of one
1.1704 + if (y.equals(ONE))
1.1705 + return this;
1.1706 +
1.1707 + // Special case for base of zero
1.1708 + if (signum==0)
1.1709 + return ZERO;
1.1710 +
1.1711 + int[] base = mag.clone();
1.1712 + int[] exp = y.mag;
1.1713 + int[] mod = z.mag;
1.1714 + int modLen = mod.length;
1.1715 +
1.1716 + // Select an appropriate window size
1.1717 + int wbits = 0;
1.1718 + int ebits = bitLength(exp, exp.length);
1.1719 + // if exponent is 65537 (0x10001), use minimum window size
1.1720 + if ((ebits != 17) || (exp[0] != 65537)) {
1.1721 + while (ebits > bnExpModThreshTable[wbits]) {
1.1722 + wbits++;
1.1723 + }
1.1724 + }
1.1725 +
1.1726 + // Calculate appropriate table size
1.1727 + int tblmask = 1 << wbits;
1.1728 +
1.1729 + // Allocate table for precomputed odd powers of base in Montgomery form
1.1730 + int[][] table = new int[tblmask][];
1.1731 + for (int i=0; i<tblmask; i++)
1.1732 + table[i] = new int[modLen];
1.1733 +
1.1734 + // Compute the modular inverse
1.1735 + int inv = -MutableBigInteger.inverseMod32(mod[modLen-1]);
1.1736 +
1.1737 + // Convert base to Montgomery form
1.1738 + int[] a = leftShift(base, base.length, modLen << 5);
1.1739 +
1.1740 + MutableBigInteger q = new MutableBigInteger(),
1.1741 + a2 = new MutableBigInteger(a),
1.1742 + b2 = new MutableBigInteger(mod);
1.1743 +
1.1744 + MutableBigInteger r= a2.divide(b2, q);
1.1745 + table[0] = r.toIntArray();
1.1746 +
1.1747 + // Pad table[0] with leading zeros so its length is at least modLen
1.1748 + if (table[0].length < modLen) {
1.1749 + int offset = modLen - table[0].length;
1.1750 + int[] t2 = new int[modLen];
1.1751 + for (int i=0; i<table[0].length; i++)
1.1752 + t2[i+offset] = table[0][i];
1.1753 + table[0] = t2;
1.1754 + }
1.1755 +
1.1756 + // Set b to the square of the base
1.1757 + int[] b = squareToLen(table[0], modLen, null);
1.1758 + b = montReduce(b, mod, modLen, inv);
1.1759 +
1.1760 + // Set t to high half of b
1.1761 + int[] t = new int[modLen];
1.1762 + for(int i=0; i<modLen; i++)
1.1763 + t[i] = b[i];
1.1764 +
1.1765 + // Fill in the table with odd powers of the base
1.1766 + for (int i=1; i<tblmask; i++) {
1.1767 + int[] prod = multiplyToLen(t, modLen, table[i-1], modLen, null);
1.1768 + table[i] = montReduce(prod, mod, modLen, inv);
1.1769 + }
1.1770 +
1.1771 + // Pre load the window that slides over the exponent
1.1772 + int bitpos = 1 << ((ebits-1) & (32-1));
1.1773 +
1.1774 + int buf = 0;
1.1775 + int elen = exp.length;
1.1776 + int eIndex = 0;
1.1777 + for (int i = 0; i <= wbits; i++) {
1.1778 + buf = (buf << 1) | (((exp[eIndex] & bitpos) != 0)?1:0);
1.1779 + bitpos >>>= 1;
1.1780 + if (bitpos == 0) {
1.1781 + eIndex++;
1.1782 + bitpos = 1 << (32-1);
1.1783 + elen--;
1.1784 + }
1.1785 + }
1.1786 +
1.1787 + int multpos = ebits;
1.1788 +
1.1789 + // The first iteration, which is hoisted out of the main loop
1.1790 + ebits--;
1.1791 + boolean isone = true;
1.1792 +
1.1793 + multpos = ebits - wbits;
1.1794 + while ((buf & 1) == 0) {
1.1795 + buf >>>= 1;
1.1796 + multpos++;
1.1797 + }
1.1798 +
1.1799 + int[] mult = table[buf >>> 1];
1.1800 +
1.1801 + buf = 0;
1.1802 + if (multpos == ebits)
1.1803 + isone = false;
1.1804 +
1.1805 + // The main loop
1.1806 + while(true) {
1.1807 + ebits--;
1.1808 + // Advance the window
1.1809 + buf <<= 1;
1.1810 +
1.1811 + if (elen != 0) {
1.1812 + buf |= ((exp[eIndex] & bitpos) != 0) ? 1 : 0;
1.1813 + bitpos >>>= 1;
1.1814 + if (bitpos == 0) {
1.1815 + eIndex++;
1.1816 + bitpos = 1 << (32-1);
1.1817 + elen--;
1.1818 + }
1.1819 + }
1.1820 +
1.1821 + // Examine the window for pending multiplies
1.1822 + if ((buf & tblmask) != 0) {
1.1823 + multpos = ebits - wbits;
1.1824 + while ((buf & 1) == 0) {
1.1825 + buf >>>= 1;
1.1826 + multpos++;
1.1827 + }
1.1828 + mult = table[buf >>> 1];
1.1829 + buf = 0;
1.1830 + }
1.1831 +
1.1832 + // Perform multiply
1.1833 + if (ebits == multpos) {
1.1834 + if (isone) {
1.1835 + b = mult.clone();
1.1836 + isone = false;
1.1837 + } else {
1.1838 + t = b;
1.1839 + a = multiplyToLen(t, modLen, mult, modLen, a);
1.1840 + a = montReduce(a, mod, modLen, inv);
1.1841 + t = a; a = b; b = t;
1.1842 + }
1.1843 + }
1.1844 +
1.1845 + // Check if done
1.1846 + if (ebits == 0)
1.1847 + break;
1.1848 +
1.1849 + // Square the input
1.1850 + if (!isone) {
1.1851 + t = b;
1.1852 + a = squareToLen(t, modLen, a);
1.1853 + a = montReduce(a, mod, modLen, inv);
1.1854 + t = a; a = b; b = t;
1.1855 + }
1.1856 + }
1.1857 +
1.1858 + // Convert result out of Montgomery form and return
1.1859 + int[] t2 = new int[2*modLen];
1.1860 + for(int i=0; i<modLen; i++)
1.1861 + t2[i+modLen] = b[i];
1.1862 +
1.1863 + b = montReduce(t2, mod, modLen, inv);
1.1864 +
1.1865 + t2 = new int[modLen];
1.1866 + for(int i=0; i<modLen; i++)
1.1867 + t2[i] = b[i];
1.1868 +
1.1869 + return new BigInteger(1, t2);
1.1870 + }
1.1871 +
1.1872 + /**
1.1873 + * Montgomery reduce n, modulo mod. This reduces modulo mod and divides
1.1874 + * by 2^(32*mlen). Adapted from Colin Plumb's C library.
1.1875 + */
1.1876 + private static int[] montReduce(int[] n, int[] mod, int mlen, int inv) {
1.1877 + int c=0;
1.1878 + int len = mlen;
1.1879 + int offset=0;
1.1880 +
1.1881 + do {
1.1882 + int nEnd = n[n.length-1-offset];
1.1883 + int carry = mulAdd(n, mod, offset, mlen, inv * nEnd);
1.1884 + c += addOne(n, offset, mlen, carry);
1.1885 + offset++;
1.1886 + } while(--len > 0);
1.1887 +
1.1888 + while(c>0)
1.1889 + c += subN(n, mod, mlen);
1.1890 +
1.1891 + while (intArrayCmpToLen(n, mod, mlen) >= 0)
1.1892 + subN(n, mod, mlen);
1.1893 +
1.1894 + return n;
1.1895 + }
1.1896 +
1.1897 +
1.1898 + /*
1.1899 + * Returns -1, 0 or +1 as big-endian unsigned int array arg1 is less than,
1.1900 + * equal to, or greater than arg2 up to length len.
1.1901 + */
1.1902 + private static int intArrayCmpToLen(int[] arg1, int[] arg2, int len) {
1.1903 + for (int i=0; i<len; i++) {
1.1904 + long b1 = arg1[i] & LONG_MASK;
1.1905 + long b2 = arg2[i] & LONG_MASK;
1.1906 + if (b1 < b2)
1.1907 + return -1;
1.1908 + if (b1 > b2)
1.1909 + return 1;
1.1910 + }
1.1911 + return 0;
1.1912 + }
1.1913 +
1.1914 + /**
1.1915 + * Subtracts two numbers of same length, returning borrow.
1.1916 + */
1.1917 + private static int subN(int[] a, int[] b, int len) {
1.1918 + long sum = 0;
1.1919 +
1.1920 + while(--len >= 0) {
1.1921 + sum = (a[len] & LONG_MASK) -
1.1922 + (b[len] & LONG_MASK) + (sum >> 32);
1.1923 + a[len] = (int)sum;
1.1924 + }
1.1925 +
1.1926 + return (int)(sum >> 32);
1.1927 + }
1.1928 +
1.1929 + /**
1.1930 + * Multiply an array by one word k and add to result, return the carry
1.1931 + */
1.1932 + static int mulAdd(int[] out, int[] in, int offset, int len, int k) {
1.1933 + long kLong = k & LONG_MASK;
1.1934 + long carry = 0;
1.1935 +
1.1936 + offset = out.length-offset - 1;
1.1937 + for (int j=len-1; j >= 0; j--) {
1.1938 + long product = (in[j] & LONG_MASK) * kLong +
1.1939 + (out[offset] & LONG_MASK) + carry;
1.1940 + out[offset--] = (int)product;
1.1941 + carry = product >>> 32;
1.1942 + }
1.1943 + return (int)carry;
1.1944 + }
1.1945 +
1.1946 + /**
1.1947 + * Add one word to the number a mlen words into a. Return the resulting
1.1948 + * carry.
1.1949 + */
1.1950 + static int addOne(int[] a, int offset, int mlen, int carry) {
1.1951 + offset = a.length-1-mlen-offset;
1.1952 + long t = (a[offset] & LONG_MASK) + (carry & LONG_MASK);
1.1953 +
1.1954 + a[offset] = (int)t;
1.1955 + if ((t >>> 32) == 0)
1.1956 + return 0;
1.1957 + while (--mlen >= 0) {
1.1958 + if (--offset < 0) { // Carry out of number
1.1959 + return 1;
1.1960 + } else {
1.1961 + a[offset]++;
1.1962 + if (a[offset] != 0)
1.1963 + return 0;
1.1964 + }
1.1965 + }
1.1966 + return 1;
1.1967 + }
1.1968 +
1.1969 + /**
1.1970 + * Returns a BigInteger whose value is (this ** exponent) mod (2**p)
1.1971 + */
1.1972 + private BigInteger modPow2(BigInteger exponent, int p) {
1.1973 + /*
1.1974 + * Perform exponentiation using repeated squaring trick, chopping off
1.1975 + * high order bits as indicated by modulus.
1.1976 + */
1.1977 + BigInteger result = valueOf(1);
1.1978 + BigInteger baseToPow2 = this.mod2(p);
1.1979 + int expOffset = 0;
1.1980 +
1.1981 + int limit = exponent.bitLength();
1.1982 +
1.1983 + if (this.testBit(0))
1.1984 + limit = (p-1) < limit ? (p-1) : limit;
1.1985 +
1.1986 + while (expOffset < limit) {
1.1987 + if (exponent.testBit(expOffset))
1.1988 + result = result.multiply(baseToPow2).mod2(p);
1.1989 + expOffset++;
1.1990 + if (expOffset < limit)
1.1991 + baseToPow2 = baseToPow2.square().mod2(p);
1.1992 + }
1.1993 +
1.1994 + return result;
1.1995 + }
1.1996 +
1.1997 + /**
1.1998 + * Returns a BigInteger whose value is this mod(2**p).
1.1999 + * Assumes that this {@code BigInteger >= 0} and {@code p > 0}.
1.2000 + */
1.2001 + private BigInteger mod2(int p) {
1.2002 + if (bitLength() <= p)
1.2003 + return this;
1.2004 +
1.2005 + // Copy remaining ints of mag
1.2006 + int numInts = (p + 31) >>> 5;
1.2007 + int[] mag = new int[numInts];
1.2008 + for (int i=0; i<numInts; i++)
1.2009 + mag[i] = this.mag[i + (this.mag.length - numInts)];
1.2010 +
1.2011 + // Mask out any excess bits
1.2012 + int excessBits = (numInts << 5) - p;
1.2013 + mag[0] &= (1L << (32-excessBits)) - 1;
1.2014 +
1.2015 + return (mag[0]==0 ? new BigInteger(1, mag) : new BigInteger(mag, 1));
1.2016 + }
1.2017 +
1.2018 + /**
1.2019 + * Returns a BigInteger whose value is {@code (this}<sup>-1</sup> {@code mod m)}.
1.2020 + *
1.2021 + * @param m the modulus.
1.2022 + * @return {@code this}<sup>-1</sup> {@code mod m}.
1.2023 + * @throws ArithmeticException {@code m} ≤ 0, or this BigInteger
1.2024 + * has no multiplicative inverse mod m (that is, this BigInteger
1.2025 + * is not <i>relatively prime</i> to m).
1.2026 + */
1.2027 + public BigInteger modInverse(BigInteger m) {
1.2028 + if (m.signum != 1)
1.2029 + throw new ArithmeticException("BigInteger: modulus not positive");
1.2030 +
1.2031 + if (m.equals(ONE))
1.2032 + return ZERO;
1.2033 +
1.2034 + // Calculate (this mod m)
1.2035 + BigInteger modVal = this;
1.2036 + if (signum < 0 || (this.compareMagnitude(m) >= 0))
1.2037 + modVal = this.mod(m);
1.2038 +
1.2039 + if (modVal.equals(ONE))
1.2040 + return ONE;
1.2041 +
1.2042 + MutableBigInteger a = new MutableBigInteger(modVal);
1.2043 + MutableBigInteger b = new MutableBigInteger(m);
1.2044 +
1.2045 + MutableBigInteger result = a.mutableModInverse(b);
1.2046 + return result.toBigInteger(1);
1.2047 + }
1.2048 +
1.2049 + // Shift Operations
1.2050 +
1.2051 + /**
1.2052 + * Returns a BigInteger whose value is {@code (this << n)}.
1.2053 + * The shift distance, {@code n}, may be negative, in which case
1.2054 + * this method performs a right shift.
1.2055 + * (Computes <tt>floor(this * 2<sup>n</sup>)</tt>.)
1.2056 + *
1.2057 + * @param n shift distance, in bits.
1.2058 + * @return {@code this << n}
1.2059 + * @throws ArithmeticException if the shift distance is {@code
1.2060 + * Integer.MIN_VALUE}.
1.2061 + * @see #shiftRight
1.2062 + */
1.2063 + public BigInteger shiftLeft(int n) {
1.2064 + if (signum == 0)
1.2065 + return ZERO;
1.2066 + if (n==0)
1.2067 + return this;
1.2068 + if (n<0) {
1.2069 + if (n == Integer.MIN_VALUE) {
1.2070 + throw new ArithmeticException("Shift distance of Integer.MIN_VALUE not supported.");
1.2071 + } else {
1.2072 + return shiftRight(-n);
1.2073 + }
1.2074 + }
1.2075 +
1.2076 + int nInts = n >>> 5;
1.2077 + int nBits = n & 0x1f;
1.2078 + int magLen = mag.length;
1.2079 + int newMag[] = null;
1.2080 +
1.2081 + if (nBits == 0) {
1.2082 + newMag = new int[magLen + nInts];
1.2083 + for (int i=0; i<magLen; i++)
1.2084 + newMag[i] = mag[i];
1.2085 + } else {
1.2086 + int i = 0;
1.2087 + int nBits2 = 32 - nBits;
1.2088 + int highBits = mag[0] >>> nBits2;
1.2089 + if (highBits != 0) {
1.2090 + newMag = new int[magLen + nInts + 1];
1.2091 + newMag[i++] = highBits;
1.2092 + } else {
1.2093 + newMag = new int[magLen + nInts];
1.2094 + }
1.2095 + int j=0;
1.2096 + while (j < magLen-1)
1.2097 + newMag[i++] = mag[j++] << nBits | mag[j] >>> nBits2;
1.2098 + newMag[i] = mag[j] << nBits;
1.2099 + }
1.2100 +
1.2101 + return new BigInteger(newMag, signum);
1.2102 + }
1.2103 +
1.2104 + /**
1.2105 + * Returns a BigInteger whose value is {@code (this >> n)}. Sign
1.2106 + * extension is performed. The shift distance, {@code n}, may be
1.2107 + * negative, in which case this method performs a left shift.
1.2108 + * (Computes <tt>floor(this / 2<sup>n</sup>)</tt>.)
1.2109 + *
1.2110 + * @param n shift distance, in bits.
1.2111 + * @return {@code this >> n}
1.2112 + * @throws ArithmeticException if the shift distance is {@code
1.2113 + * Integer.MIN_VALUE}.
1.2114 + * @see #shiftLeft
1.2115 + */
1.2116 + public BigInteger shiftRight(int n) {
1.2117 + if (n==0)
1.2118 + return this;
1.2119 + if (n<0) {
1.2120 + if (n == Integer.MIN_VALUE) {
1.2121 + throw new ArithmeticException("Shift distance of Integer.MIN_VALUE not supported.");
1.2122 + } else {
1.2123 + return shiftLeft(-n);
1.2124 + }
1.2125 + }
1.2126 +
1.2127 + int nInts = n >>> 5;
1.2128 + int nBits = n & 0x1f;
1.2129 + int magLen = mag.length;
1.2130 + int newMag[] = null;
1.2131 +
1.2132 + // Special case: entire contents shifted off the end
1.2133 + if (nInts >= magLen)
1.2134 + return (signum >= 0 ? ZERO : negConst[1]);
1.2135 +
1.2136 + if (nBits == 0) {
1.2137 + int newMagLen = magLen - nInts;
1.2138 + newMag = new int[newMagLen];
1.2139 + for (int i=0; i<newMagLen; i++)
1.2140 + newMag[i] = mag[i];
1.2141 + } else {
1.2142 + int i = 0;
1.2143 + int highBits = mag[0] >>> nBits;
1.2144 + if (highBits != 0) {
1.2145 + newMag = new int[magLen - nInts];
1.2146 + newMag[i++] = highBits;
1.2147 + } else {
1.2148 + newMag = new int[magLen - nInts -1];
1.2149 + }
1.2150 +
1.2151 + int nBits2 = 32 - nBits;
1.2152 + int j=0;
1.2153 + while (j < magLen - nInts - 1)
1.2154 + newMag[i++] = (mag[j++] << nBits2) | (mag[j] >>> nBits);
1.2155 + }
1.2156 +
1.2157 + if (signum < 0) {
1.2158 + // Find out whether any one-bits were shifted off the end.
1.2159 + boolean onesLost = false;
1.2160 + for (int i=magLen-1, j=magLen-nInts; i>=j && !onesLost; i--)
1.2161 + onesLost = (mag[i] != 0);
1.2162 + if (!onesLost && nBits != 0)
1.2163 + onesLost = (mag[magLen - nInts - 1] << (32 - nBits) != 0);
1.2164 +
1.2165 + if (onesLost)
1.2166 + newMag = javaIncrement(newMag);
1.2167 + }
1.2168 +
1.2169 + return new BigInteger(newMag, signum);
1.2170 + }
1.2171 +
1.2172 + int[] javaIncrement(int[] val) {
1.2173 + int lastSum = 0;
1.2174 + for (int i=val.length-1; i >= 0 && lastSum == 0; i--)
1.2175 + lastSum = (val[i] += 1);
1.2176 + if (lastSum == 0) {
1.2177 + val = new int[val.length+1];
1.2178 + val[0] = 1;
1.2179 + }
1.2180 + return val;
1.2181 + }
1.2182 +
1.2183 + // Bitwise Operations
1.2184 +
1.2185 + /**
1.2186 + * Returns a BigInteger whose value is {@code (this & val)}. (This
1.2187 + * method returns a negative BigInteger if and only if this and val are
1.2188 + * both negative.)
1.2189 + *
1.2190 + * @param val value to be AND'ed with this BigInteger.
1.2191 + * @return {@code this & val}
1.2192 + */
1.2193 + public BigInteger and(BigInteger val) {
1.2194 + int[] result = new int[Math.max(intLength(), val.intLength())];
1.2195 + for (int i=0; i<result.length; i++)
1.2196 + result[i] = (getInt(result.length-i-1)
1.2197 + & val.getInt(result.length-i-1));
1.2198 +
1.2199 + return valueOf(result);
1.2200 + }
1.2201 +
1.2202 + /**
1.2203 + * Returns a BigInteger whose value is {@code (this | val)}. (This method
1.2204 + * returns a negative BigInteger if and only if either this or val is
1.2205 + * negative.)
1.2206 + *
1.2207 + * @param val value to be OR'ed with this BigInteger.
1.2208 + * @return {@code this | val}
1.2209 + */
1.2210 + public BigInteger or(BigInteger val) {
1.2211 + int[] result = new int[Math.max(intLength(), val.intLength())];
1.2212 + for (int i=0; i<result.length; i++)
1.2213 + result[i] = (getInt(result.length-i-1)
1.2214 + | val.getInt(result.length-i-1));
1.2215 +
1.2216 + return valueOf(result);
1.2217 + }
1.2218 +
1.2219 + /**
1.2220 + * Returns a BigInteger whose value is {@code (this ^ val)}. (This method
1.2221 + * returns a negative BigInteger if and only if exactly one of this and
1.2222 + * val are negative.)
1.2223 + *
1.2224 + * @param val value to be XOR'ed with this BigInteger.
1.2225 + * @return {@code this ^ val}
1.2226 + */
1.2227 + public BigInteger xor(BigInteger val) {
1.2228 + int[] result = new int[Math.max(intLength(), val.intLength())];
1.2229 + for (int i=0; i<result.length; i++)
1.2230 + result[i] = (getInt(result.length-i-1)
1.2231 + ^ val.getInt(result.length-i-1));
1.2232 +
1.2233 + return valueOf(result);
1.2234 + }
1.2235 +
1.2236 + /**
1.2237 + * Returns a BigInteger whose value is {@code (~this)}. (This method
1.2238 + * returns a negative value if and only if this BigInteger is
1.2239 + * non-negative.)
1.2240 + *
1.2241 + * @return {@code ~this}
1.2242 + */
1.2243 + public BigInteger not() {
1.2244 + int[] result = new int[intLength()];
1.2245 + for (int i=0; i<result.length; i++)
1.2246 + result[i] = ~getInt(result.length-i-1);
1.2247 +
1.2248 + return valueOf(result);
1.2249 + }
1.2250 +
1.2251 + /**
1.2252 + * Returns a BigInteger whose value is {@code (this & ~val)}. This
1.2253 + * method, which is equivalent to {@code and(val.not())}, is provided as
1.2254 + * a convenience for masking operations. (This method returns a negative
1.2255 + * BigInteger if and only if {@code this} is negative and {@code val} is
1.2256 + * positive.)
1.2257 + *
1.2258 + * @param val value to be complemented and AND'ed with this BigInteger.
1.2259 + * @return {@code this & ~val}
1.2260 + */
1.2261 + public BigInteger andNot(BigInteger val) {
1.2262 + int[] result = new int[Math.max(intLength(), val.intLength())];
1.2263 + for (int i=0; i<result.length; i++)
1.2264 + result[i] = (getInt(result.length-i-1)
1.2265 + & ~val.getInt(result.length-i-1));
1.2266 +
1.2267 + return valueOf(result);
1.2268 + }
1.2269 +
1.2270 +
1.2271 + // Single Bit Operations
1.2272 +
1.2273 + /**
1.2274 + * Returns {@code true} if and only if the designated bit is set.
1.2275 + * (Computes {@code ((this & (1<<n)) != 0)}.)
1.2276 + *
1.2277 + * @param n index of bit to test.
1.2278 + * @return {@code true} if and only if the designated bit is set.
1.2279 + * @throws ArithmeticException {@code n} is negative.
1.2280 + */
1.2281 + public boolean testBit(int n) {
1.2282 + if (n<0)
1.2283 + throw new ArithmeticException("Negative bit address");
1.2284 +
1.2285 + return (getInt(n >>> 5) & (1 << (n & 31))) != 0;
1.2286 + }
1.2287 +
1.2288 + /**
1.2289 + * Returns a BigInteger whose value is equivalent to this BigInteger
1.2290 + * with the designated bit set. (Computes {@code (this | (1<<n))}.)
1.2291 + *
1.2292 + * @param n index of bit to set.
1.2293 + * @return {@code this | (1<<n)}
1.2294 + * @throws ArithmeticException {@code n} is negative.
1.2295 + */
1.2296 + public BigInteger setBit(int n) {
1.2297 + if (n<0)
1.2298 + throw new ArithmeticException("Negative bit address");
1.2299 +
1.2300 + int intNum = n >>> 5;
1.2301 + int[] result = new int[Math.max(intLength(), intNum+2)];
1.2302 +
1.2303 + for (int i=0; i<result.length; i++)
1.2304 + result[result.length-i-1] = getInt(i);
1.2305 +
1.2306 + result[result.length-intNum-1] |= (1 << (n & 31));
1.2307 +
1.2308 + return valueOf(result);
1.2309 + }
1.2310 +
1.2311 + /**
1.2312 + * Returns a BigInteger whose value is equivalent to this BigInteger
1.2313 + * with the designated bit cleared.
1.2314 + * (Computes {@code (this & ~(1<<n))}.)
1.2315 + *
1.2316 + * @param n index of bit to clear.
1.2317 + * @return {@code this & ~(1<<n)}
1.2318 + * @throws ArithmeticException {@code n} is negative.
1.2319 + */
1.2320 + public BigInteger clearBit(int n) {
1.2321 + if (n<0)
1.2322 + throw new ArithmeticException("Negative bit address");
1.2323 +
1.2324 + int intNum = n >>> 5;
1.2325 + int[] result = new int[Math.max(intLength(), ((n + 1) >>> 5) + 1)];
1.2326 +
1.2327 + for (int i=0; i<result.length; i++)
1.2328 + result[result.length-i-1] = getInt(i);
1.2329 +
1.2330 + result[result.length-intNum-1] &= ~(1 << (n & 31));
1.2331 +
1.2332 + return valueOf(result);
1.2333 + }
1.2334 +
1.2335 + /**
1.2336 + * Returns a BigInteger whose value is equivalent to this BigInteger
1.2337 + * with the designated bit flipped.
1.2338 + * (Computes {@code (this ^ (1<<n))}.)
1.2339 + *
1.2340 + * @param n index of bit to flip.
1.2341 + * @return {@code this ^ (1<<n)}
1.2342 + * @throws ArithmeticException {@code n} is negative.
1.2343 + */
1.2344 + public BigInteger flipBit(int n) {
1.2345 + if (n<0)
1.2346 + throw new ArithmeticException("Negative bit address");
1.2347 +
1.2348 + int intNum = n >>> 5;
1.2349 + int[] result = new int[Math.max(intLength(), intNum+2)];
1.2350 +
1.2351 + for (int i=0; i<result.length; i++)
1.2352 + result[result.length-i-1] = getInt(i);
1.2353 +
1.2354 + result[result.length-intNum-1] ^= (1 << (n & 31));
1.2355 +
1.2356 + return valueOf(result);
1.2357 + }
1.2358 +
1.2359 + /**
1.2360 + * Returns the index of the rightmost (lowest-order) one bit in this
1.2361 + * BigInteger (the number of zero bits to the right of the rightmost
1.2362 + * one bit). Returns -1 if this BigInteger contains no one bits.
1.2363 + * (Computes {@code (this==0? -1 : log2(this & -this))}.)
1.2364 + *
1.2365 + * @return index of the rightmost one bit in this BigInteger.
1.2366 + */
1.2367 + public int getLowestSetBit() {
1.2368 + @SuppressWarnings("deprecation") int lsb = lowestSetBit - 2;
1.2369 + if (lsb == -2) { // lowestSetBit not initialized yet
1.2370 + lsb = 0;
1.2371 + if (signum == 0) {
1.2372 + lsb -= 1;
1.2373 + } else {
1.2374 + // Search for lowest order nonzero int
1.2375 + int i,b;
1.2376 + for (i=0; (b = getInt(i))==0; i++)
1.2377 + ;
1.2378 + lsb += (i << 5) + Integer.numberOfTrailingZeros(b);
1.2379 + }
1.2380 + lowestSetBit = lsb + 2;
1.2381 + }
1.2382 + return lsb;
1.2383 + }
1.2384 +
1.2385 +
1.2386 + // Miscellaneous Bit Operations
1.2387 +
1.2388 + /**
1.2389 + * Returns the number of bits in the minimal two's-complement
1.2390 + * representation of this BigInteger, <i>excluding</i> a sign bit.
1.2391 + * For positive BigIntegers, this is equivalent to the number of bits in
1.2392 + * the ordinary binary representation. (Computes
1.2393 + * {@code (ceil(log2(this < 0 ? -this : this+1)))}.)
1.2394 + *
1.2395 + * @return number of bits in the minimal two's-complement
1.2396 + * representation of this BigInteger, <i>excluding</i> a sign bit.
1.2397 + */
1.2398 + public int bitLength() {
1.2399 + @SuppressWarnings("deprecation") int n = bitLength - 1;
1.2400 + if (n == -1) { // bitLength not initialized yet
1.2401 + int[] m = mag;
1.2402 + int len = m.length;
1.2403 + if (len == 0) {
1.2404 + n = 0; // offset by one to initialize
1.2405 + } else {
1.2406 + // Calculate the bit length of the magnitude
1.2407 + int magBitLength = ((len - 1) << 5) + bitLengthForInt(mag[0]);
1.2408 + if (signum < 0) {
1.2409 + // Check if magnitude is a power of two
1.2410 + boolean pow2 = (Integer.bitCount(mag[0]) == 1);
1.2411 + for(int i=1; i< len && pow2; i++)
1.2412 + pow2 = (mag[i] == 0);
1.2413 +
1.2414 + n = (pow2 ? magBitLength -1 : magBitLength);
1.2415 + } else {
1.2416 + n = magBitLength;
1.2417 + }
1.2418 + }
1.2419 + bitLength = n + 1;
1.2420 + }
1.2421 + return n;
1.2422 + }
1.2423 +
1.2424 + /**
1.2425 + * Returns the number of bits in the two's complement representation
1.2426 + * of this BigInteger that differ from its sign bit. This method is
1.2427 + * useful when implementing bit-vector style sets atop BigIntegers.
1.2428 + *
1.2429 + * @return number of bits in the two's complement representation
1.2430 + * of this BigInteger that differ from its sign bit.
1.2431 + */
1.2432 + public int bitCount() {
1.2433 + @SuppressWarnings("deprecation") int bc = bitCount - 1;
1.2434 + if (bc == -1) { // bitCount not initialized yet
1.2435 + bc = 0; // offset by one to initialize
1.2436 + // Count the bits in the magnitude
1.2437 + for (int i=0; i<mag.length; i++)
1.2438 + bc += Integer.bitCount(mag[i]);
1.2439 + if (signum < 0) {
1.2440 + // Count the trailing zeros in the magnitude
1.2441 + int magTrailingZeroCount = 0, j;
1.2442 + for (j=mag.length-1; mag[j]==0; j--)
1.2443 + magTrailingZeroCount += 32;
1.2444 + magTrailingZeroCount += Integer.numberOfTrailingZeros(mag[j]);
1.2445 + bc += magTrailingZeroCount - 1;
1.2446 + }
1.2447 + bitCount = bc + 1;
1.2448 + }
1.2449 + return bc;
1.2450 + }
1.2451 +
1.2452 + // Primality Testing
1.2453 +
1.2454 + /**
1.2455 + * Returns {@code true} if this BigInteger is probably prime,
1.2456 + * {@code false} if it's definitely composite. If
1.2457 + * {@code certainty} is ≤ 0, {@code true} is
1.2458 + * returned.
1.2459 + *
1.2460 + * @param certainty a measure of the uncertainty that the caller is
1.2461 + * willing to tolerate: if the call returns {@code true}
1.2462 + * the probability that this BigInteger is prime exceeds
1.2463 + * (1 - 1/2<sup>{@code certainty}</sup>). The execution time of
1.2464 + * this method is proportional to the value of this parameter.
1.2465 + * @return {@code true} if this BigInteger is probably prime,
1.2466 + * {@code false} if it's definitely composite.
1.2467 + */
1.2468 + public boolean isProbablePrime(int certainty) {
1.2469 + if (certainty <= 0)
1.2470 + return true;
1.2471 + BigInteger w = this.abs();
1.2472 + if (w.equals(TWO))
1.2473 + return true;
1.2474 + if (!w.testBit(0) || w.equals(ONE))
1.2475 + return false;
1.2476 +
1.2477 + return w.primeToCertainty(certainty, null);
1.2478 + }
1.2479 +
1.2480 + // Comparison Operations
1.2481 +
1.2482 + /**
1.2483 + * Compares this BigInteger with the specified BigInteger. This
1.2484 + * method is provided in preference to individual methods for each
1.2485 + * of the six boolean comparison operators ({@literal <}, ==,
1.2486 + * {@literal >}, {@literal >=}, !=, {@literal <=}). The suggested
1.2487 + * idiom for performing these comparisons is: {@code
1.2488 + * (x.compareTo(y)} <<i>op</i>> {@code 0)}, where
1.2489 + * <<i>op</i>> is one of the six comparison operators.
1.2490 + *
1.2491 + * @param val BigInteger to which this BigInteger is to be compared.
1.2492 + * @return -1, 0 or 1 as this BigInteger is numerically less than, equal
1.2493 + * to, or greater than {@code val}.
1.2494 + */
1.2495 + public int compareTo(BigInteger val) {
1.2496 + if (signum == val.signum) {
1.2497 + switch (signum) {
1.2498 + case 1:
1.2499 + return compareMagnitude(val);
1.2500 + case -1:
1.2501 + return val.compareMagnitude(this);
1.2502 + default:
1.2503 + return 0;
1.2504 + }
1.2505 + }
1.2506 + return signum > val.signum ? 1 : -1;
1.2507 + }
1.2508 +
1.2509 + /**
1.2510 + * Compares the magnitude array of this BigInteger with the specified
1.2511 + * BigInteger's. This is the version of compareTo ignoring sign.
1.2512 + *
1.2513 + * @param val BigInteger whose magnitude array to be compared.
1.2514 + * @return -1, 0 or 1 as this magnitude array is less than, equal to or
1.2515 + * greater than the magnitude aray for the specified BigInteger's.
1.2516 + */
1.2517 + final int compareMagnitude(BigInteger val) {
1.2518 + int[] m1 = mag;
1.2519 + int len1 = m1.length;
1.2520 + int[] m2 = val.mag;
1.2521 + int len2 = m2.length;
1.2522 + if (len1 < len2)
1.2523 + return -1;
1.2524 + if (len1 > len2)
1.2525 + return 1;
1.2526 + for (int i = 0; i < len1; i++) {
1.2527 + int a = m1[i];
1.2528 + int b = m2[i];
1.2529 + if (a != b)
1.2530 + return ((a & LONG_MASK) < (b & LONG_MASK)) ? -1 : 1;
1.2531 + }
1.2532 + return 0;
1.2533 + }
1.2534 +
1.2535 + /**
1.2536 + * Compares this BigInteger with the specified Object for equality.
1.2537 + *
1.2538 + * @param x Object to which this BigInteger is to be compared.
1.2539 + * @return {@code true} if and only if the specified Object is a
1.2540 + * BigInteger whose value is numerically equal to this BigInteger.
1.2541 + */
1.2542 + public boolean equals(Object x) {
1.2543 + // This test is just an optimization, which may or may not help
1.2544 + if (x == this)
1.2545 + return true;
1.2546 +
1.2547 + if (!(x instanceof BigInteger))
1.2548 + return false;
1.2549 +
1.2550 + BigInteger xInt = (BigInteger) x;
1.2551 + if (xInt.signum != signum)
1.2552 + return false;
1.2553 +
1.2554 + int[] m = mag;
1.2555 + int len = m.length;
1.2556 + int[] xm = xInt.mag;
1.2557 + if (len != xm.length)
1.2558 + return false;
1.2559 +
1.2560 + for (int i = 0; i < len; i++)
1.2561 + if (xm[i] != m[i])
1.2562 + return false;
1.2563 +
1.2564 + return true;
1.2565 + }
1.2566 +
1.2567 + /**
1.2568 + * Returns the minimum of this BigInteger and {@code val}.
1.2569 + *
1.2570 + * @param val value with which the minimum is to be computed.
1.2571 + * @return the BigInteger whose value is the lesser of this BigInteger and
1.2572 + * {@code val}. If they are equal, either may be returned.
1.2573 + */
1.2574 + public BigInteger min(BigInteger val) {
1.2575 + return (compareTo(val)<0 ? this : val);
1.2576 + }
1.2577 +
1.2578 + /**
1.2579 + * Returns the maximum of this BigInteger and {@code val}.
1.2580 + *
1.2581 + * @param val value with which the maximum is to be computed.
1.2582 + * @return the BigInteger whose value is the greater of this and
1.2583 + * {@code val}. If they are equal, either may be returned.
1.2584 + */
1.2585 + public BigInteger max(BigInteger val) {
1.2586 + return (compareTo(val)>0 ? this : val);
1.2587 + }
1.2588 +
1.2589 +
1.2590 + // Hash Function
1.2591 +
1.2592 + /**
1.2593 + * Returns the hash code for this BigInteger.
1.2594 + *
1.2595 + * @return hash code for this BigInteger.
1.2596 + */
1.2597 + public int hashCode() {
1.2598 + int hashCode = 0;
1.2599 +
1.2600 + for (int i=0; i<mag.length; i++)
1.2601 + hashCode = (int)(31*hashCode + (mag[i] & LONG_MASK));
1.2602 +
1.2603 + return hashCode * signum;
1.2604 + }
1.2605 +
1.2606 + /**
1.2607 + * Returns the String representation of this BigInteger in the
1.2608 + * given radix. If the radix is outside the range from {@link
1.2609 + * Character#MIN_RADIX} to {@link Character#MAX_RADIX} inclusive,
1.2610 + * it will default to 10 (as is the case for
1.2611 + * {@code Integer.toString}). The digit-to-character mapping
1.2612 + * provided by {@code Character.forDigit} is used, and a minus
1.2613 + * sign is prepended if appropriate. (This representation is
1.2614 + * compatible with the {@link #BigInteger(String, int) (String,
1.2615 + * int)} constructor.)
1.2616 + *
1.2617 + * @param radix radix of the String representation.
1.2618 + * @return String representation of this BigInteger in the given radix.
1.2619 + * @see Integer#toString
1.2620 + * @see Character#forDigit
1.2621 + * @see #BigInteger(java.lang.String, int)
1.2622 + */
1.2623 + public String toString(int radix) {
1.2624 + if (signum == 0)
1.2625 + return "0";
1.2626 + if (radix < Character.MIN_RADIX || radix > Character.MAX_RADIX)
1.2627 + radix = 10;
1.2628 +
1.2629 + // Compute upper bound on number of digit groups and allocate space
1.2630 + int maxNumDigitGroups = (4*mag.length + 6)/7;
1.2631 + String digitGroup[] = new String[maxNumDigitGroups];
1.2632 +
1.2633 + // Translate number to string, a digit group at a time
1.2634 + BigInteger tmp = this.abs();
1.2635 + int numGroups = 0;
1.2636 + while (tmp.signum != 0) {
1.2637 + BigInteger d = longRadix[radix];
1.2638 +
1.2639 + MutableBigInteger q = new MutableBigInteger(),
1.2640 + a = new MutableBigInteger(tmp.mag),
1.2641 + b = new MutableBigInteger(d.mag);
1.2642 + MutableBigInteger r = a.divide(b, q);
1.2643 + BigInteger q2 = q.toBigInteger(tmp.signum * d.signum);
1.2644 + BigInteger r2 = r.toBigInteger(tmp.signum * d.signum);
1.2645 +
1.2646 + digitGroup[numGroups++] = Long.toString(r2.longValue(), radix);
1.2647 + tmp = q2;
1.2648 + }
1.2649 +
1.2650 + // Put sign (if any) and first digit group into result buffer
1.2651 + StringBuilder buf = new StringBuilder(numGroups*digitsPerLong[radix]+1);
1.2652 + if (signum<0)
1.2653 + buf.append('-');
1.2654 + buf.append(digitGroup[numGroups-1]);
1.2655 +
1.2656 + // Append remaining digit groups padded with leading zeros
1.2657 + for (int i=numGroups-2; i>=0; i--) {
1.2658 + // Prepend (any) leading zeros for this digit group
1.2659 + int numLeadingZeros = digitsPerLong[radix]-digitGroup[i].length();
1.2660 + if (numLeadingZeros != 0)
1.2661 + buf.append(zeros[numLeadingZeros]);
1.2662 + buf.append(digitGroup[i]);
1.2663 + }
1.2664 + return buf.toString();
1.2665 + }
1.2666 +
1.2667 + /* zero[i] is a string of i consecutive zeros. */
1.2668 + private static String zeros[] = new String[64];
1.2669 + static {
1.2670 + zeros[63] =
1.2671 + "000000000000000000000000000000000000000000000000000000000000000";
1.2672 + for (int i=0; i<63; i++)
1.2673 + zeros[i] = zeros[63].substring(0, i);
1.2674 + }
1.2675 +
1.2676 + /**
1.2677 + * Returns the decimal String representation of this BigInteger.
1.2678 + * The digit-to-character mapping provided by
1.2679 + * {@code Character.forDigit} is used, and a minus sign is
1.2680 + * prepended if appropriate. (This representation is compatible
1.2681 + * with the {@link #BigInteger(String) (String)} constructor, and
1.2682 + * allows for String concatenation with Java's + operator.)
1.2683 + *
1.2684 + * @return decimal String representation of this BigInteger.
1.2685 + * @see Character#forDigit
1.2686 + * @see #BigInteger(java.lang.String)
1.2687 + */
1.2688 + public String toString() {
1.2689 + return toString(10);
1.2690 + }
1.2691 +
1.2692 + /**
1.2693 + * Returns a byte array containing the two's-complement
1.2694 + * representation of this BigInteger. The byte array will be in
1.2695 + * <i>big-endian</i> byte-order: the most significant byte is in
1.2696 + * the zeroth element. The array will contain the minimum number
1.2697 + * of bytes required to represent this BigInteger, including at
1.2698 + * least one sign bit, which is {@code (ceil((this.bitLength() +
1.2699 + * 1)/8))}. (This representation is compatible with the
1.2700 + * {@link #BigInteger(byte[]) (byte[])} constructor.)
1.2701 + *
1.2702 + * @return a byte array containing the two's-complement representation of
1.2703 + * this BigInteger.
1.2704 + * @see #BigInteger(byte[])
1.2705 + */
1.2706 + public byte[] toByteArray() {
1.2707 + int byteLen = bitLength()/8 + 1;
1.2708 + byte[] byteArray = new byte[byteLen];
1.2709 +
1.2710 + for (int i=byteLen-1, bytesCopied=4, nextInt=0, intIndex=0; i>=0; i--) {
1.2711 + if (bytesCopied == 4) {
1.2712 + nextInt = getInt(intIndex++);
1.2713 + bytesCopied = 1;
1.2714 + } else {
1.2715 + nextInt >>>= 8;
1.2716 + bytesCopied++;
1.2717 + }
1.2718 + byteArray[i] = (byte)nextInt;
1.2719 + }
1.2720 + return byteArray;
1.2721 + }
1.2722 +
1.2723 + /**
1.2724 + * Converts this BigInteger to an {@code int}. This
1.2725 + * conversion is analogous to a
1.2726 + * <i>narrowing primitive conversion</i> from {@code long} to
1.2727 + * {@code int} as defined in section 5.1.3 of
1.2728 + * <cite>The Java™ Language Specification</cite>:
1.2729 + * if this BigInteger is too big to fit in an
1.2730 + * {@code int}, only the low-order 32 bits are returned.
1.2731 + * Note that this conversion can lose information about the
1.2732 + * overall magnitude of the BigInteger value as well as return a
1.2733 + * result with the opposite sign.
1.2734 + *
1.2735 + * @return this BigInteger converted to an {@code int}.
1.2736 + */
1.2737 + public int intValue() {
1.2738 + int result = 0;
1.2739 + result = getInt(0);
1.2740 + return result;
1.2741 + }
1.2742 +
1.2743 + /**
1.2744 + * Converts this BigInteger to a {@code long}. This
1.2745 + * conversion is analogous to a
1.2746 + * <i>narrowing primitive conversion</i> from {@code long} to
1.2747 + * {@code int} as defined in section 5.1.3 of
1.2748 + * <cite>The Java™ Language Specification</cite>:
1.2749 + * if this BigInteger is too big to fit in a
1.2750 + * {@code long}, only the low-order 64 bits are returned.
1.2751 + * Note that this conversion can lose information about the
1.2752 + * overall magnitude of the BigInteger value as well as return a
1.2753 + * result with the opposite sign.
1.2754 + *
1.2755 + * @return this BigInteger converted to a {@code long}.
1.2756 + */
1.2757 + public long longValue() {
1.2758 + long result = 0;
1.2759 +
1.2760 + for (int i=1; i>=0; i--)
1.2761 + result = (result << 32) + (getInt(i) & LONG_MASK);
1.2762 + return result;
1.2763 + }
1.2764 +
1.2765 + /**
1.2766 + * Converts this BigInteger to a {@code float}. This
1.2767 + * conversion is similar to the
1.2768 + * <i>narrowing primitive conversion</i> from {@code double} to
1.2769 + * {@code float} as defined in section 5.1.3 of
1.2770 + * <cite>The Java™ Language Specification</cite>:
1.2771 + * if this BigInteger has too great a magnitude
1.2772 + * to represent as a {@code float}, it will be converted to
1.2773 + * {@link Float#NEGATIVE_INFINITY} or {@link
1.2774 + * Float#POSITIVE_INFINITY} as appropriate. Note that even when
1.2775 + * the return value is finite, this conversion can lose
1.2776 + * information about the precision of the BigInteger value.
1.2777 + *
1.2778 + * @return this BigInteger converted to a {@code float}.
1.2779 + */
1.2780 + public float floatValue() {
1.2781 + // Somewhat inefficient, but guaranteed to work.
1.2782 + return Float.parseFloat(this.toString());
1.2783 + }
1.2784 +
1.2785 + /**
1.2786 + * Converts this BigInteger to a {@code double}. This
1.2787 + * conversion is similar to the
1.2788 + * <i>narrowing primitive conversion</i> from {@code double} to
1.2789 + * {@code float} as defined in section 5.1.3 of
1.2790 + * <cite>The Java™ Language Specification</cite>:
1.2791 + * if this BigInteger has too great a magnitude
1.2792 + * to represent as a {@code double}, it will be converted to
1.2793 + * {@link Double#NEGATIVE_INFINITY} or {@link
1.2794 + * Double#POSITIVE_INFINITY} as appropriate. Note that even when
1.2795 + * the return value is finite, this conversion can lose
1.2796 + * information about the precision of the BigInteger value.
1.2797 + *
1.2798 + * @return this BigInteger converted to a {@code double}.
1.2799 + */
1.2800 + public double doubleValue() {
1.2801 + // Somewhat inefficient, but guaranteed to work.
1.2802 + return Double.parseDouble(this.toString());
1.2803 + }
1.2804 +
1.2805 + /**
1.2806 + * Returns a copy of the input array stripped of any leading zero bytes.
1.2807 + */
1.2808 + private static int[] stripLeadingZeroInts(int val[]) {
1.2809 + int vlen = val.length;
1.2810 + int keep;
1.2811 +
1.2812 + // Find first nonzero byte
1.2813 + for (keep = 0; keep < vlen && val[keep] == 0; keep++)
1.2814 + ;
1.2815 + return java.util.Arrays.copyOfRange(val, keep, vlen);
1.2816 + }
1.2817 +
1.2818 + /**
1.2819 + * Returns the input array stripped of any leading zero bytes.
1.2820 + * Since the source is trusted the copying may be skipped.
1.2821 + */
1.2822 + private static int[] trustedStripLeadingZeroInts(int val[]) {
1.2823 + int vlen = val.length;
1.2824 + int keep;
1.2825 +
1.2826 + // Find first nonzero byte
1.2827 + for (keep = 0; keep < vlen && val[keep] == 0; keep++)
1.2828 + ;
1.2829 + return keep == 0 ? val : java.util.Arrays.copyOfRange(val, keep, vlen);
1.2830 + }
1.2831 +
1.2832 + /**
1.2833 + * Returns a copy of the input array stripped of any leading zero bytes.
1.2834 + */
1.2835 + private static int[] stripLeadingZeroBytes(byte a[]) {
1.2836 + int byteLength = a.length;
1.2837 + int keep;
1.2838 +
1.2839 + // Find first nonzero byte
1.2840 + for (keep = 0; keep < byteLength && a[keep]==0; keep++)
1.2841 + ;
1.2842 +
1.2843 + // Allocate new array and copy relevant part of input array
1.2844 + int intLength = ((byteLength - keep) + 3) >>> 2;
1.2845 + int[] result = new int[intLength];
1.2846 + int b = byteLength - 1;
1.2847 + for (int i = intLength-1; i >= 0; i--) {
1.2848 + result[i] = a[b--] & 0xff;
1.2849 + int bytesRemaining = b - keep + 1;
1.2850 + int bytesToTransfer = Math.min(3, bytesRemaining);
1.2851 + for (int j=8; j <= (bytesToTransfer << 3); j += 8)
1.2852 + result[i] |= ((a[b--] & 0xff) << j);
1.2853 + }
1.2854 + return result;
1.2855 + }
1.2856 +
1.2857 + /**
1.2858 + * Takes an array a representing a negative 2's-complement number and
1.2859 + * returns the minimal (no leading zero bytes) unsigned whose value is -a.
1.2860 + */
1.2861 + private static int[] makePositive(byte a[]) {
1.2862 + int keep, k;
1.2863 + int byteLength = a.length;
1.2864 +
1.2865 + // Find first non-sign (0xff) byte of input
1.2866 + for (keep=0; keep<byteLength && a[keep]==-1; keep++)
1.2867 + ;
1.2868 +
1.2869 +
1.2870 + /* Allocate output array. If all non-sign bytes are 0x00, we must
1.2871 + * allocate space for one extra output byte. */
1.2872 + for (k=keep; k<byteLength && a[k]==0; k++)
1.2873 + ;
1.2874 +
1.2875 + int extraByte = (k==byteLength) ? 1 : 0;
1.2876 + int intLength = ((byteLength - keep + extraByte) + 3)/4;
1.2877 + int result[] = new int[intLength];
1.2878 +
1.2879 + /* Copy one's complement of input into output, leaving extra
1.2880 + * byte (if it exists) == 0x00 */
1.2881 + int b = byteLength - 1;
1.2882 + for (int i = intLength-1; i >= 0; i--) {
1.2883 + result[i] = a[b--] & 0xff;
1.2884 + int numBytesToTransfer = Math.min(3, b-keep+1);
1.2885 + if (numBytesToTransfer < 0)
1.2886 + numBytesToTransfer = 0;
1.2887 + for (int j=8; j <= 8*numBytesToTransfer; j += 8)
1.2888 + result[i] |= ((a[b--] & 0xff) << j);
1.2889 +
1.2890 + // Mask indicates which bits must be complemented
1.2891 + int mask = -1 >>> (8*(3-numBytesToTransfer));
1.2892 + result[i] = ~result[i] & mask;
1.2893 + }
1.2894 +
1.2895 + // Add one to one's complement to generate two's complement
1.2896 + for (int i=result.length-1; i>=0; i--) {
1.2897 + result[i] = (int)((result[i] & LONG_MASK) + 1);
1.2898 + if (result[i] != 0)
1.2899 + break;
1.2900 + }
1.2901 +
1.2902 + return result;
1.2903 + }
1.2904 +
1.2905 + /**
1.2906 + * Takes an array a representing a negative 2's-complement number and
1.2907 + * returns the minimal (no leading zero ints) unsigned whose value is -a.
1.2908 + */
1.2909 + private static int[] makePositive(int a[]) {
1.2910 + int keep, j;
1.2911 +
1.2912 + // Find first non-sign (0xffffffff) int of input
1.2913 + for (keep=0; keep<a.length && a[keep]==-1; keep++)
1.2914 + ;
1.2915 +
1.2916 + /* Allocate output array. If all non-sign ints are 0x00, we must
1.2917 + * allocate space for one extra output int. */
1.2918 + for (j=keep; j<a.length && a[j]==0; j++)
1.2919 + ;
1.2920 + int extraInt = (j==a.length ? 1 : 0);
1.2921 + int result[] = new int[a.length - keep + extraInt];
1.2922 +
1.2923 + /* Copy one's complement of input into output, leaving extra
1.2924 + * int (if it exists) == 0x00 */
1.2925 + for (int i = keep; i<a.length; i++)
1.2926 + result[i - keep + extraInt] = ~a[i];
1.2927 +
1.2928 + // Add one to one's complement to generate two's complement
1.2929 + for (int i=result.length-1; ++result[i]==0; i--)
1.2930 + ;
1.2931 +
1.2932 + return result;
1.2933 + }
1.2934 +
1.2935 + /*
1.2936 + * The following two arrays are used for fast String conversions. Both
1.2937 + * are indexed by radix. The first is the number of digits of the given
1.2938 + * radix that can fit in a Java long without "going negative", i.e., the
1.2939 + * highest integer n such that radix**n < 2**63. The second is the
1.2940 + * "long radix" that tears each number into "long digits", each of which
1.2941 + * consists of the number of digits in the corresponding element in
1.2942 + * digitsPerLong (longRadix[i] = i**digitPerLong[i]). Both arrays have
1.2943 + * nonsense values in their 0 and 1 elements, as radixes 0 and 1 are not
1.2944 + * used.
1.2945 + */
1.2946 + private static int digitsPerLong[] = {0, 0,
1.2947 + 62, 39, 31, 27, 24, 22, 20, 19, 18, 18, 17, 17, 16, 16, 15, 15, 15, 14,
1.2948 + 14, 14, 14, 13, 13, 13, 13, 13, 13, 12, 12, 12, 12, 12, 12, 12, 12};
1.2949 +
1.2950 + private static BigInteger longRadix[] = {null, null,
1.2951 + valueOf(0x4000000000000000L), valueOf(0x383d9170b85ff80bL),
1.2952 + valueOf(0x4000000000000000L), valueOf(0x6765c793fa10079dL),
1.2953 + valueOf(0x41c21cb8e1000000L), valueOf(0x3642798750226111L),
1.2954 + valueOf(0x1000000000000000L), valueOf(0x12bf307ae81ffd59L),
1.2955 + valueOf( 0xde0b6b3a7640000L), valueOf(0x4d28cb56c33fa539L),
1.2956 + valueOf(0x1eca170c00000000L), valueOf(0x780c7372621bd74dL),
1.2957 + valueOf(0x1e39a5057d810000L), valueOf(0x5b27ac993df97701L),
1.2958 + valueOf(0x1000000000000000L), valueOf(0x27b95e997e21d9f1L),
1.2959 + valueOf(0x5da0e1e53c5c8000L), valueOf( 0xb16a458ef403f19L),
1.2960 + valueOf(0x16bcc41e90000000L), valueOf(0x2d04b7fdd9c0ef49L),
1.2961 + valueOf(0x5658597bcaa24000L), valueOf( 0x6feb266931a75b7L),
1.2962 + valueOf( 0xc29e98000000000L), valueOf(0x14adf4b7320334b9L),
1.2963 + valueOf(0x226ed36478bfa000L), valueOf(0x383d9170b85ff80bL),
1.2964 + valueOf(0x5a3c23e39c000000L), valueOf( 0x4e900abb53e6b71L),
1.2965 + valueOf( 0x7600ec618141000L), valueOf( 0xaee5720ee830681L),
1.2966 + valueOf(0x1000000000000000L), valueOf(0x172588ad4f5f0981L),
1.2967 + valueOf(0x211e44f7d02c1000L), valueOf(0x2ee56725f06e5c71L),
1.2968 + valueOf(0x41c21cb8e1000000L)};
1.2969 +
1.2970 + /*
1.2971 + * These two arrays are the integer analogue of above.
1.2972 + */
1.2973 + private static int digitsPerInt[] = {0, 0, 30, 19, 15, 13, 11,
1.2974 + 11, 10, 9, 9, 8, 8, 8, 8, 7, 7, 7, 7, 7, 7, 7, 6, 6, 6, 6,
1.2975 + 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 5};
1.2976 +
1.2977 + private static int intRadix[] = {0, 0,
1.2978 + 0x40000000, 0x4546b3db, 0x40000000, 0x48c27395, 0x159fd800,
1.2979 + 0x75db9c97, 0x40000000, 0x17179149, 0x3b9aca00, 0xcc6db61,
1.2980 + 0x19a10000, 0x309f1021, 0x57f6c100, 0xa2f1b6f, 0x10000000,
1.2981 + 0x18754571, 0x247dbc80, 0x3547667b, 0x4c4b4000, 0x6b5a6e1d,
1.2982 + 0x6c20a40, 0x8d2d931, 0xb640000, 0xe8d4a51, 0x1269ae40,
1.2983 + 0x17179149, 0x1cb91000, 0x23744899, 0x2b73a840, 0x34e63b41,
1.2984 + 0x40000000, 0x4cfa3cc1, 0x5c13d840, 0x6d91b519, 0x39aa400
1.2985 + };
1.2986 +
1.2987 + /**
1.2988 + * These routines provide access to the two's complement representation
1.2989 + * of BigIntegers.
1.2990 + */
1.2991 +
1.2992 + /**
1.2993 + * Returns the length of the two's complement representation in ints,
1.2994 + * including space for at least one sign bit.
1.2995 + */
1.2996 + private int intLength() {
1.2997 + return (bitLength() >>> 5) + 1;
1.2998 + }
1.2999 +
1.3000 + /* Returns sign bit */
1.3001 + private int signBit() {
1.3002 + return signum < 0 ? 1 : 0;
1.3003 + }
1.3004 +
1.3005 + /* Returns an int of sign bits */
1.3006 + private int signInt() {
1.3007 + return signum < 0 ? -1 : 0;
1.3008 + }
1.3009 +
1.3010 + /**
1.3011 + * Returns the specified int of the little-endian two's complement
1.3012 + * representation (int 0 is the least significant). The int number can
1.3013 + * be arbitrarily high (values are logically preceded by infinitely many
1.3014 + * sign ints).
1.3015 + */
1.3016 + private int getInt(int n) {
1.3017 + if (n < 0)
1.3018 + return 0;
1.3019 + if (n >= mag.length)
1.3020 + return signInt();
1.3021 +
1.3022 + int magInt = mag[mag.length-n-1];
1.3023 +
1.3024 + return (signum >= 0 ? magInt :
1.3025 + (n <= firstNonzeroIntNum() ? -magInt : ~magInt));
1.3026 + }
1.3027 +
1.3028 + /**
1.3029 + * Returns the index of the int that contains the first nonzero int in the
1.3030 + * little-endian binary representation of the magnitude (int 0 is the
1.3031 + * least significant). If the magnitude is zero, return value is undefined.
1.3032 + */
1.3033 + private int firstNonzeroIntNum() {
1.3034 + int fn = firstNonzeroIntNum - 2;
1.3035 + if (fn == -2) { // firstNonzeroIntNum not initialized yet
1.3036 + fn = 0;
1.3037 +
1.3038 + // Search for the first nonzero int
1.3039 + int i;
1.3040 + int mlen = mag.length;
1.3041 + for (i = mlen - 1; i >= 0 && mag[i] == 0; i--)
1.3042 + ;
1.3043 + fn = mlen - i - 1;
1.3044 + firstNonzeroIntNum = fn + 2; // offset by two to initialize
1.3045 + }
1.3046 + return fn;
1.3047 + }
1.3048 +
1.3049 + /** use serialVersionUID from JDK 1.1. for interoperability */
1.3050 + private static final long serialVersionUID = -8287574255936472291L;
1.3051 +
1.3052 + /**
1.3053 + * Serializable fields for BigInteger.
1.3054 + *
1.3055 + * @serialField signum int
1.3056 + * signum of this BigInteger.
1.3057 + * @serialField magnitude int[]
1.3058 + * magnitude array of this BigInteger.
1.3059 + * @serialField bitCount int
1.3060 + * number of bits in this BigInteger
1.3061 + * @serialField bitLength int
1.3062 + * the number of bits in the minimal two's-complement
1.3063 + * representation of this BigInteger
1.3064 + * @serialField lowestSetBit int
1.3065 + * lowest set bit in the twos complement representation
1.3066 + */
1.3067 + private static final ObjectStreamField[] serialPersistentFields = {
1.3068 + new ObjectStreamField("signum", Integer.TYPE),
1.3069 + new ObjectStreamField("magnitude", byte[].class),
1.3070 + new ObjectStreamField("bitCount", Integer.TYPE),
1.3071 + new ObjectStreamField("bitLength", Integer.TYPE),
1.3072 + new ObjectStreamField("firstNonzeroByteNum", Integer.TYPE),
1.3073 + new ObjectStreamField("lowestSetBit", Integer.TYPE)
1.3074 + };
1.3075 +
1.3076 + /**
1.3077 + * Reconstitute the {@code BigInteger} instance from a stream (that is,
1.3078 + * deserialize it). The magnitude is read in as an array of bytes
1.3079 + * for historical reasons, but it is converted to an array of ints
1.3080 + * and the byte array is discarded.
1.3081 + * Note:
1.3082 + * The current convention is to initialize the cache fields, bitCount,
1.3083 + * bitLength and lowestSetBit, to 0 rather than some other marker value.
1.3084 + * Therefore, no explicit action to set these fields needs to be taken in
1.3085 + * readObject because those fields already have a 0 value be default since
1.3086 + * defaultReadObject is not being used.
1.3087 + */
1.3088 + private void readObject(java.io.ObjectInputStream s)
1.3089 + throws java.io.IOException, ClassNotFoundException {
1.3090 + /*
1.3091 + * In order to maintain compatibility with previous serialized forms,
1.3092 + * the magnitude of a BigInteger is serialized as an array of bytes.
1.3093 + * The magnitude field is used as a temporary store for the byte array
1.3094 + * that is deserialized. The cached computation fields should be
1.3095 + * transient but are serialized for compatibility reasons.
1.3096 + */
1.3097 +
1.3098 + // prepare to read the alternate persistent fields
1.3099 + ObjectInputStream.GetField fields = s.readFields();
1.3100 +
1.3101 + // Read the alternate persistent fields that we care about
1.3102 + int sign = fields.get("signum", -2);
1.3103 + byte[] magnitude = (byte[])fields.get("magnitude", null);
1.3104 +
1.3105 + // Validate signum
1.3106 + if (sign < -1 || sign > 1) {
1.3107 + String message = "BigInteger: Invalid signum value";
1.3108 + if (fields.defaulted("signum"))
1.3109 + message = "BigInteger: Signum not present in stream";
1.3110 + throw new java.io.StreamCorruptedException(message);
1.3111 + }
1.3112 + if ((magnitude.length == 0) != (sign == 0)) {
1.3113 + String message = "BigInteger: signum-magnitude mismatch";
1.3114 + if (fields.defaulted("magnitude"))
1.3115 + message = "BigInteger: Magnitude not present in stream";
1.3116 + throw new java.io.StreamCorruptedException(message);
1.3117 + }
1.3118 +
1.3119 + // Commit final fields via Unsafe
1.3120 + unsafe.putIntVolatile(this, signumOffset, sign);
1.3121 +
1.3122 + // Calculate mag field from magnitude and discard magnitude
1.3123 + unsafe.putObjectVolatile(this, magOffset,
1.3124 + stripLeadingZeroBytes(magnitude));
1.3125 + }
1.3126 +
1.3127 + // Support for resetting final fields while deserializing
1.3128 + private static final sun.misc.Unsafe unsafe = sun.misc.Unsafe.getUnsafe();
1.3129 + private static final long signumOffset;
1.3130 + private static final long magOffset;
1.3131 + static {
1.3132 + try {
1.3133 + signumOffset = unsafe.objectFieldOffset
1.3134 + (BigInteger.class.getDeclaredField("signum"));
1.3135 + magOffset = unsafe.objectFieldOffset
1.3136 + (BigInteger.class.getDeclaredField("mag"));
1.3137 + } catch (Exception ex) {
1.3138 + throw new Error(ex);
1.3139 + }
1.3140 + }
1.3141 +
1.3142 + /**
1.3143 + * Save the {@code BigInteger} instance to a stream.
1.3144 + * The magnitude of a BigInteger is serialized as a byte array for
1.3145 + * historical reasons.
1.3146 + *
1.3147 + * @serialData two necessary fields are written as well as obsolete
1.3148 + * fields for compatibility with older versions.
1.3149 + */
1.3150 + private void writeObject(ObjectOutputStream s) throws IOException {
1.3151 + // set the values of the Serializable fields
1.3152 + ObjectOutputStream.PutField fields = s.putFields();
1.3153 + fields.put("signum", signum);
1.3154 + fields.put("magnitude", magSerializedForm());
1.3155 + // The values written for cached fields are compatible with older
1.3156 + // versions, but are ignored in readObject so don't otherwise matter.
1.3157 + fields.put("bitCount", -1);
1.3158 + fields.put("bitLength", -1);
1.3159 + fields.put("lowestSetBit", -2);
1.3160 + fields.put("firstNonzeroByteNum", -2);
1.3161 +
1.3162 + // save them
1.3163 + s.writeFields();
1.3164 +}
1.3165 +
1.3166 + /**
1.3167 + * Returns the mag array as an array of bytes.
1.3168 + */
1.3169 + private byte[] magSerializedForm() {
1.3170 + int len = mag.length;
1.3171 +
1.3172 + int bitLen = (len == 0 ? 0 : ((len - 1) << 5) + bitLengthForInt(mag[0]));
1.3173 + int byteLen = (bitLen + 7) >>> 3;
1.3174 + byte[] result = new byte[byteLen];
1.3175 +
1.3176 + for (int i = byteLen - 1, bytesCopied = 4, intIndex = len - 1, nextInt = 0;
1.3177 + i>=0; i--) {
1.3178 + if (bytesCopied == 4) {
1.3179 + nextInt = mag[intIndex--];
1.3180 + bytesCopied = 1;
1.3181 + } else {
1.3182 + nextInt >>>= 8;
1.3183 + bytesCopied++;
1.3184 + }
1.3185 + result[i] = (byte)nextInt;
1.3186 + }
1.3187 + return result;
1.3188 + }
1.3189 +}