自学教程：C++ GET_FLOAT_WORD函数代码示例

floatacosf(float x){	float z,p,q,r,w,s,c,df;	int32_t hx,ix;	GET_FLOAT_WORD(hx,x);	ix = hx&0x7fffffff;	if(ix==0x3f800000) {		/* |x|==1 */	    if(hx>0) return 0.0;	/* acos(1) = 0  */	    else return pi+(float)2.0*pio2_lo;	/* acos(-1)= pi */	} else if(ix>0x3f800000) {	/* |x| >= 1 */	    return (x-x)/(x-x);		/* acos(|x|>1) is NaN */	}	if(ix<0x3f000000) {	/* |x| < 0.5 */	    if(ix<=0x23000000) return pio2_hi+pio2_lo;/*if|x|<2**-57*/	    z = x*x;	    p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));	    q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));	    r = p/q;	    return pio2_hi - (x - (pio2_lo-x*r));	} else  if (hx<0) {		/* x < -0.5 */	    z = (one+x)*(float)0.5;	    p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));	    q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));	    s = sqrtf(z);	    r = p/q;	    w = r*s-pio2_lo;	    return pi - (float)2.0*(s+w);	} else {			/* x > 0.5 */	    int32_t idf;	    z = (one-x)*(float)0.5;	    s = sqrtf(z);	    df = s;	    GET_FLOAT_WORD(idf,df);	    SET_FLOAT_WORD(df,idf&0xfffff000);	    c  = (z-df*df)/(s+df);	    p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));	    q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));	    r = p/q;	    w = r*s+c;	    return (float)2.0*(df+w);	}}

floatfrexpf(float x, int *eptr){	int32_t hx,ix;	GET_FLOAT_WORD(hx,x);	ix = 0x7fffffff&hx;	*eptr = 0;	if(ix>=0x7f800000||(ix==0)) return x;	/* 0,inf,nan */	if (ix<0x00800000) {		/* subnormal */	    x *= two25;	    GET_FLOAT_WORD(hx,x);	    ix = hx&0x7fffffff;	    *eptr = -25;	}	*eptr += (ix>>23)-126;	hx = (hx&0x807fffff)|0x3f000000;	*(int*)(void*)&x = hx;	return x;}

float __erff(float x){	int32_t hx,ix,i;	float R,S,P,Q,s,y,z,r;	GET_FLOAT_WORD(hx,x);	ix = hx&0x7fffffff;	if(ix>=0x7f800000) {		/* erf(nan)=nan */	    i = ((u_int32_t)hx>>31)<<1;	    return (float)(1-i)+one/x;	/* erf(+-inf)=+-1 */	}

float__ieee754_y0f(float x){	float z, s,c,ss,cc,u,v;	int32_t hx,ix;	GET_FLOAT_WORD(hx,x);        ix = 0x7fffffff&hx;    /* Y0(NaN) is NaN, y0(-inf) is Nan, y0(inf) is 0  */	if(ix>=0x7f800000) return  one/(x+x*x);        if(ix==0) return -one/zero;        if(hx<0) return zero/zero;        if(ix >= 0x40000000) {  /* |x| >= 2.0 */        /* y0(x) = sqrt(2/(pi*x))*(p0(x)*sin(x0)+q0(x)*cos(x0))         * where x0 = x-pi/4         *      Better formula:         *              cos(x0) = cos(x)cos(pi/4)+sin(x)sin(pi/4)         *                      =  1/sqrt(2) * (sin(x) + cos(x))         *              sin(x0) = sin(x)cos(3pi/4)-cos(x)sin(3pi/4)         *                      =  1/sqrt(2) * (sin(x) - cos(x))         * To avoid cancellation, use         *              sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))         * to compute the worse one.         */                s = sinf(x);                c = cosf(x);                ss = s-c;                cc = s+c;	/*	 * j0(x) = 1/sqrt(pi) * (P(0,x)*cc - Q(0,x)*ss) / sqrt(x)	 * y0(x) = 1/sqrt(pi) * (P(0,x)*ss + Q(0,x)*cc) / sqrt(x)	 */                if(ix<0x7f000000) {  /* make sure x+x not overflow */                    z = -cosf(x+x);                    if ((s*c)<zero) cc = z/ss;                    else            ss = z/cc;                }#ifdef DEAD_CODE                if(ix>0x80000000) z = (invsqrtpi*ss)/sqrtf(x);                else#endif		{                    u = pzerof(x); v = qzerof(x);                    z = invsqrtpi*(u*ss+v*cc)/sqrtf(x);                }                return z;	}	if(ix<=0x32000000) {	/* x < 2**-27 */	    return(u00 + tpi*__ieee754_logf(x));	}	z = x*x;	u = u00+z*(u01+z*(u02+z*(u03+z*(u04+z*(u05+z*u06)))));	v = one+z*(v01+z*(v02+z*(v03+z*v04)));	return(u/v + tpi*(__ieee754_j0f(x)*__ieee754_logf(x)));}

float cosf(float x) {  double y;  uint32_t ix;  unsigned n, sign;  GET_FLOAT_WORD(ix, x);  sign = ix >> 31;  ix &= 0x7fffffff;  if (ix <= 0x3f490fda) {  /* |x| ~<= pi/4 */    if (ix < 0x39800000) { /* |x| < 2**-12 */      /* raise inexact if x != 0 */      FORCE_EVAL(x + 0x1p120f);      return 1.0f;    }    return __cosdf(x);  }  if (ix <= 0x407b53d1) { /* |x| ~<= 5*pi/4 */    if (ix > 0x4016cbe3)  /* |x|  ~> 3*pi/4 */      return -__cosdf(sign ? x + c2pio2 : x - c2pio2);    else {      if (sign)        return __sindf(x + c1pio2);      else        return __sindf(c1pio2 - x);    }  }  if (ix <= 0x40e231d5) { /* |x| ~<= 9*pi/4 */    if (ix > 0x40afeddf)  /* |x| ~> 7*pi/4 */      return __cosdf(sign ? x + c4pio2 : x - c4pio2);    else {      if (sign)        return __sindf(-x - c3pio2);      else        return __sindf(x - c3pio2);    }  }  /* cos(Inf or NaN) is NaN */  if (ix >= 0x7f800000)    return x - x;  /* general argument reduction needed */  n = __rem_pio2f(x, &y);  switch (n & 3) {    case 0:      return __cosdf(y);    case 1:      return __sindf(-y);    case 2:      return -__cosdf(y);    default:      return __sindf(y);  }}

float remainderf(float x, float p){	int32_t hx,hp;	uint32_t sx;	float p_half;	GET_FLOAT_WORD(hx, x);	GET_FLOAT_WORD(hp, p);	sx = hx & 0x80000000;	hp &= 0x7fffffff;	hx &= 0x7fffffff;	/* purge off exception values */	if (hp == 0 || hx >= 0x7f800000 || hp > 0x7f800000)  /* p = 0, x not finite, p is NaN */		return (x*p)/(x*p);	if (hp <= 0x7effffff)		x = fmodf(x, p + p);  /* now x < 2p */	if (hx - hp == 0)		return 0.0f*x;	x = fabsf(x);	p = fabsf(p);	if (hp < 0x01000000) {		if (x + x > p) {			x -= p;			if (x + x >= p)				x -= p;		}	} else {		p_half = 0.5f*p;		if (x > p_half) {			x -= p;			if (x >= p_half)				x -= p;		}	}	GET_FLOAT_WORD(hx, x);	if ((hx & 0x7fffffff) == 0)		hx = 0;	SET_FLOAT_WORD(x, hx ^ sx);	return x;}

/* * On the SPU, single precision floating point returns only FP_NORMAL and * FP_ZERO, since FP_NAN, FP_INFINITE, and FP_SUBNORMAL are not * supported, base on the common f_fpclassify.c. */int__fpclassifyf (float x){  __uint32_t w;  GET_FLOAT_WORD(w,x);  if (w == 0x00000000 || w == 0x80000000)    return FP_ZERO;  return FP_NORMAL;}

float atanf(float x){	float_t w,s1,s2,z;	uint32_t ix,sign;	int id;	GET_FLOAT_WORD(ix, x);	sign = ix>>31;	ix &= 0x7fffffff;	if (ix >= 0x4c800000) {  /* if |x| >= 2**26 */		if (isnan(x))			return x;		z = atanhi[3] + 0x1p-120f;		return sign ? -z : z;	}	if (ix < 0x3ee00000) {   /* |x| < 0.4375 */		if (ix < 0x39800000) {  /* |x| < 2**-12 */			if (ix < 0x00800000)				/* raise underflow for subnormal x */				FORCE_EVAL(x*x);			return x;		}		id = -1;	} else {		x = fabsf(x);		if (ix < 0x3f980000) {  /* |x| < 1.1875 */			if (ix < 0x3f300000) {  /*  7/16 <= |x| < 11/16 */				id = 0;				x = (2.0f*x - 1.0f)/(2.0f + x);			} else {                /* 11/16 <= |x| < 19/16 */				id = 1;				x = (x - 1.0f)/(x + 1.0f);			}		} else {			if (ix < 0x401c0000) {  /* |x| < 2.4375 */				id = 2;				x = (x - 1.5f)/(1.0f + 1.5f*x);			} else {                /* 2.4375 <= |x| < 2**26 */				id = 3;				x = -1.0f/x;			}		}	}	/* end of argument reduction */	z = x*x;	w = z*z;	/* break sum from i=0 to 10 aT[i]z**(i+1) into odd and even poly */	s1 = z*(aT[0]+w*(aT[2]+w*aT[4]));	s2 = w*(aT[1]+w*aT[3]);	if (id < 0)		return x - x*(s1+s2);	z = atanhi[id] - ((x*(s1+s2) - atanlo[id]) - x);	return sign ? -z : z;}

int isposf(float x){	__int32_t wx;	GET_FLOAT_WORD(wx, x);	if (wx & 0x80000000)		return (0);	else		return (1);}

float nextafterf(float x, float y) {	int32_t hx, hy, ix, iy;	GET_FLOAT_WORD(hx, x);	GET_FLOAT_WORD(hy, y);	ix = hx & 0x7fffffff;		/* |x| */	iy = hy & 0x7fffffff;		/* |y| */	/* x is nan or y is nan? */	if ((ix > 0x7f800000) || (iy > 0x7f800000))		return x + y;	if (x == y)		return y;	if (ix == 0) { /* x == 0? */		SET_FLOAT_WORD(x, (hy & 0x80000000) | 1);		return x;	}	if (hx >= 0) { /* x > 0 */		if (hx > hy) { /* x > y: x -= ulp */			hx -= 1;		} else { /* x < y: x += ulp */			hx += 1;		}	} else { /* x < 0 */		if (hy >= 0 || hx > hy) { /* x < y: x -= ulp */			hx -= 1;		} else { /* x > y: x += ulp */			hx += 1;		}	}	hy = hx & 0x7f800000;	if (hy >= 0x7f800000) {		x = x + x; /* overflow */		return x; /* overflow */	}	SET_FLOAT_WORD(x, hx);	return x;}

float__ieee754_remainderf(float x, float p){	int32_t hx,hp;	u_int32_t sx;	float p_half;	GET_FLOAT_WORD(hx,x);	GET_FLOAT_WORD(hp,p);	sx = hx&0x80000000;	hp &= 0x7fffffff;	hx &= 0x7fffffff;    /* purge off exception values */	if((hp==0)||			 	/* p = 0 */	  (hx>=0x7f800000)||			/* x not finite */	  ((hp>0x7f800000)))			/* p is NaN */	    return nan_mix_op(x, p, *)/nan_mix_op(x, p, *);	if (hp<=0x7effffff) x = __ieee754_fmodf(x,p+p);	/* now x < 2p */	if ((hx-hp)==0) return zero*x;	x  = fabsf(x);	p  = fabsf(p);	if (hp<0x01000000) {	    if(x+x>p) {		x-=p;		if(x+x>=p) x -= p;	    }	} else {	    p_half = (float)0.5*p;	    if(x>p_half) {		x-=p;		if(x>=p_half) x -= p;	    }	}	GET_FLOAT_WORD(hx,x);	if ((hx&0x7fffffff)==0) hx = 0;	SET_FLOAT_WORD(x,hx^sx);	return x;}

float atanf(float x){	float w,s1,s2,z;	int32_t ix,hx,id;	GET_FLOAT_WORD(hx, x);	ix = hx & 0x7fffffff;	if (ix >= 0x4c800000) {  /* if |x| >= 2**26 */		if (ix > 0x7f800000)  /* NaN */			return x+x;		if (hx > 0)			return  atanhi[3] + *(volatile float *)&atanlo[3];		else			return -atanhi[3] - *(volatile float *)&atanlo[3];	}	if (ix < 0x3ee00000) {   /* |x| < 0.4375 */		if (ix < 0x39800000) {  /* |x| < 2**-12 */			/* raise inexact */			if(huge+x>1.0f)				return x;		}		id = -1;	} else {		x = fabsf(x);		if (ix < 0x3f980000) {  /* |x| < 1.1875 */			if (ix < 0x3f300000) {  /*  7/16 <= |x| < 11/16 */				id = 0;				x = (2.0f*x - 1.0f)/(2.0f + x);			} else {                /* 11/16 <= |x| < 19/16 */				id = 1;				x = (x - 1.0f)/(x + 1.0f);			}		} else {			if (ix < 0x401c0000) {  /* |x| < 2.4375 */				id = 2;				x = (x - 1.5f)/(1.0f + 1.5f*x);			} else {                /* 2.4375 <= |x| < 2**26 */				id = 3;				x = -1.0f/x;			}		}	}	/* end of argument reduction */	z = x*x;	w = z*z;	/* break sum from i=0 to 10 aT[i]z**(i+1) into odd and even poly */	s1 = z*(aT[0]+w*(aT[2]+w*aT[4]));	s2 = w*(aT[1]+w*aT[3]);	if (id < 0)		return x - x*(s1+s2);	z = atanhi[id] - ((x*(s1+s2) - atanlo[id]) - x);	return hx < 0 ? -z : z;}

floatcosf(float x){    float y[2];    int32_t n, hx, ix;    GET_FLOAT_WORD(hx,x);    ix = hx & 0x7fffffff;    if(ix <= 0x3f490fda) {		/* |x| ~<= pi/4 */        if(ix<0x39800000)		/* |x| < 2**-12 */            if(((int)x)==0) return 1.0;	/* 1 with inexact if x != 0 */        return __kernel_cosdf(x);    }    if(ix<=0x407b53d1) {		/* |x| ~<= 5*pi/4 */        if(ix>0x4016cbe3)		/* |x|  ~> 3*pi/4 */            return -__kernel_cosdf(x + (hx > 0 ? -c2pio2 : c2pio2));        else {            if(hx>0)                return __kernel_sindf(c1pio2 - x);            else                return __kernel_sindf(x + c1pio2);        }    }    if(ix<=0x40e231d5) {		/* |x| ~<= 9*pi/4 */        if(ix>0x40afeddf)		/* |x|  ~> 7*pi/4 */            return __kernel_cosdf(x + (hx > 0 ? -c4pio2 : c4pio2));        else {            if(hx>0)                return __kernel_sindf(x - c3pio2);            else                return __kernel_sindf(-c3pio2 - x);        }    }    /* cos(Inf or NaN) is NaN */    else if (ix>=0x7f800000) return x-x;    /* general argument reduction needed */    else {        n = __ieee754_rem_pio2f(x,y);        switch(n&3) {        case 0:            return  __kernel_cosdf((double)y[0]+y[1]);        case 1:            return  __kernel_sindf(-(double)y[0]-y[1]);        case 2:            return -__kernel_cosdf((double)y[0]+y[1]);        default:            return  __kernel_sindf((double)y[0]+y[1]);        }    }}

float j0f(float x){	float z, s,c,ss,cc,r,u,v;	int32_t hx,ix;	GET_FLOAT_WORD(hx, x);	ix = hx & 0x7fffffff;	if (ix >= 0x7f800000)		return 1.0f/(x*x);	x = fabsf(x);	if (ix >= 0x40000000) {  /* |x| >= 2.0 */		s = sinf(x);		c = cosf(x);		ss = s-c;		cc = s+c;		if (ix < 0x7f000000) {  /* make sure x+x does not overflow */			z = -cosf(x+x);			if (s*c < 0.0f)				cc = z/ss;			else				ss = z/cc;		}		/*		 * j0(x) = 1/sqrt(pi) * (P(0,x)*cc - Q(0,x)*ss) / sqrt(x)		 * y0(x) = 1/sqrt(pi) * (P(0,x)*ss + Q(0,x)*cc) / sqrt(x)		 */		if (ix > 0x80000000)			z = (invsqrtpi*cc)/sqrtf(x);		else {			u = pzerof(x);			v = qzerof(x);			z = invsqrtpi*(u*cc-v*ss)/sqrtf(x);		}		return z;	}	if (ix < 0x39000000) {  /* |x| < 2**-13 */		/* raise inexact if x != 0 */		if (huge+x > 1.0f) {			if (ix < 0x32000000)  /* |x| < 2**-27 */				return 1.0f;			return 1.0f - 0.25f*x*x;		}	}	z = x*x;	r =  z*(R02+z*(R03+z*(R04+z*R05)));	s =  1.0f+z*(S01+z*(S02+z*(S03+z*S04)));	if(ix < 0x3F800000) {   /* |x| < 1.00 */		return 1.0f + z*(-0.25f + (r/s));	} else {		u = 0.5f*x;		return (1.0f+u)*(1.0f-u) + z*(r/s);	}}

floatlogbf(float x){        int32_t ix;        GET_FLOAT_WORD(ix,x);        ix &= 0x7fffffff;                       /* high |x| */        if(ix==0) return (float)-1.0/fabsf(x);        if(ix>=0x7f800000) return x*x;        if((ix>>=23)==0)                        /* IEEE 754 logb */                return -126.0;        else                return (float) (ix-127);

float y1f(float x){	float z,s,c,ss,cc,u,v;	int32_t hx,ix;	GET_FLOAT_WORD(hx, x);	ix = 0x7fffffff & hx;	/* if Y1(NaN) is NaN, Y1(-inf) is NaN, Y1(inf) is 0 */	if (ix >= 0x7f800000)		return 1.0f/(x+x*x);	if (ix == 0)		return -1.0f/0.0f;	if (hx < 0)		return 0.0f/0.0f;	if (ix >= 0x40000000) {  /* |x| >= 2.0 */		s = sinf(x);		c = cosf(x);		ss = -s-c;		cc = s-c;		if (ix < 0x7f000000) {  /* make sure x+x not overflow */			z = cosf(x+x);			if (s*c > 0.0f)				cc = z/ss;			else				ss = z/cc;		}		/* y1(x) = sqrt(2/(pi*x))*(p1(x)*sin(x0)+q1(x)*cos(x0))		 * where x0 = x-3pi/4		 *      Better formula:		 *              cos(x0) = cos(x)cos(3pi/4)+sin(x)sin(3pi/4)		 *                      =  1/sqrt(2) * (sin(x) - cos(x))		 *              sin(x0) = sin(x)cos(3pi/4)-cos(x)sin(3pi/4)		 *                      = -1/sqrt(2) * (cos(x) + sin(x))		 * To avoid cancellation, use		 *              sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))		 * to compute the worse one.		 */		if (ix > 0x48000000)			z = (invsqrtpi*ss)/sqrtf(x);		else {			u = ponef(x);			v = qonef(x);			z = invsqrtpi*(u*ss+v*cc)/sqrtf(x);		}		return z;	}	if (ix <= 0x24800000)  /* x < 2**-54 */		return -tpi/x;	z = x*x;	u = U0[0]+z*(U0[1]+z*(U0[2]+z*(U0[3]+z*U0[4])));	v = 1.0f+z*(V0[0]+z*(V0[1]+z*(V0[2]+z*(V0[3]+z*V0[4]))));	return x*(u/v) + tpi*(j1f(x)*logf(x)-1.0f/x);}

float ynf(int n, float x) {    uint32_t ix, ib;    int nm1, sign, i;    float a, b, temp;    GET_FLOAT_WORD(ix, x);    sign = ix >> 31;    ix &= 0x7fffffff;    if (ix > 0x7f800000) /* nan */        return x;    if (sign && ix != 0) /* x < 0 */        return 0 / 0.0f;    if (ix == 0x7f800000)        return 0.0f;    if (n == 0)        return y0f(x);    if (n < 0) {        nm1 = -(n + 1);        sign = n & 1;    } else {        nm1 = n - 1;        sign = 0;    }    if (nm1 == 0)        return sign ? -y1f(x) : y1f(x);    a = y0f(x);    b = y1f(x);    /* quit if b is -inf */    GET_FLOAT_WORD(ib, b);    for (i = 0; i < nm1 && ib != 0xff800000;) {        i++;        temp = b;        b = (2.0f * i / x) * b - a;        GET_FLOAT_WORD(ib, b);        a = temp;    }    return sign ? -b : b;}

float __ieee754_asinf(float x){	float t,w,p,q,c,r,s;	int32_t hx,ix;	GET_FLOAT_WORD(hx,x);	ix = hx&0x7fffffff;	if(ix==0x3f800000) {		/* asin(1)=+-pi/2 with inexact */	    return x*pio2_hi+x*pio2_lo;	} else if(ix> 0x3f800000) {	/* |x|>= 1 */	    return (x-x)/(x-x);		/* asin(|x|>1) is NaN */	} else if (ix<0x3f000000) {	/* |x|<0.5 */	    if(ix<0x32000000) {		/* if |x| < 2**-27 */		if(huge+x>one) return x;/* return x with inexact if x!=0*/	    } else {		t = x*x;		w = t * (p0 + t * (p1 + t * (p2 + t * (p3 + t * p4))));		return x+x*w;	    }	}	/* 1> |x|>= 0.5 */	w = one-fabsf(x);	t = w*0.5f;	p = t * (p0 + t * (p1 + t * (p2 + t * (p3 + t * p4))));	s = __ieee754_sqrtf(t);	if(ix>=0x3F79999A) {	/* if |x| > 0.975 */	    t = pio2_hi-(2.0f*(s+s*p)-pio2_lo);	} else {	    int32_t iw;	    w  = s;	    GET_FLOAT_WORD(iw,w);	    SET_FLOAT_WORD(w,iw&0xfffff000);	    c  = (t-w*w)/(s+w);	    r  = p;	    p  = 2.0f*s*r-(pio2_lo-2.0f*c);	    q  = pio4_hi-2.0f*w;	    t  = pio4_hi-(p-q);	}	if(hx>0) return t; else return -t;}

float__ieee754_log10f(float x){	float f,hfsq,hi,lo,r,y;	int32_t i,k,hx;	GET_FLOAT_WORD(hx,x);	k=0;	if (hx < 0x00800000) {			/* x < 2**-126  */	    if ((hx&0x7fffffff)==0)		return -two25/vzero;		/* log(+-0)=-inf */	    if (hx<0) return (x-x)/zero;	/* log(-#) = NaN */	    k -= 25; x *= two25; /* subnormal number, scale up x */	    GET_FLOAT_WORD(hx,x);	}	if (hx >= 0x7f800000) return x+x;	if (hx == 0x3f800000)	    return zero;			/* log(1) = +0 */	k += (hx>>23)-127;	hx &= 0x007fffff;	i = (hx+(0x4afb0d))&0x800000;	SET_FLOAT_WORD(x,hx|(i^0x3f800000));	/* normalize x or x/2 */	k += (i>>23);	y = (float)k;	f = x - (float)1.0;	hfsq = (float)0.5*f*f;	r = k_log1pf(f);	/* See e_log2f.c and e_log2.c for details. */	if (sizeof(float_t) > sizeof(float))		return (r - hfsq + f) * ((float_t)ivln10lo + ivln10hi) +		    y * ((float_t)log10_2lo + log10_2hi);	hi = f - hfsq;	GET_FLOAT_WORD(hx,hi);	SET_FLOAT_WORD(hi,hx&0xfffff000);	lo = (f - hi) - hfsq + r;	return y*log10_2lo + (lo+hi)*ivln10lo + lo*ivln10hi + hi*ivln10hi +	    y*log10_2hi;}

float__ieee754_y1f(float x){	float z, s,c,ss,cc,u,v;	int32_t hx,ix;	GET_FLOAT_WORD(hx,x);	ix = 0x7fffffff&hx;    /* if Y1(NaN) is NaN, Y1(-inf) is NaN, Y1(inf) is 0 */	if(__builtin_expect(ix>=0x7f800000, 0)) return  one/(x+x*x);	if(__builtin_expect(ix==0, 0))		return -HUGE_VALF+x;  /* -inf and overflow exception.  */	if(__builtin_expect(hx<0, 0)) return zero/(zero*x);	if(ix >= 0x40000000) {  /* |x| >= 2.0 */		SET_RESTORE_ROUNDF (FE_TONEAREST);		__sincosf (x, &s, &c);		ss = -s-c;		cc = s-c;		if(ix<0x7f000000) {  /* make sure x+x not overflow */		    z = __cosf(x+x);		    if ((s*c)>zero) cc = z/ss;		    else            ss = z/cc;		}	/* y1(x) = sqrt(2/(pi*x))*(p1(x)*sin(x0)+q1(x)*cos(x0))	 * where x0 = x-3pi/4	 *      Better formula:	 *              cos(x0) = cos(x)cos(3pi/4)+sin(x)sin(3pi/4)	 *                      =  1/sqrt(2) * (sin(x) - cos(x))	 *              sin(x0) = sin(x)cos(3pi/4)-cos(x)sin(3pi/4)	 *                      = -1/sqrt(2) * (cos(x) + sin(x))	 * To avoid cancellation, use	 *              sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))	 * to compute the worse one.	 */		if(ix>0x48000000) z = (invsqrtpi*ss)/__ieee754_sqrtf(x);		else {		    u = ponef(x); v = qonef(x);		    z = invsqrtpi*(u*ss+v*cc)/__ieee754_sqrtf(x);		}		return z;	}	if(__builtin_expect(ix<=0x33000000, 0)) {    /* x < 2**-25 */	    z = -tpi / x;	    if (__isinff (z))		__set_errno (ERANGE);	    return z;	}	z = x*x;	u = U0[0]+z*(U0[1]+z*(U0[2]+z*(U0[3]+z*U0[4])));	v = one+z*(V0[0]+z*(V0[1]+z*(V0[2]+z*(V0[3]+z*V0[4]))));	return(x*(u/v) + tpi*(__ieee754_j1f(x)*__ieee754_logf(x)-one/x));}

float expf(float x){	float_t hi, lo, c, xx, y;	int k, sign;	uint32_t hx;	GET_FLOAT_WORD(hx, x);	sign = hx >> 31;   /* sign bit of x */	hx &= 0x7fffffff;  /* high word of |x| */	/* special cases */	if (hx >= 0x42aeac50) {  /* if |x| >= -87.33655f or NaN */		if (hx >= 0x42b17218 && !sign) {  /* x >= 88.722839f */			/* overflow */			x *= 0x1p127f;			return x;		}		if (sign) {			/* underflow */			FORCE_EVAL(-0x1p-149f/x);			if (hx >= 0x42cff1b5)  /* x <= -103.972084f */				return 0;		}	}	/* argument reduction */	if (hx > 0x3eb17218) {  /* if |x| > 0.5 ln2 */		if (hx > 0x3f851592)  /* if |x| > 1.5 ln2 */			k = invln2*x + half[sign];		else			k = 1 - sign - sign;		hi = x - k*ln2hi;  /* k*ln2hi is exact here */		lo = k*ln2lo;		x = hi - lo;	} else if (hx > 0x39000000) {  /* |x| > 2**-14 */		k = 0;		hi = x;		lo = 0;	} else {		/* raise inexact */		FORCE_EVAL(0x1p127f + x);		return 1 + x;	}	/* x is now in primary range */	xx = x*x;	c = x - xx*(P1+xx*P2);	y = 1 + (x*c/(2-c) - lo + hi);	if (k == 0)		return y;	return scalbnf(y, k);}

floatsinf(float x){    double y;    int32_t n, hx, ix;    GET_FLOAT_WORD(hx,x);    ix = hx & 0x7fffffff;    if(ix <= 0x3f490fda) {		/* |x| ~<= pi/4 */        if(ix<0x39800000)		/* |x| < 2**-12 */            if(((int)x)==0) return x;	/* x with inexact if x != 0 */        return __kernel_sindf(x);    }    if(ix<=0x407b53d1) {		/* |x| ~<= 5*pi/4 */        if(ix<=0x4016cbe3) {	/* |x| ~<= 3pi/4 */            if(hx>0)                return __kernel_cosdf(x - s1pio2);            else                return -__kernel_cosdf(x + s1pio2);        } else            return __kernel_sindf((hx > 0 ? s2pio2 : -s2pio2) - x);    }    if(ix<=0x40e231d5) {		/* |x| ~<= 9*pi/4 */        if(ix<=0x40afeddf) {	/* |x| ~<= 7*pi/4 */            if(hx>0)                return -__kernel_cosdf(x - s3pio2);            else                return __kernel_cosdf(x + s3pio2);        } else            return __kernel_sindf(x + (hx > 0 ? -s4pio2 : s4pio2));    }    /* sin(Inf or NaN) is NaN */    else if (ix>=0x7f800000) return x-x;    /* general argument reduction needed */    else {        n = __ieee754_rem_pio2f(x,&y);        switch(n&3) {        case 0:            return  __kernel_sindf(y);        case 1:            return  __kernel_cosdf(y);        case 2:            return  __kernel_sindf(-y);        default:            return -__kernel_cosdf(y);        }    }}

float sinf(float x){	double y;	uint32_t ix;	int n, sign;	GET_FLOAT_WORD(ix, x);	sign = ix >> 31;	ix &= 0x7fffffff;	if (ix <= 0x3f490fda) {  /* |x| ~<= pi/4 */		if (ix < 0x39800000) {  /* |x| < 2**-12 */			/* raise inexact if x!=0 and underflow if subnormal */			FORCE_EVAL(ix < 0x00800000 ? x/0x1p120f : x+0x1p120f);			return x;		}		return __sindf(x);	}	if (ix <= 0x407b53d1) {  /* |x| ~<= 5*pi/4 */		if (ix <= 0x4016cbe3) {  /* |x| ~<= 3pi/4 */			if (sign)				return -__cosdf(x + s1pio2);			else				return __cosdf(x - s1pio2);		}		return __sindf(sign ? -(x + s2pio2) : -(x - s2pio2));	}	if (ix <= 0x40e231d5) {  /* |x| ~<= 9*pi/4 */		if (ix <= 0x40afeddf) {  /* |x| ~<= 7*pi/4 */			if (sign)				return __cosdf(x + s3pio2);			else				return -__cosdf(x - s3pio2);		}		return __sindf(sign ? x + s4pio2 : x - s4pio2);	}	/* sin(Inf or NaN) is NaN */	if (ix >= 0x7f800000)		return x - x;	/* general argument reduction needed */	n = __rem_pio2f(x, &y);	switch (n&3) {	case 0: return  __sindf(y);	case 1: return  __cosdf(y);	case 2: return  __sindf(-y);	default:		return -__cosdf(y);	}}

floatexpm1f(float x){    float y,hi,lo,c,t,e,hxs,hfx,r1,twopk;    int32_t k,xsb;    u_int32_t hx;    GET_FLOAT_WORD(hx,x);    xsb = hx&0x80000000;		/* sign bit of x */    hx &= 0x7fffffff;		/* high word of |x| */    /* filter out huge and non-finite argument */    if(hx >= 0x4195b844) {			/* if |x|>=27*ln2 */        if(hx >= 0x42b17218) {		/* if |x|>=88.721... */            if(hx>0x7f800000)                return x+x; 	 /* NaN */            if(hx==0x7f800000)                return (xsb==0)? x:-1.0;/* exp(+-inf)={inf,-1} */            if(x > o_threshold) return huge*huge; /* overflow */        }        if(xsb!=0) { /* x < -27*ln2, return -1.0 with inexact */            if(x+tiny<(float)0.0)	/* raise inexact */                return tiny-one;	/* return -1 */        }    }    /* argument reduction */    if(hx > 0x3eb17218) {		/* if  |x| > 0.5 ln2 */        if(hx < 0x3F851592) {	/* and |x| < 1.5 ln2 */            if(xsb==0)            {                hi = x - ln2_hi;                lo =  ln2_lo;                k =  1;            }            else            {                hi = x + ln2_hi;                lo = -ln2_lo;                k = -1;            }        } else {            k  = invln2*x+((xsb==0)?(float)0.5:(float)-0.5);            t  = k;            hi = x - t*ln2_hi;	/* t*ln2_hi is exact here */            lo = t*ln2_lo;        }        STRICT_ASSIGN(float, x, hi - lo);        c  = (hi-x)-lo;    }    else if(hx < 0x33000000) {  	/* when |x|<2**-25, return x */

float__ieee754_powf(float x, float y){	float z,ax,z_h,z_l,p_h,p_l;	float y1,t1,t2,r,s,t,u,v,w;	int32_t i,j,k,yisint,n;	int32_t hx,hy,ix,iy,is;	GET_FLOAT_WORD(hx,x);	GET_FLOAT_WORD(hy,y);	ix = hx&0x7fffffff;  iy = hy&0x7fffffff;    /* y==zero: x**0 = 1 */	if(iy==0) return one;    /* x==+-1 */	if(x == 1.0) return one;	if(x == -1.0 && isinf(y)) return one;    /* +-NaN return x+y */	if(__builtin_expect(ix > 0x7f800000 ||			    iy > 0x7f800000, 0))		return x+y;    /* determine if y is an odd int when x < 0     * yisint = 0	... y is not an integer     * yisint = 1	... y is an odd int     * yisint = 2	... y is an even int     */	yisint  = 0;	if(hx<0) {	    if(iy>=0x4b800000) yisint = 2; /* even integer y */	    else if(iy>=0x3f800000) {		k = (iy>>23)-0x7f;	   /* exponent */		j = iy>>(23-k);		if((j<<(23-k))==iy) yisint = 2-(j&1);	    }	}

float powf(float x, float y){	float z,ax,z_h,z_l,p_h,p_l;	float y1,t1,t2,r,s,sn,t,u,v,w;	int32_t i,j,k,yisint,n;	int32_t hx,hy,ix,iy,is;	GET_FLOAT_WORD(hx, x);	GET_FLOAT_WORD(hy, y);	ix = hx & 0x7fffffff;	iy = hy & 0x7fffffff;	/* x**0 = 1, even if x is NaN */	if (iy == 0)		return 1.0f;	/* 1**y = 1, even if y is NaN */	if (hx == 0x3f800000)		return 1.0f;	/* NaN if either arg is NaN */	if (ix > 0x7f800000 || iy > 0x7f800000)		return x + y;	/* determine if y is an odd int when x < 0	 * yisint = 0       ... y is not an integer	 * yisint = 1       ... y is an odd int	 * yisint = 2       ... y is an even int	 */	yisint  = 0;	if (hx < 0) {		if (iy >= 0x4b800000)			yisint = 2; /* even integer y */		else if (iy >= 0x3f800000) {			k = (iy>>23) - 0x7f;         /* exponent */			j = iy>>(23-k);			if ((j<<(23-k)) == iy)				yisint = 2 - (j & 1);		}	}

float __hide_kernel_sinf(float x, float y, int iy){	float z,r,v;	int ix;	GET_FLOAT_WORD(ix,x);	ix &= 0x7fffffff;			/* high word of x */	if(ix<0x32000000)			/* |x| < 2**-27 */	   {if((int)x==0) return x;}		/* generate inexact */	z	=  x*x;	v	=  z*x;	r	=  S2+z*(S3+z*(S4+z*(S5+z*S6)));	if(iy==0) return x+v*(S1+z*r);	else      return x-((z*(half*y-v*r)-y)-v*S1);}

voidSINCOSF_FUNC (float x, float *sinx, float *cosx){  int32_t ix;  /* High word of x. */  GET_FLOAT_WORD (ix, x);  /* |x| ~< pi/4 */  ix &= 0x7fffffff;  if (ix <= 0x3f490fd8)    {      *sinx = __kernel_sinf (x, 0.0, 0);      *cosx = __kernel_cosf (x, 0.0);    }  else if (ix>=0x7f800000)    {      /* sin(Inf or NaN) is NaN */      *sinx = *cosx = x - x;      if (ix == 0x7f800000)	__set_errno (EDOM);    }  else    {      /* Argument reduction needed.  */      float y[2];      int n;      n = __ieee754_rem_pio2f (x, y);      switch (n & 3)	{	case 0:	  *sinx = __kernel_sinf (y[0], y[1], 1);	  *cosx = __kernel_cosf (y[0], y[1]);	  break;	case 1:	  *sinx = __kernel_cosf (y[0], y[1]);	  *cosx = -__kernel_sinf (y[0], y[1], 1);	  break;	case 2:	  *sinx = -__kernel_sinf (y[0], y[1], 1);	  *cosx = -__kernel_cosf (y[0], y[1]);	  break;	default:	  *sinx = -__kernel_cosf (y[0], y[1]);	  *cosx = __kernel_sinf (y[0], y[1], 1);	  break;	}    }}