[gmx-users] Box dimension size errors in MARTINI soft core simulation
Justin Lemkul
jalemkul at vt.edu
Tue Jul 30 16:44:31 CEST 2013
On 7/30/13 10:17 AM, Scott Pendley wrote:
> Justin,
>
> Thank you for the continued help. I have copied the text of the mdp file
> below. It is a little long, so fair warning. The system is a MARTINI
Perfect, thanks.
> coarse grain representation of 1,2-dipalmitoyl-phosphatidyl-glycerole (LHG)
> surrounded by 700 octanol molecules. Being a neutral MARTINI
> representation, all represented beads or atoms have a charge of 0.0.
> Octanol parameters were taken from MARTINI v 2.0 lipid parameters taken
> from the official website. 1,2-dipalmitoyl-phosphatidyl-glycerole
> parameters were created in house (and seem to work well in aqueous
> simulations). From your note yesterday, it sounded like it may have been a
> rough day yesterday. I hope today treats you better.
>
I appreciate the sentiment, not sure what I've said to give the bad impression.
Yesterday was great. Hopefully nothing I have said has come across as snarky :)
Back to the issue at hand. A couple of things I notice:
> nstxtcout = 1000
Given the frequent .xtc output, you should have a pretty good visualization of
the trajectory. Does it provide any clues?
> Pcoupl = parrinello-rahman
You said you did prior NVT and NPT without free energy options. What barostat
did you use in the previous NPT? Something forgiving like Berendsen?
> Pcoupltype = semiisotropic
I don't understand the use of semi-isotropic coupling here. That's intended for
layered systems like membranes. Sounds to me like you've got a solute in a
uniform solvent, and thus isotropic pressure coupling is appropriate.
> tau_p = 12.0 12.0 ;parrinello-rahman is more stable
> with larger tau-p, DdJ, 20130422
Very long tau_p, indeed. Normally a value no larger than 5.0 should be
necessary. That's not necessarily the source of the problem here, but your
pressure could be running away from you.
> ; and set the free energy parameters
> free-energy = yes
Does turning the free energy code off result in a stable simulation? Your first
message seemed to indicate this would be the case, but I just want to clarify.
> couple-moltype = LHG
> init-lambda = 0.4
> ; these 'soft-core' parameters make sure we never get overlapping
> ; charges as lambda goes to 0
> sc-power = 1
> sc-sigma = 0.3
> sc-alpha = 1.0
> ; we still want the molecule to interact with itself at lambda=0
> couple-intramol = no
> couple-lambda0 = vdw-q
> couple-lambda1 = none
I don't think there's anything wrong with this approach, but I generally do not
(de)couple vdW and charges simultaneously. Might be worthwhile to try a simple
vdw -> none transformation to test.
> foreign-lambda = 0 0.05 0.1 0.2 0.3 0.4 0.5 0.6 0.65 0.7 0.75 0.8
> 0.85 0.9 0.95 1
>
Nothing necessarily wrong here, though you'll get a lot of output that likely
won't improve your result and will probably just slow the simulations down.
What do results at lambda = 0 give you? In theory, this should be a fully
interacting system and would help track down any bugs in the free energy code
itself. If lambda = 0 works, then try a non-zero value.
You mentioned before that the LJ terms are affected over time, but based on this
.mdp file, they shouldn't since you're not doing any sort of slow-growth method.
Which terms are affected? By how much? Does the time at which you see any
major change correspond to any observable weirdness in the trajectory?
-Justin
>
> On Mon, Jul 29, 2013 at 4:39 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>
>>
>>
>> On 7/29/13 2:12 PM, Scott Pendley wrote:
>>
>>> Thank you, Justin. I am using gromacs version 4.5.5 and have attached
>>> .mdp
>>> file. I followed your advise and pointer to trouble shooting the system
>>>
>>
>> The mailing list does not accept attachments. Please either copy and
>> paste its contents or provide a link from which it can be downloaded.
>>
>>
>> and decomposed the energy to find potential sources for the problem. The
>>> simulation ran for 978 ps, the total energy changed from approximately
>>> -13000 kj/mol to -15000 kj/mol and plateau'ed out around 900 ps without a
>>> big energy blow-out. The majority of the change in energy was from the LJ
>>> contribution which is the expected source with the disappearing atoms of
>>> my
>>> ligand. Of more concern is the change in box coordinates/shape. The
>>> simulation cell starts with dimensions 9.45nm x 9.45nm x 9.45nm and ends
>>> at
>>> 2.8 x 2.8 x 45.8 nm which does seem to give rise to the error that I
>>> originally noted. While some change in volume and cell dimensions is
>>> expected, these changes seems to be a focused solely on expanding along
>>> the
>>> Z-coordinate. Is there any way to constrain the ratio of x:y:z
>>> coordinates
>>> during npt simulations or do you know of something that I may be missing
>>> in
>>> my simulations?
>>>
>>>
>> The magnitude of change you're seeing suggests extreme instability and
>> thus there is no real "hack" to simply "make it work." Isotropic pressure
>> coupling is normally sufficient for this purpose, but something seems to be
>> going very wrong in your system. The complete .mdp file, and a description
>> of what your ligand is (or at least how big) would be helpful. You've got
>> a very large box for what one would normally consider a "ligand."
>>
>> -Justin
>>
>> Thank you,
>>>
>>> Scott
>>>
>>>
>>> On Wed, Jul 24, 2013 at 4:59 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>>>
>>>
>>>>
>>>> On 7/24/13 4:33 PM, Scott Pendley wrote:
>>>>
>>>> I am fairly new to gromacs and I am trying to run a thermodynamic
>>>>> integration simulation of a ligand disappearing in a box of octanol at a
>>>>> single set lambda point. I have previous successful nvt and npt runs of
>>>>> this system. When I have added the free energy portions to the input
>>>>> file,
>>>>> I get the following error:
>>>>>
>>>>> Fatal error:
>>>>> One of the box vectors has become shorter than twice the cut-off length
>>>>> or
>>>>> box_yy-|box_zy| or box_zz has become smaller than the cut-off.
>>>>> For more information and tips for troubleshooting, please check the
>>>>> GROMACS
>>>>> website at http://www.gromacs.org/****Documentation/Errors<http://www.gromacs.org/**Documentation/Errors>
>>>>> <http://**www.gromacs.org/Documentation/**Errors<http://www.gromacs.org/Documentation/Errors>
>>>>>>
>>>>>
>>>>>
>>>>> This seems unusual. The box dimensions are 9.45 nm x 9.45 x 9.45 nm so
>>>>> that is fairly large even accounting for some shrinkage with a
>>>>> disappearing
>>>>> ligand.
>>>>>
>>>>>
>>>>> The available information suggests the system has become unstable and
>>>> is
>>>> imploding. See general troubleshooting information at
>>>> http://www.gromacs.org/****Documentation/Terminology/**<http://www.gromacs.org/**Documentation/Terminology/**>
>>>> Blowing_Up#Diagnosing_an_****Unstable_System<http://www.**
>>>> gromacs.org/Documentation/**Terminology/Blowing_Up#**
>>>> Diagnosing_an_Unstable_System<http://www.gromacs.org/Documentation/Terminology/Blowing_Up#Diagnosing_an_Unstable_System>
>>>>>
>>>>
>>>> .
>>>>
>>>>
>>>> Cutoffs in the input file are set as follows: rlist =1.4, rvdw=1.2,
>>>>
>>>>> rcoulomb=1.2. Doubling any of them would still be less than 3 nm which
>>>>> is
>>>>> significantly smaller than the box size. Is there anything I am missing
>>>>> or
>>>>> any suggestions that others can give me?
>>>>>
>>>>>
>>>>> I wouldn't mess with the cutoffs; they're an essential part of the
>>>> force
>>>> field. For further diagnostics, please consider the points above and
>>>> provide your .mdp file and Gromacs version.
>>>>
>>>> -Justin
>>>>
>>>> --
>>>> ==============================****====================
>>>>
>>>>
>>>> Justin A. Lemkul, Ph.D.
>>>> Postdoctoral Fellow
>>>>
>>>> Department of Pharmaceutical Sciences
>>>> School of Pharmacy
>>>> Health Sciences Facility II, Room 601
>>>> University of Maryland, Baltimore
>>>> 20 Penn St.
>>>> Baltimore, MD 21201
>>>>
>>>> jalemkul at outerbanks.umaryland.****edu <jalemkul at outerbanks.**
>>>> umaryland.edu <jalemkul at outerbanks.umaryland.edu>> | (410)
>>>> 706-7441
>>>>
>>>> ==============================****====================
>>>>
>>>> --
>>>> gmx-users mailing list gmx-users at gromacs.org
>>>> http://lists.gromacs.org/****mailman/listinfo/gmx-users<http://lists.gromacs.org/**mailman/listinfo/gmx-users>
>>>> <htt**p://lists.gromacs.org/mailman/**listinfo/gmx-users<http://lists.gromacs.org/mailman/listinfo/gmx-users>
>>>>>
>>>> * Please search the archive at http://www.gromacs.org/**
>>>> Support/Mailing_Lists/Search<h**ttp://www.gromacs.org/Support/**
>>>> Mailing_Lists/Search<http://www.gromacs.org/Support/Mailing_Lists/Search>>before
>>>> posting!
>>>>
>>>> * Please don't post (un)subscribe requests to the list. Use the www
>>>> interface or send it to gmx-users-request at gromacs.org.
>>>> * Can't post? Read http://www.gromacs.org/****Support/Mailing_Lists<http://www.gromacs.org/**Support/Mailing_Lists>
>>>> <http://**www.gromacs.org/Support/**Mailing_Lists<http://www.gromacs.org/Support/Mailing_Lists>
>>>>>
>>>>
>>>>
>>>>
>>>>
>> --
>> ==============================**====================
>>
>> Justin A. Lemkul, Ph.D.
>> Postdoctoral Fellow
>>
>> Department of Pharmaceutical Sciences
>> School of Pharmacy
>> Health Sciences Facility II, Room 601
>> University of Maryland, Baltimore
>> 20 Penn St.
>> Baltimore, MD 21201
>>
>> jalemkul at outerbanks.umaryland.**edu <jalemkul at outerbanks.umaryland.edu> | (410)
>> 706-7441
>>
>> ==============================**====================
>> --
>> gmx-users mailing list gmx-users at gromacs.org
>> http://lists.gromacs.org/**mailman/listinfo/gmx-users<http://lists.gromacs.org/mailman/listinfo/gmx-users>
>> * Please search the archive at http://www.gromacs.org/**
>> Support/Mailing_Lists/Search<http://www.gromacs.org/Support/Mailing_Lists/Search>before posting!
>> * Please don't post (un)subscribe requests to the list. Use the www
>> interface or send it to gmx-users-request at gromacs.org.
>> * Can't post? Read http://www.gromacs.org/**Support/Mailing_Lists<http://www.gromacs.org/Support/Mailing_Lists>
>>
--
==================================================
Justin A. Lemkul, Ph.D.
Postdoctoral Fellow
Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 601
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201
jalemkul at outerbanks.umaryland.edu | (410) 706-7441
==================================================
More information about the gromacs.org_gmx-users
mailing list