in commons-numbers-fraction/src/main/java/org/apache/commons/numbers/fraction/BigFraction.java [1135:1273]
private long toFloatingPointBits(int exponentLength, int significandLength) {
// Assume the following conditions:
//assert exponentLength >= 1;
//assert exponentLength <= 32;
//assert significandLength >= 1;
//assert significandLength <= 63 - exponentLength;
if (isZero()) {
return 0L;
}
final long sign = (numerator.signum() * denominator.signum()) == -1 ? 1L : 0L;
final BigInteger positiveNumerator = numerator.abs();
final BigInteger positiveDenominator = denominator.abs();
/*
* The most significant 1-bit of a non-zero number is not explicitly
* stored in the significand of an IEEE 754 normalized binary
* floating-point number, so we need to round the value of this fraction
* to (significandLength + 1) bits. In order to do this, we calculate
* the most significant (significandLength + 2) bits, and then, based on
* the least significant of those bits, find out whether we need to
* round up or down.
*
* First, we'll remove all powers of 2 from the denominator because they
* are not relevant for the significand of the prospective binary
* floating-point value.
*/
final int denRightShift = positiveDenominator.getLowestSetBit();
final BigInteger divisor = positiveDenominator.shiftRight(denRightShift);
/*
* Now, we're going to calculate the (significandLength + 2) most
* significant bits of the fraction's value using integer division. To
* guarantee that the quotient of the division has at least
* (significandLength + 2) bits, the bit length of the dividend must
* exceed that of the divisor by at least that amount.
*
* If the denominator has prime factors other than 2, i.e. if the
* divisor was not reduced to 1, an excess of exactly
* (significandLength + 2) bits is sufficient, because the knowledge
* that the fractional part of the precise quotient's binary
* representation does not terminate is enough information to resolve
* cases where the most significant (significandLength + 2) bits alone
* are not conclusive.
*
* Otherwise, the quotient must be calculated exactly and the bit length
* of the numerator can only be reduced as long as no precision is lost
* in the process (meaning it can have powers of 2 removed, like the
* denominator).
*/
int numRightShift = positiveNumerator.bitLength() - divisor.bitLength() - (significandLength + 2);
if (numRightShift > 0 &&
divisor.equals(BigInteger.ONE)) {
numRightShift = Math.min(numRightShift, positiveNumerator.getLowestSetBit());
}
final BigInteger dividend = positiveNumerator.shiftRight(numRightShift);
final BigInteger quotient = dividend.divide(divisor);
int quotRightShift = quotient.bitLength() - (significandLength + 1);
long significand = roundAndRightShift(
quotient,
quotRightShift,
!divisor.equals(BigInteger.ONE)
).longValue();
/*
* If the significand had to be rounded up, this could have caused the
* bit length of the significand to increase by one.
*/
if ((significand & (1L << (significandLength + 1))) != 0) {
significand >>= 1;
quotRightShift++;
}
/*
* Now comes the exponent. The absolute value of this fraction based on
* the current local variables is:
*
* significand * 2^(numRightShift - denRightShift + quotRightShift)
*
* To get the unbiased exponent for the floating-point value, we need to
* add (significandLength) to the above exponent, because all but the
* most significant bit of the significand will be treated as a
* fractional part. To convert the unbiased exponent to a biased
* exponent, we also need to add the exponent bias.
*/
final int exponentBias = (1 << (exponentLength - 1)) - 1;
long exponent = (long) numRightShift - denRightShift + quotRightShift + significandLength + exponentBias;
final long maxExponent = (1L << exponentLength) - 1L; //special exponent for infinities and NaN
if (exponent >= maxExponent) { //infinity
exponent = maxExponent;
significand = 0L;
} else if (exponent > 0) { //normalized number
significand &= -1L >>> (64 - significandLength); //remove implicit leading 1-bit
} else { //smaller than the smallest normalized number
/*
* We need to round the quotient to fewer than
* (significandLength + 1) bits. This must be done with the original
* quotient and not with the current significand, because the loss
* of precision in the previous rounding might cause a rounding of
* the current significand's value to produce a different result
* than a rounding of the original quotient.
*
* So we find out how many high-order bits from the quotient we can
* transfer into the significand. The absolute value of the fraction
* is:
*
* quotient * 2^(numRightShift - denRightShift)
*
* To get the significand, we need to right shift the quotient so
* that the above exponent becomes (1 - exponentBias - significandLength)
* (the unbiased exponent of a subnormal floating-point number is
* defined as equivalent to the minimum unbiased exponent of a
* normalized floating-point number, and (- significandLength)
* because the significand will be treated as the fractional part).
*/
significand = roundAndRightShift(quotient,
(1 - exponentBias - significandLength) - (numRightShift - denRightShift),
!divisor.equals(BigInteger.ONE)).longValue();
exponent = 0L;
/*
* Note: It is possible that an otherwise subnormal number will
* round up to the smallest normal number. However, this special
* case does not need to be treated separately, because the
* overflowing highest-order bit of the significand will then simply
* become the lowest-order bit of the exponent, increasing the
* exponent from 0 to 1 and thus establishing the implicity of the
* leading 1-bit.
*/
}
return (sign << (significandLength + exponentLength)) |
(exponent << significandLength) |
significand;
}