Numpy passing input array as `out` argument to ufunc

Question

Is it generally safe to provide the input array as the optional out argument to a ufunc in numpy, provided the type is correct? For example, I have verified that the following works:

>>> import numpy as np
>>> arr = np.array([1.2, 3.4, 4.5])
>>> np.floor(arr, arr)
array([ 1.,  3.,  4.])

The array type must be either compatible or identical with the output (which is a float for numpy.floor()), or this happens:

>>> arr2 = np.array([1, 3, 4], dtype = np.uint8)
>>> np.floor(arr2, arr2)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
TypeError: ufunc 'floor' output (typecode 'e') could not be coerced to provided output parameter (typecode 'B') according to the casting rule ''same_kind''

So given that an array of proper type, is it generally safe to apply ufuncs in-place? Or is floor() an exceptional case? The documentation does not make it clear, and neither do the following two threads that have tangential bearing on the question:

1. Numpy modify array in place?
2. Numpy Ceil and Floor "out" Argument

EDIT:

As a first order guess, I would assume it is often, but not always safe, based on the tutorial at http://docs.scipy.org/doc/numpy/user/c-info.ufunc-tutorial.html. There does not appear to be any restriction on using the output array as a temporary holder for intermediate results during the computation. While something like floor() and ciel() may not require temporary storage, more complex functions might. That being said, the entire existing library may be written with that in mind.

Riyal

It's not technically a ufunc, but using the out parameter in np.dot in this way with 2D arrays can produce incorrect results. – Pskov

That is almost the counterexample I was looking for, but not quite :) – Riyal

The ufunc docs mention using add(G, C, G) as an optimization of G = G + C, in the Tip under "Math operations". I'd say it's safe. (On the other hand, calling ufuncs with input and output overlapping but not identical will cause problems.) – Fragmental

True, but that would go under the category of overlapping but not identical inputs. I was just wondering if all the operations are safe to do in place given that all other conditions are met. – Riyal

See here regarding temporary storage. – Carthusian

Interesting. So that strongly implies that the output buffer is never used for intermediate calculations. – Riyal

Answer 1

This is an old question, but there is an updated answer:

Yes, it is safe. In the Numpy documentation, we see that as of v1.13:

Operations where ufunc input and output operands have memory overlap are defined to be the same as for equivalent operations where there is no memory overlap. Operations affected make temporary copies as needed to eliminate data dependency. As detecting these cases is computationally expensive, a heuristic is used, which may in rare cases result in needless temporary copies. For operations where the data dependency is simple enough for the heuristic to analyze, temporary copies will not be made even if the arrays overlap, if it can be deduced copies are not necessary. As an example, np.add(a, b, out=a) will not involve copies.

Rojo

Thanks a lot for this answer, i've been looking all over the internet for it. – Gonorrhea

Answer 2

The out parameter of a numpy function is the array where the result is written. The main advantage of using out is avoiding the allocation of new memory where it is not necessary.

Is it safe to use write the output of a function on the same array passed as input? There is no general answer, it depends on what the function is doing.

Two examples

Here are two examples of ufunc-like functions:

In [1]: def plus_one(x, out=None):
   ...:     if out is None:
   ...:         out = np.zeros_like(x)
   ...: 
   ...:     for i in range(x.size):
   ...:         out[i] = x[i] + 1
   ...:     return out
   ...: 

In [2]: x = np.arange(5)

In [3]: x
Out[3]: array([0, 1, 2, 3, 4])

In [4]: y = plus_one(x)

In [5]: y
Out[5]: array([1, 2, 3, 4, 5])

In [6]: z = plus_one(x, x)

In [7]: z
Out[7]: array([1, 2, 3, 4, 5])

Function shift_one:

In [11]: def shift_one(x, out=None):
    ...:     if out is None:
    ...:         out = np.zeros_like(x)
    ...: 
    ...:     n = x.size
    ...:     for i in range(n):
    ...:         out[(i+1) % n] = x[i]
    ...:     return out
    ...: 

In [12]: x = np.arange(5)

In [13]: x
Out[13]: array([0, 1, 2, 3, 4])

In [14]: y = shift_one(x)

In [15]: y
Out[15]: array([4, 0, 1, 2, 3])

In [16]: z = shift_one(x, x)

In [17]: z
Out[17]: array([0, 0, 0, 0, 0])

For the function plus_one there is no problem: the expected result is obtained when the parameters x and out are the same array. But the function shift_one gives a surprising result when the parameters x and out are the same array because the array

Discussion

For function of the form out[i] := some_operation(x[i]), such as plus_one above but also the functions floor, ceil, sin, cos, tan, log, conj, etc, as far as I know it is safe to write the result in the input using parameter out.

It is also safe for functions taking two input parameters of the form ``out[i] := some_operation(x[i], y[i]) such as the numpy function add, multiply, subtract.

For the other functions, it is case-by-case. As illustrated bellow, the matrix multiplication is not safe:

In [18]: a = np.arange(4).reshape((2,2))

In [19]: a
Out[19]: 
array([[0, 1],
       [2, 3]])

In [20]: b = (np.arange(4) % 2).reshape((2,2))

In [21]: b
Out[21]: 
array([[0, 1],
       [0, 1]], dtype=int32)

In [22]: c = np.dot(a, b)

In [23]: c
Out[23]: 
array([[0, 1],
       [0, 5]])

In [24]: d = np.dot(a, b, out=a)

In [25]: d
Out[25]: 
array([[0, 1],
       [0, 3]])

Last remark: if the implementation is multithreaded, the result of an unsafe function may even be non-deterministic because it depends on the order on which the array elements are processed.

Antiphonal

With numpy version 1.21.6 (and probably earlier versions), the behavior of np.dot(a, b, out=a) raises a ValueError. Ostensibly this is to prevent such a logic error – Coronet

@KyleMeador at least, in version 1.23.4 ValueError is not raised. Additionally, the example of nbedou seems to work correctly. Maybe it's currently safe to use np.dot in-place. – Peseta