001/*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements.  See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License.  You may obtain a copy of the License at
008 *
009 *      http://www.apache.org/licenses/LICENSE-2.0
010 *
011 * Unless required by applicable law or agreed to in writing, software
012 * distributed under the License is distributed on an "AS IS" BASIS,
013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014 * See the License for the specific language governing permissions and
015 * limitations under the License.
016 */
017package org.apache.commons.numbers.fraction;
018
019import java.io.Serializable;
020import java.math.BigDecimal;
021import java.math.BigInteger;
022import java.math.RoundingMode;
023import java.util.Objects;
024import org.apache.commons.numbers.core.NativeOperators;
025
026/**
027 * Representation of a rational number using arbitrary precision.
028 *
029 * <p>The number is expressed as the quotient {@code p/q} of two {@link BigInteger}s,
030 * a numerator {@code p} and a non-zero denominator {@code q}.
031 *
032 * <p>This class is immutable.
033 *
034 * <a href="https://en.wikipedia.org/wiki/Rational_number">Rational number</a>
035 */
036public final class BigFraction
037    extends Number
038    implements Comparable<BigFraction>,
039               NativeOperators<BigFraction>,
040               Serializable {
041    /** A fraction representing "0". */
042    public static final BigFraction ZERO = new BigFraction(BigInteger.ZERO);
043
044    /** A fraction representing "1". */
045    public static final BigFraction ONE = new BigFraction(BigInteger.ONE);
046
047    /** Serializable version identifier. */
048    private static final long serialVersionUID = 20190701L;
049
050    /** The default iterations used for convergence. */
051    private static final int DEFAULT_MAX_ITERATIONS = 100;
052
053    /** Message for non-finite input double argument to factory constructors. */
054    private static final String NOT_FINITE = "Not finite: ";
055
056    /** The overflow limit for conversion from a double (2^31). */
057    private static final long OVERFLOW = 1L << 31;
058
059    /** The numerator of this fraction reduced to lowest terms. */
060    private final BigInteger numerator;
061
062    /** The denominator of this fraction reduced to lowest terms. */
063    private final BigInteger denominator;
064
065    /**
066     * Private constructor: Instances are created using factory methods.
067     *
068     * <p>This constructor should only be invoked when the fraction is known
069     * to be non-zero; otherwise use {@link #ZERO}. This avoids creating
070     * the zero representation {@code 0 / -1}.
071     *
072     * @param num Numerator, must not be {@code null}.
073     * @param den Denominator, must not be {@code null}.
074     * @throws ArithmeticException if the denominator is zero.
075     */
076    private BigFraction(BigInteger num, BigInteger den) {
077        if (den.signum() == 0) {
078            throw new FractionException(FractionException.ERROR_ZERO_DENOMINATOR);
079        }
080
081        // reduce numerator and denominator by greatest common denominator
082        final BigInteger gcd = num.gcd(den);
083        if (BigInteger.ONE.compareTo(gcd) < 0) {
084            numerator = num.divide(gcd);
085            denominator = den.divide(gcd);
086        } else {
087            numerator = num;
088            denominator = den;
089        }
090    }
091
092    /**
093     * Private constructor: Instances are created using factory methods.
094     *
095     * <p>This sets the denominator to 1.
096     *
097     * @param num Numerator (must not be null).
098     */
099    private BigFraction(BigInteger num) {
100        numerator = num;
101        denominator = BigInteger.ONE;
102    }
103
104    /**
105     * Create a fraction given the double value and either the maximum
106     * error allowed or the maximum number of denominator digits.
107     *
108     * <p>
109     * NOTE: This method is called with:
110     * </p>
111     * <ul>
112     *  <li>EITHER a valid epsilon value and the maxDenominator set to
113     *      Integer.MAX_VALUE (that way the maxDenominator has no effect)
114     *  <li>OR a valid maxDenominator value and the epsilon value set to
115     *      zero (that way epsilon only has effect if there is an exact
116     *      match before the maxDenominator value is reached).
117     * </ul>
118     * <p>
119     * It has been done this way so that the same code can be reused for
120     * both scenarios. However this could be confusing to users if it
121     * were part of the public API and this method should therefore remain
122     * PRIVATE.
123     * </p>
124     *
125     * <p>
126     * See JIRA issue ticket MATH-181 for more details:
127     *     https://issues.apache.org/jira/browse/MATH-181
128     * </p>
129     *
130     * @param value Value to convert to a fraction.
131     * @param epsilon Maximum error allowed.
132     * The resulting fraction is within {@code epsilon} of {@code value},
133     * in absolute terms.
134     * @param maxDenominator Maximum denominator value allowed.
135     * @param maxIterations Maximum number of convergents.
136     * @return a new instance.
137     * @throws IllegalArgumentException if the given {@code value} is NaN or infinite.
138     * @throws ArithmeticException if the continued fraction failed to converge.
139     */
140    private static BigFraction from(final double value,
141                                    final double epsilon,
142                                    final int maxDenominator,
143                                    final int maxIterations) {
144        if (!Double.isFinite(value)) {
145            throw new IllegalArgumentException(NOT_FINITE + value);
146        }
147        if (value == 0) {
148            return ZERO;
149        }
150
151        // Remove sign, this is restored at the end.
152        // (Assumes the value is not zero and thus signum(value) is not zero).
153        final double absValue = Math.abs(value);
154        double r0 = absValue;
155        long a0 = (long) Math.floor(r0);
156        if (a0 > OVERFLOW) {
157            throw new FractionException(FractionException.ERROR_CONVERSION_OVERFLOW, value, a0, 1);
158        }
159
160        // check for (almost) integer arguments, which should not go to iterations.
161        if (r0 - a0 <= epsilon) {
162            // Restore the sign.
163            if (value < 0) {
164                a0 = -a0;
165            }
166            return new BigFraction(BigInteger.valueOf(a0));
167        }
168
169        // Support 2^31 as maximum denominator.
170        // This is negative as an integer so convert to long.
171        final long maxDen = Math.abs((long) maxDenominator);
172
173        long p0 = 1;
174        long q0 = 0;
175        long p1 = a0;
176        long q1 = 1;
177
178        long p2;
179        long q2;
180
181        int n = 0;
182        boolean stop = false;
183        do {
184            ++n;
185            final double r1 = 1.0 / (r0 - a0);
186            final long a1 = (long) Math.floor(r1);
187            p2 = (a1 * p1) + p0;
188            q2 = (a1 * q1) + q0;
189            if (Long.compareUnsigned(p2, OVERFLOW) > 0 ||
190                Long.compareUnsigned(q2, OVERFLOW) > 0) {
191                // In maxDenominator mode, fall-back to the previous valid fraction.
192                if (epsilon == 0) {
193                    p2 = p1;
194                    q2 = q1;
195                    break;
196                }
197                throw new FractionException(FractionException.ERROR_CONVERSION_OVERFLOW, value, p2, q2);
198            }
199
200            final double convergent = (double) p2 / (double) q2;
201            if (n < maxIterations &&
202                Math.abs(convergent - absValue) > epsilon &&
203                q2 < maxDen) {
204                p0 = p1;
205                p1 = p2;
206                q0 = q1;
207                q1 = q2;
208                a0 = a1;
209                r0 = r1;
210            } else {
211                stop = true;
212            }
213        } while (!stop);
214
215        if (n >= maxIterations) {
216            throw new FractionException(FractionException.ERROR_CONVERSION, value, maxIterations);
217        }
218
219        // Use p2 / q2 or p1 / q1 if q2 has grown too large in maxDenominator mode
220        long num;
221        long den;
222        if (q2 <= maxDen) {
223            num = p2;
224            den = q2;
225        } else {
226            num = p1;
227            den = q1;
228        }
229
230        // Restore the sign.
231        if (Math.signum(num) * Math.signum(den) != Math.signum(value)) {
232            num = -num;
233        }
234
235        return new BigFraction(BigInteger.valueOf(num),
236                               BigInteger.valueOf(den));
237    }
238
239    /**
240     * Create a fraction given the double value.
241     * <p>
242     * This factory method behaves <em>differently</em> to the method
243     * {@link #from(double, double, int)}. It converts the double value
244     * exactly, considering its internal bits representation. This works for all
245     * values except NaN and infinities and does not requires any loop or
246     * convergence threshold.
247     * </p>
248     * <p>
249     * Since this conversion is exact and since double numbers are sometimes
250     * approximated, the fraction created may seem strange in some cases. For example,
251     * calling {@code from(1.0 / 3.0)} does <em>not</em> create
252     * the fraction \( \frac{1}{3} \), but the fraction \( \frac{6004799503160661}{18014398509481984} \)
253     * because the double number passed to the method is not exactly \( \frac{1}{3} \)
254     * (which cannot be represented exactly in IEEE754).
255     * </p>
256     *
257     * @param value Value to convert to a fraction.
258     * @throws IllegalArgumentException if the given {@code value} is NaN or infinite.
259     * @return a new instance.
260     *
261     * @see #from(double,double,int)
262     */
263    public static BigFraction from(final double value) {
264        if (!Double.isFinite(value)) {
265            throw new IllegalArgumentException(NOT_FINITE + value);
266        }
267        if (value == 0) {
268            return ZERO;
269        }
270
271        final long bits = Double.doubleToLongBits(value);
272        final long sign = bits & 0x8000000000000000L;
273        final long exponent = bits & 0x7ff0000000000000L;
274        final long mantissa = bits & 0x000fffffffffffffL;
275
276        // Compute m and k such that value = m * 2^k
277        long m;
278        int k;
279
280        if (exponent == 0) {
281            // Subnormal number, the effective exponent bias is 1022, not 1023.
282            // Note: mantissa is never zero as that case has been eliminated.
283            m = mantissa;
284            k = -1074;
285        } else {
286            // Normalized number: Add the implicit most significant bit.
287            m = mantissa | 0x0010000000000000L;
288            k = ((int) (exponent >> 52)) - 1075; // Exponent bias is 1023.
289        }
290        if (sign != 0) {
291            m = -m;
292        }
293        while ((m & 0x001ffffffffffffeL) != 0 &&
294               (m & 0x1) == 0) {
295            m >>= 1;
296            ++k;
297        }
298
299        return k < 0 ?
300            new BigFraction(BigInteger.valueOf(m),
301                            BigInteger.ZERO.flipBit(-k)) :
302            new BigFraction(BigInteger.valueOf(m).multiply(BigInteger.ZERO.flipBit(k)),
303                            BigInteger.ONE);
304    }
305
306    /**
307     * Create a fraction given the double value and maximum error allowed.
308     * <p>
309     * References:
310     * <ul>
311     * <li><a href="http://mathworld.wolfram.com/ContinuedFraction.html">
312     * Continued Fraction</a> equations (11) and (22)-(26)</li>
313     * </ul>
314     *
315     * @param value Value to convert to a fraction.
316     * @param epsilon Maximum error allowed. The resulting fraction is within
317     * {@code epsilon} of {@code value}, in absolute terms.
318     * @param maxIterations Maximum number of convergents.
319     * @throws IllegalArgumentException if the given {@code value} is NaN or infinite;
320     * {@code epsilon} is not positive; or {@code maxIterations < 1}.
321     * @throws ArithmeticException if the continued fraction failed to converge.
322     * @return a new instance.
323     */
324    public static BigFraction from(final double value,
325                                   final double epsilon,
326                                   final int maxIterations) {
327        if (maxIterations < 1) {
328            throw new IllegalArgumentException("Max iterations must be strictly positive: " + maxIterations);
329        }
330        if (epsilon >= 0) {
331            return from(value, epsilon, Integer.MIN_VALUE, maxIterations);
332        }
333        throw new IllegalArgumentException("Epsilon must be positive: " + maxIterations);
334    }
335
336    /**
337     * Create a fraction given the double value and maximum denominator.
338     *
339     * <p>
340     * References:
341     * <ul>
342     * <li><a href="http://mathworld.wolfram.com/ContinuedFraction.html">
343     * Continued Fraction</a> equations (11) and (22)-(26)</li>
344     * </ul>
345     *
346     * <p>Note: The magnitude of the {@code maxDenominator} is used allowing use of
347     * {@link Integer#MIN_VALUE} for a supported maximum denominator of 2<sup>31</sup>.
348     *
349     * @param value Value to convert to a fraction.
350     * @param maxDenominator Maximum allowed value for denominator.
351     * @throws IllegalArgumentException if the given {@code value} is NaN or infinite
352     * or {@code maxDenominator} is zero.
353     * @throws ArithmeticException if the continued fraction failed to converge.
354     * @return a new instance.
355     */
356    public static BigFraction from(final double value,
357                                   final int maxDenominator) {
358        if (maxDenominator == 0) {
359            // Re-use the zero denominator message
360            throw new IllegalArgumentException(FractionException.ERROR_ZERO_DENOMINATOR);
361        }
362        return from(value, 0, maxDenominator, DEFAULT_MAX_ITERATIONS);
363    }
364
365    /**
366     * Create a fraction given the numerator. The denominator is {@code 1}.
367     *
368     * @param num
369     *            the numerator.
370     * @return a new instance.
371     */
372    public static BigFraction of(final int num) {
373        if (num == 0) {
374            return ZERO;
375        }
376        return new BigFraction(BigInteger.valueOf(num));
377    }
378
379    /**
380     * Create a fraction given the numerator. The denominator is {@code 1}.
381     *
382     * @param num Numerator.
383     * @return a new instance.
384     */
385    public static BigFraction of(final long num) {
386        if (num == 0) {
387            return ZERO;
388        }
389        return new BigFraction(BigInteger.valueOf(num));
390    }
391
392    /**
393     * Create a fraction given the numerator. The denominator is {@code 1}.
394     *
395     * @param num Numerator.
396     * @return a new instance.
397     * @throws NullPointerException if numerator is null.
398     */
399    public static BigFraction of(final BigInteger num) {
400        Objects.requireNonNull(num, "numerator");
401        if (num.signum() == 0) {
402            return ZERO;
403        }
404        return new BigFraction(num);
405    }
406
407    /**
408     * Create a fraction given the numerator and denominator.
409     * The fraction is reduced to lowest terms.
410     *
411     * @param num Numerator.
412     * @param den Denominator.
413     * @return a new instance.
414     * @throws ArithmeticException if {@code den} is zero.
415     */
416    public static BigFraction of(final int num, final int den) {
417        if (num == 0) {
418            return ZERO;
419        }
420        return new BigFraction(BigInteger.valueOf(num), BigInteger.valueOf(den));
421    }
422
423    /**
424     * Create a fraction given the numerator and denominator.
425     * The fraction is reduced to lowest terms.
426     *
427     * @param num Numerator.
428     * @param den Denominator.
429     * @return a new instance.
430     * @throws ArithmeticException if {@code den} is zero.
431     */
432    public static BigFraction of(final long num, final long den) {
433        if (num == 0) {
434            return ZERO;
435        }
436        return new BigFraction(BigInteger.valueOf(num), BigInteger.valueOf(den));
437    }
438
439    /**
440     * Create a fraction given the numerator and denominator.
441     * The fraction is reduced to lowest terms.
442     *
443     * @param num Numerator.
444     * @param den Denominator.
445     * @return a new instance.
446     * @throws NullPointerException if numerator or denominator are null.
447     * @throws ArithmeticException if the denominator is zero.
448     */
449    public static BigFraction of(final BigInteger num, final BigInteger den) {
450        if (num.signum() == 0) {
451            return ZERO;
452        }
453        return new BigFraction(num, den);
454    }
455
456    /**
457     * Returns a {@code BigFraction} instance representing the specified string {@code s}.
458     *
459     * <p>If {@code s} is {@code null}, then a {@code NullPointerException} is thrown.
460     *
461     * <p>The string must be in a format compatible with that produced by
462     * {@link #toString() BigFraction.toString()}.
463     * The format expects an integer optionally followed by a {@code '/'} character and
464     * and second integer. Leading and trailing spaces are allowed around each numeric part.
465     * Each numeric part is parsed using {@link BigInteger#BigInteger(String)}. The parts
466     * are interpreted as the numerator and optional denominator of the fraction. If absent
467     * the denominator is assumed to be "1".
468     *
469     * <p>Examples of valid strings and the equivalent {@code BigFraction} are shown below:
470     *
471     * <pre>
472     * "0"                 = BigFraction.of(0)
473     * "42"                = BigFraction.of(42)
474     * "0 / 1"             = BigFraction.of(0, 1)
475     * "1 / 3"             = BigFraction.of(1, 3)
476     * "-4 / 13"           = BigFraction.of(-4, 13)</pre>
477     *
478     * <p>Note: The fraction is returned in reduced form and the numerator and denominator
479     * may not match the values in the input string. For this reason the result of
480     * {@code BigFraction.parse(s).toString().equals(s)} may not be {@code true}.
481     *
482     * @param s String representation.
483     * @return an instance.
484     * @throws NullPointerException if the string is null.
485     * @throws NumberFormatException if the string does not contain a parsable fraction.
486     * @see BigInteger#BigInteger(String)
487     * @see #toString()
488     */
489    public static BigFraction parse(String s) {
490        final String stripped = s.replace(",", "");
491        final int slashLoc = stripped.indexOf('/');
492        // if no slash, parse as single number
493        if (slashLoc == -1) {
494            return of(new BigInteger(stripped.trim()));
495        }
496        final BigInteger num = new BigInteger(stripped.substring(0, slashLoc).trim());
497        final BigInteger denom = new BigInteger(stripped.substring(slashLoc + 1).trim());
498        return of(num, denom);
499    }
500
501    @Override
502    public BigFraction zero() {
503        return ZERO;
504    }
505
506    /** {@inheritDoc} */
507    @Override
508    public boolean isZero() {
509        return numerator.signum() == 0;
510    }
511
512    @Override
513    public BigFraction one() {
514        return ONE;
515    }
516
517    /** {@inheritDoc} */
518    @Override
519    public boolean isOne() {
520        return numerator.equals(denominator);
521    }
522
523    /**
524     * Access the numerator as a {@code BigInteger}.
525     *
526     * @return the numerator as a {@code BigInteger}.
527     */
528    public BigInteger getNumerator() {
529        return numerator;
530    }
531
532    /**
533     * Access the numerator as an {@code int}.
534     *
535     * @return the numerator as an {@code int}.
536     */
537    public int getNumeratorAsInt() {
538        return numerator.intValue();
539    }
540
541    /**
542     * Access the numerator as a {@code long}.
543     *
544     * @return the numerator as a {@code long}.
545     */
546    public long getNumeratorAsLong() {
547        return numerator.longValue();
548    }
549
550    /**
551     * Access the denominator as a {@code BigInteger}.
552     *
553     * @return the denominator as a {@code BigInteger}.
554     */
555    public BigInteger getDenominator() {
556        return denominator;
557    }
558
559    /**
560     * Access the denominator as an {@code int}.
561     *
562     * @return the denominator as an {@code int}.
563     */
564    public int getDenominatorAsInt() {
565        return denominator.intValue();
566    }
567
568    /**
569     * Access the denominator as a {@code long}.
570     *
571     * @return the denominator as a {@code long}.
572     */
573    public long getDenominatorAsLong() {
574        return denominator.longValue();
575    }
576
577    /**
578     * Retrieves the sign of this fraction.
579     *
580     * @return -1 if the value is strictly negative, 1 if it is strictly
581     * positive, 0 if it is 0.
582     */
583    public int signum() {
584        return numerator.signum() * denominator.signum();
585    }
586
587    /**
588     * Returns the absolute value of this fraction.
589     *
590     * @return the absolute value.
591     */
592    public BigFraction abs() {
593        return signum() >= 0 ?
594            this :
595            negate();
596    }
597
598    @Override
599    public BigFraction negate() {
600        return new BigFraction(numerator.negate(), denominator);
601    }
602
603    /**
604     * {@inheritDoc}
605     *
606     * <p>Raises an exception if the fraction is equal to zero.
607     *
608     * @throws ArithmeticException if the current numerator is {@code zero}
609     */
610    @Override
611    public BigFraction reciprocal() {
612        return new BigFraction(denominator, numerator);
613    }
614
615    /**
616     * Returns the {@code double} value closest to this fraction.
617     *
618     * @return the fraction as a {@code double}.
619     */
620    @Override
621    public double doubleValue() {
622        return Double.longBitsToDouble(toFloatingPointBits(11, 52));
623    }
624
625    /**
626     * Returns the {@code float} value closest to this fraction.
627     *
628     * @return the fraction as a {@code double}.
629     */
630    @Override
631    public float floatValue() {
632        return Float.intBitsToFloat((int) toFloatingPointBits(8, 23));
633    }
634
635    /**
636     * Returns the whole number part of the fraction.
637     *
638     * @return the largest {@code int} value that is not larger than this fraction.
639     */
640    @Override
641    public int intValue() {
642        return numerator.divide(denominator).intValue();
643    }
644
645    /**
646     * Returns the whole number part of the fraction.
647     *
648     * @return the largest {@code long} value that is not larger than this fraction.
649     */
650    @Override
651    public long longValue() {
652        return numerator.divide(denominator).longValue();
653    }
654
655    /**
656     * Returns the {@code BigDecimal} representation of this fraction.
657     * This calculates the fraction as numerator divided by denominator.
658     *
659     * @return the fraction as a {@code BigDecimal}.
660     * @throws ArithmeticException
661     *             if the exact quotient does not have a terminating decimal
662     *             expansion.
663     * @see BigDecimal
664     */
665    public BigDecimal bigDecimalValue() {
666        return new BigDecimal(numerator).divide(new BigDecimal(denominator));
667    }
668
669    /**
670     * Returns the {@code BigDecimal} representation of this fraction.
671     * This calculates the fraction as numerator divided by denominator
672     * following the passed rounding mode.
673     *
674     * @param roundingMode Rounding mode to apply.
675     * @return the fraction as a {@code BigDecimal}.
676     * @see BigDecimal
677     */
678    public BigDecimal bigDecimalValue(RoundingMode roundingMode) {
679        return new BigDecimal(numerator).divide(new BigDecimal(denominator), roundingMode);
680    }
681
682    /**
683     * Returns the {@code BigDecimal} representation of this fraction.
684     * This calculates the fraction as numerator divided by denominator
685     * following the passed scale and rounding mode.
686     *
687     * @param scale
688     *            scale of the {@code BigDecimal} quotient to be returned.
689     *            see {@link BigDecimal} for more information.
690     * @param roundingMode Rounding mode to apply.
691     * @return the fraction as a {@code BigDecimal}.
692     * @throws ArithmeticException if {@code roundingMode} == {@link RoundingMode#UNNECESSARY} and
693     *      the specified scale is insufficient to represent the result of the division exactly.
694     * @see BigDecimal
695     */
696    public BigDecimal bigDecimalValue(final int scale, RoundingMode roundingMode) {
697        return new BigDecimal(numerator).divide(new BigDecimal(denominator), scale, roundingMode);
698    }
699
700    /**
701     * Adds the specified {@code value} to this fraction, returning
702     * the result in reduced form.
703     *
704     * @param value Value to add.
705     * @return {@code this + value}.
706     */
707    public BigFraction add(final int value) {
708        return add(BigInteger.valueOf(value));
709    }
710
711    /**
712     * Adds the specified {@code value} to this fraction, returning
713     * the result in reduced form.
714     *
715     * @param value Value to add.
716     * @return {@code this + value}.
717     */
718    public BigFraction add(final long value) {
719        return add(BigInteger.valueOf(value));
720    }
721
722    /**
723     * Adds the specified {@code value} to this fraction, returning
724     * the result in reduced form.
725     *
726     * @param value Value to add.
727     * @return {@code this + value}.
728     */
729    public BigFraction add(final BigInteger value) {
730        if (value.signum() == 0) {
731            return this;
732        }
733        if (isZero()) {
734            return of(value);
735        }
736
737        return of(numerator.add(denominator.multiply(value)), denominator);
738    }
739
740    /**
741     * Adds the specified {@code value} to this fraction, returning
742     * the result in reduced form.
743     *
744     * @param value Value to add.
745     * @return {@code this + value}.
746     */
747    @Override
748    public BigFraction add(final BigFraction value) {
749        if (value.isZero()) {
750            return this;
751        }
752        if (isZero()) {
753            return value;
754        }
755
756        final BigInteger num;
757        final BigInteger den;
758
759        if (denominator.equals(value.denominator)) {
760            num = numerator.add(value.numerator);
761            den = denominator;
762        } else {
763            num = (numerator.multiply(value.denominator)).add((value.numerator).multiply(denominator));
764            den = denominator.multiply(value.denominator);
765        }
766
767        if (num.signum() == 0) {
768            return ZERO;
769        }
770
771        return new BigFraction(num, den);
772    }
773
774    /**
775     * Subtracts the specified {@code value} from this fraction, returning
776     * the result in reduced form.
777     *
778     * @param value Value to subtract.
779     * @return {@code this - value}.
780     */
781    public BigFraction subtract(final int value) {
782        return subtract(BigInteger.valueOf(value));
783    }
784
785    /**
786     * Subtracts the specified {@code value} from this fraction, returning
787     * the result in reduced form.
788     *
789     * @param value Value to subtract.
790     * @return {@code this - value}.
791     */
792    public BigFraction subtract(final long value) {
793        return subtract(BigInteger.valueOf(value));
794    }
795
796    /**
797     * Subtracts the specified {@code value} from this fraction, returning
798     * the result in reduced form.
799     *
800     * @param value Value to subtract.
801     * @return {@code this - value}.
802     */
803    public BigFraction subtract(final BigInteger value) {
804        if (value.signum() == 0) {
805            return this;
806        }
807        if (isZero()) {
808            return of(value.negate());
809        }
810
811        return of(numerator.subtract(denominator.multiply(value)), denominator);
812    }
813
814    /**
815     * Subtracts the specified {@code value} from this fraction, returning
816     * the result in reduced form.
817     *
818     * @param value Value to subtract.
819     * @return {@code this - value}.
820     */
821    @Override
822    public BigFraction subtract(final BigFraction value) {
823        if (value.isZero()) {
824            return this;
825        }
826        if (isZero()) {
827            return value.negate();
828        }
829
830        final BigInteger num;
831        final BigInteger den;
832        if (denominator.equals(value.denominator)) {
833            num = numerator.subtract(value.numerator);
834            den = denominator;
835        } else {
836            num = (numerator.multiply(value.denominator)).subtract((value.numerator).multiply(denominator));
837            den = denominator.multiply(value.denominator);
838        }
839
840        if (num.signum() == 0) {
841            return ZERO;
842        }
843
844        return new BigFraction(num, den);
845    }
846
847    /**
848     * Multiply this fraction by the passed {@code value}, returning
849     * the result in reduced form.
850     *
851     * @param value Value to multiply by.
852     * @return {@code this * value}.
853     */
854    @Override
855    public BigFraction multiply(final int value) {
856        if (value == 0 || isZero()) {
857            return ZERO;
858        }
859
860        return multiply(BigInteger.valueOf(value));
861    }
862
863    /**
864     * Multiply this fraction by the passed {@code value}, returning
865     * the result in reduced form.
866     *
867     * @param value Value to multiply by.
868     * @return {@code this * value}.
869     */
870    public BigFraction multiply(final long value) {
871        if (value == 0 || isZero()) {
872            return ZERO;
873        }
874
875        return multiply(BigInteger.valueOf(value));
876    }
877
878    /**
879     * Multiply this fraction by the passed {@code value}, returning
880     * the result in reduced form.
881     *
882     * @param value Value to multiply by.
883     * @return {@code this * value}.
884     */
885    public BigFraction multiply(final BigInteger value) {
886        if (value.signum() == 0 || isZero()) {
887            return ZERO;
888        }
889        return new BigFraction(value.multiply(numerator), denominator);
890    }
891
892    /**
893     * Multiply this fraction by the passed {@code value}, returning
894     * the result in reduced form.
895     *
896     * @param value Value to multiply by.
897     * @return {@code this * value}.
898     */
899    @Override
900    public BigFraction multiply(final BigFraction value) {
901        if (value.isZero() || isZero()) {
902            return ZERO;
903        }
904        return new BigFraction(numerator.multiply(value.numerator),
905                               denominator.multiply(value.denominator));
906    }
907
908    /**
909     * Divide this fraction by the passed {@code value}, returning
910     * the result in reduced form.
911     *
912     * @param value Value to divide by
913     * @return {@code this / value}.
914     * @throws ArithmeticException if the value to divide by is zero
915     */
916    public BigFraction divide(final int value) {
917        return divide(BigInteger.valueOf(value));
918    }
919
920    /**
921     * Divide this fraction by the passed {@code value}, returning
922     * the result in reduced form.
923     *
924     * @param value Value to divide by
925     * @return {@code this / value}.
926     * @throws ArithmeticException if the value to divide by is zero
927     */
928    public BigFraction divide(final long value) {
929        return divide(BigInteger.valueOf(value));
930    }
931
932    /**
933     * Divide this fraction by the passed {@code value}, returning
934     * the result in reduced form.
935     *
936     * @param value Value to divide by
937     * @return {@code this / value}.
938     * @throws ArithmeticException if the value to divide by is zero
939     */
940    public BigFraction divide(final BigInteger value) {
941        if (value.signum() == 0) {
942            throw new FractionException(FractionException.ERROR_DIVIDE_BY_ZERO);
943        }
944        if (isZero()) {
945            return ZERO;
946        }
947        return new BigFraction(numerator, denominator.multiply(value));
948    }
949
950    /**
951     * Divide this fraction by the passed {@code value}, returning
952     * the result in reduced form.
953     *
954     * @param value Value to divide by
955     * @return {@code this / value}.
956     * @throws ArithmeticException if the value to divide by is zero
957     */
958    @Override
959    public BigFraction divide(final BigFraction value) {
960        if (value.isZero()) {
961            throw new FractionException(FractionException.ERROR_DIVIDE_BY_ZERO);
962        }
963        if (isZero()) {
964            return ZERO;
965        }
966        // Multiply by reciprocal
967        return new BigFraction(numerator.multiply(value.denominator),
968                               denominator.multiply(value.numerator));
969    }
970
971    /**
972     * Returns a {@code BigFraction} whose value is
973     * <code>this<sup>exponent</sup></code>, returning the result in reduced form.
974     *
975     * @param exponent exponent to which this {@code BigFraction} is to be raised.
976     * @return <code>this<sup>exponent</sup></code>.
977     * @throws ArithmeticException if the intermediate result would overflow.
978     */
979    @Override
980    public BigFraction pow(final int exponent) {
981        if (exponent == 1) {
982            return this;
983        }
984        if (exponent == 0) {
985            return ONE;
986        }
987        if (isZero()) {
988            if (exponent < 0) {
989                throw new FractionException(FractionException.ERROR_ZERO_DENOMINATOR);
990            }
991            return ZERO;
992        }
993        if (exponent > 0) {
994            return new BigFraction(numerator.pow(exponent),
995                                   denominator.pow(exponent));
996        }
997        if (exponent == -1) {
998            return this.reciprocal();
999        }
1000        if (exponent == Integer.MIN_VALUE) {
1001            // MIN_VALUE can't be negated
1002            return new BigFraction(denominator.pow(Integer.MAX_VALUE).multiply(denominator),
1003                                   numerator.pow(Integer.MAX_VALUE).multiply(numerator));
1004        }
1005        // Note: Raise the BigIntegers to the power and then reduce.
1006        // The supported range for BigInteger is currently
1007        // +/-2^(Integer.MAX_VALUE) exclusive thus larger
1008        // exponents (long, BigInteger) are currently not supported.
1009        return new BigFraction(denominator.pow(-exponent),
1010                               numerator.pow(-exponent));
1011    }
1012
1013    /**
1014     * Returns the {@code String} representing this fraction.
1015     * Uses:
1016     * <ul>
1017     *  <li>{@code "0"} if {@code numerator} is zero.
1018     *  <li>{@code "numerator"} if {@code denominator} is one.
1019     *  <li>{@code "numerator / denominator"} for all other cases.
1020     * </ul>
1021     *
1022     * @return a string representation of the fraction.
1023     */
1024    @Override
1025    public String toString() {
1026        final String str;
1027        if (isZero()) {
1028            str = "0";
1029        } else if (BigInteger.ONE.equals(denominator)) {
1030            str = numerator.toString();
1031        } else {
1032            str = numerator + " / " + denominator;
1033        }
1034        return str;
1035    }
1036
1037    /**
1038     * Compares this object with the specified object for order using the signed magnitude.
1039     *
1040     * @param other {@inheritDoc}
1041     * @return {@inheritDoc}
1042     */
1043    @Override
1044    public int compareTo(final BigFraction other) {
1045        final int lhsSigNum = signum();
1046        final int rhsSigNum = other.signum();
1047
1048        if (lhsSigNum != rhsSigNum) {
1049            return (lhsSigNum > rhsSigNum) ? 1 : -1;
1050        }
1051        // Same sign.
1052        // Avoid a multiply if both fractions are zero
1053        if (lhsSigNum == 0) {
1054            return 0;
1055        }
1056        // Compare absolute magnitude
1057        final BigInteger nOd = numerator.abs().multiply(other.denominator.abs());
1058        final BigInteger dOn = denominator.abs().multiply(other.numerator.abs());
1059        return nOd.compareTo(dOn);
1060    }
1061
1062    /**
1063     * Test for equality with another object. If the other object is a {@code Fraction} then a
1064     * comparison is made of the sign and magnitude; otherwise {@code false} is returned.
1065     *
1066     * @param other {@inheritDoc}
1067     * @return {@inheritDoc}
1068     */
1069    @Override
1070    public boolean equals(final Object other) {
1071        if (this == other) {
1072            return true;
1073        }
1074
1075        if (other instanceof BigFraction) {
1076            // Since fractions are always in lowest terms, numerators and
1077            // denominators can be compared directly for equality.
1078            final BigFraction rhs = (BigFraction) other;
1079            if (signum() == rhs.signum()) {
1080                return numerator.abs().equals(rhs.numerator.abs()) &&
1081                       denominator.abs().equals(rhs.denominator.abs());
1082            }
1083        }
1084
1085        return false;
1086    }
1087
1088    @Override
1089    public int hashCode() {
1090        // Incorporate the sign and absolute values of the numerator and denominator.
1091        // Equivalent to:
1092        // int hash = 1;
1093        // hash = 31 * hash + numerator.abs().hashCode();
1094        // hash = 31 * hash + denominator.abs().hashCode();
1095        // hash = hash * signum()
1096        // Note: BigInteger.hashCode() * BigInteger.signum() == BigInteger.abs().hashCode().
1097        final int numS = numerator.signum();
1098        final int denS = denominator.signum();
1099        return (31 * (31 + numerator.hashCode() * numS) + denominator.hashCode() * denS) * numS * denS;
1100    }
1101
1102    /**
1103     * Calculates the sign bit, the biased exponent and the significand for a
1104     * binary floating-point representation of this {@code BigFraction}
1105     * according to the IEEE 754 standard, and encodes these values into a {@code long}
1106     * variable. The representative bits are arranged adjacent to each other and
1107     * placed at the low-order end of the returned {@code long} value, with the
1108     * least significant bits used for the significand, the next more
1109     * significant bits for the exponent, and next more significant bit for the
1110     * sign.
1111     *
1112     * <p>Warning: The arguments are not validated.
1113     *
1114     * @param exponentLength the number of bits allowed for the exponent; must be
1115     *                       between 1 and 32 (inclusive), and must not be greater
1116     *                       than {@code 63 - significandLength}
1117     * @param significandLength the number of bits allowed for the significand
1118     *                          (excluding the implicit leading 1-bit in
1119     *                          normalized numbers, e.g. 52 for a double-precision
1120     *                          floating-point number); must be between 1 and
1121     *                          {@code 63 - exponentLength} (inclusive)
1122     * @return the bits of an IEEE 754 binary floating-point representation of
1123     * this fraction encoded in a {@code long}, as described above.
1124     */
1125    private long toFloatingPointBits(int exponentLength, int significandLength) {
1126        // Assume the following conditions:
1127        //assert exponentLength >= 1;
1128        //assert exponentLength <= 32;
1129        //assert significandLength >= 1;
1130        //assert significandLength <= 63 - exponentLength;
1131
1132        if (isZero()) {
1133            return 0L;
1134        }
1135
1136        final long sign = (numerator.signum() * denominator.signum()) == -1 ? 1L : 0L;
1137        final BigInteger positiveNumerator = numerator.abs();
1138        final BigInteger positiveDenominator = denominator.abs();
1139
1140        /*
1141         * The most significant 1-bit of a non-zero number is not explicitly
1142         * stored in the significand of an IEEE 754 normalized binary
1143         * floating-point number, so we need to round the value of this fraction
1144         * to (significandLength + 1) bits. In order to do this, we calculate
1145         * the most significant (significandLength + 2) bits, and then, based on
1146         * the least significant of those bits, find out whether we need to
1147         * round up or down.
1148         *
1149         * First, we'll remove all powers of 2 from the denominator because they
1150         * are not relevant for the significand of the prospective binary
1151         * floating-point value.
1152         */
1153        final int denRightShift = positiveDenominator.getLowestSetBit();
1154        final BigInteger divisor = positiveDenominator.shiftRight(denRightShift);
1155
1156        /*
1157         * Now, we're going to calculate the (significandLength + 2) most
1158         * significant bits of the fraction's value using integer division. To
1159         * guarantee that the quotient of the division has at least
1160         * (significandLength + 2) bits, the bit length of the dividend must
1161         * exceed that of the divisor by at least that amount.
1162         *
1163         * If the denominator has prime factors other than 2, i.e. if the
1164         * divisor was not reduced to 1, an excess of exactly
1165         * (significandLength + 2) bits is sufficient, because the knowledge
1166         * that the fractional part of the precise quotient's binary
1167         * representation does not terminate is enough information to resolve
1168         * cases where the most significant (significandLength + 2) bits alone
1169         * are not conclusive.
1170         *
1171         * Otherwise, the quotient must be calculated exactly and the bit length
1172         * of the numerator can only be reduced as long as no precision is lost
1173         * in the process (meaning it can have powers of 2 removed, like the
1174         * denominator).
1175         */
1176        int numRightShift = positiveNumerator.bitLength() - divisor.bitLength() - (significandLength + 2);
1177        if (numRightShift > 0 &&
1178            divisor.equals(BigInteger.ONE)) {
1179            numRightShift = Math.min(numRightShift, positiveNumerator.getLowestSetBit());
1180        }
1181        final BigInteger dividend = positiveNumerator.shiftRight(numRightShift);
1182
1183        final BigInteger quotient = dividend.divide(divisor);
1184
1185        int quotRightShift = quotient.bitLength() - (significandLength + 1);
1186        long significand = roundAndRightShift(
1187                quotient,
1188                quotRightShift,
1189                !divisor.equals(BigInteger.ONE)
1190        ).longValue();
1191
1192        /*
1193         * If the significand had to be rounded up, this could have caused the
1194         * bit length of the significand to increase by one.
1195         */
1196        if ((significand & (1L << (significandLength + 1))) != 0) {
1197            significand >>= 1;
1198            quotRightShift++;
1199        }
1200
1201        /*
1202         * Now comes the exponent. The absolute value of this fraction based on
1203         * the current local variables is:
1204         *
1205         * significand * 2^(numRightShift - denRightShift + quotRightShift)
1206         *
1207         * To get the unbiased exponent for the floating-point value, we need to
1208         * add (significandLength) to the above exponent, because all but the
1209         * most significant bit of the significand will be treated as a
1210         * fractional part. To convert the unbiased exponent to a biased
1211         * exponent, we also need to add the exponent bias.
1212         */
1213        final int exponentBias = (1 << (exponentLength - 1)) - 1;
1214        long exponent = (long) numRightShift - denRightShift + quotRightShift + significandLength + exponentBias;
1215        final long maxExponent = (1L << exponentLength) - 1L; //special exponent for infinities and NaN
1216
1217        if (exponent >= maxExponent) { //infinity
1218            exponent = maxExponent;
1219            significand = 0L;
1220        } else if (exponent > 0) { //normalized number
1221            significand &= -1L >>> (64 - significandLength); //remove implicit leading 1-bit
1222        } else { //smaller than the smallest normalized number
1223            /*
1224             * We need to round the quotient to fewer than
1225             * (significandLength + 1) bits. This must be done with the original
1226             * quotient and not with the current significand, because the loss
1227             * of precision in the previous rounding might cause a rounding of
1228             * the current significand's value to produce a different result
1229             * than a rounding of the original quotient.
1230             *
1231             * So we find out how many high-order bits from the quotient we can
1232             * transfer into the significand. The absolute value of the fraction
1233             * is:
1234             *
1235             * quotient * 2^(numRightShift - denRightShift)
1236             *
1237             * To get the significand, we need to right shift the quotient so
1238             * that the above exponent becomes (1 - exponentBias - significandLength)
1239             * (the unbiased exponent of a subnormal floating-point number is
1240             * defined as equivalent to the minimum unbiased exponent of a
1241             * normalized floating-point number, and (- significandLength)
1242             * because the significand will be treated as the fractional part).
1243             */
1244            significand = roundAndRightShift(quotient,
1245                                             (1 - exponentBias - significandLength) - (numRightShift - denRightShift),
1246                                             !divisor.equals(BigInteger.ONE)).longValue();
1247            exponent = 0L;
1248
1249            /*
1250             * Note: It is possible that an otherwise subnormal number will
1251             * round up to the smallest normal number. However, this special
1252             * case does not need to be treated separately, because the
1253             * overflowing highest-order bit of the significand will then simply
1254             * become the lowest-order bit of the exponent, increasing the
1255             * exponent from 0 to 1 and thus establishing the implicity of the
1256             * leading 1-bit.
1257             */
1258        }
1259
1260        return (sign << (significandLength + exponentLength)) |
1261            (exponent << significandLength) |
1262            significand;
1263    }
1264
1265    /**
1266     * Rounds an integer to the specified power of two (i.e. the minimum number of
1267     * low-order bits that must be zero) and performs a right-shift by this
1268     * amount. The rounding mode applied is round to nearest, with ties rounding
1269     * to even (meaning the prospective least significant bit must be zero). The
1270     * number can optionally be treated as though it contained at
1271     * least one 0-bit and one 1-bit in its fractional part, to influence the result in cases
1272     * that would otherwise be a tie.
1273     * @param value the number to round and right-shift
1274     * @param bits the power of two to which to round; must be positive
1275     * @param hasFractionalBits whether the number should be treated as though
1276     *                          it contained a non-zero fractional part
1277     * @return a {@code BigInteger} as described above
1278     */
1279    private static BigInteger roundAndRightShift(BigInteger value, int bits, boolean hasFractionalBits) {
1280        // Assumptions:
1281        // assert bits > 0
1282
1283        BigInteger result = value.shiftRight(bits);
1284        if (value.testBit(bits - 1) &&
1285            (hasFractionalBits ||
1286             value.getLowestSetBit() < bits - 1 ||
1287             value.testBit(bits))) {
1288            result = result.add(BigInteger.ONE); //round up
1289        }
1290
1291        return result;
1292    }
1293}