for-loop + numpy的矢量化版本。其中

背景信息：我有大量（N）三維粒子。對於具有某些性質的所有粒子對 [i，j]，我計算幾何因子c [i，j]。然後，我想總結一個固定的i的所有對[i，j]的貢獻，並將其稱爲c [i]（並對所有粒子i重複此過程）。（N，N）維陣列C在位置[i，j]處的相關信息以及其他地方的很多零是相當公平的（相對於N^2） numpy很快，但在內存使用方面也非常低效。所以現在我只是將C [i，j]存儲在一維數組中的相關對和粒子對中。這可能是最好的例子：說，我有兩個由粒子（3,5）和（3,10）組成的對。示意性地，我的變量則是這樣的（意重複計算）：for-loop + numpy的矢量化版本。其中

i = [3,3,5,10] #list of particles i that form a pair 
j = [5,10,3,3] #corresponding particles j (not used in the later example) 
cij = [c35,c310,-c35,-c310] #(with actual numbers in reality)

現在，它真的可以歸結爲尋找一種有效的量化的方式來改寫爲循環如下：

part_list=np.arange(N) 
for a in part_list: 
    cond = np.where(i == a) 
    ci[a] = np.sum(cij[cond])

其他解決方案我有想過，但想避免：

a）並行化的for循環：不可行b/c這是嵌入在一個外部循環，我想在某些點並行化。 b）在C中編寫for循環並將其封裝到Python中：對於這個（希望）相當簡單的問題來說似乎是一種矯枉過正的行爲。

來源

2014-12-03 user4319496

你可以將其表述爲稀疏矩陣乘法（它應該並行化）。然後將cij存儲在該稀疏矩陣中。不需要查找，這將防止矢量化。你看過[scipy.sparse]（http://docs.scipy.org/doc/scipy-0.14.0/reference/sparse.html）嗎？ – smci 2014-12-03 11:15:18

將信息添加到以前的評論。稀疏矩陣不存儲零，所以你不會有內存問題。 – 2014-12-03 11:25:05

從乍一看，這將工作。目前我會堅持使用numpy解決方案，但非常感謝您的建議！ – user4319496 2014-12-03 16:37:58

你可以用np.bincount得到你想要的。如果你的顆粒，其中按順序從0編號，你可以簡單地做：

ci = np.bincount(i, weights= cij)

要看看這是什麼一樣：

>>> i = [3, 3, 5, 10] 
>>> j = [5, 10, 3, 3] 
>>> cij = [0.1, 0.2, -0.1, -0.2] 
>>> np.bincount(i, weights= cij) 
array([ 0. , 0. , 0. , 0.3, 0. , -0.1, 0. , 0. , 0. , 0. , -0.2])

如果你不希望所有這些額外的零，你可以做什麼像：

>>> unq_i, inv_i = np.unique(i, return_inverse=True) 
>>> unq_ci = np.bincount(inv_i, weights=cij) 
>>> unq_i 
array([ 3, 5, 10]) 
>>> unq_ci 
array([ 0.3, -0.1, -0.2])

而且你可以在以後通過執行類似分配這些獨特的價值觀：

ci[unq_i] = unq_ci

來源

2014-12-03 13:01:53 Jaime

太棒了，零點的第一部分（幾乎）正是我所需要的！我只需要通過在最後一個非零元素之後附加零來確保ci的維數爲N（但這可能不是我的問題中顯而易見的）。非常感謝！ – user4319496 2014-12-03 16:34:43

for-loop + numpy的矢量化版本。其中

回答

相關問題