[gmx-users] Box size/type confusion for bilayer system

[Justin Lemkul]

On 10/18/12 3:04 PM, klexa wrote:
> On 10/18/2012 11:53 AM, Justin Lemkul wrote:
>>
>>
>> On 10/18/12 2:43 PM, klexa wrote:
>>> Hi Gromacs users,
>>>
>>> I think I am a bit confused about the proper way to handle boxes that are not
>>> standard cubes. I'm trying to run a membrane simulation where a cyclic
>>> undecapeptide is inserted into the membrane and I want the water layer to be
>>> sufficiently thick that if it were pulled, the peptide could be fully solvated
>>> by the water. To avoid having an enormous box of membrane and water, I have an
>>> orthorhombic box containing my peptide and bilayer. It minimizes alright with
>>> Gromacs, but when I go to equilibrate it it fails because it's too skewed to be
>>> a triclinic box. I've tried modifying the box with editconf and converting it to
>>> a rhombic dodecahedron, sort of like the manual suggests for a membrane system.
>>> I'm not sure that even that is sensible since it seems like I would be losing
>>> content that way, yet nothing is clipped, and I did this after using trjconv to
>>> remove any periodicity from my prior simulation of this system (in Desmond) but
>>> doing so gives me a starting potential energy of NaN for the new system that I
>>> obviously cannot work around. Is what I am trying to do even possible? If it is,
>>> it seems like there is probably a better way than the way I chose, so if you
>>> have any suggestions, I would be greatly appreciative.
>>>
>>
>> I have never produced a membrane system with a hexagonal cross-section like
>> the manual describes.  The most straightforward approach in my mind is simply
>> a rectangular box.  It will save you a ton of headaches.
>>
>
> Okay, yes, it does seem much simpler. But if I can indeed just use a rectangular
> box like 7.7 7.7 10.5, why does Gromacs fail with the "triclinic too skewed" error?
>
> { -1.75e+25     0                   -0
> -0                    -1.75e+25      -0
> -0                      -0                -2.49e+25}
>

The values shown here indicate that the unit cell has gone completely haywire. 
They clearly bear no resemblance to the values you hope to use.

> Maybe it's just related to this force field mixing, but otherwise, if I should
> be able to proceed with a rectangular box, does that need to be specified
> somewhere outside of when I use genbox to solvate my system with -box 7.7 7.7 10.5?
>

I find it better to use editconf -box to build the system and simultaneously 
position the molecules within the unit cell with -center.  Then follow with 
solvation.

>
>>> I'm trying to run this simulation with AMBER FF99SB parameters for the peptide,
>>> Tieleman's lipid parameters for POPC, and SPCE waters, so just as a sanity
>>> check, is it reasonable to consider a system like that?
>>>
>>
>> I don't know how this would even run.  The AMBER protein force field and
>> Berger lipid paramters use different combination rules, and I have never seen
>> a demonstration that one can use them together.  It is most straightforward to
>> use a Gromos force field or OPLS-AA with modifications to account for the
>> changes in combination rules.
>
> Great, I'm glad to hear it as I was skeptical too. I can tell you that after 2
> months of trying to get a system like working, it still hasn't succeeded in any
> form, so the odds are not in its favor.

There is a tutorial available explaining the force field logic, if you're 
interested.

http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/membrane_protein/index.html

-Justin

-- 
========================================

Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

========================================