Particle simulation has been widely used in studying plasmas. The technique follows the motion of a large assembly of charged particles in their self-consistent electric and magnetic fields. Plasmons, collective oscillations of the free electrons in conducting media such as metals, are connected to plasmas by very similar physics, in particular, the notion of collective charge oscillations. In many cases of interest, plasmons are theoretically characterized by solving the classical Maxwell’s equations, where the electromagnetic responses can be described by bulk permittivity. That approach pays more attention to fields rather than motion of electrons. In this work, however, we apply the particle simulation method to model the kinetics of plasmons, by updating both particle position and momentum (Newton–Lorentz equation) and electromagnetic fields (Ampere and Faraday laws) that are connected by current. Particle simulation of plasmons can offer insights and information that supplement those gained by traditional experimental and theoretical approaches. Specifically, we present two case studies to show its capabilities of modeling single-electron excitation of plasmons, tracing instantaneous movements of electrons to elucidate the physical dynamics of plasmons, and revealing electron spill-out effects of ultrasmall nanoparticles approaching the quantum limit. These preliminary demonstrations open the door to realistic particle simulations of plasmons.
Plasmons, collective oscillations of a free electron gas, are usually theoretically investigated by electromagnetic field simulations through solving Maxwell’s equations. The motions of electrons are indirectly represented by the bulk permittivity of metal. Such treatment overlooks real-time motions of free electrons, especially those electrons in the vicinity of the metal surface. For example, as the size of metal nanostructures shrinks below 10 nm, the motions of the surface electrons start to play an essential role, leading to electron spill-out and quantum tunnelling effects , . Electromagnetic field simulations often introduce a fictitious layer with modified permittivity to represent these surface electrons for various effects, such as non-local responses ,  and charge transfer plasmons , . Although computationally efficient, these treatments might oversimplify the underlying physics and sometimes rely on fitting parameters to describe certain effects. Electromagnetic field simulations also encounter difficulties in modelling the dynamic interaction between a single-electron source and the sample in a typical electron energy loss spectroscopy (EELS) experiment . For example, the electron source is usually modeled as an infinitely long line-like electron beam with zero radius using the boundary element method (BEM) , or an infinitely long broadband current density with a radius of 1 nm using the finite element method (FEM) , which are both in the frequency domain and do not capture the dynamics when a single electron is passing by the sample. Even if the finite-difference time-domain (FDTD) method is applied, the electron source is approximated by a series of optical dipoles , which might not encapsulate all the physics of the electrons in motion at half the speed of light. A microscopic computational tool that traces the freely moving electrons with relatively low computational cost is thus desired.
In fact, plasmons share similar physics and equations of motion as plasmas, one of the four fundamental states of matter. A plasma consists of a gas of charged ions and free electrons with motion–dynamics that are easily affected by the external electromagnetic fields. Among many modelling methods, particle simulation, or the particle-in-cell (PIC) method , , has become the most powerful simulation tool in plasma community since its development more than half a century ago , . The PIC method is a kinetic simulation method, which tracks the motion of a large assembly of charged particles in their self-consistent electric and magnetic fields, and therefore works particularly well for ultrafast (fs-scale) and ultrasmall (nm-scale) problems. It introduces the concept of a macroparticle  (also referred to as superparticle, or finite-size particle), a computational particle that represents many real electrons or ions. This key feature makes particle simulation efficient and successful in plasma studies, especially for practical plasma systems with many degrees of freedom, for example, experimental astrophysics , plasma-based accelerators , , , inertial confinement fusion , attosecond pulse generation , and THz radiation sources , , , .
Inspired by plasma modelling that properly tracks the motions of free electrons, we start our journey with the application of particle simulation to model plasmons. After carefully examining the essential differences between field and particle simulations, we identify two classic plasmon problems to demonstrate the applicability of particle simulation. In the first example, we simulate single-electron excitation of plasmons on a gold (Au) nanoribbon in a typical EELS experimental setup ; and in particular, we analyze the electron dynamics. In the second example, we investigate the electron spill-out effects in sodium (Na) nanowires  by making use of the soft boundary condition of particle simulation. Finally, we list the limitations and possible extensions of our current model. Going further, we are interested in applying particle simulation to solve problems such as: how many free electrons indeed participate in the plasmon oscillation, whether surface electrons escape from the metal surface during a plasmon oscillation, i. e., how plasmonic hot-electrons are generated, nonlinear behavior, as well as stochastic and transport phenomena in plasmonics.
2 Particle-in-cell (PIC) simulation
Particle simulation of plasmons is substantially different from classical electromagnetic field simulation methods, as illustrated in Figure 1. The link between them can be established by the free electron model (or plasma model)  that was developed in 1927 by Arnold Sommerfeld, where a gas of free electrons of number density n moves against a background of positive ion cores. The details of the lattice potential and electron interactions are not taken into account. Instead, one simply assumes that some aspects of the band structure are incorporated into the effective optical mass m of each electron. The electrons oscillate in response to the applied electromagnetic field, and their motion is damped via collisions occurring with a characteristic collision frequency γ. The dielectric function of the free electron gas is therefore described by :
where is the plasma frequency of the free electron gas. In contrast to field simulations (Figure 1A), particle simulation does not rely on a dielectric function to describe metals. Instead, the number density n is used to randomly distribute free electrons (represented by the yellow spheres) within the metal region as shown in Figure 1B. This implementation could provide an intuitive physical picture of plasmons.
Particle simulation (or the PIC method) simultaneously solves both the Maxwell’s equations:
and the Newton–Lorentz equation:
Here, , , , and are the relativistic factor, mass, charge, and velocity of the particle α (e. g., free electrons, valence electrons, or ions). Similar to the conventional FDTD field simulation method, the electromagnetic fields E, B, and current density J in the PIC method are defined on discretized mesh nodes in space, as shown in Figure 1. On the other hand, and refer to the local electromagnetic fields at the particle position, which are interpolated from neighboring mesh nodes in Figure 1B.
Every particle moves and constitutes a local current density (note: is relevant to the concept of macroparticle in the PIC method and refers to the number of elementary particles each macroparticle represents), which is interpolated and summed up at the neighboring mesh nodes , :
Here, W (x) represents a shape function of the macroparticles , which are similarly used in the calculation of and by interpolating E and B from neighbouring mesh nodes to each particle. This Ji will be taken into the electromagnetic field solver in Eq. (2) to update the electric field Ei and magnetic field Bi at each mesh node. Then the updated fields are used to re-calculate the position and momentum of each particle. The iterative process continues, and we can record particles and fields at every instant. In this way, the dynamics of plasmon oscillations can be modelled. In this work, we employ EPOCH  and SMILEI  particle-in-cell codes.
Compared to field simulation, the extra component accessible with particle simulation is the ability to trace the particles, which leads to the second key feature, i. e., soft boundary conditions and auto-inclusion of non-local effects. As shown in Figure 1B, the electromagnetic fields are defined at fixed mesh nodes in PIC, similar to the conventional FDTD method, leading to a hard-wall boundary for the solutions of the fields (represented by the cloud shape around the medium). However, the positions of the particles are not necessarily at the mesh nodes, and the particles are allowed to move anywhere including escaping from the boundary around the medium as shown in Figure 1B. In other words, particle simulation has an intrinsic soft boundary condition for particle motion, which allows us to explore the electron spill-out in quantum plasmonics ,  or other optical charge transport problems , , . On the other hand, the motion of the electrons is related not only to the field applied at a local position, but also depends on the fields at other positions. This defines the non-local effect , . By self-consistently solving the fields and tracing the particles, our PIC method automatically takes into consideration the non-local effects.
3 Results and discussions
In this section, we will elaborate two case studies by particle simulations. The first example, single-electron excitation of plasmons in EELS, will evidently illustrate free electron tracing in particle simulation. It enables us not only to numerically excite plasmons by a single fast electron, but also to efficiently capture the dynamics starting from the very first oscillating event to the steady-state plasmon oscillation. In the second example, leveraging on the implementation of the soft boundary condition, we revisit the famous problem of studying electron spill-out effects of quantum-sized nanoparticles. The results will be benchmarked to hydrodynamic models.
3.1 Single-electron excitation in EELS
In this first case study, we focus on the excitation of plasmons in EELS experiments, where a sample is exposed to a beam of electrons with a known, narrow range of kinetic energies. To model the excitation process, the beam of electrons is often assumed as a current source term in conventional electromagnetic field simulations , . Microscopically, this ignores the kinetic nature of the focused electron beam. In scanning transmission electron microscopy (STEM), only < 1% of the electrons emitted by the electron gun are eventually focused onto the sample, with most of them being filtered out to form a monochromatic and focused electron beam. We could optimistically assume that the electron beam generates a current of 0.1 nA on the sample (note: more realistic values could be only a tenth of this number), which is equivalent to 1.5 × 108 electrons per second. Considering a typical lifetime of 100 fs for a surface plasmon resonance , only 1.5 × 10−5 electrons on average will fly past the sample before the plasmons dampen-out entirely. We will need an unattainable electron flux 105 times higher to get more than one electron on average to arrive at the sample within the plasmon decay time. Therefore, in an EELS experiment, we almost exclusively generate independent single-electron excitation events that occur sequentially. Here, with the PIC method, we are able to realistically simulate single-electron excitation of plasmons.
We study both bright and dark plasmon modes in an Au nanoribbon with dimensions of 4 × 840 × 20 nm by using single-electron excitation in a setup similar to our previous EELS experiment . As schematically depicted in Figure 2A, dark and bright plasmon modes have zero and nonzero net dipole moments, respectively. Owing to their vanishing dipole moment, dark modes have longer lifetimes and also suffer from suppressed radiative losses , , making them more efficient for storing electromagnetic energy than bright modes. These attractive properties have made them ideal for applications such as nano-scale sensing , and lossless nano-scale waveguides . While bright modes can directly be excited by illuminating the nanoribbon with linearly polarized plane waves (as we later will show), dark plasmonic modes require more intricate excitation schemes  such as spatially nonuniform fields , evanescent excitation , oblique illumination , or spatial phase shaping .
In our PIC simulation of the EELS setup, the electrons within the Au nanoribbon were initialized with a Maxwell-Boltzmann energy distribution, with mean energy of 0.0375 eV (i. e., corresponding to room temperature T = 300 K with Boltzmann constant kB), and a density of n = 5.9 × 1028 m−3. As the Au ions are much more massive than the electrons, they do not contribute significantly to the electron dynamics on femtosecond time scales. This allows us to assume that they are fixed in space only as an ionic background to keep the charge neutrality throughout our PIC simulation. We used a single macroparticle of total charge −500e (e is the elementary charge) to excite the Au nanoribbon in order to obtain clearer spectra than an actual single-electron excitation. The macroparticle has a velocity of vx = c/2, where c is the speed of light in free space, and passes 1 nm away from the edge of the nanoribbon.
Figure 2B shows the computed spectra of the plasmon magnetic field Bx, calculated by taking the Fourier transform of the time-domain intensity profile over the entire simulation duration of 100 fs with a time interval of 0.15 fs for excitation locations a (red) and d (blue) using our PIC simulation. The excitation locations and geometry of the nanoribbon are depicted in the inset panel. As expected, we observe that excitation at the edge (location a) results in the formation of both bright and dark modes (indicated by yellow and gray strips respectively) due to the asymmetry of the excitation location with respect to the nanoribbon. Conversely, the excitation at the center (location d) only results in dark modes being formed due to the inherently symmetrical excitation geometry. By plotting the spatial intensity of the plasmon modes, we are able to identify whether each peak corresponds to a dark or bright mode by checking whether the number of high-intensity regions across the Au nanoribbon is odd (dark) or even (bright). The color maps in Figure 2C show the intensity as a function of space for the dark (left panels, location d) and bright (right panels, location a) plasmon modes corresponding to several of the spectral peaks shown in Figure 2B. The cyan arrows in Figure 2C indicate the location of excitation.
To further analyze the electron dynamics within the Au nanoribbon during the plasmonic oscillations, we irradiate the Au nanoribbon with a continuous electromagnetic plane wave centered at 0.78 eV excitation energy (shown in Figure 3A), which corresponds to one of the bright plasmon modes (see Figure 2C). The resulting plasmon field Ex, as in Figure 3B, shows the gradual build up of the mode in about 15 fs (i. e., from 73 to 88 fs), until a steady-state oscillation is reached. Here, the time to establish the mode is about 3 times the plasmon mode cycle.
One benefit of PIC simulations of plasmons is the ability to track individual particles; this readily enables the study of electron dynamics of plasmons. We show this in Figure 3C (scatter-plot) and 3D (phase-plot), which depict the steady-state dynamics of the electrons in the Au nanoribbon at five typical times within a single plasmon cycle (cyan dashes in Figure 3A; white dashes in Figure 3B). The temporary local electron densities (moving upwards and downwards) vividly exhibit that electron bunching occurs to establish the plasmon mode, at a frequency of every half cycle. This bunching is strongest while the plasmon field is weakest, as seen at two instances (ii) and (iv) in the figure, suggesting a quarter-cycle time-delay between electron bunching and the plasmon field. It should be noted that the plasmon field data plotted in Figure 3B was collected 5 nm away from the edge of the nanoribbon (at the side of x < 0). For visibility, only electrons with kinetic energy greater than four times the mean value were plotted in Figure 3C, D. The electrons in the Au nanoribbon were initialized with the same parameters as in Figure 2.
The results presented in this section indicate that with regards to plasmon simulation, the PIC method enables us to not only to reproduce results of conventional field simulations, but also allows us a direct view of the actual electron oscillations which underlie the formation of plasmons. In addition to modelling single-electron excitations in EELS experiments, particle simulation could be potentially used to model other electronic excitations, for example, tunneling-electrons excitation of plasmons .
3.2 Revealing electron spill-out effects
In the second case study, we leverage on the soft boundary conditions of particle simulation to study the spill-out effect of surface electrons in mesoscopic-scale metal nanoparticles. The term “mesoscopic” refers to deep nanoscale, which lies between the granular microscopic (electronic-scale) and continuous macroscopic (wavelength-scale) regimes . Solving mesoscopic electromagnetic problems has long relied on the well-developed theoretical approaches in either macroscopic or microscopic regimes. However, electromagnetic field-simulation based macroscopic approaches often deal with individual aspects of the omissions, such as nonlocality  and free-electron spill-out . Microscopic first-principles approaches, e. g., time-dependent density functional theory , are severely constrained by computational demands, thus impractical for multiscale problems. Recently, a general and unified theoretical framework has been developed by introducing a set of mesoscopic boundary conditions to the conventional macroscopic Maxwell’s equations . Motivated by the same objective, our PIC method is an alternative approach starting from a microscopic viewpoint, but is much more computationally efficient than other first-principles approaches by emphasizing only the motions of free electrons.
In an attempt to use a plasma model to simulate plasmons in a solid-state material, there are a few things to note. In plasma simulations, the electrons are separated from the ions and totally free in space, and they are considered all equivalent following classical Maxwell-Boltzmann statistics. But in solid-state physics, the electrons follow Fermi–Dirac statistics, in which each quantum state can only be occupied by a single electron. More importantly, these Fermi electrons are only free to move within the material, as they are still loosely bound to nuclei, therefore energy is required for them to escape from the material. In PIC simulation, an external potential-wall can be set up at the material boundary to represent such an energy barrier. Depending on the height of this potential-wall, the electrons can be partially or totally confined within the material. With the freedom to set the potential-wall or leave it open, we designate such boundary as a soft boundary.
As an example, we apply the PIC method to study the electron spill-out in a sodium (Na) nanowire, which is commonly used to differentiate the standard hydrodynamic model with a hard-wall boundary condition (HW-HDM)  and the self-consistent hydrodynamic method (SC-HDM) . The major difference between these two models is that the electrons in HW-HDM are strictly confined in the metallic structure, while SC-HDM generalizes standard hydrodynamic theory to include free-electron spill-out at metal–dielectric interfaces. As a result, SC-HDM predicts red-shifted resonances of infinitely long Na nanowires with 1/R scaling (R: nanowire radius), whereas HW-HDM incorrectly predicts a blue-shift . To elucidate the effect of electron spill-out, we implement either a potential-wall or an open boundary condition in the same PIC framework, corresponding to HW-HDM and SC-HDM, respectively. The potential wall is created by applying an external static electric-field of 117 V/nm over a 2 nm shell surrounding the nanowire. This number is chosen merely to represent a very large potential wall, which is equivalent to 100 times the work function of Na over 2 nm thickness.
As evidently indicated in Figure 4A, on the same basis of Fermi–Dirac statistics, the open boundary PIC simulations predict a red-shift on the 1/R scaling (red line), whereas the potential-wall PIC simulations have an opposite blue-shift prediction (black line). Clearly, it is the spill-out of electrons that has caused the red-shift, but only for Na nanowires of radius R smaller than a critical value 5 nm. Beyond this critical point (i. e., R > Rc), a substantial deviation from the 1/R dependence is observed, which can be attributed to retardation . The results agreed reasonably well with the SC-HDM model (dashed line) . Unlike these results from Fermi–Dirac statistics, Maxwell-Boltzmann statistics, however, results in no significant spectral shift below the retardation limit (gray line). This is probably due to the underestimated kinetic energy of electrons from Maxwell-Boltzmann statistics.
Besides predicting the resonant surface plasmon energies, we are also able to compute directly from the PIC model the space distribution of free electron densities, similar to other microscopic theoretical approaches. Taking a particular Na nanowire with R = 2 nm as an example, in Figure 4B, we plot the steady-state free electron density around the edge of the nanowire, with open and potential-wall boundary conditions, after 10 fs from the start of the simulation (note: this steady state is self-evolved and developed from an initialized uniform distribution of electrons). From Figure 4B, the spill-out of the free electrons is more prominently present in the simulation with open boundary. Nevertheless, it is questionable that quite a significant number of electrons are still observable outside the material boundary for the potential-wall boundary condition. This is due to the interpolation algorithm in PIC simulation. To better distinguish the electron spill-out for the two boundary conditions, we also plot the induced charge density (i. e., the difference between ion and electron density)  in the inset of Figure 4B. For open boundary (red curve), “+” (or “−”) charge accumulation appears just inside (or outside) the edge of the Na nanowire, indicating the electron spill-out. In contrast, “−” charge accumulation occurs only within the nanowire for the potential-wall boundary (black curve), representing electron spill-in.
More vividly, for the same nanowire, we show the real-time electron movements in the x–y plane near the upper edge of the nanowire at four typical times within a single plasmon cycle T0 in Figure 4C. With red (blue) dots representing electrons moving upwards (downwards), it is pictorially evident that electron bunching occurs periodically, forming the plasmon dipole mode (see t = 0 with most red dots and t = 0.5T0 with most blue dots), for both potential-wall and open boundary conditions. However, the potential-wall boundary restricts the spill-out of the electrons in the vicinity of the Na surface, whereas the open boundary indeed allows some electrons to get out, because their movements are purely determined by the local electromagnetic fields.
The results presented in this section give a demonstration of an important feature of particle simulation, i. e., the soft boundary condition, where we have the flexibility to set the boundary to be open or a potential-wall. This allows us to capture the microscopic dynamics of electrons at the material interface (in this example, the electrons are spilling out of Na into the vacuum) so as to study the effects of such dynamics. We could further modify the settings of our potential-wall boundary condition to include more physics or even some quantum effects. By doing so, the potential of particle simulation is unlimited, e. g., an extension to Ag nanowires by including inter-band transitions (see below), bimetallic nanoplasmonic structures, electron spill-out and charge transport between nearly touching nanoparticles.
4 Limitations and extensions
4.1 Inter-band transitions
Over a wide frequency range, the optical properties of metals can be explained by the free electron model . To be more precise, the free-electron description adequately describes the optical response of metals only for photon energies below the threshold of transitions between electron bands. It breaks down when inter-band transitions occur. For alkali metals, this range extends up to the ultraviolet, while for noble metals inter-band transitions occur at visible and near-infrared frequencies, limiting the validity of this approach. Above their respective band edge thresholds, photons are very efficient in inducing inter-band transitions, where electrons from the filled band below the Fermi surface are excited to higher bands. The main consequence of these processes concerning plasmons is an increased damping and competition between the two excitations at visible frequencies. At the moment, the contribution from valence electrons is neglected in PIC simulations, while the ions (composed of the nucleus, core and valence electrons) only keep the charge neutrality in the metal. To take into account the screening effect of the valence electrons governed by the inter-band transitions, another species of particles with larger effective mass can be introduced into our PIC simulation to represent valence electrons. To simplify, as these heavier valence electrons are moving closely around the nucleus, they effectively cause a highly polarized environment, which can be modelled by a background  in future PIC simulations.
4.2 Damping via collisions
The electrons oscillate in response to the applied electromagnetic field, and their motion is damped via collisions occurring with a characteristic collision frequency . Here, τ is known as the relaxation time of the free electron gas, which is typically on the order of 10−14 s at room temperature, corresponding to γ = 100 THz . The results shown in Section 3.1 are based on the assumption that there are no collisions among the electrons. When the collision effect is included (about 4× longer computation time), plasmon resonance frequencies barely change, but a decrease in the amplitudes of the spectral peaks would be observed. The collisions between electrons should therefore be included in PIC simulations for studies concerning damping in plasmons.
4.3 Fermi-electrons in metals
To model solid metals using plasma-physics based PIC, that is normally used to model quantum plasmas , the quantum effects become important when the temperature T is lower than the Fermi temperature TF, hence the relevant statistical distribution changes from Maxwell–Boltzmann to Fermi–Dirac. Besides, energy (i. e., work function) is required to free the Fermi-electrons from the metal. It is possible to implement such an energy barrier by setting up a potential-wall at the surface of the metal, as illustrated in our electron spill-out study. So far, we have only tried a very basic implementation of a sharp potential-wall by adding an external static electric-field over an extremely short distance. This electric field will simply push away those electrons approaching the wall. To further improve the potential-wall implementation, we could consider a more accurate profile of the potential barrier and the quantum tunneling probability for each electron , . Depending on their kinetic energies, some electrons will be reflected back to the metal, while some will be allowed to transmit through the wall. This idea of a quantum potential-wall could be implemented in future.
4.4 Shape function
Non-locality in plasmonic structures is the result of electron-electron repulsion and consequently a wave-vector dependent response of the system. In PIC simulation, both factors are implicitly taken into account by charge deposition and field interpolation at mesh nodes through the shape function defined in Eq. (4) . This opens up an opportunity to use shape function as a fitting parameter for PIC simulation of non-local effects in plasmons. For example, different orders of shape function could be implemented to account for the range of electron-electron repulsion, which could affect the slope of the size-dependent curve in Figure 4A in a similar manner as the β factor of the hydrodynamic model . The shape function used in our current work is based on second order interpolation, where the function extends over a distance equal to three times the PIC cell size , which is set to be 0.1 nm. Higher orders will be more accurate but computationally more expensive . This topic can be further explored in future work to thoroughly interpret the PIC results, together with benchmarked data such as experiments or SC-HDM results.
5 Conclusions and outlook
We have demonstrated a new theoretical framework to model plasmons – particle simulation, which not only solves for electromagnetic fields, but also traces the motions of electrons. In doing so, our numerical experiment by particle simulation is able to provide a microscopic (electronic-scale) view of plasmons. In the first example of studying plasmons inside an Au nanoribbon, we have employed a single-electron source, mimicking realistic EELS experiments, to efficiently excite various plasmon modes in one go. Different excitation scenarios can be handled straightforwardly by our particle simulations. To analyze the electron dynamics for a particular plasmon mode, we have the flexibility to switch to plane-wave excitation at a central frequency and extract information of interest. For example, we have vividly illustrated the spatial bunching of free electrons with their corresponding velocities during a cycle of a plasmon oscillation, noticing that the plasmon field is a quarter-cycle delayed with respect to the electron bunching. Other diagnosis techniques could be applied to reveal different aspects of plasmons.
The other uniqueness of particle simulation lies in the soft boundary conditions. We have made use of this property to demonstrate the possibility of modelling electron spill-out effects in Na nanowires, benchmarked to the well-developed hydrodynamic models. It should be emphasized that particle simulation is able to differentiate between free and valence electrons, and can turn on and off an external potential-wall at the material boundary. This operational flexibility allows us to differentiate alkali and noble metals, and to identify regimes where electron spill-out or retardation effects dominate. Moving forward, our PIC model could be further developed by taking into consideration the inter-band transitions, the electron-electron collisions, and a more comprehensive description of Fermi-electrons in metals and beyond. This is still a largely unexplored field, and is thus likely to become an active area of future research.
Funding source: Science and Engineering Research Council
Award Identifier / Grant number: Young Individual Research Grants (YIRG No. A1784c0
Funding source: National Research Foundation Singapore
Award Identifier / Grant number: NRF2016-NRF-ANR002NRF2017-NRFNSFC002-015
Funding source: MOE Singapore
Award Identifier / Grant number: T2 grant (2018-T2-1-007) and MOE PhD RSSTier 1 (R-284-000-179-133).
Funding source: Agency for Science, Technology and Research
Funding source: National Research Foundation Singapore
Author contribution: All the authors have accepted responsibility for the entire content of this submitted manuscript and approved submission.
Employment or leadership: None declared.
Honorarium: None declared.
Conflict of interest statement: The authors declare no conflicts of interest regarding this article.
The Institute of High Performance Computing (IHPC) acknowledges financial support from A*STAR SERC Young Individual Research Grants (YIRG No. A1784c0020), the National Research Foundation Singapore (NRF2017-NRF-NSFC002-015 and NRF2016-NRF-ANR002). J. Z. J. L. and L. K. A. would like to acknowledge the support of the MOE T2 grant (2018-T2-1-007), USA ONRG grant (N62909-19-1-2047) and MOE PhD RSS. M. B. acknowledges support from the Singapore Ministry of Educationʼs Academic Research Fund Tier 1 (R-284-000-179-133). The computational work for this article was partially performed on resources of the National Supercomputing Centre, Singapore (https://www.nscc.sg).
 M. S. Tame, K. R. McEnery, Ş. K. Özdemir, J. Lee, S. A. Maier, and M. S. Kim, “Quantum plasmonics,” Nat. Phys., vol. 9, no. 6, pp. 329–340, 2013, https://doi.org/10.1038/nphys2615.Search in Google Scholar
 D. Xu, X. Xiong, L. Wu, et al., “Quantum plasmonics: new opportunity in fundamental and applied photonics,” Adv. Opt. Photonics, vol. 10, no. 4, pp. 703–756, 2018, https://doi.org/10.1364/AOP.10.000703.Search in Google Scholar
 J. M. McMahon, S. K. Gray, and G. C. Schatz, “Nonlocal optical response of metal nanostructures with arbitrary shape,” Phys. Rev. Lett., vol. 103, no. 9, 2009, Art no. 97403, https://doi.org/10.1103/PhysRevLett.103.097403.Search in Google Scholar
 A. I. Fernández-Domínguez, A. Wiener, F. J. García-Vidal, S. A. Maier, and J. B. Pendry, “Transformation-optics description of nonlocal effects in plasmonic nanostructures,” Phys. Rev. Lett., vol. 108, no. 10, 2012, Art no. 106802, https://doi.org/10.1103/PhysRevLett.108.106802.Search in Google Scholar
 R. Esteban, A. G. Borisov, P. Nordlander, and J. Aizpurua, “Bridging quantum and classical plasmonics with a quantum-corrected model,” Nat. Commun., vol. 3, 2012, Art no. 825, https://doi.org/10.1038/ncomms1806.Search in Google Scholar
 L. Wu, S. F. Tan, M. Bosman, et al., “Charge transfser plasmon resonances across silver-molecule-silver junctions: estimating the terahertz conductance of molecules at near-infrared frequencies,” RSC Adv., vol. 6, no. 75, pp. 70884–70894, 2016, https://doi.org/10.1039/C6RA16826D.Search in Google Scholar
 P. Das, T. Kumar Chini, and P. James, “Probing higher order surface plasmon modes on individual truncated tetrahedral gold nanoparticle using cathodoluminescence imaging and spectroscopy combined with fdtd simulations,” J. Phys. Chem. C, vol. 116, no. 29, pp. 15610–15619, 2012, https://doi.org/10.1021/jp3047533.Search in Google Scholar
 Z. Fan, M. Bosman, X. Huang, et al., “Stabilization of 4H hexagonal phase in gold nanoribbons,” Nat. Commun., vol. 6, 2015, Art no. 7684, https://doi.org/10.1038/ncomms8684.Search in Google Scholar
 C. K. Birdsall, and A. B. London, Plasma Physics via Computer Simulation, Boca Raton, Taylor and Francis, 2004.Search in Google Scholar
 A. B. Langdon, and K. B. Charles, “Theory of plasma simulation using finite-size particles,” Phys. Fluids, vol. 13, no. 8, pp. 2115–2122, 1970, https://doi.org/10.1063/1.1693209.Search in Google Scholar
 S. Mondal, V. Narayanan, W. J. Ding, et al., “Direct observation of turbulent magnetic fields in hot, dense laser produced plasmas,” Proc. Natl. Acad. Sci. U.S.A., vol. 109, no. 21, pp. 8011–8015, 2012, https://doi.org/10.1073/pnas.1200753109.Search in Google Scholar
 E. Esarey, C. B. Schroeder, and W. P. Leemans, “Physics of laser-driven plasma-based electron accelerators” Rev. Modern Phys., vol. 81:1229–1285, 2009, https://doi.org/10.1103/RevModPhys.81.1229.Search in Google Scholar
 A. Macchi, M. Borghesi, and M. Passoni, “Ion acceleration by superintense laser-plasma interaction,” Rev. Modern Phys., vol. 85, pp. 751–793, 2013, https://doi.org/10.1103/RevModPhys.85.751.Search in Google Scholar
 X. Q. Yan, H. C. Wu, Z. M. Sheng, J. E. Chen, and J. Meyer-ter Vehn, “Self-organizing GeV, nanocoulomb, collimated proton beam from laser foil interaction at 7 × 1021 W/cm2,” Phys. Rev. Lett., vol. 103, 2009, Art no. 135001, https://doi.org/10.1103/PhysRevLett.103.135001.Search in Google Scholar
 K. Mima, A. E. Dangor, R. G. Evans, et al., “Fast heating of ultrahigh-density plasma as a step towards laser fusion ignition,” Nature, vol. 412, pp. 798–802, 2001, https://doi.org/10.1038/35090525.Search in Google Scholar
 C. Thaury, and F. Quéré, “High-order harmonic and attosecond pulse generation on plasma mirrors: basic mechanisms,” J. Phys. B, vol. 43, 2010, Art no. 213001, https://doi.org/10.1088/0953-4075/43/21/213001.Search in Google Scholar
 G. Liao, Y. Li, H. Liu, et al., “Multimillijoule coherent terahertz bursts from picosecond laser-irradiated metal foils,” Proc. Natl. Acad. Sci., vol. 116, no. 10, pp. 3994–3999, 2019, https://doi.org/10.1073/pnas.1815256116.Search in Google Scholar
 W. J. Ding, Z. M. Sheng, and W. S. Koh, “High-field half-cycle terahertz radiation from relativistic laser interaction with thin solid targets,” Appl. Phys. Lett., 103, no. 20, 2013, Art no. 204107, https://doi.org/10.1063/1.4831684.Search in Google Scholar
 Z.-M. Sheng, K. Mima, J. Zhang, and H. Sanuki, “Emission of electromagnetic pulses from laser wakefields through linear mode conversion,” Phys. Rev. Lett., vol. 94, 2005, Art no. 095003, https://doi.org/10.1103/PhysRevLett.94.095003.Search in Google Scholar
 S. Mondal, Q. Wei, W. J. Ding, et al., “Aligned copper nanorod arrays for highly efficient generation of intense ultra-broadband thz pulses,” Sci. Rep., vol. 7, 2017, Art no. 40058, https://doi.org/10.1038/srep40058.Search in Google Scholar
 G. Toscano, J. Straubel, K. Alexander, et al., “Resonance shifts and spill-out effects in self-consistent hydrodynamic nanoplasmonics,” Nat. Commun., vol. 6, 2015, Art no. 7132, https://doi.org/10.1038/ncomms8132.Search in Google Scholar
 T. Umeda, Y. Omura, T. Tominaga, and H. Matsumoto, “A new charge conservation method in electromagnetic particle-in-cell simulations,” Comput. Phys. Commun., vol. 156, no. 1, pp. 73–85, 2003, https://doi.org/10.1016/S0010-4655(03)00437-5.Search in Google Scholar
 T.Z. Esirkepov, “Exact charge conservation scheme for particle-in-cell simulation with an arbitrary form-factor,” Comput. Phys. Commun., vol. 135, no. 2, pp. 144–153, 2001, https://doi.org/10.1016/S0010-4655(00)00228-9.Search in Google Scholar
 E. Cormier-Michel, B. A. Shadwick, C. G. R. Geddes, E. Esarey, C. B. Schroeder, and W. P. Leemans, “Unphysical kinetic effects in particle-in-cell modeling of laser wakefield accelerators,” Phys. Rev. E, vol. 78, 2008, Art no. 016404, https://doi.org/10.1103/PhysRevE.78.016404.Search in Google Scholar
 T. D. Arber, K. Bennett, C. S. Brady, et al., “Contemporary particle-in-cell approach to laser-plasma modelling,” Plasma Phys. Control. Fusion, vol. 57, no. 11, 2015, Art no. 113001, https://doi.org/10.1088/0741-3335/57/11/113001.Search in Google Scholar
 J. Derouillat, A. Beck, F. Pérez, et al., “Smilei: a collaborative, opensource, multi-purpose particle-in-cell code for plasma simulation,” Comput. Phys. Commun., vol. 222, pp. 351–373, 2018, https://doi.org/10.1016/j.cpc.2017.09.024.Search in Google Scholar
 L. Wu, H. Duan, P. Bai, M. Bosman, K. W. Y. Joel, and E. Li, “Fowler-Nordheim tunneling induced charge transfer plasmons between nearly touching nanoparticles,” ACS Nano, vol. 7, no. 1, pp. 707–716, 2012, https://doi.org/10.1021/nn304970v.Search in Google Scholar
 Z. Zhu, S. Joshi, S. Grover, and G. Moddel, “Graphene geometric diodes for terahertz rectennas,” J. Phys. D Appl. Phys., vol. 46, no. 18, 2013, Art no. 185101, https://doi.org/10.1088/0022-3727/46/18/185101.Search in Google Scholar
 S. F. Tan, L. Wu, J. K. W. Yang, P. Bai, M. Bosman, and C. A. Nijhuis, “Quantum plasmon resonances controlled by molecular tunnel junctions,” Science, vol. 343, no. 6178, pp. 1496–1499, 2014, https://doi.org/10.1126/science.1248797.Search in Google Scholar
 S. Raza, I. B. Sergey, M. Wubs, and N. A. Mortensen, “Nonlocal optical response in metallic nanostructures,” J. Phys. Condensed Matter, vol. 27, no. 18, 2015, Art no. 183204, https://doi.org/10.1088/0953-8984/27/18/183204.Search in Google Scholar
 J. A. Scholl, A. L. Koh, and J. A. Dionne, “Quantum plasmon resonances of individual metallic nanoparticles,” Nature, vol. 483, no. 7390, pp. 421–427, 2012, https://doi.org/10.1038/nature10904.Search in Google Scholar
 M. Bosman, E. Ye, S. F. Tan, et al., “Surface plasmon damping quantified with an electron nanoprobe,” Sci. Rep., vol. 3, 2013, Art no. 1312, https://doi.org/10.1038/srep01312.Search in Google Scholar
 N. J. Halas, S. Lal, W.-S. Chang, S. Link, and P. Nordlander, “Plasmons in strongly coupled metallic nanostructures,” Chem. Rev., vol. 111, no. 6, pp. 3913–3961, 2011, https://doi.org/10.1021/cr200061k.Search in Google Scholar
 D. E. Gómez, Z. Q. Teo, M. Altissimo, T. J. Davis, S. Earl, and A. Roberts, “The dark side of plasmonics,” Nano Lett., vol. 13, no. 8, pp. 3722–3728, 2013, https://doi.org/10.1021/nl401656e.Search in Google Scholar
 S. Zhang, and H. Xu, “Tunable dark plasmons in a metallic nanocube dimer: toward ultimate sensitivity nanoplasmonic sensors,” Nanoscale, vol. 8, no. 28, pp. 13722–13729, 2016, https://doi.org/10.1039/C6NR03806A.Search in Google Scholar
 D. SolisJr, B. Willingham, L. N. Scott, et al., “Electromagnetic energy transport in nanoparticle chains via dark plasmon modes,” Nano Lett., vol. 12, no. 3, pp. 1349–1353, 2012, https://doi.org/10.1021/nl2039327.Search in Google Scholar
 J.-S. Huang, J. Kern, P. Geisler, et al., “Mode imaging and selection in strongly coupled nanoantennas,” Nano Lett., vol. 10, no. 6, pp. 2105–2110, 2010, https://doi.org/10.1021/nl100614p.Search in Google Scholar
 S.-C. Yang, H. Kobori, C.-L. He, et al., “Plasmon hybridization in individual gold nanocrystal dimers: direct observation of bright and dark modes,” Nano Lett., vol. 10, no. 2, pp. 632–637, 2010, https://doi.org/10.1021/nl903693v.Search in Google Scholar
 W. Zhou, and T. W. Odom, “Tunable subradiant lattice plasmons by out-ofplane dipolar interactions,” Nat. Nanotechnol., vol. 6, no. 7, pp. 423–427, 2011, https://doi.org/10.1038/nnano.2011.72.Search in Google Scholar
 G. Volpe, S. Cherukulappurath, R. J. Parramon, G. Molina-Terriza, and R. Quidant, “Controlling the optical near field of nanoantennas with spatial phase-shaped beams,” Nano Lett., vol. 9, no. 10, pp. 3608–3611, 2009, https://doi.org/10.1021/nl901821s.Search in Google Scholar
 W. Du, T. Wang, H.-S. Chu, et al., “On-chip molecular electronic plasmon sources based on self-assembled monolayer tunnel junctions,” Nat. Photonics, vol. 10, no. 4, pp. 274–280, 2016, https://doi.org/10.1038/nphoton.2016.43.Search in Google Scholar
 Y. Yang, D. Zhu, W. Yan, et al., “A general theoretical and experimental framework for nanoscale electromagnetism,” Nature, vol. 576:248–252, 2019, https://doi.org/10.1038/s41586-019-1803-1.Search in Google Scholar
 J. Zuloaga, E. Prodan, and P. Nordlander, “Quantum description of the plasmon resonances of a nanoparticle dimer,” Nano Lett., vol. 9, no. 2, pp. 887–891, 2009, https://doi.org/10.1021/nl803811g.Search in Google Scholar
 V. T. Tatiana, P. Nordlander, J. Aizpurua, and A. G. Borisov, “Quantum effects and nonlocality in strongly coupled plasmonic nanowire dimers,” Opt. Express, vol. 21, no. 22, pp. 27306–27325, 2013, https://doi.org/10.1364/OE.21.027306.Search in Google Scholar
© 2020 Wen Jun Ding et al., published by De Gruyter, Berlin/Boston
This work is licensed under the Creative Commons Attribution 4.0 International License.