Open Access
Add to BibSonomyAbstract
The ability of the Finite-Difference Time-Domain method to model a perfect lens made of a slab of homogeneous left-handed material (LHM) is investigated. It is shown that because of the frequency dispersive nature of the medium and the time discretization, an inherent mismatch in the constitutive parameters exists between the slab and its surrounding medium. This mismatch in the real part of the permittivity and permeability is found to have the same order of magnitude as the losses typically used in numerical simulations. Hence, when the LHM slab is lossless, this mismatch is shown to be the main factor contributing to the image resolution loss of the slab.
©2005 Optical Society of America
1 Introduction
Since the conceptual introduction of the perfect lens imaging [1], the study of subwavelength imaging via both numerical simulations and experiments has been an active research topic. As shown in Ref. [1], a perfect lens can be theoretically achieved with a left-handed material (LHM) slab perfectly matched to the surrounding right-handed material (RHM), e.g. an LHM slab with its relative permittivity and permeability both equal to -1 located in vacuum. The theoretical concept of a flat LHM slab lens was verified using the Finite-Different Time-Domain (FDTD) method in Ref. [2]. Using a similar technique, the subwavelength resolution imaging was demonstrated in Ref. [3] and the effect of finite size slab was studied in Ref. [4]. Later, the Pseudospectral Time-Domain method (PSTD) was used for the subwavelength imaging study with one line source [5] and with two line sources [6]. In all these numerical studies, the refractive index of the LHM slabs are assumed to have a real part of exactly -1 while the imaginary part is considered to be the limiting factor for the image resolution. However, it was observed in simulations that the image resolution cannot be improved by reducing the material losses [3], which seems contradictory to the perfect lens theory. Several possible reasons, including the perturbation of the refractive index of the LHM, were mentioned in Ref. [3] without being studied in details.
In this paper, the limitations of the FDTD to model a perfect LHM lens is studied, with a special emphasis on the maximum resolution that can be simulated. In order to do so, we first describe a different field averaging scheme implemented in our FDTD method. Second, an analytical formula for the constitutive parameters of the LHM material is derived and illustrates the mismatch due to the time discretization. Finally, we show that for lossless LHM slabs, the image resolution limit is mainly determined by the mismatch in the real part of the constitutive parameters, which is an inherent limiting factor in FDTD simulations.
2 FDTD simulation
In our study, a 2D FDTD scheme is used in which we assume the non-zero fields to be Ey , Hx and Hz . The electric fields are setup at the grid centers while the magnetic fields are setup at the grid edges. A sinusoidal line source Ey (x,z,t) = δ(x)δ(z) sin(ωot)f(t), where f(t) is a smooth ramp function with a rise time of 30 wave periods, is located at (0,0) and the LHM slab is located at z=d/2 with a thickness of d in ẑ and a length of L (from - L/2 to L/2) in x̂. The grid size is chosen initially as λ/100 where λ is the free space wavelength at the line source frequency ωo . The size of the simulation domain is 1024 by 220 grids in x̂ and ẑ, respectively. In order to avoid the ambiguity of studying subwavelength imaging from periodic line sources, a finite LHM slab is used instead of a infinite slab which requires periodic boundary conditions. Hence PML are used at all the boundaries to absorb the radiated fields. The length of the slab L is chosen to be 10λ for all the simulation results in this paper, and has been verified to be long enough to minimize edge effects.
The permeability μ and permittivity ε of the frequency dispersive LHM are both taken to obey a Drude mode, which is implemented using the Auxiliary Differential Equation (ADE) method [2] [7], where the electric and magnetic polarization currents are introduced to achieve the dispersive constitutive parameters. However, instead of setting up both discretized electric and magnetic polarization current densities at the center of grid as in Ref. [2], we keep the current densities to be aligned with the corresponding field quantities and implement a different field averaging scheme based on the integral form of the Maxwell’s equations. At the LHM and air boundary, the field averaging is done by applying the integral equation:
where B̅ and M̅ are the magnetic flux and the magnetic polarization current density at the edges of the boundaries. For example, if the boundary is along ẑ, Eq. (1) is discretized in the finite-difference scheme as
where I and K are the grid locations in x̂ and ẑ, n is the time step, ∆x is the grid size in x and ∆t is the time step size. Reorganizing the terms, it is straightforward to obtain the update equation for Hz as
while the update equations for Ey and Hx are left unchanged, and the update equations at other boundaries can be obtained similarly. The implementation of the above field averaging technique avoids the averaging of the polarization current on every grid inside the LHM slab and also clearly defines the material boundary.
3 LHM material implementation
As mentioned above, the ADE method introduces auxiliary polarization current densities to describe the dispersive constitutive parameters of the LHM slab. Typically, the frequency dispersive permittivity is implemented using the following partial differential equations:
where εr is the relative permittivity, J̅e is the electric polarization current, ωpe and Γe are the electron plasma frequency and collision freqency respectively. In this study, we choose the electric plasma frequency to be the same as the magnetic one so the relative permittivity and permeability have a same value. After discretizing Eq. (5) and substituting Je from Eq. (5) into Eq. (4), we obtain the numerical εr as:
It is clear that Eq. (6) approaches the Drude model εr = 1 - / in the limit of ∆t → 0, which gives a value of -1 when ωpe = √2ωo . However, for a finite ∆t used in an actual simulation, εr presents a slight deviation from exactly -1 at the same ωpe . As an example, the value of εr from Eq. (6) is about -1.0003 for a typical grid size of λ/100, which is also the value of the refractive index since we choose here a magnetic plasma frequency identical to the electric one. This small perturbation does not affect the propagating waves significantly. However, the resolution of a subwavelength imaging system is critically dependent on the reconstruction of the evanescent waves, or part of the evanescent wave spectrum, by the LHM slab. This reconstruction is in turn critically dependent on the slab’s constitutive parameters [8], and the slight mismatch of 0.03% in the real part has an important impact on the resolution of the image as we shall see hereafter. It is important to note that this small perturbation in the real part of the constitutive parameters has often been overlooked and the imaginary parts with a value in the same order are typically considered to be the main contributor for limiting the image resolution.
4 Numerical examples
To further illustrate the influence of the perturbed material properties on the image resolution in simulations, a lossless LHM slab with a thickness d = 0.2λ is used to image a line source located at d/2 away from its interface. The analytical calculation is carried out using the method outlined in Ref. [9] for an infinite slab with μr and εr equal to -1.000297. The calculated and simulated electric field magnitudes at the image plane (located at z = 2d) are compared and the results are shown in Fig. 1. It can be seen that both results have a similar peak profile, which indicates that the image resolution computed from the simulation is in agreement with the one predicted analytically. The small differences between the results away from the main peak can be attributed to the difference between the infinite slab used in the analytical calculation and the finite slab used in the simulation. However, these differences have a negligible effect on the image resolution. In addition, it is clear that because of the introduced mismatch, both results differ from the electric field magnitude of a perfect image, which has an infinitely small width of the peak.
The imaging ability of the LHM slab can also be quantified by looking at the spectrum representation of the fields at the image plane. This is shown in Fig. 2 where the normalized spectrum of the simulated electric field of Fig. 1 is compared with the spectrum obtained analytically (by plotting the transmission coefficient of the LHM slab from the source plane to the image plane). The two spectra extend up to around 7ko before the cut-off, which suggest that the resolution of the image is about λo /7. In addition, it can be seen that the two spectra agree well both qualitatively and quantitatively, which confirms that the analytical calculations using the mismatched values predict the image resolution obtained from the simulation. The sharp peak in the analytical spectrum due to the singularity of the pole from the LHM slab’s transmission coefficient [9], [10], [11] is well predicted in the simulated spectrum.
It is known that the image resolution is determined by the maximum transverse wave number (maximum kx ) that can be restored at the image plane, and it was pointed out in Ref. [9] that the image resolution is related to the location of the pole in the transverse wavenumber spectrum. From Fig. 2, it is shown that the maximum restored transverse wave numbers are close to the location of the singularity, which again confirms the observation in Ref. [9]. In addition, Fig. 2 shows the same comparison study for a slab thickness of d = 0.1λ. Again, the analytical and simulated spectra are in good agreement with the image spectrum which extends to about 14ko . The comparison between the results from these two slab thicknesses show that the image resolution is approximately improved by a factor two when the LHM slab thickness is reduced by one half, which is consistent with the theoretical prediction [8][9].
The subwavelength resolution of the LHM slab can be directly visualized by considering the case of two line sources separated by 0.2λ and located 0. 1λ away from the slab interface. The time-averaged Poynting power densities at the image plane from the simulations and analytical calculations are compared in Fig. 3, where the two peaks are clearly recognizable and the features agree well between simulation and the analytical calculation. The results from the above study suggest that the simulated resolution in FDTD is limited mainly by the mismatch between the real part of constitutive parameters of the LHM slabs and the ones of the surrounding vacuum.
According to Eq. (6), the accuracy of εr is related to the time step size rather than the grid size. To illustrate this, the LHM slab used in Fig. 1 is simulated again using a different grid size while keeping the same time step size. The results are shown in Fig. 4 where two grid sizes are used, namely λ/200 and λ/100, but with an identical time step size. The resolution limits from these two results are almost identical and both are improved slightly (from 7ko to about 7.5ko ) from the resolution shown in Fig. 2 due to the reduction of time step size.
Fig. 1. Comparison of |Ey | at the image plane from the FDTD simulation and the analytical calculation for an LHM slab with a thickness of 0.2λ. The grid size used in the simulation is λ/100. The analytical calculation considers an LHM slab of μr = εr = -1.000297 in vacuum.
Fig. 2. Comparison of Ey spectra at the image plane from the FDTD simulations and the analytical calculations for two slab configurations: one with a thickness of 0.2λ and the other with a thickness of 0.1λ. Both slabs are simulated with a same grid size of λ/100.
5 Conclusions
A study of the LHM slab’s constitutive parameters in the FDTD implementation is presented. It is shown that there exist a mismatch between the slab and its surrounding medium due to the time discretization rather than space discretization. By comparing the simulation results and analytical calculations, we demonstrate that the simulated image resolution of an LHM perfect lens is mainly limited by this mismatch. In other applications such as the simulation study of surface plaritons at LHM/RHM interfaces where the matching condition is required, the understanding of this limitation in FDTD can also be very important.
Fig. 3. Comparison of time averaged Poynting power densities < Sz > at the image plane from the FDTD simulation and the analytical calculation for the two line source imaging. The LHM slab is the same as in Fig. 1. The line sources are separated by 0.2λ.
Fig. 4. Comparison of Ey spectrum at the image plane from the FDTD simulations using different grid sizes. The time step size is calculated from Courant’s criterion based on the grid size of λ/200 and is adopted for the simulation using the grid size of λ/100 as well. The LHM slab is the same as in Fig. 1.
Acknowledgments
This work was sponsored by the Department of the Air Force under Air Force Contract F19628-00-C-0002, DARPA under Contract N00014-03-1-0716, and the ONR under Contract N00014-01-1-0713. Opinions, interpretations, conclusions and recommendations are those of the author and are not necessarily endorsed by the United States Government.
References and links
1 . J. Pendry , “ Negative Refraction Makes a Perfect Lens ,” Phys. Rev. Lett. 85 , 3966 – 3969 ( 2000 ). [CrossRef] [PubMed]
2 . R. W. Ziolkowski and E. Heyman , “ Wave Propagation in Media Having Negative Permittivity and Permeability ,” Phys. Rev. E 64 , 056625 ( 2001 ). [CrossRef]
3 . S.A. Cummer , “ Simulated causal subwavelength focusing by a negative refractive index slab ,” Appl. Phys. Lett. 82 , 1503 – 1505 ( 2003 ). [CrossRef]
4 . L. Chen , S. He , and L. Shen , “ Finite-Size Effects of a Left-Handed Material Slab on the image Quality ,” Phys. Rev. Lett. 92 , 107404 ( 2004 ). [CrossRef] [PubMed]
5 . M.W. Feise , J.B. Schneider , and P.J. Bevelacqua , “ Finite-Difference and Pseudospectral Time-Domain Methods Applied to Backward-Wave .Metamaterials ,” IEEE Trans. Antennas. Propag. 52 , 2955 – 2962 ( 2004 ). [CrossRef]
6 . M.W. Feise and Y.S. Kivshar , “ Sub-wavelength imaging with a left-handed material flat lens ,” Phys. Lett. A 334 , 326 – 323 ( 2005 ). [CrossRef]
7 . A. Taflove , Computational Electrodynamics - The Finite-Difference Time-Domain Method , ( Artech House , 1995 ).
8 . D.R. Smith , D. Schurig , M. Rosenbluth , and S. Schultz , “ Limitations on subdiffraction imaging with a negative index slab ,” Appl. Phys. Lett. 82 , 1506 – 1508 ( 2003 ). [CrossRef]
9 . J. Lu , T. M. Grzegorczyk , B.-I. Wu , J. Pacheco , M. Chen , and J. A. Kong , “ Effect of poles on subwavelength focusing by an LHM slab ,” Microwave Opt. Technol. Lett. 45 , 49 ( 2005 ). [CrossRef]
10 . J. A. Kong , Electromagnetic Wave Theory , ( EMW , 2000 ).
11 . J. A. Kong , “ Electromagnetic interactions with stratified negative isotropic media ,” Progress In Electromagnetic Research 35 , 1 – 52 ( 2001 ). [CrossRef]
