## Abstract

We propose a method to obtain the resonance frequencies of coupled optical modes for a stack of two periodically corrugated slabs. The method is based on the modes in each slab, which are derived by the Fourier modal method in combination with the optical scattering matrix theory. We then use the resonant mode approximation of the scattering matrices to develop a linear eigenvalue problem with dimensions equal to the number of resonant modes. Its solutions are the resonance frequencies of the coupled system and exhibit a good agreement with exact solutions. We demonstrate the capabilities of this method for pairs of planar waveguides and gratings of one-dimensional wires.

©2010 Optical Society of America

## 1. Introduction

Nanooptics has drawn a lot of attention in the recent years. Especially the field of plasmonics and its application for constructing metamaterials has a great potential. However, the numerical calculations of the optical properties are often very time-consuming and, therefore, simple models are desired for estimating these properties.

Some attempts have been made so far to describe the coupling of nanostructures [1–3], but they require fitting parameters, because not all quantities are directly accessible from Maxwell’s equations. In this paper we present an approach for deriving the coupling in the vertical stack of two periodic layers based on the resonances in each layer. Our approach follows directly from Maxwell’s equations without any fitting parameters by using the optical scattering matrix theory [4, 5] and the resonant mode approximation of the scattering matrix [6, 7]. Our master equation will be a linear eigenvalue problem with dimensions equal to the total number of resonances in both layers. Due to the consideration of optical near fields, this method is even valid for the strong coupling regime.

In this paper, we consider the scattering matrix derived in a layer by layer calculation that is also known as the Fourier modal method (FMM) or rigorous coupled wave analysis (RCWA). First, the fields in each layer are obtained by an expansion in Fourier Bloch modes. The numerical accuracy of the solutions is limited by the truncation of the involved Fourier transforms and has been improved by various tools such as factorization rules [8] and adaptive spatial resolution [9]. Then, the scattering matrix of the whole structure is constructed iteratively from the solutions in each layer under consideration of the boundary conditions.

We will first explain the derivation of optical eigenmodes from the scattering matrix and the resonant mode approximation. Then, we will demonstrate how to obtain the coupling equations from the resonant mode approximation of two scattering matrices. The resulting equations can be used not only in the Fourier modal method, but also in other approaches, e.g., in the multiple-scattering method [10]. Finally, we will compare the results of the coupling equations with the full calculation of stacks of waveguides as well as gratings of dielectric and metallic wires.

## 2. Resonant mode approximation

The scattering matrix provides the optical output ∣*O*〉 of a stack of layers for known input ∣*I*〉 as follows:

where ∣*u*
_{b,t}〉 and ∣*d*
_{b,t}〉 define the amplitudes of the upwards and downwards propagating or
evanescently decaying part of the solutions in the bottom and top layer, respectively [see
Fig. 1(a)].

Optical eigenmodes are solutions at a complex frequency Ω–*i*Γ that can oscillate in time with frequency Ω and decay rate Γ without external excitation [5]. Hence, we must set ∣*I*〉 = 0, invert Eq. (1), and solve the problem S^{-1}(*ω*,**k**)∣*O*〉 = 0 with respect to complex frequency *ω* for a real in-plane momentum **k**. Because we fix **k** to a specific value, we henceforth suppress the dependency of the scattering matrix on this variable in our notations.

An efficient way for obtaining the optical resonances from the inverse scattering matrix is to linearize the problem:

Equation (2) defines a generalized linear eigenvalue problem for *δ* = *ω* -*ω*
_{0} that has to be solved iteratively until at least one eigenvalue *δ* ≈ 0 is found in order to determine an optical eigenmode. We can then define a diagonal matrix Δ containing all the eigenvalues and another matrix X containing all the eigenvectors as columns. Thus, the solution of Eq. (2) can be written as

By replacing the first term in the brackets in Eq. (2) by Eq. (3) and inverting, we obtain an approximated form of the scattering matrix that is similar to the Breit-Wigner formula known in quantum mechanical scattering theory (e.g., see chapt. 124 and 142 in [11]):

where

and 𝕀 denotes the identity matrix. The sum on the right side of Eq. (4) involves all terms with eigenvalues *δ _{n}* that can be considered as eigenmodes at

*ω*=

_{n}*ω*

_{0}+

*δ*with input resonance vectors ∣

_{n}*I*〉 and output resonance vectors ∣

_{n}*O*〉. The output vector ∣

_{n}*O*〉 is the eigenvector of Eq. (2) corresponding to the eigenvalue

_{n}*δ*and, hence, a column vector in X, whereas 〈

_{n}*I*∣ is a row vector in the matrix Y. It should be mentioned that the eigenvectors in X are generally not orthogonal, which means that 〈On∣ is not defined as the Hermitean conjugate of vector ∣

_{n}*O*〉, but as a row vector of X

_{n}^{-1}. Similarly, ∣

*I*〉 is a column vector in Y

_{n}^{-1}.

The background term S̃ can be either derived by subtracting the resonance terms in Eq. (4) from exact results or, as in the examples provided below, calculated straight-forward from Eq. (4). S̃ is assumed to be energy independent, as long as the other resonances are far from the region of interest *ω* ≈ *ω*
_{0} (so that the absolute values of the corresponding eigenvalues *δ* are much larger than ∣*ω* -*ω*
_{0}∣). Another limitation for the approximation of S̃ being energy independent is the opening of new diffraction orders.

## 3. Resonant mode coupling

We want to derive now the resonances of a stacked system based on the approximated scattering matrices S^{a} and S^{b} of each part:

The indices 1 and 3 refer to the superstrate and substrate of the coupled system, respectively [see Fig. 1(a)]. In order to obtain the coupled resonances, the vectors ∣*d*
_{2}〉 and ∣*u*
_{2}〉 between the structures have to be matched. Similar to the single system all inputs of the coupled system are zero, namely ∣*d*
_{1}〉 and ∣*u*
_{3}〉. Hence, we obtain the simplified equations:

We can introduce now some abbreviations for Eqs. (8) and (9):

Employing Eq. (10) in Eq. (11) and vice versa, we obtain with the relation (𝕀-AB)^{-1} A=A(𝕀-BA)^{-1}:

where

The terms in front of the resonant output vectors represent the multiple reflection between the two systems, because (𝕀-AB)^{-1} = ∑^{∞}
_{n=0}(AB)^{n}.

If we replace ∣*u*
_{2}〉 in the definition of *α _{n}* [Eq. (10)] by Eq. (13) and ∣

*d*

_{2}〉 in the definition of

*β*[Eq. (11)] by Eq. (12), we derive the linear system of equations

_{n}$$\left(\omega -{\omega}_{n}^{\text{b}}\right){\beta}_{n}=\sum _{m=1}^{N}\u3008{I}_{d,n}^{\text{b}}\mid {\mathbf{O}}_{d,m}^{a}\u3009{\alpha}_{m}+\sum _{m=1}^{M}\u3008{I}_{d,n}^{b}\mid {\tilde{S}}_{du}^{a}{\mid \mathbf{O}}_{u,m}^{b}\u3009{\beta}_{m}.$$

These equations provide non-trivial solutions for *α _{n}* and

*β*only for particular frequencies that define the resonances of the coupled system.

_{n}In general, these equations are very similar to the linear combination of atomic orbitals in quantum mechanics. The scalar products in the optical coupling equations can be also identified as overlap integrals of the resonant fields, but are modified by the background scattering matrices and the multiple reflections.

A great advantage of the above Eq. (15) is that once the solution of the resonances of each scattering matrix is derived, we can easily add a spacer layer to the approximated scattering matrices for calculating different distances (see Appendix D in Ref. [5]). Furthermore, we can introduce lateral shifts between the single systems without large computational effort, but we restrict ourselves in the following examples to the distance dependence of the mode coupling for identical layers. It should be mentioned that the propagation in homogeneous media is not energy independent. Hence, additional energy dependence has to be introduced in Eq. (15), but then the eigenvalue problem becomes nonlinear. In some cases (e.g., in the examples below) the energy dependence can be neglected and it is possible to consider the propagation between the layers in a zero order approximation to be constant.

## 4. Results

A simple test system is two coupled dielectric slab waveguides. The calculation of waveguide modes with one harmonic is exact in the Fourier modal method and, therefore, ideal for checking. We have chosen a system, where each waveguide is surrounded by air and consists of a dielectric medium with a refractive index of *n* = 2.5 and a thickness of 50nm [see Fig. 1(a)]. The bare waveguide resonance in TM-polarization is located at 1039.9THz for an in-plane wavenumber *k* = 28*μ*m^{-1}. The results of the coupling equation (15) are shown in Fig. 2(a) for a varying spacer thickness from 10nm to 200nm and agree very well with the analytical results.

In the second test we cut trenches of 100nm width with a period of 300nm in each waveguide [Fig. 1(b)] in order to study the interaction of dielectric gratings. Such gratings support optically active and inactive modes [5]. The latter merely decay exponentially in the direction normal to the slab with no interaction to propagating far field channels. The results of an optically inactive TM-mode (bare resonance of the single grating at 930.6THz) are displayed for normal incidence (*k* = 0) in Fig. 2(b) and still exhibit a very good agreement between the results of Eq. (15) and the numerical calculation of the coupled system.

Figure 2(c), 2(d) presents the results of an optically active mode in TM-polarization (bare
resonance at 900.8 - 2.5*i*THz) for (c) the real part and (d) the imaginary part of the coupled resonance. Even the oscillating behavior due to the weak coupling of the resonances with Fabry-Perot modes for distances larger than 200nm is reproduced quite well. However, strong coupling with Fabry-Perot modes would require an extension of the current equations and is beyond the scope of this paper. Another limitation of the method is the opening of additional diffraction orders that is not taken into account in the resonant mode approximation. In addition, also the derivation of the eigenmodes has some uncertainty in the grating system due to the truncation of the numerical calculation (31 harmonics with adaptive spatial resolution).

As the third test, we calculated the optical properties of stacked metallic gratings consisting of gold wires with height 15nm, width 100nm, and period 300nm [see Fig. 1(c)] at normal incidence. Such pairs of metallic wires are of high interest in the field of plasmonics and metamaterials for constructing antisymmetric resonances with a strong magnetic field pointing in the direction of the wires [12]. The surrounding medium is quartz (*n* = 1.46) and the gold is described by a simple analytical model [13]. Figure 3(a) shows the comparison of the resonance position of a localized TM-mode (bare resonance at 418.6-42.9*i*THz) derived by our model and the numerical solution of the coupled system for increasing distances. Still, the agreement is quite good for a truncation order of 31 harmonics with adaptive spatial resolution, but some deviation can be seen for larger distances. One reason maybe the energy dependence of the propagation between the layers. We neglected this dependence and approximated it by the propagation at the energy of the bare resonance. Hence, the modes propagate for larger distances longer with an incorrect propagation constant in our calculations. Figure 3(b) shows the imaginary part of the resonances as the decay rate (or half linewidth) of the resonance.

## 5. Conclusion

We have shown how to obtain simple coupling equations for the optical modes of two stacked periodic structures. These equations require only the calculation of the modes in each structure, but then provide the energy position without any further fitting procedures in dependence of the distance. The results are in good agreement with the full calculation of the optical properties.

## Acknowledgement

We acknowledge support from BMBF (13N 9155, 13N 10146), DFG (FOR 557, FOR 730), DFH/UFA, ANR Chair of Excellence Program, the RAS, and the Russian Foundation for Basic Research. Furthermore, we thank Sven Hein for providing the POV-Ray pictures.

## References

**1. **J. Petschulat, C. Menzel, A. Chipouline, C. Rockstuhl, A. Tünnermann, F. Lederer, and T. Pertsch, “Multipole approach to metamaterials,” Phys. Rev. A **78**, 438111 (2008). [CrossRef]

**2. **N. Liu, H. Liu, S. Zhu, and H. Giessen, “Stereometamaterials,” Nat. Photon. **3**, 157–162 (2009). [CrossRef]

**3. **N. Liu, L. Langguth, T. Weiss, J. Kästel, M. Fleischhauer, T. Pfau, and H. Giessen, “Plasmonic electromagnetically induced transparency at the Drude damping limit,” Nat. Mater. **8**, 758–762 (2009). [CrossRef] [PubMed]

**4. **D. M. Whittaker and I. S. Culshaw, “Scattering matrix treatment of patterned multilayer photonic structures,” Phys. Rev. B **60**, 2610–2618 (1999). [CrossRef]

**5. **S. G. Tikhodeev, A. L. Yablonskii, E. A. Muljarov, N. A. Gippius, and T. Ishihara, “Quasiguided modes and optical properties of photonic crystal slabs,” Phys. Rev. B **66**, 451021 (2002). [CrossRef]

**6. **N. A. Gippius, S. G. Tikhodeev, and T. Ishihara, “Optical properties of photonic crystal slabs with an asymmetrical unit cell,” Phys. Rev. B **72**, 451381 (2005). [CrossRef]

**7. **N. A. Gippius and S. G. Tikhodeev, “Application of the scattering matrix method for calculating the optical properties of metamaterials,” Phys.-Usp. **52**, 967–971 (2009). [CrossRef]

**8. **L. Li, “New formulation of the Fourier modal method for crossed surface-relief gratings,” J. Opt. Soc. Am. A **14**, 2758–2767 (1997). [CrossRef]

**9. **G. Granet and J. P. Plumey, “Parametric formulation of the Fourier modal method for crossed surface-relief gratings,” J. Opt. A, Pure Appl. Opt. **4**, 145–149 (2002). [CrossRef]

**10. **G. Gantzounis and N. Stefanou, “Layer-multiple-scattering method for photonic crystals of nonspherical particles,” Phys. Rev. B **73**, 035115 (2006). [CrossRef]

**11. **L. Landau and E. M. Lifshitz, *Quantum Mechanics. Non-Relativistic Theory* (Pergamon, 1966).

**12. **V. Podolskiy, A. Sarychev, and V. Shalaev, “Plasmon modes and negative refraction in metal nanowire composites,” Opt. Express **11**, 735–745 (2003). [CrossRef] [PubMed]

**13. **P. G. Etchegoin, E. C. Le Ru, and M. Meyer, “An analytic model for the optical properties of gold,” J. Chem. Phys. **125**, 164705-1-3 (2006). [CrossRef] [PubMed]