楔子

上一篇文章我們探討了 GIL 的原理，以及如何釋放 GIL 實現并行，做法是將函數聲明為 nogil，然后使用 with nogil 上下文管理器即可。在使用上非常簡單，但如果我們想讓循環(huán)也能夠并行執(zhí)行，那么該方式就不太方便了，為此 Cython 提供了一個 prange 函數，專門用于循環(huán)的并行執(zhí)行。

這個 prange 的特殊功能是 Cython 獨一無二的，并且 prange 只能與 for 循環(huán)搭配使用，不能獨立存在。

Cython 使用 OpenMP API 實現 prange，用于多平臺共享內存的處理。但 OpenMP 需要 C 或者 C++ 編譯器支持，并且編譯時需要指定特定的編譯參數來啟動。例如：當我們使用 gcc 時，必須在編譯和鏈接二進制文件的時候指定一個 -fopenmp，以確保啟用 OpenMP。

許多編譯器均支持 OpenMP ，包括免費的和商業(yè)的。但 Clang/LLVM 則是一個最顯著的例外，它只在一個單獨的分支中得到了初步的支持，而為它完全實現的 OpenMP 還在開發(fā)當中。

而使用 prange，需要從 cython.parallel 中進行導入。但是在這之前，我們先來看一個例子：

import?numpy?as?np
from?cython?cimport?boundscheck,?wraparound

cdef?inline?double?norm2(double?complex?z)?nogil:
????"""
????接收一個復數?z,?計算它的模的平方
????由于 norm2 要被下面的?escape?函數多次調用
????這里通過?inline?聲明成內聯函數
????:param?z:?
????:return:?
????"""
????return?z.real?*?z.real?+?z.imag?*?z.imag


cdef?int?escape(double?complex?z,
????????????????double?complex?c,
????????????????double?z_max,
????????????????int?n_max)?nogil:
????"""
????這個函數具體做什么,?不是我們的重點
????我們不需要關心
????"""
????cdef:
????????int?i?=?0
????????double?z_max2?=?z_max?*?z_max
????while?norm2(z)?<?z_max2?and?i?<?n_max:
????????z?=?z?*?z?+?c
????????i?+=?1
????return?i


@boundscheck(False)
@wraparound(False)
def?calc_julia(int?resolution,
???????????????double?complex?c,
???????????????double?bound=1.5,
???????????????double?z_max=4.0,
???????????????int?n_max=1000):
????"""
????我們將要在?Python?中調用的函數
????"""
????cdef:
????????double?step?=?2.0?*?bound?/?resolution
????????int?i,?j
????????double?complex?z
????????double?real,?imag
????????int[:,?::?1]?counts
????counts?=?np.zeros((resolution?+?1,?resolution?+?1),?dtype="int32")

????for?i?in?range(resolution?+?1):
????????real?=?-bound?+?i?*?step
????????for?j?in?range(resolution?+?1):
????????????imag?=?-bound?+?j?*?step
????????????z?=?real?+?imag?*?1j
????????????counts[i,?j]?=?escape(z,?c,?z_max,?n_max)

????return?np.array(counts,?copy=False)

我們手動編譯一下，然后調用 calc_julia 函數，這個函數做什么不需要關心，我們只需要將注意力放在那兩層 for 循環(huán)（準確的說是外層循環(huán)）上即可，這里我們采用手動編譯的形式。

import?cython_test
import?numpy?as?np
import?matplotlib.pyplot?as?plt
arr?=?cython_test.calc_julia(1000,?0.322?+?0.05j)
plt.imshow(np.log(arr))
plt.show()

那么 calc_julia 這個函數耗時多少呢？我們來測試一下：

使用 prange

對于上面的代碼來說，外層循環(huán)里面的邏輯是彼此獨立的，即當前循環(huán)不依賴上一層循環(huán)的結果，因此這非常適合并行執(zhí)行。所以 prange 便閃亮登場了，我們只需要做簡單的修改即可：

import?numpy?as np
from?cython?cimport boundscheck,?wraparound
from?cython.parallel?cimport prange

cdef?inline?double?norm2(double?complex?z)?nogil:
????return?z.real?*?z.real?+?z.imag?*?z.imag


cdef?int?escape(double?complex?z,
????????????????double?complex?c,
????????????????double?z_max,
????????????????int?n_max)?nogil:
????cdef:
????????int?i?=?0
????????double?z_max2?=?z_max?*?z_max
????while?norm2(z)?<?z_max2?and?i?<?n_max:
????????z?=?z?*?z?+?c
????????i?+=?1
????return?i


@boundscheck(False)
@wraparound(False)
def?calc_julia(int?resolution,
???????????????double?complex?c,
???????????????double?bound=1.5,
???????????????double?z_max=4.0,
???????????????int?n_max=1000):
????cdef:
????????double?step?=?2.0?*?bound?/?resolution
????????int?i,?j
????????double?complex?z
????????double?real,?imag
????????int[:,?::?1]?counts
????counts?=?np.zeros((resolution?+?1,?resolution?+?1),?dtype="int32")
    #?只需要將外層的 range 換成 prange
????for?i?in?prange(resolution?+?1, nogil=True):
????????real?=?-bound?+?i?*?step
????????for?j?in?range(resolution?+?1):
????????????imag?=?-bound?+?j?*?step
????????????z?=?real?+?imag?*?1j
????????????counts[i,?j]?=?escape(z,?c,?z_max,?n_max)

????return?np.array(counts,?copy=False)

我們只需要將外層循環(huán)的 range 換成 prange 即可，里面指定 nogil=True，便可實現并行的效果，至于這個函數的其它參數以及用法后面會說。而且一旦使用了 prange，那么在編譯的時候，必須啟用 OpenMP，下面看一下編譯腳本。

from?distutils.core?import?setup,?Extension
from?Cython.Build?import?cythonize

ext?=?[Extension("cython_test",
?????????????????sources=["cython_test.pyx"],
?????????????????extra_compile_args=["-fopenmp"],
?????????????????extra_link_args=["-fopenmp"])]

setup(ext_modules=cythonize(ext,?language_level=3))

編譯測試一下：

我們看到效率大概是提升了兩倍，因為我 Windows 上使用的不是 gcc，所以這里是在 CentOS 上演示的。而我的 CentOS 服務器只有兩個核，因此效率提升大概兩倍左右。

所以只是做了一些非常簡單的修改，便可帶來如此巨大的性能提升，簡直妙啊。prange 是要搭配 for 循環(huán)來使用的，如果 for 循環(huán)內部的邏輯彼此獨立，即第二層循環(huán)不依賴第一層循環(huán)的某些結果，那么不妨使用 prange 吧。

注意還沒完，我們還能做得更好，下面就來看看 prange 里面的其它的參數，這樣我們能更好利用 prange 的并行特性。

prange 的其它參數

prange 函數的原型如下：

#?第一個參數?self?我們不需要管
#?prange?實際上是類?CythonDotParallel?的成員函數
#?因為?Cython?內部執(zhí)行了下面這行邏輯
#?sys.modules['cython.parallel']?=?CythonDotParallel()
#?所以它將一個實例對象變成了一個模塊

def?prange(self,?start=0,?stop=None,?step=1,?
???????????nogil=False,?schedule=None,?
???????????chunksize=None,?num_threads=None):

我們先來看前三個參數，start、stop、step。

• prange(3): 相當于 start=0、stop=3；
• prange(1, 3): 相當于 start=1、stop=3；
• prange(1, 3, 2): 相當于 start=1、stop=3、step=2；

類似于 range，同樣不包含結尾 stop。

然后是第四個參數 nogil，它默認是 False，但事實上我們必須將其設置為 True，否則會報出編譯錯誤。

然后剩下的三個參數，如果我們不指定的話，那么 Cython 編譯器采取的策略是將整個循環(huán)分成多個大小相同的連續(xù)塊，然后給每一個可用線程一個塊。然而這個策略實際上并不是最好的，因為每一層循環(huán)用的時間不一定一樣，如果一個線程很快就完成了，那么不就造成資源上的浪費了嗎？

我們修改一下，將 schedule 指定為 "static"，chunksize 指定為 1：

for?i?in?prange(resolution?+?1,?nogil=True,?
????????????????schedule="static",?chunksize=1):

其它地方不變，只是加兩個參數，然后重新測試一下。

我們看到效率上是差不多的，原因是我的機器只有兩個核，如果核數再多一些的話，那么速度就會明顯地提升。

下面來解釋一下剩余的三個參數的含義，首先是 schedule，它有以下幾個選項：

"static"

整個循環(huán)在編譯時會以一種固定的方式分配給多個線程，如果 chunksize 沒有指定，那么會分成 num_threads 個連續(xù)塊，一個線程一個塊。如果指定了 chunksize，那么每一塊會以輪詢調度算法（Round Robin）交給線程進行處理，適用于任務均勻分布的情況。

"dynamic"

線程在運行時動態(tài)地向調度器申請下一個塊，chunksize 默認為 1，當任務負載不均時，動態(tài)調度是最佳的選擇。

"guided"

塊是動態(tài)分布的，就像 dynamic 一樣，但這與 dynamic 還不同，chunksize 的比例不是固定的，而是和 剩余迭代次數 / 線程數 成比例關系。

"runtime"

不常用。

控制 schedule 和 chunksize 可以方便地探索不同的并行執(zhí)行策略、以及工作負載分配，通常指定 schedule 為 "static"，加上設置一個合適的 chunksize 是最好的選擇。而 dynamic 和 guided 適用于動態(tài)變化的執(zhí)行上下文，但會導致運行時開銷。

當然還有最后一個參數 num_threads，很明顯不需要解釋，就是使用的線程數量。如果不指定，那么 prange 會使用盡可能多的線程。所以我們只是做了一點修改，便可以帶來巨大的性能提升，這種性能提升與 Cython 在純 Python 上帶來的性能提升成倍增關系。

在reductions操作上使用prange

我們經常會循環(huán)遍歷數組計算它們的累和、累積等等，這種數據量減少的操作我們稱之為 reduction 操作。而 prange 對這樣的操作也是支持并行執(zhí)行的，我們舉個例子：

from?cython?cimport boundscheck,?wraparound

@boundscheck(False)
@wraparound(False)
def calc_julia(int?[:,?::?1]?counts,
???????????????int?val):
????cdef:
????????int?total?=?0
????????int?i,?j,?M,?N
????N,?M?=?counts.shape[:?2]
????for?i?in?range(M):
????????for?j?in?range(N):
????????????if?counts[i,?j]?==?val:
????????????????total?+=?1

????return?total?/?float(counts.size)

顯然我們是希望計算一個數組中值 val 的元素的個數，下面測試一下：

如果改成 prange 的話，會有什么效果呢？代碼的其它部分不變，只需要導入 prange，然后將 range(M) 改成 prange(M, nogil=True) 即可。

速度比原來快了兩倍多，還是很可觀的，如果你的 CPU 是多核的，那么效率提升會更明顯。

這里我們沒有使用 schedule 和 chunksize 參數，你也可以加上去。當然啦，如果占用內存過大的話，它可能無法像預期的一樣顯著地提升性能，因為 prange 的優(yōu)化重點是在 CPU 上面。

但是可能有人會有疑問，多個線程同時對 total 變量進行自增操作，這么做不會造成沖突嗎？答案是不會的，因為加法是可交換的，即無論是 a + b 還是 b + a，結果都是相同的。Cython（通過 OpenMP）生成線程代碼，每個線程計算循環(huán)子集的和，然后所有線程再將各自的和匯總在一起。

如果是交給 Numpy 來做的話，那么等價于如下：

np.sum(counts?==?val)?/?float(counts.size)

但是效率如何呢？我們來對比一下：

我們采用并行計算用的是 6.13 毫秒，Numpy 用的是 20 毫秒，看樣子是我們贏了，并且 CPU 核心數越多，差距越明顯，這便是并行計算的威力。當然對于這種算法來說，還是直接交給 Numpy 吧，畢竟人家都幫你封裝好了，一個函數調用就可以解決了。

因此有效利用計算機硬件資源確實是最直接的辦法。

并行編程的局限性

雖然 Cython 的 prange 容易使用，但其實還是有局限性的，當然這個局限性和 Cython無關，因為理想化的并行擴展本身就是一個難以實現的事情。我們舉個例子：

def filter(nrows, ncols):
????for?i?in?range(nrows):
????????for?j?in?range(ncols):
????????????b[i,?j]?=?(a[i,?j]?+?a[i?-?1,?j]?+?a[i?+?1,?j]?+
???????????????????????a[i,?j?-?1]?+?a[i,?j?+?1])?/?5.0)

假設我們要做一個過濾器，計算每一個點加上它周圍的四個點的平均值。但如果這里將外層的 range 換成 prange，那么它的整體性能不會明顯提升。因為內層循環(huán)訪問的是不連續(xù)的數組元素，由于缺乏數據本地性，CPU 的緩存無法生效，反而導致 prange 變慢。

那么我們什么時候使用 prange 呢？遵循以下法則即可：

• 1. prange 能夠很好的利用 CPU 并行操作, 這一點我們已經說過了；
• 2. 非本地讀寫的那些和內存綁定的操作很難提高速度；
• 3. 用較少的線程更容易實現加速, 因為對于 CPU 密集而言, 即便指定了超越核心數的線程也是沒有意義的；
• 4. 使用優(yōu)化的線程并行庫是將 CPU 所有核心都用于常規(guī)計算的最佳方式；

當然，其實我們在開發(fā)的時候是可以隨時使用 prange 的，只要循環(huán)體不和 Python 對象進行交互即可。

小結

Cython 允許我們繞過全局解釋器鎖，只要我們把和 Python 無關的代碼分離出來即可。對于那些不需要和 Python 交互的 C 代碼，可以輕松地使用 prange 實現基于線程的并行。

在其它語言中，基于線程的并行很容易出錯，并且難以正確處理。而 Cython 的 prange 則不需要我們在這方面費心，能夠輕松地處理很多性能瓶頸。

到此這篇關于使用Cython中prange函數實現for循環(huán)的并行的文章就介紹到這了,更多相關Cython prange for循環(huán)并行內容請搜索腳本之家以前的文章或繼續(xù)瀏覽下面的相關文章希望大家以后多多支持腳本之家！

最新評論