Ray–Wave Correspondence in Anisotropic Mesoscopic Billiards

TL;DR

The paper explores how ray-wave correspondence applies to anisotropic billiards, using examples like bilayer graphene and optical cavities to bridge dynamics and wave behavior.

Contribution

The paper introduces anisotropic billiards as a new system class and applies ray-wave correspondence to both optical and electronic systems.

Findings

Anisotropic dispersion in bilayer graphene and birefringent optical cavities can be analyzed using ray-tracing with non-spherical index ellipsoids.

Transformation optics reveals that resonances in anisotropic systems match those of isotropic elliptical cavities.

Differences in ray-wave correspondence exist between optical and electronic anisotropic systems.

Abstract

Mesoscopic billiard systems for electrons and light, realized as quantum dots or optical microcavities, have enriched the fields of quantum chaos and nonlinear dynamics not only by enlarging the class of model systems, but also by providing access to their experimental realization. Here, we add yet another system class, two-dimensional billiards with anisotropies. One example is the anisotropic dispersion relation relevant in bilayer graphene known as trigonal warping, and another is the birefringent closed optical disk cavity. We demonstrate that the established concept of ray–wave correspondence also provides useful insight for anisotropic billiard systems. First, we approach the dynamics of the anisotropic disk from the ray-tracing side that takes the anisotropy in momentum space into account, based on the non-spherical index ellipsoid. Second, we use transformation optics to solve the wave problem and find the resonances to be those of the isotropic elliptical cavity. We illustrate ray–wave correspondence and mark differences in the description of optical and electronic anisotropic systems.

Figures

Full text

1. Introduction

Symmetries and their breaking have been an inspiring principle throughout physics, with the development of the standard model of particle physics and the gauge symmetries used to derive quantum field theories as prominent examples [1]. They rely on the Noether theorem [2], in which Emmy Noether related continuous symmetries to conserved quantities—and vice versa—more than 100 years ago. Furthermore, the ongoing progress in material sciences and mesoscopic physics [3] constantly introduces novel physical systems with various properties to explore and manipulate. The fabrication of two-dimensional electron gas enabled the investigation of ballistic electronic quantum dot billiards and boosted the field of quantum chaos [4,5], as did optical microcavities [6,7]. The variety of these systems increases the pool of properties that can be studied. For example, bilayer graphene (BLG) systems can provide an anisotropic dispersion relation in momentum space [8,9], which was found to tremendously affect the internal billiards dynamics, even for circular confinement [10]. Magnetic anisotropies play a role in the dynamics of nanomagnets subject to magnetoelectric interaction [11,12] and straintronics [13]. Anisotropic optical systems and metamaterials are receiving increasing interest [14,15,16,17,18].

The theoretical description of smaller and smaller billiard systems calls for their quantum description, that remains a challenge for mesoscopic systems with typical sizes on the micrometer scale—several orders of magnitude larger than typical quantum systems, yet corresponding to system sizes of only several wavelengths across. Semiclassical methods [4,5,19] are useful and provide deeper insight, with the concept of classical–quantum, or ray–wave, correspondence being particularly straightforward to apply. This fundamental principle provides a qualitative understanding of the dynamics, especially in ballistic mesoscopic billiards [20,21,22,23,24], based on the description of their classical counterparts. Moreover, it allows for a (semi)quantitative prediction of their emission properties not only for optical microresonators [25], but also for electronic systems like single-layer graphene [26].

How does a loss of symmetry in a system caused by the presence of anisotropies change its dynamics? For instance, let us consider whispering gallery (WG) orbits that directly reflect the isotropy of space and, according to Noether’s theorem, the conservation of angular momentum in rotational invariant cavities. How do these trajectories change in the presence of anisotropies, i.e., in the presence of preferred propagation directions that break the rotational invariance of the system.

Here, we want to explore the chances of a semiclassical description of anisotropic mesoscopic billiards, with BLG billiards and the birefringent optical disk as examples, for an electronic and an optical system, respectively. We will first briefly review the particle tracing in BLG systems [10] and outline its generalization to anisotropic optical systems. While a wave description of condensed-matter systems typically requires tight-binding modeling, Maxwell’s equation allows for a description of anisotropic dielectric systems using transformation optics (TO). We illustrate TO for the birefringent disk cavity and are able to establish the expected ray–wave correspondence. However, this does not hold for the BLG case where the presence of two so-called valleys associated with the presence of two (non-identical) Dirac points changes the situation. We end this paper with a summary and a future outlook for research in this field.

2. Trajectory Tracing in Electronic and Optical Anisotropic Billiards

The description of the particle dynamics in generic billiards is a well-established field of nonlinear dynamics [27] and is well understood both in real space—showing the trajectories—and phase space, where, typically, a Poincaré surface of section (PSOS) is taken at the billiard’s boundary to reduce the number of phase-space variables from four (as we consider two-dimensional systems) to two, namely the boundary coordinate s and $[eqn]$ (the sine of the angle of incidence). In isotropic circular billiards, the PSOS consists of points along a line of constant $[eqn]$ , directly expressing the conserved angular momentum.

2.1. Anisotropic Bilayer Graphene Billiards

The bilayer of graphene billiards has to change in the presence of anisotropies, even when the circular billiard geometry is kept: when space is anisotropic, Noether’s theorem [2] states that angular momentum cannot be conserved any more. For BLG, such an anisotropy results from its band structure (obtained from a $[eqn]$ tight-binding Hamiltonian), which is found to not be rotationally invariant around its band minima and maxima (the so-called K points) for certain yet generic parameter choices. The intersection curve of the band at a certain energy (the Fermi energy) defines the dispersion relation or Fermi line, which then differs from a circle and is thus anisotropic. Its shape is known as trigonal warping, a deformed circle where three corner-like features emerge (see Figure 1(a.2,b.2)), and is related by mirror symmetry at the two K points (or valleys), $[eqn]$ and $[eqn]$ , that result from pseudospin (the presence of two atoms per unit cell). For more details on graphene physics, we would like to refer the interested reader to Refs. [8,9,10], and leave others with the information that trigonal warping of the energy bands in BLG is an expression of anisotropy.

The PSOS for such an anisotropic BLG billiard is shown in Figure 1. This example shows the PSOS for a BLG cavity with a onigiri geometry, which somewhat resembles the shape of the BLG Fermi line and is given in polar coordinates as $[eqn]$ ). Here, $[eqn]$ is the geometric shape parameter and $[eqn]$ gives the angle of rotation between the underlying BLG lattice and the billiard’s axial symmetry axis.

First of all, we have to notice that the PSOS is not spanned by s and $[eqn]$ , but rather, by s and $[eqn]$ (with $[eqn]$ > 0 indicating a mathematically positive sense of rotation): as the angle of incidence cannot be equal to the angle of reflection in general anymore, we have to find another quantity that is conserved upon reflection at the boundary. Indeed, the momentum component parallel to the interface is conserved upon reflection and is therefore a suitable phase space variable. The curly upper and lower edges of the PSOS are caused by the non-circular shape of the Fermi line. The maximum possible value of the $[eqn]$ component changes depending on the orientation of the billiard boundary to the Fermi line (which is fixed to the BLG lattice). It is striking to see in Figure 1(a.1,b.2) that the PSOS is not symmetric with respect to $[eqn]$ anymore (for a single K point), in contrast to what applies for isotropic billiards. At this point, we have to note that the PSOS corresponding to the two K points (with Fermi lines related by mirror symmetry) are not identical (Figure 1(a.1,b.1)). Notice that symmetry with respect to $[eqn]$ = 0 is restored when both K points are considered.

The PSOS in Figure 1 contains intensity information (the reflection coefficient is derived from the BLG Hamilton operator and taken into account at each reflection [10]), with light colors indicating high intensity. It is striking that there are two chains of three islands that correspond to triangular orbits built from the three preferred propagation directions present in a trigonally warped Fermi line. Their shapes illustrate that the Fermi lines and PSOS corresponding to the two different K points are related by mirror symmetry. In turn, WG orbits are lost as they turn into chaotic trajectories.

2.2. Anisotropy and Polarization

We now turn to optical systems that are also governed by wave equations. An optical or electromagnetic wave with frequency $[eqn]$ is described by Maxwell’s equations. It is given by its wave vector $[eqn]$ , the electric field $[eqn]$ , the electric displacement field $[eqn]$ , the magnetic field strength $[eqn]$ , and the magnetic induction $[eqn]$ with the vacuum permittivity (permeability) $[eqn]$ ( $[eqn]$ ), the permittivity tensor $[eqn]$ , and the permeability tensor $[eqn]$ . The Poynting vector $[eqn]$ describes the flow direction of the energy flux. The Maxwell’s equations in anisotropic materials that we are interested in here connect the following quantities:

$$[eqn]$$ $$[eqn]$$

and we will use them whenever we argue in terms of the wave picture.

For now, we want to generalize the ray-tracing algorithm developed for BLG billiards [10] to anisotropic optical systems. To this end we have to use insight from the wave picture as an ingredient for the ray description, very similar to what we discussed for BLG. We consider a uniaxial crystal such as calcite where the non-cubic lattice or a layered structure is the origin of an anisotropy known as birefringence. We assume a refractive index $[eqn]$ in one direction, and indices $[eqn]$ along the remaining two Cartesian coordinates. An example of the resulting index ellipsoid is shown in Figure 2a for a generic wave vector $[eqn]$ . According to Maxwell’s equations, the $[eqn]$ field oscillates in a plane perpendicular to $[eqn]$ , indicated in red. The magnetic field $[eqn]$ is perpendicular to the plane spanned by $[eqn]$ and $[eqn]$ [29,30]. The boundary of the intersectional plane, which results from the intersection between the $[eqn]$ - $[eqn]$ plane and the ellipsoid, forms an ellipse, and its two semi-axes define the refractive indices for the two polarizations.

To explain the interplay between anisotropy in a two-dimensional optical system and the polarization state of light, let us consider a simpler case where the resonator plane lies in the ( $[eqn]$ ) plane. We then refer to polarization where the electric (magnetic) field $[eqn]$ ( $[eqn]$ ) oscillates in the z direction as TM (TE) polarization while the magnetic (electric) field is transverse to the z axis and lies in the resonator plane. If now, the lattice direction with the (larger) refractive index $[eqn]$ is aligned along the z axis (Figure 2b,c), it follows from Maxwell’s equations that the TM modes behave like those of an isotropic disk with refractive index $[eqn]$ . The TE modes are those with $[eqn]$ and thus remain unchanged compared to the isotropic case; they are referred to as ordinary rays. For these, the magnetic field oscillates in the resonator plane and does not experience any anisotropy. This is in agreement with the argumentation based on the index ellipsoid illustrated in Figure 2b,c. Here, all intersectional planes perpendicular to $[eqn]$ and through the center of the ellipsoid possess the same elliptical cross section. For TM polarization, the electric field $[eqn]$ is aligned with the z axis and experiences a refractive index equal to $[eqn]$ , the semi-axis in the z direction of the intersectional ellipse. In contrast, for TE polarization, the electric field $[eqn]$ lies in the $[eqn]$ plane, where it will experience the refractive index $[eqn]$ for any $[eqn]$ .

The more interesting case is the one where the larger refractive index $[eqn]$ is, without loss of generality, pointing in the y direction (Figure 2d,e). Again, there are no changes for TM polarization that behaves as in the isotropic case, but now corresponding to a refractive index $[eqn]$ , where the anisotropy in the resonator ( $[eqn]$ ) plane is not accessible to an electric field oscillating perpendicular to it. In the line of reasoning based on the index ellipsoid, we argue that the intersectional plane will always possess a semi-axis of size $[eqn]$ in the z direction. The one case remaining, a TE polarization when $[eqn]$ lies in the resonator plane, is the most interesting one, with qualitatively new behavior that differs from the isotropic disk. We will focus on this case in the following section.

2.3. From the Index Ellipsoid to the Dispersion Relation in Uniaxial Optical Billiards

The dispersion relation, i.e., the relation between the frequency $[eqn]$ of light and its wave vector $[eqn]$ , that we need to generalize the trajectory tracing algorithm developed for BLG billiards to the optical case can be directly deduced from the index ellipsoid. We recall the generic equation for an ellipsoid in $[eqn]$ space with semi-axis $[eqn]$ and $[eqn]$ to be $[eqn]$ . If the refractive index in the y direction is largest, we thus obtain an index ellipsoid elongated in the y direction in real space.

While $[eqn]$ and $[eqn]$ ; $[eqn]$ and $[eqn]$ ; and $[eqn]$ and $[eqn]$ , respectively, point in the same direction in an isotropic material (Figure 3a), this is not the case in anisotropic media (cf. Figure 3b). The isotropic condition $[eqn]$ is replaced by $[eqn]$ and $[eqn]$ . The vectors $[eqn]$ , $[eqn]$ , and $[eqn]$ are all in one plane, and $[eqn]$ and $[eqn]$ are perpendicular to it. We refer the reader to Ref. [29] for details. Maxwell’s equations yield the relation

$$[eqn]$$

providing the desired relation between $[eqn]$ and $[eqn]$ , the dispersion relation [31], in the form of an eigenvalue problem for $[eqn]$ . The existence of non-trivial solutions requires the determinant yielding the characteristic polynomial to vanish. The relation derived from this condition is referred to as the dispersion relation [29].

In the special case of a uniaxial crystal with $[eqn]$ in, say, the y direction that we consider here, the situation simplifies as the expression for the determinant factorizes. As each factor can become zero, we find two solutions: one dispersion relation is spherical (and involves only $[eqn]$ ) and represents the ordinary wave, while the other possibility is a spheroid associated with the extraordinary wave we are particularly interested in here. It corresponds to a dispersion ellipse that is elongated in the $[eqn]$ direction in momentum space, with the length of the semi-axes being directly related to the refractive indices such that the semi-axis in the $[eqn]$ direction is proportional to $[eqn]$ , while the one in the $[eqn]$ direction is governed by $[eqn]$ . Thus, the dispersion directly reflects the anisotropy of the system, and we introduce the anisotropy parameter a, $[eqn]$ ,

$$[eqn]$$

We add a qualitative argument to make the rotation of the ellipse (from an index ellipsoid elongated in the y direction to a dispersion relation elongated in the $[eqn]$ direction) required when going from the index ellipsoid to momentum space more plausible. In a medium with refractive index n, the wavelength (a measure in real space $[eqn]$ ) decreases by n in comparison to a vacuum, while the wave number (related to momentum space ( $[eqn]$ ) increases by n. This reciprocal relation between real and momentum space (we mention a similar reciprocity when constructing the reciprocal lattice, where large reciprocal lattice vectors correspond to small distances between the associated lattice planes and vice versa) suggests a dispersion relation elongated in the $[eqn]$ direction with the equation of an ellipse

$$[eqn]$$

with the anisotropy parameter a defined in Equation (4). Its consequences for the reflection at the wall of an anisotropic billiard are illustrated in Figure 3b) in comparison to the inverted isotropic situation (circular dispersion relation, now with an elliptical real-space geometry) in Figure 3a).

The PSOSs resulting from ray tracing for the situations depicted in Figure 3a,b, i.e., for a closed isotropic ellipse and an anisotropic disk, both with the same eccentricity parameter corresponding to $[eqn]$ , are shown in Figure 3c,d for comparison. They resemble each other [32,33], up to the curly boundary limits arising from the non-circular dispersion relation, as discussed above. The striking characteristics of the PSOS of the anisotropic disk billiards (cf. also Figure 4) are two islands around central elliptic fixed points, in analogy to the PSOS of the isotropic elliptic billiards, where they represent so-called bouncing-ball orbits that go back and forth between the flat sides of the ellipse. In our anisotropic case, the bouncing-ball trajectories illustrated in Figure 4(7–11) arise from the existence of preferred propagation directions in momentum space, corresponding to flat sides of the ellipse of (now) the dispersion relation. The group velocities in these regions, given by the normal to the dispersion relation, combine to form bouncing-ball orbits (similar to the triangular orbits seen in BLG billiards above (cf. also Figure 3)) along the diameter. In the case of a circular dispersion relation, diameter orbits exist in all directions but are not stable (rather, they are marginally stable). For the anisotropic case, there is also a pair of unstable fixed points in between the stable (bouncing-ball) ones, indicating a diameter orbit formed by the group velocities at the points of highest curvature of the elliptic dispersion relation.

In addition, there are WG-type trajectories that do not change the sense of their propagation direction and travel along the billiard boundary (cf. Figure 4(1–5)). Their trace in phase space depends on their sense of rotation and lies outside the island region towards the upper (lower) PSOS edge for counterclockwise (clockwise) propagation. These regions (of WG and bouncing-ball orbits, respectively) are disjoined by the separatrix, and the trajectory in cf. Figure 4(6) shows a sample trajectory very close to the separatrix (but slightly on the bouncing-ball side).

The evolution of the PSOS of the anisotropic disk with an increasing anisotropy parameter a is shown in Figure 5. Starting from the isotropic disk, $[eqn]$ , with a PSOS consisting of WG orbits represented by straight lines, the bouncing-ball island signature becomes the dominating structure in phase space as a increases. Larger a correspond to higher anisotropy and the stronger formation of preferred propagation directions. This is consistent with a larger portion of bouncing-ball orbits, as these are composed of these preferred propagation directions. The red line denotes the separatrix that separates these two regions of qualitatively different trajectory dynamics.

3. Wave and Transformation Optics for the Birefringent Disk

We start our considerations with Maxwell’s equations, which we have to solve for the electric field $[eqn]$ and the magnetic field $[eqn]$ for an anisotropic medium with permittivity tensor $[eqn]$ and permeability tensor $[eqn]$ .

$$[eqn]$$ $$[eqn]$$

While their solution is straightforward for some specific anisotropic solutions, as argued above, when discussing this in terms of the index ellipsoid, their generic solution states a mathematical challenge [34,35]. In physics, transformation optics provides a well-established method to overcome these mathematical problems.

The central idea of transformation optics [16,36] is the constancy of the optical path length, which is the product of the refractive index times the geometric path length. If there is, for example, an anisotropy in the refractive index in one direction, one may rescale this refractive index by ’counterscaling’ the respective spatial direction (such that it fits the other ones and, thus, yields an isotropic permittivity tensor $[eqn]$ ). It is essential to note that stretching or compressing the coordinate system in a given direction modifies the permeability tensor $[eqn]$ even though its components are uniform in the original system.

Here, we assume an anisotropy in the y direction and focus on TE polarized light, as before. The corresponding coordinate transformation from the old coordinates $[eqn]$ to the new ones $[eqn]$ thus reads (cf. [18])

$$[eqn]$$

with

$$[eqn]$$

While requiring constancy of the optical path length readily yields the stretching factor $[eqn]$ , it can also be derived via the permittivity tensor. Our starting point is the anisotropic situation

$$[eqn]$$

and the following ansatz for the transformation matrix $[eqn]$ :

$$[eqn]$$

The transformed quantities (denoted by ′) in the transformed coordinates $[eqn]$ follow from $[eqn]$ as $[eqn]$ for vectors, while matrices fulfill the relation [18]

$$[eqn]$$

We thus obtain

$$[eqn]$$

yielding a transformed isotropic permittivity tensor $[eqn]$ when choosing $[eqn]$ . Note that $[eqn]$ is now anisotropic such that $[eqn]$ and $[eqn]$ have changed their roles under the transformation. The dispersion relation remains unchanged as it is symmetric when interchanging $[eqn]$ and $[eqn]$ [31].

The diagonal elements of $[eqn]$ are then all equal to $[eqn]$ , corresponding to a material with a mixed refractive index $[eqn]$ . However, the transformed equation for the z component of the magnetic field $[eqn]$ reads

$$[eqn]$$

resembling the equation of an isotropic medium with the extraordinary refractive index $[eqn]$ , in contrast to the entries of $[eqn]$ . This can be understood by looking at the permeability tensor $[eqn]$ . For TE polarization, $[eqn]$ equals $[eqn]$ , and only the z component of $[eqn]$ plays a role, $[eqn]$ . Multiplied by the isotropic $[eqn]$ , it indeed yields a refractive index of $[eqn]$ in Equation (14).

Transforming the boundary of the disk accordingly (i.e., $[eqn]$ ), we arrive at an ellipse

$$[eqn]$$

with (unchanged) semi-axis lengths of one in the x (or u) direction, and a flattened semi-axis of $[eqn]$ in the y (v) direction, to constitute the boundary of the transformed, isotropic cavity. We point out that the values of the semi-axes are the same as for the ellipse representing the dispersion relation in k-space, Equation (5), which was used for the ray simulation, as it should be. The preferred group velocities, perpendicular to the flat sides of the dispersion ellipse, correspond to the bouncing-ball orbits of the corresponding isotropic ellipse (representing the transformed disk), yielding a coherent view of the dynamics occurring in the anisotropic disk. It also provides some insight into the mathematical challenges when studying the anisotropic instead of the isotropic disk. The mathematical complexity is induced by going from a geometry with one center point (circle) for the isotropic disk to a (transformed) geometry with two focal points (ellipse). This can be seen as a bifurcation. Notice, however, that the limit is well defined when taking the opposite direction and going from an ellipse to an (isotropic) disc.

The solution of closed isotropic elliptical billiards is well known [37,38] and yields eigenmodes in the form of Mathieu functions (instead of Bessel functions, as for a closed isotropic disk), with quantum numbers m for the angular motion and $[eqn]$ in the radial direction, and the requirement of continuity across the focal points of the ellipse. These solutions have to be mapped back to the original disk. By doing so, we arrive at the desired solutions for the anisotropic disk, which are shown in Figure 6 in qualitative comparison to their ray trajectory counterparts and in Figure 7 together with their isotropic ’elliptical’ counterparts. The more the (isotropic) ellipse geometry deviates from the circle, the more anisotropic is the disk cavity that it describes according to transformation optics.

4. Discussion: Ray–Particle-Wave Correspondence for Optical and Electronic Systems

What can we learn from the description of two different mesoscopic systems, an electronic and a photonic one? The issues are as follows. For the optical system, transformation optics relies on the invariance of Maxwell’s equations under coordinate transformations like the one we used here and has proven very useful in many cases. In contrast, a $[eqn]$ Hamiltonian is needed to describe bilayer graphene; moreover, electrons are charged spinful particles. To highlight the differences that are of interest here, the particle tracing point of view is useful. In the PSOS in Figure 1(a.1,b.1), we had to distinguish two different surfaces of section, depending on the valley index of the BLG billiards ( $[eqn]$ and $[eqn]$ ). It is straightforward to differentiate these cases when no intervalley scattering is present, as may be the case for low carrier densities and for external fields that lift the level of degeneracy. However, in general, both K valleys have to be taken into account. This is fulfilled in any tight-binding calculation [9] that is based on the full Hamiltonian of the system. Conversely, taking the transformation optics point of view and starting from a circular BLG billiard with a trigonally warped Fermi line, it is not guaranteed to find a non-circular BLG billiard geometry with a circular Fermi line that represents the same dynamics in a PSOS. In fact, it is impossible when only one valley ( $[eqn]$ or $[eqn]$ ) is considered. The reason is that a PSOS corresponding to a triangularly warped Fermi line is not symmetric with respect to $[eqn]$ [10,28], while it is symmetric for the associated isotropic but deformed billiard geometry. The consequence that changing the sense of rotation for an individual trajectory yields different dynamics in BLG billiards for a given pseudospin marks a distinct difference compared to the optical case.

This investigation of anisotropic BLG billiards provided the tools to study birefringent optical systems by extending ray optics to the anisotropic situation. This approach completes the wave/transformation optics point of view and establishes ray–wave correspondence for anisotropic billiards for light. To this end, we express the anisotropy of a uniaxial optical system in terms of the index ellipsoid as usual, and assume the y axis of a two-dimensional optical system in the $[eqn]$ plane to be an optical axis with a refractive index $[eqn]$ . We first used a ray-tracing approach motivated by BLG billiard studies and derived the dispersion relation (a quantity freely available in solid state physics)for the optical case, arriving at an ellipse in momentum space with the elongation in the $[eqn]$ direction. The latter becomes plausible when employing reciprocity arguments relating real and momentum space. In contrast, transformation optics does not require expedition to momentum space, though it arrives at the very same elliptical geometry now in real space in the case considered here. The reason is that we started with an anisotropic disk—its isotropic counterpart assigned by transformation optics is the isotropic elliptical cavity, which has the same shape as the dispersion relation. This occurs by chance in the sense that if we had started with a non-circular geometry, the transformation optics would yield a another transformed cavity shape, while the dispersion relation would remain unchanged. We note that in the case of BLG, changing the geometry of the cavity and the Fermi line has different effects on the overall system dynamics [39].

For non-circular cavities, the ray–particle tracing approach proves to be particularly useful as it treats any given anisotropy that can be combined with an arbitrary billiard geometry using the very same algorithm, and ray–wave correspondence is applied to generalize the results to describe the wave counterpart. When working on the wave side, transformation optics is ideally suited to transform anisotropic wave problems to isotropic ones that are much easier to solve. The knowledge gained from the (numerically cheap) tracing simulations can provide useful guidance in the interpretation of results.

We focused on uniaxial and closed optical systems; in particular, we assumed the optical axis to lie in the resonator plane, representing the most interesting yet simple case. Generalization to the generic case of an arbitrary pointing optical axis that does not coincide with the resonator plane, as illustrated in Figure 2a, is conceptually straightforward though tedious [29]. Here, we limited our consideration to the case of closed mesoscopic cavities. The important topic of realistic open systems will be discussed elsewhere.