python向量化和c哪个快_在python中向量化6 for循环累积和
編輯3:
最終(我認為)版本,更清晰,更快速地融入max9111’s answer的創(chuàng)意.
import numpy as np
from numba import as nb
@nb.njit()
def func1_jit(a, b, c, d):
# Precompute
exp_min = 5 - (a + b + c + d)
exp_max = b
exp = 2. ** np.arange(exp_min, exp_max + 1)
fact_e = np.empty((a + b - 2))
fact_e[0] = 1
for ei in range(1, len(fact_e)):
fact_e[ei] = ei * fact_e[ei - 1]
# Loops
B = 0
for ai in range(0, a):
for bi in range(0, b):
for ci in range(0, c):
for di in range(0, d):
for ei in range(0, ai + bi):
for fi in range(0, ci + di):
B += exp[ei - fi - ai - ci - di + 1 - exp_min] * (ei * ei - 2 * (ei * fi) - 7 * di) * fact_e[ei]
return B
這已經(jīng)比以前的任何選項都快得多,但我們?nèi)匀粵]有利用多個CPU.一種方法是在函數(shù)本身內(nèi),例如并行化外循環(huán).這會在每次調(diào)用時增加一些開銷來創(chuàng)建線程,因此對于小輸入實際上有點慢,但對于更大的值應(yīng)該明顯更快:
import numpy as np
from numba import as nb
@nb.njit(parallel=True)
def func1_par(a, b, c, d):
# Precompute
exp_min = 5 - (a + b + c + d)
exp_max = b
exp = 2. ** np.arange(exp_min, exp_max + 1)
fact_e = np.empty((a + b - 2))
fact_e[0] = 1
for ei in range(1, len(fact_e)):
fact_e[ei] = ei * fact_e[ei - 1]
# Loops
B = np.empty((a,))
for ai in nb.prange(0, a):
Bi = 0
for bi in range(0, b):
for ci in range(0, c):
for di in range(0, d):
for ei in range(0, ai + bi):
for fi in range(0, ci + di):
Bi += exp[ei - fi - ai - ci - di + 1 - exp_min] * (ei * ei - 2 * (ei * fi) - 7 * di) * fact_e[ei]
B[ai] = Bi
return np.sum(B)
或者,如果您有許多要評估函數(shù)的點,也可以在該級別進行并行化.這里a_arr,b_arr,c_arr和d_arr是要評估函數(shù)的值的向量:
from numba import as nb
@nb.njit(parallel=True)
def func1_arr(a_arr, b_arr, c_arr, d_arr):
B_arr = np.empty((len(a_arr),))
for i in nb.prange(len(B_arr)):
B_arr[i] = func1_jit(a_arr[i], b_arr[i], c_arr[i], d_arr[i])
return B_arr
最佳配置取決于您的輸入,使用模式,硬件等,因此您可以根據(jù)您的情況組合不同的想法.
編輯2:
實際上,忘記我之前說過的話.最好的是JIT編譯算法,但是以更有效的方式.首先計算昂貴的部分(我采用指數(shù)和階乘),然后將其傳遞給編譯的循環(huán)函數(shù):
import numpy as np
from numba import njit
def func1(a, b, c, d):
exp_min = 5 - (a + b + c + d)
exp_max = b
exp = 2. ** np.arange(exp_min, exp_max + 1)
ee = np.arange(a + b - 2)
fact_e = scipy.special.factorial(ee)
return func1_inner(a, b, c, d, exp_min, exp, fact_e)
@njit()
def func1_inner(a, b, c, d, exp_min, exp, fact_e):
B = 0
for ai in range(0, a):
for bi in range(0, b):
for ci in range(0, c):
for di in range(0, d):
for ei in range(0, ai + bi):
for fi in range(0, ci + di):
B += exp[ei - fi - ai - ci - di + 1 - exp_min] * (ei * ei - 2 * (ei * fi) - 7 * di) * fact_e[ei]
return B
在我的實驗中,這是迄今為止最快的選項,并且只占用很少的額外內(nèi)存(只有預(yù)先計算的值,輸入上的大小為線性).
a, b, c, d = 4, 6, 3, 4
# The original function
%timeit func1_orig(a, b, c, d)
# 2.07 ms ± 33.7 μs per loop (mean ± std. dev. of 7 runs, 100 loops each)
# The grid-evaluated function
%timeit func1_grid(a, b, c, d)
# 256 μs ± 25 μs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
# The precompuation + JIT-compiled function
%timeit func1_jit(a, b, c, d)
# 19.6 μs ± 3.25 μs per loop (mean ± std. dev. of 7 runs, 100000 loops each)
那么總是可以選擇網(wǎng)格評估整個事情:
import numpy as np
import scipy.special
def func1(a, b, c, d):
ai, bi, ci, di, ei, fi = np.ogrid[:a, :b, :c, :d, :a + b - 2, :c + d - 2]
# Compute
B = (2.) ** (ei - fi - ai - ci - di + 1) * (ei ** 2 - 2 * (ei * fi) - 7 * di) * scipy.special.factorial(ei)
# Mask out of range elements for last two inner loops
m = (ei < ai + bi) & (fi < ci + di)
return np.sum(B * m)
print(func1(4, 6, 3, 4))
# 21769947.844726562
顯而易見,隨著參數(shù)的增加,其內(nèi)存成本將快速增長.代碼實際上執(zhí)行的計算比必要的多,因為兩個內(nèi)部循環(huán)具有不同的迭代次數(shù),因此(在此方法中)您必須使用最大的,然后刪除您不需要的.希望是矢量化將彌補這一點.一個小的IPython基準:
a, b, c, d = 4, 6, 3, 4
# func1_orig is the original loop-based version
%timeit func1_orig(a, b, c, d)
# 2.9 ms ± 110 μs per loop (mean ± std. dev. of 7 runs, 100 loops each)
# func1 here is the vectorized version
%timeit func1(a, b, c, d)
# 210 μs ± 6.34 μs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
編輯:
請注意,以前的方法也不是一個全有或全無的事情.您可以選擇僅對某些循環(huán)進行網(wǎng)格評估.例如,兩個最里面的循環(huán)可以像這樣矢量化:
def func1(a, b, c, d):
B = 0
e = np.arange(a + b - 2).reshape((-1, 1))
f = np.arange(c + d - 2)
for ai in range(0, a):
for bi in range(0, b):
ei = e[:ai + bi]
for ci in range(0, c):
for di in range(0, d):
fi = f[:ci + di]
B += np.sum((2.) ** (ei - fi - ai - ci - di + 1) * (ei ** 2 - 2 * (ei * fi) - 7 * di) * scipy.special.factorial(ei))
return B
這仍然有循環(huán),但它確實避免了額外的計算,并且內(nèi)存要求低得多.哪一個最好取決于我猜的輸入大小.在我的測試中,使用原始值(4,6,3,4),這甚至比原始函數(shù)慢;此外,對于這種情況,似乎在每個循環(huán)上為ei和fi創(chuàng)建新數(shù)組比在預(yù)先創(chuàng)建的循環(huán)上操作更快.但是,如果將輸入乘以4(14,24,12,16),那么這比原始(約x5)快得多,盡管仍然比完全矢量化的(約x3)慢.另一方面,我可以計算輸入的值,用十(40,60,30,40)來縮放這個(在~5分鐘內(nèi))而不是前一個因為內(nèi)存(我沒有測試如何)它需要與原始功能一起使用).使用@ numba.jit有點幫助,雖然不是很大(由于階乘函數(shù)不能使用nopython).您可以嘗試使用或多或少的循環(huán)向量化,具體取決于輸入的大小.
總結(jié)
以上是生活随笔為你收集整理的python向量化和c哪个快_在python中向量化6 for循环累积和的全部內(nèi)容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: @value 静态变量_Java中的变量
- 下一篇: python中for和while可以有e