The trajectories started from the |e⟩ state for dynamics outside the cavity and from the third polaritonic state |E 2 ⟩ for the polariton dynamics when coupling the molecule inside the cavity. For all dynamics, the nuclear motion was propagated using the velocity-Verlet algorithm with a time step of dt = 0.2 fs. In each nuclear time step, 100 steps are used for the quantum subsystem propagation using Eq. 27a and Eq. 27b. All of the energy, gradient and derivative couplings of the molecules are computed on the fly. All of the population dynamics results are obtained by averaging over 500 trajectories, using the expression of Eq. 44 (for adiabatic-Fock states) and Eq. 46 (for polariton states). The decoherence correction algorithm proposed by Granucci et al. was employed in the GFSH dynamics and the decoherence parameter α was set to 0.1.