[gmx-users] How to resolve Inconsistent moves, Graphene sheet simulation with Gromacs 5.1.4

[Alex]

So, I checked your mdp file again... It is my understanding that you are 
running your md.mdp without prior energy minimization and using 
constraints=all-bonds and pressure coupling. This is a very "lucky" 
combination and internal lattice strain will cause a mini-explosion. I 
am actually surprised you're not getting a bunch of LINCS errors here -- 
which means that your initial coordinates are quite good.

In any case, for graphene and all other crystals, you need EM first, 
then maybe even a semiisotropic equilibration run at a very low 
temperature (something like 5K), only after that run your room-temp 
relaxation MD. In all of these, it should be with constraints=none and 
Berendsen barostat for any initial relaxations.

Good luck!

Alex


On 4/22/2018 4:48 AM, Edwine Tendong wrote:
> Greetings Alex,
>
> Thanks for your prompt reply and suggestions. I used the parameters you
> suggested in my itp file. Here is my itp file
> https://www.dropbox.com/s/mz5dwk7z0maq0hp/forcefield.itp?dl=0
>
> I prepared my graphene sheet in VMD using the C-C bond length = 0.142nm and
> made sure the box size respected periodicity (My sheet is rectangular:
> 2.5560000   2.214000 ). Here is a visualisation of the box and its periodic
> images in vmd
> https://www.dropbox.com/s/x4xizjppln35cle/graphene_sheet.png?dl=0
>
> I then generated the topology with x2top version 5.1.4 with the command
> ''gmx x2top -f gra.gro -o gra.top -ff graphene -name GRA -pbc -noparam''
> There are  324 proper dihedrals,    0 impropers,  648 angles     972
> pairs,   324 bonds and   216 atoms. The connectivity looks good. x2top
> automatically assigns the dihedral type 1 so I manually change to type 3 to
> match the supplied dihedral energies you suggested.
> Then doing the run  with the command
> gmx grompp -f md.mdp -c gra.gro -p gra.top -o md_01.tpr
> gmx mdrun -v -deffnm md_01
>
> Reading file md_01.tpr, VERSION 5.1.4 (single precision)
> Changing nstlist from 10 to 40, rlist from 1 to 1
>
> Using 1 MPI thread
> Using 8 OpenMP threads
>
> starting mdrun 'GRA'
> 100000 steps,    100.0 ps.
> step 99900, remaining wall clock time:     0 s
> Writing final coordinates.
> There were 44 inconsistent shifts. Check your topology
> step 100000, remaining wall clock time:     0 s
>
> I have visualised the trajectory and the atoms just vibrate about
> equilibrium as expected. Also removing dihedral does not help.
> I will be very grateful for any further suggestions as I need this to work
> properly.
>
> Thanks
>
> On Sun, Apr 22, 2018 at 5:09 AM, Alex <nedomacho at gmail.com> wrote:
>
>> Hi,
>>
>> Your mdp file contains all the necessary things for this simulation to
>> proceed normally, i.e. COM motion is removed, periodic molecules are set
>> and full pbc is all there. Visualize the trajectory and see what is going
>> on. "I corrected the dihedral function to 3" makes no sense, sorry. The
>> type of dihedrals must correspond to what it was wherever you found the
>> dihedral energy term.
>>
>> For graphene, the symmetry is incredibly simple, let's try to use my
>> parameters. I happen to think that the graphene situation for Gromacs is
>> hopeless, because so many people seem to get their parameters from very
>> strange places, including out-of-the-box stuff from FFs like OPLS-AA.
>>
>> 1. Remove all dihedrals from your topology.
>> 2. Make sure your box size corresponds exactly to the lattice constant.
>> For graphene, it is a=2.46 A and the bond length is a/sqrt(3). In Gromacs,
>> for the bulk bonds I use (CJ is my bulk graphene carbon): CJ    CJ
>> 1    0.14200   420420.0
>> 3. The angle that works for me is: CJ     CJ     CJ      1 120.000
>> 659.346
>> 4. Test if this works without dihedrals.
>> 5. Try adding dihedrals.
>>
>> Parameters for dihedral:
>> CJ     CJ     CJ     CJ      3     17.30770   0.00000 -17.30770 0.00000
>> 0.00000   0.00000
>> This absolutely works, provided you have set up your non-bonded parameters
>> and nrexcl properly.
>>
>> All parameters from http://iopscience.iop.org/arti
>> cle/10.1088/0957-4484/25/48/485701/meta (shameless self-plug). Not a
>> great paper, but can be useful. The dihedral parameters above need to be
>> consistent with the number of dihedrals. If I recall correctly, the paper
>> lists ONE dihedral per bond, so count your dihedrals. It is not expected to
>> change your results much (unless you're specifically interested in
>> ripples), but still try to double-check what the paper says and how your
>> topology is set up.
>>
>> Alex
>>
>>
>> On 4/21/2018 5:15 PM, Edwine Tendong wrote:
>>
>>> Dear Gromacs users,
>>>    I am new to gromacs and I am trying to simulate a graphene sheet in
>>> vacuum
>>> using gromacs 5.1.4. My sheet contains 216 atoms and I have followed the
>>> tips at
>>>    http://www.gromacs.org/Documentation/How-tos/Carbon_Nanotube
>>> to construct my topology. I made sure my box was periodic​​ and I used the
>>> following X2top command:
>>> gmx x2top -f gra.gro -o gra.top -ff graphene -name GRA -pbc -noparam
>>>
>>> This generates 324 bonds, 648 angles as it should. I corrected the
>>> dihedral
>>> function to 3, added periodic_molecules = yes and pbc = xyz to my mdp
>>> file.
>>>
>>> My  simulation runs fine but mdrun reports "There were 44 inconsistent
>>> shifts. Check your topology"  at the end of the run. I am not using
>>> pressure coupling.
>>>
>>> Here are my topology  mdp and coordinates:
>>> https://www.dropbox.com/s/qgdmgeo3uyrs099/gra.top?dl=0
>>> https://www.dropbox.com/s/782cwo5uf2a7ozg/gra.gro?dl=0
>>> https://www.dropbox.com/s/ymazo1p90kt9zu4/md.mdp?dl=0
>>>
>>> I have been unable to determine the source of this error even after going
>>> through many posts in the mailing list. If someone can take a look at
>>> these
>>> simulation files and point me in the right direction then I will be very
>>> grateful.
>>> Thanks a lot in advance.
>>>
>>>
>>>
>> --
>> Gromacs Users mailing list
>>
>> * Please search the archive at http://www.gromacs.org/Support
>> /Mailing_Lists/GMX-Users_List before posting!
>>
>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>
>> * For (un)subscribe requests visit
>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
>> send a mail to gmx-users-request at gromacs.org.
>
>
>