Question: I have a sphere of radius $R\_{s}$, and I would like to pick random points in its volume with uniform probability. How can I do so while preventing any sort of clustering around poles or the center of the sphere?

---

Since I'm unable to answer my own question, here's another solution:

Using the strategy suggested by [Wolfram MathWorld](http://mathworld.wolfram.com/SpherePointPicking.html) for picking points on the surface of a sphere: Let $\theta$ be randomly distributed real numbers over the interval $[0,2\pi]$, let $\phi=\arccos(2v−1)$ where $v$ is a random real number over the interval $[0,1]$, and let $r=R\_s (\mathrm{rand}(0,1))^\frac13$. Converting from spherical coordinates, a random point in $(x,y,z)$ inside the sphere would therefore be: $((r\cos(\theta)\sin(\phi)),(r\sin(\theta)\sin(\phi)),(r\cos(\phi)))$.

A quick test with a few thousand points in the unit sphere appears to show no clustering. However, I'd appreciate any feedback if someone sees a problem with this approach.

Answer: 
Nate and Kevin already answered the two I knew. Recalling [this](https://math.stackexchange.com/a/73358/312) and [this](http://mathcentral.uregina.ca/QQ/database/QQ.09.99/wilkie1.html), I think that another way to generate a uniform distribution over the sphere surface would be to generate a uniform distribution over the vertical cylinder enclosing the sphere, and then project horizontally. 

That is, generate $z \sim U[-R,R]$, $\theta \sim U[0,2\pi]$, and then $x=\sqrt{R^2-z^2} \cos(\theta)$, $y=\sqrt{R^2-z^2} \sin(\theta)$. This (if I'm not mistaken) gives a uniform distribution over the sphere surface. Then, apply Nate's recipe to get a uniform distribution over the sphere volume.

This method is a little simpler (and more efficient) than the accepted answer, though it's not generalizable to other dimensions.
Much simpler way would be to pair surfaces at different distances from the center.

The surface of a sphere of radius R has area 4 π R^2

So just pair up various surfaces that add up to the same area. 

For example:

Surface of sphere with r = 0 paired with surface of sphere with r = R to have a total area of 4π R^2

Surface of sphere with r = R/2 paired with surface of sphere with r = sqrt(3)R/2 to have a total area of 4π R^2

To get the random point, do the following steps:

1) For choosing theta (it can be uniformly chosen over its entire range, since all values of theta are equally likely):

Choose a random "theta" in range [0, 2π]

2) For choosing phi (Use same strategy as for areas, but this time it will apply by pairing circumference of circle for a phi value 2π R cos(phi) with its complement to produce total length of 2π R):

Choose a random "a" in range [0, π/2]

Choose a random "b" in range [0, 1]

Choose a random "c" to be either 1 or -1 (up or down)

```
if b <= cos(a) then
    phi = c a
else
    phi = c inv_cos(1 - cos(a))

```

3) For choosing r (Pair up every possible radius value with its complement, in such a way that the sum of the surface areas of the spheres will be equal to 4π R^2):

Choose a random "d" in range [0, 1]

Choose a random "e" in range [0, 1]

One r value is "d R" the complementary r value is "sqrt(1 - d^2) R".

Total area = 4π (d R)^2 + 4π (sqrt(1 - d^2) R)^2

Total area = 4π (d^2 R^2 + (1 - d^2) R^2)

Total area = 4π (d^2 R^2 + R^2 - d^2 R^2)

Total area = 4π R^2 (as we wanted)

```
if e <= d^2 then
    r = d R
else
    r = sqrt(1 - d^2) R

```

Your random point uniformly chosen inside the sphere in radial coordinates is:

(r, theta, phi)