Signatures of Plexcitonic States in Molecular Electroluminescence

We develop a quantum master equation (QME) approach to investigate the electroluminesence (EL) of molecules confined between metallic electrodes and coupled to quantum plasmonic modes. Within our general state-based framework, we describe electronic tunneling, vibrational damping, environmental dephasing, and the quantum coherent dynamics of coupled quantum electromagnetic field modes. As an example, we calculate the STM-induced spontaneous emission of a tetraphenylporphyrin (TPP) molecule coupled to a nanocavity plasmon. In the weak molecular exciton-plasmon coupling regime we find excellent agreement with experiments, including above-threshold hot luminescence, an effect not described by previous semiclassical calculations. In the strong coupling regime, we analyze the spectral features indicative of the formation of plexcitonic states.


Introduction
The electroluminesence (EL) of individual quantum emitters coupled to metallic electrodes has been investigated extensively since the first scanning tunneling microscope induced luminescence (STML) experiments were performed 1,2 .Through precise control of an STM probe's position, both the emitter-probe coupling and the resonant frequency of collective motion of the metallic electrons (plasmons) confined in the nanocavity formed between the probe and substrate can be tuned.By adjusting the plasmon frequency, the coupling between particle-hole excitations (excitons) on the quantum emitter and the plasmons can be controlled, leading, for instance, to the observation of plasmon enhanced photon emission 3,4 .This unprecedented control makes STML systems ideal both for exploring fundamental aspects of the nonequilibrium electro-optical response of quantum emitters and as a testbed to develop quantum-enhanced device technology, e.g.those related to biological sensing 5 , photovoltaic energy conversion 6 , or non-classical light generation [7][8][9][10][11] .
We focus on molecular emitters in particular since they can be engineered with the specific emission profiles, dipole moments, wavelengths, and symmetries necessary to harness uniquely quantum resources which may be useful in the development of novel opto-electronic devices [12][13][14] .When a voltage bias causes the source and drain electrodes' chemical potentials to align with unoccupied and occupied molecular states, respectively, a tunnel current and subsequent molecular exciton are produced.If the exciton decays radiatively, the resulting EL encodes the specific electronic and vibrational state of the molecule.Molecular vibrational states have been observed in STML experiments of single porphryin molecules [15][16][17][18][19] and fullerene C 60 and C 70 clusters 20,21 .
Dong et al. observed molecular hot-luminescence (HL) from excited vibrational modes in tetraphenylporphyrin (TPP) molecules weakly coupled to metallic electrodes 18 .Their data are a direct observation of the strong dependence of the EL on the resonant frequency of the localized nanocavity plasmons.In addition, their report of a violation of Kasha's rule, which states that the lowest vibrational transitions should dominate the molecular flourescence, indicates a strong enhancement of the spontaneous emission rate (i.e. a Purcell enhancement 22 ) caused by the formation of the nanocavity 23,24 .Interestingly, above threshold HL (i.e.eV < hω) was also observed in TPP junctions. 18In the weak coupling limit, this effect doesn't appear to be described using a classical plasmonic field, 25 although it may be explained when higher-order electron-plasmon scattering processes are included 26 .
In the study of quantum electrodynamics, Purcell enhancement is a signature of the weak coupling regime between coupled quantum emitters and optical modes.As the coupling strength is increased there is a transition into the strong coupling regime, where energy transfers coherently between the emitter and field modes, giving rise to an observable Rabi splitting between the joint emitter-field states.Systems operating in the strong coupling regime allow for the observation of quantum effects, including single-atom lasing, single photon generation, and all-optical single photon switching 10,11,[27][28][29] .
In this article, we develop a state-based quantum master equation (QME) approach to investigate the EL of molecules in both the weak and strong plexcitonic coupling regimes.We first derive an effective multi-state Jaynes-Cummings model for the molecule and quantum plasmon modes, and use the QME framework to describe finite tunneling currents, radiative and non-radiative exciton decay paths, vibrational damping, and finite plasmon lifetimes.Although similar methods have been used to investigate plasmon-enhanced EL and transport-induced EL in STML systems before 24,25,41,42 , we extend these works to describe the quantum optical regime including a full quantum many-body description of the molecule, plasmon modes, electrodes and their couplings.As a first application, we simulate the EL of a TPP molecule coupled to a single quantum plasmon mode for several voltages, plexcitonic coupling strengths, and detunings.

Theoretical Model
We consider open quantum systems composed of a molecule coupled to electromagnetic field modes, metallic leads (e.g. the substrate and STM probe), and vibrational modes subject to applied voltages and temperature gradients.A schematic of the STM-based experiments we consider is shown in Fig 1 .The Hamiltonian corresponding to this system may be partitioned as where H mol is the molecular Hamiltonian, and the independent lead, vibrational, and electromagnetic baths are described by respectively, where c k annihilates an electron in lead mode k with dispersion ε kσ and spin σ , b l annihilates a vibrational excitation (phonon) in mode l with energy hΩ l , and a j annihilates a photon in mode j with energy hω j .
The tunnel coupling, vibrational mode couplings, and electromagnetic field couplings are given by respectively, where V nk is the tunneling matrix element, d nσ annihilates an electron in molecular state n with spin σ , and W nk is the vibrational coupling between molecular orbital n and lead mode k.The current density J(r) couples to the vector potential A(r), which for quantized electromagnetic modes is given by 43 where A = ∑ j A j , η(r) is the product of the polarization vector and a function describing the spatial profile of the field, V is the effective mode volume, and hω j is the energy of mode j.

Quantum Master Equation
In general, the system described by Eq. 1 cannot be solved exactly.To proceed, we utilize a state-based quantum master equation (QME) approach, where the quantum dynamics of the joint molecule and plasmon system are treated exactly, while the other macroscopic degrees of freedom are traced over using a coarse-graining procedure.Within the QME framework, the Liouville equation for the reduced density matrix of the system is given by where ρ is the density operator, and the Liouvillian superoperators L tun , L damp , and L deph describe the nonhermitian evolution of the system due to quantum tunneling, damping of the populations and coherences of states, and pure dephasing, respectively.
Once the density matrix is determined, expectation values of observables may be calculated using O = Tr {ρ(t)O(t)}.
After coarse graining, our free-system Hamiltonian H 0 is composed of three terms: the molecular Hamiltonian, quantum plasmon modes, and their couplings.The Hamiltonian of each term is given by respectively, where nm is the one-body portion of the molecular Hamiltonian, which may be renormalized by classical electrostatic 44 (e.g.image charge) or vibrational effects induced by the electrodes; U is the Coulomb integral, bl annihilates a (renormalized) phonon in mode l with energy h Ωl ; λ is the electron-phonon coupling; and ã j annihilates a plasmon in mode j with energy h ω j .
For systems we consider, the dipole approximation of the electromagnetic coupling is sufficient.In this approximation the plexcitonic coupling parameter is given by 45 hg j nm = where j is the plasmon mode index, n and m are level indices, V is the mode volume, µ nm = −e n | r| m is the transition dipole matrix element, ω j is the mode's angular frequency, ε 0 is the permittivity of free space, and u j (x 0 ) is the mode function evaluated at the emitter's position x 0 .
Although we have expressed the molecular Hamiltonian in terms of electron and phonon operators, Eqs.8-10 are essentially a multi-state Jaynes-Cummings model where the state energies and matrix elements can be found using a variety of methods (e.g. via exact diagonalization, density functional theory, etc.).

Quantum Transport
Electron transport involves the addition and removal of charges which maintain a degree of phase coherence as they traverse a junction.We consider systems in which the molecules and electrodes are deliberately decoupled (e.g. via the growth of insulating layers on the metallic substrate 15 , or by depositing several molecular monolayers 17 ) such that the individual molecule's emission is not quenched by interactions with the metallic electrodes 46,47 .In this regime, the coherence time of electrons on the molecule are short compared to tunneling time, allowing us to neglect the excitation of coherent superposition states and instead describe the transport as a simple kinetic process. 48,49 Following an expansion of the Liouville equation for the time evolution of the density matrix to second-order in H tun , the master equation for tunneling is given by 24, 41, 50 where R i→ j is the charging rate between the N-particle state i and the N+1-particle state j, R j→i is the discharge rate between states j and i, and σ i j is a matrix in the free system's state space with element (i, j) = 1 and all other elements equal to 0. The electronic tunneling rates are given by where F i j are the Franck-Condon factors (i.e. the overlap between nuclear wave functions), and The effective tunnel coupling between electrode α and the molecule is given by 51 Γα where the bare electron tunneling rate matrix Γ α nm (E) = 2π/h∑ kσ ∈α V nk V * mk δ (E − ε kσ ) is dressed by the many-body renormalization factors 51 [C (i, j)] nσ ,mσ As indicated by Eq. 14, both the relative phase and magnitude of the many-body factors influence the effective tunneling rates and therefore the transport and optical response of systems with multiple states.In addition to the many-body wave function normalization, where the total resonance width of a molecular state is reduced by a factor of 1/N (N is the number of atomic orbitals), strong correlations can also lead to an exponential suppression of these terms.

Damping and Dephasing
When a free system interacts with the environment, an initially excited state can decay via a number of irreversible damping processes.We account for these loss mechanisms with the composite Liouvillian operator L damp = L rad + L cav + L vib , which describes radiative decay processes, the finite lifetime of the nanocavity plasmons, and vibrational relaxation processes, respectively.Assuming Markovian baths, L rad can be expressed as a Lindblad master equation 25, 41, 43, 50 where γ j→i rad is the radiative coupling rate between electronic levels i and j.The finite plasmon lifetime is included via the phenomenological decay rate κ j and master equation where ã j annihilates a plasmon in mode j.The observed optical gap of TPP is 1.89eV, with a vibrational state spacing of 0.16eV. 18Based on comparison with experiment, the source and drain tunnel couplings are set to hΓ S = hΓ D =16.4µeV, while the vibrational damping γ vib , and radiative decay γ r , are consistent with a vibrational lifetime of 2ps and a radiative lifetime of 2ns. 18We set the plasmonic decay to hκ = 10meV and plexcitonic coupling g is taken as an adjustable parameter.
The intraband vibrational damping may be described by 50 where γ vib is the vibrational coupling rate, and P(Ω) = e −hΩ/k B T /Z , with the partition function Z = ∑ k e −hΩ k vib /k B T .The states labeled i and j belong to the same electronic manifold.Although γ vib is typically several orders of magnitude larger than the radiative relaxation rate, γ rad may be enhanced (e.g. by placing an emitter in a cavity) to exceed the vibrational decay rate, resulting, for instance, in the observation of HL from excited vibrational states 18,24,25 .
In addition to relaxation, where excitations are transferred from the system to modes of the environment, pure dephasing is also possible, in which the populations are unaffected but their coherence is reduced.The Liouvillian term describing dephasing is given by 52 where γ deph is the dephasing rate.The relationship between κ, the total state coupling Γ total , and the plexcitonic coupling g distinguishes the strong (g ≫ Γ total , κ) and weak (g ≪ Γ total , κ) coupling regimes.According to Eq. 11 the strong coupling regime may be achieved by increasing the transition energy, increasing the molecular dipole moment, or decreasing the effective mode volume accessible to the plasmons (e.g. by patterning the substrate to exhibit a reduced density of modes).

Results and Discussion
As a first test of our theory, we investigate the STML of a tetraphenylporphyrin (TPP) molecule coupled to a single quantum plasmon mode.In general, the spontaneous emission spectrum may be found from the density matrix though the use of the quantum regression theorem 53,54 .Since our model for the TPP molecule is constructed ad hoc from experimental data, we consider it to be a sum of Lorentzians given by 24,25,41 5/10  hg 00 =10meV, hκ=10meV, ∆=0meV hg 00 =20meV, hκ=10meV, ∆=0meV hg 00 =30meV, hκ=10meV, ∆=0meV The calculated STML of a TPP molecule coupled to a single plasmon mode is shown as a function of photon energy for four plexcitonic coupling strengths hg 00 =0meV, 10meV, 20meV, and 30meV at three bias voltages V b =1.8V, 1.9V, and 2.5V.In the decoupled cavity limit (hg 00 = 0meV), the molecular vibrational spectrum is recovered.As g 11 is increased plexcitonic states form giving rise to the characteristic split peaks in the EL.Although all radiative transitions have non-zero coupling to the plasmonic field, the resonant Q(0, 0) transition dominates and the plexcitonic splitting is ∼ 2hg 00 .The coherent mixing of multiple states gives rise to the asymmetric peak structure and the above threshold hot HL observed when V b =1.8V; even for junctions operating in the weak coupling regime.For this system, the strong coupling boundary is hg 00 >∼15meV.Calculations are for junctions operating at T =80K, to be consistent with measured STML spectra. 18ere γ j→i rad is the radiative decay rate between states j and i, ρ ss is the steady-state solution of Eq. 7, and τ ji is the full width at half maximum (FWHM) of the EL.We assume that the plexcitonic EL linewidths are well approximated by the zero-detuned single atom vacuum Rabi splitting linewidths τ ji = (2τ 0 ji + κ)/2, where τ 0 ji =0.05eV is the FWHM extracted from experiments operating in the weak coupling limit 18,41 .
A schematic of the TPP molecule and the energy level energy diagram for the TPP system are shown in Fig. 2. The Q-band energy gap of TPP (i.e. the Q(0, 0) transition, S 1 (0) → S 0 (0)) is 1.89eV, with a vibrational level spacing of 0.16eV. 18,55 e consider the TPP molecule coupled to a single quantum plasmon mode with energy h ω = 1.89eV -∆, where ∆ is the detuning parameter.The source and drain tunnel coupling rates (Γ S = Γ D = 0.4 × 10 10 s −1 ), vibrational decay rates (γ vib = 0.5 × 10 12 s −1 ), and radiative decay rates (γ j→i rad = F ji /2ns) are established by comparison with experiment. 18,25,41 Te transition rates between vibrational levels and the transition dipole moments are scaled by the appropriate Franck-Condon factors F i j , which are found using the harmonic approximation 56 with a Huang-Rhys parameter S = 0.61 (See Supporting Information). 25e assume that the radiative lifetime of all states are equal (i.e.γ HL rad = γ rad ), and calculations were performed using a modified version of QuTIP. 57or a resonant plasmon, the plasmonic decay rate κ may be expressed in terms of the mode energy ω and quality factor Q as κ = ω/Q.Using reports for other metal-insulator-metal nanostructures, 58 we find that values of Q ≈ 100 are reasonable.Given that a plasmonic mode's lifetime (and the cavity quality factor) can vary significantly for different junction designs, probe positions, and substrate materials, we consider Q=189 in our calculations, such that hκ=10meV.
The calculated STML of a TPP junction for four plexcitonic couplings and three bias voltages are shown as a function of photon energy in Fig. 3.In the decoupled cavity limit (top left panel, hg 00 =0meV), the experimentally observed spontaneous emission peak structure 18 is recovered, where the three peaks corresponding to the Q(0, 0), Q(0, 1), and Q(0, 2) transitions of the TPP molecule.As g 00 is increased, the molecular exciton and cavity plasmon states mix, forming plexcitonic states separated in energy by ∼ 2hg i j .In the weak coupling regime the radiative transition rate is enhanced via the Purcell effect, making HL possible, while in the strong coupling regime energy exchanges coherently between plasmonic and molecular states resulting in a characteristic split peak structure of the EL.The boundary between weak and strong coupling regimes is defined by |g 00 |/|κ t | = 0.25, where for the TPP junction κ t ≈ κ + τ 0 i j = 60meV and therefore hg 00 =∼15meV.The STML of a TPP junction operating in the weak coupling regime is shown in the top right panel of Fig. 3, where a HL peak at 2.05eV is visible for all bias voltages.Although the peak was observed experimentally at 1.8V, where eV b is less than the excitation energy of the molecule, 18 it was not seen at this voltage in previous calculations using classical plasmonic fields. 25This implies that Purcell enhancement alone can't explain the measured HL.
The strong coupling STML of the TPP junction is shown in the lower two panels of Fig. 3, where the Q(0, 0) and Q(0, 1) peaks have split into two peaks.The asymmetry of these plexcitonic peaks stems from the influence of multiple detuned resonances, where the plexcitonic coupling between levels is reduced by the appropriate Franck-Condon factor.The Q(0, 2) and HL peaks are also split but can't be identified with κ = 10meV due the reduced effective couplings.
Our simulations show that HL is suppressed when the off-diagonal coupling terms g i = j are reduced or, as expected, when the vibrational relaxation rate is increased.This suggests that the HL peaks are a consequence of the (weak) coherent dynamics between off-resonant states.Although the above-threshold emission has been explained in terms of higher-order many-body processes 26,59,60 , our calculations support an additional physical explanation in which the tunnel current pumps energy into cavity modes via the nascent plexcitonic states.As shown in the lower panels of Fig. 3, when the plexcitonic coupling is increased the detuning between states is reduced and, consistent with our argument, the above-threshold emission is enhanced.
Finally, we consider the influence of the detuning between molecular transitions and cavity plasmon resonance energies on the EL.Physically, detuning can be controlled by adjusting the STM probe's height above the substrate.As evidenced by the peak at 2.05eV shown in the left panel of Fig. 4, Purcell enhancement in the weak coupling regime results in HL for all detunings.When the cavity is blue or red detuned relative to the Q(0, 0) transition of TPP, the STML spectral weight is shifted towards higher or lower energy peaks, respectively.Since the TPP junction supports a finite current, energy is constantly (albeit weakly) pumped into off resonant cavity plasmon modes, giving the observed shift to the spectrum.The EL is increased by detuning since detuning reduces the molecule's effective coupling to non-radiative plasmonic decay paths.
In the strong coupling regime, shown in the right panel of the same figure, the blue and red detuned plasmon modes again shift the STML spectral weight up or down in energy, respectively.However, in this regime the strongest peaks are split into distinguishable plexcitonic resonances.In addition to the characteristic split peak EL, plexcitonic states and the onset of the strong coupling regime can also be identified via this distinct spectral weight shift with detuning.

Conclusions
We develop a QME approach to investigate the STML of molecules coupled to quantized electromagnetic modes.Within our method we include the effects of electronic tunneling, vibrational damping, and environmental dephasing, and can describe both weak and strong plexcitonic coupling regimes.Our approach extends existing methods and includes a full quantum description of the coherent state dynamics.Our method is valid in both single-particle and many-body representations, allowing future studies to balance computational effort with chemical accuracy.
Motivated by the observation of HL in the STML of TPP, 18 and the argument that it was a consequence of STM-induced Purcell enhancement [23][24][25] , we calculated the EL of a TPP molecule coupled to a single quantum plasmon mode.In the weak coupling regime, we recover the experimentally observed spectra, including the above-threshold HL.Using a fully quantum plasmon theory, we conclude that the low-bias HL peak may be a consequence of the weakly coherent energy exchange dynamics.Finally, we identify several signatures of the formation of plexcitonic states: a split peak structure of the EL, and the shifted spectral weight as the plasmon resonance is tuned.
Although the strong coupling regime has not yet been observed in STML systems, molecular systems with coupling strengths of hundreds of meV have been fabricated [32][33][34][35][36][37][38][39][40] .For the TPP system investigated here, we find that it is physically plausible to achieve the strong coupling regime if the nanocavity losses are reduced slightly, e.g.via careful material selection or patterning of the substrate to reduce the effective plasmon mode volume (See Supporting Information).If realized, the ability to measure the spatial distribution of the electro-optical response of molecules operating in the strong coupling regime would be invaluable in the development of myriad quantum information applications, and would herald a new phase in the study of QED and molecular dynamics.

Figure 1 .
Figure 1.A schematic diagram of an STML experiment including the STM probe, molecule, and substrate.Molecular excitons generated by a finite tunnel current couple to nanocavity plasmons and may decay radiatively as EL.
the Fermi-Dirac distribution for lead α with temperature T α and chemical potential µ α .The chemical potential of the source and drain leads are given by µ S = E f − eαV b and µ D = E f − e(1 − α)V b , respectively, with Fermi energy E f , electron charge magnitude e, voltage symmetry α, and bias voltage V b .We assume a symmetric potential drop, where α = 0.5.

Figure 2 .
Figure 2. The energy level diagram for the charge neutral manifold of tetraphenylporphyrin (TPP) molecule and coupled quantum plasmon modes.The relevant (Q-band) electronic states S 0 and S 1 are shown with their associated vibrational levels.The observed optical gap of TPP is 1.89eV, with a vibrational state spacing of 0.16eV.18Based on comparison with experiment, the source and drain tunnel couplings are set to hΓ S = hΓ D =16.4µeV, while the vibrational damping γ vib , and radiative decay γ r , are consistent with a vibrational lifetime of 2ps and a radiative lifetime of 2ns.18We set the plasmonic decay to hκ = 10meV and plexcitonic coupling g is taken as an adjustable parameter.