## 1.

## Introduction

The micro-optical resonator is a structure that plays a key role in many photonic devices and can take on several geometries, such as defect sites in photonic crystals,^{1}^{,}^{2} ring waveguides,^{3} laser cavities,^{4}^{,}^{5} microspheres,^{6} capillaries,^{7} disks,^{8} and of particular interest here, the bottle resonator.^{9} The most common bottle resonator configuration consists of a radially enlarged region on a solid or hollow glass fiber. The radial and axial profile of the enlarged region determines the optical properties of the modes that may circulate in the azimuthal direction.^{10}^{,}^{11} The resulting bottles are usually quite large compared to the wavelength at which they are designed to operate as they are constructed in glass. We examine theoretically the modal properties of a silicon-based bottle resonator configuration in air. The high-index contrast that exists permits a smaller bottle radius and is on a size scale compatible with silicon photonic devices. For this article, we have chosen the silicon bottle enlargement to have a Gaussian profile in the axial direction. The structure is analyzed using the cylindrical coordinate representation of the full vector wave equation cast as an eigenvalue problem by making use of Fourier–Bessel (FB) basis functions to represent the fields and inverse dielectric profile. The resonator frequencies (wavelengths) and field profiles are computed and presented for various azimuthal mode orders.

There are a number of numerical techniques available for the calculation of the frequencies and field profiles of resonator states. In this article, we are concerned with whispering-gallery modes in a cylindrically symmetric dielectric structure. If the dielectric profile is such that the separation of variable techniques can be applied, Maxwell’s wave equation may be directly solved and field continuity applied at the various interfaces.^{12} In the case involving more general dielectric profiles, one usually resorts to either the finite-difference-time-domain (FDTD) or finite-element-method (FEM) techniques.^{13}^{,}^{14} The dielectric structure is decomposed into a three-dimensional grid and the modes are obtained. The computer resources and computation time using these techniques scale with the complexity of structure under analysis. In a number of instances, the dielectric profile examined is such that simplifications in the analysis equations can be performed resulting in a numerical approach that provides the states in a more time efficient manner. The technique presented here is general and can be placed on the same level of complexity as FDTD and FEM. As will be discussed, the present cylindrical symmetry permits a significant reduction in computer resources and computation time.

In the next section, the key mathematical steps leading to the eigen-matrix for the resonators expressed in cylindrical coordinates and using a FB basis function set are provided. Section 3 presents the geometry of the micro-optic silicon bottle resonator and examines the expansion coefficient space related to FB basis. In Sec. 4, the eigenvalue state space and corresponding field profiles, determined through the eigenvectors, are represented. The dominant features of the modes are provided. In Sec. 5, the hollow core silicon bottle resonator is examined as a design modification to the solid core bottle. Appendix A provides a set of equations required to populate the eigen-matrix.

## 2.

## General Mode Solver Development

The mode solver suitable for obtaining the frequencies and field profiles of the states supported by a cylindrically symmetric dielectric structure, “bottle configuration of interest here,” is based on solving Maxwell’s wave equations using an eigenvalue formulation. When the curl equations are combined for a charge-free, current-free, and nonmagnetic medium with a time variation expressed as ${e}^{-j\omega t}$, the equations are

## (1)

$$\nabla \times (\frac{1}{{\epsilon}_{r}}\nabla \times \overrightarrow{H})={\left(\frac{\omega}{c}\right)}^{2}\overrightarrow{H},\phantom{\rule[-0.0ex]{2em}{0.0ex}}\frac{1}{{\epsilon}_{r}}\nabla \times (\nabla \times \overrightarrow{E})={\left(\frac{\omega}{c}\right)}^{2}\overrightarrow{E}.$$In Eq. (1), $c$ is the free-space speed of light and ${\epsilon}_{r}$ is the relative dielectric constant. As for the familiar plane-wave analysis for photonic crystals,^{15} a set of basic functions is selected from which a series expansion of the inverse dielectric and field components are proposed. In cylindrical space, the basic functions selected are

## (2)

$${J}_{o}\left({\rho}_{p}\frac{r}{R}\right){e}^{jm\phi}{e}^{j{G}_{n}z}={F}_{o}(\text{\hspace{0.17em}}).$$For the radial coordinate, only the lowest order Bessel function is needed as a series expansion in ${J}_{o}(\text{\hspace{0.17em}})$ can represent any function over a finite interval $0\le r\le R$.^{16} The radial expansion index, $p$, selects the $p$’th zero of the Bessel function, ${\rho}_{p}$. Within the interval for $r$, the lowest order Bessel functions are orthogonal. In the azimuthal direction, the complex exponential function is employed as it is orthogonal over the $2\pi $ interval with $m$ as an integer, $-\infty \le m\le +\infty $. In the axial direction, a complex exponential is also selected. Here, it is assumed that the dielectric is periodic in $z$ with a period $T$. The ${G}_{n}$ in the exponential is equivalent to reciprocal lattice “vectors” (Ref. 17) for the axial direction (thus scalar) and take on the form $n(2\pi /T)$ with $-\infty \le n\le +\infty $. The axial basis function contributions are orthogonal over the periodic interval $T$. As the development to follow will involve rather lengthy expressions, the basis function notation is simplified to ${F}_{o}(\text{\hspace{0.17em}})$, which requires three indices $(p,m,n)$ and within the bracket, a parameter indicating if the basic function is related to the inverse dielectric $[\mathrm{\Omega}=(1/{\epsilon}_{r})]$ or field component $(r\to {H}_{r},\phi \to {H}_{\phi},z\to {H}_{z})$. The lower subscript indicates the order of the Bessel function. Radial derivatives give rise to the first-order Bessels and the associated basis function notation is ${F}_{1}(\text{\hspace{0.17em}})$. With these basis functions, the inverse dielectric expansion is

## (3)

$$\mathrm{\Omega}=\sum _{{p}_{\mathrm{\Omega}}{m}_{\mathrm{\Omega}}{n}_{\mathrm{\Omega}}}{\kappa}_{{p}_{\mathrm{\Omega}}{m}_{\mathrm{\Omega}}{n}_{\mathrm{\Omega}}}^{\mathrm{\Omega}}{J}_{o}\left({\rho}_{{p}_{\mathrm{\Omega}}}\frac{r}{R}\right){e}^{j{m}_{\mathrm{\Omega}}\phi}{e}^{j{G}_{{n}_{\mathrm{\Omega}}}z}=\sum _{\mathrm{\Omega}}{\kappa}^{\mathrm{\Omega}}{F}_{o}(\mathrm{\Omega}),$$## (4)

$${H}_{i}=\sum _{i}{\kappa}^{i}{J}_{o}\left({\rho}_{{p}_{i}}\frac{r}{R}\right){e}^{j{m}_{i}\phi}{e}^{j{G}_{{n}_{1}}z}=\sum _{i}{\kappa}^{i}{F}_{o}(i).$$Here, $i$ is used to represent the field component and the summation is taken over the indices related to that particular field component. The next step is to expand the double curl of Eq. (1) and form three separate equations, one for each of the field components specified by the RHS field. The derivative expressions related to $({H}_{r},{H}_{\phi},{H}_{z})$ are

## (5)

$$\frac{\partial \mathrm{\Omega}}{\partial \phi}(\frac{{H}_{\phi}}{{r}^{2}}+\frac{1}{r}\frac{\partial {H}_{\phi}}{\partial r}-\frac{1}{{r}^{2}}\frac{\partial {H}_{r}}{\partial \phi})+\frac{\partial \mathrm{\Omega}}{\partial z}(\frac{\partial {H}_{z}}{\partial r}-\frac{\partial {H}_{r}}{\partial z})\phantom{\rule{0ex}{0ex}}+\mathrm{\Omega}(\frac{{\partial}^{2}{H}_{z}}{\partial z\partial r}-\frac{{\partial}^{2}{H}_{r}}{\partial {z}^{2}})+\frac{\mathrm{\Omega}}{r}\left[\frac{1}{r}\right(\frac{\partial {H}_{\phi}}{\partial \phi}-\frac{{\partial}^{2}{H}_{r}}{\partial {\phi}^{2}})+\frac{{\partial}^{2}{H}_{\phi}}{\partial \phi \partial r}]\phantom{\rule{0ex}{0ex}}={\left(\frac{\omega}{c}\right)}^{2}{H}_{r},$$## (6)

$$\frac{\partial \mathrm{\Omega}}{\partial r}(-\frac{{H}_{\phi}}{r}-\frac{\partial {H}_{\phi}}{\partial r}+\frac{1}{r}\frac{\partial {H}_{r}}{\partial \phi})+\frac{\partial \mathrm{\Omega}}{\partial z}(\frac{1}{r}\frac{\partial {H}_{z}}{\partial \phi}-\frac{\partial {H}_{\phi}}{\partial z})\phantom{\rule{0ex}{0ex}}-\mathrm{\Omega}(\frac{{\partial}^{2}{H}_{\phi}}{\partial {r}^{2}}+\frac{{\partial}^{2}{H}_{\phi}}{\partial {z}^{2}})+\frac{\mathrm{\Omega}}{r}(\frac{1}{r}{h}_{\phi}-\frac{\partial {H}_{\phi}}{\partial r}-\frac{1}{r}\frac{\partial {H}_{r}}{\partial \phi}+\frac{{\partial}^{2}{H}_{r}}{\partial \phi \partial r}+\frac{{\partial}^{2}{H}_{z}}{\partial \phi \partial z})\phantom{\rule{0ex}{0ex}}={\left(\frac{\omega}{c}\right)}^{2}{H}_{\phi},$$## (7)

$$\frac{\partial \mathrm{\Omega}}{\partial r}(\frac{\partial {H}_{r}}{\partial z}-\frac{\partial {H}_{z}}{\partial r})+\frac{1}{r}\frac{\partial \mathrm{\Omega}}{\partial \phi}(-\frac{1}{r}\frac{\partial {H}_{z}}{\partial \phi}+\frac{\partial {H}_{\phi}}{\partial z})\phantom{\rule{0ex}{0ex}}+\mathrm{\Omega}[-\frac{{\partial}^{2}{H}_{z}}{\partial {r}^{2}}+\frac{{\partial}^{2}{H}_{r}}{\partial r\partial z}+\frac{1}{r}(-\frac{\partial {H}_{z}}{\partial r}+\frac{\partial {H}_{r}}{\partial z}+\frac{{\partial}^{2}{H}_{\phi}}{\partial \phi \partial z})\phantom{\rule{0ex}{0ex}}-\frac{1}{{r}^{2}}\frac{{\partial}^{2}{H}_{\phi}}{\partial {\phi}^{2}}]={\left(\frac{\omega}{c}\right)}^{2}{H}_{z}.$$When the derivatives of the series expansions for the field components and inverse dielectric are performed, the following equations are obtained:

## (8)

$$\sum _{r\mathrm{\Omega}}{\kappa}^{r}{\kappa}^{\mathrm{\Omega}}(\frac{{m}_{r}{m}_{\mathrm{\Omega}}+{m}_{r}^{2}}{{r}^{2}}+{G}_{{n}_{r}}{G}_{{n}_{\mathrm{\Omega}}}+{G}_{{n}_{r}}^{2}){F}_{o}(r){F}_{o}(\mathrm{\Omega})+\sum _{\phi \mathrm{\Omega}}{\kappa}^{\phi}{\kappa}^{\mathrm{\Omega}}\left[j\right(\frac{{m}_{\phi}+{m}_{\mathrm{\Omega}}}{{r}^{2}}){F}_{o}(\phi )-j(\frac{{\rho}_{{p}_{\phi}}}{rR})({m}_{\phi}+{m}_{\mathrm{\Omega}}){F}_{1}(\phi )]{F}_{o}(\mathrm{\Omega})\phantom{\rule{0ex}{0ex}}+\sum _{z\mathrm{\Omega}}{\kappa}^{z}{\kappa}^{\mathrm{\Omega}}(-j\frac{{\rho}_{{p}_{z}}}{r})({G}_{{n}_{\mathrm{\Omega}}}+{G}_{{n}_{z}}){F}_{1}(z){F}_{o}(\mathrm{\Omega})={\left(\frac{\omega}{c}\right)}^{2}{H}_{r}.$$## (9)

$$\sum _{r\mathrm{\Omega}}{\kappa}^{r}{\kappa}^{\mathrm{\Omega}}\left\{\right[(-j\frac{{m}_{r}}{{r}^{2}}){F}_{o}(r)+(-j\frac{{m}_{r}{\rho}_{{p}_{r}}}{rR}){F}_{1}(r)]{F}_{o}(\mathrm{\Omega})+[(-j\frac{{m}_{r}{\rho}_{{p}_{\mathrm{\Omega}}}}{rR}){F}_{o}(r)]{F}_{1}(\mathrm{\Omega})\}\phantom{\rule{0ex}{0ex}}+\sum _{\phi \mathrm{\Omega}}{\kappa}^{\phi}{\kappa}^{\mathrm{\Omega}}\{[({G}_{{n}_{\phi}}{G}_{{n}_{\mathrm{\Omega}}}+{G}_{{n}_{\phi}}^{2}+\frac{1}{{r}^{2}}+\frac{{\rho}_{{p}_{\phi}}^{2}}{{R}^{2}}){F}_{o}(\phi )+(\frac{{\rho}_{{p}_{\phi}}}{rR}-\frac{{\rho}_{{p}_{\phi}}}{rR}){F}_{1}(\phi )]{F}_{o}(\mathrm{\Omega})+\left[\right(\frac{{\rho}_{{p}_{\phi}}}{rR}){F}_{o}(\phi )-(\frac{{\rho}_{{p}_{\phi}}{\rho}_{{p}_{\phi}}}{{R}^{2}}){F}_{1}(\phi )]{F}_{1}(\mathrm{\Omega})\}\phantom{\rule{0ex}{0ex}}+\sum _{z\mathrm{\Omega}}{\kappa}^{z}{\kappa}^{\mathrm{\Omega}}(-\frac{{m}_{z}{G}_{{n}_{\mathrm{\Omega}}}}{r}-\frac{{m}_{z}{G}_{{n}_{z}}}{r}){F}_{o}(z){F}_{o}(\mathrm{\Omega})=\left(\frac{\omega}{c}\right){H}_{\phi},$$## (10)

$$\sum _{r\mathrm{\Omega}}{\kappa}^{r}{\kappa}^{\mathrm{\Omega}}\{\left[\right(j\frac{{G}_{{n}_{r}}}{r}){F}_{o}(r)+(-j\frac{{G}_{{n}_{r}}{\rho}_{{p}_{r}}}{R}){F}_{1}(r)]{F}_{o}(\mathrm{\Omega})+\left[\right(-j\frac{{G}_{{n}_{r}}{\rho}_{{p}_{\mathrm{\Omega}}}}{R}){F}_{o}(r)]{F}_{1}(\mathrm{\Omega})\}\phantom{\rule{0ex}{0ex}}+\sum _{\phi \mathrm{\Omega}}{\kappa}^{\phi}{\kappa}^{\mathrm{\Omega}}[(-\frac{{m}_{\mathrm{\Omega}}{G}_{{n}_{\phi}}}{r}-\frac{{m}_{\phi}{G}_{{n}_{\phi}}}{r}){F}_{o}(\phi ){F}_{o}(\mathrm{\Omega})]+\sum _{z\mathrm{\Omega}}{\kappa}^{z}{\kappa}^{\mathrm{\Omega}}\{[({\left(\frac{{\rho}_{{p}_{z}}}{R}\right)}^{2}+\frac{{m}_{z}{m}_{\mathrm{\Omega}}}{{r}^{2}}-\frac{{m}_{z}^{2}}{{r}^{2}}){F}_{o}(z)]{F}_{o}(\mathrm{\Omega})\phantom{\rule{0ex}{0ex}}+\left[\right(-\frac{{\rho}_{{p}_{z}}{\rho}_{{p}_{\mathrm{\Omega}}}}{{R}^{2}}){F}_{1}(z)]{F}_{1}(\mathrm{\Omega})\}=\left(\frac{\omega}{c}\right){H}_{z}.$$The set of Eqs. (8)–(10) can be cast into an eigenvalue problem by making use of the orthogonality properties of the basic functions and integrating over the cylindrical volume defining the dielectric profile. The resulting eigen-matrix is cast into the following general form:

## (11)

$$\left[\begin{array}{ccc}{R}_{r}& {\phi}_{r}& {Z}_{r}\\ {R}_{\phi}& {\phi}_{\phi}& {Z}_{\phi}\\ {R}_{z}& {\phi}_{z}& {Z}_{z}\end{array}\right]\left[\begin{array}{c}R\\ \mathrm{\Phi}\\ Z\end{array}\right]={\left(\frac{\omega}{c}\right)}^{2}I\left[\begin{array}{c}R\\ \mathrm{\Phi}\\ Z\end{array}\right].$$The column matrix $[\begin{array}{c}R\\ \mathrm{\Phi}\\ Z\end{array}]$ represents the expansion coefficients collected by field component and structured based on the summation index sequence. The square matrix is segmented into three row blocks and three column blocks. The first row is generated from Eq. (8), the second row from Eq. (9), and the third row from Eq. (10). The expressions required to generate the matrix elements within each block are provided in Appendix A. A fair bit of structuring is required in order to keep the square matrix elements linked to the expansion coefficients of the column matrix. $I$ is the identity matrix. As will be shown in the next sections, symmetry arguments can be used to greatly reduce the order of the matrix in Eq. (11) and can be solved on a desktop PC. Note that expressions provided in Appendix A used to produce the matrix elements in Eq. (11) have no simplifications or assumption introduced other than those indicated prior to Eq. (1) in their derivation; therefore, they can be applied to the determination of the steady states of any nonmagnetic, charge, and current free-dielectric structure with overall cylindrical symmetry. The derivation can be slightly modified to include an axial propagation contribution, ${k}_{z}$, by replacing ${e}^{j{G}_{{n}_{i}}z}$ with ${e}^{j({G}_{{n}_{i}}+{k}_{z})z}$ in the field expressions [Eq. (4)] and replacing ${G}_{{n}_{i}}$ by $({G}_{{n}_{i}}+{k}_{z})$ in the matrix element generating expression in Appendix A.

The eigenvalues obtained may be real or complex depending on which type of state they are related to, ${(\frac{\omega}{c})}^{2}=\phantom{\rule{0ex}{0ex}}{W}_{\mathrm{real}}+j{W}_{\mathrm{imag}}$. The real and imaginary parts of the state frequencies can be obtained through the following:

## (12)

$${\omega}_{\mathrm{real}}=\frac{c}{\sqrt{2}}\sqrt{{W}_{\mathrm{real}}+\sqrt{{W}_{\mathrm{real}}^{2}+{W}_{\mathrm{imag}}^{2}}},\phantom{\rule{0ex}{0ex}}{\omega}_{\mathrm{imag}}=\frac{c}{\sqrt{2}}\sqrt{-{W}_{\mathrm{real}}+\sqrt{{W}_{\mathrm{real}}^{2}+{W}_{\mathrm{imag}}^{2}}}.$$The corresponding free space wavelength is determined from

The eigenvectors represent the expansion coefficients of the field profile and may also contain real and complex contributions. These are obtained at the same time as the eigenvalues using, for instance, the $\mathrm{eig}()$ function in MATLAB©. Examination of the eigenvector coefficients and determination of the dominant field contribution provide insight into the state’s properties.

## 3.

## Silicon Bottle Dielectric Profile

The silicon bottle resonator configuration is shown in Fig. 1. The structure can be considered as the contact of a uniform dielectric silicon cylinder of height $T=3\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mu \text{m}$ and radius $D=1.75\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mu \text{m}$ and a silicon Gaussian profile extension centered on the $z$ axis with an amplitude of 0.5 *μ*m (L-D) and $\sigma $ of 0.33 *μ*m. The structure presents a cylindrical symmetry and can be efficiently characterized using the $(r,\phi ,z)$ coordinate system. For computation purposes, the structure is taken as periodic in the $z$ direction with a height sufficiently large such that the resonator bottle modes are well isolated from the axial plane periodic interfaces at $+T/2$ and $-T/2$. The external free-space region extends to a radius of $R=3\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mu \text{m}$ and is sufficient to isolate the bottle resonator modes from the radial computation domain edge. In a previous publication dealing with photonic crystals decomposed using a FB basis functions for the $(r,\theta )$ polar coordinates, we have shown that the rotational symmetry of the dielectric imposes restrictions on the nonzero expansion coefficients for the inverse dielectric.^{18} Only the base symmetry of the dielectric and its integer multiples produce nonzero expansion coefficients. The dielectric structure in Fig. 1 is uniform with respect the azimuthal coordinate $\phi $ and can thus be regarded as having infinite rotational symmetry. This property implies that the decomposition of the dielectric should be independent of the azimuthal coordinate and imposes the condition that the only nonzero inverse dielectric expansion coefficients in Eq. (3) are obtained when ${m}_{\mathrm{\Omega}}=0\xb7\infty =0$. The next index pair would be ${m}_{\mathrm{\Omega}}=\pm 1\xb7\infty $, which is not attained in a finite series expansion.

As indicated above, the expansion coefficients of the inverse dielectric are required in the eigen-matrix. Figure 2 shows a plot of the absolute value for the real part of the expansion coefficients for up to 100 Bessel Functions $(1\le {p}_{\mathrm{\Omega}}\le 100)$ and 101 axial exponential $(-50\le {G}_{{n}_{\mathrm{\Omega}}}\le 50)$ with ${m}_{\mathrm{\Omega}}=0$. As expected, the dominant expansion coefficients are observed for the zero axial order and lowest Bessel order since the structure is primarily a right cylinder. All nonzero expansion coefficients with axial order other than zero define the bottle dielectric region. The range of the indices and number of basis functions needed to accurately represent the inverse dielectric is determined by examining how well the inverse dielectric profile is reproduced through Eq. (3). The use of 100 Bessel and 101 azimuthal functions reproduces the inverse dielectric bottle region with $<1\%$ difference when compared to the original inverse dielectric profile.

## 4.

## Eigen-State Determination

## 4.1.

### Matrix Properties/Order Reduction

Determination of the eigen-state requires that the elements of the matrix be computed using the expressions given in Appendix A. The dielectric decomposition has nonzero expansion coefficients only when ${m}_{\mathrm{\Omega}}=0$ for the structure considered with infinite rotational symmetry. This plays a significant role in simplifying the eigen-matrix and reducing the order of the matrix that needs to be solved. The expressions in the eigen-matrix produce a zero matrix element unless ${m}_{\mathrm{\Omega}}={m}_{\mathrm{LHS}}-{m}_{\mathrm{RHS}}$. If all possible azimuthal mode orders are considered and the resulting matrix is examined, the column and rows for the nonzero elements of one azimuthal mode order are padded with zeros produced from the other azimuthal mode orders. As a result, the matrices for each azimuthal order are independent, which can be extracted into separate submatrices and solved independently. This significantly reduces the order of the matrix to which the eigenvalues and eigenvectors are to be determined. In addition, these submatrices may have the order further reduced by one-half, recognizing that the matrix contribution produced for $({m}_{\mathrm{LHS}},{m}_{\mathrm{RHS}})$ both positives is separable from the matrix contribution with $({m}_{\mathrm{LHS}},{m}_{\mathrm{RHS}})$ both negatives. Both of these submatrices $(+,+)$, $(-,-)$ produce the same eigenvalue space but the eigenvector space is the complex conjugate of each other. The full eigenvector is produced by combining both. As an example of matrix order reduction, using 100 Bessel functions and 101 axial exponentials and requiring the properties, the 50 lowest azimuthal order modes result in a full matrix of order 3,060,300 while through matrix order reduction, the individual azimuthal mode properties can be determined from a matrix with order 30,300, which can be accommodated on a desktop PC. Note: The plane-wave analysis technique applied to this type of dielectric structure would require a unit cell well separated such that the steady states are isolated between cells. A square unit cell would be defined in Fig. 1 rather than a right cylinder. If 100 exponentials were to be used in the expansions for each of the coordinate directions, the matrix would have an order of 3,000,000. Clearly, the FB approach in cylindrical coordinates provides a much more resource efficient and time-managed approach to the steady-state determination process for the cylindrical symmetric dielectric profiles. Using a basis function set related to the symmetry of the dielectric to reduce the order of the eigenmatrix is suggested in Refs. 19 and 20.

## 4.2.

### Mode Solutions by Azimuthal Mode Number

We first consider obtaining the eigen-states for a mode profile that has a rotational symmetry in the azimuthal direction equal to 22. For this mode order, the only possible nonzero field expansion coefficients occur when ${m}_{i}=\pm 22$. The eigen-matrix of the $H$ field for this mode order can be generated using Eq. (11) and the expressions are provided in Appendix A. Sufficiently converged eigenvalues and field profiles are obtained using 75 Bessel terms and 39 axial exponentials for the bottle structure examined here. (Note: The convergence is based on examining the variation in the eigenvalues and eigenvectors of the desired bottle states.) The matrix solution returns 8775 eigenvalues that extend from very small to very large values, real and complex. Within these are contained the desired bottle states. It is expected that the bottle states should have a wavelength located within the scale of the dielectric structure. The eigenvalues [converted to free-space wavelength using Eq. (13)] are plotted in Fig. 3 for the 1 to 3 *μ*m range. The eigenvalues for which the eigenvector expansion coefficient space is dominated by the ${H}_{z}$ field component are plotted in the upper trace. The lower trace is plotted for the remaining eigenvalues dominated by the $({H}_{r},{H}_{\phi})$ field components and are related to the other set of hybrid states and in fact are decoupled from ${H}_{z}$ for a uniform axial dielectric. Solving the E field equation in Eq. (1) and examining the eigenvectors dominated by the ${E}_{z}$ field component expansion coefficients provide the second set of hybrid states. The remainder of the presentation focusses on states dominated by ${H}_{z}$. For the wavelength range depicted in Fig. 3(a), the wavelength states of 2.23 *μ*m and larger (not shown) have eigen frequencies with a large imaginary contribution. These form the set of unphysical states and are similar to physically unrealistic states of slabs and other waveguides.^{21} The imaginary frequency of the states with wavelengths of 1.71 *μ*m and lower are of the order of ${10}^{-9}$ and are effectively zero (numerically). These form the eigen-state of the resonator with several of the states localized in the bottle region and other states extending the full axial length of the dielectric structure. The field profiles, reconstructed using the eigenvectors, are examined below and facilitate the distinction between the two types. It is worth noting that the field expression provided in Eq. (4) makes use of the lowest order Bessel function for the entire radial extent. The radial series expansion is suitable for modeling the field in the high dielectric region and in the low dielectric region even if the low dielectric field region were to be of an “evanescent” nature.

Figure 4 shows the intensity of the ${H}_{z}$ field for the states at wavelengths 1.707, 1.557, 1.507, and 1.426, plotted in the $(r,z)$ plane, top profiles, and in the $(r,\phi )$ plane, lower profiles. The solid line represents the interface between the high and low dielectric regions. Here, the modes are labeled as ${R}_{0}{Z}_{0}$, ${R}_{0}{Z}_{1}$, ${R}_{1}{Z}_{0}$, and ${R}_{1}{Z}_{1}$, which are based on the number of zero crossing of the field in the radial and axial directions. As the $RZ$ null number increases, the fields are less confined to the bottle region typical of optical waveguide solutions in slab and ridge geometries. For the traces in the azimuthal plane, the states are plotted using the axial value that corresponds to the maximum in the field value. Azimuthal mode order of 22 is confirmed as there are 44 high intensity regions per $2\pi $ rotation. In general, the excitation of bottle states is accomplished using an external waveguide (fiber) or a focused laser beam. Knowledge of the mode profile and polarization of the bottle states is crucial in designing a large overlap integral between the incident and supported fields. The FB eigen-approach presented here provides a means of accurately calculating the modes frequency, wavelength, and field profile. A parameter highly desired in resonator design is the quality factor, $Q$. The dielectric structure under analysis has no loss mechanisms presented. Thus, the eigen-frequencies of the bottle states are real and indicate an infinite $Q$. However, in a real application of the bottle resonator, input and output power coupling mechanisms must be included and the $Q$ of the resonator would be finite. An estimation of the $Q$ factor should be possible by taking the field profile of the bottle state of interest and calculating its overlap with the power bleeding environment (waveguide, lossy medium, …).

Six of the modes with wavelengths below 1.426 *μ*m are plotted in Fig. 5. They extend the full axial length and present either an even symmetry (top) or odd symmetry (lower) with respect to the structures’ center and upper/lower axial periodic boundaries. These modes resemble, to a large extent, the guided modes of a ridge waveguides resting on a low dielectric substrate. To a first approximation, whispering-gallery states have been analyzed using the Helmholtz equation with periodic azimuthal boundary conditions.^{22} It is observed that the modes inner caustic, ${R}_{c}$, decreases in relation to a decrease in mode wavelength. This is a manifestation of the conservation of angular momentum for these modes, all with the same azimuthal symmetry order ${m}_{i}$ and angular momentum $L={m}_{i}h=\phantom{\rule{0ex}{0ex}}{R}_{c}\hslash k$. Although these states are available and would produce field profiles contained primarily in the high dielectric region, they would not normally be significantly excited using a tapered fiber. A more efficient way to couple to these states would be, for instance, to adopt an end-fire coupling approach by cutting out a portion of the high dielectric region.

Applying the matrix order reduction technique, the matrices for azimuthal orders 5 to 30 were produced and the state space was determined. From the eigenvectors, the field profiles were computed and those states which are highly contained in the bottle region have their wavelengths plotted in Fig. 6. The lines between data points link states with the same general mode profiles. As the azimuthal order is reduced, the wavelength associated within a mode class increases as well as the field’s extent into the cylindrical dielectric region. The “cut-off” of the modes results when the state is converted into a mode that extends the full axial height of the dielectric structure. For example, Fig. 7 shows the mode profile for the fundamental state for several azimuthal orders. As the azimuthal order is decreased (wavelength increased), the state extends further into the cylindrical dielectric region. At azimuthal order of 4, the fundamental state extends the full axial range and is no longer considered as confined to the bottle dielectric region.

The solid light blue line (top most line in figure) in Fig. 6 is produced by considering the fundamental mode in the geometrical regime, ${R}_{0}{Z}_{0}$ geometrical. Based on the circumference, $C$, of the bottle resonator and azimuthal mode order, the free-space wavelengths that are plotted where calculated from $\lambda =C{n}_{Si}/{m}_{i}$. Good agreement between geometrical and wave optics is observed when the wavelength is small as is the general trend in the traces shown in Fig. 6. The ratio of the ${R}_{0}{Z}_{0}$ FB-computed wavelength to the ${R}_{0}{Z}_{0}$ geometrical wavelength multiplied by the index of refraction of silicon gives the ${R}_{0}{Z}_{0}$ mode’s effective index.

The effect of reducing the amplitude of the Gaussian profile making up the bottle resonator region is shown in Fig. 8 for the mode with azimuthal order of 22. As the amplitude decreases, the modes become less and less confined and extend further into the solid dielectric region. At sufficiently low amplitude, the mode may be cut-off and extend the full axial direction. The figure further shows that by adjusting the amplitude, the wavelength of the fundamental mode can be tuned to that of 1.55 *μ*m. A similar type of curve can be generated for the other azimuthal mode orders and through a combination of selecting the amplitude and $\sigma $ of the Gaussian profile, the wavelength and number of bottle modes can be adjusted.

The multimode nature of the bottle confined states permits the specification of two types of free-spectral ranges;^{23}24.25.^{–}^{26} wavelength spacing of similar states with adjacent azimuthal order; and wavelength spacing of adjacent modes within an azimuthal order. Figure 9 shows the free-spectral ranges computed for the fundamental mode, ${R}_{0}{Z}_{0}$, plotted on a log–log trace. The solid black line is the free-spectral range computed using geometrical optics for the ${R}_{0}{Z}_{0}$ mode. For the Gaussian bottle parameters chosen, the FB analysis indicates that both free-spectral ranges can be several hundred nanometers.

## 5.

## Hollow Silicon Bottle Resonator

The demonstration of hollow core rolled up semiconductor structures can be traced to the early work reported by Prinz et al.^{27}^{,}^{28} The initial structures were designed in the nanometer scale diameter regime and produced by allowing the strain mismatch in the freed layers to provide a means of rolling the sheet into a hollow structure. The technique has been extended into the micron-scale optical domain by increasing the diameter of the hollow tubes and through the utilization of semiconductor materials common employed in integrated optic devices.^{29}^{,}^{30} A great deal of effort is directed toward the determination of the modal properties that may circulate the circumference of the hollow tube.^{31}^{,}^{32} In this section of the article, the dependence of the wavelength of the fundamental mode, azimuthal order 22, with respect to the inner radius of the air-filled hollow silicon bottle resonator is examined. The hollow internal region provides an access port for fluid flow^{33}^{,}^{34} and for mechanical positioning.^{35} The Gaussian shape bottle resonator region is retained in the analysis and the presence of an inner radius equal to 1.75 *μ*m completely removes the supporting cylindrical region. Figure 10 is produced by modifying the dielectric profile, recomputing the inverse dielectric expansion coefficients and then computing the eigenvalues for the fundamental mode of azimuthal order 22. The wavelength of the mode is plotted versus hollow core radius. The inner air region starts to influence the mode’s wavelength (and profile) when the radius of the hollow region approaches and exceeds that of the original mode’s inner caustic. The trace was produced using 75 Bessel and 39 axial indices. The fluctuations on the trace (amplified here due to the small $y$ axis range) represent $<1\%$ of the wavelength value.

## 6.

## Conclusion

In this work, we have examined the steady states that are available in a micro-optic silicon optical bottle resonator. The structures are examined using a novel FB basis function set to expand the inverse dielectric and field components in order to cast Maxwell’s wave equations into an eigenvalue formulation. The expressions are provided in a general form and may be applied to solving for steady states in a large number of cylindrically symmetric dielectric configurations. Examination of the eigenvalue space and dominant field component of the associated eigenvector facilitates the determination of bottle confined states and additional extended states. For the parameters of the Gaussian bottle resonator presented, four bottle confined eigen-states are observed when the azimuthal mode order is 22. As the azimuthal order is reduced, the modal wavelengths increase and the extention of the field profile into the high dielectric region increases. A cut-off is defined when the bottle state transforms from one confined to the Gaussian region to one that extends the full axial direction. A hollow variation of the bottle resonator is examined and it is shown that the fundamental mode is unaltered by the air region provided the modal field remains confined to the high dielectric region. The analysis technique developed here along with the micro-optic bottle resonator presented in silicon indicates a rich variety of confined and extended states. Variations on the structure parameters promise interesting micro-optical components and devices on the Silicon platform.

## Appendices

## Appendix A:

### Generating Expressions for the Matrix Elements of Eq. (11)

When the orthogonality integration is applied to Eqs. (8)–(10), the resulting expressions written below are obtained. The indices $({p}_{i}^{*},{m}_{i}^{*},{n}_{i}^{*})$ relate to the application of the complex conjugate of a basis function on the right and left hand side prior to integration of the cylindrical domain.

## (14)

$${R}_{r}=\frac{2}{{R}^{2}{J}_{1}^{2}({\rho}_{{p}_{{r}^{*}}})}\{\sum _{\mathrm{\Omega}}{\kappa}^{r}{\kappa}^{\mathrm{\Omega}}([{m}_{r}{m}_{\mathrm{\Omega}}+{m}_{r}^{2}]{S}_{-1}({p}_{r},{p}_{\mathrm{\Omega}},{p}_{{r}^{*}})+[{G}_{{n}_{r}}{G}_{{n}_{\mathrm{\Omega}}}+{G}_{{n}_{r}}^{2}]{R}^{2}{S}_{1}({p}_{r},{p}_{\mathrm{\Omega}},{p}_{{r}^{*}}))\Vert \begin{array}{c}{m}_{r}+{m}_{\mathrm{\Omega}}={m}_{{r}^{*}}\\ {n}_{r}+{n}_{\mathrm{\Omega}}={n}_{{r}^{*}}\end{array}\Vert \}{|}_{{p}_{r}{m}_{r}{n}_{r}}$$## (15)

$${\phi}_{r}=\frac{2j}{{R}^{2}{J}_{1}^{2}({\rho}_{{p}_{{r}^{*}}})}\{\sum _{\mathrm{\Omega}}{\kappa}^{\phi}{\kappa}^{\mathrm{\Omega}}({m}_{\mathrm{\Omega}}+{m}_{\phi})[{S}_{-1}({p}_{\phi},{p}_{\mathrm{\Omega}},{p}_{{r}^{*}})+{\rho}_{{p}_{\phi}}{T}_{o}({p}_{\phi},{p}_{\mathrm{\Omega}},{p}_{{r}^{*}})]\Vert \begin{array}{c}{m}_{\phi}+{m}_{\mathrm{\Omega}}={m}_{{r}^{*}}\\ {n}_{\phi}+{n}_{\mathrm{\Omega}}={n}_{{r}^{*}}\end{array}\Vert \}{|}_{{p}_{\phi}{m}_{\phi}{n}_{\phi}}$$## (16)

$${Z}_{r}=\frac{-2j}{R{J}_{1}^{2}({\rho}_{{p}_{{r}^{*}}})}\{\sum _{\mathrm{\Omega}}{\kappa}^{z}{\kappa}^{\mathrm{\Omega}}[({\rho}_{{p}_{z}})[{G}_{{n}_{\mathrm{\Omega}}}+{G}_{{n}_{z}}]{T}_{1}({p}_{z},{p}_{\mathrm{\Omega}},{p}_{{r}^{*}})]\Vert \begin{array}{c}{m}_{z}+{m}_{\mathrm{\Omega}}={m}_{{r}^{*}}\\ {n}_{z}+{n}_{\mathrm{\Omega}}={n}_{{r}^{*}}\end{array}\Vert \}{|}_{{p}_{z}{m}_{z}{n}_{z}}$$## (17)

$${R}_{\phi}=\frac{-2j}{{R}^{2}{J}_{1}^{2}({\rho}_{{p}_{{\phi}^{*}}})}\{\sum _{\mathrm{\Omega}}{\kappa}^{r}{\kappa}^{\mathrm{\Omega}}{m}_{r}[{S}_{-1}({p}_{r},{p}_{\mathrm{\Omega}},{p}_{{\phi}^{*}})+{\rho}_{{p}_{r}}{T}_{o}({p}_{r},{p}_{\mathrm{\Omega}},{p}_{{\phi}^{*}})+{\rho}_{{p}_{\mathrm{\Omega}}}{U}_{o}({p}_{r},{p}_{\mathrm{\Omega}},{p}_{{\phi}^{*}})]\Vert \begin{array}{c}{m}_{r}+{m}_{\mathrm{\Omega}}={m}_{{\phi}^{*}}\\ {n}_{r}+{n}_{\mathrm{\Omega}}={n}_{{\phi}^{*}}\end{array}\Vert \}{|}_{{p}_{r}{m}_{r}{n}_{r}}$$## (18)

$${\phi}_{\phi}=\frac{2}{{J}_{1}^{2}({\rho}_{{p}_{{\phi}^{*}}})}\left\{\sum _{\mathrm{\Omega}}{\kappa}^{\phi}{\kappa}^{\mathrm{\Omega}}\left[\begin{array}{l}({G}_{{n}_{\phi}}{G}_{{n}_{\mathrm{\Omega}}}+{G}_{{n}_{\phi}}^{2}+\frac{{\rho}_{{p}_{\phi}}^{2}}{{R}^{2}}){S}_{1}({p}_{\phi},{p}_{\mathrm{\Omega}},{p}_{{\phi}^{*}})+\frac{1}{{R}^{2}}{S}_{-1}({p}_{\phi},{p}_{\mathrm{\Omega}},{p}_{{\phi}^{*}})\\ +\left(\frac{{\rho}_{{p}_{\mathrm{\Omega}}}}{{R}^{2}}\right){U}_{o}({p}_{\phi},{p}_{\mathrm{\Omega}},{p}_{{\phi}^{*}})-(\frac{{\rho}_{{p}_{\phi}}{\rho}_{{p}_{\mathrm{\Omega}}}}{{R}^{2}}{V}_{1}({p}_{\phi},{p}_{\mathrm{\Omega}},{p}_{{\phi}^{*}}))\end{array}\right]\Vert \begin{array}{c}{m}_{\phi}+{m}_{\mathrm{\Omega}}={m}_{{\phi}^{*}}\\ {n}_{\phi}+{n}_{\mathrm{\Omega}}={n}_{{\phi}^{*}}\end{array}\Vert \right\}{|}_{{p}_{\phi}{m}_{\phi}{n}_{\phi}}$$## (19)

$${Z}_{\phi}=\frac{-2}{R{J}_{1}^{2}({\rho}_{{p}_{{\phi}^{*}}})}\{\sum _{\mathrm{\Omega}}{\kappa}^{z}{\kappa}^{\mathrm{\Omega}}{m}_{z}([{G}_{{n}_{\mathrm{\Omega}}}+{G}_{{n}_{z}}]{S}_{o}({p}_{\phi},{p}_{\mathrm{\Omega}},{p}_{{\phi}^{*}}))\Vert \begin{array}{c}{m}_{z}+{m}_{\mathrm{\Omega}}={m}_{{\phi}^{*}}\\ {n}_{z}+{n}_{\mathrm{\Omega}}={n}_{{\phi}^{*}}\end{array}\Vert \}{|}_{{p}_{z}{m}_{z}{n}_{z}}$$## (20)

$${R}_{z}=\frac{2j}{R{J}_{1}^{2}({\rho}_{{p}_{{z}^{*}}})}\left\{\sum _{\mathrm{\Omega}}{\kappa}^{r}{\kappa}^{\mathrm{\Omega}}\right[\begin{array}{l}{G}_{{n}_{r}}{S}_{o}({p}_{r},{p}_{\mathrm{\Omega}},{p}_{{z}^{*}})-{G}_{{n}_{r}}{\rho}_{{p}_{r}}{T}_{1}({p}_{r},{p}_{\mathrm{\Omega}},{p}_{{z}^{*}})\\ -{G}_{{n}_{r}}{\rho}_{{p}_{\mathrm{\Omega}}}{U}_{1}({p}_{r},{p}_{\mathrm{\Omega}},{p}_{{z}^{*}})\end{array}]\Vert \begin{array}{c}{m}_{r}+{m}_{\mathrm{\Omega}}={m}_{{z}^{*}}\\ {n}_{r}+{n}_{\mathrm{\Omega}}={n}_{{z}^{*}}\end{array}\Vert \}{|}_{{p}_{r}{m}_{r}{n}_{r}}$$## (21)

$${\phi}_{z}=\frac{-2}{R{J}_{1}^{2}({\rho}_{{p}_{{z}^{*}}})}\{\sum _{\mathrm{\Omega}}{\kappa}^{\phi}{\kappa}^{\mathrm{\Omega}}{G}_{{n}_{\phi}}([{m}_{\mathrm{\Omega}}+{m}_{\phi}]{S}_{o}({p}_{\phi},{p}_{\mathrm{\Omega}},{p}_{{z}^{*}}))\Vert \begin{array}{c}{m}_{\phi}+{m}_{\mathrm{\Omega}}={m}_{{z}^{*}}\\ {n}_{\phi}+{n}_{\mathrm{\Omega}}={n}_{{z}^{*}}\end{array}\Vert \}{|}_{{p}_{\phi}{m}_{\phi}{n}_{\phi}}$$## (22)

$${Z}_{z}=\frac{2}{{J}_{1}^{2}({\rho}_{{p}_{{z}^{*}}})}\left\{\sum _{\mathrm{\Omega}}{\kappa}^{z}{\kappa}^{\mathrm{\Omega}}\right[\begin{array}{l}{\left(\frac{{\rho}_{{p}_{z}}}{R}\right)}^{2}{S}_{1}({p}_{z},{p}_{\mathrm{\Omega}},{p}_{{z}^{*}})+\left(\frac{{m}_{z}{m}_{\mathrm{\Omega}}+{m}_{z}^{2}}{{R}^{2}}\right){S}_{-1}({p}_{z},{p}_{\mathrm{\Omega}},{p}_{{z}^{*}})\\ -(\frac{{\rho}_{{p}_{z}}{\rho}_{{p}_{\mathrm{\Omega}}}}{{R}^{2}}{V}_{1}({p}_{z},{p}_{\mathrm{\Omega}},{p}_{{z}^{*}}))\end{array}\left]\Vert \begin{array}{c}{m}_{z}+{m}_{\mathrm{\Omega}}={m}_{{z}^{*}}\\ {n}_{z}+{n}_{\mathrm{\Omega}}={n}_{{z}^{*}}\end{array}\Vert \right\}{|}_{{p}_{z}{m}_{z}{n}_{z}}$$For each of the nine matrix element expressions, the conditions on the indices to produce nonzero matrix elements are given within the $\Vert \text{\hspace{0.17em}}\Vert $ bracket. These conditions arise from the orthogonality properties of the exponential portions of the basis functions. The ${|}_{{p}_{i}{m}_{i}{n}_{i}}$ indicates the matrix element based on the LHS field indices. The orthogonality condition of the Bessel function integrals on the LHS is not met and produces finite values. The resulting integrals involve the product of three Bessel functions, indicated as $(S,T,U,V)$ in Eqs. (14) to (22) and are given below:

## (23)

$${S}_{1}({p}_{1},{p}_{2},{p}_{3})={\int}_{0}^{1}{J}_{o}({\rho}_{{p}_{1}}x){J}_{o}({\rho}_{{p}_{2}}x){J}_{o}({\rho}_{{p}_{3}}x)x\mathrm{d}x$$## (24)

$${S}_{0}({p}_{1},{p}_{2},{p}_{3})={\int}_{0}^{1}{J}_{o}({\rho}_{{p}_{1}}x){J}_{o}({\rho}_{{p}_{2}}x){J}_{o}({\rho}_{{p}_{3}}x)\mathrm{d}x$$## (25)

$${S}_{-1}({p}_{1},{p}_{2},{p}_{3})={\int}_{0}^{1}{J}_{o}({\rho}_{{p}_{1}}x){J}_{o}({\rho}_{{p}_{2}}x){J}_{o}({\rho}_{{p}_{3}}x)\frac{1}{x}\mathrm{d}x$$## (26)

$${T}_{1}({p}_{1},{p}_{2},{p}_{3})={\int}_{0}^{1}{J}_{1}({\rho}_{{p}_{1}}x){J}_{o}({\rho}_{{p}_{2}}x){J}_{o}({\rho}_{{p}_{3}}x)x\mathrm{d}x$$## (27)

$${T}_{o}({p}_{1},{p}_{2},{p}_{3})={\int}_{0}^{1}{J}_{1}({\rho}_{{p}_{1}}x){J}_{o}({\rho}_{{p}_{2}}x){J}_{o}({\rho}_{{p}_{3}}x)\mathrm{d}x$$## (28)

$${U}_{o}({p}_{1},{p}_{2},{p}_{3})={\int}_{0}^{1}{J}_{o}({\rho}_{{p}_{1}}x){J}_{1}({\rho}_{{p}_{2}}x){J}_{o}({\rho}_{{p}_{3}}x)\mathrm{d}x$$## (29)

$${V}_{1}({p}_{1},{p}_{2},{p}_{3})={\int}_{0}^{1}{J}_{1}({\rho}_{{p}_{1}}x){J}_{1}({\rho}_{{p}_{2}}x){J}_{o}({\rho}_{{p}_{3}}x)x\mathrm{d}x$$The integrals are evaluated over the normalized range $x=\frac{r}{R}$ and involve various combinations of the zero and first-order Bessel functions as well as different powers of $x$. These integrals, written in the form provided, can be evaluated numerically and the results stored in look-up tables. Since the integrals are independent of the dielectric profile and particular field components required, they only need to be calculated once and can be used in the steady state determination using the expressions presented here for any dielectric profile where the radial direction is expanded in the lowest order Bessel series.

## Appendix B:

### Eigenvector Space for Fundamental Mode of Azimuthal Order 22

Figure 11 shows a plot for the real and imaginary parts of the expansion coefficients for the ${H}_{z}$ field component for the fundamental azimuthal mode of order 22 and wavelength at 1.707 *μ*m. The expansion coefficients are collected in terms of axial index ${n}_{z}$ as labeled along the $x$ axis of the plot. Within each of these are collected the Bessel indices in increasing order up to ${p}_{z}=75$. Such a figure shows that the eigenvector space is well converged.

## Acknowledgments

This work is supported by a grant provided by Natural Sciences and Engineering Research Council of Canada (NSERC).

## References

## Biography

**Robert C. Gauthier** obtained his PhD in physics in 1988, Dalhousie University. He is a registered member of the professional engineers of Ontario and is currently a faculty member in the Department of Electronics at Carleton University. His research interests include theoretical studies photonic crystals, quasi-crystals, and resonators. Theoretical efforts have resulted in the development of 3-D eigen-mode solvers using special basis functions. Past research interests include laser trapping, fiber optic sensors, and integrated optics.