Rigorous solution for optical diffraction of a sub-wavelength real-metal slit

Yann Gravel; Yunlong Sheng

doi:10.1364/OE.20.002149

1. Introduction

Light scattering by nano-sized structures such as slit, groove or hole engraved in metal surfaces is a fundamental problem of the optical diffraction. With promising properties of the surface plasmon polariton (SPP) and its possible applications to nanophotonics devices, interest in the physical understanding of this problem has been greatly increased [1,2]. After the discovery the extraordinary optical transmission of the 2D metallic nano-hole arrays [3], the theoretical explanation based on the tunneling of the surface plasmons at the front and rear interfaces of the metal layer, excited by both the free electrons in the metal and the array of perforated nano-holes [4,5], has been challenged by a debate on the role of the surface plasmons [6]. Lezec et al. proposed that the interference of the evanescent waves launched by the diffraction of a slit with the wave incident on an adjacent slit may enhance and suppress the transmission [7], which was proved by the numerical demonstration [8]. Composite diffracted evanescent waves (CDEW) have been considered in the model [9]. On the same line, Schouten et al. measured the total transmission energy of two-slits as a function of the slit separation [10]. Gay et al. [11,12] launched surface waves from a sub-wavelength groove towards an adjacent sub-wavelength slit, and measured the energy transmitted through the slit as a function of the groove-slit separation in order to prove the interference model. In this experiment, they also found a strong and fast decaying field in addition to the long-range surface wave mode within the immediate vicinity of the sub-wavelength groove. Lalanne et al. [13] explained this extra field as an additional term to the typical SPP mode related to a branch-cut integral in the solution of the Helmholtz equation for the field on the surface radiated by an electrical dipole, which is induced at the slit by the incident field. They referred this term to as creeping wave [13] and later as quasi-cylindrical wave [14].

In fact, the radiation of an electrical dipole in the vicinity of a metallic surface has been studied for long time. The integral solution of the Helmholtz equation in the frequency domain has been investigated by Sommerfeld in 1909 [15] and other authors in order to study the field radiated by Hertzian dipoles over a lossy Earth, referred to as Zenneck waves and Norton trapped surface waves [16–18]. Numerical evaluation of the Sommerfeld branch-cut integral in [13,14] showed that the creeping wave drops along the interface at a damping rate varying from 1/x^1/2 for distances from the slit less than one wavelength [13] to 1/x^3/2 for distances greater than a wavelength [14]. Asymptotic solutions of the Sommerfeld integral using a modified steepest descent method [18] also found the damping rate of 1/x^1/2 at near-distances [19] and 1/x^3/2 at far-distances [20]. We recently obtained rigorous closed-form solution of the Sommerfeld branch-cut integral for a sub-wavelength metal slit for the field at the metal/dielectric interface. The solution includes a normal SPP mode and an extra term, which we referred to as the transient SPP and is the SPP of wave-vector k_sp modulated by complex-valued envelope described by the exponential integral [21]. Similar result was obtained based on empirical modification of the CDEW model [22]. The transient SPP has a damping rate varying from –ln(x) to 1/x^1/2 and to 1/x, in the distance range from the slit x < λ/20, λ/20 < x < λ and x > λ, respectively.

However, the physical nature of the extra quasi-cylindrical wave is not fully understood. This fundamental problem is recently involved in a model of the surface wave scattering dynamic [23] in the theory of the extraordinary optical transmission of the metallic nanohole arrays [24]. In fact, most studies on the Sommerfeld branch-cut integral evaluate the field only at the surface of the conductor substrate [13–19]. This simplifies greatly the calculation, but fails to produce a complete picture of the scattered field in the dielectric space.

In this paper, we show a rigorous closed-form solution of the Sommerfeld branch-cut integral for the field scattered by a metal slit in the entire half-space over the metal substrate. The transient SPP describes the value of this field at the metal surface, while the field in the 2D space is described as a convolution between the cylindrical harmonic of first order and the transient SPP along the surface. The solution shows clearly that the so-called creeping or quasi-cylindrical field is not an extra surface wave, but is the asymptotic behavior of the scattered field in the 2D space in the region close to the slit, both at the interface and in the fee space. Moreover, the amplitude of the 2D scattered field may be enhanced significantly by surface plasmon resonance. This new theoretical understanding will be helpful in the analysis and design of nano-sized photonics devices.

2. Physical model and solution

Consider a thick planar real metal medium of dielectric function ε_m perforated by a sub-wavelength slit, which is infinitely long in the z-axis. Above the metal surface is a semi-infinite dielectric medium of permittivity ε_d. The problem is two-dimensional (2D) in the x-y space, as shown in Fig. 1 . The incident wave along the –y direction is monochromatic of Transverse Magnetic (TM) polarization with the magnetic field H parallel to the z-axis. Its electric components in the x-y plane induce displacements of the conduction electrons. The discontinuity of the surface at the slit prevents electrical currents from crossing the slit. Hence, the electrical charges accumulate at the two up-edges of the slit, resulting in an oscillating horizontal electric dipole, as have been observed in both experiments [25] and numerical simulations [8].

Fig. 1 Geometry of model for the diffraction of a metal sub-wavelength slit.

Download Full Size | PDF

The dipole acts as a source of radiation, launching the Surface Plasmons Polariton (SPP) mode along the metal/dielectric interface and radiating homogeneous waves in the dielectric half space. Thus, the problem of the diffraction by the metal slit is modeled through solving the 2D field of a horizontal electric dipole on the metallic surface. As a circulation of electric charges results in a magnetic dipole, the electric dipole can also be represented by a circulation of hypothetical magnetic charges, which takes the form of an equivalent magnetic current line along the z axis and replaces the horizontal electric dipole [13,21]. Moreover, for a sub-wavelength slit it is a good approximation to consider the equivalent magnetic current line source to have an infinitesimal spatial extent [13] and to describe it as a 2D delta function of the form: $J_{m}^{*} = J_{m}^{*} δ (x) δ (y) \hat{z}$ , where $J_{m}^{*}$ is the complex amplitude of the magnetic current density. The Maxwell’s equations in such formalism lead to a single differential equation to solve for the magnetic component instead of two coupled equations for two electric components, resulting in a simplification of the problem. For time harmonic fields described by exp(-iωt) and of TM polarization, the Maxwell’s equations in the Fourier space lead to the Helmholtz equation in the spatial frequency domain as [21]

(\partial_{y}^{2} + γ^{2}) \tilde{H} (k_{x}, y) = - i ω ε_{d} J_{m}^{*} δ (y)

where

\tilde{H} (k_{x}, y)

is the Fourier spectrum of the magnetic field

H (x, y)

with respect to x, k_x is the x-component of the wavevector, γ² = ε_dk₀²-k_x² with k₀ = ω/c₀ and c₀ is the speed of light in vacuum. After solving Eq. (1) in the dielectric medium and applying the boundary conditions at the metal/dielectric interface, the inverse Fourier transform leads to a solution of the field in the x-y space in the form of the Sommerfeld integral:

H (x, y) = - \frac{ω J_{m}^{*}}{2 π} \int_{- \infty}^{\infty} \frac{e^{i k_{x} x} e^{i γ y}}{\frac{1}{ε_{d}} \sqrt{ε_{d} k_{0}^{2} - k_{x}^{2}} + \frac{1}{ε_{m}} \sqrt{ε_{m} k_{0}^{2} - k_{x}^{2}}} d k_{x}

which has two branch-cut points in the complex plane at k_x = ε_m^1/2k₀ and k_x = ε_d^1/2k₀, and a pole at the location determined by setting the denominator of the integrant in Eq. (2) to zero, that yields

k_{x} = k_{0} \sqrt{\frac{ε_{d} ε_{m}}{ε_{d} + ε_{m}}} \equiv k_{s p}

where k_sp is the SPP propagation constant. Sommerfeld and others have computed this type of integrals at the interface with

y = 0

using the steepest descent approach in investigation of the Zenneck and Norton surface waves [15–18]. In general, the solutions for the field of an electric dipole close to or on a metal surface take complex forms [18, 26].

2.1 Scattered field in the 2D space

To find a closed-form solution of the field in the entire half space with $y \geq 0$ we consider. the integrant in Eq. (2) as a product of the term ${\tilde{H}}_{b} (k_{x}, y) = \exp (i γ y)$ with $y \neq 0$ and the remaining term, denoted by ${\tilde{H}}_{a} (k)$ , such that the field $H (x, y)$ is a Fourier transform of the product of ${\tilde{H}}_{a} (k)$ and ${\tilde{H}}_{b} (k, y)$ which can solved by the 1D convolution of the inverse Fourier transforms of ${\tilde{H}}_{a} (k)$ and ${\tilde{H}}_{b} (k, y)$ , as $H (x, y) = H_{a} (x) * H_{b} (x, y)$ . In fact, $H_{a} (x)$ corresponds to the field at the interface $y = 0$ , which has been solved by two consequent branch-cut integrals with parameterization of the branch-cut paths in Ref. [21] as

\begin{array}{l} H_{a} (x) = H (x, y = 0) \\ = i ω J_{m}^{*} ε_{d} ε_{m} [\frac{{(ε_{d} ε_{m})}^{\frac{1}{2}}}{ε_{d}^{2} + ε_{m}^{2}} e^{i k_{s p} x} + \frac{1}{2 π (ε_{d} + ε_{m})} e^{i k_{s p} x} E_{1} (i k_{s p} x)] \end{array}

where the first term is the contribution of the pole to the residue in the Cauchy integral and represents the SPP mode at the interface. The second term is the contribution of the two branch points to the residue, referred to as the transient SPP [21,22], which is the SPP, exp(ik_spx) with the complex-valued wavenumber k_sp, modulated by a complex-valued envelope described by the exponential integral of complex argument E₁(ik_spx), being the case

n = 1

of a general set of the functions defined by

E_{n} (z) = \int_{1}^{\infty} \frac{e^{- z t}}{t^{n}} d t

The inverse Fourier transform of ${\tilde{H}}_{b} (k_{x}, y) = \exp (i γ y)$ is performed analytically (see Appendix A) as

H_{b} (x, y) = \frac{i}{2} \frac{k_{d} y}{\sqrt{x^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{x^{2} + y^{2}}) = \frac{i k_{d}}{2} \sin (θ) H_{1}^{(1)} (k_{d} r)

where

k_{d} \equiv \sqrt{ε_{d}} k_{0}

, (r,θ) are the polar coordinates and

H_{1}^{(1)}

is the Hankel function of first kind and first order. Thus, H_b(x,y) represents a cylindrical harmonic of first order [27], which is a 2D solution of the scalar Helmholtz equation with a dipole source in the cylindrical coordinate system. Recall that the zero order cylindrical harmonic is the basic cylindrical wave

H_{0}^{(1)} (k_{d} r)

representing the outward propagating wave emitted by a line source.

Finally, according to Eqs. (4) and (6), the convolution $H (x, y) = H_{a} (x) * H_{b} (x, y)$ for the total magnetic field in the 2D dielectric media above the metallic surface is expressed as a sum of two convolutions: $H (x, y) = H_{s p} (x, y) + H_{s} (x, y)$ (see Appendix B). The first convolution $H_{s p} (x, y)$ is easily solved analytically (see Appendix C):

H_{s p} (x, y) = \frac{2 h_{s p}}{i k_{d}} e^{i k_{s p} x - i \sqrt{k_{d}^{2} - k_{s p}^{2}} y}

where the constant h_sp depends on

J_{m}^{*}

, k₀, ε_d and ε_m as

h_{s p}^{} = \frac{- ω k_{d} J_{m}^{*} {(ε_{d} ε_{m})}^{3 / 2}}{2 (ε_{d}^{2} + ε_{m}^{2})}

$H_{s p} (x, y)$ clearly represents a conventional SPP mode propagating along the interface with the wave-number k_sp and decaying exponentially away from $y = 0$ , with its evanescence in the direction normal to the interface as exp( $- i \sqrt{k_{d}^{2} - k_{s p}^{2}} y$ ).

The term of interest is $H_{s} (x, y)$ and can be expressed as (see Appendix B)

H_{s} (x, y) \equiv h_{s} \int_{- \infty}^{\infty} \frac{y}{\sqrt{{(x - x')}^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{{(x - x')}^{2} + y^{2}}) e^{i k_{s p} x'} E_{1} (i k_{s p} x') d x' ​

where the constant h_s depends on

J_{m}^{*}

, k₀, ε_d and ε_m as

h_{s} = \frac{- ω k_{d} J_{m}^{*} ε_{d} ε_{m}}{4 π (ε_{d} + ε_{m})}

Equation (9) shows clearly that the field scattered into the 2D space above the metallic surface is a 1D convolution along the x-axis of the cylindrical harmonic of first order, $\sin (θ) H_{1}^{(1)} (k_{d} r)$ with the transient SPP [21,22], $\exp (i k_{s p} x) E_{1} (i k_{s p} x)$ . As a comparison, recall that in the case of a perfect conductor slit, which does not support the SPP mode, the diffracted field would be a summation of cylindrical wavelets, or the cylindrical harmonics of zero order, emitted from hypothetical sources inside the slit, according to the Huygens principle. As the width of the slit becomes negligible with respect to the wavelength, the scattered field takes the form of a cylindrical harmonic of order zero, or a typical cylindrical wave, as shown by the simulation in Ref. [28]. However, in the case of a real metal slit, according to Eq. (9) the scattered field is a continuous sum of cylindrical harmonics of first order, emitted from a set of hypothetic dipolar sources with their amplitude and phase following a spatial distribution as the transient SPP function, which extends from the inside to the outside of the slit along the metal/dielectric interface.

2.2 Field at the interface

The cylindrical harmonic of first order, $\sin (θ) H_{1}^{(1)} (k_{d} r)$ , is equal to zero at the interface y = 0 and $θ = 0$ , except at the origin $x = 0$ where it is singular. In fact, at the interface $y = 0$ , ${\tilde{H}}_{b} (k_{x}, 0)$ = exp(iγy)|_{y = 0} = 1, so that the Fourier transform H_b(x,0) is a delta function and its convolution with the transient SPP gives back the transient SPP. This is consistent with previous rigorous solution of the field at the interface [21].

The exponential integral of complex argument, k_sp, is complex valued. Its phase, amplitude, real and imaginary parts are shown in Fig. 2 . The latters can be calculated respectively by the cosine and sine integral. Our calculation shows that the exponential integral of complex argument representing the envelope of the transient SPP has a negative phase, which cancels the positive phase of the SPP propagation, so that the phase of the transient SPP tends to a negative constant of -π/2 within the transient distance range 0 < x < λ from the origin.

Fig. 2 Exponential integral function E₁(ik_spx): (a): Real, imaginary parts and amplitude;(b): Phase. Constitutive parameters have been taken for an air-silver interface at λ=1 μm.

Download Full Size | PDF

In fact, as shown in Ref. [21], the transient SPP is given by

e^{i k_{s p} x} E_{1} (i k_{s p} x) = \int_{0}^{\infty} \frac{e^{i k x}}{k - k_{s p}} d k

which is the inverse Fourier transform of a simple pole centered at k_sp, a continuous sum of all the homogeneous and evanescent plane waves propagating along the x-direction with the complex amplitude 1/(k-k_sp). In this summation the evanescent components are dominant as can be seen in the spectrum of the transient SPP shown in Fig. 3 , where the evanescent components with the frequencies higher than k₀ have greater amplitudes than that have the propagating components with the frequencies lower than k₀, and the spectrum is asymmetric with respect to k/k₀ = 1, so that the scattered field at the surface is damped at a greater rate than the typically SPP mode with a propagation constant k_sp greater than k₀.

Fig. 3 Spectrum of the transient SPP for an air-silver interface at λ=1 μm with the pole at k_sp slightly higher than k₀.

Download Full Size | PDF

The transient SPP is first launched from the slit as a typical cylindrical wave with its amplitude dropping as characterized by –ln(x) and then 1/x^1/2, corresponding to the field naturally radiated by a punctual source. Then, the transient SPP along the surface suffers from an additional drop in amplitude due to the loss of energy through the excitation of the long-range SPP mode.

2.3 Asymptotic behavior near the source

For a small argument 0 < k_dr < $\sqrt{2}$ , the Hankel function is approximated as

H_{1}^{(1)} (k_{d} r) = \frac{2}{i π k_{d} r}

so that

\sin (θ) H_{1}^{(1)} (k_{d} r)

tends to infinity at the origin

r = 0

. Within the distance of x’<λ the series representation of the exponential integral E₁(ik_spx’) can be used, as for an argument z we have

E_{1} (z) = - γ - \ln (z) + \sum_{n = 1}^{\infty} \frac{{(- 1)}^{n + 1} z^{n}}{n n!}

where γ is the Euler constant. On substituting the asymptotic form of

H_{1}^{(1)} (k_{d} r)

and the series representation of E₁(ik_spx), Eqs. (12) and (13), into Eq. (9) we obtain the 2D scattered filed in the region close to the slit as (see Appendix D)

H_{s} (x, y) \approx \frac{2 h_{s}}{i k_{d}} e^{i k_{s p} (x - i y)} E_{1} (i k_{s p} (x - i y))

which includes a new exponential integral of complex argument E₁(ik_sp(x-iy)). When both x and y, or the distance to the origin in the 2D space, tend to zero the amplitude of H_s(x,y) tends to infinity. Figure 4 shows the scattered field in the near-slit region obtained by numerical plot of Eq. (14). In fact, the amplitude drops with the distance from the origin first as a cylindrical wave and then with an increasing damping rate, which is slightly higher in the y-direction, as shown in Fig. 4b, so that in the region close to the slit the field behaves as a flattened cylindrical wave, as shown in Fig. 4c. Such behavior explains the previous experimental results on the extra field [11,12], which are in fact not due to an extra surface wave, but are simply related to the asymptotic behavior of the 2D scattered field near the slit, as expressed in Eqs. (9) and (14). Hence, the rapid drop of the scattered field should be observed not only at the metallic surface, but also in the 2D dielectric region near the slit. Also, the field damping rate along the surface obtained from our theoretical results is in agreement with previous works based on numerical solution [13,14], asymptotic steepest descent solution [19] and rigorous solution [21].

Fig. 4 Field scattered by a metal nano-slit in the near-slit region at λ=1 μm. (a) and (b): Amplitude profile along respectively the x and y axis (solid lines) compared to a typical cylindrical wave (dashed line). (c): 2D plot of amplitude in arbitrary units and (d): the phase in radians. The black semicircle is the upper limit of the approximation (x² + y²)^1/2<< λ/(2^1/2π). For all figures, the complex permittivity of silver was calculated using the Lorentz-Drude dispersive model.

Download Full Size | PDF

2.4 Surface plasmon resonance

The scattered field described in Eq. (9) can have strong resonant enhancement by several orders of magnitude due to the coefficient h_s, which depends on the permittivities of the metal and dielectric, and the incident wavelength. If Re{ε_m}≈-ε_d and Im{ε_m} is small, then from Eq. (10), the amplitude of h_s is strongly enhanced. Figure 5 shows the amplitude of h_s as a function of the wavelength for Ag and Au slit, with the index of refraction in the dielectric medium over the metal slit n = 1, 1.5, 2, 2.5 and 3. The dielectric function of the metal is computed by Lorentz-Drude dispersive model. Strong enhancements, which are more than 3000 times for Ag and 1000 times for Au are clearly shown at the resonant wavelengths for both Ag and Au slits in the visible spectrum. In the case of Ag, the strongest enhancement is observed for a dielectric refractive index of n = 3 at λ = 528 nm, corresponding to a complex permittivity of ε_m = −8.9594 + 0.7958i. In the case of Au, the strongest enhancement has less magnitude than for Ag, and is observed for a dielectric refractive index of n = 3 at λ = 623 nm corresponding to a complex permittivity of ε_m = −9.2504 + 1.9880i. Hence, the enhancement frequencies and magnitudes clearly depend on the metal considered. Furthermore, for a given metallic surface, the resonant wavelength is red shifted and the resonance is stronger, with the increase of the dielectric refractive index.

Fig. 5 Amplitude of the coefficient h_s of the scattered field in Eq. (9) in arbitrary units as a function of the incident wavelength for (a): Ag; and (b): Au. The index of refraction in the dielectric medium over the metal slit is n = 1 (black), 1.5 (blue), 2 (red), 2.5 (green), 3 (orange).

Download Full Size | PDF

2.5 Field far from the source

In the region where the distance to the slit is much larger than one wavelength, the width of the transient SPP peak at $x = 0$ , as shown in Fig. 2a, is negligible. Hence, exp(ik_spx)E₁(ik_spx) may be approximated by a delta function and the scattered field expressed in Eq. (9) is simply the first order cylindrical harmonic with the propagation constant k_d:

H_{s} (x, y) \approx - h_{s} \frac{2 π}{k_{s p}} \frac{y}{\sqrt{x^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{x^{2} + y^{2}})

which differs from the basic cylindrical wave by the presence of the factor sin(θ), corresponding to a diffraction shadow for small θ. Along the y-axis, the field drops as 1/y^1/2 as the typical behavior of the Hankel function, but not as an evanescent surface wave. In the diffraction shadow region close to the interface where x >> y, the field takes an approximate form:

H_{s} (x, y) \approx - h_{s} \frac{2 π}{k_{s p}} \frac{y}{x} H_{1}^{(1)} (k_{d} x)

As the Hankel function drops in amplitude as 1/x^1/2, $H_{s} (x, y)$ will drop as 1/x^3/2. This result is in agreement with the work of Nikitin et al. [20] which evaluated the field in the diffraction shadow using the asymptotic steepest descent method. Similar result has also been obtained by Lalanne et al. [13] with numerical integration methods. Note that when the substrate is a perfect conductor and no SPP mode is present in the solution, the scattered field is a cylindrical wave with a damping rate of 1/x^1/2 away from the slit. When the substrate is real metal, the scattered field behaves as a flattened cylindrical wave only for the distances x<λ/10 as shown in Fig. 4a, and is damped more rapidly as the wave propagates along the interface. Clearly, this behavior is different from the perfect conductor case and can therefore be linked to the real metal and to the pumping of energy to the SPP mode. This result is in agreement with the numerical evaluation by Lalanne et al. [13,14] and the asymptotic steepest descent solution by Ung et al. [19]

3. Numerical solution

To support the analytical result, a Finite Difference Time Domain (FDTD) simulation of the field radiated from a punctual magnetic line source located on a metallo-dielectric interface of a perfect conductor (Fig. 6a ) and silver (Fig. 6b) has been performed. Snapshots of the real part of the magnetic intensity field are shown. The mesh size was 0.01 μm and time step was 0.02238 fs, at the wavelength λ = 500 nm. Results for the perfect conductor, which does not support the SPP mode, show no surface wave and the scattered field in the half space above the surface takes the form of the zero order cylindrical harmonic. As discussed in Section 2.1, this result is expected because the zero order cylindrical harmonic corresponds to the 2D Green’s function of Helmholtz equation. The scattered field behaves as the Hankel function $H_{0}^{(1)} (k_{d} r)$ which has a damping rate along radial axis as a logarithmic drop in the region close to the source, and changes to $1 / \sqrt{r}$ for distance greater than a wavelength from the source. In the case of a silver slit, the SPP mode is clearly visible along the interface as a consequence of the real metal permittivity allowing the metal surface to support SPP modes. However, the scattered field in dielectric region takes the form of the first order cylindrical harmonic, $\sin (θ) H_{1}^{(1)} (k_{d} r)$ , as predicted by Eq. (6). The diffraction shadow is clearly present.

Fig. 6 FDTD simulations (Optiwave^TM) for real part of the magnetic intensity at λ = 500 nm from a punctual magnetic line source located on a metallo-dielectric interface with (a): Perfect electric conductor and (b): Silver.

Download Full Size | PDF

In the direction normal to the surface along the y-axis and sinθ = 1 according to Eq. (6), the field behaves as $H_{1}^{(1)} (k_{d} y)$ , which drops as 1/y^1/2 for y larger than a wavelength, in agreement with our results presented in Section 2.5, and similarly to that in perfect conductor case, which drops as $H_{0}^{(1)} (k_{d} y)$ . In fact, the difference between $H_{0}^{(1)} (k_{d} y)$ and $H_{1}^{(1)} (k_{d} y)$ is when the argument tends to zero where the damping rates are –ln(y) and 1/y, respectively. In the diffraction shadow, the SPP mode is dominant along the interface. The mismatch between the wavelength of the scattered field in the dielectric, k₀, and that of the SPP, k_sp, can be seen clearly from the distortion of the wave fronts, again in agreement with our rigorous closed solution.

4. Conclusion

We have solved the Sommerfeld branch-cut integral in the geometry of 1D slit and found a rigorous closed-form solution in the entire half-space over the slit for the optical diffraction of a metallic nano-slit at the first-surface. Apart from the launched SPP mode, the scattered field propagating in the 2D space is a 1D convolution of the first order cylindrical harmonic with the transient SPP, which is a fast dropping surface component dominant within one wavelength from the slit and is described by the exponential integral. Our results show that the experimentally discovered extra field described as the creeping wave or quasi-cylindrical wave is in fact not an extra evanescent surface wave, but the asymptotic behavior of the 2D space solution at the proximity of the slit. Our rigorous model is consistent with numerical FDTD simulations, as well as previous results based on asymptotic steepest descent solutions. Furthermore, our solution predicts a strong resonant enhancement of the scattered field at the proximity of the slit, depending on the materials and wavelength. This new elementary physical knowledge on the light diffraction by a sub-wavelength metallic nano-slit may be used in the future study on the physical property of the similar quasi-cylindrical waves, which have been newly introduced by the numerical calculation on a 1D array of holes [23]. Our rigorous closed-form solution is limited to the 1D slit geometry. The next challenge would be to analyze the diffraction of the 2D array of nanoholes of non-rectangular forms.

Appendix A

The inverse Fourier transform Eq. (10) is expressed as

H_{b} (x, y) = \frac{1}{2 π} \int_{- \infty}^{\infty} e^{i \sqrt{ε_{d} k_{0}^{2} - k^{2}} y} e^{i k x} d k

Using Euler’s formula, we rewrite Eq. (A1) as

H_{b} (x, y) = \frac{1}{2 π} \int_{- \infty}^{\infty} e^{i \sqrt{ε_{d} k_{0}^{2} - k^{2}} y} (\cos (k x) + i \sin (k x)) d k

The first exponential function is even with respect to k, so that

H_{b} (x, y) = \frac{1}{π} \int_{0}^{\infty} \cos (k x) e^{- \sqrt{k^{2} - ε_{d} k_{0}^{2}} y} d k

In this form, the following identity (see [29], Section 3.914) can be used to solve the integral:

\int_{0}^{\infty} \cos (b x) e^{- β \sqrt{γ^{2} + x^{2}}} d x = \frac{β γ}{\sqrt{β^{2} + b^{2}}} K_{1} (γ \sqrt{β^{2} + b^{2}})

where K₁ is the modified Bessel function of the second kind. Hence, Eq. (A3) becomes

H_{b} (x, y) = \frac{1}{π} \frac{y (\pm i k_{0} \sqrt{ε_{d}})}{\sqrt{x^{2} + y^{2}}} K_{1} (\pm i \sqrt{ε_{d}} k_{0} \sqrt{x^{2} + y^{2}})

The K₁ function with an imaginary argument can be expressed in term of the Hankel function:

H_{b} (x, y) = \frac{i}{2} \frac{k_{d} y}{\sqrt{x^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{x^{2} + y^{2}})

where we replaced

\sqrt{ε_{d}}

k₀ by k_d. The negative sign has been chosen in order to retrieve the Hankel function of the first kind representing waves out-coming from the source. We then obtain the result presented in Eq. (9).

Appendix B

Expressing Eq. (1) is separated as a Fourier transform of a product of two functions:

H (x, y) = \frac{1}{2 π} \int_{- \infty}^{\infty} {\tilde{H}}_{a} (k) {\tilde{H}}_{b} (k, y) e^{i k x} d k

with

{\tilde{H}}_{a} (k) = - ω J_{m}^{*} \frac{1}{\frac{1}{ε_{d}} \sqrt{ε_{d} k_{0}^{2} - k^{2}} + \frac{1}{ε_{m}} \sqrt{ε_{m} k_{0}^{2} - k^{2}}}

{\tilde{H}}_{b} (k, y) = e^{i \sqrt{ε_{d} k_{0}^{2} - k^{2}} y}

The total magnetic field can thus be obtained using the convolution theorem:

H (x, y) = H_{a} (x) * H_{b} (x, y)

Since $H_{a} (x)$ and $H_{b} (x, y)$ are given by respectively Eqs. (4) and (7), the convolution can be written explicitly:

\begin{array}{l} H (x, y) = h_{s p}^{} \int_{- \infty}^{\infty} \frac{y}{\sqrt{{(x - x')}^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{{(x - x')}^{2} + y^{2}}) e^{i k_{s p} x'} d x' \\ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ + h_{s} \int_{- \infty}^{\infty} \frac{y}{\sqrt{{(x - x')}^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{{(x - x')}^{2} + y^{2}}) e^{i k_{s p} x'} E_{1} (i k_{s p} x') d x' \end{array}

where the convolution has been distributed on the sum of H_a and the following constants have been defined:

h_{s p}^{} = \frac{- ω k_{d} J_{m}^{*} {(ε_{d} ε_{m})}^{\frac{3}{2}}}{2 (ε_{d}^{2} + ε_{m}^{2})}

h_{s} = \frac{- ω k_{d} J_{m}^{*} ε_{d} ε_{m}}{4 π (ε_{d} + ε_{m})}

The following form is also valid by virtue of the commutability of the convolution:

\begin{array}{l} H (x, y) = h_{s p}^{} \int_{- \infty}^{\infty} \frac{y}{\sqrt{{(x')}^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{{(x')}^{2} + y^{2}}) e^{i k_{s p} (x - x')} d x' \\ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ + h_{s} \int_{- \infty}^{\infty} \frac{y}{\sqrt{{(x')}^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{{(x')}^{2} + y^{2}}) e^{i k_{s p} (x - x')} E_{1} (i k_{s p} (x - x')) d x' \end{array}

The first term is labeled $H_{s p} (x, y)$ for the SPP mode and the second $H_{s} (x, y)$ for the scattered field such that

H (x, y) = H_{s p} (x, y) + H_{s} (x, y)

Appendix C

The first term of the convolution given by Eq. (B5) can be solved analytically:

H_{s p} (x, y) = h_{s p} \int_{- \infty}^{\infty} \frac{y}{\sqrt{{(x')}^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{{(x')}^{2} + y^{2}}) e^{i k_{s p} (x - x')} d x' ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​

We use Euler formula to rewrite Eq. (C1):

\begin{array}{l} H_{s p} (x, y) = h_{s p} y e^{i k_{s p} x} \\ \int_{- \infty}^{\infty} \frac{1}{\sqrt{x'^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{x'^{2} + y^{2}}) (\cos (k_{s p} x') - i \sin (k_{s p} x')) d x' ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ \end{array}

All the functions of x’ in the integrant are even, except sin(k_spx’), which is odd. We also separate the Hankel function in its Bessel components:

\begin{array}{l} H_{s p} (x, y) = 2 h_{s p} y e^{i k_{s p} x} \\ [\int_{0}^{\infty} \frac{J_{1} (k_{d} \sqrt{x'^{2} + y^{2}})}{\sqrt{x'^{2} + y^{2}}} \cos (k_{s p} x') d x' + i \int_{0}^{\infty} \frac{Y_{1} (k_{d} \sqrt{x'^{2} + y^{2}})}{\sqrt{x'^{2} + y^{2}}} \cos (k_{s p} x') d x'] ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ \end{array}

The following identities (see [23], Section 6.726) can be used to solve the two integrals:

\begin{array}{l} \int_{0}^{\infty} {(x^{2} + b^{2})}^{\frac{- ν}{2}} ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ J_{ν} (a \sqrt{x^{2} + b^{2}}) ​ ​ ​ \cos (c x) d x \\ = \sqrt{\frac{π}{2}} a^{- ν} b^{- ν + \frac{1}{2}} {(a^{2} - c^{2})}^{\frac{ν}{2} - \frac{1}{4}} J_{ν - \frac{1}{2}} (b \sqrt{a^{2} - c^{2}}) \end{array}

\begin{array}{l} \int_{0}^{\infty} {(x^{2} + b^{2})}^{\frac{- ν}{2}} ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ Y_{ν} (a \sqrt{x^{2} + b^{2}}) ​ ​ ​ \cos (c x) d x \\ = \sqrt{\frac{π}{2}} a^{- ν} b^{- ν + \frac{1}{2}} {(a^{2} - c^{2})}^{\frac{ν}{2} - \frac{1}{4}} Y_{ν - \frac{1}{2}} (b \sqrt{a^{2} - c^{2}}) \end{array}

Using Eqs. (C4) and (C5), Eq. (C3) becomes

\begin{array}{l} H_{s p} (x, y) = 2 h_{s p} y e^{i k_{s p} x} \\ [\sqrt{\frac{π}{2}} \frac{{(k_{d}^{2} - k_{s p}^{2})}^{\frac{1}{4}}}{k_{d} \sqrt{y}} J_{\frac{1}{2}} (\pm y \sqrt{k_{d}^{2} - k_{s p}^{2}}) + i \sqrt{\frac{π}{2}} \frac{{(k_{d}^{2} - k_{s p}^{2})}^{\frac{1}{4}}}{k_{d} \sqrt{y}} Y_{\frac{1}{2}} (\pm y \sqrt{k_{d}^{2} - k_{s p}^{2}})] ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ \end{array}

Simplifying the equation and combining the two Bessel Functions in a Hankel function gives

H_{s p} (x, y) = h_{s p} e^{i k_{s p} x} \sqrt{2 π} \frac{{(k_{d}^{2} - k_{s p}^{2})}^{\frac{1}{4}}}{k_{d}} \sqrt{y} H_{\frac{1}{2}}^{(1)} (\pm y \sqrt{k_{d}^{2} - k_{s}^{2}})

The Hankel function can be expressed in an analytical form with the following identity (see [29], Section 8.469):

H_{\frac{1}{2}}^{(1)} (z) = \sqrt{\frac{2}{π z}} ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ \frac{e^{i z}}{i}

Equation (C7) becomes

H_{s p} (x, y) = h_{s p} \frac{- 2 i}{k_{d}} e^{i k_{s p} x} e^{\pm i y \sqrt{k_{d}^{2} - k_{s p}^{2}}} ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​

Selecting the minus sign to satisfy the radiation condition gives the result presented in Eq. (C4).

Appendix D

To obtain a solution for the scattered field in the near-field region, which is when the distance to the slit is smaller than one wavelength, we take the convolution integral for H_s(x,y) in the form shown in Eq. (B8):

H (x, y) = ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ ​ h_{s} \int_{- \infty}^{\infty} \frac{y}{\sqrt{{(x')}^{2} + y^{2}}} H_{1}^{(1)} (k_{d} \sqrt{{(x')}^{2} + y^{2}}) e^{i k_{s p} (x - x')} E_{1} (i k_{s p} (x - x')) d x'

For a small arguments z, the Hankel function can be approximated as

H_{1}^{(1)} (z) \approx \frac{2}{i π z}

which is valid for 0<z<<2^1/2 and z = k_d[(x-x’)² + y²]^1/2where x’ is smaller than a wavelength for E₁(ik_spx)≠0. The exponential integral can be expressed as

E_{1} (z) = - γ - \ln (z) + \sum_{n = 1}^{\infty} \frac{{(- 1)}^{n + 1} z^{n}}{n n!}

where γ is the Euler constant. Hence, in the region with the distances to the slit smaller than the wavelength, Eq. (D1) can be approximated as

H_{s} (x, y) ≃ \frac{2 h_{s} y}{i π k_{d}} e^{i k_{s p} x} \int_{- \infty}^{\infty} \frac{1}{{(x')}^{2} + y^{2}} ​ e^{- i k_{s p} x'} [\begin{array}{l} - γ - \ln (i k_{s p} (x - x')) \\ + \sum_{n = 1}^{\infty} \frac{{(- 1)}^{n + 1} {(i k_{s p} (x - x'))}^{n}}{n n!} \end{array}] d x'

Equation (D4) can be solved using complex contour integration. The poles of the integrant are at x’=±iy. Closing the contour with a half-circular path in the upper half space encloses the positive pole. The integral is then computed with Cauchy’s theorem trough evaluation of the residue around the pole x’=iy, resulting in the following solution satisfying the radiation condition:

H_{s} (x, y) \approx \frac{2 h_{s}}{i k_{d}} e^{i k_{s p} (x - i y)} [- γ - \ln (i k_{s p} (x - i y)) + \sum_{n = 1}^{\infty} \frac{{(- 1)}^{n + 1} {(i k_{s p} (x - i y))}^{n}}{n n!}]

which according to Eq. (D3) can be expressed in terms of a new exponential integral function:

H_{s} (x, y) \approx \frac{2 h_{s}}{i k_{d}} e^{i k_{s p} (x - i y)} E_{1} (i k_{s p} (x - i y))

References and links

1. W. L. Barnes, A. Dereux, and T. W. Ebbesen, “Surface plasmon subwavelength optics,” Nature 424(6950), 824–830 (2003). [CrossRef] [PubMed]

2. S. A. Maier, Plasmonics: Fundamentals and Applications (Springer, Berlin, 2007).

3. T. W. Ebbesen, H. J. Lezec, H. F. Ghaemi, T. Thio, and P. A. Wolff, “Extraordinary optical transmission through sub-wavelength hole arrays,” Nature 391(6668), 667–669 (1998). [CrossRef]

4. L. Matrin-Moreno, F. J. Gracia-Vidal, H. J. Lezec, K. M. Pellerin, T. Thio, J. B. Pendry, and T. W. Ebbesen, “Theory of extraordinary optical transmission through sub-wavelength hole arrays,” Phys. Rev. Lett. 86, 1112–1117 (2001).

5. J. B. Pendry, L. Martín-Moreno, and F. J. Garcia-Vidal, “Mimicking surface plasmons with structured surfaces,” Science 305(5685), 847–848 (2004). [CrossRef] [PubMed]

6. Q. Cao and P. Lalanne, “Negative role of surface plasmons in the transmission of metallic gratings with very narrow slits,” Phys. Rev. Lett. 88(5), 057403–057407 (2002). [CrossRef] [PubMed]

7. H. J. Lezec and T. Thio, “Diffracted evanescent wave model for enhanced and suppressed optical transmission through subwavelength hole arrays,” Opt. Express 12(16), 3629–3651 (2004). [CrossRef] [PubMed]

8. B. Ung and Y. Sheng, “Interference of surface waves in a metallic nanoslit,” Opt. Express 15(3), 1182–1190 (2007). [CrossRef] [PubMed]

9. M. W. Kowarz, “Homogeneous and evanescent contributions in scalar near-field diffraction,” Appl. Opt. 34(17), 3055–3063 (1995). [CrossRef] [PubMed]

10. H. F. Schouten, N. Kuzmin, G. Dubois, T. D. Visser, G. Gbur, P. F. A. Alkemade, H. Blok, G. W. Hooft, D. Lenstra, and E. R. Eliel, “Plasmon-assisted two-slit transmission: Young’s experiment revisited,” Phys. Rev. Lett. 94(5), 053901 (2005). [CrossRef] [PubMed]

11. G. Gay, O. Alloschery, B. Viaris De Lesegno, C. O'Dwyer, J. Weiner, and H. J. Lezec, “The optical response of nanostructured surfaces and the composite diffracted evanescent wave model,” Nat. Phys. 2(4), 262–267 (2006). [CrossRef]

12. G. Gay, O. Alloschery, J. Weiner, H. J. Lezec, C. O’Dwyer, M. Sukharev, and T. Seideman, “Surface quality and surface waves on subwavelength-structured silver films,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 75(1), 016612 (2007). [CrossRef] [PubMed]

13. P. Lalanne and J. P. Hugonin, “Interaction between optical nano-objects at metallo-dielectric interfaces,” Nat. Phys. 2(8), 551–556 (2006). [CrossRef]

14. P. Lalanne, J. P. Hugonin, H. T. Liu, and B. Wang, “A microscopic view of the electromagnetic properties of sub-λ metallic surfaces,” Surf. Sci. Rep. 64(10), 453–469 (2009). [CrossRef]

15. A. N. Sommerfeld, “Propagation of waves in wireless telegraphy,” Ann. Phys. (Leipzig) 28, 665–737 (1909).

16. J. Zenneck, “Propagation of plane EM waves along a plane conducting surface,” Ann. Phys. (Leipzig) 23, 846–866 (1907).

17. K. A. Norton, “The propagation of radio waves over the surface of the earth and in the upper atmosphere,” Proc. IRE 24(10), 1367–1387 (1936). [CrossRef]

18. R. E. Collin, “Hertzian dipole radiating over a lossy earth or sea: Some early and late 20^th-century controversies,” IEEE Antennas Propagat. Mag. 46(2), 64–79 (2004). [CrossRef]

19. B. Ung and Y. Sheng, “Optical surface waves over metallo-dielectric nanostructures: Sommerfeld integrals revisited,” Opt. Express 16(12), 9073–9086 (2008). [CrossRef] [PubMed]

20. A. Y. Nikitin, S. G. Rodrigo, F. J. Garcia-Vidal, and L. Martín-Moreno, “In the diffraction shadow: Norton waves versus surface plasmon polaritons in the optical region,” New J. Phys. 11(12), 123020 (2009). [CrossRef]

21. Y. Gravel and Y. Sheng, “Rigorous solution for the transient Surface Plasmon Polariton launched by subwavelength slit scattering,” Opt. Express 16(26), 21903–21913 (2008). [CrossRef] [PubMed]

22. G. Lévêque, O. J. F. Martin, and J. Weiner, “Transient behaviour of surface plasmon polaritons scattered at a sub-wavelength groove,” Phys. Rev. B 76(15), 155418 (2007). [CrossRef]

23. H. Liu and P. Lalanne, “Microscopic theory of the extraordinary optical transmission,” Nature 452(7188), 728–731 (2008). [CrossRef] [PubMed]

24. T. W. Ebbesen, H. J. Lezec, H. F. Ghaemi, T. Thio, and P. A. Wolff, “Extraordinary optical transmission through sub-wavelength hole arrays,” Nature 391(6668), 667–669 (1998). [CrossRef]

25. L. Yin, V. K. Vlasko-Vlasov, J. Pearson, J. M. Hiller, J. Hua, U. Welp, D. E. Brown, and C. W. Kimball, “Subwavelength focusing and guiding of surface plasmons,” Nano Lett. 5(7), 1399–1402 (2005). [CrossRef] [PubMed]

26. R. W. P. King, M. Owens, and T. T. Wu, Lateral Electromagnetic Waves: Theory and Applications to Communications, Geophysical Exploration and Remote Sensing (Springer-Verlag, New York, 1992).

27. R. F. Harrington, Time-Harmonic Electromagnetic Fields (McGraw-Hill, 1961).

28. L. Chen, J. T. Robinson, and M. Lipson, “Role of radiation and surface plasmon polaritons in the optical interactions between a nano-slit and a nano-groove on a metal surface,” Opt. Express 14(26), 12629–12636 (2006). [CrossRef] [PubMed]

29. I. S. Gradshteyn and I. m. Ryzhik, Table of Integrals, Series and Products, 7th ed. (Academic, 2007).

Rigorous solution for optical diffraction of a sub-wavelength real-metal slit

Abstract

1. Introduction

2. Physical model and solution

2.1 Scattered field in the 2D space

2.2 Field at the interface

2.3 Asymptotic behavior near the source

2.4 Surface plasmon resonance

2.5 Field far from the source

3. Numerical solution

4. Conclusion

Appendix A

Appendix B

Appendix C

Appendix D

References and links

Cited By

Figures (6)

Equations (46)

Optics Express