|
@@ -397,3 +397,186 @@ static float64 roundAndPackFloat64( struct roundingData *roundData, flag zSign,
|
|
|
roundIncrement = 0x200;
|
|
|
if ( ! roundNearestEven ) {
|
|
|
if ( roundingMode == float_round_to_zero ) {
|
|
|
+ roundIncrement = 0;
|
|
|
+ }
|
|
|
+ else {
|
|
|
+ roundIncrement = 0x3FF;
|
|
|
+ if ( zSign ) {
|
|
|
+ if ( roundingMode == float_round_up ) roundIncrement = 0;
|
|
|
+ }
|
|
|
+ else {
|
|
|
+ if ( roundingMode == float_round_down ) roundIncrement = 0;
|
|
|
+ }
|
|
|
+ }
|
|
|
+ }
|
|
|
+ roundBits = zSig & 0x3FF;
|
|
|
+ if ( 0x7FD <= (bits16) zExp ) {
|
|
|
+ if ( ( 0x7FD < zExp )
|
|
|
+ || ( ( zExp == 0x7FD )
|
|
|
+ && ( (sbits64) ( zSig + roundIncrement ) < 0 ) )
|
|
|
+ ) {
|
|
|
+ //register int lr = __builtin_return_address(0);
|
|
|
+ //printk("roundAndPackFloat64 called from 0x%08x\n",lr);
|
|
|
+ roundData->exception |= float_flag_overflow | float_flag_inexact;
|
|
|
+ return packFloat64( zSign, 0x7FF, 0 ) - ( roundIncrement == 0 );
|
|
|
+ }
|
|
|
+ if ( zExp < 0 ) {
|
|
|
+ isTiny =
|
|
|
+ ( float_detect_tininess == float_tininess_before_rounding )
|
|
|
+ || ( zExp < -1 )
|
|
|
+ || ( zSig + roundIncrement < LIT64( 0x8000000000000000 ) );
|
|
|
+ shift64RightJamming( zSig, - zExp, &zSig );
|
|
|
+ zExp = 0;
|
|
|
+ roundBits = zSig & 0x3FF;
|
|
|
+ if ( isTiny && roundBits ) roundData->exception |= float_flag_underflow;
|
|
|
+ }
|
|
|
+ }
|
|
|
+ if ( roundBits ) roundData->exception |= float_flag_inexact;
|
|
|
+ zSig = ( zSig + roundIncrement )>>10;
|
|
|
+ zSig &= ~ ( ( ( roundBits ^ 0x200 ) == 0 ) & roundNearestEven );
|
|
|
+ if ( zSig == 0 ) zExp = 0;
|
|
|
+ return packFloat64( zSign, zExp, zSig );
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Takes an abstract floating-point value having sign `zSign', exponent `zExp',
|
|
|
+and significand `zSig', and returns the proper double-precision floating-
|
|
|
+point value corresponding to the abstract input. This routine is just like
|
|
|
+`roundAndPackFloat64' except that `zSig' does not have to be normalized in
|
|
|
+any way. In all cases, `zExp' must be 1 less than the ``true'' floating-
|
|
|
+point exponent.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+static float64
|
|
|
+ normalizeRoundAndPackFloat64( struct roundingData *roundData, flag zSign, int16 zExp, bits64 zSig )
|
|
|
+{
|
|
|
+ int8 shiftCount;
|
|
|
+
|
|
|
+ shiftCount = countLeadingZeros64( zSig ) - 1;
|
|
|
+ return roundAndPackFloat64( roundData, zSign, zExp - shiftCount, zSig<<shiftCount );
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+#ifdef FLOATX80
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns the fraction bits of the extended double-precision floating-point
|
|
|
+value `a'.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+INLINE bits64 extractFloatx80Frac( floatx80 a )
|
|
|
+{
|
|
|
+
|
|
|
+ return a.low;
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns the exponent bits of the extended double-precision floating-point
|
|
|
+value `a'.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+INLINE int32 extractFloatx80Exp( floatx80 a )
|
|
|
+{
|
|
|
+
|
|
|
+ return a.high & 0x7FFF;
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns the sign bit of the extended double-precision floating-point value
|
|
|
+`a'.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+INLINE flag extractFloatx80Sign( floatx80 a )
|
|
|
+{
|
|
|
+
|
|
|
+ return a.high>>15;
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Normalizes the subnormal extended double-precision floating-point value
|
|
|
+represented by the denormalized significand `aSig'. The normalized exponent
|
|
|
+and significand are stored at the locations pointed to by `zExpPtr' and
|
|
|
+`zSigPtr', respectively.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+static void
|
|
|
+ normalizeFloatx80Subnormal( bits64 aSig, int32 *zExpPtr, bits64 *zSigPtr )
|
|
|
+{
|
|
|
+ int8 shiftCount;
|
|
|
+
|
|
|
+ shiftCount = countLeadingZeros64( aSig );
|
|
|
+ *zSigPtr = aSig<<shiftCount;
|
|
|
+ *zExpPtr = 1 - shiftCount;
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Packs the sign `zSign', exponent `zExp', and significand `zSig' into an
|
|
|
+extended double-precision floating-point value, returning the result.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+INLINE floatx80 packFloatx80( flag zSign, int32 zExp, bits64 zSig )
|
|
|
+{
|
|
|
+ floatx80 z;
|
|
|
+
|
|
|
+ z.low = zSig;
|
|
|
+ z.high = ( ( (bits16) zSign )<<15 ) + zExp;
|
|
|
+ z.__padding = 0;
|
|
|
+ return z;
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Takes an abstract floating-point value having sign `zSign', exponent `zExp',
|
|
|
+and extended significand formed by the concatenation of `zSig0' and `zSig1',
|
|
|
+and returns the proper extended double-precision floating-point value
|
|
|
+corresponding to the abstract input. Ordinarily, the abstract value is
|
|
|
+rounded and packed into the extended double-precision format, with the
|
|
|
+inexact exception raised if the abstract input cannot be represented
|
|
|
+exactly. If the abstract value is too large, however, the overflow and
|
|
|
+inexact exceptions are raised and an infinity or maximal finite value is
|
|
|
+returned. If the abstract value is too small, the input value is rounded to
|
|
|
+a subnormal number, and the underflow and inexact exceptions are raised if
|
|
|
+the abstract input cannot be represented exactly as a subnormal extended
|
|
|
+double-precision floating-point number.
|
|
|
+ If `roundingPrecision' is 32 or 64, the result is rounded to the same
|
|
|
+number of bits as single or double precision, respectively. Otherwise, the
|
|
|
+result is rounded to the full precision of the extended double-precision
|
|
|
+format.
|
|
|
+ The input significand must be normalized or smaller. If the input
|
|
|
+significand is not normalized, `zExp' must be 0; in that case, the result
|
|
|
+returned is a subnormal number, and it must not require rounding. The
|
|
|
+handling of underflow and overflow follows the IEC/IEEE Standard for Binary
|
|
|
+Floating-point Arithmetic.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+static floatx80
|
|
|
+ roundAndPackFloatx80(
|
|
|
+ struct roundingData *roundData, flag zSign, int32 zExp, bits64 zSig0, bits64 zSig1
|
|
|
+ )
|
|
|
+{
|
|
|
+ int8 roundingMode, roundingPrecision;
|
|
|
+ flag roundNearestEven, increment, isTiny;
|
|
|
+ int64 roundIncrement, roundMask, roundBits;
|
|
|
+
|
|
|
+ roundingMode = roundData->mode;
|
|
|
+ roundingPrecision = roundData->precision;
|
|
|
+ roundNearestEven = ( roundingMode == float_round_nearest_even );
|
|
|
+ if ( roundingPrecision == 80 ) goto precision80;
|
|
|
+ if ( roundingPrecision == 64 ) {
|
|
|
+ roundIncrement = LIT64( 0x0000000000000400 );
|
|
|
+ roundMask = LIT64( 0x00000000000007FF );
|
|
|
+ }
|
|
|
+ else if ( roundingPrecision == 32 ) {
|
|
|
+ roundIncrement = LIT64( 0x0000008000000000 );
|