Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Efficient Bloch mode calculation of periodic systems with arbitrary geometry and open boundary conditions in the complex wavevector domain

Open Access Open Access

Abstract

We show how existing iterative methods can be used to efficiently and accurately calculate Bloch periodic solutions of Maxwell’s equations in arbitrary geometries. This is carried out in the complex-wavevector domain using a commercial frequency-domain finite-element solver that is available to the general user. The method is capable of dealing with leaky Bloch mode solutions, and is extremely efficient even for 3D geometries with non-trivial material distributions. We perform independent finite-difference time-domain simulations of Maxwell’s equations to confirm our results. This comparison demonstrates that the iterative mode finder is more accurate, since it provides the true solutions in the complex-wavevector domain and removes the need for additional signal processing and fitting. Due to its efficiency, generality and reliability, this technique is well suited for complex and novel design tasks in integrated photonics, and also for a wider range of photonics problems.

© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

1. Introduction

Design of photonic components often involves using geometries with periodic modulation of the material index [1]. Well-known examples are grating couplers, photonic crystals, metamaterials, and metasurfaces [26]. The study of such periodic geometries can significantly benefit from Bloch mode analysis, in order to improve numerical efficiency and provide a better physical insight into device behavior. This is due to the fact that Bloch modes are only calculated using a single unit-cell of the periodic device and they form the modal basis for light propagation through the media. For example, the dispersion diagram of photonic crystals, effective propagation constant of subwavelength waveguides, and diffraction angle of surface grating couplers, can all be extracted from the system Bloch modes. Since such devices (and all planar devices in the context of integrated photonics) are only periodic in one or two dimensions, their Bloch modes are naturally subject to open boundary conditions [7] in other directions. Such boundary conditions are commonly implemented in simulations via a perfectly matched layer (PML) [8], in order to accurately handle the energy flow leaving the computation box. The Bloch modes in such situations can be intrinsically leaky. Mathematically, leaky Bloch modes can be described using complex-wavevectors for which the imaginary part describes the losses through light scattering (note that material losses, can also be included in such a picture).

There are different approaches for calculating Bloch modes. The exact method is to solve an eigenvalue problem in the complex-wavevector domain [9,10]. While such an approach can provide the full dispersion diagram, it can be resource intensive, particularly for the purpose of optimization where many simulation runs are required, and is not easily available to all researchers/designers. There are also techniques such as guided mode expansion [11] and Fourier mode expansion [12,13] methods that can be very efficient in simulating Maxwell’s equations. However, among other technical limitations, these methods require considerable expertise to assure their reliable implementation [14]. Moreover, dealing with arbitrary shaped geometries and 3D configurations can become expensive to the point that they cannot be justified against standard Maxwell solvers such as finite-difference time-domain (FDTD). Indeed, commercially available FDTD solvers are very popular and can be used to work around solving eigenvalue problems in the complex-wavevector domain. Such an approach relies on identifying system resonances through a time-domain simulation on a periodic unit-cell where Bloch boundary condition with real-valued wavevector is enforced. This method commonly involves using signal processing techniques such as harmonic inversion [15] to estimate the degree of the losses associated with the identified resonances, from which one can infer the complex-wavevector, $k_{\textrm {Bloch}}=k_{\textrm {Bloch}}^{\textrm {R}}+i\,k_{\textrm {Bloch}}^{\textrm {I}}$ associated with the Bloch mode. However, such FDTD implementations are also resource intensive, particularly for 3D simulations (a common concern when using FDTD). Moreover, multiple simulations over different wavevector values are needed to interpolate for the target $k_{\textrm {Bloch}}$ which is not known in advance. There are also many situations (such as the example case studies presented here) where the designer only needs to know one particular Bloch mode at a given frequency and not the entire band diagram of the system. Therefore, a rigorous and efficient method of directly calculating the true solutions of the Bloch eigenvalue problem subjected to open boundary conditions would be very useful in many optical design and simulation applications.

In this work, we present an iterative technique that runs frequency-domain finite-element method (FEM) calculations at a given frequency to find the target complex $k_{\textrm {Bloch}}$ of the system subject to both Bloch-periodic and open boundary conditions. This is inspired by the success of the mode analysis techniques previously developed for efficient simulation of resonator systems [16,17] in the complex frequency domain, even for complex 3D geometries [18]. Once the complex $k_{\textrm {Bloch}}$ is available, we show that for various devices of practical interest in integrated photonics, only simple analytical equations are needed to arrive at important physical quantities. These include the effective Bloch mode index, $n_{\textrm {Bloch}}$, describing light propagation through periodic devices such as (surface/edge) grating couplers, as well as the scattering strength $\gamma$ and the diffracting angle $\theta$ of the surface grating couplers. The approach is very efficient and usually converges in less than 10 iterations (including 3 iterations for initialization) in all of our case studies. We have tested our technique on grating couplers (using 2D simulations) where the entire process took less than 30 seconds, as well as subwavelength waveguides (using 3D simulations) for less than 3 minutes on a desktop computer. We additionally tested our approach on a more complex 3D geometry that was introduced in Ref. [19] for the purpose of on-chip shaping of optical beams and confirmed convergence in only 7 minutes, despite the increased demand on meshing. To evaluate the accuracy, reliability and efficiency of our method, we perform completely independent FDTD calculations to compare against the iterative method results. We implemented our method using the commercial FEM solver from COMSOL [20] and therefore it can be easily adopted by others in the community.

2. Method

Bloch modes provide an elegant way of solving Maxwell’s differential equations in the presence of a periodic modulation of the material properties. Mathematically, the Bloch modes $\mathbf {E}_{\textrm {Bloch}}\left (\mathbf {r}\right )$ and $\mathbf {H}_{\textrm {Bloch}}\left (\mathbf {r}\right )$ are solutions to source-free Maxwell equations in periodic systems,

$$\nabla\times\mathbf{E}_{\textrm{Bloch}}\left(\mathbf{r}\right) = i\omega_{k_{\textrm{Bloch}}}\mu(\omega_{k_{\textrm{Bloch}}})\,\mathbf{H}_{\textrm{Bloch}}\left(\mathbf{r}\right)$$
$$\nabla\times\mathbf{H}_{\textrm{Bloch}}\left(\mathbf{r}\right) ={-}i\omega_{k_{\textrm{Bloch}}}\varepsilon(\omega_{k_{\textrm{Bloch}}})\,\mathbf{E}_{\textrm{Bloch}}\left(\mathbf{r}\right),$$
having the property
$$\mathbf{E}_{\textrm{Bloch}}\left(\mathbf{r}\right) = e^{ik_{\textrm{Bloch}}x}\,\mathbf{e}_{\textrm{Bloch}}\left(\mathbf{r}\right)$$
$$\mathbf{H}_{\textrm{Bloch}}\left(\mathbf{r}\right) = e^{ik_{\textrm{Bloch}}x}\,\mathbf{h}_{\textrm{Bloch}}\left(\mathbf{r}\right).$$

Here, $\mathbf {E}_{\textrm {Bloch}}\left (\mathbf {r}\right )$ and $\mathbf {H}_{\textrm {Bloch}}\left (\mathbf {r}\right )$ are the electric and magnetic field distributions of the Bloch mode over space, $\mathbf {r}$, with the periodic components $\mathbf {e}_{\textrm {Bloch}}\left (\mathbf {r}\right )$ and $\mathbf {h}_{\textrm {Bloch}}\left (\mathbf {r}\right )$ along the direction $x$, $k_{\textrm {Bloch}}$ is the Bloch wavevector along the periodic direction, and $\omega _k$ is the wavevector-dependent angular frequency. The material properties are described using frequency dependent $\varepsilon (\omega _k)$ and $\mu (\omega _k)$ functions.

As discussed in the introduction, when solving for the Bloch modes of a given system, the correct exercise is to impose an open boundary condition along dimensions that are not periodic. Therefore, in general, Bloch modes can be solutions to a non-Hermitian Maxwell problem with complex eigenvalues, where the corresponding eigensolutions are the so-called "leaky" modes of the system. In the past, iterative techniques have been successful and efficient in calculating leaky modes of various resonator systems in the complex frequency domain [16,18,21]. The iteration procedure starts with an initial guess for the complex eigenvalue to eventually find its actual value. Inspired by the method presented in Ref. [16] for the analysis of resonant systems in the complex frequency domain, to guide an iterative predictor in the complex-wavevector domain, we use the equation,

$$E = \dfrac{c_0 + c_1\,k}{k-k_{\textrm{Bloch}}},$$
with two unknown constants $c_0$ and $c_1$, as well as the unknown target $k_{\textrm {Bloch}}$ of the Bloch mode. Note that this equation serves as an approximate representation of the electric field $E$ (a certain component in practice) at a given location, that is calculated from solving Maxwell’s equations subjected to a dipole excitation [16]. Although in general the choice of the excitation and observation locations is not fundamentally important (this has been numerically tested by choosing multiple observation points in each case study), regions with expected high field values are favorable. Since there are three unknowns, we initially solve Maxwell’s equations for three times at $k$, $k+\delta$ and $k-\delta$ where $k$ is our initial guess for the complex wavevector of the Bloch mode and $\delta$ is an arbitrary shift parameter. Each solution generates a corresponding $E$ value to provide an instance of Eq. 5. Now, we solve the set of three equations to get our first estimates on $c_0,\,c_1$ and $k_{\textrm {Bloch}}$. Next, we run a new simulation at the estimated $k=k_{\textrm {Bloch}}$ that was just obtained. From this new set of four equations, only the three with the strongest field responses ($E$) are kept for the next iteration. This procedure is repeated until the predicted values converge and we arrive at the target $k_{\textrm {Bloch}}$.

It is important to note that, as discussed in Ref. [16], in order to obtain the system modes correctly, one needs to solve the Maxwell equation twice. Once in the presence of the full geometry where the light scattering behavior is captured; and once again in the presence of a homogeneous background where the source contributions are captured to be then subtracted. However, this is practically needed to obtain the modal field distribution, while the iterative procedure itself can be successfully guided only using the simulations in the presence of the full geometry. This makes sense as the target solutions are the poles of the Green function for the entire problem (see Ref. [16] and Ref. [17] for more details). Therefore one can be more efficient (almost by a factor of 2) by running the additional background simulations needed for calculating the proper system modes only once when the iteration process has converged.

Once the $k_{\textrm {Bloch}}$ value that satisfies Eqs. 14 is known, various quantities relevant to the design procedure can be analytically obtained. This is an instantaneous post-processing stage employing simple algebra, as we will see below. For example, the effective index associated with the Bloch mode, $n_{\textrm {Bloch}}$, is readily available using:

$$n_{\textrm{Bloch}} ={\pm}\,k_{\textrm{Bloch}}^{\textrm{R}}/k_0,$$
where $k_0=2\pi /\lambda$ is the propagation wavevector in vacuum at operation wavelength of $\lambda$. Note that $\pm$ correspond to when Bloch wavevector is in the same/opposite direction of the group velocity. Additionally, the imaginary part of athe Bloch wavevector dictates the exponential amplitude drop for the Bloch mode. Therefore, the field intensity decay (scattering strength) through one unit-cell of propagation is quantified via
$$\gamma = \textrm{exp}(2k_{\textrm{Bloch}}^{\textrm{I}}/\Lambda),$$
with a negative imaginary part for the Bloch wavevector. Moreover, in the event of diffracting behavior such as in surface grating couplers, momentum conservation arguments can be used to obtain the diffracting angle for the device. Following the convention shown in Fig. 1(a) and using geometrical identities, one can write (for the zeroth order diffraction):
$$\theta = \sin^{{-}1}(n_{\textrm{Bloch}}/n_{\textrm{cl}}),$$
where $n_{\textrm {cl}}$ is the refractive index of the cladding that in general can be different from that of air. Note that the wavevector associated with the free propagation into the cladding medium is given by $k_{\textrm {cl}}=n_{\textrm {cl}}k_0$.

 figure: Fig. 1.

Fig. 1. (a) Side view of an example surface grating coupler with periodicity $\Lambda$ that can redirect light from an input waveguide on the integrated chip to a fiber located above. (b) Side view of a subwavelength waveguide that offers flexible effective index engineering. (c) Top view of an integrated Bragg deflector coupler for on-chip optical beam shaping. We show a periodic example of the device studied in Ref. [19].

Download Full Size | PDF

3. Case studies

For our case studies, we first consider two important classes of devices commonly used in integrated photonics. Schematically shown in Fig. 1(a) is a surface grating coupler that is used to interface optical fibers and integrated photonic chips. Also in Fig. 1(b) is an optical waveguide with subwavelength patterning to engineer its effective index. These devices are typically formed by periodic structuring of the waveguide core medium and therefore can be efficiently studied and modeled using Bloch mode theory. To demonstrate the power of our technique in dealing with arbitrary-shaped 3D geometries, we additionally consider a periodic version of the distributed Bragg deflector coupler studied in Ref. [19] as shown in Fig. 1(c). The proposed unit-cell in 3D is non-trivial and redirects the light from an input waveguide by 90 degree toward an on-chip output waveguide with significantly different mode size.

From the mode theory perspective, the subwavelength waveguides might be considered (theoretically) lossless devices where in general their Bloch modes do not couple to the states above the light-line. On the other hand, surface grating couplers that are designed to redirect light from the in-plane integrated chip to the above environment and then (sometimes) a fiber, are necessarily leaky devices. This is also the case for the deflector coupler where the light is redirected on the same chip rather than to free space. The iterative complex-wavevector approach presented here can be reliably applied to both of these scenarios and extract useful parameters, as we described before. Indeed, as we will see, the challenging scenario of the perfectly vertical diffracting behavior of the surface grating couplers (near zero wavevector) can be also accurately studied by the iterative leaky Bloch mode calculator.

3.1 Surface grating coupler, an open system example

We begin with the surface grating coupler example. First, we present the convergence and performance of the iterative method in finding the complex-wavevector Bloch modes. Then we study three example grating couplers with different diffraction angles including the particular case of perfectly vertical diffraction at $k_{\textrm {Bloch}}\approx 0$ to emphasize the reliability of the technique.

For the grating coupler we choose a geometry that is non-trivial and has shown success in achieving perfectly vertical operation with high fiber-coupling efficiency and low back-reflection [22,23]. As shown in Fig. 1(a), this geometry uses a pillar section (for controlling back-reflection) and a L-shaped step (to redirect the light upward more effectively). This can be done using a two-step etch process, for which we consider a shallow etch of 50%, half-way through the waveguide height. We consider a standard 220 nm silicon-on-insulator (SOI) platform with a $3\,\mu \textrm{m}$ buried oxide layer and a semi-infinite oxide cladding on top. The refractive index of silicon and oxide used in these simulations are respectively $n_{\textrm{Si}}=3.479$ and $n_{{\textrm{SiO}}_2}=1.448$. Our first case study has segment parameters of $L_i=[78,\,84,\,117,\,249,\,166]$ nm and performs perfectly vertical with a near-zero $k_{\textrm{Bloch}}$. In this example we run 2D simulations of Maxwell’s equation as we have verified that they agree very well with the 3D simulations.

In Fig. 2, the general iterative workflow for this example grating coupler is presented. As can be seen, we start with an initial guess of $k_{\textrm{Bloch}}\Lambda /\pi =0.2-0.2\,i$ and request a relative tolerance of $10^{-4}$ on the absolute value of the wavevector as the stop criterion. Convergence is achieved within only six iterations after the initial three steps that determine the starting estimates on $c_0,\,c_1$ and $k_{\textrm{Bloch}}$. This is remarkably efficient as each iteration takes approximately 3 seconds on a desktop computer (a complete study using this technique took approximately 30 seconds in comparison to 5 minutes for the corresponding FDTD implementation that we will discuss later). The final estimated value of the complex-wavevector in this case is $k_{\textrm{Bloch}}\Lambda /\pi =-0.002-0.032\,i$. Through using Eq. (5), this leads to an estimated value of $\theta =-0.08^{\circ }$ for the grating diffraction angle which is considered perfectly vertical.

 figure: Fig. 2.

Fig. 2. Iterative calculation of the complex-wavevector Bloch modes for the surface grating coupler shown in Fig. 1(a). The bottom plot show the real and imaginary parts of the estimated $k_{\textrm{Bloch}}$ for each iteration. Note that the first three steps shown here correspond to three closely spaced solutions of Maxwell’s equations that are used to initialize the iteration procedure, as described in Section 2. After only the six subsequent iterations following initialization, a relative convergence of $10^{-4}$ is obtained. The colormaps on top show the field profile associated with the corresponding Bloch mode estimation at each iteration. In this case, we place the dipole excitation at (81, 0) nm and the observation point at (320, 0) nm, when the origin is placed in the middle of the $L_1$ segment and half way through the waveguide height.

Download Full Size | PDF

Next, we demonstrate the technique reliability and ease of use for devices with different operational characteristics. To this end, we started with the previously discussed grating coupler and made segment adjustments such that the grating diffraction angle changes. This can be done by adjusting the aspect ratio of the $L_4/L_5$ as the L-shaped component largely determines the upward diffraction of light from the grating. We simply make $L_4$ smaller/larger to have a set of three devices, for which the results of the analysis are summarized in Table 1. For all of these calculations the same level of tolerance as before, i.e. $10^{-4}$, was demanded, and convergence was achieved after less than nine iterations (including the three initialization steps). Note that the imaginary part of the Bloch wavevector $k_{\textrm{Bloch}}^{\textrm {I}}$ in all three cases is negative, as it must be to generate a decaying field profile along the direction of propagation. The real part $k_{\textrm {Bloch}}^{\textrm {R}}$ on the other hand, can be either positive or negative. These represent different band structures for different devices, that here translate into different diffracting angles (see Table 1).

Tables Icon

Table 1. Comparison of the Bloch mode analysis results for three surface grating couplers with different diffracting angles $\theta$, using the iterative complex-wavevector method and real-wavevector FDTD method. The diffracted beam angles $\theta$ obtained by independent fiber coupling simulations in FDTD are also given in the final row.

To evaluate the accuracy of the various quantities obtained with the iterative calculations, we carry out a completely independent set of analyses in the time domain using FDTD solver from Ansys/Lumerical [24]. This includes a real-wavevector band dispersion analysis on a single unit-cell, as well as a scattering experiment where the light injected via a slab waveguide propagates through the full grating coupler with 25 periods and scatters up toward the fiber. Included in Table 1 is a comparison of the real part of the wavevector $k_{\textrm {Bloch}}^{\textrm {R}}$ obtained from the FDTD band calculation that shows a very good agreement up to the second decimal place with the iterative complex-wavevector. The relative discrepancy is more noticeable for $k_{\textrm {Bloch}}\approx 0$. It is important to note that the iterative analysis results are expected to be closer to the true values, not only because it works in the correct complex-wavevector domain but also because it does not rely on any major post-processing steps. This is confirmed by the independent fiber coupling simulations also carried out using FDTD. In these simulations, the fundamental TE mode of an input waveguide is injected into the grating, and the resulting upward diffracted field above the structure is used to perform a mode overlap analysis with a $10.4\,\mu \textrm {m}$ MFD pick-up fiber. We vary the fiber angle to find the angle at which the coupling efficiency is maximized. This $\theta$ is also included in Table 1 along with the estimated diffraction angle obtained from both the iterative complex-wavevector and real-wavevector FDTD studies. As seen, the complex-wavevector analysis results consistently have the best agreement with the FDTD fiber coupling simulations. Additionally, in Figs. 3(a)–3(c), we plot the field profile of the input waveguide mode propagating through the grating coupler calculated with FDTD and compare it with the calculated Bloch mode profile using iterative complex-wavevector method in Fig. 3(d)–3(f). The grating geometry is also overlayed to show the excellent agreement between the Bloch mode study and the field propagation study. Note that the Bloch mode calculation provides the actual solution to the source-free Maxwell’s equations in wavevector domain, compared to the FDTD where propagation effects subjected to the choice of excitation can distort the resulted field profile.

 figure: Fig. 3.

Fig. 3. (a)-(c) Scattered electric field profile for grating couplers discussed in Table 1 calculated with FDTD for the entire structure subjected to an injected TE mode. The geometry is overlayed for each device to aid visual comparison with Bloch mode results. (d)-(f) Bloch modes calculated using iterative complex-wavevector analysis on a unit-cell, corresponding to the extended light propagation shown in (a)-(c).

Download Full Size | PDF

3.2 Periodically patterned waveguides

Waveguides with subwavelength patterning can be of great interest in the field of integrated photonics for a variety of applications. In these structures diffraction is suppressed because of very short periodicity, and the patterning is used to manipulate the local effective index and dispersion of the device. For example, extremely efficient waveguide edge couplers have been demonstrated in standard SOI platforms [25]. Non-trivial geometries have also been used to engineering band dispersion [26]. Such devices can be easily treated using the proposed iterative Bloch mode approach in a matter of minutes, even using fully 3D computations.

Our example subwavelength waveguide geometry is schematically shown in Fig. 1(b) and is studied using full 3D simulations. As in the case of the grating couplers, we consider standard SOI platform where the 220 nm thick waveguide core with refractive index of $n_{\textrm{Si}}=3.479$ is embedded in a glass background with $n_{{\textrm{SiO}}_2}=1.448$. The waveguide width is $w=450\,\textrm{nm}$, the pitch of the subwavelength patterning is $\Lambda =400\,\textrm{nm}$ and the duty cycle (DC), that is the ratio of the non-etched silicon over the pitch, is varied in our case study.

We study three different values of the duty cycle: 0.7, 0.5 and 0.3. These are respectively labeled as A, B and C and the results are summarized in Table 2. The imaginary parts of the wavevectors are considered to be zero if they fall below the estimated numerical accuracy of the simulations used here that is $10^{-4}$. If necessary the accuracy can be improved by decreasing the convergence tolerance and allocating more computational resources. In all cases the solver finds $k^{\textrm{I}}_{\textrm {Bloch}}\approx 0$ consistent with the expectation that the corresponding Bloch modes are theoretically lossless in the subwavelength regime. Table 2 also includes the results of an independent 3D-FDTD band dispersion analysis for comparison. For the lossless Bloch modes, the real-wavevector FDTD method has a better agreement with the iterative complex-wavevector method than seen for the previous grating coupler examples. In the same table, we also provide the effective Bloch mode index, $n_{\textrm {Bloch}}$, using both the iterative and the FDTD method. As expected, for smaller DC values, the effective mode index decreases and becomes closer to the refractive index of the background medium.

Tables Icon

Table 2. Comparison of the Bloch mode analysis for three subwavelength grating waveguides using the iterative complex-wavevector method and real-wavevector FDTD analysis. The SWG pitch of $\Lambda =400\,\textrm {nm}$ and duty cycle values of $\textrm {DC}=0.7,\,0.5,\,0.3$ are studied.

In terms of time efficiency, in all of these three cases, each iteration took approximately 20 seconds and a total number of eight iterations was sufficient to satisfy the $10^{-4}$ convergence tolerance. Therefore the total study time per device is no more than 3 minutes, in comparison to 1 hour/device for the FDTD analysis. Figure 4 depicts the electric field profiles of the Bloch modes obtained using the iterative method for three different DC values as labeled. Note that, as the DC decreases, the field profile of the Bloch mode is less confined within the high index material, as expected.

 figure: Fig. 4.

Fig. 4. Electric field distribution ($E_z$) of the Bloch mode for three different waveguide with subwavelength patterning calculated using the iterative complex-wavevector method. The DC values are 0.7, 0.5 and 0.3 from left to right, in decreasing order. The corresponding $k_{\textrm {Bloch}}$ and $n_{\textrm {Bloch}}$ values are listed in Table 2. For these calculations, we place th dipole at (0, 0, 50) nm and the observation point at (150, 110, 0) nm, when the origin is in the middle of the SWG material.

Download Full Size | PDF

3.3 Bragg deflector coupler

To complete our demonstration of the efficiency and generality of the iterative Bloch mode calculations, we study a periodic version of the distributed Bragg deflector coupler schematically shown in Fig. 1(c). This device consists of a single mode waveguide connected to a Bragg deflector that diffracts light at right angles to the incident propagation direction in order to form a shaped slab waveguide output beam. The tapered sections provide adiabatic index matching to the slab waveguide. This device has been proposed for on-chip optical beam shaping [19]. This device uses a non-trivial geometry and accurate modelling requires a complex mesh structure and the use of 3D simulations.

The material characteristics of this silicon photonic device are the same as of the previous examples. One unit-cell of the device is characterized by several parameters as shown in Fig. 5(a), showing the actual values that we use here. With so many parameters to vary and the need to perform 3D simulations, optimization of such a device can significantly benefit from efficient and accurate complex-wavevector Bloch mode calculations. In our FEM implementation of this device, we use PML termination to properly capture the outgoing power flow at the output port of the device.

 figure: Fig. 5.

Fig. 5. (a) Schematic top-view of the unit-cell geometry for the Bragg deflector coupler studied in Ref. [19] showing the structural dimensions that we used to carry out our iterative Bloch mode calculations. (b) Real part of the electric field profile ($E_z$) of the calculated Bloch mode within the unit-cell of the device. We also highlight the direction of the input and the output for this coupler. In this case, we place the dipole excitation at (0, 0, 50) nm and the observation point at (150, 110, 0) nm, when the origin is in between the two large etch areas and 300 nm away from the sharp ends along the y axis.

Download Full Size | PDF

For this Bragg deflector coupler, using our iterative process, it took 9 iterations to arrive at the converged wavevector values of $k_{\textrm {Bloch}}\Lambda /\pi =0.779 - 0.025\,i$. Note that, similar to the case of surface grating coupler, here, the imaginary part of the wavevector can be used to accurately estimate the scattering strength of the unit-cell under study that is an essential parameter for the design process [19]. In Fig. 5(b), we also plot the electric field profile of the calculated Bloch mode where the transition from the input waveguide mode to the output waveguide mode is clear. Given the involved geometry and the fully 3D nature of the simulations, this is a highly challenging task for existing numerical approaches but it took only approximately 7 minutes for the total of 9 iterations with the proposed method.

4. Conclusion

We have presented an efficient, easy to implement and reliable method for calculating Bloch modes of periodic systems, for both leaky and lossless cases. The method uses a commercially available frequency-domain Maxwell solver to iteratively look for the target Bloch mode at a given frequency in the complex-wavevector domain. Three different types of devices have been show-cased: (a) surface grating couplers, (b) subwavelength optical waveguides, and (c) Bragg deflector on-chip couplers. The method provides reliable numbers for various quantities that are of major importance in the design of integrated photonics components, such as the effective modal index, diffraction angle, and scattering strength. It also out-performs the time-domain methods not only in terms of simulation efficiency, but also in accuracy. In contrast to the standard time-domain methods, our technique avoids any numerical processing as well as fitting exercises that can only provide estimates of the real and imaginary parts of the eigenvalues for the Bloch-periodic Maxwell problem. We have presented case studies both in 2D and 3D scenarios and have imposed no constraints on the geometry of the objects. Therefore, no matter the complexity of the device under study, we foresee great benefit in using this method by designers and researchers in the photonics community.

Disclosures

The authors declare no conflicts of interest.

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. Y. D’Mello, O. Reshef, S. Bernal, E. El-fiky, Y. Wang, M. Jacques, and D. V. Plant, “Integration of periodic, sub-wavelength structures in silicon-on-insulator photonic device design,” IET Optoelectron. 14(3), 125–135 (2020). [CrossRef]  

2. L. Cheng, S. Mao, Z. Li, Y. Han, and H. Y. Fu, “Grating couplers on silicon photonics: Design principles, emerging trends and practical issues,” Micromachines 11(7), 666 (2020). [CrossRef]  

3. J. D. Joannopoulos, S. G. Johnson, J. N. Winn, and R. D. Meade, Photonic Crystals: Molding the Flow of Light (Princeton University, 2008), 2nd ed.

4. S. Jahani and Z. Jacob, “All-dielectric metamaterials,” Nat. Nanotechnol. 11(1), 23–36 (2016). [CrossRef]  

5. M. Kadic, G. W. Milton, M. van Hecke, and M. Wegener, “3d metamaterials,” Nat. Rev. Phys. 1(3), 198–210 (2019). [CrossRef]  

6. H.-T. Chen, A. J. Taylor, and N. Yu, “A review of metasurfaces: physics and applications,” Rep. Prog. Phys. 79(7), 076401 (2016). [CrossRef]  

7. P. T. Kristensen, R.-C. Ge, and S. Hughes, “Normalization of quasinormal modes in leaky optical cavities and plasmonic resonators,” Phys. Rev. A 92(5), 053810 (2015). [CrossRef]  

8. J.-P. Berenger, “A perfectly matched layer for the absorption of electromagnetic waves,” J. Comput. Phys. 114(2), 185–200 (1994). [CrossRef]  

9. V. Laude, Y. Achaoui, S. Benchabane, and A. Khelif, “Evanescent Bloch waves and the complex band structure of phononic crystals,” Phys. Rev. B 80(9), 092301 (2009). [CrossRef]  

10. J. Notaros and M. A. Popovic, “Finite-difference complex-wavevector band structure solver for analysis and design of periodic radiative microphotonic structures,” Opt. Lett. 40(6), 1053–1056 (2015). [CrossRef]  

11. L. C. Andreani and D. Gerace, “Photonic-crystal slabs with a triangular lattice of triangular holes investigated using a guided-mode expansion method,” Phys. Rev. B 73(23), 235114 (2006). [CrossRef]  

12. Q. Cao, P. Lalanne, and J.-P. Hugonin, “Stable and efficient Bloch-mode computational method for one-dimensional grating waveguides,” J. Opt. Soc. Am. A 19(2), 335–338 (2002). [CrossRef]  

13. A. Ortega-Moñux, I. Molina-Fernández, and J. G. Wangüemert Pérez, “3D-Scalar Fourier Eigenvector Expansion Method (Fourier-EEM) for analyzing optical waveguide discontinuities,” Opt. Quantum Electron. 37(1-3), 213–228 (2005). [CrossRef]  

14. D. F. G. Gallagher and T. P. Felici, “Eigenmode expansion methods for simulation of optical propagation in photonics: pros and cons,” in Integrated Optics: Devices, Materials, and Technologies VII, vol. 4987Y. S. Sidorin and A. Tervonen, eds., International Society for Optics and Photonics (SPIE, 2003), pp. 69–82.

15. Harminv: http://ab-initio.mit.edu/wiki/index.php?title=Harminv.

16. Q. Bai, M. Perrin, C. Sauvan, J.-P. Hugonin, and P. Lalanne, “Efficient and intuitive method for the analysis of light scattering by a resonant nanostructure,” Opt. Express 21(22), 27371–27382 (2013). [CrossRef]  

17. C. Sauvan, J. P. Hugonin, I. S. Maksymov, and P. Lalanne, “Theory of the Spontaneous Optical Emission of Nanosize Photonic and Plasmon Resonators,” Phys. Rev. Lett. 110(23), 237401 (2013). [CrossRef]  

18. M. Kamandar Dezfouli, R. Gordon, and S. Hughes, “Modal theory of modified spontaneous emission of a quantum emitter in a hybrid plasmonic photonic-crystal cavity system,” Phys. Rev. A 95(1), 013846 (2017). [CrossRef]  

19. A. Hadij-ElHouati, P. Cheben, A. Ortega-Moñux, J. G. Wangüemert-Pérez, R. Halir, J. H. Schmid, and I. Molina-Fernández, “Distributed bragg deflector coupler for on-chip shaping of optical beams,” Opt. Express 27(23), 33180–33193 (2019). [CrossRef]  

20. COMSOL Multiphysics: https://www.comsol.com.

21. M. Kamandar Dezfouli, C. Tserkezis, N. A. Mortensen, and S. Hughes, “Nonlocal quasinormal modes for arbitrarily shaped three-dimensional plasmonic resonators,” Optica 4(12), 1503–1509 (2017). [CrossRef]  

22. T. Watanabe, M. Ayata, U. Koch, Y. Fedoryshyn, and J. Leuthold, “Perpendicular Grating Coupler Based on a Blazed Antiback-Reflection Structure,” J. Lightwave Technol. 35(21), 4663–4669 (2017). [CrossRef]  

23. D. Melati, Y. Grinberg, M. Kamandar Dezfouli, S. Janz, P. Cheben, J. H. Schmid, A. Sánchez-Postigo, and D.-X. Xu, “Mapping the global design space of nanophotonic components using machine learning pattern recognition,” Nat. Commun. 10(1), 4775 (2019). [CrossRef]  

24. Lumerical: https://www.lumerical.com.

25. P. Cheben, D.-X. Xu, S. Janz, and A. Densmore, “Subwavelength waveguide grating for mode conversion and light coupling in integrated optics,” Opt. Express 14(11), 4695–4702 (2006). [CrossRef]  

26. P. Cheben, J. Čtyroký, J. H. Schmid, S. Wang, J. Lapointe, J. J. G. Wangüemert-Pérez, I. Molina-Fernández, A. Ortega-Moñux, R. Halir, D. Melati, D. Xu, S. Janz, and M. Dado, “Bragg filter bandwidth engineering in subwavelength grating metamaterial waveguides,” Opt. Lett. 44(4), 1043–1046 (2019). [CrossRef]  

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.

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (5)

Fig. 1.
Fig. 1. (a) Side view of an example surface grating coupler with periodicity $\Lambda$ that can redirect light from an input waveguide on the integrated chip to a fiber located above. (b) Side view of a subwavelength waveguide that offers flexible effective index engineering. (c) Top view of an integrated Bragg deflector coupler for on-chip optical beam shaping. We show a periodic example of the device studied in Ref. [19].
Fig. 2.
Fig. 2. Iterative calculation of the complex-wavevector Bloch modes for the surface grating coupler shown in Fig. 1(a). The bottom plot show the real and imaginary parts of the estimated $k_{\textrm{Bloch}}$ for each iteration. Note that the first three steps shown here correspond to three closely spaced solutions of Maxwell’s equations that are used to initialize the iteration procedure, as described in Section 2. After only the six subsequent iterations following initialization, a relative convergence of $10^{-4}$ is obtained. The colormaps on top show the field profile associated with the corresponding Bloch mode estimation at each iteration. In this case, we place the dipole excitation at (81, 0) nm and the observation point at (320, 0) nm, when the origin is placed in the middle of the $L_1$ segment and half way through the waveguide height.
Fig. 3.
Fig. 3. (a)-(c) Scattered electric field profile for grating couplers discussed in Table 1 calculated with FDTD for the entire structure subjected to an injected TE mode. The geometry is overlayed for each device to aid visual comparison with Bloch mode results. (d)-(f) Bloch modes calculated using iterative complex-wavevector analysis on a unit-cell, corresponding to the extended light propagation shown in (a)-(c).
Fig. 4.
Fig. 4. Electric field distribution ($E_z$) of the Bloch mode for three different waveguide with subwavelength patterning calculated using the iterative complex-wavevector method. The DC values are 0.7, 0.5 and 0.3 from left to right, in decreasing order. The corresponding $k_{\textrm {Bloch}}$ and $n_{\textrm {Bloch}}$ values are listed in Table 2. For these calculations, we place th dipole at (0, 0, 50) nm and the observation point at (150, 110, 0) nm, when the origin is in the middle of the SWG material.
Fig. 5.
Fig. 5. (a) Schematic top-view of the unit-cell geometry for the Bragg deflector coupler studied in Ref. [19] showing the structural dimensions that we used to carry out our iterative Bloch mode calculations. (b) Real part of the electric field profile ($E_z$) of the calculated Bloch mode within the unit-cell of the device. We also highlight the direction of the input and the output for this coupler. In this case, we place the dipole excitation at (0, 0, 50) nm and the observation point at (150, 110, 0) nm, when the origin is in between the two large etch areas and 300 nm away from the sharp ends along the y axis.

Tables (2)

Tables Icon

Table 1. Comparison of the Bloch mode analysis results for three surface grating couplers with different diffracting angles θ , using the iterative complex-wavevector method and real-wavevector FDTD method. The diffracted beam angles θ obtained by independent fiber coupling simulations in FDTD are also given in the final row.

Tables Icon

Table 2. Comparison of the Bloch mode analysis for three subwavelength grating waveguides using the iterative complex-wavevector method and real-wavevector FDTD analysis. The SWG pitch of Λ = 400 nm and duty cycle values of DC = 0.7 , 0.5 , 0.3 are studied.

Equations (8)

Equations on this page are rendered with MathJax. Learn more.

× E Bloch ( r ) = i ω k Bloch μ ( ω k Bloch ) H Bloch ( r )
× H Bloch ( r ) = i ω k Bloch ε ( ω k Bloch ) E Bloch ( r ) ,
E Bloch ( r ) = e i k Bloch x e Bloch ( r )
H Bloch ( r ) = e i k Bloch x h Bloch ( r ) .
E = c 0 + c 1 k k k Bloch ,
n Bloch = ± k Bloch R / k 0 ,
γ = exp ( 2 k Bloch I / Λ ) ,
θ = sin 1 ( n Bloch / n cl ) ,
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.