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 answerSo 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

~~two~~three 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*