From: nsz Date: Wed, 7 Mar 2012 09:12:47 +0000 (+0100) Subject: enable complex in libm.h X-Git-Url: http://nsz.repo.hu/git/?p=libm;a=commitdiff_plain;h=5ed902d4e232ceac6db6c14f52c6c972ae7cf029 enable complex in libm.h --- diff --git a/src/internal/libm.h b/src/internal/libm.h index 13a4378..fecc50a 100644 --- a/src/internal/libm.h +++ b/src/internal/libm.h @@ -16,26 +16,11 @@ #include #include #include -#if 0 #include -#endif // FIXME #include "ldhack.h" -/* - * The original fdlibm code used statements like: - * n0 = ((*(int*)&one)>>29)^1; * index of high word * - * ix0 = *(n0+(int*)&x); * high word of x * - * ix1 = *((1-n0)+(int*)&x); * low word of x * - * to dig two 32 bit words out of the 64 bit IEEE floating point - * value. That is non-ANSI, and, moreover, the gcc instruction - * scheduler gets it wrong. We instead use the following macros. - * Unlike the original code, we determine the endianness at compile - * time, not at run time; I don't see much benefit to selecting - * endianness at run time. - */ - union fshape { float value; uint32_t bits; @@ -167,7 +152,7 @@ long double __p1evll(long double, long double *, int); */ void _scan_nan(uint32_t *__words, int __num_words, const char *__s); -// TODO: not needed when -fexcess-precision=standard is supported (>=gcc4.5) +// FIXME: not needed when -fexcess-precision=standard is supported (>=gcc4.5) /* * Attempt to get strict C99 semantics for assignment with non-C99 compilers. */ @@ -181,66 +166,32 @@ void _scan_nan(uint32_t *__words, int __num_words, const char *__s); #endif -#if 0 -/* - * C99 specifies that complex numbers have the same representation as - * an array of two elements, where the first element is the real part - * and the second element is the imaginary part. - */ -typedef union { - float complex f; - float a[2]; -} float_complex; -typedef union { - double complex f; - double a[2]; -} double_complex; -typedef union { - long double complex f; - long double a[2]; -} long_double_complex; -#define REALPART(z) ((z).a[0]) -#define IMAGPART(z) ((z).a[1]) +/* complex */ -/* - * Inline functions that can be used to construct complex values. - * - * The C99 standard intends x+I*y to be used for this, but x+I*y is - * currently unusable in general since gcc introduces many overflow, - * underflow, sign and efficiency bugs by rewriting I*y as - * (0.0+I)*(y+0.0*I) and laboriously computing the full complex product. - * In particular, I*Inf is corrupted to NaN+I*Inf, and I*-0 is corrupted - * to -0.0+I*0.0. - */ -static inline float complex -cpackf(float x, float y) -{ - float_complex z; - - REALPART(z) = x; - IMAGPART(z) = y; - return (z.f); -} - -static inline double complex -cpack(double x, double y) -{ - double_complex z; - - REALPART(z) = x; - IMAGPART(z) = y; - return (z.f); -} - -static inline long double complex -cpackl(long double x, long double y) -{ - long_double_complex z; - - REALPART(z) = x; - IMAGPART(z) = y; - return (z.f); -} -#endif +union dcomplex { + double complex z; + double a[2]; +}; +union fcomplex { + float complex z; + float a[2]; +}; +union lcomplex { + long double complex z; + long double a[2]; +}; + +// FIXME: move to complex.h ? +#define creal(z) ((double)(z)) +#define crealf(z) ((float)(z)) +#define creall(z) ((long double)(z)) +#define cimag(z) ((union dcomplex){(z)}.a[1]) +#define cimagf(z) ((union fcomplex){(z)}.a[1]) +#define cimagl(z) ((union lcomplex){(z)}.a[1]) + +/* x + y*I is not supported properly by gcc */ +#define cpack(x,y) ((union dcomplex){.a={(x),(y)}}.z) +#define cpackf(x,y) ((union fcomplex){.a={(x),(y)}}.z) +#define cpackl(x,y) ((union lcomplex){.a={(x),(y)}}.z) #endif