Accelerated simplex-laminated and vectorized function

Question

I use sympy to create some functions for numerical calculations. So I bind the expression to vectorize it, to use it with numpy arrays. Here is an example:

import numpy as np import sympy as sp def numpy_function(): x, y, z = np.mgrid[0:1:40*1j, 0:1:40*1j, 0:1:40*1j] T = (1 - np.cos(2*np.pi*x))*(1 - np.cos(2*np.pi*y))*np.sin(np.pi*z)*0.1 return T def sympy_function(): x, y, z = sp.Symbol("x"), sp.Symbol("y"), sp.Symbol("z") T = (1 - sp.cos(2*sp.pi*x))*(1 - sp.cos(2*sp.pi*y))*sp.sin(sp.pi*z)*0.1 lambda_function = np.vectorize(sp.lambdify((x, y, z), T, "numpy")) x, y, z = np.mgrid[0:1:40*1j, 0:1:40*1j, 0:1:40*1j] T = lambda_function(x,y,z) return T 

The problem between the sympy version and the pure numpy version is speed, i.e.

 In [3]: timeit test.numpy_function() 100 loops, best of 3: 11.9 ms per loop 

against.

 In [4]: timeit test.sympy_function() 1 loops, best of 3: 634 ms per loop 

So, is there a way to get closer to the speed of the numpy version? I think np.vectorize is rather slow, but somehow part of my code does not work without it. Thank you for any suggestions.

EDIT : Therefore, I found the reason why the vectorization function is needed, i.e.

 In [35]: y = np.arange(10) In [36]: f = sp.lambdify(x,sin(x),"numpy") In [37]: f(y) Out[37]: array([ 0. , 0.84147098, 0.90929743, 0.14112001, -0.7568025 , -0.95892427, -0.2794155 , 0.6569866 , 0.98935825, 0.41211849]) 

this seems to work fine:

 In [38]: y = np.arange(10) In [39]: f = sp.lambdify(x,1,"numpy") In [40]: f(y) Out[40]: 1 

So, for a simple expression like 1 , this function does not return an array. Is there a way to fix this and isn't this some kind of mistake or at least an inconsistent design?

Answer 1

lambdify returns a single value for constants because no numpy functions are involved. This is due to how lambdify works (see fooobar.com/questions/1202407 / ... ).

But this, as a rule, is not a problem, because a constant will be automatically translated into the correct form in any operation in which you use it with an array. On the other hand, if you were explicitly working with an array of the same constant, that would be much less efficient, because you could calculate the same operations multiple times.

Answer 2

Using np.vectorize() in this case is like a loop over the first size x , y and z , and therefore it becomes slower. You do not need np.vectorize() IF to tell lambdify() use the NumPy functions, what exactly are you doing. Then using:

 def sympy_function(): x, y, z = sp.Symbol("x"), sp.Symbol("y"), sp.Symbol("z") T = (1 - sp.cos(2*sp.pi*x))*(1 - sp.cos(2*sp.pi*y))*sp.sin(sp.pi*z)*0.1 lambda_function = sp.lambdify((x, y, z), T, "numpy") x, y, z = np.mgrid[0:1:40*1j, 0:1:40*1j, 0:1:40*1j] T = lambda_function(x,y,z) return T 

makes performance comparable:

 In [26]: np.allclose(numpy_function(), sympy_function()) Out[26]: True In [27]: timeit numpy_function() 100 loops, best of 3: 4.08 ms per loop In [28]: timeit sympy_function() 100 loops, best of 3: 5.52 ms per loop