1. Introduction

Localized surface plasmons (LSPs) are collective excitations of the free electrons in metal nanoparticles (NPs). At resonance, a resulting near-field enhancement enables light concentration well below the diffraction limit [1]. The LSP resonance frequency depends not only on the geometry of the NPs [2], but also on the dielectric environment [3] and the NP material [4]. Manipulation of these parameters allows tuning of the plasmon frequency from ultraviolet to mid-infrared [5, 6]. The unique properties of LSPs make plasmonic devices, based on arrays of metal NPs, relevant in a variety of applications. Examples include light concentrators [1,6], plasmonic wave guides [7], nano-scale polarization control [8], imaging below the diffraction limit [9], negative-refractive-index materials [10] and nanosensors [11].

A key element in the plasmonic-device design process is optical modeling. If both the NP dimensions and separations are larger than 10 nm, quantum effects can be neglected [12], and the problem is reduced to solving Maxwell’s equations. While the simple case of a spherical NP can be solved analytically using Mie theory [13], the general case of an arbitrary NP geometry can only be solved numerically. For sheets of periodically arranged NPs, an ensemble can be modeled by considering a unit cell domain to which boundary conditions that enforce the appropriate translational symmetry are applied. Non-periodic structures are intrinsically more diverse and can generally only be modeled by considering the full multi-scattering problem. A broad variety of non-periodic sheet-like structures have been studied including random planar composites [14], large disordered clusters [15] and deterministic aperiodic nanostrutures [16] with ensemble properties ranging from chaotic to tuneable. For spheres, the multi-scattering problem can be solved analytically using generalized Mie theory [17]. Non-spherical NPs can be treated using the T-matrix method [18], but convergence remains poor for non-spheroids [19]. NPs of arbitrary geometry in inhomogeneous environments can be simulated using the boundary element method (BEM) [20], the finite-difference time-domain method (FDTD) [21], or the frequency-domain finite-element method (FEM) [22]. All of these methods discretize Maxwell’s equations on a mesh, BEM across the involved surfaces and FEM/FDTD throughout the simulation volume. As the volume mesh must be fine enough to resolve the electric-field oscillation, the available memory on a modern computer limits the domain size for FEM/FDTD to the micrometer scale for visible wavelengths. In the case of BEM, the resulting system of equations is smaller, but typically dense, which causes the memory consumption to increase rapidly with increasing number and/or complexity of the NPs. The memory-imposed restrictions are usually not a problem for a unit cell representation, but it severely limits the applicability of these methods to NP ensembles. While small NP arrays have been simulated using FDTD [23], large arrays cannot be treated with currently available computational resources.

At sufficiently large inter-particle spacings, particle-particle interactions become negligible. In this quasi-single-particle case, the ensemble response can be predicted from single-particle simulations. A such one-particle model was employed by [4] and [6] which consider similarly sized nanodisk ensembles with center-to-center spacings of several diameters. While in the former case a very good fit with experimental data was obtained, important features in the extinction spectrum remained unaccounted for in the latter. In this paper a two-particle model is proposed which relaxes the assumption of no particle-particle interactions in the simplest possible way. While being sufficiently simple in concept to provide physical insight and in computational complexity to run on a modern personal computer, it captures the main deviations in the experiment [6] versus the one-particle model.

2. Methods

The field distribution in the vicinity of metal NPs can be obtained by solving Maxwell’s equations. For non-magnetic, isotropic, linear media in the absence of external charges and currents, the wave equation for a time-harmonic field is

In accordance with the samples considered in [6], the example structure is chosen as truncated gold nanocones of 50 nm height placed on a 100 nm TiO₂:Er³⁺ thin film on a SiO₂ substrate. The geometries for the one- and two-particle models are shown in Fig. 1. To emulate that the NPs are embedded in infinite space, all boundaries are truncated by perfectly matched layers (PMLs).

 

Fig. 1 Schematic illustration of the xz-plane at y = 0 for the (a) one- and (b) two-particle models. Origo is placed at the air/film interface at (a) the NP center and (b) midway between the NPs. The bottom NP diameter D and the two-particle center-to-center distance d are indicated. The NPs are modeled as truncated cones with top diameter 0.9D.

Download Full Size | PDF

The background field is normal incident along the positive z-axis, i.e. E_B lies in the xy-plane. Relative to the line connecting the NPs in the two-particle geometry, polarization states with the electric field parallel (‖) and orthogonal (⊥) can be defined. With the NPs located along the x-axis the ‖ and ⊥ polarization states correspond to x- and y-polarization, respectively.

3. Simulation results

The NP diameter was selected as D = 270 nm, one of the diameters for which a peak splitting was observed experimentally by [6]. Calculated field distributions for the two-particle geometry for d = 2D are shown in Fig. 2(a) and 2(b). While the near-field component is most intense at the NP edges along the axis of polarization, the far-field component is most pronounced orthogonal to it. For dimer-like two-NP systems with small interparticle spacings (d < 1.5D) the near field dominates, and the maximum coupling is obtained for ‖ polarization [26]. However, for d ≥ 2D as considered in this work, the far-field prevails and the strongest coupling is observed for ⊥ polarization. To emphasize the difference between the two models, all calculations in this section are shown for ⊥ polarization rather than a polarization average.

 

Fig. 2 (a,b) Simulated field distributions at z = 0 for (a) ‖ and (b) ⊥ polarization, D = 270 nm, d = 2D and λ = 1290 nm (single particle resonance). To allow visualization of the far-field component, the color scale has been truncated at 2.5 even though the local field enhancement in (b) exceeds 20 within 1 nm of the NP edges. (c) Extinction cross section σ_ext in units of the geometrical cross section σ_geo as a function of λ for one- (dashed lines) and two-particle (solid lines) models for different d. The curves for which a peak sharpening/splitting is observed are shown in blue/red. (d,e) Simulated phase distributions. Shadings indicate particle positions for which peak sharpening (white) and splitting (black) is observed. Shown are Δϕ in (d) the xy-plane at the NP center (z = −25 nm) for λ = 1290 nm and (e) along the line (y, z) = (0, −25 nm) for different values of λ. The two plots (d,e) coincide along the dashed, white line.

Download Full Size | PDF

Figure 2(c) shows calculated extinction spectra for selected values of d along with a single-particle calculation. Alongside, the phase difference between the background field and the scattered field emitted by a single particle, Δϕ, is plotted in Fig. 2(d) and 2(e). Intuitively, the strongest coupling is expected when the second particle is located at points of in-phase addition [27], i.e, along the Δϕ = 0 lines. When d is small, a blueshift in the extinction cross section (Fig. 2(c)) is observed, but as d increases the blueshift decreases. At d = 580 nm the two-particle peak aligns with the single-particle peak resulting in a peak sharpening. From Fig. 2(d) we see that the position d = 580 nm (white shading) agrees well with the first Δϕ = 0 line. Subsequently, the peak is redshifted and around d = 900 nm a new, blueshifted peak is formed while the amplitude of the redshifted peak decreases. At d = 1100 nm, the amplitudes of the two peaks become similar in magnitude and they appear as a split peak. Figure 2(e) illustrates that while the Δϕ = 0 condition is not satisfied at the single-particle resonance, it is so at higher/lower λ-values (intersections of Δϕ = 0 lines with the black, shaded column) which explains the observed peak splitting. Upon increasing d further, the peak movements repeat continuously. Spectra for the first two/three distances for which a peak sharpening/splitting is observed are included in Fig. 2(c). At large values of d, the two-particle model curves converge towards the single-particle curve apart from small wiggles.

4. Comparison to experimental results

To assess the applicability of the two-particle model for random NP arrays, the simulation results are compared with previously published experimental data [6]. The data are grouped in three series denoted S4k, S6k and S8k according to their NP densities of N = 4000, 6000, and 8000 NPs per 100 × 100 μm². The NPs where fabricated by electron-beam lithography with the pseudo-random distributions constructed using a random number generator. The resulting mean distance to the four nearest neighbors, a measure analog to the nearest-neighbor distance in a square lattice, appears normally distributed with a mean value very close the characteristic particle spacing $1/\sqrt{N}$ and a standard deviation around 10% of the mean value. For each sample, the extinction cross section was measured as a function of λ, as exemplified in Fig. 3(a). In all three cases, the main resonance peak is formed by two subpeaks. Together with the experimental data, simulated spectra calculated from the one- and two-particle models are shown. In the two-particle model d was chosen as $1/\sqrt{N}$. As the orientation of the electric field varies randomly for the nearest-neighbor pairs across the sample, polarization averaged calculations were performed.

 

Fig. 3 Comparison of one- (dashed line) and two-particle models (solid lines) with experimental data (symbols, connected by a thin line in (a)). (a) Extinction spectra for a selected diameter in each series, 345 nm for S4k (left), 315 nm for S6k (center) and 270 nm for S8k (right). (b) Fitted peak position(s) as a function of inverse diameter D for the S4k, S6k and S8k series. The selected diameters in (a) are marked by grey, vertical shadings. The inserts show SEM images of a sample from each series.

Download Full Size | PDF

While the one-particle model reproduces the overall peak location fairly well, the absence of additional features clearly demonstrates it being too simplistic. The two-particle model, on the other hand, is able to predict the main features in the spectrum including the positions of the subpeaks. However, a disagreement in scale remains. A possible explanation could be additional light trapping in the TiO₂:Er³⁺ thin film due to presence of metal NPs on the surface. The peak location(s) were estimated by fitting a (double) Gaussian function to the spectra for both experiment and calculation. Applying this procedure for all samples, Fig. 3(b) was obtained. For each data series, the peak location(s) are shown as a function of the inverse of the NP diameter. The one-particle model predicts a single curve independent of particle density, which approaches a straight line for large diameters [4]. Overall it fits the data reasonably well, but the peak splitting is not captured. On contrary, the two-particle model includes the necessary symmetry breaking to reproduce the peak splitting, thus providing a much better fit with the experimental data.

When NPs are deposited directly on a substrate, particle-particle couplings are strongly suppressed due to the lack of phase matching between waves propagating along each side of the interface. As reported by [28–30], the suppression can be lifted by applying an index-matching liquid as superstrate. Two-particle model simulations (not shown) successfully reproduce these observations. The stronger couplings observed in this work compared to [4], where similar NP distributions are considered, are due the presence of the TiO₂ thin film. As the scattered field is guided along the film, a significant particle-particle coupling is observed even in the absence of an index-matched superstrate. The reproduction of the guiding effect, in terms of a simulated non-negligible particle-particle coupling, demonstrates the applicability of the two-particle model to NP ensembles in inhomogeneous environments that can otherwise be challenging to simulate.

5. Summary and discussion

Particle-particle interactions have been investigated in a two-particle system for center-to-center spacings larger than two diameters. It was found that the coupling is dominated by the far field, and that it depends strongly on the dielectric environment. The applicability of a two-particle model for investigating particle-particle couplings in large, sparse arrays of randomly distributed nanoparticles in inhomogeneous dielectric environments was assessed by comparison to experimental data. While the agreement in absolute extinction cross section is not perfect, the two-particle model succeeds in predicting the main spectral features including the observed peak splitting. At the same time, the computational complexity remains low enough that calculations can be carried out on a personal computer.

By tuning the inter-particle spacing d, the particle-particle interactions can be exploited to selectively sharpen (or broaden) the plasmonic resonance of the ensemble. The effect is analogous to diffraction in periodic arrays, but weaker and less sensitive to incidence angle. Randomly distributed nanoparticles can be manufactured inexpensively on industrial scale by self-assembly processes, while periodic arrays typically require expensive, low-throughput processing. The combination of cost-effective assembly and low angular sensitivity make plasmonic devices based on random arrangements of nanoparticles favorable for many real life applications, e.g. solar cells. The two-particle model provides a valuable tool for tailoring the particle-particle interactions in such devices.