如何在python中加速随机数组的生成

Question

我想使用 Monte-Carlo 积分方法求出 n 球体的体积并将其与解析结果进行比较。我想为 10^6 点和 10^9 点执行此操作，虽然它在某种程度上适用于 10^6 点（对于 n = 2 (circle) 、 n = 3 (sphere) 和 n = 12 大约需要一分钟）， 10 ^ 9分非常慢。

MC 方法的简短解释：要找到半径 r = 1 的 n 球体的体积，我想象一个简单已知的体积（比如一个边长为 2*r 的 n 立方体），它完全包含 n -领域。然后我从 n 立方体中的均匀分布点采样并检查该点是否位于球体中。我计算了 n 球内所有如此生成的点。 V_sphere/V_cube的比值可以近似为N_inside/N_total，因此，V_sphere = V_cube * N_inside/N_total

这里是函数：

def hyp_sphere_mc(d,samples):    

    inside = 0                     #number of points inside sphere                                             
    sum = 0                        #sum of squared components

    for j in range(0,samples):        

        x2 = np.random.uniform(0,1,d)**2     #generate random point in d-dimensions                           
        sum = np.sum(x2)                     #sum its components

        if(sum < 1.0):                                              
            inside += 1                      #count points inside sphere

    V = ((2)**d)*(float(inside)/samples)     #V = V_cube * N_inside/N_total           

    V_true = float(math.pi**(float(d)/2))/math.gamma(float(d)/2 + 1) #analytical result 

    ERR = (float(abs(V_true-V))/V_true)*100        #relative Error

    print "Numerical:", V, "\t" , "Exact: ", V_true, "\t", "Error: ", ERR

我想问题是，对于每次迭代，我都会生成一个新的随机数组，这需要很多时间，尤其是如果我有 10^9 次迭代。有什么办法可以加快速度吗？

Answer 1

您可以将循环替换为以下内容：

inside = np.sum(np.sum(np.random.rand(samples,d)**2,1)<1)

使用 numpy 时，应尽量避免循环。这个想法是，您可以一次在矩阵中生成所有样本，然后对所有后续操作进行矢量化。