问如何在Python中用蒙特卡罗方法计算10维球体的体积？

EN

你在日常工作中遇到了多个问题。

您收到的错误消息来自您的线路

r = 2.0*np.random.uniform(0,1,iter) - 1.0

函数调用np.random.uniform(0,1,iter)不会创建单个随机数。相反，与大多数numpy函数一样，它返回一个数组--在本例中是您声明的长度的向量(在本例中为iter)。因此，r也是一个数组，因为它使用这个数组，而total0也是一个数组，原因与此相同。

然后，您将尝试评估total0 < 1。但是左边是一个数组，右边是一个标量，所以numpy不喜欢这种比较。我可以更详细地了解错误消息的含义，但这是基本的想法。

解决方案是将r视为一个向量--实际上，将其视为具有所需边长2的球体中的随机点。您犯了一个拼写错误，并且使用iter而不是dim作为随机向量的大小，但是如果您进行了这样的更改，您将获得您想要的点。您可以使用numpy快速获得其长度，并查看该点是否位于以原点为中心的半径球面内。

以下是一些已更正的代码。我删除了一个循环--在这个循环中，您试图构建一个大小正确的向量，但现在我们让numpy一次性构建它。我还用更具描述性的名称替换了您的变量名称，并做了一些其他更改。

import numpy as np

def nSphereVolume(dim, iterations):
    count_in_sphere = 0

    for count_loops in range(iterations):
        point = np.random.uniform(-1.0, 1.0, dim)
        distance = np.linalg.norm(point)
        if distance < 1.0:
            count_in_sphere += 1

    return np.power(2.0, dim) * (count_in_sphere / iterations)

print(nSphereVolume(10, 100000))

请注意，仅仅1000次蒙特卡洛迭代就会给出非常差的精度。您将需要更多的迭代才能获得有用的答案，因此我将重复次数改为100,000。该例程现在速度较慢，但在2.5附近给出了更一致的答案。这与的theoretical answer非常一致。

np.pi**(10 // 2) / math.factorial(10 // 2)

它的计算结果为2.550164039877345。

(这是对注释的响应，用于解释返回值np.power(2.0, dim) * (count_in_sphere / iterations。)

此例程生成随机点，其中每个坐标都在间隔[-1.0, 1.0)中。这意味着这些点均匀分布在维度dim的超立方体中。超立方体的体积是它的边的乘积。每条边都有length 2.0，因此我们可以使用2.0 power dim或np.power(2.0, dim)来计算超立方体的体积。

我们生成了iterations点，其中的count_in_sphere在半径为1的超球面中以原点为中心。因此，超立方体中也在超球面中的点的分数或比例是count_in_sphere / iterations。这些点均匀地分布在超立方体中，因此我们可以估计超球体的体积占超立方体体积的比例与这些集合中随机点的比例相同。从高中的比例来看，我们得到了

volume_of_hypersphere / volume_of_hypercube = points_in_hypersphere / points_in_hypercube

意识到方程式仅仅是近似值。将方程式的两边乘以volume_of_hypercube，我们得到

volume_of_hypersphere = volume_of_hypercube * points_in_hypersphere / points_in_hypercube

代入，我们得到

volume_of_hypersphere = np.power(2.0, dim) * count_in_sphere / iterations

它是函数返回的值。同样，它只是一个近似值，但概率论告诉我们，我们生成的随机点越多，近似值越好。