|
@@ -1431,3 +1431,94 @@ float32 float32_rem( struct roundingData *roundData, float32 a, float32 b )
|
|
|
aSig -= bSig;
|
|
|
} while ( 0 <= (sbits32) aSig );
|
|
|
sigMean = aSig + alternateASig;
|
|
|
+ if ( ( sigMean < 0 ) || ( ( sigMean == 0 ) && ( q & 1 ) ) ) {
|
|
|
+ aSig = alternateASig;
|
|
|
+ }
|
|
|
+ zSign = ( (sbits32) aSig < 0 );
|
|
|
+ if ( zSign ) aSig = - aSig;
|
|
|
+ return normalizeRoundAndPackFloat32( roundData, aSign ^ zSign, bExp, aSig );
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns the square root of the single-precision floating-point value `a'.
|
|
|
+The operation is performed according to the IEC/IEEE Standard for Binary
|
|
|
+Floating-point Arithmetic.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+float32 float32_sqrt( struct roundingData *roundData, float32 a )
|
|
|
+{
|
|
|
+ flag aSign;
|
|
|
+ int16 aExp, zExp;
|
|
|
+ bits32 aSig, zSig;
|
|
|
+ bits64 rem, term;
|
|
|
+
|
|
|
+ aSig = extractFloat32Frac( a );
|
|
|
+ aExp = extractFloat32Exp( a );
|
|
|
+ aSign = extractFloat32Sign( a );
|
|
|
+ if ( aExp == 0xFF ) {
|
|
|
+ if ( aSig ) return propagateFloat32NaN( a, 0 );
|
|
|
+ if ( ! aSign ) return a;
|
|
|
+ roundData->exception |= float_flag_invalid;
|
|
|
+ return float32_default_nan;
|
|
|
+ }
|
|
|
+ if ( aSign ) {
|
|
|
+ if ( ( aExp | aSig ) == 0 ) return a;
|
|
|
+ roundData->exception |= float_flag_invalid;
|
|
|
+ return float32_default_nan;
|
|
|
+ }
|
|
|
+ if ( aExp == 0 ) {
|
|
|
+ if ( aSig == 0 ) return 0;
|
|
|
+ normalizeFloat32Subnormal( aSig, &aExp, &aSig );
|
|
|
+ }
|
|
|
+ zExp = ( ( aExp - 0x7F )>>1 ) + 0x7E;
|
|
|
+ aSig = ( aSig | 0x00800000 )<<8;
|
|
|
+ zSig = estimateSqrt32( aExp, aSig ) + 2;
|
|
|
+ if ( ( zSig & 0x7F ) <= 5 ) {
|
|
|
+ if ( zSig < 2 ) {
|
|
|
+ zSig = 0xFFFFFFFF;
|
|
|
+ }
|
|
|
+ else {
|
|
|
+ aSig >>= aExp & 1;
|
|
|
+ term = ( (bits64) zSig ) * zSig;
|
|
|
+ rem = ( ( (bits64) aSig )<<32 ) - term;
|
|
|
+ while ( (sbits64) rem < 0 ) {
|
|
|
+ --zSig;
|
|
|
+ rem += ( ( (bits64) zSig )<<1 ) | 1;
|
|
|
+ }
|
|
|
+ zSig |= ( rem != 0 );
|
|
|
+ }
|
|
|
+ }
|
|
|
+ shift32RightJamming( zSig, 1, &zSig );
|
|
|
+ return roundAndPackFloat32( roundData, 0, zExp, zSig );
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns 1 if the single-precision floating-point value `a' is equal to the
|
|
|
+corresponding value `b', and 0 otherwise. The comparison is performed
|
|
|
+according to the IEC/IEEE Standard for Binary Floating-point Arithmetic.
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+*/
|
|
|
+flag float32_eq( float32 a, float32 b )
|
|
|
+{
|
|
|
+
|
|
|
+ if ( ( ( extractFloat32Exp( a ) == 0xFF ) && extractFloat32Frac( a ) )
|
|
|
+ || ( ( extractFloat32Exp( b ) == 0xFF ) && extractFloat32Frac( b ) )
|
|
|
+ ) {
|
|
|
+ if ( float32_is_signaling_nan( a ) || float32_is_signaling_nan( b ) ) {
|
|
|
+ float_raise( float_flag_invalid );
|
|
|
+ }
|
|
|
+ return 0;
|
|
|
+ }
|
|
|
+ return ( a == b ) || ( (bits32) ( ( a | b )<<1 ) == 0 );
|
|
|
+
|
|
|
+}
|
|
|
+
|
|
|
+/*
|
|
|
+-------------------------------------------------------------------------------
|
|
|
+Returns 1 if the single-precision floating-point value `a' is less than or
|
|
|
+equal to the corresponding value `b', and 0 otherwise. The comparison is
|
|
|
+performed according to the IEC/IEEE Standard for Binary Floating-point
|