#ifdef __APPLE_CC__
#if __APPLE_CC__ > 930
#include "math.h"
#include "fenv.h"
#include "fp_private.h"
#include "fenv_private.h"
#define POWER_NAN "37"
#define SET_INVALID 0x01000000
static const double ln2 = 6.931471805599452862e-1;
static const double ln2Tail = 2.319046813846299558e-17;
static const double kRint = 6.7553994410557440000e+15; static const double maxExp = 7.0978271289338397000e+02; static const double min2Exp = -7.4513321910194122000e+02; static const double denormal = 2.9387358770557188000e-39; static const double twoTo128 = 3.402823669209384635e+38; static const double twoTo52 = 4503599627370496.0; static const double oneOverLn2 = 1.4426950408889633000e+00; static const hexdouble huge = HEXDOUBLE(0x7ff00000, 0x00000000);
static const double c2 = -0.5000000000000000042725;
static const double c3 = 0.3333333333296328456505;
static const double c4 = -0.2499999999949632195759;
static const double c5 = 0.2000014541571685319141;
static const double c6 = -0.1666681511199968257150;
static const double d2 = -0.5000000000000000000000; static const double d3 = 0.3333333333333360968212; static const double d4 = -0.2500000000000056849516; static const double d5 = 0.1999999996377879912068; static const double d6 = -0.1666666662009609743117; static const double d7 = 0.1428690115662276259345; static const double d8 = -0.1250122079043555957999;
static const double cc4 = 0.5000000000000000000;
static const double cc3 = 0.16666666666663809522820;
static const double cc2 = 0.04166666666666111110939;
static const double cc1 = 0.008333338095239329810170;
static const double cc0 = 0.001388889583333492938381;
extern unsigned long logTable[];
extern unsigned long expTable[];
struct logTableEntry
{
double X;
double G;
hexdouble F;
};
struct expTableEntry
{
double x;
double f;
};
#ifdef notdef
static double PowerInner ( double, double );
static double _NearbyInt ( double x );
double pow ( double base, double exponent )
{
register long int isExpOddInt;
double absBase, result;
hexdouble u, v, OldEnvironment;
register unsigned long int expExp;
v.d = exponent;
u.d = base;
if ( ( ( v.i.hi & 0x000fffff ) | v.i.lo ) == 0 )
{ expExp = v.i.hi & 0xfff00000ul;
if ( expExp == 0x40000000 )
return base*base; else if ( exponent == 0.0 )
return 1.0; else if ( expExp == 0x3ff00000 )
return base; else if ( expExp == 0xbff00000 )
return 1.0/base; }
if ( ( v.i.hi & 0x7ff00000ul ) < 0x7ff00000ul )
{ if ( ( u.i.hi & 0x7ff00000ul ) < 0x7ff00000ul )
{ if ( base > 0.0 )
return PowerInner( base, exponent ); else if ( base < 0.0 )
{
if ( _NearbyInt ( exponent ) != exponent ) {
FEGETENVD( OldEnvironment.d );
OldEnvironment.i.lo |= SET_INVALID;
result = nan ( POWER_NAN );
FESETENVD( OldEnvironment.d );
return result;
}
result = PowerInner( -base, exponent );
if ( _NearbyInt ( 0.5 * exponent ) != 0.5 * exponent ) result = - result;
return ( result );
}
else
{ isExpOddInt = ( ( _NearbyInt( exponent ) == exponent ) ?
( _NearbyInt( 0.5 * exponent ) != 0.5 * exponent ) : 0 );
if ( exponent > 0.0 )
return ( ( isExpOddInt ) ? base : 0.0 );
else return ( ( isExpOddInt ) ? 1.0/base : 1.0/__FABS( base ) );
}
}
else if (base != base)
return base;
else
{ if ( base > 0 )
return ( exponent > 0 ) ? huge.d : +0.0;
else
{
isExpOddInt = ( ( _NearbyInt ( exponent ) == exponent ) ?
( _NearbyInt ( 0.5 * exponent ) != 0.5 * exponent ) : 0 );
return ( exponent > 0 ) ?
( isExpOddInt ? -huge.d : +huge.d ) : ( isExpOddInt ? -0.0 : +0.0 );
}
}
}
else if ( exponent != exponent )
return base + exponent;
else
{ if ( base != base )
return base;
absBase = __FABS( base );
if ( ( exponent > 0 && absBase > 1 ) || ( exponent < 0 && absBase < 1 ) )
return huge.d;
else if ( ( exponent > 0 && absBase < 1 ) || ( exponent < 0 && absBase > 1 ) )
return +0.0;
else return 1.0;
}
}
static double _NearbyInt ( double x )
{
hexdouble xInHex;
double OldEnvironment;
xInHex.d = x;
if ( ( xInHex.i.hi & 0x7ffffffful ) < 0x43300000ul )
{ FEGETENVD( OldEnvironment ); FESETENVD( 0.0 );
if ( xInHex.i.hi & 0x80000000ul ) x = ( ( x - twoTo52 ) + twoTo52 );
else
x = ( ( x + twoTo52 ) - twoTo52 );
FESETENVD( OldEnvironment );
return x;
}
else return ( x );
}
static double PowerInner ( double base, double exponent )
{
hexdouble u, v, OldEnvironment, scale, exp;
register long int i, n;
register unsigned long int head;
register double z, zLow, high, low, zSquared, temp1, temp2, temp3, result, tail,
d, x, xLow, y, yLow, power;
struct logTableEntry *logTablePointer = ( struct logTableEntry * ) logTable;
struct expTableEntry *expTablePointer = ( struct expTableEntry * ) expTable + 177;
u.d = base;
head = u.i.hi;
if ( head >= 0x000ffffful ) d = ( u.i.hi >> 20 ) - 1023.0;
else
{ u.d = base * twoTo128; d = ( u.i.hi >> 20 ) - 1151.0;
}
i = u.i.hi & 0x000fffff;
u.i.hi = i | 0x3ff00000;
if ( ( ( u.i.hi & 0x000fffff ) | u.i.lo ) == 0 ) {
z = d * exponent;
if ( z == _NearbyInt( z ) )
{
if (z > 2097)
n = 2098;
else if (z < -2098)
n = -2099;
else
n = z;
return scalbn( 1.0, n );
}
}
if ( u.i.hi & 0x00080000 )
{ n = 1;
u.d *= 0.5;
i = ( i >> 13 ) & 0x3f; }
else
{ n = 0;
i = ( i >> 12 ) + 64; }
FEGETENVD( OldEnvironment.d ); FESETENVD( 0.0 );
z = ( u.d - logTablePointer[i].X ) * logTablePointer[i].G;
zLow = ( ( u.d - logTablePointer[i].X ) - z * logTablePointer[i].X ) * logTablePointer[i].G;
zSquared = z * z;
if ( n == 0 )
{ v.d = base;
if ( ( v.i.hi == 0x3ff00000 ) && ( v.i.lo == 0x0 ) )
{ FESETENVD( OldEnvironment.d );
return 1.0;
}
temp1 = d8 * zSquared + d6;
temp2 = d7 * zSquared + d5;
temp1 = temp1 * zSquared + d4;
temp2 = temp2 * zSquared + d3;
temp1 = temp1 * z + temp2;
temp2 = temp1 * z + d2;
temp3 = z + logTablePointer[i].F.d;
low = logTablePointer[i].F.d - temp3 + z + ( zLow - z * zLow );
result = ( temp2 * zSquared + low ) + temp3;
tail = temp3 - result + ( temp2 * zSquared + low );
}
else if ( logTablePointer[i].F.i.hi != 0 )
{ low = ln2Tail + zLow;
high = z + logTablePointer[i].F.d;
zLow = logTablePointer[i].F.d - high + z;
temp1 = ln2 + low;
low = ( ln2 - temp1 ) + low;
temp3 = c6 * zSquared + c4;
temp2 = c5 * zSquared + c3;
temp3 = temp3 * zSquared;
temp2 = ( temp2 * z + temp3 ) + c2;
temp3 = high + temp1;
temp1 = ( temp1 - temp3 ) + high;
result = ( ( ( temp2 * zSquared + low ) + temp1 ) + zLow ) + temp3;
tail = temp3 - result + zLow + temp1 + ( temp2 * zSquared + low );
}
else
{ low = ln2Tail + zLow;
temp3 = ln2 + low;
low = ( ln2 - temp3 ) + low;
temp1 = d8 * zSquared + d6;
temp2 = d7 * zSquared + d5;
temp1 = temp1 * zSquared + d4;
temp2 = temp2 * zSquared + d3;
temp1 = temp1 * z + temp2;
temp2 = temp1 * z + d2;
result = ( ( temp2 * zSquared + low ) + z ) + temp3;
tail = temp3 - result + z + ( temp2 * zSquared + low );
}
temp1 = d * ln2;
temp2 = __FMSUB(d, ln2, temp1); z = temp1 + result;
zLow = temp1 - z + result + __FMADD(d, ln2Tail, tail + temp2); temp3 = exponent * zLow;
x = __FMADD(exponent, z, temp3); xLow = __FMSUB(exponent, z, x) + temp3;
if ( base == _NearbyInt(base) && exponent > 0 && exponent == _NearbyInt(exponent))
{
}
else
OldEnvironment.i.lo |= FE_INEXACT;
u.d = x;
if ( ( u.i.hi & 0x7ff00000ul ) < 0x40800000ul )
{ scale.d = 0.0;
exp.d = x * oneOverLn2 + kRint;
scale.i.hi = ( exp.i.lo + 1023 ) << 20;
exp.d -= kRint;
y = x - ln2 * exp.d;
yLow = ln2Tail * exp.d;
u.d = 512.0 * y + kRint;
i = u.i.lo;
d = y - expTablePointer[i].x;
z = d - yLow;
zLow = d - z - yLow + xLow;
zSquared = z * z;
temp1 = cc0 * zSquared + cc2;
temp2 = cc1 * zSquared + cc3;
temp1 = temp1 * zSquared + cc4;
temp2 = temp1 + temp2 * z;
temp1 = temp2 * zSquared + z;
result = scale.d * expTablePointer[i].f;
#if 1
temp2 = __FMUL( result , ( temp1 + ( zLow + zLow * temp1 ) )); #else
temp2 = result * ( temp1 + ( zLow + zLow * temp1 ) );
#endif
result += temp2;
FESETENVD( OldEnvironment.d );
return result;
}
if ( ( x <= maxExp ) && ( x > min2Exp ) )
{
scale.d = 0.0;
exp.d = x * oneOverLn2 + kRint;
if ( x >= 512.0 )
{
power = 2.0;
scale.i.hi = ( exp.i.lo + 1023 - 1 ) << 20;
}
else
{
power = denormal;
scale.i.hi = ( exp.i.lo + 1023 + 128 ) << 20;
}
exp.d -= kRint;
y = x - ln2 * exp.d;
yLow = ln2Tail * exp.d;
u.d = 512.0 * y + kRint;
i = u.i.lo;
d = y - expTablePointer[i].x;
z = d - yLow;
zLow = d - z - yLow + xLow;
zSquared = z * z;
temp1 = cc0 * zSquared + cc2;
temp2 = cc1 * zSquared + cc3;
temp1 = temp1 * zSquared + cc4;
temp2 = temp1 + temp2 * z;
temp1 = temp2 * zSquared + z;
result = scale.d * expTablePointer[i].f;
temp2 = __FMUL( result , ( temp1 + ( zLow + zLow * temp1 ) ) ); u.d = ( result + temp2 ) * power;
}
else if ( x > maxExp )
u.d = huge.d;
else
u.d = +0.0;
i = u.i.hi & 0x7ff00000;
if ( i == 0x7ff00000 )
OldEnvironment.i.lo |= (FE_OVERFLOW | FE_INEXACT); else if ( i == 0x0 )
OldEnvironment.i.lo |= (FE_UNDERFLOW | FE_INEXACT);
FESETENVD( OldEnvironment.d );
return u.d;
}
#else
static double PowerInner ( double, double, hexdouble );
static double _NearbyInt ( double x );
double pow ( double base, double exponent )
{
register long int isExpOddInt;
double absBase, result;
hexdouble u, v, OldEnvironment;
register unsigned long int expExp;
register unsigned long GPR_mant, GPR_exp;
register double FPR_z, FPR_half, FPR_one, FPR_t;
v.d = exponent; u.d = base;
GPR_mant = 0x000ffffful; GPR_exp = 0x7ff00000ul;
FPR_z = 0.0; FPR_one = 1.0;
__ENSURE( FPR_z, FPR_z, FPR_one );
if ( ( ( v.i.hi & GPR_mant ) | v.i.lo ) == 0 )
{ expExp = v.i.hi & 0xfff00000ul;
if ( expExp == 0x40000000 )
return base * base; else if ( exponent == FPR_z )
return FPR_one; else if ( expExp == 0x3ff00000 )
return base; else if ( expExp == 0xbff00000 )
return FPR_one/base; }
if ( ( v.i.hi & GPR_exp ) < GPR_exp )
{ if ( ( u.i.hi & GPR_exp ) < GPR_exp )
{ if ( base > FPR_z )
return PowerInner( base, exponent, u ); else if ( base < FPR_z )
{
if ( _NearbyInt ( exponent ) != exponent ) {
FEGETENVD( OldEnvironment.d );
OldEnvironment.i.lo |= SET_INVALID;
result = nan ( POWER_NAN );
FESETENVD( OldEnvironment.d );
return result;
}
FPR_half = 0.5;
u.i.hi ^= 0x80000000;
result = PowerInner( -base, exponent, u );
FPR_t = FPR_half * exponent; if ( _NearbyInt ( FPR_t ) != FPR_t ) result = - result;
return ( result );
}
else
{ if ( _NearbyInt ( exponent ) == exponent ) {
FPR_half = 0.5;
FPR_t = FPR_half * exponent; isExpOddInt = ( _NearbyInt ( FPR_t ) != FPR_t );
}
else
isExpOddInt = 0;
if ( exponent > FPR_z )
return ( ( isExpOddInt ) ? base : FPR_z );
else return ( ( isExpOddInt ) ? FPR_one/base : FPR_one/__FABS( base ) );
}
}
else if (base != base)
return base;
else
{ if ( base > FPR_z )
return ( exponent > FPR_z ) ? huge.d : FPR_z;
else
{
if ( _NearbyInt ( exponent ) == exponent ) {
FPR_half = 0.5;
FPR_t = FPR_half * exponent; isExpOddInt = ( _NearbyInt ( FPR_t ) != FPR_t );
}
else
isExpOddInt = 0;
return ( exponent > 0 ) ?
( isExpOddInt ? -huge.d : +huge.d ) : ( isExpOddInt ? -FPR_z : FPR_z );
}
}
}
else if ( exponent != exponent )
return base + exponent;
else
{ if ( base != base )
return base;
absBase = __FABS( base );
if ( ( exponent > FPR_z && absBase > FPR_one ) || ( exponent < FPR_z && absBase < FPR_one ) )
return huge.d;
else if ( ( exponent > FPR_z && absBase < FPR_one ) || ( exponent < FPR_z && absBase > FPR_one ) )
return FPR_z;
else return FPR_one;
}
}
static double _NearbyInt ( double x )
{
register double result, FPR_env, FPR_absx, FPR_Two52, FPR_z;
FPR_absx = __FABS( x );
FPR_z = 0.0; FPR_Two52 = twoTo52;
FEGETENVD( FPR_env );
__ENSURE( FPR_z, FPR_Two52, FPR_z );
if ( FPR_absx < FPR_Two52 ) {
FESETENVD( FPR_z );
if ( x < FPR_z )
result = ( ( x - FPR_Two52 ) + FPR_Two52 );
else
result = ( ( x + FPR_Two52 ) - FPR_Two52 );
FESETENVD( FPR_env );
return result;
}
else return ( x );
}
static inline double _NearbyIntDfltEnv ( double x )
{
register double result, FPR_absx, FPR_Two52, FPR_z;
FPR_absx = __FABS( x );
FPR_z = 0.0; FPR_Two52 = twoTo52;
__ENSURE( FPR_z, FPR_Two52, FPR_z );
if ( FPR_absx < FPR_Two52 ) {
if ( x < FPR_z )
result = ( ( x - FPR_Two52 ) + FPR_Two52 );
else
result = ( ( x + FPR_Two52 ) - FPR_Two52 );
return result;
}
else return ( x );
}
static const double kMinNormal = 2.2250738585072014e-308; static const double kMaxNormal = 1.7976931348623157e308;
static const hexdouble kConvULDouble = HEXDOUBLE(0x43300000, 0x00000000);
static double PowerInner ( double base, double exponent, hexdouble u )
{
hexdouble dInHex, scale, exp;
register long int i, n;
register double z, zLow, high, low, zSquared, temp1, temp2, temp3, result, tail,
d, x, xLow, y, yLow, power;
struct logTableEntry *logTablePointer = ( struct logTableEntry * ) logTable;
struct expTableEntry *expTablePointer = ( struct expTableEntry * ) expTable + 177;
register double FPR_env, FPR_z, FPR_half, FPR_one, FPR_512, FPR_t, FPR_G, FPR_X, FPR_ud;
register double FPR_q, FPR_kConvDouble;
register struct logTableEntry *pT;
register struct expTableEntry *pE;
FPR_z = 0.0; FPR_one = 1.0;
FPR_half = 0.5; FPR_512 = 512.0;
FPR_kConvDouble = kConvULDouble.d; dInHex.i.hi = 0x43300000;
if ( base == FPR_one )
return FPR_one;
if ( u.i.hi >= 0x000ffffful ) {
FPR_q = 1023.0;
}
else
{ u.d = __FMUL( base, twoTo128 ); FPR_q = 1151.0;
}
dInHex.i.lo = u.i.hi >> 20;
i = u.i.hi & 0x000fffff;
u.i.hi = i | 0x3ff00000;
FEGETENVD( FPR_env); __ENSURE( FPR_z, FPR_half, FPR_512 ); __ENSURE( FPR_z, FPR_half, FPR_one );
FESETENVD( FPR_z );
if ( 0 == (( u.i.hi & 0x000fffff ) | u.i.lo) ) {
d = dInHex.d; d -= FPR_kConvDouble; d -= FPR_q;
z = __FMUL( d, exponent );
if ( z == _NearbyIntDfltEnv( z ) )
{
if (z > 2097.0)
n = 2098;
else if (z < -2098.0)
n = -2099;
else
n = z;
FESETENVD( FPR_env );
return scalbn( FPR_one, n ); }
}
if ( u.i.hi & 0x00080000 )
{ n = 1;
i = ( i >> 13 ) & 0x3f; FPR_ud = __FMUL( u.d, FPR_half );
u.d = FPR_ud;
}
else
{ n = 0;
i = ( i >> 12 ) + 64; FPR_ud = u.d;
}
pT = &(logTablePointer[i]);
FPR_X = pT->X; FPR_G = pT->G;
FPR_t = FPR_ud - FPR_X;
z = __FMUL( FPR_t, FPR_G );
zSquared = __FMUL( z, z ); FPR_t = __FNMSUB( z, FPR_X, FPR_t );
zLow = __FMUL( FPR_t, FPR_G );
if ( n == 0 )
{ register double FPR_Fd, FPR_d2, FPR_d3, FPR_d4, FPR_d5, FPR_d6, FPR_d7, FPR_d8;
FPR_d8 = d8; FPR_d6 = d6;
FPR_d7 = d7; FPR_d5 = d5;
FPR_Fd = pT->F.d;
temp1 = __FMADD( FPR_d8, zSquared, FPR_d6); temp2 = __FMADD( FPR_d7, zSquared, FPR_d5);
FPR_d4 = d4;
FPR_d3 = d3;
temp3 = z + FPR_Fd; FPR_t = __FNMSUB( z, zLow, zLow );
temp1 = __FMADD( temp1, zSquared, FPR_d4); temp2 = __FMADD( temp2, zSquared, FPR_d3);
FPR_d2 = d2;
temp1 = __FMADD( temp1, z, temp2); low = FPR_Fd - temp3 + z + FPR_t;
temp2 = __FMADD( temp1, z, FPR_d2);
FPR_t = __FMADD( temp2, zSquared, low );
result = FPR_t + temp3;
tail = temp3 - result + FPR_t;
}
else if ( pT->F.i.hi != 0 )
{ register double FPR_Fd, FPR_ln2, FPR_ln2Tail, FPR_c2, FPR_c3, FPR_c4, FPR_c5, FPR_c6;
FPR_c6 = c6; FPR_c4 = c4;
FPR_c5 = c5; FPR_c3 = c3;
temp3 = __FMADD( FPR_c6, zSquared, FPR_c4 ); temp2 = __FMADD( FPR_c5, zSquared, FPR_c3);
FPR_ln2Tail = ln2Tail;
FPR_ln2 = ln2; FPR_Fd = pT->F.d;
low = FPR_ln2Tail + zLow;
high = z + FPR_Fd; temp3 = __FMUL( temp3, zSquared );
FPR_c2 = c2;
zLow = FPR_Fd - high + z; temp1 = FPR_ln2 + low;
temp2 = __FMADD( temp2, z, temp3 ) + FPR_c2; low = ( FPR_ln2 - temp1 ) + low;
temp3 = high + temp1;
temp1 = ( temp1 - temp3 ) + high;
FPR_t = __FMADD( temp2, zSquared, low );
result = ( ( FPR_t + temp1 ) + zLow ) + temp3;
tail = temp3 - result + zLow + temp1 + FPR_t;
}
else
{ register double FPR_ln2, FPR_ln2Tail, FPR_d2, FPR_d3, FPR_d4, FPR_d5, FPR_d6, FPR_d7, FPR_d8;
FPR_d8 = d8; FPR_d6 = d6;
FPR_d7 = d7; FPR_d5 = d5;
temp1 = __FMADD( FPR_d8, zSquared, FPR_d6); temp2 = __FMADD( FPR_d7, zSquared, FPR_d5);
FPR_ln2Tail = ln2Tail;
FPR_ln2 = ln2;
low = FPR_ln2Tail + zLow; __ORI_NOOP;
FPR_d4 = d4; FPR_d3 = d3;
temp1 = __FMADD( temp1, zSquared, FPR_d4); temp2 = __FMADD( temp2, zSquared, FPR_d3);
FPR_d2 = d2;
temp1 = __FMADD( temp1, z, temp2); temp3 = FPR_ln2 + low;
temp2 = __FMADD( temp1, z, FPR_d2); low = ( FPR_ln2 - temp3 ) + low;
FPR_t = __FMADD( temp2, zSquared, low );
result = ( FPR_t + z ) + temp3;
tail = temp3 - result + z + FPR_t;
}
{
register double FPR_ln2, FPR_ln2Tail;
d = dInHex.d; d -= FPR_kConvDouble; d -= FPR_q;
FPR_ln2 = ln2; FPR_ln2Tail = ln2Tail;
temp1 = __FMUL( d, FPR_ln2 );
temp2 = __FMSUB( d, FPR_ln2, temp1 ); z = temp1 + result;
FPR_t = tail + temp2;
zLow = temp1 - z + result + __FMADD(d, FPR_ln2Tail, FPR_t);
temp3 = exponent * zLow;
x = __FMADD( exponent, z, temp3 );
xLow = __FMSUB( exponent, z, x ) + temp3;
}
if ( __FABS( x ) < FPR_512 )
{ register double FPR_ln2, FPR_ln2Tail, FPR_ln2Inv, FPR_kRint, FPR_f;
register double FPR_cc0, FPR_cc1, FPR_cc2, FPR_cc3, FPR_cc4;
FPR_ln2 = ln2; FPR_ln2Tail = ln2Tail;
FPR_ln2Inv = oneOverLn2; FPR_kRint = kRint;
FPR_t = __FMADD( x, FPR_ln2Inv, FPR_kRint );
exp.d = FPR_t;
FPR_t -= FPR_kRint;
y = __FNMSUB( FPR_ln2, FPR_t, x); yLow = __FMUL( FPR_ln2Tail, FPR_t );
u.d = __FMADD( FPR_512, y, FPR_kRint );
scale.i.lo = 0;
scale.i.hi = ( exp.i.lo + 1023 ) << 20;
FPR_cc0 = cc0; FPR_cc2 = cc2;
FPR_cc1 = cc1; FPR_cc3 = cc3;
FPR_cc4 = cc4;
i = u.i.lo;
pE = &(expTablePointer[i]);
d = y - pE->x;
z = d - yLow;
zSquared = __FMUL( z, z ); zLow = d - z - yLow + xLow;
FPR_t = scale.d; FPR_f = pE->f;
temp1 = __FMADD( FPR_cc0, zSquared, FPR_cc2 ); temp2 = __FMADD( FPR_cc1, zSquared, FPR_cc3 );
temp1 = __FMADD( temp1, zSquared, FPR_cc4 );
temp2 = __FMADD( temp2, z, temp1 );
temp1 = __FMADD( temp2, zSquared, z ); result = __FMUL( FPR_t, FPR_f );
temp2 = __FMUL( result , ( temp1 + ( zLow + zLow * temp1 ) ));
result += temp2;
if ( exponent > FPR_z && base == _NearbyIntDfltEnv(base) && exponent == _NearbyIntDfltEnv(exponent))
{
FESETENVD( FPR_env );
}
else
{
FESETENVD( FPR_env );
__PROG_INEXACT ( FPR_ln2 );
}
return result;
}
if ( ( x <= maxExp ) && ( x > min2Exp ) )
{
register double FPR_ln2, FPR_ln2Tail, FPR_ln2Inv, FPR_kRint, FPR_f;
register double FPR_cc0, FPR_cc1, FPR_cc2, FPR_cc3, FPR_cc4;
FPR_ln2 = ln2; FPR_ln2Tail = ln2Tail;
FPR_ln2Inv = oneOverLn2; FPR_kRint = kRint;
FPR_t = __FMADD( x, FPR_ln2Inv, FPR_kRint );
exp.d = FPR_t;
if ( x >= FPR_512 )
{
power = 2.0;
scale.i.lo = 0;
scale.i.hi = ( exp.i.lo + 1023 - 1 ) << 20;
}
else
{
power = denormal;
scale.i.lo = 0;
scale.i.hi = ( exp.i.lo + 1023 + 128 ) << 20;
}
FPR_t -= FPR_kRint;
exp.d = FPR_t;
y = __FNMSUB( FPR_ln2, FPR_t, x); yLow = __FMUL( FPR_ln2Tail, FPR_t );
u.d = __FMADD( FPR_512, y, FPR_kRint );
FPR_cc0 = cc0; FPR_cc2 = cc2;
FPR_cc1 = cc1; FPR_cc3 = cc3;
FPR_cc4 = cc4;
i = u.i.lo;
pE = &(expTablePointer[i]);
d = y - pE->x;
z = d - yLow;
zSquared = __FMUL( z, z ); zLow = d - z - yLow + xLow;
FPR_t = scale.d; FPR_f = pE->f;
temp1 = __FMADD( FPR_cc0, zSquared, FPR_cc2 ); temp2 = __FMADD( FPR_cc1, zSquared, FPR_cc3 );
temp1 = __FMADD( temp1, zSquared, FPR_cc4 );
temp2 = __FMADD( temp2, z, temp1 );
temp1 = __FMADD( temp2, zSquared, z );
result = __FMUL( FPR_t, FPR_f );
temp2 = __FMUL( result , ( temp1 + ( zLow + zLow * temp1 ) ) );
result = ( result + temp2 ) * power;
}
else if ( x > maxExp )
result = huge.d;
else
result = FPR_z;
FESETENVD( FPR_env );
if ( __FABS( result ) == huge.d )
__PROG_OF_INEXACT( kMaxNormal ); else if ( result == FPR_z )
__PROG_UF_INEXACT( kMinNormal ); else
__PROG_INEXACT ( ln2 );
return result;
}
#endif
#ifdef notdef
float powf ( float x, float y )
{
return (float)pow ( x, y );
}
#endif
#else
#error Version gcc-932 or higher required. Compilation terminated.
#endif
#endif