1 #include "libm.h"
2 
3 #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
modfl(long double x,long double * iptr)4 long double modfl(long double x, long double* iptr) {
5     double d;
6     long double r;
7 
8     r = modf(x, &d);
9     *iptr = d;
10     return r;
11 }
12 #elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
13 
14 static const long double toint = 1 / LDBL_EPSILON;
15 
modfl(long double x,long double * iptr)16 long double modfl(long double x, long double* iptr) {
17     union ldshape u = {x};
18     int e = (u.i.se & 0x7fff) - 0x3fff;
19     int s = u.i.se >> 15;
20     long double absx;
21     long double y;
22 
23     /* no fractional part */
24     if (e >= LDBL_MANT_DIG - 1) {
25         *iptr = x;
26         if (isnan(x))
27             return x;
28         return s ? -0.0 : 0.0;
29     }
30 
31     /* no integral part*/
32     if (e < 0) {
33         *iptr = s ? -0.0 : 0.0;
34         return x;
35     }
36 
37     /* raises spurious inexact */
38     absx = s ? -x : x;
39     y = absx + toint - toint - absx;
40     if (y == 0) {
41         *iptr = x;
42         return s ? -0.0 : 0.0;
43     }
44     if (y > 0)
45         y -= 1;
46     if (s)
47         y = -y;
48     *iptr = x + y;
49     return -y;
50 }
51 #endif
52