Covalent bond breaks after energy minimization

GROMACS version:VERSION 5.1.1
GROMACS modification: Yes/No
Here post your question

I’m a beginner to md simulations and just started learning. I’ve desinged a silver nanoparticle and packed butanethiol as a capping agent using Packmol. When I try to do energy minimization, covalent bonds of the resulting structure are broken( visualized using vmd) . I also solvated my nanoparticle - capping agent system in cyclohexane and water (separately) and when energy minimization is done again I see the covalent bonds are broken. I cannot figure out where it has gone wrong? Does anyone have an idea about what I should do?

This is the minim.mdp file I used.

; minim.mdp - used as input into grompp to generate em.tpr
; Parameters describing what to do, when to stop and what to save
integrator = steep ; Algorithm (steep = steepest descent minimization)
emtol = 1000.0 ; Stop minimization when the maximum force < 1000.0 kJ/mol/nm
emstep = 0.01 ; Minimization step size
nsteps = 50000 ; Maximum number of (minimization) steps to perform

; Parameters describing how to find the neighbors of each atom and how to calculate the interactions
nstlist = 1 ; Frequency to update the neighbor list and long range forces
cutoff-scheme = Verlet ; Buffered neighbor searching
ns_type = grid ; Method to determine neighbor list (simple, grid)
coulombtype = PME ; Treatment of long range electrostatic interactions
rcoulomb = 1.2 ; Short-range electrostatic cut-off
rvdw = 1.2 ; Short-range Van der Waals cut-off
pbc = xyz ; Periodic Boundary Conditions in all 3 dimensions

Thank you.

Covalent bonds cannot break or form in a molecular mechanics calculation. What you are seeing is distortion of your structure that your visualization software can’t handle. It means the structured is strained and you should not consider it stable.

Thank you sir!

I have another question. I want make a box using cyclohexane as the solvent. The steps I followed was,

  1. obtaining .pdb file for cyclohexane and use pdb2gmx to get .gro file
  2. using the command ‘gmx editconf -f chx.gro -o chx_newbox.gro -c -d 1.0 -bt cubic’ define a box
  3. then use ’ gmx insert-molecules -ci chx.gro -nmol 200 -box 3 3 3 -o chx_box.gro’ to make a box
  4. Use ’ gmx solvate -cp NP10_newbox.gro -cs chx_box.gro -o NP10_chx.gro ’ to place a nanoparticle - with a capping agent inside the cyclohexane box.

When visualized I can see the Nanoparticle is properly placed inside box. But when I try to energy minimize it, I get an error saying, " number of coordinates in coordinate file (NP10_chx.gro, 6979) does not match topology (, 319)"

In website what they have suggested is,

" Choose a box density and box size. The size does not have to be that of your eventual simulation box - a 1nm cube is probably fine. Generate a single molecule of the solvent. Work out how much volume a single molecule would have in the box of your chosen density and size. Use [editconf] to place a box of that size around your single molecule. Then use [editconf] to move the molecule a little bit off center. Then use genconf -rot to replicate that box into a large one of the right size and density. Then equilibrate thoroughly to remove the residual ordering of the molecules, using NVT and periodic boundary conditions. Now you have a box you can pass to genbox -cs, which will replicate it to fit the size of the actual simulation box."

But it’s not clear to me. Do you think the steps I followed are not correct? If what is not your suggestion?

I truly appreciate your help. Thank you very much!

You should equilibrate the cyclohexane box before trying to use it.

Thank you sir! It was of great help.