Differential method for modelling dielectric-loaded surface plasmon polariton waveguides

This paper demonstrates the efficiency of the differential method, a conventional grating theory, to investigate dielectric loaded surface plasmon polariton waveguides (DLSPPWs), known to be a potential solution for optical interconnects. The method is used to obtain the mode effective indices (both real and imaginary parts) and the mode profiles. The results obtained with the differential method are found to be in good agreement with those provided by the effective index method or finite elements. The versatility of the differential method is demonstrated by considering complex configurations such as trapezoidal waveguides or DLSPPWs lying on a finite width metal stripe. © 2008 Optical Society of America OCIS codes: (130.2790) Guided waves; (240.6680) Surface plasmons; (050.1960) Diffraction theory. References and links 1. J.-H. Yeh and R. K. Kostuk , “Substrate-mode holograms used in optical interconnects: design issues,” Appl. Opt. 34, 3152-3164 (1995). 2. S. H. Song, S. Park, C. H. Oh, P. S. Kim, M. H. Cho and Y. S. Kim, “Gradient-index planar optics for optical interconnections,” Opt. Lett. 23, 1025-1027 (1998). 3. S. Kawai, “Handbook of optical interconnects,” Marcel Dekker Inc., New-York 2005. 4. J.-C.Weeber, Y. Lacroute, and A. Dereux, “Optical near-field distribution of surface plasmon waveguide modes,” Phys. Rev. B 68, 115401 (2003). 5. T. Holmgaard and S. I. Bozhevolnyi, “Theoretical analysis of dielectric-loaded surface plasmon-polariton waveguides,” Phys. Rev. B 75, 245405 (2007). 6. A. V. Krasavin and A. V. Zayats, “Passive photonic elements based on dielectric-loaded surface plasmon polariton waveguides,” Appl. Phys. Lett. 90, 211101 (2007). 7. B. Steinberger, A. Hohenau, H. Ditlbacher, A. L. Stepanov, A. Drezet, F. R. Aussenegg and J. R. Krenn, “Dielectric stripes on gold as surface plasmon waveguides,” Appl. Phys. Lett. 88, 094104 (2006). 8. B. Steinberger, A. Hohenau, H. Ditlbacher, F. R. Aussenegg, A. Leitner, and J. R. Krenn, “Dielectric stripes on gold as surface plasmon waveguides: Bends and directional couplers,” Appl. Phys. Lett. 91, 081111 (2007). 9. A. Hohenau, J. R. Krenn, A. L. Stepanov, A. Drezet, H. Ditlbacher, B. Steinberger, A. Leitner, and F. R. Aussenegg, “Dielectric optical elements for surface plasmons,” Opt. Lett. 30, 893-895 (2005). 10. R. Kiyan, C. Reinhardt, S. Passinger, A. L. Stepanov, A. Hohenau, J. R. Krenn, and B. N. Chichkov, “Rapid prototyping of optical components for surface plasmon polaritons,” Opt. Express 15, 4205-4215 (2007). #93728 $15.00 USD Received 11 Mar 2008; revised 14 Jun 2008; accepted 15 Jun 2008; published 17 Oct 2008 (C) 2008 OSA 27 October 2008 / Vol. 16, No. 22 / OPTICS EXPRESS 17599 11. S. Massenot, J. Grandidier, A. Bouhelier, G. Colas des Francs, L. Markey, J.-C. Weeber, A. Dereux, J. Renger, M. U. Gonzàlez and R. Quidant, “Polymer-metal waveguides characterization by Fourier plane leakage radiation microscopy,” Appl. Phys. Lett. 91, 243102 (2007). 12. E. Amenogiannis, E. N. Glytsis and T. K. Gaylord, “Determination of guided and leaky modes in lossless and lossy planar multilayer optical waveguides: Reflection pole method andWavevector density method,” J. Lightwave Technol. 17, 929-941 (1999). 13. K. Kawano and T. Kitoh, “Introduction to optical waveguide analysis,” Wiley, New-York (2001). 14. M. Nevière and E. Popov, “Light Propagation in Periodic Media, Differential Theory and Design,” Marcel Dekker Inc., New-York (2003). 15. L. Li, “Formulation and comparison of two recursive matrix algorithms for modelling layered diffraction gratings,” J. Opt. Soc. Am. A 13, 1024-1035 (1996). 16. P. Dawson, F. de Fornel and J.-P. Goudonnet , “Imaging of surface plasmon propagating and edge interaction using a photon scanning tunnelling microscope,” Phys. Rev. Lett. 72, 2927 (1994). 17. S.-D. Wu, T. K. Gaylord, E. N. Glytsis, and Y.-M. Wu , “Three-dimensional converging-diverging Gaussian beam diffraction by a volume grating,” J. Opt. Soc. Am. A 22, 1293-1303 (2005). 18. P. Berini, “Plasmon-polariton waves guided by thin lossy metal films of finite widths: Bound modes of asymmetric structures,” Phys. Rev. B 63, 125417 (2001). 19. R. Zia , J. A. schuller and M. L. Brongersma, “Near-field characterization of guided polariton propagation and cutoff in surface plasmons waveguides,” Phys. Rev. B 74, 165415 (2006).


Introduction
Optical interconnects are promising for very short distance communications applications in their ability to support higher bandwidths than conventional metal wires or stripes used for carrying data from board to board, chip to chip or component to component within a single chip. Several solutions have been studied including free-space (substrate-mode holograms, grin substrates for example) [1,2], plastic optical fibers or integrated optics [3]. In a logic of massive integration, the latter solution seems the more suitable especially if the refractive index contrast of the components is high. This can be obtained with Silicon On Insulator (SOI) technology or metal based components for example. If SOI is a more mature technology, metal based integrated components keep a potential practical interest since they can combine both the electronic and optical properties of metals by using surface plasmon polaritons (SPP).
If SPP are generally propagating along a metal / dielectric interface, it is preferable to confine them in the lateral direction by using metal stripes or dielectric loaded surface plasmon polariton waveguides (DLSPPWs) on a metal film [4,5]. The latter are the more promising to confine SPPs in submicron structures hence their interest for large scale integration. The fact that DLSPP waveguides are intrinsically lossy due to the presence of the metal is not necessary an inconvenient if propagation lengths involved are not too high in order to keep a sufficient signal to noise ratio. It is well known that surface plasmon based devices are strongly polarization dependant and intrinsically not well suited for telecommunications applications. If low polarization dependance is a crucial parameter for transport and access networks due to the presence of polarization dependant losses, it is a criteria of lower relevance for on-board optical interconnects hence the possibility of using surface plasmon polaritons to carry the information.
DLSPPWs have recently rise a great interest [5][6][7][8][9][10][11], in this context, it is highly desirable to dispose of suitable tools for the design of the waveguides. In particular, it is of great importance to assess the fundamental properties of these waveguides such as the complex effective index or the cross-section mode profile from which the mode confinement, the propagation length (such waveguides are intrinsically lossy) or the optimum coupling conditions can be determined. The goal of this paper is to demonstrate the suitability of a diffraction grating theory known as the differential method to study DLSPPWs. In this respect, the differential method can be an easy-to-implement alternative to more conventional method such as FDTD or finite elements. In section 2, the principle of the differential method is reminded and the modifications allowing the study of DLSPPWs are presented. Next a comparison with results obtained for DLSPPWs with other numerical methods, effective index (EIM) and finite elements (FEM), is given in section 3. Once the reliability of the differential method has been established for modelling DLSPPWs, the properties of waveguides having triangularly shaped profile or lying on a finitewidth metal stripe are investigated in section 4.

Principle
The differential method belongs to rigorous vector diffraction theories. It allows the study of diffraction grating and more generally periodic media [14]. Several works during the last decade [15] improved the performance of this method especially in terms of convergence rate and numerical stability. The aim of this section is to remind briefly its principle.
The geometry of the problem to be solved is shown in Fig. 1. It is comprised of a modulated layer (e.g. the grating) and two homogeneous semi-infinite media (noted I and III in Fig. 1). Λ is the grating period, θ the incident angle, δ the azimuthal angle (δ = 0 corresponds to a conical diffraction regime), Ψ is the polarisation angle in the local system of axis of the field and t the grating thickness.
The method relies on the fact that the electromagnetic field is pseudo-periodic in the modulated layer. This field can Fourier-expanded and injected into Maxwell equations. From this, the field in the two homogeneous media on both sides of the grating can be written as plane-waves series called Rayleigh expansions. For example, the more general form of the E y component of the field is given in Eq. 1: (k 0 being the vacuum wavevector, k 0 = 2π/λ ).
A m and B m are the Rayleigh coefficients and correspond respectively to the amplitudes of the downgoing and upgoing waves of the structure. N is the number of Fourier harmonics retained for the calculation. The output of the differential method is the set of A m and B m coefficients in media I and III from which the diffraction efficiencies of the grating can be computed.
The computation of A m and B m coefficients requires the propagation of the electromagnetic field through the modulated region. Writing Maxwell-Faraday and Maxwell-Ampère equations leads to a linear system of coupled differential equations for the transverse components of the field E x , E y , H x and H y . By introducing the Fourier expansions of each field component (using the correct factorization rules [14]), one obtain a system that can be noted: T , the size of this system is 4(2N+1) (subscript m corresponds to the m th Fourier harmonic of the field component). Two different cases can be considered: if the profile of the grating is constant along the vertical coordinate (lamellar grating for example), the system 2 has constant coefficients and it can be solved by calculating the eigenvalues and eigenvectors of matrix M. In the other case, M = M(z) and the system is solved by a standard routine for ordinary differential equations. An example is given in section 3.2.1 where a trapezoidal profile is considered.
Practically the computation of the unknowns of the problem (the Rayleigh coefficients) can be achieved by means of a shooting method leading to the scattering matrix of the modulated zone. At this stage it is worth to note that the use of a scattering matrix instead of a transfer matrix to describe the propagation of the electromagnetic field through the modulated zone is necessary to prevent numerical contaminations [15]. Once the scattering matrix of the grating is known, the electromagnetic field can be computed from the Rayleigh expansions in each homogeneous medium surrounding the structure. In addition, it is also possible to consider modelling of multilayer structures comprising either homogeneous (with arbitrary dielectric permittivity) or periodically modulated layers (having the same period as the initial grating).

Modelling DLSPPWs with the differential method
As reminded in the previous section, the differential method is dedicated to the study of periodic objects. However this method can also be used to investigate isolated structures provided that the grating period is large enough in order to avoid coupling effects between isolated elements. A typical DLSPPW structure is shown in Fig. 2(a). A dielectric waveguide is loaded on a metal layer. A TM polarized incident light field with a plane of incidence parallel to the longitudinal axis of the waveguide (δ = 90°) is used to excite the guided modes (see Fig. 2(b), the incident wave comes now from the bottom of the structure for convenience contrary to Fig. 1). By adjusting the angle of incidence θ , the k y component of the incident wavevector can be chosen such that it matches the phase constant of the considered mode allowing a selective excitation of either the modes guided within the dielectric stripe or the surface plasmon at the gold / air interface.
Because of the presence of the metal, the modes supported by DLSPPW are intrinsically lossy. For a full characterization of these modes, the complex effective index (n e f f = ν + ıα) and the mode profile inside and outside the waveguide must be computed.
The real part ν of the effective index can be obtained from the computation of the reflectivity of the structure versus θ . Indeed, as expected for an ATR (Attenuated Total Reflection) setup, the excitation of a guided mode manifests itself as an absorption dip into the reflectivity angular spectrum. The real part of the effective index is then computed from the angular position θ of the dip by ν = n 1 sin θ sin δ (δ = 90°in this case).
For the computation of α, the most straightforward solution consists in searching for the poles of the scattering matrix of the structure. However, from a numerical point of view this technique requires the use of optimisation routines searching for complex roots that are rather time consuming. A less elaborate, but numerically more efficient solution to compute α is to excite the structure with a finite size beam and to determine the 1/e damping distance L spp of the field intensity guided mode [16]. In this way, α is obtained according the relation L spp = λ /4πα. The excitation of the guided mode with a finite-size beam is performed by using a plane wave decomposition of a gaussian beam (see [17] for example). Expression 3 gives the amplitude of the Gaussian beam and its angular spectrum (plane wave) decomposition. For the numerical simulations, this decomposition is first restricted to a an interval of values for k x and k y and then converted to a discrete sum. Each plane wave of the decomposition is associated to a couple of angles (θ , δ ) (see [17] for example).
Then, the total field is reconstructed and the near field intensity above the waveguide (in the xy plane) can be computed. An exponential fit along the longitudinal axis of the waveguide with the formula I = I 0 exp − y L spp provides directly the damping distance (hence α). To simplify the problem and to reduce the computation time, it is possible to use a two-dimensional Gaussian beam instead of the expression given by Eq. 3 (the beam is infinite along the x-axis).
Concerning the mode profiles, they are calculated in the case of a single plane wave excitation with an angle of incidence θ = arcsin ℜ(n e f f ) n 1 . The profiles are obtained by performing an additional integration of Maxwell equations by taking as initial values the field components obtained from solving the two-point boundary value problem. In addition to the mode profile in the cross-section, it is possible to obtain the confinement factor Γ of the considered guided mode (Γ is defined as the ration between the mode power in the waveguide and the whole space).

Rectangular profile
In this section, the fundamental mode supported by the DLSPPW structure shown in Fig. 2 is studied. It corresponds to the plasmon mode confined along the metal / dielectric interface [5] contrary to conventional waveguides where the energy is usually located near the center of the guide. Without dielectric, only the metal / air plasmon mode exists with an effective index ν =1.0038 (for λ = 1.55 μm on a 100 nm gold film, value obtained with the reflection pole method [12]) and this effective index increases with the presence of the dielectric waveguide. The procedure described in the previous section has been applied on a DLSPPW structure whose properties have already been established by another numerical method. Concerning the simulation parameters, a grating period Λ = 8 μm has been used in order to prevent a significant crosstalk between adjacent waveguides. This period being large compared to the incident wavelength (in the near-infrared range), the number of Fourier harmonics in the field expansion must be increased accordingly to ensure the convergence of the field profile. A total number of Fourier harmonics (2N+1) of 201 has been found to satisfy this convergence criterion and has been used for all the calculations presented. Following the notations of Fig.  2, the structural parameters retained are w=t=600 nm, h=100 nm, λ = 1550 nm, n Glass = 1.6, n w = 1.535, n m = 0.55 + 11.5i and n Air = 1.0 [5]. For this structure, the finite elements method gives n e f f = 1.291 + 2.78 × 10 −3 i (corresponding to a propagation length of 44.4 μm). Results are given in Fig. 3.
In Fig. 3(a), a dip in the reflectivity curve is observed, it is not very pronounced due to the fact that the coupling conditions in ATR configuration are not optimal (gold thickness of 100 nm and λ = 1.55 μm). With the aim of improving the visibility of the dip, the derivative with respect of θ of the reflectivity has been computed in Fig. 3(b) The dip occurs for θ 2 = 53.8 o which corresponds to an effective index value of ν = 1.291, which is the same that the one obtained with FEM for the fundamental guided mode supported by the DLSPPW. Concerning the propagation length, the lossy character of this waveguide is clearly shown in Fig. 3(c) and 3(d). From these data a value of L spp = 44.0 μm is obtained, one more time in good agreement with FEM.
To check the validity of the differential method, several profiles have been studied especially by changing the width of the dielectric waveguide (all the other parameters being unchanged compared the reference structure). The complex effective indices have been determined and compared to those obtained with FEM and the effective index method EIM [5], which is known to be an approximate method rather easy to implement [13]. The results are given in table 1 for both the real part of the effective index and the damping distance. Table 1. Comparison of the results obtained by the differential method (DM) to those provided by the Effective Index Method (EIM) and Finite Elements (FEM) for waveguides of different widths (t=600 nm, n w = 1.535, n m = 0.55 + 11.5ı and λ = 1.55 μm) Good agreements are observed between the different methods, which validates the use of the differential method for studying DLSPPWs. The differences in the propagation lengths obtained with the differential method and FEM could be due to the quality of the exponential fitting, whereas the approximate nature of the EIM is at the origin of the discrepancies.
The intensity distributions of the guided modes in the cross-section of different waveguides are shown in Fig. 4. For these four waveguides, the fundamental mode is clearly confined at the metal / dielectric interface demonstrating the plasmonic nature of this guided mode. The confinement factor Γ has been calculated for each waveguide and the best value is obtained for the 600 nm × 600 nm waveguide (Γ=72.2%). The high field intensity observed at the top of the waveguide especially for the case 300 nm × 300 nm is due to the discontinuity of the E z component of the electric field and that the mode tends to behave as the plasmon mode of the gold / air interface.

Trapezoidal profiles
Practically, DLSPPWs are realized by UV or e-beam lithography and their profiles are generally not perfectly rectangular due to proximity effects or exposition parameters. It is then important to know the impact of the realistic waveguide profile on the propagation properties. The case of trapezoidal waveguides is considered first since it has been observed that such profiles can be obtained for specific exposure parameters in photolithography. It can be anticipated that the effective index of the modes supported by a trapezoidal DLSPPW will be different from that of a rectangular waveguide, however the impact of the waveguide geometrical profile change onto the damping distance of the guided mode is not straightforward and desserves a carefull examination.
Keeping the same notations as in Fig. 2, the trapezoidal profile of the waveguides is defined by two parameters w 1 and w 2 corresponding to the width respectively at the top and the bottom of the waveguide (Fig. 5). In all the subsequent calculations, w 2 is kept at 600 nm whereas w 1 is changed from 600 nm to 510 nm.
The real effective index and the damping distances of the fundamental mode supported by t=300 nm -w=600 nm t=600 nm -w=600 nm t=300 nm -w=300 nm t=600 nm -w=300 nm  each waveguide are given in table 2. For decreasing values of w 1 , a decrease of the real part of the effective index of the guided mode is observed. This can be understood from the reduction of the waveguide cross-section surface when decreasing w 1 . One can also observe that the change of waveguide geometry does not impact significantly the damping distance L spp leading to the conclusion that tolerance exists on the fabrication parameters of the DLSPPW. This conclusion is supported by the results obtained in the case of a triangularly shaped crosssection waveguide obtained with w 1 =0 nm and w 2 =600 nm. In this case, the real part of the effective index of the mode is, as expected, significantly lower than for the rectangular profile (n eff =1.161) whereas the damping distance is only reduced by 5% for the triangular waveguide (L spp =42 μm) compared to the rectangular one. This rather weak impact of parameter w 1 onto the attenuation constant of the guided modes can be understood from the fact that the modes we consider originate from a SPP resonance and thus exhibit a field confinement at the surface of the metal. It is then likely that the properties of those modes are more affected by the metal / dielectric interface conditions rather than by a change in the volume of the waveguide.

Rectangular waveguide on a finite-width metal stripe
In this last paragraph, a situation that is of practical interest in the perspective of DLSPPW integration is investigated. A rectangular DLSPPW is lying onto a finite-width metal stripe instead of an infinitely extended thin metal film (see Fig. 6(a)). It has been already shown that plasmonic waveguides consisting of metal stripe have a cutoff width [18,19]. The goal of this calculation is then to establish the minimum width of the metal stripe that does not degrade the propagation length of a given DLSPPW. To that aim the 600 nm × 600 nm waveguide investigated in section 3.1 is considered. It is deposited on a metal stripe of different widths. From a numerical point of view, this situation is modelled with the differential method as a stack of two lamellar gratings having the same period. The propagation lengths obtained for metal stripe width ranging from 0.8 μm to 4 μm are plotted in Fig. 6(b) and it is shown that a metal stripe width of 3 μm is necessary to recover a propagation length similar to that of an infinite thin gold film (noted L ∞ SPP in Fig. 6(b)). This result can be attributed to the larger damping distance of the SPP modes supported by metal stripes compared to an infinite thin film and demonstrates unambiguously the plasmonic nature of the fundamental mode supported by the DLSPPW.

Conclusion
In summary, The suitability of the differential method for studying dielectric-loaded surface plasmon polariton waveguides has been demonstrated. This numerically effective method provides results in excellent agreement with those obtained with more sophisticated techniques, such as the finite element method for example. The versatility of the differential method has been used to investigate DLSPPWs with realistic profiles. In particular, it has been shown that the properties of DLSPPWs are not significantly modified by changing the geometrical profile of te waveguides, provided that these changes do not impact the metal / dielectric interface. In addition, it has been shown that when deposited on a metal stripe, the damping distance of the DLSPPW modes depends drastically on the stripe width. This property that proves unambiguously the plasmonic nature of the DLSPPW modes must be taken into account for a practical implementation of such waveguides.