(__ieee754_sqrt): Remove unused variables b and n.
[kopensolaris-gnu/glibc.git] / sysdeps / ieee754 / ldbl-96 / s_modfl.c
1 /* s_modfl.c -- long double version of s_modf.c.
2  * Conversion to long double by Ulrich Drepper,
3  * Cygnus Support, drepper@cygnus.com.
4  */
5
6 /*
7  * ====================================================
8  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
9  *
10  * Developed at SunPro, a Sun Microsystems, Inc. business.
11  * Permission to use, copy, modify, and distribute this
12  * software is freely granted, provided that this notice
13  * is preserved.
14  * ====================================================
15  */
16
17 #if defined(LIBM_SCCS) && !defined(lint)
18 static char rcsid[] = "$NetBSD: $";
19 #endif
20
21 /*
22  * modfl(long double x, long double *iptr)
23  * return fraction part of x, and return x's integral part in *iptr.
24  * Method:
25  *      Bit twiddling.
26  *
27  * Exception:
28  *      No exception.
29  */
30
31 #include "math.h"
32 #include "math_private.h"
33
34 #ifdef __STDC__
35 static const long double one = 1.0;
36 #else
37 static long double one = 1.0;
38 #endif
39
40 #ifdef __STDC__
41         long double __modfl(long double x, long double *iptr)
42 #else
43         long double __modfl(x, iptr)
44         long double x,*iptr;
45 #endif
46 {
47         int32_t i0,i1,j0;
48         u_int32_t i,se;
49         GET_LDOUBLE_WORDS(se,i0,i1,x);
50         j0 = (se&0x7fff)-0x3fff;        /* exponent of x */
51         if(j0<32) {                     /* integer part in high x */
52             if(j0<0) {                  /* |x|<1 */
53                 SET_LDOUBLE_WORDS(*iptr,se&0x8000,0,0); /* *iptr = +-0 */
54                 return x;
55             } else {
56                 i = (0x7fffffff)>>j0;
57                 if(((i0&i)|i1)==0) {            /* x is integral */
58                     *iptr = x;
59                     SET_LDOUBLE_WORDS(x,se&0x8000,0,0); /* return +-0 */
60                     return x;
61                 } else {
62                     SET_LDOUBLE_WORDS(*iptr,se,i0&(~i),0);
63                     return x - *iptr;
64                 }
65             }
66         } else if (j0>63) {             /* no fraction part */
67             *iptr = x*one;
68             /* We must handle NaNs separately.  */
69             if (j0 == 0x4000 && ((i0 & 0x7fffffff) | i1))
70               return x*one;
71             SET_LDOUBLE_WORDS(x,se&0x8000,0,0); /* return +-0 */
72             return x;
73         } else {                        /* fraction part in low x */
74             i = ((u_int32_t)(0x7fffffff))>>(j0-32);
75             if((i1&i)==0) {             /* x is integral */
76                 *iptr = x;
77                 SET_LDOUBLE_WORDS(x,se&0x8000,0,0);     /* return +-0 */
78                 return x;
79             } else {
80                 SET_LDOUBLE_WORDS(*iptr,se,i0,i1&(~i));
81                 return x - *iptr;
82             }
83         }
84 }
85 weak_alias (__modfl, modfl)