QM/MM¶

In QM/MM methods, the energy is partitioned into quantum mechanical (QM) and molecular mechanical (MM) contributions [79] corresponding to different parts of the system. Thus the QM/MM approach combines the accuracy of QM with the speed of MM. Typically, the QM subsystem is a small region where the interesting chemistry takes place, while the MM part corresponds to some larger environment, represented by simple point charges, surrounding the QM region. In general, there are two ways of expressing the total energy in QM/MM, “subtractive” and “additive” [80]. In InQuanto the additive approach is adopted, which means the total energy takes the form

()¶\[E_{\text{tot}} = E_{\text{QM}} + E_{\text{QM/MM}} + E_{\text{MM}}\]

where \(E_{\text{QM}}\) represents the contributions internal to the QM subsystem (calculated using high-level methods which include a two-body interaction or approximation thereof), \(E_{\text{MM}}\) comes from interactions between particles in the MM environment (a low cost classical term like the Coulomb interaction between point charges), and \(E_{\text{QM/MM}}\) represents the interactions between the QM and MM regions (also a classical interaction). The latter involves a modification of the one-body component of the electronic Hamiltonian. The modified one-body part of the Hamiltonian becomes

()¶\[h_{i,j} \rightarrow h_{i,j}^{\text{QM/MM}} = \int{d\mathbf{x}\phi_{i}^{*}(\mathbf{x})\bigg(\frac{-\nabla^{2}}{2}} - \sum_{A}\frac{Z_{A}}{|\mathbf{R}_{A} - \mathbf{r}|} - \sum_M \frac{q_M}{|\mathbf{R}_{M} - \mathbf{r}|} \bigg)\phi_{j}(\mathbf{x})\]

where \(\mathbf{x} = \{\mathbf{r}, s\}\) with \(\mathbf{r}\) (\(s\)) the electronic position (spin), and the modified classical Coulomb interaction for nuclei becomes

()¶\[V_{\text{nuc}} \rightarrow V_{\text{nuc}}^{\text{QM/MM}} = \sum_{A, B}\frac{Z_{A}Z_{B}}{|\mathbf{R}_{A} - \mathbf{R}_{B}|} + \sum_{A, M} \frac{Z_{A} q_M}{|\mathbf{R}_{A} - \mathbf{R}_{M}|}\]

where \(M\) labels the MM point charges, with charge \(q_M\) and position \(\mathbf{R}_{M}\).

In inquanto-pyscf the specialized driver class: ChemistryDriverPySCFMolecularRHFQMMMCOSMO, uses the QM/MM scheme of PySCF to modify the one-body integrals and nuclear Coulomb interaction. With this driver, QM/MM embedding is performed as follows; we consider an H₂ molecule surrounded by ten random point charges:

from inquanto.extensions.pyscf import ChemistryDriverPySCFMolecularRHFQMMMCOSMO
# define the MM region
mm_charges = 0.1  # this will set all MM particles to have charge = +0.1 in units of electron charge.
# lists and numpy arrays are also accepted to specify charges individually.
mm_geometry = [  # pregenerated random geometry using nunmpy.random.seed(1)
    [-0.49786797, 1.32194696, -2.99931375],
    [-1.18600456, -2.11946466, -2.44596843],
    [-1.88243873, -0.92663564, -0.61939515],
    [0.2329004, -0.48483291, 1.111317],
    [-1.7732865, 2.26870462, -2.83567444],
    [1.02280506, -0.49617119, 0.35213897],
    [-2.15767837, -1.81139107, 1.80446741],
    [2.80956945, -1.11945493, 1.15393569],
    [2.25833491, 2.36763998, -2.48973473],
    [-2.7656713, -1.98101748, 2.26885502],
]
# e_mm_coulomb should be 0.068638356203663

# Execute the chemistry drivers.
driver_qmmm = ChemistryDriverPySCFMolecularRHFQMMMCOSMO(
    zmatrix="""
        H
        H  1 0.7122
    """,
    basis="sto3g",
    mm_charges=mm_charges,
    mm_geometry=mm_geometry,
    do_mm_coulomb=True,
    do_qmmm=True,
)

The Hamiltonian returned by get_system() will now contain the modified one-body terms due to the QM/MM interaction. To calculate \(E_{\text{tot}}\) as above including the Coulomb interaction between MM point charges, the optional boolean keyword do_mm_coulomb must be set to True, which stores the MM Coulomb interaction energy \(E_{\text{MM}}\) as a driver property: e_mm_coulomb. However, this does not add \(E_{\text{MM}}\) to the total energy, which must be done manually to recover \(E_{\text{QM}} + E_{\text{QM/MM}} + E_{\text{MM}}\).

QM/MM can also be used with ROHF and UHF methods with the corresponding drivers ChemistryDriverPySCFMolecularROHFQMMMCOSMO and ChemistryDriverPySCFMolecularUHFQMMMCOSMO. In addition, QM/MM can be combined with other embedding schemes.