# Generate points in circle uniformly

## Introduction

This chapter’s problem concerns uniformly generating a (potentially large) number of random points in a circle of a certain radius. Despite its simplicity the problem poses some unexpected challenges. We will discuss the best approach to this problem as well as one solution that many candidates provide which, whilst initially appearing correct actually fails in one crucial aspect (spoiler: it does not distribute points uniformly).

## Problem statement

Write a function that, given a circle of radius $$r$$ and centered at $$(x,y)$$ where $$r,x,y \in \mathcal{R}$$ returns a uniformly distributed point in the circle.

## Clarification Questions

What exactly does it mean for the point to be uniformly distributed?

It means that every point of the circle has the same probability of being picked/generated by the function.

## Discussion

Before discussing solutions it is worth mentioning that the fact that the circle is centered at $$(x,y)$$ makes very little difference and we can continue our discussion as if it were centered at $$(0,0)$$. This is the case because all the points we generate can then be translated to $$(x,y)$$ by simply adding $$x$$ and $$y$$ to the $$x$$-coordinate and $$y$$-coordinate of the generated point.

### Polar Coordinates - The wrong approach

Let’s start by discussing an intuituve, but ultimately incorrect, approach. One might think that in order to pick a point in the circle it is sufficient to

1. Pick a random angle $$\theta \in [0, 2\pi[$$

2. Pick a random radius $$\overline{r} \in [0,r]$$

3. Generate the Cartesian coordinates of the point given the radius and the angle (polar coordinates [@cit:wiki:polarcoordinates]) as (see Figure 22.1): $\begin{gathered} x=\overline{r}\sin(\theta) \\ y=\overline{r}\cos(\theta) \end{gathered}$

[fig:random_points_in_cirle:polar_coordinates] Unfortunately, despite its appealing simplicity, this approach is wrong as it fails to produce points that are distributed uniformly in the circle. Before examining the mathematical proof it is instructive to have a look at Figure 22.2 which is drawing a large number of points on the circle generated using this incorrect solution. As you can see, the points are not generated uniformly as their density is higher towards the center. The bottom line is, do not use this solution in an interview. A possible matlab implementation of this buggy approach is shown at [list:random_points_in_circle:buggy]

Listing 1: Non-uniform random point in a circle generation using Matlab

function [px, py] = buggy_random_point(radius, x,y)
theta = rand()*2*pi;
px = r * sin(theta);
py = r *cos(theta);
endfunction

[fig:random_points_in_cirle:buggy] ### Loop approach

A good way to ensure that the point density is uniform across on the surface of the circle is to pick a point randomly in an enclosing square and make sure that we discard all the points that lie outside the circle. In other words, we keep asking for a random point $$(p_x=\text{rand()}, p_y=\text{rand()})$$ in the enclosing square until the following is true: $$p_x^2 + p_y^2 \leq r$$. In this way we are guaranteed to generate uniformly distributed points because we pick those points from a set of points that are already uniformly distributed in a square, and we exclude those which are not inside the circle. This method is also known as the exclusion method. Figure 22.3 depicts a large number of points generated with this method.

The downside here is that we might need to generate a number of points in the square before getting lucky and picking one lying in the circle. We need to make on average $$\approx 1.2732$$ tries before getting a point in the circle. This number is the ratio between the are of enclosing square and the area of the enclosed circle i.e. $$\frac{(2r)^2}{\pi r^2} = \frac{4}{\pi}$$.

[fig:random_points_in_cirle:loop] A Matlab implementation of this approach is shown in Listing [list:random_points_in_circle:loop].

Listing 2: Random point in a circle generation using the exclusion method.

function [px, py, t] = random_point_loop(radius, x,y)
px = 100;
py = 100;
t=0;
while (px*px + py*py > 1)
signx = 2*randi([0,1])-1;
signy = 2*randi([0,1])-1;
t=t+1;
endwhile

endfunction

### Polar Coordinates - The right approach

In order for the points to be distributed uniformly it is necessary that the average distance between the points be the same regardless of how far they lie from the center of the circle. This means that, looking at the points generated on a circumference of radius $$2$$, there have to be twice as many points as the the number of points on a circumference of radius $$1$$. A circumference that is twice as long translates to needing twice as many points to maintain the same density. Another intuitive way to understand why simply picking a random angle and a random radius is not enough would be to think about having to distribute $$10$$ points at random on a circle of radius $$1$$ and $$2$$. It is clear that the circumference of radius $$2$$ would look emptier than the one with radius $$1$$ simply because there is more circumference to be filled but a constant amount of points.

The fundamental problem with the appraoch described in Section 22.3.1 is that the area of the circle is not uniformly covered. The random radius cuts through the area of the circle and this is the only parameter that affects how the points are going to be distributed across the full area of the circle. Therefore we should focus our attention how we can pick a better radius by making sure that larger radii are picked more often to accommodate for the larger area they define. In other words, we need to ensure that our random function for picking the radius takes the area of our circle into account.

Consider the area $$A$$ of a circle of radius $$r$$ i.e. $$A = \pi r^2$$. We can rearrange the formula so that $$r = \sqrt{\frac{A}{\pi}}$$. What this formula is really telling us is that the radius is proportional to $$\sqrt{A}$$. Now we have a way of choosing the radius that depends on the area of the circle. We can simply pick an area at random and then calculate the radius accordingly. This will make sure that the radius is picked taking into consideration the area of the circle. Figure 22.4 shows many points generated using this method. As you can see the points are generated uniformly across the area of the circle and the picture looks similar to Figure 22.3.

A C++ implementation of this method is shown in Listing [list:random_points_in_circle:sqrtcpp]. Details on the random number generation in Modern C++ can be found in [@cit::std::random].

Listing 3: C++ implementation of the function for generating a random point in a circle described in Section \ref{random_points_in_circle:sec:polar_sqrt}

auto generate_random_point_in_circle()
{
static std::uniform_real_distribution<double> dist_angle(0, 2 * M_PI);
const auto theta = dist_angle(rnd);
const auto x     = r * cos(theta);
const auto y     = r * sin(theta);
return std::make_pair{x + x_center, y + y_center};
}

[fig:random_points_in_cirle:polar_sqrt] A Matlab implementation of this approach is shown in Listing [list:random_points_in_circle:sqrt].

Listing 4: Random point in a circle generation using polar coordinates and the $\approx \sqrt{A}$ dependency of the radius on the area of the circle.

function [px, py] = random_sqrt_area(radius, x,y)
r = sqrt(area/pi);
theta = rand()*2*pi;
px = r * sin(theta);
py = r *cos(theta);
endfunction

### Conclusion

For both the viable methods for generating random points withing a circle which we have discussed the time and space complexity is constant although the one presented in Section 22.3.3 will probably have better performance when included in a hot path i.e. in a loop for the generation of many random points.

All the code used to generate the Figures in this chapter is shown in Listing [list:random_points_in_circle:drivercode].

Listing 5: Matlab driver code for the generation of all figures in Chapter \ref{ch:random_points_in_circle}

function draw_points(n)
clf(1);
% n is the total number of points
px = zeros(1,n);
py = zeros(1,n);

tries = 0;
for i =0:n
%[x,y] = buggy_random_point(1,0,0);
% [x,y,t] = random_point_loop(1,0,0);
[x,y] = random_sqrt_area(1,0,0);
% tries = tries + t;
px(i+1) = x;
py(i+1) = y;
endfor

average = tries/n

% Plot a circle.
angles = linspace(0, 2*pi, n);
xCenter = 0;
yCenter = 0;
cx = radius * cos(angles) + xCenter;
cy = radius * sin(angles) + yCenter;
% Plot circle.
plot(cx, cy, 'b-', 'LineWidth', 2);
% Plot center.
hold on;
plot(xCenter, yCenter, 'k+', 'LineWidth', 2, 'MarkerSize', 16);
grid on;
axis equal;
xlabel('X', 'FontSize', 16);
ylabel('Y', 'FontSize', 16);

% Plot random points.
plot(px, py, 'r.', 'MarkerSize', 1);

rectangle('Position',[-1 -1 2 2], 'LineWidth',3,  'EdgeColor' , [0 .5 .5])

endfunction