1
2 /* @(#)s_nextafter.c 5.1 93/09/24 */
3 /*
4 * ====================================================
5 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
6 *
7 * Developed at SunPro, a Sun Microsystems, Inc. business.
8 * Permission to use, copy, modify, and distribute this
9 * software is freely granted, provided that this notice
10 * is preserved.
11 * ====================================================
12 */
13
14 /*
15 FUNCTION
16 <<nextafter>>, <<nextafterf>>---get next number
17
18 INDEX
19 nextafter
20 INDEX
21 nextafterf
22
23 SYNOPSIS
24 #include <math.h>
25 double nextafter(double <[val]>, double <[dir]>);
26 float nextafterf(float <[val]>, float <[dir]>);
27
28 DESCRIPTION
29 <<nextafter>> returns the double-precision floating-point number
30 closest to <[val]> in the direction toward <[dir]>. <<nextafterf>>
31 performs the same operation in single precision. For example,
32 <<nextafter(0.0,1.0)>> returns the smallest positive number which is
33 representable in double precision.
34
35 RETURNS
36 Returns the next closest number to <[val]> in the direction toward
37 <[dir]>.
38
39 PORTABILITY
40 Neither <<nextafter>> nor <<nextafterf>> is required by ANSI C
41 or by the System V Interface Definition (Issue 2).
42 */
43
44 /* IEEE functions
45 * nextafter(x,y)
46 * return the next machine floating-point number of x in the
47 * direction toward y.
48 * Special cases:
49 */
50
51 #include "fdlibm.h"
52
53 #ifdef _NEED_FLOAT64
54
55 __float64
nextafter64(__float64 x,__float64 y)56 nextafter64(__float64 x, __float64 y)
57 {
58 __int32_t hx,hy,ix,iy;
59 __uint32_t lx,ly;
60
61 EXTRACT_WORDS(hx,lx,x);
62 EXTRACT_WORDS(hy,ly,y);
63 ix = hx&0x7fffffff; /* |x| */
64 iy = hy&0x7fffffff; /* |y| */
65
66 if(((ix>=0x7ff00000)&&((ix-0x7ff00000)|lx)!=0) || /* x is nan */
67 ((iy>=0x7ff00000)&&((iy-0x7ff00000)|ly)!=0)) /* y is nan */
68 return x+y;
69 if(x==y) return y; /* x=y, return y (follow y sign for 0) */
70 if((ix|lx)==0) { /* x == 0 */
71 INSERT_WORDS(x,hy&0x80000000,1); /* return +-minsubnormal */
72 force_eval_float64(x*x);
73 return x;
74 }
75 if(hx>=0) { /* x > 0 */
76 if(hx>hy||((hx==hy)&&(lx>ly))) { /* x > y, x -= ulp */
77 if(lx==0) hx -= 1;
78 lx -= 1;
79 } else { /* x < y, x += ulp */
80 lx += 1;
81 if(lx==0) hx += 1;
82 }
83 } else { /* x < 0 */
84 if(hy>=0||hx>hy||((hx==hy)&&(lx>ly))){/* x < y, x -= ulp */
85 if(lx==0) hx -= 1;
86 lx -= 1;
87 } else { /* x > y, x += ulp */
88 lx += 1;
89 if(lx==0) hx += 1;
90 }
91 }
92 hy = hx&0x7ff00000;
93 if(hy>=0x7ff00000)
94 return __math_oflow(hx<0); /* overflow */
95 INSERT_WORDS(x,hx,lx);
96 if(hy<0x00100000) /* underflow */
97 return __math_denorm(x);
98 return (x);
99 }
100
101 _MATH_ALIAS_d_dd(nextafter)
102
103 #endif /* _NEED_FLOAT64 */
104
