001/*
002 * Copyright (C) 2011 The Guava Authors
003 *
004 * Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except
005 * in compliance with the License. You may obtain a copy of the License at
006 *
007 * http://www.apache.org/licenses/LICENSE-2.0
008 *
009 * Unless required by applicable law or agreed to in writing, software distributed under the License
010 * is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express
011 * or implied. See the License for the specific language governing permissions and limitations under
012 * the License.
013 */
014
015package com.google.common.math;
016
017import static com.google.common.base.Preconditions.checkArgument;
018import static com.google.common.base.Preconditions.checkNotNull;
019import static com.google.common.math.MathPreconditions.checkNoOverflow;
020import static com.google.common.math.MathPreconditions.checkNonNegative;
021import static com.google.common.math.MathPreconditions.checkPositive;
022import static com.google.common.math.MathPreconditions.checkRoundingUnnecessary;
023import static java.lang.Math.abs;
024import static java.lang.Math.min;
025import static java.math.RoundingMode.HALF_EVEN;
026import static java.math.RoundingMode.HALF_UP;
027
028import com.google.common.annotations.GwtCompatible;
029import com.google.common.annotations.GwtIncompatible;
030import com.google.common.annotations.VisibleForTesting;
031import com.google.common.primitives.Longs;
032import com.google.common.primitives.UnsignedLongs;
033import java.math.BigInteger;
034import java.math.RoundingMode;
035
036/**
037 * A class for arithmetic on values of type {@code long}. Where possible, methods are defined and
038 * named analogously to their {@code BigInteger} counterparts.
039 *
040 * <p>The implementations of many methods in this class are based on material from Henry S. Warren,
041 * Jr.'s <i>Hacker's Delight</i>, (Addison Wesley, 2002).
042 *
043 * <p>Similar functionality for {@code int} and for {@link BigInteger} can be found in {@link
044 * IntMath} and {@link BigIntegerMath} respectively. For other common operations on {@code long}
045 * values, see {@link com.google.common.primitives.Longs}.
046 *
047 * @author Louis Wasserman
048 * @since 11.0
049 */
050@GwtCompatible(emulated = true)
051@ElementTypesAreNonnullByDefault
052public final class LongMath {
053  // NOTE: Whenever both tests are cheap and functional, it's faster to use &, | instead of &&, ||
054
055  @VisibleForTesting static final long MAX_SIGNED_POWER_OF_TWO = 1L << (Long.SIZE - 2);
056
057  /**
058   * Returns the smallest power of two greater than or equal to {@code x}. This is equivalent to
059   * {@code checkedPow(2, log2(x, CEILING))}.
060   *
061   * @throws IllegalArgumentException if {@code x <= 0}
062   * @throws ArithmeticException of the next-higher power of two is not representable as a {@code
063   *     long}, i.e. when {@code x > 2^62}
064   * @since 20.0
065   */
066  public static long ceilingPowerOfTwo(long x) {
067    checkPositive("x", x);
068    if (x > MAX_SIGNED_POWER_OF_TWO) {
069      throw new ArithmeticException("ceilingPowerOfTwo(" + x + ") is not representable as a long");
070    }
071    return 1L << -Long.numberOfLeadingZeros(x - 1);
072  }
073
074  /**
075   * Returns the largest power of two less than or equal to {@code x}. This is equivalent to {@code
076   * checkedPow(2, log2(x, FLOOR))}.
077   *
078   * @throws IllegalArgumentException if {@code x <= 0}
079   * @since 20.0
080   */
081  public static long floorPowerOfTwo(long x) {
082    checkPositive("x", x);
083
084    // Long.highestOneBit was buggy on GWT.  We've fixed it, but I'm not certain when the fix will
085    // be released.
086    return 1L << ((Long.SIZE - 1) - Long.numberOfLeadingZeros(x));
087  }
088
089  /**
090   * Returns {@code true} if {@code x} represents a power of two.
091   *
092   * <p>This differs from {@code Long.bitCount(x) == 1}, because {@code
093   * Long.bitCount(Long.MIN_VALUE) == 1}, but {@link Long#MIN_VALUE} is not a power of two.
094   */
095  @SuppressWarnings("ShortCircuitBoolean")
096  public static boolean isPowerOfTwo(long x) {
097    return x > 0 & (x & (x - 1)) == 0;
098  }
099
100  /**
101   * Returns 1 if {@code x < y} as unsigned longs, and 0 otherwise. Assumes that x - y fits into a
102   * signed long. The implementation is branch-free, and benchmarks suggest it is measurably faster
103   * than the straightforward ternary expression.
104   */
105  @VisibleForTesting
106  static int lessThanBranchFree(long x, long y) {
107    // Returns the sign bit of x - y.
108    return (int) (~~(x - y) >>> (Long.SIZE - 1));
109  }
110
111  /**
112   * Returns the base-2 logarithm of {@code x}, rounded according to the specified rounding mode.
113   *
114   * @throws IllegalArgumentException if {@code x <= 0}
115   * @throws ArithmeticException if {@code mode} is {@link RoundingMode#UNNECESSARY} and {@code x}
116   *     is not a power of two
117   */
118  @SuppressWarnings("fallthrough")
119  // TODO(kevinb): remove after this warning is disabled globally
120  public static int log2(long x, RoundingMode mode) {
121    checkPositive("x", x);
122    switch (mode) {
123      case UNNECESSARY:
124        checkRoundingUnnecessary(isPowerOfTwo(x));
125        // fall through
126      case DOWN:
127      case FLOOR:
128        return (Long.SIZE - 1) - Long.numberOfLeadingZeros(x);
129
130      case UP:
131      case CEILING:
132        return Long.SIZE - Long.numberOfLeadingZeros(x - 1);
133
134      case HALF_DOWN:
135      case HALF_UP:
136      case HALF_EVEN:
137        // Since sqrt(2) is irrational, log2(x) - logFloor cannot be exactly 0.5
138        int leadingZeros = Long.numberOfLeadingZeros(x);
139        long cmp = MAX_POWER_OF_SQRT2_UNSIGNED >>> leadingZeros;
140        // floor(2^(logFloor + 0.5))
141        int logFloor = (Long.SIZE - 1) - leadingZeros;
142        return logFloor + lessThanBranchFree(cmp, x);
143    }
144    throw new AssertionError("impossible");
145  }
146
147  /** The biggest half power of two that fits into an unsigned long */
148  @VisibleForTesting static final long MAX_POWER_OF_SQRT2_UNSIGNED = 0xB504F333F9DE6484L;
149
150  /**
151   * Returns the base-10 logarithm of {@code x}, rounded according to the specified rounding mode.
152   *
153   * @throws IllegalArgumentException if {@code x <= 0}
154   * @throws ArithmeticException if {@code mode} is {@link RoundingMode#UNNECESSARY} and {@code x}
155   *     is not a power of ten
156   */
157  @GwtIncompatible // TODO
158  @SuppressWarnings("fallthrough")
159  // TODO(kevinb): remove after this warning is disabled globally
160  public static int log10(long x, RoundingMode mode) {
161    checkPositive("x", x);
162    int logFloor = log10Floor(x);
163    long floorPow = powersOf10[logFloor];
164    switch (mode) {
165      case UNNECESSARY:
166        checkRoundingUnnecessary(x == floorPow);
167        // fall through
168      case FLOOR:
169      case DOWN:
170        return logFloor;
171      case CEILING:
172      case UP:
173        return logFloor + lessThanBranchFree(floorPow, x);
174      case HALF_DOWN:
175      case HALF_UP:
176      case HALF_EVEN:
177        // sqrt(10) is irrational, so log10(x)-logFloor is never exactly 0.5
178        return logFloor + lessThanBranchFree(halfPowersOf10[logFloor], x);
179    }
180    throw new AssertionError();
181  }
182
183  @GwtIncompatible // TODO
184  static int log10Floor(long x) {
185    /*
186     * Based on Hacker's Delight Fig. 11-5, the two-table-lookup, branch-free implementation.
187     *
188     * The key idea is that based on the number of leading zeros (equivalently, floor(log2(x))), we
189     * can narrow the possible floor(log10(x)) values to two. For example, if floor(log2(x)) is 6,
190     * then 64 <= x < 128, so floor(log10(x)) is either 1 or 2.
191     */
192    int y = maxLog10ForLeadingZeros[Long.numberOfLeadingZeros(x)];
193    /*
194     * y is the higher of the two possible values of floor(log10(x)). If x < 10^y, then we want the
195     * lower of the two possible values, or y - 1, otherwise, we want y.
196     */
197    return y - lessThanBranchFree(x, powersOf10[y]);
198  }
199
200  // maxLog10ForLeadingZeros[i] == floor(log10(2^(Long.SIZE - i)))
201  @VisibleForTesting
202  static final byte[] maxLog10ForLeadingZeros = {
203    19, 18, 18, 18, 18, 17, 17, 17, 16, 16, 16, 15, 15, 15, 15, 14, 14, 14, 13, 13, 13, 12, 12, 12,
204    12, 11, 11, 11, 10, 10, 10, 9, 9, 9, 9, 8, 8, 8, 7, 7, 7, 6, 6, 6, 6, 5, 5, 5, 4, 4, 4, 3, 3, 3,
205    3, 2, 2, 2, 1, 1, 1, 0, 0, 0
206  };
207
208  @GwtIncompatible // TODO
209  @VisibleForTesting
210  static final long[] powersOf10 = {
211    1L,
212    10L,
213    100L,
214    1000L,
215    10000L,
216    100000L,
217    1000000L,
218    10000000L,
219    100000000L,
220    1000000000L,
221    10000000000L,
222    100000000000L,
223    1000000000000L,
224    10000000000000L,
225    100000000000000L,
226    1000000000000000L,
227    10000000000000000L,
228    100000000000000000L,
229    1000000000000000000L
230  };
231
232  // halfPowersOf10[i] = largest long less than 10^(i + 0.5)
233  @GwtIncompatible // TODO
234  @VisibleForTesting
235  static final long[] halfPowersOf10 = {
236    3L,
237    31L,
238    316L,
239    3162L,
240    31622L,
241    316227L,
242    3162277L,
243    31622776L,
244    316227766L,
245    3162277660L,
246    31622776601L,
247    316227766016L,
248    3162277660168L,
249    31622776601683L,
250    316227766016837L,
251    3162277660168379L,
252    31622776601683793L,
253    316227766016837933L,
254    3162277660168379331L
255  };
256
257  /**
258   * Returns {@code b} to the {@code k}th power. Even if the result overflows, it will be equal to
259   * {@code BigInteger.valueOf(b).pow(k).longValue()}. This implementation runs in {@code O(log k)}
260   * time.
261   *
262   * @throws IllegalArgumentException if {@code k < 0}
263   */
264  @GwtIncompatible // TODO
265  public static long pow(long b, int k) {
266    checkNonNegative("exponent", k);
267    if (-2 <= b && b <= 2) {
268      switch ((int) b) {
269        case 0:
270          return (k == 0) ? 1 : 0;
271        case 1:
272          return 1;
273        case (-1):
274          return ((k & 1) == 0) ? 1 : -1;
275        case 2:
276          return (k < Long.SIZE) ? 1L << k : 0;
277        case (-2):
278          if (k < Long.SIZE) {
279            return ((k & 1) == 0) ? 1L << k : -(1L << k);
280          } else {
281            return 0;
282          }
283        default:
284          throw new AssertionError();
285      }
286    }
287    for (long accum = 1; ; k >>= 1) {
288      switch (k) {
289        case 0:
290          return accum;
291        case 1:
292          return accum * b;
293        default:
294          accum *= ((k & 1) == 0) ? 1 : b;
295          b *= b;
296      }
297    }
298  }
299
300  /**
301   * Returns the square root of {@code x}, rounded with the specified rounding mode.
302   *
303   * @throws IllegalArgumentException if {@code x < 0}
304   * @throws ArithmeticException if {@code mode} is {@link RoundingMode#UNNECESSARY} and {@code
305   *     sqrt(x)} is not an integer
306   */
307  @GwtIncompatible // TODO
308  public static long sqrt(long x, RoundingMode mode) {
309    checkNonNegative("x", x);
310    if (fitsInInt(x)) {
311      return IntMath.sqrt((int) x, mode);
312    }
313    /*
314     * Let k be the true value of floor(sqrt(x)), so that
315     *
316     *            k * k <= x          <  (k + 1) * (k + 1)
317     * (double) (k * k) <= (double) x <= (double) ((k + 1) * (k + 1))
318     *          since casting to double is nondecreasing.
319     *          Note that the right-hand inequality is no longer strict.
320     * Math.sqrt(k * k) <= Math.sqrt(x) <= Math.sqrt((k + 1) * (k + 1))
321     *          since Math.sqrt is monotonic.
322     * (long) Math.sqrt(k * k) <= (long) Math.sqrt(x) <= (long) Math.sqrt((k + 1) * (k + 1))
323     *          since casting to long is monotonic
324     * k <= (long) Math.sqrt(x) <= k + 1
325     *          since (long) Math.sqrt(k * k) == k, as checked exhaustively in
326     *          {@link LongMathTest#testSqrtOfPerfectSquareAsDoubleIsPerfect}
327     */
328    long guess = (long) Math.sqrt((double) x);
329    // Note: guess is always <= FLOOR_SQRT_MAX_LONG.
330    long guessSquared = guess * guess;
331    // Note (2013-2-26): benchmarks indicate that, inscrutably enough, using if statements is
332    // faster here than using lessThanBranchFree.
333    switch (mode) {
334      case UNNECESSARY:
335        checkRoundingUnnecessary(guessSquared == x);
336        return guess;
337      case FLOOR:
338      case DOWN:
339        if (x < guessSquared) {
340          return guess - 1;
341        }
342        return guess;
343      case CEILING:
344      case UP:
345        if (x > guessSquared) {
346          return guess + 1;
347        }
348        return guess;
349      case HALF_DOWN:
350      case HALF_UP:
351      case HALF_EVEN:
352        long sqrtFloor = guess - ((x < guessSquared) ? 1 : 0);
353        long halfSquare = sqrtFloor * sqrtFloor + sqrtFloor;
354        /*
355         * We wish to test whether or not x <= (sqrtFloor + 0.5)^2 = halfSquare + 0.25. Since both x
356         * and halfSquare are integers, this is equivalent to testing whether or not x <=
357         * halfSquare. (We have to deal with overflow, though.)
358         *
359         * If we treat halfSquare as an unsigned long, we know that
360         *            sqrtFloor^2 <= x < (sqrtFloor + 1)^2
361         * halfSquare - sqrtFloor <= x < halfSquare + sqrtFloor + 1
362         * so |x - halfSquare| <= sqrtFloor.  Therefore, it's safe to treat x - halfSquare as a
363         * signed long, so lessThanBranchFree is safe for use.
364         */
365        return sqrtFloor + lessThanBranchFree(halfSquare, x);
366    }
367    throw new AssertionError();
368  }
369
370  /**
371   * Returns the result of dividing {@code p} by {@code q}, rounding using the specified {@code
372   * RoundingMode}.
373   *
374   * @throws ArithmeticException if {@code q == 0}, or if {@code mode == UNNECESSARY} and {@code a}
375   *     is not an integer multiple of {@code b}
376   */
377  @GwtIncompatible // TODO
378  @SuppressWarnings("fallthrough")
379  public static long divide(long p, long q, RoundingMode mode) {
380    checkNotNull(mode);
381    long div = p / q; // throws if q == 0
382    long rem = p - q * div; // equals p % q
383
384    if (rem == 0) {
385      return div;
386    }
387
388    /*
389     * Normal Java division rounds towards 0, consistently with RoundingMode.DOWN. We just have to
390     * deal with the cases where rounding towards 0 is wrong, which typically depends on the sign of
391     * p / q.
392     *
393     * signum is 1 if p and q are both nonnegative or both negative, and -1 otherwise.
394     */
395    int signum = 1 | (int) ((p ^ q) >> (Long.SIZE - 1));
396    boolean increment;
397    switch (mode) {
398      case UNNECESSARY:
399        checkRoundingUnnecessary(rem == 0);
400        // fall through
401      case DOWN:
402        increment = false;
403        break;
404      case UP:
405        increment = true;
406        break;
407      case CEILING:
408        increment = signum > 0;
409        break;
410      case FLOOR:
411        increment = signum < 0;
412        break;
413      case HALF_EVEN:
414      case HALF_DOWN:
415      case HALF_UP:
416        long absRem = abs(rem);
417        long cmpRemToHalfDivisor = absRem - (abs(q) - absRem);
418        // subtracting two nonnegative longs can't overflow
419        // cmpRemToHalfDivisor has the same sign as compare(abs(rem), abs(q) / 2).
420        if (cmpRemToHalfDivisor == 0) { // exactly on the half mark
421          increment = (mode == HALF_UP || (mode == HALF_EVEN && (div & 1) != 0));
422        } else {
423          increment = cmpRemToHalfDivisor > 0; // closer to the UP value
424        }
425        break;
426      default:
427        throw new AssertionError();
428    }
429    return increment ? div + signum : div;
430  }
431
432  /**
433   * Returns {@code x mod m}, a non-negative value less than {@code m}. This differs from {@code x %
434   * m}, which might be negative.
435   *
436   * <p>For example:
437   *
438   * <pre>{@code
439   * mod(7, 4) == 3
440   * mod(-7, 4) == 1
441   * mod(-1, 4) == 3
442   * mod(-8, 4) == 0
443   * mod(8, 4) == 0
444   * }</pre>
445   *
446   * @throws ArithmeticException if {@code m <= 0}
447   * @see <a href="http://docs.oracle.com/javase/specs/jls/se7/html/jls-15.html#jls-15.17.3">
448   *     Remainder Operator</a>
449   */
450  @GwtIncompatible // TODO
451  public static int mod(long x, int m) {
452    // Cast is safe because the result is guaranteed in the range [0, m)
453    return (int) mod(x, (long) m);
454  }
455
456  /**
457   * Returns {@code x mod m}, a non-negative value less than {@code m}. This differs from {@code x %
458   * m}, which might be negative.
459   *
460   * <p>For example:
461   *
462   * <pre>{@code
463   * mod(7, 4) == 3
464   * mod(-7, 4) == 1
465   * mod(-1, 4) == 3
466   * mod(-8, 4) == 0
467   * mod(8, 4) == 0
468   * }</pre>
469   *
470   * @throws ArithmeticException if {@code m <= 0}
471   * @see <a href="http://docs.oracle.com/javase/specs/jls/se7/html/jls-15.html#jls-15.17.3">
472   *     Remainder Operator</a>
473   */
474  @GwtIncompatible // TODO
475  public static long mod(long x, long m) {
476    if (m <= 0) {
477      throw new ArithmeticException("Modulus must be positive");
478    }
479    long result = x % m;
480    return (result >= 0) ? result : result + m;
481  }
482
483  /**
484   * Returns the greatest common divisor of {@code a, b}. Returns {@code 0} if {@code a == 0 && b ==
485   * 0}.
486   *
487   * @throws IllegalArgumentException if {@code a < 0} or {@code b < 0}
488   */
489  public static long gcd(long a, long b) {
490    /*
491     * The reason we require both arguments to be >= 0 is because otherwise, what do you return on
492     * gcd(0, Long.MIN_VALUE)? BigInteger.gcd would return positive 2^63, but positive 2^63 isn't an
493     * int.
494     */
495    checkNonNegative("a", a);
496    checkNonNegative("b", b);
497    if (a == 0) {
498      // 0 % b == 0, so b divides a, but the converse doesn't hold.
499      // BigInteger.gcd is consistent with this decision.
500      return b;
501    } else if (b == 0) {
502      return a; // similar logic
503    }
504    /*
505     * Uses the binary GCD algorithm; see http://en.wikipedia.org/wiki/Binary_GCD_algorithm. This is
506     * >60% faster than the Euclidean algorithm in benchmarks.
507     */
508    int aTwos = Long.numberOfTrailingZeros(a);
509    a >>= aTwos; // divide out all 2s
510    int bTwos = Long.numberOfTrailingZeros(b);
511    b >>= bTwos; // divide out all 2s
512    while (a != b) { // both a, b are odd
513      // The key to the binary GCD algorithm is as follows:
514      // Both a and b are odd. Assume a > b; then gcd(a - b, b) = gcd(a, b).
515      // But in gcd(a - b, b), a - b is even and b is odd, so we can divide out powers of two.
516
517      // We bend over backwards to avoid branching, adapting a technique from
518      // http://graphics.stanford.edu/~seander/bithacks.html#IntegerMinOrMax
519
520      long delta = a - b; // can't overflow, since a and b are nonnegative
521
522      long minDeltaOrZero = delta & (delta >> (Long.SIZE - 1));
523      // equivalent to Math.min(delta, 0)
524
525      a = delta - minDeltaOrZero - minDeltaOrZero; // sets a to Math.abs(a - b)
526      // a is now nonnegative and even
527
528      b += minDeltaOrZero; // sets b to min(old a, b)
529      a >>= Long.numberOfTrailingZeros(a); // divide out all 2s, since 2 doesn't divide b
530    }
531    return a << min(aTwos, bTwos);
532  }
533
534  /**
535   * Returns the sum of {@code a} and {@code b}, provided it does not overflow.
536   *
537   * @throws ArithmeticException if {@code a + b} overflows in signed {@code long} arithmetic
538   */
539  @SuppressWarnings("ShortCircuitBoolean")
540  public static long checkedAdd(long a, long b) {
541    long result = a + b;
542    checkNoOverflow((a ^ b) < 0 | (a ^ result) >= 0, "checkedAdd", a, b);
543    return result;
544  }
545
546  /**
547   * Returns the difference of {@code a} and {@code b}, provided it does not overflow.
548   *
549   * @throws ArithmeticException if {@code a - b} overflows in signed {@code long} arithmetic
550   */
551  @GwtIncompatible // TODO
552  @SuppressWarnings("ShortCircuitBoolean")
553  public static long checkedSubtract(long a, long b) {
554    long result = a - b;
555    checkNoOverflow((a ^ b) >= 0 | (a ^ result) >= 0, "checkedSubtract", a, b);
556    return result;
557  }
558
559  /**
560   * Returns the product of {@code a} and {@code b}, provided it does not overflow.
561   *
562   * @throws ArithmeticException if {@code a * b} overflows in signed {@code long} arithmetic
563   */
564  @SuppressWarnings("ShortCircuitBoolean")
565  public static long checkedMultiply(long a, long b) {
566    // Hacker's Delight, Section 2-12
567    int leadingZeros =
568        Long.numberOfLeadingZeros(a)
569            + Long.numberOfLeadingZeros(~a)
570            + Long.numberOfLeadingZeros(b)
571            + Long.numberOfLeadingZeros(~b);
572    /*
573     * If leadingZeros > Long.SIZE + 1 it's definitely fine, if it's < Long.SIZE it's definitely
574     * bad. We do the leadingZeros check to avoid the division below if at all possible.
575     *
576     * Otherwise, if b == Long.MIN_VALUE, then the only allowed values of a are 0 and 1. We take
577     * care of all a < 0 with their own check, because in particular, the case a == -1 will
578     * incorrectly pass the division check below.
579     *
580     * In all other cases, we check that either a is 0 or the result is consistent with division.
581     */
582    if (leadingZeros > Long.SIZE + 1) {
583      return a * b;
584    }
585    checkNoOverflow(leadingZeros >= Long.SIZE, "checkedMultiply", a, b);
586    checkNoOverflow(a >= 0 | b != Long.MIN_VALUE, "checkedMultiply", a, b);
587    long result = a * b;
588    checkNoOverflow(a == 0 || result / a == b, "checkedMultiply", a, b);
589    return result;
590  }
591
592  /**
593   * Returns the {@code b} to the {@code k}th power, provided it does not overflow.
594   *
595   * @throws ArithmeticException if {@code b} to the {@code k}th power overflows in signed {@code
596   *     long} arithmetic
597   */
598  @GwtIncompatible // TODO
599  @SuppressWarnings("ShortCircuitBoolean")
600  public static long checkedPow(long b, int k) {
601    checkNonNegative("exponent", k);
602    if (b >= -2 & b <= 2) {
603      switch ((int) b) {
604        case 0:
605          return (k == 0) ? 1 : 0;
606        case 1:
607          return 1;
608        case (-1):
609          return ((k & 1) == 0) ? 1 : -1;
610        case 2:
611          checkNoOverflow(k < Long.SIZE - 1, "checkedPow", b, k);
612          return 1L << k;
613        case (-2):
614          checkNoOverflow(k < Long.SIZE, "checkedPow", b, k);
615          return ((k & 1) == 0) ? (1L << k) : (-1L << k);
616        default:
617          throw new AssertionError();
618      }
619    }
620    long accum = 1;
621    while (true) {
622      switch (k) {
623        case 0:
624          return accum;
625        case 1:
626          return checkedMultiply(accum, b);
627        default:
628          if ((k & 1) != 0) {
629            accum = checkedMultiply(accum, b);
630          }
631          k >>= 1;
632          if (k > 0) {
633            checkNoOverflow(
634                -FLOOR_SQRT_MAX_LONG <= b && b <= FLOOR_SQRT_MAX_LONG, "checkedPow", b, k);
635            b *= b;
636          }
637      }
638    }
639  }
640
641  /**
642   * Returns the sum of {@code a} and {@code b} unless it would overflow or underflow in which case
643   * {@code Long.MAX_VALUE} or {@code Long.MIN_VALUE} is returned, respectively.
644   *
645   * @since 20.0
646   */
647  @SuppressWarnings("ShortCircuitBoolean")
648  public static long saturatedAdd(long a, long b) {
649    long naiveSum = a + b;
650    if ((a ^ b) < 0 | (a ^ naiveSum) >= 0) {
651      // If a and b have different signs or a has the same sign as the result then there was no
652      // overflow, return.
653      return naiveSum;
654    }
655    // we did over/under flow, if the sign is negative we should return MAX otherwise MIN
656    return Long.MAX_VALUE + ((naiveSum >>> (Long.SIZE - 1)) ^ 1);
657  }
658
659  /**
660   * Returns the difference of {@code a} and {@code b} unless it would overflow or underflow in
661   * which case {@code Long.MAX_VALUE} or {@code Long.MIN_VALUE} is returned, respectively.
662   *
663   * @since 20.0
664   */
665  @SuppressWarnings("ShortCircuitBoolean")
666  public static long saturatedSubtract(long a, long b) {
667    long naiveDifference = a - b;
668    if ((a ^ b) >= 0 | (a ^ naiveDifference) >= 0) {
669      // If a and b have the same signs or a has the same sign as the result then there was no
670      // overflow, return.
671      return naiveDifference;
672    }
673    // we did over/under flow
674    return Long.MAX_VALUE + ((naiveDifference >>> (Long.SIZE - 1)) ^ 1);
675  }
676
677  /**
678   * Returns the product of {@code a} and {@code b} unless it would overflow or underflow in which
679   * case {@code Long.MAX_VALUE} or {@code Long.MIN_VALUE} is returned, respectively.
680   *
681   * @since 20.0
682   */
683  @SuppressWarnings("ShortCircuitBoolean")
684  public static long saturatedMultiply(long a, long b) {
685    // see checkedMultiply for explanation
686    int leadingZeros =
687        Long.numberOfLeadingZeros(a)
688            + Long.numberOfLeadingZeros(~a)
689            + Long.numberOfLeadingZeros(b)
690            + Long.numberOfLeadingZeros(~b);
691    if (leadingZeros > Long.SIZE + 1) {
692      return a * b;
693    }
694    // the return value if we will overflow (which we calculate by overflowing a long :) )
695    long limit = Long.MAX_VALUE + ((a ^ b) >>> (Long.SIZE - 1));
696    if (leadingZeros < Long.SIZE | (a < 0 & b == Long.MIN_VALUE)) {
697      // overflow
698      return limit;
699    }
700    long result = a * b;
701    if (a == 0 || result / a == b) {
702      return result;
703    }
704    return limit;
705  }
706
707  /**
708   * Returns the {@code b} to the {@code k}th power, unless it would overflow or underflow in which
709   * case {@code Long.MAX_VALUE} or {@code Long.MIN_VALUE} is returned, respectively.
710   *
711   * @since 20.0
712   */
713  @SuppressWarnings("ShortCircuitBoolean")
714  public static long saturatedPow(long b, int k) {
715    checkNonNegative("exponent", k);
716    if (b >= -2 & b <= 2) {
717      switch ((int) b) {
718        case 0:
719          return (k == 0) ? 1 : 0;
720        case 1:
721          return 1;
722        case (-1):
723          return ((k & 1) == 0) ? 1 : -1;
724        case 2:
725          if (k >= Long.SIZE - 1) {
726            return Long.MAX_VALUE;
727          }
728          return 1L << k;
729        case (-2):
730          if (k >= Long.SIZE) {
731            return Long.MAX_VALUE + (k & 1);
732          }
733          return ((k & 1) == 0) ? (1L << k) : (-1L << k);
734        default:
735          throw new AssertionError();
736      }
737    }
738    long accum = 1;
739    // if b is negative and k is odd then the limit is MIN otherwise the limit is MAX
740    long limit = Long.MAX_VALUE + ((b >>> (Long.SIZE - 1)) & (k & 1));
741    while (true) {
742      switch (k) {
743        case 0:
744          return accum;
745        case 1:
746          return saturatedMultiply(accum, b);
747        default:
748          if ((k & 1) != 0) {
749            accum = saturatedMultiply(accum, b);
750          }
751          k >>= 1;
752          if (k > 0) {
753            if (-FLOOR_SQRT_MAX_LONG > b | b > FLOOR_SQRT_MAX_LONG) {
754              return limit;
755            }
756            b *= b;
757          }
758      }
759    }
760  }
761
762  @VisibleForTesting static final long FLOOR_SQRT_MAX_LONG = 3037000499L;
763
764  /**
765   * Returns {@code n!}, that is, the product of the first {@code n} positive integers, {@code 1} if
766   * {@code n == 0}, or {@link Long#MAX_VALUE} if the result does not fit in a {@code long}.
767   *
768   * @throws IllegalArgumentException if {@code n < 0}
769   */
770  @GwtIncompatible // TODO
771  public static long factorial(int n) {
772    checkNonNegative("n", n);
773    return (n < factorials.length) ? factorials[n] : Long.MAX_VALUE;
774  }
775
776  static final long[] factorials = {
777    1L,
778    1L,
779    1L * 2,
780    1L * 2 * 3,
781    1L * 2 * 3 * 4,
782    1L * 2 * 3 * 4 * 5,
783    1L * 2 * 3 * 4 * 5 * 6,
784    1L * 2 * 3 * 4 * 5 * 6 * 7,
785    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8,
786    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9,
787    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10,
788    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11,
789    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12,
790    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13,
791    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14,
792    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15,
793    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16,
794    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17,
795    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17 * 18,
796    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17 * 18 * 19,
797    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17 * 18 * 19 * 20
798  };
799
800  /**
801   * Returns {@code n} choose {@code k}, also known as the binomial coefficient of {@code n} and
802   * {@code k}, or {@link Long#MAX_VALUE} if the result does not fit in a {@code long}.
803   *
804   * @throws IllegalArgumentException if {@code n < 0}, {@code k < 0}, or {@code k > n}
805   */
806  public static long binomial(int n, int k) {
807    checkNonNegative("n", n);
808    checkNonNegative("k", k);
809    checkArgument(k <= n, "k (%s) > n (%s)", k, n);
810    if (k > (n >> 1)) {
811      k = n - k;
812    }
813    switch (k) {
814      case 0:
815        return 1;
816      case 1:
817        return n;
818      default:
819        if (n < factorials.length) {
820          return factorials[n] / (factorials[k] * factorials[n - k]);
821        } else if (k >= biggestBinomials.length || n > biggestBinomials[k]) {
822          return Long.MAX_VALUE;
823        } else if (k < biggestSimpleBinomials.length && n <= biggestSimpleBinomials[k]) {
824          // guaranteed not to overflow
825          long result = n--;
826          for (int i = 2; i <= k; n--, i++) {
827            result *= n;
828            result /= i;
829          }
830          return result;
831        } else {
832          int nBits = LongMath.log2(n, RoundingMode.CEILING);
833
834          long result = 1;
835          long numerator = n--;
836          long denominator = 1;
837
838          int numeratorBits = nBits;
839          // This is an upper bound on log2(numerator, ceiling).
840
841          /*
842           * We want to do this in long math for speed, but want to avoid overflow. We adapt the
843           * technique previously used by BigIntegerMath: maintain separate numerator and
844           * denominator accumulators, multiplying the fraction into result when near overflow.
845           */
846          for (int i = 2; i <= k; i++, n--) {
847            if (numeratorBits + nBits < Long.SIZE - 1) {
848              // It's definitely safe to multiply into numerator and denominator.
849              numerator *= n;
850              denominator *= i;
851              numeratorBits += nBits;
852            } else {
853              // It might not be safe to multiply into numerator and denominator,
854              // so multiply (numerator / denominator) into result.
855              result = multiplyFraction(result, numerator, denominator);
856              numerator = n;
857              denominator = i;
858              numeratorBits = nBits;
859            }
860          }
861          return multiplyFraction(result, numerator, denominator);
862        }
863    }
864  }
865
866  /** Returns (x * numerator / denominator), which is assumed to come out to an integral value. */
867  static long multiplyFraction(long x, long numerator, long denominator) {
868    if (x == 1) {
869      return numerator / denominator;
870    }
871    long commonDivisor = gcd(x, denominator);
872    x /= commonDivisor;
873    denominator /= commonDivisor;
874    // We know gcd(x, denominator) = 1, and x * numerator / denominator is exact,
875    // so denominator must be a divisor of numerator.
876    return x * (numerator / denominator);
877  }
878
879  /*
880   * binomial(biggestBinomials[k], k) fits in a long, but not binomial(biggestBinomials[k] + 1, k).
881   */
882  static final int[] biggestBinomials = {
883    Integer.MAX_VALUE,
884    Integer.MAX_VALUE,
885    Integer.MAX_VALUE,
886    3810779,
887    121977,
888    16175,
889    4337,
890    1733,
891    887,
892    534,
893    361,
894    265,
895    206,
896    169,
897    143,
898    125,
899    111,
900    101,
901    94,
902    88,
903    83,
904    79,
905    76,
906    74,
907    72,
908    70,
909    69,
910    68,
911    67,
912    67,
913    66,
914    66,
915    66,
916    66
917  };
918
919  /*
920   * binomial(biggestSimpleBinomials[k], k) doesn't need to use the slower GCD-based impl, but
921   * binomial(biggestSimpleBinomials[k] + 1, k) does.
922   */
923  @VisibleForTesting
924  static final int[] biggestSimpleBinomials = {
925    Integer.MAX_VALUE,
926    Integer.MAX_VALUE,
927    Integer.MAX_VALUE,
928    2642246,
929    86251,
930    11724,
931    3218,
932    1313,
933    684,
934    419,
935    287,
936    214,
937    169,
938    139,
939    119,
940    105,
941    95,
942    87,
943    81,
944    76,
945    73,
946    70,
947    68,
948    66,
949    64,
950    63,
951    62,
952    62,
953    61,
954    61,
955    61
956  };
957  // These values were generated by using checkedMultiply to see when the simple multiply/divide
958  // algorithm would lead to an overflow.
959
960  static boolean fitsInInt(long x) {
961    return (int) x == x;
962  }
963
964  /**
965   * Returns the arithmetic mean of {@code x} and {@code y}, rounded toward negative infinity. This
966   * method is resilient to overflow.
967   *
968   * @since 14.0
969   */
970  public static long mean(long x, long y) {
971    // Efficient method for computing the arithmetic mean.
972    // The alternative (x + y) / 2 fails for large values.
973    // The alternative (x + y) >>> 1 fails for negative values.
974    return (x & y) + ((x ^ y) >> 1);
975  }
976
977  /*
978   * This bitmask is used as an optimization for cheaply testing for divisibility by 2, 3, or 5.
979   * Each bit is set to 1 for all remainders that indicate divisibility by 2, 3, or 5, so
980   * 1, 7, 11, 13, 17, 19, 23, 29 are set to 0. 30 and up don't matter because they won't be hit.
981   */
982  private static final int SIEVE_30 =
983      ~((1 << 1) | (1 << 7) | (1 << 11) | (1 << 13) | (1 << 17) | (1 << 19) | (1 << 23)
984          | (1 << 29));
985
986  /**
987   * Returns {@code true} if {@code n} is a <a
988   * href="http://mathworld.wolfram.com/PrimeNumber.html">prime number</a>: an integer <i>greater
989   * than one</i> that cannot be factored into a product of <i>smaller</i> positive integers.
990   * Returns {@code false} if {@code n} is zero, one, or a composite number (one which <i>can</i> be
991   * factored into smaller positive integers).
992   *
993   * <p>To test larger numbers, use {@link BigInteger#isProbablePrime}.
994   *
995   * @throws IllegalArgumentException if {@code n} is negative
996   * @since 20.0
997   */
998  @GwtIncompatible // TODO
999  public static boolean isPrime(long n) {
1000    if (n < 2) {
1001      checkNonNegative("n", n);
1002      return false;
1003    }
1004    if (n < 66) {
1005      // Encode all primes less than 66 into mask without 0 and 1.
1006      long mask =
1007          (1L << (2 - 2))
1008              | (1L << (3 - 2))
1009              | (1L << (5 - 2))
1010              | (1L << (7 - 2))
1011              | (1L << (11 - 2))
1012              | (1L << (13 - 2))
1013              | (1L << (17 - 2))
1014              | (1L << (19 - 2))
1015              | (1L << (23 - 2))
1016              | (1L << (29 - 2))
1017              | (1L << (31 - 2))
1018              | (1L << (37 - 2))
1019              | (1L << (41 - 2))
1020              | (1L << (43 - 2))
1021              | (1L << (47 - 2))
1022              | (1L << (53 - 2))
1023              | (1L << (59 - 2))
1024              | (1L << (61 - 2));
1025      // Look up n within the mask.
1026      return ((mask >> ((int) n - 2)) & 1) != 0;
1027    }
1028
1029    if ((SIEVE_30 & (1 << (n % 30))) != 0) {
1030      return false;
1031    }
1032    if (n % 7 == 0 || n % 11 == 0 || n % 13 == 0) {
1033      return false;
1034    }
1035    if (n < 17 * 17) {
1036      return true;
1037    }
1038
1039    for (long[] baseSet : millerRabinBaseSets) {
1040      if (n <= baseSet[0]) {
1041        for (int i = 1; i < baseSet.length; i++) {
1042          if (!MillerRabinTester.test(baseSet[i], n)) {
1043            return false;
1044          }
1045        }
1046        return true;
1047      }
1048    }
1049    throw new AssertionError();
1050  }
1051
1052  /*
1053   * If n <= millerRabinBases[i][0], then testing n against bases millerRabinBases[i][1..] suffices
1054   * to prove its primality. Values from miller-rabin.appspot.com.
1055   *
1056   * NOTE: We could get slightly better bases that would be treated as unsigned, but benchmarks
1057   * showed negligible performance improvements.
1058   */
1059  private static final long[][] millerRabinBaseSets = {
1060    {291830, 126401071349994536L},
1061    {885594168, 725270293939359937L, 3569819667048198375L},
1062    {273919523040L, 15, 7363882082L, 992620450144556L},
1063    {47636622961200L, 2, 2570940, 211991001, 3749873356L},
1064    {
1065      7999252175582850L,
1066      2,
1067      4130806001517L,
1068      149795463772692060L,
1069      186635894390467037L,
1070      3967304179347715805L
1071    },
1072    {
1073      585226005592931976L,
1074      2,
1075      123635709730000L,
1076      9233062284813009L,
1077      43835965440333360L,
1078      761179012939631437L,
1079      1263739024124850375L
1080    },
1081    {Long.MAX_VALUE, 2, 325, 9375, 28178, 450775, 9780504, 1795265022}
1082  };
1083
1084  private enum MillerRabinTester {
1085    /** Works for inputs ≤ FLOOR_SQRT_MAX_LONG. */
1086    SMALL {
1087      @Override
1088      long mulMod(long a, long b, long m) {
1089        /*
1090         * lowasser, 2015-Feb-12: Benchmarks suggest that changing this to UnsignedLongs.remainder
1091         * and increasing the threshold to 2^32 doesn't pay for itself, and adding another enum
1092         * constant hurts performance further -- I suspect because bimorphic implementation is a
1093         * sweet spot for the JVM.
1094         */
1095        return (a * b) % m;
1096      }
1097
1098      @Override
1099      long squareMod(long a, long m) {
1100        return (a * a) % m;
1101      }
1102    },
1103    /** Works for all nonnegative signed longs. */
1104    LARGE {
1105      /** Returns (a + b) mod m. Precondition: {@code 0 <= a}, {@code b < m < 2^63}. */
1106      private long plusMod(long a, long b, long m) {
1107        return (a >= m - b) ? (a + b - m) : (a + b);
1108      }
1109
1110      /** Returns (a * 2^32) mod m. a may be any unsigned long. */
1111      private long times2ToThe32Mod(long a, long m) {
1112        int remainingPowersOf2 = 32;
1113        do {
1114          int shift = min(remainingPowersOf2, Long.numberOfLeadingZeros(a));
1115          // shift is either the number of powers of 2 left to multiply a by, or the biggest shift
1116          // possible while keeping a in an unsigned long.
1117          a = UnsignedLongs.remainder(a << shift, m);
1118          remainingPowersOf2 -= shift;
1119        } while (remainingPowersOf2 > 0);
1120        return a;
1121      }
1122
1123      @Override
1124      long mulMod(long a, long b, long m) {
1125        long aHi = a >>> 32; // < 2^31
1126        long bHi = b >>> 32; // < 2^31
1127        long aLo = a & 0xFFFFFFFFL; // < 2^32
1128        long bLo = b & 0xFFFFFFFFL; // < 2^32
1129
1130        /*
1131         * a * b == aHi * bHi * 2^64 + (aHi * bLo + aLo * bHi) * 2^32 + aLo * bLo.
1132         *       == (aHi * bHi * 2^32 + aHi * bLo + aLo * bHi) * 2^32 + aLo * bLo
1133         *
1134         * We carry out this computation in modular arithmetic. Since times2ToThe32Mod accepts any
1135         * unsigned long, we don't have to do a mod on every operation, only when intermediate
1136         * results can exceed 2^63.
1137         */
1138        long result = times2ToThe32Mod(aHi * bHi /* < 2^62 */, m); // < m < 2^63
1139        result += aHi * bLo; // aHi * bLo < 2^63, result < 2^64
1140        if (result < 0) {
1141          result = UnsignedLongs.remainder(result, m);
1142        }
1143        // result < 2^63 again
1144        result += aLo * bHi; // aLo * bHi < 2^63, result < 2^64
1145        result = times2ToThe32Mod(result, m); // result < m < 2^63
1146        return plusMod(result, UnsignedLongs.remainder(aLo * bLo /* < 2^64 */, m), m);
1147      }
1148
1149      @Override
1150      long squareMod(long a, long m) {
1151        long aHi = a >>> 32; // < 2^31
1152        long aLo = a & 0xFFFFFFFFL; // < 2^32
1153
1154        /*
1155         * a^2 == aHi^2 * 2^64 + aHi * aLo * 2^33 + aLo^2
1156         *     == (aHi^2 * 2^32 + aHi * aLo * 2) * 2^32 + aLo^2
1157         * We carry out this computation in modular arithmetic.  Since times2ToThe32Mod accepts any
1158         * unsigned long, we don't have to do a mod on every operation, only when intermediate
1159         * results can exceed 2^63.
1160         */
1161        long result = times2ToThe32Mod(aHi * aHi /* < 2^62 */, m); // < m < 2^63
1162        long hiLo = aHi * aLo * 2;
1163        if (hiLo < 0) {
1164          hiLo = UnsignedLongs.remainder(hiLo, m);
1165        }
1166        // hiLo < 2^63
1167        result += hiLo; // result < 2^64
1168        result = times2ToThe32Mod(result, m); // result < m < 2^63
1169        return plusMod(result, UnsignedLongs.remainder(aLo * aLo /* < 2^64 */, m), m);
1170      }
1171    };
1172
1173    static boolean test(long base, long n) {
1174      // Since base will be considered % n, it's okay if base > FLOOR_SQRT_MAX_LONG,
1175      // so long as n <= FLOOR_SQRT_MAX_LONG.
1176      return ((n <= FLOOR_SQRT_MAX_LONG) ? SMALL : LARGE).testWitness(base, n);
1177    }
1178
1179    /** Returns a * b mod m. */
1180    abstract long mulMod(long a, long b, long m);
1181
1182    /** Returns a^2 mod m. */
1183    abstract long squareMod(long a, long m);
1184
1185    /** Returns a^p mod m. */
1186    private long powMod(long a, long p, long m) {
1187      long res = 1;
1188      for (; p != 0; p >>= 1) {
1189        if ((p & 1) != 0) {
1190          res = mulMod(res, a, m);
1191        }
1192        a = squareMod(a, m);
1193      }
1194      return res;
1195    }
1196
1197    /** Returns true if n is a strong probable prime relative to the specified base. */
1198    private boolean testWitness(long base, long n) {
1199      int r = Long.numberOfTrailingZeros(n - 1);
1200      long d = (n - 1) >> r;
1201      base %= n;
1202      if (base == 0) {
1203        return true;
1204      }
1205      // Calculate a := base^d mod n.
1206      long a = powMod(base, d, n);
1207      // n passes this test if
1208      //    base^d = 1 (mod n)
1209      // or base^(2^j * d) = -1 (mod n) for some 0 <= j < r.
1210      if (a == 1) {
1211        return true;
1212      }
1213      int j = 0;
1214      while (a != n - 1) {
1215        if (++j == r) {
1216          return false;
1217        }
1218        a = squareMod(a, n);
1219      }
1220      return true;
1221    }
1222  }
1223
1224  /**
1225   * Returns {@code x}, rounded to a {@code double} with the specified rounding mode. If {@code x}
1226   * is precisely representable as a {@code double}, its {@code double} value will be returned;
1227   * otherwise, the rounding will choose between the two nearest representable values with {@code
1228   * mode}.
1229   *
1230   * <p>For the case of {@link RoundingMode#HALF_EVEN}, this implementation uses the IEEE 754
1231   * default rounding mode: if the two nearest representable values are equally near, the one with
1232   * the least significant bit zero is chosen. (In such cases, both of the nearest representable
1233   * values are even integers; this method returns the one that is a multiple of a greater power of
1234   * two.)
1235   *
1236   * @throws ArithmeticException if {@code mode} is {@link RoundingMode#UNNECESSARY} and {@code x}
1237   *     is not precisely representable as a {@code double}
1238   * @since 30.0
1239   */
1240  @SuppressWarnings("deprecation")
1241  @GwtIncompatible
1242  public static double roundToDouble(long x, RoundingMode mode) {
1243    // Logic adapted from ToDoubleRounder.
1244    double roundArbitrarily = (double) x;
1245    long roundArbitrarilyAsLong = (long) roundArbitrarily;
1246    int cmpXToRoundArbitrarily;
1247
1248    if (roundArbitrarilyAsLong == Long.MAX_VALUE) {
1249      /*
1250       * For most values, the conversion from roundArbitrarily to roundArbitrarilyAsLong is
1251       * lossless. In that case we can compare x to roundArbitrarily using Longs.compare(x,
1252       * roundArbitrarilyAsLong). The exception is for values where the conversion to double rounds
1253       * up to give roundArbitrarily equal to 2^63, so the conversion back to long overflows and
1254       * roundArbitrarilyAsLong is Long.MAX_VALUE. (This is the only way this condition can occur as
1255       * otherwise the conversion back to long pads with zero bits.) In this case we know that
1256       * roundArbitrarily > x. (This is important when x == Long.MAX_VALUE ==
1257       * roundArbitrarilyAsLong.)
1258       */
1259      cmpXToRoundArbitrarily = -1;
1260    } else {
1261      cmpXToRoundArbitrarily = Longs.compare(x, roundArbitrarilyAsLong);
1262    }
1263
1264    switch (mode) {
1265      case UNNECESSARY:
1266        checkRoundingUnnecessary(cmpXToRoundArbitrarily == 0);
1267        return roundArbitrarily;
1268      case FLOOR:
1269        return (cmpXToRoundArbitrarily >= 0)
1270            ? roundArbitrarily
1271            : DoubleUtils.nextDown(roundArbitrarily);
1272      case CEILING:
1273        return (cmpXToRoundArbitrarily <= 0) ? roundArbitrarily : Math.nextUp(roundArbitrarily);
1274      case DOWN:
1275        if (x >= 0) {
1276          return (cmpXToRoundArbitrarily >= 0)
1277              ? roundArbitrarily
1278              : DoubleUtils.nextDown(roundArbitrarily);
1279        } else {
1280          return (cmpXToRoundArbitrarily <= 0) ? roundArbitrarily : Math.nextUp(roundArbitrarily);
1281        }
1282      case UP:
1283        if (x >= 0) {
1284          return (cmpXToRoundArbitrarily <= 0) ? roundArbitrarily : Math.nextUp(roundArbitrarily);
1285        } else {
1286          return (cmpXToRoundArbitrarily >= 0)
1287              ? roundArbitrarily
1288              : DoubleUtils.nextDown(roundArbitrarily);
1289        }
1290      case HALF_DOWN:
1291      case HALF_UP:
1292      case HALF_EVEN:
1293        {
1294          long roundFloor;
1295          double roundFloorAsDouble;
1296          long roundCeiling;
1297          double roundCeilingAsDouble;
1298
1299          if (cmpXToRoundArbitrarily >= 0) {
1300            roundFloorAsDouble = roundArbitrarily;
1301            roundFloor = roundArbitrarilyAsLong;
1302            roundCeilingAsDouble = Math.nextUp(roundArbitrarily);
1303            roundCeiling = (long) Math.ceil(roundCeilingAsDouble);
1304          } else {
1305            roundCeilingAsDouble = roundArbitrarily;
1306            roundCeiling = roundArbitrarilyAsLong;
1307            roundFloorAsDouble = DoubleUtils.nextDown(roundArbitrarily);
1308            roundFloor = (long) Math.floor(roundFloorAsDouble);
1309          }
1310
1311          long deltaToFloor = x - roundFloor;
1312          long deltaToCeiling = roundCeiling - x;
1313
1314          if (roundCeiling == Long.MAX_VALUE) {
1315            // correct for Long.MAX_VALUE as discussed above: roundCeilingAsDouble must be 2^63, but
1316            // roundCeiling is 2^63-1.
1317            deltaToCeiling++;
1318          }
1319
1320          int diff = Longs.compare(deltaToFloor, deltaToCeiling);
1321          if (diff < 0) { // closer to floor
1322            return roundFloorAsDouble;
1323          } else if (diff > 0) { // closer to ceiling
1324            return roundCeilingAsDouble;
1325          }
1326          // halfway between the representable values; do the half-whatever logic
1327          switch (mode) {
1328            case HALF_EVEN:
1329              return ((DoubleUtils.getSignificand(roundFloorAsDouble) & 1L) == 0)
1330                  ? roundFloorAsDouble
1331                  : roundCeilingAsDouble;
1332            case HALF_DOWN:
1333              return (x >= 0) ? roundFloorAsDouble : roundCeilingAsDouble;
1334            case HALF_UP:
1335              return (x >= 0) ? roundCeilingAsDouble : roundFloorAsDouble;
1336            default:
1337              throw new AssertionError("impossible");
1338          }
1339        }
1340    }
1341    throw new AssertionError("impossible");
1342  }
1343
1344  private LongMath() {}
1345}