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