LINK bonds
When performing a QM/MM simulation, one often wants to partition the system such that some covalent bonds cross the QM/MM boundary. These so-called 'link' bonds demand special attention in any QM/MM implementation. The link bonds are a critical aspect of the QM/MM method used here and a good understanding of the concepts is essential. In this section we describe how they are treated in ADF and we introduce the nomenclature that is used throughout the manual.
Figure 1-1 Example of QM/MM partitioning and details of the naming conventions adopted in ADF.
Figure 1-1a depicts a simple molecular system that has been divided into QM and MM regions by the dotted polygon. In this example there is only one link bond, or one covalent bond that traverses the QM-MM boundary (that cross the dotted polygon). When a covalent bond traverses the QM-MM boundary, the electronic system of the QM region must in some way be truncated across this bond. Several methods of dealing with this problem have been proposed in the literature. By far the most commonly adopted method, which was originally introduced by Singh and Kollman [6], involves capping the QM system with a 'dummy atom' or what we call a 'capping atom' (We use capping atom to avoid confusion with the dummy atoms used in ADF). Since the pioneering work of Singh and Kollman [6] many variations of the basic capping atom approach have evolved. In ADF, one of the adopted approaches is the one developed by Maseras and Morokuma which has been given the name the 'Integrated Molecular Orbital and Molecular Mechanics' or the IMOMM method by the authors. The key feature of the capping atom approach is that the electronic structure calculation is performed on what is referred to as the 'QM model' system where the MM region is removed and replaced with capping dummy atoms (often hydrogen but not necessarily so). Figure 1-1b depicts the QM model system corresponding to the example system presented in Figure 1-1a. The capping atom satisfies the valence requirements of the QM region and allows for a standard electronic structure calculation to be performed on the QM fragment. It is important to realize that the capping atom is not part of the real system, but is simply an atom that is introduced to truncate the electronic system of the QM region. This is why it is often referred to as the dummy atom. For every 'link' bond there are three atoms of importance, the capping atom and the two atoms that are part of the 'real' link bond - one from the QM region and one from the MM region. Figure 1-1c illustrates the three atoms involved in the link bond. From this point on, we will refer to the MM atom that is part of the 'real' link bond as the link atom; it is labeled 'LI' in Figure 1-1c. Although both the QM and MM atoms that are part of the link bond could be considered 'link' atoms, we designate only the MM atom as the link atom, because it has a special place in the ADF QM/MM input. It is this atom that is replaced by the capping atom in the electronic structure calculation of the QM model system.
Although the capping atom approach is convenient from the standpoint of the electronic structure calculation, the 'extra' capping atoms complicates the situation, as they do not exist in the real system, see Figure 1-1c. For each link bond, there are potentially three extra nuclear degrees of freedom (corresponding to the Cartesian coordinates of the capping atom) that are not present in the real system. In the IMOMM/ADF implementation [1] we alleviate the problem by removing the MM atom that is part of the 'real' link bond as a free variable. Instead we define its position in terms of the QM atom it is bonded to and the capping atom that replaces it in the QM model system. More specifically, the MM link atom is constrained to lie along the bond vector of the capping atom bond, via the simple relationship expressed in equation 1.1 and depicted in Figure 1-1d. Here, XMM, Xcap and XQM refer to the Cartesian coordinates of the subscripted atoms and a is a constant defined as the ratio of the real link bond length to that of the length of the capping bond.
XLI = a Xcap + (1 - a) XQM (1.1)
For each link bond, there is a unique a parameter that is held constant throughout the simulation. Since the capping atom is often at a shorter distance than the real MM atom, alpha is usually greater than unity. For example, when a Hydrogen capping atom is used to cap a C-C single bond, a is around 1.38. Note that the energy depends on the value chosen for this parameter a, and comparisons between two systems are meaningful only when the values for them are chosen equal.
Although the position of the MM atom is not an independent variable (or free degree of freedom), the bond length of the link bond can change during a geometry optimization. If the capping bond in the model QM system stretches or contracts, so does the link bond in the full system. Note that any forces exerted on the LI atom are projected onto the connected QM atom and onto the capping atom. For more details see references [1, 2].
In the AddRemove model [3] the position dependence is reversed. Instead of the real classical LI atom following the artificial capping atom, it is the position of the capping atom that is constrained:
Xcap = XQM + Req ULI-QM (1.2)
Here Req is the equilibrium distance (from the force field !) for the QM-cap bond, while ULI-QM is the unit vector pointing from QM to LI.
This means also that the artificial degrees of freedom of the capping atom are being removed ! Moreover, after every QM geometry optimization cycle, the interactions of the artificial capping atoms with the real QM atoms are corrected for by subtracting the corresponding MM interactions [3]. In this way, the capping atoms are Added and Removed. It has been shown [3] that this choice gives a much better description for the bonds involved in the QM/MM boundary region than the IMOMM/ADF scheme.
