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