Surface-wave propagation through a metal gap with ... - Springer Link

J. Appl. Math. & Computing Vol. 25(2007), No. 1 - 2, pp. 315 - 327 Website: http://jamc.net

SURFACE-WAVE PROPAGATION THROUGH A METAL GAP WITH THE DIELECTRIC CORE SUBDIVIDED INTO MULTIPLE THIN FILMS

Jinsik Mok∗ and Hyoung-In Lee

Abstract. Mathematical aspects of the electromagnetic surface-wave propaga-

tion are examined for the dielectric core consisting of multiple sub-layers, which are embedded in the gap between the two bounding cladding metals. For this purpose, the linear problem with a partial differential wave equation is formulated into a nonlinear eigenvalue problem. The resulting eigenvalue is found to exist only for a certain combination of the material densities and the number of the multiple sub-layers. The implications of several limiting cases are discussed in terms of electromagnetic characteristics. AMS Mathematics Subject Classification : : 35L05, 35Q60, 35P20, 78A48, Key words and phrases : Wave equation, equations of electromagnetic theory and optics, asymptotic distribution of eigenvalues and eigenfunctions for PDE, composite media, random media

1. Introduction Though the partial differential equations (PDE) describing the electromagneticwave propagation are in general linear, they defy analytical solutions except in a very few cases. In multi-dimensional geometry, ordinary differential equations (ODE) result when the separation of variables is applied to the PDEs [3]. The penalty of formulating ODEs from PDEs is the introduction of additional separation parameters, one of which serves as an eigenvalue. The resulting ODEs are solved by applying the appropriate boundary conditions, leading to a nonlinear eigenvalue problem. This eigenvalue depends in turn strongly on the characteristic length and the material densities, and in the present problem, on the spatial periodicity of the sub-layers. The existence of solutions is of primary interest in the present paper, but the uniqueness remains to be answered in the future. Received February 5, 2007. ∗ Corresponding author. c 2007 Korean Society for Computational & Applied Mathematics.

315

316

Jinsik Mok and Hyoung-In Lee

The sinusoidal volume waves have been of traditional interest to both theory and application of the wave equation. In comparison, the hyperbolic surface waves currently receive increasing attentions, since they could divert incoming waves into directions favorable to many applications such as in flat-panel displays [2]. The surface waves are established in optical systems, if a pair of positive and negative dielectric materials exists in contact with each other. Two simple configurations are possible to meet this condition: one is the metal film situated between two dielectric materials, and the other is the metal gap where a dielectric core is located between two surrounding cladding bulk metals as shown in Fig. 1 [4, 5]. The alteration of wave-propagation characteristics is oftentimes achieved by the presence of the metal-dielectric interfaces, although the presence of metal may incur energy loss to some degrees. The effect of the dielectric core of periodically stratified structure on the optical performance in the metal gap is of primary interest in this paper. The problem is formulated as an eigenvalue problem, where the propagation constant in the longitudinal direction normal to the depth serves as an eigenvalue. The entire problem domain is first divided into two regions: one embedded core region and the other surrounding cladding region. The inner core region of fixed depth is then subdivided into multiple sub-layers, where a pair consisting of one sublayer of one data and the other adjoining sub-layer of another data is repeatedly stacked over one another. The present structure is in contrast to the traditional structure, where the periodic part of the structure is of infinite or semi-infinite extent [6, 7, 8]. The convergence of the eigenvalue is investigated, if exist, as the number of these sub-layers is increased. The non-trivial convergence property thus found is a well-known optical characteristics from physical viewpoint on one hand, but it spurs further investigations from the mathematical viewpoint on the other hand. 2. Wave equations Figure 1 shows schematically the problem configuration along with the necessary data. The Cartesian coordinate x refers to the longitudinal direction of the main wave propagation or in-plane direction, whereas z indicates the depth direction into the multiple sub-layers. The relative dielectric constants of the material, k , are prescribed data, where the subscript k indicates different parts of the layers. Its square root, (k )1/2 for k > 0 is the refractive index, denoting the slowness of the electromagnetic wave propagation with respect to c, the light speed in vacuum. For simplicity, only a real value is ascribed to k > 0. Decoupled from the full set of the Maxwell equations, the transverse magnetic field f (t, x, z) is then governed by the following scalar wave equation in time t and (x, z) [1]: ∂2f ∂2f ∂2f + 2 = 2 2 (1) 2 ∂x ∂z c ∂t Although the dependence of f on the transverse direction y is missing, the magnetic field itself is described by f only in this direction, hence being called

Surface-wave propagation through a metal gap

317

“ transverse-magnetic (TM)”. An infinite domain {x : −∞ < x < ∞} is considered in the longitudinal direction. The core region in the domain 0 ≤ z ≤ h is occupied by a total of 2J dielectric media as illustrated in Fig. 1, where pairs of two layers of different properties, u > 0 and l > 0, are stacked by one over another. The whole configuration is symmetric with respect to the bottom plane located at z = 0. On top of the multiple dielectric sub-layers for h ≤ z, there is a semi-infinite cladding metal of negative dielectric constant, m < 0. As a whole, the dielectric core of a total of 4J sub-layers is embedded within the two surrounding cladding metals from above and below.

Figure 1. Schematic illustration of the problem domain. The semi-infinite metallic cladding surrounds the dielectric core of depth 2δ. The structure is symmetric with respect to the plane at z = 0. The dielectric core is further subdivided into a total of 2J sub-layers. The linear Eq. (1) allows a separation of variables through the normal-mode √ expression f = F (z)exp[i(ωt − βx)] with i = −1, where ω denotes the timefrequency and β the propagation constant. The depth profile F (z) is in general complex-valued, but it can be confined to being real-valued as long as both ω and β are kept real-valued [2, 9]. The wave equation is then reduced to the following ordinary differential equation: d2 F/dz 2 = (β 2 − k02 )F , where k0 = ω/c is the real wave number in vacuum. A set of dimensionless variables are defined such that Z − Z0 = k0 z, δ = k0 h, and b = β/k0 . The constant Z0 is chosen for convenience of the ensuing mathematical manipulations. With B ≡ b2 , d2 F = (B − )F. dZ 2

(2)

318


Since the dimensionless phase speed is defined by VP = /(cβ) = 1/b, the longitudinal wave gets faster with decreasing b. 3. Transfer matrices A boundary-value problem is now formulated and solved for F (Z). A twolayer system of equal depth ∆J ≡ δ/(2J) > 0, is considered as displayed in the upper right corner of Figure 1. Different data are specified such that = u for 0 ≤ Z ≤ ∆J and = l for −∆J ≤ Z ≤ 0, where the subscripts u and l imply “ upper” and “ lower”, respectively. Figure 1 shows that each pair consisting of the two adjoining layers is numbered by the subscript j, raging from j = 0 for the bottommost pair to j = J for the topmost pair. For the given data k , the decay function (DF) and modified decay function (MDF) are defined as follows: Dk (B; k ) ≡

p B − k ,

√ Uk (B; k ) ≡

B − k k

(3)

In this paper, only a restricted range of B ≥ Bmin ≡ max(u , l ) is treated, thus ensuring hyperbolic depth profiles. Hence, the solutions to Eq. (2) are obtained in the corresponding domains of the j-th pair as follows. Fu,j (Z) = Au,j cosh(Du Z) + Bu,j sinh(Du Z) for 0 ≤ Z ≤ ∆J h i h i Fl,j (Z) = Al,j cosh Dl (Z+∆J ) +Bl,j sinh Dl (Z+∆J ) for −∆J ≤ Z ≤ 0 (4) The hyperbolic depth profiles in Eq. (4) change monotonically away from the local interface at Z = 0, that is measured from the absolute location at z = h(2j − 1)/(2J). The solutions in Eq. (4) are associated with surface waves. In contrast, the sinusoidal profiles occur in the other range B ≤ min(u , l ), being associated with volume waves. In general, the difference |u − l | is taken here to be non-negligible, or has a finite index contrast. Formally expressed, |u − l |/Bmin = O(1), where O(1) implies the magnitude of order unity. In electromagnetic-wave theory, the boundary condition across the interface at Z = 0 dictates that both function and its depth-wise derivative modified by k are continuous such that Fu,j (0) = Fl,j (0) and −1 u dFu,j /dZ(0) = −1 dF /dZ(0)[1]. They are given in the vector-matrix notation as follows. l,j l   bl,J bl,J S C Al Au (5) = Teu←l,J , Teu←l,J ≡  Ul b Ul b  Bu j Bl j Uu Sl,J Uu Cl,J bk,J ≡ cosh(Dk ∆J ), C

Sbk,J ≡ sinh(Dk ∆J )

(6)

The transfer matrix Teu←l,J denotes the influence of the one layer of l on the lower side onto the other layer of u on the upper side [6]. For all values of j, Teu←l,J remains the same. Interchanging the role of the layer of u with that of the layer of l , another transfer matrix Tel←u,J is obtained, which denotes the


319

influence of one layer of u on the lower side onto the other layer of l on the upper side. Next, a three-layer system is considered, where three layers with equal depth ∆J are given different data: u , l , and u . The “ u-l-u unit-cell transfer matrix” Teu←l←u,J ≡ Teu←l,J Tel←u,J relates one layer of u to next close layer of u through Au Au = Teu←l←u,J . (7) Bu j Bu j−1 Component-wise, the composite 2-by-2 matrix Teu←l←u,J ≡ Teul,2 ≡ (Tijul,2 ) with i, j = 0, 1 is expressed as follows:    bu,J bl,J C Sbl,J Sbu,J C Teu←l←u,J =  Ul b Ul b   Uu b Uu b  , Sl,J Cl,J Su,J Cu,J Uu Uu Ul Ul   bu,J + Uu Sbl,J Sbu,J C bl,J Sbu,J + Uu Sbl,J C bu,J bl,J C C   Ul Ul = (8) Ul b b Ul b b  , b b b b Cl,J Su,J + Sl,J Cu,J Cl,J Su,J + Sl,J Su,J Uu Uu The components of Teu←l←u,J are found to satisfy the uni-modularity condition: ul,2 ul,2 ul,2 ul,2 T00 T11 −T01 T10 = 1. Accounting for the influence of the bottommost layer of u on the topmost layer of u , the total composite transfer matrix is therefore equal to (Teu←l←u,J )J−1 , because the total number of unit cells involved is J − 1. Teul,J−1 ≡ (Teu←l←u,J )J−1 ,

Au Bu

j=J

= Teul,J−1

Au Bu

.

(9)

j=l

Component-wise, Teul,J−1 ≡ (Tijul,J−1 ). As with Teu←l←u,J , the components of Teul,J−1 are found to satisfy the uni-modularity condition: ul,J−1 ul,J−1 ul,J−1 ul,J−1 T00 T11 − T01 T10 = 1.

However, ul,J−1 ul,J−1 ul,J−1 ul,J−1 MJ ≡ T00 T11 − T01 T10

is found not to remain unity as B → ∞, when it is evaluated on finite-precision computers. Instead, lim MJ (B) shows a highly oscillatory behavior for large B→∞

J.

It will turn out essential to examine the two limiting behaviors of Teul,J−1 with respect to B: lim Teul,J−1 and lim Teul,J−1 . For convenience, only B→Bmin

B→∞

the case u > l is chosen. Consequently, B ≥ Bmin = max(u , l ) = u . The Bmin = u is now investigated. behavior of Teu←l←u,J as B → √ √ It is found in lim Dl = u − l (finite), this limit that lim Du = u − u = 0, B→Bmin

B→Bmin

320

and


lim

B→Bmin

Dm =

√

u − m . Furthermore,

lim

B→Bmin

bu,J = 1 , C

lim

B→Bmin

Sbu,J =

Du ∆J from Eq. (6) and p by linearization p for small Du ∆J . In consequence, lim (−Um )/Uu = u B − m /(−m B − u ) = ∞ as B → Bmin = u .

B→Bmin

It is found from a more detailed analysis not presented here that an analytic expression of limB→Bmin Teul,J−1 is very complicated, and therefore it may be better just to numerically evaluate the value of Teul,J−1 , say, at B = (1.0001)Bmin . Attention is thereby turned to a little restricted case of the limit as B → Bmin = u . In other words, simultaneously with the limit as B → Bmin = u , the limit of large sub-layer number, i.e., ∆J → 0 as J → ∞ is considered. From Eq. (8), the following double limit is obtained. " # 1 0 u lim lim Teu←l←u,J = Ul D ∆ (10) 1 + (Dl ∆J )2 . l J J→∞ B→Bmin Uu l The point here is to keep the large term Dl ∆J (Ul /Uu ) and small term (Dl ∆J )2 momentarily, and to carry out further manipulation. With Ψ11 ≡ 1 + (u /l )(Dl ∆J )2 , the induction argument given in Appendix A shows that " # 1 0 ul,J−1 √ e δ U l lim lim T = (11) u − l 1 . J→∞ B→Bmin Uu 2 Meanwhile, the parameter φ ≡ u /l is introduced for the analysis of lim Teul,J−1 . In the limit as B → ∞, lim Uu /Ul = φ−1 < 1, and

B→∞

B→∞

√ bk,J = lim Sbk,J = lim 1 exp (Dk ∆J ) = 1 exp (∆J B) lim C B→∞ B→∞ B→∞ 2 2 √ because of lim Dk = B. Therefore, the transfer matrix in Eq. (8) becomes B→∞

asymptotically √ 1 + φ−1 1 lim Teu←l←u,J = exp 2∆J B 1+φ B→∞ 4

1 + φ−1 . 1+φ

(12)

As detailed in Appendix B, the corresponding composite matrix Teul,J−1 behaves in this limit according to √ 1 FJ−1 (φ−1 ) FJ−1 (φ−1 ) ul,J−1 e lim T = J−1 exp (2(J − 1)∆J B) . (13) FJ−1 (φ) FJ−1 (φ) B→∞ 4 4. Eigenvalue problem A closed form of an eigenvalue problem is now derived by applying appropriate boundary conditions. The symmetry condition at z = 0 leads to the


321

symmetric depth profile Fl,1 (Z) = cosh(Dl,J Z) for the bottommost sub-layer of l , when the local coordinate is attached to the bottom of this sub-layer, i.e, at the symmetry plane. Here, the magnitude of unity is arbitrarily prescribed to the constant multiplying the function cosh(Dl,J Z), because of the homogeneity of the eigenvalue problem. The interface conditions with the first layer of u bl,J and Uu Bu,1 = Ul Sbl,J , resulting in the following transfer read then as Au,1 = C matrix.   bl,J C Au (14) =  Ul b . Bu j=1 Sl,J Uu At the other end in the cladding metal region, the solution is taken to be Fm (Z) = Am exp (−Dm Z) when the local coordinate is attached to the metalcore interface at z = h. Because m < 0, it is always true that Dm ≡ (B − m )1/2 > 0, thereby satisfying the far-field radiation condition Fm (Z) → 0 as Z → ∞. Incidentally, Dm > Du and Dm > Dl . The application of the interface conditions with the topmost dielectric sub-layer of u produces bu,J + Bu,J Sbu,J and − Um Am = Uu Au,J Sbu,J + Uu Bu,J C bu,J . Am = Au,J C Here Um ≡ Sm /m , in accordance with Eq. (3). The corresponding transfer matrix denoting the influence of the topmost sub-layer of u on the surrounding metallic cladding of m is then expressed as usual:

Au 1 e , A = Tm←u,J −1 m Bu j=J

Tem←u,J



bu,J C  ≡ Uu b Su,J Um

 Sbu,J Uu b  Cu,J Um

(15)

The combination of Eqs. (9), (14), and (15) gives rise to the total transfer matrix, which expressed the influence of the bottommost sub-layer of l just above the symmetry plane on the surrounding cladding metal of m :     bl,J bu,J C C Sbu,J 1 (16) A =  Uu b Uu b  Teul,J−1  Ul b  −1 m Su,J Cu,J Sl,J Um Um Uu As a result, the two equations in Eq. (16) provides a solvability condition, or more conventionally called, dispersion relation for the single unknown Am : V00 e V ≡ V10

V01 V11

GJ (B) ≡ −

b C ≡ bu,J Su,J

Sbu,J eul,J−1 bu,J T C

bl,J + Ul V11 Sbl,J Um Uu V10 C − bl,J + Ul V01 Sbl,J Uu Uu V00 C

(17)

(18)

Normally, the square of the propagation constant, B = b2 , serves as the eigenvalue of the nonlinear dispersion Eq. (18). The subscript J in GJ emphasizes

322


its strong dependence on the number of multiple sub-layers. When the set of the real material properties (u , l , m ) and the depth δ are explicitly stated, the dispersion relation then reads as GJ (B; u , l , m , δ). For a true eigenvalue B, GJ (B) = 0 exactly. The function GJ (B) is usually called the residual function, because |GJ (B)| takes its minimum value for an approximate eigenvalue B. 5. Existence of a solution Being nonlinear, Eq. (18) defies a general theory for the existence and/or uniqueness of solutions. A direct numerical evaluation of Eq. (18) over a vast parameter space (u , l , m , δ) requires excessive human efforts on the side. Here, some aspects of the existence are investigated. The idea is to examine the sign change of GJ (B) over the entire domain {B : Bmin ≤ B