Galperin’s billiard method of computing pi

For the lack of original pi-related ideas, I describe a billiard experiment invented and popularized by Gregory Galperin (Григорий Гальперин). I read about it in a short note he published in 2001 in the Russian journal “Mathematical Education” (“Математическое Просвещение”). Although the table of contents is in Russian, it’s not hard to find the article with “pi” in the title. It’s clear from the article that the observation precedes its publication by some years; the author gave talks on this subject for a while before committing it to print.

Place two billiard balls (red and blue) on the half-line [0,\infty): the red ball is at rest, the blue ball is to the right of it and is moving the the left. At x=0 there is an elastic wall.

1-dimensional billiard table

Assume perfectly elastic collisions: the kinetic energy is preserved. How many collisions will happen depends on the masses of the ball, or rather their ratio. Suppose they have equal mass (ratio=1). Then we’ll see the following:

  1. Blue hits Red and stops; Red begins to move to the left
  2. Red hits the wall and bounces back
  3. Red hits Blue and stops; Blue begins to move to the right

That’s all, 3 collisions in total. Here is a Scilab-generated illustration: I placed red at x=1 and blue at x=1.5, with initial velocity -1. The time axis is vertical.

Mass ratio 1 leads to 3 collisions

Number 3 is also the first digit of \pi. Coincidence? No.

Consider the blue/red mass ratio of 100. Now we really need scilab. Count the collisions as you scroll down:

Mass ratio 100 leads to 31 collisions

Yes, that’s exactly 31 collisions. The mass ratio of 10000 would yield 314 collisions, and so on.

How does this work? The trick is to consider the two-dimensional configuration space in which the axes are scaled proportional to the square roots of masses. A point on the configuration plane describes the position of two balls at once, but it also behaves just like a billiard ball itself, bouncing off the lines (red=0) and (red=blue). The scaling by square roots of masses is needed to make the (red=blue) bounce work correctly: the angle of incidence is equal to the angle of reflection. (See Update 2 below.)

Configuration space

The opening angle is \alpha=\tan^{-1}\sqrt{M_r/M_b}. Now we can straighten the billiard trajectory by repeated reflection: this will help us count the number of collisions.

Straightening a billiard trajectory by reflection

It remains to count how many angles of size \alpha=\tan^{-1}\sqrt{M_r/M_b} fit inside of angle \pi. By taking \sqrt{M_r/M_b}=10^{-n}, we should get 10^n \pi collisions, rounded to an integer one way or the other. It seems that the rounding is always down, that is, we get exactly the first n digits of \pi. I once tried to prove that this always happens, but unsuccessfully. The error of approximation \tan^{-1} 10^{-n}\approx 10^{-n} needs to be compared to the fractional part of 10^n \pi, which looks like tricky business.

Here’s the scilab source. The plots in this post were obtained with galperin(1) and galperin(100).

function galperin(m)
    h=0.001; steps=4000;  x1=zeros(steps); x2=zeros(steps);
    x1(1)=1.5;  v1=-1;  x2(1)=1;  v2=0; collision=0; j=1;
    while j<steps,
        select collision
        case 0 then
            newx1=x1(j)+h*v1; newx2=x2(j)+h*v2;
            if (newx2<0) then collision=1;
            elseif (newx1<newx2) then collision=2;
                x1(j+1)=newx1; x2(j+1)=newx2; j=j+1;
        case 1 then
            v2=-v2; x1(j+1)=x1(j); x2(j+1)=x2(j); collision=0; j=j+1;
        case 2 then
            newv1=(v1*(m-1)+2*v2)/(m+1);  newv2=(v2*(1-m)+2*m*v1)/(m+1);
            v1=newv1; v2=newv2; x1(j+1)=x1(j); x2(j+1)=x2(j);
            collision=0; j=j+1;
    t=linspace(0, h*steps, steps); plot(t, x1); plot(t, x2, 'red');

Update: this post was referred to from Math StackExchange, with a request for a more intuitive explanation of how the scaling by \sqrt{M} works. Let V_c be the velocity vector of the point in the configuration space. Thanks to the scaling, we have the following:

  • The magnitude of V_c is equal to the total kinetic energy of the system, up to a constant factor.
  • The projection of V_c onto the collision line x_r=x_b is equal to the total momentum of the system, up to a constant factor.

When the balls collide with each other, both the energy and the momentum are preserved. It follows that the new vector V_c will have the same magnitude and the same projection onto the collision line. By trigonometry, this implies the equality of the incident and reflected angles.

3 thoughts on “Galperin’s billiard method of computing pi

  1. I tried implementing this code myself (in MATLAB, as usual, though the approach is the same), and the most I can get is the first four digits before I start running into a bunch of computation errors. In particular, with a mass ratio of 100,000,000:1, the change in velocity imparted by the collisions with the smaller billiard ball is never enough to turn the larger ball around (I think this may even be a floating point error, where the change in velocity is rounded to 0). I wonder if there is any reasonable way to implement this and get, say, 10 digits.

  2. I think that this is a wonderful post! I want to confirm that I got similar results by running a simulation using Mathematica. A truly fascinating result!

Comments are closed.