AMBER Archive (2008)

Subject: Re: AMBER: pmemd iwrap trouble

From: Lars Skjærven (lars.skjarven_at_biomed.uib.no)
Date: Thu Aug 28 2008 - 05:38:02 CDT

Hi there,

Thanks Bob. forgive me for still squawking =)

Thomas C: you mentioned that sander/pmemd imaged in triclinic space,
and that what you would observe after "iwrap=1" is the triclinic box,
instead of the "familiar" octahedron shape(?). so that "iwrap=1"
should be equivalent to the ptraj command "image center", without the
"famililar" keyword(?). Forgive me if I misunderstand this, but what I
observe is that "iwrap=1" is equivalent to the "image center
familiar". I do get the octahedron shape after "iwrap=1" and with a
large part of the protein sticking out, and not the triclinic shape.

Bob: I tried with one step of "iwrap=1" with dt = 0.00001, which
results in similar energies compared to what I had with dt=0.001. now
the energies are -2.1750x10E6 and -2.1749x10E6. I guess that does not
look too bad..

Good energies or not, two of my subunits do still not fit into the
rest of the complex properly. Total rmsd of 70Å when I do not reimage
the wrapped restart file, and 1.6Å when I do reimage the wrapped
restart file (with ptraj).

You mentioned that my subunit complex could be spanning the primary
unit cell, and that two subunits are therefore individually
wrapped(?). This makes sense, but should it not be possible to reimage
this with ptraj and not get an rmsd of 1,6Å? I can clearly see that
two subunits are shifted.

Yes, I am using the latest ptraj version.

Thanks again. Sorry for being such a pain in the neck.. ;-)

Lars

On Wed, Aug 27, 2008 at 3:26 PM, Carlos Simmerling
<carlos.simmerling_at_gmail.com> wrote:
> Hi Bob,
> thanks for looking into this so thoroughly!
> carlos
>
> On Tue, Aug 26, 2008 at 11:17 PM, Robert Duke <rduke_at_email.unc.edu> wrote:
>> Okay gang, I have reviewed the pmemd vs. sander code as it applies to
>> wrapping in general, and wrapping the truncated octahedron in particular.
>> The code, as I said, is basically derived from code that existed in amber 6
>> sander. Since then, there has been 1 bugfix, applied to both sander and
>> pmemd. In addition, in sander 6 there was conditionally compiled code that
>> opted to use the fortran 90 sign() intrinsic instead of a conditional (foo
>> .gt. bar) to determine whether or not to add or subtract a value to
>> coordinates in wrapping molecules for the truncated octahedron. Since amber
>> 6 the code that uses sign() has become the default for sander but not pmemd.
>> This is probably a bit faster, but insignificantly so, given the frequency
>> of wrapping. The other "feature" you get from using sign() is slightly
>> enhanced precision, whereby this function distinguishes between +0.d0 and
>> -0.d0. Using the conditional as pmemd does, all values that evaluate to
>> 0.d0 to within the precision of storage wrap in one direction; for sander
>> using sign(), half the 0.d0 values go one way and the other half go the
>> other way. This is basically an insignificant difference in a rare event,
>> but would have potential to cause sander and pmemd to do a slightly
>> different wrap. I don't think you can do something equivalent in C, and the
>> f90 spec states that sign() is the only place where this distinction is made
>> in f90. So for all intents and purposes, sander and pmemd do the same
>> thing, and it is correct... I do not know what the correct machinations to
>> get what you want in ptraj are, but Tom C. is obviously the authority, so I
>> would do what he says. I tend to just do single step wraps in pmemd to
>> prevent restart overflow, and I never worry about imaging truncated
>> octahedra (I rarely use them). There was an indication from Tom C. that a
>> recent release of ptraj was broken, so that is potentially a source of some
>> grief here. I think we have beat this to death by now, and I hope everyone
>> has what they need (if not, keep squawking); I am reasonably confident that
>> pmemd is doing the "right" thing.
>> Regards - Bob Duke
>>
>> ----- Original Message ----- From: "Robert Duke" <rduke_at_email.unc.edu>
>> To: <amber_at_scripps.edu>
>> Sent: Tuesday, August 26, 2008 3:49 PM
>> Subject: Re: AMBER: pmemd iwrap trouble
>>
>>
>>> You might try wrapping in 1 step of 0.01 fsec (ie., dt = .00001). That
>>> way the before and after delta in all atom positions should be reduced
>>> 100-fold, and you should see a smaller energy delta. You would expect some
>>> difference, as you are really minimizing two structures that differ by more
>>> than just the wrapping operation - they also differ by the distance moved by
>>> all atoms in the dynamics step. Also, using the 32 angstom water buffer
>>> sounds pretty excessive. Is your protein complex perhaps not a good
>>> candidate for fitting in a truncated octahedron? (so a globular protein
>>> approximating a sphere is the best candidate; elongated structures may fit
>>> better in an orthogonal unit cell). Also note my earlier comments about
>>> no_intermolecular_bonds.
>>> Regards - Bob
>>>
>>> ----- Original Message ----- From: "Lars Skjærven"
>>> <lars.skjarven_at_biomed.uib.no>
>>> To: <amber_at_scripps.edu>
>>> Sent: Tuesday, August 26, 2008 2:36 PM
>>> Subject: Re: AMBER: pmemd iwrap trouble
>>>
>>>
>>>> Hi Tom,
>>>> I performed a 1 step energy minimization of the restart files before
>>>> and after the wrapping in pmemd10. The energies are indeed similar,
>>>> but not identical: -2.1748E+06 and -2.1750E+06.
>>>>
>>>> When I do "image center familiar" in ptraj of the restart files before
>>>> and after pmemd wrapping, and then an rmsd in ptraj of the two new
>>>> restart files, it gives an rmsd of 1.6Å. visualizing in vmd shows that
>>>> most of the protein aligns, but two subunits are translated slightly
>>>> with respect to the 10 other subunits.
>>>>
>>>> Ptraj scripts used:
>>>> ---
>>>> trajin b4_wrapping.rst
>>>> trajin after_wrapping.rst
>>>> center :1-6456
>>>> image center familiar
>>>> trajout reimaged.rst restart
>>>> ----
>>>>
>>>> and then,
>>>> ---
>>>> trajin reimaged.rst.1
>>>> trajin reimaged.rst.2
>>>> rms first mass out rmsd.dat @N,CA,C
>>>> ----
>>>>
>>>> Yes, I am still concerned about the validity of my simulation.. :-)
>>>>
>>>> Do I have any other choice to tackle this overflow issue? Can I use
>>>> ptraj to reimage the water back in the box (which works fine), and
>>>> make a new restart file (which I assume do not contain velocities,
>>>> right?), and copy the velocities from the old restart file (which are
>>>> about to overflow)? then, use this new and edited restart file to
>>>> continue my simulation?
>>>>
>>>> When I built my system I had trouble with "solvateoct 10Å" completely
>>>> solvating my protein, instead I had to use "solvateoct 32Å"
>>>> (http://structbio.vanderbilt.edu/archives/amber-archive/2008/1109.php).
>>>> using only 10Å left big chunks of my protein out of the waterbox. for
>>>> any other proteins besides this big monster protein I have worked
>>>> with, "solvateoct 10Å" works fine. can my problem with wrapping have
>>>> anything to do with this issue? probably not.. ?
>>>>
>>>> Thanks for all help on this issue.
>>>>
>>>> Best regards,
>>>> Lars Skjaerven
>>>> University of Bergen, Norway
>>>>
>>>>
>>>> On Tue, Aug 26, 2008 at 6:59 PM, Robert Duke <rduke_at_email.unc.edu> wrote:
>>>>>
>>>>> Hi Tom,
>>>>> Thanks, that was my understanding that the imaging in sander and pmemd
>>>>> was
>>>>> the same, and "equivalent" to a truncated octahedron, but that the image
>>>>> would not look correct in something like vmd. I was hoping to get Tom
>>>>> D. to
>>>>> comment on this since he did the code. I don't use truncated
>>>>> octahedrons
>>>>> that much myself, but have had 0 problems with other imaging, aside from
>>>>> problems with not being quite sure what to do in ptraj sometimes, and
>>>>> also
>>>>> ptraj understandably can't keep up with box size changes associated with
>>>>> a
>>>>> constant pressure simulation (so wrapping a whole trajectory rather than
>>>>> a
>>>>> frame really blows up; one might expect that...). Anyway, thanks for
>>>>> the
>>>>> input. I am going to confirm that my code is identical to sander in
>>>>> amber
>>>>> 10 with regard to wrapping, but that is what I expect. There IS a known
>>>>> problem with the amber 9 release of pmemd, for which a patch was
>>>>> released,
>>>>> so some folks with an unpatched version of 9 may see a few molecules not
>>>>> wrap correctly (so do the patch...).
>>>>> Regards - Bob
>>>>>
>>>>> ----- Original Message ----- From: "Thomas Cheatham III" <tec3_at_utah.edu>
>>>>> To: <amber_at_scripps.edu>
>>>>> Sent: Tuesday, August 26, 2008 12:48 PM
>>>>> Subject: Re: AMBER: pmemd iwrap trouble
>>>>>
>>>>>
>>>>>>
>>>>>>> I get exactly the same results from pmemd9, pmemd10, and sander10 in
>>>>>>> regard the wrapping (iwrap=1); 2 of my subunits are displaced
>>>>>>> (regardless of ioutfm=1 or 0). However, by using the following ptraj
>>>>>>> script (instead of iwrap) it results in a truncated octahedral fitted
>>>>>>> nicely around the protein complex (as expected):
>>>>>>>
>>>>>>> center :1-5332
>>>>>>> image center familiar
>>>>>>> trajout reimaged.rst restart
>>>>>>>
>>>>>>> However, I see no trouble when I use solvatebox instead of solvateoct
>>>>>>> (for both sander and pmemd). Unfortunately I cant go back using
>>>>>>> solvatebox at this time, after 2 months of simulations..
>>>>>>
>>>>>> I've been following this thread and think that it is just a
>>>>>> misunderstanding and that nothing is wrong with the imaging, it just
>>>>>> doesn't look like you would expect. sander/pmemd image in triclinic
>>>>>> space
>>>>>> so for a protein this looks like a slanted box with the protein
>>>>>> potentially sticking out...
>>>>>>
>>>>>> __PPP__
>>>>>> / PPP /
>>>>>> /___PPP/
>>>>>>
>>>>>> ...rather than the more "familiar" truncated octahedron shape. To
>>>>>> verify
>>>>>> this, look at the ptraj imaged trajectory without the familiar keyword;
>>>>>> this will look like what you are seeing from sander/pmemd and they are
>>>>>> equivalent, albeit different ways of looking at things.
>>>>>>
>>>>>> There is also a way in pmemd to generate images in alternative unit
>>>>>> cells,
>>>>>>
>>>>>> image xoffset 1.0 yoffset 0.0 zoffset 0.0
>>>>>>
>>>>>> which you can do for each direction and then "see" the whole packing in
>>>>>> the unit cell.
>>>>>>
>>>>>> If you are still concerned, you could also try a single point energy
>>>>>> minimization on both the "familiar" and alternative restrt file to
>>>>>> verify
>>>>>> the energies are similar.
>>>>>>
>>>>>> imaging in ptraj is not broken (now) to the best of my knowledge; there
>>>>>> was a brief problem with AmberTools 1.0 that lead to errors in the box
>>>>>> size, but this was working in earlier versions and is working in
>>>>>> current
>>>>>> versions.
>>>>>>
>>>>>> --tom
>>>>>>
>>>>>> -----------------------------------------------------------------------
>>>>>> The AMBER Mail Reflector
>>>>>> To post, send mail to amber_at_scripps.edu
>>>>>> To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
>>>>>> to majordomo_at_scripps.edu
>>>>>>
>>>>>
>>>>> -----------------------------------------------------------------------
>>>>> The AMBER Mail Reflector
>>>>> To post, send mail to amber_at_scripps.edu
>>>>> To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
>>>>> to majordomo_at_scripps.edu
>>>>>
>>>> -----------------------------------------------------------------------
>>>> The AMBER Mail Reflector
>>>> To post, send mail to amber_at_scripps.edu
>>>> To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
>>>> to majordomo_at_scripps.edu
>>>>
>>>
>>> -----------------------------------------------------------------------
>>> The AMBER Mail Reflector
>>> To post, send mail to amber_at_scripps.edu
>>> To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
>>> to majordomo_at_scripps.edu
>>>
>>
>> -----------------------------------------------------------------------
>> The AMBER Mail Reflector
>> To post, send mail to amber_at_scripps.edu
>> To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
>> to majordomo_at_scripps.edu
>>
>
>
>
> --
> ===================================================================
> Carlos L. Simmerling, Ph.D.
> Associate Professor Phone: (631) 632-1336
> Center for Structural Biology Fax: (631) 632-1555
> CMM Bldg, Room G80
> Stony Brook University E-mail: carlos.simmerling_at_gmail.com
> Stony Brook, NY 11794-5115 Web: http://comp.chem.sunysb.edu
> ===================================================================
> -----------------------------------------------------------------------
> The AMBER Mail Reflector
> To post, send mail to amber_at_scripps.edu
> To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
> to majordomo_at_scripps.edu
> 

-- 
mvh Lars Skjærven
Institutt for Biomedisin
Universitetet i Bergen
-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber_at_scripps.edu
To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
      to majordomo_at_scripps.edu