[SciPy-User] Parallelizing a for loop

Anne Archibald peridot.faceted@gmail....
Sat Dec 26 13:49:41 CST 2009


2009/12/24 Scott Askey <scotta_2002@yahoo.com>:
> Anne:
>
> A pointer to the tools in numpy to apply array operation to my function would be appreciated.
> I found found the scipy.vectorize but it did not seem applicable to my 42 input

That's a little difficult, since vectorizing operations is the basic
function of numpy. (You're right, the vectorize function isn't really
relevant here; it doesn't give any speed benefit, just a convenience
one.)

For a start, pretty much any calculation that you can write as a
simple formula is automatically vectorized:

In [2]: x = np.linspace(0,1,10)

In [3]: y = x**3 - x**2 + x - 1

In [4]: y
Out[4]:
array([-1.        , -0.89986283, -0.81618656, -0.74074074, -0.66529492,
       -0.58161866, -0.48148148, -0.35665295, -0.19890261,  0.        ])

If your coefficients vary, that works too:

In [5]: a = np.linspace(0,1,10)

In [6]: b = np.linspace(-1,1,10)

In [7]: a*x**3+b
Out[7]:
array([-1.        , -0.77762536, -0.5531169 , -0.32098765, -0.07209267,
        0.20637098,  0.5308642 ,  0.92150587,  1.40207285,  2.        ])

There's lots more, but the point is, since each of these operations
runs in only a handful of python-level instructions, and all the
element-wise operations are done in tight C loops, the calculations
run vastly faster when you write them this way.

> As best I could tell I seemed that the time required to calculate two elements (2 x dof in residual 500 times) was the same as calculating  1 x dof residual 1000 times in a for loop.
>
> Chuck:
>
> It is a time marching structural dynamics 1-d beam problem.  The 100 or so parameters are the initial conditions, geometric and material properties.  The 42 dof are the configuration of the  element at the end of the time step.  Given the array of (42 * elements) initial conditions the residuals for each element may be independently calculated.
>
> Currently I calculate the 42*elements residual by running  "for i in range(element)" and getting  42 residuals in each iteration.
>
> My element residual function is at most at worst cubic polynomials.
> The element wise residual and fprime functions are lambdafiedsympy functions.

I strongly recommend rewriting these so that they are expressed as
(say) a single cubic with arrays as coefficients. You may have to sue
some cleverness to reexpress the sympy functions in a convenient way,
but you will almost certainly find the speed gain worth it.

Anne


More information about the SciPy-User mailing list