11.3 Stand-Alone QM/MM Calculations

11.3.4 QM/MM Examples

(June 30, 2021)

Features of this job:

  • Geometry optimization using ONIOM mechanical embedding.

  • MM region (water 1) described using OPLSAA.

  • QM region (water 2) described using PBE0/6-31G*.

  • $molecule input section contains user-defined MM bonds. A zero is used as a placeholder if there are no more connections.

Example 11.22  ONIOM optimization of water dimer.

$molecule
   0 1
   O   -0.790909    1.149780    0.907453  186  2  3  0  0
   H   -1.628044    1.245320    1.376372  187  1  0  0  0
   H   -0.669346    1.913705    0.331002  187  1  0  0  0
   O    1.178001   -0.686227    0.841306  186  5  6  0  0
   H    0.870001   -1.337091    1.468215  187  4  0  0  0
   H    0.472696   -0.008397    0.851892  187  4  0  0  0
$end

$rem
   JOBTYPE           opt
   METHOD            pbe0
   BASIS             6-31G*
   QM_MM_INTERFACE   oniom
   FORCE_FIELD       oplsaa
   USER_CONNECT      true
   MOLDEN_FORMAT     true
$end

$qm_atoms
   4 5 6
$end

View output

Features of this job:

  • Janus electronic embedding with a YingYang link atom (the glycosidic carbon at the C1 position of the deoxyribose).

  • MM region (deoxyribose) is described using AMBER99.

  • QM region (adenine) is described using HF/6-31G*.

  • The first 5 electronically excited states are computed with CIS. MM energy interactions between a QM atom and an MM atom (e.g., van der Waals interactions, as well as angles involving a single QM atom) are assumed to be the same in the excited states as in the ground state.

  • $molecule input section contains user-defined MM bonds.

  • Gaussian-blurred charges are used on all MM atoms, with a width set to 1.5 Å.

Example 11.23  Excited-state single-point QM/MM calculation on deoxyadenosine.

$molecule
   0 1
   O     0.000000     0.000000     0.000000   1244      2      9      0      0
   C     0.000000     0.000000     1.440000   1118      1      3     10     11
   C     1.427423     0.000000     1.962363   1121      2      4      6     12
   O     1.924453    -1.372676     1.980293   1123      3      5      0      0
   C     2.866758    -1.556753     0.934073   1124      4      7     13     18
   C     2.435730     0.816736     1.151710   1126      3      7      8     14
   C     2.832568    -0.159062     0.042099   1128      5      6     15     16
   O     3.554295     1.211441     1.932365   1249      6     17      0      0
   H    -0.918053     0.000000    -0.280677   1245      1      0      0      0
   H    -0.520597    -0.885828     1.803849   1119      2      0      0      0
   H    -0.520597     0.885828     1.803849   1120      2      0      0      0
   H     1.435560     0.337148     2.998879   1122      3      0      0      0
   H     3.838325    -1.808062     1.359516   1125      5      0      0      0
   H     1.936098     1.681209     0.714498   1127      6      0      0      0
   H     2.031585    -0.217259    -0.694882   1129      7      0      0      0
   H     3.838626     0.075227    -0.305832   1130      7      0      0      0
   H     4.214443     1.727289     1.463640   1250      8      0      0      0
   N     2.474231    -2.760890     0.168322   1132      5     19     27      0
   C     1.538394    -2.869204    -0.826353   1136     18     20     28      0
   N     1.421481    -4.070993    -1.308051   1135     19     21      0      0
   C     2.344666    -4.815233    -0.582836   1134     20     22     27      0
   C     2.704630    -6.167666    -0.619591   1140     21     23     24      0
   N     2.152150    -7.057611    -1.455273   1142     22     29     30      0
   N     3.660941    -6.579606     0.239638   1139     22     25      0      0
   C     4.205243    -5.691308     1.066416   1138     24     26     31      0
   N     3.949915    -4.402308     1.191662   1137     25     27      0      0
   C     2.991769    -4.014545     0.323275   1133     18     21     26      0
   H     0.951862    -2.033257    -1.177884   1145     19      0      0      0
   H     2.449361    -8.012246    -1.436882   1143     23      0      0      0
   H     1.442640    -6.767115    -2.097307   1144     23      0      0      0
   H     4.963977    -6.079842     1.729564   1141     25      0      0      0
$end

$rem
   METHOD             cis
   BASIS              6-31G*
   QM_MM_INTERFACE    janus
   USER_CONNECT       true
   FORCE_FIELD        amber99
   GAUSSIAN_BLUR      true
   GAUSS_BLUR_WIDTH   15000
   CIS_N_ROOTS        5
   CIS_TRIPLETS       false
   MOLDEN_FORMAT      true
   PRINT_ORBITALS     true
$end

$qm_atoms
   18 19 20 21 22 23 24 25 26 27 28 29 30 31
$end

View output

Features of this job:

  • An MM-only calculation. BASIS and EXCHANGE need to be defined, in order to prevent a crash, but no electronic structure calculation is actually performed.

  • All atom types and MM interactions are defined in $force_field_params using the CHARMM27 force field. Atomic charges, equilibrium bond distances, and equilibrium angles have been extracted from a HF/6-31G* calculation, but the force constants and van der Waals parameters are fictitious values invented for this example.

  • Molecular dynamics is propagated for 10 steps within a microcanonical ensemble (NVE), which is the only ensemble available at present. Initial velocities are sampled from a Boltzmann distribution at 400 K.

Example 11.24  MM molecular dynamics with user-defined MM parameters.

$molecule
  -2 1
   C    0.803090    0.000000    0.000000  -1  2  3  6  0
   C   -0.803090    0.000000    0.000000  -1  1  4  5  0
   H    1.386121    0.930755    0.000000  -2  1  0  0  0
   H   -1.386121   -0.930755    0.000000  -2  2  0  0  0
   H   -1.386121    0.930755    0.000000  -2  2  0  0  0
   H    1.386121   -0.930755    0.000000  -2  1  0  0  0
$end

$rem
   METHOD            hf
   BASIS             sto-3g
   QM_MM_INTERFACE   MM
   FORCE_FIELD       charmm27
   USER_CONNECT      true
   JOBTYPE           aimd
   TIME_STEP         42
   AIMD_STEPS        10
   AIMD_INIT_VELOC   thermal
   AIMD_TEMP         400
$end

$force_field_params
   NumAtomTypes 2
   AtomType -1  -0.687157  2.0000  0.1100
   AtomType -2  -0.156422  1.3200  0.0220
   Bond     -1  -1  250.00  1.606180
   Bond     -1  -2  300.00  1.098286
   Angle    -2  -1  -2  50.00  115.870
   Angle    -2  -1  -1  80.00  122.065
   Torsion  -2  -1  -1  -2  2.500  180.0  2
$end

View output

Further examples of QM/MM calculations can be found in the $QC/samples directory, including a QM/MM/PCM example, QMMMPCM_crambin.in. This calculation consists of a protein molecule (crambin) described using a force field, but with one tyrosine side chain described using electronic structure theory. The entire QM/MM system is placed within an implicit solvent model, of the sort described in Section 11.2.2.