Position restrains - why is it heavy atoms only?

GROMACS version: 2020.4
GROMACS modification: Yes/No
Hi all

I’m trying to define position restrains prior to equilibration. I’m curious if there’s a reason why the default is to position restrain heavy atoms/omitting hydrogens?

Sam

Bonds to H atoms are almost always constrained. The combination of constraints and restraints may in some cases oppose each other and lead to instability. Another point is hydrogen bonding. You want to allow H atoms to reorient to engage in favorable H-bonds to solvent or other entities that may be stabilizing but perhaps we’re absent from the initial coordinates.

On an even higher level, the purpose of restraints is to prevent distortion of the structure after randomization of velocities in conjunction with unoptimized solvent. H atoms don’t dictate the overall folded structure of the molecule and it is sufficient to restrain the heavier atoms, which do more directly govern the structure.

Yes, this all makes perfect sense.

One last thing, the numbers/order of numbers present in a position restraint file should match the numbers/order of the atoms you wish to restrain in the different [ moleculetype ] in the .top file? I have a cofactor position restraint file whereby I initially restrained all atoms and then I removed the numbers of atoms in the restraint file that the .top file told me were hydrogens. Hopefully this would give me just cofactor heavy atom position restraints but was skeptical if this is a valid method.

Thank you in advance

It doesn’t really matter how you create the file as long as you specify the atoms to which you want to apply restraints.