numpy · mhvk · Apr 26, 2018 · Apr 24, 2018
diff --git a/doc/release/1.15.0-notes.rst b/doc/release/1.15.0-notes.rst
@@ -158,6 +158,12 @@ Added experimental support for the 64-bit RISC-V architecture.
 Improvements
 ============
 
+``np.ufunc.reduce`` and related functions now accept an initial value
+---------------------------------------------------------------------
+``np.ufunc.reduce``, ``np.sum``, ``np.prod``, ``np.min`` and ``np.max`` all
+now accept an ``initial`` keyword argument that specifies the value to start
+the reduction with.
+
 ``np.flip`` can operate over multiple axes
 ------------------------------------------
 ``np.flip`` now accepts None, or tuples of int, in its ``axis`` argument. If

diff --git a/numpy/add_newdocs.py b/numpy/add_newdocs.py
@@ -5850,7 +5850,7 @@ def luf(lamdaexpr, *args, **kwargs):
 
 add_newdoc('numpy.core', 'ufunc', ('reduce',
     """
-    reduce(a, axis=0, dtype=None, out=None, keepdims=False)
+    reduce(a, axis=0, dtype=None, out=None, keepdims=False, initial)
 
     Reduces `a`'s dimension by one, by applying ufunc along one axis.
 
@@ -5906,6 +5906,14 @@ def luf(lamdaexpr, *args, **kwargs):
         the result will broadcast correctly against the original `arr`.
 
         .. versionadded:: 1.7.0
+    initial : scalar, optional
+        The value with which to start the reduction.
+        If the ufunc has no identity or the dtype is object, this defaults
+        to None - otherwise it defaults to ufunc.identity.
+        If ``None`` is given, the first element of the reduction is used,
+        and an error is thrown if the reduction is empty.
+
+        .. versionadded:: 1.15.0
 
     Returns
     -------
@@ -5937,7 +5945,24 @@ def luf(lamdaexpr, *args, **kwargs):
     >>> np.add.reduce(X, 2)
     array([[ 1,  5],
            [ 9, 13]])
-
+
+    You can use the ``initial`` keyword argument to initialize the reduction with a
+    different value.
+
+    >>> np.add.reduce([10], initial=5)
+    15
+    >>> np.add.reduce(np.ones((2, 2, 2)), axis=(0, 2), initializer=10)
+    array([14., 14.])
+
+    Allows reductions of empty arrays where they would normally fail, i.e.
+    for ufuncs without an identity.
+
+    >>> np.minimum.reduce([], initial=np.inf)
+    inf
+    >>> np.minimum.reduce([])
+    Traceback (most recent call last):
+        ...
+    ValueError: zero-size array to reduction operation minimum which has no identity
     """))
 
 add_newdoc('numpy.core', 'ufunc', ('accumulate',

diff --git a/numpy/core/_methods.py b/numpy/core/_methods.py
@@ -11,6 +11,7 @@
 from numpy.core import umath as um
 from numpy.core.numeric import asanyarray
 from numpy.core import numerictypes as nt
+from numpy._globals import _NoValue
 
 # save those O(100) nanoseconds!
 umr_maximum = um.maximum.reduce
@@ -22,17 +23,21 @@
 
 # avoid keyword arguments to speed up parsing, saves about 15%-20% for very
 # small reductions
-def _amax(a, axis=None, out=None, keepdims=False):
-    return umr_maximum(a, axis, None, out, keepdims)
+def _amax(a, axis=None, out=None, keepdims=False,
+          initial=_NoValue):
+    return umr_maximum(a, axis, None, out, keepdims, initial)
 
-def _amin(a, axis=None, out=None, keepdims=False):
-    return umr_minimum(a, axis, None, out, keepdims)
+def _amin(a, axis=None, out=None, keepdims=False,
+          initial=_NoValue):
+    return umr_minimum(a, axis, None, out, keepdims, initial)
 
-def _sum(a, axis=None, dtype=None, out=None, keepdims=False):
-    return umr_sum(a, axis, dtype, out, keepdims)
+def _sum(a, axis=None, dtype=None, out=None, keepdims=False,
+         initial=_NoValue):
+    return umr_sum(a, axis, dtype, out, keepdims, initial)
 
-def _prod(a, axis=None, dtype=None, out=None, keepdims=False):
-    return umr_prod(a, axis, dtype, out, keepdims)
+def _prod(a, axis=None, dtype=None, out=None, keepdims=False,
+          initial=_NoValue):
+    return umr_prod(a, axis, dtype, out, keepdims, initial)
 
 def _any(a, axis=None, dtype=None, out=None, keepdims=False):
     return umr_any(a, axis, dtype, out, keepdims)

diff --git a/numpy/core/fromnumeric.py b/numpy/core/fromnumeric.py
@@ -1812,7 +1812,7 @@ def clip(a, a_min, a_max, out=None):
     return _wrapfunc(a, 'clip', a_min, a_max, out=out)
 
 
-def sum(a, axis=None, dtype=None, out=None, keepdims=np._NoValue):
+def sum(a, axis=None, dtype=None, out=None, keepdims=np._NoValue, initial=np._NoValue):
     """
     Sum of array elements over a given axis.
 
@@ -1851,6 +1851
8000
,10 @@ def sum(a, axis=None, dtype=None, out=None, keepdims=np._NoValue):
         `ndarray`, however any non-default value will be.  If the
         sub-class' method does not implement `keepdims` any
         exceptions will be raised.
+    initial : scalar, optional
+        Starting value for the sum. See `~numpy.ufunc.reduce` for details.
+
+        .. versionadded:: 1.15.0
 
     Returns
     -------
@@ -1898,6 +1902,10 @@ def sum(a, axis=None, dtype=None, out=None, keepdims=np._NoValue):
     >>> np.ones(128, dtype=np.int8).sum(dtype=np.int8)
     -128
 
+    You can also start the sum with a value other than zero:
+
+    >>> np.sum([10], initial=5)
+    15
     """
     if isinstance(a, _gentype):
         # 2018-02-25, 1.15.0
@@ -1912,7 +1920,8 @@ def sum(a, axis=None, dtype=None, out=None, keepdims=np._NoValue):
             return out
         return res
 
-    return _wrapreduction(a, np.add, 'sum', axis, dtype, out, keepdims=keepdims)
+    return _wrapreduction(a, np.add, 'sum', axis, dtype, out, keepdims=keepdims,
+                          initial=initial)
 
 
 def any(a, axis=None, out=None, keepdims=np._NoValue):
@@ -2209,7 +2218,7 @@ def ptp(a, axis=None, out=None, keepdims=np._NoValue):
     return _methods._ptp(a, axis=axis, out=out, **kwargs)
 
 
-def amax(a, axis=None, out=None, keepdims=np._NoValue):
+def amax(a, axis=None, out=None, keepdims=np._NoValue, initial=np._NoValue):
     """
     Return the maximum of an array or maximum along an axis.
 
@@ -2241,6 +2250,13 @@ def amax(a, axis=None, out=None, keepdims=np._NoValue):
         sub-class' method does not implement `keepdims` any
         exceptions will be raised.
 
+    initial : scalar, optional
+        The minimum value of an output element. Must be present to allow
+        computation on empty slice. See `~numpy.ufunc.reduce` for details.
+
+        .. versionadded:: 1.15.0
+
+
     Returns
     -------
     amax : ndarray or scalar
@@ -2293,11 +2309,26 @@ def amax(a, axis=None, out=None, keepdims=np._NoValue):
     >>> np.nanmax(b)
     4.0
 
+    You can use an initial value to compute the maximum of an empty slice, or
+    to initialize it to a different value:
+
+    >>> np.max([[-50], [10]], axis=-1, initial=0)
+    array([ 0, 10])
+
+    Notice that the initial value is used as one of the elements for which the
+    maximum is determined, unlike for the default argument Python's max
+    function, which is only used for empty iterables.
+
+    >>> np.max([5], initial=6)
+    6
+    >>> max([5], default=6)
+    5
     """
-    return _wrapreduction(a, np.maximum, 'max', axis, None, out, keepdims=keepdims)
+    return _wrapreduction(a, np.maximum, 'max', axis, None, out, keepdims=keepdims,
+                          initial=initial)
 
 
-def amin(a, axis=None, out=None, keepdims=np._NoValue):
+def amin(a, axis=None, out=None, keepdims=np._NoValue, initial=np._NoValue):
     """
     Return the minimum of an array or minimum along an axis.
 
@@ -2329,6 +2360,12 @@ def amin(a, axis=None, out=None, keepdims=np._NoValue):
         sub-class' method does not implement `keepdims` any
         exceptions will be raised.
 
+    initial : scalar, optional
+        The maximum value of an output element. Must be present to allow
+        computation on empty slice. See `~numpy.ufunc.reduce` for details.
+
+        .. versionadded:: 1.15.0
+
     Returns
     -------
     amin : ndarray or scalar
@@ -2381,8 +2418,22 @@ def amin(a, axis=None, out=None, keepdims=np._NoValue):
     >>> np.nanmin(b)
     0.0
 
+    >>> np.min([[-50], [10]], axis=-1, initial=0)
+    array([-50,   0])
+
+    Notice that the initial value is used as one of the elements for which the
+    minimum is determined, unlike for the default argument Python's max
+    function, which is only used for empty iterables.
+
+    Notice that this isn't the same as Python's ``default`` argument.
+
+    >>> np.min([6], initial=5)
+    5
+    >>> min([6], default=5)
+    6
     """
-    return _wrapreduction(a, np.minimum, 'min', axis, None, out, keepdims=keepdims)
+    return _wrapreduction(a, np.minimum, 'min', axis, None, out, keepdims=keepdims,
+                          initial=initial)
 
 
 def alen(a):
@@ -2418,7 +2469,7 @@ def alen(a):
         return len(array(a, ndmin=1))
 
 
-def prod(a, axis=None, dtype=None, out=None, keepdims=np._NoValue):
+def prod(a, axis=None, dtype=None, out=None, keepdims=np._NoValue, initial=np._NoValue):
     """
     Return the product of array elements over a given axis.
 
@@ -2458,6 +2509,10 @@ def prod(a, axis=None, dtype=None, out=None, keepdims=np._NoValue):
         `ndarray`, however any non-default value will be.  If the
         sub-class' method does not implement `keepdims` any
         exceptions will be raised.
+    initial : scalar, optional
+        The starting value for this product. See `~numpy.ufunc.reduce` for details.
+
+        .. versionadded:: 1.15.0
 
     Returns
     -------
@@ -2515,8 +2570,13 @@ def prod(a, axis=None, dtype=None, out=None, keepdims=np._NoValue):
     >>> np.prod(x).dtype == int
     True
 
+    You can also start the product with a value other than one:
+
+    >>> np.prod([1, 2], initial=5)
+    10
     """
-    return _wrapreduction(a, np.multiply, 'prod', axis, dtype, out, keepdims=keepdims)
+    return _wrapreduction(a, np.multiply, 'prod', axis, dtype, out, keepdims=keepdims,
+                          initial=initial)
 
 
 def cumprod(a, axis=None, dtype=None, out=None):

diff --git a/numpy/core/src/umath/override.c b/numpy/core/src/umath/override.c
@@ -123,11 +123,16 @@ normalize_reduce_args(PyUFuncObject *ufunc, PyObject *args,
     npy_intp nargs = PyTuple_GET_SIZE(args);
     npy_intp i;
     PyObject *obj;
-    static char *kwlist[] = {"array", "axis", "dtype", "out", "keepdims"};
+    static PyObject *NoValue = NULL;
+    static char *kwlist[] = {"array", "axis", "dtype", "out", "keepdims",
+        "initial"};
+
+    npy_cache_import("numpy", "_NoValue", &NoValue);
+    if (NoValue == NULL) return -1;
 
-    if (nargs < 1 || nargs > 5) {
+    if (nargs < 1 || nargs > 6) {
         PyErr_Format(PyExc_TypeError,
-                     "ufunc.reduce() takes from 1 to 5 positional "
+                     "ufunc.reduce() takes from 1 to 6 positional "
                      "arguments but %"NPY_INTP_FMT" were given", nargs);
         return -1;
     }
@@ -151,6 +156,10 @@ normalize_reduce_args(PyUFuncObject *ufunc, PyObject *args,
             }
             obj = PyTuple_GetSlice(args, 3, 4);
         }
+        /* Remove initial=np._NoValue */
+        if (i == 5 && obj == NoValue) {
+            continue;
+        }
         PyDict_SetItemString(*normal_kwds, kwlist[i], obj);
         if (i == 3) {
             Py_DECREF(obj);