vectorize is indeterminate

Question

I'm trying to vectorize a simple function in numpy and getting inconsistent behavior. I expect my code to return 0 for values < 0.5 and the unchanged value otherwise. Strangely, different runs of the script from the command line yield varying results: sometimes it works correctly, and sometimes I get all 0's. It doesn't matter which of the three lines I use for the case when d <= T. It does seem to be correlated with whether the first value to be returned is 0. Any ideas? Thanks.

import numpy as np

def my_func(d, T=0.5):
    if d > T:   return d
    #if d <= T:  return 0
    else:  return 0
    #return 0

N = 4
A = np.random.uniform(size=N**2)
A.shape = (N,N)
print A
f = np.vectorize(my_func)
print f(A)

$ python x.py
[[ 0.86913815  0.96833127  0.54539153  0.46184594]
 [ 0.46550903  0.24645558  0.26988519  0.0959257 ]
 [ 0.73356391  0.69363161  0.57222389  0.98214089]
 [ 0.15789303  0.06803493  0.01601389  0.04735725]]
[[ 0.86913815  0.96833127  0.54539153  0.        ]
 [ 0.          0.          0.          0.        ]
 [ 0.73356391  0.69363161  0.57222389  0.98214089]
 [ 0.          0.          0.          0.        ]]
$ python x.py
[[ 0.37127366  0.77935622  0.74392301  0.92626644]
 [ 0.61639086  0.32584431  0.12345342  0.17392298]
 [ 0.03679475  0.00536863  0.60936931  0.12761859]
 [ 0.49091897  0.21261635  0.37063752  0.23578082]]
[[0 0 0 0]
 [0 0 0 0]
 [0 0 0 0]
 [0 0 0 0]]

Answer 1

If this really is the problem you want to solve, then there's a much better solution:

A[A<=0.5] = 0.0

The problem with your code, however, is that if the condition passes, you are returning the integer 0, not the float 0.0. From the documentation:

The data type of the output of vectorized is determined by calling the function with the first element of the input. This can be avoided by specifying the otypes argument.

So when the very first entry is <0.5, it tries to create an integer, not float, array. You should change return 0 to

return 0.0

Alternately, if you don't want to touch my_func, you can use

f = np.vectorize(my_func, otypes=[np.float])