[gmx-users] problems with pull code on Martini CG

[Berk Hess]

Hi,

You should not use direction-periodic, you don't want to pull more than half the box size, I presume.

Your pull rate is extremely high 50 nm/ns. I would start by lowering this.

Cheers,

Berk

________________________________
From: gromacs.org_gmx-users-bounces at maillist.sys.kth.se <gromacs.org_gmx-users-bounces at maillist.sys.kth.se> on behalf of Bennett Addison <baddison at sdsu.edu>
Sent: Thursday, April 11, 2019 12:50 AM
To: gmx-users at gromacs.org
Subject: [gmx-users] problems with pull code on Martini CG

Hello gromacs users,

I am fairly new to MD, trying hard to troubleshoot this issue myself before
posting. But I am at a loss.

I am trying to simulate pulling on a large protein, using Martini CG
simulations and adding the pull code. I am following tutorials for umbrella
sampling (all atom) and for martini protein simulations. Basically, the
exact same setup procedure works when the box size is small, but when I
scale up, things get screwy.

For example: 200 AA protein, placed 50 x 12 x 12 box and pull along X (see
pull code below: using direction-periodic geometry, pull along vector -1.0
0.0 0.0)

gmx editconf -f CG.pdb -o newbox.gro -box 50 12 12 -center 7 6 6

Followed by minimization in vac, solvate, minimize, equilibrate, then
create groups with gmx make_ndx, assign Chain_A as anchored residue, assign
Chain_B as residue to pull, make sure position restraints are correct, and
finally setup and run md_pull.mdp. This simulation works great. I see a
nice unfolding using martini CG model.
Pull code:

; Pull code

pull                    = yes

pull_ncoords            = 1         ; only one reaction coordinate

pull_ngroups            = 2         ; two groups defining one reaction
coordinate

pull_group1_name        = Chain_A

pull_group2_name        = Chain_B

pull_coord1_type        = umbrella  ; harmonic potential

pull_coord1_geometry = direction-periodic ;

pull_coord1_vec = -1.0 0.0 0.0

pull_coord1_groups      = 1 2

pull_coord1_start       = yes       ; define initial COM distance > 0

pull_coord1_rate        = 0.05      ; 0.05 nm per ps = 5 nm per ns

pull_coord1_k           = 1000      ; kJ mol^-1 nm^-2



Now when I try to replicate the exact same setup but using a larger box (100nm
20nm 20nm) everything fails. When I visualize the trajectory my protein
does not appear to be pulled, but the position output pullx.xvg shows
increase in distance, pullf.xvg shows no changes in force, and if I plot
distance between Chain_A and Chain_B using gmx distance, I see an expected
distance increase.


DOES NOT WORK: gmx editconf -f CG.pdb -o new box.gro -box 100 20 20 -center
8 10 10 (pullx looks ok, but pullf is flat, protein does not extend)

WORKS: gmx editconf -f CG.pdb -o new box.gro -box 50 12 12 -center 5 6 6
(trajectory looks great, pullf and pullx as expected)


Every other parameter / setup procedure is the same, just box size is
different. I need a larger box on my bigger protein (200 AAs) because it
will fully extend past 50nm.


I have tried with two configurations, same issues

- gromacs version 5.1.5 on a computing cluster, single node, AMD with 32
cores

- gromacs version 2018.3 on local machine with intel xeon E5 processors, 40
cores with GPU acceleration


Does anybody know why this may be happening, and how I can get around it?


Thank you!

-Bennett


________________________
Dr. Bennett Addison
Manager, SDSU NMR Facility
baddison at sdsu.edu
206-235-5415 <//206-235-5415> (cell)
--
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.