Distributing points evenly in a circle

As pointed out earlier, in fact its question is mathematics.

The source of the problem

The problem is that you are taking a uniform distribution of Rays, i.e., the probability that a point is at any distance Rfrom the center does not depend on R. The intuitive way to see that this does not produce a uniform distribution in the area of the circle is as follows: if the probability of a point being at a distance r from the center it is the same for any R, you expect to find on average the same number of points at a distance R from the Center, for all R.

Furthermore, since the probability is independent of the angle, the points at a distance r must be equally distributed over a radius circumference r. It turns out that for R large the circumferences are getting larger, and since you are distributing the same number of points over them, it is wait for the points to be more spaced from each other!

In order for the distribution to be uniform in terms of area, you need to place more points over the larger circumferences, or in other words, the probability of placing a certain point at a distance R should grow with R. But grow how? Since the size of the circumferences grows linearly with the radius (=2πr ), this probability is expected to grow as well linearly with r . To see this more rigorously we need more math.

Ray distribution

It occurs that the area elements in spherical coordinates depend on the radius, becoming larger for R larger:

dA = r dr dθ,

And for the distribution to be indeed uniform, you want the probability that there is a point in a given area element to be constant for any element. This probability is, by the definition of probability density, P (r, θ) dA , that is, P (r, θ) R dr dθ.

Integrating this for R between 0 and R and θ between 0 and 2π you have the total probability of the point falling somewhere in the circle, which should be 1. It's not hard to see that this fixes

P(r,θ) = 1 / (πR²),

What is obvious: in a homogeneous distribution over a given total area, the probability associated with an area element must be the ratio between the area of the element and the total area. See that this is exactly what gives the product P (r, θ) da :

P(r,θ) dA = r dr dθ / (πR²).

Inhomogeneous variable

Well, so now you can interpret this result a little differently. Instead of the variables r and θ , we changed to two random variables x (ranging from 0 to R) and y (ranging from 0 to 2π) that now have a probability density

P(x,y) = x / (πR²)

So that the integral in addition to x going from 0 to R and y going from 0 to 2π let 1:

∫∫ P(x,y) dx dy = 1.

See that now your probability density P (x, y) does not actually depend on y but is proportional to x ! This is the result I mention at the beginning: you are actually looking for a random variable x between 0 and R but whose distribution is not uniform, actually being proportional to x.

To understand square root in rand(1, N) we still need to advance a little further.

A small technicality: the variables x and y are independent, since we can rewrite P (x,y) as a product

P(x,y) = f(x)g(y).

It is easy to see that the proper choice is

g(y) = 1/2π

So that g g (y) dy = 1 for the integral between 0 and 2π, and

f(x) = 2x/R²

Such that f f (X) dx = 1 for x between 0 and R . Notice that f (x) g (y) actually retrieves P (x, y) .

Since g (y) does not depend on y , we see that the angle is actually homogeneously distributed(which explains why you can use 2*pi*rand (1,N) inside the exponential).

The challenge then now is to fetch a variable x distributed according to the probability density f(x) above.

Transformation method

A method for doing this from evenly distributed variables it is the so-called method of the transformation .

For variables x distributed from 0 to R according to a distribution f (X) , the fraction of numbers falling between 0 and x , which we will call F (X) , is given by

F(x) = ∫f(x') dx',

With the integral performed between 0 and x .

It turns out that for numbers evenly distributed between 0 and 1, the fraction of numbers that fall between 0 and G (for a given G between 0 and 1) is exactly g . See: the fraction of numbers between 0 and 1/2 is 1/2, between 0 and 1/3 is 1/3, and so on.

Then the fraction of numbers between 0 and a given F(X) for numbers evenly distributed between 0 and 1 is f(x), which is exactly the fraction of numbers between 0 and x in the distribution f(x) defined between 0 and R.

So what we want is to map the range between 0 and F (x) of numbers r distributed homogeneously in the range from 0 to x of the new distribution. And this is quite easy to do: the extremes of the intervals must match! Thus the number x of the non-uniform distribution should be obtained when the number

r = F(x)

Is obtained in the uniform distribution. Conclusion: given the evenly distributed numbers r, just invert this relationship - finding the inverse function of F(X) - to get the x. Come on apply this idea to our F(x) above.

We have f (x) = 2x / R2, so that, remembering that the integral is performed between 0 and x,

F(x) = ∫f(x') dx' = ∫2x'/R² dx' = x²/R².

Reversing this relationship is easy:

r = F(x) => r = x²/R² => x = R √r

What (finally!) shows that a variable x distributed between 0 and R with probability density f (x) = 2x / R2 can be obtained simply by multiplying R by the square root of a variable r homogeneously distributed between 0 and 1. This is exactly what you are doing with

R * sqrt(rand(1,N))

In your second line of code. This 'empirical' solution therefore actually provides exactly homogeneously distributed numbers over a circle of Radius R .