我有以下代码,涉及总结一些嵌套的for循环。
有什么方法可以加快执行这段代码的速度?我感兴趣的不仅仅是一种加速它的方法,我真的很想看到一系列方法,例如: 'Pure Python',Numpy,Scipy,Cython等。
这样,对于类似的,但(更多)更复杂的代码我必须编写,我可以选择一个加速选项,可以很好地权衡执行速度与实现的复杂性。还有什么可以让我免于编写C ++代码而导致我失去生存意愿。
def f(a,b,c,d):
return a+b+c+d
x = [0.04691008, 0.23076534, 0.5, 0.76923466, 0.95308992]
w = [0.11846344, 0.23931434, 0.28444444, 0.23931434, 0.11846344]
numQuadNodes = 5
def tensorGauss(func):
sum = 0;
for i in range(0,numQuadNodes):
for j in range(0,numQuadNodes):
for k in range(0,numQuadNodes):
for l in range(0,numQuadNodes):
sum += w[i]*w[j]*w[k]*w[l]*func(x[l],x[k],x[j],x[i])
return sum
print(tensorGauss(f))
编辑 - 更真实的代码正如你所看到的,tensorGauss已经比nquad快得多(在我的机器上0.07秒对20.86秒),但我真的想有一些方法让tensorGauss再次更快,因为我必须计算一吨张量高斯评价!
import numpy as np
import numpy.linalg as LA
from scipy.integrate import nquad
import time
##################################################
# Triangle vertices
##################################################
v_a_1 = np.array([[4,0,0]]).T
v_a_2 = np.array([[5,1,0]]).T
v_a_3 = np.array([[4,2,0]]).T
v_b_1 = np.array([[4,0,0]]).T
v_b_2 = np.array([[5,-1,0]]).T
v_b_3 = np.array([[4,-2,0]]).T
##################################################
# g_tau
##################################################
def g_tau():
J_tau = v_a_2-v_a_1
J_tau = np.append(J_tau, v_a_3-v_a_2,axis=1)
G = np.dot(J_tau.T,J_tau)
return np.sqrt(LA.det(G))
##################################################
# g_t
##################################################
def g_t():
J_t = v_b_2-v_b_1
J_t = np.append(J_t, v_b_3-v_b_2,axis=1)
G = np.dot(J_t.T,J_t)
return np.sqrt(LA.det(G))
##################################################
# chi_tau
##################################################
def chi_tau(x):
return v_a_1 + (v_a_2-v_a_1)*x[0] + (v_a_3-v_a_2)*x[1]
##################################################
# chi_t
##################################################
def chi_t(y):
return v_b_1 + (v_b_2-v_b_1)*y[0] + (v_b_3-v_b_2)*y[1]
##################################################
# k_
##################################################
def k_(x,y):
return LA.norm(x+y)
##################################################
# k
##################################################
def k(x,y):
return k_(chi_tau(x),chi_t(y))*g_tau()*g_t()
start=time.time()
##################################################
# tensorGauss
##################################################
x = [0.04691008, 0.23076534, 0.5, 0.76923466, 0.95308992]
w = [0.11846344, 0.23931434, 0.28444444, 0.23931434, 0.11846344]
numQuadNodes = 5
def f(z, y, x, w):
a_1_1 = z;
a_1_2 = z * w;
a_2_1 = z * x;
a_2_2 = z * x * y;
a_1 = np.array([a_1_1,a_1_2]).T
a_2 = np.array([a_2_1,a_2_2]).T
res = k(a_1,a_2)
a_1_1 = z * x;
a_1_2 = z * x * y;
a_2_1 = z;
a_2_2 = z * w;
a_1 = np.array([a_1_1,a_1_2]).T
a_2 = np.array([a_2_1,a_2_2]).T
res += k(a_1,a_2)
a_1_1 = z * y;
a_1_2 = z * w;
a_2_1 = z * x;
a_2_2 = z;
a_1 = np.array([a_1_1,a_1_2]).T
a_2 = np.array([a_2_1,a_2_2]).T
res += k(a_1,a_2)
return res
def tensorGauss(func):
sum = 0;
for i in range(0,numQuadNodes):
for j in range(0,numQuadNodes):
for k in range(0,numQuadNodes):
for l in range(0,numQuadNodes):
sum += w[i]*w[j]*w[k]*w[l]*func(x[l],x[k],x[j],x[i])
return sum
start=time.time()
tensorGauss_res = tensorGauss(f)
end=time.time()
tensorGauss_time = end-start
start=time.time()
[nquad_res, err] = nquad(f, [[0,1], [0,1], [0,1], [0,1]])
end=time.time()
nquad_time = end-start
print(f'tensor-gauss: {tensorGauss_res}')
print(f'nquad: {nquad_res}')
print('\n')
print(f'tensor-gauss time: {tensorGauss_time}')
print(f'nquad time: {nquad_time}')
我重写你的张量Gauss()函数如下:
def tensorGauss(func):
w_gen = np.meshgrid(w,w,w,w,indexing='ij')
x_gen = np.meshgrid(x,x,x,x,indexing='ij')
sum = np.sum(w_gen[0] * w_gen[1] * w_gen[2] * w_gen[3] *
f(x_gen[3], x_gen[2], x_gen[1], x_gen[0]))
return sum
它打印了2.0
的结果,而不是你发布的简化的1.9999999999999971
打印的tensorGauss
的值(使用简化的f()
函数的值)。
但是,一些免责声明:
g_tau
这样的名字等等)。我这样说是因为,这个解决方案假设你的f
是一个矢量化函数,如果传递数组而不是标量,它将在元素方面起作用。我可以看到这个假设适用于你的虚拟f()
,但是,我不知道它是否也适用于你真正的f()
f
以及您的实际数据大小进行检查