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

Question

我想用蒙特卡洛积分法求出n球体的体积，并与解析结果进行比较。我想对 10^6 点和 10^9 点执行此操作，虽然它在某种程度上适用于 10^6 点（n = 2（圆）、n = 3（球）和 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 时应尽量避免循环。这个想法是你可以在一个矩阵中一次生成所有样本，然后向量化所有后续操作。

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

How to speed up random array generation in python

python-2.7

math

optimization

numerical-methods

montecarlo