|
@@ -2774,3 +2774,192 @@ static floatx80 addFloatx80Sigs( struct roundingData *roundData, floatx80 a, flo
|
|
|
else {
|
|
|
if ( aExp == 0x7FFF ) {
|
|
|
if ( (bits64) ( ( aSig | bSig )<<1 ) ) {
|
|
|
+ return propagateFloatx80NaN( a, b );
|
|
|
+ }
|
|
|
+ return a;
|
|
|
+ }
|
|
|
+ zSig1 = 0;
|
|
|
+ zSig0 = aSig + bSig;
|
|
|
+ if ( aExp == 0 ) {
|
|
|
+ normalizeFloatx80Subnormal( zSig0, &zExp, &zSig0 );
|
|
|
+ goto roundAndPack;
|
|
|
+ }
|
|
|
+ zExp = aExp;
|
|
|
+ goto shiftRight1;
|
|
|
+ }
|
|
|
+
|
|
|
+ zSig0 = aSig + bSig;
|
|
|
+
|
|
|
+ if ( (sbits64) zSig0 < 0 ) goto roundAndPack;
|
|
|
+ shiftRight1:
|
|
|
+ shift64ExtraRightJamming( zSig0, zSig1, 1, &zSig0, &zSig1 );
|
|
|
+ zSig0 |= LIT64( 0x8000000000000000 );
|
|
|
+ ++zExp;
|
|
|
+ roundAndPack:
|
|
|
+ return
|
|
|
+ roundAndPackFloatx80(
|
|
|
+ roundData, zSign, zExp, zSig0, zSig1 );
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns the result of subtracting the absolute values of the extended
|
|
|
+double-precision floating-point values `a' and `b'. If `zSign' is true,
|
|
|
+the difference is negated before being returned. `zSign' is ignored if the
|
|
|
+result is a NaN. The subtraction is performed according to the IEC/IEEE
|
|
|
+Standard for Binary Floating-point Arithmetic.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+static floatx80 subFloatx80Sigs( struct roundingData *roundData, floatx80 a, floatx80 b, flag zSign )
|
|
|
+{
|
|
|
+ int32 aExp, bExp, zExp;
|
|
|
+ bits64 aSig, bSig, zSig0, zSig1;
|
|
|
+ int32 expDiff;
|
|
|
+ floatx80 z;
|
|
|
+
|
|
|
+ aSig = extractFloatx80Frac( a );
|
|
|
+ aExp = extractFloatx80Exp( a );
|
|
|
+ bSig = extractFloatx80Frac( b );
|
|
|
+ bExp = extractFloatx80Exp( b );
|
|
|
+ expDiff = aExp - bExp;
|
|
|
+ if ( 0 < expDiff ) goto aExpBigger;
|
|
|
+ if ( expDiff < 0 ) goto bExpBigger;
|
|
|
+ if ( aExp == 0x7FFF ) {
|
|
|
+ if ( (bits64) ( ( aSig | bSig )<<1 ) ) {
|
|
|
+ return propagateFloatx80NaN( a, b );
|
|
|
+ }
|
|
|
+ roundData->exception |= float_flag_invalid;
|
|
|
+ z.low = floatx80_default_nan_low;
|
|
|
+ z.high = floatx80_default_nan_high;
|
|
|
+ z.__padding = 0;
|
|
|
+ return z;
|
|
|
+ }
|
|
|
+ if ( aExp == 0 ) {
|
|
|
+ aExp = 1;
|
|
|
+ bExp = 1;
|
|
|
+ }
|
|
|
+ zSig1 = 0;
|
|
|
+ if ( bSig < aSig ) goto aBigger;
|
|
|
+ if ( aSig < bSig ) goto bBigger;
|
|
|
+ return packFloatx80( roundData->mode == float_round_down, 0, 0 );
|
|
|
+ bExpBigger:
|
|
|
+ if ( bExp == 0x7FFF ) {
|
|
|
+ if ( (bits64) ( bSig<<1 ) ) return propagateFloatx80NaN( a, b );
|
|
|
+ return packFloatx80( zSign ^ 1, 0x7FFF, LIT64( 0x8000000000000000 ) );
|
|
|
+ }
|
|
|
+ if ( aExp == 0 ) ++expDiff;
|
|
|
+ shift128RightJamming( aSig, 0, - expDiff, &aSig, &zSig1 );
|
|
|
+ bBigger:
|
|
|
+ sub128( bSig, 0, aSig, zSig1, &zSig0, &zSig1 );
|
|
|
+ zExp = bExp;
|
|
|
+ zSign ^= 1;
|
|
|
+ goto normalizeRoundAndPack;
|
|
|
+ aExpBigger:
|
|
|
+ if ( aExp == 0x7FFF ) {
|
|
|
+ if ( (bits64) ( aSig<<1 ) ) return propagateFloatx80NaN( a, b );
|
|
|
+ return a;
|
|
|
+ }
|
|
|
+ if ( bExp == 0 ) --expDiff;
|
|
|
+ shift128RightJamming( bSig, 0, expDiff, &bSig, &zSig1 );
|
|
|
+ aBigger:
|
|
|
+ sub128( aSig, 0, bSig, zSig1, &zSig0, &zSig1 );
|
|
|
+ zExp = aExp;
|
|
|
+ normalizeRoundAndPack:
|
|
|
+ return
|
|
|
+ normalizeRoundAndPackFloatx80(
|
|
|
+ roundData, zSign, zExp, zSig0, zSig1 );
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns the result of adding the extended double-precision floating-point
|
|
|
+values `a' and `b'. The operation is performed according to the IEC/IEEE
|
|
|
+Standard for Binary Floating-point Arithmetic.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+floatx80 floatx80_add( struct roundingData *roundData, floatx80 a, floatx80 b )
|
|
|
+{
|
|
|
+ flag aSign, bSign;
|
|
|
+
|
|
|
+ aSign = extractFloatx80Sign( a );
|
|
|
+ bSign = extractFloatx80Sign( b );
|
|
|
+ if ( aSign == bSign ) {
|
|
|
+ return addFloatx80Sigs( roundData, a, b, aSign );
|
|
|
+ }
|
|
|
+ else {
|
|
|
+ return subFloatx80Sigs( roundData, a, b, aSign );
|
|
|
+ }
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns the result of subtracting the extended double-precision floating-
|
|
|
+point values `a' and `b'. The operation is performed according to the
|
|
|
+IEC/IEEE Standard for Binary Floating-point Arithmetic.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+floatx80 floatx80_sub( struct roundingData *roundData, floatx80 a, floatx80 b )
|
|
|
+{
|
|
|
+ flag aSign, bSign;
|
|
|
+
|
|
|
+ aSign = extractFloatx80Sign( a );
|
|
|
+ bSign = extractFloatx80Sign( b );
|
|
|
+ if ( aSign == bSign ) {
|
|
|
+ return subFloatx80Sigs( roundData, a, b, aSign );
|
|
|
+ }
|
|
|
+ else {
|
|
|
+ return addFloatx80Sigs( roundData, a, b, aSign );
|
|
|
+ }
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns the result of multiplying the extended double-precision floating-
|
|
|
+point values `a' and `b'. The operation is performed according to the
|
|
|
+IEC/IEEE Standard for Binary Floating-point Arithmetic.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+floatx80 floatx80_mul( struct roundingData *roundData, floatx80 a, floatx80 b )
|
|
|
+{
|
|
|
+ flag aSign, bSign, zSign;
|
|
|
+ int32 aExp, bExp, zExp;
|
|
|
+ bits64 aSig, bSig, zSig0, zSig1;
|
|
|
+ floatx80 z;
|
|
|
+
|
|
|
+ aSig = extractFloatx80Frac( a );
|
|
|
+ aExp = extractFloatx80Exp( a );
|
|
|
+ aSign = extractFloatx80Sign( a );
|
|
|
+ bSig = extractFloatx80Frac( b );
|
|
|
+ bExp = extractFloatx80Exp( b );
|
|
|
+ bSign = extractFloatx80Sign( b );
|
|
|
+ zSign = aSign ^ bSign;
|
|
|
+ if ( aExp == 0x7FFF ) {
|
|
|
+ if ( (bits64) ( aSig<<1 )
|
|
|
+ || ( ( bExp == 0x7FFF ) && (bits64) ( bSig<<1 ) ) ) {
|
|
|
+ return propagateFloatx80NaN( a, b );
|
|
|
+ }
|
|
|
+ if ( ( bExp | bSig ) == 0 ) goto invalid;
|
|
|
+ return packFloatx80( zSign, 0x7FFF, LIT64( 0x8000000000000000 ) );
|
|
|
+ }
|
|
|
+ if ( bExp == 0x7FFF ) {
|
|
|
+ if ( (bits64) ( bSig<<1 ) ) return propagateFloatx80NaN( a, b );
|
|
|
+ if ( ( aExp | aSig ) == 0 ) {
|
|
|
+ invalid:
|
|
|
+ roundData->exception |= float_flag_invalid;
|
|
|
+ z.low = floatx80_default_nan_low;
|
|
|
+ z.high = floatx80_default_nan_high;
|
|
|
+ z.__padding = 0;
|
|
|
+ return z;
|
|
|
+ }
|
|
|
+ return packFloatx80( zSign, 0x7FFF, LIT64( 0x8000000000000000 ) );
|
|
|
+ }
|
|
|
+ if ( aExp == 0 ) {
|
|
|
+ if ( aSig == 0 ) return packFloatx80( zSign, 0, 0 );
|
|
|
+ normalizeFloatx80Subnormal( aSig, &aExp, &aSig );
|
|
|
+ }
|
|
|
+ if ( bExp == 0 ) {
|
|
|
+ if ( bSig == 0 ) return packFloatx80( zSign, 0, 0 );
|