[Numpy-discussion] How to remove fortran-like loops with numpy?

Anne Archibald peridot.faceted@gmail....
Mon Jun 8 17:38:46 CDT 2009


2009/6/8 Robert Kern <robert.kern@gmail.com>:
> On Mon, Jun 8, 2009 at 17:04, David Goldsmith<d_l_goldsmith@yahoo.com> wrote:
>>
>> I look forward to an instructive reply: the "Pythonic" way to do it would be to take advantage of the facts that Numpy is "pre-vectorized" and uses broadcasting, but so far I haven't been able to figure out (though I haven't yet really buckled down and tried real hard) how to broadcast a conditionally-terminated iteration where the number of iterations will vary among the array elements.  Hopefully someone else already has. :-)
>
> You can't, really. What you can do is just keep iterating with the
> whole data set and ignore the parts that have already converged. Here
> is an example:

Well, yes and no. This is only worth doing if the number of problem
points that require many iterations is small - not the case here
without some sort of periodicity detection - but you can keep an array
of not-yet-converged points, which you iterate. When some converge,
you store them in a results array (with fancy indexing) and remove
them from your still-converging array.

It's also worth remembering that the overhead of for loops is large
but not enormous, so you can often remove only the inner for loop, in
this case perhaps iterating over the image a line at a time.

Anne

>
> z = np.zeros((201,201), dtype=complex)
> Y, X = np.mgrid[1:-1:-201j, -1.5:0.5:201j]
> c = np.empty_like(z)
> c.real = X
> c.imag = Y
> N = np.zeros(z.shape, dtype=int)
>
> while ((N<30) | (abs(z)<2)).all():
>    N += abs(z) < 2
>    z = z ** 2 + c
>
> N[N>=30] = 0
>
> --
> Robert Kern
>
> "I have come to believe that the whole world is an enigma, a harmless
> enigma that is made terrible by our own mad attempt to interpret it as
> though it had an underlying truth."
>  -- Umberto Eco
> _______________________________________________
> Numpy-discussion mailing list
> Numpy-discussion@scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>


More information about the Numpy-discussion mailing list