Seite 1 von 3
Mandelbrot - Optimierung
Verfasst: Freitag 23. April 2010, 20:32
von Huy
Hallo,
Ich habe vor einigen Wochen die Aufgabe bekommen, ein Programm zu programmieren, welches dazu in der Lage ist, die Mandelbrot-Menge darzustellen. Da ich mit Python nicht eine extra Klasse für komplexe Zahlen oder mühsame Umrechnereien brauchte, entschied ich mich für diese Sprache. Nun bin ich soweit, dass die Mandelbrot-Menge bereits ziemlich sauber dargestellt wird, nur dauert die Rechnerei noch unglaublich lange und das Programm selber stockt auch ziemlich stark.
Code: Alles auswählen
#!/usr/bin/python
# Filename: mandelbrot.py
import matplotlib.pyplot as plt
def mandelbrot():
global c, c0
c = c ** 2 + c0
plt.xlabel(r'$re(z)$')
plt.ylabel(r'$i \cdot im(z)$')
plt.title('The Mandelbrot Set')
start = complex(-3,-3) # starting point
re = start.real # store real and imaginary part in c (part to calculate with) and c0 (fix during calculations)
im = start.imag
c = complex(re,im)
c0 = complex(re,im)
stepx = complex(0.01,0)
stepy = complex(0,0.01)
steps = 601
iterations = 100
for h in range(0,steps):
for k in range(0,steps): # do this for all points in the selected area
if abs(c) >= 2: # check divergence
start += stepx
else:
for i in range(0,iterations):
if abs(c ** 2 + c0) >= 2: # check divergence of value
break
elif i < iterations - 1:
mandelbrot()
elif i == iterations - 1:
plt.plot(round(c0.real,2), round(c0.imag,2), 'k,') # print(c0)
start += stepx
re = start.real
im = start.imag
c = complex(re,im)
c0 = complex(re,im)
start = start + stepy - complex(steps * stepx.real,0)
plt.axis([-3,2,-2,2])
plt.grid(True)
plt.show()
Nun frage ich mich, wie ich das Script bzw. die Rechnereien optimieren kann, damit es nicht mehr so stockt und die Rechenzeit deutlich verringert wird? Ich habe vor Python noch nie irgendetwas programmiert und kenne mich demzufolge ziemlich schlecht mit Programmieren aus... Gibt es vllt. einen "besseren" Plotter für Python als matplot, der schneller arbeitet? Wie würdet ihr eine Zoomfunktion einbauen, die mir erlaubt, unendlich weit hineinzuzoomen? Dabei müsste ich halt einfach immer wieder die Rechnung für alle Punkte im ausgewählten Bereich durchführen, aber wie mache ich sowas?
Danke im Voraus für Hilfe...
MfG
Verfasst: Freitag 23. April 2010, 21:07
von Gabelmensch
Schneller Plotten geht auf jedenfall mit Gnuplot, vermutlich kann Gnuplot von sich selbst aus auch Mandelbrote darstellen.
Re: Optimierung des Scripts
Verfasst: Freitag 23. April 2010, 21:11
von numerix
Vorweg: Ich habe matplotlib noch nie benutzt, habe es auch nicht installiert und kann deinen Code darum nicht testen.
Dennoch ein paar Gedanken:
Deine innere Schleife wird 36 Mio. mal durchlaufen. Falls du die Anzahl dieser Iterationen reduzieren kannst, solltet du das tun.
Auf jeden Fall solltest du bei der inneren Schleife ansetzen:
Der Funktionsaufruf ist - so wie die Funktion gemacht ist - nicht nur hässlich, sondern auch total überflüssig und zudem teurer, als wenn du die eine Zeile direkt in die Schleife setzt.
Die letzte Abfrage mit elif ist teuer und überflüssig. Es genügt doch, diesen Schritt durchzuführen, wenn die Iterationen durchgelaufen sind.
Und: Falls du mit Python 2.x arbeitest, solltest du xrange statt range verwenden, um nicht unnötigerweise immer wieder Listen zu erzeugen.
Verfasst: Freitag 23. April 2010, 21:21
von Huy
@Gabelmensch: Gibt es irgendwo eine richtige Dokumentation von gnuplot.py? Habe es jetzt installiert, aber keine Ahnung, wie man es benutzt...
@numerix: Die Iterationen sollten schon um die 100 liegen, um saubere Resultate liefern zu können. Die Zeile setze ich also direkt in die Schleife. Bezüglich der letzten elif-Abfrage: Wenn ich diesen Schritt durchführe, wenn die Iterationen durchgelaufen sind, wird dann nicht auch im Falle "abs(c ** 2 + c0) >= 2" geplottet? Das möchte ich dann ja nicht... Oder stehe ich gerade aufm Schlauch?
Und ja, ich arbeite mit Python 2.6. Danke für die Tipps bisher.
MfG
Verfasst: Freitag 23. April 2010, 21:23
von gkuhl
Verfasst: Freitag 23. April 2010, 21:47
von numerix
Huy hat geschrieben:Bezüglich der letzten elif-Abfrage: Wenn ich diesen Schritt durchführe, wenn die Iterationen durchgelaufen sind, wird dann nicht auch im Falle "abs(c ** 2 + c0) >= 2" geplottet?
Nicht, wenn du die for-Schleife um einen else-Teil ergänzt. Der wird nur ausgeführt, wenn die Schleife nicht mittels break vorzeitig verlassen wurde.
Re: Optimierung des Scripts
Verfasst: Freitag 23. April 2010, 22:52
von cofi
numerix hat geschrieben:Und: Falls du mit Python 2.x arbeitest, solltest du xrange statt range verwenden, um nicht unnötigerweise immer wieder Listen zu erzeugen.
Code: Alles auswählen
Python 2.6.5 (r265:79063, Mar 18 2010, 23:38:15)
[GCC 4.4.3] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> xrange(10**10)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
OverflowError: long int too large to convert to int
In Python 2.x muss man da leider vorsichtig sein, wenn man mit sehr grossen Zahlen arbeitet. Aber ein einfacher eigener Generator ist da auch kein Problem.
Verfasst: Samstag 24. April 2010, 08:37
von BlackJack
@cofi: Was soll die überflüssige Bemerkung denn jetzt? Der OP verwendet 2.x und `range()` da würde eine Länge von ``10**10`` auch nicht funktionieren -- alleine die Zeiger auf die Objekte würden da 37 Gigabyte benötigen. Das ist ein "Problem" was hier überhaupt nicht existiert.
@Huy: Pack die ganze Berechnung mal in eine Funktion. Lokale Namen werden schneller nachgeschlagen.
Verfasst: Samstag 24. April 2010, 10:25
von Huy
Was meinst du mit der ganzen Berechnung konkret? Alle drei for-Schleifen? Die innerste for-Schleife?
MfG
Verfasst: Samstag 24. April 2010, 11:20
von Gabelmensch
Schau mal hier:
http://rosettacode.org/wiki/Mandelbrot_set
gnuplot.py habe ich nur mal zum testen verwendet, wenn ich mal gnuplot in Py nutze, nehme ich meist subprocess.
Verfasst: Samstag 24. April 2010, 12:19
von BlackJack
@Huy: Ich meine "alles". Auf Modulebene sollten nach Möglichkeit sowieso nur ``import``\e und Definitionen von Funktionen, Klassen, und Konstanten stehen. Das macht Module einfacher wiederverwendbar und testbar.
Verfasst: Samstag 24. April 2010, 16:54
von Huy
Sorry, aber ich glaube, ich stehe gerade aufm Schlauch.
Code: Alles auswählen
#!/usr/bin/python
# Filename: mandelbrot.py
import matplotlib.pyplot as plt
plt.xlabel(r'$re(z)$')
plt.ylabel(r'$i \cdot im(z)$')
plt.title('The Mandelbrot Set')
start = complex(-3,-3) # starting point
re = start.real # store real and imaginary part in c (part to calculate with) and c0 (fix during calculations)
im = start.imag
c = complex(re,im)
c0 = complex(re,im)
stepx = complex(0.01,0)
stepy = complex(0,0.01)
steps = 601
iterations = 100
def mandelbrot():
global start, re, im, c, c0, stepx, stepy, steps, iterations
for h in xrange(0,steps):
for k in xrange(0,steps): # do this for all points in the selected area
if abs(c) >= 2: # check divergence
start += stepx
else:
for i in xrange(0,iterations):
if abs(c ** 2 + c0) >= 2: # check divergence of value
break
else:
c = c ** 2 + c0
else:
plt.plot(round(c0.real,2), round(c0.imag,2), 'k,') # print(c0)
start += stepx
re = start.real
im = start.imag
c = complex(re,im)
c0 = complex(re,im)
start = start + stepy - complex(steps * stepx.real,0)
mandelbrot()
plt.axis([-3,2,-2,2])
plt.grid(True)
plt.show()
Warum sollte das hier denn schneller berechnet werden so? Oder wie soll ich die Funktion gestalten? Sorry für die dummen Fragen, aber habe echt kein Plan von...
MfG
Verfasst: Samstag 24. April 2010, 18:27
von EyDu
Weil der Zugriff auf lokale Namen schneller geht als auf globale.
Dein "global" in Zeile 24 solltest du übrigens loswerden. Funktionen sollten so wenige Seiteneffekte haben wie möglich. Du möchtest Parameter und Rückgabewerte verwenden.
Verfasst: Samstag 24. April 2010, 20:44
von Huy
Code: Alles auswählen
#!/usr/bin/python
# Filename: mandelbrot.py
import matplotlib.pyplot as plt
def mandelbrot(steps, start, stepx, stepy, iterations):
c = start
c0 = start
for h in xrange(0,steps):
for k in xrange(0,steps): # do this for all points in the selected area
if abs(c) >= 2: # check divergence
start += stepx
else:
for i in xrange(0,iterations):
if abs(c ** 2 + c0) >= 2: # check divergence of value
break
else:
c = c ** 2 + c0
else:
plt.plot(round(c0.real,2), round(c0.imag,2), 'k,') # print(c0)
start += stepx
c = start
c0 = start
start = start + stepy - complex(steps * stepx.real,0)
mandelbrot(601, complex(-3,-3), complex(0.01,0), complex(0,0.01), 100)
plt.xlabel(r'$re(z)$')
plt.ylabel(r'$i \cdot im(z)$')
plt.title('The Mandelbrot Set')
plt.axis([-3,2,-2,2])
plt.grid(True)
plt.show()
So? Oder habe ich das immer noch falsch verstanden?
MfG
Verfasst: Samstag 24. April 2010, 22:21
von numerix
Huy hat geschrieben:So? Oder habe ich das immer noch falsch verstanden?
Hast du richtig verstanden. Die innere Schleife lässt sich allerdings noch weiter optimieren. Statt if .. else kannst du die Berechnung von c als ersten Schritt in der Schleife durchführen und danach nur mittels abs(c) das Abbruchkriterium prüfen. Wenn du dann noch multiplizierst statt zu potenzieren wird die Berechnung fast doppelt so schnell. Wenn du dann noch psyco hinzunimmst, ist die Berechnung gegenüber dem von dir gezeigten Code etwa um den Faktor 6 schneller. Das eigentliche Plotten ist da nicht enthalten, weil ich - wie gesagt - matplotlib nicht installiert habe.
Verfasst: Samstag 24. April 2010, 22:36
von Huy
Das Potenzieren und die if-else-Abfrage habe ich jetzt entsprechend verändert. Habe noch versucht psyco zu installieren, allerdings rausgefunden, dass das auf 64bit-System nicht läuft, bzw. es nicht unbedingt deutlich schneller ist als sonst, da das Script auf nem 64bit-System ohnehin ähnlich schnell läuft wie auf nem 32bit-System + psyco.
Gibt es noch irgendetwas zu optimieren? Wie kann man jetzt die Zoomfunktion verändern, bzw. wie verändere ich, was matplotlib.pyplot "macht"? Oder ist das im Endeffekt zu kompliziert?
MfG
Verfasst: Sonntag 25. April 2010, 11:02
von HerrHagen
Matplotlib ist schnell genug. Das es so langsam läuft, hat IMO andere Ursachen.
Dein Hauptproblem ist das du matplotlib mit plt.plot jeden Punkt einzeln übergibst. So ist die Verwendung dieser Funktion einfach nicht gedacht. Die Punkte werden in vielfacher Ausführung in einer Datenstruktur gehalten die für diese Aufgabe nicht geeignet ist. Du erzeugst bei jedem Aufruf von plt.plot ein neues matplotlib-Objekt. Schau dir noch mal das folgende Beispiel an:
http://www.scipy.org/Tentative_NumPy_Tu ... et_Example
Hier wird das Ergebnis komplett in einer sehr effizienten Datenstruktur gehalten, nämlich einem numpy-array (entspricht c-Array + Zusatzinformationen) und dann als ganzes matplotlib übergeben. Entsprechend schnell liegt das Ergebnis vor.
Verfasst: Sonntag 25. April 2010, 11:34
von Huy
Ich verstehe leider das Beispiel auf scipy nicht wirklich; also was genau gerechnet wird etc. Habe mal versucht, die Daten in einem numpy-Array festzuhalten, allerdings bin ich nicht ganz sicher, ob ich das richtig mache. Ich definiere zuerst ein Array mit empty() und füge dann immer Daten hinzu. Allerdings muss ich beim empty() Befehl die Dimensionen des Arrays festlegen, wobei ich ja keine Ahnung habe, wie lange es wird... Ist das schon die falsche Richtung? Wie muss ich das angehen?
MfG
Verfasst: Sonntag 25. April 2010, 11:43
von BlackJack
@Huy: Natürlich weisst Du wie gross das Array werden muss -- Du weist doch wieviele Punkte Du in X- und Y-Richtung berechnest.
Verfasst: Sonntag 25. April 2010, 11:47
von Huy
Ehm, nein. Bzw. ich weiss wieviele Punkte ich berechne, aber ich weiss ja nicht, wieviele Punkte nicht divergieren. Ob ein Punkt divergiert oder nicht, finde ich ja erst mittels for-Schleife heraus. Oder?
MfG