Core design

The BMO package is intended to provide optical simulation capabilites with as much "out-of-the-box" comfort as possible, meaning that users should not have to worry about e.g. providing the correct optical sequence and exact alignment of objects. The following five design principles are core assumptions of the underlying API:

1. Optical interactions are decoupled from the underlying geometry representation

2. Optical elements are closed volumes or must mimic as such (exceptions apply, e.g. coatings)

3. Elements should be easily moveable and have working interactions for most angles of incidence

4. Without additional knowledge, tracing is performed non-sequentially

5. With additional knowledge, tracing is performed sequentially

Intersect-Interact-Repeat-Loop

The first two principles will be elaborated upon in more detail in the Geometry representation section. For the latter two design decisions, the following high-level solver schematic can be used to abstract the steps that are performed when calling solve_system! with an input system and beam:

This scheme is loosely referred to as the Intersect-Interact-Repeat-Loop and consists of the following steps:

1. Calculate the closest Intersection between a ray/beam and the objects within the system

2. Calculate the optical Interaction that occurs at the surface or within the volume of the element

3. Attach the next part of the ray chain or beam tree

4. Use the new information to repeat 1.

Once this procedure has been completed, the alignment of the system or other time-dependent optical properties (e.g. the phase of a Gaussian beamlet) can be updated. A new call of solve_system! solves the beam again from its start. This is described in more detail in the section: Tracing systems.

The next sections will focus on the Intersection and Interaction steps.

Intersections

Calculating intersections between straight lines, i.e. rays, and surfaces is a central challenge for every geometrical optics simulator. This must be done with high numerical precision, since many optical effects are sensitive on the order of the wavelength of the light under consideration with respect to position and direction [26, p. 265 ff]. In order to define this mathematically or algorithmically, many different methods exist [27]. The first question is, how is the geometry of the problem defined. This topic is treated in the Geometry representation section. The second question concerns then the algorithm or equation that allows to calculate the point of intersection between a ray and the surface of the element. This function is called intersect3d and is, at its core, defined for each shape and ray:

BeamletOptics.intersect3d Method

intersect3d(shape::AbstractShape, ::AbstractRay)

Defines the intersection between an AbstractShape and an AbstractRay, must return an Intersection or nothing. The default behavior for concrete shapes and rays is to indicate no intersection, that is nothing, which will inform the tracing algorithm to stop. Refer to the Intersection documentation for more information on the return type value.

source

Regardless of the underlying concrete implementation, each call of intersect3d must return nothing or the following type:

BeamletOptics.Intersection Type

Intersection{T}

Stores data calculated by the intersect3d method.

Fields:

• object: a Nullable reference to the AbstractObject that has been hit (optional but recommended)

• shape: a Nullable reference to the AbstractShape of the object that has been hit (optional but recommended)

• t: length of the ray parametrization in [m]

• n: normal vector at the point of intersection

source

Since an optical element can consist of multiple joint shapes, the return type must store which specific part of the object was hit.

Bounding sphere pretest

Since the tracing is non-sequential, every ray is tested against every shape of the system, and most of these tests are misses. A shape can therefore provide a sphere that encloses it. A ray that does not pass through this sphere is a miss without a call of intersect3d. This is optional: a shape without a sphere is tested as described above.

BeamletOptics.bounding_sphere_of Method

bounding_sphere_of(shape::AbstractShape)
bounding_sphere_of(object::AbstractObject)
bounding_sphere_of(group::AbstractObjectGroup)

Computes a sphere that encloses the shape, object or group in its current pose and returns it as a SingleBoundingSphere or MultiBoundingSphere, or returns NoBoundingSphere() if there is none.

For a shape, the method is an optional part of the AbstractShape interface. With the default NoBoundingSphere() every ray is tested via intersect3d. With a sphere, a solve tests a ray against the sphere first and skips the shape if it misses, which pays off for shapes with a costly intersect3d. The sphere is computed once per solve_system! and not stored in the shape, hence the method may take some time, e.g. look at all vertices of a mesh.

The sphere must enclose every point of the shape, for all its parameters: a hit outside of the sphere is lost without a warning. It should also be tight, since a ray that hits the sphere is tested as before. A safety margin is not needed, the solver adds its own tolerance. Use render_bounding_sphere! to look at the result.

# a cylinder of the radius `r` and the height `h`, with the origin at the center of its base
function BeamletOptics.bounding_sphere_of(c::MyCylinder)
    return BeamletOptics.SingleBoundingSphere(c, Point3(0, c.h / 2, 0), sqrt(c.r^2 + (c.h / 2)^2))
end

For an object, the result follows from its shapes, see AbstractShapeTrait: the sphere of the shape of a SingleShape object, and for a MultiShape object or an object group the MultiBoundingSphere that encloses the spheres of all its parts, or NoBoundingSphere() if one part has none.

A solve does not call this method per ray. It computes the spheres once into a table and looks them up via bounding_sphere_of(table, x).

source

All Signed Distance Functions (SDFs) and Meshes have a bounding sphere. An object of several shapes and an object group have the sphere around the spheres of their parts, such that a ray that misses it is tested against none of them.

The spheres are not stored in the shapes or objects. solve_system! computes them once per call, for the poses that the objects have at that time, via bounding_sphere_of(x) into a table, which the intersection code reads per ray via bounding_sphere_of(table, x). Outside of a solve there is no table, i.e. a direct call of intersect3d(object, ray) tests the object without a sphere.

BeamletOptics.bounding_sphere_of Method

bounding_sphere_of(table::BoundingSphereTable, x)
bounding_sphere_of(::Nothing, x)

Returns the bounding sphere that the table holds for the shape or object x, or NoBoundingSphere() if it holds none. Nothing is computed, in contrast to bounding_sphere_of(x), which fills the table, see bounding_spheres.

The intersection code calls this method per ray with the table of the solve that is running in the current task, current_bounding_spheres(). Outside of a solve this table is nothing and the result is always NoBoundingSphere(), i.e. intersect3d(object, ray) then tests the object without a sphere.

source

BeamletOptics.SingleBoundingSphere Type

SingleBoundingSphere{T} <: AbstractBoundingSphere

A sphere that encloses one shape, see bounding_sphere_of.

Fields

• pos: the center in world coordinates [m]

• radius: the radius [m]

source

BeamletOptics.MultiBoundingSphere Type

MultiBoundingSphere{T} <: AbstractBoundingSphere

A sphere that encloses several bounding spheres, i.e. all parts of a MultiShape object or all objects of an object group, see bounding_sphere_of. A ray that misses it hits none of the parts, hence a solve tests it before the parts. It is tested against a ray like a SingleBoundingSphere and differs only in what it states: the sphere of several things.

Fields

• pos: the center in world coordinates [m]

• radius: the radius [m]

source

BeamletOptics.NoBoundingSphere Type

NoBoundingSphere <: AbstractBoundingSphere

States that a shape or object has no bounding sphere. It is then tested against every ray via intersect3d, see bounding_sphere_of.

It has no center and no radius: position and radius throw an ArgumentError. Code that needs them handles NoBoundingSphere by a method of its own, or tests for it first.

source

The spheres of a system can be shown via render_bounding_sphere! or the show_bounding_sphere keyword of render!.

Interactions

Optical interactions are performed after the point of intersection has been determined. The interact3d interface allows users to implement algorithms that calculate or try to mimic optical effects. The fidelity of the algorithm is effectively only limited by the amount of information that can be passed into the interact3d interface. The method is defined as follows:

BeamletOptics.interact3d Method

interact3d(::AbstractSystem, object::AbstractObject, ::AbstractBeam, ::AbstractRay)

Defines the optical interaction between an incoming/outgoing beam/ray of light and an optical element, must return an AbstractInteraction or nothing. The default behavior is that no interaction occurs, i.e. return of nothing, which should stop the system tracing procedure. Refer to the AbstractInteraction typedocs for more information on the return type value.

A MultiShape object receives this call for whichever of its parts was hit. Its method finds that part, e.g. via shape(intersection(ray)), and hands the beam to it, see CubeBeamsplitter. Parts that do not interact as one belong in an ObjectGroup or another AbstractObjectGroup instead, where each part interacts on its own.

source

As with the intersect3d method, a predefined return type must be provided in order to make the solve_system! interface work. The AbstractInteraction is used in order to create "building blocks" from which the output beam is constructed.

BeamletOptics.AbstractInteraction Type

AbstractInteraction

Describes how an AbstractBeam and an AbstractObject interact with each other. This data type stores information from the interact3d function and provides it to the solver. The solver can use this data to extend the AbstractBeam.

Implementation reqs.

Subtypes of AbstractInteraction must implement the following:

Fields

• hint: a nullable Hint for the solver (optional but recommended)

Beam data

It is required that concrete implementations of this type provide some form of data on how to extend the beam. For instance, refer to theBeamInteraction and GaussianBeamletInteraction.

source

An object whose interact3d method stores data in the object, e.g. a detector its hits, also implements initialize!, which discards this data before a new solve.

The interact3d return type limits the interface to only accepting one new beam segment per interaction at the moment. The developer needs to take into account that after e.g. a lens surface air-to-glass interaction, the solver "forgets" that the next logical step is to immediately test against the lens again, since the most likely step will be the refraction at the glass-to-air surface. In order to alleviate this issue, the Hint type can be used.

Re-emitting components

Some components do not continue a beam where it ended, but hand it to another solver and emit the result elsewhere, e.g. a fiber that takes the light at its input facet and releases it at its output facet. Their interact3d method samples the incoming beamlet segment as a field, lets the other solver compute the field at the exit and starts new beamlets there. The exchange format is the PlaneField of OpticsBase.jl (see OpticsBase interoperability):

function BeamletOptics.interact3d(::AbstractSystem, fiber::MyFiber, agb::AstigmaticGaussianBeamlet, id::Int)
    f = PlaneField(agb, id; size, spacing, origin, axes)      # the segment that hit the fiber, on its input plane
    g = other_solver(f)                                       # PlaneField at the exit, with the full phase
    BeamletOptics.relaunch!(agb, [GaussianModeDecomposition(g)])
    return nothing
end

The new beamlets are attached with BeamletOptics.relaunch! instead of children!. Children added by children! continue the parent: its length, optical path and reference plane are counted on. A re-emitted beamlet already carries the full phase in its amplitude and starts at another place, so it must count its optical path from its own start.

BeamletOptics.relaunch! Function

relaunch!(beam::B, new_beams::AbstractVector{B}) where {B <: AbstractBeam}

Attaches new_beams as children of beam that start their own optical reference: they are traced and removed together with beam, but their length, optical path and phase are counted from their own start instead of being continued from beam. The new beams carry the full phase in their amplitude.

Use it in the interact3d of a component that re-emits light elsewhere, e.g. the output end of a fiber or a solver coupled through a field: the light leaves at another place than where beam ended, so the geometric continuation of children! (parent length, reference plane of a beamlet) does not apply.

The new_beams replace any children that beam already has. A beam that the component keeps and attaches again is solved from its first ray, like every beam.

source

OpticsBase.PlaneField Method

OpticsBase.PlaneField(beam::BeamletOptics.AstigmaticGaussianBeamlet, id::Integer;
    size::NTuple{2,Int}, spacing::NTuple{2,Real}, origin, axes, progress = false)
OpticsBase.PlaneField(beam::BeamletOptics.GaussianBeamlet, id::Integer; kwargs...)

Samples the field of segment id of a traced beamlet beam on a plane, as PlaneField(detector; ...) does for the beamlets that hit a detector. Use it where there is no detector, e.g. in the interact3d of a component that hands the incoming beamlet to another solver: the segment that just hit the component is id.

size and spacing are as for PlaneField(detector; ...); origin (in [m]) and axes (columns u, v, n) place the plane and have no default. The field is that of the segment's straight continuation through the plane (no further interaction), including the optical path from the source (parents included), and follows the conventions of PlaneField(detector; ...): no projection factor, H from the beamlet direction, a scalar GaussianBeamlet polarized along u, R = Inf, n the refractive index of the segment.

An AstigmaticGaussianBeamlet can be sampled on any of its segments, also within a medium. The scalar GaussianBeamlet model can be sampled on any of its segments as well, but in a medium of index n its amplitude is √n too large.

source

Hints

As mentioned in the previous section, the BeamletOptics.Hint interface allows developers to manipulate the non-sequential solver algorithm into testing against a specific component and shape during the next cycle of the Intersect-Interact-Repeat-Loop. This interface has very high priority during intersection testing.

BeamletOptics.Hint Type

Hint

A Hint can be passed as part of an AbstractInteraction and will inform the tracing algorithm about which AbstractObject in the AbstractSystem will be hit next.

Info

The Hint does not need to result in a guaranteed Intersection. However, if the hinted shape is intersected, it will be immediately assumed as the correct global intersection.

Fields

• object: the object that might or will be intersected next

• shape: the underlying shape that will be intersected next, i.e. shape(object), relevant for multi-shape objects

source

The main reason for this is the intersection ambiguity encountered at interfaces between air-tight component interfaces, e.g. Plate beamsplitters or cemented Doublet lenses. This is caused by the fact that for a ray with a starting point at this interface, technically both shapes are being "touched" at the same time. Additional program logic considering the ray direction of propagation can not always resolve this ambiguity. Therefore the task of providing additional information to the solver via the Hint interface is placed as a burden on the developer.

Tip

Use of the Hint interface is primarily intended for multishape objects with joint surfaces.

Tracing logic

Tracing systems

In the initial state, is is assumed that the problem consists of objects <: AbstractObject (in a system) and a beam <: AbstractBeam with a defined starting position and direction. No additional information is provided, and the specific path of the beam is not known beforehand. Consequently, brute force tracing of the optical system is required, involving testing against each individual element to determine the trajectory of the beam.

BeamletOptics.trace_system! Function

trace_system!(system::AbstractSystem, beam::Beam{T}; r_max = get_default_r_max()) where {T <: Real}

Trace a Beam through an optical system. Maximum number of tracing steps can be capped by r_max.

Tracing logic

The intersection of the last ray of the beam with any objects contained within the system is tested. If an object is hit, the optical interaction is calculated. If no interaction occurs or no further objects are hit, the tracing procedure is stopped.

Arguments

• system:: The optical system through which the Beam is traced.

• beam: The Beam object to be traced.

• r_max: Maximum number of tracing iterations.

source

trace_system!(system::System, gauss::GaussianBeamlet{T}; r_max = get_default_r_max()) where {T <: Real}

Trace a GaussianBeamlet through an optical system. Maximum number of tracing steps can be capped by r_max.

Tracing logic

The chief, waist and divergence beams are traced step-by-step through the system. For each intersection after a tracing_step!, the intersections are compared. If all rays hit the same target, the optical interaction is analyzed, else the tracing stops.

Arguments

• system: The optical system through which the GaussianBeamlet is traced.

• gauss: The GaussianBeamlet object to be traced.

• r_max: Maximum number of tracing iterations.

source

trace_system!(system, agb::AstigmaticGaussianBeamlet; r_max = get_default_r_max(), check_invariant = true, threshold = get_invariant_threshold())

Trace an AstigmaticGaussianBeamlet through an optical system. All 9 component beams (chief + 8 parabasal) are traced in lockstep: the chief ray is traced first for each intersection, followed by the auxiliary rays. All rays must hit the same shape; otherwise tracing stops.

source

This non-sequential mode is comparatively safe in determining the "true" beam path, but will scale suboptimally in time-complexity with the amount of optical elements.

Object order

Unlike with classic, surface-based ray tracers, the order in which objects are listed in the System object vector/tuple is not considered for the purpose of tracing.

CPU and GPU support

Parallizing the execution of a solve_system! call on the CPU is straight-forward for systems that do not feature objects which can be mutated during runtime, e.g. detectors like the Detector. For each beam or ray the solution is independent and the solver can run on multiple threads. Special consideration needs to be taken when implementing mutable elements as mentioned above, since multiple threads might be able to access the underlying memory, leading to race conditions. Specifically, this means ensuring e.g. atomic write and read access.

With respect to GPU acceleration, this is not the case. Currently, all available implementations of solve_system! are highly branching algorithms which can not be implemented in a parallized way easily. This will most likley require a specific new subtype of the BeamletOptics.AbstractSystem with determinable sequential properties. This is not a development goal as of the writing of this section.