TODO update
[libm] / src / math / nextafterf.c
1 /* origin: FreeBSD /usr/src/lib/msun/src/s_nextafterf.c */
2 /*
3  * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
4  */
5 /*
6  * ====================================================
7  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
8  *
9  * Developed at SunPro, a Sun Microsystems, Inc. business.
10  * Permission to use, copy, modify, and distribute this
11  * software is freely granted, provided that this notice
12  * is preserved.
13  * ====================================================
14  */
15
16 #include "libm.h"
17
18 float nextafterf(float x, float y)
19 {
20         volatile float t;
21         int32_t hx,hy,ix,iy;
22
23         GET_FLOAT_WORD(hx, x);
24         GET_FLOAT_WORD(hy, y);
25         ix = hx & 0x7fffffff;  /* |x| */
26         iy = hy & 0x7fffffff;  /* |y| */
27
28         if (ix > 0x7f800000 || /* x is nan */
29             iy > 0x7f800000)   /* y is nan */
30                 return x+y;
31         if (x == y)            /* x == y */
32                 return y;
33         if (ix == 0) {         /* x == 0 */
34                 SET_FLOAT_WORD(x, (hy&0x80000000)|1); /* return +-minsubnormal */
35                 /* raise underflow flag */
36                 t = x*x;
37                 if (t == x)
38                         return t;
39                 return x;
40         }
41         if (hx >= 0) {         /* x > 0 */
42                 if (hx > hy) {             /* x > y, x -= ulp */
43                         hx--;
44                 } else {                   /* x < y, x += ulp */
45                         hx++;
46                 }
47         } else {               /* x < 0 */
48                 if (hy >= 0 || hx > hy) {  /* x < y, x -= ulp */
49                         hx--;
50                 } else {                   /* x > y, x += ulp */
51                         hx++;
52                 }
53         }
54         hy = hx & 0x7f800000;
55         if (hy >= 0x7f800000)  /* overflow */
56                 return x+x;
57         if (hy < 0x00800000) { /* underflow */
58                 /* raise underflow flag */
59                 t = x*x;
60                 if (t != x) {
61                         SET_FLOAT_WORD(y, hx);
62                         return y;
63                 }
64         }
65         SET_FLOAT_WORD(x, hx);
66         return x;
67 }