-/*
- * 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)