Energy minimization and Stereochemistry

GROMACS version: 2019.1
GROMACS modification: No
Here post your question

I used opls-forcefield to simulate a peptide 30 times (new coordinates from previous run and random velocities). Surprisingly, when these 30 structures were energy minimized, they did not have similar potential energy. For example, say structure #4 was 40 KJ/mol higher in energy than #5. Two main differences were: More H-bonds and less deviation from peptide bond stereochemistry in the lower energy structure.
Does that imply I should use any external program for proper energy minimization, specifically proper stereochemistry?

Observation: Energy minimized structures from gromacs have significant deviations from allowed stereochemistry (dihedral angles around peptide bond). When I use PyMol ‘clean’ function to energy minimize structures obtained from gromacs, all strains are removed and shows expected stereochemistry.

Question: Does that imply I should use an external program for proper energy minimization, specifically proper stereochemistry?

When you do energy minimization, energy of entire system is minimized. Energy is consist of both bonded (bond, angle, dihedral) and nonbonded(LJ and electrostatic) interactions. Due to this manybody effects, you can not expect all of bonds, angles and dihedrals are in minima in their respective potential energy surfaces. Please note, dihedral has swallow potential energy surface compared to bond and angle.

I don’t see any reason to use external program for energy minimization.

When, you will run actual molecular dynamics ( energy minimization is not molecular dynamics), molecule will evolve according to underlying force field. Looking at a single frame for bond, angle and dihedral is less conclusive. Only thing that matters, is their distributions.

Pymol uses universal force field (MMFF94) for cleaning small molecules. Unless you are using same Force field for gromacs, you can not expect same structure. If you really think structure is distorted, that is due to either efficacy of force field or mistakes while setting topology file.

Understood thoroughly! Thank you for sharing your knowledge @Masrul