123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535 |
- #include <limits.h>
- #include <math.h>
- #define SET_INVALID 0x01000000UL
- typedef union
- {
- struct {
- #if defined(__BIG_ENDIAN__)
- unsigned long int hi;
- unsigned long int lo;
- #else
- unsigned long int lo;
- unsigned long int hi;
- #endif
- } words;
- double dbl;
- } DblInHex;
- static const unsigned long int signMask = 0x80000000ul;
- static const double twoTo52 = 4503599627370496.0;
- static const double doubleToLong = 4503603922337792.0;
- static const DblInHex Huge = {{ 0x7FF00000, 0x00000000 }};
- static const DblInHex TOWARDZERO = {{ 0x00000000, 0x00000001 }};
-
- double nearbyint ( double x )
- {
- double y;
- double OldEnvironment;
-
- y = twoTo52;
-
- asm ("mffs %0" : "=f" (OldEnvironment));
- if ( fabs ( x ) >= y )
- return x;
- if ( x < 0 ) y = -y;
- y = ( x + y ) - y;
- if ( y == 0.0 )
- y = copysign ( y, x );
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment ));
- return ( y );
- }
-
- long int rinttol ( double x )
- {
- register double y;
- DblInHex argument, OldEnvironment;
- unsigned long int xHead;
- register long int target;
-
- argument.dbl = x;
- target = ( argument.words.hi < signMask );
- xHead = argument.words.hi & 0x7ffffffful;
-
- if ( target )
- {
- if ( xHead < 0x41dffffful )
- {
- y = ( x + twoTo52 ) - twoTo52;
- argument.dbl = y + doubleToLong;
- return ( ( long ) argument.words.lo );
- }
-
- asm ("mffs %0" : "=f" (OldEnvironment.dbl));
-
- if ( xHead > 0x41dffffful )
- {
- OldEnvironment.words.lo |= SET_INVALID;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- return ( LONG_MAX );
- }
-
- y = ( x + twoTo52 ) - twoTo52;
- if ( y > ( double ) LONG_MAX )
- {
- OldEnvironment.words.lo |= SET_INVALID;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- return ( LONG_MAX );
- }
- argument.dbl = y + doubleToLong;
- return ( ( long ) argument.words.lo );
- }
-
- if ( xHead < 0x41e00000ul )
- {
- y = ( x - twoTo52 ) + twoTo52;
- argument.dbl = y + doubleToLong;
- return ( ( long ) argument.words.lo );
- }
-
- asm ("mffs %0" : "=f" (OldEnvironment.dbl));
-
- if ( xHead > 0x41e00000ul )
- {
- OldEnvironment.words.lo |= SET_INVALID;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- return ( LONG_MIN );
- }
-
- y = ( x - twoTo52 ) + twoTo52;
- if ( y < ( double ) LONG_MIN )
- {
- OldEnvironment.words.lo |= SET_INVALID;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- return ( LONG_MIN );
- }
- argument.dbl = y + doubleToLong;
- return ( ( long ) argument.words.lo );
- }
-
- double round ( double x )
- {
- DblInHex argument, OldEnvironment;
- register double y, z;
- register unsigned long int xHead;
- register long int target;
-
- argument.dbl = x;
- xHead = argument.words.hi & 0x7fffffffUL;
- target = ( argument.words.hi < signMask );
-
- if ( xHead < 0x43300000ul )
- {
- if ( xHead < 0x3ff00000ul )
- {
- asm ("mffs %0" : "=f" (OldEnvironment.dbl));
- if ( xHead < 0x3fe00000ul )
- {
- if ( ( xHead | argument.words.lo ) != 0ul )
- OldEnvironment.words.lo |= 0x02000000ul;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- if ( target )
- return ( 0.0 );
- else
- return ( -0.0 );
- }
- OldEnvironment.words.lo |= 0x02000000ul;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- if ( target )
- return ( 1.0 );
- else
- return ( -1.0 );
- }
- if ( target )
- {
- y = ( x + twoTo52 ) - twoTo52;
- if ( y == x )
- return ( x );
- z = x + 0.5;
- y = ( z + twoTo52 ) - twoTo52;
- if ( y > z )
- return ( y - 1.0 );
- else
- return ( y );
- }
-
- else
- {
- y = ( x - twoTo52 ) + twoTo52;
- if ( y == x )
- return ( x );
- z = x - 0.5;
- y = ( z - twoTo52 ) + twoTo52;
- if ( y < z )
- return ( y + 1.0 );
- else
- return ( y );
- }
- }
- return ( x );
- }
- long int roundtol ( double x )
- {
- register double y, z;
- DblInHex argument, OldEnvironment;
- register unsigned long int xhi;
- register long int target;
- const DblInHex kTZ = {{ 0x0, 0x1 }};
- const DblInHex kUP = {{ 0x0, 0x2 }};
-
- argument.dbl = x;
- xhi = argument.words.hi & 0x7ffffffful;
- target = ( argument.words.hi < signMask );
-
- if ( xhi > 0x41e00000ul )
- {
- asm ("mffs %0" : "=f" (OldEnvironment.dbl));
- OldEnvironment.words.lo |= SET_INVALID;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- if ( target )
- return ( LONG_MAX );
- else
- return ( LONG_MIN );
- }
-
- if ( target )
- {
- if ( x < 2147483647.5 )
- {
- y = ( x + doubleToLong ) - doubleToLong;
- if ( y != x )
- {
- asm ("mffs %0" : "=f" (OldEnvironment.dbl));
- asm ("mtfsf 255,%0" : : "f" ( kTZ.dbl ));
- z = x + 0.5;
- argument.dbl = z + doubleToLong;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- return ( ( long ) argument.words.lo );
- }
-
- argument.dbl = y + doubleToLong;
- return ( ( long ) argument.words.lo );
- }
- asm ("mffs %0" : "=f" (OldEnvironment.dbl));
- OldEnvironment.words.lo |= SET_INVALID;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- return ( LONG_MAX );
- }
- if ( x > -2147483648.5 )
- {
- y = ( x + doubleToLong ) - doubleToLong;
- if ( y != x )
- {
- asm ("mffs %0" : "=f" (OldEnvironment.dbl));
- asm ("mtfsf 255,%0" : : "f" ( kUP.dbl ));
- z = x - 0.5;
- argument.dbl = z + doubleToLong;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- return ( ( long ) argument.words.lo );
- }
-
- argument.dbl = y + doubleToLong;
- return ( ( long ) argument.words.lo );
- }
- asm ("mffs %0" : "=f" (OldEnvironment.dbl));
- OldEnvironment.words.lo |= SET_INVALID;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- return ( LONG_MIN );
- }
-
- double trunc ( double x )
- {
- DblInHex argument,OldEnvironment;
- register double y;
- register unsigned long int xhi;
- register long int target;
-
- argument.dbl = x;
- xhi = argument.words.hi & 0x7fffffffUL;
- target = ( argument.words.hi < signMask );
-
- if ( xhi < 0x43300000ul )
- {
- if ( xhi < 0x3ff00000ul )
- {
- if ( ( xhi | argument.words.lo ) != 0ul )
- {
- asm ("mffs %0" : "=f" (OldEnvironment.dbl));
- OldEnvironment.words.lo |= 0x02000000ul;
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment.dbl ));
- }
- if ( target )
- return ( 0.0 );
- else
- return ( -0.0 );
- }
- if ( target )
- {
- y = ( x + twoTo52 ) - twoTo52;
- if ( y > x )
- return ( y - 1.0 );
- else
- return ( y );
- }
-
- else
- {
- y = ( x - twoTo52 ) + twoTo52;
- if ( y < x )
- return ( y + 1.0 );
- else
- return ( y );
- }
- }
- return ( x );
- }
- double modf ( double x, double *iptr )
- {
- register double OldEnvironment, xtrunc;
- register unsigned long int xHead, signBit;
- DblInHex argument;
-
- argument.dbl = x;
- xHead = argument.words.hi & 0x7ffffffful;
- signBit = ( argument.words.hi & 0x80000000ul );
- if (xHead == 0x7ff81fe0)
- signBit = signBit | 0;
-
- if ( xHead < 0x43300000ul )
- {
- if ( xHead < 0x3ff00000ul )
- {
- argument.words.hi = signBit;
- argument.words.lo = 0ul;
- *iptr = argument.dbl;
- return ( x );
- }
- asm ("mffs %0" : "=f" (OldEnvironment));
-
- asm ("mtfsf 255,%0" : : "f" ( TOWARDZERO.dbl ));
- if ( signBit == 0ul )
- xtrunc = ( x + twoTo52 ) - twoTo52;
- else
- xtrunc = ( x - twoTo52 ) + twoTo52;
-
- asm ("mtfsf 255,%0" : : "f" ( OldEnvironment ));
- *iptr = xtrunc;
- if ( x != xtrunc )
- return ( x - xtrunc );
- else
- {
- argument.words.hi = signBit;
- argument.words.lo = 0ul;
- return ( argument.dbl );
- }
- }
-
- *iptr = x;
- if ( x != x )
- return x;
- else
- {
- argument.words.hi = signBit;
- argument.words.lo = 0ul;
- return ( argument.dbl );
- }
- }
|