1. Introduction

Recently, the work on metamaterial [1–4] and “perfect lenses” [5] revived Veselago’s lefthanded materials [6] and triggered intensive discussion on negative refraction. Meanwhile, negative refraction [7] and imaging by negative refraction [8,9] were also demonstrated using 2D photonic crystals (PhCs) by dispersion-engineering. Subsequently, we experimentally achieved subwavelength imaging by a 3D PhC flat lens exhibiting full 3D negative refraction [10]. Although significant experimental progress has been made in observing and employing this novel phenomenon, several basic questions are still unclear — for example, what is the relation between metamaterial and negative-refraction PhCs, and specifically whether negative-refraction PhCs have negative effective permittivity and negative permeability. More recently, Decoopman, et al., proposed to calculate the impedance of a negative-refraction PhC by observing the reflection into its matched medium, and numerically demonstrated that the value of the impedance of a PhC strongly depends on its boundary — especially the way of truncation of the elementary cells located on this boundary [11]. In their work, they assumed that the relations among effective permittivity, permeability, impedance, and refractive index, ${\eta }_{\mathrm{eff}}=\sqrt{\frac{{\mu }_{\mathrm{eff}}}{{\epsilon }_{\mathrm{eff}}}}$ and $n_{\mathrm{eff}}=\sqrt{{\mu }_{\mathrm{eff}}}\sqrt{{\epsilon }_{\mathrm{eff}}}$, are still valid. Their work provides a very novel approach to calculate effective impedance, but several issues are still open to discuss. First, effective permittivity and permeability as bulk parameters for material should not depend on the boundary. If they only describe the surface properties of the PhC, then the effective parameters inside the PhC are still unknown. Second, negative-refraction PhCs work by the unique photonic dispersion, which however does not provide any clue for negative permittivity and negative permeability. In fact, these two parameters in a PhC are known and accurately given as periodic functions, i.e. μ_r(r) = μ_r(r + R _lmn and ε_r = ε_r (r + R), where R _lmn represents a unit lattice translation vector. Moreover, in their work the determination of the effective parameters relies on numerically simulated results and each simulation provides the parameters for only one boundary. To address these issues, we present a method based on averaging eigen-fields to calculate the effective permittivity, permeability, and surface impedance on all possible boundaries of negative-refraction PhCs.

2. Calculation of effective permittivity and permeability

We consider a 2D triangular-lattice PhC formed by patterning air holes in a low loss material with dielectric constant ε_r = 20. The lattice constant is a and the diameter of the air holes is 2r = 0.7a, as shown in Fig. 1(a). Subwavelength imaging by this PhC was reported in our recent work [9]. For transverse-magnetic-like (TM-like, electric field is polarized along the z-axis) modes, the electromagnetic fields have three independent components, E_z, H_x, and H_y. The coordinate system as illustrated in Fig. 1(a) will be used throughout this paper. Using μ_r(r) = 1 in Maxwell equations, we obtain the governing function for E_z

According to Bloch’s Theorem,

where E ^(k) _z(r) = E ^(k) _z(r + R _lmn) is a function with the same periodicity as the PhC. Plugging ε_r(r) = ε_r (r + R _lmn) into the governing function and using the plane-wave method [12], we calculate the corresponding frequencies for a series of wave vectors and depict them as the dispersion surfaces [Fig. 1(b)], where the wave vectors are normalized to k _a/2π, and the frequency to ωa/2πc (c is the speed of light in vacuum). To avoid clutter, we only depict the first two bands and the light cone. As can be seen from Fig. 1(b), the second band is curved downward. As a result, the group velocity, v _g = ∇_k ω(k), is opposite to the phase velocity and negative refraction ensues. In this work, we focused on the intersection of the second band and the light cone because it was suggested that the corresponding effective index n _eff = -1. In order to see that region clearly, we project the second band in the vicinity of the intersection into a series of equi-frequency contours. As shown in Fig. 1(c), these contours are inward-growing and centered at the origin (k_x = 0, k_y = 0). In particular, we can see that the intersection occurs at ω = 0.236. Since the contour of ω = 0.236 is approximately circular, we can treat it as an isotropic homogenous medium. Using the plane wave method, we can also calculate the eigen-field, E ^(k) _z(r), for a given wave vector. Note that the Bloch wave, E_z(k,r,t) = exp(j k ∙ r - jωt)E(k)_z, can be understood as the eigen-field E(k)_z(r) packaged in a plane wave, exp(j k ∙ r - jωt). The corresponding eigen-fields H_x and H_y can be calculated through

 

Fig. 1. (a) The illustration of the PhC and coordinate system for the study. (b) The photonic dispersion diagram for the PhC and the light in this frequency range. (c) The equi-frequency dispersion contours for the second band and the light cone.

Download Full Size | PDF

For simplicity, we consider the case where light propagates along the x-axis. We first tried (k_x = 0.236, k_y = 0) and found that the power propagates opposite to the x-axis, which is predicted by the dispersion of the PhC. So we switched the wave vector to (k_x = -0.236, k_y = 0) and repeated the computation. Figures 2 (a, b) show the amplitude and phase distributions of the corresponding eigen-field, E(k)_z(r), respectively: the amplitude distribution has the same periodicity as ε_r(r) = ε_r (r + R _lmn) and the phase of E(k)_z(r) increases along the x-axis. For comparison, the plane wave, exp(j k ∙ r - jωt), has longer period along the x-axis than that of the PhC, and its phase decreases along the x-axis. See Fig. 2(c). In other words, the phase growth of the package and the carrier waves is along opposite directions, i.e. the carrier wave propagates backward while the package wave propagates forward. The phase distribution of the package wave determines the power flow direction because S(av)_x(r) = -0.5Re{E_z(k,r,t)[H_y(k,r,t)]^*} = -0.5Re{E^(k)_z(r)[H^(k)_y(r]^*}. As a result, the power flows forward. Compared to H ^(k) _y(r),H ^(k) _x(r) is very small and can be neglected. When they evolve with time, the Bloch wave moves forwards, while the carrier wave moves backwards. The supplemental movies illustrate Re{exp(j k ∙ r - jωt)E ^(k) _z(r)} and Re{exp(j k ∙ r - jωt)}, respectively. Different phase growth directions of the carrier and package explain why the group velocity is opposite to the phase velocity.

To obtain the effective permittivity and effective permeability, we need to treat the PhC as a homogenous medium. In a homogenous medium, the solution of Maxwell equations for a single wave vector is plane waves. Therefore, we need to approximate the Bloch waves as plane waves, or equivalently approximate the package waves as constants. A direct approach is to take the average value of the package.

where the averaging region is chosen as unit cells. In our case, we chose the hexagonal cell located in - √3a/4 < x < √3a/4 and -0.5a < y < 0.5a. To avoid positive and negative values cancelling each other, like the case in sin-function, as well as to keep the “handness” of the electromagnetic wave, the points r are grouped into region I with Re{H ^(k) _y(r)} > 0 and region II with Re{H ^(k) _y(r)} < 0, and the averaging is applied in these two regions separately. In region I, we obtained E ^(2Dav) _z = 0.0202 + 3×10^-5, j, H ^(2Dav) _y = 0.0861 + 6×10^-5; in region II, we obtained E ^(2Dav) _z = -0.0241 - 8×10^-5, H ^(2Dav) _y = -0.0822 - 2×10^-4. In both regions, the imaginary part of average values is negligible. If the electric and magnetic fields of a plane wave are both known, we can in turn calculate the effective parameters

In region I, we obtain ε^(I) _eff = - 4.270 and μ^(I) _eff = -0.235. In region II, we obtain ε^(II) _eff = -3.405 and μ^(II) _eff = -0.294. The final effective parameters can be estimated to be ${\epsilon }_{\text{eff}}=\sqrt{{\epsilon }_{\text{eff}}^{(\text{I)}}}\sqrt{{\epsilon }_{\text{eff}}^{\text{II}}}=-3.813$ and ${\mu }_{\text{eff}}=\sqrt{{\mu }_{\text{eff}}^{(\text{I)}}}\sqrt{{\mu }_{\text{eff}}^{\text{II}}}=-0.263$ (If we divide the region by Re{E ^(k) _z (r)} > 0 and Re{E ^(k) _z(r)} < 0, very close result can be obtained because the signs of electric and magnetic fields ensure S ^(av) _x(r) > 0 at every point). Consequently, the PhC simultaneously exhibits negative effective permittivity and negative effective permeability. Note that $n_{\text{eff}}=\sqrt{{\epsilon }_{\text{eff}}}\sqrt{{\mu }_{\text{eff}}}=-1$ is valid and ε_eff μ_eff = (k/ω)² is always satisfied. In addition, these effective parameters are derived from the eigen-fields using effective medium theory [2]. In this approach, $n_{\text{eff}}=\sqrt{{\epsilon }_{\text{eff}}}\sqrt{{\mu }_{\text{eff}}}=-1$ is a derived result, rather than a presumption as used in Ref. [11].

 

Fig. 2. (a) The amplitude distribution of the eigen-field, E ^(k) _z(r). (b) The phase distribution of the eigen-field,E ^(k) _z(r). [Media 1] (c) The phase distribution of the plane wave, exp(j k ∙ r - jωt). [Media 2] (d) The impedance with respect to the location of PhC surface.

Download Full Size | PDF

3. Calculation of effective surface impedance

To obtain the effective impedance, we average the field distribution on the surface according to

Note that in the spatial frequency domain, $E_z^{\left(k\right)}(x_0;k_y)=\frac{1}{L_y}{\int }_{\mathrm{Unit}\mathrm{Cell}}{E}_z^{\left(k\right)}(x_0,y)\exp \left(-j{k}_{y}y\right)\mathrm{dy}$ and $E_z^{\left(k\right)}(x_0;k_y=0)=\frac{1}{L_y}{\int }_{\mathrm{Unit}\mathrm{Cell}}{E}_z^{\left(k\right)}(x_0,y)\mathrm{dy}=E_z^{\left(1\mathrm{Dav}\right)}$. The convolution with the Fourier transform of the plane wave, exp(j k ∙ r), will shift spatial frequency to k. As a result, the averaging process results in the component with the spatial frequency matching that of the plane wave the best. In this fashion, the effective impedance can be calculated using ${\eta }_{\mathrm{eff}}=-\frac{E_z^{\left(1\mathrm{Dav}\right)}}{H_z^{\left(1\mathrm{Dav}\right)}}.$

 

Fig. 3. The simulations (amplitudes) of propagation of light from (a, c) the matched medium and (b, d) the air into the PhC. In (a, b), the PhC surface is located at x = 0; in (c, d), the PhC surface is located at x = √3a/4.

Download Full Size | PDF

We shifted the termination surface along the x-axis and obtained different values, which mean the impedance is location-dependent. So we refer to the effective impedance as effective surface impedance. The imaginary part of surface impedance is negligible. As can be seen from Fig. 2(d), the surface has the maximum impedance when it goes through the center of holes (x = 0), while the surface has the minimum impedance when it is located at the half way between two neighboring hole columns (x = √3a/4). The effective impedance is not constant, so ${\eta }_{\mathrm{eff}}=\sqrt{\frac{{\mu }_{\mathrm{eff}}}{{\epsilon }_{\mathrm{eff}}}}$ is not valid. The reason lies in: 1D-averaging is taken to calculate the surface impedance, while 2D-averaging is used to calculate the effective permittivity and permeability.

To validate the results, we apply the effective impedance into 2D finite-difference time-domain simulations. In each simulation, a collimated light beam (with width 6λ) is normally incident into the PhC from a matched material and the mesh size is 0.004λ (The smaller the mesh size is, the more accurate result can be obtained). The PhC has a specific boundary and the matched material has the corresponding impedance calculated from field-averaging. The edges of the simulated region are paved by perfectly-matched layers (PML). Good transmission requires impedance match, but the refractive index of matched material can be arbitrary — or “ ε_m = ±μ_eff and μ_m = ±μ_eff” are sufficient but not necessary. There are many possibilities to match the surface impedance of the photonic crystal only if $\sqrt{\frac{{\mu }_m}{{\epsilon }_m}}={\eta }_{\mathrm{PhC}}$, but for convenience we assume the permeability of the matched medium, μ_m = 1. The selection of permeability is different from that presented in Ref. [12], where the constraint, $\sqrt{{\epsilon }_m}\sqrt{{\mu }_m}=-1$, was applied.

In particular, when the boundary goes through the center of holes (x = 0), we apply the calculated effective impedance 0.283. As can be seen in Fig. 3(a), very high transmission is achieved. The back-refraction is less than 0.1%. In contrast, the back-reflection increases to 31.2% if the matched medium is switched to the air, as shown in Fig. 3(b). When the boundary goes along the middle line between two hole columns (x = √3a/4), we apply the calculated effective impedance 0.221. As can be seen in Fig. 3(c), very high transmission is achieved. The back-refraction is less than 0.1%. In contrast, the back-reflection increases to 40.7% if the matched medium is switched to the air, as shown in Fig. 3(d).

In addition, different wave vectors have different eigen-fields and these effective parameters depend on eigen-fields. As a result, these parameters are also different when light is incident at different angles. We calculated the eigen-fields for all wave vectors along the equi-frequency dispersion contour corresponding to ω = 0.236 as shown in Fig. 1(c). Figures 4(a,b) show the amplitude and phase distributions, respectively, when wave vector is at 45° to the x-axis. The phase distribution of E ^(k) _z(r) is also opposite to that of the carrier as shown in Fig. 4(c). The surface impedance strongly depends on incident angles. Figure 4(d) plots the effective impedance for the truncation surfaces located at x = 0 and x = √3a/4 changing with respect to incident angle, respectively.

 

Fig. 4. (a) The amplitude distribution of the eigen-field, E ^(k) _z(r)(ω_n = 0.236 and θ = 45° ). (b) The phase distribution of the eigen-field, E ^(k) _z(r)(ω_n = 0.236 and θ = 45° ). (c) The phase distribution of the plane wave, exp(j k ∙ r - jωt). (d) The impedance with respect to the location of PhC surface (ω_n = 0.236 and θ = 45°).

Download Full Size | PDF

4. Summary

To summarize, a plane wave method was employed to solve the eigen-fields of a negative refraction PhC. By averaging the eigen-fields, negative effective permittivity and permeability were obtained and calculated. Effective surface impedance was shown to change with truncation surfaces and verified by the finite-difference time-domain simulations. The power propagation mechanism was also investigated by observing the time-evolution of eigen-fields. The application of this method can be extended to analyze the impedance of PhCs for other applications, such as self-collimation, super-prism, photonic bandgap waveguide, and slow light.