Merge branch 'master' into force-tuple

11 files changed, 212 insertions, 78 deletions

diff --git a/doc/source/reference/arrays.classes.rst b/doc/source/reference/arrays.classes.rst
index 2719f9239..f17cb932a 100644
--- a/doc/source/reference/arrays.classes.rst
+++ b/doc/source/reference/arrays.classes.rst
@@ -215,6 +215,13 @@ Matrix objects
 .. index::
    single: matrix
 
+.. note::
+   It is strongly advised *not* to use the matrix subclass.  As described
+   below, it makes writing functions that deal consistently with matrices
+   and regular arrays very difficult. Currently, they are mainly used for
+   interacting with ``scipy.sparse``. We hope to provide an alternative
+   for this use, however, and eventually remove the ``matrix`` subclass.
+
 :class:`matrix` objects inherit from the ndarray and therefore, they
 have the same attributes and methods of ndarrays. There are six
 important differences of matrix objects, however, that may lead to
diff --git a/doc/source/reference/arrays.nditer.rst b/doc/source/reference/arrays.nditer.rst
index 76f5991cf..239f4296b 100644
--- a/doc/source/reference/arrays.nditer.rst
+++ b/doc/source/reference/arrays.nditer.rst
@@ -78,20 +78,28 @@ order='C' for C order and order='F' for Fortran order.
     ...
     0 3 1 4 2 5
 
+.. _nditer-context-manager:
+
 Modifying Array Values
 ----------------------
 
-By default, the :class:`nditer` treats the input array as a read-only
-object. To modify the array elements, you must specify either read-write
-or write-only mode. This is controlled with per-operand flags.
+By default, the :class:`nditer` treats the input operand as a read-only
+object. To be able to modify the array elements, you must specify either
+read-write or write-only mode using the `'readwrite'` or `'writeonly'`
+per-operand flags.
+
+The nditer will then yield writeable buffer arrays which you may modify. However,
+because  the nditer must copy this buffer data back to the original array once
+iteration is finished, you must signal when the iteration is ended, by one of two
+methods. You may either:
 
-Regular assignment in Python simply changes a reference in the local or
-global variable dictionary instead of modifying an existing variable in
-place.  This means that simply assigning to `x` will not place the value
-into the element of the array, but rather switch `x` from being an array
-element reference to being a reference to the value you assigned. To
-actually modify the element of the array, `x` should be indexed with
-the ellipsis.
+ - used the nditer as a context manager using the `with` statement, and
+   the temporary data will be written back when the context is exited.
+ - call the iterator's `close` method once finished iterating, which will trigger
+   the write-back.
+
+The nditer can no longer be iterated once either `close` is called or its
+context is exited.
 
 .. admonition:: Example
 
@@ -99,8 +107,9 @@ the ellipsis.
     >>> a
     array([[0, 1, 2],
            [3, 4, 5]])
-    >>> for x in np.nditer(a, op_flags=['readwrite']):
-    ...     x[...] = 2 * x
+    >>> with np.nditer(a, op_flags=['readwrite']) as it:
+    ...    for x in it:
+    ...        x[...] = 2 * x
     ...
     >>> a
     array([[ 0,  2,  4],
@@ -178,9 +187,10 @@ construct in order to be more readable.
     0 <(0, 0)> 1 <(0, 1)> 2 <(0, 2)> 3 <(1, 0)> 4 <(1, 1)> 5 <(1, 2)>
 
     >>> it = np.nditer(a, flags=['multi_index'], op_flags=['writeonly'])
-    >>> while not it.finished:
-    ...     it[0] = it.multi_index[1] - it.multi_index[0]
-    ...     it.iternext()
+    >>> with it:
+    ....    while not it.finished:
+    ...         it[0] = it.multi_index[1] - it.multi_index[0]
+    ...         it.iternext()
     ...
     >>> a
     array([[ 0,  1,  2],
@@ -385,10 +395,10 @@ parameter support.
 .. admonition:: Example
 
     >>> def square(a):
-    ...     it = np.nditer([a, None])
-    ...     for x, y in it:
-    ...          y[...] = x*x
-    ...     return it.operands[1]
+    ...     with np.nditer([a, None]) as it:
+    ...         for x, y in it:
+    ...             y[...] = x*x
+    ...         return it.operands[1]
     ...
     >>> square([1,2,3])
     array([1, 4, 9])
@@ -426,9 +436,10 @@ reasons.
     ...             flags = ['external_loop', 'buffered'],
     ...             op_flags = [['readonly'],
     ...                         ['writeonly', 'allocate', 'no_broadcast']])
-    ...     for x, y in it:
-    ...         y[...] = x*x
-    ...     return it.operands[1]
+    ...     with it:
+    ...         for x, y in it:
+    ...             y[...] = x*x
+    ...         return it.operands[1]
     ...
 
     >>> square([1,2,3])
@@ -480,10 +491,12 @@ Everything to do with the outer product is handled by the iterator setup.
     >>> b = np.arange(8).reshape(2,4)
     >>> it = np.nditer([a, b, None], flags=['external_loop'],
     ...             op_axes=[[0, -1, -1], [-1, 0, 1], None])
-    >>> for x, y, z in it:
-    ...     z[...] = x*y
+    >>> with it:
+    ...     for x, y, z in it:
+    ...         z[...] = x*y
+    ...     result = it.operands[2]  # same as z
     ...
-    >>> it.operands[2]
+    >>> result
     array([[[ 0,  0,  0,  0],
             [ 0,  0,  0,  0]],
            [[ 0,  1,  2,  3],
@@ -491,6 +504,9 @@ Everything to do with the outer product is handled by the iterator setup.
            [[ 0,  2,  4,  6],
             [ 8, 10, 12, 14]]])
 
+Note that once the iterator is closed we can not access :func:`operands <nditer.operands>`
+and must use a reference created inside the context manager.
+
 Reduction Iteration
 -------------------
 
@@ -505,9 +521,10 @@ For a simple example, consider taking the sum of all elements in an array.
 
     >>> a = np.arange(24).reshape(2,3,4)
     >>> b = np.array(0)
-    >>> for x, y in np.nditer([a, b], flags=['reduce_ok', 'external_loop'],
-    ...                     op_flags=[['readonly'], ['readwrite']]):
-    ...     y[...] += x
+    >>> with np.nditer([a, b], flags=['reduce_ok', 'external_loop'],
+    ...                     op_flags=[['readonly'], ['readwrite']]) as it:
+    ...     for x,y in it:
+    ...         y[...] += x
     ...
     >>> b
     array(276)
@@ -525,11 +542,13 @@ sums along the last axis of `a`.
     >>> it = np.nditer([a, None], flags=['reduce_ok', 'external_loop'],
     ...             op_flags=[['readonly'], ['readwrite', 'allocate']],
     ...             op_axes=[None, [0,1,-1]])
-    >>> it.operands[1][...] = 0
-    >>> for x, y in it:
-    ...     y[...] += x
+    >>> with it:
+    ...     it.operands[1][...] = 0
+    ...     for x, y in it:
+    ...         y[...] += x
+    ...     result = it.operands[1]
     ...
-    >>> it.operands[1]
+    >>> result
     array([[ 6, 22, 38],
            [54, 70, 86]])
     >>> np.sum(a, axis=2)
@@ -558,12 +577,14 @@ buffering.
     ...                                  'buffered', 'delay_bufalloc'],
     ...             op_flags=[['readonly'], ['readwrite', 'allocate']],
     ...             op_axes=[None, [0,1,-1]])
-    >>> it.operands[1][...] = 0
-    >>> it.reset()
-    >>> for x, y in it:
-    ...     y[...] += x
+    >>> with it:
+    ...     it.operands[1][...] = 0
+    ...     it.reset()
+    ...     for x, y in it:
+    ...         y[...] += x
+    ...     result = it.operands[1]
     ...
-    >>> it.operands[1]
+    >>> result
     array([[ 6, 22, 38],
            [54, 70, 86]])
 
@@ -609,11 +630,12 @@ Here's how this looks.
     ...                 op_flags=[['readonly'], ['readwrite', 'allocate']],
     ...                 op_axes=[None, axeslist],
     ...                 op_dtypes=['float64', 'float64'])
-    ...     it.operands[1][...] = 0
-    ...     it.reset()
-    ...     for x, y in it:
-    ...         y[...] += x*x
-    ...     return it.operands[1]
+    ...     with it:
+    ...         it.operands[1][...] = 0
+    ...         it.reset()
+    ...         for x, y in it:
+    ...             y[...] += x*x
+    ...         return it.operands[1]
     ...
     >>> a = np.arange(6).reshape(2,3)
     >>> sum_squares_py(a)
@@ -661,16 +683,17 @@ Here's the listing of sum_squares.pyx::
                     op_flags=[['readonly'], ['readwrite', 'allocate']],
                     op_axes=[None, axeslist],
                     op_dtypes=['float64', 'float64'])
-        it.operands[1][...] = 0
-        it.reset()
-        for xarr, yarr in it:
-            x = xarr
-            y = yarr
-            size = x.shape[0]
-            for i in range(size):
-               value = x[i]
-               y[i] = y[i] + value * value
-        return it.operands[1]
+        with it:
+            it.operands[1][...] = 0
+            it.reset()
+            for xarr, yarr in it:
+                x = xarr
+                y = yarr
+                size = x.shape[0]
+                for i in range(size):
+                   value = x[i]
+                   y[i] = y[i] + value * value
+            return it.operands[1]
 
 On this machine, building the .pyx file into a module looked like the
 following, but you may have to find some Cython tutorials to tell you
diff --git a/doc/source/reference/c-api.array.rst b/doc/source/reference/c-api.array.rst
index ad7c725a8..5ea7bfcfc 100644
--- a/doc/source/reference/c-api.array.rst
+++ b/doc/source/reference/c-api.array.rst
@@ -1360,7 +1360,7 @@ Special functions for NPY_OBJECT
 .. c:function:: int PyArray_SetWritebackIfCopyBase(PyArrayObject* arr, PyArrayObject* base)
 
     Precondition: ``arr`` is a copy of ``base`` (though possibly with different
-    strides, ordering, etc.) Sets the :c:data:`NPY_ARRAY_WRITEBACKIFCOPY` flag 
+    strides, ordering, etc.) Sets the :c:data:`NPY_ARRAY_WRITEBACKIFCOPY` flag
     and ``arr->base``, and set ``base`` to READONLY. Call
     :c:func:`PyArray_ResolveWritebackIfCopy` before calling
     `Py_DECREF`` in order copy any changes back to ``base`` and
@@ -3260,12 +3260,14 @@ Memory management
 .. c:function:: int PyArray_ResolveWritebackIfCopy(PyArrayObject* obj)
 
     If ``obj.flags`` has :c:data:`NPY_ARRAY_WRITEBACKIFCOPY` or (deprecated)
-    :c:data:`NPY_ARRAY_UPDATEIFCOPY`, this function copies ``obj->data`` to
-    `obj->base->data`, clears the flags, `DECREF` s `obj->base` and makes it
-    writeable, and sets ``obj->base`` to NULL. This is the opposite of 
+    :c:data:`NPY_ARRAY_UPDATEIFCOPY`, this function clears the flags, `DECREF` s
+    `obj->base` and makes it writeable, and sets ``obj->base`` to NULL. It then
+    copies ``obj->data`` to `obj->base->data`, and returns the error state of
+    the copy operation. This is the opposite of
     :c:func:`PyArray_SetWritebackIfCopyBase`. Usually this is called once
     you are finished with ``obj``, just before ``Py_DECREF(obj)``. It may be called
-    multiple times, or with ``NULL`` input.
+    multiple times, or with ``NULL`` input. See also
+    :c:func:`PyArray_DiscardWritebackIfCopy`.
 
     Returns 0 if nothing was done, -1 on error, and 1 if action was taken.
 
@@ -3487,12 +3489,14 @@ Miscellaneous Macros
 
 .. c:function:: PyArray_DiscardWritebackIfCopy(PyObject* obj)
 
-    Reset the :c:data:`NPY_ARRAY_WRITEBACKIFCOPY` and deprecated
-    :c:data:`NPY_ARRAY_UPDATEIFCOPY` flag. Resets the
-    :c:data:`NPY_ARRAY_WRITEABLE` flag on the base object. It also
-    discards pending changes to the base object. This is
-    useful for recovering from an error condition when
-    writeback semantics are used.
+    If ``obj.flags`` has :c:data:`NPY_ARRAY_WRITEBACKIFCOPY` or (deprecated)
+    :c:data:`NPY_ARRAY_UPDATEIFCOPY`, this function clears the flags, `DECREF` s
+    `obj->base` and makes it writeable, and sets ``obj->base`` to NULL. In
+    contrast to :c:func:`PyArray_DiscardWritebackIfCopy` it makes no attempt
+    to copy the data from `obj->base` This undoes
+    :c:func:`PyArray_SetWritebackIfCopyBase`. Usually this is called after an
+    error when you are finished with ``obj``, just before ``Py_DECREF(obj)``.
+    It may be called multiple times, or with ``NULL`` input.
 
 .. c:function:: PyArray_XDECREF_ERR(PyObject* obj)
 
diff --git a/doc/source/reference/c-api.coremath.rst b/doc/source/reference/c-api.coremath.rst
index d3f7fcf75..ad92235da 100644
--- a/doc/source/reference/c-api.coremath.rst
+++ b/doc/source/reference/c-api.coremath.rst
@@ -183,14 +183,46 @@ Those can be useful for precise floating point comparison.
     * NPY_FPE_UNDERFLOW
     * NPY_FPE_INVALID
 
+    Note that :c:func:`npy_get_floatstatus_barrier` is preferable as it prevents
+    agressive compiler optimizations reordering the call relative to
+    the code setting the status, which could lead to incorrect results.
+
     .. versionadded:: 1.9.0
 
+.. c:function:: int npy_get_floatstatus_barrier(char*)
+
+    Get floating point status. A pointer to a local variable is passed in to
+    prevent aggresive compiler optimizations from reodering this function call
+    relative to the code setting the status, which could lead to incorrect
+    results.
+
+    Returns a bitmask with following possible flags:
+
+    * NPY_FPE_DIVIDEBYZERO
+    * NPY_FPE_OVERFLOW
+    * NPY_FPE_UNDERFLOW
+    * NPY_FPE_INVALID
+
+    .. versionadded:: 1.15.0
+
 .. c:function:: int npy_clear_floatstatus()
 
     Clears the floating point status. Returns the previous status mask.
 
+    Note that :c:func:`npy_clear_floatstatus_barrier` is preferable as it
+    prevents agressive compiler optimizations reordering the call relative to
+    the code setting the status, which could lead to incorrect results.
+
     .. versionadded:: 1.9.0
 
+.. c:function:: int npy_clear_floatstatus_barrier(char*)
+
+    Clears the floating point status. A pointer to a local variable is passed in to
+    prevent aggresive compiler optimizations from reodering this function call.
+    Returns the previous status mask.
+
+    .. versionadded:: 1.15.0
+n
 Complex functions
 ~~~~~~~~~~~~~~~~~
 
diff --git a/doc/source/reference/c-api.iterator.rst b/doc/source/reference/c-api.iterator.rst
index 4c59bce51..392dcb730 100644
--- a/doc/source/reference/c-api.iterator.rst
+++ b/doc/source/reference/c-api.iterator.rst
@@ -110,6 +110,7 @@ number of non-zero elements in an array.
             /* Increment the iterator to the next inner loop */
         } while(iternext(iter));
 
+        NpyIter_Close(iter) /* best practice, not strictly required in this case */
         NpyIter_Deallocate(iter);
 
         return nonzero_count;
@@ -194,6 +195,7 @@ is used to control the memory layout of the allocated result, typically
         ret = NpyIter_GetOperandArray(iter)[1];
         Py_INCREF(ret);
 
+        NpyIter_Close(iter);
         if (NpyIter_Deallocate(iter) != NPY_SUCCEED) {
             Py_DECREF(ret);
             return NULL;
@@ -490,7 +492,10 @@ Construction and Destruction
 
             Indicate how the user of the iterator will read or write
             to ``op[i]``.  Exactly one of these flags must be specified
-            per operand.
+            per operand. Using ``NPY_ITER_READWRITE`` or ``NPY_ITER_WRITEONLY``
+            for a user-provided operand may trigger `WRITEBACKIFCOPY``
+            semantics. The data will be written back to the original array
+            when ``NpyIter_Close`` is called.
 
         .. c:var:: NPY_ITER_COPY
 
@@ -502,12 +507,12 @@ Construction and Destruction
 
             Triggers :c:data:`NPY_ITER_COPY`, and when an array operand
             is flagged for writing and is copied, causes the data
-            in a copy to be copied back to ``op[i]`` when the iterator
-            is destroyed.
+            in a copy to be copied back to ``op[i]`` when ``NpyIter_Close`` is
+            called.
 
             If the operand is flagged as write-only and a copy is needed,
             an uninitialized temporary array will be created and then copied
-            to back to ``op[i]`` on destruction, instead of doing
+            to back to ``op[i]`` on calling ``NpyIter_Close``, instead of doing
             the unnecessary copy operation.
 
         .. c:var:: NPY_ITER_NBO
@@ -704,6 +709,10 @@ Construction and Destruction
     the functions will pass back errors through it instead of setting
     a Python exception.
 
+    :c:func:`NpyIter_Deallocate` must be called for each copy. One call to
+    :c:func:`NpyIter_Close` is sufficient to trigger writeback resolution for
+    all copies since they share buffers.
+
 .. c:function:: int NpyIter_RemoveAxis(NpyIter* iter, int axis)``
 
     Removes an axis from iteration.  This requires that
@@ -754,10 +763,23 @@ Construction and Destruction
 
     Returns ``NPY_SUCCEED`` or ``NPY_FAIL``.
 
+.. c:function:: int NpyIter_Close(NpyIter* iter)
+
+    Resolves any needed writeback resolution. Should be called before
+    :c:func::`NpyIter_Deallocate`. After this call it is not safe to use the operands.
+    When using :c:func:`NpyIter_Copy`, only one call to :c:func:`NpyIter_Close`
+    is sufficient to resolve any writebacks, since the copies share buffers.
+
+    Returns ``0`` or ``-1`` if unsuccessful.
+
 .. c:function:: int NpyIter_Deallocate(NpyIter* iter)
 
-    Deallocates the iterator object.  This additionally frees any
-    copies made, triggering UPDATEIFCOPY behavior where necessary.
+    Deallocates the iterator object.  
+
+    :c:func:`NpyIter_Close` should be called before this. If not, and if
+    writeback is needed, it will be performed at this point in order to maintain
+    backward-compatibility with older code, and a deprecation warning will be
+    emitted. Old code should be updated to call `NpyIter_Close` beforehand.
 
     Returns ``NPY_SUCCEED`` or ``NPY_FAIL``.
 
diff --git a/doc/source/reference/constants.rst b/doc/source/reference/constants.rst
new file mode 100644
index 000000000..46de7552a
--- /dev/null
+++ b/doc/source/reference/constants.rst
@@ -0,0 +1,5 @@
+*********
+Constants
+*********
+
+.. automodule:: numpy.doc.constants
diff --git a/doc/source/reference/index.rst b/doc/source/reference/index.rst
index 4f246096d..2140c57f7 100644
--- a/doc/source/reference/index.rst
+++ b/doc/source/reference/index.rst
@@ -19,6 +19,7 @@ For learning how to use NumPy, see also :ref:`user`.
    :maxdepth: 2
 
    arrays
+   constants
    ufuncs
    routines
    distutils
@@ -34,10 +35,4 @@ Large parts of this manual originate from Travis E. Oliphant's book
 `Guide to NumPy <https://archive.org/details/NumPyBook>`__ (which generously
 entered Public Domain in August 2008). The reference documentation for many of
 the functions are written by numerous contributors and developers of
-NumPy, both prior to and during the
-`NumPy Documentation Marathon
-<http://scipy.org/Developer_Zone/DocMarathon2008>`__.
-
-Please help to improve NumPy's documentation! Instructions on how to
-join the ongoing documentation marathon can be found
-`on the scipy.org website <http://scipy.org/Developer_Zone/DocMarathon2008>`__
+NumPy.
\ No newline at end of file
diff --git a/doc/source/reference/routines.io.rst b/doc/source/reference/routines.io.rst
index 573498792..55489951f 100644
--- a/doc/source/reference/routines.io.rst
+++ b/doc/source/reference/routines.io.rst
@@ -14,7 +14,7 @@ NumPy binary files (NPY, NPZ)
    savez_compressed
 
 The format of these binary file types is documented in
-http://numpy.github.io/neps/npy-format.html
+:py:mod:`numpy.lib.format`
 
 Text files
 ----------
@@ -78,3 +78,11 @@ Data sources
    :toctree: generated/
 
    DataSource
+
+Binary Format Description
+-------------------------
+.. autosummary::
+   :template: autosummary/minimal_module.rst
+   :toctree: generated/
+
+    lib.format
diff --git a/doc/source/reference/routines.statistics.rst b/doc/source/reference/routines.statistics.rst
index d359541aa..e287fe9c8 100644
--- a/doc/source/reference/routines.statistics.rst
+++ b/doc/source/reference/routines.statistics.rst
@@ -17,6 +17,8 @@ Order statistics
    ptp
    percentile
    nanpercentile
+   quantile
+   nanquantile
 
 Averages and variances
 ----------------------
diff --git a/doc/source/reference/swig.testing.rst b/doc/source/reference/swig.testing.rst
index 13642a52e..594df952e 100644
--- a/doc/source/reference/swig.testing.rst
+++ b/doc/source/reference/swig.testing.rst
@@ -22,7 +22,7 @@ typemaps are working as expected.
 Testing Organization
 --------------------
 
-There are three indepedent testing frameworks supported, for one-,
+There are three independent testing frameworks supported, for one-,
 two-, and three-dimensional arrays respectively.  For one-dimensional
 arrays, there are two C++ files, a header and a source, named::
 
diff --git a/doc/source/reference/ufuncs.rst b/doc/source/reference/ufuncs.rst
index 3711f660f..68d85ccf5 100644
--- a/doc/source/reference/ufuncs.rst
+++ b/doc/source/reference/ufuncs.rst
@@ -327,6 +327,13 @@ advanced usage and will not typically be used.
     multiple outputs is deprecated, and will raise a warning in numpy 1.10,
     and an error in a future release.
 
+    If 'out' is None (the default), a uninitialized return array is created.
+    The output array is then filled with the results of the ufunc in the places
+    that the broadcast 'where' is True. If 'where' is the scalar True (the
+    default), then this corresponds to the entire output being filled.
+    Note that outputs not explicitly filled are left with their
+    uninitialized values.
+
 *where*
 
     .. versionadded:: 1.7
@@ -336,6 +343,35 @@ advanced usage and will not typically be used.
     of False indicate to leave the value in the output alone. This argument
     cannot be used for generalized ufuncs as those take non-scalar input.
 
+    Note that if an uninitialized return array is created, values of False
+    will leave those values **uninitialized**.
+
+*axes*
+
+    .. versionadded:: 1.15
+
+    A list of tuples with indices of axes a generalized ufunc should operate
+    on. For instance, for a signature of ``(i,j),(j,k)->(i,k)`` appropriate
+    for matrix multiplication, the base elements are two-dimensional matrices
+    and these are taken to be stored in the two last axes of each argument.
+    The corresponding axes keyword would be ``[(-2, -1), (-2, -1), (-2, -1)]``.
+    For simplicity, for generalized ufuncs that operate on 1-dimensional arrays
+    (vectors), a single integer is accepted instead of a single-element tuple,
+    and for generalized ufuncs for which all outputs are scalars, the output
+    tuples can be omitted.
+
+*keepdims*
+
+    .. versionadded:: 1.15
+
+    If this is set to `True`, axes which are reduced over will be left in the
+    result as a dimension with size one, so that the result will broadcast
+    correctly against the inputs. This option can only be used for generalized
+    ufuncs that operate on inputs that all have the same number of core
+    dimensions and with outputs that have no core dimensions , i.e., with
+    signatures like ``(i),(i)->()`` or ``(m,m)->()``. If used, the location of
+    the dimensions in the output can be controlled with ``axes``.
+
 *casting*
 
     .. versionadded:: 1.6