Hydrogen Capping:

An alternative way would be to use the CAPPING flag in order to introduce additional (dummy) hydrogen atoms to saturate the dangling bonds. These capping hydrogen atoms have to be hidden from the MM hamiltonian so the Gromos INPUT and TOPOLOGY files have to be modified for subsequent runs, in addition to adding an explicit EXCLUSION LIST to the cpmd input. The whole procedure is a bit complicated so here is a short protocol of the required steps.

1a

Set up a normal QM/MM run as for using link atoms with the additional keyword CAPPING and instead of the link-atom potential use a hydrogen potential with the additional flag ADD_H. Note that you have to provide the correct number of hydrogens, but no atom index number.

1b

Run a short MD (a couple of steps) and use the resulting CRD_FIN.grm file (under a different name) in the COORDINATES section.

2a

Modify the Gromos input file to match the new coordinate file.

• Increase the number of atoms per solute molecule in the SUBMOLECULES section.
• increase the total number of atoms in the FORCE section.

2b

Modify the Gromos topology file to match the new coordinate file.

• Add a DUM atom type at the end of the ATOMTYPENAME section if not already present (and increase NRATT accordingly).
• if you have added a new atom type, you have to add the corresponding entries in the LJPARAMETERS section as well. Since the capping hydrogen atoms should be invisible from the MM hamiltonian all Lennard-Jones parameters are set to 0.0 for those new entries. This section is a triangular matrix, so you have to add NRATT lines (and increase NRATT2 accordingly).
• Add new residues named DUM (one for each capping hydrogen) to the RESNAME section (and increase NRAA2).
• In the SOLUTEATOM section you have to increase NRP and add the dummy hydrogens at the end of the solute. The structure of the entry is: <atom nr> <residue nr> <atom type name> <vdw type index> <mass> <charge> 1 0 and a single '0' on the next line. Use a mass of 1.008 and a charge of 0.000.

2c

Modify the CPMD input file.

• Make sure that the TOPOLOGY, INPUT, and COORDINATES keywords in the &QMMM section match the newly created or modified files.
• Add the capping hydrogens to the &ATOMS section as normal QM atoms, but add the DUMMY flag to the pseudopotential line.
• Build an EXCLUSION LIST entry that lists for each capping atom the respective QM/MM atoms pairs that should be excluded from the electrostatic coupling (all other MM interactions are set to zero already in the topology file). For consistency only full charge groups should be excluded. In the supplementary material should be a script genexcl.tcl which can help you in building that list (it needs the modified Gromos coordinate and topology file as well as the QMMM_ORDER file as input).
• Update the ARRAYSIZES ... END ARRAYSIZES entry to match the new topology.

With the three modified files you should be able to run a regular QM/MM run. Note, that you may have to update the exclusion list occasionally, depending on your system and that you should pick the bond(s) to cut very carefully.