1 #include "libm.h"
2 
3 #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
nexttoward(double x,long double y)4 double nexttoward(double x, long double y)
5 {
6 	return nextafter(x, y);
7 }
8 #else
nexttoward(double x,long double y)9 double nexttoward(double x, long double y)
10 {
11 	union {double f; uint64_t i;} ux = {x};
12 	int e;
13 
14 	if (isnan(x) || isnan(y))
15 		return x + y;
16 	if (x == y)
17 		return y;
18 	if (x == 0) {
19 		ux.i = 1;
20 		if (signbit(y))
21 			ux.i |= 1ULL<<63;
22 	} else if (x < y) {
23 		if (signbit(x))
24 			ux.i--;
25 		else
26 			ux.i++;
27 	} else {
28 		if (signbit(x))
29 			ux.i++;
30 		else
31 			ux.i--;
32 	}
33 	e = ux.i>>52 & 0x7ff;
34 	/* raise overflow if ux.f is infinite and x is finite */
35 	if (e == 0x7ff)
36 		FORCE_EVAL(x+x);
37 	/* raise underflow if ux.f is subnormal or zero */
38 	if (e == 0)
39 		FORCE_EVAL(x*x + ux.f*ux.f);
40 	return ux.f;
41 }
42 #endif
43