Theor Chem Acc (2015) 134:132
1 3
one hand, and toward the delocalized region on the other
hand. These latter hybrids contribute to the basis of the
delocalized central subsystem. This computational scheme
was shown to be appropriate to study various properties of
extended systems like rotational barriers, protonation energies and conformational energy differences [ 1 – 3 ]. Furthermore, following the idea of Warshel and Levitt [ 4 ] strictly
localized orbitals and predefi ned hybrids provided means
for separating the covalently bound subsystems in several
QM/MM schemes. In the local self-consistent fi eld (LSCF)
method [ 5 ] an atom on the boundary has a single hybrid
orbital that forms a strictly localized molecular orbital with
an atom in the QM subsystem, while other electrons of the
atom are not treated explicitly; the interaction of this atom
with atoms in the MM subsystem is described by a classical force fi eld. The strictly localized molecular orbital is
not included in the self-consistent fi eld procedure, rather
its coeffi cients are kept frozen at some predefi ned values. A
similar separation with frozen localized orbitals but using
a different orbital equation and specifi c parametrization is
proposed by Friesner et al. [ 6 , 7 ]. Another related approach,
the generalized hybrid orbital (GHO) method [ 8 , 9 ] puts
hybrid orbitals on the QM boundary atom and includes the
hybrid directed toward the QM subsystem in the SCF procedure, while other hybrids are kept fi xed.
Other ways to separate covalently bound QM and MM
subsystems have also been proposed, and most notable
among them is the link atom approach. It cuts the bonds at
the QM–MM interface and saturates the dangling bonds of
the QM subsystem typically by hydrogen atoms. The principal advantage of the link atom approach is that its computational realization can be reduced to almost standard QM
and MM calculations that greatly facilitates its implementation into existing codes. On the other hand, the link atom
approach introduces extra atoms not present in the original
system, and this creates artifacts that requires special considerations in the computation. Nevertheless, early comparisons of the link atom and the LSCF methods concluded
that both can give valuable results when applied with precaution [ 10 ].
Early QM/MM methods were introduced at a semiempirical level, and later they were generalized to ab initio and
DFT schemes. Such a generalization is typically straightforward for the link atom approach but poses diffi culties for
methods using frozen localized orbitals. The principal diffi culty is that these methods assume orthogonality between
the frozen and the optimized orbitals, and this is not automatically guaranteed at ab initio level. It is not straightforward to set up orbital equations that include the effect of
frozen orbitals and give solutions that are orthogonal to
these orbitals. A possible approach to circumvent this diffi culty is to orthogonalize the basis functions of the orbitals
to be optimized to the frozen orbitals [ 9 , 11 ]. Alternative
propositions include the neglect of the overlap between the
frozen and optimized orbitals [ 9 ] or the explicit inclusion of
the orbital overlap [ 12 ]. The orthogonality requirement can
also be in included in the orbital equation. This was realized in Ref. [ 6 ] in the framework of a specifi c model with
frozen localized orbitals. We have recently proposed [ 13 ]
the application of the Huzinaga equation [ 14 ] to ensure
the orthogonality of the optimized and frozen orbitals.
Although the Huzinaga equation was originally proposed to
calculate valence orbitals orthogonal to frozen core orbitals
[ 15 ], it was shown [ 13 ] that the equation is well suited to
calculate QM/MM wave functions with optimized orbitals
orthogonal to the strictly localized frozen orbitals separating the covalently bound QM and MM subsystems. Having
established the validity of this approach, it is highly desirable to implement it in a computer code that, fi rst, extends
its applicability with various techniques to a wide range of
problems and, second, makes it conveniently available for
the scientifi c community.
We selected the AMBER [ 16 , 17 ] molecular mechanics
and the MRCC [ 18 ] quantum chemical codes to develop a
versatile and user-friendly QM/MM program which is
freely available for academic purposes. The choice of
AMBER is motivated by its capabilities that include molecular mechanics, molecular dynamics, and techniques that
allow the effi cient calculation of free energies. Since the
QM/MM program is organized in a way that the MM code
drives the calculation and the input instructions are primarily based on the MM code, the widespread use of AMBER
ensures that the resulting QM/MM code can be used with
reasonable extra effort for those familiar with AMBER . It is
also to be noted that although force fi elds in neither AMBER ,
nor the AMBER – MRCC code restrict the use of force fi elds to
the AMBER force fi eld, its application in a QM/MM code is
expected to be advantageous owing to the charge derivation
scheme that guaranties a sensible electrostatic potential in
the QM region. The MRCC code was selected for the QM
region owing to its highly effi cient implementation of local
correlation techniques that are benefi cial for calculating
reaction mechanisms, a fi eld with primary importance for
QM/MM applications.
It has to be noted that the link atom method and the frozen localized orbital method require signifi cantly different
approaches in their implementation into existing computer
codes. The link atom method basically performs standard MM and QM calculations, and it is the generation of
the subsystems and the communication between the QM
and MM programs that require special attention and coding. This is the reason, while interfaces connecting various MM and QM programs [ 19 – 24 ] can be developed, and
only small modifi cations in the original codes are necessary
to be introduced. The situation is different for the frozen
localized orbital method where specifi c orbital equations
134
Reprinted from the journal
1 3
one hand, and toward the delocalized region on the other
hand. These latter hybrids contribute to the basis of the
delocalized central subsystem. This computational scheme
was shown to be appropriate to study various properties of
extended systems like rotational barriers, protonation energies and conformational energy differences [ 1 – 3 ]. Furthermore, following the idea of Warshel and Levitt [ 4 ] strictly
localized orbitals and predefi ned hybrids provided means
for separating the covalently bound subsystems in several
QM/MM schemes. In the local self-consistent fi eld (LSCF)
method [ 5 ] an atom on the boundary has a single hybrid
orbital that forms a strictly localized molecular orbital with
an atom in the QM subsystem, while other electrons of the
atom are not treated explicitly; the interaction of this atom
with atoms in the MM subsystem is described by a classical force fi eld. The strictly localized molecular orbital is
not included in the self-consistent fi eld procedure, rather
its coeffi cients are kept frozen at some predefi ned values. A
similar separation with frozen localized orbitals but using
a different orbital equation and specifi c parametrization is
proposed by Friesner et al. [ 6 , 7 ]. Another related approach,
the generalized hybrid orbital (GHO) method [ 8 , 9 ] puts
hybrid orbitals on the QM boundary atom and includes the
hybrid directed toward the QM subsystem in the SCF procedure, while other hybrids are kept fi xed.
Other ways to separate covalently bound QM and MM
subsystems have also been proposed, and most notable
among them is the link atom approach. It cuts the bonds at
the QM–MM interface and saturates the dangling bonds of
the QM subsystem typically by hydrogen atoms. The principal advantage of the link atom approach is that its computational realization can be reduced to almost standard QM
and MM calculations that greatly facilitates its implementation into existing codes. On the other hand, the link atom
approach introduces extra atoms not present in the original
system, and this creates artifacts that requires special considerations in the computation. Nevertheless, early comparisons of the link atom and the LSCF methods concluded
that both can give valuable results when applied with precaution [ 10 ].
Early QM/MM methods were introduced at a semiempirical level, and later they were generalized to ab initio and
DFT schemes. Such a generalization is typically straightforward for the link atom approach but poses diffi culties for
methods using frozen localized orbitals. The principal diffi culty is that these methods assume orthogonality between
the frozen and the optimized orbitals, and this is not automatically guaranteed at ab initio level. It is not straightforward to set up orbital equations that include the effect of
frozen orbitals and give solutions that are orthogonal to
these orbitals. A possible approach to circumvent this diffi culty is to orthogonalize the basis functions of the orbitals
to be optimized to the frozen orbitals [ 9 , 11 ]. Alternative
propositions include the neglect of the overlap between the
frozen and optimized orbitals [ 9 ] or the explicit inclusion of
the orbital overlap [ 12 ]. The orthogonality requirement can
also be in included in the orbital equation. This was realized in Ref. [ 6 ] in the framework of a specifi c model with
frozen localized orbitals. We have recently proposed [ 13 ]
the application of the Huzinaga equation [ 14 ] to ensure
the orthogonality of the optimized and frozen orbitals.
Although the Huzinaga equation was originally proposed to
calculate valence orbitals orthogonal to frozen core orbitals
[ 15 ], it was shown [ 13 ] that the equation is well suited to
calculate QM/MM wave functions with optimized orbitals
orthogonal to the strictly localized frozen orbitals separating the covalently bound QM and MM subsystems. Having
established the validity of this approach, it is highly desirable to implement it in a computer code that, fi rst, extends
its applicability with various techniques to a wide range of
problems and, second, makes it conveniently available for
the scientifi c community.
We selected the AMBER [ 16 , 17 ] molecular mechanics
and the MRCC [ 18 ] quantum chemical codes to develop a
versatile and user-friendly QM/MM program which is
freely available for academic purposes. The choice of
AMBER is motivated by its capabilities that include molecular mechanics, molecular dynamics, and techniques that
allow the effi cient calculation of free energies. Since the
QM/MM program is organized in a way that the MM code
drives the calculation and the input instructions are primarily based on the MM code, the widespread use of AMBER
ensures that the resulting QM/MM code can be used with
reasonable extra effort for those familiar with AMBER . It is
also to be noted that although force fi elds in neither AMBER ,
nor the AMBER – MRCC code restrict the use of force fi elds to
the AMBER force fi eld, its application in a QM/MM code is
expected to be advantageous owing to the charge derivation
scheme that guaranties a sensible electrostatic potential in
the QM region. The MRCC code was selected for the QM
region owing to its highly effi cient implementation of local
correlation techniques that are benefi cial for calculating
reaction mechanisms, a fi eld with primary importance for
QM/MM applications.
It has to be noted that the link atom method and the frozen localized orbital method require signifi cantly different
approaches in their implementation into existing computer
codes. The link atom method basically performs standard MM and QM calculations, and it is the generation of
the subsystems and the communication between the QM
and MM programs that require special attention and coding. This is the reason, while interfaces connecting various MM and QM programs [ 19 – 24 ] can be developed, and
only small modifi cations in the original codes are necessary
to be introduced. The situation is different for the frozen
localized orbital method where specifi c orbital equations
134
Reprinted from the journal