The equilibrated coordinates were then converted to the Drude model in CHARMM as described below. Finally, Drude systems were equilibrated for 1 ns. An isotropic Monte Carlo barostat was used to maintain the pressure at 1 bar by attempting box scaling every 25 integration steps. Position restraints were applied to heavy atoms during equilibration with a force constant of 800 kJ mol -1 nm 2 . For the solvated and in vacuo systems of both additive and Drude models, unbiased simulations were performed for 100 ns, saving coordinates every 10 ps. While the production simulations of solvated systems were carried out under an NPT ensemble, the in vacuo simulations employed an NVT ensemble. In all systems, a weak position restraining force of 50 kJ mol -1 nm 2 was applied to the bottom layer of the crystal film to facilitate post-processing reimaging and analyses.