1.2.2 Mandelbrot集

接下來，我們介紹一個非常經典的并行計算的例子，并且將在本書中多次用到這個例子—— 一個生成Mandelbrot集圖形的算法。首先，讓我們來定義Mandelbrot集。

對于給定的復數c，當時，我們可以定義這樣一個遞歸序列，其中，而時，該序列可以表示為。如果n增加到無窮大，仍然以2為界，就說c是Mandelbrot集的元素。

回想一下，我們可以將復數表示在二維笛卡兒平面上，其中x軸表示實數分量，y軸表示虛數分量。因此，我們可以很容易地用一個令人印象深刻的（也是眾所周知的）圖形來可視化Mandelbrot集。這里，我們將在復笛卡兒平面上用較淺的陰影表示Mandelbrot集的元素，用較深的陰影表示不屬于Mandelbrot集的元素，具體如圖1-1所示。

圖1-1

現在，讓我們考慮一下如何用Python代碼來生成Mandelbrot集。因為我們不可能檢查每一個復數是否屬于Mandelbrot集，所以必須選擇檢查的范圍。同時，我們必須確定需要檢查每個范圍（由width、height確定）中的多少個點，還要在檢查|z_n|時，確定n的最大值（max_iters）。接下來，我們將編寫一個函數，用它來生成Mandelbrot集的圖形——就本例來說，我們通過連續迭代圖形中的每一個點來完成該任務。

首先，我們需要導入NumPy庫，它是本書中經常用到的一個數值運算庫。具體來說，這里的功能是通過simple_mandelbrot函數實現的。我們先用NumPy的linspace函數生成一個充當離散復平面的網格（下面的代碼應該是相當簡單的）：

import numpy as np

def simple_mandelbrot(width, height, real_low, real_high, imag_low,
imag_high, max_iters):
     real_vals = np.linspace(real_low, real_high, width)
     imag_vals = np.linspace(imag_low, imag_high, height)
     # we will represent members as 1, non-members as 0.
     mandelbrot_graph = np.ones((height,width), dtype=np.float32)
     for x in range(width):
         for y in range(height):
             c = np.complex64( real_vals[x] + imag_vals[y] * 1j )
             z = np.complex64(0)
             for i in range(max_iters):
                 z = z**2 + c
                 if(np.abs(z) > 2):
                     mandelbrot_graph[y,x] = 0
                     break
     return mandelbrot_graph

我們要添加一些代碼，以便將Mandelbrot集的圖形轉儲到一個PNG格式的文件中。為此，首先需要導入相應的庫：

from time import time
import matplotlib
# the following will prevent the figure from popping up
matplotlib.use('Agg')
from matplotlib import pyplot as plt

現在，讓我們添加一些代碼來生成Mandelbrot集，將其圖形轉儲到一個文件中，并使用time函數對兩個操作進行計時：

if __name__ == '__main__':
     t1 = time()
     mandel = simple_mandelbrot(512,512,-2,2,-2,2,256, 2)
     t2 = time()
     mandel_time = t2 - t1
     t1 = time()
     fig = plt.figure(1)
     plt.imshow(mandel, extent=(-2, 2, -2, 2))
     plt.savefig('mandelbrot.jpg', dpi=fig.dpi)
     t2 = time()
     dump_time = t2 - t1
     print 'It took {} seconds to calculate the Mandelbrot
graph.'.format(mandel_time)
     print 'It took {} seconds to dump the image.'.format(dump_time)

現在，讓我們運行這個程序（該程序也可以從本書配套資源的文件夾1下面的Mandelbrot0.py文件中找到），如圖1-2所示。

圖1-2

可以看到，生成Mandelbrot集耗時約14.62秒，轉儲圖形耗時約0.11秒。我們是以逐點方式生成Mandelbrot集的。不同點的坐標值之間并不存在依賴關系，因此生成Mandelbrot集實際上就是一個可并行操作。相比之下，轉儲圖形的代碼則無法并行化。

現在，讓我們用阿姆達爾定律來分析一下。就本例來說，如果將相應的代碼并行化，我們可以得到多大的加速比？如上所述，該程序的兩個部分總共用時約14.73秒。我們可以并行化生成Mandelbrot集的代碼，也就是說，可并行化的那部分代碼的執行時間占比為p=14.62/14.73≈0.99。因此，這個程序可并行化比例約為99%！

那么，我們可以獲得多大的加速比呢？筆者目前用的是一臺擁有640個內核的GTX 1050的筆記本電腦，使用阿姆達爾定律時，N將是640。所以，計算加速比的公式為

這個加速比的值無疑是非常好的，足以表明將算法并行化到GPU上的努力是非常值得的。記住，阿姆達爾定律只是給出了一個非常粗略的估計值！將計算任務轉移到GPU上時，我們還要考慮其他因素，比如在CPU和GPU之間發送和接收數據所需的額外時間，或者轉移到GPU上的算法只能部分并行，等等。