Multi-scale Molecular Dynamics¶

The design of the ForceJob class allows for flexible extension of its behavior, while at the same time keeping the client code unaware of its nature: it can either act as a simple wrapper for ADF programs, or it can be a more complex orchestrating class, combining simpler ForceJob classes to implement multi-scale strategies. One application of this extensible design can be found in the QMMMForceJob object, which combines a QM and an MM method in an IMOMM-type scheme (mechanical embedding). The QMMMForceJob object is assigned two other ForceJob objects, the first representing the high-resolution calculation (QM), while the other represents the low resolution (MM). Both ForceJob objects contain an MDMolecule object for the full molecular system. The selection of the QM-region is handled by the QMMMForceJob, which contains the information about the part of the molecule that constitutes the QM region. When forces are requested from the QMMMForceJob, the following behavior is orchestrated: first, a MM force calculation is performed on the full system; then, the QM-region is selected, a QM calculation is executed solely for that region, and energy and forces are added to those from the full system MM calculation. Finally, an MM calculation is computed for the small QM-region, and the energy and forces are subtracted, yielding the final result, returned to the invoker. In symbols:

EQM/MM(Full) = EMM(Full) + EQM(QMRegion) – EMM(QMRegion)

The QMMMForceJob handles periodic boundary conditions if the low-level (MM) method supports this feature (i.e. NAMD). Whether the periodic interaction of the QM region with itself is handled at high or low resolution depends on the method used for the QM calculation. An example of QM/MM MD calculations can be found in the examples directory examples/scmlib/qmmm_dftbUFF_2h2o. The QMMMForceJob allows the use of link atoms when the QM boundary cuts through covalent bonds. However, this feature comes at the price of an increased script complexity. An example of a QMMM link-atom MD simulation is provided in the examples directory, under examples/scmlib/qmmm_linkatom_dftbNAMD_glutamate.

For more complex multi-scale calculations the HybridForceJob class can be used. This class allows the combination of a large set of different ForceJobs, each of them describing either the same, or different molecular systems. Each ForceJob can either involve a calculation on the full MDMolecule object it contains, or restricted to a specified region of the corresponding molecule. The forces from each contributing ForceJob can either be added or subtracted from the total force according to user preference, as specified at construction of the HybridForceJob object.

In order to perform QM/MM simulations on chemical reactivity in solution, it is important that the description of the solvent molecules can change on the fly, as the molecules move towards or away from the reactive region. To facilitate this, an AdaptiveQMMMForceJob class is available to provide adaptive QM/MM simulations using several available schemes, as described by Bulo et al.[4] and P. Fleurat-Lessard et al.[8] In these schemes, the description of the diffusing molecules changes gradually from QM to MM and vice versa, based on the distance of those molecules to a predefined reactive site. Various schemes are available for assigning the QM and MM character of the molecules. The class contains a QMMMForceJob object, as well as a partitioning object that assigns the partial QM and MM character to the molecules. An examples python script for such an adaptive QM/MM simulation, using the DAS [4] method, is provided in the examples directory, examples/scmlib/adqmmm_mopacscmUFF_h2o.