Unified Treatment for Scattering and Photoluminescence Properties of Strongly Coupled Metallic Nanoparticle Chains based on a Coupling Classic Harmonic Oscillator Model

Growing interest in coupled systems, especially strongly coupled metallic nanostructures, has been realized by numerous investigators, because the ability and potential of them in optical devices makes it possible to tune the electric field at subwavelength scale. Particularly, localized surface plasmon resonance (LSPR) of metallic nanostructures, e.g., gold nanorods and nanospheres etc., are widely investigated. Due to the amazing optical properties of LSPR, e.g., strongly enhancement and highly confinement of electric field, it has been applied in numerous technological applications such as biosensing Lu et al. (2012); Zhang et al. (2018); Wu et al. (2018); Qiu, Ng, and Wu (2018), optical recording Zijlstra, Chon, and Gu (2009); Taylor, Kim, and Chon (2012), optical waveguiding Catalano et al. (2018); Diez et al. (2018); Ebrahimpouri et al. (2018), nonlinear optics Mi et al. (2019); Wurtz et al. (2011); Kauranen and Zayats (2012), and nano-optical devices Bhuyan et al. (2018); Atabaki et al. (2018); Cheben et al. (2018).

Multiple nanoparticles including dimers result in new LSPR properties when they are strongly coupled with each other, compared with individual nanoparticle. This phenomenon has been widely investigated both experimentally and theoretically by researchers. S. Biswas $et\ al$. employ self-assembled gold nanorod heterodimers to achieve plasmon-induced optical transparency (PIT) uniquely in the visible wavelengths with a slowdown factor of 10 and extreme dispersion Biswas et al. (2013). L.-J. Black $et\ al$. tune the polarization conversion of coupled gold nanorod dimer by controlling the gap width. They conclude that much higher conversion efficiencies may be obtained in “kissing” antennas Black et al. (2014). C.-Y. Tsai $et\ al$. investigate the optical properties of gold nanoring dimers and find that LSPR peak of them strongly depends on the polarization direction and gap distance, i.e., longitudinal and transverse polarizations result in the red-shift and blue-shift, respectively, when the nanorings approach each other Tsai et al. (2012). P. Mulvaney $et\ al$. investigate the scattering properties of linear chains of strongly coupled gold nanoparticles with particle number from 1 to 6 by employing self-assembly techniques. They obtain the red-shift of LSPR peak when particle number increases and conclude that a maximum resonance wavelength should be at a chain length of 10–12 particles Barrow et al. (2011). A. M. Soehartono $et\ al$. employ tall nanosquare dimer for biosensing applications, where the resonance mode and hybrid modes are investigated and the longitudinal mode is found to be red-shifts when increasing the height Soehartono et al. (2020). L. S. Slaughter $et\ al$. investigate LSPR properties of symmetry breaking gold nanorod dimers with different configurations. They find that symmetry breaking has the strongest effect on the collective plasmon modes Slaughter et al. (2010). S. Yoo $et\ al$. report the novel strategy for the synthesis of complex 3-dimensional (3D) nanostructures and apply them to surface enhanced Raman scattering (SERS). The coupling phenomenon is illustrated to be directly related to the intrananogap and interior volume size Yoo et al. (2020). All these obtained coupling phenomena have pointed to the fact that strong coupling indicates large splitting of LSPR peak. Explanations based on simulations and models have been done by these researchers. Particularly, the the models usually provides coupled LSPR peak positions. However, the spectra shapes of these coupled systems are seldom estimated in theory. In a previous work, we investigate the spectra of metallic nanoparticles with individual mode without any coupling Cheng et al. (2018). In another previous work, we only consider the dimer case for the coupled system. Multi-particle coupled system was not realized due to its complexity Cheng and Sun .

In this study, we investigate a metallic multi-particle system which is arranged as a linear chain. We develop a coupling classic harmonic oscillator model for multi-particle to study the optical properties. Taking 1 to 6 particles as examples, we carefully analyze the scattering and PL properties varying with particle number $n$, coupling strength $g$ and effective free electron number $N$. The theoretical results illustrate that, for scattering spectra, LSPR peak red-shifts and full width at half maximum (FWHM) decreases as $n$ or/and $g$ increases; for PL spectra, the splitting of the modes gets larger with larger $g$, and the amplitudes are dependent on the excitation wavelength. Comparisons between our model and experimental data from Mulvaney’s work Barrow et al. (2011) demonstrate that our model is practical and accuracy. This work would give rise to understanding complex coupling phenomena more deeply.

Now we consider a metallic nanostructure chain consisted by $n$ nanoparticles arranged in a line. The schematic is shown in Fig. 1. The $x$-polarized incident light excites the system with frequency $\omega_{ex}$ and amplitude $E_{0}$. Due to the fact that the coupling strength $g$ decreases rapidly with the increase of $r_{0}$, we only consider the interaction between the neighboring particles. For the simplest case, we assume that all the particles are the same. Therefore, the dynamical equations of these oscillators should be in this form:

		$\displaystyle\ddot{x}_{1}+2\beta_{0}\dot{x}_{1}+\omega_{0}^{2}x_{1}-\gamma\dot{x}_{2}-g^{2}x_{2}=K_{0}\mathrm{exp}(-\mathrm{i}\omega_{ex}t),$		(2)
		$\displaystyle......$
		$\displaystyle\ddot{x}_{j}+2\beta_{0}\dot{x}_{j}+\omega_{0}^{2}x_{j}$
		$\displaystyle-\gamma\dot{x}_{j-1}-g^{2}x_{j-1}-\gamma\dot{x}_{j+1}-g^{2}x_{j+1}=K_{0}\mathrm{exp}(-\mathrm{i}\omega_{ex}t),$
		$\displaystyle(\mathrm{for}\ 2\leq j\leq n-1)$
		$\displaystyle......$
		$\displaystyle\ddot{x}_{n}+2\beta_{0}\dot{x}_{n}+\omega_{0}^{2}x_{n-1}-\gamma\dot{x}_{n-1}-g^{2}x_{n-1}=K_{0}\mathrm{exp}(-\mathrm{i}\omega_{ex}t).$

Here, $x_{j}(t)$, $\dot{x}_{j}(t)$ and $\ddot{x}_{j}(t)$ are the displacement relative to equilibrium position, velocity and acceleration of $j$-th oscillator, respectively, $K_{0}=-eE_{0}/m_{e}$, and $e$ is the charge of electron.

Firstly, we deal with the scattering properties. Assume that $x_{j}(t)=A_{j}\mathrm{exp}(\alpha t)$ ($j=1,\ 2,\ ...,\ n$) are the solutions of Eq. (2). Define $B=-\gamma\alpha-g^{2}$, $C=\alpha^{2}+2\beta_{0}\alpha+\omega_{0}^{2}$, and

$$D=\begin{pmatrix}C&B&&&&&\\ B&C&B&&&&\\ &B&C&&&&\\ &&&\ddots&&&\\ &&&&C&B\\ &&&&B&C&B\\ &&&&&B&C\end{pmatrix}_{n\times n},$$

(3)

where $D_{j_{1}j_{2}}=C$ for $j_{1}=j_{2}$, $D_{j_{1}j_{2}}=B$ for $|j_{1}-j_{2}|=1$, and $D_{j_{1}j_{2}}=0$ for other cases. Also define:

	$\displaystyle X=(A_{1},\ ...,\ A_{j},\ ...,\ A_{n})^{T},$		(4)
	$\displaystyle K=K_{0}(1,\ ...,\ 1,\ ...,\ 1)^{T}.$

After substituting $x_{j}(t)$ into Eq. (2), the amplitudes should satisfy the following equation:

$$DX=K,$$

(5)

where the solution of $\alpha$ is $\alpha=-\mathrm{i}\omega_{ex}$. The solutions of the amplitudes are:

$$X=D^{-1}K,$$

(6)

where $D^{-1}$ is the inverse of matrix $D$. The total far field emission of the electric field is proportional to the accelerations of the oscillators Griffiths (2013); Cheng et al. (2018); Cheng and Sun :

	$\displaystyle\ddot{x}(t)=\sum_{j=1}^{n}\ddot{x}_{j}(t)$	$\displaystyle=(-\mathrm{i}\omega_{ex})^{2}\left(\sum_{j=1}^{n}A_{j}\right)\mathrm{exp}({-\mathrm{i}\omega_{ex}t})$		(7)
		$\displaystyle=-\omega_{ex}^{2}A^{(n)}\ \mathrm{exp}({-\mathrm{i}\omega_{ex}t}).$

Actually, $A^{(n)}$ satisfy:

$$A^{(n)}=\sum_{j=1}^{n}A_{j}=\sum_{j=1}^{n}X_{j}=K_{0}\sum_{j_{1}=1}^{n}\sum_{j_{2}=1}^{n}D_{j_{1}j_{2}}^{-1}\ .$$

(8)

We give out the solutions of the amplitudes $A^{(n)}$ for $n=1-\-6$ as examples:

		$\displaystyle A^{(1)}=\frac{1}{C}K_{0},\ \ \ \ \ \ \ \ \ \ \ \ A^{(2)}=\frac{2}{B+C}K_{0},$		(9)
		$\displaystyle A^{(3)}=\frac{4B-3C}{2B^{2}-C^{2}}K_{0},\ \ A^{(4)}=\frac{2(B-2C)}{B^{2}-BC-C^{2}}K_{0},$
		$\displaystyle A^{(5)}=\frac{B^{2}+8BC-5C^{2}}{3B^{2}C-C^{3}}K_{0},$
		$\displaystyle A^{(6)}=\frac{2(2B^{2}+2BC-3C^{2})}{B^{3}+2B^{2}C-BC^{2}-C^{3}}K_{0}\ .$

Define $A^{{}^{\prime}(n)}=-\omega_{ex}^{2}A^{(n)}$, and employ the Fourier transform in Ref. [ Cheng et al. (2018); Cheng and Sun ], thus obtaining the elastic emission spectrum:

$$I_{ela}(\omega)=\left|A^{{}^{\prime}(n)}\right|^{2}\sqrt{2\pi}\delta(\omega-\omega_{ex}).$$

(10)

Therefore, the white light scattering spectra should be given as:

$$I_{sca}(\omega)=I_{ela}(\omega_{ex}\to\omega)=\sqrt{2\pi}\left|A^{{}^{\prime}(n)}(\omega_{ex}\to\omega)\right|^{2}.$$

(11)

Secondly, we deal with the PL properties. As our previous work demonstrates Cheng et al. (2018), PL term origins from the general solutions of the homogeneous linear equations:

$$DX=0.$$

(12)

The necessary and sufficient conditions for the existence of nontrivial solutions of Eq. (12) is that the determinant of $D$ is zero:

$$\mathrm{det}(D)=\frac{z_{+}^{n+1}-z_{-}^{n+1}}{z_{+}-z_{-}}=0,$$

(13)

where $z_{\pm}=\frac{1}{2}(C\pm\sqrt{C^{2}-4B^{2}})$. Eq. (13) determines the solutions of $\alpha$. Obviously, there are $2n$ solutions for Eq. (13). We can rewrite the solutions as $\alpha_{k\pm}=-\beta_{k}\pm\mathrm{i}\omega_{k}$, $k=1,\ 2,\ ...,\ n$. Due to the large difference from the excitation frequency, we omit the solutions of “$\alpha_{k+}$” as our previous work does Cheng et al. (2018); Cheng and Sun . Therefore, the total solutions for Eq. (2) can be assumed and solved by:

		$\displaystyle x_{j}(t)=S_{j}\mathrm{exp}(-\mathrm{i}\omega_{ex}t)+\sum_{k=1}^{n}P_{jk}\mathrm{exp}(\alpha_{k-}t),$		(14)
		$\displaystyle\mathrm{initial\ conditions:}$
		$\displaystyle x_{j}(0)=0,\ \dot{x}_{j}(0)=0,\ \ddot{x}_{j}(0)=K_{0},\ ...,\ \frac{\mathrm{d}^{n}x_{j}(0)}{\mathrm{d}t^{n}},$
		$\displaystyle\mathrm{for}\ j=1,\ 2,\ ...,\ n.$

Here, $S_{j}$ is the amplitude of the elastic term (scattering) of the $j$-th oscillator; $P_{jk}$ is the amplitude of the $k$-th inelastic term (PL) of the $j$-th oscillator. After considering the PL term for the solutions, the total far field emission of the electric field can be written as:

	$\displaystyle\ddot{x}(t)=\sum_{j=1}^{n}\ddot{x}_{j}(t)=$	$\displaystyle(-\mathrm{i}\omega_{ex})^{2}\sum_{j=1}^{n}S_{j}\mathrm{exp}(-\mathrm{i}\omega_{ex}t)$		(15)
	$\displaystyle+$	$\displaystyle\sum_{j=1}^{n}\sum_{k=1}^{n}\alpha_{k-}^{2}P_{jk}\mathrm{exp}(\alpha_{k-}t).$

Define $P_{k}^{\prime}=\sum_{j=1}^{n}\alpha_{k-}^{2}P_{jk}$, employing Fourier transform and Fermi-Dirac distribution Cheng et al. (2018); Cheng and Sun , the total PL spectrum of the system can be written as:

	$\displaystyle I_{PL}(\omega)=$	$\displaystyle\sum_{k=1}^{n}\left|P_{k}^{\prime}\right|^{2}\frac{1-\mathrm{exp}(-2\beta_{k}t_{0})}{2\beta_{k}t_{0}}\frac{\beta_{k}}{(\omega-\omega_{k})^{2}+\beta_{k}^{2}}$		(16)
		$\displaystyle\times\frac{1}{1+\mathrm{exp}[(\hbar\omega-\hbar\omega_{f})/k_{B}T]}.$

Here, $t_{0}$ is the effective interaction time between the excitation light and the oscillators, $\hbar$ is the reduced Planck constant, $\omega_{f}$ is the so-called chemical potential and $k_{B}$ is Boltzmann’s constant, and $T$ is the temperature.

For $n=1$, the solution degenerates into the one of an individual particle, which has been discussed in Ref. [ Cheng et al. (2018) ]. For $n=2$, the solutions degenerate into the dimer case, which has been discussed in Ref. [ Cheng and Sun ].

First of all, scattering properties are analyzed.

Fig. 2a-c show the normalized scattering spectra of the chains with different particle number $n$ varying with effective free electrons number $N$. The coupling strength is $g=1.3$ eV. The primary LSPR peak red-shifts as $n$ increases. When $n\geq 3$, there exist other peaks at blue side of the primary peaks with small amplitudes for each case. Fig. 2d shows the LSPR peak positions as a function of $N$ for different cases. It indicates that for small $n$, the peaks almost stay unchanged as the increase of $N$. However, for large $n$, the peak decrease at about $N=10^{6}$. Notice that in Fig. 2, we keep $g$ unchanged, and Eq. (1) can be rewritten as $\gamma^{3}=g^{4}/\kappa$. When $N$ is small, $\gamma$ is much less then $g$, thus influencing the interaction parts little. However, when $N$ increases, $\kappa$ decreases, resulting in the increase of $\gamma$, which influences the interaction parts greatly if $\gamma$ is comparable with or even larger than $g$. Hence, the peak positions are influenced by $N$. Fig. 2e shows FWHM of different cases as a function of $N$. For $n\geq 2$, FWHM decreases first and then increases as $N$ increases. The minimums of FWHM occur at about $N=10^{5}$, resulting in narrow shapes of the spectra as shown in Fig. 2c.

Fig. 3a-b show the normalized scattering spectra of the chains with different $n$ varying with coupling strength $g$. Here, $N=1.5\times 10^{4}$. The same as Fig. 2a-c, the primary LSPR peak red-shifts as $n$ increases. Obviously, the red-shift of $g=1.3$ eV is larger than the one of $g=1.0$ eV. In Fig. 3c, LSPR wavelength increases as $n$ increases (red-shift), and it also increases as $g$ increases (red-shift). That is, the larger $g$ or/and $n$ is, the larger red-shift is. In Fig. 3d, FWHM decreases as $n$ increases, and it also decreases as $g$ increases, resulting in narrower shapes of the spectra with larger $g$ or/and $n$ in unit of “eV”.

In Fig. 4, we show the comparison between the scattering spectra of these chains calculated from this model and the experimental ones from P. Mulvaney $et\ al$.Barrow et al. (2011). Here, the intrinsic frequency $\omega_{0}$ and damping coefficient $\beta_{0}$ of each identical nanoparticle are $\omega_{0}=2.204$ eV and $\beta_{0}=0.156$ eV, resulting in the resonant wavelength at 560 nm of an individual nanoparticle. The coupling strength and effective free electrons number are $g=1.317$ eV and $N=1.5\times 10^{4}$. Fig. 4a shows the scattering spectra of these chains calculated from this model. As $n$ increases, the scattering intensity increases and meanwhile the peak red-shifts. The inset of Fig. 4a illustrates the peaks with small values (for $n\geq 3$). The intensities of these small peaks increase and they also red-shift as $n$ increases. In Fig. 4b, we compare the normalized scattering spectra of our model with the experiments of P. Mulvaney Barrow et al. (2011). In general, they agree well with each other. Notice that the FWHMs of $n=2$ and $n=6$ agree not very well. For $n=2$, the dimer case, the FWHM of our model is smaller than the one of the experiment. For $n=6$ the result reverses, and the peak of the model is a bit blue-shifted compared to the experiment. The main reason may be the nonuniformity of the samples. Of course, other reasons should be considered. For example, in our model, we omit the interaction parts of the non-neighboring particles, which indeed influences the spectra.

The next, PL properties are analyzed.

The new generated PL modes are shown in Fig. 5, taking $n=2-\-5$ as examples. For a given particle number $n$, there are $n$ modes of PL. Particularly, when $n$ is even, the modes spit, half ($n/2$) of which blue shift, called the blue branches, the other half ($n/2$) of which red shift, called the red branches. When $n$ is odd, however, there exists one mode with frequency and damping unchanged and equaling to $\omega_{0}$ and $\beta_{0}$, respectively, and the rest modes split the same as even case. More special, Mode 1 and Mode 2 of $n=2$ behave the same as Mode 2 and Mode 4 of $n=5$. There exist cut-off coupling strengthes $g_{cut}$ for the red branches, at which the frequency decreases to zero. For Mode 1 of all the cases, $g_{cut}$ decreases as $n$ increases. Notice that in Eq. (16), the shape of each mode is a Lorentzian curve with FWHM $2\beta_{j}$. Hence, the increase (decrease) of $\beta_{j}$ in Fig. 5 indicates the increase (decrease) of FWHM of Mode $j$.

In Fig. 6, PL spectra calculated from Eq. (16) are shown for two different excitation wavelengths, i.e., 532 nm and 633 nm. In Fig. 6a, it reveals the number of modes as 1, 2, 3, 4 and 3 for $n=$1, 2, 3, 4 and 5, respectively. There is only 3 modes in the spectrum of $n=5$, because one mode’s wavelength is larger than 1000 nm so that it is not shown, the other mode’s wavelength is smaller than 532 nm so that Fermi-Dirac distribution makes it disappear. In Fig. 6b, due to the fact that some modes’ wavelengths are smaller than 633 nm, they disappear. The numbers of the modes illustrated here are 1, 2, 2, 3 and 3 for $n=$ 1, 2, 3, 4 and 5, respectively. Notice that, for $n=2$ and $n=5$, there are two modes behaving similarly at around 624 nm and 811 nm, respectively. This phenomenon has been revealed in Fig. 5.

Furthermore, we should emphasize that, in this work, the excitation electric field felt by these particles is treated as the same, i.e., the amplitudes and phases are the same. This approximation is reasonable when dealing with small $n$, e.g., $n\leq 6$. In Mulvaney’s experiments Barrow et al. (2011), the diameter of one particle is 64 nm with gap of 1 nm, resulting in the length of the chain for $n=6$ of about 390 nm (omitting the bend of the chain). For larger $n$, the length of the chain would be so large that it is larger than the excitation wavelength. In experiments, usually, the light is focused through the lens and the size of the facula is near the excitation wavelength due to the diffraction limit. In other words, a large $n$ breaks the identity of the right side in Eq. (2), i.e., the electric field felt by the particles is not the same and it depends on the position of the facula. Therefore, if employing our model to calculate the scattering and PL spectra of a long chain, one just needs to rewrite Eq. (2) according to the situation, and then derives the rest equations as this work does.

The authors declare no conflicts of interest.

The data that support the findings of this study are available from the corresponding author upon reasonable request.