Open Access
Add to BibSonomy
Abstract
A prime objective of modeling optical fibers is capturing mode confinement losses correctly. This paper demonstrates that specific modeling choices, especially regarding the outer fiber cladding regions and the placement of the computational boundary, have significant impacts on the calculated mode losses. This sensitivity of the computed mode losses is especially high for microstructure fibers that do not guide light by total internal reflection. Our results illustrate that one can obtain disparate mode confinement loss profiles for the same optical fiber design simply by moving the boundary to a new material region. We conclude with new recommendations for how to better model these losses.
© 2023 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement
1. Introduction
The primary concern of this effort is exploring the effect of modeling choices on the computed mode confinement losses (CLs) of microstructured optical fibers. Many microstructured fibers guide light by interference and/or anti-resonant reflection effects [1,2]. Such guidance mechanisms are often inherently lossy, and this loss is modeled by finding so-called leaky modes [3]. These modes satisfy an eigenvalue problem derived from Maxwell equations on the fiber’s transverse geometry. The confinement loss can then be found from the imaginary component of the associated eigenvalue. The loss value, as a function of input wavelength, gives the fiber’s spectral loss profile. This work shows that these loss profiles are highly sensitive to various modeling choices, especially those regarding the material in which the boundary of the computational domain is placed.
It seems to be common practice to place this boundary in the innermost cladding region of the fiber, effectively allowing the cladding material to extend to infinity [4, p. 62ff]. The assumption underlying this choice is that any material changes or structures beyond this inner cladding region have little to no effect on the core modes of the given fiber. This reasoning may be well-justified for optical fibers that guide light by total internal reflection, such as step-index or graded-index fiber designs, where the guided mode profiles exponentially decay in the radial direction outside of the guiding region. In contrast, leaky modes do not decay exponentially in the innermost cladding [5]; in fact, they grow exponentially radially far from the guiding region. Therefore, for leaky modes, unlike guided modes, exponential decay cannot be used to justify a model with infinite cladding.
Premature termination of the computational domain can have significant effects on the loss profile of a given fiber. Bragg fibers, which consist of concentric rings of dielectric material, referred to as “anti-resonant layers” in [6], have resonant wavelengths at which their losses spike [7]. The location of these wavelengths depends on the number and thickness of these concentric rings. During modeling, if the computational boundary truncates what would otherwise be a full anti-resonant layer, an entire set of loss peaks may be removed from the spectral profile.
To elucidate the possibility of such effects in modern optical waveguides, we model two anti-resonant fibers (ARFs) from the literature [8,9] and compare their spectral loss profiles for different choices of boundary placement. We will see that placing the computational boundary beyond the cladding (in a lower index ambient material) results in mode loss profiles with extreme variability, with losses that vary by several orders of magnitude over small changes in wavelength. In contrast, the loss profiles of these same fibers are smooth and less variable over wide bandwidths of wavelengths when the boundary of the computational domain is placed in a high-index material such as the cladding. It is then observed that including a lossy material layer in the mode solver model, such as the polymer jacket that most fibers possess, leads to intermediate solutions between the highly variable and smooth loss profiles, suggesting that knowledge of the material loss properties of fiber coatings may be important for accurate model predictions.
Such modeling is of particular importance for microstructured fibers. The microstructures of the anti-resonant fibers consists of glass capillary tubes attached at regular angular intervals to the inner wall of a glass cladding that encompasses the guiding (e.g. core) region (see Fig. 1). Many use mode modeling to predict optical fiber performances, including how the mode confinement losses change, under various geometry modifications [8, Figs. 6, 8–9]. In this work, a numerical study that varies “embedding distance” of the capillary tubes into their supporting cladding, i.e. the portion of the tube wall that overlaps the cladding (see Fig. 1 inset), is used to demonstrate that the computed mode losses are highly dependent on the outer material configurations of the model. This embedding distance is not entirely physical because the amount of glass material is not conserved. The main point is that the geometry is being altered by small distances, near or less than the optical field wavelength ($\lambda$), and yet, as our study shall show, the outer material conditions, including the placement of the computational boundary, that are many tens of wavelengths away from these geometry alterations strongly impact the confinement loss of the fundamental mode.
Outline of the remainder of this paper is as follows. In Section 2, we leverage semi-analytic solutions for Bragg fibers to establish the simulation effects of outer material configuration on the core modes, and to verify the accuracy of our numerical mode solver. Then, in Section 3, two anti-resonant fibers from the literature [8,9] are examined, showing that their spectral mode loss curves are highly sensitive to different choices of computational boundary placement. Next, in Section 4, we repeat these studies incorporating an additional lossy dielectric (polymer) layer into the model. Section 5. displays the sensitivity of mode CL to the capillary tube embedding distance for all choices of outer material configurations. Exploiting the understanding gained from these preceding studies, it is demonstrated in Section 6. that informed meshing practices can accelerate the convergence of the numerical results. The implications of this study and our recommendations for modeling microstructured optical fibers are discussed in Section 7. Finally, Appendix A. describes the eigenproblem that is solved to produce the mode profiles and their corresponding propagation constants and loss values. This includes the description of our numerical approach, which formulates the mode solving into a generalized eigenvalue problem. We utilize the FEAST algorithm, which provides multiple propitious qualities for mode finding and achieving mathematically converged results.
2. Bragg fiber and numerical verification
In this section, we establish the impact of changing the simulated outer fiber configuration on the fundamental mode (FM) CL in Bragg fibers. Note that Bragg fibers have been used as reference fibers for understanding the guiding mechanism in anti-resonant fibers (ARFs) [2], and provide analytical spectral mode loss profiles, against which one can ensure the computational validity of their numerical mode solving method.
It was shown in [10] that a waveguide with a low index core region could support lossless modes provided that the cladding consists of periodic high and low index dielectric materials. The Bragg fiber is an optical fiber design that implements this insight. It is radially symmetric (azimuthally invariant) and formed of concentric rings (annuli) of homogeneous material, with the final ambient material around the fiber extending to infinity. The necessity of terminating the infinitely periodic cladding of the theoretical fiber introduces some loss. The analytic calculations used here for finding this loss are based on the matrix methodology established in the 1978 publication [11].
This work utilizes the terminology for Bragg fiber designs developed in [6], where $N_i$ denotes a Bragg fiber with $i$ number of complete (or finite-width) anti-resonant layers. Thus, the $N_0$ fiber is a hollow-core surrounded by an infinite high-index material, and the $N_1$ Bragg fiber design is a glass ring suspended in low index material (see Fig. 2), et cetera. Figure 3 depicts the 1-D cross sections of the refractive index profile for the first four Bragg fiber designs.
Fig. 2. $N_1$ Bragg fiber fundamental mode (FM) vector component profiles and their 1-D horizontal cross-sections.
In our studies, the low index material is air with refractive index of $n_{\textrm{air}} = 1.00028$ and the high-index material is fused silica glass with a refractive index of $n_{\textrm{glass}} = 1.44$. The core radius is $R_{\text {core}} = 40.775 \, \mu\textrm{m}$ and thickness of the glass ring is $t = 10 \, \mu\textrm{m}$. For the $N_2$ Bragg fiber, the thickness of the complete air layer beyond the glass ring was set equal to the core radius. The $N_3$ fiber had the same inner measurements as the $N_2$ fiber, but with a second complete layer of glass of the same thickness as the inner glass layer after the second ring of air. Our numerical studies, not shown here, demonstrated that chromatic dispersion does not affect any of our conclusions to this effort, and since including chromatic dispersion necessitates an extra numerical solve in relation (1) without adversely affecting this effort, we omit it from this work. Using these fiber parameters, Fig. 2 depicts the FM profile of the $N_1$ Bragg fiber. Note how the glass ring supports high-frequency ripples/oscillations in the mode profile vector components, which illustrates a non-negligible effect that this entire finite cladding layer (anti-resonant layer) has on the mode, an effect that would be lost if the computational boundary were to be placed in this layer (material). This effect will appear again in when considering realistic ARF designs.
Figure 4 plots the FM spectral loss profiles for several Bragg fibers with different numbers of complete anti-resonant layers. Numerical results (green dots) are in excellent agreement with the curves, verifying that our numerical method aligns with the semi-analytic Bragg fiber theory predictions.
Fig. 4. Spectral loss profiles of Bragg fibers; the numerically computed results (dots) are in good agreement with the analytical results (curves; from Bragg fiber theory).
The $N_0$ fiber, consisting of a low index core surrounded by an ambient high-index material (in which the mode eigenproblem’s computational boundary starts $10 \, \mu\textrm{m}$ away from the low-high index transition), does not exhibit (resonant) spikes in the CL profile, unlike for the $N_1$ and $N_2$, and $N_3$ Bragg fiber designs. This is in accordance with the Bragg fiber theory predictions, which indicate that the presence of complete anti-resonant layers are needed for these resonant loss peaks. These loss resonances (spikes) occur approximately at
Note that one effectively models the $N_0$ Bragg fiber by beginning with the $N_1$ or $N_2$ fiber designs, but placing the computational boundary inside of the first (innermost) high-index anti-resonant layer. This choice of boundary placement erases any information about the outer material layers of the $N_1$ or $N_2$ designs and thus models the $N_0$ design. Figure 4 indicates that the FM loss profile for the $N_0$ Bragg fiber is smooth over a broad set of operational wavelengths, and quite distinct from the Bragg designs with complete anti-resonant layers. This indicates that capturing the outer layers of a fiber design (not guided by total internal reflection) is essential to computing the mode CLs correctly.
The results above were for a general Bragg fiber with dimensions chosen to clarify the structure of the loss profiles and modes. Now we will consider $N_2$ and $N_3$ Bragg fibers that are comparable to the more complex 6- and 8-tube ARFs (recall Fig. 1), chosen as examples for detailed study in later sections. This comparison approximates the portions of the ARF capillaries that surround the core region as if they form a circular high-index ring around the core, and simultaneously neglects the remaining portions of the capillaries such that there would be an air layer separation from the glass cladding that supported the capillaries. If the computational boundary is set in this outer high-index cladding, then the ARF approximately appears as a $N_2$ Bragg fiber; but, if the boundary is placed in an air material outside of the cladding, then the ARF approximately imitates a $N_3$ Bragg fiber.
For the reference $N_2$ Bragg fiber, the thickness of the inner glass ring is chosen to be the same thickness as the capillary tube in the ARF. The width of the subsequent air layer is chosen based on the optimization process in [13], and corresponds approximately to the diameter of the capillary tubes. For the 6-tube ARF, the reference core radius is $15 \, \mu\textrm{m}$, the ring thickness is $0.42 \, \mu\textrm{m}$, and the air layer width is $12.48 \, \mu\textrm{m}$. For the 8-tube fiber, these measurements are $59.5 \, \mu\textrm{m}$, $6\, \mu\textrm{m},$ and $25.5 \, \mu\textrm{m}$ respectively. Their loss profiles are drawn as red curves in both plots of Fig. 5.
For the reference $N_3$ Bragg fiber that compares with the 8-tube ARF, we used the published scanning electron micrograph image [9, Fig. 3] to estimate the thickness of the outer glass cladding. Likewise, when comparing to the 6-tube ARF, we estimated this thickness from the images in [8, Fig. 3]. The spectral loss profiles for these reference $N_3$ Bragg fibers are drawn as blue curves in both plots of Fig. 5. These extra resonant peaks compared to the $N_2$ design CL profiles correspond to values predicted by (1) when $t$ is the thickness of the cladding.
Thus, we have seen the effects of changing the material configuration at the boundary of the computational domain, and verified that our numerical method can capture the results of these changes. In all cases, if the boundary is placed in what would otherwise be a complete anti-resonant layer, a set of loss spikes is eliminated from the mode spectral loss profile. In the next section, we will further explore the two aforementioned ARF designs.
3. ARF loss profiles in two outer configurations
In this section, we extend our previous conclusions about the sensitivity of the mode CL profiles on the material configuration at the computational boundary from Bragg fiber designs to real ARF designs found in the literature [8,9]. We have already seen in Section 2. that the number of resonant peaks in the loss profiles of the reference Bragg fibers depends on how the outer layer is modeled. Such differences will become more pronounced in a full-featured numerical model of the ARF.
Since the ARFs, or more generally microstructure fibers, do not perfectly correspond to Bragg fiber designs, we will not use the $N_{i}$ fiber design notation for these fibers. Instead, for notational convenience, we will set the high configuration ($C_{\text {high}}$) to be the case where the microstructure fiber is immersed in an infinite ambient high-index material (e.g., polymer coating or glass). For the model, this is equivalent to saying that the mode solver’s outer computational boundary is placed in a high-index material. Similarly, the low configuration ($C_{\text {low}}$) will refer to the case where the given microstructure fiber is immersed in an infinite ambient low index material (e.g., air). This notation holds even when the fiber may have more than one cladding layer such as having soft and hard polymer coatings around a glass cladding that supports the capillary tubes of the ARF designs. Note that when $i$ is an even number, the $N_{i}$ Bragg fibers can be considered to be of the $C_{\text {high}}$ configuration, and to be of the $C_{\text {low}}$ configuration when $i$ is odd.
Figure 6 shows the geometry and FM profile of the 8 capillary tube ARF presented in [9], which is in the $C_{\text {low}}$ configuration; i.e. immersed in air. In the fiber schematic (a), the blue regions are all comprised of air, whereas the red regions are comprised of glass. The black lines delineate regions of different materials, and/or are used for regulating the discretization of the domain (i.e. mesh size control), and/or distinguish the outer computational boundary layer (outer ring in blue encompassing the fiber and all other layers/materials). We implement the perfectly matched later (PML) to achieve the outer computational boundary that effectively extends that material to infinity (i.e. prevents/suppresses electromagnetic field reflections back into the fiber guiding region).
Fig. 6. Depictions of the fiber geometry (a), and the corresponding FM profile (magnitude of the transverse vector components) (b), of the 8 capillary tube anti-resonant hollow-core fiber presented in [9].
The geometry of the microstructured core region of the ARFs considered here can be determined by 5 parameters: (1) number of capillary tubes, $N_{\textrm{cap}}$, (2) core radius, $R_{\textrm{core}}$, (3) inner capillary radius, $R_{\textrm{cap}}$, (4) capillary thickness, $t_{\textrm{cap}}$, and (5) embedding distance, $e_{\textrm{cap}}$. The thicknesses of the remaining fiber material layers completes the model. For the 6- and 8-tube ARFs, we used the parameters delineated in Table 1.
Table 1. ARF parameter values used for the modeling experiments.
Figure 7 portrays the FM loss profiles of the 8-tube fiber. Clearly, the differences in the computed CL values for the $C_{\text {high}}$ and $C_{\text {low}}$ configurations are dramatic and disconcerting in regards to the objective of predicting the expected losses in real optical fiber waveguides. The $C_{\text {low}}$ loss profile has profound variation, on a scale difficult to resolve even with the 800 samples taken over the given wavelength interval. The total computed loss variation is nearly 6 orders of magnitude. In stark contrast, the $C_{\text {high}}$ loss profile is quite smooth on a granular level, and its total loss variation is a little less than 2 orders of magnitude.
Fig. 7. Spectral loss profile of the 8-tube ARF using the $C_{\text {high}}$ (red curve) and $C_{\text {low}}$ (blue curve) outer boundary configurations.
Figure 8 examines the FM losses over a smaller wavelength band for the same 8-tube ARF in the $C_{\text {low}}$ configuration as was seen in Fig. 7. Here the goal is to demonstrate that numerical convergence in the mode solver can be achieved even at specific wavelengths, apparently associated with resonant loss spikes, where the relative errors may remain higher than is desired as the degrees of freedom in the model are increased. The lower plot of Fig. 8 shows that under one further increase of the finite element polynomial order ($p$ described in §A.3), the relative error in the computed FM CL value is less than $1{\% }$ for almost all wavelengths at the chosen model resolution. Relative error is calculated as the absolute value of the difference between the CL values found for the subsequent polynomial order divided by the CL value found using the lower order. Locations in this wavelength interval, where error exceeded $1{\% }$, are associated with resonant loss spikes, but can be, and, in this case, have been, resolved by using a higher polynomial order (more degrees of freedom) in the mode solver.
Fig. 8. Wavelength interval of the 8-tube ARF’s FM spectral loss profile, indicating the numerical stability of the mode solver model. The green curve in the lower plot provides the relative percentage error in the CL value when the polynomial order is increased from 5 to 6. The red portions of the curve depict the lower relative error achieved by using even higher order elements.
This study is repeated for the 6-tube ARF, again exploring the $C_{\text {high}}$ and $C_{\text {low}}$ configurations; see Fig. 9. The message is the same: the FM CL profiles from the two computational boundary configurations are vastly different from one another. In comparison to the results for the 8-tube fiber, more structure appears in the $C_{\text {low}}$ spectral loss profile. As in that fiber, there is a great deal of variation, but here there are bandwidths of intense variation, interspersed with regions of relative smoothness. The lower portion of Fig. 9 is an inset focused on one of the high variance regions. Sampling more wavelengths for high resolution, the loss profile structure is well-resolved, depicting clear resonant peaks and anti-resonant troughs.
In order to ensure that the regions of high variability in the loss profiles within Fig. 9 are not just numerical noise, we assess a yet smaller bandwidth of wavelengths, as seen in plot (a) of Fig. 10, which illustrates that the FM CL profile is smooth over fine changes in wavelength. Furthermore, at a specific resonant loss peak (point A), and at a specific anti-resonant loss trough (point B), the FM’s Poynting vector longitudinal components are respectively depicted in plots (b) and (c), which have been normalized to have the same maximum value. The colorbar spans several orders of magnitude. As one might expect, the mode field profile at the loss peak (point A) displays a greater amount of its mode profile leaking into the glass cladding and beyond. Whereas the FM profile at the loss trough (point B) is much more contained within the core and capillary tube areas. The reader may wish to compare these FM Poynting vector longitudinal components to that which was seen in Poletti’s paper [8, Fig. 3], and further consider why the nested capillary design, also discussed in that article, better mitigates mode losses than this typical ARF design.
Fig. 10. Wavelength sub-interval (from Fig. 9) of the 6-tube ARF’s FM CL spectral profile (plot a), and the corresponding longitudinal components to the FM’s Poynting vector at a local loss peak (A; plot b) and at a local loss minimum (B; plot c). Both fields are scaled to same maximum value and plots have identical color profiles.
These studies illustrate that varying the outer material layer configuration has a similar but more pronounced effect on the ARFs’ loss profiles than for the Bragg fiber designs. The reader may find it interesting to compare the levels of CL variation that are computed over very fine changes in the wavelength in this effort to those observed experimentally: [8, Figs. 2 & 6], [14, Figs. 2, 3, & 5], [15, Fig. 4], likely measured over coarser wavelength changes. In the next section we analyze how the inclusion of a layer of lossy material affects these results.
4. Inclusion of lossy polymer
The vast difference between the FM loss profiles of the ARF $C_{\text {high}}$ and $C_{\text {low}}$ configurations naturally prompts this question: which constitutes the better modeling approach? The infinite high-index cladding of the $C_{\text {high}}$ configuration is obviously not realizable in practice, and their loss profiles are smoother than experimentally observed results (see [8,14,15]). Similarly, the $C_{\text {low}}$ configuration, immersed in infinite air, does not represent a fiber in any real life scenario. Since optical fibers usually have coatings that are typically many times thicker than their core diameters and these polymers are somewhat lossy, it seems prudent to understand how such a jacketing can affect the FM CL spectral profiles. Therefore, we will first examine the addition of a lossy material layer for the Bragg fiber design, where semi-analytic solutions are available, and then will investigate a similar lossy material coating for the ARF designs.
In [18, p. 63], it is noted that “in most practical cases cladding modes will be suppressed by a lossy coating on the outside of the fiber.” Such lossy coatings are present for most optical fiber applications to provide damage protection and stress relief. Even small distortions to an uncoated fiber can cause significant performance problems, hence a good deal of research has been done on optimizing optical jackets for the minimization of transmission loss [19–21].
Material losses can be accounted for with the addition of a complex component to the refractive index [22] called the extinction coefficient ($k$), where $k \equiv \text {imag}(n) = c \alpha _{\text {loss}} / (2 \omega )$ and $\alpha _{\text {loss}} \geq 0$ is the power loss value [1/m]. Measured loss values for polymer coatings can be difficult to find, sometimes because the information is kept proprietary. Figure 11 plots the dispersion of the extinction coefficient values over the near-to-short infrared regime for poly vinyl-chloride (PVC) and poly(methyl methacrylate) (PMMA), extracted from [16] and [17].
For the purpose of our current studies, we will use polymer thicknesses that are about the size of the core radius and extinction coefficients within $10^{-4} \leq k \leq 10^{-1}$. For reference, the corresponding power loss value at a wavelength of $\lambda = 3.25 \, \mu\textrm{m}$ is $\alpha _{\text {loss}} \approx 386$ 1/m or $1680$ dB/m for $k = 10^{-4}$, and grows by factors of $10$ as $k$ is increased to $10^{-3}$, $10^{-2}$, and $10^{-1}$. Our goal here is to draw attention to how the given fiber spectral loss profiles are controlled by this extinction coefficient ($k$), and could be used to match with experimental mode CL measurements if known/available. When using realistic polymer layer thicknesses in the loss analysis, it would be necessary to balance increased computational load with the level of desired detail sought in the spectral loss profile.
4.1 Bragg fiber with a polymer layer
We begin by extending the Bragg model to include a lossy polymer layer. The Bragg fiber still has semi-analytic solutions in this setting. This allows us to quickly explore the ramifications of including this layer in our model and to verify that our numerical method can still capture the loss values in this scenario.
First, consider the $N_1$ Bragg fiber design with an added lossy polymer layer directly outside the glass ring, yet before the low index ambient material. We set the thickness of the polymer equal to the core radius, while keeping all remaining dimensions as stated in Section 2. We fix the real part of the polymer refractive index to be slightly higher than that of glass, and we vary its imaginary component: $n_{\text {polymer}} = 1.5 + k \hat {\imath }$, where $\hat {\imath }$ is the imaginary unit.
Figure 12 depicts the differences between the FM profile longitudinal components for the polymer layer with (a) and without (b) material loss.
Fig. 12. Longitudinal component of the FM for a polymer coated $N_1$ Bragg fiber, with (b) and without (a) material loss.
Notice that plot (a) for the lossless polymer is very comparable to Fig. 2 plots (a)) and (c), which illustrated the FM longitudinal component for the glass/air-only $N_1$ Bragg fiber. Therefore, the presence of the lossless polymer layer is approximately the same as just extending the glass anti-resonance layer to a greater thickness, even though, in this example, the real-valued polymer refractive index is slightly higher than that of the glass. Alternatively, the mode profile in plot (b) is seen attenuating in the lossy polymer region, as one might expect to happen to an out-going wave passing through this region.
The FM CL profiles for the $N_1$ Bragg fiber jacketed with a polymer layer are portrayed in Fig. 13 for several extinction coefficient values. Since, when $k = 0$, this fiber is very similar to the glass/air-only $N_1$ Bragg fiber with a thicker glass anti-resonant layer, the loss profile exhibits resonant loss spikes, much like what was observed in Fig. 4. The computational results demonstrate that polymer losses ($k > 0$) dampen the magnitude of these resonant loss spikes. The other limiting case, as $k$ increases, seems to be that of the $N_0$ Bragg fiber of the same core radius. It appears that increasing material loss in the polymer makes both the glass annulus and its surrounding polymer layer progressively less resonant, until the loss profile becomes smooth. This is effectively as if the glass ring is given an infinite thickness; i.e. it converges to a $N_0$ Bragg fiber.
Fig. 13. Analytically determined spectral loss profiles (solid lines) for the $N_1$ Bragg fiber with a lossy polymer jacket. The green dots show the numerically computed losses for the $k = 0.002$ case.
The pattern from the $N_1$ case carries over to the $N_3$ configuration when lossy polymer is included there. In Fig. 14 we see the loss profile of the polymer coated $N_3$ Bragg fiber. As before, the measurements of the fiber are identical to those of section 2, but with an additional polymer layer added between the final layer of glass and the ambient low index material, with thickness equal to the core radius. Since this fiber contains resonant layers of different widths, it displays loss peaks on multiple scales. The wavelength region displayed is between the sixth and seventh loss spikes associated with resonances from the inner glass ring. The spikes visible in the figure are associated with resonances of the outer glass ring / polymer regions. As in the previous case, the limiting profile as $k$ increases is that of the Bragg fiber with one less complete resonant layer.
Fig. 14. Analytically determined spectral loss profiles for the polymer coated $N_3$ Bragg fiber. Green dots show the numerically computed losses for the $k = 0.0025$ case.
4.2 ARF with a polymer layer
With this intuition gathered from the Bragg fiber augmented with a lossy polymer coating, we repeat our previous FM CL computations for the 8-tube ARF, likewise encompassing the glass cladding with a lossy polymer material. As with the Bragg fiber of the previous subsection, the polymer layer is followed by a low index material; a $C_{\text {low}}$-type configuration. The polymer thickness is set to $15.3 \, \mu\textrm{m}$, as is the width of the buffer layer of air between the polymer and computational boundary layer. The PML layer is $22.95 \, \mu\textrm{m}$ thick (see Table 1).
The computed 8-tube ARF spectral FM CL profile results are plotted in Fig. 15 for polymer refractive indices of $1.5 + \hat {\imath } k$, where $k = 10^{-1}, 10^{-2}, 10^{-3}, 10^{-4}$. Also, the loss profile for the ARF in the $C_{\text {high}}$ configuration is displayed on the same plot since it seems to be the limiting case for high values of $k$. As $k$ decreases, the variation in the spectral profile increases, approaching the level of variation seen in the $C_{\text {low}}$ configuration results of Fig. 7.
Therefore, the inclusion of a lossy polymer layer smooths out the variations in the fiber’s spectral loss profile. These observations suggest that careful modeling of polymer layers will likely be necessary when validating numerically predicted mode loss profiles against measured values. Conversely, not precisely knowing the material properties of the polymer coatings may also limit experimental validation. Next, we further exemplify how profound of an effect that the model’s outer material configuration has on the computed mode CL when an interior fiber geometry is varied.
5. Sensitivity to geometric parameters in $C_{\text {high}}$ and $C_{\text {low}}$ configurations
It is natural to investigate how susceptible the computed losses of a microstructured fiber are to expected geometrical variations of the design. To this end, we study the effects of varying the radial embedding distance of capillary tubes into their supporting glass cladding as was illustrated in Fig. 1. In this study, we allow the tubes to embed some distance $e_{\textrm{cap}}$ into the cladding, varying the distance between 0 and the thickness of the capillary tube. We keep the core region and capillary tube radii constant, which implies that the radius to the start of the glass support cladding shifts to achieve the desired embedding.
The upper plot of Fig. 16 shows how the FM CL values ($\lambda = 1.8 \, \mu\textrm{m}$) change as the embedding fraction ($e_{\textrm{cap}} / t_{\textrm{cap}}$) is varied from 0 to 1 for the 6-tube ARF design. The results strongly depend on the choice of outer material layer configurations. In the $C_{\text {high}}$ configuration, where the glass cladding extends to infinity, the CL values exhibit essentially no sensitivity to changes in embedding fraction. Whereas for a $C_{\text {low}}$ configuration, where the ambient material is air, one observes that the CL values vary by about three orders of magnitude, which would have drastic impacts on optical field propagation in such a waveguide. The lower plot of Fig. 16 illustrates the numerical results for the case of a lossy polymer layer (for various extinction coefficients $k$) between the glass cladding and the infinite layer of air. Again, one observes loss profiles that act as a gradient between the two limiting cases as was demonstrated for the polymer coated Bragg fiber design in subsection 4.1.
Fig. 16. Embedding sensitivity of the FM CLs of the 6-tube ARF operating at a wavelength of $\lambda = 1.8 \, \mu\textrm{m}$.
These results suggest that if the fiber coating has a small enough extinction coefficient at the operating optical wavelength, then the embedding distance of the capillary tubes may play a large role in determining the mode loss. For higher polymer coating losses, perturbations in geometry will not affect the computed mode loss. However, if the loss in the polymer constitutes an absorption (rather than scattering), then one may find sufficient levels of thermal heating to induce polymer degradation/destruction [22]. Overall, this is a demonstration that emphasizes the fact that the outer material layer configurations must be treated correctly in the model to assess the robustness of a fiber design to the expected tolerances in geometry. Potentially, this outer layer treatment issue may explain previously unknown sources of loss and/or polymer heating in hollow-core/anti-resonant fiber designs.
6. Accelerating convergence with finer mesh in high-index regions
It has been reported [8,23] that numerical modeling of ARFs requires extremely fine mesh sizes for accurate loss predictions. In [8], the use of finite element sizes of, at most, $\lambda / 4$ in the air and $\lambda / 6$ in the glass regions were found to be necessary to achieve model convergence when the element order was quadratic ($p = 2$). Here we report that convergence of the FM CL values for these fibers may be accelerated by reducing the element size in the glass portions of the fiber alone. While there is broad agreement that mesh sizes must be small enough to capture the thin curvilinear capillary structures, what may be non-intuitive is that the homogeneous glass cladding region also needs a fine mesh, as we shall reason presently.
Careful examination of the computed mode profiles shows that in the cladding region, a fine radial wave is present, even at wavelengths where the fiber is in anti-resonance. This is in line with the high frequency mode profile oscillations seen in the high index regions of the Bragg fiber designs (see Figs. 2 and 12). Such features are also visible in [1, Fig. 4] and [24, Fig. 2], where their properties are directly linked to transmission minima in the fiber’s spectral loss profile. When using a coarse mesh in this region, it is not possible to capture this fine wave pattern accurately with low order elements. Decreasing the mesh size in this region allows these high-frequency oscillations to be better approximated with elements of the same order.
To clearly portray this phenomenon, consider the zoomed-in view of a mode profile in Fig. 17, near a capillary structure, computed on a series of four meshes from coarse (a) to fine (d), using the same finite element order. To obtain the finer meshes, we decreased the mesh size only in the glass cladding region, with the element size changing in low index regions only as necessary to maintain geometrical conformity of elements at the interfaces. One can see that the mode profile oscillations in the cladding are already well-approximated in the third mesh (c). The fourth mesh (d) confirms that these oscillations are not a numerical artifact and remain unchanged under further refinements. Next, we investigate the possibility that capturing this high-frequency mode profile structure in the cladding is critical for efficiently achieving numerical convergence of the modeled mode CL values.
Fig. 17. Longitudinal component of the FM for the $C_{\text {low}}$ 6-tube ARF in the cladding region, found with order 3 elements, on four distinct meshes for the glass cladding region. Note that $h$ represents the finite element mesh size, and $\lambda$ denotes the optical field wavelength.
To test this hypothesis, we computed the mode CL value starting from two meshes, a “naive” mesh with large elements in the cladding, and an “informed” mesh with fine elements in the cladding, both depicted in Fig. 18. The maximum diameter of elements in low index regions was kept near $\lambda / 0.4$. The FM CL values are computed for a sequence of increasing element orders, from linear ($p = 1$), to quadratic ($p = 2$), and higher. Specifically, the corresponding element order, as described in A.3, ranges between $p = 0$ and $p = 14$ on the mesh that is shown in plot (b), and between $p = 0$ and $p = 12$ on the mesh that is drawn in (c). The vertical dashed lines of plot (a) delineate the point where the model transitions from its preasymptotic regime to the regime of converged mode CL values.
Fig. 18. Comparison of numerical convergence for the naive and informed cladding mesh cases. Note that, in plot (a), the computational mode solver model is said to be in the preasymptotic regime to the left of the vertical dashed lines, and in the converged regime to the right of those vertical lines.
For both meshes, numerical convergence is achieved as the order of the elements used is increased. However, for the informed mesh, with a smaller mesh size in the high-index cladding, the CL convergence occurs at a lower number of degrees of freedom, such that there is little need to go beyond element order 6 or 7. For the naive mesh, variations in computed CL values persist even as one nears $10^{6}$ degrees of freedom. In other words, using a naive mesh, one faces a large preasymptotic regime, which must be crossed before trustworthy mode CL values are found. The presence of such a large preasymptotic regime for the ARF designs was first reported in [23]. In [23, Fig. 8], a plot similar to plot (a) of Fig. 18 was obtained using results from a scalar mode solver, while this effort leverages the hybrid vector mode solver described in Appendix A. Therefore, by using an informed mesh, one can obtain converged mode CL values with fewer degrees of freedom, accelerating the mode solving process.
7. Discussion
In this work it has been established that confinement loss values for leaky modes of anti-resonant fibers depend strongly on the choice of material configurations at and before the computational boundary. Allowing the cladding material to extend to the model boundary leads to smooth loss profiles, while terminating the cladding thickness and immersing the geometry in low index material leads to highly varying loss profiles. The inclusion of a lossy polymer layer led to profiles between these two extremes, with the variability of the profile decreasing as material loss is increased.
Because of this strong dependence, we recommend that numerical studies involving optical fiber modes provide details of any included outer material layers and of the implementation of the computational boundary, as it is now clear that this information is essential for reproducibility. This information further allows for others to make comparisons with different modeling choices and/or later experimental measurements.
The establishment of this dependence also raises the possibility that modelers may be able to more accurately predict experimental results by including complete representations of the coating materials into their optical fiber models. Of course, model benchmarking/validation would need to be accomplished while considering other known loss mechanisms such as surface scattering, macro- and micro-bending, geometry variations, and/or relevant material scatterings and absorptions [25]. It is our hope that this work inspires experiments for better, more complete, model validation. Such experiments would require detailed information on the wavelength dependent material losses of included coating materials. Ideally, polymer coating manufacturers would release material loss values over relevant optical wavelength ranges for their products. More public domain information on these materials could also lead to further studies on energy conversions (such as optical field absorption to thermal energy) within polymers and how that affects mode losses.
When initially designing an optical fiber, information on coatings may not be present. In this case, we recommend that modelers investigate losses using both the $C_{\text {high}}$ and $C_{\text {low}}$ configurations. The loss profile found using both configurations contain pertinent information about the fiber design. This is attested to by the fact that both are a distinct limiting case for the resulting spectral loss profile as the extinction coefficient was varied in the polymer layer studies of Section 4.
We conclude by offering a few general recommendations gathered from our numerical experience with modern anti-resonant fibers. First, one must recognize the importance of not presuming a continuous, nor differentiable, refractive index profile, particularly for microstructured fibers (which are not weakly guiding), having large index contrasts. The importance of choosing the correct mode solving approach is emphasized at the beginning of Appendix A. Furthermore, as it has been observed that the mode profiles of ARFs exhibit high-frequency oscillations in high-index cladding regions (anti-resonant layers), we encourage modelers to use meshes and discretization orders capable of resolving them accurately (such as the informed meshes we presented in Section 6). Finally, we strongly recommend that researchers ensure numerical convergence of computed mode losses, especially as this may require higher degrees of freedom than one might initially expect due to the above-mentioned high-frequency oscillations.
A. Appendix on the eigenproblem and the numerical method
Simulation techniques for hybrid vector modes of optical fibers are based on formulations that fall into two categories. The first uses derivatives of the refractive index in the equations: examples include [26–28]. The method we describe below falls into the second category of formulations without index derivatives [29–32]. Note that the first category is unsuitable for application to modern microstructured fibers with discontinuous jumps in refractive index: such refractive index functions are not differentiable (not even continuous).
A.1. Constrained eigenproblem for hybrid vector modes
Consider the time-harmonic Maxwell system for the electric field $\hat E$ and the magnetic field $\hat H$, assuming that time variations are of the form $e^{-\hat {\imath } \omega t}$, namely
We are interested in Maxwell solutions that propagate along the fiber’s longitudinal direction. Accordingly, we seek solutions of the form $\hat E(x, y, z) = E(x, y) e^{\hat {\imath } \beta z}$ and $\hat H(x, y, z) = H(x, y) e^{\hat {\imath } \beta z},$ where $\beta$ is the propagation constant. The fields $E, H$ are hybrid modes, i.e., vector fields with both transverse and longitudinal components, all of which depend only on the transverse coordinates $x, y$. Using the unit coordinate vectors, $e_x, e_y,$ and $e_z$, in the $x,y$ and $z$ directions, respectively, we may decompose a hybrid field $u$ into $u(x,y) = u_x(x,y) e_x + u_y(x,y) e_y + u_z(x,y) e_z.$ Denote its transverse projection by $u_\tau = u_x e_x + u_y e_y$. Let $J = \left [\begin {smallmatrix} \phantom {-}0 & 1 \\ -1 & 0 \end {smallmatrix}\right ]$ denote the operator that rotates vectors in the $xy$-plane clockwise by 90 degrees about the positive $z$-axis. For transverse fields, we define the 2D (scalar) curl and (vector) rot by $\operatorname {curl} u_\tau = \operatorname {div} (J u_\tau ) = \partial _x u_y - \partial _y u_x,$ and $\operatorname {rot} \phi = J (\operatorname {grad} \phi ) = (\partial _y\phi ) e_x - (\partial _x \phi ) e_y.$ Using the identity $\nabla \times \left (u\, e^{\hat {\imath }\beta z}\right ) = \left ( e^{\hat {\imath }\beta z} \operatorname {curl} u_\tau \right )e_z + e^{\hat {\imath }\beta z}\left ( \operatorname {rot} u_z - \hat {\imath }\beta J u_\tau \right ),$ the system (2) simplifies into the following coupled system for the two transverse vector fields $E_\tau, H_\tau,$ the two longitudinal components $E_z, H_z,$ and the associated propagation constant $\beta$.
Next, we eliminate the variables $H_z$ and $H_\tau$ from (4) and simplify to obtain the following system for $E_\tau$ and the scaled longitudinal component $\varphi = \hat {\imath } \beta E_z$,
A.2. Nondimensionalization
Considering an optical field in the infrared wavelength regime and practical fiber geometries, one finds propagation constants $\beta$ of the order of millions, with a length scale in micrometers. Hence nondimensionalization is necessary for preserving digits of accuracy in computations. Outside some finite radius in the transverse plane, the refractive index $n(x, y)$ is isotropic and equal to a constant $n_0$.
Define the index well $V(x, y) = k^2 (n_0^2 - n_\tau (x, y)^2).$ Next, using a length scale $L$ appropriate for the fiber cross section, we introduce the nondimensional transverse coordinates $\tilde x = x/L,$ and $\tilde y = y/L,$ and set
Then (6) implies the following nondimensional system
A.3. Discretization using Nédélec and Lagrange finite elements
For guided modes, we use the equations of (7) on a finite domain (a large enough circle) $\varOmega$ in the $xy$-plane, together with the boundary condition $E_\tau = \varphi = 0$ on $\partial \varOmega$. This condition is appropriate when $\partial \varOmega$ is far enough away from the compact support of $V$ for the guided mode components to have exponentially decayed to machine zero. Let $\varOmega$ be meshed by a geometrically conforming finite element mesh of triangles, denoted by $\varOmega _h$. On a triangle $K$, let $P_p(K)$ denote the space of polynomials in two variables of degree at most $p$. The degree $p$ and the maximal element diameter $h = \max _{K \in \varOmega _h} \mathrm {diam}(K)$ together can be used to compute the dimension of the discrete spaces. The Lagrange and the Nédélec spaces on $\varOmega _h$ are defined, respectively, for any $p\ge 0$ by $\mathcal {V}^{h} = \{ \psi : \psi |_K \in P_{p+1}(K) \text { for all } K \in \varOmega _h, \psi\,\textrm{is continuous} \},$ and $\mathcal {N}^{h} = \{ v \in \mathring {H}(\operatorname {curl},\varOmega ): v|_K \in P_p(K)^2 + \left [\begin {smallmatrix} \phantom {-}y \\ -x \end {smallmatrix}\right ] P_p(K) \text { for all } K \in \varOmega _h, v\,\textrm{is tangentially continuous}\}.$ Adjacent to curved material interfaces, we use curved triangles, in which case, the space within such an element is revised as usual to equal the pullback of the above-indicated polynomial spaces from a unit reference triangle.
Using both the boundary conditions $E_\tau =0$ and $\varphi = 0$ after multiplying (7) by test functions and integrating by parts, we obtain a “mixed finite element discretization” [33] that computes approximations of $E_\tau$ and $\varphi$, denoted by $E^{h}_\tau \in \mathcal {N}^{h}$ and $\varphi ^{h}\in \mathcal {V}^{h}$, respectively, by solving
The system (8) leads to a generalized matrix eigenproblem with $Z_h^2$ as the eigenvalue once a basis is used. To show it, let $\{\psi _i: i=1, \ldots, n_\mathcal {V}\}$ denote a basis for $\mathcal {V}^{h}$ and let $\{v_k: k = 1, \ldots, n_\mathcal {N}\}$ denote a basis for $\mathcal {N}^{h}$. Then using the matrices
For computation of leaky modes, the procedure is similar. The only modification needed is that standard artificial material coefficients dictated by a cylindrical PML [33] (set, e.g., in the outermost ring of plot (a) of Fig. 6), are incorporated in (7) (and consequently in (8) and (10)). For assembling these matrices, we use the open source finite element library NGSolve [34] which has facilities to robustly compute them for high values of $p$. Our implementation of the method is built as an add-on to this software.
A.4. Solving the eigenproblem via the FEAST algorithm
There are many algorithmic options for solving the constrained eigenproblem (10). In this subsection, we show how an eigensolver named “FEAST” [35–37] can be adapted for solving such eigenproblems. Letting
In our application, the right hand side matrix $\mathcal {B}$ is not invertible. The nullspace of $\mathcal {B}$ generates an eigenspace corresponding to an infinite eigenvalue. If the subspaces $Y_\ell$ and $\tilde {Y}_\ell$ produced by iteration (14) contain components in the nullspace of $\mathcal {B}$, computation of eigenvalue approximations from $Y_\ell, \tilde {Y}_\ell$ can numerically fail. This problem can be avoided in two ways. (i) The first is to periodically filter out the nullspace of $\mathcal {B}$ while conducting the algorithm, as described in [23], 16-22 of Algorithm 1]. (ii) A second alternative is to compute eigenvalue approximations from $Y_\ell, \tilde {Y}_\ell$ in a nonstandard way, but leveraging the constraint equation, as described next.
Writing the constrained eigenvalue Eq. (10) as a system of two equations,
we eliminate $b$. Since the second equation yields $b = -D^{-1} B a$, the first equation impliesSuppose that at the $\ell$th iteration in (14), we have a basis for the eigenspace approximations $Y_\ell$ and $\tilde {Y}_\ell$, listed as columns of $n \times m$ matrices $X$ and $\tilde {X}$, respectively. Block partitioning them as in (11), extract the top $n_\mathcal {N} \times m$ blocks and denote them by $x, \tilde {x} \in \mathbb {C}^{n_\mathcal {N} \times m}$. As the algorithm converges, each column of $x$ and $\tilde {x}$ approximates the $a$-block of a generalized left or right eigenvector in (15). It is possible to compute Ritz values using $x$ and $\tilde x$, ignoring the remaining blocks. Since the exact eigenvalue solves (15), an approximate Ritz value can be obtained using the small dense $m \times m$ matrices
By means of standard dense linear algebra routines, we then compute the Ritz values $\varLambda _\ell = \mathrm {diag}(\lambda _1, \ldots, \lambda _m)$, $w, \tilde w \in \mathbb {C}^{m \times m}$ satisfying $\tilde w^* T_x w = \varLambda _\ell$ and $\tilde w^* M_x w = I$. Unlike $\mathcal {B}$, the matrix $M$ is invertible. Hence whenever $x$ has linearly independent columns (readily ensured in practice), $M_x$ is invertible and the small dense eigenproblem is robustly solvable.
Funding
National Science Foundation (DMS-1912779, DMS-2136228); Air Force Office of Scientific Research (FA9550-23-1-0103).
Acknowledgments
P. Vandenberge gratefully acknowledges being selected for an NSF RTG Graduate Fellowship and for an internship in the AFRL Scholars Program, both of which facilitated this work.
Disclosures
The authors declare no conflicts of interest.
Disclaimers. This article is approved for public release; distribution unlimited. Public Affairs release approval $\# \text {AFRL-}2023\text {-}1168$. The views expressed are those of the authors and do not necessarily reflect the official policy or position of the Department of the Air Force, the Department of Defense, or the U.S. government. Portions of this work were published in Cornell University’s arXiv (10 March 2023; arXiv.2303.05713v1) and presented at the 2023 Laser Systems Modeling & Simulation Workshop in Albuquerque, NM USA (15 March 2023) by the same authors as this paper.
Data Availability
Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.
References
1. N. M. Litchinitser, A. K. Abeeluck, C. Headley, and B. J. Eggleton, “Antiresonant reflecting photonic crystal optical waveguides,” Opt. Lett. 27(18), 1592 (2002). [CrossRef]
2. C. Wei, R. Joseph Weiblen, C. R. Menyuk, and J. Hu, “Negative curvature fibers,” Adv. Opt. Photonics 9(3), 504 (2017). [CrossRef]
3. J. Hu and C. R. Menyuk, “Understanding leaky modes: Slab waveguide revisited,” Adv. Opt. Photonics 1(1), 58 (2009). [CrossRef]
4. D. Marcuse, Theory of Dielectric Optical Waveguides (Elsevier, 2013).
5. W. Belardi and J. C. Knight, “Hollow antiresonant fibers with reduced attenuation,” Opt. Lett. 39(7), 1853–1856 (2014). [CrossRef]
6. D. Bird, “Attenuation of model hollow-core, anti-resonant fibres,” Opt. Express 25(19), 23215–23237 (2017). [CrossRef]
7. M. A. Duguay, Y. Kokubun, T. L. Koch, and L. Pfeiffer, “Antiresonant reflecting optical waveguides in SiO2-Si multilayer structures,” Appl. Phys. Lett. 49(1), 13–15 (1986). [CrossRef]
8. F. Poletti, “Nested antiresonant nodeless hollow core fiber,” Opt. Express 22(20), 23807–23828 (2014). [CrossRef]
9. A. N. Kolyadin, A. F. Kosolapov, A. D. Pryamikov, A. S. Biriukov, V. G. Plotnichenko, and E. M. Dianov, “Light transmission in negative curvature hollow core fiber in extremely high material loss region,” Opt. Express 21(8), 9514–9519 (2013). [CrossRef]
10. P. Yeh and A. Yariv, “Bragg reflection waveguides,” Opt. Commun. 19(3), 427–430 (1976). [CrossRef]
11. P. Yeh, A. Yariv, and E. Marom, “Theory of Bragg Fiber,” J. Opt. Soc. Am. 68(9), 1196–1201 (1978). [CrossRef]
12. J.-L. Archambault, R. J. Black, S. Lacroix, and J. Bures, “Loss calculations for antiresonant waveguides,” J. Lightwave Technol. 11(3), 416–423 (1993). [CrossRef]
13. F. Poletti, J. R. Hayes, and D. J. Richardson, “Optimising the performances of hollow antiresonant fibres,” in 2011 37th European Conference and Exhibition on Optical Communication, (2011), pp. 1–3.
14. J. R. Hayes, S. R. Sandoghchi, T. D. Bradley, Z. Liu, R. Slavík, M. A. Gouveia, N. V. Wheeler, G. Jasion, Y. Chen, E. N. Fokoua, M. N. Petrovich, D. J. Richardson, and F. Poletti, “Antiresonant Hollow Core Fiber With an Octave Spanning Bandwidth for Short Haul Data Communications,” J. Lightwave Technol. 35(3), 437–442 (2017). [CrossRef]
15. K. S. R. Shaha, A. Khaleque, and M. S. Hosen, “Wideband Low Loss Hollow Core Fiber With Nested Hybrid Cladding Elements,” J. Lightwave Technol. 39(20), 6585–6591 (2021). [CrossRef]
16. X. Zhang, J. Qiu, X. Li, J. Zhao, and L. Liu, “Complex refractive indices measurements of polymers in visible and near-infrared bands,” Appl. Opt. 59(8), 2337–2344 (2020). [CrossRef]
17. M. N. Polyanskiy, “RefractiveIndex.Info,” https://refractiveindex.info.
18. D. Marcuse, Light Transmission Optics2nd Edition (Van Nostrand Reinhold Company Inc., 1982).
19. D. Gloge, “Optical-fiber packaging and its influence on fiber straightness and loss,” The Bell Syst. Tech. J. 54(2), 245–262 (1975). [CrossRef]
20. S.-T. Shiue, “Design of double-coated optical fibers to minimize hydrostatic pressure induced microbending losses,” IEEE Photonics Technol. Lett. 4(7), 746–748 (1992). [CrossRef]
21. F. Cocchini, “The lateral rigidity of double-coated optical fibers,” J. Lightwave Technol. 13(8), 1706–1710 (1995). [CrossRef]
22. A. W. Snyder and J. Love, Optical Waveguide Theory (Springer US, 2024).
23. J. Gopalakrishnan, B. Q. Parker, and P. VandenBerge, “Computing leaky modes of optical fibers using a FEAST algorithm for polynomial eigenproblems,” Wave Motion 108, 102826 (2022). [CrossRef]
24. N. M. Litchinitser, S. C. Dunn, B. Usner, B. J. Eggleton, T. P. White, R. C. McPhedran, and C. M. de Sterke, “Resonances in microstructured optical waveguides,” Opt. Express 11(10), 1243–1251 (2003). [CrossRef]
25. E. N. Fokoua, S. A. Mousavi, G. T. Jasion, D. J. Richardson, and F. Poletti, “Loss in hollow-core optical fibers: mechanisms, scaling rules, and limits,” Adv. Opt. Photonics 15(1), 1–85 (2023). [CrossRef]
26. A. W. Snyder and W. R. Young, “Modes of optical waveguides,” J. Opt. Soc. Am. 68(3), 297–309 (1978). [CrossRef]
27. T. Monro, D. Richardson, N. Broderick, and P. Bennett, “Modeling large air fraction holey optical fibers,” J. Lightwave Technol. 18(1), 50–56 (2000). [CrossRef]
28. H. P. Uranus and H. J. W. M. Hoekstra, “Modelling of microstructured waveguides using a finite-element-based vectorial mode solver with transparent boundary conditions,” Opt. Express 12(12), 2795 (2004). [CrossRef]
29. M. Koshiba and K. Inoue, “Simple and efficient finite-element analysis of microwave and optical waveguides,” IEEE Trans. Microwave Theory Tech. 40(2), 371–377 (1992). [CrossRef]
30. M. Koshiba and Y. Tsuji, “Curvilinear hybrid edge/nodal elements with triangular shape for guided-wave problems,” J. Lightwave Technol. 18(5), 737–743 (2000). [CrossRef]
31. J.-F. Lee, D.-K. Sun, and Z. Cendes, “Full-wave analysis of dielectric waveguides using tangential vector finite elements,” IEEE Trans. Microwave Theory Tech. 39(8), 1262–1271 (1991). [CrossRef]
32. L. Vardapetyan and L. Demkowicz, “Full-Wave Analysis of Dielectric Waveguides at a Given Frequency,” Math. Comp. 72(241), 105–130 (2003). [CrossRef]
33. P. Monk, Finite Element Methods for Maxwell’s Equations (Clarendon Press, 2003).
34. J. Schoberl, “NGSolve,” ngsolve.org.
35. J. Gopalakrishnan, L. Grubišić, and J. Ovall, “Spectral discretization errors in filtered subspace iteration,” Math. Comp. 89(321), 203–228 (2020). [CrossRef]
36. S. Güttel, E. Polizzi, P. T. P. Tang, and G. Viaud, “Zolotarev quadrature rules and load balancing for the FEAST eigensolver,” SIAM J. Sci. Comput. 37(4), A2100–A2122 (2015). [CrossRef]
37. J. Kestyn, E. Polizzi, and P. T. Peter Tang, “Feast Eigensolver for Non-Hermitian Problems,” SIAM J. Sci. Comput. 38(5), S772–S799 (2016). [CrossRef]
