Calculation of Linear and Non-linear Electric Response Properties of Systems in Aqueous Solution: A Polarizable Quantum/Classical Approach with Quantum Repulsion Effects

We present a computational study of polarizabilities and hyperpolarizabilities of organic molecules in aqueous solutions, focusing on solute–water interactions and the way they affect a molecule’s linear and non-linear electric response properties. We employ a polarizable quantum mechanics/molecular mechanics (QM/MM) computational model that treats the solute at the QM level while the solvent is treated classically using a force field that includes polarizable charges and dipoles, which dynamically respond to the solute’s quantum-mechanical electron density. Quantum confinement effects are also treated by means of a recently implemented method that endows solvent molecules with a parametric electron density, which exerts Pauli repulsion forces upon the solute. By applying the method to a set of aromatic molecules in solution we show that, for both polarizabilities and first hyperpolarizabilities, observed solution values are the result of a delicate balance between electrostatics, hydrogen-bonding, and non-electrostatic solute solvent interactions.


S0.1 Convergence of QM/MM calculations
The quality of QM/MM results rests on the assumption that the entire solute-solvent configuration space has been sampled by the MD and the number of extracted snapshots is sufficient. Final results are obtained by averaging the computed values across the whole set of snapshots. In Fig. S1 we report the average of each property (static polarizability α(0; 0), dynamic polarizability α(−ω; ω), and dynamic hyperpolarizability β(−2ω; ω, ω)) as a function of the number of snapshots, obtained at the CAM-B3LYP/FQ b level of theory.
We can see that for most molecules convergence is achieved within the first 100 snapshots.
The only exception to this seems to be molecule 1, for which we observe further oscillation past the 100 snapshot mark, and for which 200 snapshots are necessary to achieve convergence. Nonetheless, oscillations around the average values remain quite small and most of S1 the sampling is still achieved within the first 100 snapshots. These conclusions hold for all three properties analyzed, both with and without repulsion effects.  Figure S1: Average values of α(0; 0), α(−ω; ω) and β(−2ω; ω, ω) as functions of the number of snapshots exctracted from MD for all the molecules considered using the CAM-B3LYP/FQ b method. On the left, snapshots without repulsion, while, on the right, snapshots with repulsion.

S0.2 Tabulated QM/MM polarizabilities and hyperpolarizabilities
We report here tabulated numerical values corresponding to those presented in figures within the main manuscript. QM/FQ and QM/FQFµ data also come with an error bar given by the variance of the computed property across the snapshots extracted from the classical MD. Gas-phase values, by contrast, do not have error bars because they are just single-point calculations peformed on an optimized molecular geometry. Table S1 shows the static and dynamic polarizability of the six molecular systems studied, both in water in vacuo, with and without repulsion effects, with three different functionals. Table S2 shows

S0.3 Plotted density matrix derivatives
We report here a list of the electric field derivatives of the density matrix for a randomly selected snapshot of the molecule 1 (see main article text for more details) with and without the inclusion of repulsion effects. The density matrices are in the molecular orbital basis and shown in figures S2 and S3. It can be observed that derivatives along the x direction are less affected by repulsion compared to the other two direction. In addition, off-diagonal derivatives are also less influenced by repulsion.
P (1) rep−norep (x) P (1) rep−norep (y) Figure S2: Difference between the density matrix first derivatives with and without Pauli repulsion of a randomly snapshot of molecule 1 extracted from MD simulation. Derivatives are taken with respect to the components of the electric field indicated in parenthesis.