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

Bound states in the continuum in waveguide arrays within a symmetry classification scheme

Open Access Open Access

Abstract

We study a photonic implementation of a modified Fano-Anderson model – a waveguide array with two additional waveguides and by using the coupled mode theory we calculate its spectral and scattering properties. We classify eigenmodes according to vertical symmetry of the structure given by self-coupling coefficients of the additional waveguides and establish the conditions for bound states in the continuum (BIC) existence. The main predictions drawn from the theoretical model are verified by rigorous full-wave simulations of realistic structures. We use the Weierstrass factorization theorem and interpret the scattering spectra of the systems with broken symmetry in terms of the eigenmodes. The Fano resonance related with excitation of quasi-BIC is explained as arising from the interference between this mode and another leaky mode.

© 2022 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement

1. Introduction

A counterintuitive concept of bound states in the continuum (BICs) was initially proposed by von Neumann and Wigner [1] as a mathematical construct represented by special states which lie in the continuum but remain localized without any radiation as a consequence of confining potential. By inspecting various scenarios of light scattering and their relation to location of the poles and zeros of the scattering matrix, BICs were identified as the limiting case of general Hermitian scatterer, when the pole and zero coalesce on the real axis with mutual coherent destruction [2]. BICs have been found in a wide range of different material platforms and have been experimentally observed in electromagnetic, acoustic and water waves [35].

In this paper we investigate a discrete photonic system supporting a symmetry-protected BIC. The structure, which can be regarded as generalization of a famous finite system in which an optical BIC was confirmed experimentally for the first time [6], consists of an infinite array of single-mode optical waveguides with two additional waveguides above and below the array (Fig. 1). Our configuration can be viewed as an extension of a photonic implementation of the Fano-Anderson model [7] with a single additional waveguide [8]. We note that a similar structure consisting of a discrete linear chain, which in addition to Fano defects contains one $\delta$-like defect which gives rise to the Fano-Feshbach resonance, was considered in [9].

 figure: Fig. 1.

Fig. 1. A coupled system consisting of a periodic array with two additional waveguides above and below the array.

Download Full Size | PDF

To explore physics and to provide a theoretical analysis of the formation of the symmetry-protected BIC, we determine spectrum of eigenmodes as well as scattering properties of the waveguide array and present the theoretical framework describing the symmetry classification scheme. Namely, we establish the conditions, given by relations among the self-coupling coefficients associated with side-coupled waveguides $\delta _{\pm 1}$ and the change of the propagation constant $\varepsilon$ that define regimes in which the governing coupled dynamic equations support the different eigenmode spectra. By solving the equations corresponding with regimes with preserved and broken-vertical symmetry we unveiled (i) yet undiscovered spectral properties of the non-Hermitian system, which reveal then existence of the EPs fundamentally different from those occurring in the $\cal {PT}$-symmetric systems, and (ii) the role of the leaky and bound states in interpretation of the scattering properties that has been elucidated on the basis of the Weierstrass factorization technique. To investigate the limits of the model we carried out rigorous simulations of realistic structures and verified the most important predictions arising from the classification scheme. In addition, we examined applicability of the model in respect to nonlocal coupling by evaluating the structures with appropriately selected parameters.

The eigenvalues $\varepsilon$ observed in the model reveal new and interesting features, such as the level repulsion and exceptional points (EPs), that have been found by tracing their position within the complex plane with increasing difference between $\delta _1$ and $\delta _{-1}$. The knowledge of a set of the discrete singular points associated with the scattering matrix, provides a widely applicable framework to interpret and to engineer the resonant properties of the nanostructures [10,11]. All the modes supported by the structure manifest themselves in the scattering spectra as the resonances that can be identified by means of the generalized Weierstrass factorization approach [10,11]. We found that the resonance associated with quasi-BIC state possesses a Fano line-shape arising from the interference between quasi-BIC and another leaky mode.

In practice, BIC manifests itself as quasi-BIC with finite $Q$-factor limited by material absorption, finiteness of structure, geometrical disorder, and fabrication imperfection [5]. In our model we assume only radiative losses, which can be compensated when the structure supports a symmetry protected BIC, which is completely decoupled from the radiating waves. The losses increase with the increasing asymmetry strength $\delta$, while the contribution associated with the inherent non-radiative losses is neglected. In practice (full-wave simulation or experiment), this assumption has to be verified and contributions of both radiative and non-radiative losses require a rigorous assessment.

For the sake of compactness of the present paper we refer reader interested in more subtle details of the complex band structure and the scattering properties to our companion paper [12], which contains additional results which illustrate the nature of the regimes inspected and interpret the mechanisms behind the salient features found in the scattering spectra. We also provide in Ref. [12] a more detailed description of the Weierstrass factorization technique and prove that the transmittance near the resonance can be rewritten into the form of the Fano formula, where the shape parameter $f$ can be expressed in terms of the poles associated with the interacting modes.

2. Model

To describe the coupled structure (Fig. 1) we use the standard approach based on the coupled mode theory (CMT) [13,14]. The total field is represented with slowly varying modal amplitudes $\psi _m(z)$ ($m\in Z$) and $\varphi _{\mp 1}(z)$ on the $m$th waveguide in the array and on the upper and lower additional waveguide, respectively. The evolution of the amplitudes is described by a set of coupled equations

$$ i\psi_m' = C\left(\psi_{m-1} + \psi_{m+1}\right),\; m\neq~0 $$
$$i\psi_0' = \delta_{0} \psi_{0}+C\left(\psi_{{-}1} + \psi_{1}\right) + C_1\varphi_{{-}1}+C_2\varphi_{1}$$
$$i\varphi_{{-}1}' = \delta_{{-}1} \varphi_{{-}1} + C_1\psi_{0} $$
$$i\varphi_{1}' = \delta_{1} \varphi_{1} + C_2\psi_{0},$$
where the prime stands for the derivative according to the $z$, $C$ is the coupling coefficient in the array, $C_{1,2}$ are the coupling coefficients between the lower and upper additional waveguide and the $0$th waveguide in the array, and $\delta _{0}$ and $\delta _{\mp 1}$ are self-coupling coefficients in the $0$th waveguide in the array and in the upper and lower additional waveguide, respectively. The self-coupling coefficients may include the effect of nearest-neighbors due to the geometry as well as the possible external perturbations of $0$-th waveguide and the additional waveguides. We assume that all the parameters which enter Eqs. (1)–(4) are real.

The unperturbed array supports propagation of Bloch waves $\psi _m=A\exp \left (-i k_xam-i\varepsilon z\right )$, where $k_x$ is the Bloch wavenumber, $a$ is the period of the lattice and the “energy” (change of the propagation constant $\beta$) $\varepsilon$ fulfils the well-known dispersion relation

$$\varepsilon = 2C\cos\left( k_x a\right).$$

3. Scattering problem

The additional waveguides cause scattering of the Bloch waves. To evaluate the effect, we follow the standard approach, see, e.g., [8,15], and represent the field in the left part of the periodic array in terms of incident and reflected waves

$$\psi_m=e^{{-}i k_xam-i\varepsilon z} + r e^{i k_xam-i\varepsilon z}, \; m < 0$$
in the right part as the transmitted wave
$$\psi_m=t e^{{-}i k_xam-i\varepsilon z}, \; m > 0$$
and the fields in the additional waveguides as
$$\varphi_{{\pm} 1}=B_{{\pm} 1} e^{{-}i\varepsilon z}.$$
After substituting Eqs. (6)–(8) into Eqs. (1)–(4) and few straightforward steps one finds that Eqs. (6), (7) are valid also for $m=0$ and obtains the analytical expressions for the amplitude transmittance $t$ and reflectance $r$ as
$$t = r+1= \frac{2C\sin(k_xa)}{2C\sin(k_xa)+i \mu}$$
where $\mu$ reads
$$\mu = \delta_0 + \frac{C_1^{2}}{\varepsilon-\delta_{{-}1}} + \frac{C_2^{2}}{\varepsilon-\delta_{1}}.$$

Figure 2 demonstrates typical behavior of the intensity transmittance $T = |t|^{2}$ evaluated for various structural parameters. Figures 2(a) and (b) show results for symmetric structures where the self-coupling coefficients in both additional waveguides are the same $\delta _1 = \delta _{-1}$. The most dominant feature is vanishing $T$ at single point $\varepsilon = \delta _{\pm 1}$ (for fixed $k_x$, $0<k_xa<\pi$) due to diverging parameter $|\mu |\rightarrow \infty$ in Eq. (9), see also Eq. (10). This resonance corresponds to existence of the symmetry protected BIC and, for varying $k_x$, is indeed located at positions given by the dispersion relation (5) [see dashed lines in Figs. 2(a) and (b)]. Note also, that the zero transmittance at $\varepsilon = \delta _{\pm 1}$ is double degenerate and identical with energies of localized states of the additional waveguides, i.e., when $\varepsilon$ is calculated from Eqs. (3) and (4) for $C_{1,2}=0$.

 figure: Fig. 2.

Fig. 2. The transmittance $T$ of the waveguide array. (a), (b) Systems with the symmetry $\delta _{-1}= \delta _1$; (c), (d) asymmetric systems defined as $\delta =-\delta _{-1}=\delta _1$; $k_x$ is the Bloch wavenumber, $C_{1,2}$ are shown in the boxes, $\delta _0=0$. The dashed line in (a) and (b) displays the dispersion relation (5) with $\delta _{\pm 1} = \varepsilon$.

Download Full Size | PDF

By breaking the vertical symmetry, the degeneracy is lifted as is illustrated in Figs. 2(c), (d) and Fig. 3 for asymmetric structures defined as $\delta =-\delta _{-1}=\delta _1$, where the parameter $\delta$ describes an asymmetry strength. We observe two different zeros with energies equal to those of states of isolated additional waveguides, $\varepsilon = \pm \delta$. It will be shown further, that the resonances correspond to the excitation of two leaky modes of the composed structure, one of them being quasi-BIC, and the related spectral features may be explained by interference of the two leaky modes that leads to the Fano profile.

 figure: Fig. 3.

Fig. 3. The transmittance $T$ vs. the Bloch wavenumber $k_x$ for systems with preserved $\delta =0$ (dashed line) and broken $\delta /C=0.3$ (solid line) symmetry for the coupling coefficients indicated and $\delta _0=0$.

Download Full Size | PDF

Importantly, the zero positions observed in Figs. 2 and 3 do not depend on the coupling coefficients $C_{1,2}$. On the other hand, $C_{1,2}$ have profound influence on the shape of the resonances: decreasing of $C_1/C$ and $C_2/C$ (weaker coupling) leads to narrower dip in the transmittance $T$, while asymmetry in $C_{1,2}$ strongly affects asymmetric behavior of the two dips [compare Figs. 3(a) and (b)].

4. Modes

The spectral features described above can be explained in terms of the modes of the composed structure. The modal field in the linear chain has the form

$$\psi_m=A e^{- i k_xa|m|-i\varepsilon z}, $$
while that in the additional waveguides is given by Eq. (8). By substituting Eqs. (11) and (8) into Eqs. (2)–(4) we obtain the following nonlinear eigenvalue problem
$$ \left[\delta_0-\varepsilon+2C\exp\left(- ik_x a\right)\right]A+C_1B_{{-}1}+C_2B_1=0 $$
$$C_1A+ \left(\delta_{{-}1}-\varepsilon\right)B_{{-}1}=0 $$
$$C_2A+ \left(\delta_{1}-\varepsilon\right)B_{1}=0.$$
The relation between the self-coupling coefficients $\delta _{\pm 1}$ and change in the propagation constant $\varepsilon$ plays the key role in a symmetry classification scheme, which defines the regimes with different families of the eigenmodes.

A. We first consider the case when $\varepsilon =\delta _1$ or $\varepsilon =\delta _{-1}$. According to Eqs. (12)–(14) a nontrivvial solution is not supported unless the structure possesses the vertical symmetry $\delta _{-1}=\delta _1$ when one obtains a doubly-degenerate bound state $\varepsilon = \delta _{\pm 1}$, which yields a vertically antisymmetric eigenfunction

$$C_1B_{{-}1}={-}{C_2}B_1, \,\, A = 0$$
and represents the symmetry protected BIC provided $-2C<\delta _{\pm 1}<2C$.

B. Now we examine the cases $\varepsilon \neq \delta _{1}$ and $\varepsilon \neq \delta _{-1}$. Equations (13)–(14) yield nontrivial solutions for the amplitudes associated with the attached waveguides in the form

$$B_{{\mp} 1}=\frac{C_{1,2}}{\varepsilon - \delta_{{\mp} 1}}A .$$
By substituting the both expressions into Eq. (12) one obtains the eigenvalue equation
$$2C\sin(k_xa)+ i\mu = 0.$$
It is important to keep in mind, that none of the solutions of Eq. (17) is a true BIC as we have $A\neq 0$, however, for the vertically asymmetric structures ($\delta _{-1}\neq \delta _1$), one of them may represent a quasi-BIC as will be shown below.

B.1 Symmetric states. The family of the symmetric states can be divided to several groups according to the vertical symmetry of the structure.

B.1a We first consider the structure with the trivial vertical symmetry, $\delta _{-1}= \delta _1 = \delta _0=0$, for which Eq. (17) provides 4 eigenvalues

$$\varepsilon ={\pm} \sqrt{2C^{2} \pm 2\sqrt{C^{4}+C_{\rm av}^{4}}},$$
where all combinations of the signs are allowed and $C_{\rm av}^{2} = (C_1^{2}+C_2^{2})/2$. The modes (in contrast to antisymmetric modes supported by the same vertically symmetric structures) exhibit symmetry
$$C_2B_{{-}1}={C_1}B_1.$$

The full spectrum of the trivially symmetric structures, shown in Figs. 4(a) and (b), consists of symmetry protected doubly-degenerate BIC within the continuum, two bound modes (B1, B2) with real energy, which bifurcate from the lower and upper edges of the continuum when the parameter $C_{\rm av}$ is increased, and two modes (L, L’) with pure imaginary energy. We note that Eq. (18) illustrates the general feature of Eq. (17): the energy eigenvalues are either real or occur in complex conjugate pairs. The states with $\mathrm{Im}(\varepsilon )<0$ (solid lines in Figs. 4 and 5) fulfil the conventional definition of leaky modes and satisfy outgoing boundary conditions $0<\mathrm{Re}(k_xa)\le \pi$, the states with $\mathrm{Im}(\varepsilon )>0$ (dashed lines in Figs. 4 and 5) correspond to ingoing boundary conditions $-\pi <\mathrm{Re}(k_xa)<0$ [12].

 figure: Fig. 4.

Fig. 4. Spectrum of structures with (a), (b) trivial, $\delta _{\pm 1}=\delta _0=0$, and (c), (d) nontrivial, $\delta _{-1}=\delta _{1}$, vertical symmetry. (a) The real and (b) imaginary part of the energy $\varepsilon$ as a function of $C_{\rm av}/C$. (c) The real and (d) imaginary part of $\varepsilon$ as a function of $\delta _{\pm 1}/C$ for $\delta _0=0$ and $C_1=C_2=C$. The modes in (c) and (d) are labeled according to their asymptotics when $\delta _{\pm 1} \rightarrow 0$. Dashed lines indicate the states with $\mathrm{Im}(\varepsilon )>0$. The gray area in (a) and (c) indicates the continuum.

Download Full Size | PDF

 figure: Fig. 5.

Fig. 5. Spectrum of structures with broken vertical symmetry $\delta =-\delta _{-1}=\delta _1$: the energy $\varepsilon$ vs. the asymmetry strength $\delta$. (a) and (b) $C_1=C_2=C$; (c) and (d) $C_1/C=0.1$, $C_2/C=0.4$; $\delta _0=0$. The gray area in (a) and (c) indicates the continuum.

Download Full Size | PDF

B.1b For the structure belonging to the class characterized by the general vertical symmetry $\delta _{-1}= \delta _1 \neq 0$, Eq. (17) yields four solutions, the modes are symmetric according to Eq. (19). The spectrum shown in Figs. 4(c) and (d) also consists of the bound modes (B1, B2, BIC); the latter one enters the gap at $\delta _{\pm 1}/C =2$ when $\delta _{\pm 1}/C$ is increased. The L, L’ modes initially exhibit complex conjugate energies $\varepsilon$, the real component of their energy $\mathrm{Re}(\varepsilon )$ also enters the gap at $\delta _{\pm 1}/C = 2.5$. The imaginary component of their energy $\mathrm{Im}(\varepsilon )$ remains nonzero in the range $2.5 < \delta /C < 3.1$ until the (L, L’) doublet reaches the EP beyond which the level repulsion effect occurs. Here the level repulsion arises from the coupling between two degenerate modes with positive energy [16]. To this end it is important to keep in mind, that this EP which occurs at $\delta /C = 3.8$ differs from those known in PT-symmetric structures and indicates a more general non-Hermitian - Hermitian transition [17]. We will show that breaking the vertical symmetry leads to an intriguing behavior of the complex band structure – see Fig. 5 – which reveals attraction and repulsion of the levels in certain ranges of the asymmetry strength $\delta /C$.

B.2 Structures with broken vertical symmetry (asymmetric states). Breaking the vertical symmetry, $\delta _{-1}\neq \delta _1$, yields coupling of the true BIC antisymmetric state to the continuum and the doubly degenerate mode is split into quasi-BIC (qBIC) and qBIC’ mode as shown in Fig. 5, which demonstrates a rich spectral behavior of the structure. Equation (17) yields six solutions, the amplitudes $B_{\pm 1}$ instead of Eq. (19) satisfy the condition

$$C_2B_{{-}1}(\delta_{{-}1}-\varepsilon)= C_1B_1(\delta_1-\varepsilon).$$

The spectra for the structures with $C_1=C_2$, see Figs. 5(a) and (b), reveal in comparison with those shown in Fig. 4 the EP which occurs when the asymmetry strength $\delta /C = 0.5$. Likewise in the case of system with nontrivial vertical symmetry, the EP indicates a more general non-Hermitian – Hermitian transition [17]. In contrast to the results shown in Figs. 4(c) and (d) the complex band structure shown in Fig. 5 exhibits rather peculiar behavior. Specifically, the L, L’ and qBIC, qBIC’ modes form doublets in the range $0 < \delta /C < 0.5$ with complex conjugate values of energy $\varepsilon$ and vanishing real component, $\mathrm{Re}(\varepsilon )=0$ and coalesce at the EP at $\delta /C = 0.5$. In the range between $0.5 < \delta /C < 3.1$, the common value of $\mathrm{Re}(\varepsilon )$ for the L, L’ (or qBIC, qBIC’) doublet increases (decreases) and approaches the upper (lower) edge of the continuum, whereas $\mathrm{Im}(\varepsilon )$ for the L, qBIC and L’, qBIC’ pairs coalesce. Due to the identical value of $\mathrm{Im}(\varepsilon )$, the distinction between L and qBIC (or L’ and qBIC’) in this range is arbitrary, in Figs. 5(a) and (b) we used the notation corresponding to the limit $C_1 \rightarrow C_2$ when $C_1 < C_2$. Beyond another EP at $\delta /C = 3.1$ the L, L’ and qBIC, qBIC’ doublets split into the singlets and transform into the bound states. In fact, we observe here hybridized modes which attract each other in the range $0.5 < \delta /C < 3.1$ when one of interacting mode has negative energy, whereas the interaction of two modes with positive energy leads to level repulsion [16] for $\delta /C > 3.1$.

The structures with nonequal coupling coefficients $C_1\neq C_2$ exhibit, when increasing $\delta /C$ from zero, thresholdless behavior, see Figs. 5(c) and (d). For $C_1<C_2$, the common value of $\mathrm{Re}(\varepsilon )$ for the L, L’ (or qBIC, qBIC’) doublets starts to increase (decrease) with $\delta /C$ and enters the gap above (below) the continuum. The decrease of $\mathrm{Re}(\varepsilon )$ for qBIC is described by the approximation leading to Eq. (20) in [12]. Indeed, for $C_1>C_2$, the previously described trends are the opposite, such structures will be considered in Fig. 6(a). With further increasing $\delta /C$ the doublets L, L’ and qBIC, qBIC’ reach EPs at which they also transform into the bound states with $\mathrm{Im}(\varepsilon )=0$.

 figure: Fig. 6.

Fig. 6. (a) Real part of the effective mode index Re$(N_{\rm eff})$ vs. the asymmetry $\Delta n$ for the structures with $a=1.6\,\mu$m, $d_1=1.7\,\mu$m, and $d_2=2.0\,\mu$m. Thick lines: COMSOL simulation of B1, B2, qBIC a L modes. Thin lines: CMT calculation (the colors correspond with those used in Figs. 4 and 5). The gray area indicates the continuum as obtained from CMT. The black dash-dotted horizontal lines mark edges of the continuum as obtained from COMSOL. (b) The loss $b$ of qBIC mode obtained with COMSOL vs. the asymmetry $\Delta n$ for the structures with $a=1.6\,\mu$m, $d_1$ specified in the graph, and $d_2=a$. Detailed description of the structures used in (a) and (b) is in the text.

Download Full Size | PDF

5. Validation of CMT results

CMT is an approximate technique, which is valid only for structures consisting of weakly coupled waveguides. Moreover, in our CMT formulation, we neglected the coupling between the additional waveguides and the waveguides with $m=-1$ and $m=1$ in the array (nonlocal coupling). Naturally, one may wonder whether the presented classification still holds for realistic structures, namely, in the cases when the mentioned approximations are not adequate. To answer the question we carried out a series of rigorous COMSOL simulations and compared their results with CMT predictions. For the calculations we assumed waveguide arrays being fabricated by direct 3D laser writing combined with infiltration, often employed in experiments with coupled systems [18,19]. We consider propagation of light dominantly polarized in $x$-direction (see Fig. 1) at wavelength $\lambda =750\,$nm. All the waveguides possess a circular cross section with radius $r=0.5\,\mu$m. Refractive indices are $n_{\rm core}=1.59$ for waveguides in the array, $n_{\mp 1}=n_{\rm core} \mp \Delta n$ for upper and lower additional waveguide, respectively, and $n_{\rm back}=1.59$ for the background medium. Here we introduced the asymmetry parameter $\Delta n$ that plays similar role as the asymmetry strength $\delta$. The waveguides are single mode (for the additional waveguides in the range of $\Delta n$ considered here) with the propagation constant $\beta _0 = 1.554296k_0$, $k_0=2\pi /\lambda$ (for $\Delta n=0$). Centre-to-centre distance between the upper/lower additional waveguide and waveguide array is $d_{1,2}$ (Fig. 1).

The CMT parameters were determined as follows. The energy is defined as $\varepsilon =\beta - \beta _0 - 2\delta _{\rm arr}$, where $\beta \equiv N_{\rm eff} k_0$ is the propagation constant of a mode of the composed structure, $N_{\rm eff}$ is the effective mode index and $\delta _{\rm arr}$ is the self-coupling coefficient in the array. $C$ and $\delta _{\rm arr}$ were calculated from the propagation constants of supermodes of a structure consisting of two waveguides at distance $a$ (for details of the numerical procedure see, e.g., [20]); similarly we obtained $C_{1,2}$. $\delta _{\mp 1}$ was calculated from the propagation constant of a mode of an individual additional waveguide with specified $n_{\mp 1}$. In the CMT calculations we neglected possible dependencies of $\delta _{\mp 1}$ on $d_{1,2}$ and $C_{1,2}$ on $n_{\mp 1}$; however, we numerically verified that these effects have marginal impact on the presented results.

For COMSOL simulations we considered finite structures consisting of 59 waveguides in the array and the outgoing boundary conditions were modelled using the perfectly matched layers (PML). We fixed the period of the array as $a=1.6\,\mu$m, resulting in the coupling coefficient $C=0.0176517\,\mu {\rm m}^{-1}$ and the coupling length about $90\,\mu$m. Note that structures with lengths up to 1 mm can be fabricated [19]. In the first example, we simulated structures with $d_1=1.7\,\mu$m and $d_2=2.0\,\mu$m ($C_1=0.801C$, $C_2=0.434C$). Figure 6(a) shows evolution of the spectrum when the asymmetry parameter $\Delta n$ increases in the range corresponding approximately to $0<\delta /C<4$. The results obtained from the rigorous COMSOL simulations indicated by thick lines in Fig. 6(a) obviously confirm the main spectral features, namely the existence of B1, B2, qBIC a L modes (qBIC’ and L’ are not accessible due to the used outgoing boundary conditions and for the same reason one cannot observe EPs) and exhibit good qualitative agreement in trends predicted with CMT (thin lines; for the reference we also included qBIC’ and L’ modes). However, Fig. 6(a) also reveals certain differences between the techniques, mostly due to the limits of CMT for the selected distances $a$, and $d_{1,2}$. For example, the edges of the continuum determined with COMSOL (the black dash-dotted horizontal lines; their positions were determined from simulation of one period of the array) do not exactly agree with CMT (gray area). We verified numerically, that the agreement cannot be improved by considering the next-nearest-neighbor coupling in CMT, in fact, the accuracy of CMT becomes even worse (similar observation was done in [14]). On the other hand, we confirmed that most of the differences disappear with increasing $a$, $d_{1,2}$ (weaker coupling). There is an exception from the latter rule concerning modes near band edges: for example BIC and L are always embedded in the continuum, even for strong values of $\Delta n$ ($\Delta n>0.008$ for qBIC, $\Delta n>0.01$ for L), in contrast to the CMT predictions. Finally, CMT describes infinite structures. The discontinuities observed in the dependencies obtained with COMSOL for qBIC and L modes arise from the discretization of the continuum due to the finite extent of the simulated structures.

Considering imaginary part of $N_{\rm eff}$, rigorous simulation provides only very crude approximation and therefore cannot be used as a reference. This is due to PML that are not suitable for systems where the refractive index changes in the direction perpendicular to the boundary in the PML region [21]. However, Im($N_{\rm eff}$) can still be used for finding of structures supporting BIC [i.e., modes with negligible Im($N_{\rm eff}$)]. This is used in Fig. 6(b) where we show loss of qBIC mode calculated with COMSOL as a function of asymmetry parameter $\Delta n$ in narrower interval around zero for structures with various $d_1$, specified in the annotations of the individual curves, and $d_2=a=1.6\,\mu$m. Corresponding values of $C_1/C$ are 2.31, 1.50, 0.985, 0.434 and 0.196 ($d_1$ in increasing order). These parameters are clearly beyond the assumptions of our CMT model and include systems with strong geometrical asymmetry in the vertical direction (i.e., in $C_{1,2}$). Even in this regime BIC is always observed, however, in contrast to CMT predictions, its position has to be carefully tuned by (generally) non-zero value of $\Delta n$. The shift of $\Delta n$ indicates the possible effect of the nonlocal coupling. The curves illustrate that when $d_1$ is further increased, the position of BIC becomes stable and depends only on $d_2$. We also confirmed, that in the limit of weak coupling (with $a$ and $d_{1,2}$ being sufficiently large) BIC always occurs at $\Delta n\approx 0$, regardless of actual values of $d_{1,2}$.

6. Factorization approach

The scattering properties of an arbitrary structure can be characterized by the scattering matrix, which, in frame of our CMT model, has two eigenvalues [12] $S_o=r-t=-1$ and

$$S_e=r+t=2t-1=\frac{2C\sin(k_xa)-i \mu}{2C\sin(k_xa)+i \mu}.$$
The dependence $S_e(q)$, $q\equiv k_xa$, can be expressed through the Weierstrass factorization theorem [10,11] which, for our structure, takes the simple form [12]
$$S_e(q) = \prod_n {\frac{q - p_n^{{\ast}}}{q - p_n}},$$
where $p_n$ are poles of $S_e(q)$. It follows from Eq. (21) that the poles $p_n$ are solutions of the eigenvalue problem (17) in terms of $k_xa$. Therefore, the knowledge of eigenmodes is sufficient to restore $S_e(q)$ through Eq. (22) and to evaluate the scattering spectra as $r=(S_e-1)/2$ and $t=(S_e+1)/2$, where, in the reflectance $R=|r|^{2}$, each pole creates a peak of Lorentzian shape [10]. This way, we can interpret observed spectral features in terms of the eigenmodes. The procedure is illustrated in Fig. 7 where, in addition to the reflectance (solid line), we also display contributions of the individual modes (shaded areas) to the total spectrum. The used CMT parameters correspond to the structure from Fig. 6(a) with $\Delta n \;{\dot =}\; 0.0035$. To demonstrate the effect of poles corresponding with all 6 previously described eigenmodes we extended the range of $k_x$ (B2’ refers to the pole B2 shifted by $-2\pi$). Clearly, B1 and B2 are responsible for peaks at band edges whereas L and qBIC for the peaks in the continuum. Even L’ and qBIC’ (or other poles such as B2’) may affect the spectrum in the physical range $0\le k_xa \le \pi$ provided that the overlap with their tails is sufficiently large.

 figure: Fig. 7.

Fig. 7. The reflectance $R$ vs. the Bloch wavenumber $k_x$ for the system with the parameters $C_1 = 0.801C$, $C_2 = 0.434C$, $\delta =-\delta _{-1}=\delta _{1}=C$, and $\delta _{0}=0$. The shaded areas correspond to the calculations when only the indicated pole is used in the expansion given by Eq. (22).

Download Full Size | PDF

The destructive interference mainly between L and qBIC modes leads to the asymmetric resonance associated with qBIC mode and minima $R=0$. This feature corresponds to the asymmetric dip and maxima $T=1$ in the transmittance, such as in Fig. 3(b). In fact, it can be shown [12] that superposition of L and qBIC modes, playing the roles of the continuum and a narrow discrete state with poles $p_c$ and $p_s$, respectively, leads to the Fano formula

$$T(\Omega)=\left(\frac{1}{1+f^{2}}\right)\frac {\left(\Omega +f\right)^{2}}{\Omega^{2}+1},$$
where
$$\Omega = \frac{q-\mathrm{Re}(p_s)}{\mathrm{Im}(p_s)}, \;\; f \approx \frac{\mathrm{Im}(p_c)}{\mathrm{Re}(p_c)-\mathrm{Re}(p_s)}.$$

The approximation used in deriving of Fano formula (23) manifests itself in the deviation from the the transmittance calculated from the scattering matrix – Fig. (8), in particular in the structures with symmetric coupling coefficients – see Fig. 8(a)

 figure: Fig. 8.

Fig. 8. Comparison of the Fano formula with the exact solution: The transmittance $T$ vs. the Bloch wavenumber $k_x$ for systems with $\delta /C=0.3$ and the indicated coupling coefficients.

Download Full Size | PDF

7. Conclusion

In conclusion, we studied systematically spectral and scattering properties of a photonic analog of an extended Fano-Anderson tight-binding model. We classified the modes according to the relation between the self-coupling coefficients and the change of the propagation constant. The classification of the modes predicted by theoretical model was verified by means of full-wave rigorous numerical simulations of realistic photonic structures with the parameters beyond the limits of the approximations used in the model. We confirmed the main findings of the model, namely, the occurrence of B1, B2, qBIC and L modes, the trends in their behavior and the existence of BIC. BIC is found even in the systems with strong asymmetry in positions of the additional waveguides as well as in the regime when nonlocal coupling may become non-negligible.

The resonant features in the scattering spectra of the structures with broken symmetry can be interpreted through a generalized Weierstrass factorization. In particular, the Fano resonance associated with quasi-BIC arises from the interference between this state and another leaky mode. The transmittance near the quasi-BIC resonance can be rewritten into the form of the Fano formula, where the shape parameter $f$ can be expressed in terms of the poles of the interacting modes. Our paper provides the theoretical framework which describes transformation of the symmetry protected BIC into a leaky mode and allows to interpret the resonant properties of the more complex systems. We envisage extending our approach to be applied to investigation of lossy structure which offer a possibility to explore their topological properties through the engineering of zeros of the transmission matrix, which have shown to subjected to strong topological constraints [22] and may prove to be useful in various fields of optics as well as in cold matter and quantum confined systems.

Funding

Grantová Agentura České Republiky (1900062S); The Faculty of Mechanical Engineering of Brno University of Technology (FSI-S-20-6353).

Acknowledgments

We acknowledge financial support by the Czech Science Foundation (CSF) through Project No. 1900062S. JP acknowledges the financial support provided by the Faculty of Mechanical Engineering of Brno University of Technology through Project No. FSI-S-20-6353.

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. J. von Neumann and E. Wigner, “Über merkwürdige diskrete Eigenwerte,” Phys. Z. 30, 465–467 (1929).

2. A. Krasnok, D. Baranov, H. Li, M.-A. Miri, F. Monticone, and A. Alu, “Anomalies in light scattering,” Adv. Opt. Photonics 11(4), 892–951 (2019). [CrossRef]  

3. C. W Hsu, B. Zhen, A. D. Stone, J. D. Joannopoulos, and M. Soljačić, “Bound states in the continuum,” Nat. Rev. Mater. 1(9), 16048 (2016). [CrossRef]  

4. S. I. Azzam and A. V. Kildishev, “Photonic Bound States in the Continuum: From Basics to Applications,” Adv. Opt. Mater. 9(1), 2001469 (2021). [CrossRef]  

5. K. Koshelev, A. Bogdanov, and Y. Kivshar, “Engineering with the Bound States in the Continuum,” Opt. Photonics News 31(1), 38 (2020). [CrossRef]  

6. Y. Plotnik, O. Peleg, F. Dreisow, M. Heinrich, S. Nolte, A. Szameit, and M. Segev, “Experimental Observation of Opical Bound States in the Continuum,” Phys. Rev. Lett. 107(18), 183901 (2011). [CrossRef]  

7. G. D. Mahan, Many-Particle Physics (Kluwer Academic/Plenum Press, 2000), p. 44; pp. 207–217.

8. A. E. Miroshnichenko, S. F. Mingaleev, S. Flach, and Y. S. Kivshar, “Nonlinear Fano resonance and bistable wave transmission,” Phys. Rev. E 71(3), 036626 (2005). [CrossRef]  

9. B. P. Nguyen, K. Kim, F. Rotermund, and H. Lim, “Transmission resonance induced by a δ-like defect in the Fano-Anderson model with two Fano defects,” Phys. Status Solidi B 249(9), 1765–1770 (2012). [CrossRef]  

10. V. Grigoriev, A. Tahri, S. Varault, B. Rolly, B. Stout, J. Wenger, and N. Bonod, “Optimization of resonant effects in nanostructures via Weierstarss factorization,” Phys. Rev. A 88(1), 011803(R) (2013). [CrossRef]  

11. V. Grigoriev, S. Varault, G. Boudarham, B. Stout, J. Wenger, and N. Bonod, “Singular analysis of the Fano resonances in plasmonic nanostructures,” Phys. Rev. A 88(6), 063805 (2013). [CrossRef]  

12. J. Petráček and V. Kuzmiak, “Effect of symmetry breaking on bound states in the continuum in waveguide arrays,” Phys. Rev. A 105(6), 063505 (2022). [CrossRef]  

13. A. Hardy and W. Streifer, “General coupled-mode theory in non-Hermitian waveguides,” J. Lightwave Technol. 3(5), 1135–1146 (1985). [CrossRef]  

14. X. Shi, X. Chen, B. A. Malomed, N. C. Panoiu, and F. Ye, “Anderson localization at the subwavelength scale for surface plasmon polaritons in disordered arrays of metallic nonowires,” Phys. Rev. B 89(19), 195428 (2014). [CrossRef]  

15. S. Weimann, Y. Xu, R. Keil, A. E. Miroshnichenko, A. Tünnermann, S. Nolte, A. A. Sukhorukov, A. Szameit, and Y. S. Kivshar, “Compact Surface Fano States Embedded in the Continuum of Waveguide arrays,” Phys. Rev. Lett. 111(24), 240403 (2013). [CrossRef]  

16. N. R. Bernier, L. D. Tóth, A. K. Feofanov, and T. J. Kippenberg, “Level attraction a a microweve optomechanical circuit,” Phys. Rev. A 98(2), 023841 (2018). [CrossRef]  

17. W. D. Hess, “Physics of exceptional points,” J. Phys. A: Math. Theor. 45(44), 444016 (2012). [CrossRef]  

18. Z. Fedorova, C. Jörg, C. Dauer, F. Letscher, M. Fleischhauer, S. Eggert, S. Linden, and G. von Freymann, “Limits of topological protection under local periodic driving,” Light: Sci. Appl. 8(1), 63 (2019). [CrossRef]  

19. J. Schulz, S. Vaidya, and C. Jörg, “Topological photonics in 3D micro-printed systems,” APL Photonics 6(8), 080901 (2021). [CrossRef]  

20. J. Petráček and V. Kuzmiak, “Transverse Anderson localization of channel plasmon polaritons,” Phys. Rev. A 98(2), 023806 (2018). [CrossRef]  

21. A. F. Oskooi, L. Zhang, Y. Avniel, and S. G. Johnson, “The failure of perfectly matched layers, and towards their redemption by adiabatic absorbers,” Opt. Express 16(15), 11376 (2008). [CrossRef]  

22. Y. Kang and A. Z. Genack, “Transmission zeroes with topological symmetry in complex systems,” Phys. Rev. B 103(10), L100201 (2021). [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 (8)

Fig. 1.
Fig. 1. A coupled system consisting of a periodic array with two additional waveguides above and below the array.
Fig. 2.
Fig. 2. The transmittance $T$ of the waveguide array. (a), (b) Systems with the symmetry $\delta _{-1}= \delta _1$; (c), (d) asymmetric systems defined as $\delta =-\delta _{-1}=\delta _1$; $k_x$ is the Bloch wavenumber, $C_{1,2}$ are shown in the boxes, $\delta _0=0$. The dashed line in (a) and (b) displays the dispersion relation (5) with $\delta _{\pm 1} = \varepsilon$.
Fig. 3.
Fig. 3. The transmittance $T$ vs. the Bloch wavenumber $k_x$ for systems with preserved $\delta =0$ (dashed line) and broken $\delta /C=0.3$ (solid line) symmetry for the coupling coefficients indicated and $\delta _0=0$.
Fig. 4.
Fig. 4. Spectrum of structures with (a), (b) trivial, $\delta _{\pm 1}=\delta _0=0$, and (c), (d) nontrivial, $\delta _{-1}=\delta _{1}$, vertical symmetry. (a) The real and (b) imaginary part of the energy $\varepsilon$ as a function of $C_{\rm av}/C$. (c) The real and (d) imaginary part of $\varepsilon$ as a function of $\delta _{\pm 1}/C$ for $\delta _0=0$ and $C_1=C_2=C$. The modes in (c) and (d) are labeled according to their asymptotics when $\delta _{\pm 1} \rightarrow 0$. Dashed lines indicate the states with $\mathrm{Im}(\varepsilon )>0$. The gray area in (a) and (c) indicates the continuum.
Fig. 5.
Fig. 5. Spectrum of structures with broken vertical symmetry $\delta =-\delta _{-1}=\delta _1$: the energy $\varepsilon$ vs. the asymmetry strength $\delta$. (a) and (b) $C_1=C_2=C$; (c) and (d) $C_1/C=0.1$, $C_2/C=0.4$; $\delta _0=0$. The gray area in (a) and (c) indicates the continuum.
Fig. 6.
Fig. 6. (a) Real part of the effective mode index Re$(N_{\rm eff})$ vs. the asymmetry $\Delta n$ for the structures with $a=1.6\,\mu$m, $d_1=1.7\,\mu$m, and $d_2=2.0\,\mu$m. Thick lines: COMSOL simulation of B1, B2, qBIC a L modes. Thin lines: CMT calculation (the colors correspond with those used in Figs. 4 and 5). The gray area indicates the continuum as obtained from CMT. The black dash-dotted horizontal lines mark edges of the continuum as obtained from COMSOL. (b) The loss $b$ of qBIC mode obtained with COMSOL vs. the asymmetry $\Delta n$ for the structures with $a=1.6\,\mu$m, $d_1$ specified in the graph, and $d_2=a$. Detailed description of the structures used in (a) and (b) is in the text.
Fig. 7.
Fig. 7. The reflectance $R$ vs. the Bloch wavenumber $k_x$ for the system with the parameters $C_1 = 0.801C$, $C_2 = 0.434C$, $\delta =-\delta _{-1}=\delta _{1}=C$, and $\delta _{0}=0$. The shaded areas correspond to the calculations when only the indicated pole is used in the expansion given by Eq. (22).
Fig. 8.
Fig. 8. Comparison of the Fano formula with the exact solution: The transmittance $T$ vs. the Bloch wavenumber $k_x$ for systems with $\delta /C=0.3$ and the indicated coupling coefficients.

Equations (24)

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

i ψ m = C ( ψ m 1 + ψ m + 1 ) , m   0
i ψ 0 = δ 0 ψ 0 + C ( ψ 1 + ψ 1 ) + C 1 φ 1 + C 2 φ 1
i φ 1 = δ 1 φ 1 + C 1 ψ 0
i φ 1 = δ 1 φ 1 + C 2 ψ 0 ,
ε = 2 C cos ( k x a ) .
ψ m = e i k x a m i ε z + r e i k x a m i ε z , m < 0
ψ m = t e i k x a m i ε z , m > 0
φ ± 1 = B ± 1 e i ε z .
t = r + 1 = 2 C sin ( k x a ) 2 C sin ( k x a ) + i μ
μ = δ 0 + C 1 2 ε δ 1 + C 2 2 ε δ 1 .
ψ m = A e i k x a | m | i ε z ,
[ δ 0 ε + 2 C exp ( i k x a ) ] A + C 1 B 1 + C 2 B 1 = 0
C 1 A + ( δ 1 ε ) B 1 = 0
C 2 A + ( δ 1 ε ) B 1 = 0.
C 1 B 1 = C 2 B 1 , A = 0
B 1 = C 1 , 2 ε δ 1 A .
2 C sin ( k x a ) + i μ = 0.
ε = ± 2 C 2 ± 2 C 4 + C a v 4 ,
C 2 B 1 = C 1 B 1 .
C 2 B 1 ( δ 1 ε ) = C 1 B 1 ( δ 1 ε ) .
S e = r + t = 2 t 1 = 2 C sin ( k x a ) i μ 2 C sin ( k x a ) + i μ .
S e ( q ) = n q p n q p n ,
T ( Ω ) = ( 1 1 + f 2 ) ( Ω + f ) 2 Ω 2 + 1 ,
Ω = q R e ( p s ) I m ( p s ) , f I m ( p c ) R e ( p c ) R e ( p s ) .
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.