語境:
我有 3 個 3D 陣列(“precursor陣列”),我正在使用逆距離加權方法進行上采樣。為此,我計算了一個 3Dweights陣列,我在陣列的每個點上的 for 回圈中使用該precursor陣列。
我的陣列的每個 2D 切片weights都用于計算一個partial陣列。一旦我生成了所有 28 個,它們被求和以給出一個最終host陣列。
我想并行化這個 for 回圈以減少我的計算時間。我試過這樣做,但我無法正確更新我的host陣列。
問題:
我怎樣才能并行化我的主要功能(我的代碼的最后一部分)?
編輯:或者有沒有一種方法可以“切片”我i的 for 回圈(例如,一個核心在 i = 0 到 5 之間運行,一個核心在 i = 6 到 9 上運行)?
概括:
- 3 個
precursor陣列(溫度、降水、雪):10x4x7(10是時間維度) - 1 個
weight陣列(寬):28x1101x2101 - 28x3
partial陣列:1101x2101 - 3 個
host陣列(temp、prec、eprec):1101x2101
這是我的代碼(可運行,因為它在該MAIN ALGORITHM PARALLEL部分之外,請參閱最后的MAIN ALGORITHM NOT PARALLEL部分以了解我的代碼的非并行版本):
import numpy as np
import multiprocessing as mp
import time
#%% ------ Create data ------ ###
temperatures = np.random.rand(10,4,7)*100
precipitation = np.random.rand(10,4,7)
snow = np.random.rand(10,4,7)
# Array of altitudes to "adjust" the temperatures
alt = np.random.rand(4,7)*1000
#%% ------ Functions to run in parallel ------ ###
# This function upsamples the precursor arrays and creates the partial arrays
def interpolator(i, k, mx, my):
T = ((temperatures[i,mx,my]-272.15) (-alt[mx, my] * -6/1000)) * w[k,:,:]
P = (precipitation[i,mx,my])*w[k,:,:]
S = (snow[i,mx,my])*w[k,:,:]
return(T, P, S)
# We add each partial array to each other to create the host array
def get_results(results):
global temp, prec, Eprec
temp = results[0]
prec = results[1]
Eprec = results[2]
#%% ------ IDW Interpolation ------ ###
# We create a weight matrix that we use to upsample our temperatures, precipitations and snow matrices
# This part is not that important, it works well as it is
MX,MY = np.shape(temperatures[0])
N = 300
T = np.zeros([N*MX 1, N*MY 1])
# create NxM inverse distance weight matrices based on Gaussian interpolation
x = np.arange(0,N*MX 1)
y = np.arange(0,N*MY 1)
X,Y = np.meshgrid(x,y)
k = 0
w = np.zeros([MX*MY,N*MX 1,N*MY 1])
for mx in range(MX):
for my in range(MY):
# Gaussian
add_point = np.exp(-((mx*N-X.T)**2 (my*N-Y.T)**2)/N**2)
w[k,:,:] = add_point
k = 1
sum_weights = np.sum(w, axis=0)
for k in range(MX*MY):
w[k,:,:] /= sum_weights
#%% ------ MAIN ALGORITHM PARALLEL ------ ###
if __name__ == '__main__':
# Create an empty array to use as a template
dummy = np.zeros((w.shape[1], w.shape[2]))
# Start a timer
ts = time.time()
# Iterate over the time dimension
for i in range(temperatures.shape[0]):
# Initialize the host arrays
temp = dummy.copy()
prec = dummy.copy()
Eprec = dummy.copy()
# Create the pool based on my amount of cores
pool = mp.Pool(mp.cpu_count())
# Loop through every weight slice, for every cell of the temperatures, precipitations and snow arrays
for k in range(0,w.shape[0]):
for mx in range(MX):
for my in range(MY):
# Upsample the temperatures, precipitations and snow arrays by adding the contribution of each weight slice
pool.apply_async(interpolator, args = (i, k, mx, my), callback = get_results)
pool.close()
pool.join()
# Print the time spent on the loop
print("Time spent: ", time.time()-ts)
#%% ------ MAIN ALGORITHM NOT PARALLEL ------ ###
if __name__ == '__main__':
# Create an empty array to use as a template
dummy = np.zeros((w.shape[1], w.shape[2]))
ts = time.time()
for i in range(temperatures.shape[0]):
# Create empty host arrays
temp = dummy.copy()
prec = dummy.copy()
Eprec = dummy.copy()
k = 0
for mx in range(MX):
for my in range(MY):
get_results(interpolator(i, k, mx, my))
k = 1
print("Time spent:", time.time()-ts)
uj5u.com熱心網友回復:
多處理的問題在于它創建了許多新行程,它們在 main 之前(即 before )執行代碼。if __name__ == '__main__'這會導致初始化非常緩慢(因為所有行程都這樣做)并且大量 RAM 被無用使用。您當然應該在 main 中移動所有內容,或者如果可能的話在函式中移動(這通常會導致更快的執行,并且無論如何都是一種很好的軟體工程實踐,尤其是對于并行代碼)。即使這樣,多處理還有另一個大問題:行程間通信很慢。一種解決方案是使用通過使用Numba 或 Cython實作的多執行緒方法(您可以使用它們禁用 GIL,而不是基本的 CPython 執行緒)。事實上,它們通常比多處理更易于使用。但是,您應該更加小心,因為并行訪問不受保護,并且資料競爭可能出現在偽造的并行代碼中。
在您的情況下,計算主要是memory-bound。這意味著多處理毫無用處。事實上,并行性在這里幾乎沒有用處,除非您在具有高吞吐量的計算服務器上運行此代碼。事實上,記憶體是一種共享資源,使用更多的計算核心并沒有多大幫助,因為一個核心幾乎可以使普通 PC 上的記憶體帶寬飽和(而計算服務器上需要很少的核心)。
加速記憶體系結代碼的關鍵是避免創建臨時陣列并使用快取友好的演算法。在您的情況下T,P和S被填充只是為了稍后讀取,以便更新temp,prec和Eprec陣列。這個臨時步驟在這里非常昂貴且必要(尤其是填充陣列)。洗掉它會增加算術強度,從而導致代碼在順序上肯定會更快,并且可以更好地在多個內核上擴展。我的機器上就是這種情況。
這是使用 Numba 來并行化代碼的代碼示例:
import numba as nb
# get_results interpolator
@nb.njit('void(float64[:,::1], float64[:,::1], float64[:,::1], float64[:,:,::1], int_, int_, int_, int_)', parallel=True)
def interpolate_and_get_results(temp, prec, Eprec, w, i, k, mx, my):
factor1 = ((temperatures[i,mx,my]-272.15) (-alt[mx, my] * -6/1000))
factor2 = precipitation[i,mx,my]
factor3 = snow[i,mx,my]
for i in nb.prange(w.shape[1]):
for j in range(w.shape[2]):
val = w[k, i, j]
temp[i, j] = factor1 * val
prec[i, j] = factor2 * val
Eprec[i, j] = factor3 * val
# Example of usage:
interpolate_and_get_results(temp, prec, Eprec, w, i, k, mx, my)
請注意,其中的字串nb.njit稱為簽名,并為 JIT 指定型別,以便它可以急切地編譯它。
這段代碼在我的 6 核機器上快了 4.6 倍(而在沒有合并的情況下它幾乎快不了get_results和interpolator)。事實上,它在順序上快了 3.8 倍,所以執行緒并沒有多大幫助,因為計算仍然是記憶體系結的。實際上,與記憶體讀/寫相比,乘加的成本可以忽略不計。
轉載請註明出處,本文鏈接:https://www.uj5u.com/houduan/481501.html
上一篇:撰寫一個for回圈,創建5個隨機生成的長度分別為100,200,500,800和1000的資料幀,并列印每個資料幀
