1. Introduction

In optical resonators of high finesse, the sensitivity to boundary conditions grows with decreasing size. For electrodynamic mode calculations, it is therefore important not only to accurately account for geometric shape of bounding interfaces, but also for the way in which fields penetrate these interfaces. For dielectric mirrors, the resulting phase retardation depends on angle of incidence, and this can make even simple resonator geometries hard to treat if the cavity modes contain a large spread of incident angles in their plane-wave spectrum.

Paraxial optics largely circumvents this problem because it relies, by definition, on small angular spread[1]. However, limitations of paraxial theory have long been recognized both in the context of macroscopic and microscopic high-finesse cavities when the relevant numerical aperture becomes sufficiently large; for recent examples see, e.g., [2] and [4]. In non-paraxial resonator geometries, Gaussian beams become less suitable as a starting point, so that one retreats to more general approximation methods, such as iterative application of the Rayleigh-Sommerfeld diffraction integral [5]. It can also be attempted to connect diffraction theory back to the Gaussian beams, thereby arriving at non-paraxial correction terms[6].

In this paper, we explore an alternative approach which can be connected smoothly to the limit of paraxial optics, preserves significantly more physical insight than numerical diffraction theory, and enables us to investigate corrections to paraxial results in analytic or semi-analytic form. As a specific aplication, we address mixed boundary conditions that can be characterized by angle-dependent reflectivities. r(θ). Here, θ is the angle of incidence of a ray with respect to the surface normal, and r(θ) = exp(iφ) with a real phase φ= φ(θ) if there is no absorption or leakage at the mirror. It will be shown that even in the paraxial limit, the Fourier-relationship between position and angle makes it necessary to take the angular spread of Gaussian cavity modes into account when modeling the interaction with dielectric multilayer mirrors.

 

Fig. 1. (a) Sketch of the cavity geometry. The vertical axis is z, and the radial distance from the z axis is r. (b) Example for a paraxial mode (gray scale: intensity). The height function of the top mirror is z _M, and its maximum is at z = h. The Born-Oppenheimer potential V _BO, superimposed as a sketch, is responsible for the transverse confinement of the mode. Here, h = 3μm and kh = 37.5; the field is calculated using the method described in this paper.

Download Full Size | PDF

In particular, the relationship between the Gouy phase and transverse mode splittings in a paraxial dome-shaped resonator is modified when dielectric Bragg mirrors introduce phase shifts that depend on the angle of incidence. Gaussian-beam parameters such as the Gouy phase can be calculated purely based on ray optics[1], despite the fact that the transverse field pattern is typically not in the classical limit of large quantum numbers. A reflectivity that depends on angle of incidence is easily introduced in the ray picture, but not as easily incorporated into the paraxial theory, precisely because the transverse beam structure requires a wave description[3].

To describe the solution to this problem, we restrict attention to the scalar wave equation. Although the method is quite general, we specialize to a resonator geometry shown schematically in Fig. 1. We have previously studied this plano-concave geometry using numerical computation of the vectorial electromagnetic modes [4, 7, 8], motivated by recent experimental progress [9, 10, 11, 12].

2. Adiabatic separation of variables

The central ingredient in our approach is a “local-mode” approximation, known in molecular physics a the Born-Oppenheimer (BO) method[13]. It relies on an adiabaticity assumption whose limits of validity differ from that of the paraxial method. We will show, however, that the two methods not only have an intersecting regime of validity, but that the BO approach indeed extends beyond paraxiality.

As shown in Fig. 1, the z axis is assumed to point perpendicular to the plane of the flat mirror at the base of the cavity. In the wave equation

the central step is to assume that the x,y variation is of long wavelength, and use the envelope ansatz ψ = η(z;x,y)χ(x ,y) in the wave equation. In the lowest (adiabatic) approximation, we neglect the dependence of η on x,y, allowing us to separate off the equation for η.

For simplicity, the top mirror is assumed to be an ideal metal, with reflectivity r ≡ - 1, while the bottom mirror has reflectivity r = exp(iφ). The height of the top mirror at x,y is z _M(x,y); e.g., for a spherical mirror as in Fig. 1(a) we have z _M(x,y) = h + √R ² - x ² - y ² - R, where R is the radius of curvature and h the maximum cavity length, as reached along the z axis. In the paraxial approximation, this function is replaced by the quadratic expansion for small x,y, z _M(x,y) ≈ h - r ²/(2R), where r ² = x ² + y ². The solutions for η are of the form η (z;x,y) = e ^-ik_zz + e ^iφ e ^ik_zz At a given (x,y), the condition fixing k_z is then

where the integer v is the longitudinal quantum number.

The transverse field satisfies the analogue of a two-dimensional Schrödinger equation, with the units chosen such that h̄ ² ≡ 2 and mass M = 1:

where we have defined

and

which we shall call the Born-Oppenheimer potential. The shape of this potential V ^(ν) _BO depends on v as defined in Eq. (1). The quantity K ² plays the role of the energy eigenvalue, except that K as defined in Eq. (3) is a true constant only if the reflection phase φ is a constant We will make this assumption φ = const in the following sections, and relax it again in Section 5.

3. Gaussian beams from the Born-Oppenheimer method

3.1. Transverse harmonic oscillator and mode waist

In this section we show that the results of Gaussian-beam theory can be recovered with our approach if V ^(ν) _BO in Eq. (4) is a harmonic potential. The main goal is to show how the Gouy phase and spreading angles of paraxial theory arise in our approach.

We arrive at a harmonic BO potential using z _M(r) ≈ h - r ²/(2R), cf. Fig. 1 (b):

This represents an axially symmetric cavity, for which the preferred basis of solutions are the Laguerre-Gauss functions,

where L_p ^∣ℓ∣(u) is the associated Laguerre polynomial. The integer ℓ is the orbital angular momentum, and p counts the number of radial nodes in the transverse wave field.

The length w_v appearing in this equation is related to Ω_v of Eq. (5) by

This is just the fundamental mode waist radius at the focus. It agrees with the Gaussian-beam result W _paraxial = √2 (h (R - h) /k ²)^1/4 when h ≪ R. The eigenvalues, Eq. (3), are now

 

Fig. 2. Wave number k for four longitudinal modes, each with their lowest three transverse levels, as a function of cavity height h, as given by Eq. (10). In varying h, the longitudinal phase, vπ+π/2 -φ/2, is adjusted such that a fundamental transverse mode exists at k = 12.5(μm)^-1; this fixes the lowest mode in the plots to be a horizontal line. The grouping of curves into triplets corresponds to clustering of modes with the same longitudinal quantum number v, i.e., same phase vπ+π/2 -φ/2. The transverse levels in each cluster are enumerated by the principal quantum number, N = 0, 1, 2. For mirror radius of curvature R = 45μm, the transverse mode spacings are seen to be much smaller than the longitudinal spacing, so that the displayed clusters do not intersect. For smaller radius of curvature, R = 10μm, the spacing of k ² between transverse and longitudinal modes becomes comparable near h ≈ 17μm when R = 10μm. Note that for R = 45μm, the Born-Oppenheimer method is valid over the entire range of h shown.

Download Full Size | PDF

where N = 2p + ∣ℓ∣. For a given N, the allowed angular momenta are ℓ = -N, -N + 2, …N - 2, N, with the highest value ∣ℓ∣ = N occurring when p = 0. The radial wave equation for Ψ contains an effective radial potential, V _eff(r) = V ^(v) _BO(r) + ℓ ²/r ²; it depends both on v and on ℓ , but not on the principal quantum number N.

In analogy to free-space optics, we can aso introduce the paraxiality parameter θ_p, which measures the far-field divergence angle of the wave if continued to large z. For a Gaussian beam, this angle is diffractively related to the beam waist, i.e., θ_p and w_v satisfy an uncertainty relation:

3.2. Gouy phase

From Eq. (8), we can obtain the transverse mode splitting by solving the definition of K in Eq. (3) for the mode wave number k in terms of the mode indices N and v:

It is worth pointing out that within the adiabatic approach, wavenumbers always appear as squared quantities. It is therefore k ² that becomes quantized, not k. In paraxial theory, on the other hand, one quantizes k directly by imposing the round-trip single-valuedness of the field; the transverse excitation is viewed as giving rise to an additional phase shift in each round trip, called the Gouy phase, which must be added on to k. We now show how these two approaches are connected.

Taking the square root of Eq. (10) and assuming the transverse “energy” K _N;v ² to be much smaller than the longitudinal contribution k ² _z (r = 0), we get

The validity of this approximation can be verified by inspection of the spectra shown in Fig. 2. The plots use Eq. (10), where we fixed the wave number window to around k = 12.5μm and vary the cavity length from h ≈ 2μm to h ≈ 18μm. Comparing Eq. (11) with a Fabry-Pérot spectrum (plane-parallel mirrors), we see that the last term is the Gouy phase contribution to the wave number,

This result agrees with conventional resonator theory where the Gouy phase θ_G for a dome geometry is given by [1] cosθ_G = √1-h/R, whence θ_G → √h/R for R ≫ h. This reduces to k _G as defined above when divided by h to turn the angle θ_G into a wavenumber.

The Gouy phase describes a phase shift of the resonator wave solution which is purely determined by geometric optics. Formally, this phase is the eigenphase of the round-trip monodromy matrix describing the stable ray orbit on which the paraxial cavity mode is built [14, 15], and hence it measures the winding number of the off-axis rays corresponding to the transverse mode profile; cf. also [2] and references therein.

3.3. Range of validity

The calculations in this section have shown that the BO method is an alternative to the paraxial approximation for calculating all the relevant properties of a Gaussian beam. The two approaches become identical in the limit h≪R. This is also the limit in which the transverse mode splittings of size k _G become small compared to the longitudinal free spectral range, π/h. The validity of the BO approximation relies precisely on this separation of frequency scales between neighboring longitudinal (v) and transverse (N) mode indices [13]: it insures the adiabaticity that prevents couplings between different v. It is also worth noting here that the magnitude of v does not enter the conditions of validity in our approach, and therefore we can describe modes that have either small or large values of the dimensionless wave number kh, as is illustrated in Fig. 2.

Equation (10) reduces to the standard paraxial resonator formula, Eq. (11), only in the limit K ² _N;v ≪ k ² _z(r = 0) which requires kh ≫ hk _G = θ_G; but this dependence on the absolute size of kh is a limitation of the square root expansion leading to Eq. (11), and not of the the BO treatment. In Fig. 2, this condition is in fact valid for all data shown.

In the Helmholtz equation, any cavity with axial symmetry conserves angular momentum ℓ, and therefore it is straightforward to incorporate anharmonic potentials into V _BO ^(v)(r) without losing the effective-potential picture: in Fig. 3, the curves representing the effective radial potential will change shape, but the transverse field is still determined by a one-dimensional radial potential well with discrete bound states. This is one of the main advantages of the BO method in cavity mode calculations: the quadratic approximation loses the special status which it occupies in the paraxial approximation. In particular, for an axisymmetric cavity, the breaking of the degeneracies (8) is necessarily non-paraxial, but still falls within the adiabatic regime on which the BO approach relies.

In the absence of axial symmetry, V _BO ^(v)(x,y) ≠ V _BO ^(v)(r), we can still use Eq. (4), except that the transverse wave equation is not separable in cylinder coordinates. Aberrations such as astigmatism[17] can therefore be treated with the BO method as well.

When judging the validity of the adiabaticity assumption using Eq. (11), one has to keep in mind that this result is derived by making the approximate separation of variables in cylinder (or Cartesian) coordinates. Even when Eq. (11) predicts the breakdown of spectral separation between longitudinal and transverse modes, there may exist other coordinate choices (e.g., parabolic cylinder coordinates [16]) in which this separation is restored. This greatly increases the range of geometries to which the BO approach can be applied, even allowing near-confocal cavities with h ≈ R. Examples for BO approximations in different two-dimensional coordinate systems have been given in [19].

 

Fig. 3. Effective radial potential V _eff(r) for the Laguerre-Gauss modes (solid lines, labeled by orbital angular momentum ℓ = 0, 1, 2). Shown as horizontal lines are three values of the effective energy K ² for mode orders N = 0, 1, 2. Laguerre-Gaussians are superimposed with the same line style and at the same ordinate as the K ² to which they correspond, and labeled by their ℓ-value.

Download Full Size | PDF

4. The local plane wave approximation

4.1. Comparing Bessel-wave expansions

The eigenfunctions of a harmonic oscillator have the same form in position space as in momentum space. Formally, this is due to the symmetric appearance of the non-dimensionalized positions q and momenta p in the Hamiltonian H = p ² + q ². In this notation, if ψ(q) is the position representation of an eigenstate, then ψ(p) is its momentum representation, which means ψ(q) is equal to its own Fourier transform.

Applying the same reasoning to the Laguerre-Gauss solutions in polar coordinates, Eq. (6), one finds that the radial part R_p,ℓ,v(r) of Ψ_p,ℓ,v(r,ו) is formally identical to its own Fourier-Bessel transform[18], which we call g_p,ℓ,v(u) This is defined to be the rotationally symmetrized form of the Fourier transform, expressed as the integral over the dimensionless radial coordinate ρ = r/w_v:

(its inverse is identical in form). If we define the radial function without regard for normalization constants,

then according to the above argument its Fourier-Bessel transform must read

We now relate this general property of harmonic oscillators to the angular spectrum of the paraxial Laguerre-Gauss beams; that is the symmetrized plane-wave decomposition of the field, in terms of the Bessel waves, J_ℓ(kr sinθ)exp{ikz cosθ θ iℓϕ} (in cylinder coordinates). We define the plane-wave spectrum g(θ) as the expansion coefficients in the integral decomposition

Here, the left-hand side is the upward propagating beam contained in the BO mode. Even without assuming the BO form, it can be shown within Gaussian-beam theory that g(θ) for a Laguerre-Gauss mode is given by[20]

when the approximation sin θ ≈ θ is valid. In the angular distribution, the spreading angle θ_p defined in Eq. (9) thus takes over the role of the mode waist w_v via the substitution r → w_v θ/θ_p. We now want to derive Eq. (16) with the BO approach by reducing Eq. (15) to Eq. (13).

4.2. Local angle of incidence from local wave number

To get from the integral over θ in Eq. (15) to the form of Eq. (13), we drop the ϕ dependence and then proceed as follows to eliminate z: integrate both sides over z, from the bottom of the resonator at z = 0 to the top at z = z _M(r), with the factor exp(-ik_z(r)). This yields

The paraxial limit allows us to simplify the integral: first we can assume that the angular spectrum g(θ) is narrowly confined to θ ≈ θ, so that the upper integration limit, θ = π/2, can be replaced by ∞, and we can also replace sin θ ≈ θ wherever the integrand is appreciable.

A second assumption motivated by paraxiality would be to set

allowing us to use lim_x→0sin(x)/x = 1 in Eq. (17) and thus leading to

Comparing this result with Eq. (13), we identify g(θ) as a rescaled version of the Fourier-Bessel transform S_p,ℓ,v(u), i.e, g_p,ℓ,v(θ) = C′ S_p,ℓ,v(θ).

By reproducing this result which is known from Gaussian-beam theory, the assumption made in Eq. (18) is therefore justified a posteriori. Equation (18) is important for what follows, and hence it is worth discussing its implications further. Accordingly, there is then a special angle θ(r) for every radial position r, given by

Physically, this represents a local plane wave approximation because Eq. (15) identifies θ(r) as the angle of the Bessel wave that we would obtain if the cavity had the same height h = z _M(r) for all r.

Next, we bring Eq. (20) into a more explicit form. Since most of the “energy” of the modes resides in the longitudinal motion, we can use k_z ≈ k to approximate

using Eqs. (4), (9) and (10). The Bessel-wave limit, θ = const, is recovered from this when w_v → ∞. Although this limit makes θ_p go to zero, we can get finite θ(r) by simultaneously choosing N large. The term containing N in Eq. (21) originates from the transverse mode splitting k _G(N + 1) of Eq. (11), which depends on the Gouy phase θ_G.

Clearly, this expression allows only a limited range of r and θ. The largest allowed radius is the outer classical turning point of the BO potential V _BO ^(v)(r),

Here, we have also noted the connection to the root-mean-square width, r _rms, of a two-dimensional harmonic oscillator in the N-th excite state.

4.3. Semiclassical discussion

As illustrated in Fig. 3, there is in fact a well-defined interval of allowed r, delimited by the classical turning points, r _min and r _max, in the effective radial potential, V _eff(r). The calculation of rmin and r _max amounts to reading off the intersections of the potential curves (solid lines) in Fig. 3 with the lines of K ² = K _N ² = const (the three horizontal lines in the graph). For the lowest orders N of transverse excitation, we can expect the turning points to be washed out in the wave solutions, due to tunneling into the effective-potential barrier. This is illustrated by the transverse Laguerre-Gauss mode functions superimposed as dashed lines in Fig. 3. One can verify that in all cases the classical turning points give a reasonable estimate of the interval within which the respective mode functions have large amplitude, and that only exponential-type decay is observed in these functions when r ∉ [r _min, r _max].

The outer turning points, r _max, of V _eff(r) are less than or equal to r _outer. Specifically, r _max → r _outer for N ≫ ℓ . The points r _max are seen in Fig. 3 to coalesce at a single value that depends on N, but not on ℓ, when N ≫ ℓ; this is because the shapes of the effective potentials for different ℓ approach each other on the large-r side of the graph. The angle θ(r) is largest when r takes on its smallest value, r _min. Both r _min and r _max depend on the angular momentum of the mode, but r _min has a stronger ℓ dependence due to the effect of the centrifugal barrier ℓ ²/r ² in V _eff(r).

The transverse modes do not behave classically, but their splittings are much smaller than the longitudinal free spectral range. The formal analogue of this in molecular applications of BO theory is the dense spectrum of the nuclear degrees of freedom in the effective potential created by the electronic motion. The inner and outer turning points of the LG modes with angular momentum ℓ and principal transverse quantum number N are given by

For ℓ = N = 0, we get r _min = 0 and r _max = w_v, the conventional mode waist radius. Note also that we always have

i.e., the root-mean-square radius defined in Eq. (22) is bracketed by the classical turning points.

Equation (21) predicts a monotonic decrease of θ with increasing r. This is physically significant: e.g., for N = ℓ = 0 we find that there is a finite angle of incidence θ = θ_p at r = 0. In the BO approach, the smallest angle of incidence with respect to the normal occurs at large radial distance from the axis. This is what one also expects from the ray picture.

Inserting r _min/max into θ(r) from Eq. (21), the extremal plane-wave angles of incidence are

This can equivalently be deduced from the duality between θ and r, cf. Eq. (16). The significant content of θ_max/min is that they depend on the mode indices N and ℓ, whereas θ_p is independent of these transverse quantum numbers. For fixed value of the fundamental mode waist, w_n, both r _max and θ_max increase when the transverse order N increases, in contrast to the reciprocal relation between w_n and θ_p.

5. Angle-dependent reflectivity

5.1. Including variable reflection phase shifts

In this section, we will pursue the interplay between position and Bessel-wave representation by making use of the local relation for θ(r), Eq. (21). In a realistic cavity we want to allow the mirror reflectivity to depend both on position and on the angle of incidence. Here, we introduce a non-constant reflectivity only for the planar bottom mirror, and neglect absorption and leakage because we want to focus on the closed-resonator mode structure. If the planar interface is moreover translationally invariant, its reflectivity can be written as r(θ) = exp(iφ(θ)), i.e., it is a function of incident angle θ but not of position r.

The longitudinal wave number is determined by the mirror shape and reflectivities according to Eq. (1), where θ enters through the reflection phase in the form[21]

Here, ε need not be small, but will be assumed to be a constant (the factor of two in front is introduced merely for convenience later). In fact, for multilayer mirrors, ε in general also depends on the wavelength; however, in what follows we shall neglect that dependence. It will be interesting for a complete analysis to allow either positive or negative sign for ε, even though in practice (see, e.g., [4]) for dielectric multilayer stacks one always finds ε < 0, and often ε < −1. The case ε > 0 can be realized in particular for meta-materials with negative index of refraction, and in view of recent advances in scaling such materials to optical wavelengths[22] it seems appropriate to consider this a physically realizable and novel scenario. In order to distinguish the constant and variable parts of the phase entering the longitudinal quantization condition, we abbreviate

Now insert the expression for the local incident angle, Eq. (20), into the expansion for the phase shift, Eq. (25). With θ² ≈ 1-k ² _z/k ², Eq. (1) becomes

which determines k_z(r) implicitly. Because the angle of incidence θ enters squared and not linearly, the square-root singularity of Eq. (21) plays no role here, and hence the classical turning points do not limit the interval of r on which k_z(r) is defined. This is important to note because it means that Eq. (21) will allow us derive a modified BO potential for the radial wave equation, which is defined for r ∈ [0, ∞[. Quartic or higher powers of θ in the expansion for φ(θ) can be taken into account in straightforward generalization of this approach, but we will focus on the leading-order effects which directly modify the paraxial results.

5.2. Modified transverse wave equation

Paraxial theory is, in our framework, the limit of a quadratic BO potential. Equation (26) is a quadratic equation for k ²(r), and the solution to quadratic order in r can be written as

where we define the angle γ through

for small ε this is the only physically acceptable solution. The parameter γ has been introduced only to make the notation more compact; it may be real or imaginary.

For ε≠0, Eq. (27) will make k_z = k_z(r,k), where the additional k dependence in k_z is not present for φ = const. In other words, the eigenvalue k, which is determined by the radial equation, appears in the longitudinal quantization. This is a new coupling that also makes the effective potential k-dependent: V _BO ^(v)(r,k) ≈ Ω² _v,ε (k) r ²/2 with a k-dependent oscillator frequency

Laguerre-Gaussians are still the solutions to the new radial equation, analogous to Eq. (6), but with some modifications. We know from Eq. (8) that Laguerre-Gaussians yield the eigenvalues K ² _p,ℓ;v;ε(k) = √2Ω_v,ε(k) (N+1) = K _N;v;ε ²(k), where we use the radial (p) and azimuthal (ℓ) quantum numbers as mode labels, and N = 2p + ∣ℓ∣. The solution to our problem consists in finding k such that K _N;v;ε ²(k) = k ² - k ²(0), i.e.,

where Eq. (27) was used to rewrite k ²(0). Because Ω² _v,ε(k) and k ²(0) depend on k ² through γ, this is an implicit equation for k ², but even without solving it we can see that the wavenumber k will depend only on the principal quantum number N, not on ℓ or p separately. Numerical solutions of Eq. (29) are shown in Fig. 4(a). They first of all confirm that despite this more complicated situation the BO separation of spectral scales is still valid. In fact, the global appearance of the spectrum is almost indistinguishable from the ε = 0 case, Fig. 2(a), except at small cavity lengths h.

5.3. Results and discussion

Short cavities, h ≪ R, are a domain where our approximations are especially good, so we will look at cavities of this type more closely, to illustrate the method we have introduced in the foregoing discussion. The angle-dependent reflectivity giving rise to Eq. (29) will not break the transverse-mode degeneracy of Eq. (8), and therefore we can label the multiplets of identical eigenvalues k _p,ℓ,v;ε = k _N,v;ε. For ε ≠ 0 the k ² _N,v;ε for N = 0, 1,2… are no longer strictly equidistant, as they are for the standard harmonic oscillator. This is intriguing because we have therefore broken one spectral property of the harmonic oscillator (its equidistant spectrum) while strictly preserving another (its degeneracies). Quantitatively, however, the deviations from equidistant transverse mode spacing are negligible for most cavity parameters.

 

Fig. 4. Wave number spectra: (a), same as in Fig. 2(a), but for angle dependent reflection phase, with ε = -5. The Born-Oppenheimer separation of transverse and longitudinal splittings is still valid, but deviations from Fig. 2(a) are visible at small h. In (b), the phase coefficient ε is varied while keeping the resonator geometry fixed at h = 2μm and R = 45μm. In both plots, we adjust the constant part of the longitudinal phase following the same prescription as in Fig. 2. The three lowest transverse modes from part (a) shift and become more degenerate as ε increases.

Download Full Size | PDF

Now we turn to the interplay between Gouy phase and boundary reflection phase. In Fig. 4(b), we inspect the level spacings on the scale of the transverse mode splitting, for a single longitudinal manifold. The most striking observation is that the spectral spacing of adjacent k _N,v,ε is no longer equal to the purely geometric quantity k _G of Eq. (12), despite that fact that there is no splitting of the harmonic-oscillator degeneracies. In other words, the Gouy phase θ_G cannot be read off directly from the mode splittings. The spectrum now shows a spacing between k _N,v,ε that is determined not only by k _G, but in addition by the angle-dependence of the mirror reflection phase (ε). Recall that a constant reflection phase shift φ₀ would not affect the paraxial mode spacing, only cause a global spectral shift.

The relations derived above allow us to disentangle these effects, which are due to the interplay of the k-dependencies in the transverse eigenvalues K ² _N;n;ε(k) and the longitudinal wavenumber k _z ²(0), entering in k ² = K ² _N;v;ε + k _z ²(0). Figure 5 shows the ε-dependence of (a) the oscillator frequency Ω_v,ε entering K ² _N;v;ε, and (b) the wave number k _z ²(0). Increasing ε will reduce k ² _z(0) but increase Ω_v,ε, but in Ω_v,ε this is outweighed by the growth in cos^-1/2〻 which in fact diverges for ε(α - ε) → (hk/2)², cf. Eq. (28). The meaning of Fig. 5 (a) is that increasing ε causes the transverse modes to be confined more tightly. In other words, the mode waist Ω_v,ε given by Eq. (7) will expand for ε < 0 and shrink for ε > 0. In a usual paraxial resonator without angle-dependent reflection phase shifts, one would expect larger waist radius to come with smaller transverse mode splitting. However, in our case the opposite happens: according to Fig. 5 (b), value of k _z ²(0) increases for ε < 0, and it does so in a mode-dependent way, with larger slope for the higher-order transverse modes. As a result, the total spectral separation between the actual mode wave numbers k ² = K ² _N;v;ε + k ² _z (0) increases for ε < 0.

One can ask whether this outcome of the balance between Figs. 5 (a) and (b) is a numerical coincidence, or whether there is reason to expect this counter-intuitive spectral phenomenon more generally. The latter is in fact true, as one can argue both mathematically and physically. To understand the functional behavior of Fig. 4 (b), it is useful to analyze the ε > 0 region where the spectral coalescence occurs. Because k ²(0) and Ω_v,ε both determine and are determined by k, the limit ε (α - ε) → (hk/2)² produces a fixed point to which all k of a given longitudinal mode number v (and hence phase α) flow. The curve k(ε) near coalescence then approaches k(ε) ≈ 2√ε(α - ε)/h, which is independent of the transverse index N.

 

Fig. 5. Competition between (a) oscillator frequency, Ω_v, and (b) the effective “zero of energy” of the harmonic Born-Oppenheimer potential, k _z ²(0), as a function of the quadratic reflection phase coefficient ε. The cavity parameters are the same as in Fig. 4.

Download Full Size | PDF

To gain more physical insight into Fig. 4 (b), one can begin by asking why modes with the same longitudinal quantum number v should be affected differently by the angle-dependent reflectivity (measured by ε), given the fact that the spreading angle θ_p of Eq. (9) is the same for all modes with the same v, assuming standard paraxial theory. In particular, θ_p does not depend on the transverse mode index N, but Fig. 4 (b) shows that modes with different N have different slope in k versus ε.

Qualitative insight into this question can be gained by returning to the analysis of Section 4. There, we point out that θ_p does not measure the details of the angular spectrum, and consequently introduce the angles θ_max, θ_min and θ_rms. The root-mean-square angle is defined through the corresponding harmonic-oscillator length, Eq. (22), which depends on the quantum number N: θ_rms = √(N + 1)/2θ_p. As was noted below Eq. (24), both r _max and θ_max increase with N, and this also holds for θ_rms. Now assume that θ_rms characterizes the plane-wave content of the mode, as it interacts with the mirrors. In paraxial theory, one obtains the modal wave number by quantizing the sum of the dynamical phase along the geometric path length plus the phases due to reflections (φ) and transverse excitation (the Gouy phase), cf. Eq. (11). When we introduce ε < 0, the reflection phase retardation φ(θ_rms) is decreased without changing the other phase contributions. Writing Eq. (25) as φ(ε) = φ₀ + 2εθ_rms ² , this amounts to a shortening of the “effective cavity length” which is more pronounced for larger N because such modes have larger Ω_rms. As a result, the quantized wave number is increased more rapidly for modes with larger N (hence larger θ_rms), as ε is made more negative. While this basic interpretation explains the trend in the left half of Fig. 4 (b), it is not designed to explain the dramatic coalescence in the vicinity of the square-root singularity for ε > 0.

6. Conclusion and outlook

Position and angle of incidence are dual in the sense of Fourier theory, with an uncertainty relation that permits both of them to be sharply peaked only in the ray limit. This limit is not valid in the small cavities we consider here, but one may still hope that Gaussian-beam optics should be sufficient to describe the cavity modes. As we have pointed out, this is no longer true at the level of transverse mode splittings when the boundary conditions involve angle-dependent reflection phase shifts.

To address this problem, a Born-Oppenheimer calculation was carried out, compared to paraxial theory, and then augmented by a local-plane wave approximation. In a straightforward BO treatment for mirrors with constant reflectivity, the transverse mode profile is governed by a fictitious potential which allows us to illuminate the connection between paraxial Gaussian beams and the theory of quantum harmonic oscillators in a novel way. Investigations exploiting this connection [23, 24, 17] within the paraxial context can thereby pushed into size regimes where Gaussian beam theory itself reaches its limitations.

The BO method does not place special emphasis on harmonic approximations, and extending the present work to non-paraxial mirror shapes is thus tantamount to considering anharmonic BO potentials. This also makes it possible to apply our method to the inverse problem of designing specific mirror shapes to create desired transverse mode profiles[25]. The numerical parameters we have chosen for the present study closely match those of Ref. [11]. There, deviations from Gaussian-beam theory were observed for plano-convex microresonators fabricated using template spheres, and were attributed to non-sphericity of the curved mirror. To asses the interplay between mirror boundary conditions and non-ideal mirror shape in these resonators, it appears promising to combine the BO approach with independent measurements of mirror cross sections and reflectivities[10].

Because of the new type of BO potential we arrive at, not all of the paraxial transverse mode splitting is due to the geometry of the cavity and the associated ray orbit; some of it is caused by the dependence of reflection phase shift on angle of incidence, i.e. it originates outside of ray optics. As we have seen, the geometric and boundary-induced contributions to the mode splittings can be disentangled to a certain extent, but in the BO picture they both manifest themselves by changing the effective harmonic-oscillator frequency for the transverse modes. If we re-define the term “Gouy phase” to account for the additional phase shift due to angle-dependent reflectivities, the result is no longer strictly a geometric constant, but instead depends on the modes under consideration.

The parameter ε, which measures the sensitivity of the boundary phase shift to incident angle, affects this modified Gouy phase significantly, especially at short cavity length. In characterizing this dependence, we also considered values of ε corresponding to metamaterials with negative refractive index, and found a spectral coalescence where the transverse mode splitting decreases while the spatial mode confinement becomes stronger.

Our treatment of ε as an arbitrary parameter allowed us to give expressions for the BO spectrum in which the specific mirror design can be chosen freely. This is of great potential interest for the design of VCSELs with large numerical aperture, where both dielectric layer structure [26] and mirror shape[27] have previously been exploited for transverse mode control. The present work constitutes a novel approach to combine these two design aspects. It has recently been observed that oval lateral side walls in broad-area VCSELS can lead to transverse mode structure that bears no resemblance to Gaussian-beam theory at all[28]. To a first approximation, this can be related to ray chaos in the transverse plane, and our method provides a way to couple this phenomenon to the longitudinal degree of freedom.

It is possible, by way of generalizing the present calculation, to allow variable reflectivity for the top mirror as well; one only needs to know the local normal vector n(x,y) of the curved mirror in order to define the angle of incidence entering the reflectivity. An additional extension of our work will be the inclusion of polarization-dependent reflectivities, in order to apply the BO method to the phenomenon of photonic spin-orbit coupling in paraxial dome resonators[4].