陣列迭代#
注意
陣列支援迭代器協定,並且可以像 Python 列表一樣進行迭代。請參閱快速入門指南中的索引、切片和迭代章節,以了解基本用法和範例。本文件的其餘部分介紹了 nditer
物件,並涵蓋更進階的用法。
迭代器物件 nditer
,在 NumPy 1.6 中引入,提供了許多彈性的方式,可以系統性地訪問一個或多個陣列的所有元素。本頁介紹了一些在 Python 中使用此物件對陣列進行運算的基本方法,然後總結了如何在 Cython 中加速內部迴圈。由於 nditer
的 Python 介面是 C 陣列迭代器 API 的相對直接的映射,這些想法也將有助於在 C 或 C++ 中進行陣列迭代。
單一陣列迭代#
使用 nditer
可以完成的最基本任務是訪問陣列的每個元素。每個元素都使用標準 Python 迭代器介面逐一提供。
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3)
>>> for x in np.nditer(a):
... print(x, end=' ')
...
0 1 2 3 4 5
對於此迭代,需要注意的重要一點是,順序的選擇是為了匹配陣列的記憶體佈局,而不是使用標準的 C 或 Fortran 順序。這樣做是為了提高訪問效率,反映了預設情況下,人們只想訪問每個元素,而不必擔心特定的順序。我們可以通過迭代先前陣列的轉置來看到這一點,並與以 C 順序複製該轉置進行比較。
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3)
>>> for x in np.nditer(a.T):
... print(x, end=' ')
...
0 1 2 3 4 5
>>> for x in np.nditer(a.T.copy(order='C')):
... print(x, end=' ')
...
0 3 1 4 2 5
a 和 a.T 的元素都以相同的順序遍歷,即它們在記憶體中儲存的順序,而 a.T.copy(order=’C’) 的元素則以不同的順序訪問,因為它們已被放入不同的記憶體佈局中。
控制迭代順序#
有時,以特定順序訪問陣列的元素非常重要,而與元素在記憶體中的佈局無關。nditer
物件提供了一個 order 參數來控制迭代的這個方面。預設值,具有上述行為,是 order='K' 以保持現有的順序。這可以用 order='C' 來覆蓋為 C 順序,以及 order='F' 為 Fortran 順序。
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3)
>>> for x in np.nditer(a, order='F'):
... print(x, end=' ')
...
0 3 1 4 2 5
>>> for x in np.nditer(a.T, order='C'):
... print(x, end=' ')
...
0 3 1 4 2 5
修改陣列值#
預設情況下,nditer
將輸入運算元視為唯讀物件。為了能夠修改陣列元素,您必須使用每個運算元的 ‘readwrite’ 或 ‘writeonly’ 標誌指定讀寫或僅寫模式。
然後,nditer 將產生可寫入的緩衝區陣列,您可以修改它們。但是,由於 nditer 必須在迭代完成後將此緩衝區資料複製回原始陣列,因此您必須通過兩種方法之一發出迭代結束的信號。您可以選擇
一旦調用 close
或退出其上下文,nditer 就無法再被迭代。
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3)
>>> a
array([[0, 1, 2],
[3, 4, 5]])
>>> with np.nditer(a, op_flags=['readwrite']) as it:
... for x in it:
... x[...] = 2 * x
...
>>> a
array([[ 0, 2, 4],
[ 6, 8, 10]])
如果您正在編寫需要支援舊版本 numpy 的程式碼,請注意在 1.15 之前,nditer
不是上下文管理器,也沒有 close
方法。相反,它依賴析構函數來啟動緩衝區的寫回。
使用外部迴圈#
在到目前為止的所有範例中,a 的元素都是由迭代器一次提供一個,因為所有的迴圈邏輯都在迭代器內部。雖然這很簡單方便,但效率不高。更好的方法是將一維最內層迴圈移動到您的程式碼中,在迭代器外部。這樣,NumPy 的向量化運算可以用於訪問的元素的大塊資料上。
nditer
將嘗試為內部迴圈提供盡可能大的塊。通過強制 'C' 和 'F' 順序,我們得到不同的外部迴圈大小。此模式通過指定迭代器標誌啟用。
請注意,在保持原生記憶體順序的預設情況下,迭代器能夠提供單個一維塊,而當強制 Fortran 順序時,它必須提供三個每個包含兩個元素的塊。
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3)
>>> for x in np.nditer(a, flags=['external_loop']):
... print(x, end=' ')
...
[0 1 2 3 4 5]
>>> for x in np.nditer(a, flags=['external_loop'], order='F'):
... print(x, end=' ')
...
[0 3] [1 4] [2 5]
追蹤索引或多重索引#
在迭代期間,您可能希望在計算中使用目前元素的索引。例如,您可能希望以記憶體順序訪問陣列的元素,但使用 C 順序、Fortran 順序或多維索引來查找不同陣列中的值。
索引由迭代器物件本身追蹤,並且可以通過 index 或 multi_index 屬性訪問,具體取決於請求的內容。下面的範例顯示了演示索引進展的輸出
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3)
>>> it = np.nditer(a, flags=['f_index'])
>>> for x in it:
... print("%d <%d>" % (x, it.index), end=' ')
...
0 <0> 1 <2> 2 <4> 3 <1> 4 <3> 5 <5>
>>> it = np.nditer(a, flags=['multi_index'])
>>> for x in it:
... print("%d <%s>" % (x, it.multi_index), end=' ')
...
0 <(0, 0)> 1 <(0, 1)> 2 <(0, 2)> 3 <(1, 0)> 4 <(1, 1)> 5 <(1, 2)>
>>> with np.nditer(a, flags=['multi_index'], op_flags=['writeonly']) as it:
... for x in it:
... x[...] = it.multi_index[1] - it.multi_index[0]
...
>>> a
array([[ 0, 1, 2],
[-1, 0, 1]])
追蹤索引或多重索引與使用外部迴圈不相容,因為它需要每個元素不同的索引值。如果您嘗試組合這些標誌,nditer
物件將引發例外。
範例
>>> import numpy as np
>>> a = np.zeros((2,3))
>>> it = np.nditer(a, flags=['c_index', 'external_loop'])
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
ValueError: Iterator flag EXTERNAL_LOOP cannot be used if an index or multi-index is being tracked
替代迴圈和元素訪問#
為了使其屬性在迭代期間更易於訪問,nditer
有一個用於迭代的替代語法,它顯式地與迭代器物件本身一起工作。使用這種迴圈結構,可以通過索引到迭代器來訪問當前值。其他屬性,例如追蹤的索引保持不變。下面的範例產生與上一節中相同的結果。
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3)
>>> it = np.nditer(a, flags=['f_index'])
>>> while not it.finished:
... print("%d <%d>" % (it[0], it.index), end=' ')
... is_not_finished = it.iternext()
...
0 <0> 1 <2> 2 <4> 3 <1> 4 <3> 5 <5>
>>> it = np.nditer(a, flags=['multi_index'])
>>> while not it.finished:
... print("%d <%s>" % (it[0], it.multi_index), end=' ')
... is_not_finished = it.iternext()
...
0 <(0, 0)> 1 <(0, 1)> 2 <(0, 2)> 3 <(1, 0)> 4 <(1, 1)> 5 <(1, 2)>
>>> with np.nditer(a, flags=['multi_index'], op_flags=['writeonly']) as it:
... while not it.finished:
... it[0] = it.multi_index[1] - it.multi_index[0]
... is_not_finished = it.iternext()
...
>>> a
array([[ 0, 1, 2],
[-1, 0, 1]])
緩衝陣列元素#
當強制迭代順序時,我們觀察到外部迴圈選項可能會以較小的塊提供元素,因為無法以恆定步幅以適當的順序訪問元素。在編寫 C 程式碼時,這通常沒問題,但是在純 Python 程式碼中,這可能會導致效能顯著下降。
通過啟用緩衝模式,迭代器提供給內部迴圈的塊可以變得更大,從而顯著減少 Python 直譯器的開銷。在強制 Fortran 迭代順序的範例中,啟用緩衝時,內部迴圈可以一次看到所有元素。
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3)
>>> for x in np.nditer(a, flags=['external_loop'], order='F'):
... print(x, end=' ')
...
[0 3] [1 4] [2 5]
>>> for x in np.nditer(a, flags=['external_loop','buffered'], order='F'):
... print(x, end=' ')
...
[0 3 1 4 2 5]
以特定資料類型迭代#
有時,有必要將陣列視為與其儲存的資料類型不同的資料類型。例如,即使操作的陣列是 32 位浮點數,也可能希望對 64 位浮點數執行所有計算。除非編寫低階 C 程式碼,否則通常最好讓迭代器處理複製或緩衝,而不是在內部迴圈中自己轉換資料類型。
有兩種機制可以實現這一點,臨時複製和緩衝模式。使用臨時複製,將使用新的資料類型製作整個陣列的副本,然後在副本中完成迭代。通過一種模式允許寫入訪問,該模式在所有迭代完成後更新原始陣列。臨時複製的主要缺點是臨時副本可能會消耗大量記憶體,特別是如果迭代資料類型的 itemsize 比原始資料類型大。
緩衝模式減輕了記憶體使用問題,並且比製作臨時副本更具快取友善性。除了特殊情況外,在迭代器外部一次需要整個陣列的情況下,建議使用緩衝而不是臨時複製。在 NumPy 中,ufuncs 和其他函數使用緩衝來支援靈活的輸入,並具有最小的記憶體開銷。
在我們的範例中,我們將使用複數資料類型處理輸入陣列,以便我們可以對負數取平方根。如果不啟用複製或緩衝模式,如果資料類型不完全匹配,迭代器將引發例外。
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3) - 3
>>> for x in np.nditer(a, op_dtypes=['complex128']):
... print(np.sqrt(x), end=' ')
...
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
TypeError: Iterator operand required copying or buffering, but neither copying nor buffering was enabled
在複製模式下,'copy' 被指定為每個運算元的標誌。這樣做是為了以每個運算元的方式提供控制。緩衝模式被指定為迭代器標誌。
範例
>>> import numpy as np
>>> a = np.arange(6).reshape(2,3) - 3
>>> for x in np.nditer(a, op_flags=['readonly','copy'],
... op_dtypes=['complex128']):
... print(np.sqrt(x), end=' ')
...
1.7320508075688772j 1.4142135623730951j 1j 0j (1+0j) (1.4142135623730951+0j)
>>> for x in np.nditer(a, flags=['buffered'], op_dtypes=['complex128']):
... print(np.sqrt(x), end=' ')
...
1.7320508075688772j 1.4142135623730951j 1j 0j (1+0j) (1.4142135623730951+0j)
迭代器使用 NumPy 的類型轉換規則來確定是否允許特定的轉換。預設情況下,它強制執行 'safe' 類型轉換。這意味著,例如,如果您嘗試將 64 位浮點數陣列視為 32 位浮點數陣列,它將引發例外。在許多情況下,規則 'same_kind' 是最合理的規則,因為它允許從 64 位浮點數轉換為 32 位浮點數,但不允許從浮點數轉換為整數或從複數轉換為浮點數。
範例
>>> import numpy as np
>>> a = np.arange(6.)
>>> for x in np.nditer(a, flags=['buffered'], op_dtypes=['float32']):
... print(x, end=' ')
...
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
TypeError: Iterator operand 0 dtype could not be cast from dtype('float64') to dtype('float32') according to the rule 'safe'
>>> for x in np.nditer(a, flags=['buffered'], op_dtypes=['float32'],
... casting='same_kind'):
... print(x, end=' ')
...
0.0 1.0 2.0 3.0 4.0 5.0
>>> for x in np.nditer(a, flags=['buffered'], op_dtypes=['int32'], casting='same_kind'):
... print(x, end=' ')
...
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
TypeError: Iterator operand 0 dtype could not be cast from dtype('float64') to dtype('int32') according to the rule 'same_kind'
需要注意的一件事是在使用讀寫或僅寫運算元時轉換回原始資料類型。常見的情況是以 64 位浮點數實現內部迴圈,並使用 'same_kind' 類型轉換以允許處理其他浮點類型。雖然在唯讀模式下,可以提供整數陣列,但讀寫模式會引發例外,因為轉換回陣列會違反類型轉換規則。
範例
>>> import numpy as np
>>> a = np.arange(6)
>>> for x in np.nditer(a, flags=['buffered'], op_flags=['readwrite'],
... op_dtypes=['float64'], casting='same_kind'):
... x[...] = x / 2.0
...
Traceback (most recent call last):
File "<stdin>", line 2, in <module>
TypeError: Iterator requested dtype could not be cast from dtype('float64') to dtype('int64'), the operand 0 dtype, according to the rule 'same_kind'
廣播陣列迭代#
NumPy 有一套規則來處理形狀不同的陣列,當函數接受多個按元素組合的運算元時,就會應用這些規則。這稱為廣播。nditer
物件可以為您應用這些規則,當您需要編寫這樣的函數時。
作為範例,我們印出將一維和二維陣列一起廣播的結果。
範例
>>> import numpy as np
>>> a = np.arange(3)
>>> b = np.arange(6).reshape(2,3)
>>> for x, y in np.nditer([a,b]):
... print("%d:%d" % (x,y), end=' ')
...
0:0 1:1 2:2 0:3 1:4 2:5
當發生廣播錯誤時,迭代器會引發例外,其中包含輸入形狀以幫助診斷問題。
範例
>>> import numpy as np
>>> a = np.arange(2)
>>> b = np.arange(6).reshape(2,3)
>>> for x, y in np.nditer([a,b]):
... print("%d:%d" % (x,y), end=' ')
...
Traceback (most recent call last):
...
ValueError: operands could not be broadcast together with shapes (2,) (2,3)
迭代器分配的輸出陣列#
NumPy 函數中的常見情況是根據輸入的廣播分配輸出,並且額外有一個名為 'out' 的可選參數,當提供該參數時,結果將放置在其中。nditer
物件提供了一種方便的慣用語,可以非常輕鬆地支援此機制。
我們將通過建立一個 square
函數來展示其工作原理,該函數將其輸入平方。讓我們從一個最小的函數定義開始,不包括 'out' 參數支援。
範例
>>> import numpy as np
>>> def square(a):
... with np.nditer([a, None]) as it:
... for x, y in it:
... y[...] = x*x
... return it.operands[1]
...
>>> square([1,2,3])
array([1, 4, 9])
預設情況下,nditer
對於作為 None 傳入的運算元使用標誌 'allocate' 和 'writeonly'。這意味著我們只需向迭代器提供兩個運算元,它就會處理其餘部分。
當添加 'out' 參數時,我們必須顯式提供這些標誌,因為如果有人將陣列作為 'out' 傳入,迭代器將預設為 'readonly',並且我們的內部迴圈將失敗。'readonly' 是輸入陣列的預設值的原因是為了防止對意外觸發縮減運算感到困惑。如果預設值是 'readwrite',任何廣播運算也會觸發縮減,這是本文檔後面涵蓋的主題。
在此同時,我們也引入 'no_broadcast' 標誌,它將阻止輸出被廣播。這很重要,因為我們每個輸出只需要一個輸入值。聚合多個輸入值是一個縮減運算,需要特殊處理。由於必須在迭代器標誌中顯式啟用縮減,它已經會引發錯誤,但是禁用廣播產生的錯誤訊息對於終端使用者來說更容易理解。要了解如何將平方函數推廣到縮減,請查看關於 Cython 的章節中的平方和函數。
為了完整性,我們還將添加 'external_loop' 和 'buffered' 標誌,因為這些是您通常出於效能原因而需要的。
範例
>>> import numpy as np
>>> def square(a, out=None):
... it = np.nditer([a, out],
... flags = ['external_loop', 'buffered'],
... op_flags = [['readonly'],
... ['writeonly', 'allocate', 'no_broadcast']])
... with it:
... for x, y in it:
... y[...] = x*x
... return it.operands[1]
...
>>> square([1,2,3])
array([1, 4, 9])
>>> b = np.zeros((3,))
>>> square([1,2,3], out=b)
array([1., 4., 9.])
>>> b
array([1., 4., 9.])
>>> square(np.arange(6).reshape(2,3), out=b)
Traceback (most recent call last):
...
ValueError: non-broadcastable output operand with shape (3,) doesn't
match the broadcast shape (2,3)
外積迭代#
任何二元運算都可以像 outer
中那樣以外積方式擴展到陣列運算,並且 nditer
物件提供了一種通過顯式映射運算元的軸來完成此操作的方法。也可以使用 newaxis
索引來完成此操作,但我們將向您展示如何直接使用 nditer op_axes 參數來完成此操作,而無需中間視圖。
我們將執行一個簡單的外積,將第一個運算元的維度放在第二個運算元的維度之前。op_axes 參數需要每個運算元一個軸列表,並提供從迭代器的軸到運算元軸的映射。
假設第一個運算元是一維的,第二個運算元是二維的。迭代器將有三個維度,因此 op_axes 將有兩個 3 元素列表。第一個列表挑選出第一個運算元的一個軸,對於其餘的迭代器軸為 -1,最終結果為 [0, -1, -1]。第二個列表挑選出第二個運算元的兩個軸,但不應與第一個運算元中挑選出的軸重疊。其列表為 [-1, 0, 1]。輸出運算元以標準方式映射到迭代器軸,因此我們可以提供 None 而不是構造另一個列表。
內部迴圈中的運算是簡單的乘法。與外積相關的所有事情都由迭代器設定處理。
範例
>>> import numpy as np
>>> a = np.arange(3)
>>> b = np.arange(8).reshape(2,4)
>>> it = np.nditer([a, b, None], flags=['external_loop'],
... op_axes=[[0, -1, -1], [-1, 0, 1], None])
>>> with it:
... for x, y, z in it:
... z[...] = x*y
... result = it.operands[2] # same as z
...
>>> result
array([[[ 0, 0, 0, 0],
[ 0, 0, 0, 0]],
[[ 0, 1, 2, 3],
[ 4, 5, 6, 7]],
[[ 0, 2, 4, 6],
[ 8, 10, 12, 14]]])
請注意,一旦迭代器關閉,我們就無法訪問 operands
,並且必須使用在上下文管理器內部建立的參考。
縮減迭代#
每當可寫入運算元的元素少於完整的迭代空間時,該運算元正在進行縮減。nditer
物件要求任何縮減運算元都標記為讀寫,並且僅當提供 'reduce_ok' 作為迭代器標誌時才允許縮減。
對於一個簡單的範例,請考慮取陣列中所有元素的總和。
範例
>>> import numpy as np
>>> a = np.arange(24).reshape(2,3,4)
>>> b = np.array(0)
>>> with np.nditer([a, b], flags=['reduce_ok'],
... op_flags=[['readonly'], ['readwrite']]) as it:
... for x,y in it:
... y[...] += x
...
>>> b
array(276)
>>> np.sum(a)
276
當組合縮減和分配的運算元時,事情會稍微複雜一些。在迭代開始之前,任何縮減運算元都必須初始化為其起始值。以下是我們如何做到這一點,沿著 a 的最後一個軸取總和。
範例
>>> import numpy as np
>>> a = np.arange(24).reshape(2,3,4)
>>> it = np.nditer([a, None], flags=['reduce_ok'],
... op_flags=[['readonly'], ['readwrite', 'allocate']],
... op_axes=[None, [0,1,-1]])
>>> with it:
... it.operands[1][...] = 0
... for x, y in it:
... y[...] += x
... result = it.operands[1]
...
>>> result
array([[ 6, 22, 38],
[54, 70, 86]])
>>> np.sum(a, axis=2)
array([[ 6, 22, 38],
[54, 70, 86]])
要進行緩衝縮減,需要在設定期間進行另一次調整。通常,迭代器構造涉及將第一個資料緩衝區從可讀陣列複製到緩衝區。任何縮減運算元都是可讀的,因此可以讀取到緩衝區中。不幸的是,在此緩衝操作完成後初始化運算元將不會反映在迭代開始的緩衝區中,並且會產生垃圾結果。
迭代器標誌 “delay_bufalloc” 的存在是為了允許迭代器分配的縮減運算元與緩衝一起存在。設定此標誌後,迭代器將使其緩衝區保持未初始化狀態,直到收到重置,之後它將準備好進行常規迭代。以下是如果我們也啟用緩衝,則先前範例的外觀。
範例
>>> import numpy as np
>>> a = np.arange(24).reshape(2,3,4)
>>> it = np.nditer([a, None], flags=['reduce_ok',
... 'buffered', 'delay_bufalloc'],
... op_flags=[['readonly'], ['readwrite', 'allocate']],
... op_axes=[None, [0,1,-1]])
>>> with it:
... it.operands[1][...] = 0
... it.reset()
... for x, y in it:
... y[...] += x
... result = it.operands[1]
...
>>> result
array([[ 6, 22, 38],
[54, 70, 86]])
將內部迴圈放入 Cython 中#
那些希望從其低階運算中獲得真正良好效能的人應強烈考慮直接使用 C 中提供的迭代 API,但對於那些不熟悉 C 或 C++ 的人來說,Cython 是一個具有合理效能權衡的好中間地帶。對於 nditer
物件,這意味著讓迭代器處理廣播、dtype 轉換和緩衝,同時將內部迴圈交給 Cython。
對於我們的範例,我們將建立一個平方和函數。首先,讓我們用簡單的 Python 實現此函數。我們想要支援類似於 numpy sum
函數的 'axis' 參數,因此我們需要為 op_axes 參數構造一個列表。以下是它的外觀。
範例
>>> def axis_to_axeslist(axis, ndim):
... if axis is None:
... return [-1] * ndim
... else:
... if type(axis) is not tuple:
... axis = (axis,)
... axeslist = [1] * ndim
... for i in axis:
... axeslist[i] = -1
... ax = 0
... for i in range(ndim):
... if axeslist[i] != -1:
... axeslist[i] = ax
... ax += 1
... return axeslist
...
>>> def sum_squares_py(arr, axis=None, out=None):
... axeslist = axis_to_axeslist(axis, arr.ndim)
... it = np.nditer([arr, out], flags=['reduce_ok',
... 'buffered', 'delay_bufalloc'],
... op_flags=[['readonly'], ['readwrite', 'allocate']],
... op_axes=[None, axeslist],
... op_dtypes=['float64', 'float64'])
... with it:
... it.operands[1][...] = 0
... it.reset()
... for x, y in it:
... y[...] += x*x
... return it.operands[1]
...
>>> a = np.arange(6).reshape(2,3)
>>> sum_squares_py(a)
array(55.)
>>> sum_squares_py(a, axis=-1)
array([ 5., 50.])
為了對此函數進行 Cython 化,我們用專門針對 float64 dtype 的 Cython 程式碼替換內部迴圈 (y[…] += x*x)。啟用 'external_loop' 標誌後,提供給內部迴圈的陣列將始終是一維的,因此幾乎不需要進行檢查。
這是 sum_squares.pyx 的列表
import numpy as np
cimport numpy as np
cimport cython
def axis_to_axeslist(axis, ndim):
if axis is None:
return [-1] * ndim
else:
if type(axis) is not tuple:
axis = (axis,)
axeslist = [1] * ndim
for i in axis:
axeslist[i] = -1
ax = 0
for i in range(ndim):
if axeslist[i] != -1:
axeslist[i] = ax
ax += 1
return axeslist
@cython.boundscheck(False)
def sum_squares_cy(arr, axis=None, out=None):
cdef np.ndarray[double] x
cdef np.ndarray[double] y
cdef int size
cdef double value
axeslist = axis_to_axeslist(axis, arr.ndim)
it = np.nditer([arr, out], flags=['reduce_ok', 'external_loop',
'buffered', 'delay_bufalloc'],
op_flags=[['readonly'], ['readwrite', 'allocate']],
op_axes=[None, axeslist],
op_dtypes=['float64', 'float64'])
with it:
it.operands[1][...] = 0
it.reset()
for xarr, yarr in it:
x = xarr
y = yarr
size = x.shape[0]
for i in range(size):
value = x[i]
y[i] = y[i] + value * value
return it.operands[1]
在這台機器上,將 .pyx 檔案建置為模組看起來像下面這樣,但您可能需要找到一些 Cython 教學來告訴您系統配置的具體細節。
$ cython sum_squares.pyx
$ gcc -shared -pthread -fPIC -fwrapv -O2 -Wall -I/usr/include/python2.7 -fno-strict-aliasing -o sum_squares.so sum_squares.c
從 Python 直譯器執行此操作會產生與我們的原生 Python/NumPy 程式碼相同的答案。
範例
>>> from sum_squares import sum_squares_cy
>>> a = np.arange(6).reshape(2,3)
>>> sum_squares_cy(a)
array(55.0)
>>> sum_squares_cy(a, axis=-1)
array([ 5., 50.])
在 IPython 中進行少量計時表明,Cython 內部迴圈減少的開銷和記憶體分配比直接的 Python 程式碼和使用 NumPy 內建 sum 函數的表達式都提供了非常好的加速。
>>> a = np.random.rand(1000,1000)
>>> timeit sum_squares_py(a, axis=-1)
10 loops, best of 3: 37.1 ms per loop
>>> timeit np.sum(a*a, axis=-1)
10 loops, best of 3: 20.9 ms per loop
>>> timeit sum_squares_cy(a, axis=-1)
100 loops, best of 3: 11.8 ms per loop
>>> np.all(sum_squares_cy(a, axis=-1) == np.sum(a*a, axis=-1))
True
>>> np.all(sum_squares_py(a, axis=-1) == np.sum(a*a, axis=-1))
True