/* -*- c -*- */

/*
 * vim:syntax=c
 */

/*
 *****************************************************************************
 **                            INCLUDES                                     **
 *****************************************************************************
 */

/*
 * _UMATHMODULE IS needed in __ufunc_api.h, included from numpy/ufuncobject.h.
 * This is a mess and it would be nice to fix it. It has nothing to do with
 * __ufunc_api.c
 */
#define _UMATHMODULE

#include "Python.h"
#include "numpy/noprefix.h"
#include "numpy/ufuncobject.h"
#include "abstract.h"
#include "config.h"

/*
 * Looks like some versions of Python.h do naughty things, so math.h needs
 * to come after.
 */
#include <math.h>

/*
 *****************************************************************************
 **                    INCLUDE GENERATED CODE                               **
 *****************************************************************************
 */
#include "umath_funcs_c99.inc"
#include "umath_funcs.inc"
#include "umath_loops.inc"
#include "umath_ufunc_object.inc"
#include "__umath_generated.c"
#include "__ufunc_api.c"


/*
 *****************************************************************************
 **                            SETUP UFUNCS                                 **
 *****************************************************************************
 */

/* Less automated additions to the ufuncs */

static PyUFuncGenericFunction frexp_functions[] = {
#ifdef HAVE_FREXPF
    FLOAT_frexp,
#endif
    DOUBLE_frexp
#ifdef HAVE_FREXPL
    ,LONGDOUBLE_frexp
#endif
};

static void * blank3_data[] = { (void *)NULL, (void *)NULL, (void *)NULL};
static char frexp_signatures[] = {
#ifdef HAVE_FREXPF
    PyArray_FLOAT, PyArray_FLOAT, PyArray_INT,
#endif
    PyArray_DOUBLE, PyArray_DOUBLE, PyArray_INT
#ifdef HAVE_FREXPL
    ,PyArray_LONGDOUBLE, PyArray_LONGDOUBLE, PyArray_INT
#endif
};

static PyUFuncGenericFunction ldexp_functions[] = {
#ifdef HAVE_LDEXPF
    FLOAT_ldexp,
#endif
    DOUBLE_ldexp
#ifdef HAVE_LDEXPL
    ,LONGDOUBLE_ldexp
#endif
};

static char ldexp_signatures[] = {
#ifdef HAVE_LDEXPF
    PyArray_FLOAT, PyArray_INT, PyArray_FLOAT,
#endif
    PyArray_DOUBLE, PyArray_INT, PyArray_DOUBLE
#ifdef HAVE_LDEXPL
    ,PyArray_LONGDOUBLE, PyArray_INT, PyArray_LONGDOUBLE
#endif
};

static void
InitOtherOperators(PyObject *dictionary) {
    PyObject *f;
    int num=1;

#ifdef HAVE_FREXPL
    num += 1;
#endif
#ifdef HAVE_FREXPF
    num += 1;
#endif
    f = PyUFunc_FromFuncAndData(frexp_functions, blank3_data,
                                frexp_signatures, num,
                                1, 2, PyUFunc_None, "frexp",
                                "Split the number, x, into a normalized"\
                                " fraction (y1) and exponent (y2)",0);
    PyDict_SetItemString(dictionary, "frexp", f);
    Py_DECREF(f);

    num = 1;
#ifdef HAVE_LDEXPL
    num += 1;
#endif
#ifdef HAVE_LDEXPF
    num += 1;
#endif
    f = PyUFunc_FromFuncAndData(ldexp_functions, blank3_data, ldexp_signatures, num,
                                2, 1, PyUFunc_None, "ldexp",
                                "Compute y = x1 * 2**x2.",0);
    PyDict_SetItemString(dictionary, "ldexp", f);
    Py_DECREF(f);
    return;
}

/* Setup +inf and +0 */

static double
pinf_init(void)
{
    double mul = 1e10;
    double tmp = 0.0;
    double pinf;

    pinf = mul;
    for (;;) {
        pinf *= mul;
        if (pinf == tmp) break;
        tmp = pinf;
    }
    return pinf;
}

static double
pzero_init(void)
{
    double div = 1e10;
    double tmp = 0.0;
    double pinf;

    pinf = div;
    for (;;) {
        pinf /= div;
        if (pinf == tmp) break;
        tmp = pinf;
    }
    return pinf;
}

/* Setup the umath module */
/* Remove for time being, it is declared in __ufunc_api.h */
/*static PyTypeObject PyUFunc_Type;*/

static struct PyMethodDef methods[] = {
    {"frompyfunc", (PyCFunction) ufunc_frompyfunc,
     METH_VARARGS | METH_KEYWORDS, doc_frompyfunc},
    {"seterrobj", (PyCFunction) ufunc_seterr,
     METH_VARARGS, NULL},
    {"geterrobj", (PyCFunction) ufunc_geterr,
     METH_VARARGS, NULL},
    {NULL, NULL, 0, NULL}                /* sentinel */
};

PyMODINIT_FUNC initumath(void) {
    PyObject *m, *d, *s, *s2, *c_api;
    double pinf, pzero, mynan;
    int UFUNC_FLOATING_POINT_SUPPORT = 1;

#ifdef NO_UFUNC_FLOATING_POINT_SUPPORT
    UFUNC_FLOATING_POINT_SUPPORT = 0;
#endif
    /* Create the module and add the functions */
    m = Py_InitModule("umath", methods);

    /* Import the array */
    if (_import_array() < 0) {
        if (!PyErr_Occurred()) {
            PyErr_SetString(PyExc_ImportError,
                            "umath failed: Could not import array core.");
        }
        return;
    }

    /* Initialize the types */
    if (PyType_Ready(&PyUFunc_Type) < 0)
        return;

    /* Add some symbolic constants to the module */
    d = PyModule_GetDict(m);

    c_api = PyCObject_FromVoidPtr((void *)PyUFunc_API, NULL);
    if (PyErr_Occurred()) goto err;
    PyDict_SetItemString(d, "_UFUNC_API", c_api);
    Py_DECREF(c_api);
    if (PyErr_Occurred()) goto err;

    s = PyString_FromString("0.4.0");
    PyDict_SetItemString(d, "__version__", s);
    Py_DECREF(s);

    /* Load the ufunc operators into the array module's namespace */
    InitOperators(d);

    InitOtherOperators(d);

    PyDict_SetItemString(d, "pi", s = PyFloat_FromDouble(M_PI));
    Py_DECREF(s);
    PyDict_SetItemString(d, "e", s = PyFloat_FromDouble(exp(1.0)));
    Py_DECREF(s);

#define ADDCONST(str) PyModule_AddIntConstant(m, #str, UFUNC_##str)
#define ADDSCONST(str) PyModule_AddStringConstant(m, "UFUNC_" #str, UFUNC_##str)

    ADDCONST(ERR_IGNORE);
    ADDCONST(ERR_WARN);
    ADDCONST(ERR_CALL);
    ADDCONST(ERR_RAISE);
    ADDCONST(ERR_PRINT);
    ADDCONST(ERR_LOG);
    ADDCONST(ERR_DEFAULT);
    ADDCONST(ERR_DEFAULT2);

    ADDCONST(SHIFT_DIVIDEBYZERO);
    ADDCONST(SHIFT_OVERFLOW);
    ADDCONST(SHIFT_UNDERFLOW);
    ADDCONST(SHIFT_INVALID);

    ADDCONST(FPE_DIVIDEBYZERO);
    ADDCONST(FPE_OVERFLOW);
    ADDCONST(FPE_UNDERFLOW);
    ADDCONST(FPE_INVALID);

    ADDCONST(FLOATING_POINT_SUPPORT);

    ADDSCONST(PYVALS_NAME);

#undef ADDCONST
#undef ADDSCONST
    PyModule_AddIntConstant(m, "UFUNC_BUFSIZE_DEFAULT", (long)PyArray_BUFSIZE);

    pinf = pinf_init();
    pzero = pzero_init();
    mynan = pinf / pinf;

    PyModule_AddObject(m, "PINF", PyFloat_FromDouble(pinf));
    PyModule_AddObject(m, "NINF", PyFloat_FromDouble(-pinf));
    PyModule_AddObject(m, "PZERO", PyFloat_FromDouble(pzero));
    PyModule_AddObject(m, "NZERO", PyFloat_FromDouble(-pzero));
    PyModule_AddObject(m, "NAN", PyFloat_FromDouble(mynan));

    s = PyDict_GetItemString(d, "conjugate");
    s2 = PyDict_GetItemString(d, "remainder");
    /* Setup the array object's numerical structures with appropriate
       ufuncs in d*/
    PyArray_SetNumericOps(d);

    PyDict_SetItemString(d, "conj", s);
    PyDict_SetItemString(d, "mod", s2);

    return;

 err:
    /* Check for errors */
    if (!PyErr_Occurred()) {
        PyErr_SetString(PyExc_RuntimeError,
                        "cannot load umath module.");
    }
    return;
}