楔子

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

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

Cython 使用 OpenMP API 實(shí)現(xiàn) prange，用于多平臺(tái)共享內(nèi)存的處理。但 OpenMP 需要 C 或者 C++ 編譯器支持，并且編譯時(shí)需要指定特定的編譯參數(shù)來(lái)啟動(dòng)。例如：當(dāng)我們使用 gcc 時(shí)，必須在編譯和鏈接二進(jìn)制文件的時(shí)候指定一個(gè) -fopenmp，以確保啟用 OpenMP。

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

而使用 prange，需要從 cython.parallel 中進(jìn)行導(dǎo)入。但是在這之前，我們先來(lái)看一個(gè)例子：

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

cdef?inline?double?norm2(double?complex?z)?nogil:
????"""
????接收一個(gè)復(fù)數(shù)?z,?計(jì)算它的模的平方
????由于 norm2 要被下面的?escape?函數(shù)多次調(diào)用
????這里通過(guò)?inline?聲明成內(nèi)聯(lián)函數(shù)
????: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:
????"""
????這個(gè)函數(shù)具體做什么,?不是我們的重點(diǎn)
????我們不需要關(guān)心
????"""
????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?中調(diào)用的函數(shù)
????"""
????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)

我們手動(dòng)編譯一下，然后調(diào)用 calc_julia 函數(shù)，這個(gè)函數(shù)做什么不需要關(guān)心，我們只需要將注意力放在那兩層 for 循環(huán)（準(zhǔn)確的說(shuō)是外層循環(huán)）上即可，這里我們采用手動(dòng)編譯的形式。

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 這個(gè)函數(shù)耗時(shí)多少呢？我們來(lái)測(cè)試一下：

使用 prange

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

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，便可實(shí)現(xiàn)并行的效果，至于這個(gè)函數(shù)的其它參數(shù)以及用法后面會(huì)說(shuō)。而且一旦使用了 prange，那么在編譯的時(shí)候，必須啟用 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))

編譯測(cè)試一下：

我們看到效率大概是提升了兩倍，因?yàn)槲?Windows 上使用的不是 gcc，所以這里是在 CentOS 上演示的。而我的 CentOS 服務(wù)器只有兩個(gè)核，因此效率提升大概兩倍左右。

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

注意還沒(méi)完，我們還能做得更好，下面就來(lái)看看 prange 里面的其它的參數(shù)，這樣我們能更好利用 prange 的并行特性。

prange 的其它參數(shù)

prange 函數(shù)的原型如下：

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

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

我們先來(lái)看前三個(gè)參數(shù)，start、stop、step。

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

類似于 range，同樣不包含結(jié)尾 stop。

然后是第四個(gè)參數(shù) nogil，它默認(rèn)是 False，但事實(shí)上我們必須將其設(shè)置為 True，否則會(huì)報(bào)出編譯錯(cuò)誤。

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

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

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

其它地方不變，只是加兩個(gè)參數(shù)，然后重新測(cè)試一下。

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

下面來(lái)解釋一下剩余的三個(gè)參數(shù)的含義，首先是 schedule，它有以下幾個(gè)選項(xiàng)：

"static"

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

"dynamic"

線程在運(yùn)行時(shí)動(dòng)態(tài)地向調(diào)度器申請(qǐng)下一個(gè)塊，chunksize 默認(rèn)為 1，當(dāng)任務(wù)負(fù)載不均時(shí)，動(dòng)態(tài)調(diào)度是最佳的選擇。

"guided"

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

"runtime"

不常用。

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

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

在reductions操作上使用prange

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

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)

顯然我們是希望計(jì)算一個(gè)數(shù)組中值 val 的元素的個(gè)數(shù)，下面測(cè)試一下：

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

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

這里我們沒(méi)有使用 schedule 和 chunksize 參數(shù)，你也可以加上去。當(dāng)然啦，如果占用內(nèi)存過(guò)大的話，它可能無(wú)法像預(yù)期的一樣顯著地提升性能，因?yàn)?prange 的優(yōu)化重點(diǎn)是在 CPU 上面。

但是可能有人會(huì)有疑問(wèn)，多個(gè)線程同時(shí)對(duì) total 變量進(jìn)行自增操作，這么做不會(huì)造成沖突嗎？答案是不會(huì)的，因?yàn)榧臃ㄊ强山粨Q的，即無(wú)論是 a + b 還是 b + a，結(jié)果都是相同的。Cython（通過(guò) OpenMP）生成線程代碼，每個(gè)線程計(jì)算循環(huán)子集的和，然后所有線程再將各自的和匯總在一起。

如果是交給 Numpy 來(lái)做的話，那么等價(jià)于如下：

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

但是效率如何呢？我們來(lái)對(duì)比一下：

我們采用并行計(jì)算用的是 6.13 毫秒，Numpy 用的是 20 毫秒，看樣子是我們贏了，并且 CPU 核心數(shù)越多，差距越明顯，這便是并行計(jì)算的威力。當(dāng)然對(duì)于這種算法來(lái)說(shuō)，還是直接交給 Numpy 吧，畢竟人家都幫你封裝好了，一個(gè)函數(shù)調(diào)用就可以解決了。

因此有效利用計(jì)算機(jī)硬件資源確實(shí)是最直接的辦法。

并行編程的局限性

雖然 Cython 的 prange 容易使用，但其實(shí)還是有局限性的，當(dāng)然這個(gè)局限性和 Cython無(wú)關(guān)，因?yàn)槔硐牖牟⑿袛U(kuò)展本身就是一個(gè)難以實(shí)現(xiàn)的事情。我們舉個(gè)例子：

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)

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

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

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

當(dāng)然，其實(shí)我們?cè)陂_發(fā)的時(shí)候是可以隨時(shí)使用 prange 的，只要循環(huán)體不和 Python 對(duì)象進(jìn)行交互即可。

小結(jié)

Cython 允許我們繞過(guò)全局解釋器鎖，只要我們把和 Python 無(wú)關(guān)的代碼分離出來(lái)即可。對(duì)于那些不需要和 Python 交互的 C 代碼，可以輕松地使用 prange 實(shí)現(xiàn)基于線程的并行。

在其它語(yǔ)言中，基于線程的并行很容易出錯(cuò)，并且難以正確處理。而 Cython 的 prange 則不需要我們?cè)谶@方面費(fèi)心，能夠輕松地處理很多性能瓶頸。

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

最新評(píng)論