Replace Cactus complex number type with C/C++ complex numbers

Map CCTK_COMPLEX to "double complex" in C, and "complex<double>" in
C++. (It is already mapped to "double complex" in Fortran.)

Update type definitions.

Re-implement Cactus complex number math functions by calling the
respective C functions.

Update thorn that access real and imaginary parts of complex numbers
to use standard-conforming methods instead.


git-svn-id: http://svn.cactuscode.org/flesh/trunk@4979 17b73243-c579-4c4c-a9d2-2d5706c11dac

4 files changed, 185 insertions, 438 deletions

diff --git a/src/include/cctk_Complex.h b/src/include/cctk_Complex.h
index 977896d3..96eaf17f 100644
--- a/src/include/cctk_Complex.h
+++ b/src/include/cctk_Complex.h
@@ -12,92 +12,32 @@
 #define _CCTK_COMPLEX_H_
 
 #ifdef __cplusplus
-#  define EXTERN_C_BEGIN extern "C" {
-#  define EXTERN_C_END   }
-#else
-#  define EXTERN_C_BEGIN
-#  define EXTERN_C_END
-#endif
-
-
-#ifndef DEFINE_CCTK_COMPLEX_EXTERN_FUNCTIONS
-#  define PREFIX static inline
-#else
-#  define PREFIX
-#endif
-
-
-#ifdef __cplusplus
-  /* declare C++ overloaded operators */
-#  include <complex>
-#  include <ostream>
-#  define DECLARE_CMPLX_CXX_OPERATORS(CCTK_Cmplx, cctk_real, cctk_complex) \
-PREFIX cctk_complex operator+  (cctk_complex const & a);                \
-PREFIX cctk_complex operator-  (cctk_complex const & a);                \
-PREFIX cctk_complex conj       (cctk_complex const & a);                \
-PREFIX cctk_real    abs        (cctk_complex const & a);                \
-PREFIX cctk_real    arg        (cctk_complex const & a);                \
-PREFIX cctk_real    norm       (cctk_complex const & a);                \
-PREFIX cctk_complex operator+  (cctk_complex const & a, cctk_complex const & b); \
-PREFIX cctk_complex operator+  (cctk_real    const & a, cctk_complex const & b); \
-PREFIX cctk_complex operator+  (cctk_complex const & a, cctk_real    const & b); \
-PREFIX cctk_complex operator-  (cctk_complex const & a, cctk_complex const & b); \
-PREFIX cctk_complex operator-  (cctk_real    const & a, cctk_complex const & b); \
-PREFIX cctk_complex operator-  (cctk_complex const & a, cctk_real    const & b); \
-PREFIX cctk_complex operator*  (cctk_complex const & a, cctk_complex const & b); \
-PREFIX cctk_complex operator*  (cctk_real    const & a, cctk_complex const & b); \
-PREFIX cctk_complex operator*  (cctk_complex const & a, cctk_real    const & b); \
-PREFIX cctk_complex operator/  (cctk_complex const & a, cctk_complex const & b); \
-PREFIX cctk_complex operator/  (cctk_real    const & a, cctk_complex const & b); \
-PREFIX cctk_complex operator/  (cctk_complex const & a, cctk_real    const & b); \
-PREFIX cctk_complex operator+= (cctk_complex       & a, cctk_complex const & b); \
-PREFIX cctk_complex operator+= (cctk_complex       & a, cctk_real    const & b); \
-PREFIX cctk_complex operator-= (cctk_complex       & a, cctk_complex const & b); \
-PREFIX cctk_complex operator-= (cctk_complex       & a, cctk_real    const & b); \
-PREFIX cctk_complex operator*= (cctk_complex       & a, cctk_complex const & b); \
-PREFIX cctk_complex operator*= (cctk_complex       & a, cctk_real    const & b); \
-PREFIX cctk_complex operator/= (cctk_complex       & a, cctk_complex const & b); \
-PREFIX cctk_complex operator/= (cctk_complex       & a, cctk_real    const & b); \
-PREFIX cctk_complex pow        (cctk_complex const & a, int                  i); \
-PREFIX cctk_complex pow        (cctk_complex const & a, cctk_real    const & r); \
-PREFIX cctk_complex pow        (cctk_complex const & a, cctk_complex const & b); \
-PREFIX cctk_complex sin        (cctk_complex const & a);                \
-PREFIX cctk_complex cos        (cctk_complex const & a);                \
-PREFIX cctk_complex exp        (cctk_complex const & a);                \
-PREFIX cctk_complex log        (cctk_complex const & a);                \
-PREFIX cctk_complex sqrt       (cctk_complex const & a);                \
-PREFIX std::ostream & operator << (std::ostream & os, cctk_complex const & a);
-#else
-  /* declare no C++ overloaded operators */
-#  define DECLARE_CMPLX_CXX_OPERATORS(CCTK_Cmplx, cctk_real, cctk_complex)
+extern "C" {
 #endif
 
 
 /* macro to declare a set of complex functions for a given precision */
 #define DECLARE_CMPLX_FUNCTIONS(CCTK_Cmplx, cctk_real, cctk_complex)    \
-EXTERN_C_BEGIN                                                          \
-PREFIX cctk_complex CCTK_Cmplx        (cctk_real Re, cctk_real Im);     \
-PREFIX cctk_real    CCTK_Cmplx##Real  (cctk_complex a);                 \
-PREFIX cctk_real    CCTK_Cmplx##Imag  (cctk_complex a);                 \
-PREFIX cctk_complex CCTK_Cmplx##Neg   (cctk_complex a);                 \
-PREFIX cctk_complex CCTK_Cmplx##Conjg (cctk_complex a);                 \
-PREFIX cctk_real    CCTK_Cmplx##Abs   (cctk_complex a);                 \
-PREFIX cctk_real    CCTK_Cmplx##Arg   (cctk_complex a);                 \
-PREFIX cctk_real    CCTK_Cmplx##Norm  (cctk_complex a);                 \
-PREFIX cctk_complex CCTK_Cmplx##Add   (cctk_complex a, cctk_complex b); \
-PREFIX cctk_complex CCTK_Cmplx##Sub   (cctk_complex a, cctk_complex b); \
-PREFIX cctk_complex CCTK_Cmplx##Mul   (cctk_complex a, cctk_complex b); \
-PREFIX cctk_complex CCTK_Cmplx##Div   (cctk_complex a, cctk_complex b); \
-PREFIX cctk_complex CCTK_Cmplx##CPow  (cctk_complex a, cctk_complex b); \
-PREFIX cctk_complex CCTK_Cmplx##Sin   (cctk_complex a);                 \
-PREFIX cctk_complex CCTK_Cmplx##Cos   (cctk_complex a);                 \
-PREFIX cctk_complex CCTK_Cmplx##Exp   (cctk_complex a);                 \
-PREFIX cctk_complex CCTK_Cmplx##Log   (cctk_complex a);                 \
-PREFIX cctk_complex CCTK_Cmplx##Sqrt  (cctk_complex a);                 \
-PREFIX cctk_complex CCTK_Cmplx##Pow   (cctk_complex a, cctk_real b);    \
-PREFIX cctk_complex CCTK_Cmplx##IPow  (cctk_complex a, int b);          \
-EXTERN_C_END                                                            \
-DECLARE_CMPLX_CXX_OPERATORS(CCTK_Cmplx, cctk_real, cctk_complex)
+cctk_complex CCTK_Cmplx        (cctk_real Re, cctk_real Im);            \
+cctk_real    CCTK_Cmplx##Real  (cctk_complex a);                        \
+cctk_real    CCTK_Cmplx##Imag  (cctk_complex a);                        \
+cctk_complex CCTK_Cmplx##Neg   (cctk_complex a);                        \
+cctk_complex CCTK_Cmplx##Conjg (cctk_complex a);                        \
+cctk_real    CCTK_Cmplx##Abs   (cctk_complex a);                        \
+cctk_real    CCTK_Cmplx##Arg   (cctk_complex a);                        \
+cctk_real    CCTK_Cmplx##Norm  (cctk_complex a);                        \
+cctk_complex CCTK_Cmplx##Add   (cctk_complex a, cctk_complex b);        \
+cctk_complex CCTK_Cmplx##Sub   (cctk_complex a, cctk_complex b);        \
+cctk_complex CCTK_Cmplx##Mul   (cctk_complex a, cctk_complex b);        \
+cctk_complex CCTK_Cmplx##Div   (cctk_complex a, cctk_complex b);        \
+cctk_complex CCTK_Cmplx##CPow  (cctk_complex a, cctk_complex b);        \
+cctk_complex CCTK_Cmplx##Sin   (cctk_complex a);                        \
+cctk_complex CCTK_Cmplx##Cos   (cctk_complex a);                        \
+cctk_complex CCTK_Cmplx##Exp   (cctk_complex a);                        \
+cctk_complex CCTK_Cmplx##Log   (cctk_complex a);                        \
+cctk_complex CCTK_Cmplx##Sqrt  (cctk_complex a);                        \
+cctk_complex CCTK_Cmplx##Pow   (cctk_complex a, cctk_real b);           \
+cctk_complex CCTK_Cmplx##IPow  (cctk_complex a, int b);
 
 
 /* declare complex functions for all available precisions */
@@ -113,8 +53,6 @@ DECLARE_CMPLX_FUNCTIONS (CCTK_Cmplx16, CCTK_REAL8, CCTK_COMPLEX16)
 DECLARE_CMPLX_FUNCTIONS (CCTK_Cmplx32, CCTK_REAL16, CCTK_COMPLEX32)
 #endif
 
-#undef PREFIX
-
 
 /* declare the default precision complex functions as #define'd macros */
 #ifdef CCTK_REAL_PRECISION_4
@@ -182,10 +120,8 @@ DECLARE_CMPLX_FUNCTIONS (CCTK_Cmplx32, CCTK_REAL16, CCTK_COMPLEX32)
 #define CCTK_CmplxIPow  CCTK_Cmplx32IPow
 #endif
 
-#ifndef DEFINE_CCTK_COMPLEX_EXTERN_FUNCTIONS
-#  define DEFINE_CCTK_COMPLEX_INLINE_FUNCTIONS
-#  include "../main/Complex.c"
-#  undef DEFINE_CCTK_COMPLEX_INLINE_FUNCTIONS
+#ifdef __cplusplus
+}
 #endif
 
 #endif /* _CCTK_COMPLEX_H_ */
diff --git a/src/include/cctk_Types.h b/src/include/cctk_Types.h
index e3fb55a6..eeda7d77 100644
--- a/src/include/cctk_Types.h
+++ b/src/include/cctk_Types.h
@@ -34,45 +34,37 @@ typedef const char * CCTK_STRING;
 #define HAVE_CCTK_CHAR 1
 #define HAVE_CCTK_STRING 1
 
-/* Structures for complex types */
+/* Declarations for complex types */
 
-#ifdef HAVE_CCTK_REAL16
-#define HAVE_CCTK_COMPLEX32 1
-typedef struct CCTK_COMPLEX32
-{
-  CCTK_REAL16 Re;
-  CCTK_REAL16 Im;
 #ifdef __cplusplus
-  CCTK_REAL16 real() const { return Re; }
-  CCTK_REAL16 imag() const { return Im; }
+#  include <complex>
 #endif
-} CCTK_COMPLEX32;
+
+#ifdef HAVE_CCTK_REAL16
+#  define HAVE_CCTK_COMPLEX32 1
+#  ifdef __cplusplus
+typedef std::complex<CCTK_REAL16> CCTK_COMPLEX32;
+#  else
+typedef long double _Complex CCTK_COMPLEX32;
+#  endif
 #endif
 
 #ifdef HAVE_CCTK_REAL8
 #define HAVE_CCTK_COMPLEX16 1
-typedef struct CCTK_COMPLEX16
-{
-  CCTK_REAL8 Re;
-  CCTK_REAL8 Im;
-#ifdef __cplusplus
-  CCTK_REAL8 real() const { return Re; }
-  CCTK_REAL8 imag() const { return Im; }
-#endif
-} CCTK_COMPLEX16;
+#  ifdef __cplusplus
+typedef std::complex<CCTK_REAL8> CCTK_COMPLEX16;
+#  else
+typedef double _Complex CCTK_COMPLEX16;
+#  endif
 #endif
 
 #ifdef HAVE_CCTK_REAL4
-#define HAVE_CCTK_COMPLEX8 1
-typedef struct CCTK_COMPLEX8
-{
-  CCTK_REAL4 Re;
-  CCTK_REAL4 Im;
-#ifdef __cplusplus
-  CCTK_REAL4 real() const { return Re; }
-  CCTK_REAL4 imag() const { return Im; }
-#endif
-} CCTK_COMPLEX8;
+#  define HAVE_CCTK_COMPLEX8 1
+#  ifdef __cplusplus
+typedef std::complex<CCTK_REAL4> CCTK_COMPLEX8;
+#  else
+typedef float _Complex CCTK_COMPLEX8;
+#  endif
 #endif
 
 /* Small positive integer type */
@@ -262,7 +254,7 @@ typedef unsigned char CCTK_BYTE;
 
 /* We have __attribute__((unused)), so use it */
 #define CCTK_DECLARE_INIT(typ,nam,val)          \
-  typ nam __attribute__((unused)) = (val);
+  typ nam CCTK_ATTRIBUTE_UNUSED = (val);
 
 #else
 
diff --git a/src/main/Complex.c b/src/main/Complex.c
index 47e3a888..84bcf7c6 100644
--- a/src/main/Complex.c
+++ b/src/main/Complex.c
@@ -8,6 +8,7 @@
    @version   $Id$
  @@*/
 
+#include <complex.h>
 #include <math.h>
 
 #include "cctk_Flesh.h"
@@ -72,16 +73,10 @@ CCTK_FILEVERSION(main_Complex_c);
                The complex number
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX(CCTK_Cmplx, cctk_real, cctk_complex)                \
-cctk_complex CCTK_Cmplx (cctk_real Re, cctk_real Im)                          \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result.Re = Re;                                                             \
-  result.Im = Im;                                                             \
-                                                                              \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX(CCTK_Cmplx, cctk_real, cctk_complex, type)    \
+cctk_complex CCTK_Cmplx (cctk_real Re, cctk_real Im)                    \
+{                                                                       \
+  return Re + _Complex_I * Im;                                          \
 }
 
 
@@ -93,7 +88,7 @@ cctk_complex CCTK_Cmplx (cctk_real Re, cctk_real Im)                          \
                Returns the real part of a complex number.
    @enddesc
 
-   @var        complex_number
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -104,10 +99,10 @@ cctk_complex CCTK_Cmplx (cctk_real Re, cctk_real Im)                          \
                The real part
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_REAL(CCTK_Cmplx, cctk_real, cctk_complex)           \
-cctk_real CCTK_Cmplx##Real (cctk_complex complex_number)                      \
-{                                                                             \
-  return (complex_number.Re);                                                 \
+#define DEFINE_CCTK_CMPLX_REAL(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_real CCTK_Cmplx##Real (cctk_complex x)                             \
+{                                                                       \
+  return creal##type(x);                                                \
 }
 
 
@@ -119,7 +114,7 @@ cctk_real CCTK_Cmplx##Real (cctk_complex complex_number)                      \
                Returns the imaginary part of a complex number.
    @enddesc
 
-   @var        complex_number
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -130,10 +125,10 @@ cctk_real CCTK_Cmplx##Real (cctk_complex complex_number)                      \
                The imaginary part
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_IMAG(CCTK_Cmplx, cctk_real, cctk_complex)           \
-cctk_real CCTK_Cmplx##Imag (cctk_complex complex_number)                      \
-{                                                                             \
-  return (complex_number.Im);                                                 \
+#define DEFINE_CCTK_CMPLX_IMAG(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_real CCTK_Cmplx##Imag (cctk_complex x)                             \
+{                                                                       \
+  return cimag##type(x);                                                \
 }
 
 
@@ -145,7 +140,7 @@ cctk_real CCTK_Cmplx##Imag (cctk_complex complex_number)                      \
                Returns the negative of a complex number.
    @enddesc
  
-   @var       in
+   @var       x
    @vdesc     The complex number
    @vtype     CCTK_COMPLEX
    @vio       in
@@ -156,15 +151,10 @@ cctk_real CCTK_Cmplx##Imag (cctk_complex complex_number)                      \
                The negative
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_NEG(CCTK_Cmplx, cctk_real, cctk_complex)            \
-cctk_complex CCTK_Cmplx##Neg (cctk_complex complex_number)                    \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result.Re = -complex_number.Re;                                             \
-  result.Im = -complex_number.Im;                                             \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_NEG(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Neg (cctk_complex x)                           \
+{                                                                       \
+  return -x;                                                            \
 }
 
 
@@ -176,7 +166,7 @@ cctk_complex CCTK_Cmplx##Neg (cctk_complex complex_number)                    \
                Returns the complex conjugate of a complex number.
    @enddesc
  
-   @var       in
+   @var       x
    @vdesc     The complex number
    @vtype     CCTK_COMPLEX
    @vio       in
@@ -187,15 +177,10 @@ cctk_complex CCTK_Cmplx##Neg (cctk_complex complex_number)                    \
                The complex conjugate
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_CONJG(CCTK_Cmplx, cctk_real, cctk_complex)          \
-cctk_complex CCTK_Cmplx##Conjg (cctk_complex complex_number)                  \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result.Re =  complex_number.Re;                                             \
-  result.Im = -complex_number.Im;                                             \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_CONJG(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Conjg (cctk_complex x)                         \
+{                                                                       \
+  return conj##type(x);                                                 \
 }
 
 
@@ -207,7 +192,7 @@ cctk_complex CCTK_Cmplx##Conjg (cctk_complex complex_number)                  \
                Return the absolute value of a complex number.
    @enddesc
  
-   @var        in
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -218,10 +203,10 @@ cctk_complex CCTK_Cmplx##Conjg (cctk_complex complex_number)                  \
                The absolute value of the complex number
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_ABS(CCTK_Cmplx, cctk_real, cctk_complex, type)      \
-cctk_real CCTK_Cmplx##Abs (cctk_complex complex_number)                       \
-{                                                                             \
-  return ((cctk_real)hypot##type (complex_number.Re, complex_number.Im));     \
+#define DEFINE_CCTK_CMPLX_ABS(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_real CCTK_Cmplx##Abs (cctk_complex x)                              \
+{                                                                       \
+  return cabs##type(x);                                                 \
 }
 
 
@@ -233,7 +218,7 @@ cctk_real CCTK_Cmplx##Abs (cctk_complex complex_number)                       \
                Return the absolute value squared of a complex number.
    @enddesc
  
-   @var        in
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -244,10 +229,10 @@ cctk_real CCTK_Cmplx##Abs (cctk_complex complex_number)                       \
                The absolute value squared of the complex number
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_NORM(CCTK_Cmplx, cctk_real, cctk_complex, type)     \
-cctk_real CCTK_Cmplx##Norm (cctk_complex complex_number)                      \
-{                                                                             \
-  return (pow##type (complex_number.Re, 2) + pow##type (complex_number.Im, 2));\
+#define DEFINE_CCTK_CMPLX_NORM(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_real CCTK_Cmplx##Norm (cctk_complex x)                             \
+{                                                                       \
+  return creal##type(x)*creal##type(x) + cimag##type(x)*cimag##type(x); \
 }
 
 
@@ -259,7 +244,7 @@ cctk_real CCTK_Cmplx##Norm (cctk_complex complex_number)                      \
                Return the argument of a complex number.
    @enddesc
  
-   @var        in
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -270,10 +255,10 @@ cctk_real CCTK_Cmplx##Norm (cctk_complex complex_number)                      \
                The argument of the complex number
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_ARG(CCTK_Cmplx, cctk_real, cctk_complex, type)      \
-cctk_real CCTK_Cmplx##Arg (cctk_complex complex_number)                       \
-{                                                                             \
-  return (atan2##type (complex_number.Im, complex_number.Re));                \
+#define DEFINE_CCTK_CMPLX_ARG(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_real CCTK_Cmplx##Arg (cctk_complex x)                              \
+{                                                                       \
+  return carg##type(x);                                                 \
 }
 
 
@@ -301,15 +286,10 @@ cctk_real CCTK_Cmplx##Arg (cctk_complex complex_number)                       \
                The sum of a and b
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_ADD(CCTK_Cmplx, cctk_real, cctk_complex)            \
-cctk_complex CCTK_Cmplx##Add (cctk_complex a, cctk_complex b)                 \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result.Re = a.Re + b.Re;                                                    \
-  result.Im = a.Im + b.Im;                                                    \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_ADD(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Add (cctk_complex a, cctk_complex b)           \
+{                                                                       \
+  return a+b;                                                           \
 }
 
 
@@ -337,15 +317,10 @@ cctk_complex CCTK_Cmplx##Add (cctk_complex a, cctk_complex b)                 \
                The difference
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_SUB(CCTK_Cmplx, cctk_real, cctk_complex)            \
-cctk_complex CCTK_Cmplx##Sub (cctk_complex a, cctk_complex b)                 \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result.Re = a.Re - b.Re;                                                    \
-  result.Im = a.Im - b.Im;                                                    \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_SUB(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Sub (cctk_complex a, cctk_complex b)           \
+{                                                                       \
+  return a-b;                                                           \
 }
 
 
@@ -373,15 +348,10 @@ cctk_complex CCTK_Cmplx##Sub (cctk_complex a, cctk_complex b)                 \
                The product
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_MUL(CCTK_Cmplx, cctk_real, cctk_complex)            \
-cctk_complex CCTK_Cmplx##Mul (cctk_complex a, cctk_complex b)                 \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result.Re = a.Re*b.Re - a.Im*b.Im;                                          \
-  result.Im = a.Im*b.Re + a.Re*b.Im;                                          \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_MUL(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Mul (cctk_complex a, cctk_complex b)           \
+{                                                                       \
+  return a*b;                                                           \
 }
 
 
@@ -391,10 +361,6 @@ cctk_complex CCTK_Cmplx##Mul (cctk_complex a, cctk_complex b)                 \
    @author     Gabrielle Allen
    @desc
                Divide two complex numbers.
-               Make sure that the intermediate values do not overflow.
-               See e.g. the C99 compliant implementation in libgcc
-               (which ships with gcc).
-               This implementation here does not handle nan and inf.
    @enddesc
  
    @var        a
@@ -413,24 +379,10 @@ cctk_complex CCTK_Cmplx##Mul (cctk_complex a, cctk_complex b)                 \
                The quotient
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_DIV(CCTK_Cmplx, cctk_real, cctk_complex, type)      \
-cctk_complex CCTK_Cmplx##Div (cctk_complex a, cctk_complex b)                 \
-{                                                                             \
-  cctk_real afact, bfact, factor;                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  afact = (cctk_real)(fabs##type(a.Re) + fabs##type(a.Im));                   \
-  a.Re /= afact;                                                              \
-  a.Im /= afact;                                                              \
-  bfact = (cctk_real)(fabs##type(b.Re) + fabs##type(b.Im));                   \
-  b.Re /= bfact;                                                              \
-  b.Im /= bfact;                                                              \
-  factor = b.Re*b.Re + b.Im*b.Im;                                             \
-  result.Re = afact / bfact * (a.Re*b.Re + a.Im*b.Im) / factor;               \
-  result.Im = afact / bfact * (a.Im*b.Re - a.Re*b.Im) / factor;               \
-                                                                              \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_DIV(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Div (cctk_complex a, cctk_complex b)           \
+{                                                                       \
+  return a/b;                                                           \
 }
 
 
@@ -440,8 +392,6 @@ cctk_complex CCTK_Cmplx##Div (cctk_complex a, cctk_complex b)                 \
    @author     Erik Schnetter
    @desc
                Raises a complex number to a given power
-               This algorithm was taken from glibc 2.3.5,
-               file sysdeps/generic/s_cpow.c.
    @enddesc
  
    @var        a
@@ -460,10 +410,10 @@ cctk_complex CCTK_Cmplx##Div (cctk_complex a, cctk_complex b)                 \
                a to the power of b
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_CPOW(CCTK_Cmplx, cctk_real, cctk_complex)           \
-cctk_complex CCTK_Cmplx##CPow (cctk_complex a, cctk_complex b)                \
-{                                                                             \
-  return CCTK_Cmplx##Exp (CCTK_Cmplx##Mul (b, CCTK_Cmplx##Log (a)));          \
+#define DEFINE_CCTK_CMPLX_CPOW(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##CPow (cctk_complex a, cctk_complex b)          \
+{                                                                       \
+  return cpow##type(a, b);                                              \
 }
 
 
@@ -475,7 +425,7 @@ cctk_complex CCTK_Cmplx##CPow (cctk_complex a, cctk_complex b)                \
                Returns the sine of a complex number.
    @enddesc
  
-   @var        complex_number
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -486,16 +436,10 @@ cctk_complex CCTK_Cmplx##CPow (cctk_complex a, cctk_complex b)                \
                The sine
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_SIN(CCTK_Cmplx, cctk_real, cctk_complex, type)      \
-cctk_complex CCTK_Cmplx##Sin (cctk_complex complex_number)                    \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result.Re = sin##type (complex_number.Re) * cosh##type (complex_number.Im); \
-  result.Im = cos##type (complex_number.Re) * sinh##type (complex_number.Im); \
-                                                                              \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_SIN(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Sin (cctk_complex x)                           \
+{                                                                       \
+  return csin##type(x);                                                 \
 }
 
 
@@ -507,7 +451,7 @@ cctk_complex CCTK_Cmplx##Sin (cctk_complex complex_number)                    \
                Returns the cosine of a complex number.
    @enddesc
  
-   @var        complex_number
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -518,16 +462,10 @@ cctk_complex CCTK_Cmplx##Sin (cctk_complex complex_number)                    \
                The cosine
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_COS(CCTK_Cmplx, cctk_real, cctk_complex, type)      \
-cctk_complex CCTK_Cmplx##Cos (cctk_complex complex_number)                    \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result.Re = cos##type (complex_number.Re) * cosh##type (complex_number.Im); \
-  result.Im = sin##type (complex_number.Re) * sinh##type (complex_number.Im); \
-                                                                              \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_COS(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Cos (cctk_complex x)                           \
+{                                                                       \
+  return ccos##type(x);                                                 \
 }
 
 
@@ -539,7 +477,7 @@ cctk_complex CCTK_Cmplx##Cos (cctk_complex complex_number)                    \
                Returns the exponential of a complex number.
    @enddesc
  
-   @var        complex_number
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -550,19 +488,10 @@ cctk_complex CCTK_Cmplx##Cos (cctk_complex complex_number)                    \
                The exponential
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_EXP(CCTK_Cmplx, cctk_real, cctk_complex, type)      \
-cctk_complex CCTK_Cmplx##Exp (cctk_complex complex_number)                    \
-{                                                                             \
-  cctk_real rho, theta;                                                       \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  rho = (cctk_real)exp##type (complex_number.Re);                             \
-  theta = complex_number.Im;                                                  \
-  result.Re = rho * (cctk_real)cos##type (theta);                             \
-  result.Im = rho * (cctk_real)sin##type (theta);                             \
-                                                                              \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_EXP(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Exp (cctk_complex x)                           \
+{                                                                       \
+  return cexp##type(x);                                                 \
 }
 
 
@@ -572,11 +501,9 @@ cctk_complex CCTK_Cmplx##Exp (cctk_complex complex_number)                    \
    @author     Erik Schnetter
    @desc
                Returns the logarithm of a complex number.
-               This algorithm was taken from glibc 2.3.5,
-               file sysdeps/generic/s_clog.c.
    @enddesc
  
-   @var        complex_number
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -587,16 +514,10 @@ cctk_complex CCTK_Cmplx##Exp (cctk_complex complex_number)                    \
                The logarithm
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_LOG(CCTK_Cmplx, cctk_real, cctk_complex, type)      \
-cctk_complex CCTK_Cmplx##Log (cctk_complex complex_number)                    \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result.Re = log##type (hypot##type (complex_number.Re, complex_number.Im)); \
-  result.Im = atan2##type (complex_number.Im, complex_number.Re);             \
-                                                                              \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_LOG(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Log (cctk_complex x)                           \
+{                                                                       \
+  return clog##type(x);                                                 \
 }
 
 
@@ -606,15 +527,9 @@ cctk_complex CCTK_Cmplx##Log (cctk_complex complex_number)                    \
    @author     Thomas Radke
    @desc
                Returns the square root of a complex number.
-               This algorithm was taken from glibc 2.3.5,
-               file sysdeps/generic/s_csqrt.c.
    @enddesc
-   @history
-               This code bears a resemblance to that in the GSL.  That
-	       code contains references.
-   @endhistory
 
-   @var        complex_number
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -625,31 +540,10 @@ cctk_complex CCTK_Cmplx##Log (cctk_complex complex_number)                    \
                The square root
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_SQRT(CCTK_Cmplx, cctk_real, cctk_complex, type)     \
-cctk_complex CCTK_Cmplx##Sqrt (cctk_complex complex_number)                   \
-{                                                                             \
-  cctk_real d, r, s;                                                          \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  d = hypot##type (complex_number.Re, complex_number.Im);                     \
-  /* Use the identity   2  Re res  Im res = Im x                              \
-     to avoid cancellation error in  d +/- Re x.  */                          \
-  if (complex_number.Re > 0)                                                  \
-  {                                                                           \
-    r = sqrt##type (0.5##type * d + 0.5##type * complex_number.Re);           \
-    s = (0.5##type * complex_number.Im) / r;                                  \
-  }                                                                           \
-  else                                                                        \
-  {                                                                           \
-    s = sqrt##type (0.5##type * d - 0.5##type * complex_number.Re);           \
-    r = fabs##type ((0.5##type * complex_number.Im) / s);                     \
-  }                                                                           \
-                                                                              \
-  result.Re = r;                                                              \
-  result.Im = copysign##type (s, complex_number.Im);                          \
-                                                                              \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_SQRT(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Sqrt (cctk_complex x)                          \
+{                                                                       \
+  return csqrt##type(x);                                                \
 }
 
 
@@ -659,11 +553,9 @@ cctk_complex CCTK_Cmplx##Sqrt (cctk_complex complex_number)                   \
    @author     Yaakoub Y El Khamra
    @desc
                Raises a complex number to a given power
-               This algorithm was taken from glibc 2.3.5,
-               file sysdeps/generic/s_cpow.c.
    @enddesc
  
-   @var        cx_number
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -679,18 +571,10 @@ cctk_complex CCTK_Cmplx##Sqrt (cctk_complex complex_number)                   \
                The power of the complex number
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_POW(CCTK_Cmplx, cctk_real, cctk_complex)            \
-cctk_complex CCTK_Cmplx##Pow (cctk_complex complex_number, cctk_real w)       \
-{                                                                             \
-  cctk_complex result;                                                        \
-                                                                              \
-                                                                              \
-  result = CCTK_Cmplx##Log (complex_number);                                  \
-  result.Re *= w;                                                             \
-  result.Im *= w;                                                             \
-  result = CCTK_Cmplx##Exp (result);                                          \
-                                                                              \
-  return (result);                                                            \
+#define DEFINE_CCTK_CMPLX_POW(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##Pow (cctk_complex x, cctk_real w)              \
+{                                                                       \
+  return cpow##type(x, w);                                              \
 }
 
 
@@ -700,11 +584,9 @@ cctk_complex CCTK_Cmplx##Pow (cctk_complex complex_number, cctk_real w)       \
    @author     Erik Schnetter
    @desc
                Raises a complex number to a given power
-               This algorithm was taken from glibc 2.3.5,
-               file sysdeps/generic/s_cpow.c.
    @enddesc
  
-   @var        cx_number
+   @var        x
    @vdesc      The complex number
    @vtype      CCTK_COMPLEX
    @vio        in
@@ -720,132 +602,68 @@ cctk_complex CCTK_Cmplx##Pow (cctk_complex complex_number, cctk_real w)       \
                The power of the complex number
    @endreturndesc
 @@*/
-#define DEFINE_CCTK_CMPLX_IPOW(CCTK_Cmplx, cctk_real, cctk_complex)     \
-cctk_complex CCTK_Cmplx##IPow (cctk_complex complex_number, int w)      \
+#define DEFINE_CCTK_CMPLX_IPOW(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+cctk_complex CCTK_Cmplx##IPow (cctk_complex x, int w)                   \
 {                                                                       \
-  cctk_complex result;                                                  \
-                                                                        \
-                                                                        \
-  result = CCTK_Cmplx (1, 0);                                           \
-  if (w < 0)                                                            \
+  int w0 = w;                                                           \
+  cctk_complex result = 1;                                              \
+  while (w)                                                             \
   {                                                                     \
-    result = CCTK_Cmplx##Div (result,                                   \
-                              CCTK_Cmplx##IPow (complex_number, - w));  \
-  }                                                                     \
-  else                                                                  \
-  {                                                                     \
-    while (w)                                                           \
+    if (w % 2)                                                          \
     {                                                                   \
-      if (w % 2)                                                        \
-      {                                                                 \
-        result = CCTK_Cmplx##Mul (result, complex_number);              \
-      }                                                                 \
-      w /= 2;                                                           \
-      complex_number = CCTK_Cmplx##Mul (complex_number, complex_number); \
+      result *= x;                                                      \
     }                                                                   \
+    w /= 2;                                                             \
+    x *= x;                                                             \
+  }                                                                     \
+  if (w0 < 0)                                                           \
+  {                                                                     \
+    result = 1 / result;                                                \
   }                                                                     \
-                                                                        \
-  return (result);                                                      \
+  return result;                                                        \
 }
 
 
 
-#ifdef DEFINE_CCTK_COMPLEX_INLINE_FUNCTIONS
-#  define PREFIX static inline
-#else
-#  define PREFIX
-#endif
-
-
-#ifdef __cplusplus
-  /* define C++ overloaded operators */
-#  define DEFINE_CMPLX_CXX_OPERATORS(CCTK_Cmplx, cctk_real, cctk_complex) \
-PREFIX cctk_complex operator+  (cctk_complex const & a)                         { return a; } \
-PREFIX cctk_complex operator-  (cctk_complex const & a)                         { return CCTK_Cmplx##Neg (a); } \
-PREFIX cctk_complex conj       (cctk_complex const & a)                         { return CCTK_Cmplx##Conjg (a); } \
-PREFIX cctk_real    abs        (cctk_complex const & a)                         { return CCTK_Cmplx##Abs (a); } \
-PREFIX cctk_real    arg        (cctk_complex const & a)                         { return CCTK_Cmplx##Arg (a); } \
-PREFIX cctk_real    norm       (cctk_complex const & a)                         { return CCTK_Cmplx##Norm (a); } \
-PREFIX cctk_complex operator+  (cctk_complex const & a, cctk_complex const & b) { return CCTK_Cmplx##Add (a, b); } \
-PREFIX cctk_complex operator+  (cctk_real    const & a, cctk_complex const & b) { return CCTK_Cmplx##Add (CCTK_Cmplx (a, 0), b); } \
-PREFIX cctk_complex operator+  (cctk_complex const & a, cctk_real    const & b) { return CCTK_Cmplx##Add (a, CCTK_Cmplx (b, 0)); } \
-PREFIX cctk_complex operator-  (cctk_complex const & a, cctk_complex const & b) { return CCTK_Cmplx##Sub (a, b); } \
-PREFIX cctk_complex operator-  (cctk_real    const & a, cctk_complex const & b) { return CCTK_Cmplx##Sub (CCTK_Cmplx (a, 0), b); } \
-PREFIX cctk_complex operator-  (cctk_complex const & a, cctk_real    const & b) { return CCTK_Cmplx##Sub (a, CCTK_Cmplx (b, 0)); } \
-PREFIX cctk_complex operator*  (cctk_complex const & a, cctk_complex const & b) { return CCTK_Cmplx##Mul (a, b); } \
-PREFIX cctk_complex operator*  (cctk_real    const & a, cctk_complex const & b) { return CCTK_Cmplx##Mul (CCTK_Cmplx (a, 0), b); } \
-PREFIX cctk_complex operator*  (cctk_complex const & a, cctk_real    const & b) { return CCTK_Cmplx##Mul (a, CCTK_Cmplx (b, 0)); } \
-PREFIX cctk_complex operator/  (cctk_complex const & a, cctk_complex const & b) { return CCTK_Cmplx##Div (a, b); } \
-PREFIX cctk_complex operator/  (cctk_real    const & a, cctk_complex const & b) { return CCTK_Cmplx##Div (CCTK_Cmplx (a, 0), b); } \
-PREFIX cctk_complex operator/  (cctk_complex const & a, cctk_real    const & b) { return CCTK_Cmplx##Div (a, CCTK_Cmplx (b, 0)); } \
-PREFIX cctk_complex operator+= (cctk_complex       & a, cctk_complex const & b) { return a = CCTK_Cmplx##Add (a, b); } \
-PREFIX cctk_complex operator+= (cctk_complex       & a, cctk_real    const & b) { return a = CCTK_Cmplx##Add (a, CCTK_Cmplx (b, 0)); } \
-PREFIX cctk_complex operator-= (cctk_complex       & a, cctk_complex const & b) { return a = CCTK_Cmplx##Sub (a, b); } \
-PREFIX cctk_complex operator-= (cctk_complex       & a, cctk_real    const & b) { return a = CCTK_Cmplx##Sub (a, CCTK_Cmplx (b, 0)); } \
-PREFIX cctk_complex operator*= (cctk_complex       & a, cctk_complex const & b) { return a = CCTK_Cmplx##Mul (a, b); } \
-PREFIX cctk_complex operator*= (cctk_complex       & a, cctk_real    const & b) { return a = CCTK_Cmplx##Mul (a, CCTK_Cmplx (b, 0)); } \
-PREFIX cctk_complex operator/= (cctk_complex       & a, cctk_complex const & b) { return a = CCTK_Cmplx##Div (a, b); } \
-PREFIX cctk_complex operator/= (cctk_complex       & a, cctk_real    const & b) { return a = CCTK_Cmplx##Div (a, CCTK_Cmplx (b, 0)); } \
-PREFIX cctk_complex pow        (cctk_complex const & a, int                  i) { return CCTK_Cmplx##IPow (a, i); } \
-PREFIX cctk_complex pow        (cctk_complex const & a, cctk_real    const & r) { return CCTK_Cmplx##Pow (a, r); } \
-PREFIX cctk_complex pow        (cctk_complex const & a, cctk_complex const & b) { return CCTK_Cmplx##CPow (a, b); } \
-PREFIX cctk_complex sin        (cctk_complex const & a)                         { return CCTK_Cmplx##Sin (a); } \
-PREFIX cctk_complex cos        (cctk_complex const & a)                         { return CCTK_Cmplx##Cos (a); } \
-PREFIX cctk_complex exp        (cctk_complex const & a)                         { return CCTK_Cmplx##Exp (a); } \
-PREFIX cctk_complex log        (cctk_complex const & a)                         { return CCTK_Cmplx##Log (a); } \
-PREFIX cctk_complex sqrt       (cctk_complex const & a)                         { return CCTK_Cmplx##Sqrt (a); } \
-PREFIX std::ostream & operator << (std::ostream & os, cctk_complex const & a) { return os << CCTK_Cmplx##Real (a) << " " << CCTK_Cmplx##Imag (a); }
-#else
-  /* define no C++ overloaded operators */
-#  define DEFINE_CMPLX_CXX_OPERATORS(CCTK_Cmplx, cctk_real, cctk_complex)
-#endif
-
-
 /* macro to define a set of complex functions for a given precision */
-#define DEFINE_CMPLX_FUNCTIONS(CCTK_Cmplx, cctk_real, cctk_complex, type)     \
-PREFIX DEFINE_CCTK_CMPLX       (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_REAL  (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_IMAG  (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_NEG   (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_CONJG (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_ABS   (CCTK_Cmplx, cctk_real, cctk_complex, type)    \
-PREFIX DEFINE_CCTK_CMPLX_NORM  (CCTK_Cmplx, cctk_real, cctk_complex, type)    \
-PREFIX DEFINE_CCTK_CMPLX_ARG   (CCTK_Cmplx, cctk_real, cctk_complex, type)    \
-PREFIX DEFINE_CCTK_CMPLX_ADD   (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_SUB   (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_MUL   (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_DIV   (CCTK_Cmplx, cctk_real, cctk_complex, type)    \
-PREFIX DEFINE_CCTK_CMPLX_CPOW  (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_SIN   (CCTK_Cmplx, cctk_real, cctk_complex, type)    \
-PREFIX DEFINE_CCTK_CMPLX_COS   (CCTK_Cmplx, cctk_real, cctk_complex, type)    \
-PREFIX DEFINE_CCTK_CMPLX_EXP   (CCTK_Cmplx, cctk_real, cctk_complex, type)    \
-PREFIX DEFINE_CCTK_CMPLX_LOG   (CCTK_Cmplx, cctk_real, cctk_complex, type)    \
-PREFIX DEFINE_CCTK_CMPLX_SQRT  (CCTK_Cmplx, cctk_real, cctk_complex, type)    \
-PREFIX DEFINE_CCTK_CMPLX_POW   (CCTK_Cmplx, cctk_real, cctk_complex)    \
-PREFIX DEFINE_CCTK_CMPLX_IPOW  (CCTK_Cmplx, cctk_real, cctk_complex)    \
-DEFINE_CMPLX_CXX_OPERATORS(CCTK_Cmplx, cctk_real, cctk_complex)
+#define DEFINE_CMPLX_FUNCTIONS(CCTK_Cmplx, cctk_real, cctk_complex, type) \
+DEFINE_CCTK_CMPLX       (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_REAL  (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_IMAG  (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_NEG   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_CONJG (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_ABS   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_NORM  (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_ARG   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_ADD   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_SUB   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_MUL   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_DIV   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_CPOW  (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_SIN   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_COS   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_EXP   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_LOG   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_SQRT  (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_POW   (CCTK_Cmplx, cctk_real, cctk_complex, type)     \
+DEFINE_CCTK_CMPLX_IPOW  (CCTK_Cmplx, cctk_real, cctk_complex, type)
 
 
 /* define complex functions for all available precisions */
-/* NOTE: The code in this file will function correctly with C float
- * and double types, but not with long double.
-*/
 #ifdef HAVE_CCTK_REAL4
-  #define _REAL_TYPE   f
-  DEFINE_CMPLX_FUNCTIONS (CCTK_Cmplx8, CCTK_REAL4, CCTK_COMPLEX8, _REAL_TYPE)
-  #undef _REAL_TYPE
+  #define KIND_SUFFIX f
+  DEFINE_CMPLX_FUNCTIONS (CCTK_Cmplx8, CCTK_REAL4, CCTK_COMPLEX8, KIND_SUFFIX)
+  #undef KIND_SUFFIX
 #endif
 
 #ifdef HAVE_CCTK_REAL8
-  #define _REAL_TYPE
-  DEFINE_CMPLX_FUNCTIONS (CCTK_Cmplx16, CCTK_REAL8, CCTK_COMPLEX16, _REAL_TYPE)
-  #undef _REAL_TYPE
+  #define KIND_SUFFIX
+  DEFINE_CMPLX_FUNCTIONS (CCTK_Cmplx16, CCTK_REAL8, CCTK_COMPLEX16, KIND_SUFFIX)
+  #undef KIND_SUFFIX
 #endif
 
 #ifdef HAVE_CCTK_REAL16
-  #define _REAL_TYPE   l
-  DEFINE_CMPLX_FUNCTIONS (CCTK_Cmplx32, CCTK_REAL16, CCTK_COMPLEX32, _REAL_TYPE)
-  #undef _REAL_TYPE
+  #define KIND_SUFFIX l
+  DEFINE_CMPLX_FUNCTIONS (CCTK_Cmplx32, CCTK_REAL16, CCTK_COMPLEX32, KIND_SUFFIX)
+  #undef KIND_SUFFIX
 #endif
-
-#undef PREFIX
diff --git a/src/util/Table.c b/src/util/Table.c
index 1c9a3db8..9a613a41 100644
--- a/src/util/Table.c
+++ b/src/util/Table.c
@@ -94,6 +94,7 @@ typedef int bool;
 #endif
 
 #include "cctk_Types.h"
+#include "cctk_Complex.h"
 #include "cctk_Constants.h"
 #include "cctk_Groups.h"
 #include "cctk_Flesh.h"
@@ -5591,8 +5592,8 @@ int Util_TablePrint(FILE *stream, int handle)
           for (i = 0 ; i < tep->N_elements ; ++i)
           {
             fprintf(stream, "\t(%g,%g)",
-                     (double) value_ptr_complex[i].Re,
-                     (double) value_ptr_complex[i].Im);
+                    (double) CCTK_CmplxReal(value_ptr_complex[i]),
+                    (double) CCTK_CmplxImag(value_ptr_complex[i]));
           }
             }
           break;
@@ -5729,8 +5730,8 @@ int Util_TablePrintPretty(FILE *stream, int handle)
             fprintf(stream, " ");
           }
           fprintf(stream, "(%#.17g,%#.17g)",
-                  (double)value_ptr_complex[i].Re,
-                  (double)value_ptr_complex[i].Im);
+                  (double)CCTK_CmplxReal(value_ptr_complex[i]),
+                  (double)CCTK_CmplxImag(value_ptr_complex[i]));
         }
         if (tep->N_elements != 1)
         {