BeginPackage["UVW`Billiard`"]
NextBounce::usage = "
NextBounce[currentposition,r] computes the next bouncing point of a ball
inside a square table with a circular obstacle of radius r. The initial
conditions are given in the list currentposition that has four real
coordinates, respectively the absissa, ordinate, incoming angle of the ball,
and the current time. The result is returned as another list of four elements,
the abscissa and ordinate of the new bouncing point, the new direction of
the ball and the current time incremented by the running time between the
two bounces."
Trajectory::usage ="
Trajectory[shootingangle,tmax,r] draws a square table with a centered
circular obstacle of radius r. Draws inside this support the
trajectory of a ball starting at the bottom left corner with initial
direction shootingangle, up to time tmax. "
Bundle::usage ="
Bundle[listofangles,tmax,r] draws a square table with a centered circular
obstacle of radius r. Draws inside the trajectories of balls starting at the
bottom left corner with initial directions read in listofangles, up to
time tmax. "
Differences::usage = "
Differences[angle,dangle,tmax,r] simulates two trajectories of balls in a
square table with a centered circular obstacle with radius r. Both
trajectories start from the bottom left corner. They are followed up to time
tmax. The shooting angle of the first trajectory is angle, its difference
with the second shooting angle is dangle. The function represents three
consecutive graphics. The first one is the billiard table with the two
trajectories . The second one is the evolution of the absolute difference
of angles as a function of time. The third one is the norm of the difference
of positions as a function of time. "
Begin ["`Private`"]
NextBounce[currentposition_List,r_]:=
Block[
{xold,yold,aold,cosaold,sinaold,epsi=10^(-6),
tnorth,teast,tsouth,twest,tcircle1,tcircle2,delta,beta,
trun,xnew,ynew,anew},
(* Old coordinates ----------------------------*)
xold = currentposition[[1]];
yold = currentposition[[2]];
aold = currentposition[[3]];
cosaold = Cos[aold];
sinaold = Sin[aold];
(* Hitting time of the four edges -------------*)
tnorth = (+1-yold)/sinaold; If[tnorth0,
(
delta = Sqrt[delta];
tcircle1 = -(xold*cosaold+yold*sinaold) + delta;
If[tcircle11-epsi),
If[(aold1
]
]
Bundle[listofangles_List,tmax_,r_]:=
Block[
{nbtraj,running=0.,traj,next,i,g1,g2},
nbtraj = Length[listofangles];
traj=Table[{{-1.,-1.}},{nbtraj}];
Do[(
next = {-1.,-1.,N[listofangles[[i]]],0.};
running = 0.;
While[running1];
];
Differences[angle_,dangle_,tmax_,r_]:=
Block[
{listofangles,running=0.,traj,next,i,g1,g2,
angles,diffangles,aold,diffpos,i1,i2},
listofangles = N[{angle,angle+dangle}];
traj=Table[{{-1.,-1.,listofangles[[i]],0.}},{i,1,2}];
Do[(
next = {-1.,-1.,listofangles[[i]],0.};
running = 0.;
While[running1];
(* Representation of angles differences -------*)
Show[Graphics[{
Thickness[0.005],
Line[diffangles]}],
Axes->True,
Frame->True,
AxesLabel->{"Time","Difference of angles"}
];
(* Representation of positions differences ----*)
Show[Graphics[{
Thickness[0.005],
Line[diffpos]}],
Axes->True,
Frame->True,
AxesLabel->{"Time","Difference of positions"}
]
];
End[]
EndPackage[]