Re: [AMBER] SHAKE failure

From: Pengfei Li <ambermailpengfei.gmail.com>
Date: Sat, 23 Sep 2017 16:34:31 -0500

Hi Lizelle,

I guess you met the problem mentioned in the end of the webpage: http://ambermd.org/tutorials/advanced/tutorial20/mcpbpy.htm <http://ambermd.org/tutorials/advanced/tutorial20/mcpbpy.htm>. To keep it short, this is because that the hydrogen atoms in water has zero VDW interaction with surroundings. However, since they have positive charges, when they have a 1-4 interaction any atoms with negative charges, the structural collapse would happen. This is because there is only 1-4 electrostatic attraction but no 1-4 VDW repulsion under this situation, this is an error can be meet in the bonded model of metal ion. To solve the problem, you can consider to use the restrained nonbonded model: http://ambermd.org/tutorials/advanced/tutorial20/rnb.htm <http://ambermd.org/tutorials/advanced/tutorial20/rnb.htm>.

Kind regards,
Pengfei

> On Sep 22, 2017, at 1:40 PM, Lizelle Lubbe <LBBLIZ002.myuct.ac.za> wrote:
>
> Sorry just a slight correction... During minimization 2 restraintmask :WAR & !.H=
>
> Lizelle
> ________________________________________
> From: Lizelle Lubbe <LBBLIZ002.myuct.ac.za>
> Sent: 22 September 2017 08:37:05 PM
> To: amber.ambermd.org; Daniel Roe; ross.cgl.ucsf.edu
> Subject: Re: [AMBER] SHAKE failure
>
> Dear Daniel and Bill,
>
> Thanks so much for your replies. Daniel, I have repeated the equilibration with printing at each step and saw that just as in the minimizations when having the water unrestrained, the oxygen moves sideways and causes an overlap of the H with OE1 of the Zn-bound glutamate. The initial pdb input for leap and leap output look fine and cpptraj checks out with no overlap messages. The overlap only occurs when allowing the water to move freely. I suspect that there is something wrong in my water parameters somewhere but don't know where to start searching. I've attached some images of the equilibration printing every step and some energy output below.
> Bill, the input files of all steps are given below. Also the lib file of the water molecule. I obtained this from another group and perhaps don't understand the content very well. I replaced the coordinates under !entry.'resname'.unit.positions table with my Zn-bound water coordinates from the pdb. Should I rather have retained the original values matching the ones for WAT in solvents.lib? I also replaced the Zn-bound histidine and glutamate coordinates with my own and am wondering if this was a mistake and results in this strange behaviour?
>
> original water lib file from publisher's of Zn parameters:
>
> !!index array str
> "WAR"
> !entry.WAR.unit.atoms table str name str type int typex int resx int flags int seq int elmnt dbl chg
> "O" "OK" 0 1 131072 3 8 -0.751392
> "H1" "HW" 0 1 131072 1 1 0.396784
> "H2" "HW" 0 1 131072 2 1 0.424564
> !entry.WAR.unit.atomspertinfo table str pname str ptype int ptypex int pelmnt dbl pchg
> "O" "OK" 0 -1 0.0
> "H1" "HW" 0 -1 0.0
> "H2" "HW" 0 -1 0.0
> !entry.WAR.unit.boundbox array dbl
> -1.000000
> 0.0
> 0.0
> 0.0
> 0.0
> !entry.WAR.unit.childsequence single int
> 2
> !entry.WAR.unit.connect array int
> 0
> 0
> !entry.WAR.unit.connectivity table int atom1x int atom2x int flags
> 3 1 1
> 2 1 1
> 2 3 1
> !entry.WAR.unit.hierarchy table str abovetype int abovex str belowtype int belowx
> "U" 0 "R" 1
> "R" 1 "A" 1
> "R" 1 "A" 3
> "R" 1 "A" 2
> !entry.WAR.unit.name single str
> "TP3"
> !entry.WAR.unit.positions table dbl x dbl y dbl z
> 0.0 0.0 0.0
> 0.957200 0.0 0.0
> -0.239988 0.926627 0.0
> !entry.WAR.unit.residueconnect table int c1x int c2x int c3x int c4x int c5x int c6x
> 0 0 0 0 0 0
> !entry.WAR.unit.residues table str name int seq int childseq int startatomx str restype int imagingx
> "WAR" 1 4 1 "w" 1
> !entry.WAR.unit.residuesPdbSequenceNumber array int
> 0
> !entry.WAR.unit.solventcap array dbl
> -1.000000
> 0.0
> 0.0
> 0.0
> 0.0
> !entry.WAR.unit.velocities table dbl x dbl y dbl z
> 0.0 0.0 0.0
> 0.0 0.0 0.0
> 0.0 0.0 0.0
>
> minimization 1 of solvent:
>
> # Control section
> &cntrl
> imin=1,
> irest = 0,
> ntmin = 1, maxcyc = 5000, ncyc = 2500, dx0 = 0.01, drms = 0.0001,
> cut = 10.0,
> ntb = 1,
> ntp = 0,
> ntwx = 50,
> ntpr = 50,
> ntwr = 5000,
> ioutfm = 1,
> iwrap = 1,
> ntr = 1, restraintmask='!:WAT,Na+', restraint_wt=50.0
>
> /
>
> minimization 2 of system (WAR = Zn-bound water):
>
> # Control section
> &cntrl
> imin=1,
> irest = 0,
> ntmin = 1, maxcyc = 20000, ncyc = 10000, dx0 = 0.01, drms = 0.0001,
> cut = 10.0,
> ntb = 1,
> ntp = 0,
> ntwx = 100,
> ntpr = 100,
> ntwr = 20000,
> ioutfm = 1,
> iwrap = 1,
> ntr = 1, restraintmask=':WAR', restraint_wt=50.0
>
> /
>
> heating solvent:
>
> # Control section
> &cntrl
> imin = 0,
> irest = 0,
> ntx = 1,
> nstlim = 300000, dt = 0.002,
> ntc = 2, ntf = 2,
> tempi = 10.0, temp0 = 300.0,
> ntwx = 5000, ntpr = 5000, ntwe = 5000,
> ntb = 1, ntp = 0,
> ntt = 3, gamma_ln = 2, ig=-1,
> cut = 10.0,
> ntxo = 2,
> ioutfm = 1,
> iwrap = 1,
> nmropt = 1,
> ntr = 1, restraintmask='!:WAT,Na+', restraint_wt=50.0
> /
> &wt type='TEMP0', istep1=0, istep2=100000, value1=10.0, value2=100.0/
> &wt type='TEMP0', istep1=100001, istep2=200000, value1=100.0, value2=200.0/
> &wt type='TEMP0', istep1=200001, istep2=300000, value1=200.0, value2=300.0/
> &wt type='END'/
>
> equilibration 1 of solvent:
>
> # Control section
> &cntrl
> imin = 0,
> irest = 1,
> ntx = 5,
> ntb = 2, ntp = 1, pres0 = 1.0, taup = 2.0,
> cut = 10.0,
> ntc = 2, ntf = 2,
> temp0 = 300.0,
> ntt = 3, gamma_ln = 2, ig = -1,
> nstlim = 100000, dt = 0.002,
> ntwx = 1000, ntpr = 1000, ntwe = 1000,
> ntwr = 100000,
> ioutfm=1,
> iwrap=1,
> ntr = 1, restraintmask='!:WAT,Na+', restraint_wt=10.0,
> /
>
> equilibration 2 of system:
>
> # Control section
> &cntrl
> imin = 0,
> irest = 1,
> ntx = 5,
> ntb = 2, ntp = 1, pres0 = 1.0, taup = 2.0,
> cut = 10.0,
> ntc = 2, ntf = 2,
> temp0 = 300.0,
> ntt = 3, gamma_ln = 2, ig = -1,
> nstlim = 100, dt = 0.002,
> ntwx = 1, ntpr = 1, ntwe = 1,
> ntwr = 100,
> ioutfm=1,
> iwrap=1,
> ntr = 0,
> /
>
> Selected output of eq2:
>
> NSTEP = 1 TIME(PS) = 800.002 TEMP(K) = 299.47 PRESS = 111.4
> Etot = -457615.5261 EKtot = 115240.5321 EPtot = -572856.0582
> BOND = 2615.1688 ANGLE = 6117.2451 DIHED = 7678.7143
> 1-4 NB = 3205.5098 1-4 EEL = 46188.6133 VDWAALS = 79065.9051
> EELEC = -717727.2146 EHBOND = 0.0000 RESTRAINT = 0.0000
> EKCMT = 53294.0126 VIRIAL = 48629.1977 VOLUME = 1939919.9138
> Density = 0.9935
> Ewald error estimate: 0.4355E-04
> ------------------------------------------------------------------------------
> NSTEP = 20 TIME(PS) = 800.040 TEMP(K) = 297.96 PRESS = 23.7
> Etot = -457190.3038 EKtot = 114659.3744 EPtot = -571849.6782
> BOND = 2699.5720 ANGLE = 6578.6660 DIHED = 7858.9118
> 1-4 NB = 3248.6770 1-4 EEL = 46041.9874 VDWAALS = 79365.4920
> EELEC = -717642.9844 EHBOND = 0.0000 RESTRAINT = 0.0000
> EKCMT = 53460.8164 VIRIAL = 52469.7550 VOLUME = 1939949.3169
> Density = 0.9935
> Ewald error estimate: 0.5955E-04
> ------------------------------------------------------------------------------
> NSTEP = 40 TIME(PS) = 800.080 TEMP(K) = 299.86 PRESS = -9.9
> Etot = -455438.8761 EKtot = 115389.4418 EPtot = -570828.3179
> BOND = 3420.4131 ANGLE = 6819.4033 DIHED = 7909.4072
> 1-4 NB = 3488.2059 1-4 EEL = 45783.8401 VDWAALS = 78987.9064
> EELEC = -717237.4939 EHBOND = 0.0000 RESTRAINT = 0.0000
> EKCMT = 53403.0882 VIRIAL = 53819.6086 VOLUME = 1939986.1076
> Density = 0.9935
> Ewald error estimate: 0.4069E-04
> ------------------------------------------------------------------------------
>
> vlimit exceeded for step 40; vmax = 22.4968
>
> Coordinate resetting (SHAKE) cannot be accomplished,
> deviation is too large
> NITER, NIT, LL, I and J are : 0 1 6057 12494 12496
>
> Thanks for the assistance
>
> Lizelle Lubbe
>
> PhD (Medical biochemistry) candidate
> Department of Integrative Biomedical Sciences
> University of Cape Town
> ________________________________________
> From: Bill Ross <ross.cgl.ucsf.edu>
> Sent: 22 September 2017 04:18 PM
> To: amber.ambermd.org
> Subject: Re: [AMBER] SHAKE failure
>
> Please copy/paste your min.in and md.in files for minimization and heating.
>
> Bill
>
>
> On 9/22/17 7:08 AM, Daniel Roe wrote:
>> On Fri, Sep 22, 2017 at 8:30 AM, Lizelle Lubbe <LBBLIZ002.myuct.ac.za> wrote:
>>> I've prepared my glycosylated metalloprotein using tleap and attempted minimization, heating and equilibration. Upon equilibration (having ntc = ntf = 2) I almost instantly get the following error:
>>>
>>> vlimit exceeded for step 61; vmax = 27.0336
>>> vlimit exceeded for step 62; vmax = 536.4434
>>>
>>> Coordinate resetting (SHAKE) cannot be accomplished,
>>> deviation is too large
>>> NITER, NIT, LL, I and J are : 0 5 6058 12494 12495
>>>
>>> Note: This is usually a symptom of some deeper
>>> problem with the energetics of the system.
>>>
>>> Atom 12494 is the O and atom 12495 the H1 of a Zn-coordinated water molecule (described using published forcefield parameters). Since the O and H1 are explicitly bonded I do not understand why SHAKE would want to constrain these bonds. Am I reading the error correctly? In parmed.py printBonds, printAngles, printDetails and printDihedrals for this water molecule all seem fine.
>> When you specify ntc = 2, *all* bonds to hydrogen are constrained. The
>> vlimit errors indicate your system is blowing up. It may be worthwhile
>> to re-run the MD and print out energies and trajectory frames every
>> step so you can see what is happening in more detail (especially
>> because it happens so quickly the file sizes should remain
>> reasonable). What may be happening is that the H of the water ends up
>> having a bad overlap with the Zn and shoots off into oblivion or
>> something. One thing you could try to go is check for bad overlaps in
>> your initial structure using the 'check' command from cpptraj:
>> https://github.com/Amber-MD/cpptraj. You may need more minimization
>> before running your MD phase.
>>
>> -Dan
>>
>>> Could anyone please help me to solve this issue?
>>>
>>> Lizelle Lubbe
>>>
>>> PhD (Medical biochemistry) candidate
>>> Department of Integrative Biomedical Sciences
>>> University of Cape Town
>>> Disclaimer - University of Cape Town This e-mail is subject to UCT policies and e-mail disclaimer published on our website at http://www.uct.ac.za/about/policies/emaildisclaimer/ or obtainable from +27 21 650 9111. If this e-mail is not related to the business of UCT, it is sent by the sender in an individual capacity. Please report security incidents or abuse via csirt.uct.ac.za
>>>
>>> _______________________________________________
>>> AMBER mailing list
>>> AMBER.ambermd.org
>>> http://lists.ambermd.org/mailman/listinfo/amber
>>
>>
>
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
> Disclaimer - University of Cape Town This e-mail is subject to UCT policies and e-mail disclaimer published on our website at http://www.uct.ac.za/about/policies/emaildisclaimer/ or obtainable from +27 21 650 9111. If this e-mail is not related to the business of UCT, it is sent by the sender in an individual capacity. Please report security incidents or abuse via csirt.uct.ac.za
> Disclaimer - University of Cape Town This e-mail is subject to UCT policies and e-mail disclaimer published on our website at http://www.uct.ac.za/about/policies/emaildisclaimer/ or obtainable from +27 21 650 9111. If this e-mail is not related to the business of UCT, it is sent by the sender in an individual capacity. Please report security incidents or abuse via csirt.uct.ac.za
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber

_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Sat Sep 23 2017 - 15:00:03 PDT
Custom Search