Competition between collective and individual conical intersection dynamics in an optical cavity

Light-induced nonadiabatic phenomena arise when molecules or molecular ensembles are exposed to resonant external electromagnetic fields. The latter can either be classical laser or quantized cavity radiation fields, which can couple to either the electronic, nuclear or rotational degrees of freedom of the molecule. In the case of quantized radiation fields, the light–matter coupling results in the formation of two new hybrid light–matter states, namely the upper and lower ‘polaritons’. Light-induced avoided crossings and light-induced conical intersections (CIs) between polaritons exist as a function of the vibrational and rotational coordinates of single molecules. For ensembles of N molecules, the N − 1 dark states between the two optically active polaritons feature, additionally, so-called collective CIs, involving the coordinates of more than one molecule to form. Here, we study the competition between intramolecular and collective light-induced nonadiabatic phenomena by comparing the escape rate from the Franck–Condon region of a single molecule and of a molecular ensemble coupled to a cavity mode. In situations where the polaritonic gap would be large and the dark-state decay channels could not be reached effectively, the presence of a seam of light-induced CI between the polaritons facilitates again the participation of the dark manifold, resulting in a cooperative effect that determines the overall non-radiative decay rate from the upper into the lower polaritonic states.


Introduction
Conical intersections (CIs) are degeneracy points between multidimensional potential energy surfaces (PES) which may be present under natural conditions in polyatomic molecules [1][2][3][4][5]. They play an essential role in the nonadiabatic photochemical and photophysical processes-such as ultrafast radiationless relaxation, photodissociation, photofragmentation or isomerization of polyatomic molecules, to name a few-providing highly efficient pathways for a significant energy exchange between the nuclei and electrons [4,[6][7][8][9]. In such a situation the well-known Born-Oppenheimer approximation [10] loses its validity.
Exposing the molecules to resonant laser light, a new phenomenon arises. This phenomenon is a CI induced by light (LICI, light-induced CI), which can be formed even in the case of diatomic molecules [11,12]. By varying the laser parameters (intensity and frequency), it is possible to modify the position and structure of the LICIs, and thus manipulate the dynamics so as to investigate how molecules behave near a LICI. Similarly to the case of natural degeneracies, LICIs can also give rise to a variety of nonadiabatic phenomena and several theoretical/experimental works have demonstrated that the LICIs provide considerable effect on different molecular properties [13][14][15][16][17][18][19][20][21].
LICIs can also be formed if classical laser light is replaced by the quantized electromagnetic field of an optical or plasmonic nanocavity [22][23][24][25]. In this situation, the confined photonic mode of the cavity can resonantly couple the electronic states of the molecule resulting in hybrid light-matter states carrying both photonic and excitonic properties, the so-called upper and lower polaritonic states. LICIs emerge in this case when, as a function of the vibrational or rotational coordinates of the molecule, the transition dipole moment between the cavity-coupled electronic states vanishes [22][23][24][25]. For diatomic molecules, the coupling coordinate of the LICI corresponds to the angle between the molecular axis and the light polarization. As a function of this angle the two polaritonic states either split into the hybridized light-matter states, or become degenerate at the angle for which the transition dipole vanishes. It is widely accepted that hybrid light-matter polaritonic states play an important role in steering and modifying various aspects of the reactivity and ultrafast dynamics of molecular systems, thus defining the novel research field of polaritonic chemistry .
So far, we have introduced CIs in single molecules, either natural or LICIs. Feist and coworkers pointed out that molecular ensembles coherently coupled to a cavity mode result in so-called collective CIs (CCI) in the dark-states manifold [36,38]. Assuming N molecules in the same configuration N − 1 dark states are formed. These dark states are superpositions of molecular excitons that are not coupled to the photonic mode. Thus their energy corresponds to the bare molecular exciton. Upon considering the nuclear motion, these N − 1 dark states become degenerate at the point of the CCI and the degeneracy is lifted along the displacement of any nuclear degree of freedom. Therefore, at least three molecules are needed so as to form a CCI. Interestingly, CCIs can be understood as special case of the well-known dynamical Jahn-Teller effect [38], which in isolated molecules is responsible for ultrafast non-radiative electronic relaxation processes. Previous studies have shown how the presence of CCIs can substantially modify the collective photo-dissociation dynamics of diatomic molecules as a function of the ensemble size [38], as well as affect the non-radiative decay through natural CIs in polyatomic molecules [46,54].
In diatomic molecules with a rotational degree of freedom, both individual and collective (for ensembles) LICIs appear as a function of the rotation angle. Moreover, CCI between dark states are present for N > 2 molecular ensembles. These different kinds of nonadiabatic coupling mechanisms can, separately, strongly affect the photodissociation dynamics in a dissociative excited electronic state [38,44]. In this paper, we compare the effect of individual and collective LICIs and of dark-states CCIs on the process of ultrafast photodissociation. In particular, we examine how the escape rate from the Franck-Condon (FC) region depend upon the size of the molecular ensemble, ranging from the N = 1 case, where only an individual LICI is present, to N > 2, where the collective LICI and dark-states CCI both play a role simultaneously. This work generalizes previous theoretical investigations where the photodissociation dynamics was investigated for CCI without rotational degrees of freedom [38], and in the vicinity of individual quantum LICIs [44]. In the model applied we neglect the permanent dipole as we expect a very small effect due to the large detuning of the cavity with pure rotational excitations.
We note here that neither natural nonadiabatic effects, nor cavity losses are treated in the paper as our primary aim is to understand how the interplay between the individual and collective light-induced nonadiabatic effects modifies the molecular dynamics. Possible experimental realizations of these effects may involve plasmonic nanocavities, with lifetimes of about 10-12 femtoseconds, or Fabry-Perot cavities with longer lifetimes. As shown in [63], lifetimes shorter than the time for significant nuclear displacements leave the molecular dynamics mostly unmodified, and this needs to be considered in possible experimental realizations. Theoretically, the interplay between LICI and CCI for lossy cavities shall be considered separately.
The structure of the paper is as follows. Section 2 provides the description of the (i) cavity-coupled molecular model, (ii) the applied numerical method for the solution of the time-dependent Schrödinger equation and (iii) the computed quantities. Section 3 presents and discusses the results of the quantum dynamics simulations, while conclusions and outlook are found in section 4.

The Hamiltonian
We consider an ensemble of molecules interacting with the quantized electromagnetic mode of a cavity. The density of molecules inside the cavity is considered to be low such that the direct interaction of the molecules is negligible. In such a situation the molecules are coupled to the cavity mode and to the external electric field, but not to each other. The total Hamiltonian of such a cavity-coupled molecular ensemble can be given as follows [32,37] (Atomic units used throughout the paper except otherwise stated.) where the individual Hamiltonian of the κth molecule,Ĥ (κ) mol is given as the sum of the kinetic energy operators for the nuclei and the electrons plus the Coulombic interaction terms for the electron-electron repulsion, electron-nuclei attraction and nuclei-nuclei repulsion, respectivelŷ In equation (1), the cavity Hamiltonian has the form [32,37] whereâ † andâ are the photon creation and annihilation operators, respectively, while ω c is the angular frequency of the cavity mode, ε c is the polarization vector, and g is the molecule-cavity coupling strength and it is defined in terms of ω c , g = α * ω c (here the α factor has the physical dimension of inverse dipole moment, but since atomic units are applied throughout, it is not written out explicitly). The total dipole operator of the ensemble is given as the sum of the individual molecular dipoles,ˆ D = N κ=1 μ (κ) . The last term in equation (1) describes the interaction of an external electric field with the molecules in the dipole approximationĤ where the actual form of the electric field is considered as In equation (5), E 0 is the electric field amplitude, ε L is the polarization unit vector, ω L is the angular frequency and f (t) = e −2ln(2)(t−t 0 ) 2 /τ 2 is the envelope function with the τ pulse duration parameter. For the sake of simplicity, we will assume that the external field (pulsed laser) is not coupled with the cavity. A typical pulse duration of τ = 20 fs will be applied to pump the molecules and initiate the dynamics of the hybrid cavity-molecular ensemble system. For an explicit demonstration, the widely-studied sodium iodine (NaI) molecule will be considered throughout this work. This molecule has been the subject of investigation on cavity-induced nonadiabaticity in former works authored by one of us [38]. For consistency and to avoid the interplay between the natural and light-induced nonadiabatic phenomena, we apply the same model of NaI (for further details, see reference [37]). The original model [38] corresponds to an ensemble of vibrating molecules aligned with the cavity mode (termed as aligned). Individual light-induced nonadiabatic effects, though, require extending the previous model by incorporating the dynamical rotation of each molecule in the ensemble. This extended model describes rotating-vibrating molecules (termed as 2D) and also vibrating molecules averaged over different fixed molecule-cavity orientations (termed as 1D). In the extended model, the inclusion of molecular rotation implies the addition of the angular momentum operator term in the kinetic energy operator of the nuclei, where M κ is the reduced mass, R κ is the internuclear coordinate, whileL κ is the angular momentum operator of the κth molecule.

Numerical propagation
As in previous works [38], the time-dependent Schrödinger equation of the hybrid cavity-ensemble system is solved by the very efficient multi-configuration time-dependent Hartree (MCTDH) method [64][65][66]. The total time-dependent MCTDH wave function of the hybrid system with N molecules has the form where n l and n p are the number of single-particle functions for the lth molecule and for the photonic mode, respectively. Mode combination is applied to combine the electronic and nuclear (vibrational and rotational) DOFs of each individual molecule into a single mode. N s denotes the number of electronic states, while N p is the maximal number of photons inside the cavity in the primitive basis representation. φ (l) s l ,j l (t) is the nuclear wave packet of the lth molecule in the s l electronic state corresponding to the j l configuration space index. B P,j p (t) is the expansion coefficient for the |P photonic state related to the j p configuration space index. The vibrational DOF was discretized in a Fourier basis, the rotational coordinate was described with Legendre polynomials, while for the photonic mode a harmonic oscillator basis was utilized. Further details on the MCTDH treatment of the ensemble-cavity system can be found in reference [37].
The actual values of the MCTDH wave function parameters varied depending on the specific problem. Typical values were: n l = 10-40, n p = 5-15, N s = 2, N p = 41. The R grid was discretized on 2048 points, while the θ grid included 121 points. These basis set parameters ensured numerically converged wave packet propagations.

Calculated quantities
The nonadiabatic relaxation dynamics of the photo-excited hybrid system is traced through the rate of vanishing the wavefunction from the FC region corresponding to a single molecule. For that purpose, single-molecule dissociation and excitation probabilities have been calculated as described below. The single-molecule excitation probability was calculated according to [38] where N is the number of coupled molecules and P tot gs is the probability that the ensemble of N molecules in the cavity remains in the absolute ground state after the photo-pumping process. P tot gs is calculated as the absolute square of the autocorrelation function some time (t * ) after the pump pulse is over, The probability that any of the molecules has dissociated independently from the other members of the ensemble can be calculated as escape probability (EP) using the propagated wave packet. This amount can be defined as the loss of overlap with the initial wavefunction. For convenience, the EP is calculated by setting a CAP with the only requirement that it is zero on the FC region: Here Θ(R) is the Heaviside step function and R d = 10 a.u. is the starting point of the escape region. Finally, our quantity of interest is the single-molecule EP relative to the single-molecule excitation probability, P esc /P sm ex . We note that the ω L central angular frequency of the pump laser was tuned such that the collective Rabi splitting, Ω R = 2gμ √ N of the lower and upper polaritonic branch (UPB) is compensated and the P sm ex single-molecule excitation probability remain more-or-less constant for the different ensemble sizes.

Results and discussion
Let us start by considering a single NaI molecule and then an ensemble of molecules placed into a cavity, aligned with the direction of the electromagnetic mode. We assume that a direct interaction between the molecules is negligible, rather the molecules interact with the quantized electromagnetic mode of the cavity. Initially, all the molecules and also the cavity are in the absolute ground state. The angular frequency of the cavity mode is chosen such that it resonantly couples (ω c = 3.81 eV) the molecular ground electronic state with the first singlet excited electronic state, which is dissociative. A single laser pulse of E 0 = 0.001 a.u. field strength and τ = 20 fs duration coherently and resonantly excites the hybrid cavity-matter system into the UPB (see figure 1), thus above all CIs. Since the collective Rabi splitting depends upon the number of molecules coupled to the cavity, the ω L central angular frequency of the laser is adjusted in each case such that the single-molecule excitation probability (P sm ex ) to the UPB remains roughly constant for the different ensembles.
By solving the TDSE of the hybrid cavity-matter system, characterized by the Hamiltonian in equation (1), one can study the cavity-induced nonadiabatic relaxation dynamics of the molecular ensemble. The results for aligned molecules are presented in figure 2 for different regimes of the g cavity-matter coupling strength. In the weak field limit, the excited molecules almost completely dissociate on PES of the lower polaritonic branch (LPB) following an ultrafast decay from the UPB, as shown by the final values of the single-molecule escape probabilities in figure 2(a). For low g values, the final EP per molecule is almost independent from the number of coupled molecules. The small amount of population that is prevented from dissociation is due to the shallow well on the excited state of NaI (see figure 1(a)). At the strongest considered couplings, on the other hand, the rate of escape from the FC region (specifically the EP at the final simulation time) is much smaller irrespective of the size of the ensemble. In the intermediate regime, the final EP substantially depends on the size of the molecular ensemble and on the coupling strength. The UPB survival time is largest when a single molecule (N = 1) is coupled to the cavity. The UPB-LPB avoided crossing prevents an efficient non-radiative decay in this case. For N > 1, instead, the final EP (and rate) increases with the number of molecules for all coupling strengths. This is due to the new decay channels through the dark states, which are absent for a single molecule. As the coupling strength increases for a specific ensemble size, however, the EP decreases due to the increasing gap between the bright UPB and LPB.
The time-dependent single-molecule escape probabilities for three different g (g = 0.005ω c , 0.01ω c and 0.015ω c ) values are shown in panels (b)-(d) of figure 2. For g = 0.005ω c , due to the small gap between UPB and LPB, the molecules quickly escape for all N after the pumping process ( figure 2(b)). At stronger coupling in the intermediate regime (g = 0.01ω c ), the trapping of the population in the UPB becomes more prominent. The decay rate, and thus the time needed for the complete escape, decreases in direct correlation with increasing ensemble size (figure 2(c)). As already mentioned, this is directly related to the increasing number of non-radiative decay channels through the dark states as the number of molecules increases [38]. Further increasing the cavity coupling strength, the trapping of the molecular-photonic wave packet in the UPB starts to dominate the escape process ( figure 2(d)). At this coupling strength, an inversion of the decay rate with the increasing number of molecules can be observed between N = 2 and N = 3 molecules in figure 2(d). As the number of molecules increases, the number channels through dark states increases as well, thus increasing the decay rate. However, this trend is opposed by the increasing gap between the UPB and LPB as N increases, which eventually reverts the situation and slows down the nonadiabatic decay again. For smaller couplings (e.g. g = 0.005ω c or g = 0.01ω c ) this inverted regime would appear for larger ensemble sizes than considered in figures 2(b) and (c).
So far, we have considered 1D molecules aligned with the polarization direction of the cavity mode and discussed their collective nonadiabatic effects. These nonadiabatic dynamics are caused by the modulation of the gap between ground and excited electronic states as a function of the vibrational DOF of each molecule in the ensemble [38]. As already mentioned, a different kind of nonadiabatic effects occur in rotating diatomics coupled to either classical [20] or quantum light [44]. Here, the nonadiabatic effects are related to the dependence of the transition dipole moment with the rotation angle of the molecules and can be very strong already for a single monomer. For a fixed molecular orientation one obtains a LIAC (1D molecules) as a function of the internuclear distance, whereas the consideration of the θ rotation as a dynamical DOF results in a LICI between the UPB and LPB (2D molecules). Specifically, the UPB-LPB LICI occurs for both a single molecule, and when all N molecules are perpendicular to the polarization direction of the cavity mode. A clear signature of the LICI for a single molecule coupled to the cavity is seen in figure 3(a). Here, we consider the stronger coupling, g = 0.015ω c . At this coupling strength, a single, fixed molecule, (N = 1, black dashed line) has a low escape rate due to the large UPB-LPB gap. When considering the rotational degree of freedom in the 2D model, which modulates the UPB-LPB gap and constitutes the coupling coordinate of the LICI, the molecule efficiently escape from the FC region (blue line with squares). In contrast, the 1D case (red line with circles; average over molecular orientations), where only a LIAC is formed, features a somewhat slower escape process. The difference between the blue and red curves represents the dynamical contribution of the LICI compared to the static average over all possible molecular orientations [44].
The addition of further molecules into the cavity (figure 3(b)) greatly enhances the escape rate, i.e. the non-radiative decay from the UPB. For N = 2 we see that the escape rate of 2D and 1D (averaged) molecules is faster than for two fixed molecules. The escape process is also faster than for N = 1 2D (and 1D) molecules (cf figure 3(a)). Hence, the dynamics is determined by a combination of both types of nonadiabatic effects. The rotational DOF determines the light-matter coupling and thus the UPB-LPB gap. This way, the molecules have access to configurations with a smaller polaritonic gap, which is required for an effective non-radiative transition. This effect is also active for N = 1. With N > 1, though, the decay channels through the dark states are present and further enhance the nonadiabatic relaxation towards the LPB. Only this combined effect can explain the faster rate of the blue and red curves in figure 3(b) compared to the black curve in figure 3(b) and the curves in figure 3(a). When considering one further molecule, N = 3 in figure 3(c), the combined effect results in slightly faster decay. More importantly, the presence of more decay channels through dark states washes out the difference between the dynamic and static treatments of the rotational coordinate, i.e. the difference between considering collective LICIs or rotationally averaged LIACs. In any case, allowing different orientations of the diatomic molecules in the cavity, either dynamically or in an averaged way, leads to significantly faster non-radiative decay from the UPB and escape rate than found for fixed, aligned molecules.
Finally, we directly compare the non-radiative decay for a single 2D molecule forming an individual LICI, with N = 3 fixed molecules forming a CCI between dark states. In the weak (g = 0.005ω c ), and even in the intermediate (g = 0.01ω c ) coupling regime, the short-time dynamics via the collective and individual  1, 2, 3, 4), calculated according to equations (8) and (9) for aligned molecules (black dashed lines), rotating-vibrating molecules (blue lines with squares) and vibrating molecules averaged over different fixed orientations of the cavity mode and the molecular dipole (red lines with dots). The trapping of the molecular-photonic wave packet is efficient for aligned molecules. On the other hand, a remarkably faster escape process occurs when different molecule-cavity orientations are allowed. When a single molecule is coupled to the cavity (panel (a)), the strong nonadiabaticity induced by the cavity allows for a significantly fast decay via the created LICI in 2D (blue line with squares). Upon increasing the number of coupled molecules, the 1D vs 2D differences are washed out as the collective nonadiabatic effects start to dominate over the individual ones. The considered cavity-matter coupling strength is g = 0.015ω c (ω c = 3.81 eV).  (8) and (9) for N = 3 aligned molecules (green lines) and for a single rotating-vibrating molecule (blue lines with square symbols). Three different cavity-matter coupling strengths are considered: g = 0.005ω c (dashed lines), g = 0.01ω c (solid lines), and g = 0.015ω c (dotted lines). The cavity mode frequency is resonant at the FC point of NaI (ω c = 3.81 eV).
CI is quantitatively similar in terms of single-molecule escape probabilities (cf figure 4). The differences become more noticeable in the strong coupling regime (g = 0.015ω c ). Here, the fixed molecules coupled to the cavity result in a large UPB-LPB gap that prevents the non-radiative decay from the UPB. Only the addition of more molecules could lead to a faster decay, as previously discussed. Comparatively, the single 2D molecule still has access to the LICI through the rotational coordinate and therefore it can still efficiently decay even if the light-matter coupling g for the aligned geometry is large. As already demonstrated, the fastest non-radiative decay and subsequent photodissociation occurs when both types of nonadiabatic effects operate simultaneously.

Conclusions
In summary, we have investigated how the escape dynamics of ensembles of diatomic molecules coupled to a cavity mode and initially excited to the upper polaritonic state depend upon the size of the molecular ensemble. For a single molecule, only an individual LICI between the upper and lower polaritons is present. Besides a LICI, decay channels through dark states are also active for N > 1 and, for N > 2, the dark-states manifold features CCIs. Both types of nonadiabatic coupling mechanisms, namely collective LICIs and dark-states CCI, cooperatively determine the decay rate from the UPB to the LPB and thus the escape rate on the LPB PES. While only considering the role of CCIs in fixed molecules, we showed how, by increasing N, the number of the decay channels through dark states increases as well, thus leading to an increased non-radiative decay rate. However, this trend is opposed by the increasing gap between the UPB and LPB as N increases, which eventually reverts the situation and slows down the nonadiabatic decay rate again. When molecules are allowed to rotate, the gap between the UPB and LPB can be modulated by the orientation of the molecules. This provides a mechanism to efficiently access the dark-states manifold again, even if the polaritonic gap would be large for aligned molecules. The overall nonadiabatic decay rate when all effects are considered is larger than for multiple fixed and aligned molecules, and is also larger than for a single 2D rotator, demonstrating the involvement of a cooperative effect between both kinds of light-induced nonadiabatic effects.
We hope that our findings will stimulate laser as well as photochemical cavity experiments. Due to tremendous advances in light generation technology recent development in pump-probe optical spectroscopy might be able to follow the early time escape dynamics from the FC region [28,67,68].