## Abstract

A comprehensive study of the plasmonic thin-film solar cell with the periodic strip structure is presented in this paper. The finite-difference frequency-domain method is employed to discretize the inhomogeneous wave function for modeling the solar cell. In particular, the hybrid absorbing boundary condition and the one-sided difference scheme are adopted. The parameter extraction methods for the zeroth-order reflectance and the absorbed power density are also discussed, which is important for testing and optimizing the solar cell design. For the numerical results, the physics of the absorption peaks of the amorphous silicon thin-film solar cell are explained by electromagnetic theory; these peaks correspond to the waveguide mode, Floquet mode, surface plasmon resonance, and the constructively interference between adjacent metal strips. The work is therefore important for the theoretical study and optimized design of the plasmonic thin-film solar cell.

© 2010 Optical Society of America

## 1. Introduction

Solar cells (SCs) [1, 2], which can provide renewable and clean energy by converting sunlight to electrical power, have attracted much attention in the past few years. Despite the growing importance, we need to reduce the cost of the SCs and increase the energy conversion efficiency before they can successfully replace fossil fuel for electrical power generation.

A careful optical design of the device structure is crucial for optimizing the performance particularly for thin film SCs due to the typically weak light absorption of the materials such as amorphous silicon (A-Si) [3]. Figure 1 shows the schematic diagram of a solar cell structure with the textured back reflector (BR) and the antireflection (AR) coatings. As a light-trapping configuration, the textured BRs [4, 5, 6] have been proposed for extending the optical path length. However, the configuration may suffer from the back surface recombination loss. Moreover, the textured AR coatings [7] have been employed to reduce the reflection of light at the top surface of the SCs. But the AR coatings cannot offer the sufficient light concentration. Therefore, better light trapping and concentration schemes are desirable for further improving the external quantum efficiency of SCs. Surface plasmon resonances (SPRs) are collective oscillations of the free electrons that are confined to surfaces and interact strongly with light resulting in a polariton. SPR usually occurs at the interface between a dielectric with the positive dielectric constant ε_{r}^{d} and a metal with the negative dielectric constant ε_{r}^{m}. Meanwhile, SPR, which is the eigenstate of the Maxwell’s equations for the dielectric-metal structure, only exist when Re(-ε_{r}^{m}) > ε_{r}^{d} is satisfied [8][9]. To excite the SPR by light, a coupling technique providing the wavevector mismatch is required. For SCs, the metallic periodic nanostructures and the sub-wavelength scatterers have been used to excite the SPR [10, 11, 12, 13, 14, 15]. Some of the unique features of the plasmonic thin-film SC are the broadband and wideangle absorption enhancements. Both theoretical and experimental results have demonstrated that the strong-peak absorption enhancements appear at specific wavelengths. As a result, the short-circuit current or the open-circuit voltage is also increased.

In order to describe the propagation and scattering of sunlight within the SC and to optimize the process, Maxwell’s equations have to be rigorously solved. Some of the well-adopted full-wave solvers for optical simulation are the finite-difference time-domain (FDTD) [16, 17, 18] and the FDFD methods [19, 20]. Due to the fact that most optical materials are dispersive, the recursive convolution (RC) method [21] or the piecewise linear recursive convolution (PLRC) method [22] should be used for the FDTD method. For noble metals, such as Ag and Au, the complex dielectric constants have to been described by a large number of summation terms in Lorentz-Drude model [23]. Hence, these recursive convolution implementations will consume a lot of CPU time. However, for the FDFD method, one can use an experimentally tabulated dielectric constants of the dispersive materials directly. Furthermore, compared with the FDFD method, the FDTD method is not easy to treat the periodic boundary conditions particularly for the oblique incidence case [16, 24, 25]. In most places, unless a solar panel is mounted on an expensive tracking system, most of the time light is incident on the array obliquely. Hence, the ability of FDFD method to handle the case of oblique incidence is clearly an important advantage over the FDTD method for modeling the SCs. In addition, the FDFD method can use a variety of flexible difference techniques, such as one-sided difference scheme, to improve the accuracy. But the FDTD method will suffer from instability problem if these techniques are adopted [26].

For the FDFD method, most papers adopted the staggered grid, because it can satisfy the divergence-free condition automatically. However, the averaged material model [19, 20, 27] cannot accurately resolve the local field at the high-contrast interface between dielectric and metal. Besides, few papers have modeled SCs with periodic structures by using the FDFD method. Consequently, it is highly desirable to establish a rigorous and efficient method for modeling SCs, in particular, the periodically patterned structures.

In this paper, we use the FDFD method to discretize the inhomogeneous wave function for modeling the plasmonic thin-film SC. In order to accurately treat the dielectric-metal interface, the flexible one-sided difference scheme is introduced. While the perfectly matched layer (PML) cannot work very well under periodic boundary condition [28], the hybrid absorbing boundary condition (ABC) is proposed here to reduce the spurious numerical reflections from the outmost boundary of the PML. Moreover, we propose the rigorous phase and attenuation constants conditions of the SPR for general lossy materials. In addition, we will unveil the origins of the absorption peaks of SC structures using electromagnetic theory.

## 2. Theoretical modeling

A two-dimensional (2-D) plasmonic thin-film silicon SC structure is shown in Fig. 2. Since the s-polarized incident light cannot excite the SPR [9], we mainly consider the p-polarized light with the electromagnetic components of *H _{z}, E_{x}*, and

*E*. The exp(

_{y}*j*

_{0}

*ω*) time convention is used. For the SCs, all the materials are non-magnetic (i.e.

_{t}*μ*= 1).

_{r}#### 2.1. Finite-difference equation

For the 2-D isotropic and inhomogeneousmedia with the complex dielectric constant of *ε _{r}*(

*x,y*), the wave function of the total field

*H*is given by [29]

_{z}^{t}where *k*_{0} is the wave number of free space. Figure 3 shows the general geometry for the inhomogeneous material treatment. Using the second-order central differences, we have

$$-\frac{{H}_{z}^{t}\left(i,j\right)-{H}_{z}^{t}\left(i-1,j\right)}{{\epsilon}_{r}\left(\frac{i-1}{2,j}\right){\mathrm{\Delta}}_{x}})+O\left({\mathrm{\Delta}}_{x}^{2}\right)$$

where ∆_{x} is the spatial step along the *x* direction. For the p-polarized incident light, the following averaging techniques can be adopted for the dielectric constants, i.e.

where the subscript 1, 2, 3, and 4 denote the small rectangular regions as shown in Fig. 3. Using the notations of Φ_{1} = *H _{z}^{t}*(

*i, j*− 1), Φ

_{2}=

*H*(

_{z}^{t}*i*− 1,

*j*), Φ

_{3}=

*H*(

_{z}^{t}*i, j*), Φ

_{4}=

*H*(

_{z}^{t}*i*+ 1,

*j*), and Φ

_{5}=

*H*(

_{z}_{t}*i, j*+ 1), the continuous inhomogeneous wave equation Eq. (1) can be discretized into FDFD equation as

where

It should be noted that the FDFD form of the inhomogeneous wave equation can be converted to that of the homogeneous wave equation if the four rectangular regions have the same dielectric constant. For the p-polarized incident plane wave, it has been reported that the averaging techniques of Eq. (3) and Eq. (4) are effective for the inhomogeneous material treatment [30]. Here, we show that the averaging techniques can also be rigorously derived by the integral form of the wave equation and the boundary conditions (See Appendix A). Since *H _{z}^{inc}* is known (the incident light) and the FDFD equation of

*H*is expressed as Eq. (5), the FDFD equation for the scattered-field

_{z}^{t}*H*can be derived by

_{z}^{s}#### 2.2. Boundary conditions

As shown in Fig. 2, the ABCs along the *y* directions are used to avoid the spurious reflections of waves at the top and bottom boundaries of computational domain. The complex-coordinate PML [31, 32] is an efficient and accurate ABC which has the form of

where

where *ε*_{0} is the permittivity of free space, and ω is the angular frequency of the incident light. The polynomial variation of the conductivities σ is employed, i.e.

where *L* is the layer number of the PML, *Q* is the order of the polynomial, and *C* is a constant. For reducing the spurious numerical reflections, the optimized settings are set to *L* = 8, *Q* = 3.7, and *C* = 0.02. The proper discretization form [33] for the coordinate-stretched term of Eq. (12) is given by

$$\frac{{H}_{z}^{s}\left(i,j\right)-{H}_{z}^{s}\left(i,j-1\right)}{{s}_{y}\left(j-1/2\right){\mathrm{\Delta}}_{y}}]$$

At the outermost boundary of the PML, the Mur ABC [34] replacing the traditional perfectly electric conductor truncation condition is employed to further reduce the spurious numerical reflections. Taking the top plane *y* = 0 as an example, the second-order Mur ABC can be written as

and its discretized form is given by

where

$${f}_{2}={f}_{3}=1-\mathrm{exp}\left({j}_{0}{k}_{0}{\mathrm{\Delta}}_{y}\right)$$

$${f}_{4}=2{k}_{0}^{2}{\mathrm{\Delta}}_{x}^{2}$$

Particularly, we introduce the exponential difference strategy [35] to improve accuracy.

Regarding the periodic structure, the periodic boundary conditions along the *x* directions need to be implemented. Based on the Floquet theorem, we have

$${H}_{z}^{s}\left(x,y\right)={H}_{z}^{s}\left(x+P,y\right)\mathrm{exp}({j}_{0}{k}_{0}\mathrm{cos}\theta \xb7P)$$

where *P* is the periodicity and *θ* is the incident angle with respect to *x* direction. Compared with the FDFD method, the second equation of Eq. (20) is not easy to be treated by the FDTD method, particularly for the oblique incidence case because one does not know the scattered-field values at a future time in periodic device structures. It should be noted that the periodic boundary conditions should also be implemented at the regions of the hybrid ABC.

Although the FDFD form Eq. (5) can generally treat the dielectric-dielectric and dielectric-metal interfaces, the accuracy will degrade at the interfaces. Thus the one-sided difference scheme is used to rectify the problem. For the horizontal interface (*y* = *y _{h}*) between the media 1 and the media 2, the boundary condition for the scattered magnetic field is

Using the high-order-accurate one-sided differences, we get

$$\frac{\partial}{\partial y}{H}_{z}^{s2}{\mid}_{x=i{\mathrm{\Delta}}_{x}}\approx \frac{-1.5{H}_{z}^{s2}\left(i,j\right)+2{H}_{z}^{s2}\left(i,j+1\right)-0.5{H}_{z}^{s2}\left(i,j+2\right)}{{\mathrm{\Delta}}_{y}}$$

The one-sided difference scheme is flexible and can have higher-order accuracy. For the vertical interfaces, one can also use the one-sided differences.

#### 2.3. Parameter extraction

In this section, we will discuss the details to extract the important parameters including the absorption power density and the zeroth-order reflectance and transmittance. By applying the above FDFD equations and the boundary conditions to all free *N* nodes in the solution region, a sparse matrix equation is formed because only the nearest adjacent nodes affect the value of *H _{z}^{s}* at each node. Hence, the scattered magnetic field can be solved by the iterative methods with the memory and computational complexity of

*O*(

*N*). After this, the total electric field components are calculated by

The electron-hole pair generation depends on the photon energy absorbed by the absorbing material per unit time per unit area, i.e.

where η is the power density, *S _{a}* denotes the region of the absorbing material, ∆

_{Sa}is the area of

*S*, and σ

_{a}_{a}= −ωε

_{0}Im(

*ε*) is the conductivity of the absorbing material.

_{ra}For the SC with the periodic structure, the zeroth-order reflectance and transmittance are the important parameters for optimizing the SC structures and comparing the theoretical results to the experimental ones. Physically, the dips of the reflectance correspond to the absorption peaks. The Floquet modes (space harmonics) for the two-dimensional periodic structure are given by

where

According to the orthogonal properties of the Floquet modes, we get

and

where *A* is the amplitude of the incident light, and *y _{r}* and

*y*are the virtual boundaries for computing the zeroth-order reflectance

_{t}*R*and the zeroth-order transmittance

_{p}*T*, respectively.

_{p}As a testing structure, we consider the periodic dielectric strip model in free space. The incident plane wave is given by

where θ = 85° and *k*_{0} = 2π. The dielectric constant of each strip is taken as *ε _{r}* = 4 − 0.1

*j*The periodicity is

*P*= 0.6

*m*, the thickness of each strip is

*d*

_{3}= 0.5

*m*, and the distance between the adjacent dielectric strips is

*d*= 0.3

_{s}*m*. The spatial steps are set to ∆

_{x}= ∆

_{y}= 0.01

*m*. The zeroth-order reflectance and transmittance are calculated by the FDFD method and the rigorous coupled-wave algorithm [36]. As shown in Fig. 4, the two approaches agree with each other very well.

## 3. Simulation results

We start with a simple semi-infinite structure of A-Si/Au which can allow us to intrusively study SPRs. For the semi-infinite structure, the corresponding *H _{z}* fields in the A-Si and Au layers can be assumed as

where *k, β*, and *α* are the propagation, phase, and attenuation constants, respectively. For the p-polarized plane wave, the reflection coefficient for the upgoing wave in the A-Si layer reflected by the Au layer is given by

The poles of Eq. (34) are determined by

Using the facts that

we can derive the *x̃*-directed propagation constant is

The propagation constant *k _{x}* is also the momentum of the SPR. It is well known that SPR will exist if the condition Re(−

*ε*

_{r}^{Au}) >

*ε*

_{r}

^{Si}is satisfied. But the condition is based on the assumption that

*ε*

_{r}^{Si}and

*ε*

_{r}^{Au}are predominately real. For the real situation, the loss of them cannot be ignored. As a result, the condition is not accurate.

Here we will re-study the problem through comprehensively taking into account the complex dielectric constant and defining the attenuation and phase conditions for the formation of SPR. The *ỹ*-directed propagation constants in the A-Si and Au layers are the double-value functions of *k _{x}*. Considering that the SPR is a surface wave decayed away from the dielectric-metal interface (

*ỹ*= 0), we have

In other words, the surface plasmon poles should be on the proper Riemann sheets of both *k _{y}*

^{Si}and

*k*

_{y}^{Au}, i.e. Im(

*k*

_{y}^{Si}) < 0 and Im(

*k*

_{y}^{Au}) < 0. As shown in Fig. 5, from solving for the root of Eq. (35), the attenuation constant of the field at A-Si layer change sign when the incident wavelength goes through 560 nm. Meanwhile, the locations of the poles of Eq. (34) are changed from the improper Riemann sheets to the proper Riemann sheets. According to the Drude model, the metallic dielectric function is approximated as

where *ω _{p}* is the plasma frequency of Au. Hence, the SPRs exist at long wavelength range where the metal is opaque. This also agrees with the fact that the SPRs exist when the incident wavelength is larger than 560nm. Thus the eigenstates of the Maxwell’s equations for the semi-infinite A-Si/Au structure are the SPRs if the conditions of Eq. (39) are satisfied. Besides the attenuation constants conditions, the SPRs should satisfy the phase constants conditions also, i.e.

The phase constants conditions agree with the oscillation property of the SPRs. Figure 6 shows the contour plot of the eigenstate for *E _{x}* field at 735nm, where the maximum phase constant

*β*of the SPR is achieved by the resonance condition

_{x}*ε*

_{r}^{Si}+

*ε*

_{r}^{Au}≈ 0. It is interesting to note that the field profile looks very symmetric because

*k*becomes very large at 735nm and thus

_{x}*k*

_{y}^{Si}in Eq. (36) and

*k*

_{y}^{Au}in Eq. (37) are comparable to each other. At the wavelength, the maximum amplitudes of the attenuation constants (i.e. the minimum attenuation lengths along

*ỹ*directions) are achieved as well.

For the SC with periodic metal nanopatterns as shown in Fig. 2, the absorbing material is A-Si and the substrate is glass (SiO_{2}). The complex dielectric constants of the materials (Au, A-Si, etc) are taken from [37, 38]. The geometric parameters of the device are set as *d*_{1} = 25nm, *d*_{2} = 120nm, *d*_{3} = 40nm, *d*_{4} = 30nm, *d _{s}* = 100nm, and

*P*= 200nm. The

*y*-directed incident field is the p-polarized plane wave with the amplitude of 1 and the frequency spectrum from 400nm to 800nm. The spatial step is set to ∆

_{x}= ∆

_{y}= 0.5nm. Fig. 7 shows the absorbed power density of the A-Si layer. Using the planar Au layer, the non-strip (planar) structure is also modeled. For the non-strip structure,

*d*

_{2}= 140nm is adopted for achieving the same A-Si area while other parameters are unchanged.

The A-Si bulk material has insufficient absorption from 650 nm to 800 nm due to its lower conductivity at the long wavelength region as shown in the inset of Fig. 7. For both strip and non-strip structures, the absorption η as described by Eq. (25) is significantly enhanced at the long wavelength region as shown in Fig. 7. The enhancement is due to the substantial increase of |**E**|^{2}. Moreover, the strip structure shows even stronger absorption due to the excited SPR and the constructive interference between strips. In fact, the efficient absorption can enhance the external quantum efficiency.

More importantly, our results can show the origins of the absorption peaks which are explained by electromagnetic theory. Along ±*y* directions, we consider the structures as the multilayered medium. The waveguide modes can be approximately found by computing the generalized reflection coefficient *R̃*_{i,i+1} of the medium between the i^{th} layer and the (*i* + 1)^{th} layer where (*d*_{i+1} - *d _{i}*) is the thickness of the (

*i*+ 1) layer. Considering that the excitation is

*y*-directed plane-wave,

*k*

_{i+1,y}=

*k*

_{i+1}=

*k*

_{0}√

*ε*

_{r,i+1}. The waveguide modes can be obtained from the local minima of the generalized reflection coefficient as shown in Fig. 8. As shown in Fig. 7 and Fig. 8, the waveguide modes contribute to all the absorption peaks (A and B) of the non-strip structure. In this case, the interface between the A-Si and Au layers can be considered as a good mirror for trapping the light in the non-strip structure. This is the reason why the absorption enhancement still happens in the planar structure. However, for the non-strip structure, the SPR cannot be excited due to the momentum mismatch. For the strip structure, Fig. 9 shows

*H*field distributions at the absorption peaks. There are two different multilayered media in the strip structure as shown in Fig. 2. At the regions where the Au strips are present, the medium is Air/ITO/A-Si/Au/SiO2/Air. In the other region, the medium is Air/ITO/A-Si/SiO2/Air. The waveguide modes of the former medium contribute to the absorption peaks 1, 2, and 6 (See Fig. 7). The absorption peaks 2 and 4 are mainly due to the waveguide modes of the medium without the Au strip.

_{z}^{t}Since the SC has the periodic structure, the Floquet modes also enhance the absorption. The *H _{z}^{t}* field distributions for the peaks 3 and 4 show the resonant states (Floquet modes) along

*x*directions. The periodic boundaries (PBC in Fig. 2) behave like a magnetic wall or electric wall for the peaks of 3 and 4 respectively. The concentrated field at the absorption layer is due to the Floquet modes that can be observed clearly in Fig. 9(c) and Fig. 9(d).

The SPR is successfully excited by the sub-scatterer strip at the wavelength 745nm as shown in Fig. 5 and Fig. 7. According to the mode conversion theory, the sub-wavelength strip can excite the evanescent wave components, which may provide the momentum mismatch ∆_{k} with the continuous spectrum up to 2π/*d*_{3}. The field profile is shown in Fig. 9(e). The *x*-directed boundaries of the strips achieve better field concentration than the *y*-directed boundaries. The
phenomenon agrees well with the semi-infinite model in which the momentum of the SPR is approximately calculated by Eq. (35)-Eq. (38).

The interference pattern between the adjacent metal strips is shown in Fig. 9(d), Fig. 9(e), and Fig. 9(f). The surface waves decayed away from the y-directed boundaries of the strips interfere with each other. The constructive interference leads to the broadband absorption especially for the red light and even extends to the infrared region. From 700nm to 800nm, the absorption peaks of 4, 5, and 6 are not obvious due to the low conductivities of the A-Si. Thus, a good broadband absorbing material with the larger conductivity is very important for improving the absorption properties of SCs. It should be noted that besides Fig. 9(e), there are some concentrated fields at the dielectric-metal interfaces in Fig. 9(d) and Fig. 9(f). They are also SPRs. The SPRs are weak because there are no obvious momentum dips as shown in Fig. 5.

## 4. Conclusion

Using the FDFD method, the rigorous and efficient model of the thin-film plasmonic SC with the periodic strip has been developed. Compared with the PML, the hybrid ABC works better especially for the periodic structure. Moreover, the material discontinuities are accurately treated by the inhomogeneous wave equation and the one-sided difference scheme. We also proposed the phase and attenuation constants conditions of the SPR for lossy material systems.

By using the semi-infinite dielectric-metal structure, our results show that SPRs will exist if the vertical phase and attenuation constants are negative in both dielectric and metal layers. The sub-wavelength scatterers can excite the evanescent waves, which provide the size-dependent continuous spectrum components. Hence, the subwavelength scatterers can excite the SPRs for the broadband light enhancement.

Our results also show that the origins of the absorption peaks for the periodic strip SC structure can be explained by the waveguide mode, the Floquet mode, the SPR, or the constructive interference between strips. By changing the geometric parameters, the locations of peaks can be modified for optimizing the performance of solar cells.

## A. Integral form of the wave equation

Using the Gauss’s law, the contour integral form for the wave equation is given by

where *S _{m}* denotes the four small rectangles enclosing the center square point as shown in Fig. 3 and ∂

*H*/∂

_{z}^{t}*n*denotes the derivatives of

_{m}*H*normal to the contours ∂

_{z}^{t}*S*. Applying the central differences to Eq. (A-1) yields

_{m}$$-{\int}_{{\mathrm{\Gamma}}_{\mathrm{12}}}\frac{\partial \mathrm{\Phi}}{\partial {n}_{1}}\mathit{dl}-{\int}_{{\mathrm{\Gamma}}_{\mathrm{12}}}\frac{\partial \mathrm{\Phi}}{\partial {n}_{1}}\mathit{dl}$$

$$-{\int}_{{\mathrm{\Gamma}}_{\mathrm{21}}}\frac{\partial \mathrm{\Phi}}{\partial {n}_{2}}\mathit{dl}-{\int}_{{\mathrm{\Gamma}}_{\mathrm{23}}}\frac{\partial \mathrm{\Phi}}{\partial {n}_{2}}\mathit{dl}$$

$$-{\int}_{{\mathrm{\Gamma}}_{\mathrm{32}}}\frac{\partial \mathrm{\Phi}}{\partial {n}_{3}}\mathit{dl}-{\int}_{{\mathrm{\Gamma}}_{\mathrm{34}}}\frac{\partial \mathrm{\Phi}}{\partial {n}_{3}}\mathit{dl}$$

$$-{\int}_{{\mathrm{\Gamma}}_{\mathrm{43}}}\frac{\partial \mathrm{\Phi}}{\partial {n}_{4}}\mathit{dl}-{\int}_{{\mathrm{\Gamma}}_{\mathrm{41}}}\frac{\partial \mathrm{\Phi}}{\partial {n}_{4}}\mathit{dl}$$

According to the boundary conditions at the interfaces

we get the FDFD form of the inhomogeneous wave function, which is the same as Eq. (5).

## Acknowledgement

We would like to acknowledge the support of the grants (#712108) from the Research Grant Council of the Hong Kong Special Administrative Region, China. The first author expresses his sincere gratitude to Mr. Kwan Wing Keung, the University of Hong Kong for his technical support.

## References and links

**1. **J. Nelson, The Physics of Solar Cells (Imperial College Press, London, 2003)

**2. **P. Würfel, Physics of Solar Cells: From Principles to New Concepts (Wiley-VCH, Berlin, 2004).

**3. **K. L. Chopra, P. D. Paulson, and V. Dutta, “Thin-film solar cells: An overview,” Prog. Photovoltaics **12**, 69–92 (2004). [CrossRef]

**4. **L. Zeng, Y. Yi, C. Hong, J. Liu, N. Feng, X. Duan, L. C. Kimerling, and B. A. Alamariu, “Efficiency enhancement in Si solar cells by textured photonic crystal back reflector,” Appl. Phys. Lett . **89**, 111111 (2006). [CrossRef]

**5. **C. Haase and H. Stiebig, “Thin-film silicon solar cells with efficient periodic light trapping texture,” Appl. Phys. Lett . **91**, 061116 (2007). [CrossRef]

**6. **H. Sai, H. Fujiwara, M. Kondo, and Y. Kanamori, “Enhancement of light trapping in thin-film hydrogenated microcrystalline Si solar cells using back reflectors with self-ordered dimple pattern,” Appl. Phys. Lett . **93**, 143501 (2008). [CrossRef]

**7. **W. Zhou, M. Tao, L. Chen, and H. Yang, “Microstructured surface design for omnidirectional antireflection coatings on solar cells,” J. Appl. Phys . **102**, 103105 (2007). [CrossRef]

**8. **R. Zia, M. D. Selker, P. B. Catrysse, and M. L. Brongersma, “Geometries and materials for subwavelength surface plasmon modes,” J. Opt. Soc. Am. A **21**, 2442–2446 (2004). [CrossRef]

**9. **J. Homola, S. S. Yee, and G. Gauglitz, “Surface plasmon resonance sensors: Review,” Sens. Actuator B . **54**, 3–15 (1999). [CrossRef]

**10. **H. Ditlbacher, J. R. Krenn, N. Felidj, B. Lamprecht, G. Schider, M. Salerno, A. Leitner, and F. R. Aussenegg, “Fluorescence imaging of surface plasmon fields,” Appl. Phys. Lett . **80**, 404–406 (2002). [CrossRef]

**11. **H. Raether, Surface Plasmons on Smooth and Rough Surfaces and on Gratings (Springer-Verlag, Berlin, 1988).

**12. **K. Kato, H. Tsuruta, T. Ebe, K. Shinbo, F. Kaneko, and T. Wakamatsu, “Enhancement of optical absorption and photocurrents in solar cells of merocyanine Langmuir-Blodgett films utilizing surface plasmon excitations,” Mater. Sci. Eng. C-Biomimetic Supramol. Syst . **22**, 251–256 (2002)
[CrossRef]

**13. **K. Tvingstedt, N. K. Persson, O. Inganas, A. Rahachou, and I. V. Zozoulenko, “Surface plasmon increase absorption in polymer photovoltaic cells,” Appl. Phys. Lett . **91**, 113514 (2007). [CrossRef]

**14. **V. E. Ferry, L. A. Sweatlock, D. Pacifici, and H. A. Atwater, “Plasmonic nanostructure design for efficient light coupling into solar cells,” Nano Lett . **8**, 4391–4397 (2008) [CrossRef]

**15. **R. A. Pala, J. White, E. Barnard, J. Liu, and M. L. Brongersma, “Design of plasmonic thin-film solar cells with broadband absorption enhancements,” Adv. Mater . **21**, 3504–3509 (2009) [CrossRef]

**16. **M. Qiu and S. L. He, “A nonorthogonal finite-difference time-domain method for computing the band structure of a two-dimensional photonic crystal with dielectric and metallic inclusions,” J. Appl. Phys . **87**, 8268–8275 (2000). [CrossRef]

**17. **A. Farjadpour, D. Roundy, A. Rodriguez, M. Ibanescu, P. Bermel, J. D. Joannopoulos, S. G. Johnson, and G. W. Burr, “Improving accuracy by subpixel smoothing in the finite-difference time domain,” Opt. Lett . **31**, 2972–2974 (2006). [CrossRef] [PubMed]

**18. **K. G. Ong, O. K. Varghese, G. K. Mor, K. Shankar, and C. A. Grimes, “Application of finite-difference time domain to dye-sensitized solar cells: The effect of nanotube-array negative electrode dimensions on light absorption,” Sol. Energy Mater. Sol. Cells **91**, 250–257 (2007). [CrossRef]

**19. **Z. M. Zhu and T. G. Brown, “Full-vectorial finite-difference analysis of microstructured optical fibers,” Opt. Express **10**, 853–864 (2002) [PubMed]

**20. **C. P. Yu and H. C. Chang, “Compact finite-difference frequency-domain method for the analysis of two-dimensional photonic crystals,” Opt. Express **12**, 1397–1408 (2004). [CrossRef] [PubMed]

**21. **R. Luebbers, F. P. Hunsberger, K. S. Kunz, R. B. Standler, and M. Schneider, “A frequency-dependent finite-difference time-domain formulation for dispersive materials,” IEEE Trans. Electromagn. Compat . **32**, 222–227 (1990). [CrossRef]

**22. **D. F. Kelley and R. J. Luebbers, “Piecewise linear recursive convolution for dispersive media using FDTD,” IEEE Trans. Antennas Propag . **44**, 792–797 (1996). [CrossRef]

**23. **G. Veronis and S. Fan, “Overview of Simulation Techniques for Plasmonic Devices,” in Surface Plasmon Nanophotonics,
M. L. Brongersma and P. G. Kik, eds., (Springer, Dordrecht, The Netherlands, 2007)

**24. **M. E. Veysoglu, R. T. Shin, and J. A. Kong, “A finite-difference time-domain analysis of wave scattering from periodic surfaces: Oblique-incidence case,” J. Electromagn. Waves Appl . **7**, 1595–1607 (1993). [CrossRef]

**25. **A. Taflove and S. C. Hagness, Computational Electrodynamics: the Finite-Difference Time-Domain Method, Third Edition (Artech House, Boston, 2005).

**26. **S. Zhao and G. W. Wei, “High-order FDTD methods via derivative matching for Maxwell’s equations with material interfaces,” J. Comput. Phys . **200**, 60–103 (2004). [CrossRef]

**27. **D. M. Sullivan, Electromagnetic Simulation Using the FDTD Method (Wiley-IEEE Press, New York, 2000).

**28. **A. F. Oskooi, L. Zhang, Y. Avniel, and S. G. Johnson, “The failure of perfectly matched layers and towards their redemption by adiabatic absorbers,” Opt. Express **16**, 11376–11392 (2008). [CrossRef] [PubMed]

**29. **W. C. Chew, Waves and Fields in Inhomogenous Media (Van Nostrand Reinhold, New York, 1990).

**30. **S. G. Johnson and J. D. Joannopoulos, “Block-iterative frequency-domain methods for Maxwell’s equations in a planewave basis,” Opt. Express **8**, 173–190 (2001). [CrossRef] [PubMed]

**31. **W. C. Chew and W. H. Weedon, “A 3-D perfectly matched medium from modified Maxwell’s equations with stretched coordinates,” Microw. Opt. Technol. Lett . **7**, 599–604 (1994). [CrossRef]

**32. **W. C. Chew, J. M. Jin, and E. Michielssen, “Complex coordinate stretching as a generalized absorbing boundary condition,” Microw. Opt. Technol. Lett . **15**, 363–369 (1997). [CrossRef]

**33. **C. M. Rappaport, M. Kilmer, and E. Miller, “Accuracy considerations in using the PML ABC with FDFD Helmholtz equation computation,” Int. J. Numer. Model.-Electron. Netw. Device Fields **13**, 471–482 (2000). [CrossRef]

**34. **G. Mur, “Absorbing boundary-conditions for the finite-difference approximation of the time-domain electromagnetic-field equations,” IEEE Trans. Electromagn. Compat . **23**, 377–382 (1981). [CrossRef]

**35. **J. P. Berenger, “A perfectly matched layer for the absorption of electromagnetic-waves,” J. Comput. Phys . **114**, 185–200 (1994). [CrossRef]

**36. **M. G. Moharam, E. B. Grann, D. A. Pommet, and T. K. Gaylord, “Formulation for stable and efficient implementation of the rigorous coupled-wave analysis of binary gratings,” J. Opt. Soc. Am. A-Opt. Image Sci. Vis . **12**, 1068–1076 (1995). [CrossRef]

**37. **E. D. Palik, Handbook of Optical Constants of Solids (Academic Press, London, 1998).

**38. **X. W. Chen, W. C. H. Choy, and S. L. He, “Efficient and rigorous modeling of light emission in planar multilayer organic light-emitting diodes,” J. Disp. Technol . **3**, 110–117 (2007). [CrossRef]