1 /* 2 * ==================================================== 3 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 4 * 5 * Developed at SunPro, a Sun Microsystems, Inc. business. 6 * Permission to use, copy, modify, and distribute this 7 * software is freely granted, provided that this notice 8 * is preserved. 9 * ==================================================== 10 */ 11 12 /* IEEE functions 13 * nextafter(x,y) 14 * return the next machine floating-point number of x in the 15 * direction toward y. 16 * Special cases: 17 */ 18 19 #include <float.h> 20 21 #include "fpmath.h" 22 #include "math.h" 23 #include "math_private.h" 24 25 #if LDBL_MAX_EXP != 0x4000 26 #error "Unsupported long double format" 27 #endif 28 29 long double 30 nextafterl(long double x, long double y) 31 { 32 volatile long double t; 33 union IEEEl2bits ux, uy; 34 35 if (isnan(x) || isnan(y)) 36 return x+y; /* x or y is nan */ 37 if(x==y) return y; /* x=y, return y */ 38 39 ux.e = x; 40 uy.e = y; 41 42 if(x==0.0) { 43 ux.bits.manh = 0; /* return +-minsubnormal */ 44 ux.bits.manl = 1; 45 ux.bits.sign = uy.bits.sign; 46 t = ux.e*ux.e; 47 if(t==ux.e) return t; else return ux.e; /* raise underflow flag */ 48 } 49 if(x>0.0 ^ x<y) { /* x -= ulp */ 50 if(ux.bits.manl==0) { 51 if ((ux.bits.manh&~LDBL_NBIT)==0) 52 ux.bits.exp -= 1; 53 ux.bits.manh = (ux.bits.manh - 1) | (ux.bits.manh & LDBL_NBIT); 54 } 55 ux.bits.manl -= 1; 56 } else { /* x += ulp */ 57 ux.bits.manl += 1; 58 if(ux.bits.manl==0) { 59 ux.bits.manh = (ux.bits.manh + 1) | (ux.bits.manh & LDBL_NBIT); 60 if ((ux.bits.manh&~LDBL_NBIT)==0) 61 ux.bits.exp += 1; 62 } 63 } 64 if(ux.bits.exp==0x7fff) return x+x; /* overflow */ 65 if(ux.bits.exp==0) { /* underflow */ 66 mask_nbit_l(ux); 67 t = ux.e * ux.e; 68 if(t!=ux.e) /* raise underflow flag */ 69 return ux.e; 70 } 71 return ux.e; 72 } 73 74 __strong_reference(nextafterl, nexttowardl); 75 76 #if LDBL_MANT_DIG == 113 77 __weak_reference(nextafterl, nextafterf128); 78 #endif 79