0

I have a question regarding numpy.

Let's assume I have a function

def calcSomething(A,t):

    return t*A

I want to pass A a constant numpy-array, e.g. A=np.array([1.0,2.0]) and t as equidistant points, e.g. t = np.linspace(0.0,1.0,10). Now, when I execute

x = calcSomething(A,t)

I'd like to get x[0] = [0,0], x[1] = [.1,.2],...,x[9] = [1.0,2.0]. Is there an easy way to achieve that, or do I have to implement it using for-loops and appending?

Edit: 'calcSomething(): is just a placeholder for a more complex function.

E.g. y=np.sin(t) yields y[i]=np.sin(t[i]). That's what I want to achieve for a generic function.

Thanks in advance!

J.Doe
  • 69
  • 2
  • 7
  • `np.vectorize` takes a 'scalar' function, and feeds it values from arrays. But it does not promise speed. It's just a convenience for when you can't write the function to take full advantage of compiled code. This topic comes up quite often, with people trying to avoid nested loops or trying to implement numpy vectorization in a magical, easy, way. – hpaulj Apr 11 '18 at 16:36

3 Answers3

1

Where possible, using numpy functions and operators that operate on the whole arrays, and do the necessary broadcasting for you:

In [24]: A = np.array([1.0,2.0]); t = np.linspace(0.0, 1.0, 10)
In [25]: x = t[:,None] * A[None,:]
In [26]: x.shape
Out[26]: (10, 2)
In [27]: x[:3,:]
Out[27]: 
array([[0.        , 0.        ],
       [0.11111111, 0.22222222],
       [0.22222222, 0.44444444]])
In [28]: np.sin(t)
Out[28]: 
   array([0.        , 0.11088263, 0.22039774, 0.3271947 , 0.42995636,
          0.52741539, 0.6183698 , 0.70169788, 0.77637192, 0.84147098])

If you have a function that only works with scalar inputs, you can use np.vectorize to feed it those values - from broadcasted arrays:

In [30]: def foo(A,t):
    ...:     return t*A
    ...: 
In [31]: f = np.vectorize(foo, otypes=[float])
In [32]: f(A[None,:], t[:,None])
Out[32]: 
array([[0.        , 0.        ],
       [0.11111111, 0.22222222],
       [0.22222222, 0.44444444],
       ....
       [1.        , 2.        ]])

vectorize is a convenience function; it does not promise speed. Compare its time with In[25]:

In [33]: timeit f(A[None,:], t[:,None])
41.3 µs ± 48.5 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)
In [34]: timeit x = t[:,None] * A[None,:]
5.41 µs ± 8.2 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

8x slower.

Or with the scalar math.sin function:

In [35]: import math
In [36]: g = np.vectorize(math.sin)
...
In [39]: timeit g(t)
39 µs ± 72.4 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)
In [40]: timeit np.sin(t)
1.4 µs ± 2.72 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)

Called iteratively math.sin is faster than np.sin, but np.sin is much faster when given a whole array.

The basic difference is that vectorize, like explicit loops is iterating in interpreted Python, and calling your function once for each output element. np.sin and array * is iterating, but in compiled code.

There are various ways of iterating in Python and numpy, but rarely do they give you more than a 2x speedup.

There are tools for moving calculations to compiled code, such as cython and numba. Search on those tags to get ideas on how to use them.

hpaulj
  • 221,503
  • 14
  • 230
  • 353
  • Thanks very much for your detailed answer! Problem is that I have to operate in `Abaqus-Python` on Windows which makes it very hard to use additional packages like `cython`. Unless you also know an answer for that ;) – J.Doe Apr 12 '18 at 06:56
0

You can use matrix multiplication:

>>> A = np.array([1.0,2.0])
>>> t = np.linspace(0.0,1.0,10)
>>> np.matmul(t.reshape(10, 1), A.reshape(1, 2))
array([[0.        , 0.        ],
       [0.11111111, 0.22222222],
       [0.22222222, 0.44444444],
       [0.33333333, 0.66666667],
       [0.44444444, 0.88888889],
       [0.55555556, 1.11111111],
       [0.66666667, 1.33333333],
       [0.77777778, 1.55555556],
       [0.88888889, 1.77777778],
       [1.        , 2.        ]])
ducminh
  • 1,312
  • 1
  • 8
  • 17
0

Just use None to specify a new axis:

A[None, :] * t[:, None] 
ben thorne
  • 61
  • 5
  • Thanks very much. However, the 'calcSomething' function is just a placeholder for a more complex function. – J.Doe Apr 11 '18 at 15:19