Unified treatment for photoluminescence and scattering of coupled metallic nanostructures: I. Two-body system

Photoluminescence (PL) and scattering phenomena from metallic nanostructures have been explained and understood by several point of views. One of them is based on the classic harmonic oscillator model, which describes PL of single mode. In this study, we continue to expand this classic model to a coupling case, which involves two oscillators that interact with each other together with the excitation electric field. The new generated modes due to the coupling are carefully analyzed, including their behaviors varying with the coupling coefficients in different cases. Furthermore, for practical purpose, PL spectra and white light scattering spectra of two individual metallic nanostructures are calculated as examples employing the model to verify its validity. This work would give a deeper understanding on coupling PL and scattering phenomena and is helpful to relative applications.

The origin of PL has been discussed in plenty of studies, with different explanations such as interband transitions enhanced by LSPR [16], microscopic explanation for enhanced PL from gold nanoparticles [17], classic oscillator model assisted with electron distributions for single mode emission [18], and non-equilibrium electron dynamics affecting PL of metallic nanostructures [19]. Nevertheless, the coupling PL phenomena are seldom investigated in theory. For example, Prodan et al present a molecular orbital theory to describe the coupling plasmon modes introduced by the metallic nanostructures of arbitrary shape [20]. Jain et al provide a semiempirical 'plasmon ruler equation' based on discrete dipole approximation simulation method to estimate the plasmon shifts as a function of the separation between the nanoparticles [21,22]. However, the developed models based on quantum theories are neither lack of details on the emission spectra, especially for PL, nor lack of intrinsic physical pictures. Besides, Davis and Gómez present a simple algebraic model for LSPRs including their interactions [23]. They consider that the interaction part is influenced by both the incident electric field and the one induced by other particles, and the latter is treated with dipole approximation. However, the induced electric field they use is from Coulomb's law, which is approximate, because the surface charges are moving rather than stationary. Coincidentally, a recent work of Lu et al presents the exotic coupling between plasmonic nanoparticles through geometric configurations, in which several configurations have been studied carefully and an improved method for the coupled systems is developed [24]. The plasmon-induced electric field Figure 1. Schematic of the coupling harmonic oscillator model. The electrons (grey, negative charged) oscillate collectively along x-axis near their equilibrium positions. The ions (white, positive charged) is at rest. r 0 is the distance between the two ions, while r is the distance between the two electrons. x 1 and x 2 are the displacements relative to equilibrium positions of each oscillator. The two oscillators both oscillate along x-axis when excited by the excitation light at the circular frequency of ω ex which is x-polarized. distribution in this improved method is based on the finite-difference time-domain (FDTD) simulations, which is more accurate than the one calculated from the dipole approximation. Though FDTD method is effective as an aid, it costs simulation time, which is not that convenient. Hence, a clear picture and an algebraic model for coupling PL as well as scattering spectra is required to be built up.
In this study, we present a practical model to give a deep understanding on PL and scattering from coupled metallic nanostructures, e.g., gold nanorods or nanospheres. This model is based on the classic harmonic oscillator model, considering two oscillators that interact with each other. We treat the interaction part, i.e., coupling coefficients, between the two in a non-phenomenological way. That is, the coupling coefficients are obtained from the intrinsic physics rather than just assuming as several parameters. The model shows reasonable results to explain PL and white light scattering spectra of coupled metallic nanostructures in different situations. This work would help to understand coupling PL and scattering phenomena in a classical way.

Model
Since there are plenty free electrons in the metallic nanostructure, and these electrons oscillate when excited by the external electric field, we treat the nanostructure as a resonator, the oscillators of which are the electrons. Due to the collectively oscillating, we can simplify the multiple electrons as only one electron. Consider that two metallic nanostructures treated as two oscillators are close to each other and are driven by the electric field of the excitation light. The schematic is shown in figure 1. In order to obtain the emission electric field from them, we need to find out the differential equations of them. Define x 1 (t) and x 2 (t) as the displacements relative to equilibrium positions of each oscillator, thusẋ 1 (t) andẋ 2 (t) the velocities, andẍ 1 (t) andẍ 2 (t) the accelerations. The equations should be in this form: Here, F 21 and F 12 are the interaction forces between oscillator 1 and oscillator 2, and m e is the mass of electron. E 1 and E 2 are the amplitudes of the excitation electric field at the positions of the two oscillators, and usually E 1 = E 2 = E 0 is a good approximation; q 1 and q 2 are the charge of oscillator 1 and oscillator 2, respectively; β 01 and β 02 represent the damping coefficients, and ω 01 and ω 02 represent the inherent circular frequencies. The next step is to find out the interaction parts of the equations. The electric field introduced by a moving charged particle is given by [25]: where q is the charge of the particle, ε 0 is the permittivity of vacuum, and u ≡ cr/r − v. Here, c is the velocity of light in vacuum, v and a are the velocity and the acceleration of the particle, respectively, and r is the displacement vector from the particle to field point. In our one-dimension case, when considering the electric field introduced by one oscillator acting on the other oscillator, the second part of equation (2) is zero due to the fact that u and a are parallel. Besides, we notice that the charged particle that moves is the electron while the positive ion is assumed to be at rest. Hence, the interacted electric field at one oscillator should be contributed to both positive charged ion and negative charged electron of the other oscillator. Therefore, the electric field can be written as: Here, we use the conditions v/c 1 and x/r 0 1 for approximation. Notice that E 21 is the electric field in oscillator 1 introduced by one pair of the electrons and ions in oscillator 2, and E 12 is the electric field in oscillator 2 introduced by one pair of the electrons and ions in oscillator 1. Hence, the interaction forces should be written as F 21 = −N 2 eE 21 and F 12 = −N 1 eE 12 , where N 1 and N 2 are the effective numbers of free electrons in oscillator 1 and oscillator 2, respectively. We define the coupling coefficients as: thus equation (1) can be written as: For simplicity, we assume N 1 = N 2 = N, thus γ 21 = γ 12 = γ, g 21 = g 12 = g, and if we define it results in a simple form for γ and g: Firstly, we consider the situation without coupling, i.e., γ = g = 0. In such conditions, equation (5) is degenerated into the simple form:ẍ The general solutions are: Here, ω c1 = ω 2 01 − β 2 01 and ω c2 = ω 2 02 − β 2 02 represent the resonant circular frequencies, respectively, which are different from the inherent ones (ω 01 , ω 02 ). The coefficients that represent amplitudes are omitted for the moment, which can be obtained with the initial conditions. The details of this kind of individual oscillator has been discussed carefully in our previous work [18].
Secondly, we start to consider the coupling situation without excitation light, i.e., C 1 = C 2 = 0. The equations are:ẍ To solve equation (9), we can assume that x 1 (t) = A exp(αt) and x 2 (t) = B exp(αt), and substitute them back into equation (9), thus obtaining: Obviously, to obtain non-zero solutions, α should satisfy: Notice that equation (11) has analytic solutions for α, marked as α 1 , α 2 , α 3 and α 4 . However, the expressions are so complex that we would not write in the text. Instead, to illustrate the physical significance of α, we rewrite it in this form: Here, ω 1 and ω 2 are the new generated resonant circular frequencies when the two oscillators couple. We can call them mode 1 and mode 2, respectively. In a more special case, i.e., β 01 = β 02 = β 0 , ω 01 = ω 02 = ω 0 , the solutions of equation (11) are expressed easily: Thirdly, notice that the particular solutions for equation (5) are x 1 (t) = exp(−iω ex t) and x 2 (t) = exp(−iω ex t) (amplitudes are omitted). Therefore, combining these particular solutions and the general ones (equation (8)), we obtain the total solutions of equation (5) in a symmetric form: where We emphasize here that equation (13) is just a special case for ω 1 and ω 2 , and the general case for them should satisfy equation (12). The initial conditions are where we assume that C 1 = C 2 = C 0 due to the subwavelength scale of the system. Hence, these coefficients are obtained as: This results in the fact that x 1 (t) = x 2 (t) = x(t).
At last, we deal with the far field radiation. For simplicity, we consider the electric field at the position d, where d is perpendicular to x-axis, and d = |d| is the distance between field point and the center of the two oscillators. The assumption of d r 0 is reasonable for far field radiation. Hence, the first part of equation (2) is ignored compared with the second part, thus giving the electric field introduced by oscillator 1 and oscillator 2 as: where D = Ne 2πε 0 c 2 d , and E far is x-polarized. The emission intensity in the frequency domain, i.e., emission spectrum, can be evaluated by [18]: where Re Q is the real part of Q , and Q = 1 t 0 t 0 0 Q dt is the time average of quantity Q. The calculated result is: where A j = A j Ω 2 j D for j = 1, 2, 3. Here, we ignore the cross terms in the calculation because the time average is zero when ω 1 = ω 2 .
As our previous work explains [18], the emission spectrum is separated into two parts, one is the inelastic part (I inela ) which corresponds to PL spectrum, and the other is the elastic part (I ela ) which corresponds to white light scattering spectrum. Rewrite equation (18) as: Therefore, the PL spectrum is given by equation (19a), i.e., while the white light scattering spectrum is given from equation (19b) as long as ω ex is substitute by ω: To show the coupling modes for PL more clearly and to understand PL phenomenon more easily, we do not consider the electron distributions here as before [18], which contributes mostly to the anti-Stokes part of PL spectra, unless otherwise specified, though this model would be more accuracy for PL when assisted with the electron distributions.

Results and discussions
After obtaining these formulas, we would analyze in details to understand them more deeply.
Start from the coupling coefficients, g and γ. Figure 2(a) shows g and γ varying with the distance r 0 , calculated from equation (6). It implicates that the coupling coefficients decrease with the increase of r 0 , and γ is smaller than g. When r 0 is small enough, the coupling coefficients get large. Since these two coefficients are both related to r 0 , we take one of them, i.e., g, as the coupling strength in the rest of this work. Figures 2(b)-(d) show the new generated resonant circular frequencies (ω 1 , ω 2 ) and damping coefficients (β 1 , β 2 ) in different cases of the coupled oscillators as a function of the coupling strength g, calculated from equations (11) and (12). The simplest one ( figure 2(b)) is that the two oscillators are the same. The two new modes split when coupling, and the splitting increases with the increase of g. Here, we generally call the increasing ω 'blue branch', and the decreasing ω 'red branch'. On the other hand, the two damping coefficients also splits, and one increases (corresponding blue branch), the other decreases (red branch). Notice that there is a cut-off coupling strength for the red branch at around g cut ≈ ω 0 . In figure 2(c), the situation is almost the same, i.e., the difference of ω and β between the two branches increase with the increase of g, and there also exists g cut . The difference between figures 2(b) and (c) is that, to obtain the same level of splitting, the former needs a smaller g than the latter does. That is, the former gets a better coupling efficient than the latter does. In figure 2(d), due to the fact that ω c1 = ω 2 01 − β 2 01 = 0.954ω 0 and ω c2 = ω 2 02 − β 2 02 = 0.995ω 0 , it results in ω c1 < ω c2 with a small difference. The difference of ω between the two branches (ω 1 − ω 2 ) increases from negative value to zero and then increase to positive value as the increase of g. On the other hand, the difference of β between the Furthermore, the effective free electrons number affects the splitting as shown in figure 3, giving three values of N as examples. Firstly, we find out the behaviors of ω 1 and ω 2 as g increases for each figure. In figure 3(a), ω 1 increases and then decreases, while ω 2 decreases, indicating that ω 1 has a maximal value. In figure 3(b), ω 1 decreases and then increases and finally decreases, while ω 2 increases and then decreases, indicating that both ω 1 and ω 2 have maximal values. In figure 3(c), ω 1 decreases, while ω 2 increases slightly and then decreases. In figure 3(d), ω 1 decreases, while ω 2 increases and then decreases, the curves of which almost coincide with each other at the range around g = 0.3ω 0 to g = 0.5ω 0 . In figure 3(e), the behaviors are similar to the ones in figure 3(c). In figure 3(f), the behaviors are similar to the ones in figure 3(d), but the two curves cross rather than coincide. Secondly, we find out the similar behaviors for these parameters in a general view. In all cases, there are cut-off coupling strengths for both modes, writing as g cut1 and g cut2 , at which ω 1 = 0 and ω 2 = 0, respectively. The differences are, for smaller N (10 8 or 10 9 ), g cut1 > g cut2 , while for larger N (10 10 ), g cut1 < g cut2 . As g increases, the splitting of damping coefficients β 1 and β 2 gets larger. Furthermore, another interesting result is that there is a point g 0 (shown with black cross circle) at which β 2 (g 0 ) = 0 for each case, and g 0 < g cut2 . This is different from the one in figure 2 where g 0 > g cut2 . When g < g 0 , mode 2 behaves normally. However, when g > g 0 , β 2 < 0 indicates that this is an exponentially increasing mode, which should be removed from the total solutions (equation (14)), resulting in the absence of mode 2. The most special case is when g = g 0 (or g → g − 0 ), which corresponds to a lossless (or low loss) mode. In frequency domain this mode would result in a narrow spectrum. However, the effective free electrons number that satisfy this condition is so large that it is almost impossible for a metallic nanostructure. Therefore, in the rest of this work, we only consider the number at the order of magnitudes of N = 10 6 unless otherwise specified.
In equation (14), A 1 , A 2 and A 3 represent the amplitudes of the three corresponding modes of x(t). When considering the far field, one should use the amplitudes ofẍ(t), i.e., A 1 , A 2 and A 3 . Obviously, the frequency of the excitation light plays a significant role in the amplitudes. Figure 4 shows these amplitudes as a function of ω ex . In the first case ( figure 4(a)), i.e., two same oscillators, the coupled resonant circular frequencies (relative to ω 0 ) are calculated as (ω 1 − ω 0 )/ω 0 = 0.11 and (ω 2 − ω 0 )/ω 0 = −0.13. We find that to obtain the maximum intensities of mode 1 and mode 2 of the emission field, the circular frequency of the excitation light ω ex should be close the corresponding circular resonant frequencies. For mode 3, there are two peaks when varying ω ex , which correspond to around ω 1 and ω 2 , respectively. In the second case ( figure 4(b)), i.e., two different oscillators (ω 1 > ω 2 ), the coupled resonant circular frequencies (relative to ω 0 ) are calculated as 0.22 and −0.24, respectively, which, however, corresponds to a weak coupling due to the frequency splitting is small. This result has been identified in figure 2. Also, the intensities of mode 1 and mode 2 for far field reach their maximums when ω ex is close to the resonant circular frequencies for each of them, and the two corresponding peaks appear for mode 3.
For practical purpose, we consider two metallic nanostructures, e.g., gold nanorods or nanospheres, as the two oscillators, each of which has an individual resonant mode. Figure 5 shows the coupling PL spectra for different coupling strengths at two different excitation wavelengths, calculated from equation (20). With the increase of g, the splitting of the two modes of PL increases, and the total emission intensities decrease. The decrease of the intensities origin from equations (15a) and (15b). Take equation (15a) as an example to explain. The amplitudes depend not only on |ω 1 − ω ex | (this has been discussed in figure 4), but also on |ω 1 − ω 2 |. When g increases, |ω 1 − ω 2 | increases, resulting in the decrease of the amplitude of mode 1. So Figure 5. PL spectra of the two coupled oscillators at g = 0.8 eV (black), 1.0 eV (orange) and 1.2 eV (blue), respectively, calculated from equation (20). The excitation light is at the wavelength of λ ex = 532 nm (a) and λ ex = 633 nm (b). Here, λ c1 = 550 nm and λ c2 = 650 nm represent the resonant wavelengths for each oscillator (before coupling), respectively; β 01 = β 02 = 0.247 eV. Vertical dashed lines stand for the position of 532 nm (green) and 633 nm (red), respectively. does mode 2. Therefore, the PL intensities decrease as g increases. Besides, when excited by 532 nm laser, mode 1 is close to it, resulting in a larger intensity than the one of mode 2. While excited by 633 nm laser, mode 2 is close to it, resulting in a larger intensity than the one of mode 1. This is consistent with the results in figure 4. Here, unit 'eV' and unit 'Hz' for g satisfy the following relationship: where is the reduced Planck constant. So does the damping coefficient β. Figure 6 shows the coupling white light scattering spectra for different coupling strengths in different cases, calculated from equation (21). In figure 6(a), i.e., two oscillators with different resonant wavelengths, with the increase of g, the splitting of the two modes increases, which behaves the same as PL does. However, the scattering intensities stay in the same level which is different from PL spectra. In figure 6(b), i.e., two same oscillators, with the increase of g, mode 2 red-shifts, while mode 1 is hardly to be obtained. Also, the intensities stay in the same level. This behavior agrees well with the experiments [26][27][28]. For example, figure 7 shows the calculated scattering and PL spectra from our model and the experimental ones from Shen et al [28]. In figure 7(a), the scattering spectra of the experiment and the theory agree well with each other. The peak of the coupled one red-shifts from the single one. In figure 7(b), the PL spectra of the experiment and the theory agree a little well. Here, we use Fermi-Dirac distribution which is employed in a previous work to modulate the spectra [18]. The peaks do not agree very well, because more physical processes are not considered here, such as electron-phonon and electron-electron interactions [19]. These processes are useful for single mode of PL spectra. However, the coupled cases are more complex so that more detailed processes should be considered and discussed, which is another topic and is worth studying in the future. Normalized PL spectra. Black circles and red squares stand for the experimental data of single and coupled oscillators, respectively, copied from Shen et al [28]. Black solid curves and red dashed curves stand for the theory data of single and coupled oscillators, respectively, calculated from this model combined with Fermi-Dirac distribution. Here, λ c1 = λ c2 = 610 nm, β 01 = β 02 = 0.187 eV, N = 3.7 × 10 5 , λ ex = 532 nm, g = 0.8 eV. The temperatures in PL spectra are T = 400 K for single case, and T = 1400 K for coupled case.

Conclusions
In summary, we develop a coupling classic harmonic oscillator model to explain the coupling PL spectra as well as the white light scattering spectra from two coupled metallic nanostructures. Each nanostructure is treated as a classic charged oscillator with its own single mode. The coupling coefficients are obtained from the electric interactions between the charges, and are proportional to the velocity and the acceleration of the oscillator, respectively. The behaviors of the two new generated modes due to the coupling are different under different conditions. In general, they split and the splitting gets large as the coupling strength g increases at the beginning. Meanwhile, tuning effective free electron number N, when g gets large enough, there exist cut-off coupling strengths for both modes, and a maximum frequency for one of the modes. Besides, PL spectra and white light scattering spectra are calculated from the model, and their behaviors varying with the coupling strength agree well with the experimental ones of other researchers' work. It is worth noting that this coupling model could be expanded to other wavebands dealing with two coupled single-mode resonators.