This problem arose in a different context at work, but I have translated it to pizza.
Suppose you have a circular pizza of radius $R$. Upon this disc, $n$ pepperoni will be distributed completely randomly. All pepperoni have the same radius $r$.
A pepperoni is “free” if it does not overlap any other pepperoni.
You are free to choose $n$.
Suppose you choose a small $n$. The chance that any given pepperoni is free are very large. But $n$ is small so the total number of free pepperoni is small. Suppose you choose a large $n$. The chance that any given pepperoni is free are small. But there are a lot of them.
Clearly, for a given $R$ and $r$, there is some optimal $n$ that maximizes the expected number of free pepperoni. How to find this optimum?
Edit: picking the answer
So it looks like leonbloy’s answer given the best approximation in the cases I’ve looked at:
r/R n* by simulation n_free (sim) (R/2r)^2 0.1581 12 4.5 10 0.1 29 10.4 25 0.01 2550 929.7 2500(There’s only a few hundred trials in the r=0.01 sim, so 2550 might not be super accurate.)
So I’m going to pick it for the answer. I’d like to thank everyone for their contributions, this has been a great learning experience.Here are a few pictures of a simulation for r/R = 0.1581, n=12:
Edit after three answers posted:
I wrote a little simulation. I’ll paste the code below so it can be checked (edit: it’s been fixed to correctly pick points randomly on a unit disc). I’ve looked at
twothree cases so far. First case, r = 0.1581, R = 1, which is roughly p = 0.1 by mzp’s notation. At these parameters I got n* = 12 (free pepperoni = 4.52). Arthur’s expression did not appear to be maximized here. leonbloy’s answer would give 10. I also did r = 0.1, R = 1. I got n* = 29 (free pepperoni = 10.38) in this case. Arthur’s expression was not maximized here and leonbloy’s answer would give 25. Finally for r = 0.01 I get roughly n*=2400 as shown here:Here’s my (ugly) code, now edited to properly pick random points on a disc:
from __future__ import division import numpy as np # the radius of the pizza is fixed at 1 r = 0.1 # the radius of the pepperoni n_to_try = [1,5,10,20,25,27,28,29,30,31,32,33,35] # the number of pepperoni trials = 10000# the number of trials (each trial randomly places n pepperoni) def one_trial(): # place the pepperoni pepperoni_coords = [] for i in range(n): theta = np.random.rand()*np.pi*2 # a number between 0 and 2*pi a = np.random.rand() # a number between 0 and 1 coord_x = np.sqrt(a) * np.cos(theta) # see http://mathworld.wolfram.com/DiskPointPicking.html coord_y = np.sqrt(a) * np.sin(theta) pepperoni_coords.append((coord_x, coord_y)) # how many pepperoni are free? num_free_pepperoni = 0 for i in range(n): # for each pepperoni pepperoni_coords_copy = pepperoni_coords[:] # copy the list so the orig is not changed this_pepperoni = pepperoni_coords_copy.pop(i) coord_x_1 = this_pepperoni[0] coord_y_1 = this_pepperoni[1] this_pepperoni_free = True for pep in pepperoni_coords_copy: # check it against every other pepperoni coord_x_2 = pep[0] coord_y_2 = pep[1] distance = np.sqrt((coord_x_1 - coord_x_2)**2 + (coord_y_1 - coord_y_2)**2) if distance < 2*r: this_pepperoni_free = False break if this_pepperoni_free: num_free_pepperoni += 1 return num_free_pepperoni for n in n_to_try: results = [] for i in range(trials): results.append(one_trial()) x = np.average(results) print "For pizza radius 1, pepperoni radius", r, ", and number of pepperoni", n, ":" print "Over", trials, "trials, the average number of free pepperoni was", x print "Arthur's quantity:", x* ((((1-r)/1)**(x-1) - (r/1)) / ((1-r) / 1))
Answer
Updated: see below (Update 3)
Here’s another approximation. Consider the center of the disks (pepperoni) as an homogeneous point process of density $\lambda = n/A$, where $A=\pi R^2$ is the pizza surface. Let $D$ be nearest neighbour distance from a given center. Then
$$P(D\le d) = 1- \exp(-\lambda \pi d^2)=1- \exp(-n \,d^2/R^2) \tag{1}$$
A pepperoni is free if $D > 2r$. Let $E$ be the expected number of free peperoni.
Then $$E= n\, P(D>2r) = n \exp (-n \,4 \, r^2/R^2) = n \exp(-n \, p)$$
where $p=(2r)^2/R^2$ (same notation as mzp’s answer).
The maximum is attained for (ceil or floor of) $n^*=1/p=(R/2r)^2\tag{2}$
Update 1: Formula $(1)$ could be corrected for the border effects, the area near the border would be computed as the intersection of two circles. It looks quite cumbersome, though.
Update 2: In the above, I assumed that the center of the pepperoni could be placed anywhere inside the pizza circle. If that’s not the case, if the pepperoni must be fully inside the pizza, then $R$ should be replaced by the “effective radius” $R’ = R-r$
Update 3: The Poisson approach is really not necessary. Here’s an exact solution
Let $$t = \frac{R}{2r}$$
(Equivalently, think of $t$ as the pizza radius, and assume a pepperoni of radius $1/2$). Assume $t>1$. Let $g(x)$ be the area of a unit circle, at a distance $x$ from the origin, intersected with the circle of radius $t$. Then
$$g(x)=\begin{cases}\pi & {\rm if}& 0\le x \le t-1\\
h(x) & {\rm if}& t-1<x \le t \\
0 & {\rm elsewhere}
\end{cases}
\tag{3}$$
where
$$h(x)=\cos^{-1}\left(\frac{x^2+1-t^2}{2x}\right)+t^2
\cos^{-1}\left(\frac{x^2+t^2-1}{2xt}\right) -\frac{1}{2}\sqrt{[(x+t)^2-1][1-(t-x)^2]} \tag{4}$$
Here’s a graph of $g(x)/\pi$ for $t=5$
Let the random variable $e_i$ be $1$ if the $i$ pepperoni is free, $0$ otherwise.
Then
$$E[e_i \mid x] = \left(1- \frac{g(x)}{\pi \, t^2}\right)^{n-1} \tag{5}$$
(remaining pepperoni fall in the free area). And
$$E[e_i] =E[E(e_i \mid x)]= \int_{0}^t \frac{2}{t^2} x \left(1- \frac{g(x)}{\pi \, t^2}\right)^{n-1} dx = \\
=\left(1- \frac{1}{t^2}\right)^{n-1} \left(1- \frac{1}{t}\right)^2
+\frac{2}{t^2} \int_{t-1}^t x \left(1- \frac{h(x)}{\pi\, t^2}\right)^{n-1} dx
\tag{6}$$
The objective function (expected number of free pepperoni) is then given by:
$$J(n)=n E[e_i ] \tag{7} $$
This is exact… but (almost?) intractable. However, it can be evaluated numerically [**].
We can also take as approximation
$$g(x)=\pi$$ for $0\le x < t$ (neglect border effects) and then it gets simple:
$$E[e_i ] =E[e_i \mid x]= \left(1- \frac{1}{t^2}\right)^{n-1}$$
$$J(n)= n \left(1- \frac{1}{t^2}\right)^{n-1} \tag{8}$$
To maximize, we can write
$$\frac{J(n+1)}{J(n)}= \frac{n+1}{n} \left(1- \frac{1}{t^2}\right)=1 $$
which gives
$$ n^{*}= t^2-1 = \frac{R^2}{4 r^2}-1 \tag{9}$$
quite similar to $(2)$.
Notice that as $t \to \infty$, $J(n^{*})/n^{*} \to e^{-1}$, i.e. the proportion of free pepperoni (when using the optimal number) is around $36.7\%$. Also, the “total pepperoni area” is 1/4 of the pizza.
[**] Some Maxima code to evaluate (numerically) the exact solution $(7)$:
h(x,t) := acos((x^2+1-t^2)/(2*x))+t^2*acos((x^2-1+t^2)/(2*x*t))
-sqrt(((x+t)^2-1)*(1-(t-x)^2))/2 $
j(n,t) := n * ( (1-1/t)^2*(1-1/t^2)^(n-1)
+ (2/t^2) * quad_qag(x * (1-h(x,t)/(%pi*t^2))^(n-1),x,t-1,t ,3)[1]) $
j0(n,t) := n *(1-1/t^2)^(n-1)$
tt : 1/(2*0.1581) $
j(11,tt);
4.521719308511862
j(12,tt);
4.522913706608645
j(13,tt);
4.494540361913981
tt : 1/(2*0.1) $
j(27,tt);
10.37509984083333
j(28,tt);
10.37692859747294
j(29,tt);
10.36601271146961
fpprintprec: 4$
nn : makelist(n, n, 2, round(tt^2*1.4))$
jnn : makelist(j(n,tt),n,nn) $
j0nn : makelist(j0(n,tt),n,nn) $
plot2d([[discrete,nn,jnn],[discrete,nn,j0nn]],
[legend,"exact","approx"],[style,[linespoints,1,2]],[xlabel,"n"], [ylabel,"j(n)"],
[title,concat("t=",tt, " r/R=",1/(2*tt))],[y,2,tt^2*0.5]);
The first result agrees with the OP simulation, the second is close: I got $n^{*}=28$ instead of $29$. For the third case, I get $n^{*}=2529$ ($j(n)=929.1865331$)
Attribution
Source : Link , Question Author : Ben S. , Answer Author : leonbloy
