RE: Constraint failure in RATTLE algorithm

From: Douglas Houston (DouglasR.Houston_at_ed.ac.uk)
Date: Fri Jun 27 2014 - 09:29:42 CDT

Hi Tristan,

Yes, you've got it. The CMPS patch is the CHARMM36 CMAP potential for
the end of the staple - every amino acid residue in the CHARMM36
topology file (top_all36_prot.rtf) seems to have one so I assumed it
was appropriate.

The reason I did it this way was because I had example patches I could
modify - I wouldn't know how to write a patch that linked the the two
copies together at their terminal carbons as you suggest (note the
double bond in the middle).

One thing that I was unsure about was the IC tables - each amino acid
residue has one of these but I didn't create one for my OL1 residue. I
thought using "regenerate angles dihedrals" in psfgen accomplished the
same thing?

Here are the relevant entries in my modified .rtf topology file:

RESI OL1 0.00
GROUP
ATOM N NH1 -0.47 ! |
           |
ATOM HN H 0.31 ! HN-N
           NL-HNL
ATOM CA CT1 0.16 ! HB1 | HH1 HG1 HD1 HE1 HP1 HZL1 HDL1
HGL1 HHL1 | HBL1
                         ! \ | | | | | | | |
  | | | /
GROUP
!HB2-CB-CA--CH--CG--CD--CE--CP--CZL==CEL--CDL--CGL--CHL--CAL--CBL-HBL2
ATOM CB CT3 -0.27 ! / | | | | | | | |
| | | \
ATOM HB1 HA3 0.09 ! HB3 | HH2 HG2 HD2 HE2 HP2 HEL1 HDL2
HGL2 HHL2 | HBL3
ATOM HB2 HA3 0.09 ! O=C
           CL=OL
ATOM HB3 HA3 0.09 ! |
           |
GROUP !
ATOM C C 0.51
ATOM O O -0.51
GROUP
ATOM CH CT2 -0.18 !
ATOM HH1 HA2 0.09 !
ATOM HH2 HA2 0.09 !
GROUP
ATOM CG CT2 -0.18 !
ATOM HG1 HA2 0.09 !
ATOM HG2 HA2 0.09 !
GROUP
ATOM CD CT2 -0.18 !
ATOM HD1 HA2 0.09 !
ATOM HD2 HA2 0.09 !
GROUP
ATOM CE CT2 -0.18 !
ATOM HE1 HA2 0.09 !
ATOM HE2 HA2 0.09 !
GROUP
ATOM CP CT2 -0.18 !
ATOM HP1 HA2 0.09 !
ATOM HP2 HA2 0.09 !
GROUP
ATOM CZL CEL1 -0.15 !
ATOM HZL1 HEL1 0.15 !
GROUP
ATOM CEL CEL1 -0.15 !
ATOM HEL1 HEL1 0.15 !
GROUP
ATOM CDL CT2 -0.18 !
ATOM HDL1 HA2 0.09 !
ATOM HDL2 HA2 0.09 !
GROUP
ATOM CGL CT2 -0.18 !
ATOM HGL1 HA2 0.09 !
ATOM HGL2 HA2 0.09 !
GROUP
ATOM CHL CT2 -0.18 !
ATOM HHL1 HA2 0.09 !
ATOM HHL2 HA2 0.09 !
GROUP
ATOM NL NH1 -0.47 ! |
           |
ATOM HNL H 0.31 ! HN-N
           NL-HNL
ATOM CAL CT1 0.16 ! HB1 | HH1 HG1 HD1 HE1 HP1 HZL1 HDL1
HGL1 HHL1 | HBL1
                         ! \ | | | | | | | |
  | | | /
GROUP
!HB2-CB-CA--CH--CG--CD--CE--CP--CZL==CEL--CDL--CGL--CHL--CAL--CBL-HBL2
ATOM CBL CT3 -0.27 ! / | | | | | | | |
| | | \
ATOM HBL1 HA3 0.09 ! HB3 | HH2 HG2 HD2 HE2 HP2 HEL1 HDL2
HGL2 HHL2 | HBL3
ATOM HBL2 HA3 0.09 ! O=C
           CL=OL
ATOM HBL3 HA3 0.09 ! |
           |
GROUP !
ATOM CL C 0.51
ATOM OL O -0.51
BOND CH CA CG CH CD CG CE CD CP CE
BOND CZL CP CDL CEL CGL CDL CHL CGL
BOND CAL CHL CZL HZL1 CEL HEL1 CDL HDL1 CDL HDL2
BOND CGL HGL1 CGL HGL2 CHL HHL1 CHL HHL2
BOND CH HH1 CH HH2 CG HG1 CG HG2 CD HD1
BOND CD HD2 CE HE1 CE HE2
BOND CP HP1 CP HP2
BOND CB CA N HN N CA
BOND CBL CAL NL HNL NL CAL
BOND C CA C +N CB HB1 CB HB2 CB HB3
BOND CL CAL CBL HBL1 CBL HBL2 CBL HBL3
DOUBLE O C OL CL
DOUBLE CZL CEL
IMPR N -C CA HN C CA +N O
CMAP -C N CA C N CA C +N
DONOR HN N
DONOR HNL NL
ACCEPTOR O C
ACCEPTOR OL CL

PRES CMPS 0.00 ! CMAP potential for end of staple
                        ! 1 refers to previous (N terminal)
                        ! 2 refers to staple
                        ! 3 refers to next (C terminal)
                        ! use in a patch statement
CMAP 1C 2NL 2CAL 2CL 2NL 2CAL 2CL 3N

PRES STCN 0.00 ! linkage for joining staple (atom names end in
L) with backbone - pre-staple
                        ! 1 refers to previous (N terminal)
                        ! 2 refers to next (C terminal)
                        ! use in a patch statement
DELETE IMPR 2N 1C 2CA 2HN !Improper specified by IMPR N -C CA HN
DELETE IC 1C 2CA *2N 2HN!Specified by IC -C CA *N HN
BOND 1C 2NL
IMPR 2NL 1C 2CAL 2HNL 1C 1CA 2NL 1O

PRES STNC 0.00 ! linkage for joining staple (atom names end in
L) with backbone - post-staple
                        ! 1 refers to previous (N terminal)
                        ! 2 refers to next (C terminal)
                        ! use in a patch statement
DELETE IMPR 2N 1CL 2CA 2HN !Improper specified by IMPR N -C CA HN
DELETE IC 1C 2CA *2N 2HN!Specified by IC -C CA *N HN
BOND 1CL 2N
IMPR 2N 1CL 2CA 2HN 1CL 1CAL 2N 1OL

cheers,
Doug

Quoting Tristan Croll <tristan.croll_at_qut.edu.au> on Sat, 21 Jun 2014
23:43:43 +0000:

> Let me see if I understand how you've put this together...
>
> You've got your special "linker" residue OL1 which is a 10-carbon
> chain with a peptide backbone on each end, and a methyl group in
> place of the alpha hydrogen. It's effectively symmetrical, but
> you've got one end defined as the "backbone" and inserted as resid
> 3038, and are patching the other end in between 3045 and 3047. So I
> take it...
>
> "patch DELB" deletes the existing amide bond between 3045 and 3047;
> "patch STCN" makes a new amide bond between 3045 and 3038; patch
> STNC makes the bond between 3038 and 3047... what does CMPS do?
>
> I'll let someone more knowledgeable than me comment on the behaviour
> you're seeing... but may I suggest an easier way to get to the same
> end goal? Why not chop your OL1 residue in half? Define it as simply
> a (C-alpha) methylated amino acid with a linear 5-carbon sidechain,
> and insert a copy at positions 3038 and 3046. Then it's one very
> simple patch to link the two copies together at their terminal
> carbons. The fewer (and simpler) patches you introduce, the less
> likely you are to introduce problems - and those problems you do
> introduce will be easier to troubleshoot.
>
>
> Cheers,
>
> Tristan
>
> -----Original Message-----
> From: owner-namd-l_at_ks.uiuc.edu [mailto:owner-namd-l_at_ks.uiuc.edu] On
> Behalf Of Douglas Houston
> Sent: Friday, 20 June 2014 11:32 PM
> To: Aron Broom
> Cc: Namd Mailing List
> Subject: Re: namd-l: Constraint failure in RATTLE algorithm
>
> Hi Aron,
>
> Thanks, this allows me to reproduce the problem. The following files
> should be enough for anyone to visualise the instability:
>
> https://drive.google.com/file/d/0B3IbJv_UNHlWN3p0bTR6NjFER2M/edit?usp=sharing
> https://drive.google.com/file/d/0B3IbJv_UNHlWYXYwYkZCTmkzc2M/edit?usp=sharing
> https://drive.google.com/file/d/0B3IbJv_UNHlWQmlfTjhoNjZUZWc/edit?usp=sharing
>
> Apart from the problem atoms flying around, I notice the whole system
> is "pulsing" - is this normal or indicative of a problem with e.g.
> pressure?
>
> cheres,
> Doug
>
>
> Quoting Aron Broom <broomsday_at_gmail.com> on Thu, 19 Jun 2014 11:37:59 -0400:
>
>> You'd just need to save that frame's coordinates using VMD, and then use
>> those as an input for a "new" simulation. But of course don't do any
>> minimization.
>>
>>
>> On Thu, Jun 19, 2014 at 10:51 AM, Douglas Houston <DouglasR.Houston_at_ed.ac.uk
>>> wrote:
>>
>>> Hi Aron,
>>>
>>> I had considered starting from just before where the distortion sets in
>>> but wasn't sure how to do that?
>>>
>>> As I understand it the "firsttimestep" keyword doesn't specify the frame
>>> to start on, merely where to start the numbering. Is there another keyword
>>> I need to use to get the simulation to start from a particular frame?
>>>
>>> cheers,
>>> Doug
>>>
>>>
>>>
>>>
>>> Quoting Aron Broom <broomsday_at_gmail.com> on Thu, 19 Jun 2014 10:15:49
>>> -0400:
>>>
>>> You're right Normal, I did also mean DCDFreq.
>>>>
>>>> Douglas,
>>>>
>>>> What happens if you take your "penultimate frame" and restart from that,
>>>> does it fail again very quickly?
>>>>
>>>> I agree that if it's taking any more than a few ps to crash, minimization
>>>> likely isn't the problem. It sounds like there is a very sensitive
>>>> problem
>>>> that only appears under specific circumstances, maybe something with a
>>>> dihedral.
>>>>
>>>> ~Aron
>>>>
>>>>
>>>> On Thu, Jun 19, 2014 at 6:26 AM, Douglas Houston <
>>>> DouglasR.Houston_at_ed.ac.uk>
>>>> wrote:
>>>>
>>>> Hi all,
>>>>>
>>>>> Thanks for all the suggestions. A few more details:
>>>>>
>>>>> This is not particularly reproducible. The error may occur after 5ns or
>>>>> 50ns. I have been running 20ns simulations and usually they finish, but
>>>>> now
>>>>> I want to do a 200ns simulation.
>>>>>
>>>>> When I visualise the trajectory just before the failure I see that the
>>>>> residue the atom belongs to (3045) is heavily distorted in the final
>>>>> frame;
>>>>> it looks OK in the penultimate frame. I will rerun with smaller dcdfreq
>>>>> as
>>>>> Norman suggested. I will also increase the timestep from 0.5 to 2 fs as
>>>>> this will hopefully throw the error sooner (my current simulation failed
>>>>> after 3 days of running).
>>>>>
>>>>> I do 10,000 minimization steps to start. The fact that the system is
>>>>> stable for up to 50ns suggests to me that further initial minimization
>>>>> won't help ... ?
>>>>>
>>>>> I understand that I may well have have a problem with bad parameters, or
>>>>> a
>>>>> bad/incomplete topology. Residue 3045 has a number of patches applied to
>>>>> it.
>>>>>
>>>>> The following lists the psfgen commands I have been using to generate the
>>>>> topology:
>>>>>
>>>>> package require psfgen
>>>>> topology top_all36_prot_ole.rtf
>>>>> topology top_all36_lipid.rtf
>>>>> topology toppar_water_ions.str
>>>>> pdbalias residue HIS HSE
>>>>> pdbalias atom ILE CD1 CD
>>>>> segment A {pdb 5CstapleW2A.pdb
>>>>> first ACE
>>>>> last CT2}
>>>>> patch DELB A:3045 A:3047
>>>>> patch STCN A:3045 A:3038
>>>>> patch STNC A:3038 A:3047
>>>>> patch CMPS A:3045 A:3038 A:3047
>>>>> coordpdb 5CstapleW2A.pdb A
>>>>> regenerate angles dihedrals
>>>>> guesscoord
>>>>> writepdb 5CstapleW2A_psfgen.pdb
>>>>> writepsf 5CstapleW2A_psfgen.psf
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> Quoting Axel Kohlmeyer <akohlmey_at_gmail.com> on Wed, 18 Jun 2014 14:23:05
>>>>> -0400:
>>>>>
>>>>> On Wed, Jun 18, 2014 at 1:29 PM, Douglas Houston
>>>>>
>>>>>> <DouglasR.Houston_at_ed.ac.uk> wrote:
>>>>>>
>>>>>> Hi all,
>>>>>>>
>>>>>>> I keep encountering the following fatal error:
>>>>>>>
>>>>>>> ERROR: Constraint failure in RATTLE algorithm for atom 189!
>>>>>>> ERROR: Constraint failure; simulation has become unstable.
>>>>>>> ERROR: Exiting prematurely; see error messages above.
>>>>>>>
>>>>>>> The atom itself varies. I have searched previous messages and tried
>>>>>>> the
>>>>>>> suggestions (smaller timestep, minimization overkill, etc.) to no
>>>>>>> avail.
>>>>>>>
>>>>>>> What else could I try to get my simulation to finish? I have attached
>>>>>>> my
>>>>>>>
>>>>>>>
>>>>>> you need to look at this from a different perspective. first you need
>>>>>> to find out the reason, not try to suppress it.
>>>>>>
>>>>>> how reproducible is this failure? how soon does this happen after the
>>>>>> start of your simulation. have you visualized your simulation around
>>>>>> the time of the failure and seen where exactly an atom experiences
>>>>>> (too) large forces. you may have a problem with bad parameters, or a
>>>>>> bad/incomplete topology (= .psf) file. or you are very very far away
>>>>>> from equilibrium and may need to do multiple iterations of
>>>>>> minimization and relaxation. and so on and so on. there are many ways,
>>>>>> but no simple general solution that works always.
>>>>>>
>>>>>> axel.
>>>>>>
>>>>>>
>>>>>> .conf file so you can see my system and the simulation parameters I am
>>>>>>
>>>>>>> specifying.
>>>>>>>
>>>>>>>
>>>>>>
>>>>>>
>>>>>> cheers,
>>>>>>
>>>>>>> Doug
>>>>>>>
>>>>>>> _____________________________________________________
>>>>>>> Dr. Douglas R. Houston
>>>>>>> Lecturer
>>>>>>> Institute of Structural and Molecular Biology
>>>>>>> Room 3.23, Michael Swann Building
>>>>>>> King's Buildings
>>>>>>> University of Edinburgh
>>>>>>> Edinburgh, EH9 3JR, UK
>>>>>>> Tel. 0131 650 7358
>>>>>>> http://tinyurl.com/douglasrhouston
>>>>>>>
>>>>>>>
>>>>>>> --
>>>>>>> The University of Edinburgh is a charitable body, registered in
>>>>>>> Scotland, with registration number SC005336.
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>
>>>>>> --
>>>>>> Dr. Axel Kohlmeyer akohlmey_at_gmail.com http://goo.gl/1wk0
>>>>>> College of Science & Technology, Temple University, Philadelphia PA, USA
>>>>>> International Centre for Theoretical Physics, Trieste. Italy.
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>> _____________________________________________________
>>>>> Dr. Douglas R. Houston
>>>>> Lecturer
>>>>> Institute of Structural and Molecular Biology
>>>>> Room 3.23, Michael Swann Building
>>>>> King's Buildings
>>>>> University of Edinburgh
>>>>> Edinburgh, EH9 3JR, UK
>>>>> Tel. 0131 650 7358
>>>>> http://tinyurl.com/douglasrhouston
>>>>>
>>>>> --
>>>>> The University of Edinburgh is a charitable body, registered in
>>>>> Scotland, with registration number SC005336.
>>>>>
>>>>>
>>>>>
>>>>>
>>>>
>>>> --
>>>> Aron Broom M.Sc
>>>> PhD Student
>>>> Department of Chemistry
>>>> University of Waterloo
>>>>
>>>>
>>>
>>>
>>>
>>> _____________________________________________________
>>> Dr. Douglas R. Houston
>>> Lecturer
>>> Institute of Structural and Molecular Biology
>>> Room 3.23, Michael Swann Building
>>> King's Buildings
>>> University of Edinburgh
>>> Edinburgh, EH9 3JR, UK
>>> Tel. 0131 650 7358
>>> http://tinyurl.com/douglasrhouston
>>>
>>> --
>>> The University of Edinburgh is a charitable body, registered in
>>> Scotland, with registration number SC005336.
>>>
>>>
>>>
>>
>>
>> --
>> Aron Broom M.Sc
>> PhD Student
>> Department of Chemistry
>> University of Waterloo
>>
>
>
>
>
> _____________________________________________________
> Dr. Douglas R. Houston
> Lecturer
> Institute of Structural and Molecular Biology
> Room 3.23, Michael Swann Building
> King's Buildings
> University of Edinburgh
> Edinburgh, EH9 3JR, UK
> Tel. 0131 650 7358
> http://tinyurl.com/douglasrhouston
>
>
> --
> The University of Edinburgh is a charitable body, registered in
> Scotland, with registration number SC005336.
>
>
>
>

_____________________________________________________
Dr. Douglas R. Houston
Lecturer
Institute of Structural and Molecular Biology
Room 3.23, Michael Swann Building
King's Buildings
University of Edinburgh
Edinburgh, EH9 3JR, UK
Tel. 0131 650 7358
http://tinyurl.com/douglasrhouston

-- 
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.

This archive was generated by hypermail 2.1.6 : Thu Dec 31 2015 - 23:20:56 CST