diff options
-rw-r--r-- | docs/probe_math_exception.c | 41 | ||||
-rw-r--r-- | libm/math_private.h | 55 | ||||
-rw-r--r-- | libm/s_nextafterf.c | 15 |
3 files changed, 100 insertions, 11 deletions
diff --git a/docs/probe_math_exception.c b/docs/probe_math_exception.c new file mode 100644 index 000000000..dbeccc5cc --- /dev/null +++ b/docs/probe_math_exception.c @@ -0,0 +1,41 @@ +/* Small test program for probing how various math functions + * with specific operands set floating point exceptions + */ + +#define _ISOC99_SOURCE 1 +#define _GNU_SOURCE 1 + +#include <math.h> +#include <fenv.h> +#include <stdio.h> + +int main(int argc, char **argv) +{ + float infF = HUGE_VALF * 2; + + 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 + +#define PREX(ex) do { if (fetestexcept(ex)) printf(#ex); } while(0) +#ifdef FE_INEXACT + PREX(FE_INEXACT); +#endif +#ifdef FE_DIVBYZERO + PREX(FE_DIVBYZERO); +#endif +#ifdef FE_UNDERFLOW + PREX(FE_UNDERFLOW); +#endif +#ifdef FE_OVERFLOW + PREX(FE_OVERFLOW); +#endif +#ifdef FE_INVALID + PREX(FE_INVALID); +#endif + if (fetestexcept(FE_ALL_EXCEPT)) + printf("\n"); + printf("done\n"); + return 0; +} diff --git a/libm/math_private.h b/libm/math_private.h index bdd0aba48..0601f2c3d 100644 --- a/libm/math_private.h +++ b/libm/math_private.h @@ -187,4 +187,59 @@ extern double __kernel_cos (double,double) attribute_hidden; 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_force_eval(x): force expression x to be evaluated and put 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). + */ +#if defined(__i386__) +#define math_opt_barrier(x) ({ \ + __typeof(x) __x = (x); \ + /* "t": load x into top-of-stack fpreg */ \ + __asm ("" : "=t" (__x) : "0" (__x)); \ + __x; \ +}) +#define math_force_eval(x) do { \ + if (sizeof(x) <= sizeof(double)) \ + /* "m": store x into a memory location */ \ + __asm __volatile ("" : : "m" (x)); \ + else /* long double */ \ + /* "f": load x into (any) fpreg */ \ + __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)) \ + /* "x": load into XMM SSE register */ \ + __asm ("" : "=x" (__x) : "0" (__x)); \ + else /* long double */ \ + /* "t": load x into top-of-stack fpreg */ \ + __asm ("" : "=t" (__x) : "0" (__x)); \ + __x; \ +}) +#define math_force_eval(x) do { \ + if (sizeof(x) <= sizeof(double)) \ + /* "x": load into XMM SSE register */ \ + __asm __volatile ("" : : "x" (x)); \ + else /* long double */ \ + /* "f": load x into (any) fpreg */ \ + __asm __volatile ("" : : "f" (x)); \ +} while (0) +#endif + +/* Default implementation forces 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 + + #endif /* _MATH_PRIVATE_H_ */ diff --git a/libm/s_nextafterf.c b/libm/s_nextafterf.c index 8dee00ff7..5fc44e31c 100644 --- a/libm/s_nextafterf.c +++ b/libm/s_nextafterf.c @@ -16,11 +16,6 @@ #include "math.h" #include "math_private.h" -#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 - float nextafterf(float x, float y) { int32_t hx, hy, ix, iy; @@ -38,12 +33,12 @@ float nextafterf(float x, float y) return y; if (ix == 0) { /* x == 0? */ - float u; +// glibc 2.4 does not seem to set underflow? +// float u; /* return +-minsubnormal */ SET_FLOAT_WORD(x, (hy & 0x80000000) | 1); - u = math_opt_barrier(x); - u = u * u; - math_force_eval(u); /* raise underflow flag */ +// u = x * x; /* raise underflow flag */ +// math_force_eval(u); return x; } @@ -63,8 +58,6 @@ float nextafterf(float x, float y) hy = hx & 0x7f800000; if (hy >= 0x7f800000) { x = x + x; /* overflow */ -//?? if (FLT_EVAL_METHOD != 0) -// asm ("" : "+m"(x)); return x; /* overflow */ } if (hy < 0x00800000) { |