Repository: pydata/sparse

Issue:
In-place operations
I was wondering about operators such as `operator.iadd`, etc. There are a few ways we can go about this:

1. Support them by mutating the object in-place and invalidating the cache, i.e. performing the operation and then making `self` a copy of the returned object.
2. Support them only when the sparsity structure is the same, and modify `data` in-place.
3. Don't support them at all.

If we want to maintain compatibility with Numpy code (I hope to make `COO` a mostly drop-in replacement for `ndarray` with a few exceptions at some point), I would go with 1, with a warning in the docs that in-place isn't really "in-place".

If we want to do our own thing... Then we have options 2 and 3.

Discussion:
codecov-io: # [Codecov](https://codecov.io/gh/pydata/sparse/pull/146?src=pr&el=h1) Report
> Merging [#146](https://codecov.io/gh/pydata/sparse/pull/146?src=pr&el=desc) into [master](https://codecov.io/gh/pydata/sparse/commit/8f2a9aebe595762eace6bc48531119462f979e21?src=pr&el=desc) will **increase** coverage by `0.65%`.
> The diff coverage is `100%`.

[![Impacted file tree graph](https://codecov.io/gh/pydata/sparse/pull/146/graphs/tree.svg?token=H212u0Uxxw&height=150&src=pr&width=650)](https://codecov.io/gh/pydata/sparse/pull/146?src=pr&el=tree)

```diff
@@            Coverage Diff             @@
##           master     #146      +/-   ##
==========================================
+ Coverage   95.81%   96.46%   +0.65%     
==========================================
  Files          10       10              
  Lines        1195     1189       -6     
==========================================
+ Hits         1145     1147       +2     
+ Misses         50       42       -8
```


| [Impacted Files](https://codecov.io/gh/pydata/sparse/pull/146?src=pr&el=tree) | Coverage Δ | |
|---|---|---|
| [sparse/coo/umath.py](https://codecov.io/gh/pydata/sparse/pull/146/diff?src=pr&el=tree#diff-c3BhcnNlL2Nvby91bWF0aC5weQ==) | `97.11% <ø> (+0.44%)` | :arrow_up: |
| [sparse/coo/core.py](https://codecov.io/gh/pydata/sparse/pull/146/diff?src=pr&el=tree#diff-c3BhcnNlL2Nvby9jb3JlLnB5) | `94.65% <100%> (+0.98%)` | :arrow_up: |
| [sparse/coo/indexing.py](https://codecov.io/gh/pydata/sparse/pull/146/diff?src=pr&el=tree#diff-c3BhcnNlL2Nvby9pbmRleGluZy5weQ==) | `100% <0%> (ø)` | :arrow_up: |
| [sparse/coo/common.py](https://codecov.io/gh/pydata/sparse/pull/146/diff?src=pr&el=tree#diff-c3BhcnNlL2Nvby9jb21tb24ucHk=) | `97.04% <0%> (+0.92%)` | :arrow_up: |
| [sparse/dok.py](https://codecov.io/gh/pydata/sparse/pull/146/diff?src=pr&el=tree#diff-c3BhcnNlL2Rvay5weQ==) | `95.32% <0%> (+1.74%)` | :arrow_up: |

------

[Continue to review full report at Codecov](https://codecov.io/gh/pydata/sparse/pull/146?src=pr&el=continue).
> **Legend** - [Click here to learn more](https://docs.codecov.io/docs/codecov-delta)
> `Δ = absolute <relative> (impact)`, `ø = not affected`, `? = missing data`
> Powered by [Codecov](https://codecov.io/gh/pydata/sparse/pull/146?src=pr&el=footer). Last update [8f2a9ae...09edda0](https://codecov.io/gh/pydata/sparse/pull/146?src=pr&el=lastupdated). Read the [comment docs](https://docs.codecov.io/docs/pull-request-comments).

mrocklin: This looks pretty slick :)

Might want to test something like the following as well:

```python
np.sin(x, out=x)
```
hameerabbasi: Docs updated, warning added to docs that "in-place ops aren't really in-place".

Supported `out` kwarg for all ops now, including `round` and `COO.sum`, etc.
hameerabbasi: cc @mrocklin I think this is ready for a final review. :-)
hameerabbasi: Removed some unnecessary code for `astype` -- It doesn't have an `out` argument. It was what was causing coverage to fail.

Patch:
diff --git a/.codecov.yml b/.codecov.yml
index a385d75..6cfd0cf 100644
--- a/.codecov.yml
+++ b/.codecov.yml
@@ -12,7 +12,6 @@ flags:
   project:
     enabled: yes
     target: 95%
-    threshold: 1%
     if_no_uploads: error
     if_not_found: success
     if_ci_failed: error
@@ -21,7 +20,6 @@ flags:
     default:
       enabled: yes
       target: 95%
-      threshold: 1%
       if_no_uploads: error
       if_not_found: success
       if_ci_failed: error
diff --git a/.coveragerc b/.coveragerc
index d838d37..6371c70 100644
--- a/.coveragerc
+++ b/.coveragerc
@@ -5,3 +5,9 @@ source=
 omit=
     sparse/_version.py
     sparse/tests/*
+
+[report]
+exclude_lines =
+    pragma: no cover
+    return NotImplemented
+    raise NotImplementedError
diff --git a/docs/operations.rst b/docs/operations.rst
index a6c1a10..fa83e92 100644
--- a/docs/operations.rst
+++ b/docs/operations.rst
@@ -20,13 +20,12 @@ results for both Numpy arrays, COO arrays, or a mix of the two:
 
    np.log(X.dot(beta.T) + 1)
 
-However some operations are not supported, like inplace operations,
-operations that implicitly cause dense structures,
-or numpy functions that are not yet implemented for sparse arrays
+However some operations are not supported, like operations that
+implicitly cause dense structures, or numpy functions that are not
+yet implemented for sparse arrays.
 
 .. code-block:: python
 
-   x += y     # inplace operations not supported
    x + 1      # operations that produce dense results not supported
    np.svd(x)  # sparse svd not implemented
 
@@ -34,7 +33,7 @@ or numpy functions that are not yet implemented for sparse arrays
 This page describes those valid operations, and their limitations.
 
 :obj:`elemwise`
-~~~~~~~~~~~~~~~~~~~
+~~~~~~~~~~~~~~~
 This function allows you to apply any arbitrary broadcasting function to any number of arguments
 where the arguments can be :obj:`SparseArray` objects or :obj:`scipy.sparse.spmatrix` objects.
 For example, the following will add two arrays:
@@ -155,7 +154,9 @@ be a :obj:`scipy.sparse.spmatrix` The following operators are supported:
    * :obj:`operator.lshift` (:code:`x << y`)
    * :obj:`operator.rshift` (:code:`x >> y`)
 
-.. note:: In-place operators are not supported at this time.
+.. warning::
+   While in-place operations are supported for compatibility with Numpy,
+   they are not truly in-place, and will effectively calculate the result separately.
 
 .. _operations-elemwise:
 
diff --git a/sparse/coo/core.py b/sparse/coo/core.py
index 331e21f..33c427e 100644
--- a/sparse/coo/core.py
+++ b/sparse/coo/core.py
@@ -676,8 +676,7 @@ class COO(SparseArray, NDArrayOperatorsMixin):
         >>> s.sum()
         25
         """
-        assert out is None
-        return self.reduce(np.add, axis=axis, keepdims=keepdims, dtype=dtype)
+        return np.add.reduce(self, out=out, axis=axis, keepdims=keepdims, dtype=dtype)
 
     def max(self, axis=None, keepdims=False, out=None):
         """
@@ -738,8 +737,7 @@ class COO(SparseArray, NDArrayOperatorsMixin):
         >>> s.max()
         8
         """
-        assert out is None
-        return self.reduce(np.maximum, axis=axis, keepdims=keepdims)
+        return np.maximum.reduce(self, out=out, axis=axis, keepdims=keepdims)
 
     def min(self, axis=None, keepdims=False, out=None):
         """
@@ -800,8 +798,7 @@ class COO(SparseArray, NDArrayOperatorsMixin):
         >>> s.min()
         0
         """
-        assert out is None
-        return self.reduce(np.minimum, axis=axis, keepdims=keepdims)
+        return np.minimum.reduce(self, out=out, axis=axis, keepdims=keepdims)
 
     def prod(self, axis=None, keepdims=False, dtype=None, out=None):
         """
@@ -867,8 +864,7 @@ class COO(SparseArray, NDArrayOperatorsMixin):
         >>> s.prod()
         0
         """
-        assert out is None
-        return self.reduce(np.multiply, axis=axis, keepdims=keepdims, dtype=dtype)
+        return np.multiply.reduce(self, out=out, axis=axis, keepdims=keepdims, dtype=dtype)
 
     def transpose(self, axes=None):
         """
@@ -1039,16 +1035,27 @@ class COO(SparseArray, NDArrayOperatorsMixin):
 
     def __array_ufunc__(self, ufunc, method, *inputs, **kwargs):
         out = kwargs.pop('out', None)
-        if out is not None:
+        if out is not None and not all(isinstance(x, COO) for x in out):
             return NotImplemented
 
         if method == '__call__':
-            return elemwise(ufunc, *inputs, **kwargs)
+            result = elemwise(ufunc, *inputs, **kwargs)
         elif method == 'reduce':
-            return COO._reduce(ufunc, *inputs, **kwargs)
+            result = COO._reduce(ufunc, *inputs, **kwargs)
         else:
             return NotImplemented
 
+        if out is not None:
+            (out,) = out
+            if out.shape != result.shape:
+                raise ValueError('non-broadcastable output operand with shape %s'
+                                 'doesn\'t match the broadcast shape %s' % (out.shape, result.shape))
+
+            out._make_shallow_copy_of(result)
+            return out
+
+        return result
+
     def __array__(self, dtype=None, **kwargs):
         x = self.todense()
         if dtype and x.dtype != dtype:
@@ -1366,10 +1373,11 @@ class COO(SparseArray, NDArrayOperatorsMixin):
         The :code:`out` parameter is provided just for compatibility with Numpy and isn't
         actually supported.
         """
-        assert out is None
-        return elemwise(np.round, self, decimals=decimals)
+        if out is not None and not isinstance(out, tuple):
+            out = (out,)
+        return self.__array_ufunc__(np.round, '__call__', self, decimals=decimals, out=out)
 
-    def astype(self, dtype, out=None):
+    def astype(self, dtype):
         """
         Copy of the array, cast to a specified type.
 
@@ -1385,8 +1393,7 @@ class COO(SparseArray, NDArrayOperatorsMixin):
         The :code:`out` parameter is provided just for compatibility with Numpy and isn't
         actually supported.
         """
-        assert out is None
-        return elemwise(np.ndarray.astype, self, dtype=dtype)
+        return self.__array_ufunc__(np.ndarray.astype, '__call__', self, dtype=dtype)
 
     def maybe_densify(self, max_size=1000, min_density=0.25):
         """
diff --git a/sparse/coo/umath.py b/sparse/coo/umath.py
index 45de789..d5ebf58 100644
--- a/sparse/coo/umath.py
+++ b/sparse/coo/umath.py
@@ -64,8 +64,7 @@ def elemwise(func, *args, **kwargs):
         elif isinstance(arg, SparseArray) and not isinstance(arg, COO):
             args[i] = COO(arg)
         elif not isinstance(arg, COO):
-            raise ValueError("Performing this operation would produce "
-                             "a dense result: %s" % str(func))
+            return NotImplemented
 
     # Filter out scalars as they are 'baked' into the function.
     func = PositinalArgumentPartial(func, pos, posargs)


Test patch:
diff --git a/sparse/tests/test_coo.py b/sparse/tests/test_coo.py
index 74bcce0..2163892 100644
--- a/sparse/tests/test_coo.py
+++ b/sparse/tests/test_coo.py
@@ -264,6 +264,22 @@ def test_elemwise(func):
     assert_eq(func(x), fs)
 
 
+@pytest.mark.parametrize('func', [np.expm1, np.log1p, np.sin, np.tan,
+                                  np.sinh, np.tanh, np.floor, np.ceil,
+                                  np.sqrt, np.conj,
+                                  np.round, np.rint, np.conjugate,
+                                  np.conj, lambda x, out: x.round(decimals=2, out=out)])
+def test_elemwise_inplace(func):
+    s = sparse.random((2, 3, 4), density=0.5)
+    x = s.todense()
+
+    func(s, out=s)
+    func(x, out=x)
+    assert isinstance(s, COO)
+
+    assert_eq(x, s)
+
+
 @pytest.mark.parametrize('func', [
     operator.mul, operator.add, operator.sub, operator.gt,
     operator.lt, operator.ne
@@ -279,6 +295,23 @@ def test_elemwise_binary(func, shape):
     assert_eq(func(xs, ys), func(x, y))
 
 
+@pytest.mark.parametrize('func', [
+    operator.imul, operator.iadd, operator.isub
+])
+@pytest.mark.parametrize('shape', [(2,), (2, 3), (2, 3, 4), (2, 3, 4, 5)])
+def test_elemwise_binary_inplace(func, shape):
+    xs = sparse.random(shape, density=0.5)
+    ys = sparse.random(shape, density=0.5)
+
+    x = xs.todense()
+    y = ys.todense()
+
+    xs = func(xs, ys)
+    x = func(x, y)
+
+    assert_eq(xs, x)
+
+
 @pytest.mark.parametrize('func', [
     lambda x, y, z: x + y + z,
     lambda x, y, z: x * y * z,
@@ -497,7 +530,7 @@ def test_ndarray_densification_fails():
     xs = sparse.random((3, 4), density=0.5)
     y = np.random.rand(3, 4)
 
-    with pytest.raises(ValueError):
+    with pytest.raises(TypeError):
         xs + y
 
 
@@ -624,6 +657,30 @@ def test_bitwise_binary(func, shape):
     assert_eq(func(xs, ys), func(x, y))
 
 
+@pytest.mark.parametrize('func', [
+    operator.iand, operator.ior, operator.ixor
+])
+@pytest.mark.parametrize('shape', [
+    (2,),
+    (2, 3),
+    (2, 3, 4),
+    (2, 3, 4, 5)
+])
+def test_bitwise_binary_inplace(func, shape):
+    # Small arrays need high density to have nnz entries
+    # Casting floats to int will result in all zeros, hence the * 100
+    xs = (sparse.random(shape, density=0.5) * 100).astype(np.int_)
+    ys = (sparse.random(shape, density=0.5) * 100).astype(np.int_)
+
+    x = xs.todense()
+    y = ys.todense()
+
+    xs = func(xs, ys)
+    x = func(x, y)
+
+    assert_eq(xs, x)
+
+
 @pytest.mark.parametrize('func', [
     operator.lshift, operator.rshift
 ])
@@ -649,7 +706,7 @@ def test_bitshift_binary(func, shape):
 
 
 @pytest.mark.parametrize('func', [
-    operator.and_
+    operator.ilshift, operator.irshift
 ])
 @pytest.mark.parametrize('shape', [
     (2,),
@@ -657,13 +714,37 @@ def test_bitshift_binary(func, shape):
     (2, 3, 4),
     (2, 3, 4, 5)
 ])
-def test_bitwise_scalar(func, shape):
+def test_bitshift_binary_inplace(func, shape):
     # Small arrays need high density to have nnz entries
     # Casting floats to int will result in all zeros, hence the * 100
     xs = (sparse.random(shape, density=0.5) * 100).astype(np.int_)
 
     # Can't merge into test_bitwise_binary because left/right shifting
     # with something >= 64 isn't defined.
+    ys = (sparse.random(shape, density=0.5) * 64).astype(np.int_)
+
+    x = xs.todense()
+    y = ys.todense()
+
+    xs = func(xs, ys)
+    x = func(x, y)
+
+    assert_eq(xs, x)
+
+
+@pytest.mark.parametrize('func', [
+    operator.and_
+])
+@pytest.mark.parametrize('shape', [
+    (2,),
+    (2, 3),
+    (2, 3, 4),
+    (2, 3, 4, 5)
+])
+def test_bitwise_scalar(func, shape):
+    # Small arrays need high density to have nnz entries
+    # Casting floats to int will result in all zeros, hence the * 100
+    xs = (sparse.random(shape, density=0.5) * 100).astype(np.int_)
     y = np.random.randint(100)
 
     x = xs.todense()
@@ -1376,6 +1457,17 @@ def test_two_arg_where():
         sparse.where(cs, xs)
 
 
+@pytest.mark.parametrize('func', [
+    operator.imul, operator.iadd, operator.isub
+])
+def test_inplace_invalid_shape(func):
+    xs = sparse.random((3, 4), density=0.5)
+    ys = sparse.random((2, 3, 4), density=0.5)
+
+    with pytest.raises(ValueError):
+        func(xs, ys)
+
+
 def test_nonzero():
     s = sparse.random((2, 3, 4), density=0.5)
     x = s.todense()