Mandelbrot simple mis en Python

Question

J'ai regardé d'autres questions à ce sujet, mais je n'arrive pas à comprendre où je vais mal. Le but est d'écrire un programme pour faire une image de l'ensemble de Mandelbrot en effectuant l'itération pour toutes les valeurs de c = x + iy sur une grille N × N couvrant la région où -2 ≤ x ≤ 2 et -2 ≤ y ≤ 2. Faire un tracé de densité dans lequel les points de la grille à l'intérieur de l'ensemble Mandelbrot sont colorés en noir et ceux à l'extérieur sont en blanc. " Lorsque l'on considère que l'on est dans l'ensemble de Mandelbrot, l'amplitude de z n'est jamais supérieure à 2 lorsqu'elle est répétée comme z '= z^2 + c. Mon code produit une grille presque entièrement noire avec une petite entaille de blanc.

from pylab import imshow,show,gray 
from numpy import zeros,linspace 

z = 0 + 0j 
n=100 

M = zeros([n,n],int) 
xvalues = linspace(-2,2,n) 
yvalues = linspace(-2,2,n) 

for x in xvalues: 
    for y in yvalues: 
     c = complex(x,y) 
     for i in range(100): 
      z = z*z + c 
      if abs(z) > 2.0: 
       M[y,x] = 1 
       break 

imshow(M,origin="lower") 
gray() 
show() 

Pour tout futurs lecteurs, voici comment mon nouveau code a fini à la recherche:

from pylab import imshow,show,gray 
from numpy import zeros,linspace 

n=1000 

M = zeros([n,n],int) 
xvalues = linspace(-2,2,n) 
yvalues = linspace(-2,2,n) 

for u,x in enumerate(xvalues): 
    for v,y in enumerate(yvalues): 
     z = 0 + 0j 
     c = complex(x,y) 
     for i in range(100): 
      z = z*z + c 
      if abs(z) > 2.0: 
       M[v,u] = 1 
       break 

imshow(M,origin="lower") 
gray() 
show() 

Answer 1

Il y a quelques problèmes avec votre code.

Tout d'abord, vous utilisez xvalues et yvalues à l'index M, mais ces indices devraient être des entiers d'index de pixels dans la plage 0 .. (n-1). Numpy convertira les flottants en xvalues et en yvalues en nombres entiers, mais les nombres qui en résulteront seront en -2..2, donc il n'y aura pas beaucoup de pixels tracés, l'image sera minuscule, et vous obtiendrez l'emballage en raison du Les index négatifs fonctionnent en Python. Un moyen simple d'obtenir les indices de pixels requis est d'utiliser la fonction Python intégrée enumerate, mais il existe peut-être un moyen de réorganiser votre code pour ce faire avec les fonctions Numpy.

L'autre problème est que vous devez réinitialiser z à zéro pour chaque pixel. Actuellement, votre code réutilise la dernière z pour le pixel précédent, et si ce pixel était dans l'ensemble Mandelbrot, alors z sera trop grand.

Voici une version réparée de votre code. Je n'ai pas de pylab, j'ai donc écrit un visualiseur bitmap simple en utilisant PIL. Vous pouvez les enregistrer dans un fichier en appelant img.save(filename) dans ma fonction show; PIL va trouver le bon format de fichier à partir de l'extension du nom de fichier.

import numpy as np 
from PIL import Image 

def show(data): 
    img = Image.frombytes('1', data.shape[::-1], np.packbits(data, 1)) 
    img.show() 

n = 100 
maxiter = 100 

M = np.zeros([n, n], np.uint8) 
xvalues = np.linspace(-2, 2, n) 
yvalues = np.linspace(-2, 2, n) 

for u, x in enumerate(xvalues): 
    for v, y in enumerate(yvalues): 
     z = 0 
     c = complex(x, y) 
     for i in range(maxiter): 
      z = z*z + c 
      if abs(z) > 2.0: 
       M[v, u] = 1 
       break 

show(M) 

est ici l'image de sortie:

---

Bien sûr, chaque fois que vous vous trouvez itérer sur indices de tableau numpy, qui est un signe que vous faites mal. Le point principal de l'utilisation de Numpy est qu'il peut effectuer des opérations sur des tableaux entiers à la fois en itérant en interne sur eux à la vitesse C; Le code ci-dessus pourrait également utiliser des listes Python au lieu de tableaux Numpy.

Voici une version qui permet à Numpy de faire la majeure partie de la boucle. Il utilise plus de RAM, mais il est environ 2,5 fois plus rapide que la version précédente, et il est un peu plus court.

Ce code utilise plusieurs tableaux 2D. c contient tous les numéros de semences complexes, nous effectuons le calcul de base Mandelbrot en z, qui est initialisé à zéro. mask est un tableau booléen qui contrôle où le calcul Mandelbrot doit être effectué.Tous ses éléments sont initialement définis sur True et à chaque itération True éléments dans mask qui correspondent aux éléments dans z qui ont échappé à l'ensemble Mandelbrot sont définis sur False. Pour tester si un point s'est échappé, nous utilisons z.real**2 + z.imag**2 > 4.0 au lieu de abs(z) > 2.0, ce qui économise un peu de temps car il évite le calcul coûteux de la racine carrée et l'appel de fonction abs.

Nous pouvons utiliser la valeur finale de mask pour tracer l'ensemble Mandelbrot, mais pour rendre les points de l'ensemble Mandelbrot noirs, nous devons inverser ses valeurs, ce que nous pouvons faire avec 1 - mask.

import numpy as np 
from PIL import Image 

def show(data): 
    img = Image.frombytes('1', data.shape[::-1], np.packbits(data, 1)) 
    img.show() 
    img.save('mset.png') 

n = 100 
maxiter = 100 

a = np.linspace(-2, 2, n) 
c = a + 1.j * a[:, None] 
z = np.zeros((n, n), np.complex128) 
mask = np.ones((n, n), np.bool) 

for i in range(maxiter): 
    mask[mask] = z[mask].real**2 + z[mask].imag**2 < 4.0 
    z[mask] = z[mask]**2 + c[mask] 

show(1 - mask)