Trouble in the alley
You're walking alone at night through a dark alley in the middle of the night[1] when a hooded figure stops you. They tell you:
"Solve the following mathematical problem right now OR I WILL SHOOT YOU DEAD.". They give you a piece of paper with the problem statement.
Three points are chosen independently and uniformly at random from the interior of a unit circle. What is the probability that their circumcircle is entirely contained within the unit circle?
They have a gun so you have no choice but to solve it. Now, because you're super stressed you have no idea how to do this in some funky clever way by appealing to symmetry or whatever. You'll just have to do it the old fashioned way and painstakingly calculate the result. Luckily the hoodlum doesn't seem to care how long it takes you to arrive at the answer. He told you as much. He provides a notebook and pen.
The Cartesian start
You are dealing with 3 points, so 6 coordinates in total: 3 s, 3 s.
Now you know you're dealing with a circular problem (heh). So even though you're in grave danger (the ruffian keeps twirling the gun on his fingers as if to emphasize just how much danger) you've got the mental wherewithal to realise there's a potential variable transformation you can do to make it easier for yourself.
Let's express everything in terms of the circumcircle center , its radius as well as the three angles the 3 points make with the center : , , .
This gives you a mapping from to .
You have
Now, you have to calculate the Jacobian since this is a non-trivial transformation, and you'll have to adjust the integral appropriately. The hoodlum reminds you the columns represent the partials with respect to the given output variable, and the rows are with respect to the same input variable. You proceed.
The good news is that this is quite a regular Jacobian. The first two columns have a nice pattern, the third is straightforward (as long as you remember the and derivatives). The only slightly tricky columns are 4th, 5th and 6th - but they too have a cool pattern to them. You end up with:
Nice, you think to yourself. Now this is a 6x6 matrix. You really don't want to be calculating the determinant... But you need it. You begin trying to calculate the determinant but the masked figure stops you. They remind you of determinant properties that will make your job much easier:
- You can factor out a number from a whole row / column and multiply the determinant by the same factor instead.
- The determinant doesn't change if you add or subtract one row from another.
- For
Block Triangular Matrices
Where:
- is a square matrix
- is a block of zeros
- and are the other chunks...
Then you know:
First, you notice the common factor of in 3 columns. Factored out, they will give you a factor of altogether.
Now you have
Now you think about the block triangular matrix thing the guy mentioned. Surely this is relevant here. You notice that you could clear the first two columns and turn the top left 2x2 section into an identity matrix. Then which would make everything else much easier.
You continue. You want to clear the in the first column, rows 3 and 5. Also the in column 2, rows 4 and 6.
You're very happy. Since the top left is now an identity matrix with and you can skip the whole submatrix in the top right. You've now boiled this all down to calculating:
is a matrix, so not the easiest to calculate by hand. But it's tractable. Let's just remember to put that factor out front at the end.
You can now subtract row 1 from row 3, and subtract row 2 from row 4. You get
Now you can expand down colum 2 instead of the dense column 1.
You go on.
You start wondering if maybe you should have stuck with the original Cartesian product. Those sines look scary. You crack on nonetheless.
And then, remembering about the
Now, you know that and you also remember the formula for . Also, the hoodlum helpfully points out that he's happy to give you a refresher on the sine of a difference in case you need it.
You're ecstatic. You love sines and cosines now. You loudly point out that
and hence
Your thug gives you a supportive nod. You don't care about the sign of so you flip the difference inside of the sines and add absolute value for good measure to remove the stupid minus sign.
At this point you notice that this kind of makes sense because three collinear points won't make any kind of circumcirle. You carry on but you make a mental note to remember this problem if you ever need to interview someone for a software role. This is trivial now after all.
Okay, the Jacobian is all sorted. You move on to the probability calculation. In general, since we're picking points from the unit circle (area ) and you do it 3 times ("independently" as the villain points out excitedly!) you have a "volume" of to scale everything by.
Very generally
and less generally
where the condition you're after is since you know you only care about points whose circumcirle fits nicely inside the unit circle. That is: is the distance from the origin to the circumcircle, and you want to make sure that .
You hit this head on: you can split this into the angular part (the sine difference bit from the determinant of the Jacobian) and the geometric part.
For the geometric part, you realise that can go from to , and for any specific the center must be within a smaller circle of radius . "Area of the circle is pi times radius squared", the gangster points out. He's right.
So, the simplifies to (the comes from that pesky factor you promised not to forget. You did forget but the brigand helpfully reminded you to include it.).
Now, the angular part is just
Now, you grab a spare bit of paper and try to tackle this triple angular integral that has cost you so much effort already.
You start scribbling on the second page to quickly work out this triple integral...
Let the integrand be .
Using the identity , you can rewrite the absolute value:
By shift invariance[2], you integrate over the differences and and multiply by the length of the interval :
For more detail
Formally, you can perform a change of variables:
The Jacobian determinant of this transformation is , so .
Substituting this into the integral:
Notice that the integrand only depends on and , not on . This allows you to separate the integrals:
The first integral simply evaluates to . The remaining double integral accounts for the relative differences and over their periods.
Now, substitute and . The limits become and the differentials become :
The term is positive when and negative when . Due to symmetry, the integral over these two triangular regions is identical. You calculate the positive region and multiply by 2.
You use the product-to-sum formula to integrate with respect to .
Using the product-to-sum formula on the inner term:
Substituting back into :
Evaluating the two integrals separately:
- (via integration by parts)
and you finally get the answer: .
You feel you're very close now. You go back to the integral expression for and you plug in the angular integral value:
You simplify and expand the polynomial.
Distribute everything. The criminal seems very happy you're at this stage of the calculation.
Just the fractions now...
You've done it! The hooligan looks happy. He says you can go. He also says that everyone always calculates it this way, and he's yet to see a clever way that doesn't involve a sextuple integration. He invites you to share a clever solution later on by emailing it to him.
You wonder if you actually got the right answer. It was a nice fraction but did you not overcount? Undercount? Did you mess up the integral? The whole angular thing had like 6 permutations to it after all... Well, the villain seemed satisfied either way.
You go home. You're very tired. You plan on not doing the "dark alley at night" thing any more.
You don't remember why you were in this situation. ↩︎
Since the integrand depends only on the differences between angles, not their absolute values. ↩︎