As a note to Alex (and the rest of the list), the coarse-grained Martini forcefield is usually run with timesteps between 20-40 fs. 15fs is already rather low. I do agree that longer equilibration at low timestep (5 or 10) might help.
Alternatively, Do you think a semiisotropic pressure coupling might be applicable in this case, since it's an infinite collagen polymer? Peter (PhD in the Martini group) On 16-12-16 00:21, Nash, Anthony wrote: > Alex and Mark, thanks for the information. I’ll drop dt down, > significantly, drop the temperature, and run it for a long time. > > Thanks for the ideas. > Anthony > > > On 15/12/2016 21:54, "gromacs.org_gmx-users-boun...@maillist.sys.kth.se on > behalf of Alex" <gromacs.org_gmx-users-boun...@maillist.sys.kth.se on > behalf of nedoma...@gmail.com> wrote: > >> Mark is right, no two ways about it. For initial equilibration and >> assessing preexisting structural strains try vacuum, _much_ smaller >> timesteps and possibly low temperatures in vacuum, only then transfer to >> solvent, etc. Algorithmically, LINCS requires convergence and you already >> are using a pretty high LINCS order... From what I see, dt = 15 fs at 310K >> looks like a cowboy mode simulation in this case. >> >> Alex >> >> On Thu, Dec 15, 2016 at 2:32 PM, Mark Abraham <mark.j.abra...@gmail.com> >> wrote: >> >>> Hi, >>> >>> If a simulation isn't stable with a small time step (as I think you are >>> saying) then moving to a larger time step is guaranteed to make that >>> worse. >>> Try an even smaller time step, for a long time, and see what happens. Or >>> take a subset of your protein and see what happens. Or simulate in vacuo >>> for a while. Your topology could be unsuited to your starting structure, >>> e.g. some part is under a lot of tension that gets released at some >>> point >>> and no finite time step can in practice deal with the velocity of the >>> recoil... >>> >>> Mark >>> >>> On Thu, 15 Dec 2016 23:06 Nash, Anthony <a.n...@ucl.ac.uk> wrote: >>> >>>> Hi all, >>>> >>>> I¹m hoping for some help. I¹m very sorry, this is a bit of a long one. >>>> >>>> I¹ve been struggling for almost a month trying to run a CG >>> representation >>>> of our all-atom model of a collagen protein (3 polypeptide chains in a >>>> protein). Our original AMBER all-atom model had been successful >>> modelling >>>> using MD. I went on to use the latest version of Martinize.py with the >>>> latest version of the MARTINI forcefield fields. >>>> >>>> After a little tweaking (the way AMBER names histidine residues), I >>>> successful converted the molecule (approx 3100 amino acids) into a CG >>>> representation. I successfully energy min the protein in vacuum to a >>>> threshold of 500, and in solvent to a threshold of 750 using steepest >>>> descent. Looking for a system at an energy min of a threshold around >>> 300 >>> I >>>> begin to see LINCS warnings. Observing the initial structure, there is >>>> nothing obviously wrong with the bond network (both protein and >>> polarized >>>> CG water). >>>> >>>> I take the system that energy mins at 750 (protein-water mix, with no >>>> fault reported), and went straight to NPT, 20fs step. Blew up. After a >>> bit >>>> of chatting with the MARTINI community, I¹ve started with an NVT >>> ensemble, >>>> beginning at 5s then through 10fs, 15fs, and 20fs. I only run for 1000 >>>> steps before switching. Keeping any of the simulations running for >>> longer >>>> throws lincs warnings followed by a segmentation fault from the >>> warning: >>>> "3 particles communicated to PME rank 7 are more than 2/3 times the >>>> cut-off out of the domain decomposition cell of their charge group in >>>> dimension x." >>>> >>>> Observing the trajectories of any of the extended simulations shows >>> the >>>> protein snapping like a rope, and always at the same place. I have >>> watched >>>> every trajectory at this point, using numerous energy min start >>> points, >>> to >>>> try and understand why it is blowing up. I can¹t see any obvious >>> reason. >>> I >>>> was told to consider how the temperature is changing. Below is an >>> example >>>> of the temperature and pressure from an NPT of 20fs step continued >>> from >>>> the very short 20fs step NVT simulation (hoping that perhaps CG >>> without >>>> pressure just doesn¹t behave happily; I was wrong). >>>> >>>> >>>> TEMP: >>>> Š >>>> 6.630000 311.000336 >>>> 6.645000 311.371643 >>>> 6.660000 311.724213 >>>> 6.675000 313.878693 >>>> 6.690000 681558.937500 >>>> >>>> >>>> PRESSURE: >>>> Š >>>> 6.630000 3.559879 >>>> 6.645000 3.901433 >>>> 6.660000 3.589078 >>>> 6.675000 4.158611 >>>> 6.690000 81762.437500 >>>> >>>> The final LINCS warning from this same run: >>>> >>>> Step 300, time 4.5 (ps) LINCS WARNING >>>> relative constraint deviation after LINCS: >>>> rms 0.000035, max 0.003386 (between atoms 2125 and 2126) >>>> bonds that rotated more than 45 degrees: >>>> atom 1 atom 2 angle previous, current, constraint length >>>> 2125 2126 68.3 0.2781 0.2691 0.2700 >>>> 2125 2127 45.9 0.2789 0.2701 0.2700 >>>> >>>> >>>> At this stage the structure ruptures as described above. >>>> >>>> >>>> My NVT settings (with NPT included to save space) are: >>>> >>>> ----------------- >>>> title = Martini >>>> >>>> integrator = md >>>> dt = 0.015 >>>> nsteps = 1000 >>>> nstcomm = 100 >>>> ;comm-grps = >>>> >>>> nstxout = 1000 >>>> nstvout = 1000 >>>> nstfout = 0 >>>> nstlog = 1 >>>> nstenergy = 1 >>>> nstxout-compressed = 0 >>>> compressed-x-precision = 0 >>>> ;compressed-x-grps = >>>> energygrps = collagen solvent >>>> >>>> cutoff-scheme = Verlet >>>> nstlist = 20 >>>> ns_type = grid >>>> >>>> pbc = xyz >>>> verlet-buffer-tolerance = 0.005 >>>> >>>> coulombtype = PME ;reaction-field >>>> rcoulomb = 1.1 >>>> fourierspacing = 0.16 ;0.2 ;0.12 >>>> >>>> epsilon_r = 2.5 ;15 ; 2.5 (with polarizable water) >>>> epsilon_rf = 0 >>>> vdw_type = cutoff >>>> vdw-modifier = Potential-shift-verlet >>>> rvdw = 1.1 >>>> >>>> tcoupl = v-rescale ;berendsen ;v-rescale >>>> tc-grps = collagen solvent >>>> tau_t = 0.5 0.5 ;1.0 1.0 >>>> ref_t = 310 310 >>>> >>>> Pcoupl = berendsen ;parrinello-rahman >>>> Pcoupltype = isotropic >>>> tau_p = 12.0 ; parrinello-rahman is more stable >>> with >>>> larger tau-p, DdJ, 20130422 >>>> compressibility = 10e-4 >>>> ref_p = 1.0 >>>> refcoord_scaling = com >>>> >>>> gen_vel = no >>>> gen_temp = 310 >>>> gen_seed = 473529 >>>> >>>> continuation = yes >>>> constraints = none >>>> constraint_algorithm = lincs >>>> lincs-warnangle = 45 >>>> lincs-order=8 >>>> lincs-iter=4 >>>> >>>> >>>> ‹‹‹‹‹‹‹‹‹‹ >>>> >>>> Every setting bar the lincs iter, order, warnangle were supplied with >>> the >>>> latest version of MARTINI. During many NVT runs I have adjusted the >>> tau-t >>>> to try and keep the thermostat from oscillating its way into infinity. >>>> >>>> I¹m curious, will an out of control thermostat break a structure, or >>> will >>>> a structure breaking (for what ever reason this structure is breaking) >>>> cause the thermostat to go out of control? >>>> >>>> My only thought thought is the initial .itp file that Martinize >>> created. >>> I >>>> informed the script that this was collagen, therefor it sets ³F² and >>> the >>>> corresponding bonded parameters. A human collagen is not perfect in >>> its >>>> helical structure. Could there be underlying forces contributed from >>> badly >>>> bonded backbone-backbone arrangements? >>>> >>>> [ atoms ] >>>> 1 N0 1 GLY BB 1 0.0000 ; F >>>> 2 N0 2 LEU BB 2 0.0000 ; F >>>> 3 C1 2 LEU SC1 3 0.0000 ; F >>>> 4 N0 3 SER BB 4 0.0000 ; F >>>> >>>> >>>> Many thanks >>>> Anthony >>>> >>>> Thanks >>>> Anthony >>>> >>>> -- >>>> 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-requ...@gromacs.org. >>>> >>> -- >>> 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-requ...@gromacs.org. >> -- >> 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-requ...@gromacs.org.
-- 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-requ...@gromacs.org.