## Abstract

A parallelized 3D FDTD (Finite-Difference Time-Domain) solver has been used to study the near-field electromagnetic intensity upon plasmonics nanostructures. The studied structures are obtained from AFM (Atomic Force Microscopy) topography measured on real disordered gold layers deposited by thermal evaporation under ultra-high vacuum. The simulation results obtained with these 3D metallic nanostructures are in good agreement with previous experimental results: the localization of the electromagnetic intensity in subwavelength areas (“hot spots”) is demonstrated; the spectral and polarization dependences of the position of these “hot spots” are also satisfactory; the enhancement factors obtained are realistic compared to the experimental ones. These results could be useful to further our understanding of the electromagnetic behavior of random metal layers.

©2012 Optical Society of America

The understanding of the interaction between light and disordered matter is of great interest for many research fields. For example, the scattering through random media have attracted considerable attraction in a double perspective: on account of fundamental issues in physics [1] and for imaging problems in complex media [2]. In particular, due to their unique optical properties, the metallic disordered rough layers on dielectric surfaces have been widely studied both theoretically and experimentally over the two last decades [3–7]. These random films are made of clusters of metal nanoparticles. Depending on the metallic filling fraction, three different regimes have been identified: (i) for low concentration, metal nano-islands are distributed on the dielectric substrate. (ii) For intermediate concentration at the vicinity of the percolation transition, growing clusters have coalesced and they exhibit self-similarity with fractal dimension. (iii) For high metal filling fraction, the deposited film is now an almost continuous metallic layer with dielectric voids. In the intermediate regime (ii) of interest, two main properties have been identified: a large absorption over the red and near-infrared spectrum [8] and strong spatial fluctuations of the near-field intensity [9]. Theoretical studies give a quite simple description of these behaviors [10]. In this specific metallic filling fraction, the multi-scale clusters exhibit a self-similarity character, meaning that parts resemble the whole in a statistical sense. In the case of these disordered metal films, the self-similarity corresponds to a wide distribution of the sizes and shapes of the clusters. Due to surface plasmon resonances [11], an anomalous spectral absorption and a spectral broadening are observed. The random character of the self-organized metal clusters is also responsible of the multiple scattering, leading to light confinement and enhancement at a sub-wavelength scale. This spatial localization of the electromagnetic field leads to “hot spots” [3] and huge fluctuations of the near-field intensity [8].

Fluctuations and enhancements of the near-field intensity above such random metal films have been reported by different groups using Near-Field Optical Microscopy (NSOM) [3,5,12,13]. More recently, hot spots have been studied using the PEEM (Photo-Emission Electron Microscopy) technique [14,15]. The fluctuations of the near-field electromagnetic energy, evidenced by the hot spots, have also been related to the local density of states (LDOS) through the statistical distribution of fluorescent-molecule lifetimes [16].

Numerical simulations, based on effective-medium theory [17], quasi-static calculations [7] or scaling theory [18,19], have also predicted anomalous absorption and localized “hot spots”. More recently, full wave numerical simulations (FDTD for Finite-Difference Time Domain) have been used to describe the near and the far-field optical responses of disordered metal-dielectric composites [20–22]. But in all these numerical studies, the metal-dielectric structure was either: a sample computer-generated random film [23]; or generated with a Kinetic Monte-Carlo algorithm [24]; or obtained from SEM (Scanning Electron Microscopy) or TEM (Transmission Electron Microscopy) images of a real structure [21]. Despite the diversity of the techniques used, the structures under study were in all cases 2D films: the roughness of the metal layer is never taken into account. Some unresolved discrepancies between experimental and numerical results still remain: for example, the measured enhancement factors are much lower than the theoretical or numerical predicted ones; the effect of the incident polarization is not described by theoretical predictions.

In this article, we show that the combination between FDTD simulations and the roughness of random metal samples (real 3D nanostructures) gives results that are in very good agreement with experimental data. We study a disordered gold layer at the vicinity of the percolation threshold. We verify some previously confirmed experimental results: the spatial distribution of the near-field intensity strongly depends on the incident photon energy; the incident polarization state plays a major role on the local energy distributions; the near-field intensity is dominated by the non-propragative contribution; the sizes of the hot spots is in good agreement with experimental measurements. We also show that the Probability Density Function (PDF) of the enhancement factors is well fitted by a log-normal function [15,25,26].

The metallic structures studied are random gold films near the percolation threshold, just below the transition’s concentration. Gold is evaporated on a glass cover-slip, under ultra-high vacuum conditions (10^{−9} Torr). The mass thicknesses of the samples are between 6.1 nm and 6.4 nm and the typical rms roughness is about 4.3 nm. An AFM image [27] of a standard layer, 512 x 512 points, is shown in Fig. 1(a)
. The filling factor is about 64%. This image and the corresponding typical profile shown in Fig. 1(b) were obtained with an ultra-sharp probe (Veeco TESP-SS), in order to get the most accurate topography.

At the vicinity of the percolation threshold, it is well-known that Euclidian and fractal clusters co-exist. The perimeter (*P*) versus the surface (*S*) of the clusters for a given sample is shown in Fig. 1(c). From the well-known relation *P* ∝ *S ^{D}*

^{/2}, it is possible to determine the dimension

*D*of the largest non-Euclidian clusters. We obtain

*D*~1.9 which is in very good agreement with the theoretical value 1.88.

The topographic data *z*(*x*,*y*) describing the surface of the layer are then transformed in an array of complex ($\tilde{n}=n+ik$) refractive index [28], for the dielectric (air) and for the metal (gold) on each point. Then the *(n, k)* array is implemented in the parallelized 3D FDTD solver (Lumerical, FDTD Solutions) in order to define the 3D structure studied, supported by a glass substrate. The ensemble (substrate + gold layer) is then illuminated with various incident wavelengths under normal incidence through the substrate.

The FDTD method or Yee algorithm [29,30] is a well-known technique based on space and time discretizations of the Maxwell curl equations

*leapfrog*scheme.

So the time and space evolutions of every field component can be written as, for example the *x* component of the electric field,

*i*,

*j*,

*k*,

*n*are integers that are running on the space and time mesh, the terms $\pm 1/2$ are due to the staggered field components, $\Delta y$ and $\Delta z$ (and non appearing Δ

*x*) are the space mesh steps and $\Delta t$ is the time mesh step.

Due to the limitation of the computer memory, the simulation cell has to be restricted to a finite region of space (500x500x2500 nm^{3}). To avoid non-physical reflections on the boundary of the simulation cell and due to the strong scattering behavior of the metallic structure which leads to a wide range of incidence of the scattered fields on the boundaries, we have used the Perfectly Matched Layers (PML) as boundary conditions with at least 216 layers. This ensures a maximum reflection of 10^{−7} over all incidences in the half-space. Due to the cubic nature of the FDTD mesh, it is not possible to avoid stairs effects on the typical random metallic surfaces described in Fig. 1. To minimize this issue as much as possible, we have defined the value of the maximum mesh step equal to 0.8 nm on the volume of the layer that gives at least 11 x 10^{9} simulation grid points to describe the metallic layer. The substrate which is a homogeneous glass coverslip, is also discretized on half of its thickness, just below the metallic nanostructure. Before performing the simulations whose results are presented in this article, we have checked that a double-sized meshing (1,6 nm) did not fundamentally modify the localization of the hot spots but, of course, with a poorest resolution, wider hot spots and weaker intensities.

In Fig. 2
, the results of our FDTD simulations are shown for a single geometry in order to compare the localization and the values of local intensities. All the calculated intensities above the gold layer are normalized with respect to the incident intensity: in other words, rather than intensity maps, these are enhancement factor's maps. For presenting the simulation results in the same manner as the experimental ones previously obtained with SNOM technique [6,8], these maps were calculated for each (*x*, *y*) point, at a distance 10 nm [3] on top of the corresponding *z*(*x*, *y*) value. This leads to representations similar to the ones obtained in “constant distance mode” for SNOM technique.

To improve clarity, in the plots of Fig. 2(b) and (c), the maximum value of the color bar was set at the mean value of the respective maximum. One has to underline that the polarization states of the normal incident plane waves are identical. From the enhancement factors (Γ) maps, large local-field enhancements ($\Gamma \simeq 14$ at maximum for $\lambda =650$nm and $\Gamma \simeq 27$ for $\lambda =1$µm) are revealed. And the positions of the highest enhancements are clearly wavelength-dependent.

Two points have to be noticed. First, the Γ values are rather low (one or two orders of magnitude less) compared to the previously simulated and published ones [4,21], but they are in a quite good agreement with the experimental results measured with aperture-SNOM technique [6,8] or with the PEEM technique [15]. And the maximum Γ values increase with the wavelength.

Secondly, we can analyze the Probability Density Functions (PDF) of the Γ distribution. The PDFs, shown in Fig. 2(d) and (e), are well fitted by the log-normal distribution given by

*K*is a constant,

*σ*is the standard deviation of the distribution lnΓ and $\u3008\mathrm{ln}\Gamma \u3009$ is the mean value of the logarithm of the enhancement factor Γ. This behavior was previously identified in [31]. It is clearly seen that, for the same topography, larger values of Γ are obtained for wavelengths in the near infrared. And, as theoretically predicted [10], the probability of occurrence for the largest Γ values is very low which means that the enhancement factor distribution displays a long tail. Our simulation results also demonstrate the evidence of a distribution’s broadening as it was recently observed in [15].

Another result of our FDTD simulations is that the locations of the hot spots are depending on the incident polarization state. In Fig. 3 , we present, for the same topography presented in Fig. 2(a), the spatial fluctuations of the enhancement factor for two incident wavelengths and for two linear orthogonal polarization states: λ = 1 µm and (a) and (b) for linear X and Y polarizations respectively; λ = 1.3 µm and (c) and (d) for linear X and Y polarizations respectively.

This behavior, previously observed in experimental results [8,14], clearly means that, in addition to a dependence of the incident photon’s energy, the hot spots positions are also due to a complex interaction between the topographic configurations (nanogaps, tips,…) and the incident polarization of the electromagnetic field.

One can observe in Fig. 2 and 3 that the sizes and the shapes of the hot spots are very different, depending both on the incident wavelength and on the local structure. Considering the hot spots with intensities larger than 90% of the maximum intensity, we can deduce that the sizes (Half-Width at Half Maximum) of the hot spots are ranging between 7 nm to 32 nm which is in very good agreement with the values determined in [3,14,32]. The shapes of the hot spots in the (*x*, *y*) plane vary from slightly to strongly stretched, depending on the gaps of coupled metal junction and/or the polarization of the incident electromagnetic field. The relationship between the polarization direction of the incident field and the orientation of the inter-particles is not as clear as previously demonstrated in [33,34], but, in our case, the metallic structure studied is much more complex than 2 or 3 metal spheres in interaction.

Extracting the different components (*x*, *y*, *z*) of the calculated intensity is possible with the FDTD method. If we consider the near-field intensity distribution at constant height, for example *z* = 12 nm in Fig. 4(a), (b), (c)
, it is clear that the main contribution to the intensity is due to ${I}_{z}={\left|{E}_{z}\right|}^{2}$, where *E _{z}* is the component of the local electromagnetic field parallel to the incident direction. It is no more true at some distance (

*z*= 50 nm)) from the metal layer (Fig. 4(d), (e) and (f)) and the main contribution is parallel to the incident field.

So it is confirmed that, with our FDTD simulations, the major contribution to the intensity in the very near-field zone is due to the non-propagative component, parallel to the incident wave-vector.

In conclusion, the introduction of surface roughness in the 3D parallel FDTD simulations performed on fractal gold films gives results in good agreement with the experimental near-field measurements. The local fluctuations in the near-field intensity are quite well reproduced and their statistical behavior is very satisfying compared to the experimentally observed ones. These results could be useful to further the understanding of the local fields. These local field fluctuations are a signature of the inhomogeneous character of the disordered metallic films, as for the fluctuations of the LDOS. So, by modifying the nature of the excitation source, the FDTD simulations could be a useful tool to understand the interactions between such random metallic layers and single nano-emitters deposited upon it.

## Acknowledgments

X. Q. wishes to acknowledge support from the ANR PNANO “NANOEC”.

## References and links

**1. **P. Sheng, *Introduction to Wave Scattering, Localization, and Mesoscopic Phenomena* (Academic Press, 1995).

**2. **P. Sebbah, *Waves and Imaging through Complex Media*, (Kluwer, 2001).

**3. **S. Grésillon, L. Aigouy, A. C. Boccara, J.-C. Rivoal, X. Quélin, C. Desmarest, P. Gadenne, V. A. Shubin, A. K. Sarychev, and V. M. Shalaev, “Experimental observation of localized optical excitations in random metal-dielectric films,” Phys. Rev. Lett. **82**(22), 4520–4523 (1999). [CrossRef]

**4. **A. K. Sarychev and V. M. Shalaev, “Electromagnetic field fluctuations and optical nonlinearities in metal-dielectric composites,” Phys. Rep. **335**(6), 275–371 (2000). [CrossRef]

**5. **S. Ducourtieux, V. A. Podolskiy, S. Grésillon, S. Buil, B. Berini, P. Gadenne, A. C. Boccara, J.-C. Rivoal, W. D. Bragg, K. Banerjee, V. P. Safonov, V. P. Drachev, Z. C. Ying, A. K. Sarychev, and V. M. Shalaev, “Near-field optical studies of semicontinuous metal films,” Phys. Rev. B **64**(16), 165403 (2001). [CrossRef]

**6. **S. I. Bozhevolnyi and V. Coello, “Statistics of local field intensity enhancements at nanostructured surfaces investigated with a near-field optical microscope,” Phys. Rev. B **64**(11), 115414 (2001). [CrossRef]

**7. **M. I. Stockman, S. V. Faleev, and D. J. Bergman, “Localization versus delocalization of surface plasmons in nanosystems: can one state have both characteristics?” Phys. Rev. Lett. **87**(16), 167401 (2001). [CrossRef] [PubMed]

**8. **S. Buil, J. Aubineau, J. Laverdant, and X. Quélin, “Local field intensity enhancements on gold semicontinuous films investigated with an aperture near-field optical microscope in collection mode,” J. Appl. Phys. **100**(6), 063530 (2006). [CrossRef]

**9. **K. Seal, A. K. Sarychev, H. Noh, D. A. Genov, A. Yamilov, V. M. Shalaev, Z. C. Ying, and H. Cao, “Near-field intensity correlations in semi-continuous metal-dielectric films,” Phys. Rev. Lett. **94**(22), 226101 (2005). [CrossRef] [PubMed]

**10. **V. M. Shalaev, *Optical Properties of Nanostructured Random Media,* (Springer, 2002).

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

**12. **D. P. Tsai, J. Kovacs, Z. Wang, M. Moskovits, V. M. Shalaev, J. S. Suh, and R. Botet, “Photon scanning tunneling microscopy images of optical excitations of fractal metal colloid clusters,” Phys. Rev. Lett. **72**(26), 4149–4152 (1994). [CrossRef] [PubMed]

**13. **S. I. Bozhevolnyi, V. A. Markel, V. Coello, W. Kim, and V. M. Shalaev, “Direct observation of localized dipolar excitations on rough nanostructured surfaces,” Phys. Rev. B **58**(17), 11441–11448 (1998). [CrossRef]

**14. **R. C. Word, J. Fitzgerald, and R. Könenkamp, “Photoelectron emission control with polarized light in plasmonics metal random structure,” Appl. Phys. Lett. **99**(4), 041106 (2011). [CrossRef]

**15. **C. Awada, G. Barbillon, F. Charra, L. Douillard, and J.-J. Greffet, “Experimental study of hot spots in gold/glass nanocomposite films by photoemission electron microscopy,” Phys. Rev. B **85**(4), 045438 (2012). [CrossRef]

**16. **V. Krachmalnicoff, E. Castanié, Y. De Wilde, and R. Carminati, “Fluctuations of the local density of states probe localized surface plasmons on disordered metal films,” Phys. Rev. Lett. **105**(18), 183901 (2010). [CrossRef] [PubMed]

**17. **P. Sheng, “Theory of the dielectric function of granular composite media,” Phys. Rev. Lett. **45**(1), 60–63 (1980). [CrossRef]

**18. **V. M. Shalaev and A. Sarychev, “Nonlinear optics of random metal-dielectric films,” Phys. Rev. B **57**(20), 13265–13288 (1998). [CrossRef]

**19. **V. M. Shalaev, *Nonlinear Optics of Random Media,* (Springer, 2000).

**20. **C. H. Chan, S. H. Lou, L. Tsang, and J. A. Kong, “Electromagnetic scattering of waves by random rough surface: a finite-difference time-domain approach,” Microw. Opt. Technol. Lett. **4**(9), 355–359 (1991). [CrossRef]

**21. **U. K. Chettiar, P. Nyga, M. D. Thoreson, A. V. Kildishev, V. P. Drachev, and V. M. Shalaev, “FDTD modeling of realistic semicontinuous metal films,” Appl. Phys. B **100**(1), 159–168 (2010). [CrossRef]

**22. **M. D. Thoreson, J. Fang, A. V. Kildishev, L. J. Prokopeva, P. Nyga, U. K. Chettiar, V. M. Shalaev, and V. P. Drachev, “Fabrication and realistic modeling of three-dimensional metal-dielectric composites,” J. Nanophotonics **5**(1), 051513 (2011). [CrossRef]

**23. **V. A. Markel, V. M. Shalaev, E. B. Stechel, W. Kim, and R. L. Armstrong, “Small-particle composites. I. linear optical properties,” Phys. Rev. B Condens. Matter **53**(5), 2425–2436 (1996). [CrossRef] [PubMed]

**24. **J. Aubineau, “Modélisation de couches minces métalliques fractales et calcul d’exaltations de champs électromagnétiques, » PhD thesis, Université de Versailles Saint-Quentin, 2005.

**25. **M. I. Stockman, L. N. Pandey, L. S. Muratov, and T. F. George, “Giant fluctuations of local optical fields in fractal clusters,” Phys. Rev. Lett. **72**(15), 2486–2489 (1994). [CrossRef] [PubMed]

**26. **J. A. Sánchez-Gil and J. V. Garcia-Ramos, “Calculations of the direct electromagnetic enhancement in surface enhanced Raman scattering on random self-affine fractal metal surfaces,” J. Chem. Phys. **108**(1), 317–325 (1998). [CrossRef]

**27. **I. Horcas, R. Fernández, J. M. Gómez-Rodríguez, J. Colchero, J. Gómez-Herrero, and A. M. Baro, “WSXM: a software for scanning probe microscopy and a tool for nanotechnology,” Rev. Sci. Instrum. **78**(1), 013705 (2007). [CrossRef] [PubMed]

**28. **P. B. Johnson and R. W. Christy, “Optical constants of the noble metals,” Phys. Rev. B **6**(12), 4370–4379 (1972). [CrossRef]

**29. **K. S. Yee, “Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media,” IEEE Trans. Antenn. Propag. **14**(3), 302–307 (1966). [CrossRef]

**30. **A. Taflove and S. C. Hagness, *Computational Electrodynamics: The Finite-Difference Time-Domain Method* (Artech, 2005).

**31. **D. A. Genov, A. K. Sarychev, and V. M. Shalaev, “Plasmon localization and local field distribution in metal-dielectric films,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. **67**(5), 056611 (2003). [CrossRef] [PubMed]

**32. **H. Cang, A. Labno, C. Lu, X. Yin, M. Liu, C. Gladden, Y. Liu, and X. Zhang, “Probing the electromagnetic field of a 15-nanometre hotspot by single molecule imaging,” Nature **469**(7330), 385–388 (2011). [CrossRef] [PubMed]

**33. **H. Xu and M. Käll, “Polarization-dependent surface-enhanced Raman spectroscopy of isolated silver nanoaggregates,” ChemPhysChem **4**(9), 1001–1005 (2003). [CrossRef] [PubMed]

**34. **H. Wei, F. Hao, Y. Huang, W. Wang, P. Nordlander, and H. Xu, “Polarization dependence of surface-enhanced Raman scattering in gold nanoparticle-nanowire systems,” Nano Lett. **8**(8), 2497–2502 (2008). [CrossRef] [PubMed]