DynamicalBilliards.Testing.billiards_testset
— Methodbilliards_testset(description, f, args...; caller = all_tests)
Wrap a testset around caller(f, args...)
which times the result and separates the test from others (println
).
DynamicalBilliards.MushroomTools
— ModuleMushroomTools
Module containing many functions helpful in simulating (perfect) mushroom billiards, see billiard_mushroom
. Contains stuff like initializing efficiently regular or chaotic particles and functions that return the corresponding chaotic or regular phase-space volumes or portions. The functions V_3D_tot
and V_3D_reg
use equations derived in ref. [1].
Made by Lukas Hupe.
References
[1] A. Barnett & T. Betcke, Chaos 17, 043125 (2007).
DynamicalBilliards.MushroomTools.V_2D_tot
— MethodV_2D_tot(l,w,r)
Return the total boundary map volume (2D) of a billiard_mushroom
parameterized by (l,w,r)
.
DynamicalBilliards.MushroomTools.V_3D_tot
— MethodV_3D_tot(l,w,r)
Return the total phasespace volume (3D) of a billiard_mushroom
parameterized by (l,w,r)
.
DynamicalBilliards.MushroomTools.g_c_2D
— Methodg_c_2D(l, w, r)
Return the chaotic phasespace portion of the boundary map (2D) of a billiard_mushroom
with stem length l
, stem width w
and cap radious r
. This result is known analytically, see MushroomTools
for references.
DynamicalBilliards.MushroomTools.g_c_3D
— Methodg_c_3D(l, w, r)
Return the chaotic phasespace portion of the full (3D) phase-space of a billiard_mushroom
with stem length l
, stem width w
and cap radious r
. This result is known analytically, see MushroomTools
for references.
DynamicalBilliards.MushroomTools.g_r_2D
— Methodg_r_2D(l, w, r)
Return the regular phasespace portion of the boundary map (2D) of a billiard_mushroom
with stem length l
, stem width w
and cap radious r
. This result is known analytically, see MushroomTools
for references.
DynamicalBilliards.MushroomTools.g_r_3D
— Methodg_r_3D(l, w, r)
Return the regular phasespace portion of the full (3D) phase-space of a billiard_mushroom
with stem length l
, stem width w
and cap radious r
. This result is known analytically, see MushroomTools
for references.
DynamicalBilliards.MushroomTools.randin_mushroom
— Methodrandin_mushroom(l, w, r [, ω])
Generate a random particle within the billiard_mushroom
parameterised by l
, w
and r
. If ω
is given the particle is magnetic instead.
This function is much more efficient than randominside
.
DynamicalBilliards.MushroomTools.randomchaotic
— Methodrandomchaotic(l, w, r)
Generate a chaotic particle, i.e. not trapped in the cap.
DynamicalBilliards.MushroomTools.randomregular
— Methodrandomregular(l, w, r)
Generate a regular particle (i.e. trapped in the cap).
DynamicalBilliards.DynamicalBilliards
— ModuleA Julia package for dynamical billiard systems in two dimensions.
It provides a flexible, easy-to-use, and intuitive framework for fast implementation of billiards of arbitrary construction.
DynamicalBilliards.AbstractParticle
— TypeAbstractParticle
Particle supertype.
DynamicalBilliards.Antidot
— TypeAntidot{T<:AbstractFloat} <: Circular{T}
Disk-like obstacle that allows propagation both inside and outside of the disk (mutable type). Used in ray-splitting billiards.
Fields:
c::SVector{2,T}
: Center.r::T
: Radius.pflag::Bool
: Flag that keeps track of where the particle is currently propagating (pflag
= propagation-flag).true
stands for outside the disk,false
for inside the disk. Defaults totrue
.name::String
: Name of the obstacle given for user convenience. Defaults to "Antidot".
DynamicalBilliards.Billiard
— MethodBilliard(obstacles...)
Construct a Billiard
from given obstacles
(tuple, vector, varargs).
For functions like boundarymap
, it is expected (if possible) that the obstacles of the billiard are sorted, such that the arc-coordinate ξ
around the billiard is increasing counter-clockwise.
ξ
is measured as:
- the distance from start point to end point in
Wall
s - the arc length measured counterclockwise from the open face in
Semicircle
s - the arc length measured counterclockwise from the rightmost point in
Circular
s
DynamicalBilliards.Circular
— TypeCircular{T<:AbstractFloat} <: Obstacle{T}
Circular obstacle supertype.
DynamicalBilliards.Disk
— TypeDisk{T<:AbstractFloat} <: Circular{T}
Disk-like obstacle with propagation allowed outside of the circle (immutable type).
Fields:
c::SVector{2,T}
: Center.r::T
: Radius.name::String
: Some name given for user convenience. Defaults to "Disk".
DynamicalBilliards.Ellipse
— TypeEllipse{T<:AbstractFloat} <: Obstacle{T}
Ellipse obstacle that also allows ray-splitting. The ellipse is always oriented on the x and y axis (although you can make whichever you want the major one).
Fields:
c::SVector{2,T}
: Center.a::T
: x semi-axis.b::T
: y semi-axis.pflag::Bool
: Flag that keeps track of where the particle is currently propagating.true
(default) is associated with being outside the ellipse.name::String
: Some name given for user convenience. Defaults to"Ellipse"
.
The ellipse equation is given by
\[\left(\frac{x - c[1]}{a} \right)^2+ \left(\frac{y - c[2]}{b}\right)^2 = 1\]
DynamicalBilliards.FiniteSplitterWall
— TypeFiniteSplitterWall{T<:AbstractFloat} <: Wall{T}
Finite wall obstacle allowing fow ray-splitting (mutable type).
Fields:
sp::SVector{2,T}
: Starting point of the Wall.ep::SVector{2,T}
: Ending point of the Wall.normal::SVector{2,T}
: Normal vector to the wall, pointing to where the particle will come from before a collision (pointing towards the inside of the billiard). The size of the vector is irrelevant since it is internally normalized.isdoor::Bool
: Flag of whether thisFiniteSplitterWall
instance is a "Door". Defaults tofalse
.pflag::Bool
: Flag that keeps track of where the particle is currently propagating (pflag
= propagation flag).true
is associated with thenormal
vector the wall is instantiated with. Defaults totrue
.name::String
: Name of the obstacle, given for user convenience. Defaults to "Finite Splitter Wall".
DynamicalBilliards.FiniteWall
— TypeFiniteWall{T<:AbstractFloat} <: Wall{T}
Wall obstacle imposing specular reflection during collision (immutable type). Slower than InfiniteWall
, meant to be used for non-convex billiards.
Giving a true
value to the field isdoor
designates this obstacle to be a Door
. This is used in escapetime
function. A Door
is a obstacle of the billiard that the particle can escape from, thus enabling calculations of escape times.
Fields:
sp::SVector{2,T}
: Starting point of the Wall.ep::SVector{2,T}
: Ending point of the Wall.normal::SVector{2,T}
: Normal vector to the wall, pointing to where the particle will come from before a collision (pointing towards the inside of the billiard). The size of the vector is irrelevant since it is internally normalized.isdoor::Bool
: Flag of whether thisFiniteWall
instance is a "Door". Defaults tofalse
.name::String
: Name of the obstacle, given for user convenience. Defaults to "Finite Wall".
DynamicalBilliards.InfiniteWall
— TypeInfiniteWall{T<:AbstractFloat} <: Wall{T}
Wall obstacle imposing specular reflection during collision (immutable type). Faster than FiniteWall
, meant to be used for convex billiards.
Fields:
sp::SVector{2,T}
: Starting point of the Wall.ep::SVector{2,T}
: Ending point of the Wall.normal::SVector{2,T}
: Normal vector to the wall, pointing to where the particle will come from before a collision (pointing towards the inside of the billiard). The size of the vector is irrelevant since it is internally normalized.name::String
: Name of the obstacle, given for user convenience. Defaults to "Wall".
DynamicalBilliards.MagneticParticle
— TypeMagneticParticle(ic::AbstractVector{T}, ω::Real) # where ic = [x0, y0, φ0]
MagneticParticle(x0, y0, φ0, ω)
MagneticParticle(pos::SVector, vel::SVector, ω)
MagneticParticle(p::AbstractParticle, ω)
Create a magnetic particle with initial conditions x0, y0, φ0
and angular velocity ω
. It propagates as a circle instead of a line, with radius 1/abs(ω)
.
The field current_cell
shows at which cell of a periodic billiard is the particle currently located.
DynamicalBilliards.Obstacle
— TypeObstacle{<:AbstractFloat}
Obstacle supertype.
DynamicalBilliards.Particle
— TypeParticle(ic::Vector{T}) #where ic = [x0, y0, φ0]
Particle(x0, y0, φ0)
Particle(pos::SVector, vel::SVector)
Create a particle with initial conditions x0, y0, φ0
. It propagates as a straight line.
The field current_cell
shows at which cell of a periodic billiard is the particle currently located.
DynamicalBilliards.PeriodicWall
— TypePeriodicWall{T<:AbstractFloat} <: Wall{T}
Wall obstacle that imposes periodic boundary conditions upon collision (immutable type).
Fields:
sp::SVector{2,T}
: Starting point of the Wall.ep::SVector{2,T}
: Ending point of the Wall.normal::SVector{2,T}
: Normal vector to the wall, pointing to where the particle will come from (to the inside the billiard). The size of the vector is important! This vector is added to a particle'spos
during collision. Therefore the size of the normal vector must be correctly associated with the size of the periodic cell.name::String
: Name of the obstacle, given for user convenience. Defaults to "Periodic wall".
DynamicalBilliards.RandomDisk
— TypeRandomDisk{T<:AbstractFloat} <: Circular{T}
Disk-like obstacle that randomly (and uniformly) reflects colliding particles. The propagation is allowed outside of the circle.
Fields:
c::SVector{2,T}
: Center.r::T
: Radius.name::String
: Some name given for user convenience. Defaults to "Random disk".
DynamicalBilliards.RandomWall
— TypeRandomWall{T<:AbstractFloat} <: Wall{T}
Wall obstacle imposing (uniformly) random reflection during collision (immutable type).
Fields:
sp::SVector{2,T}
: Starting point of the Wall.ep::SVector{2,T}
: Ending point of the Wall.normal::SVector{2,T}
: Normal vector to the wall, pointing to where the particle is expected to come from (pointing towards the inside of the billiard).name::String
: Name of the obstacle, given for user convenience. Defaults to "Random wall".
DynamicalBilliards.RaySplitter
— TypeRaySplitter(idxs, transmission, refraction [, newangular]; affect)
Return a RaySplitter
instance, used to perform raysplitting. idxs
is a Vector{Int}
with the indices of the obstacles that this RaySplitter
corresponds to.
transmission
, refraction
and newangular
are functions. Let φ
be the angle of incidence and ω
be the angular velocity and pflag
the propagation flag (before transmission). The functions have the following signatures:
transmission(φ, pflag, ω) -> T
, transmission probability.refraction(φ, pflag, ω) -> θ
, refraction angle. This angle is relative to the normal vector.newangular(ω, pflag) -> newω
, new angular velocity after transmission.
The above three functions use the same convention: the argument pflag
is the one the obstacle has before transmission. For example, if a particle is outside an Antidot
(with pflag = true
here) and is transmitted inside the Antidot
(pflag
becomes false
here), then all three functions will be given their second argument (the Boolean one) as true
!
affect
is a function, and denotes which obstacles of the billiard are affected when transmission occurs at obstacle i
(for which obstacles should the field pflag
be reversed). Defaults to idxs = (i) -> i
, i.e. only the colliding obstacle is affected. If you want many obstacles to be affected you could write idxs = (i) -> SVector(2,3,5)
, etc. Keep in mind that the only values of i
that can be passed into this function are the ones that are given in the argument idxs
!
DynamicalBilliards.Semicircle
— TypeSemicircle{T<:AbstractFloat} <: Circular{T}
Obstacle that represents half a circle. Propagation is allowed only inside the semicircle.
Fields:
c::SVector{2,T}
: Center.r::T
: Radius.facedir::SVector{2,T}
: Direction where the open face of the Semicircle is facing.name::String
: Name of the obstacle given for user convenience. Defaults to "Semicircle".
DynamicalBilliards.SplitterWall
— TypeSplitterWall{T<:AbstractFloat} <: Wall{T}
Wall obstacle imposing allowing for ray-splitting (mutable type).
Fields:
sp::SVector{2,T}
: Starting point of the Wall.ep::SVector{2,T}
: Ending point of the Wall.normal::SVector{2,T}
: Normal vector to the wall, pointing to where the particle will come from before a collision. The size of the vector is irrelevant.pflag::Bool
: Flag that keeps track of where the particle is currently propagating (pflag
= propagation flag).true
is associated with thenormal
vector the wall is instantiated with. Defaults totrue
.name::String
: Name of the obstacle, given for user convenience. Defaults to "Splitter wall".
DynamicalBilliards.Wall
— TypeWall{T<:AbstractFloat} <: Obstacle{T}
Wall obstacle supertype.
All Wall
subtypes (except PeriodicWall
) can be called as Wall(sp, ep)
, in which case the normal vector is computed automatically to point to the left of v = ep - sp
.
DynamicalBilliards.acceptable_raysplitter
— Methodacceptable_raysplitter(raysplitters, bd::Billiard)
Return true
if the given raysplitters
can be used in conjuction with given billiard bd
.
DynamicalBilliards.acos1mx
— MethodApproximate arccos(1 - x) for x very close to 0.
DynamicalBilliards.arcintervals
— Methodarcintervals(bd::Billiard) -> s
Generate a vector s
, with entries being the delimiters of the arclengths of the obstacles of the billiard. The arclength from s[i]
to s[i+1]
is the arclength spanned by the i
th obstacle.
s
is used to transform an arc-coordinate ξ
from local to global and vice-versa. A local ξ
becomes global by adding s[i]
(where i
is the index of current obstacle). A global ξ
becomes local by subtracting s[i]
.
See also boundarymap
, to_bcoords
, from_bcoords
.
DynamicalBilliards.bdplot
— Functionbdplot(x; kwargs...) → fig, ax
bdplot!(ax::Axis, x; kwargs...)
Plot an object x
from DynamicalBilliards
into a given axis (or a new figure). x
can be an obstacle, a particle, a vector of particles, or a billiard.
bdplot!(ax,::Axis, o::Obstacle; kwargs...)
Keywords are propagated to lines!
or poly!
. Functions obfill, obcolor, obls, oblw
(not exported) decide global defaults for linecolor, fillcolor, linestyle, linewidth, when plotting obstacles.
bdplot!(ax,::Axis, bd::Billiard; clean = true, kwargs...)
If clean = true
, all axis elements are removed and an equal aspect ratio is establised. Other keywords are propagated to the obstacle plots.
bdplot!(ax,::Axis, bd::Billiard, xmin, xmax, ymin, ymax; kwargs...)
This call signature plots periodic billiards: it plots bd
along its periodic vectors so that it fills the total amount of space specified by xmin, xmax, ymin, ymax
.
bdplot!(ax,::Axis, ps::Vector{<:AbstractParticle}; kwargs...)
Plot particles as a scatter plot (positions) and a quiver plot (velocities). Keywords particle_size = 5, velocity_size = 0.05
set the size of plotted particles. Keyword colors = JULIADYNAMICS_CMAP
decides the color of the particles, and can be either a colormap or a vector of colors with equal length to ps
. The rest of the keywords are propagated to the scatter plot of the particles.
DynamicalBilliards.bdplot!
— Functionbdplot!(ax, args...; kwargs...)
See bdplot
.
DynamicalBilliards.bdplot_boundarymap
— Functionbdplot_boundarymap(bmap, intervals; figkwargs = NamedTuple(), kwargs...)
Plot the output of DynamicalBilliards.boundarymap
into an axis that correctly displays information about obstacle arclengths. Return figure, axis
.
Also works for the parallelized version of boundary map.
Keyword Arguments
figkwargs = NamedTuple()
keywords propagated toFigure
.color
: The color to use for the plotted points. Can be either a single color or a vector of colors of lengthlength(bmap)
, in order to give each initial condition a different color (for parallelized version).- All other keywords propagated to
scatter!
.
DynamicalBilliards.bdplot_interactive
— Functionbdplot_interactive(bd::Billiard, ps::Vector{<:AbstractParticle}; kwargs...)
Create a new figure with bd
plotted, and in it initialize various data for animating the evolution of ps
in bd
. Return fig, phs, chs
, where fig
is a figure instance and phs, chs
can be used for making custom animations, see below.
Keywords (interactivity-related)
playback_controls = true
: If true, add controls that allow live-interaction with the figure, such as pause/play, reset, and creating new particles by clicking and dragging on the billiard plot.
Keywords (visualization-related)
dt = 0.001
: The animation always occurs in steps of timedt
. A slider can decide how many stepsdt
to evolve before updating the plots.plot_bmap = false
: If true, add a second plot with the boundary map.colors = JULIADYNAMICS_CMAP
: If a symbol (colormap name) each particle gets a color from the map. If Vector of lengthN
, each particle gets a color form the vector. If Vector with length <N
, linear interpolation across contained colors is done.tail_length = 1000
: The lasttail_length
positions of each particle are visualized.tail_width = 1
: Width of the dtail plot.fade = true
: Tail color fades away.plot_particles = true
: Besides their tails, also plot the particles as a scatter and quiver plot.particle_size = 5
: Marker size for particle scatter plot.velocity_size = 0.05
: Multiplication of particle velocity before plotted as quiver.bmap_size = 4
: Marker size of boundary map scatter plot.figure = NamedTuple()
: Keywords propagated toFigure
creation.kwargs...
: Remaining keywords are propagated to the billiard plotting.
Custom Animations
Two helper structures are defined for each particle:
ParticleHelper
: Contains quantities that are updated eachdt
step: the particle, time elapsed since last collision, total time ellapsed, tail (positions in the lasttail_length
dt
-sized steps).CollisionHelper
: Contains quantities that are only updated at collisions: index of obstacle to be collided with, time to next collision, total collisions so far, boundary map point at upcoming collision.
These two helpers are necessary to transform the simulation into real-time stepping (1 step = dt
time), instead of the traditional DynamicalBilliards.jl setup of discrete time stepping (1 step = 1 collision).
The returned phs, chs
are two observables, one having vector of ParticleHelpers
, the other having vector of CollisionHelpers
. Every plotted element is lifted from these observables.
An exported high-level function bdplot_animstep!(phs, chs, bd, dt; update, intervals)
progresses the simulation for one dt
step. Users should be using bdplot_animstep!
for custom-made animations, examples are shown in the documentation online. The only thing the update
keyword does is notify!(phs)
. You can use false
for it if you want to step for several dt
steps before updating plot elements. Notice that chs
is always notified when collisions occur irrespectively of update
. They keyword intervals
is nothing
by default, but if it is arcintervals(bd)
instead, then the boundary map field of chs
is also updated at collisions.
DynamicalBilliards.bdplot_video
— Functionbdplot_video(file::String, bd::Billiard, ps::Vector{<:AbstractParticle}; kwargs...)
Create an animation of ps
evolving in bd
and save it into file
. This function shares all visualization-related keywords with bdplot_interactive
. Other keywords are:
steps = 10
: How manydt
-steps are taken between producing a new frame.frames = 1000
: How many frames to produce in total.framerate = 60
.
DynamicalBilliards.billiard_bunimovich
— Functionbilliard_bunimovich(l=1.0, w=1.0)
Return a vector of Obstacle
s that define a Buminovich billiard, also called a stadium. The length is considered without the attached semicircles, meaning that the full length of the billiard is l + w
. The left and right edges of the stadium are Semicircle
s.
billiard_stadium
is an alias of billiard_bunimovich
.
DynamicalBilliards.billiard_hexagonal_sinai
— Functionbilliard_hexagonal_sinai(r, R, center = [0,0]; setting = "standard")
Create a sinai-like billiard, which is a hexagon of outer radius R
, containing at its center (given by center
) a disk of radius r
. The setting
keyword is passed to billiard_polygon
.
DynamicalBilliards.billiard_iris
— Functionbilliard_iris(a=0.2, b=0.4, w=1.0; setting = "standard")
Return a billiard that is a square of side w
enclosing at its center an ellipse with semi axes a
, b
.
DynamicalBilliards.billiard_logo
— Methodbilliard_logo(;h=1.0, α=0.8, r=0.18, off=0.25, T = Float64) -> bd, ray
Create the billiard used as logo of DynamicalBilliards
and return it along with the tuple of raysplitters.
DynamicalBilliards.billiard_lorentz
— Functionbilliard_lorentz(r=0.25, x=1.0, y=1.0)
Alias for billiard_sinai(r,x,y; setting = "periodic")
.
DynamicalBilliards.billiard_mushroom
— Functionbilliard_mushroom(sl = 1.0, sw = 0.2, cr = 1.0, sloc = 0.0; door = true)
Create a mushroom billiard with stem length sl
, stem width sw
and cap radius cr
. The center of the cap (which is Semicircle) is always at [0, sl]
. The center of the stem is located at sloc
.
Optionally, the bottom-most Wall
is a Door
(see escapetime
).
DynamicalBilliards.billiard_polygon
— Functionbilliard_polygon(n::Int, R, center = [0,0]; setting = "standard")
Return a vector of obstacles that defines a regular-polygonal billiard with n
sides, radius r
and given center
.
Note: R
denotes the so-called outer radius, not the inner one.
Settings
- "standard" : Specular reflection occurs during collision.
- "periodic" : The walls are
PeriodicWall
type, enforcing periodicity at the boundaries. Only available forn=4
orn=6
. - "random" : The velocity is randomized upon collision.
DynamicalBilliards.billiard_raysplitting_showcase
— Functionbilliard_raysplitting_showcase(x=2.0, y=1.0, r1=0.3, r2=0.2) -> bd, rayspl
Showcase example billiard for ray-splitting processes. A rectangle (x,y)
with a SplitterWall at x/2
and two disks at each side, with respective radii r1
, r2
.
Notice: This function returns a billiard bd
as well as a rayspl
dictionary!
DynamicalBilliards.billiard_rectangle
— Functionbilliard_rectangle(x=1.0, y=1.0; setting = "standard")
Return a vector of obstacles that defines a rectangle billiard of size (x
, y
).
Settings
- "standard" : Specular reflection occurs during collision.
- "periodic" : The walls are
PeriodicWall
type, enforcing periodicity at the boundaries - "random" : The velocity is randomized upon collision.
- "ray-splitting" : All obstacles in the billiard allow for ray-splitting.
DynamicalBilliards.billiard_sinai
— Functionbilliard_sinai(r=0.25, x=1.0, y=1.0; setting = "standard")
Return a vector of obstacles that defines a Sinai billiard of size (x
, y
) with a disk in its center, of radius r
.
In the periodic case, the system is also known as "Lorentz Gas".
Settings
- "standard" : Specular reflection occurs during collision.
- "periodic" : The walls are
PeriodicWall
type, enforcing periodicity at the boundaries - "random" : The velocity is randomized upon collision.
- "ray-splitting" : All obstacles in the billiard allow for ray-splitting.
DynamicalBilliards.billiard_vertices
— Functionbilliard_vertices(v, type = FiniteWall)
Construct a polygon billiard that connects the given vertices v
(vector of 2-vectors). The vertices should construct a billiard in a counter-clockwise orientation (i.e. the normal vector always points to the left of v[i+1] - v[i]
.).
type
decides what kind of walls to use. Ths function assumes that the first entry of v
should be connected with the last.
DynamicalBilliards.bounce!
— Methodbounce!(p::AbstractParticle, bd::Billiard) → i, t, pos, vel
"Bounce" the particle (advance for one collision) in the billiard. Takes care of finite-precision issues.
Return:
- index of the obstacle that the particle just collided with
- the time from the previous collision until the current collision
t
- position and velocity of the particle at the current collision (after the collision has been resolved!). The position is given in the unit cell of periodic billiards. Do
pos += p.current_cell
for the position in real space.
bounce!(p, bd, raysplit) → i, t, pos, vel
Ray-splitting version of bounce!
.
DynamicalBilliards.boundarymap
— Methodboundarymap(p, bd, t [,intervals]) → bmap, arclengths
Compute the boundary map of the particle p
in the billiard bd
by evolving the particle for total amount t
(either float for time or integer for collision number).
Return a vector of 2-vectors bmap
and also arclengths(bd)
. The first entry of each element of bmap
is the arc-coordinate at collisions $\xi$, while the second is the sine of incidence angle $\sin(\phi_n)$.
The measurement direction of the arclengths of the individual obstacles is dictated by their order in bd
. The sine of the angle is computed after specular reflection has taken place.
The returned values of this function can be used in conjuction with the function bdplot_boundarymap
(requires using PyPlot
) to plot the boundary map in an intuitive way.
Notice - this function only works for normal specular reflection. Random reflections or ray-splitting will give unexpected results.
See also to_bcoords
, boundarymap_portion
. See parallelize
for a parallelized version.
DynamicalBilliards.boundarymap_portion
— Methodboundarymap_portion(bd::Billiard, t, p::AbstractParticle, δξ, δφ = δξ)
Calculate the portion of the boundary map of the billiard bd
covered by the particle p
when it is evolved for time t
(float or integer). Notice that the
The boundary map is partitioned into boxes of size (δξ, δφ)
and as the particle evolves visited boxes are counted. The returned ratio is this count divided by the total boxes of size (δξ, δφ)
needed to cover the boundary map.
Important: This portion does not equate the portion the particle's orbit covers on the full, three dimensional phase space. Use the function phasespace_portion
for that!
DynamicalBilliards.cellsize
— Methodcellsize(bd)
Return the delimiters xmin, ymin, xmax, ymax
of the given obstacle/billiard.
Used in randominside()
, error checking and plotting.
DynamicalBilliards.collision
— Methodcollision(p::AbstractParticle, o::Obstacle) → t, cp
Find the collision (if any) between given particle and obstacle. Return the time until collision and the estimated collision point cp
.
Return Inf, SV(0, 0)
if the collision is not possible or if the collision happens backwards in time.
It is the duty of collision
to avoid incorrect collisions when the particle is on top of the obstacle (or very close).
DynamicalBilliards.default_normal
— Methoddefault_normal(sp, ep)
Return a vector to v = ep - sp
, pointing to the left of v
.
DynamicalBilliards.distance
— Methoddistance(p::AbstractParticle, o::Obstacle)
Return the signed distance between particle p
and obstacle o
, based on p.pos
. Positive distance corresponds to the particle being on the allowed region of the Obstacle
. E.g. for a Disk
, the distance is positive when the particle is outside of the disk, negative otherwise.
distance(p::AbstractParticle, bd::Billiard)
Return minimum distance(p, obst)
for all obst
in bd
. If the distance(p, bd)
is negative and bd
is convex, this means that the particle is outside the billiard.
WARNING : distance(p, bd)
may give negative values for non-convex billiards, or billiards that are composed of several connected sub-billiards.
All distance
functions can also be given a position (vector) instead of a particle.
DynamicalBilliards.ellipse_arclength
— Methodellipse_arclength(θ, e::Ellipse)
Return the arclength of the ellipse that spans angle θ
(in normal coordinates, not in the ellipse parameterization). Expects θ
to be in [0, 2π]
.
After properly calculating the
\[d=b\,E\bigl(\tan^{-1}(a/b\,\tan(\theta))\,\big|\,1-(a/b)^2\bigr)\]
DynamicalBilliards.escapetime
— Methodescapetime([p,] bd, t; warning = false)
Calculate the escape time of a particle p
in the billiard bd
, which is the time until colliding with any "door" in bd
. As a "door" is considered any FiniteWall
with field isdoor = true
.
If the particle evolves for more than t
(integer or float) without colliding with the Door
(i.e. escaping) the returned result is Inf
.
A warning can be thrown if the result is Inf
. Enable this using the keyword warning = true
.
If a particle is not given, a random one is picked through randominside
. See parallelize
for a parallelized version.
DynamicalBilliards.evolve!
— Methodevolve!([p::AbstractParticle,] bd::Billiard, t)
Evolve the given particle p
inside the billiard bd
. If t
is of type AbstractFloat
, evolve for as much time as t
. If however t
is of type Int
, evolve for as many collisions as t
. Return the states of the particle between collisions.
This function mutates the particle, use evolve
otherwise. If a particle is not given, a random one is picked through randominside
.
Return
ct::Vector{T}
: Collision times.poss::Vector{SVector{2,T}}
: Positions at the collisions.vels::Vector{SVector{2,T}})
: Velocities exactly after the collisions.ω
, eitherT
orVector{T}
: Angular velocity/ies (returned only for magnetic particles).
The time ct[i+1]
is the time necessary to reach state poss[i+1], vels[i+1]
starting from the state poss[i], vels[i]
. That is why ct[1]
is always 0 since poss[1], vels[1]
are the initial conditions. The angular velocity ω[i]
is the one the particle has while propagating from state poss[i], vels[i]
to i+1
.
Notice that at any point, the velocity vector vels[i]
is the one obdained after the specular reflection of the i-1
th collision.
Ray-splitting billiards
evolve!(p, bd, t, raysplitters)
To implement ray-splitting, the evolve!
function is supplemented with a fourth argument, raysplitters
which is a tuple of RaySplitter
instances. Notice that evolve
always mutates the billiard if ray-splitting is used! For more information and instructions on using ray-splitting please visit the official documentation.
DynamicalBilliards.evolve
— Methodevolve(p, args...)
Same as evolve!
but copies the particle instead.
DynamicalBilliards.extrapolate
— Methodextrapolate(particle, prevpos, prevvel, ct, dt[, ω]) → x, y, vx, vy, t
Create the timeseries that connect a particle
's previous position and velocity prevpos, prevvel
with the particle
's current position and velocity, provided that the collision time between previous and current state is ct
.
dt
is the sampling time and if the particle
is MagneticParticle
then you should provide ω
, the angular velocity that governed the free flight.
Here is how this function is used (for example)
prevpos, prevvel = p.pos + p.current_cell, p.vel
for _ in 1:5
i, ct, pos, vel = bounce!(p, bd)
x, y, x, vy, t = extrapolate(p, prevpos, prevvel, ct, 0.1)
# append x, y, ... to other vectors or do whatever with them
prevpos, prevvel = p.pos + p.current_cell, p.vel
end
DynamicalBilliards.find_cyclotron
— Methodfind_cyclotron(p::MagneticParticle)
Return the center of cyclotron motion of the particle.
DynamicalBilliards.from_bcoords
— Methodfrom_bcoords(ξ, sφ, bd::Billiard, intervals = arcintervals(bd))
Same as above, but now ξ
is considered to be the global arclength, parameterizing the entire billiard, instead of a single obstacle.
DynamicalBilliards.from_bcoords
— Methodfrom_bcoords(ξ, sφ, o::Obstacle) -> pos, vel
Convert the boundary coordinates ξ, sφ
on the obstacle to real coordinates pos, vel
.
Note that vel
always points away from the obstacle.
This function is the inverse of to_bcoords
.
DynamicalBilliards.isphysical
— Methodisphysical(raysplitter(s))
Return true
if the given raysplitters
have physically plausible properties.
Specifically, check if (φ is the incidence angle, θ the refraction angle):
- Critical angle means total reflection: If θ(φ) ≥ π/2 then Tr(φ) = 0
- Transmission probability is even function: Tr(φ) ≈ Tr(-φ) at ω = 0
- Refraction angle is odd function: θ(φ) ≈ -θ(-φ) at ω = 0
- Ray reversal is true: θ(θ(φ, pflag, ω), !pflag, ω) ≈ φ
- Magnetic conservation is true: (ωnew(ωnew(ω, pflag), !pflag) ≈ ω
DynamicalBilliards.ispinned
— Methodispinned(p::MagneticParticle, bd::Billiard)
Return true
if the particle is pinned with respect to the billiard. Pinned particles either have no valid collisions (go in circles forever) or all their valid collisions are with periodic walls, which again means that they go in cirles for ever.
DynamicalBilliards.law_of_refraction
— Functionlaw_of_refraction(n1, n2 = 1.0) -> t, r
Create transmission and refraction functions t, r
that follow Snell's law, i.e. the transmission probability is set to 1.0 except for the case of total internal reflection.
n1
is the index of refraction for the pflag = false
side of an obstacle, while n2
is the index of refraction for pflag = true
.
DynamicalBilliards.lyapunovspectrum
— Methodlyapunovspectrum([p::AbstractParticle,] bd::Billiard, t)
Returns the finite time lyapunov exponents (averaged over time t
) for a given particle in a billiard table using the method outlined in [1]. t
can be either Float
or Int
, meaning total time or total amount of collisions.
Returns zeros for pinned particles.
If a particle is not given, a random one is picked through randominside
. See parallelize
for a parallelized version.
[1] : Ch. Dellago et al, Phys. Rev. E 53 (1996)
DynamicalBilliards.meancollisiontime
— Methodmeancollisiontime([p,] bd, t) → κ
Compute the mean collision time κ
of the particle p
in the billiard bd
by evolving for total amount t
(either float for time or integer for collision number).
Collision times are counted only between obstacles that are not PeriodicWall
.
If a particle is not given, a random one is picked through randominside
. See parallelize
for a parallelized version.
DynamicalBilliards.next_collision
— Functionnext_collision(p::AbstractParticle, bd::Billiard) -> i, tmin, cp
Compute the collision
across all obstacles in bd
and find the minimum one. Return the index of colliding obstacle, the time and the collision point.
DynamicalBilliards.normalvec
— Methodnormalvec(obst::Obstacle, position)
Return the vector normal to the obstacle's boundary at the given position (which is assumed to be very close to the obstacle's boundary).
DynamicalBilliards.parallelize
— Methodparallelize(f, bd::Billiard, t, particles; partype = :threads)
Parallelize function f
across the available particles. The parallelization type can be :threads
or :pmap
, which use threads or a worker pool initialized with addprocs
before using DynamicalBilliards
.
particles
can be:
- A
Vector
of particles. - An integer
n
optionally followed by an angular velocityω
. This usesrandominside
.
The functions usable here are:
meancollisiontime
escapetime
lyapunovspectrum
(returns only the maximal exponents)boundarymap
(returns vector of vectors of 2-vectors andarcintervals
)
DynamicalBilliards.particlebeam
— Functionparticlebeam(x0, y0, φ, N, dx, ω = nothing, T = eltype(x0)) → ps
Make N
particles, all with direction φ
, starting at x0, y0
. The particles don't all have the same position, but are instead spread by up to dx
in the direction normal to φ
.
The particle element type is T
.
DynamicalBilliards.periodicity!
— Methodperiodicity!(p::AbstractParticle, w::PeriodicWall)
Perform periodicity conditions of w
on p
.
DynamicalBilliards.perturbationgrowth
— Methodperturbationgrowth([p,] bd, t) -> ts, Rs, is
Calculate the evolution of the perturbation vector Δ
along the trajectory of p
in bd
for total time t
. Δ
is initialised as [1,1,1,1]
.
If a particle is not given, a random one is picked through randominside
. Returns empty lists for pinned particles.
Description
This function safely computes the time evolution of a perturbation vector using the linearized dynamics of the system, as outlined by [1]. Because the dynamics are linear, we can safely re-normalize the perturbation vector after every collision (otherwise the perturbations grow to infinity).
Immediately before and after every collison, this function computes
- the current time.
- the element-wise ratio of Δ with its previous value
- the obstacle index of the current obstacle
and returns these in three vectors ts, Rs, is
.
To obtain the actual evolution of the perturbation vector you can use the function perturbationevolution(Rs)
which simply does
Δ = Vector{SVector{4,Float64}}(undef, length(R))
Δ[1] = R[1]
for i in 2:length(R)
Δ[i] = R[i] .* Δ[i-1]
end
[1] : Ch. Dellago et al, Phys. Rev. E 53 (1996)
DynamicalBilliards.phasespace_portion
— Methodphasespace_portion(bd::Billiard, t, p::AbstractParticle, δξ, δφ = δξ)
Calculate the portion of the phase space of the billiard bd
covered by the particle p
when it is evolved for time t
(float or integer).
This function extends boundarymap_portion
using a novel approach. For each visited box of the boundary map, bounce!
attributes a third dimension (the collision time, equal to collision distance) which expands the two dimensions of the boundary map to the three dimensions of the phase space.
The true phase space portion is then the weighted portion of boxes visited by the particle, divided by the total weighted sum of boxes. The weights of the boxes are the collision times.
DynamicalBilliards.polygon_vertices
— Functionpolygon_vertices(r, sides, center = [0, 0], θ=0) -> v
Return the vertices that define a regular polygon with sides
sides and radius r
, centered at center
.
Optionally rotate by θ
degrees.
DynamicalBilliards.project_to_line
— Methodproject_to_line(point, c, n)
Project given point
to line that contains point c
and has normal vector n
.
DynamicalBilliards.propagate!
— Methodpropagate!(p::AbstractParticle, t)
Propagate the particle p
for given time t
, changing appropriately the the p.pos
and p.vel
fields.
propagate!(p, position, t)
Do the same, but take advantage of the already calculated position
that the particle should end up at.
DynamicalBilliards.propagate_offset!
— Methodpropagate_offset!(offset::MArray{Tuple{4,4},T}, p::AbstractParticle)
Computes the linearized evolution of the offset vectors during propagation for a time interval t
DynamicalBilliards.psos
— Methodpsos(bd::Billiard, plane::InfiniteWall, t, particles)
Compute the Poincaré section of the particles
with the given plane
, by evolving each one for time t
(either integer or float) inside bd
.
The plane
can be an InfiniteWall
of any orientation, however only crossings of the plane
such that dot(velocity, normal) < 0
are allowed, with normal
the normal unit vector of the plane
.
particles
can be:
- A single particle.
- A
Vector
of particles. - An integer
n
optionally followed by an angular velocityω
.
Return the positions poss
and velocities vels
at the instances of crossing the plane
. If given more than one particle, the result is a vector of vectors of vectors.
Notice - This function can handle pinned particles. If a pinned particle can intersect with the plane
, then an intersection is returned. If however it can't then empty vectors are returned.
DynamicalBilliards.randominside
— Functionrandominside(bd::Billiard [, ω]) → p
Return a particle p
with random allowed initial conditions inside the given billiard. If supplied with a second argument the type of the returned particle is MagneticParticle
, with angular velocity ω
.
WARNING : randominside
works for any convex billiard but it does not work for non-convex billiards. (this is because it uses distance
)
DynamicalBilliards.randominside_xyφ
— Methodrandominside_xyφ(bd::Billiard) → x, y, φ
Same as randominside
but only returns position and direction.
DynamicalBilliards.raysplit_indices
— Methodraysplit_indices(bd::Billiard, raysplitters::Tuple)
Create a vector of integers. The i
th entry tells you which entry of the raysplitters
tuple is associated with the i
th obstacle of the billiard.
If the i
th entry is 0
, this means that the obstacle does not do raysplitting.
DynamicalBilliards.real_pos
— Methodreal_pos(ξ, o::Obstacle)
Converts the arclength coordinate ξ
relative to the obstacle o
into a real space position vector.
DynamicalBilliards.realangle
— Methodrealangle(p::MagneticParticle, o::Obstacle, I) -> θ
Given the intersection point I
of the trajectory of a magnetic particle p
with some obstacle o
, find the real angle that will be spanned until the particle collides with the obstacle.
The function also takes care of problems that may arise when particles are very close to the obstacle's boundaries, due to floating-point precision.
DynamicalBilliards.relocate!
— MethodDynamicalBilliards.reset_billiard!
— Methodreset_billiard!(bd)
Sets the pflag
field of all ray-splitting obstacles of a billiard table to true
.
DynamicalBilliards.resolvecollision!
— Methodresolvecollision!(p::AbstractParticle, o::Obstacle)
Resolve the collision between particle p
and obstacle o
, depending on the type of o
(do specular!
or periodicity!
).
resolvecollision!(p, o, T::Function, θ::Function, new_ω::Function)
This is the ray-splitting implementation. The three functions given are drawn from the ray-splitting dictionary that is passed directly to evolve!()
. For a calculated incidence angle φ, if T(φ) > rand(), ray-splitting occurs.
DynamicalBilliards.specular!
— Methodspecular!(p::AbstractParticle, o::Obstacle)
Perform specular reflection based on the normal vector of the Obstacle.
In the case where the given obstacle is a RandomObstacle
, the specular reflection randomizes the velocity instead (within -π/2+ε to π/2-ε of the normal vector).
DynamicalBilliards.timeseries!
— Methodtimeseries!([p::AbstractParticle,] bd::Billiard, t; dt, warning)
Evolve the given particle p
inside the billiard bd
for the condition t
and return the x, y, vx, vy timeseries and the time vector. If t
is of type AbstractFloat
, then evolve for as much time as t
. If however t
is of type Int
, evolve for as many collisions as t
. Otherwise, t
can be any function, that takes as an input t(n, τ, i, p)
and returns true
when the evolution should terminate. Here n
is the amount of obstacles collided with so far, τ
the amount time evolved so far, i
the obstacle just collided with and p
the particle (so you can access e.g. p.pos
).
This function mutates the particle, use timeseries
otherwise. If a particle is not given, a random one is picked through randominside
.
The keyword argument dt
is the time step used for interpolating the time series in between collisions. dt
is capped by the collision time, as the interpolation always stops at collisions. For straight propagation dt = Inf
, while for magnetic dt = 0.01
.
For pinned magnetic particles, timeseries!
issues a warning and returns the trajectory of the particle. If t
is integer, the trajectory is evolved for one full circle only.
Internally uses DynamicalBilliards.extrapolate
.
Ray-splitting billiards
timeseries!(p, bd, t, raysplitters; ...)
To implement ray-splitting, the timeseries!
function is supplemented with a fourth argument, raysplitters
which is a tuple of RaySplitter
instances. Notice that timeseries
always mutates the billiard if ray-splitting is used! For more information and instructions on using ray-splitting please visit the official documentation.
DynamicalBilliards.timeseries
— Methodtimeseries(p, args...; kwargs...)
Non-mutating version of DynamicalBilliards.timeseries!
DynamicalBilliards.to_bcoords
— Methodto_bcoords(pos, vel, o::Obstacle) -> ξ, sφ
Convert the real coordinates pos, vel
to boundary coordinates (also known as Birkhoff coordinates) ξ, sφ
, assuming that pos
is on the obstacle.
ξ
is the arc-coordinate, i.e. it parameterizes the arclength of the obstacle. sφ
is the sine of the angle between the velocity vector and the vector normal to the obstacle.
The arc-coordinate ξ
is measured as:
- the distance from start point to end point in
Wall
s - the arc length measured counterclockwise from the open face in
Semicircle
s - the arc length measured counterclockwise from the rightmost point in
Circular
/Ellipse
s
Notice that this function returns the local arclength. To get the global arclength parameterizing an entire billiard, simply do ξ += arcintervals(bd)[i]
if the index of obstacle o
is i
.
See also from_bcoords
, which is the inverse function.
DynamicalBilliards.totallength
— Methodtotallength(o::Obstacle)
Return the total boundary length of o
.
DynamicalBilliards.translate
— Functiontranslate(obst::Obstacle, vector)
Create a copy of the given obstacle with its position translated by vector
.
DynamicalBilliards.visited_obstacles!
— Methodvisited_obstacles!([p::AbstractParticle,] bd::Billiard, t)
Evolve the given particle p
inside the billiard bd
exactly like evolve!
. However return only:
ts::Vector{T}
: Vector of time points of when each collision occured.obst::Vector{Int}
: Vector of obstacle indices inbd
that the particle collided with at the time points ints
.
The first entries are 0.0
and 0
. Similarly with evolve!
the function does not record collisions with periodic walls.
Currently does not support raysplitting. Returns empty arrays for pinned particles.
DynamicalBilliards.visited_obstacles
— Methodvisited_obstacles(p, args...)
Same as visited_obstacles!
but copies the particle instead.