在Python中設置簡單的Mandelbrot

Question

我已經看過其他關於這個問題，但我似乎無法弄清楚我要去哪裏錯了。我們的目標是「編寫一個程序來生成Mandelbrot集的圖像，方法是對跨越-2≤x≤2和-2≤y區域的N×N網格執行c = x + iy的所有值的迭代≤2.製作一張密度圖，Mandelbrot集內的網格點顏色爲黑色，而外部的顏色爲白色。「當在Mandelbrot集合中被認爲當z'= z^2 + c迭代時，z的大小從不大於2。

我的代碼產生一個幾乎完全黑色的網格，並帶有一個白色的小切口。

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() 

對於任何未來的讀者，這是我新的代碼是如何結束了尋找：

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

有幾個與你的代碼的問題。

首先，您使用xvalues和yvalues來索引M，但這些索引應該是範圍爲0 ..（n-1）的像素索引整數。 Numpy會將xvalues和yvalues中的浮點數轉換爲整數，但結果數字將在-2..2，因此繪製的像素不會很多，圖像將很小，並且由於方式負指數在Python中起作用。

獲得所需像素索引的一種簡單方法是使用內置Python函數enumerate，但可能有一種方法可以重新組織您的代碼，以便使用Numpy函數執行此操作。

另一個問題是，您需要將z重置爲每個像素爲零。目前，您的代碼會重複使用前一個像素的最後一個z，如果該像素位於Mandelbrot集合中，那麼z將會太大。

下面是您的代碼的修復版本。我沒有pylab，所以我使用PIL編寫了一個簡單的位圖查看器。您可以通過在我的show函數中調用img.save(filename)將它們保存爲文件; PIL將從文件擴展名中找出正確的文件格式。

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) 

這裏的輸出圖像：

---

當然，無論何時你發現自己遍歷numpy的數組下標，這是你這樣做是錯誤的標誌。使用Numpy的主要觀點是它可以通過以C速度對它們進行內部迭代來一次對整個數組執行操作;上面的代碼可能會使用普通的Python列表而不是Numpy數組。

這是一個讓Numpy做大部分循環的版本。它使用更多的內存，但比前一版本快大約2.5倍，並且略短一些。

此代碼使用幾個二維數組。 c包含所有複雜的種子編號，我們在z執行核心Mandelbrot計算，它被初始化爲零。 mask是一個布爾數組，用於控制需要執行Mandelbrot計算的位置。其所有項目初始設置爲True，並且在每個迭代True項目mask中對應於z中已從Mandelbrot集合中跳出的項目設置爲False。

要測試一個點是否已經逃脫，我們使用z.real**2 + z.imag**2 > 4.0而不是abs(z) > 2.0，這可以節省一點時間，因爲它避免了昂貴的平方根計算和abs函數調用。

我們可以使用mask的最終值繪製Mandelbrot集合，但要將Mandelbrot集合中的點設置爲黑色，我們需要將其值反轉，我們可以使用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)