1*49393c00Smartynas /* @(#)s_nextafter.c 5.1 93/09/24 */
2*49393c00Smartynas /*
3*49393c00Smartynas * ====================================================
4*49393c00Smartynas * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
5*49393c00Smartynas *
6*49393c00Smartynas * Developed at SunPro, a Sun Microsystems, Inc. business.
7*49393c00Smartynas * Permission to use, copy, modify, and distribute this
8*49393c00Smartynas * software is freely granted, provided that this notice
9*49393c00Smartynas * is preserved.
10*49393c00Smartynas * ====================================================
11*49393c00Smartynas */
12*49393c00Smartynas
13*49393c00Smartynas /* IEEE functions
14*49393c00Smartynas * nexttoward(x,y)
15*49393c00Smartynas * return the next machine floating-point number of x in the
16*49393c00Smartynas * direction toward y.
17*49393c00Smartynas * Special cases:
18*49393c00Smartynas */
19*49393c00Smartynas
20*49393c00Smartynas #include <math.h>
21*49393c00Smartynas #include <float.h>
22*49393c00Smartynas
23*49393c00Smartynas #include "math_private.h"
24*49393c00Smartynas
25*49393c00Smartynas double
nexttoward(double x,long double y)26*49393c00Smartynas nexttoward(double x, long double y)
27*49393c00Smartynas {
28*49393c00Smartynas int32_t hx,ix;
29*49393c00Smartynas int64_t hy,iy;
30*49393c00Smartynas u_int32_t lx;
31*49393c00Smartynas u_int64_t ly;
32*49393c00Smartynas
33*49393c00Smartynas EXTRACT_WORDS(hx,lx,x);
34*49393c00Smartynas GET_LDOUBLE_WORDS64(hy,ly,y);
35*49393c00Smartynas ix = hx&0x7fffffff; /* |x| */
36*49393c00Smartynas iy = hy&0x7fffffffffffffffLL; /* |y| */
37*49393c00Smartynas
38*49393c00Smartynas if(((ix>=0x7ff00000)&&((ix-0x7ff00000)|lx)!=0) || /* x is nan */
39*49393c00Smartynas ((iy>=0x7fff000000000000LL)&&((iy-0x7fff000000000000LL)|ly)!=0))
40*49393c00Smartynas /* y is nan */
41*49393c00Smartynas return x+y;
42*49393c00Smartynas if((long double) x==y) return y; /* x=y, return y */
43*49393c00Smartynas if((ix|lx)==0) { /* x == 0 */
44*49393c00Smartynas volatile double u;
45*49393c00Smartynas INSERT_WORDS(x,(u_int32_t)((hy>>32)&0x80000000),1);/* return +-minsub */
46*49393c00Smartynas u = x;
47*49393c00Smartynas u = u * u; /* raise underflow flag */
48*49393c00Smartynas return x;
49*49393c00Smartynas }
50*49393c00Smartynas if(hx>=0) { /* x > 0 */
51*49393c00Smartynas if (hy<0||(ix>>20)>(iy>>48)-0x3c00
52*49393c00Smartynas || ((ix>>20)==(iy>>48)-0x3c00
53*49393c00Smartynas && (((((int64_t)hx)<<28)|(lx>>4))>(hy&0x0000ffffffffffffLL)
54*49393c00Smartynas || (((((int64_t)hx)<<28)|(lx>>4))==(hy&0x0000ffffffffffffLL)
55*49393c00Smartynas && (lx&0xf)>(ly>>60))))) { /* x > y, x -= ulp */
56*49393c00Smartynas if(lx==0) hx -= 1;
57*49393c00Smartynas lx -= 1;
58*49393c00Smartynas } else { /* x < y, x += ulp */
59*49393c00Smartynas lx += 1;
60*49393c00Smartynas if(lx==0) hx += 1;
61*49393c00Smartynas }
62*49393c00Smartynas } else { /* x < 0 */
63*49393c00Smartynas if (hy>=0||(ix>>20)>(iy>>48)-0x3c00
64*49393c00Smartynas || ((ix>>20)==(iy>>48)-0x3c00
65*49393c00Smartynas && (((((int64_t)hx)<<28)|(lx>>4))>(hy&0x0000ffffffffffffLL)
66*49393c00Smartynas || (((((int64_t)hx)<<28)|(lx>>4))==(hy&0x0000ffffffffffffLL)
67*49393c00Smartynas && (lx&0xf)>(ly>>60))))) { /* x < y, x -= ulp */
68*49393c00Smartynas if(lx==0) hx -= 1;
69*49393c00Smartynas lx -= 1;
70*49393c00Smartynas } else { /* x > y, x += ulp */
71*49393c00Smartynas lx += 1;
72*49393c00Smartynas if(lx==0) hx += 1;
73*49393c00Smartynas }
74*49393c00Smartynas }
75*49393c00Smartynas hy = hx&0x7ff00000;
76*49393c00Smartynas if(hy>=0x7ff00000) {
77*49393c00Smartynas x = x+x; /* overflow */
78*49393c00Smartynas return x;
79*49393c00Smartynas }
80*49393c00Smartynas if(hy<0x00100000) {
81*49393c00Smartynas volatile double u = x*x; /* underflow */
82*49393c00Smartynas }
83*49393c00Smartynas INSERT_WORDS(x,hx,lx);
84*49393c00Smartynas return x;
85*49393c00Smartynas }
86