# Theory for the response of hybrid light-matter modes in the presence of vibrations

###### Abstract

We construct a model describing the response of a hybrid system where the electromagnetic field — in particular, surface plasmon polaritons — couples strongly with electronic excitations of atoms or molecules. Our approach is based on the input-output theory of quantum optics, and in particular it takes into account the thermal and quantum vibrations of the molecules. The latter is taken into account by the theory analogous to that used for example in the theory of dynamical Coulomb blockade. As a result, we are able to include the effect of the molecular Stokes shift on the strongly coupled response of the system. Our model then accounts for the asymmetric emission from upper and lower polariton modes. It also allows for an accurate description of the partial decoherence of the light emission from the strongly coupled system. Our results can be readily used to connect the response of the hybrid modes to the emission and fluorescence properties of the individual molecules, and thus are relevant in understanding any utilization of such systems, like coherent light harvesting.

Recent experiments have shown the possibility of coupling the electromagnetic field and electronic excitations so strongly that the coupling energy shows up in absorption and emission spectra of such systems, suggesting formation of hybrid light-matter states, called polaritons Törmä and Barnes (2014); Yu et al. (2019). Because the coupling can change the potential energy surface of these hybrid states with respect to the bare molecules, polaritons are considered a promising paradigm for controlling photochemical reactions Hertzog et al. (2019); Ribeiro et al. (2018); Feist et al. (2017). While the avoided crossing between the excitations of the molecules and confined electromagnetic light mode as such is not a quantum effect, the detailed response of such systems can exhibit quantum interference effects. In particular, some of the present authors showed recently Baieva et al. (2017) that a quantum coherent coupled system of plasmons and a large number of molecules with random polarizations can emit only light that is polarized along the plasmon propagation direction, whereas loss of coherence leads to the presence of also the orthogonal polarization in the emitted light.

Here we construct a detailed description of the main decoherence mechanisms due to the quantum and thermal vibrations of the molecules. This hence allows for describing the effect of inhomogeneous broadening of the molecular response due to such vibrations on the strongly coupled response. We take the vibrations into account via the theory analogous to that used in Coulomb blockade Ingold and Nazarov (1992); Heikkilä (2013). This theory describes the probability of absorbing (for ) or emitting () the energy to/from the vibrations. For the specific models of harmonic vibrations, such a function can be calculated exactly, although in many cases as a result of a somewhat tedious integral. In particular, we here find how this shows up in the absorption and emission spectrum of individual molecules. The resulting fluorescence spectrum is similar to that found via other approaches Mahan (2013). Therefore, an alternative approach is to deduce an effective for the measured spectra of individual molecules. What is more, we then go on connecting this fluorescence spectrum directly to the absorption/emission spectrum of the strongly coupled system. In a certain limit of parameters, the resulting inhomogeneous broadening of the molecular absorption/emission then determines the linewidth of the strongly coupled spectrum exhibiting the avoided crossing between the hybridized plasmon–molecular excitations known as polaritons. In particular, our model explains the asymmetric emission spectra of upper and lower polaritons seen in many experiments Baieva et al. (2017); Bellessa et al. (2004); Hakala et al. (2009); Baieva et al. (2012); Koponen et al. (2013); Guebrou et al. (2012).

To be specific, we consider the plasmon–molecule system in the strong-coupling regime and describe the plasmon by a single bosonic mode of frequency and a given polarization with respect to its wavevector . A concrete example of such a plasmon is the surface plasmon polariton traveling along an interface in the -plane in the -direction with as in Fig. 1a. The plasmon interacts with identical molecules which we consider to be two-level systems with transition frequency . We denote the rising (lowering) operator of the molecules with (). As in typical experiments, we assume that the molecules are distributed to a region much larger than the wavelength of the plasmon and that their electric dipole moments point in uniformly random directions . Following the standard approach of quantum optics Walls and Milburn (2007); González-Tudela et al. (2013), the Hamiltonian of the strong-coupled system is in the rotating wave approximation ()

(1) |

The position of a molecule affects the coupling in two ways: it contains a complex phase factor due to the phase of the plasmon, and the coupling strength depends on the distance to the interface. We disregard the latter effect here for brevity. Also, the coupling strength depends on the angle between the plasmon polarization and dipole moment of a molecule. Thus, we write .

In addition to the strong-coupled system we include the vibrational modes of the molecules. For brevity, we assume that there is a single vibration mode per molecule with vibrational frequency but the generalization to multiple modes is straightforward sup (). These vibrations and their interactions are described by

(2) |

The coupling between electronic and vibrational modes is quantified with a dimensionless parameter , the Huang–Rhys factor Huang and Rhys (1950), which can be related to the Stokes shift measured in fluorescent emission.

We seek an approach to find the response of the strongly-coupled plasmon–molecule system in the presence of vibrations. Instead of transforming the Hamiltonian to the polariton basis, we use the polaron transformation to work in the eigenbasis of the molecule–vibration system. That is, we introduce new operators and . The operator is the generator of vibrations and formally identical to the displacement operator that generates coherent states Glauber (1963). These new operators reduce the molecular Hamiltonian to a diagonal form

(3) |

where . Note that the commutation relations of and are identical to those of and .

Next, we employ the input–output formalism of quantum optics Gardiner and Collett (1985). We assume that there are separate bosonic baths for each molecule, vibration and the plasmon to which the coupling is linear in and , respectively. In the Markov approximation these couplings are described by the dissipation rates , and of the molecules, vibrations, and plasmon. The dynamical equations are in this case

(4a) | ||||

(4b) | ||||

(4c) |

where is the coupling to the external optical field. These nonlinear equations can be simplified with a sequence of approximations. First, we assume the single-excitation limit which corresponds to low driving power where and . This uncouples the dynamics of the vibrational mode in Eq. (4c) from the plasmon–molecule system which allows us to use the Caldeira–Leggett model Caldeira and Leggett (1981) to obtain the vibrational dynamics. Furthermore, we denote as the total effective damping rate of the individual molecules. In addition, we suppose identical molecules and vibrations so that and . Lastly, we neglect the thermal fluctuations of the molecules and set since at room temperature and to represent a coherent driving of the plasmons. When also we may neglect the second row of Eq. (4b) sup (). With these assumptions the molecular equation (4b) simplifies to

(5) |

while the plasmon equation (4a) is unaffected.

We model a measurement on the spectral properties of the system so that the incoming light produces a reflected and transmitted field . These fields contain both the plasmon and the molecular emission but not the emission of phonons from the vibrations. We also separately include coupling to s- and p-polarized light represented by and (see Fig. 1a). Since the propagating plasmon cannot emit s-polarized light to the direction perpendicular to the interface but the molecules have no directional preference, we consider s- and p-polarized output fields separately. The output fields can be described by a general formula

(6) |

In this equation and which describes the fact that only the reflected field interferes with the input field. The is defined similarly to describe the plasmon which couples only to p-polarized modes. The constants describe the coupling of the molecule electronic states to the environmental s- and p-polarized free space modes, and thus . These fields and couplings to the system are represented schematically in Fig. 1b. The output spectral density is obtained from

(7) |

where the first argument is the frequency of the output field and is the driving frequency.

The presence of the vibrations makes the input-output equations (4) non-linear as they contain products of different dynamical fields. We solve this problem using the theory similar to the one in dynamical Coulomb blockade Ingold and Nazarov (1992) as follows. Let us first define and its Fourier transform

(8) |

We assume that does not depend on the molecule index ; this assumption can be relaxed if needed. The function normalizes to unity and is real for stationary vibrations if for any time . We can thus interpret the function as a probability distribution of transforming energy to the vibrations () or vice versa (). It is possible to find an exact expression of when the fluctuations are Gaussian sup (). The function is characterized by four parameters: vibration eigenfrequency , their linewidth , Huang–Rhys factor , and temperature of their bath. In this Letter we use an analytical limit in which vanishes which results in a discrete , depicted in Fig. 2a. In this regime, is related to the absorption function defined by Huang and Rhys Huang and Rhys (1950). However, our analytic results for the response apply also in the case of general .

Stokes shift. Before solving the full plasmon–molecule problem we illustrate how the theory can be used to describe a measurement of the Stokes shift in a molecule–vibration system described by the Hamiltonian (3). This is achieved by removing the plasmon terms from Eqs. (4) and driving the molecules incoherently where represents a random phase. The driving is scaled so that the total input power spectral density is given by . Making the same assumptions as above the equation for simplifies to

(9) |

This equation can be solved by the use of Fourier transform and convolution theorem. The spectra are found from Eq. (7) when the output fields are changed to . Here, ’reflected’ field should not be understood literally but rather as the field that contains the driving field. ’Transmitted’ field is fully from the molecular fluorescence. Since and the solution of Eq. (9) depends on we encounter a four-point correlator in the calculation of . Here, refers to the Fourier transform of . Assuming that the vibration modes are independent and identical, the correlator can be factorized into two-point correlators when . These resulting two- and four-point correlators are related to the function by

(10a) | |||

(10b) |

where is the Fourier transform of

(11) |

The function is related to and its reciprocal function as shown in sup (). The Fourier transform can be calculated formally by noting that inverting the sign of in results in .

After averaging over the random phases the resulting spectra are

(12) | |||

where is the detuning between the driving field and the (renormalized) molecular frequency, and

(13) |

with as the bare molecular response function. The term describes inelastic scattering (output field frequency different from driving frequency ) and can be written as

(14) |

The information about the spectral properties of the molecules including their vibrations is fully contained in the functions and .

The emission at is given directly by which in turn relates to . Since we set , is also delta-peaked at frequencies with an integer . The absorption spectrum is obtained from power conservation evaluated at the driving frequency and it is mostly determined by . In Fig. 2b we have plotted the emission spectrum along with the absorption spectrum using analytical expressions of and sup (). The emission and absorption maxima differ by which is exactly the Stokes shift. The central frequency between the peaks is the renormalized frequency meaning that the absorption peak is at the bare molecular frequency . This corresponds to the results describing the transient response obtained with Green functions Mahan (2013). However, in our stationary model the absorption is not a mirror image of the emission because it is possible for the emission to happen also from the excited vibrational states.

Plasmon–molecule system. To solve Eqs. (4a)–(4b) we first integrate Eq. (5) from an initial time to and neglect the initial condition which has no role in a stationary situation. We then substitute this into Eq. (4a) which leads to

(15) | ||||

At this point we use a mean-field approximation and average the equation over the fluctuating vibrations. This leads to the function since . Consequently, the elastic response of the plasmon field is given by where the response is given by

(16) |

Vibrations provide a channel of relaxation broadening the response which is associated with the real part of . The imaginary part contains information about the frequencies of the polariton modes. When the vibrations are absent, and , the usual strong-coupling response is obtained as with Rabi splitting proportional to at . When the vibrations are present, especially the upper polariton branch is perturbed as Fig. 3 shows.

From Eq. (5) can now be solved in terms of by Fourier transformation using the convolution theorem and .

Incoherent polaritonic response. Finally, let us consider a large number of identical molecules with random dipole moment directions. In this case we can replace the sums over the molecule index with an integral over a surface of a sphere . Then, the square of the Rabi splitting in Eq. (16) is . We also assume that the positions of the molecules are random so that we may replace for an ensemble average. Using these assumptions the polarization dependence shows up in the spectra as the coefficients

(17) |

Above, only the terms where contribute in the sum which results in four-point correlators as in Eq. (10b). If the assumption on the random positions of the molecules is relaxed and the molecules are driven by the plasmon coherently, we obtain also the two-point correlators of Eq. (10a) and different coefficients sup ().

The s- and p-polarized emission spectra are

(18a) | ||||

(18b) |

The main difference between the s- and p-polarized spectra is due to the fact that the plasmon can only emit p-polarized light. The molecular fluorescence is also slightly enhanced in this polarization for .

Both the s- and p-polarized emission spectra are now represented with and found in molecular fluorescence Eq. (12). Therefore, one can relate the emission of strongly coupled plasmon–molecule system to the properties of the plasmon and the molecules separately with a few parameters describing the plasmon–molecule coupling strength and their intrinsic decay rates. Furthermore, and can be determined theoretically from the molecule’s electronic states and its vibrational modes via .

Since only the terms with contribute to the inelastic emission, we see that for a given driving frequency the ratio of p- and s-polarized emission at is given by . For example, an interface of vacuum and silver in the Drude model results in and for the ratio at Yang et al. (2015); González-Tudela et al. (2013). This ratio is otherwise independent of the system as long as is finite.

In Fig. 4a we have plotted the elastic emission spectra from Eqs. (18). The s- and p-polarized emission are similar for small but for larger values the competition between the plasmon and molecule emission becomes more noticeable. In total, more p-polarized light is emitted due to the plasmon polarization. The ratio between p- and s-polarization power is controlled by the ratio between and and the detuning . The latter dependence we have plotted in Fig. 4b which shows that the ratio between p- and s-polarized light increases as a function of for the lower polariton branch, as observed in Baieva et al. (2017), and diminishes for the upper branch.

We find that the upper and lower polariton modes emit asymmetrically as in other approaches Neuman and Aizpurua (2018); del Pino et al. (2018); Herrera and Spano (2017); Zeb et al. (2017); Reitz et al. (2019) and experiments Baieva et al. (2017); Bellessa et al. (2004); Hakala et al. (2009); Baieva et al. (2012). This is caused by the asymmetry of the effective dissipation rate which is related to the molecule’s absorption spectrum. Analytical insight can be obtained when and which allows us to neglect processes that involve multiple vibrational quanta. From the response function we can find the polariton frequencies for zero detuning to be approximatively when . At these frequencies the dissipation rate in the first order of is given by

(19) |

Due to vibrations the dissipation rate of the upper polariton (upper sign, ) is larger than the dissipation rate of the lower polariton (), in agreement with femto-second pump-probe experiments on strongly coupled molecule-cavity systems Virgili et al. (2011). As a consequence the emission from the upper polariton is suprresed compared to the lower polariton, in line with experiment Chovan et al. (2008) and calculations Michetti and La Rocca (2008).

To summarize, we have constructed a model that allows describing the effect of vibrations on the strongly coupled stationary response of driven coupled light-matter modes. Depending on the case, one can either find the function describing the absorption and emission of vibrations in a given model system as in sup (), or relate the measured absorption and fluorescence of uncoupled molecules to . With small modifications, this approach can be extended also to the case of molecule-cavity systems Yu et al. (2019); Canaguier-Durand et al. (2013); Schwartz et al. (2013); Chovan et al. (2008); Virgili et al. (2011), plasmonic lattices Hakala et al. (2018) and/or higher-order correlation functions of the emitted light. Our quantum Langevin equation approach allows describing the stationary driven system, and hence it complements the often-used computational methods usually concentrating on transient response del Pino et al. (2018).

###### Acknowledgements.

We acknowledge the support from the Academy of Finland Center of Excellence program (project no. 284594) and project numbers 289947, 290677, and 317118.## References

- Törmä and Barnes (2014) P. Törmä and W. L. Barnes, Rep. Prog. Phys. 78, 013901 (2014).
- Yu et al. (2019) X. Yu, Y. Yuan, J. Xu, K.-T. Yong, J. Qu, and J. Song, Laser Photonics Rev. 13, 1800219 (2019).
- Hertzog et al. (2019) M. Hertzog, M. Wang, J. Mony, and K. Börjesson, Chem. Soc. Rev. 48, 937 (2019).
- Ribeiro et al. (2018) R. F. Ribeiro, L. A. Martínez-Martínez, M. Du, J. Campos-Gonzalez-Angulo, and J. Yuen-Zhou, Chem. Sci. 9, 6325 (2018).
- Feist et al. (2017) J. Feist, J. Galego, and F. J. Garcia-Vidal, ACS Photonics 5, 205 (2017).
- Baieva et al. (2017) S. Baieva, O. Hakamaa, G. Groenhof, T. T. Heikkilä, and J. J. Toppari, ACS Photonics 4, 28 (2017).
- Ingold and Nazarov (1992) G.-L. Ingold and Y. V. Nazarov, in Single charge tunneling (Springer, 1992) pp. 21–107.
- Heikkilä (2013) T. T. Heikkilä, The Physics of Nanoelectronics: Transport and Fluctuation Phenomena at Low Temperatures (Oxford University Press, 2013).
- Mahan (2013) G. D. Mahan, Many-Particle Physics (Springer Science & Business Media, 2013).
- Bellessa et al. (2004) J. Bellessa, C. Bonnand, J. C. Plenet, and J. Mugnier, Phys. Rev. Lett. 93, 036404 (2004).
- Hakala et al. (2009) T. K. Hakala, J. J. Toppari, A. Kuzyk, M. Pettersson, H. Tikkanen, H. Kunttu, and P. Törmä, Phys. Rev. Lett. 103, 053602 (2009).
- Baieva et al. (2012) S. V. Baieva, T. K. Hakala, and J. J. Toppari, Nanoscale Res. Lett. 7, 191 (2012).
- Koponen et al. (2013) M. A. Koponen, U. Hohenester, T. K. Hakala, and J. J. Toppari, Phys. Rev. B 88, 085425 (2013).
- Guebrou et al. (2012) S. A. Guebrou, C. Symonds, E. Homeyer, J. C. Plenet, Y. N. Gartstein, V. M. Agranovich, and J. Bellessa, Phys. Rev. Lett. 108, 066401 (2012).
- Walls and Milburn (2007) D. F. Walls and G. J. Milburn, Quantum Optics (Springer Science & Business Media, 2007).
- González-Tudela et al. (2013) A. González-Tudela, P. Huidobro, L. Martín-Moreno, C. Tejedor, and F. García-Vidal, Phys. Rev. Lett. 110, 126801 (2013).
- (17) See Supplemental Material for details on theory, approximation of Eq. (4b) and the effect of coherent driving on the spectra (Eqs. (12) and (18)). It also contains Refs. Giovannetti and Vitali (2001); Skellam (1946); Dicke (1954).
- Huang and Rhys (1950) K. Huang and A. Rhys, Proc. Royal Soc. A 204, 406 (1950).
- Glauber (1963) R. J. Glauber, Phys. Rev. 131, 2766 (1963).
- Gardiner and Collett (1985) C. W. Gardiner and M. J. Collett, Phys. Rev. A 31, 3761 (1985).
- Caldeira and Leggett (1981) A. O. Caldeira and A. J. Leggett, Phys. Rev. Lett. 46, 211 (1981).
- Yang et al. (2015) H. U. Yang, J. D’Archangel, M. L. Sundheimer, E. Tucker, G. D. Boreman, and M. B. Raschke, Phys. Rev. B 91, 235137 (2015).
- Neuman and Aizpurua (2018) T. Neuman and J. Aizpurua, Optica 5, 1247 (2018).
- del Pino et al. (2018) J. del Pino, F. A. Y. N. Schröder, A. W. Chin, J. Feist, and F. J. Garcia-Vidal, Phys. Rev. Lett. 121, 227401 (2018).
- Herrera and Spano (2017) F. Herrera and F. C. Spano, ACS Photonics 5, 65 (2017).
- Zeb et al. (2017) M. A. Zeb, P. G. Kirton, and J. Keeling, ACS Photonics 5, 249 (2017).
- Reitz et al. (2019) M. Reitz, C. Sommer, and C. Genes, (2019), arXiv:1812.08592 [quant-ph] .
- Virgili et al. (2011) T. Virgili, D. Coles, A. M. Adawi, C. Clark, P. Michetti, S. K. Rajendran, D. Brida, D. Polli, G. Cerullo, and D. G. Lidzey, PHYSICAL REVIEW B 83 (2011).
- Chovan et al. (2008) J. Chovan, I. E. Perakis, S. Ceccarelli, and D. G. Lidzey, PHYSICAL REVIEW B 78 (2008).
- Michetti and La Rocca (2008) P. Michetti and G. C. La Rocca, PHYSICAL REVIEW B 77 (2008).
- Canaguier-Durand et al. (2013) A. Canaguier-Durand, E. Devaux, J. George, Y. Pang, J. A. Hutchison, T. Schwartz, C. Genet, N. Wilhelms, J.-M. Lehn, and T. W. Ebbesen, Angew. Chem. Int. Ed. 52, 10533 (2013).
- Schwartz et al. (2013) T. Schwartz, J. A. Hutchison, J. Léonard, C. Genet, S. Haacke, and T. W. Ebbesen, ChemPhysChem 14, 125 (2013).
- Hakala et al. (2018) T. K. Hakala, A. J. Moilanen, A. I. Väkeväinen, R. Guo, J.-P. Martikainen, K. S. Daskalakis, H. T. Rekola, A. Julku, and P. Törmä, Nature Phys. 14, 739 (2018).
- Giovannetti and Vitali (2001) V. Giovannetti and D. Vitali, Phys. Rev. A 63, 023812 (2001).
- Skellam (1946) J. G. Skellam, J. Royal Stat. Soc. A 109, 296 (1946).
- Dicke (1954) R. H. Dicke, Phys. Rev. 93, 99 (1954).