QMMM

There are two QM/MM calculators native to ASE:

Explicit Interaction QMMM EIQMMM
Simple, subtrative QMMM SimpleQMMM

Explicit Interaction QMMM

In Explicit Interaction QMMM, the QM and MM regions are explicitly coupled with an electrostatic interaction term. This requires that the electrostatic potential from the classical charges of the MM subsystem is fed into the QM calculator. This is built into GPAW. More info In this paper, which should be cited if the method is used.

See also

The ASE for QM/MM Simulations tutorial, on how to use the Explicit Interaction QMMM calculator

class ase.calculators.qmmm.EIQMMM(selection, qmcalc, mmcalc, interaction, vacuum=None, embedding=None, output=None)[source]

Explicit interaction QMMM calculator.

EIQMMM object.

The energy is calculated as:

        _          _         _    _
E = E  (R  ) + E  (R  ) + E (R  , R  )
     QM  QM     MM  MM     I  QM   MM

parameters:

selection: list of int, slice object or list of bool
Selection out of all the atoms that belong to the QM part.
qmcalc: Calculator object
QM-calculator.
mmcalc: Calculator object
MM-calculator.
interaction: Interaction object
Interaction between QM and MM regions.
vacuum: float or None
Amount of vacuum to add around QM atoms. Use None if QM calculator doesn’t need a box.
embedding: Embedding object or None
Specialized embedding object. Use None in order to use the default one.
output: None, ‘-‘, str or file-descriptor.
File for logging information - default is no logging (None).

Here, you need to specify the interaction:

from ase.calculators.qmmm import EIQMMM, LJInteraction
from ase.calculators.tip3p import epsilon0, sigma0
lj = LJInteraction({'OO': (epsilon0, sigma0)})
atoms.calc = EIQMMM([0, 1, 2],
                    QMCalculator(...),
                    MMCalculator(...),
                    interaction=lj)

For Lennard-Jones type of interactions you can use:

class ase.calculators.qmmm.LJInteractions(parameters)[source]

Lennard-Jones type explicit interaction.

parameters: dict
Mapping from pair of atoms to tuple containing epsilon and sigma for that pair.

Example:

lj = LJInteractions({(‘O’, ‘O’): (eps, sigma)})

You can control how the QM part is embedded in the MM part by supplying your own embedding object when you construct the EIQMMM instance. The Embedding object will be specific to the QM calculator you want to use. The default is this one:

class ase.calculators.qmmm.Embedding(molecule_size=3, **parameters)[source]

Point-charge embedding.

Simple, subtractive QMMM calculations

This QM/MM calculator is similar to the original ONIOM model, doing simple, subtractive QM/MM between any two calculators.

class ase.calculators.qmmm.SimpleQMMM(selection, qmcalc, mmcalc1, mmcalc2, vacuum=None)[source]

Simple QMMM calculator.

SimpleQMMM object.

The energy is calculated as:

        _          _          _
E = E  (R  ) - E  (R  ) + E  (R   )
     QM  QM     MM  QM     MM  all

parameters:

selection: list of int, slice object or list of bool
Selection out of all the atoms that belong to the QM part.
qmcalc: Calculator object
QM-calculator.
mmcalc1: Calculator object
MM-calculator used for QM region.
mmcalc2: Calculator object
MM-calculator used for everything.
vacuum: float or None
Amount of vacuum to add around QM atoms. Use None if QM calculator doesn’t need a box.

This type of QMMM can combine any pair of ASE calculators:

from ase.calculators.qmmm import SimpleQMMM
atoms = ...
atoms.calc = SimpleQMMM([0, 1, 2],
                        QMCalculator(...),
                        MMCalculator(...))

where [0, 1, 2] would select the first three atoms for the QM-part.