--- Python-2.6a3/Modules/cmathmoduleORIG.c 2008-04-20 11:30:05.000000000 -0700 +++ Python-2.6a3/Modules/cmathmodule.c 2008-05-10 16:34:28.000000000 -0700 @@ -1048,6 +1048,165 @@ Checks if the real or imaginary part of z is infinite."); +/* Summation function like function msum() in this Python recipe + + + See function math_sum() in file ./mathmodule.c for more details + and for the original code. Make any changes in both places. +*/ + +static int +_cmath_sum_add(double x, double y, double *s_ptr) +{ + double s; + + errno = 0; + PyFPE_START_PROTECT("sum", return -1) + s = x + y; + PyFPE_END_PROTECT(s) + /* like math_2() in mathmodule.c */ + if (Py_IS_NAN(s)) { + if (!Py_IS_NAN(x) && !Py_IS_NAN(y)) + errno = EDOM; + else + errno = 0; + } + else if (Py_IS_INFINITY(s)) { + if (Py_IS_FINITE(x) && Py_IS_FINITE(y)) + errno = ERANGE; + else + errno = 0; + } + if (errno != 0 && math_error() == NULL) + return -1; + + *s_ptr = s; + return 0; +} + +static int +_cmath_sum_partial(double x, double **p_ptr, Py_ssize_t *n_ptr, + double *ps, Py_ssize_t *m_ptr) +{ + Py_ssize_t i, j, n = *n_ptr, m = *m_ptr; + double y, hi, lo, *p = *p_ptr; + + /* i = 0; for y in partials */ + for (i = j = 0; j < n; j++) { + y = p[j]; + if (_cmath_sum_add(x, y, &hi) < 0) + return -1; + lo = fabs(x) < fabs(y) + ? x - (hi - y) + : y - (hi - x); + if (lo != 0.0) + p[i++] = lo; + x = hi; + } + /* ps[i:] = [x] */ + n = i; + if (n >= m) { /* extend partials */ + void *v = NULL; + m += m / 2; /* plus 50% */ + if (n < m && m < (PY_SSIZE_T_MAX / sizeof(double))) { + if (p == ps) { + v = PyMem_Malloc(sizeof(double) * m); + if (v != NULL) + memcpy(v, ps, sizeof(double) * n); + } else + v = PyMem_Realloc(p, sizeof(double) * m); + } + if (v == NULL) { /* no memory or size overflow */ + PyErr_SetString(PyExc_MemoryError, "math sum error"); + return -1; + } + *p_ptr = p = (double*) v; + *m_ptr = m; + } + p[n++] = x; + *n_ptr = n; + return 0; +} + +static PyObject* +cmath_sum(PyObject *self, PyObject *args) +{ +# define CMATH_SUM_PS 128 /* static partials */ + + Py_complex c; + Py_ssize_t i, n_r = 0, m_r = CMATH_SUM_PS, + n_i = 0, m_i = CMATH_SUM_PS; + double ps_r[CMATH_SUM_PS], *p_r = ps_r, + ps_i[CMATH_SUM_PS], *p_i = ps_i; + PyObject *seq, *iter, *item, *sum = NULL; + + if (!PyArg_UnpackTuple(args, "sum", 1, 2, &seq, &sum)) + return NULL; + + iter = PyObject_GetIter(seq); + if (iter == NULL) + return NULL; + + if (sum != NULL) { + c = PyComplex_AsCComplex(sum); + sum = NULL; + if (PyErr_Occurred()) + goto cmath_sum_exit; + p_r[n_r++] = c.real; + p_i[n_i++] = c.imag; + } else + p_r[0] = p_i[0] = 0.0; + + /* for x in iterable */ + for(;;) { + item = PyIter_Next(iter); + if (item == NULL) { + if (PyErr_Occurred()) + goto cmath_sum_exit; + else + break; + } + c = PyComplex_AsCComplex(item); + Py_DECREF(item); + if (PyErr_Occurred()) + goto cmath_sum_exit; + /* for y in partials */ + if (_cmath_sum_partial(c.real, &p_r, &n_r, ps_r, &m_r) < 0 || + _cmath_sum_partial(c.imag, &p_i, &n_i, ps_i, &m_i) < 0) + goto cmath_sum_exit; + } + + /* sum(ps, 0.0) */ + c.real = p_r[0]; /* always valid */ + for (i = 1; i < n_r; i++) { + if (_cmath_sum_add(p_r[i], c.real, &c.real) < 0) + goto cmath_sum_exit; + } + c.imag = p_i[0]; /* always valid */ + for (i = 1; i < n_i; i++) { + if (_cmath_sum_add(p_i[i], c.imag, &c.imag) < 0) + goto cmath_sum_exit; + } + sum = PyComplex_FromCComplex(c); + +cmath_sum_exit: + Py_DECREF(iter); + if (p_r != ps_r) + PyMem_Free(p_r); + if (p_i != ps_i) + PyMem_Free(p_i); + return sum; + +# undef CMATH_SUM_PS +} + +PyDoc_STRVAR(cmath_sum_doc, +"sum(sequence[, start=0])\n" +"\n" +"Returns the full precision complex sum of a sequence of complex numbers plus the\n" +"value of parameter 'start'. When the sequence is empty, return complex start."); + + PyDoc_STRVAR(module_doc, "This module is always available. It provides access to mathematical\n" "functions for complex numbers."); @@ -1072,6 +1231,7 @@ {"sin", cmath_sin, METH_VARARGS, c_sin_doc}, {"sinh", cmath_sinh, METH_VARARGS, c_sinh_doc}, {"sqrt", cmath_sqrt, METH_VARARGS, c_sqrt_doc}, + {"sum", cmath_sum, METH_VARARGS, cmath_sum_doc}, {"tan", cmath_tan, METH_VARARGS, c_tan_doc}, {"tanh", cmath_tanh, METH_VARARGS, c_tanh_doc}, {NULL, NULL} /* sentinel */