The June edition of the IBM Ponder This website poses the following puzzle:

*Assume that cars have a length of two units and that they are parked along the circumference of a circle whose length is 100 units, which is marked as 100 segments, each one exactly one unit long.*

*A car can park on any two adjacent free segments (i.e., it does not need any extra maneuvering space).*

*Our question is as follows: Let's assume that we start with an empty circle. We add one car at a time, and each car parks in a random free space (aligned to a unit length), till no such place exists. What is the expected number of cars that can park along that 100-unit-long circle?*

Let C(n) be the expected number of cars that can park in a circular lot with n spaces. If n>1, then at least one car can park, and after doing so there will be n-2 contiguous spaces remaining. Let L(n) be the expected number of cars that can park in n contiguous spaces (not in a circle).

Clearly for n>1,

(1) C(n)=L(n-2)+1

Let's now find an expression for L(n). Number the stripes that mark the boundaries of the spaces from 0 to n. Because a car is two spaces long, the middle of a parked car is directly over a stripe. The first car can park on the stripes numbered 1 to n-1. Let the first car park over stripe k. That leaves two contiguous sections, with spaces k-1 and n-1-k. The sum of the expected number of cars that can park in those sections is L(k-1) + L(n-1-k). For the first driver, there is an equal probability, 1/(n-1), of choosing any one of the n-1 available spaces. Consequently,

(2) L(0) = L(1) = 0

L(n) = 1/(n-1)*sum(1 + L(k-1) + L(n-1-k), k=1..n-1)

This recurrence can be easily converted to a Maple procedure, that, with the ?remember option, is reasonably efficient and makes short work of this puzzle.

L := proc (n) local k;

option remember;

if n < 2 then

return 0;

else

1/(n-1)*add(1+L(k-1)+L(n-1-k), k = 1 .. n-1);

end if

end proc:

For the given problem,

C(100) = 1 + L(98) = 42.2332

While the puzzle only required a numerical result, we would like to analyze the general case and find a closed-form expression for L(n). A nice way to do this is with generating functions. Herbert Wilf's *GeneratingFunctionology* is a pleasant introduction packed with useful and fun examples. The following analysis references the third edition. The author has generously made the second edition available at his website.

The ordinary generating function for the sequence L(0), L(1), ... is the formal power series (FPS)

(3) lambda(x) = L(0) + L(1)*x + L(2)*x^2 + L(3)*x^3 + ...

Let's rewrite recurrence (2) using the Iverson convention whereby [boolean] = 1 if boolean is true, 0 otherwise:

(n-1)*L(n) = sum([n>1] + L(k-1) + L(n-1-k), k=1..n-1)

= (n-1)*[n>1] + 2*sum(L(k)=0..n-2)

Multiply both sides by x^n and sum over n

(4) sum((n-1)*L(n)*x^n,n) = sum((n-1)*[n>1]*x^n,n) + 2*sum(sum(L(k-1),k=0..n-2),n)

The left side of (4) can be evaluated using rule 2, section 2 of Wilf. Here we want the ordinary power series form of the sequence given by (n-1)*L(n). As explained in the book, this can be computed by substituting the operator x*D for the n in the polynomial and applying that to L(n). That is,

sum((n-1)*L(n)*x^n,n) (x.D-Id)(lambda(x)) = x*lambda'(x) - lambda(x)

The first term of the right side of (4) can be evaluated directly with Maple

sum((n-1)*x^n,n=1..infinity) = x^2/(1-x)^2

To evaluate the double summation in (4), recast the sum over k as a sum from 0 to n,

sum(L(k-1),k=1..n-1) = sum(L(k),k=0..n) - L(n-1) - L(n)

Multiplying by x^n, summing over n, and then applying rule 3, section 2 of Wilf gives an expression for the double summation

lambda(x)*(1/(1-x) - x - 1)

Plugging these results into (4) gives the differential equation

(5) deqs := x*lambda'(x)*x^n + 1/n*lambda(x)*x^n - lambda(x)*x^n

= 1/n*x^n - x^n + 2*lambda(x)*(1/(1-x)-x-1)

From the definition of lambda(x), lambda(0) = L(0) = 0 and lambda'(0) = L(1) = 0, so

ics := ( lambda(0) = 0, lambda'(0) = 0 ):

Maple's dsolve procedure returns the following solution

dsol := dsolve({deq, ics});

lambda(x) = (1-exp(-2*x))*x/(1-x)^2/2

To check that it is correct, use series to compute several coefficients and compare them with a direct computation using our procedure (above).

series(rhs(dsol),x,10);

x^2+x^3+(5/3)*x^4+2*x^5+(37/15)*x^6+(26/9)*x^7+O(x^8)

Now that we have a closed-form for the generating function, we can, with a bit of further analysis, extract an expression for L(n).

Move the x/2 term to the left side of the equation to get

2/x*lambda(x) = (1-exp(-2*x))/(1-x)^2

Compute the FPS of each of the two terms of the product on the right (that can be done using convert(., Sum, 'dummy'=n):

2/x*lambda(x) = (1 - Sum((-2*x)^n/n!,n>=0))*Sum((n+1)*x^n,n>=0)

= Sum((n+1)*x^n,n>=0) - Sum((-2*x)^n/n!,n>=0))*Sum((n+1)*x^n,n>=0)

Compute the coefficient of x^n on each side of the equation

2*L(n+1) = (n+1) - Sum((n+1-k)*(-2)^k/k!,k=0..n)

Substitute n -> n-1 and move the factor of 2 back to the right

L(n) = 1/2*(n - Sum((n-k)*(-2)^k/k!,k=0..n-1))

Maple can evaluate the sum,

L(n) = 1/2*(n-(n+2)*exp(-2)

+(-1)^n*2^((1/2)*n)*exp(-1)*WhittakerM(1-(1/2)*n,(1/2)*n+1/2,2)/(n+1)!)

While the rhs is messy, numerical evaluation reveals that the final term quickly goes to 0 as n increases. The result is

L(n) ~ Lx(n) = 1/2*(n-(n+2)*exp(-2))

Using this to compute C(100) gives

C(100) ~ 1 + Lx(98) = 42.2332

which agrees with our exact computation.

From the result, we see that the expected utilization (ratio of expected to maximum) of a linear parking lot quickly reaches

limit(L(n)/(n/2),n=infinity) = 1-exp(-2) = 0.865.

With that modest success, the natural inclination is to generalize the problem until we arrive at one we cannot solve. Alas, that doesn't take long. Suppose each car takes up three spaces. While the same approach should theoretically work, I could not get Maple to symbolically solve the resulting differential equation. Is there a way to continue without it?

A more interesting generalization is to remove the stripes and allow cars to park anywhere. Instead of a difference equation we get a differential equation, specifically, a delay-differential equation. An analytical solution looks hopeless, though for a given size lot it can be solved numerically by integrating a piece at a time. Can we find the expected utilization for large n?