我想從帶有索引的短陣列中為大陣列賦值。簡單代碼如下:
import numpy as np
def assign_x():
a = np.zeros((int(3e6), 20))
index_a = np.random.randint(int(3e6), size=(int(3e6), 20))
b = np.random.randn(1000, 20)
for i in range(20):
index_b = np.random.randint(1000, size=int(3e6))
a[index_a[:, i], i] = b[index_b, i]
return a
%timeit x = assign_x()
# 2.79 s ± 18.4 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
我嘗試了其他可能相關的方法,例如np.take和 numba's jit,但似乎上面是最快的方法。也可以使用multiprocessing. 我有profile代碼,大部分時間在下面的行,因為它運行了很多次(這里是 20 次)
a[index_a[:, i], i] = b[index_b, i]
在使用之前我有沒有機會讓它更快multiprocessing?
uj5u.com熱心網友回復:
為什么這很慢
這很慢,因為記憶體訪問模式非常低效。事實上,隨機訪問很慢,因為處理器無法預測它們。因此,它會導致代價高昂的高速快取未命中(如果陣列不適合 L1/L2 高速快取),而提前預取資料是無法避免的。問題是陣列太大以適合快取:每個 457 MiB 和index_a156 KiB。因此,訪問通常在具有較高延遲的 L2 快取中完成,而對其他兩個陣列的訪問則在 RAM 中完成。這很慢,因為當前的 DDR RAM 具有巨大的延遲abb在典型的 PC 上為 60-100 ns。更糟糕的是:這種延遲在不久的將來可能不會小得多:自過去二十年以來,RAM 延遲并沒有太大變化。這被稱為記憶墻。另請注意,當請求隨機位置的值時,現代處理器從 RAM 中獲取通常為 64 位元組的完整高速快取行(導致僅浪費 56/64=87.5% 的帶寬)。最后,生成亂數是一個非常昂貴的程序,尤其是大整數,并且np.random.randint可以針對目標平臺生成 32 位或 64 位整數。
如何改善這一點
第一個改進是在最連續的維度上更喜歡間接,這通常是最后一個維度,因為a[:,i]它比a[i,:]. 您可以轉置陣列并交換索引值。但是,Numpy 轉置函式只回傳一個視圖,并沒有真正轉置記憶體中的陣列。因此,當前需要一個顯式副本。這里最好的方法是直接生成陣列,以便訪問高效(而不是使用昂貴的轉置)。請注意,您可以使用簡單精度,因此陣列可以更好地適應快取,但會降低精度。
這是一個回傳轉置陣列的示例:
import numpy as np
def assign_x():
a = np.zeros((20, int(3e6)))
index_a = np.random.randint(int(3e6), size=(20, int(3e6)))
b = np.random.randn(20, 1000)
for i in range(20):
index_b = np.random.randint(1000, size=int(3e6))
a[i, index_a[i, :]] = b[i, index_b]
return a
%timeit x = assign_x()
可以使用Numba進一步改進代碼,以便并行運行代碼(由于 RAM 延遲,一個內核不應該足以使記憶體飽和,但許多內核可以更好地使用它,因為可以同時進行多個提取)。此外,它可以幫助避免創建大型臨時陣列。
這是一個優化的 Numba 代碼:
import numpy as np
import numba as nb
import random
@nb.njit('float64[:,:]()', parallel=True)
def assign_x():
a = np.zeros((20, int(3e6)))
b = np.random.randn(20, 1000)
for i in nb.range(20):
for j in range(3_000_000):
index_a = random.randint(0, 3_000_000)
index_b = random.randint(0, 1000)
a[i, index_a] = b[i, index_b]
return a
%timeit x = assign_x()
以下是 10 核 Skylake Xeon 處理器的結果:
Initial code: 2798 ms
Better memory access pattern: 1741 ms
With Numba: 318 ms
請注意,理論上并行化最內層回圈會更快,因為一行a更有可能適合最后一級快取。但是,這樣做會導致競爭條件,只能通過 Numba 中尚不可用的原子存盤(在 CPU 上)有效修復。
Note that the final code does not scale well because it is memory-bound. This is due to 87.5% of the memory throughput being wasted as explained before. Additionally, on many processors (like all Intel and AMD-Zen processors) the write allocate cache policy force data being read from memory for each store in this case. This makes the computation much more inefficient raising the wasted throughput to 93.7%... AFAIK, there is no way to prevent this in Python. In C/C , the write allocate issue can be fixed using low-level instructions. The rule of thumb is avoid memory random access patterns on big arrays like the plague.
轉載請註明出處,本文鏈接:https://www.uj5u.com/yidong/436247.html
