1/* 2 * inspired by glibc-2.0.6/sysdeps/libm-ieee754/s_nextafterf.c 3 * 4 * gcc -O2 -S -DOP=+ gives faddp %st(1),%st 5 * gcc -O2 -S -DOP=* gives fmulp %st(1),%st 6 * gcc -O2 -S -DOP=- gives fsubrp %st(1),%st 7 * gcc -O2 -S -DOP=/ gives fdivrp %st(1),%st 8 */ 9 10#ifndef OP 11#define OP * 12#endif 13 14typedef int int32_t __attribute__ ((__mode__ ( __SI__ ))) ; 15typedef unsigned int u_int32_t __attribute__ ((__mode__ ( __SI__ ))) ; 16 17typedef union 18{ 19 float value; 20 u_int32_t word; 21} ieee_float_shape_type; 22 23float __nextafterf(float x, float y) 24{ 25 int32_t hx,hy,ix,iy; 26 27 { 28 ieee_float_shape_type gf_u; 29 gf_u.value = x; 30 hx = gf_u.word; 31 } 32 { 33 ieee_float_shape_type gf_u; 34 gf_u.value = y; 35 hy = gf_u.word; 36 } 37 ix = hx&0x7fffffff; 38 iy = hy&0x7fffffff; 39 40 if ( ix > 0x7f800000 || iy > 0x7f800000 ) 41 return x+y; 42 if (x == y) return x; 43 if (ix == 0) 44 { 45 { 46 ieee_float_shape_type sf_u; 47 sf_u.word = (hy&0x80000000) | 1; 48 x = sf_u.value; 49 } 50 y = x*x; 51 if (y == x) return y; else return x; 52 } 53 if (hx >= 0) 54 { 55 if (hx > hy) 56 hx -= 1; 57 else 58 hx += 1; 59 } 60 else 61 { 62 if (hy >= 0 || hx > hy) 63 hx -= 1; 64 else 65 hx += 1; 66 } 67 hy = hx & 0x7f800000; 68 if (hy >= 0x7f800000) 69 return x+x; 70 if (hy < 0x00800000) 71 { 72 y = x OP x; 73 if (y != x) 74 { 75 ieee_float_shape_type sf_u; 76 sf_u.word = hx; 77 y = sf_u.value; 78 return y; 79 } 80 } 81 { 82 ieee_float_shape_type sf_u; 83 sf_u.word = hx; 84 x = sf_u.value; 85 } 86 return x; 87} 88 89 90