首先我们生成一个cdf[x]

一点到圆心的距离小于x的概率。假设圆的半径为R。

显然，如果x = 0，那么cdf[0] = 0

显然，如果x是R，则cdf[R] = 1

显然，如果x = r，则cdf[r] = (r^2)/(r^2)

这是因为圆上的每个“小区域”都有相同的被选中的概率，所以概率与问题区域成比例。距离圆心x的面积是r^2

所以cdf[x] = x^2/R^2因为两者相互抵消了

我们有cdf[x]=x^2/R^2其中x从0到R

我们解出x

R^2 cdf[x] = x^2

x = R Sqrt[ cdf[x] ]

现在我们可以用一个从0到1的随机数来替换cdf

x = R Sqrt[ RandomReal[{0,1}] ]

最后

r = R Sqrt[  RandomReal[{0,1}] ];
theta = 360 deg * RandomReal[{0,1}];
{r,theta}

我们得到极坐标 {0.601168 R, 311.915°}

半径和“靠近”该半径的点的数量之间存在线性关系，因此他需要使用半径分布，这也使得半径r附近的数据点的数量与r成正比。

如何在半径为R的圆内随机生成一个点:

r = R * sqrt(random())
theta = random() * 2 * PI

(假设random()均匀地给出0到1之间的值)

如果你想把它转换成笛卡尔坐标，你可以做到

x = centerX + r * cos(theta)
y = centerY + r * sin(theta)

为什么sqrt(随机())?

让我们看看sqrt(random())之前的数学运算。为简单起见，假设我们是在单位圆上工作，即R = 1。

点与点之间的平均距离应该是相同的，不管我们看的距离中心有多远。这意味着，例如，观察一个周长为2的圆的周长，我们应该找到的点的数量是周长为1的圆周长上点的数量的两倍。

由于圆的周长(2πr)随r线性增长，因此随机点的数量应该随r线性增长。换句话说，期望的概率密度函数(PDF)线性增长。由于PDF的面积应该等于1，最大半径是1，我们有

所以我们知道随机值的理想密度应该是什么样的。 现在:当我们只有一个0到1之间的均匀随机值时，我们如何生成这样一个随机值?

我们用了一个叫做反变换采样的技巧

从PDF中创建累积分布函数(CDF) 沿着y = x镜像 将得到的函数应用于0到1之间的统一值。

听起来复杂吗?让我插入一段带有小侧轨的引语来传达直觉:

Suppose we want to generate a random point with the following distribution:                  That is 1/5 of the points uniformly between 1 and 2, and 4/5 of the points uniformly between 2 and 3. The CDF is, as the name suggests, the cumulative version of the PDF. Intuitively: While PDF(x) describes the number of random values at x, CDF(x) describes the number of random values less than x. In this case the CDF would look like:                  To see how this is useful, imagine that we shoot bullets from left to right at uniformly distributed heights. As the bullets hit the line, they drop down to the ground:                  See how the density of the bullets on the ground correspond to our desired distribution! We're almost there! The problem is that for this function, the y axis is the output and the x axis is the input. We can only "shoot bullets from the ground straight up"! We need the inverse function! This is why we mirror the whole thing; x becomes y and y becomes x:                  We call this CDF-1. To get values according to the desired distribution, we use CDF-1(random()).

所以，回到生成随机半径值，其中PDF等于2x。

步骤1:创建CDF: 由于我们处理的是实数，CDF表示为PDF的积分。

CDF（x） = ∫ 2x = x2

步骤2:沿y = x镜像CDF:

从数学上讲，这可以归结为交换x和y并求解y:

CDF: y = x2 交换:x = y2 解:y =√x CDF-1: y =√x

步骤3:将得到的函数应用于0到1之间的统一值

CDF-1(random()) =√random()

这就是我们要推导的:-)

这样想。如果你有一个矩形，其中一个轴是半径，一个是角，你取这个矩形内半径为0的点。它们都离原点很近(在圆上很近)然而，半径R附近的点，它们都落在圆的边缘附近(也就是说，彼此相距很远)。

这可能会让你知道为什么你会有这种行为。

在这个链接上导出的因子告诉你，矩形中有多少对应的区域需要调整，以便在映射到圆后不依赖于半径。

编辑:所以他在你分享的链接中写道，“通过计算累积分布的倒数，这很容易做到，我们得到r:”。

这里的基本前提是，通过将均匀分布映射为期望概率密度函数的累积分布函数的逆函数，可以从均匀分布创建一个具有期望分布的变量。为什么?现在把它当做理所当然，但这是事实。

这是我对数学的一些直观解释。密度函数f(r)关于r必须与r本身成比例。理解这个事实是任何微积分基础书的一部分。请参阅有关极区元素的部分。其他一些海报也提到了这一点。

我们记作f(r) = C*r;

这就是大部分的工作。现在，由于f(r)应该是一个概率密度，你可以很容易地看到，通过对f(r)在区间(0,r)上积分，你可以得到C = 2/ r ^2(这是给读者的练习)。

因此，f(r) = 2*r/ r ^2

好，这就是如何得到链接中的公式。

然后，最后一部分是从(0,1)中的均匀随机变量u你必须从这个期望密度f(r)映射到累积分布函数的逆函数。要理解为什么会这样，你可能需要找到像Papoulis这样的高级概率文本(或者自己推导)。

对f(r)积分得到f(r) = r^2/ r^2

为了求出它的反函数你设u = r^2/ r^2然后解出r，得到r = r *√(u)

直观上讲，u = 0映射到r = 0。同样，u = 1应该映射到r = r。同样，它通过平方根函数，这是有意义的，与链接匹配。

让我们像阿基米德那样处理这个问题。

我们如何在三角形ABC中均匀地生成一个点，其中|AB|=|BC|?让我们把它扩展到平行四边形ABCD。在ABCD中很容易均匀地生成点。我们均匀地选择AB上的X点和BC上的Y点并选择Z使XBYZ是一个平行四边形。为了在原始三角形中得到一个均匀选择的点，我们只需将ADC中出现的任何点沿AC折叠回ABC。

现在考虑一个圆。在极限情况下，我们可以把它想象成无穷多个等腰三角形ABC, B在原点，A和C在周长上，彼此逐渐接近。我们可以从这些三角形中选择一个角。所以我们现在需要通过在ABC条上选择一点来生成到中心的距离。同样，延伸到ABCD, D现在是圆中心半径的两倍。

使用上述方法可以很容易地在ABCD中选择一个随机点。在AB上随机选一个点，在BC上随机选一个点。Ie。在[0,R]上取一对随机数字x和y，给出离中心的距离。三角形是一条细条AB和BC本质上是平行的。所以Z点到原点的距离是x+y。如果x+y >r我们向下折叠。

这是R=1的完整算法。我希望你同意这很简单。它使用三角函数，但您可以保证它需要多长时间，以及需要多少次random()调用，这与拒绝抽样不同。

t = 2*pi*random()
u = random()+random()
r = if u>1 then 2-u else u
[r*cos(t), r*sin(t)]

这里是Mathematica。

f[] := Block[{u, t, r},
  u = Random[] + Random[];
  t = Random[] 2 Pi;
  r = If[u > 1, 2 - u, u];
  {r Cos[t], r Sin[t]}
]

ListPlot[Table[f[], {10000}], AspectRatio -> Automatic]

这可能会帮助那些对选择速度算法感兴趣的人;最快的方法是(可能?)拒绝抽样。

只需在单位正方形内生成一个点，并拒绝它，直到它在圆内。如(伪代码),

def sample(r=1):
    while True:
        x = random(-1, 1)
        y = random(-1, 1)
        if x*x + y*y <= 1:
            return (x, y) * r

虽然有时它可能运行不止一次或两次(而且它不是常量时间，也不适合并行执行)，但它要快得多，因为它不使用像sin或cos这样复杂的公式。