In Sec. 2, we discuss the rigorous theoretical background of molecular cavity QED. We first review the basic theory of the molecular Hamiltonian (Sec. 2.1) and introduce the necessary formalism for molecular quantum electrodynamics (Sec. 2.2). Sec.3 focuses on solving the polaritonic eigenvalue equation using various approaches from direct diagonalization in the adiabatic-Fock basis (Sec. 3.1.1) and polarized Fock basis (Sec. 3.1.2) as well as the selfconsistent solution at the mean-field (Sec. 3.2.2) and correlated (Sec. 3.2.3) levels of theory. Using the tools put forth in the previous two sections, Sec. 4 explores the rich quantum dynamics of strongly coupled light-matter systems, laying out common approaches toward simulating exact (Sec. 4.1) and realistic (Sec. 4.2.1-4.2.3) molecules in the extended Hilbert space of the coupled electron-photon system as well as motivating the need for an accurate description and efficient calculation of the nuclear gradients (Sec. 4.3) used in the exact Pauli-Fierz Hamiltonian. Finally, in Sec. 5, we provide three directions toward a more complete picture of the molecular polariton picture in experiments, including machine learning polariton gradients (Sec. 5.1), the simulation and effects of cavity loss (Sec. 5.2), and the extension to many molecule and many cavity modes (Sec. 5.3).