diff options
| -rw-r--r-- | docs/probe_math_exception.c | 33 | ||||
| -rw-r--r-- | libc/sysdeps/linux/i386/bits/mathinline.h | 17 | ||||
| -rw-r--r-- | libm/math_private.h | 44 | ||||
| -rw-r--r-- | libm/s_finite.c | 9 | ||||
| -rw-r--r-- | libm/s_finitef.c | 9 |
5 files changed, 83 insertions, 29 deletions
diff --git a/docs/probe_math_exception.c b/docs/probe_math_exception.c index dbeccc5cc..dbc9020d4 100644 --- a/docs/probe_math_exception.c +++ b/docs/probe_math_exception.c @@ -5,20 +5,40 @@ #define _ISOC99_SOURCE 1 #define _GNU_SOURCE 1 +#include <stdint.h> #include <math.h> #include <fenv.h> #include <stdio.h> int main(int argc, char **argv) { - float infF = HUGE_VALF * 2; + float largest, small, t, inf_float; + + largest = small = 1; + while (1) { + t = largest + small; + /* optimizations may make plain "t == largest" unreliable */ + if (memcmp(&t, &largest, sizeof(float)) == 0) + break; + if (isfinite(t)) { + largest = t; + small *= 2; + continue; + } + small /= 2; + } + inf_float = largest + largest; + //printf("%.40g ", largest); + //printf("[%llx]\n", (long long) (*(uint32_t *)&largest)); feclearexcept(FE_ALL_EXCEPT); -// printf("%.40e\n", 1.0 / 0.0); // FE_DIVBYZERO -// printf("%.40e\n", nextafterf(HUGE_VALF, infF)); // no exceptions in glibc 2.4 + //t = 1.0 / 0.0; // simple test: FE_DIVBYZERO + //t = nextafterf(largest, 1); // glibc 2.8: no math exceptions raised + //t = nextafterf(largest, largest); // glibc 2.8: no math exceptions raised + //t = nextafterf(largest, inf_float); // glibc 2.8: FE_INEXACT FE_OVERFLOW -#define PREX(ex) do { if (fetestexcept(ex)) printf(#ex); } while(0) +#define PREX(ex) do { if (fetestexcept(ex)) printf(#ex " "); } while(0) #ifdef FE_INEXACT PREX(FE_INEXACT); #endif @@ -36,6 +56,9 @@ int main(int argc, char **argv) #endif if (fetestexcept(FE_ALL_EXCEPT)) printf("\n"); - printf("done\n"); + else + printf("no math exceptions raised\n"); + + printf("%.40g\n", t); return 0; } diff --git a/libc/sysdeps/linux/i386/bits/mathinline.h b/libc/sysdeps/linux/i386/bits/mathinline.h index ca7277d76..5caf73353 100644 --- a/libc/sysdeps/linux/i386/bits/mathinline.h +++ b/libc/sysdeps/linux/i386/bits/mathinline.h @@ -713,11 +713,22 @@ __inline_mathcodeNP2 (drem, __x, __y, \ __MATH_INLINE int __NTH (__finite (double __x)) { - return (__extension__ - (((((union { double __d; int __i[2]; }) {__d: __x}).__i[1] - | 0x800fffffu) + 1) >> 31)); + union { double __d; int __i[2]; } u; + u.__d = __x; + /* Finite numbers have at least one zero bit in exponent. */ + /* All other numbers will result in 0xffffffff after OR: */ + return (u.__i[1] | 0x800fffff) != 0xffffffff; } +__MATH_INLINE int +__NTH (__finitef (float __x)) +{ + union { float __d; int __i; } u; + u.__d = __x; + return (u.__i | 0x807fffff) != 0xffffffff; +} + + /* Miscellaneous functions */ # ifdef __FAST_MATH__ __inline_mathcode (__coshm1, __x, \ diff --git a/libm/math_private.h b/libm/math_private.h index 0601f2c3d..2c5a30ac2 100644 --- a/libm/math_private.h +++ b/libm/math_private.h @@ -188,14 +188,24 @@ extern double __kernel_tan (double,double,int) attribute_hidden; extern int __kernel_rem_pio2 (double*,double*,int,int,int,const int*) attribute_hidden; /* - * math_opt_barrier(x): force expression x to be evaluated and put into - * a floating point register or memory. This macro returns the value. + * math_opt_barrier(x): safely load x, even if it was manipulated + * by non-floationg point operations. This macro returns the value of x. + * This ensures compiler does not (ab)use its knowledge about x value + * and don't optimize future operations. Example: + * float x; + * SET_FLOAT_WORD(x, 0x80000001); // sets a bit pattern + * y = math_opt_barrier(x); // "compiler, do not cheat!" + * y = y * y; // compiler can't optimize, must use real multiply insn * - * math_force_eval(x): force expression x to be evaluated and put into + * math_force_eval(x): force expression x to be evaluated. + * Useful if otherwise compiler may eliminate the expression + * as unused. This macro returns no value. + * Example: "void fn(float f) { f = f * f; }" + * versus "void fn(float f) { f = f * f; math_force_eval(f); }" + * + * Currently, math_force_eval(x) stores x into * a floating point register or memory *of the appropriate size*. - * This forces floating point flags to be set correctly - * (for example, when float value is overflowing, but FPU registers - * are wide enough to "hide" this). + * There is no guarantee this will not change. */ #if defined(__i386__) #define math_opt_barrier(x) ({ \ @@ -205,19 +215,20 @@ extern int __kernel_rem_pio2 (double*,double*,int,int,int,const int*) attribu __x; \ }) #define math_force_eval(x) do { \ - if (sizeof(x) <= sizeof(double)) \ + __typeof(x) __x = (x); \ + if (sizeof(__x) <= sizeof(double)) \ /* "m": store x into a memory location */ \ - __asm __volatile ("" : : "m" (x)); \ + __asm __volatile ("" : : "m" (__x)); \ else /* long double */ \ /* "f": load x into (any) fpreg */ \ - __asm __volatile ("" : : "f" (x)); \ + __asm __volatile ("" : : "f" (__x)); \ } while (0) #endif #if defined(__x86_64__) #define math_opt_barrier(x) ({ \ __typeof(x) __x = (x); \ - if (sizeof(x) <= sizeof(double)) \ + if (sizeof(__x) <= sizeof(double)) \ /* "x": load into XMM SSE register */ \ __asm ("" : "=x" (__x) : "0" (__x)); \ else /* long double */ \ @@ -226,19 +237,22 @@ extern int __kernel_rem_pio2 (double*,double*,int,int,int,const int*) attribu __x; \ }) #define math_force_eval(x) do { \ - if (sizeof(x) <= sizeof(double)) \ + __typeof(x) __x = (x); \ + if (sizeof(__x) <= sizeof(double)) \ /* "x": load into XMM SSE register */ \ - __asm __volatile ("" : : "x" (x)); \ + __asm __volatile ("" : : "x" (__x)); \ else /* long double */ \ /* "f": load x into (any) fpreg */ \ - __asm __volatile ("" : : "f" (x)); \ + __asm __volatile ("" : : "f" (__x)); \ } while (0) #endif -/* Default implementation forces store to a memory location */ +/* Default implementations force store to a memory location */ #ifndef math_opt_barrier #define math_opt_barrier(x) ({ __typeof(x) __x = (x); __asm ("" : "+m" (__x)); __x; }) -#define math_force_eval(x) __asm __volatile ("" : : "m" (x)) +#endif +#ifndef math_force_eval +#define math_force_eval(x) do { __typeof(x) __x = (x); __asm __volatile ("" : : "m" (__x)); } while (0) #endif diff --git a/libm/s_finite.c b/libm/s_finite.c index 50fc11b50..9bbc00286 100644 --- a/libm/s_finite.c +++ b/libm/s_finite.c @@ -22,8 +22,11 @@ int __finite(double x) { - int32_t hx; - GET_HIGH_WORD(hx,x); - return (int)((u_int32_t)((hx&0x7fffffff)-0x7ff00000)>>31); + u_int32_t hx; + + GET_HIGH_WORD(hx, x); + /* Finite numbers have at least one zero bit in exponent. */ + /* All other numbers will result in 0xffffffff after OR: */ + return (hx | 0x800fffff) != 0xffffffff; } libm_hidden_def(__finite) diff --git a/libm/s_finitef.c b/libm/s_finitef.c index b7c853d4a..b427ea691 100644 --- a/libm/s_finitef.c +++ b/libm/s_finitef.c @@ -23,8 +23,11 @@ int __finitef(float x) { - int32_t ix; - GET_FLOAT_WORD(ix,x); - return (int)((u_int32_t)((ix&0x7fffffff)-0x7f800000)>>31); + u_int32_t ix; + + GET_FLOAT_WORD(ix, x); + /* Finite numbers have at least one zero bit in exponent. */ + /* All other numbers will result in 0xffffffff after OR: */ + return (ix | 0x807fffff) != 0xffffffff; } libm_hidden_def(__finitef) |