要在 3D 球中生成指定数量的随机点，一种可能性是使用拒绝采样来获得初始点分布，然后使用 Lloyd 松弛来改进采样。

具有 30、300 和 10000 点的示例

拒绝采样： 使用统一的随机数生成器生成 x,y,z 坐标，忽略球外的点，并继续直到获得所需的点数。

劳埃德松弛 是一种迭代方法，使点的分布更加规则。每次迭代的工作如下：计算点的 Voronoi 图，然后将每个点移动到其 Voronoi 单元的重心。有关更详细的说明，请参见 [1]。在我的 Geogram 软件库 [2] 中有一个实现。使用 Geogram，您需要使用四面体版本的球体。

在体积/表面上 还有一个算法版本适用于 3D 表面（球体或任意表面），称为约束质心 Voronoi Tesselation [3]，也在 Geogram [2] 中实现。

您也可以使用我的软件 Graphite [4]（围绕地理的图形用户界面）。在 Graphite 中对球体进行采样：

Scene->Create object->OGF::MeshGrob
Surface->Shapes->Create sphere
Points->Sample volume or Points->Sample surface

minitutorial：在 Graphite 中采样一个球

[1] https://en.wikipedia.org/wiki/Lloyd%27s_algorithm

[2] http://alice.loria.fr/software/geogram/doc/html/index.html

[3] 约束质心 Voronoi Tesselation，Du、Gunzburger 和 Ju，SIAM J. 计算，2003

[4] http://alice.loria.fr/software/graphite/doc/html/

有一个很好的文件：

• M. Deserno，“如何在球体表面生成均匀分布的点”

它基本上总结了在球体表面上生成点分布的两种最简单的方法：

1. 均匀随机分布（生成$N$随机的$z$和$\varphi$在$[-R, R]$和$[0,2\pi]$，分别计算$(x,y)$作为$(\sqrt{r^2-z^2}\cos{\varphi},\sqrt{r^2-z^2}\sin{\varphi})$. 然而，用足够小的$N$输出看起来不太统一。

2. 生成具有恒定间隔的纬度圈，点之间的距离几乎相同（参见上面的论文）。不会产生完全要求的点数，但会尝试使其足够接近。

我很快写了一个示例 Python 脚本：

#!/usr/bin/python
import numpy as np 

n = 5000
r = 1
generateRandom = True
generateFixed = True

if (generateRandom):
    print("Generating randomly %d points on a sphere centered at the origin" % (n))
    theta = np.random.uniform(0.0,2*np.pi,n)
    z = np.random.uniform(-1.0,1.0,n)
    for i in range (0,n):
        zp = z[i]
        xp = np.sqrt(r*r - zp*zp)* np.cos(theta[i])
        yp = np.sqrt(r*r - zp*zp)* np.sin(theta[i])
if (generateFixed):
    print("Generating fixed %d points on a sphere centered at the origin" % (n))
    alpha = 4.0*np.pi*r*r/n
    d = np.sqrt(alpha)
    m_nu = int(np.round(np.pi/d))
    d_nu = np.pi/m_nu
    d_phi = alpha/d_nu
    count = 0
    for m in range (0,m_nu):
        nu = np.pi*(m+0.5)/m_nu
        m_phi = int(np.round(2*np.pi*np.sin(nu)/d_phi))
        for n in range (0,m_phi):
            phi = 2*np.pi*n/m_phi
            xp = r*np.sin(nu)*np.cos(phi)
            yp = r*np.sin(nu)*np.sin(phi)
            zp = r*np.cos(nu)
            count = count +1