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. :-)
DG --- On Mon, 6/8/09, Juanjo Gomez Navarro <juanjo.gomeznava...@gmail.com> wrote: > From: Juanjo Gomez Navarro <juanjo.gomeznava...@gmail.com> > Subject: [Numpy-discussion] How to remove fortran-like loops with numpy? > To: numpy-discussion@scipy.org > Date: Monday, June 8, 2009, 2:52 PM > > Hi all, > > I'm new in numpy. Actually, I'm new in Python. In > order to learn a bit, I want to create a program to plot the > Mandelbrot set. This program is quite simple, and I have > already programmed it. The problem is that I come from > fortran, so I use to think in "for" loops. I know > that it is not the best way to use Python and in fact the > performance of the program is more than poor. > > > > Here is the program: > > #!/usr/bin/python > > > > import numpy as np > > import matplotlib.pyplot as plt > > > > # Some parameters > Xmin=-1.5 > > Xmax=0.5 > > Ymin=-1 > > Ymax=1 > > > Ds = 0.01 > > > # Initialization of varibles > X = np.arange(Xmin,Xmax,Ds) > > > Y = np.arange(Ymax,Ymin,-Ds) > > > > N = np.zeros((X.shape[0],Y.shape[0]),'f') > > > ############## Here are inefficient the calculations > ################ > for i in range(X.shape[0]): > for j in range(Y.shape[0]): > z= complex(0.0, 0.0) > c = complex(X[i], Y[j]) > > while N[i, j] < 30 and abs(z) < 2: > > N[i, j] += 1 > z = z**2 + c > if N[i, j] == 29: > N[i, j]=0 > #################################################### > > # And now, just for ploting... > N = N.transpose() > > fig = plt.figure() > > plt.imshow(N,cmap=plt.cm.Blues) > > > plt.title('Mandelbrot set') > > plt.xticks([]); plt.yticks([]) > > > plt.show() > > fig.savefig('test.png') > > > > As you can see, it is very simple, but it takes several > seconds running just to create a 200x200 plot. Fortran takes > same time to create a 2000x2000 plot, around 100 times > faster... So the question is, do you know how to programme > this in a Python-like fashion in order to improve seriously > the performance? > > > > Thanks in advance > -- > Juan José Gómez Navarro > > Edificio CIOyN, Campus de Espinardo, 30100 > Departamento de Física > Universidad de Murcia > Tfno. (+34) 968 398552 > > Email: juanjo.gomeznava...@gmail.com > > Web: http://ciclon.inf.um.es/Inicio.html > > > -----Inline Attachment Follows----- > > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion@scipy.org > http://mail.scipy.org/mailman/listinfo/numpy-discussion > _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion