Numerical Experiments in Billiards
Michael Turaev
1
Abstract
In this work we first talk about several physical motivations behind billiards and the rea-son why studying billiards is useful in certain areas of physics. From these, we deduce a mathematical model and talk about several methods for performing numerical compu-tations that allow us to construct and analyse the dynamics of closed smooth boundary billiards in a simulation. Finally, we will analyse how the dynamics of the elliptical billiard changes with the introduction of a perturbation.
2
Introduction
A billiard is a model in which a particle moves freely within a confined region and bounces off the walls elastically. Physically this can be interpreted in several different ways.
The goal of this work is to analyse whether billiards is an appropriate model for kinetic gas theory. In this report we will be focusing on two dimensional smooth closed convex billiard and we will check whether they follow the ergodic hypothesis, which is assumed for kinetic gas theory and for thermodynamics. The ergodic hypothesis states that in a system of fixed energy all possible states will be achieved over a long period of time.
The billiards that will be considered are the circle and the ellipse. These can be fully explained using trigonometry, the results of which will be shown. A mathematical model of the billiards will be constructed and using numerical methods will be simulated. Then a perturbation on the elliptical billiard will be considered and how the perturbation changes the dynamics of a billiard will be discussed.
3
Physical motivation
A billiard can be interpreted in different ways, each way will result in allowing us to study a different area of physics in more detail.
3.1
Geodesic flow
Consider a particle confined to geodesic flow on the surface of a volume. The billiard can be seen as the limit of flattening the volume until a 2 dimensional plate is formed. In the limit when the particle hits the boundary it will go over it and resume the trajectory from underneath the surface, but when viewed from the top the trajectory looks identical to a billiard. For example take the volume to be a sphere, in the limit we will obtain the circular billiard. As the particle moves along the geodesics it satisfies the principle of least action, so the particle can also be considered a light ray travelling inside a boundary with total internal reflection.
3.2
Fermat’s principle
Consider a light ray travelling with total internal reflection inside a boundary moving from point A to point B with exactly 1 bounce. Due to the total internal reflection the light ray does not lose energy to the surrounding and energy is conserved. If we let the system evolve we obtain a picture that looks exactly like a billiard.
3.3
Hamiltonian
Consider a particle which travels freely inside a boundary B that has the following Hamil-tonian
The boundary can be described by a function. For example the circular billiard is constructed by the B satisfying the following zero level set
f (x, y) = x2+ y2− 1. (3) As the collisions are elastic the kinetic energy of the particle has to be conserved. This is achieved by keeping the p2
i term constant, where i = x, y. This can be achieved
by changing pi into −pi, this corresponds to changing the direction. This is equivalent to
the laws of reflection.
4
Mathematical model
Now that the physical motivations have been discussed a mathematical model can built to satisfy the conditions we need.
Firstly we discuss the coordinates of billiards. We define them in 2 different ways. The first way is done by assigning normal Cartesian coordinates to the particle. These will mainly be used in the program and to draw the trajectories.
The second way is done by the parametrisation of the arclength i.e. we only use one coordinate s ∈ [0, 1). These will mainly be used for the phase map which plots the current point si on the horizontal axis and the landing position si+1 on the vertical axis.
Secondly we discuss the initial conditions. They are given as the initial position of the particle s0 and the angle θ that the trajectory makes to the boundary. This is equivalent
to knowing the initial position and the position of the first collision with the boundary s1, as each angle will give a unique position of the landing point. By the laws of reflection
if the first two initial points on the boundary are known, then the third point will be unique. For a specific billiard to fully construct the trajectory all we need to know are the numeric values of s0 and s1.
4.1
The circle
The following are the possible trajectories in the case of the circular billiard.
If the initial point is at s0 = 0 and the landing point is halfway across of the boundary
s1 = 0.5, then it will bounce back to the initial point.
Figure 1: Trajectory of a circular billiard with initial conditions at s0 = 0 and s1 = 0.5.
If the landing point s1 is a rational number of the arclength the trajectory is periodic
and the particle will always return to the initial starting point.
Figure 2: Periodic trajectory of a circular billiard with initial conditions at s0 = 0 and
If the landing point s1 is an irrational number, the trajectory is quasi periodic and
the particle will not return to the initial position, but it will get arbitrarily close to the initial starting point.
Figure 3: Trajectory of a circular billiard with initial conditions at s0 = 0 and s1 = 2.45/π.
All trajectories are tangent to a circle inside the circular boundary, this is called a caustic [1].
4.2
The ellipse
The following are the possible trajectories in the case of the elliptical billiard.
If the initial point is at s0 = 0 and the landing point is halfway the boundary at
s1 = 0.5 then the particle will bounce back to the initial position, same as in the circle.
If the first trajectory crosses a focus point of the ellipse, it will hit the second focus and converge to the line between the two foci.
Figure 5: Trajectory of an elliptical billiard with initial conditions at s0 = 0.9825 and
s1 = 0.4686.
We observe two different trajectories depending on whether the initial trajectory crosses the line that connects the two foci in between or outside the two foci. If it crosses in between the two foci then all the reflections will do the same.
If the initial trajectory crosses outside of the two foci then all the reflections will also do the same.
Figure 7: Trajectory of an elliptical billiard with initial conditions at s0 = 0.1 and s1 = 0.4.
All these trajectories are tangent to cofocal conics. In the first case it is a hyperbolic conic, and in the second case it is an elliptical conic.
4.3
The boundary
The boundary is modelled as a zero level set of a function f (x, y). In the case of the circle the function is given by
f (x, y) = x2+ y2− 1. (4) In the case of the ellipse the function is given by
f (x, y) =x a 2 +y b 2 − 1 (5)
5
Numerical methods
The billiard model will be simulated using python (see Appendix). Therefore the model needs to be constructed using numerical methods.
We start with two points in Cartesian coordinates which will correspond to the initial conditions of the billiard. Define M to be the vector pointing from the first point to the second point. It is obtained from the equation of the line joining the two points. Then we follow the vector incrementally in small steps. As we go in steps we will jump over the boundary when we get close to it. This is checked by plugging in the values for the coordinates of the initial position of the particle into the equation of the boundary, this gives us a value v1. Then we subtract the obtained value from the value we would have
received if we used the initial position, v0. If the answer is negative we are still within the
boundary, and if the answer is positive we have jumped over the boundary. When this occurs we refine the position of the particle using the Newton method on the equation of the boundary and the equation of the line until the difference of v1− v0 is of order 10−12.
reflected along the tangent of boundary and then the process is repeated to find the next landing point.
The boundary is plotted by turning the equation of the boundary into two coupled dif-ferential equations. We consider the following initial value problem
˙x ˙ y = ∇f (x, y) ||∇f (x, y)|| n , x(0) = 1, y(0) = 0, (6) where (v, w)n= (w, −v).
This is computed using the Runge-Kutta 4 method numerically. This is also used to find the arclength coordinate of the particle, which is done by integrating up to the landing point and recording the arclength.
6
Results
These are the results obtained from the simulations.
6.1
The circle
Figure 8: The trajectory of a circular billiard obtained from 100 reflections with s0 = 0
and s1 = 2.45/π.
6.2
The ellipse
The following plot corresponds to the trajectory and phase map when the initial trajectory crosses inside of the two foci.
(a) Plot of the trajectory in Cartesian space. (b) Plot of the trajectory in a phase map.
Figure 9: The trajectory (a) and the phase map (b) of an elliptical billiard obtained from 100 reflections with s0 = 0.1 and s1 = 0.8.
The following plot corresponds to the trajectory and phase map when the initial trajectory crosses outside of the two foci.
(a) Plot of the trajectory in Cartesian space. (b) Plot of the trajectory in a phase map.
Figure 10: The trajectory (a) and the phase map (b) of an elliptical billiard obtained from 100 reflections with s0 = 0.1 and s1 = 0.4.
Now we will plot many initial conditions on one phase map, so that many initial conditions can be analysed at the same time.
Figure 11: The phase map of an elliptical billiard obtained with 90 different initial con-ditions and 100 reflections each.
We can see that the ellipse also only has caustic solutions and is therefore completely non ergodic.
7
Pertubation
We will now introduce a perturbation on to the elliptical billiard in the following way f (x, y) =x
a 2
+ y2− 1 + sin(x), (7) where is small.
First we will look at how the perturbation affects the single trajectory.
(a) = 0 (b) = 0.1
(c) = 0.2 (d) = 0.3
Figure 12: How the trajectory of an elliptical billiard changes when a perturbation is applied with initial conditions s0 = 0.05 and s1 = 0.45 with increasing in steps of 0.1
Now we look at the same trajectory as above, but plotted in a phase map.
(a) = 0 (b) = 0.1
(c) = 0.2 (d) = 0.3
Figure 13: How the phase map of an elliptical billiard changes when a perturbation is applied with initial conditions s0 = 0.05 and s1 = 0.45 with increasing in steps of 0.1
We see that firstly the caustic breaks up into islands, which in turn spread out into the chaotic region. Lastly we will look at how the perturbation affects the phase map with many initial conditions.
(c) = 0.2 (d) = 0.3
Figure 15: How the phase map of an elliptical billiard changes when a perturbation is applied with 90 different initial conditions, 100 reflections each.
As we can see more caustics break down and become chaotic as we increase the per-turbation, but there always remain some caustic solutions even when the system is highly perturbed.
An example trajectory that shows this, is a trajectory that has initial points that are very close to each other. This is portrayed below
(a) Plot of the trajectory in Cartesian space. (b) Plot of the trajectory in a phase map.
We will now look at a more detailed phase map with each trajectory having its own colour.
Figure 17: The phase map of an elliptical billiard obtained with 210 different initial conditions and 1000 reflections each with = 0.5.
We can see that the caustic solutions are a singular colour, and the chaotic region have a lot of colour in them.
Near the hyperbolic caustics we can see that there are regions of single colour chaos. Which corresponds to a hyperbolic caustic breaking.
Near the elliptical caustics we can see a lot of chaos, and small islands forming near them.
8
Conclusion
We have seen that there are billiards which show no ergodicity at all, and there are billiards that partially show ergodicity. This shows that the billiard model is an incomplete model, so we need to consider adding more hypotheses to develop the model further so that it will satisfy the kinetic gas theory. As in a gas there are many particles confined in a volume, collisions could be a factor in the dynamics of the particles. This should be considered in a further model that will expand on this one. Non convex billiards should also be studied further, as they will provide more results on the dynamics of billiards.
