typedef cumo_dcomplex dtype;
typedef double rtype;
#define cT  cumo_cDComplex
#define cRT cumo_cDFloat
#define mTM cumo_mDComplexMath

#include "complex_macro.h"
#include "cublas_v2.h"
#include "cumo/cuda/cublas.h"

static inline bool c_nearly_eq(dtype x, dtype y) {
    return c_abs(c_sub(x,y)) <= (c_abs(x)+c_abs(y))*DBL_EPSILON*2;
}

#ifdef SFMT_H
/* generates a random number on [0,1)-real-interval */
inline static dtype m_rand(dtype max)
{
    dtype z;
    CUMO_REAL(z) = genrand_res53_mix() * CUMO_REAL(max);
    CUMO_IMAG(z) = genrand_res53_mix() * CUMO_IMAG(max);
    return z;
}

/* generates random numbers from the normal distribution
   using Box-Muller Transformation.
 */
inline static void m_rand_norm(dtype mu, rtype sigma, dtype *a0)
{
    rtype x1, x2, w;
    do {
	x1 = genrand_res53_mix();
	x1 = x1*2-1;
	x2 = genrand_res53_mix();
	x2 = x2*2-1;
	w = x1 * x1 + x2 * x2;
    } while (w>=1);
    w = sqrt( (-2*log(w)) / w );
    CUMO_REAL(*a0) = x1*w * sigma + CUMO_REAL(mu);
    CUMO_IMAG(*a0) = x2*w * sigma + CUMO_IMAG(mu);
}
#endif

#define M_EPSILON rb_float_new(2.2204460492503131e-16)
#define M_MIN     rb_float_new(2.2250738585072014e-308)
#define M_MAX     rb_float_new(1.7976931348623157e+308)