可以使用什么方法来加速这段涉及嵌套 for 循环的 Python 代码？

Question

我有以下代码，其中涉及对多个嵌套 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 更快，因为我将不得不计算大量的 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}')

Answer 1

我重写了你的 tensorGauss() 函数如下：

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 的结果，与您发布的简化 tensorGauss 打印的 1.9999999999999971 的值相反（使用简化的 f() 函数).

但是，一些免责声明：

这是否甚至运行不会抛出错误，将取决于您的真实代码（名称如 g_tau 等的代码）。我这样说是因为，此解决方案假定您的 f 是一个矢量化函数，如果传递数组而不是标量，它将按元素工作。我可以看到这个假设适用于你的虚拟 f()，但是，我不知道它是否也适用于你的真实 f()
最好使用您的真实代码而不是虚拟代码 f 以及您的真实数据大小

可以使用什么方法来加速这段涉及嵌套 for 循环的 Python 代码？

What methods can be used to speed up this Python code involving nested for loops?

python

numpy

cython

scipy

numba