Re: [gmx-users] RE: Gibbs Energy Calculation and charges

2013-10-30 Thread Michael Shirts
I likely won't have much time to look at it tonight, but you can see
exactly what the option is doing to the topology.  run gmxdump on the
tpr.  All of the stuff that couple-intramol does is in grompp, so the
results will show up in the detailed listings of the interactions, and
which ones have which values set for the A and B states.

On Wed, Oct 30, 2013 at 5:36 PM, Dallas Warren dallas.war...@monash.edu wrote:
 Michael, thanks for taking the time to comment and have a look.

 The real issue I am having is a bit deeper into the topic than that, my last 
 reply was just an observation on something else.  Will summarise what I have 
 been doing etc.

 I have a molecule that are calculating the Gibbs energy of hydration and 
 solvation (octanol).  In a second topology the only difference is that the 
 atomic charges have been doubled.  Considering that charges are scaled 
 linearly with lambda, the normal charge values of dH/dl from lambda 0 to 1 
 obtained should reproduce that of the double charged molecule from lambda 0.5 
 to 1.0.  Is that a correct interpretation?  Since the only difference should 
 be that charge of the atoms and over that range the charge will be identical.

 I was using couple-intramol = no and the following are the results from those 
 simulations.

 For the OE atom within the molecule, I have plotted the following graphs of 
 dH/dl versus charge of that atom for both of the topologies.
 octanol - http://ozreef.org/stuff/octanol.gif
 water - http://ozreef.org/stuff/water.gif
 mdp file - http://ozreef.org/stuff/gromacs/mdout.mdp

 The mismatch between the two topologies is the real issue that I am having.  
 I was hoping to get the two to overlap.

 My conclusion based on this is that there is actually something else being 
 changed with the topology by GROMACS when the simulations are being run.  The 
 comments in the manual allude to that, but not entirely sure what is going on.

 From the manual:

couple-intramol:

no
 All intra-molecular non-bonded interactions for moleculetype
couple-moltype are replaced by exclusions and explicit pair
interactions. In this manner the decoupled state of the molecule
corresponds to the proper vacuum state without periodicity effects.
yes
 The intra-molecular Van der Waals and Coulomb interactions are
also turned on/off. This can be useful for partitioning free-energies
of relatively large molecules, where the intra-molecular non-bonded
interactions might lead to kinetically trapped vacuum conformations.
The 1-4 pair interactions are not turned off.

 Chris Neale commented:
Ah, I see. I guess that you are using couple-intramol = no (the default in 
v4.6.3 at least). That means
that the intramolecular charge-charge interactions are always at 
full-strength (and therefore different).
I would expect that normal at lambda=0 should be the same as double at 
lambda=0.5 only for
couple-intramol = yes
If you were using couple-intramol = yes already, then I am as confused as you 
are.

 From Chris' comment, I then went back and repeated the calculation using 
 couple-intramol = yes with the results of this in the below image (plus the 
 previous when set to no for comparison)

 http://ozreef.org/stuff/gromacs/couple-intramol.png

 My rather garbled comment that you saw was my attempt to make sense of the 
 fact that this setting made things worse (and try to add something to get 
 this issue to have another pass on the list, hoping someone like yourself 
 will comment), and the value of dH/dl for c-i=yes at lambda=1 matches with 
 c-i=no at lambad=0.  This you can see with the previously linked to graphs.

 Hopefully that helps understand what I have been doing and the issues.

 Catch ya,

 Dr. Dallas Warren
 Drug Delivery, Disposition and Dynamics
 Monash Institute of Pharmaceutical Sciences, Monash University
 381 Royal Parade, Parkville VIC 3052
 dallas.war...@monash.edu
 +61 3 9903 9304
 -
 When the only tool you own is a hammer, every problem begins to resemble a 
 nail.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] RE: Gibbs Energy Calculation and charges

2013-10-29 Thread Michael Shirts
I think the grammar got a little garbled there, so I'm not sure quite
what you are claiming.

One important thing to remember; 1-4 interactions are treated as
bonded interactions right now FOR COUPLE intramol (not for lambda
dependence of the potential energy function), so whether
couple-intramol is set to yes or no does not affect these interactions
at all.  It only affects the nonbondeds with distances greater than
1-5.  At least to me, this is nonintuitive (and we're coming up with a
better scheme for 5.0), but might that explain what you are getting?

On Tue, Oct 29, 2013 at 9:44 PM, Dallas Warren dallas.war...@monash.edu wrote:
 Just want this to make another pass, just in case those in the know missed it.

 Using couple-intrmol = yes the resulting dH/dl plot actually looks like that 
 at lamba = 1 it is actually equal to couple-intramol = no with lambda = 0.

 Should that be the case?

 Catch ya,

 Dr. Dallas Warren
 Drug Delivery, Disposition and Dynamics
 Monash Institute of Pharmaceutical Sciences, Monash University
 381 Royal Parade, Parkville VIC 3052
 dallas.war...@monash.edu
 +61 3 9903 9304
 -
 When the only tool you own is a hammer, every problem begins to resemble a 
 nail.

 -Original Message-
 From: gmx-users-boun...@gromacs.org [mailto:gmx-users-
 boun...@gromacs.org] On Behalf Of Dallas Warren
 Sent: Tuesday, 22 October 2013 2:49 PM
 To: Discussion list for GROMACS users
 Subject: [gmx-users] RE: Gibbs Energy Calculation and charges

 I have completed the simulations using couple-intramol = yes

 Reminder, have a molecule that are calculating the Gibbs energy of
 hydration and solvation (octanol).  In one topology the only difference
 is the atomic charges have been doubled.  Considering that charges are
 scaled linearly with lambda, the normal charge values of dH/dl from
 lambda 0 to 1 obtained should reproduce that of the double charged
 molecule from lambda 0.5 to 1.0.

 Here is the mdout.mdp file from one of the simulations using couple-
 intramol = yes.  The difference with the simulations run using couple-
 intramol = no is only that setting, checked using the diff command.
 Copy of the mdout.mdp file can be found:

 http://ozreef.org/stuff/gromacs/mdout.mdp

 Changing to couple-intrmol = yes makes things significantly worse.

 Here are the plots of dH/dl obtained as a function of the charge of one
 of the atoms, OE (calculated based on fact that at lambda = 0 fully
 charged, lambda = 1 no charged)

 http://ozreef.org/stuff/gromacs/couple-intramol.png

 Anyone with some insight on what is going on here?

 From the manual:
 
 couple-intramol:

 no
 All intra-molecular non-bonded interactions for moleculetype
 couple-moltype are replaced by exclusions and explicit pair
 interactions. In this manner the decoupled state of the molecule
 corresponds to the proper vacuum state without periodicity effects.
 yes
 The intra-molecular Van der Waals and Coulomb interactions are
 also turned on/off. This can be useful for partitioning free-energies
 of relatively large molecules, where the intra-molecular non-bonded
 interactions might lead to kinetically trapped vacuum conformations.
 The 1-4 pair interactions are not turned off.
 

 Catch ya,

 Dr. Dallas Warren
 Drug Delivery, Disposition and Dynamics
 Monash Institute of Pharmaceutical Sciences, Monash University
 381 Royal Parade, Parkville VIC 3052
 dallas.war...@monash.edu
 +61 3 9903 9304
 -
 When the only tool you own is a hammer, every problem begins to
 resemble a nail.


  -Original Message-
  From: gmx-users-boun...@gromacs.org [mailto:gmx-users-
  boun...@gromacs.org] On Behalf Of Dallas Warren
  Sent: Thursday, 17 October 2013 3:32 PM
  To: Discussion list for GROMACS users
  Subject: [gmx-users] RE: Gibbs Energy Calculation and charges
 
  Chris,
 
  Thank you, that appears to be the issue then.
 
  Running them again now with couple-intramol = yes
 
  Will report back once that is completed with the results.
 
  Catch ya,
 
  Dr. Dallas Warren
  Drug Delivery, Disposition and Dynamics
  Monash Institute of Pharmaceutical Sciences, Monash University
  381 Royal Parade, Parkville VIC 3052
  dallas.war...@monash.edu
  +61 3 9903 9304
  -
  When the only tool you own is a hammer, every problem begins to
  resemble a nail.
 
 
   -Original Message-
   From: gmx-users-boun...@gromacs.org [mailto:gmx-users-
   boun...@gromacs.org] On Behalf Of Christopher Neale
   Sent: Thursday, 17 October 2013 2:15 PM
   To: gmx-users@gromacs.org
   Subject: [gmx-users] Gibbs Energy Calculation and charges
  
   Ah, I see. I guess that you are using couple-intramol = no (the
  default
   in v4.6.3 at least). That means
   that the intramolecular charge-charge interactions are always at
  full-
   strength (and therefore different).
  
   I would expect that 

Re: [gmx-users] lmc-stats

2013-10-28 Thread Michael Shirts
Hi, Andrew-

The choice of lmc-stats only affects the calculation of the weights,
it does not affect the calculation of the transition matrix.

There are two possibilities for the transition matrix; the estimated
one (which is just called 'Transition Matrix' -- we should probably
have a better name), and the empirical one.  The empirical transition
matrix is just transition counts. The first one is an ensemble average
of P(k|x); for neighbor moves, this is just the average of the barker
transition probabilities (plus the choice of going up or down),

Let me know if this is useful.  Not that many people have used this
code, thus there are likely ways that it can be improved for better
utility.  I generally haven't tried to calculate free energy
differences from the transition matrix, though it should give
consistent results with the other methods.

Best,

Michael Shirts
Assistant Professor
Department of Chemical Engineering
University of Virginia
michael.shi...@virginia.edu
(434)-243-1821

On Mon, Oct 28, 2013 at 11:13 PM, Andrew S. Paluch paluc...@miamioh.edu wrote:
 When performing free energy calculations using the expanded ensemble method,
 there is a barker and metropolis option for lmc-stats. Are the corresponding
 transition probabilities computed with or without the weighting factors?
 That is, are these probabilities biased or not?

 Thank you,

 Andrew

 --

 Andrew S. Paluch, PhD
 Department of Chemical, Paper, and Biomedical Engineering
 Miami University
 paluc...@miamioh.edu
 (513) 529-0784
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Replica Exchange with Solute Tempering

2013-10-26 Thread Michael Shirts
Hi, all-

Rest essentially scales the solute-solvent interactions, but maintains
the solute-solute interactions. This can be done solely with
Hamiltonian replica exchange, which is in 4.6.  It's a bit tricky,
though.  We plan on having something that does this automatically in
5.0 or 5.1, but it's not there yet.

What do you intend to do?  It could be that there's a better way to do
what you want to do with Hamiltonian replica exchange, as well.



On Sat, Oct 26, 2013 at 3:52 PM, David Osguthorpe
david.osgutho...@gmail.com wrote:

 Hi

 Ive been looking for a sample input for a REST simulation but
 so far have not been able to find one.

 Does anybody have an example to share?

 Ive seen the tutorial by Mark which mentions REST but I dont
 see this included in the archived examples.

 Am I missing something?

 Thanks for any help

 David
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] energy drift - comparison of double and single precision

2013-10-25 Thread Michael Shirts
Hi, all-

At this point, any fixes are going to be in the 5.0 version, where the
integrators will be a bit different.  If you upload your system files
to redmine.gromacs.org (not just the .mdp), then I will make sure this
gets tested there.

On Fri, Oct 25, 2013 at 10:14 AM, Guillaume Chevrot
guillaume.chev...@gmail.com wrote:
 Dear Xavier,

 2013/10/12 XAvier Periole x.peri...@rug.nl


 Could you try to reduce the nstcalcenergy flag from 100 to 10 and then one?


 FYI, I tried 10 and 1 and the energy drift is exactly the same.




  Similar flags apply to temperature and pressure and I believe might
 seriously affect energy conservation.

 XAvier.

  On Oct 12, 2013, at 0:50, Mark Abraham mark.j.abra...@gmail.com wrote:
 
  Didn't see any problem in the .mdp. -4500 kJ/mol in 10ns over (guessing)
  30K atoms is 0.015 kJ/mol/ns/atom. k_B T is at least 100 times larger
 than
  that. I bet the rest of the lysozyme model physics is not accurate to
 less
  than 1% ;-) There are some comparative numbers at
  http://dx.doi.org/10.1016/j.cpc.2013.06.003 - the two systems are rather
  different but they share the use of SETTLE.
 
  Note that using md-vv guarantees the 2007 paper is inapplicable, because
  GROMACS did not have a velocity Verlet integrator back then. Sharing the
  .log files might be informative.
 
  Mark
 
 
  On Fri, Oct 11, 2013 at 11:38 PM, Guillaume Chevrot 
  guillaume.chev...@gmail.com wrote:
 
  Hi,
 
  sorry for my last post! I re-write my e-mail (with some additional
  information) and I provide the links to my files ;-)
 
  I compared the total energy of 2 simulations:
  lysozyme in water / NVE ensemble / single precision / Gromacs 4.6.3
  lysozyme in water / NVE ensemble / double precision / Gromacs 4.6.3
 
  ... and what I found was quite ... disturbing (see the plots of the
 total
  energy: http://dx.doi.org/10.6084/m9.figshare.820153). I observe a
  constant
  drift in energy in the case of the single precision simulation.
 
  Did I do something wrong*? Any remarks are welcomed! Here is the link to
  the ‘mdout.mdp’ file (http://dx.doi.org/10.6084/m9.figshare.820154) so
 you
  can check what mdp options I used.
 
  My second question is: if I did not do something wrong, what are the
  consequences on the simulation? Can I trust the results of single
 precision
  simulations?
 
  Regards,
 
  Guillaume
 
  *PS: I am not the only one encountering this behavior. In the
 literature,
  this problem has already been mentioned:
  http://jcp.aip.org/resource/1/jcpsa6/v126/i4/p046101_s1
 
 
 
 
  2013/10/11 Mark Abraham mark.j.abra...@gmail.com
 
  On Oct 11, 2013 7:59 PM, Guillaume Chevrot 
  guillaume.chev...@gmail.com
  wrote:
 
  Hi all,
 
  I recently compared the total energy of 2 simulations:
  lysozyme in water / NVE ensemble / single precision
  lysozyme in water / NVE ensemble / double precision
 
  ... and what I found was quite ... disturbing (see the attached figure
  -
  plots of the total energy). I observe a constant drift in energy in
 the
  case of the single precision simulation.
 
  Did I do something wrong*? Any remarks are welcomed! I join the
  ‘mdout.mdp’
  file so you can check what mdp options I used.
 
  Maybe. Unfortunately we cannot configure the mailing list to allow
 people
  to send attachments to thousands of people, so you will need to do
  something like provide links to files on a sharing service.
 
 
  My second question is: if I did not do something wrong, what are the
  consequences on the simulation? Can I trust the results of single
  precision
  simulations?
 
  Yes, as you have no doubt read in the papers published by the GROMACS
  team.
 
  Regards,
 
  Guillaume
 
  *PS: I am not the only one encountering this behavior. In the
  literature,
  this problem has already been mentioned:
  http://jcp.aip.org/resource/1/jcpsa6/v126/i4/p046101_s1
 
  ... which is six years old, examining the properties of code seven
 years
  old. Life has moved on! :-) Even if you have found a problem, it is a
 big
  assumption that this is (still) the cause.
 
  Mark
 
  --
  gmx-users mailing listgmx-users@gromacs.org
  http://lists.gromacs.org/mailman/listinfo/gmx-users
  * Please search the archive at
  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-requ...@gromacs.org.
  * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
  --
  gmx-users mailing listgmx-users@gromacs.org
  http://lists.gromacs.org/mailman/listinfo/gmx-users
  * Please search the archive at
  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-requ...@gromacs.org.
  * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
  --
  gmx-users mailing listgmx-users@gromacs.org
  http://lists.gromacs.org/mailman/listinfo/gmx-users
  

Re: [gmx-users] a new GROMACS simulation tool

2013-10-22 Thread Michael Shirts
Is there a link to the documentation?  It's a little difficult to know
exactly what this supposed to be doing.  Is it a GUI interface to
gromacs?

In general, it would be great to get these sort of extensions
coordinated with the main gromacs development tree, since otherwise
they would tend to get out of sync, lowering utility to everyone.

On Tue, Oct 22, 2013 at 10:34 AM, Kevin Chen fch6...@gmail.com wrote:
 Hi Everyone,

 I'm writing to let you guys know that we have developed a web-based tool MD
 simulation tool for GROMACS.  It is a software package primarily developed
 for biological MD and offers a huge amount of possible options and settings
 for tailoring the simulations. Seamlessly integrated with newly developed
 GUI interfaces, the tool provides comprehensive setup, simulation, analysis
 and job submission tools. Most importantly, unlike other GROMACS GUI
 applications, user can actually run really simulations using the dedicated
 HPC resources. That been said, there's no proposal and installation
 required.  This tool could be a great fit for both teaching and research
 projects. Users inexperienced in MD can work along prepared workflows, while
 experts may enjoy a significant relief from the tedium of typing and
 scripting. As for now, we'd like to invite people to participate in user
 testing on this newly developed tool. Let me know if you'd like to try it
 out. We will set up an account for you.

 Best Regards,

 Kevin Chen, Ph.D.
 Information Technology at Purdue (ITaP)
 West Lafayette, IN 47907-2108

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] energy drift - comparison of double and single precision

2013-10-11 Thread Michael Shirts
Hi, Guillaume-

No one can tell you if you did anything wrong if you didn't tell us
what you did!  There are literally thousands of combinations of
options in running an NVE simulation, a substantial fraction of which
are guaranteed not to conserve energy.

If you post the files (inputs and relevant output files) to
redmine.gromacs.org, then people would be able to see if it was an
error in how you ran the simulation, or an actual bug.



On Fri, Oct 11, 2013 at 1:58 PM, Guillaume Chevrot
guillaume.chev...@gmail.com wrote:
 Hi all,

 I recently compared the total energy of 2 simulations:
 lysozyme in water / NVE ensemble / single precision
 lysozyme in water / NVE ensemble / double precision

 ... and what I found was quite ... disturbing (see the attached figure -
 plots of the total energy). I observe a constant drift in energy in the
 case of the single precision simulation.

 Did I do something wrong*? Any remarks are welcomed! I join the ‘mdout.mdp’
 file so you can check what mdp options I used.

 My second question is: if I did not do something wrong, what are the
 consequences on the simulation? Can I trust the results of single precision
 simulations?

 Regards,

 Guillaume

 *PS: I am not the only one encountering this behavior. In the literature,
 this problem has already been mentioned:
 http://jcp.aip.org/resource/1/jcpsa6/v126/i4/p046101_s1

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] question about OPLS-AA force field -required bond constraints

2013-10-10 Thread Michael Shirts
OPLS-AA was generally derived with Monte Carlo, which means that all
bonds were exactly constrained.  But read the papers!

On Thu, Oct 10, 2013 at 12:27 PM, Justin Lemkul jalem...@vt.edu wrote:


 On 10/10/13 12:17 PM, Martin, Erik W wrote:


 Hi, I'm new to both Gromacs and OPLS.  I have always used CHARMM or AMBER
 before so am a bit confused about some idiosyncrasies.

 When I go to a production run on my system (a protein in a truncated
 octahedron box of Tip4 water) I have removed all constraints.  When I run
 grompp, I get an error telling me that a bond between a HG proton and OG
 oxygen of protein has an oscillation period of 9e-3 ps and is incompatible
 with my 1 fs time step.  I am used to being able to run a completely
 unconstrained simulation with a 1 fs time step… so I'm assuming this is a
 property of OPLS-AA?


 No, it's not.  It's a property of the bonds.  Nearly all force fields (if
 not all of them) constrain X-H bonds at minimum, to remove the fastest
 motions in the system and allow for a longer time step.  Using dt on the
 same order as the fastest frequencies leads to instability.  See manual
 section 6.7.


 Could someone please let me know what the bond constraint requirements are
 for OPLS-AA and GROMOS forcefields (I'm going to need to use that next).
 I
 can't seem to find this information by searching but imagine its on a
 website
 somewhere?


 Information like this should be in the primary literature for the parameter
 sets.  I don't know of any sort of database online that tracks such things,
 and moreover there are no absolute rules AFAIK.  I generally follow the
 methods of the paper that derived the parameter set.  Unfortunately,
 sometimes they don't say.  For Gromos96 53A5/53A6, validation simulations
 were carried out with all bonds constrained.

 -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

 jalem...@outerbanks.umaryland.edu | (410) 706-7441

 ==

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Problem running free energy simulations

2013-10-03 Thread Michael Shirts
Hmm. This really isn't quite enough information to go on.  Can you
file a redmine issue, and include the files used to generate the run
that crashed (.mdp, .gro,. .top), as well as the files that show it
failing (.log)?

http://redmine.gromacs.org/

On Thu, Oct 3, 2013 at 4:32 AM, Jernej Zidar jernej.zi...@gmail.com wrote:
 Dear Michael.
   The simulations at each lambda point starts from the same structure
 that I equilibrated (NPT ensemble) for 20 nanoseconds. The system has
 ~7500 atoms in a box the size 5 nm x 5 nm x 5 nm.  The molecule of
 interest is located in the center of the unit cell.

 Thanks,
 JErnej

 Sounds like the simulation is blowing up.  How soon does it start crashing.

 Also, what configurations are you using to start your free energy
 simulations at each lambda?




 --
 Windows: Re-Boot, Linux: Be-Root.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Problem running free energy simulations

2013-10-02 Thread Michael Shirts
Sounds like the simulation is blowing up.  How soon does it start crashing.

Also, what configurations are you using to start your free energy
simulations at each lambda?

On Wed, Oct 2, 2013 at 10:31 PM, Jernej Zidar jernej.zi...@gmail.com wrote:
 Hi all,
   I'm trying to determine the free energy of solvation for a molecule
 in n-octanol. I'm separately turning off the Coulomb and Lennard-Jones
 interactions as instructed in the free energy tutorial. The
 Lennard-Jones simulations keep on crashing for most values of lambda
 with the message:
 Program mdrun, VERSION 4.6.3
 Source code file: /home/zidar/gromacs-4.6.3/src/mdlib/domdec_top.c, line: 393

 Fatal error:
 1 of the 79007 bonded interactions could not be calculated because
 some atoms involved moved further apart than the multi-body cut-off
 distance (1 nm) or the two-body cut-off distance (1 nm), see option
 -rdd, for pairs and tabulated bonds also see option -ddcheck
 For more information and tips for troubleshooting, please check the GROMACS
 website at http://www.gromacs.org/Documentation/Errors

 - - -

   I tried using -noddcheck, -pd (separate simulations) but the
 simulations still fail.

   The molecule surveyed is not charged (branched hydrocarbon) and I
 can simulate it using any ensemble without a problem. I'm using CgenFF
 without adding any new atom types so all the interactions should be
 properly accounted for.

   Any advice on how to solve this will be greatly appreciated. I was
 able to perform free energy simulations of a similar molecule in
 octanol without any major problem.

 Thanks,
 Jernej
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] segfault when running the alchemistry.org ethanol solvation free energy tutorial

2013-09-26 Thread Michael Shirts
Just to be clear, is this the expanded ensemble version of the calculation?

On Thu, Sep 26, 2013 at 5:25 PM, Mark Abraham mark.j.abra...@gmail.com wrote:
 I found the -multi version of that tutorial a bit temperamental...
 Michael Shirts suggested that double precision is more reliable for
 expanded ensemble. Hopefully he can chime in in a day or two.

 Mark

 On Thu, Sep 26, 2013 at 9:00 PM, Christopher Neale
 chris.ne...@mail.utoronto.ca wrote:
 Dear Users:

 Has anyone successfully run the free energy tutorial at 
 http://www.alchemistry.org/wiki/GROMACS_4.6_example:_Direct_ethanol_solvation_free_energy
  ?

 I just tried it and I get a segmentation fault immediately (see output at 
 the end of this post).

 I get a segfault with both 4.6.3 and 4.6.1.

 Note that if I modify the .mdp file to set free-energy = no , then the 
 simulation runs just fine. (I have, of course, set init-lambda-state in the 
 .mdp file that I downloaded from the aforementioned site and I get a 
 segfault with any value of init-lambda-state from 0 to 8).

 gpc-f103n084-$ mdrun -nt 1 -deffnm ethanol.1 -dhdl ethanol.1.dhdl.xvg
  :-)  G  R  O  M  A  C  S  (-:

   GROup of MAchos and Cynical Suckers

 :-)  VERSION 4.6.3  (-:

 Contributions from Mark Abraham, Emile Apol, Rossen Apostolov,
Herman J.C. Berendsen, Aldert van Buuren, Pär Bjelkmar,
  Rudi van Drunen, Anton Feenstra, Gerrit Groenhof, Christoph Junghans,
 Peter Kasson, Carsten Kutzner, Per Larsson, Pieter Meulenhoff,
Teemu Murtola, Szilard Pall, Sander Pronk, Roland Schulz,
 Michael Shirts, Alfons Sijbers, Peter Tieleman,

Berk Hess, David van der Spoel, and Erik Lindahl.

Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  Copyright (c) 2001-2012,2013, The GROMACS development team at
 Uppsala University  The Royal Institute of Technology, Sweden.
 check out http://www.gromacs.org for more information.

  This program is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public License
 as published by the Free Software Foundation; either version 2.1
  of the License, or (at your option) any later version.

 :-)  mdrun  (-:

 Option Filename  Type Description
 
   -s  ethanol.1.tpr  InputRun input file: tpr tpb tpa
   -o  ethanol.1.trr  Output   Full precision trajectory: trr trj cpt
   -x  ethanol.1.xtc  Output, Opt. Compressed trajectory (portable xdr format)
 -cpi  ethanol.1.cpt  Input, Opt.  Checkpoint file
 -cpo  ethanol.1.cpt  Output, Opt. Checkpoint file
   -c  ethanol.1.gro  Output   Structure file: gro g96 pdb etc.
   -e  ethanol.1.edr  Output   Energy file
   -g  ethanol.1.log  Output   Log file
 -dhdl ethanol.1.dhdl.xvg  Output, Opt! xvgr/xmgr file
 -field  ethanol.1.xvg  Output, Opt. xvgr/xmgr file
 -table  ethanol.1.xvg  Input, Opt.  xvgr/xmgr file
 -tabletf  ethanol.1.xvg  Input, Opt.  xvgr/xmgr file
 -tablep  ethanol.1.xvg  Input, Opt.  xvgr/xmgr file
 -tableb  ethanol.1.xvg  Input, Opt.  xvgr/xmgr file
 -rerun  ethanol.1.xtc  Input, Opt.  Trajectory: xtc trr trj gro g96 pdb cpt
 -tpi  ethanol.1.xvg  Output, Opt. xvgr/xmgr file
 -tpid ethanol.1.xvg  Output, Opt. xvgr/xmgr file
  -ei  ethanol.1.edi  Input, Opt.  ED sampling input
  -eo  ethanol.1.xvg  Output, Opt. xvgr/xmgr file
   -j  ethanol.1.gct  Input, Opt.  General coupling stuff
  -jo  ethanol.1.gct  Output, Opt. General coupling stuff
 -ffout  ethanol.1.xvg  Output, Opt. xvgr/xmgr file
 -devout  ethanol.1.xvg  Output, Opt. xvgr/xmgr file
 -runav  ethanol.1.xvg  Output, Opt. xvgr/xmgr file
  -px  ethanol.1.xvg  Output, Opt. xvgr/xmgr file
  -pf  ethanol.1.xvg  Output, Opt. xvgr/xmgr file
  -ro  ethanol.1.xvg  Output, Opt. xvgr/xmgr file
  -ra  ethanol.1.log  Output, Opt. Log file
  -rs  ethanol.1.log  Output, Opt. Log file
  -rt  ethanol.1.log  Output, Opt. Log file
 -mtx  ethanol.1.mtx  Output, Opt. Hessian matrix
  -dn  ethanol.1.ndx  Output, Opt. Index file
 -multidir ethanol.1  Input, Opt., Mult. Run directory
 -membed  ethanol.1.dat  Input, Opt.  Generic data file
  -mp  ethanol.1.top  Input, Opt.  Topology file
  -mn  ethanol.1.ndx  Input, Opt.  Index file

 Option   Type   Value   Description
 --
 -[no]h   bool   no  Print help info and quit
 -[no]version bool   no  Print version info and quit
 -niceint0   Set the nicelevel
 -deffnm  string ethanol.1  Set the default filename for all file options
 -xvg enum   xmgrace  xvg plot formatting: xmgrace, xmgr or none
 -[no]pd  bool   no  Use particle decompostion
 -dd  vector 0 0 0   Domain decomposition grid, 0 is optimize
 -ddorder

Re: [gmx-users] segfault when running the alchemistry.org ethanol solvation free energy tutorial

2013-09-26 Thread Michael Shirts
I thought I had just managed to solve the issue :)

If you look at the soft core parameters, there are two types listed --
one with sc-r-power = 48, and one with sc-r-power = 6.  The sc-r-power
are more stable with single precision calculations.

I have changed the files on the website to make the single precision
ones the default.  Expanded.mdp warned about the issues with
precision, but left the sc-r-power in place; the Ethanol.mdp did not
warn about the potential issue.  Now both include the more efficient
path commented out.

NOTE that this means the exact numbers of the tutorials are not quite
right anymore; the process is unchanged, as is the final answer, but
the intermediate dG's will be slightly different.

However, I need to redo them anyway to make them easier in the next
1-2 weeks, so I will update them then.

If it still fails with double, that's a different issue -- because I'm
running them fine with double.





On Thu, Sep 26, 2013 at 6:32 PM, Christopher Neale
chris.ne...@mail.utoronto.ca wrote:
 My mistake  I still get a segfault even when using double precision. (EM 
 doesn't help, nor does switching to Berendsen pressure coupling).

 Note that I can stop the segfault when running at init-lambda-state = 0 if I 
 set:

 couple-lambda0   = none
 couple-lambda1   = none

 instead of

 couple-lambda0   = vdw-q
 couple-lambda1   = none

 This looks like http://bugzilla.gromacs.org/issues/1306

 Thank you,
 Chris.

 -- original message --

 Thank you Mark.

 I actually found that it crashed wihout the -multi part (no hamiltonian 
 exchange).
 The command that I used was: mdrun -nt 1 -deffnm ethanol.1 -dhdl 
 ethanol.1.dhdl.xvg

 If I use the double precision version, there is no segfault. That's a working 
 solution, but it is worrysome.

 Thank you,
 Chris.

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] segfault when running the alchemistry.org ethanol solvation free energy tutorial

2013-09-26 Thread Michael Shirts
And with this change, single is running fine as well.

This was a known issue, but was only documented in the expanded.mdp
files, which was an oversight.  After this, I switched so the default
is less likely to cause problems.  Because of some theory improvements
developed in the group in free energy calculation pathways, the
sc-r-power=48 pathway will now be phased out anyway by 5.1.

On Thu, Sep 26, 2013 at 6:37 PM, Michael Shirts mrshi...@gmail.com wrote:
 I thought I had just managed to solve the issue :)

 If you look at the soft core parameters, there are two types listed --
 one with sc-r-power = 48, and one with sc-r-power = 6.  The sc-r-power
 are more stable with single precision calculations.

 I have changed the files on the website to make the single precision
 ones the default.  Expanded.mdp warned about the issues with
 precision, but left the sc-r-power in place; the Ethanol.mdp did not
 warn about the potential issue.  Now both include the more efficient
 path commented out.

 NOTE that this means the exact numbers of the tutorials are not quite
 right anymore; the process is unchanged, as is the final answer, but
 the intermediate dG's will be slightly different.

 However, I need to redo them anyway to make them easier in the next
 1-2 weeks, so I will update them then.

 If it still fails with double, that's a different issue -- because I'm
 running them fine with double.





 On Thu, Sep 26, 2013 at 6:32 PM, Christopher Neale
 chris.ne...@mail.utoronto.ca wrote:
 My mistake  I still get a segfault even when using double precision. (EM 
 doesn't help, nor does switching to Berendsen pressure coupling).

 Note that I can stop the segfault when running at init-lambda-state = 0 if I 
 set:

 couple-lambda0   = none
 couple-lambda1   = none

 instead of

 couple-lambda0   = vdw-q
 couple-lambda1   = none

 This looks like http://bugzilla.gromacs.org/issues/1306

 Thank you,
 Chris.

 -- original message --

 Thank you Mark.

 I actually found that it crashed wihout the -multi part (no hamiltonian 
 exchange).
 The command that I used was: mdrun -nt 1 -deffnm ethanol.1 -dhdl 
 ethanol.1.dhdl.xvg

 If I use the double precision version, there is no segfault. That's a 
 working solution, but it is worrysome.

 Thank you,
 Chris.

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] segfault when running the alchemistry.org ethanol solvation free energy tutorial

2013-09-26 Thread Michael Shirts
Sounds like we've resolved all the confusion.  Thanks for prompt help
in making this clearer and better.

On Thu, Sep 26, 2013 at 7:01 PM, Christopher Neale
chris.ne...@mail.utoronto.ca wrote:
 Agreed, the following parameters do not segfault in single or double 
 precision:
 sc-alpha = 0.5
 sc-power = 1
 sc-r-power   = 6
 Same goes for http://bugzilla.gromacs.org/issues/1306

 The following parameters give a segfault in single precision but are ok in 
 double precision
 sc-alpha = 0.001
 sc-power = 1
 sc-r-power   = 48
 Same goes for http://bugzilla.gromacs.org/issues/1306

 Sorry for the confusion earlier. I was using a compilation that I thought was 
 double precision
 but it was actually single. Recompiling in double precision gave me the 
 stability outlined above.

 Thank you for your assitance Mark and Michael.

 Chris.

 -- original message --

 Michael Shirts mrshirts at gmail.com
 Fri Sep 27 00:41:17 CEST 2013
 Previous message: [gmx-users] segfault when running the alchemistry.org 
 ethanol solvation free energy tutorial
 Next message: [gmx-users] segfault when running the alchemistry.org ethanol 
 solvation free energy tutorial
 Messages sorted by: [ date ] [ thread ] [ subject ] [ author ]
 And with this change, single is running fine as well.

 This was a known issue, but was only documented in the expanded.mdp
 files, which was an oversight.  After this, I switched so the default
 is less likely to cause problems.  Because of some theory improvements
 developed in the group in free energy calculation pathways, the
 sc-r-power=48 pathway will now be phased out anyway by 5.1.

 On Thu, Sep 26, 2013 at 6:37 PM, Michael Shirts mrshirts at gmail.com wrote:
 I thought I had just managed to solve the issue :)

 If you look at the soft core parameters, there are two types listed --
 one with sc-r-power = 48, and one with sc-r-power = 6.  The sc-r-power
 are more stable with single precision calculations.

 I have changed the files on the website to make the single precision
 ones the default.  Expanded.mdp warned about the issues with
 precision, but left the sc-r-power in place; the Ethanol.mdp did not
 warn about the potential issue.  Now both include the more efficient
 path commented out.

 NOTE that this means the exact numbers of the tutorials are not quite
 right anymore; the process is unchanged, as is the final answer, but
 the intermediate dG's will be slightly different.

 However, I need to redo them anyway to make them easier in the next
 1-2 weeks, so I will update them then.

 If it still fails with double, that's a different issue -- because I'm
 running them fine with double.


 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Need protein-ligand free energy calculation tutorial

2013-09-19 Thread Michael Shirts
There are some starter files here:

http://www.gromacs.org/Documentation/Tutorials/GROMACS_USA_Workshop_and_Conference_2013/An_introduction_to_free_energy_calculations%3a_Michael_Shirts%2c_Session_2A

Which can be used in conjunction with the Alchemistry.org instructions.

But it needs to be updated a bit more . . .

On Thu, Sep 19, 2013 at 3:38 AM, hsp85 sfil...@gmail.com wrote:
 Actually,
 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/
 there are two different manuals - one about free energy calculation and
 other about protein ligand complexes.

 This page is also can be helpfull
 http://www.alchemistry.org/wiki/Category:Free_Energy_How-to%27s

 Sergey Filkin

 --
 View this message in context: 
 http://gromacs.5086.x6.nabble.com/Need-protein-ligand-free-energy-calculation-tutorial-tp5011299p5011302.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Long range Lennard Jones

2013-08-29 Thread Michael Shirts
IPS in CHARMM involves additional calculations beyond a simple
homogeneous approximation -- roughly equivalent to PME for dispersion,
though its a bit messier.

http://www.ncbi.nlm.nih.gov/pmc/articles/PMC2723858/



On Thu, Aug 29, 2013 at 7:56 AM, Justin Lemkul jalem...@vt.edu wrote:


 On 8/29/13 1:18 AM, Gianluca Interlandi wrote:

 Justin,

 I respect your opinion on this. However, in the paper indicated below by
 BR
 Brooks they used a cutoff of 10 A on LJ when testing IPS in CHARMM:

 Title: Pressure-based long-range correction for Lennard-Jones interactions
 in
 molecular dynamics simulations: Application to alkanes and interfaces
 Author(s): Lague, P; Pastor, RW; Brooks, BR
 Source: JOURNAL OF PHYSICAL CHEMISTRY B Volume: 108 Issue: 1 Pages:
 363-368
 Published: JAN 8 2004


 I cannot say for sure whether the DispCorr implementation in Gromacs is the
 same or comparable to IPS in CHARMM.  You would have to test that, and also
 test a more complex system than simple liquids.


 There is also a paper by Piana and Shaw where different cutoffs for
 non-bonded
 are tested with CHARMM22 on Anton:

 http://www.plosone.org/article/info%3Adoi%2F10.1371%2Fjournal.pone.0039918

 They found some subtle differences, in particular for cutoffs shorter than
 9 A.
 However, Anton uses abrupt truncation (no switching) and I believe that
 the
 differences they found at cutoffs  9 A would be much smaller if they had
 used a
 finer mesh (as they show at the 8 A cutoff). I always use
 fourierspacing=1.0


 Note that Shaw's group is using their own modified force field, CHARMM22*,
 so that factors in here, as well.  I do recall that paper, though I haven't
 read it in a while, so details may be fuzzy.  Wasn't the principal point to
 test methods for long-range electrostatics and the influence of cutoffs in
 that context?  I seem to recall LJ taking a backseat there.  It is certainly
 true that short-range Coulomb cutoffs are more flexible when using Ewald
 methods.


 I agree though that it strongly depends on the system and I have always
 run
 control simulations but never found significant differences in the case of
 just
 proteins.


 That's good to know that you've tested different setups.  Are the water
 interaction energies comparable?  Have you tested model compounds or just
 full proteins?  Over how long of a time period?


 Finally, I have not tested it in gromacs, but in NAMD there is a
 performance
 gain of 25% when using the shorter cutoff. This is a factor to consider.
 When I
 asked for Teragrid supercomputing allocations back in 2006 and 2007 and I
 suggested 10/12/14 cutoff, the reviewers always complained and cut my
 requested
 time of 20% with the justification that I must use a shorter cutoff.


 That's unfortunate.  In my opinion, it is shameful that desired performance
 exceeds proven integrity of the data.  If all we're after is performance, I
 can write up a dozen papers in the next few months with completely useless
 data, but they will be produced quickly! ;)


 -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

 jalem...@outerbanks.umaryland.edu | (410) 706-7441

 ==
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] A question about velocity rescaling thermostat and velocity verlet integrator

2013-08-29 Thread Michael Shirts
The only integrators with stochastic force components are sd and bd.

vrescale has a small stochastic component, but that is for the target
kinetic energy, and is not a random force acting on each particle.

On Thu, Aug 29, 2013 at 3:15 PM, Ali Sinan Saglam
asinansag...@gmail.com wrote:
 Hi,

 I was planning to use the velocity verlet ingetrator (md-vv) with velocity
 rescaling thermostat (tcoupl = v-rescale). I know that velocity rescaling
 thermostat is a stochastic thermostat and uses a random seed and generally
 with the stochastic thermostats one would use the sd integrator that uses
 ld-seed = -1 option for a random seed.

 I was just wondering if one uses the velocity verlet integrator with a
 stochastic thermostat does the ld-seed option work as it does with the sd
  (or bd) integrator?  (the manual currently doesn't say anything about this
 http://manual.gromacs.org/online/mdp_opt.html#run)

 Best,
 Ali Sinan Saglam
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] NPT-REMD

2013-08-25 Thread Michael Shirts
Pressure should fluctuate significantly.  The estimator for the
pressure that is generally used is very noisy. The question is, do the
pressure averages over, say, 500 ps or 1 ns look about right?

On Sun, Aug 25, 2013 at 5:07 PM, Mark Abraham mark.j.abra...@gmail.com wrote:
 On Sun, Aug 25, 2013 at 6:22 PM, rahul seth
 rahul.seth.grom...@gmail.com wrote:
 Hi All,

 I have been trying to perform NPT-REMD with a Protein Substrate and water.
 I am trying to use md+nose-hoover thermostat and parrinello-rahman
 barostat. However, I am not quite sure about the tau-t and the tau-p that I
 should be using here. I paste a part of my mdp file below:
 define   =  -DPOSRES_subs
 constraints=  all-bonds
 pbc   =  xyz
 integrator   =  md
 ;ld_seed= %8i
 dt =  0.002; ps !
 nsteps   =  25   ; 50 ns
 ;
 nstcomm =  1
 nstcalcenergy   =  10
 nstxout   =  0
 nstxtcout=  5000 ;every 10 ps
 nstvout   =  0
 nstfout   =  0
 nstlog=  1000
 nstenergy   =  1000
 nstlist=  10
 ns_type =  grid
 rlist   =  1.0
 vdwtype =  cut-off
 rvdw  =  1.0
 coulombtype =  pme
 rcoulomb=  1.0
 fourierspacing  = 0.12
 pme_order   = 4
 ewald_rtol   = 1e-5
 optimize_fft =  yes
 tcoupl= nose-hoover
 nsttcouple  = 5
 tc-grps   =  subs Protein SOL
 tau_t  = 0.2 0.2 0.2
 ref_t   =  %8.3f %8.3f %8.3f
 ; Energy monitoring
 energygrps  =  subs Protein SOL
 ; Isotropic pressure coupling is now on
 Pcoupl   = parrinello-rahman
 nstpcouple  = 1
 refcoord-scaling   = all
 Pcoupltype  = isotropic
 tau_p  = 2.0
 compressibility = 4.5e-5
 ref_p= 1.01325


 Although the temperatures in this case stays close to the desired values,
 the pressures of the individual replicas fluctuate significantly and hence
 I do not think I have weird replica exchange probabilities.

 Seemingly large
 http://www.gromacs.org/Documentation/Terminology/Pressure fluctuations
 are normal - see link. You should satisfy yourself that your .mdp file
 produces ensembles that are respectable before you get involved with
 REMD, and gathering statistics for that will take longer for pressure
 than temperature.

 I have varied tau-t and tau-p from 0.2 to 1 and 1 to 5 respectively. It
 doesn't really seem to improve anything. Can anyone suggest a set of
 optimized values for these parameters? I understand these maybe context
 dependent but to give you a rough idea I have about ~18,000 atoms with 170
 for atoms 1156 for substrate and the rest for water.

 Your literature survey of similar simulations should give you an idea
 of usage in the field, but (e.g. per GROMACS manual thermostat
 section) I would think tau-t for N-H should not go below 1. Otherwise,
 first observe whether you have a real problem before making haphazard
 changes :-)

 Mark

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: NPT-REMD

2013-08-25 Thread Michael Shirts
Can you clarify - Do you mean that different replicas have different
average pressures?

WITHIN each replica, the +/- 2000 bar changing from step to step is
very common for using an atomic virial like gromacs does.

The AVERAGES of EACH replica should each be the average pressure they
are set as (+/- much less than 2000 -- maybe 10-20 bar).

You could also try md-vv + MTTK + Nose-Hoover.  I have had decent luck
with this combination and tested it for a while.

You should also try running without any exchange.  It could be that
the problem is not replica exchange, but something else.



On Sun, Aug 25, 2013 at 8:36 PM, rahul.seth.gromacs
rahul.seth.grom...@gmail.com wrote:
 Thanks everyone for your suggestions. Pressures for some of my replicas
 fluctuate between -2000 bars to +2000 bars. I do not think that can be
 termed as about right isn't it? The averages I am talking are over 1ns.



 --
 View this message in context: 
 http://gromacs.5086.x6.nabble.com/NPT-REMD-tp5010721p5010728.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Early registration period ending for 2013 GROMACS USA Workshop and Conference

2013-08-20 Thread Michael Shirts
Hi,  Rajat-

We are working on ways to make the workshop materials available
online, either streaming or posted afterwards.  The tutorials will
have online components as well. We were unable to obtain travel
funding for the conference this year, unfortunately!

 I would be willing to pay the whole workshop fee (and more) if you could 
 kindly video the whole workshop and make it available to us. Thank you.

We are not planning on charging for the material, but we should make
it easier to donate if you wanted to after we get this material out!


On Tue, Aug 20, 2013 at 12:14 AM, rajat desikan rajatdesi...@gmail.com wrote:
 Hi Prof. Shirts,
 I am sure that I speak for a lot of GROMACS users when I say that I want to
 attend the workshop, but I cannot. I would be willing to pay the whole
 workshop fee (and more) if you could kindly video the whole workshop and
 make it available to us. Thank you.


 On Tue, Aug 20, 2013 at 10:50 AM, Michael Shirts 
 michael.shi...@virginia.edu wrote:

 Dear GROMACS users-

 I'd like to remind you all about the 2013 GROMACS USA Workshop and
 Conference at the University of Virginia Sept 13th-15th.  There are
 still registration slots available, but the early registration
 deadline of Aug 22nd is coming up in just a few days; after the 22nd,
 the registration price will rise from the very low $60 for the two day
 conference to the moderately-low-but-why-pay-more $95.

 As mentioned in previous emails to the list, the workshop and
 conference will consist of two days of tutorials, discussions of
 future plans for GROMACS, face time with many of the main GROMACS
 developers, plenary software and application sessions, and an optional
 last day of development planning to which attendees are also invited
 to help influence the future directions of GROMACS.

 Please visit http://faculty.virginia.edu/gromacsworkshop for
 registration and for much more information about the workshop,
 including travel logistics. The website was also recently has been
 updated with more specifics about the program.

 For specific questions about registering or logistics after visiting
 the website, please write to
 gromacsworkshop-registrat...@virginia.edu.

 Sincerely,
 The 2013 GROMACS USA Workshop and Conference Steering Committee
 Michael Shirts (chair)
 Angel Garcia
 Berk Hess
 Yu-Shan Lin
 Erik Lindahl
 Peter Kasson
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists




 --
 Rajat Desikan (Ph.D Scholar)
 Prof. K. Ganapathy Ayappa's Lab (no 13),
 Dept. of Chemical Engineering,
 Indian Institute of Science, Bangalore
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


[gmx-users] Early registration period ending for 2013 GROMACS USA Workshop and Conference

2013-08-19 Thread Michael Shirts
Dear GROMACS users-

I'd like to remind you all about the 2013 GROMACS USA Workshop and
Conference at the University of Virginia Sept 13th-15th.  There are
still registration slots available, but the early registration
deadline of Aug 22nd is coming up in just a few days; after the 22nd,
the registration price will rise from the very low $60 for the two day
conference to the moderately-low-but-why-pay-more $95.

As mentioned in previous emails to the list, the workshop and
conference will consist of two days of tutorials, discussions of
future plans for GROMACS, face time with many of the main GROMACS
developers, plenary software and application sessions, and an optional
last day of development planning to which attendees are also invited
to help influence the future directions of GROMACS.

Please visit http://faculty.virginia.edu/gromacsworkshop for
registration and for much more information about the workshop,
including travel logistics. The website was also recently has been
updated with more specifics about the program.

For specific questions about registering or logistics after visiting
the website, please write to
gromacsworkshop-registrat...@virginia.edu.

Sincerely,
The 2013 GROMACS USA Workshop and Conference Steering Committee
Michael Shirts (chair)
Angel Garcia
Berk Hess
Yu-Shan Lin
Erik Lindahl
Peter Kasson
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] LINCS Constraints - all-bonds or h-bonds?

2013-08-15 Thread Michael Shirts
I don't go beyond 2 fs with either all- bonds or h-bonds. Things like kinetic 
energy start being subtly off.

H-bonds has less chance of failing with large numbers of constraints- less 
iteration required, especially if bond system cross parallelization boundaries.

If your molecules are  10 atoms, it probably doesn't matter either way.

Sent from my iPhone

On Aug 15, 2013, at 9:11, Barnett, James W. jbarn...@tulane.edu wrote:

 Searching through this mailing list it seems like some have stated that 
 with a 2 fs timestep (dt=0.002), constraints=h-bonds is fine in general.
 
 The questions I have are:
 
 1) What are some personal opinions on when it is ok to switch to h-bonds 
   from all-bonds for LINCS constraints? Is 2 fs and h-bonds a general 
 practice?
 
 2) Also, if typically 2 fs and h-bonds are ok, what time-step do users 
   (or you personally) generally go to with all-bonds? 
 
 I am speaking generally here of course. Thanks for your responses.
 
 -- 
 Wes Barnett | jbarn...@tulane.edu
 
 -- 
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Unkwown Keyword HILLS

2013-08-14 Thread Michael Shirts
This is

Sent from my iPhone

On Aug 14, 2013, at 1:40, Albert mailmd2...@gmail.com wrote:

 Does anybody have any idea what's the problem?
 
 I use the tutorial example and I don't know why it doesn't work.
 
 THX
 
 
 On 08/13/2013 07:19 PM, Albert wrote:
 Dear:
 
 I am trying to run plumed with gromacs plugin. Here is my plumed.dat file 
 which I defined two dihedral angels as cvs:
 
 *HILLS HEIGHT 0.3 W_STRIDE 450
 WELLTEMPERED SIMTEMP 310 BIASFACTOR 1.96
 TORSION LIST 1 4 65 344 SIGMA 0.12
 TORSION LIST 2 46 80 656 SIGMA 0.12
 
 ENDMETA*
 
 I am using plumed-1.3+gromacs-4.6.2 with command:
 
 
 *mpirun -np 24 mdrun_mpi -s md.tpr -plumed plumed.dat -g md.log -v -x md.xtc 
 -o md.trr -e md.edr*
 
 
 but it failed with messages:
 
 
 *starting mdrun 'protein'
 1 steps, 20.0 ps.
 ! PLUMED ERROR: Line 1 Unkwown Keyword HILLS
 
 ! ABORTING RUN
 ! PLUMED ERROR: Line 1 Unkwown Keyword HILLS
 
 ! ABORTING RUN
 --*
 
 thank you very much
 
 Albert
 
 -- 
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Unkwown Keyword HILLS

2013-08-14 Thread Michael Shirts
This is a plumed error, not a gromacs error. Gromacs never handles those 
keywords. 

Sent from my iPhone

On Aug 14, 2013, at 1:40, Albert mailmd2...@gmail.com wrote:

 Does anybody have any idea what's the problem?
 
 I use the tutorial example and I don't know why it doesn't work.
 
 THX
 
 
 On 08/13/2013 07:19 PM, Albert wrote:
 Dear:
 
 I am trying to run plumed with gromacs plugin. Here is my plumed.dat file 
 which I defined two dihedral angels as cvs:
 
 *HILLS HEIGHT 0.3 W_STRIDE 450
 WELLTEMPERED SIMTEMP 310 BIASFACTOR 1.96
 TORSION LIST 1 4 65 344 SIGMA 0.12
 TORSION LIST 2 46 80 656 SIGMA 0.12
 
 ENDMETA*
 
 I am using plumed-1.3+gromacs-4.6.2 with command:
 
 
 *mpirun -np 24 mdrun_mpi -s md.tpr -plumed plumed.dat -g md.log -v -x md.xtc 
 -o md.trr -e md.edr*
 
 
 but it failed with messages:
 
 
 *starting mdrun 'protein'
 1 steps, 20.0 ps.
 ! PLUMED ERROR: Line 1 Unkwown Keyword HILLS
 
 ! ABORTING RUN
 ! PLUMED ERROR: Line 1 Unkwown Keyword HILLS
 
 ! ABORTING RUN
 --*
 
 thank you very much
 
 Albert
 
 -- 
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Hamiltonian replica exchange not working in 4.6

2013-08-06 Thread Michael Shirts
Hi, Sanku-

The way to invoke Hamiltonian replica exchange has changed to be a bit
more flexible.  We should go back and make sure that this legacy way
is supported (I thought this invocation was supported, but apparently
it isn't), but what you should be able to do to get it working quickly
is include in all the .mdp files the line:

fep-lambdas  0 0.08 . . . . . . 1.0

I.e. list all of the lambdas that all of the files use.  This is
described in the documentation.

It should run correctly then.  If not, then bug me on the list again.


On Tue, Aug 6, 2013 at 6:52 PM, Sanku M msank...@yahoo.com wrote:
 Dear Gromacs user,
   I am not sure whether this is any potential bug. But I found that in 
 gromacs4.6.3, hamiltonian replica exchange simulation is not
 working. The SAME simulation works perfectly fine in gromacs 4.5.4
 But, when I try to run the simulation using gromacs4.6.3, I get following 
 error:
 The properties of the 10 systems are all the same, there is nothing to 
 exchange
 For more information and tips for troubleshooting, please check the GROMACS


 Here is the details on the error:
I have first defined A and B state by scaling the relevant parameter.. 
 Then I have  generated different tpr file by using mdp file which only differ 
 by 'init-lambda' parameter:
 For example, I have 10 different mdp files: Free energy part of two of them 
 are shown below:
 for hremd0.mdp, the relevant part for free energy part is
   ; Free energy control stuff
 free_energy  = yes
 init-lambda  = 0
 delta_lambda = 0

 for hremd1.mdp
 ; Free energy control stuff
 free_energy  = yes
 init-lambda  = .080
 delta_lambda = 0
 ..
 and so on..
 The rest of mdp portions are same for all the mdp files...



 However, with gromacs4.5.4, the same simulation works without problem. I do 
 not get the error that all the systems are identical.

 In both cases, I am using same system with OPLS forcefield..
 The way I run the simulations is:
 mpirun -np 10 mdrun_463 -multi 10 -v -s hremd -replex 250 -c log_repl


 Any help will be appreciated.

 Sanku

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Unphysical conformations in decoupled free energy simulation

2013-08-05 Thread Michael Shirts
Hi, all-

That particular redmine was in the case of couple-intramol NOT being
specified.  So it may be different.

Can you upload samples files for this test case to redmine.gromacs.org
with specific instructions?  Please upload both the 'behaving
correctly' and 'not behaving correctly' versions.

Note that in may cases, you can get the correct behavior by changing
the topologies themselves, though this is more complicated.

On Mon, Aug 5, 2013 at 10:25 AM, Justin Lemkul jalem...@vt.edu wrote:


 On 8/5/13 9:52 AM, Joerg Sauter wrote:

 On 08/05/2013 02:15 PM, Justin Lemkul wrote:



 On 8/5/13 8:12 AM, Joerg Sauter wrote:

 On 08/05/2013 01:54 PM, Justin Lemkul wrote:



 On 8/5/13 7:44 AM, Joerg Sauter wrote:

 Dear all,

 I am trying to compute the free energy of hydration for cellobiose (a
 beta
 (1-4)
 glucose dimer) using BAR in Gromacs 4.6.1. However, I encounter a
 problem.

 I find that the vacuum conformations of the molecule in a regular
 vacuum
 simulation differ from the conformations in the decoupled simulation
 in the
 free
 energy case i.e., with the additional mdp entries:

 free-energy  = yes
 init-lambda = 0
 couple-lambda0   = none
 couple-moltype   = solute
 couple-intramol  = no

 Here is a histogram of the dihedral angles of the glycosidic linkage
 in
 the vacuum case
 https://dl.dropboxusercontent.com/u/70358077/reg.pdf (it stays in the
 global
 free energy minimum)
 and this is the decoupled free energy case (same starting conformation
 in the
 global free energy minimum)
 https://dl.dropboxusercontent.com/u/70358077/fe.pdf

 I understand that Gromacs replaces the intramolecular interactions
 with
 explicit
 pair
 interactions. Therefore, I had to increase table-extension but that
 did not
 change much. I was thinking that maybe this could be a problem
 specific to
 this
 topology (a GLYCAM06h conversion from Amber), however, I do not see
 how
 this can
 occur.
 I hope someone has an idea what is going wrong.

 mdp file: https://dl.dropboxusercontent.com/u/70358077/md.mdp (same
 behaviour
 for longer runs)
 top file: https://dl.dropboxusercontent.com/u/70358077/cellob.top
 gro file: https://dl.dropboxusercontent.com/u/70358077/cellob.gro

 This is sort of a minimal example. The same problem occurs
 in a regular simulation when decoupling from water using multiple
 lambdas.


 Have you tried running with rlist=rcoulomb=rvdw = 0 rather than 99? If
 nothing
 else, I would try again with version 4.6.3.  There have been tons of
 fixes to
 the free energy code since 4.6.1, though I don't know offhand which, if
 any,
 would be affecting your results. Regardless, the current version should
 always
 be the best version ;)

 -Justin

 Hi Justin,

 I just tried it with version 4.6.3. and rlist=rcoulomb=rvdw = 0.
 Unfortunately,
 I get the same results.


 OK, good to know.  We can rule out 4.6.1-era bugs, then.  Does
 visualization
 or an analysis of energy terms stored in the .edr file provide any clues
 as to
 the source of the problem?

 -Justin

 Visualization does not show any abnormalities beyond odd confomations in
 the
 free energy coupling case.

 Regarding the energies: I am not sure whether this is normal behaviour,
 but, in
 the free energy coupling case I do not see energy values for LJ 14 or
 Coulomb 14
 interactions. I thought that by selecting couple-intramol=0 all
 intramolecular
 interactions would be replaced by pairtype (14) interactions.

 In addition, for Coulomb (SR) I get in the free energy case
 LJ (SR)-22.8084  0.0765.90879 -0.379238
 (kJ/mol)
 in comparisson to
 LJ (SR) 16.8295   0.8113.2303 4.51371
 (kJ/mol)
 in the regular vacuum simulation. So this looks like the 14 interactions
 get
 shifted to the LJ/Coulomb terms in the coupling simulation?! (However, the
 numbers do not add up)


 This is smelling buggy to me, and there have been issues with 1-4
 interactions in the free energy code before
 (http://redmine.gromacs.org/issues/1225), so there could be a lingering
 problem.  Hopefully Michael will chime in.  It might be worthwhile to open
 up a bug report on redmine.gromacs.org with all of the information you have
 provided (input files, output data, and the very nice description of the
 issue you have posted).


 -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

 jalem...@outerbanks.umaryland.edu | (410) 706-7441

 ==
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
 * Please don't post (un)subscribe requests to the list. 

Re: [gmx-users] restraint-lambdas for position restraints in hamiltonian exchange

2013-08-03 Thread Michael Shirts
Short answer is anything that has a B state parameter can be included
in in Hamiltonian exchange.

If it's pull code or explicit restraints, it's controlled by restraint lambda.

 I went through the manual and couldn't find any definite answers to the
 following questions.

 First, I wonder if the reference positions of position restraints, not just
 the force constants, of different replicas are exchanged in hamiltonian
 exchange based on restraint-lambdas? For example, if I have 1 molecule that
 has 2 structures, say, A and B and the following 2-component
 restraint-lambdas:  1.0 1.0 with replica 0 starting in A and replica 1 in B,
 would the exchange between replica 0 and 1 yield anything meaningful ?

Right now, the pull code doesn't have B state positions, just force
constants. It CAN be done using distance restraints, but then you have
to make the atoms involved part of the same molecule.  Adding B state
distances to pull code was recently requested, and will probably be
straightforward to put in soon.

 In
 other word, would the relaxation towards the other structure occur in both
 states once an exchange was accepted?

 Another question is what types of restraints can restraint-lambdas act on?
 For example, bond/angle/distance/position restraints, COM-pulling potential
 ?

Anything with a B state in gromacs topologies.
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: restraint-lambdas for position restraints in hamiltonian exchange

2013-08-03 Thread Michael Shirts
That is correct. Such a functionality wouldn't be that hard to
implement - but there are a long list of easy functionalities to be
implemented. You can submit a request to redmine.gromacs.org so that
the request is archived.

On Sat, Aug 3, 2013 at 6:12 PM, Dejun Lin dejun@gmail.com wrote:
 So I take it that in the position restraint case (not COM-pulling), where
 the reference positions are determined by the starting structure instead of
 a B-state topology, the reference positions won't be swapped ?


 2013/8/3 Michael Shirts-2 [via GROMACS] 
 ml-node+s5086n5010324...@n6.nabble.com

 Short answer is anything that has a B state parameter can be included
 in in Hamiltonian exchange.

 If it's pull code or explicit restraints, it's controlled by restraint
 lambda.

  I went through the manual and couldn't find any definite answers to the
  following questions.
 
  First, I wonder if the reference positions of position restraints, not
 just
  the force constants, of different replicas are exchanged in hamiltonian
  exchange based on restraint-lambdas? For example, if I have 1 molecule
 that
  has 2 structures, say, A and B and the following 2-component
  restraint-lambdas:  1.0 1.0 with replica 0 starting in A and replica 1
 in B,
  would the exchange between replica 0 and 1 yield anything meaningful ?

 Right now, the pull code doesn't have B state positions, just force
 constants. It CAN be done using distance restraints, but then you have
 to make the atoms involved part of the same molecule.  Adding B state
 distances to pull code was recently requested, and will probably be
 straightforward to put in soon.

  In
  other word, would the relaxation towards the other structure occur in
 both
  states once an exchange was accepted?
 
  Another question is what types of restraints can restraint-lambdas act
 on?
  For example, bond/angle/distance/position restraints, COM-pulling
 potential
  ?

 Anything with a B state in gromacs topologies.
 --
 gmx-users mailing list[hidden 
 email]http://user/SendEmail.jtp?type=nodenode=5010324i=0
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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 [hidden 
 email]http://user/SendEmail.jtp?type=nodenode=5010324i=1.

 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


 --
  If you reply to this email, your message will be added to the discussion
 below:

 http://gromacs.5086.x6.nabble.com/restraint-lambdas-for-position-restraints-in-hamiltonian-exchange-tp5010322p5010324.html
  To unsubscribe from restraint-lambdas for position restraints in
 hamiltonian exchange, click 
 herehttp://gromacs.5086.x6.nabble.com/template/NamlServlet.jtp?macro=unsubscribe_by_codenode=5010322code=ZGVqdW4ubGluQGdtYWlsLmNvbXw1MDEwMzIyfDE0MzAxNDIxNA==
 .
 NAMLhttp://gromacs.5086.x6.nabble.com/template/NamlServlet.jtp?macro=macro_viewerid=instant_html%21nabble%3Aemail.namlbase=nabble.naml.namespaces.BasicNamespace-nabble.view.web.template.NabbleNamespace-nabble.view.web.template.NodeNamespacebreadcrumbs=notify_subscribers%21nabble%3Aemail.naml-instant_emails%21nabble%3Aemail.naml-send_instant_email%21nabble%3Aemail.naml





 --
 View this message in context: 
 http://gromacs.5086.x6.nabble.com/restraint-lambdas-for-position-restraints-in-hamiltonian-exchange-tp5010322p5010325.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Expanded ensemble simulation died with fatal error: Something wrong in choosing new lambda state with a Gibbs move

2013-07-31 Thread Michael Shirts
Hi Dejun-

The basic problem is that for this particular configuration, the
current state is the only state with nonzero weight. Note that the
state with the second highest weight has weight 10^-7.  When it tries
to compare weights in single precision, it has a numerical overflow
and fails.

A few things:

1. This really should be more robust, so that it will realize it's
supposed to stay in the most likely state, since that's the only state
with nonzero weight.  I have a fix that I've been working on for
exactly this problem, but it's not quite ready yet.  Hopefully in the
next couple of days.

2. This problem is very unlikely to occur in double precision, if you
can afford the performance hit in the meantime.

3. If this is a typical average difference, exchanges will be very
unlikely.  You should probably choose your lambda intervals to be a
bit closer together at the end range.

Hopefully this will give you enough information to move forward for
the time being until a better fix is implemented.


On Wed, Jul 31, 2013 at 2:15 PM, Dejun Lin dejun@gmail.com wrote:
 Hi all,

 I'm running an expanded ensemble simulation using gromacs 4.6.3 and it
 crashed with the error:

 Fatal error:
 Something wrong in choosing new lambda state with a Gibbs move -- probably
 underflow in weight determination.
 Denominator is:   0 1.002384e+00
   idEnumerator  weights
   0 -9.1451739502e+02 0.00e+00 0.00e+00
   1 -9.128174e+02 0.00e+00 0.00e+00
   2 -8.8548516846e+02 0.00e+00 0.00e+00
   3 -8.7096899414e+02 0.00e+00 0.00e+00
   4 -8.5645288086e+02 0.00e+00 0.00e+00
   5 -8.4193676758e+02 0.00e+00 0.00e+00
   6 -8.2742059326e+02 0.00e+00 0.00e+00
   7 -8.1290447998e+02 0.00e+00 0.00e+00
   8 -7.9838836670e+02 0.00e+00 0.00e+00
   9 -7.8387219238e+02 0.00e+00 0.00e+00
  10 -7.6935607910e+02 0.00e+00 0.00e+00
  11 -7.5483990479e+02 0.00e+00 0.00e+00
  12 -7.4032379150e+02 0.00e+00 0.00e+00
  13 -7.2580767822e+02 0.00e+00 0.00e+00
  14 -7.1129150391e+02 0.00e+00 0.00e+00
  15 -6.9677539062e+02 0.00e+00 0.00e+00
  16 -6.8225927734e+02 0.00e+00 0.00e+00
  17 -6.6774316406e+02 0.00e+00 0.00e+00
  18 -6.5322698975e+02 0.00e+00 0.00e+00
  19 -6.3871087646e+02 0.00e+00 0.00e+00
  20 -6.2419470215e+02 0.00e+00 0.00e+00
  21 -6.0967858887e+02 0.00e+00 0.00e+00
  22 -5.9516247559e+02 0.00e+00 0.00e+00
  23 -5.8064630127e+02 0.00e+00 0.00e+00
  24 -5.6613018799e+02 0.00e+00 0.00e+00
  25 -5.5161407471e+02 0.00e+00 0.00e+00
  26 -5.3709790039e+02 0.00e+00 0.00e+00
  27 -5.2258178711e+02 0.00e+00 0.00e+00
  28 -5.0806564331e+02 0.00e+00 0.00e+00
  29 -4.9354953003e+02 0.00e+00 0.00e+00
  30 -4.7903335571e+02 0.00e+00 0.00e+00
  31 -4.6451724243e+02 0.00e+00 0.00e+00
  32 -4.518311e+02 0.00e+00 0.00e+00
  33 -4.3548400879e+02 0.00e+00 0.00e+00
  34 -4.2096792603e+02 0.00e+00 0.00e+00
  35 -4.0645178223e+02 0.00e+00 0.00e+00
  36 -3.9193563843e+02 0.00e+00 0.00e+00
  37 -3.8107025146e+02 0.00e+00 0.00e+00
  38 -3.6290338135e+02 0.00e+00 0.00e+00
  39 -3.4838723755e+02 0.00e+00 0.00e+00
  40 -3.3387109375e+02 0.00e+00 0.00e+00
  41 -3.1935494995e+02 0.00e+00 0.00e+00
  42 -3.0483883667e+02 0.00e+00 0.00e+00
  43 -2.9032269287e+02 0.00e+00 0.00e+00
  44 -2.7580654907e+02 0.00e+00 0.00e+00
  45 -2.6129040527e+02 0.00e+00 0.00e+00
  46 -2.4677430725e+02 0.00e+00 0.00e+00
  47 -2.3225816345e+02 0.00e+00 0.00e+00
  48 -2.1774200439e+02 0.00e+00 0.00e+00
  49 -2.0322586060e+02 0.00e+00 0.00e+00
  50 -1.8970976257e+02 0.00e+00-1.00e+00
  51 -1.7419361877e+02 0.00e+00 0.00e+00
  52 -1.5967747498e+02 0.00e+00 0.00e+00
  53 -1.4516131592e+02 0.00e+00 0.00e+00
  54 -1.3064523315e+02 0.00e+00 0.00e+00
  55 -1.1612908173e+02 0.00e+00 0.00e+00
  56 -1.0161293030e+02 7.0064923216e-45 0.00e+00
  57 -8.7096786499e+01 1.4939846888e-38 0.00e+00
  58 -7.2580688477e+01 3.0102835162e-32 0.00e+00
  59 -5.8064540863e+01 6.0658294505e-26 0.00e+00
  60 -4.3548393250e+01 1.865341e-19 0.00e+00
  61 -2.9032243729e+01 2.4629560544e-13 0.00e+00
  62 -1.5516148567e+01 

[gmx-users] Reminder about US GROMACS workshop + soliciting presenters for talks and tutorials

2013-07-26 Thread Michael Shirts
Dear gmx-users list members-

I'd like to remind you all about the 2013 GROMACS Workshop and
Conference at the University of Virginia Sept 13th-15th; registration
costs will rise from $60 to $95 on Aug 22nd. The original invitation,
with link to the website, is at the end of this email.

We are issuing a last call for contributed 20 min talks, with
decisions on the schedule made next week, and would like to invite you
to apply to present your GROMACS-related research. Email
gromacsworkshop-registrat...@virginia.edu if you wish to submit a talk
for consideration; see the preliminary schedule on the conference
webpage for more information on the design of the talks.

We still also want to give a opportunity to advanced users who may be
interested in either presenting specific tutorial subjects at the
workshop or alternatively having an online tutorial hosted in the
GROMACS website. See the list of tutorial topics on the workshop
schedule, or propose your own. We want to try to leverage the energy
and experience of the GROMACS community to enrich the tutorial
experience, rather than have the same people up there the whole
weekend!

Sincerely,

The 2013 GROMACS USA Workshop and Conference Steering Committee

Michael Shirts (chair)
Angel Garcia
Berk Hess
Yu-Shan Lin
Erik Lindahl
Peter Kasson

Original announcement:

~~~

We are pleased to announce the 2013 GROMACS USA Workshop and
Conference at the University of Virginia in Charlottesville, Virginia,
on September 13-15th. This will be the first full GROMACS workshop
held here in the U.S. The workshop and conference will consist of two
days of tutorials, discussions of future plans for GROMACS, face time
with a large percentage of the developers, plenary software and
application sessions, and a last day of development planning to which
attendees are also invited to help influence the future directions of
GROMACS.


Please visit http://faculty.virginia.edu/gromacsworkshop for
registration and for much more information about the program.

For specific questions about registering or logistics after visiting
the website, please write to
gromacsworkshop-registrat...@virginia.edu.


Michael Shirts
Assistant Professor
Department of Chemical Engineering
University of Virginia
michael.shi...@virginia.edu
(434)-243-1821
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


[gmx-users] Question on Hessian and forces not being preserved on restart with Hessian calculation?

2013-07-20 Thread Michael Shirts
When I minimize a structure, I can get down to the force max being 0.01

Low-Memory BFGS Minimizer converged to Fmax  0.01 in 6839 steps
Potential Energy  = -5.12340607768673e+03
Maximum force =  6.68907856457542e-03 on atom 3029
Norm of force =  2.19978176343026e-03

kJ/nm.  However, when I try to perform a normal mode analysis, it
complains that the maximum energy is 30 kJ/mn.


Non-cutoff electrostatics used, forcing full Hessian format.
Allocating Hessian memory...

Maximum force: 3.91984e+01
Maximum force probably not small enough to ensure that you are in an
energy well. Be aware that negative eigenvalues may occur when the
resulting matrix is diagonalized.

I've restarted from the binary .trr output (-t option) , and I'm using
double precision. Any suggestions as to why the force max is not the
asme on restarting?  Only difference is going from l-bfgs to nm.

I've seen some older comments on this, but no resolution and nothing recent.
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


[gmx-users] Re: Question on Hessian and forces not being preserved on restart with Hessian calculation?

2013-07-20 Thread Michael Shirts
To follow up -- if I try to minimize again using -t, I get the same
low forces as in the minimization in the previous step.  So it appears
to be something with what do_nm is doing, not with errors in the
output structure.

On Sat, Jul 20, 2013 at 9:37 AM, Michael Shirts mrshi...@gmail.com wrote:
 When I minimize a structure, I can get down to the force max being 0.01

 Low-Memory BFGS Minimizer converged to Fmax  0.01 in 6839 steps
 Potential Energy  = -5.12340607768673e+03
 Maximum force =  6.68907856457542e-03 on atom 3029
 Norm of force =  2.19978176343026e-03

 kJ/nm.  However, when I try to perform a normal mode analysis, it
 complains that the maximum energy is 30 kJ/mn.


 Non-cutoff electrostatics used, forcing full Hessian format.
 Allocating Hessian memory...

 Maximum force: 3.91984e+01
 Maximum force probably not small enough to ensure that you are in an
 energy well. Be aware that negative eigenvalues may occur when the
 resulting matrix is diagonalized.

 I've restarted from the binary .trr output (-t option) , and I'm using
 double precision. Any suggestions as to why the force max is not the
 asme on restarting?  Only difference is going from l-bfgs to nm.

 I've seen some older comments on this, but no resolution and nothing recent.
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


[gmx-users] Re: Question on Hessian and forces not being preserved on restart with Hessian calculation?

2013-07-20 Thread Michael Shirts
Problem partly addressed.  If I run normal modes in -nt 1, then I get
the same force as after the minimization.

I'll file a redmine.

On Sat, Jul 20, 2013 at 9:43 AM, Michael Shirts mrshi...@gmail.com wrote:
 To follow up -- if I try to minimize again using -t, I get the same
 low forces as in the minimization in the previous step.  So it appears
 to be something with what do_nm is doing, not with errors in the
 output structure.

 On Sat, Jul 20, 2013 at 9:37 AM, Michael Shirts mrshi...@gmail.com wrote:
 When I minimize a structure, I can get down to the force max being 0.01

 Low-Memory BFGS Minimizer converged to Fmax  0.01 in 6839 steps
 Potential Energy  = -5.12340607768673e+03
 Maximum force =  6.68907856457542e-03 on atom 3029
 Norm of force =  2.19978176343026e-03

 kJ/nm.  However, when I try to perform a normal mode analysis, it
 complains that the maximum energy is 30 kJ/mn.


 Non-cutoff electrostatics used, forcing full Hessian format.
 Allocating Hessian memory...

 Maximum force: 3.91984e+01
 Maximum force probably not small enough to ensure that you are in an
 energy well. Be aware that negative eigenvalues may occur when the
 resulting matrix is diagonalized.

 I've restarted from the binary .trr output (-t option) , and I'm using
 double precision. Any suggestions as to why the force max is not the
 asme on restarting?  Only difference is going from l-bfgs to nm.

 I've seen some older comments on this, but no resolution and nothing recent.
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] segfault with an otherwise stable system when I turn on FEP (complete decoupling)

2013-07-18 Thread Michael Shirts
Chris, can you post a redmine on this so I can look at the files?

Also, does it crash immediately, or after a while?

On Thu, Jul 18, 2013 at 2:45 PM, Christopher Neale
chris.ne...@mail.utoronto.ca wrote:
 Dear Users:

 I have a system with water and a drug (54 total atoms; 27 heavy atoms). The 
 system is stable when I simulate it for 1 ns. However, Once I add the 
 following options to the .mdp file, the run dies after a few ps with a 
 segfault.

 free-energy = yes
 init-lambda = 1
 couple-lambda0 = vdw-q
 couple-lambda1 = none
 couple-intramol = no
 couple-moltype = drug

 I do not get any step files or any lincs warnings. If I look at the .xtc and 
 .edr files, there is no indication of something blowing up before the 
 segfault. I have also verified that the drug runs without any problem in 
 vacuum. I get the same behaviour if I remove constraints and use a timestep 
 of 0.5 fs. The segfault is reproducible with v4.6.1 and v4.6.3. I am using 
 the charmm FF, but I converted all UB angles in my drug to type-1 angles and 
 still got the segfault. I also get the segfault with particle decomposition 
 and/or while running a single thread. I am currently using the SD integrator, 
 but I get the same segfault with md and md-vv. Couple-intramol=yes doesn't 
 resolve it, neither does using separate T-coupling groups for the water and 
 drug. Neither does turning off pressure coupling.

 Here is the .mdp file that works fine, but gives me a segfault when I add the 
 free energy stuff (above):

 constraints = all-bonds
 lincs-iter =  1
 lincs-order =  6
 constraint_algorithm =  lincs
 integrator = sd
 dt = 0.002
 tinit = 0
 nsteps = 10
 nstcomm = 1
 nstxout = 0
 nstvout = 0
 nstfout = 0
 nstxtcout = 500
 nstenergy = 500
 nstlist = 10
 nstlog=0 ; reduce log file size
 ns_type = grid
 vdwtype = cut-off
 rlist = 0.8
 rvdw = 0.8
 rcoulomb = 0.8
 coulombtype = cut-off
 tc_grps =  System
 tau_t   =  1.0
 ld_seed =  -1
 ref_t = 310
 gen_temp = 310
 gen_vel = yes
 unconstrained_start = no
 gen_seed = -1
 Pcoupl = berendsen
 pcoupltype = isotropic
 tau_p = 4
 compressibility = 4.5e-5
 ref_p = 1.0

 I do realize that some of these settings are not ideal for a production run. 
 I started with the real Charmm cutoffs + PME, etc, (which also gives the 
 segfault) but this is what I am using right now for quick testing.

 The only thing keeping me from filing a redmine issue is that if I remove my 
 drug and do the FEP on one of the water molecules (using the FEP code listed 
 above), I have no segfault. Therefore it is clearly related to the drug, 
 whose parameters I built so I may have caused the problem somehow. 
 Nevertheless, the drug runs fine in water and in vacuum without the FEP code, 
 so I can't imagine what could be causing this segfault (also, the fact that 
 it's a segfault means that I don;t get any useful info from mdrun as to what 
 might be going wrong).

 Thank you,
 Chris.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: window exchange umbrella sampling

2013-07-17 Thread Michael Shirts
4.5.7 does not support Hamiltonian exchange.  It says all properties
are the same, because all the temperatures and pressures are the same
-- it won't switch the umbrellas.

On Wed, Jul 17, 2013 at 3:30 PM, Parisa akhshi...@gmail.com wrote:
 Hi Michael,

 I think that this is an issue with the gromacs version I am using (4.6.1). I
 switched to 4.5.7 and now, I get the error:

 Fatal error:
 The properties of the ... systems are all the same, there is nothing to
 exchange

 I don't understand what I am missing here. I suppose if I have 3 replicas to
 fix a molecule at positions 1.0 and 2.0 and 3.0 while the force constant is
 being exchanged between them in a range of 300 to 200, I can define 3 mdp
 files as below:

 For replica 1 at position 1.0:

 pull_group0  = POPC
 pull_weights0=
 pull_pbcatom0= 0
 pull_group1  = DNA
 pull_weights1=
 pull_pbcatom1= 0
 pull_vec1= 0.0 0.0 0.0
 pull_init1   = 0.0 0.0 1.0
 pull_rate1   = 0
 pull_k1  = 300
 pull_kB1 = 200
 free_energy  = yes
 restraint-lambdas  =  1.0 0.5 0.0
 init-lambda-state   =  0
 nstdhdl = 10

 For replica 2 at position 2.0:

 pull_group0  = POPC
 pull_weights0=
 pull_pbcatom0= 0
 pull_group1  = DNA
 pull_weights1=
 pull_pbcatom1= 0
 pull_vec1= 0.0 0.0 0.0
 pull_init1   = 0.0 0.0 2.0
 pull_rate1   = 0
 pull_k1  = 300
 pull_kB1 = 200
 free_energy  = yes
 restraint-lambdas  =  1.0 0.5 0.0
 init-lambda-state   =  1
 nstdhdl = 10

 For replica 3 at position 3.0:

 pull_group0  = POPC
 pull_weights0=
 pull_pbcatom0= 0
 pull_group1  = DNA
 pull_weights1=
 pull_pbcatom1= 0
 pull_vec1= 0.0 0.0 0.0
 pull_init1   = 0.0 0.0 3.0
 pull_rate1   = 0
 pull_k1  = 300
 pull_kB1 = 200
 free_energy  = yes
 restraint-lambdas  =  1.0 0.5 0.0
 init-lambda-state   =  2
 nstdhdl = 10

 Could you please comment if I am missing something here? Unfortunately,
 there is not much about this posted before.

 Many thanks,

 Parisa










 --
 View this message in context: 
 http://gromacs.5086.x6.nabble.com/window-exchange-umbrella-sampling-tp5009894p5009947.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: gmx 4.6.1, Expanded ensemble: weird balancing factors

2013-07-17 Thread Michael Shirts
 It seems to have something to do with the integrator/pressure-coupling.

that is what I expected based on some preliminary testing earlier.

When
 I ran the tutorial on
 http://www.alchemistry.org/wiki/GROMACS_4.6_example:_Ethanol_solvation_with_expanded_ensemble,
 everything seems fine

OK good!

 unless I switched to sd integrator instead of md-vv

There may be a problem with sd and expanded ensemble. This MAY have
been fixed in 4.6.3, but I'll need a little bit of time to check.

 (and Pcoupl from MTTK to berendsen),

Berendsen is provably wrong. The ensemble is incorrect.  It is evil
whenever a true distribution is needed. There are warnings that should
be printed when you set up the system.  Perhaps they should be made
stronger with expanded ensemble!

in which case SHAKE didn't seem to be
 able to settle the constraints.

Expanded ensemble is rough on constraints.  The integration is
formally correct, but it is harder to converge pressure, especially
with constraints.  The larger the system, however, the more stable it
is.

 file from the tutorial. Any idea?

Well, one solution is to run your molecule with the parameters in the
.mdp!  Would that solve your problem for now?

Longer term. I will check on the SD to see if it is fixed in 4.6.3 --
I seem to re call this, but I'm not sure.   However, I am not sure
that it can be made formally correct.  I will either fix this or make
the warnings more clear -- or not permit it at all.
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: gmx 4.6.1, Expanded ensemble: weird balancing factors

2013-07-17 Thread Michael Shirts
Unfortunately, the semiiso fix for MTTK is not trivial, and will have
to wait until 5.0.

Berendsen volume + expanded ensemble would likely give very bad
results, especially since surface area fluctuations are important.

Replica exchange works with sd and parrinello-Rahman, if that helps!

I am working on checking various combinations of expanded ensemble and
integrators that can more easily be supported, but that may take a few
weeks.







On Wed, Jul 17, 2013 at 6:22 PM, Dejun Lin dejun@gmail.com wrote:
 The test I did was using gmx-4.6.3. I'm currently working on a membrane
 system, which requires semiisotropic pressure coupling. I know that MTTK in
 gmx-4.6.3 doesn't support semiiso yet so the only combination available for
 expanded ensemble is md-vv + v-rescale tcoupl + berendsen pcoupl, which I'm
 testing right now. A preliminary run gave me some sane dG values although
 it's not optimal for data collection as you pointed out. I wonder if it's a
 trivial fix that might have been done to add semiiso to MTTK?

 Thanks again for your help!


 2013/7/17 Michael Shirts mrshi...@gmail.com

  It seems to have something to do with the integrator/pressure-coupling.

 that is what I expected based on some preliminary testing earlier.

 When
  I ran the tutorial on
 
 http://www.alchemistry.org/wiki/GROMACS_4.6_example:_Ethanol_solvation_with_expanded_ensemble
 ,
  everything seems fine

 OK good!

  unless I switched to sd integrator instead of md-vv

 There may be a problem with sd and expanded ensemble. This MAY have
 been fixed in 4.6.3, but I'll need a little bit of time to check.

  (and Pcoupl from MTTK to berendsen),

 Berendsen is provably wrong. The ensemble is incorrect.  It is evil
 whenever a true distribution is needed. There are warnings that should
 be printed when you set up the system.  Perhaps they should be made
 stronger with expanded ensemble!

 in which case SHAKE didn't seem to be
  able to settle the constraints.

 Expanded ensemble is rough on constraints.  The integration is
 formally correct, but it is harder to converge pressure, especially
 with constraints.  The larger the system, however, the more stable it
 is.

  file from the tutorial. Any idea?

 Well, one solution is to run your molecule with the parameters in the
 .mdp!  Would that solve your problem for now?

 Longer term. I will check on the SD to see if it is fixed in 4.6.3 --
 I seem to re call this, but I'm not sure.   However, I am not sure
 that it can be made formally correct.  I will either fix this or make
 the warnings more clear -- or not permit it at all.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: gmx 4.6.1, Expanded ensemble: weird balancing factors

2013-07-16 Thread Michael Shirts
Hi, all-

This not a problem with W-L, but is instead something that is wrong
with a particular combination of mdp options that are not working for
expanded ensemble simulations.  W-L can equilibrate to incorrect
distributions because it decreases the weights too fast (more on that
later), but that is not the problem here. The option wl-oneovert
allows switching to a slower scaling that is more likely to converge.

The reason that it is not a W-L problem is that after the WL weights
are equilibrated, it is equilibrium sampling in the expanded ensemble.
 If the W-L equilibration was a problem, then the distribution of
states would not be flat.  They are flat, so therefore the expanded
ensemble dynamics are wrong.

I have example expanded ensemble setups.

http://www.alchemistry.org/wiki/GROMACS_4.6_example:_Ethanol_solvation_with_expanded_ensemble

These mdp files should work. Note that you should be able to swap out
any top, and it will still work.  If you get CORRECT results (and it
should take just a few ns to see) with these tops, then I will go
through and try to figure out which differences are causing the
problems.

Thanks for the patience while we test this new functionality.
Frequently when one puts new code in the hands of others it breaks in
ways that were not foreseen!

On Tue, Jul 16, 2013 at 12:53 PM, Dejun Lin dejun@gmail.com wrote:
 I have a similar question here. Can anyone tell if the Wang-Landau algorithm
 in lambda space is robust in that it doesn't depend on the choice of the
 convergence factor (weight-equil-wl-delta), the flatness of the histogram
 (wl-ratio) and/or the frequency of trial move (nstexpanded), especially in
 the case of barely overlapping energy distribution corresponding to
 different lambdas? I can imagine in the extreme case, with only 2 lambdas (
 l = 1 or 0), the gap between the 2 energy distributions might be so large
 such that for most of the time, trial moves from the center of
 distribution 1 would frequently end up in the tail of distribution 0. In
 this case, the Wang-Landau weights would be biased if there's not enough
 time for the system to relax back to equilibrium.

 Thanks,
 Dejun



 --
 View this message in context: 
 http://gromacs.5086.x6.nabble.com/gmx-4-6-1-Expanded-ensemble-weird-balancing-factors-tp5007681p5009892.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] window exchange umbrella sampling

2013-07-16 Thread Michael Shirts
You need to have different pull parameters at the end states. Right
now, pull-kB1 is not defined in your code, so there is nothing to
interpolate to: it assume pull-kB1 = pull-k1.

Longer scale -- one would want to define reference distances that
change with lambda within the same simulation, but looking over the
code, I'm not seeing anything.  One would want to add a pull_vect1B to
change the location of the restraint as a function of lambda, but this
isn't yet defined.  That should be something we look at in the future
. . .

On Tue, Jul 16, 2013 at 2:41 PM, Parisa Akhshi akhshi...@gmail.com wrote:
 I am trying to use Hamiltonian replica exchange as implemented in Gromacs.
 The idea is to use different pulling umbrella forces and then obtain a PMF
 profile. I am specifically interested in exchanging force constants between
 windows.

 To do a quick test, I successfully ran replica exchange in which the
 temperature is exchanged. In this case, I prepared for example 8 tpr files
 at 8 different temperatures and run them using

 mpirun -np 8  mdrun_mpi -s prefix_.tpr -multi 8 -replex 100

 And it ran just fine.

 For window exchange umbrella sampling, though, I am confused with the input
 preparation. I am not sure how to specify in the .mdp file how the values
 of pull force should be used to define different force constants for
 different replicas. Also, when submitting the job, how should I ask the
 program to exchange pull forces between windows? I tried using  this:

 http://lists.gromacs.org/pipermail/gmx-users/2011-April/060448.html

 But, it gives me the error:

 Fatal error:
 The properties of the ... systems are all the same, there is nothing to
 exchange

 I noticed a previous discussion on below link:

 http://gromacs.5086.x6.nabble.com/Hamiltonian-replica-exchange-umbrella-sampling-with-gmx-4-6-td5006187.html

 However, I am not sure how to use the restraint-lambda parameter exactly.
 Is there any example, ... how to use it?

 I have copied the related portion of one of the the .mdp files below for
 moving MOL2 with respect to MOL1

 Thanks for your help,

 Parisa




 .
 .
 .
 pull_group0  = MOL1
 pull_weights0=
 pull_pbcatom0= 0
 pull_group1  = MOL2
 pull_weights1=
 pull_pbcatom1= 0
 pull_vec1= 0.0 0.0 0.0
 pull_init1   = 0.0 0.0 0.0
 pull_rate1   = 0
 pull_k1  = 300

 free_energy  = yes
 restraint-lambdas  =  0.0 0.5 1.0 0.0 0.5 1.0 0.0 0.5 1.0
 .
 .
 .
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: gmx 4.6.1, Expanded ensemble: weird balancing factors

2013-07-16 Thread Michael Shirts
That's a good question.  My understanding and experience is the
difference in required overlap in replica exchange and expanded
ensemble methods is not significant.  There are probably some more
detailed studies out there. Expanded ensemble can be somewhat more
efficient (see a paper by Park and Pande) But you really do want to
have at least 30% of the states moving each time in any case.

However, you DO need less overlap to get fairly good free energy
estimates.  I.e. the overlap for BAR to work fairly well is lower than
the what you need for good mixing.  That's more anecdotal -- I'd have
to dig a bit to get you good results on that . . .

On Tue, Jul 16, 2013 at 3:19 PM, Dejun Lin dejun@gmail.com wrote:
 Hi Michael,

 Thanks for the reply. Just a quick follow-up. Do you think the overlap of
 energy histogram between different lambdas matter for lambda-dynamics in
 general? (Maybe not in this particular case or in the case of the tutorial
 you just posted.) I wonder if we need as many intermediate lambdas as in
 the case of replica exchange since the weights compensate the difference in
 energy.

 Thanks,
 Dejun


 2013/7/16 Michael Shirts mrshi...@gmail.com

 Hi, all-

 This not a problem with W-L, but is instead something that is wrong
 with a particular combination of mdp options that are not working for
 expanded ensemble simulations.  W-L can equilibrate to incorrect
 distributions because it decreases the weights too fast (more on that
 later), but that is not the problem here. The option wl-oneovert
 allows switching to a slower scaling that is more likely to converge.

 The reason that it is not a W-L problem is that after the WL weights
 are equilibrated, it is equilibrium sampling in the expanded ensemble.
  If the W-L equilibration was a problem, then the distribution of
 states would not be flat.  They are flat, so therefore the expanded
 ensemble dynamics are wrong.

 I have example expanded ensemble setups.


 http://www.alchemistry.org/wiki/GROMACS_4.6_example:_Ethanol_solvation_with_expanded_ensemble

 These mdp files should work. Note that you should be able to swap out
 any top, and it will still work.  If you get CORRECT results (and it
 should take just a few ns to see) with these tops, then I will go
 through and try to figure out which differences are causing the
 problems.

 Thanks for the patience while we test this new functionality.
 Frequently when one puts new code in the hands of others it breaks in
 ways that were not foreseen!

 On Tue, Jul 16, 2013 at 12:53 PM, Dejun Lin dejun@gmail.com wrote:
  I have a similar question here. Can anyone tell if the Wang-Landau
 algorithm
  in lambda space is robust in that it doesn't depend on the choice of the
  convergence factor (weight-equil-wl-delta), the flatness of the histogram
  (wl-ratio) and/or the frequency of trial move (nstexpanded), especially
 in
  the case of barely overlapping energy distribution corresponding to
  different lambdas? I can imagine in the extreme case, with only 2
 lambdas (
  l = 1 or 0), the gap between the 2 energy distributions might be so large
  such that for most of the time, trial moves from the center of
  distribution 1 would frequently end up in the tail of distribution 0.
 In
  this case, the Wang-Landau weights would be biased if there's not enough
  time for the system to relax back to equilibrium.
 
  Thanks,
  Dejun
 
 
 
  --
  View this message in context:
 http://gromacs.5086.x6.nabble.com/gmx-4-6-1-Expanded-ensemble-weird-balancing-factors-tp5007681p5009892.html
  Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
  --
  gmx-users mailing listgmx-users@gromacs.org
  http://lists.gromacs.org/mailman/listinfo/gmx-users
  * Please search the archive at
 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-requ...@gromacs.org.
  * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support

Re: [gmx-users] Re: window exchange umbrella sampling

2013-07-16 Thread Michael Shirts
Ah, this is a force field issue -- urey-bradley terms are not
supported free energy calculations. However, since only restraints are
changing, this warning doesn't really need to be there.

It would be relatively simple to put in a check to allow this to work,
but it might take a week or two to get around to it.

I'll file a redmine in the meantime, in case someone else can get to it first!

http://redmine.gromacs.org/issues/1302



On Tue, Jul 16, 2013 at 5:32 PM, Parisa akhshi...@gmail.com wrote:
 Thanks for quick reply. I fixed the issue with kB1. Below is what I have
 changed in the .mdp file. I use this to run it mpirun -np 3  mdrun_mpi -s
 prefix_.tpr -multi 3 -replex 6 -dhdl a.dhdl. But I get this error (I am
 using version 4.6.1):

 Fatal error:
 Function type U-B not implemented in ip_pert
 For more information and tips for troubleshooting, please check the GROMACS

 Here is the input file:
 .
 .
 .
 pull_group0  = POPC
 pull_weights0=
 pull_pbcatom0= 0
 pull_group1  = DNA
 pull_weights1=
 pull_pbcatom1= 0
 pull_vec1= 0.0 0.0 0.0
 pull_init1   = 0.0 0.0 0.0
 pull_rate1   = 0
 pull_k1  = 300
 pull_kB1 = 200

 free_energy  = yes
 restraint-lambdas  =  1.0 0.5 0.0
 init-lambda-state   =  0
 nstdhdl = 10
 .
 .
 .

 Any idea about it?

 Thanks,

 Parisa




 --
 View this message in context: 
 http://gromacs.5086.x6.nabble.com/window-exchange-umbrella-sampling-tp5009894p5009899.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Larger number of decimal places for coordinates with velocities

2013-07-05 Thread Michael Shirts
Have you checked out the -ndec option for trjconv?  If you have a high
precision format (.trr, or .xtc if they are stored with sufficient
precision) you can print out a .gro file (that gromacs can read) with
higher precision.

Gromacs can read .gro files with increased precisions in the
coordinates -- as long as the precision does not CHANGE anywhere in
the file.

On Fri, Jul 5, 2013 at 11:16 AM, C.M.Sampson cms1...@soton.ac.uk wrote:
 Dear all,

 I use Gromacs 4.5.5 with the AMBER ff99SB force field.

 I'm working on a method that requires me to run short NVE simulations
 one after the other, but randomly generating velocities at the start of
 each simulation. i.e. I would run 10 ps, use the final structure with
 new random velocities and run another 10 ps.

 I had issues with the final potential energy of the previous simulation
 not matching the first potential energy of the current simulation, but
 managed to fix that by converting the coordinates from the .xtc file to
 a .gro file.

 My problem is that when I then put the velocities into my new .gro file
 the simulation breaks after the first step and the kinetic energy is way
 too high:

 
Step   Time Lambda
   00.00.0

Energies (kJ/mol)
  Bond  Angle Ryckaert-Bell.  LJ-14 Coulomb-14
 5.92468e+031.01216e+037.68471e-016.38163e-013.95610e+00
  LJ (SR)   Coulomb (SR)   Coul. recip.  PotentialKinetic En.
  8.66770e+03   -2.27543e+035.12418e+046.45763e+042.11278e+07
Total EnergyTemperature Pressure (bar)
 2.11924e+077.81749e+051.08075e+07
 

 The temperature should be 300K and if I use a .gro file with less
 decimal places these same velocities work.

 I thought it had to be the way I had formatted my .gro file, but found
 the following:

 
 This format is fixed, ie. all columns are in a fixed position.
 Optionally (for now only yet with trjconv) you can write gro files with
 any number of decimal places, the format will then be n+5 positions with
 n decimal places (n+1 for velocities) in stead of 8 with 3 (with 4 for
 velocities). Upon reading, the precision will be inferred from the
 distance between the decimal points (which will be n+5). Columns contain
 the following information (from left to right):
   * residue number (5 positions, integer)
   * residue name (5 characters)
   * atom name (5 characters)
   * atom number (5 positions, integer)
   * position (in nm, x y z in 3 columns, each 8 positions with 3
 decimal places)
   * velocity (in nm/ps (or km/s), x y z in 3 columns, each 8
 positions with 4 decimal places)
 

 within http://manual.gromacs.org/online/gro.html .

 An extract of my .gro file can be seen below:

 
 Generated by trjconv : 2168 system t=  15.0
  2168
 1ETH C11   2.735383   2.672010   1.450194  0.2345 -0.1622
 0.2097
 1ETHH112   0.015804   2.716597   1.460588  0.8528 -0.7984
 0.6605
 1ETHH123   2.744822   2.565544   1.409227 -2.3812  2.8618
 1.8101
 

 I was wondering if there's a way to use a larger number of decimal
 places for the coordinates with velocities?

 Best Wishes

 Chris Sampson

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] a question concerning on entropy

2013-07-04 Thread Michael Shirts
No.

This is a statistical mechanical issue, not a GROMACS issue.  For
interacting systems, entropy is a quantity describing the system as a
whole, and cannot be defined for different parts of the system, at
least not in any way such that the individual components can be added
together.

I'm also not certain the md.edr file will give the entropy of the system . . .

On Thu, Jul 4, 2013 at 2:28 PM, Albert mailmd2...@gmail.com wrote:
 Hello :

  I've got a question about the the entropy. As we all know that in the
 md.edr file it will give us the entropy value of the system along the
 simulations.

 However, my system is a protein/membrane system, and I am only would like to
 make statics for the protein/water related entropy. I am just wondering, is
 it possible to do this?

 thank you very much

 Albert
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: 1-4 interactions free energy calculations

2013-06-25 Thread Michael Shirts
Hi, Sonia-

Gromacs 4.6.2 (some bug fixes vs 4.6.1) with the example files you
point out from Alchemistry.org should work well for expanded ensemble.
 David Mobley and I have been validating expanded ensemble and replica
exchange, and the files posted there now are stable for all sizes of
systems, and give reasonable results.  expanded ensemble with NPT is
somewhat unstable for smaller systems (1000 total molecules), though
is much more stable for larger systems -- the posted parameter should
work.

However, note that you have to be careful with expanded ensemble, as
the weights can sometime converge to a poor value.

I'm working on some more guidance on expanded ensemble besides the
examples on Alchemistry.org, but it may take a few weeks.  Check back
in a few weeks for some more information on Alchemistry.org.

On Tue, Jun 25, 2013 at 4:36 PM, Sonia Aguilera
sm.aguiler...@uniandes.edu.co wrote:
 Hi Justin,

 Thank you for your answer. I’m performing several tests to see what is the
 best for my system. The reason I decided not to couple the intramolecular
 interactions is because I think that the annihilation of the molecule will
 lead to very extreme configurations and that will affect my phase-space
 overlap. I know I can just increase the number of intermediate states, but I
 wanted to first test my theory. Also, I have limited resources to run.

 There is a gap between the indirect and direct calculation methods  for my
 system. If I set up couple-intramol=no, then I can´t use domain
 decomposition but particle decomposition in all simulations (charge groups
 are too big, and there is no domain decomposition… ) and I got the 1-4
 interactions warning. If I restrain the h-bonds, then the minimization
 converges but not to Fmax less than 100 in only 499 steps and 12 minutes (8
 cores machines). If I continue and run the NVT stabilization, then it takes
 17h20min to complete only 100 ps, and I still got the 1-4 interactions
 warning. I still have to run the NPT and 3 ns MD, which will take too much
 time. I think it worth it  because I think it will help with the phase-space
 overlap, and I will only have to perform 20 series of simulations (vdw and
 coulomb separately) to get the free energy change.

 On the other hand, if I set up couple-intramol=yes, then I have to run 20
 extra simulations in vacuo to counter the annihilation effect. I think that
 it will also affect the phase-space overlap. The good thing is that I can
 use domain decomposition and the minimization converges to machine precision
 (takes 7000 steps and around 40 minutes in a machine with the same
 characteristics that the previous one). Also, the same 100-ps NVT takes only
 3h20min. If I think of the overall picture, it seems that the total time I
 will spend doing dd (coupleintramol=yes) rather than pd (coupleintramol=no)
 will be less. However I’m worried about the phase-space overlap and the
 available resources I have for running since it makes difficult to just
 increase the intermediate states.

 Also, I think I can try doing the same but with the 4.6.1 version. It looks
 like I can now couple both vdw and coulomb in the same simulation
 (http://www.alchemistry.org/wiki/GROMACS_4.6_example:_Direct_ethanol_solvation_free_energy).
 I have also wanted to try expanded ensemble to improve my sampling, but I
 have not found that much documentation about the implementation on gromacs.
 It would be great if you can provide an example.

 So, do you know another way to improve my sampling and the phase-space
 overlap that can help me to solve my problem?  Also, I’m working with a
 charged molecule. Would it help in any aspect if I neutralize the system?
 Thank you very much in advance,

 Best regards,

 Sonia Aguilera
 Graduate Assistant
 Department of Chemical Engineering
 Universidad de los Andes




 --
 View this message in context: 
 http://gromacs.5086.x6.nabble.com/1-4-interactions-free-energy-calculations-tp5009364p5009378.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Enthalpy Confusion

2013-06-11 Thread Michael Shirts
If you are computing enthaply in the NPT ensemble, P is constant, and
is the applied pressure.

The pressure quantity calculated from the KE and the virial is not
the pressure.  It is a quantity that when averaged over time is equal
the pressure.  Only the average is meaningful macroscopically.

If you are computing enthalpy in another ensemble (which is possible,
though it may be harder to interpret) then you would use the average
pressure from this quantity.



On Tue, Jun 11, 2013 at 3:08 PM, Jeffery Perkins jeffery.perk...@ufv.ca wrote:
 that's what i thought, and what i tried to do, my pressure is a bit higher
 then that, we want a Lennard-Jones liquid so it's running at 1000+ bar, and
 while I agree that gromacs is giving H as Etot + pV it appears that when i
 calculate pV i get a different value from what g_energy returns for it I did
 p in Pa, V in m^3, as the whole box, to get J of energy, and then multiply
 by 6.02E23 to get J/mol of my box, and then convert down to kJ/mol to be the
 same units as g_energy.

 When you say the applied pressure you mean that p = ref_p? or the calculated
 pressure from the virial and Ke?

 Thanks for the help,

 Jeffery



 --
 View this message in context: 
 http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009058.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


[gmx-users] 2013 GROMACS USA Workshop and Conference

2013-06-11 Thread Michael Shirts
Dear GROMACS users list members:

We are pleased to announce the 2013 GROMACS USA Workshop and
Conference at the University of Virginia in Charlottesville, Virginia,
on September 13-15th.  This will be the first full GROMACS workshop
held here in the U.S.

The workshop and conference will consist of two days of tutorials,
discussions of future plans for GROMACS, face time with a large
percentage of the developers, plenary software and application
sessions, and a last day of development planning to which attendees
are also invited to help influence the future directions of GROMACS.

Please visit http://faculty.virginia.edu/gromacsworkshop for
registration and for much more information about the program.

For specific questions about registering or logistics after visiting
the website, please write to
gromacsworkshop-registrat...@virginia.edu.

Sincerely,
The 2013 GROMACS USA Workshop and Conference Steering Committee
Michael Shirts (chair)
Angel Garcia
Berk Hess
Yu-Shan Lin
Erik Lindahl
Peter Kasson
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Enthalpy Confusion

2013-06-11 Thread Michael Shirts
 or should i be doing  U+V*ref_p  = H?

More specifically, U + V*ref_p = H

H isn't really meaningful thing.  I mean, you can define something
such that H* = H, but that's not really thermodynamics.

 example system gives H = -1168 kJ/mol and i find H = -725 kJ/mol either

Interesting.  What material at what phase conditions?  For liquids,
the PV contribution should be very small.
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Free Energy Calculations in Gromacs

2013-06-10 Thread Michael Shirts
An important final point is that you can always see EXACTLY what
grompp is putting into the B state by running gmxdump on the resulting
tpr.  It's a LOT of information, but all in text all the interactions
are listed explicitly there.

On Mon, Jun 10, 2013 at 6:20 PM, Justin Lemkul jalem...@vt.edu wrote:


 On 6/10/13 2:50 PM, JW Gibbs wrote:

 Hi,

 I have been trying to perform the simulations using the amber forcefield,
 in
 which the [ pairtypes ] directive is not defined explicitly in the
 ffnonbonded.itp file, but rather the 1-4 interactions are generated as per
 the defaults section in the forcefield.itp file. I was wondering what
 happens to these 1-4 interactions at lambda=1 state?

 Suppose 1 corresponds to CA and 4 corresponds to NA at state A (lambda =
 0)
 and the CA_per and NA_per are the corresponding atomtypes at state B
 (lambda
 = 1).

 It is defined in the topology file that

 ...
 ...
 [ pairs ]

 1  4   1
 
 

 So does it mean that at state B (lambda = 1), the 1-4 nonbonded
 interactions
 will be calculated between CA_per and  NA_per?

 Although the [ pairs_nb ] section is described in the manual, I would
 appreciate if someone elaborates a little more.


 The actual answer depends on exactly how you're treating the system.  Are
 you using [pairs_nb] in your topology or just [pairs]?  The outcome will be
 different depending on which you're using.  Also note that, as the manual
 says, you don't really need to mess with [pairs_nb] at all; you can achieve
 unscaled internal interactions using couple-intramol = no in the .mdp
 file.  Without seeing the .mdp file, it's even more difficult to know what
 you're doing and what you should expect.

 -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

 

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Difference between the electrostatic treatments PME/Cut-offs and Reaction Field

2013-06-05 Thread Michael Shirts
 It should also be noted (and obvious now that I actually look into it) that 
 using dispersion correction results in both the latent heat of vapourisation 
 and density of the alkanes being over estimated (for both Cut-off and 
 Reaction Field, and by the same amount).

That may not be quite the right way to think about it.  The dispersion
correction gives a result that (for homogeneous fluids) is essentially
cutoff independent.  So it's a little truer calculation in that it
is relatively independent of your potential cutoff, which is a
nonphysical quantity that's really an artifact of the simulation.

So saying that the properties are overestimated with a dispersion
correction, you're essentially saying that the latent heat of
vaporization and the densities would also be overestimated if your ran
with longer cutoffs, and that they are only right if run with quite
short cutoffs.

On Tue, Jun 4, 2013 at 10:51 PM, Dallas Warren dallas.war...@monash.edu wrote:
 Problem solved.

 Just as well I held off reporting the issue in full until I had explored 
 everything, I would ended have look a bit stupid ;-)  But asking the initial 
 question helped direct my thinking to the reason.

 The issue was the difference in van der Waals cut off between the PME/Cut-off 
 and Reaction Field methods that I was using, 0.9/0.9/0.9 vs 0.8/1.4/1.4  The 
 difference was hidden in the results until you turned off Dispersion 
 Correction, which was what was confusing me and final led to realizing what 
 is going on.  My forehead hurt from slapping it after coming to the 
 realization ...

 It should also be noted (and obvious now that I actually look into it) that 
 using dispersion correction results in both the latent heat of vapourisation 
 and density of the alkanes being over estimated (for both Cut-off and 
 Reaction Field, and by the same amount).

 Catch ya,

 Dr. Dallas Warren
 Drug Discovery Biology
 Monash Institute of Pharmaceutical Sciences, Monash University
 381 Royal Parade, Parkville VIC 3052
 dallas.war...@monash.edu
 +61 3 9903 9304
 -
 When the only tool you own is a hammer, every problem begins to resemble a 
 nail.


 -Original Message-
 From: gmx-users-boun...@gromacs.org [mailto:gmx-users-
 boun...@gromacs.org] On Behalf Of Mark Abraham
 Sent: Thursday, 30 May 2013 9:32 PM
 To: Vitaly Chaban; Discussion list for GROMACS users
 Subject: Re: [gmx-users] Difference between the electrostatic
 treatments PME/Cut-offs and Reaction Field

 Things should be identical - any quantity computed from a zero charge
 has
 to be zero :-).

 Mark

 On Thu, May 30, 2013 at 1:26 PM, Dr. Vitaly Chaban
 vvcha...@gmail.comwrote:

  Hmmm...
 
  And what does the observed difference look like, numerically?
 
 
 
 
 
  On Thu, May 30, 2013 at 10:14 AM, Mark Abraham
 mark.j.abra...@gmail.com
  wrote:
 
   No charges = no problem. You can trivially test this yourself with
 mdrun
   -rerun ;-)  Manual 4.1.4 talks about what RF is doing.
  
   Mark
  
  
   On Thu, May 30, 2013 at 6:38 AM, Dallas Warren
 dallas.war...@monash.edu
   wrote:
  
In a system that has no charges, should we observe a difference
 between
simulations using PME/Cut-offs or Reaction Field?
   
From my understanding there should not be, since there are no
 charges
which treatment you use shouldn't' make a difference.
   
However, it does and I am trying to work out why.
   
Any suggestions on the reason?
   
What is it that Reaction Field is doing, does it influence
 anything
  other
than long range charge interactions?
   
Catch ya,
   
Dr. Dallas Warren
Drug Discovery Biology
Monash Institute of Pharmaceutical Sciences, Monash University
381 Royal Parade, Parkville VIC 3052
dallas.war...@monash.edu
+61 3 9903 9304
-
When the only tool you own is a hammer, every problem begins to
  resemble
   a
nail.
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
   
   --
   gmx-users mailing listgmx-users@gromacs.org
   http://lists.gromacs.org/mailman/listinfo/gmx-users
   * Please search the archive at
   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-requ...@gromacs.org.
   * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
  
  --
  gmx-users mailing listgmx-users@gromacs.org
  http://lists.gromacs.org/mailman/listinfo/gmx-users
  * Please search the archive at
  http://www.gromacs.org/Support/Mailing_Lists/Search 

Re: [gmx-users] Re: Nose-Hover chains for membrane protein simulation

2013-06-02 Thread Michael Shirts
And I my point was I didn't think that there was going to be a
measurable ergodicity difference between NH chains and NH.  Given the
size of the system, the chaoticness of atomic motions will likely give
you configurational sampling indistinguishable from full ergodicity.
Most of the errors that are solved by NH chains are for small toy
systems.

On Sun, Jun 2, 2013 at 2:17 AM, James Starlight jmsstarli...@gmail.com wrote:
 Michael,


 thanks for suggestions.

 the main reason of ussage N-H with chains is the assumption that simple N-H
 does not provide ergodicity of the system assuming that I'd like to sample
 all temperature acceptable conformations on the 100 ns trajectory.

 But as I understood the chain regime does not compatible with the membrane
 protein simulation due to the artifacts arising with MTTK batostat.

 James

 2013/6/1 Michael Shirts mrshi...@gmail.com

 I can't think of any instance where nose-hoover chains provides an
 advantage over nose-hoover in a large system -- all the demonstrations
 of superiority are in model systems that are not particularly chaotic.
  As the system gets more chaotic, it matters less.

 I would go with md, nose-hoover (w/o chains), and Parrinello-Rahman
 with semiisotropic scaling.

 On Sat, Jun 1, 2013 at 1:07 AM, James Starlight jmsstarli...@gmail.com
 wrote:
  Performing 5ns simulation after 2 ns equilibration in NPT ensemble (MTTK
  barostat 5ps coupling ) I've observed non-physical behaviour of my system
  with the constant drift of the protein molecule as the rigid body in the
  y-z plane
 
  Energy  Average   Err.Est.   RMSD  Tot-Drift
 
 ---
  Pressure   -760.137 --193.776218.059
  (bar)
 
 
  From manual I've noticed that MMTK doest not support *semiisotropic
  scalling.  *Doest it means that with the Nose-hover chains and md-vv I
  should use only weaker coupling during productions runs (I cant use
  Parinello;s barostat with such options too)
 
  Thanks for help
 
  James
 
 
 
  2013/5/31 James Starlight jmsstarli...@gmail.com
 
  Dear Gromacs users!
 
  I'd like to perform simulation of the membrane protein in lipid-water
  system using Nose-Hover with chains.
 
  From manual I've found that with such thermostat I should use (1) md-vv
  integrator (2) MTTK  instead of Parinello's batostat  and (3) shake
 instead
  of LINCS.
 
 
  How doest such options compatible with the simulation of membrane
 proteins
  in general ? On what other options should I pay attention during
 simulation
  of membrane protein with NH chains ?
 
 
 
  Thanks for help,
  James
 
  --
  gmx-users mailing listgmx-users@gromacs.org
  http://lists.gromacs.org/mailman/listinfo/gmx-users
  * Please search the archive at
 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-requ...@gromacs.org.
  * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Nose-Hover chains for membrane protein simulation

2013-06-01 Thread Michael Shirts
I can't think of any instance where nose-hoover chains provides an
advantage over nose-hoover in a large system -- all the demonstrations
of superiority are in model systems that are not particularly chaotic.
 As the system gets more chaotic, it matters less.

I would go with md, nose-hoover (w/o chains), and Parrinello-Rahman
with semiisotropic scaling.

On Sat, Jun 1, 2013 at 1:07 AM, James Starlight jmsstarli...@gmail.com wrote:
 Performing 5ns simulation after 2 ns equilibration in NPT ensemble (MTTK
 barostat 5ps coupling ) I've observed non-physical behaviour of my system
 with the constant drift of the protein molecule as the rigid body in the
 y-z plane

 Energy  Average   Err.Est.   RMSD  Tot-Drift
 ---
 Pressure   -760.137 --193.776218.059  (bar)


 From manual I've noticed that MMTK doest not support *semiisotropic
 scalling.  *Doest it means that with the Nose-hover chains and md-vv I
 should use only weaker coupling during productions runs (I cant use
 Parinello;s barostat with such options too)

 Thanks for help

 James



 2013/5/31 James Starlight jmsstarli...@gmail.com

 Dear Gromacs users!

 I'd like to perform simulation of the membrane protein in lipid-water
 system using Nose-Hover with chains.

 From manual I've found that with such thermostat I should use (1) md-vv
 integrator (2) MTTK  instead of Parinello's batostat  and (3) shake instead
 of LINCS.


 How doest such options compatible with the simulation of membrane proteins
 in general ? On what other options should I pay attention during simulation
 of membrane protein with NH chains ?



 Thanks for help,
 James

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] issue in replica exchange

2013-05-03 Thread Michael Shirts
Summarizing!

On Fri, May 3, 2013 at 12:31 AM, XAvier Periole x.peri...@rug.nl wrote:

 Are confirming that you reproduce the problem with gmx-4.6.1 or simply 
 summarizing in case we lose track :))

 On May 2, 2013, at 23:31, Michael Shirts mrshi...@gmail.com wrote:

 So to summarize -- the problem appears to be with particle decomposition.

 On Thu, May 2, 2013 at 4:15 PM, XAvier Periole x.peri...@rug.nl wrote:

 I'll look at the 4.6.1 version next week, I could install it but I got a 
 conflict between the environmental variable defining openMP variable but I 
 turned it off during compilation …

 You could try to run on particle decomposition to see if you get a problem 
 … it should one quite quick.

 On May 2, 2013, at 2:36 PM, Michael Shirts mrshi...@gmail.com wrote:

 Both.  So if 4.6.1 doesn't work, I want to know so we can patch it
 before 4.6.2 comes out.  If it does work, then there is probably stuff
 that can be backported.

 On Thu, May 2, 2013 at 8:32 AM, XAvier Periole x.peri...@rug.nl wrote:

 You mean working with or working on the code?

 I'll try gmx-4.6.1

 On May 2, 2013, at 2:26 PM, Michael Shirts mrshi...@gmail.com wrote:

 Quick check here -- is 4.6 behaving correctly?  I actually spent some
 time working on REMD in 4.6, and it seems to be behaving  correctly in
 my hands with temperature and pressure control.

 Thanks for any additional info on this!

 On Thu, May 2, 2013 at 8:18 AM, Mark Abraham mark.j.abra...@gmail.com 
 wrote:
 On Thu, May 2, 2013 at 12:58 PM, XAvier Periole x.peri...@rug.nl 
 wrote:


 I saw that redmine report, which could be related but it seems to 
 happen
 only for runs done outside the domain and particle decompositions.

 I'll fill up a red mine.

 Anything I could do to help speeding the fix?

 What'd be really nice is some thought on how one can demonstrate that 
 the
 implementation of the exchange matches what would be expected from the
 theory. For T-exchange under NVT, it is sufficient to rescale velocities
 and quantities derived from them by the correct factor. That includes
 various things like T-coupling history and integrator half-step 
 quantities
 (and does REMD with leap-frog make sense anyway?). For NPT, there's
 probably also some P-coupling quantities to scale, and the box to 
 exchange.
 Anything I've missed? Hopefully virial contributions don't matter either
 way?

 Perhaps a decent first step is to hack the code to do a self 
 exchange, by
 clearing the entire state and rebuilding with what would/should be 
 received
 from an exchange with a hypothethetical replica in an identical
 pre-exchange state. Only if the code can do that (i.e. mdrun -reprod
 produces a trajectory indistinguishable from a run that does not attempt
 this self exchange) is it worth considering proper state exchanges, and 
 the
 process of making the code do the former should illustrate what is 
 required
 for the latter.

 Mark
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org

Re: [gmx-users] issue in replica exchange

2013-05-02 Thread Michael Shirts
Quick check here -- is 4.6 behaving correctly?  I actually spent some
time working on REMD in 4.6, and it seems to be behaving  correctly in
my hands with temperature and pressure control.

Thanks for any additional info on this!

On Thu, May 2, 2013 at 8:18 AM, Mark Abraham mark.j.abra...@gmail.com wrote:
 On Thu, May 2, 2013 at 12:58 PM, XAvier Periole x.peri...@rug.nl wrote:


 I saw that redmine report, which could be related but it seems to happen
 only for runs done outside the domain and particle decompositions.

 I'll fill up a red mine.

 Anything I could do to help speeding the fix?


 What'd be really nice is some thought on how one can demonstrate that the
 implementation of the exchange matches what would be expected from the
 theory. For T-exchange under NVT, it is sufficient to rescale velocities
 and quantities derived from them by the correct factor. That includes
 various things like T-coupling history and integrator half-step quantities
 (and does REMD with leap-frog make sense anyway?). For NPT, there's
 probably also some P-coupling quantities to scale, and the box to exchange.
 Anything I've missed? Hopefully virial contributions don't matter either
 way?

 Perhaps a decent first step is to hack the code to do a self exchange, by
 clearing the entire state and rebuilding with what would/should be received
 from an exchange with a hypothethetical replica in an identical
 pre-exchange state. Only if the code can do that (i.e. mdrun -reprod
 produces a trajectory indistinguishable from a run that does not attempt
 this self exchange) is it worth considering proper state exchanges, and the
 process of making the code do the former should illustrate what is required
 for the latter.

 Mark
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] issue in replica exchange

2013-05-02 Thread Michael Shirts
Both.  So if 4.6.1 doesn't work, I want to know so we can patch it
before 4.6.2 comes out.  If it does work, then there is probably stuff
that can be backported.

On Thu, May 2, 2013 at 8:32 AM, XAvier Periole x.peri...@rug.nl wrote:

 You mean working with or working on the code?

 I'll try gmx-4.6.1

 On May 2, 2013, at 2:26 PM, Michael Shirts mrshi...@gmail.com wrote:

 Quick check here -- is 4.6 behaving correctly?  I actually spent some
 time working on REMD in 4.6, and it seems to be behaving  correctly in
 my hands with temperature and pressure control.

 Thanks for any additional info on this!

 On Thu, May 2, 2013 at 8:18 AM, Mark Abraham mark.j.abra...@gmail.com 
 wrote:
 On Thu, May 2, 2013 at 12:58 PM, XAvier Periole x.peri...@rug.nl wrote:


 I saw that redmine report, which could be related but it seems to happen
 only for runs done outside the domain and particle decompositions.

 I'll fill up a red mine.

 Anything I could do to help speeding the fix?


 What'd be really nice is some thought on how one can demonstrate that the
 implementation of the exchange matches what would be expected from the
 theory. For T-exchange under NVT, it is sufficient to rescale velocities
 and quantities derived from them by the correct factor. That includes
 various things like T-coupling history and integrator half-step quantities
 (and does REMD with leap-frog make sense anyway?). For NPT, there's
 probably also some P-coupling quantities to scale, and the box to exchange.
 Anything I've missed? Hopefully virial contributions don't matter either
 way?

 Perhaps a decent first step is to hack the code to do a self exchange, by
 clearing the entire state and rebuilding with what would/should be received
 from an exchange with a hypothethetical replica in an identical
 pre-exchange state. Only if the code can do that (i.e. mdrun -reprod
 produces a trajectory indistinguishable from a run that does not attempt
 this self exchange) is it worth considering proper state exchanges, and the
 process of making the code do the former should illustrate what is required
 for the latter.

 Mark
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] issue in replica exchange

2013-05-02 Thread Michael Shirts
So to summarize -- the problem appears to be with particle decomposition.

On Thu, May 2, 2013 at 4:15 PM, XAvier Periole x.peri...@rug.nl wrote:

 I'll look at the 4.6.1 version next week, I could install it but I got a 
 conflict between the environmental variable defining openMP variable but I 
 turned it off during compilation …

 You could try to run on particle decomposition to see if you get a problem … 
 it should one quite quick.

 On May 2, 2013, at 2:36 PM, Michael Shirts mrshi...@gmail.com wrote:

 Both.  So if 4.6.1 doesn't work, I want to know so we can patch it
 before 4.6.2 comes out.  If it does work, then there is probably stuff
 that can be backported.

 On Thu, May 2, 2013 at 8:32 AM, XAvier Periole x.peri...@rug.nl wrote:

 You mean working with or working on the code?

 I'll try gmx-4.6.1

 On May 2, 2013, at 2:26 PM, Michael Shirts mrshi...@gmail.com wrote:

 Quick check here -- is 4.6 behaving correctly?  I actually spent some
 time working on REMD in 4.6, and it seems to be behaving  correctly in
 my hands with temperature and pressure control.

 Thanks for any additional info on this!

 On Thu, May 2, 2013 at 8:18 AM, Mark Abraham mark.j.abra...@gmail.com 
 wrote:
 On Thu, May 2, 2013 at 12:58 PM, XAvier Periole x.peri...@rug.nl wrote:


 I saw that redmine report, which could be related but it seems to happen
 only for runs done outside the domain and particle decompositions.

 I'll fill up a red mine.

 Anything I could do to help speeding the fix?


 What'd be really nice is some thought on how one can demonstrate that the
 implementation of the exchange matches what would be expected from the
 theory. For T-exchange under NVT, it is sufficient to rescale velocities
 and quantities derived from them by the correct factor. That includes
 various things like T-coupling history and integrator half-step quantities
 (and does REMD with leap-frog make sense anyway?). For NPT, there's
 probably also some P-coupling quantities to scale, and the box to 
 exchange.
 Anything I've missed? Hopefully virial contributions don't matter either
 way?

 Perhaps a decent first step is to hack the code to do a self exchange, 
 by
 clearing the entire state and rebuilding with what would/should be 
 received
 from an exchange with a hypothethetical replica in an identical
 pre-exchange state. Only if the code can do that (i.e. mdrun -reprod
 produces a trajectory indistinguishable from a run that does not attempt
 this self exchange) is it worth considering proper state exchanges, and 
 the
 process of making the code do the former should illustrate what is 
 required
 for the latter.

 Mark
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http

Re: [gmx-users] Free Energy Calculations in Gromacs

2013-04-20 Thread Michael Shirts
You have to change atom types.  For example:

[ atomtypes ]
;name  bond_typemasscharge   ptype  sigma  epsilon
h1h1 0.  0.  A  2.47135e-01  6.56888e-02
h1_pert h1 0.  0.  A  2.47135e-01  3.56888e-02
; perturbed

The mass and charge can be zero, because they will be defined in the [
atoms ] part

[ atoms ]
;   nrtype  resnr residue  atom  cgnr  chargemasstypeB
 chargeB  MassB
1 h1  1 BLAHH1a 1 -0.0891.008h1_pert
 -0.030  1.008;  perturbed



On Sat, Apr 20, 2013 at 4:04 PM, HANNIBAL LECTER hanniballecte...@gmail.com
 wrote:

 Hi,

 I am trying to perform expanded ensemble simulations between 2 states A
 (lambda=0) and B (lambda =1) with 6 intermediate lambda values.

 In state B the Hamiltonian is rescaled, such that the epsilons of the vdW
 interactions, the charges, the bond, angle and dihedral constants are a
 multiple of the similar terms of State A.

 It's not quite clear to me (going through the archives of the gmx-users
 mailing list how to perform these calculations. One way to do this, is to
 use a single topology file in which the charges and the other terms are
 specified for both states A and B. However, it is not clear as to how
 should I scale the epsilons in the topology file. (My atoms are not mutated
 in state B. They are the same atoms with scaled epsilons.)

 Secondly, since I have the topology files for states A and B, is there a
 way I could perform the simulations (any particular way in grompp) where
 both the topology files can be preprocessed designating the end states and
 using the mdp options, the simulations corresponding to the intermediate
 lambda values performed??
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Hamiltonian replica exchange umbrella sampling with gmx 4.6

2013-04-02 Thread Michael Shirts
Hi, Dejun-

Right now, the vector of lambda parameters is simply vdw, coul, bonded,
restraint, temperature.  You can't have, say, 2 different coul vectors or
two different restraint vectors for different restraints.  But you can
change any of these components.

You define the vector manually by writing out the list.  So to do 2D in vdw
and restraints, you would define states:

vdw-lambdas=  0.0 0.0 0.0 0.5 0.5 0.5 1.0 1.0 1.0
restraint-lambdas  =  0.0 0.5 1.0 0.0 0.5 1.0 0.0 0.5 1.0

Multidimensional movement in parameter space with replica exchange can be
done with the md setting '-nex Nswaps'.  This does multiple pair swaps,
instead of just neighbor pair swaps.  I'd suggest something like Nswaps =
N^3, where N is the number of lambda states.  This is a close approximation
to Gibbs sampling (http://jcp.aip.org/resource/1/jcpsa6/v135/i19/p194110_s1),
where _all_ permutations are selected based on their Boltzmann weight.
 It's an approximation because with a finite number of swaps, you don't
quite get uncorrelated moves in state space, but it is rigorously correct
thermodynamically.  It is more powerful than randomly selecting a dimension
to move in, since it allows moves in all dimensions simultaneously with
proper weight.

Note that the number of lambdas states is hardcoded as 1024, but that can
(I think) be edited as desired without causing direct problems (other than
perhaps needing to make the .mdp lines longer).

We are working on more general multistate lambda vector formalism for 5.0.
 If you have suggestions, let me know.

I'd be happy to look over input files or give additional advice on a
specific setup.


On Tue, Apr 2, 2013 at 12:52 PM, Dejun Lin dejun@gmail.com wrote:

 Hi Michael,

 Do the codes now support walking in multidimensional parameter space? i.e.,
 a state is defined by a set of lambda parameters {l1,l2,l3,...,ln} and a MC
 move is attempted along one of the parameter, which is randomly picked.

 Thanks,
 Dejun



 --
 View this message in context:
 http://gromacs.5086.n6.nabble.com/Hamiltonian-replica-exchange-umbrella-sampling-with-gmx-4-6-tp5006187p5006842.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Hamiltonian replica exchange umbrella sampling with gmx 4.6

2013-03-08 Thread Michael Shirts
Great -- if it doesn't seem to be working the way it should after some
playing around, then submit it as a redmine issue, and I'll take a
look.

On Fri, Mar 8, 2013 at 2:51 AM, Joakim Jämbeck jamb...@me.com wrote:
 Dear Michael,

 Thank you for your reply.

 Yes, it is relatively clear now. Will try to play around with this later 
 today.

 Best,
 Joakim


 Date: Thu, 7 Mar 2013 08:55:31 -0500
 From: Michael Shirts mrshi...@gmail.com
 Subject: Re: [gmx-users] Hamiltonian replica exchange umbrella
   sampling with   gmx 4.6
 To: Discussion list for GROMACS users gmx-users@gromacs.org

 Hi, Joakim-

 Hamiltonian exchange only should work if there is a lambda coupling
 parameter that defines the potential at each state.  You need to
 define your pulling potential so that the coupling-lambda parameter
 can be used to define the different pulling location centers along
 your trajectory.  Does this make it clearer?
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Hamiltonian replica exchange umbrella sampling with gmx 4.6

2013-03-07 Thread Michael Shirts
Hi, Joakim-

Hamiltonian exchange only should work if there is a lambda coupling
parameter that defines the potential at each state.  You need to
define your pulling potential so that the coupling-lambda parameter
can be used to define the different pulling location centers along
your trajectory.  Does this make it clearer?


On Thu, Mar 7, 2013 at 7:26 AM, Joakim Jämbeck jamb...@me.com wrote:
 Dear gmx-users,

 I am currently trying to runt Hamiltonian replica exchange umbrella sampling 
 in hope to do some better sampling.

 I have generated a number of tpr-files along my reaction coordinate and they 
 all run fine in independent simulations. The issue comes when I would like 
 the replica exchange to start.

 The following line is used to initiate the exchange:

 mdrun_mpi -replex 1000 -s md -pf pullf -px pullx -multi 40 -maxh 0.5

 All replicas have the same temperature and the following error is what I face 
 seconds after submitting the job:

 The properties of the 40 systems are all the same, there is nothing to 
 exchange
 For more information and tips for troubleshooting, please check the GROMACS
 website at http://www.gromacs.org/Documentation/Errors

 I could simply change the temperature between the replicas by 0.001 K and it 
 would run I think. But that is not very elegant.

 Does anyone have any suggestions?

 Thanks in advance!

 Best regards,
 Joakim--
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] The time for the temperature and pressure coupling

2013-02-07 Thread Michael Shirts
Ah, now perhaps I see that I misread the question - it could have been
phrased more clearly.  If Erik understood it correctly, then the
answer to the question is: It depends on the integrator.  The
simulation is not constrained to a particular temperature or pressure
- rather, the dynamics are modified such that the ensemble is
consistent with the external temperature and pressure are specified.
Exactly where in the integration routine the dynamics are modified
depends on the method used.  See the manual for more of the details!

On Thu, Feb 7, 2013 at 3:59 AM, Erik Marklund er...@xray.bmc.uu.se wrote:
 Hi,

 Perhaps a side point: Temperature and pressure can not be seen as
 constraints to the system at any given instant in the sense that e.g. the
 instantaneous kinetic energy perfectly match the temperature at every time
 step just because you have a thermostat. Time and ensemble averages will,
 however, reflect the temperature and pressure coupling.

 Erik


 On Feb 6, 2013, at 11:28 PM, Bao Kai wrote:

 Dear Gromacs Team,

 I have a small question related to the scheme of the MD in Gromacs.

 When are the temperature and pressure constrains are enforced, before
 the update of the velocity and position or after?

 Thank you very much.

 Best Regards,
 Kai
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
 * Please don't post (un)subscribe requests to the list. Use thewww interface
 or send it to gmx-users-requ...@gromacs.org.

 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] The time for the temperature and pressure coupling

2013-02-06 Thread Michael Shirts
The coordinates and velocities that are printed (and that are used to
calculate the properties like energy, virial, etc) are always
consistent with the constraints.  The exact order of how things are
done often depends on the integrator.  For example, velocity scaling
can be done before or after constraints, because there are no
velocities parallel to the bond vector, so after velocity scaling, the
constraints are still obeyed.

On Wed, Feb 6, 2013 at 5:28 PM, Bao Kai paeanb...@gmail.com wrote:
 Dear Gromacs Team,

 I have a small question related to the scheme of the MD in Gromacs.

 When are the temperature and pressure constrains are enforced, before
 the update of the velocity and position or after?

 Thank you very much.

 Best Regards,
 Kai
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Expanded Ensemble and Gromacs 4.6

2013-02-05 Thread Michael Shirts
Hi, Joakim-

Expanded ensemble is still a bit experimental.  I don't immediately
see any problem that jump right out, but if you go to
http://redmine.gromacs.org/ and file a bug report, including giving
example files that cause the problem, I can take a look at it.

On Tue, Feb 5, 2013 at 6:00 AM, Joakim Jämbeck jamb...@me.com wrote:
 Hi,

 I am trying to use the Expanded Ensemble (EE) method to compute the free 
 energy of solvation of a small organic molecule.

 Basically I am playing around but I cannot get the simulations to run.

 Here are my EE and free energy settings:

 %---
 free-energy = expanded  ; Expanded ensemble
 couple-moltype  = C1X   ; Molecule to introduce
 couple-lambda0  = none  ; Go from no interactions with 
 solvent...
 couple-lambda1  = vdw-q ; ...to full interactions
 couple-intramol = no; Do not decouple internal 
 interactions
 init_lambda_state   = 0 ; Start from the first column of the 
 lambda vectors
 delta-lambda= 0 ; No increments in lambda
 nstdhdl = 100   ; Frequency for writing dH/dlambda
 coul-lambdas= 0.0 0.25 0.5 0.75 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0
 vdw-lambdas = 0.0 0.0 0.0 0.0 0.0 0.2 0.4 0.6 0.7 0.8 0.9 1.0
 dhdl-print-energy   = yes   ; Include total energy in the dhdl 
 file
 sc-alpha= 0.5   ; Soft core alpha parameter
 sc-power= 1 ; Power of lambda in soft core 
 function
 sc-r-power  = 6 ; Power of radial term in the soft 
 core function
 sc-sigma= 0.3   ; Soft core sigma
 sc-coul = no; No soft core for Coulomb
 separate-dhdl-file  = yes   ; Seperate dhdl files
 dhdl-derivatives= yes   ; Print derivatives of the Hamiltonian

 ; expanded ensemble variables
 nstexpanded = 100   ; Number of steps between attempts to 
 change the Hamiltonian
 lmc-stats   = wang-landau   ; WL algorithm to explore state space
 lmc-move= gibbs ; Decides which state to move to
 lmc-weights-equil   = number-steps  ; EE weight updating stops after...
 weight-equil-number-steps = 500 ; ...10ns of simulation

 ; Seed for Monte Carlo in lambda space
 lmc-seed= 1993
 lmc-repeats = 1
 lmc-gibbsdelta  = -1
 lmc-forced-nstart   = 0
 symmetrized-transition-matrix = no
 nst-transition-matrix   = -1
 mininum-var-min = 100
 wl-scale= 0.6
 wl-ratio= 0.8
 init-wl-delta   = 1
 wl-oneovert = yes
 %---

 Besides this I use LINCS to constrain all bonds, sd-integrator and a normal 
 cut-off scheme with a 1 fs time step.

 However, once I try to run the files on the cluster I always end up with 
 LINCS warnings and after a few seconds the program crashes due to too many 
 LINCS warnings. If I increase the time step to 2 fs I run into problems 
 SETTLE for the water. I always start from a nice structure that is taken as 
 the final snap shot from a simple 1 ns MD simulation of my system.

 What could be the problem?

 If I use free_energy = yes instead of EE things work fine.

 Did I perhaps mess up the EE settings or something?

 Thanks in advance.

 Best,
 Joakim

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] meta-dynamics in gromacs-4.6

2013-01-16 Thread Michael Shirts
I assume PLUMED will be implemented for Gromacs 4.6, as many PLUMED
developers use Gromacs.  Perhaps any PLUMED lurkers on the list can
speak up. . . .

On Wed, Jan 16, 2013 at 9:20 AM, Mark Abraham mark.j.abra...@gmail.com wrote:
 The GROMACS team has no plans for that. The usual problem here is that
 everybody would like every algorithm included, but that developers with
 time and experience are scarce :-) It's an open source project though, so
 anyone can do whatever they like. We're prepared to consider inclusions to
 the main project.

 Note that we're doing a major rework of the code base to C++ in the next
 year (or more), so people implementing new features may wish to consider
 that in their choice to write code :-)

 Mark

 On Tue, Jan 15, 2013 at 2:26 PM, James Starlight 
 jmsstarli...@gmail.comwrote:

 Dear Gromacs developers!


 There is well-known plugin plumed which can be used for implementation
 of meta-dynamics simulation if Gromacs-4.5. I wounder to know if it
 possible to include some meta-dynamics options in the new gromacs
 release ( similar to inclusion of essential dynamics sampling in
 previous gromacs versions) ?


 Thanks for attention,
 James
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Re: Is vacuum simulation NVT?

2012-12-16 Thread Michael Shirts
 Could you tell me is there any difference of different Tau_t ussage (
 inverse friction in case of Stochastic dynamics) for simulation of
 water-soluble as well as membrane-proteins ? In the first case I'm
 using tau_t 2ps that is lower than internal water friction. In the
 second case one part of protein could be in the membrane ( e.g
 helixes) but other ( e.g loops) in water media. Both of that solvents
 have different characteristic viscousity. So what Tau_t should be used
 in stochastic dynamics for such biphastic systems?

I don't think there is a right sd tau_t.  You probably want at
tau_t^-1 that is lower (tau_t longer) than the natural viscosity of
either system, such that the viscosity of the actual intermolecular
interactions dominates, not the unphysical viscosity of the
thermostat.  I don't know how much lower is good enough, though.  You
may have to experiment.
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Is vacuum simulation NVT?

2012-12-11 Thread Michael Shirts
 In the absence of PBC, you simply have an infinite system.  In a loose
 sense, that may be NVT, but V is infinite, so whether or not you can
 consider that to be constant or not is theoretical math above what I know :)

A real molecule in vacuum is usually NVE -- it is not coupled to the
environment, and thus must have conserved energy.  You can certainly
add a thermostat, and then it will be NVT, though it won't be very
much like a real isolated gas molecule.  If you try to run NPT, the
simulation will likely crash because of numerical instabilities, and
there's not much of a point, since you are essentially either 1) in
the ideal gas limit if running with no periodic boundary conditions 2)
in some sort of weird superdilute crystal that really doesn't resemble
anything real if run with a periodic boundary conditions

If you are subtracting out the center of mass motion, then V is not
infinite -- you remove the center of mass degree of freedom, and thus
you have a very different ensemble than if you include the center of
mass motion.  You would need to multiply by V to get the partition
function for an actual gas.

Note that I would strongly suggest sd as the integrator/thermostat,
since there are real issues with ergodicity in systems with only a few
degrees of freedom
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Question about conserved energy in MTTK

2012-11-28 Thread Michael Shirts
Hi, all-

I would recommend using Parrinellio-Rahman + Nose-Hoover md + at least
until 4.6.

A random-walk drift in the conserved energy is actually what MTTK
gives -- it's not as conserved as, say, energy conservation, it just
has an expectation value of zero drift over time, which means that the
RMSD will increase with time according to sqrt(dt).

But if you are seeing constant drift, something is wrong.

On Wed, Nov 28, 2012 at 7:20 AM, tarak karmakar tarak20...@gmail.com wrote:
 Hi,

 I was also facing the same problem. If you check your pressure during
 this NPT run, u can see that it got increased to a higher value. I had
 posted the same problem few days back, u can follow the thread. It
 seems MTTK is not stable enough and is not performing well in this
 context. So I have moved to Leap-frog, Nose-Hoover, Parinello-Rahman
 combination for the NPT simulation. There is one paper as well by
 Prof. Shirts in JCTC.

 Cheers,
 Tarak

 On Thu, Nov 22, 2012 at 2:08 PM, Shun Sakuraba sakuraba.s...@jaea.go.jp 
 wrote:
 Dear list,

 I am trying to use MTTK barostat in GROMACS 4.5.5.
 After analyzing the result for a while, I found that the conserved energy 
 (not total energy) of MTTK is drifting during the simulation.
 The .xvg, .edr files are uploaded at [1] and [2]. It is drifting with a 
 constant ratio of ca. -185 kJ/mol/ps.

 I cannot believe this is an expected behavior, so could anyone point out 
 where I am wrong in my simulation setup? I found similar report at [3] but 
 seems it was when 4.5 was in pre-release stage.

 Thanks in advance for your help!

 * Simulation detail
 The system consists of 1000 SPC-E water molecules, and the time step is set 
 to 0.5 fs, just in case the long timestep harms the conserevation (c.f. 
 [3]). The interaction energy is set to switching version, just in case, too. 
 Changing these parameter does not seem to improve the conservation.
 The double precision version of GROMACS is used (single precision version 
 also has the same problem).
 The system has been pre-equilibrated with Berendsen pressure coupling 
 simulation with the same pressure and temperature.

 [1] https://www.dropbox.com/s/16bgfhvavqcw1f8/mttk05.xvg
 [2] https://www.dropbox.com/s/tg8butw39rmz8qk/mttk05.edr
 [3] 
 http://lists.gromacs.org/pipermail/gmx-developers/2010-January/003979.html

 == .mdp file contents follow

 integrator = md-vv
 define =

 dt = 0.0005
 nsteps = 100 ; 500 ps

 coulombtype = PME-Switch
 vdwtype = Switch
 pbc = xyz

 rlist = 1.2
 rcoulomb = 1.0
 rcoulomb_switch = 0.9
 rvdw = 1.0
 rvdw_switch = 0.9
 nstlist = 1

 tinit = 0
 tcoupl = nose-hoover
 tc_grps = System
 tau_t = 0.5
 ref_t = 300.0
 nsttcouple = 1

 pcoupl = MTTK
 pcoupltype = isotropic
 compressibility = 4.5e-5
 ref_p = 1.01325
 tau_p = 0.5
 refcoord_scaling = no
 nstpcouple = 1

 constraints = hbonds
 constraint_algorithm = LINCS

 nstxtcout = 100
 nstlog = 100
 nstenergy = 100
 nstvout = 0
 nstxout = 1000

 --
 Shun SAKURABA, Ph.D.
 Postdoc @ Molecular Modeling  Simulation Group, Japan Atomic Energy Agency
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] pressure_coupling

2012-11-22 Thread Michael Shirts
Hi, all-

There are some issues with MTTK + constraints that are being worked
out for 4.6.  The good thing is, I have developed some sensitive tests
of the correct volume distribution (see
http://arxiv.org/abs/1208.0910) and the errors in PR are very, very
small. I would recommend using md + PR for projects with code before
4.6.

On Thu, Nov 22, 2012 at 4:27 AM, Florian Dommert
domm...@icp.uni-stuttgart.de wrote:
 -Ursprüngliche Nachricht-
 Von: gmx-users-boun...@gromacs.org [mailto:gmx-users-
 boun...@gromacs.org] Im Auftrag von tarak karmakar
 Gesendet: Donnerstag, 22. November 2012 10:15
 An: Discussion list for GROMACS users
 Betreff: Re: [gmx-users] pressure_coupling

 U r right FLorian
 I have also tried playing around the tau_p but in vain.
 Even in absence of any constraints, it is giving almost same result.
 Em thinking to move again to Leap-Frog, NH , PR.  I see people generally
 use this
 combination a lot.
 Thanks

 Tarak

 Yes, that is right. The reason might be, that it is stable, working with
 Leap-Frog and implemented in Gromacs. However actually PR does not produce
 an NpT but an isoenthalpic ensemble. It also does not conform to both
 pressure virial theorems (see appendix of the Nose paper cited in the
 Gromacs manual). For this reason it would be very very good, if MTTK would
 work in Gromacs, because it fulfills all requirements for an NpT ensemble.
 On the other hand side the deviations vanish in the thermodynamic limit, so
 if your system is large enough, there should be no significant difference.

 /Flo


 On Wed, Nov 21, 2012 at 8:08 PM, Florian Dommert domm...@icp.uni-
 stuttgart.de wrote:
 
 
  ---
  Florian Dommert
  Dipl. Phys.
 
  Institut für Computerphysik
  Universität Stuttgart
  Allmandring 3
  D-70569 Stuttgart
 
  Tel.: 0711-68563613
  Fax: 0711-68563658
 
 
  -Ursprüngliche Nachricht-
  Von: gmx-users-boun...@gromacs.org [mailto:gmx-users-
  boun...@gromacs.org] Im Auftrag von tarak karmakar
  Gesendet: Mittwoch, 21. November 2012 15:03
  An: Discussion list for GROMACS users
  Betreff: Re: [gmx-users] pressure_coupling
 
  Thanks for the information Flo.
  Before doing NPT I have already equilibrated my system by heating it
  from
  0K to
  300K in 300 ps, then the pressure has reached to 1 bar. Now while
  doing
  NPT I'm
  getting the excess pressure.
  Is there any problem with the coupling constant ? I am checking it by
  taking
  different tau_p values. Let's see.
 
 
  I don't think that playing around with the coupling constant will help
 you.
  You can set it to extreme values, but you won't see any difference.
  The coupling constant determines, how fast the system pressure should
  relax to the reference pressure. I would see a better possibility to
  play around by simulating for a longer time. Then observing the
  variation of the pressure in time, the size of the fluctuation and the
  excess pressure. Perhaps something will change, but I don't think so.
  I play around with the coupling constants but observed no change.
 
  Maybe, but this is really speculation, there is a problem with the
  combination of constraints and MTTK. Please check the archives of the
  user and developer list to obtain more information.
 
  /Flo
 
 
  On Wed, Nov 21, 2012 at 1:16 AM, Florian Dommert domm...@icp.uni-
  stuttgart.de wrote:
   -Ursprüngliche Nachricht-
   Von: gmx-users-boun...@gromacs.org [mailto:gmx-users-
   boun...@gromacs.org] Im Auftrag von Justin Lemkul
   Gesendet: Dienstag, 20. November 2012 18:33
   An: Discussion list for GROMACS users
   Betreff: Re: [gmx-users] pressure_coupling
  
  
  
   On 11/20/12 12:29 PM, tarak karmakar wrote:
Thanks Justin for the quick reply.
Is there any problem with the algorithms ??
   
I have used Velocity Verlet , Nose-Hoover and MTTK combination.
SHAKE has been used to constrains the H-covalent bonds.
tau_t = 1 ps
tau_P = 1 ps
I got the mean pressure at ~130 bar.
   
Previously with the same initial coordinates I have used
Leap-Frog, NH, Parinello-Rehman with LINCS to constrain H-covalent
 bonds.
tau_t was 0.1 ps
and tau P was 2 ps.
The I have seen the pressure fluctuating around 1 bar( as
expected) So can you please inform me from where this problem is
coming - algorithms and/ tau_t and tau_P parameters ?
   
  
   I have no personal experience with the md-vv/MTTK combination.
   The way to test if there is a bug or something is to take an
   equilibrated system (as
   suggested
   before) and continue it with the desired parameters.  If they
   deviate or incorrectly report pressure, then there's probably a
   bug.  I'm not ready
   to
   conclude that until it is tested though.
  
   -Justin
  
  
   I once tried to use the same combination of T and P coupling and
   MD-vv for a system, which could be simulated with PR at 1bar
   without problems. But I also observed this large pressure. Somehow
   I have in mind, that there was 

Re: [gmx-users] Re: pressure_coupling

2012-11-22 Thread Michael Shirts
It's in review with JCTC right now.

On Thu, Nov 22, 2012 at 2:19 PM, ABEL Stephane 175950
stephane.a...@cea.fr wrote:
 Hello,

 This is a very nice and interesting work, Michael. Thank you for the efforts 
 you made in writing this paper. I hope you will publish it.

 Best

 Stephane


 

 Hi, all-

 There are some issues with MTTK + constraints that are being worked
 out for 4.6.  The good thing is, I have developed some sensitive tests
 of the correct volume distribution (see
 http://arxiv.org/abs/1208.0910) and the errors in PR are very, very
 small. I would recommend using md + PR for projects with code before
 4.6.

 On Thu, Nov 22, 2012 at 4:27 AM, Florian Dommert
 domm...@icp.uni-stuttgart.de wrote:
 -Urspr?ngliche Nachricht-
 Von: gmx-users-boun...@gromacs.org [mailto:gmx-users-
 boun...@gromacs.org] Im Auftrag von tarak karmakar
 Gesendet: Donnerstag, 22. November 2012 10:15
 An: Discussion list for GROMACS users
 Betreff: Re: [gmx-users] pressure_coupling

 U r right FLorian
 I have also tried playing around the tau_p but in vain.
 Even in absence of any constraints, it is giving almost same result.
 Em thinking to move again to Leap-Frog, NH , PR.  I see people generally
 use this
 combination a lot.
 Thanks

 Tarak

 Yes, that is right. The reason might be, that it is stable, working with
 Leap-Frog and implemented in Gromacs. However actually PR does not produce
 an NpT but an isoenthalpic ensemble. It also does not conform to both
 pressure virial theorems (see appendix of the Nose paper cited in the
 Gromacs manual). For this reason it would be very very good, if MTTK would
 work in Gromacs, because it fulfills all requirements for an NpT ensemble.
 On the other hand side the deviations vanish in the thermodynamic limit, so
 if your system is large enough, there should be no significant difference.

 /Flo


 On Wed, Nov 21, 2012 at 8:08 PM, Florian Dommert domm...@icp.uni-
 stuttgart.de wrote:
 
 
  ---
  Florian Dommert
  Dipl. Phys.
 
  Institut f?r Computerphysik
  Universit?t Stuttgart
  Allmandring 3
  D-70569 Stuttgart
 
  Tel.: 0711-68563613
  Fax: 0711-68563658
 
 
  -Urspr?ngliche Nachricht-
  Von: gmx-users-boun...@gromacs.org [mailto:gmx-users-
  boun...@gromacs.org] Im Auftrag von tarak karmakar
  Gesendet: Mittwoch, 21. November 2012 15:03
  An: Discussion list for GROMACS users
  Betreff: Re: [gmx-users] pressure_coupling
 
  Thanks for the information Flo.
  Before doing NPT I have already equilibrated my system by heating it
  from
  0K to
  300K in 300 ps, then the pressure has reached to 1 bar. Now while
  doing
  NPT I'm
  getting the excess pressure.
  Is there any problem with the coupling constant ? I am checking it by
  taking
  different tau_p values. Let's see.
 
 
  I don't think that playing around with the coupling constant will help
 you.
  You can set it to extreme values, but you won't see any difference.
  The coupling constant determines, how fast the system pressure should
  relax to the reference pressure. I would see a better possibility to
  play around by simulating for a longer time. Then observing the
  variation of the pressure in time, the size of the fluctuation and the
  excess pressure. Perhaps something will change, but I don't think so.
  I play around with the coupling constants but observed no change.
 
  Maybe, but this is really speculation, there is a problem with the
  combination of constraints and MTTK. Please check the archives of the
  user and developer list to obtain more information.
 
  /Flo
 
 
  On Wed, Nov 21, 2012 at 1:16 AM, Florian Dommert domm...@icp.uni-
  stuttgart.de wrote:
   -Urspr?ngliche Nachricht-
   Von: gmx-users-boun...@gromacs.org [mailto:gmx-users-
   boun...@gromacs.org] Im Auftrag von Justin Lemkul
   Gesendet: Dienstag, 20. November 2012 18:33
   An: Discussion list for GROMACS users
   Betreff: Re: [gmx-users] pressure_coupling
  
  
  
   On 11/20/12 12:29 PM, tarak karmakar wrote:
Thanks Justin for the quick reply.
Is there any problem with the algorithms ??
   
I have used Velocity Verlet , Nose-Hoover and MTTK combination.
SHAKE has been used to constrains the H-covalent bonds.
tau_t = 1 ps
tau_P = 1 ps
I got the mean pressure at ~130 bar.
   
Previously with the same initial coordinates I have used
Leap-Frog, NH, Parinello-Rehman with LINCS to constrain H-covalent
 bonds.
tau_t was 0.1 ps
and tau P was 2 ps.
The I have seen the pressure fluctuating around 1 bar( as
expected) So can you please inform me from where this problem is
coming - algorithms and/ tau_t and tau_P parameters ?
   
  
   I have no personal experience with the md-vv/MTTK combination.
   The way to test if there is a bug or something is to take an
   equilibrated system (as
   suggested
   before) and continue it with the desired parameters.  If they
   deviate or incorrectly report pressure, then there's 

[gmx-users] Experiences with Gromacs scaling on US supercomputer centers?

2012-09-26 Thread Michael Shirts
Hi all,

I'd be interested to know about people's experiences with Gromacs on
US national computing centers.  Which machines have it set up to scale
the
best?  We're putting in an XSEDE request soon, and I'm trying to
figure out which resource to request.  Our system is
semi-coarse-grained, using
reaction field electrostatics, so we don't have to worry about PME.
Of course, couldn't hurt to know about PME scaling as well.   I'm
interesting
scalings with 100K - 300K atoms.

Of course, best performance will probably change with 4.6 because of
all the setup tweaks, but let's start with 4.5 scaling info!

Best,

Michael Shirts
Assistant Professor
Department of Chemical Engineering
University of Virginia
michael.shi...@virginia.edu
(434)-243-1821
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Fast exchanges for REMD

2012-09-26 Thread Michael Shirts
 However, the time value (4 in this example) is limited to 6 digits.

Sounds like this should be increased?  There's a pending change to
replica exchange, so this could be added to 4.6 without disrupting the
release timing.

On Wed, Sep 26, 2012 at 11:22 AM, Andreas Zink zinka...@googlemail.com wrote:
 Dear all,

 I could finally demux my REMD trajectories with high EAF. They look fine,
 but I'm not 100% sure about it.

 Unfortunately, there seems to be a bug in mdrun. As you migh know, the log
 files contain the exchange attempts like:

 Replica exchange at step 2000 time 4
 Repl ex  0123 x  45
 Repl pr.00   1.0

 However, the time value (4 in this example) is limited to 6 digits. It's not
 really a bug because I think GROMACS recommends EAFs of max. 1/ps. But if
 you exchange every 0.5 ps (EAF= 2/ps) you can run the simulation for max.
 9 ps only. Otherwise your log will simply chop the time after the
 decimal point, e.g.

 Replica exchange at step ... time 9.5
 ...

 Replica exchange at step ... time 10
 ...

 Replica exchange at step ... time 10
 ...

 Replica exchange at step ... time 11
 ...

 I have written a Perl script which fixes the time values, based on the given
 number of steps and stepsize. Additionally the demux.pl script was changed
 because it also chops after 6 digits. I simply changed:

 printf(NDX %-20g,$t);

 in the pr_order and pr_revorder subroutines to:

 printf(NDX %-20.2f,$t);

 and it works just fine.

 Like I said the trajectories look fine, but I'm not really sure if it's
 actually correct that way. I would be happy if anyone would like to discuss
 this :)


 Cheers,

 Andi




 Am 24.09.2012 08:49, schrieb Andreas Zink:

 Dear all,

 I've done some REMD simulations using a quite high exchange attempt
 frequency (10 attempts per ps) as proposed by Sindhikara et al. (Exchange
 Often and Properly in Replica Exchange Molecular Dynamics,J. Chem. Theory
 Comput. 2010, 6, 2804–2808 ).
 Unfortunately, I have now recognized that the demux perl script cannot
 account for an EAF which is higher than the saving frequency in the
 trajectory.

 Comments from demux.pl:
 # If your exchange was every N ps and you saved every M ps you can make
 for
 # the missing frames by setting extra to (N/M - 1). If N/M is not integer,
 # you're out of luck and you will not be able to demux your trajectories
 at all.

 In my case exchanges every 0.1 ps and saved every 5 ps

 Changing the demux.pl script, so that it writes the replica_index.xvg
 with a higher precision (time in ps) should be no problem. However, my
 question is: will this work together with trjcat? Does trjcat search for the
 timeframe given in the first column of replica_index.xvg, or does it links
 each line to one saved timeframe? If so, could I just delete the additional
 lines in replica_index.xvg?

 I would be really happy if someone could help me with this!
 Thanks!

 Andi


 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: v-rescale

2012-09-20 Thread Michael Shirts
I've done some extensive testing (paper on testing method in the
works) and vrescale gives a very accurate ensemble very well for NVT.
Parrinello-Rahman and MTTK are the only algorithms that are correct
for NPT.  Berendsen barostat is not.  Note that there is a bug with
vrescale + md-vv + that is fixed in 4.6 and (hopefully) for the next
patch of 4.5.

For equilibrating a system, Berendsen for both temperature and
pressure is the best bet.  They artificially minimize fluctuations,
which is great for equilibration, bad for data collection.

On Thu, Sep 20, 2012 at 6:05 PM, Mark Abraham mark.abra...@anu.edu.au wrote:
 On 20/09/2012 9:35 AM, Peter C. Lai wrote:

 I am not sure where the idea of using berendsen barostat with the
 v-rescale
 thermostat for equilibration came from, however. Doesn't the typical
 equilibration begin with v-rescale for temperature equilibration then
 adding parinello-rahman barostat then switching to nose-hoover for
 production
 runs (as nose-hoover chains result in the correct canonical distribution)?


 N-H does have known issues, see
 http://link.aip.org/link/doi/10.1063/1.2989802 and
 http://link.aip.org/link/doi/10.1063/1.2408420. I am not aware of any
 shortcomings of the Bussi v-rescale thermostat in GROMACS.

 Mark



 On 2012-09-19 04:24:27PM -0700, Ladasky wrote:

 Dear Sara,

 I just had a problem with my simulations that I traced to the use of the
 V-rescale temperature algorithm.  Here is my recent post:


 http://gromacs.5086.n6.nabble.com/Re-Water-molecules-cannot-be-settled-why-td4999302.html;cid=1348087067061-71#a5001121

 V-rescale may be appropriate in certain simulations, but it is apparently
 NOT appropriate when used with Berendsen pressure coupling during the
 initial equilibration.  I don't know if that is related to your problem,
 but
 it's something that I just discovered the hard way.



 --
 View this message in context:
 http://gromacs.5086.n6.nabble.com/v-rescale-tp5001066p5001122.html
 Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] BAR / g_bar problems

2012-08-24 Thread Michael Shirts
Hi, David-

Perhaps we can work with you to compare the internal g_bar
implementation with our external BAR and MBAR implementation run from
the .dhdl.xvg data output in 4.6.  It would probably be easier to run
these calculations after the optimization updates are added in the
next few days(?), but we can plan now.  Note that withe errors this
big, 100-200 ps should be enough to see what's going on, so it can be
done rapidly.  Drop me a line off the list to figure out details?

Best,
Michael

On Fri, Aug 24, 2012 at 9:22 AM, David van der Spoel
sp...@xray.bmc.uu.se wrote:
 Hi,

 we have terrible problems to get reasonable results from BAR free energy
 calculations. We basically follow Justin's tutorial for solvation of
 methanol, ethanol, diethyl-ether and neopentane in water using OPLS but get
 very strange results. The free energy curves are here:

 http://folding.bmc.uu.se/images/koko.jpg

 Molecule  Energy Exper
 Methanol -8-21
 Ethanol  -9-21
 Diethylether -11-7
 Neopentane   -10   +11

 any clue?
 --
 David van der Spoel, Ph.D., Professor of Biology
 Dept. of Cell  Molec. Biol., Uppsala University.
 Box 596, 75124 Uppsala, Sweden. Phone:  +46184714205.
 sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Only plain text messages are allowed!
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Only plain text messages are allowed!
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] When to use Dispersion Correction for Lipid Bilayers

2012-08-19 Thread Michael Shirts
Short answer -- use a dispersion correction if the force field was
parameterized to include one.  Make sure you use the cutoff that the
lipid was parameterized for as well.

Long answer -- neither on nor off is correct for a lipid bilayer.  A
dispersion corrections corrects for the fact that you are neglecting
the r^-6 term at long range.  However, it assumes an isotropic
distribution outside the cutoff.  Lipid bilayers have long range
order, so a standard dispersion correction is inappropriate.  The
lipid properties will be cutoff dependent, which is not a very good
thing.  See  PJ in't Veld, AE Ismail, GS Grest, J. Chem. Phys. 127,
144711 (2007) for a solution, using Ewald summation for the
Lennard-Jones terms.  This is in the works for Gromacs (and has been
for a while), but might still be a while because it's a little lower
on the to do list.  Once methods like this are in, it will be possible
to parameterize lipids that behave correctly.

Best,
Michael

On Sun, Aug 19, 2012 at 3:16 PM, David Ackerman da...@cornell.edu wrote:
 Hello,

 I was wondering when it is appropriate to use Dispersion Correction
 for lipid bilayers, or which setting (no, EnerPres, or Ener) is best.
 I have seen it used in some discussions, and not used in others. As
 for myself, I have simulated a DPPC bilayer. With DispCorr=EnerPres, I
 get an area per lipid of ~.615 nm^2 (whereas other reported values are
 closer to .64 nm^2) and slower diffusion than is reported. However,
 when I don't have a DispCorr term, my area per lipid becomes ~.656
 nm^2 and my diffusion more closely matches other reported values.

 So, should DispCorr be used at all for bilayer simulations, and if so,
 which type is most appropriate?

 Thank you for your help,
 David
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Only plain text messages are allowed!
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Only plain text messages are allowed!
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] lincs with mttk

2012-07-27 Thread Michael Shirts
Good question.  Short answer, no -- LINCS doesn't play well with a
velocity verlet based pressure control algorithm.

Long answer: MTTK has ended up not being robust because you need to
solve a self consistent set of equations every timestep using the
pressure estimator, which is extremely noisy, so it crashes randomly
for no good reason.

The good news is that I've done some fairly extensive testing, and
leapfrog MD + Parrinello-Rahman + Nose-Hoover or vrescale is very
close to the correct distribution -- there's really not that much
difference for any practical purpose (I have a preprint with more
details).  The volume fluctuations are almost exactly right. So I
think that for anything that doesn't require getting free energies to
kJ/mol accuracy you are OK with that combination (just don't use
Berendsen anything).

Longer term (5.0), MTTK + shake will be removed -- its just too
numerically unstable, and the self-consistent iteration part means
it's very hard to parallelize.  MTTK will be left in, but only without
constraints.  Hopefully, for 5.0 we will also have RESPA-like
functionality, so that the bondeds can be performed at a higher
frequency than the nonbondeds, which will be another useful way to get
longer times steps. We also have planned to implement a Monte Carlo
barostat, which will give exactly the correct NPT distribution for any
integrator.

Best,
Michael

On Fri, Jul 27, 2012 at 2:01 PM, Katie Maerzke maerz...@gmail.com wrote:
 Hi all -

 Is there any plan to get LINCS working with the MTTK barostat in the near 
 future?

 Thanks
 Katie

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Only plain text messages are allowed!
 * Please search the archive at 
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Only plain text messages are allowed!
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] BAR gives different result than TI

2012-07-07 Thread Michael Shirts
The implementation of BAR in gromacs is pretty hard for me to follow
because of how everything is stored noncompactly in the histogram.  In
4.6, both can be computed from the same dhdl.xvg file, so it might be
easier to track down possible bugs.

On Fri, Jun 29, 2012 at 2:24 PM, David van der Spoel
sp...@xray.bmc.uu.se wrote:
 Hi,

 we've been trying to do free energy calculations for solvation of two small
 molecules in water, n-butylamine (NBA) and diethyl-ether (DEE). For one of
 them the result with BAR (using Justin's tutorial) is significantly
 different from TI:
   BAR TI  Exper. (kJ/mol)
 NBA   -11.1   -10.8   -17.9
 DEE   -11.0   -1.8 -7.4

 Since for NBA the results are pretty close it seems that the protocol should
 be right, and the numbers seem to have converged pretty good as well.

 Any clues why these results could be so different?

 --
 David van der Spoel, Ph.D., Professor of Biology
 Dept. of Cell  Molec. Biol., Uppsala University.
 Box 596, 75124 Uppsala, Sweden. Phone:  +46184714205.
 sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se


 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 * Only plain text messages are allowed!
 * Please search the archive at
 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-requ...@gromacs.org.
 * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Only plain text messages are allowed!
* Please search the archive at 
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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] What does --enable-fahcore mean?

2012-06-26 Thread Michael Shirts
It only matters for running on Folding@Home.  For other users of
gromacs, it doesn't do anything.


On Tue, Jun 26, 2012 at 3:50 PM, Bao Kai paeanb...@gmail.com wrote:
 Hi, all,

 I am wondering what the --enable-fahcore option of configure means.  I
 got the explanation from configure --help of create a library with
 mdrun functionality, while it is not very clear to me.

 Best Regards,
 Kai
 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at 
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Configurational bias Monte Carlo

2012-06-01 Thread Michael Shirts
There's a fair amount of interest for more general support for Monte
Carlo methods in GROMACS 5.0.  However, there is no any current
support for configuration bias Monte Carlo (or any other kind of MC)
currently in the code.

On Fri, Jun 1, 2012 at 12:10 PM, Benjamin Haley bha...@purdue.edu wrote:
 Hello,

   I have found the excellent documentation for creating polymer
 chains in GROMACS.  I want to create several chains and pack
 them into a volume, adjusting torsion angles to minimize
 interactions with other chains (and self-interaction within a chain).
 One method for doing this is the configurational bias Monte Carlo
 sampling.  Can this (or something similar) be done in GROMACS?

 Thank you!
 Ben Haley
 Purdue University
 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] REMD question

2012-05-28 Thread Michael Shirts
Gromacs already supports replica exchange -- what particularly are you
implementing?

Equilibration of pressure is always a good idea -- even if you are
running NVT simulations, you want to get them to be at the equilibrium
volume for your system and temperature choice, which will require
equilibration at constant pressure.

On Mon, May 28, 2012 at 4:37 PM, Nathalia Garces natsgar...@gmail.com wrote:
 Dear Gromacs Users,



 We are implementing REMD method in Gromacs in protein folding, in your web
 page you give some steps that don´t mention any step about NPT
 stabilization.  This step is necessary to run REMD simulations?



 Thank you in advance,



 Nathalia




 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] poor performance in Hemiltonian Replica Exchange

2012-05-10 Thread Michael Shirts
4.6 will include Hamiltonian replia exchange functionality built into
the MPI version.

Currently the description of the error is very vague -- if you can
write up what exactly the numbers are, and what they should be, with
files that exactly replicate the error, then I can take a look.  But
unless I can reproduce the error you are describing out of the box,
its unlikely I will be able to find it.

Additionally, it would be easiest if the files were deposited as a
redmine bug report, so that the information is centrally located.

Best,
Michael

On Thu, May 10, 2012 at 12:23 PM, francesco oteri
francesco.ot...@gmail.com wrote:
 Dear gromacs users,

 I performed a Hemiltonian Replica Exchange (i.e. replica exchange where each
 replica has a init_lambda=0, delta_lambda=0 and init_lambda ranging
 uniformely from 0 to 1).

 Since I have only ten fixed discrete lambda, I run a Temperature Replica
 Exchange where, for each replica I generated a .top file with the parameter
 rescaled through a
 python script ( in practice I did through python the same thing gromacs is
 supposed to do with the H-REM previously described). Now gromacs complained
 because
 every replica has the same setup, so I changed the temperatures using very
 close values (300.0001K,
  300.0002K,300.0003K,300.0004K,300.0005K,300.0006K,300.0007K,300.0008K, 300.0009K)
 With this setup the simulation runs fine and I expect to have similar
 result.

 Then I compared the results observing two phenomena:

 1) In the second case exchange rate is 100%, while in the first case I have
 an exchange rate close to 30%.
 Does it rise  because the temperatures are too close?

 2) The second setup is 3x faster!
 In particular I observe an imbalance between PME and force calculation
 ranging from 10% to 60%.
 I tried to run each replia indipendently (a different mdrun instance for
 each .tpr file) but still I observe the same performance slowdown.
 I guess the free energy impairs the efficient force calculation, but I dont
 understand why.

 Can someone explain me the two observations?



 Francesco

 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] BUG: Free energy calculation

2012-03-23 Thread Michael Shirts
Hi, Sabine-

If you can go to http://redmine.gromacs.org/projects/gromacs/issues
and file a bug report (including attaching files), I can look at it.
If that's ends up not working well, you can send me the files off of
the list, but it's usually better to have things in the redmine system
so problems are documented.  I'm also working on the updates to free
energy code in 4.6, so I want to make sure this will be solved there
as well.

Best,
Michael

On Fri, Mar 23, 2012 at 7:16 AM, Sabine Reisser sabine.reis...@kit.edu wrote:
 Hi Mark,

 with FE, without PR : same error
 without FE, with PR: stable
 without FE, without PR: stable

 I've never had this error before.

 Logfile says:
 [...]
 Initializing Domain Decomposition on 2 nodes
 Dynamic load balancing: no
 Will sort the charge groups at every domain (re)decomposition
 Initial maximum inter charge-group distances:
    two-body bonded interactions: 3.341 nm, LJC Pairs NB, atoms 22476 22761
  multi-body bonded interactions: 0.649 nm, Angle, atoms 1668 1686
 Minimum cell size due to bonded interactions: 3.675 nm
 Using 0 separate PME nodes
 Optimizing the DD grid for 2 cells with a minimum initial size of 3.675 nm
 The maximum allowed number of cells is: X 1 Y 1 Z 2
 Domain decomposition grid 1 x 1 x 2, separate PME nodes 0
 PME domain decomposition: 2 x 1 x 1
 Domain decomposition nodeid 0, coordinates 0 0 0
 [...]
 Linking all bonded interactions to atoms
 There are 55376 inter charge-group exclusions,
 will use an extra communication step for exclusion forces for PME

 The initial number of communication pulses is: Z 1
 The initial domain decomposition cell size is: Z 5.09 nm

 The maximum allowed distance for charge groups involved in interactions is:
                 non-bonded interactions           1.200 nm
            two-body bonded interactions  (-rdd)   4.213 nm
          multi-body bonded interactions  (-rdd)   4.213 nm





 There it stops. In the 1-thread case, the second part is replaced by
 Initiating Steepest Descents and then writing out the energies for every
 step.

 Gromacs version is 4.5.5.

 I also attached the whole logfile.

 Cheers
 Sabine





 On 03/23/2012 11:52 AM, Mark Abraham wrote:

 On 23/03/2012 9:17 PM, Sabine Reisser wrote:


 Dear gromacs users/developers,

 when trying to couple in a peptide into a membrane with:

 ; Define position restraints for peptide
 define          = -DPOSRES

 ; couple in peptide
 free_energy     = yes
 init_lambda     = 0.05
 sc_alpha        = 0.7
 sc_power        = 1
 couple-moltype  = Protein
 couple-lambda0  = none
 couple-lambda1  = vdw-q


 grompp works fine, but mdrun (2 threads) gives me

 Making 1D domain decomposition 1 x 1 x 2
 *** glibc detected *** mdrun: realloc(): invalid next size:
 0x7f0f30305810 ***

 and breaks up.


 When running mdrun -nt 1  on only one thread, it works fine.

 Is this a known bug?


 First, is it likely not to be a problem with your setup... is your
 system stable in parallel without FE code? Without position restraints?
 What does your .log file say? What GROMACS version is it?

 Mark




 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] BUG: Free energy calculation

2012-03-23 Thread Michael Shirts
Sabine, thanks for filing it in redmine!  Having a record helps a lot.
 Can you also attach all your input files to the redmine filing?  It
can only really be debugged if the input files you used are included.

Best,
Michael
On Fri, Mar 23, 2012 at 9:11 AM, Justin A. Lemkul jalem...@vt.edu wrote:


 Sabine Reisser wrote:

 Hi Mark,

 with FE, without PR : same error
 without FE, with PR: stable
 without FE, without PR: stable

 I've never had this error before.

 Logfile says:
 [...]
 Initializing Domain Decomposition on 2 nodes
 Dynamic load balancing: no
 Will sort the charge groups at every domain (re)decomposition
 Initial maximum inter charge-group distances:
    two-body bonded interactions: 3.341 nm, LJC Pairs NB, atoms 22476 22761
  multi-body bonded interactions: 0.649 nm, Angle, atoms 1668 1686
 Minimum cell size due to bonded interactions: 3.675 nm
 Using 0 separate PME nodes
 Optimizing the DD grid for 2 cells with a minimum initial size of 3.675 nm
 The maximum allowed number of cells is: X 1 Y 1 Z 2
 Domain decomposition grid 1 x 1 x 2, separate PME nodes 0
 PME domain decomposition: 2 x 1 x 1
 Domain decomposition nodeid 0, coordinates 0 0 0
 [...]
 Linking all bonded interactions to atoms
 There are 55376 inter charge-group exclusions,
 will use an extra communication step for exclusion forces for PME

 The initial number of communication pulses is: Z 1
 The initial domain decomposition cell size is: Z 5.09 nm

 The maximum allowed distance for charge groups involved in interactions
 is:
                 non-bonded interactions           1.200 nm
            two-body bonded interactions  (-rdd)   4.213 nm
          multi-body bonded interactions  (-rdd)   4.213 nm



 These quantities look very weird to me.  They indicate interactions that are
 very far apart are influencing one another.  Can you provide a complete .mdp
 file?  It seems like some aspect of the free energy settings (perhaps
 couple-intramol?) and DD aren't getting along.  The other possibility is to
 try particle decomposition instead of DD (i.e. mdrun -pd).

 -Justin

 --
 

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

 

 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Exchange interval in REMD

2012-01-26 Thread Michael Shirts
 It also depends whether or not you use constant pressure, in which case it
 makes sense to increase the time to long enough to let the volume relax.
 I still do not understand why people do NVT REMD, because it makes all but
 one replica have a pressure that is not the ambient pressure.

If all you care about is speedup at a single temperature, then this
fact doesn't really matter all that much.

Note that if you have an incorrect barostat (for example, Berendsen),
then NPT REMD is likely to produce some very bad artifacts.
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] amber99sb in GROMACS vs amber99sb in AMBER

2011-11-21 Thread Michael Shirts
 Is anyone aware of any benchmark analysis about the implementation of the
 amber99sb force field (or any of its variants: 99sb-ildn, 99sb-nmr) in
 GROMACS and AMBER. I am interested to know in what extend the energies
 correlate and if the results agree with experimental data.

Whether the results compare agree with experimental data is irrelevant
to the correctness of the implementation -- that has to do with the
validity of AMBER99sb to begin with.  The only question is, do GROMACS
and AMBER give the same energies for the same configurations?

I actually do not know the answer, and would be interested to hear if
it's been tested.  http://ffamber.cnsm.csulb.edu/ does not appear to
have the very comprehensive tests that appear for earlier models at
the current time.  I SUSPECT there should not be a problem, since
AMBER99sb just changed a couple of torsions, so errors are unlikely
(though in theory possible).
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] How DispCorr influsnces pressure

2011-11-18 Thread Michael Shirts
 Dispcorr is only for homogeneous liquids. It should not be used for
 membranes.

More precisely -- the dispersion correction is an analytical
correction to both the pressure and energy that is rigorously correct
in the limit of the radial distribution function g(r)=1 outside the
van der Waals cutoff.  For a lipid bilayer, this not a good
approximation, since the system is not isotropic.

Since there is not right way to do membrane systems currently, the
best bet is to go with the same cutoff treatment the lipids were
parameterized with.

DispCorr = Ener is probably still a good idea, since g(r)=1 isn't
perfect, but it's better than completely ignoring the long range
energy.  DispCorr = Ener will not change the configurations sampled,
so any phase properties, pressures, etc. will be unaffected, but the
results will be somewhat less cutoff dependent.

See:
http://pubs.acs.org/doi/abs/10.1021/jp0735987
for some additional information.

Members of the Gromacs team are working behind the scenes for various
possible solutions to the problem for nonisotropic systems.  Don't
expect anything until 5.0, though.

Best,
Michael
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] FEP

2011-10-10 Thread Michael Shirts
FEP is a poorly defined term.  It can either mean 1) making small
changes 'alchemical' changes in the molecules and computing the free
energies by any method (BAR, TI, etc), or 2) it can mean specifically
computing the free energies by exponentially averaging the potential
energy differences.  Basically, using the exponential averaging
formula is a bad idea -- if you have code that only uses that method,
you can get decent results out if you are careful, but TI, BAR, and
MBAR are all better choices.

On Mon, Oct 10, 2011 at 10:58 AM, Justin A. Lemkul jalem...@vt.edu wrote:


 Steven Neumann wrote:

 Thank you guys! So, is there any tutorial in Gromacs for calculating free
 energy of ligand binding using FEP?


 TI or BAR are better methods for calculating binding free energies, I would
 think.  FEP is best for mutating between different species.

 -Justin

 Steven

 On Mon, Oct 10, 2011 at 2:02 PM, Justin A. Lemkul jalem...@vt.edu
 mailto:jalem...@vt.edu wrote:



    mohsen ramezanpour wrote:

        Hi
        Please have a look at Dr.Justin tutorial page at the following
 link:


  http://www.bevanlab.biochem.__vt.edu/Pages/Personal/justin/__gmx-tutorials/index.html

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


    This tutorial is not for FEP explicitly, but may be of some use.
     There is discussion on using the BAR algorithm for binding free
    energy calculations.

    -Justin

        Cheers


        On Mon, Oct 10, 2011 at 12:27 PM, Steven Neumann
        s.neuman...@gmail.com mailto:s.neuman...@gmail.com
        mailto:s.neuman...@gmail.com mailto:s.neuman...@gmail.com__
        wrote:

           Hi Gmx Users,
                Can you suggest some reading and some tutorial in
        calculations of
           binding free energy (ligand binding) in Gromacs? ?I want to
        use Free
           Energy Perturbation method.
                Steven
           --
           gmx-users mailing list    gmx-users@gromacs.org
        mailto:gmx-users@gromacs.org
           mailto:gmx-users@gromacs.org mailto:gmx-users@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-requ...@gromacs.org
        mailto:gmx-users-requ...@gromacs.org
           mailto:gmx-users-request@__gromacs.org
        mailto:gmx-users-requ...@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. Candidate
    ICTAS Doctoral Scholar
    MILES-IGERT Trainee
    Department of Biochemistry
    Virginia Tech
    Blacksburg, VA
    jalemkul[at]vt.edu http://vt.edu/ | (540) 231-9080
    tel:%28540%29%20231-9080
    http://www.bevanlab.biochem.__vt.edu/Pages/Personal/justin
    http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

    ==__==

    --     gmx-users mailing list    gmx-users@gromacs.org
    mailto:gmx-users@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-requ...@gromacs.org
    mailto:gmx-users-requ...@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. Candidate
 ICTAS Doctoral Scholar
 MILES-IGERT Trainee
 Department of Biochemistry
 Virginia Tech
 Blacksburg, VA
 jalemkul[at]vt.edu | (540) 231-9080
 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

 
 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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 

Re: [gmx-users] FEP

2011-10-10 Thread Michael Shirts
 So BAR is only
 referring to the mathematical code used to calculate the overall free
 energy for the FEP, correct?

Yes.  The information required is the same, with the exception that
exponential averaging requires energy differences from only one state,
whereas BAR requires energy differences from both directions.  Since
you are likely running simulations at a range of states anyway,
there's no extra cost.

 My question is, for a mutation of a -CH3 group to a -H group, is it
 better to simply run:
 [+ from (Lambda=0 ,  R-CH3, full charges and interactions -STATE A)
 -- (Lambda=1, R-CH, full charges and interactions -STATE B)]

 OR

 [1) from (Lambda=0 ,  R-CH3, STATE A : Charges and LJ Interactions: ON)
 2)  (Lambda=?, -CH3 Charges: OFF ,LJ Interactions: ON and unmutated)
 3)  (Lambda=?, R-CH3, -CH3 Charges: OFF ,LJ Interactions: OFF)
 4)  (Lambda=?, R-CH3, -CH3 Charges: OFF ,LJ interactions: ON and Mutated to 
 -H)
 5)  (Lambda=1, R-CH3, -CH3  STATE B : Charges and LJ Interactions: ON)

Currently, you'd need to run three simulations to do this (will likely
be simplified somewhat in 4.6).

Set 1: turn R-CH3 charges off in a way that preserves the total charge.
Set 2: change CH3 LJ to H LJ
Set 3: Turn R-H charges on in a way that preserves the total charge.

Note that the first and last transformations will need only one two
states -- their energies will be VERY small.

 Reason I'm asking is because when I try the first choice to do it
 STATE A to STATE B in one step, when I reach Lambda=0.85 and above on
 the NVT equilibration right after EM, I receive errors saying that
 bonds are moving way to far off their constraints which leads me to
 believe that the system is moving too far from where it was energy
 minimized.  Errors such as:

 Step 188, time 0.376 (ps)
 LINCS WARNING
 relative constraint deviation after LINCS:
 rms 0.17, max 0.000636 (between atoms 9 and 68)
 bonds that rotated more than 30 degrees:
  atom 1 atom 2  angle  previous, current, constraint length
      9     68   31.2    0.   0.1110      0.

 Step 188, time 0.376 (ps)  LINCS WARNING
 relative constraint deviation after LINCS:
 rms 0.15, max 0.000531 (between atoms 9 and 68)
 bonds that rotated more than 30 degrees:
  atom 1 atom 2  angle  previous, current, constraint length
      9     68   31.0    0.   0.1110      0.

You probably have problems with the charge + LJ terms on the hydrogen.
 Note that you could possibly use soft core for the Coloumb
interactions for hydrogen -- it's sometimes causes sampling problems,
but for H's it probably doesn't matter.  I haven't done this approach
much so, so I can't tell for sure if it will work.
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Question about Justin's Free Energy Tutorial

2011-10-07 Thread Michael Shirts
I don't think this is a question about new free energy code -- I think
this is asking about the fact if you can do a free energy calculation
by specifying the A and B variables in the topology, instead of using
the MDP coupl-moltype arguments.  This is actually the way free energy
calculations were managed before 4.0, and is still supported.

The answer is yes, though it requires some experience with what's
going on -- and I don't have time to document it all quite right now,
as it's described relatively thoroughly in the manual.

The one difference is that the coupl-moltype keywords makes it easy to
decouple molecules (turn off the intermolecular but not intramolecular
energies) from their environment, which is very hard (and sometimes
impossible) to do by changing A and B states. If you want to turn off
all interactions, then you can do it very straightforwardly by
changing the top files.

Best,
Michael

On Thu, Oct 6, 2011 at 2:42 PM, Justin A. Lemkul jalem...@vt.edu wrote:


 Fabian Casteblanco wrote:

 Hello Justin,

 I have a question about your tutorial.  If I want to mutate one small
 group of a molecule, I would have to not provide 'couple_lambda0' and
 'couple_lambda1', correct?  I would essentially have to follow sec
 5.7.4 in the Gromacs manual and I have to actually provide all state A
 variable and all state B variables.  Gromacs would calculate the new B
 state parameters for bond lengths, angles, etc, correct?  Are there
 any other major differences to account for?


 I have never attempted FEP with the new free energy code.  Try a simple test
 system first and make sure it works as expected.

 -Justin

 --
 

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

 
 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Recommended parameters for NVE simulation of SPCE water

2011-08-17 Thread Michael Shirts
Hmm.  Even now, I'm noticing problems with what I sent:

It should be;

 rcoulomb                    = 1.3.

For PME, rlist should equal rcoul.  PME-switch can improve energy
conservation, but the wrong PME-switch parameters can affect the
results too much.  PME w/o switch should be appropriate for most
'standard' usage.  In the next month or so (by 4.6 at least), the
manual will include some more guidance on accuracy in cutoffs.

Summarizing again:
 For pretty good energy conservation, I would suggest:

 rlist                           = 1.3
 coulombtype              = PME
 rcoulomb                    = 1.3
 vdw-type                    = Switch
 rvdw-switch                = 1.0
 rvdw                          = 1.1

 This should work quite well -- you might get some drift after 1-2 ns,
 but not much.  I'm working on developing suggested PME parameters
 right now for highly quantitative work, but it's not quite ready yet.


 
 Michael Shirts
 Assistant Professor
 Department of Chemical Engineering
 University of Virginia
 michael.shi...@virginia.edu
 (434)-243-1821

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Hamiltonian replica exchange?

2011-08-16 Thread Michael Shirts
Hamiltonian replica exchange is planned for 4.6, and is being beta
tested by some users.


On Tue, Aug 16, 2011 at 2:39 PM, Sanku M msank...@yahoo.com wrote:
 Hi,
    I was wondering whether hamiltonian replica exchange simulation has been
 implemented in latest version of  gromacs . Or, is there any other way of
 performing the hamiltonian replica exchange using some variants of
 lambda-dynamics ?
 Sanku
 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Recommended parameters for NVE simulation of SPCE water

2011-08-10 Thread Michael Shirts
 1.  NOTE 1 above suggests that I use vdwtype = Shift.  When I do this, do
 you recommend that I apply long range dispersion corrections for both energy
 and pressure, using DispCorr = EnerPres, or for only energy, using DispCorr
 = Ener?  Typically, for various (non-NVE) calculations, I have been using
 DispCorr = no, but I am not sure if this is a good idea.  Pages 97-98 of the
 Gromacs 4.5.4 manual seem to suggest that the energy correction due to
 DispCorr is small and usually only significant for free energy calculations
 (which I will not be doing here).  As a rule of thumb, do you typically turn
 dispersion corrections off?

For constant pressure simulations, or for reaching the constant
pressure equilibrium simulation, you should definitely include a
dispersion correction -- the density will be too large, and will be
cutoff dependent.

For constant volume simulations, the dispersion correction will be
constant.  It will thus NOT affect energy conservation, but WILL
affect average potential energy and average total energy,
significantly.

 2.  NOTE 2 above suggests that I use either coulombtype = PME-Switch or
 coulombtype = Reaction-Field-zero.  Do you have any advice or
 recommendation?

For pretty good energy conservation, I would suggest:

rlist   = 1.3
coulombtype  = PME
rcoulomb= 1.1
vdw-type= Switch
rvdw-switch= 1.0
rvdw  = 1.1

This should work quite well -- you might get some drift after 1-2 ns,
but not much.  I'm working on developing suggested PME parameters
right now for highly quantitative work, but it's not quite ready yet.



Michael Shirts
Assistant Professor
Department of Chemical Engineering
University of Virginia
michael.shi...@virginia.edu
(434)-243-1821
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] autocorrelation time of dVpot/dlambda?

2011-03-11 Thread Michael Shirts
 I am doing free energy calculation in Gromacs and want to get an error
 estimate of my results. Is it possible to compute the autocorrelation time
 of dVpot/dlambda in Gromacs using a certain length of trajectory such as 10
 ns? Thanks a lot

The amount of simulation time required to compute the autocorrelation
time of dVpot/dlambda will depend on the timescale of the system.  You
should be able to use g_analyze to analyze the energy.xvg output, and
get an autocorrelation time.  If your simulation time is 50x longer
than the autocorrelation time, then you are very likely converged.
However, if there are very long correlation times, you still might not
have captured all of the time dependent variations in  Vpot/dlambda.
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Possible free energy bug?

2011-03-10 Thread Michael Shirts
Hi, all-

Have you tried running

constraints = hbonds?

That might eliminate some of the constraint issues.  Much less likely
for LINCS to break or have DD issues if only the hbonds are
constrained.  2 fs is not that big a deal for the heteroatom bonds.

Best,
Michael

On Thu, Mar 10, 2011 at 8:04 PM, Justin A. Lemkul jalem...@vt.edu wrote:

 Hi Matt,

 Thanks for the extensive explanation and tips.  I'll work through things and
 report back.  It will take a while to get things going through (unless one
 of the early solutions works!) since I have no admin access to install new
 compilers, libraries, etc. and for some reason the only thing I can ever get
 to work in my home directory is Gromacs itself.  The joys of an aging
 cluster.

 We recently got access to gcc-4.4.5 on Linux, but we're stuck with 3.3 on OS
 X, so there's at least a bit of hope for one partition.

 Thanks again.

 -Justin

 Matthew Zwier wrote:

 Hi Justin,

 I should have specified that the segfault happened for us after we got
 similar warnings and errors (DD and/or LINCS), so the segfault may
 have been tangential.  Given that everything about your system worked
 before GROMACS 4.5, it's possible that your older compilers are
 generating code that's incompatible with the GROMACS assembly loops
 (which you are likely running with, as they are the default option on
 most mainstream processors).  The bug you mentioned in your original
 post also has my antennae twitching about bad machine code.

 If that's indeed happening, it's almost certainly some bizarre
 alignment issue, something like half of a float is getting overwritten
 on the way into or out of the assembly code, which corruption would
 trigger the results you describe.  It's also distantly possible that
 GROMACS is working fine, but your copy of FFTW or BLAS/LAPACK (more
 likely the latter) has alignment problems.  One final possibility
 (which would explain the failure on YellowDog but unfortunately not
 the failure on OS X) is that GCC is generating badly-aligned code for
 auto-vectorized Altivec loops, which is still a problem for Intel's
 SIMD instructions on 32-bit x86 architectures even with GCC 4.4.  I've
 also observed MPI gather/reduce operations to foul up alignment (or
 rigidly enforce it where badly compiled code is relying on broken
 alignment) under exceedingly rare circumstances, usually involving
 different libraries compiled with different compilers (which is
 generally a bad idea for scientific code anyway).

 Okay...so all of that said, there are a few things to try:

 1) Recompile GROMACS using -O2 instead of -O3; that'll turn off the
 automatic vectorizer (on Yellow Dog) and various other relatively
 risky optimizations (on both platforms).  CFLAGS=-O2 -march=powerpc
 in the environment AND on the configure command line would do that.
 Check your build logs to make sure it took, though, because if you
 don't do it exactly right, configure will ignore your directives and
 merrily set up GROMACS to compile with -O3, which is the most likely
 culprit for badly-aligned code.

 2) Recompile GROMACS specifying a forced alignment flag.  I have no
 experience with PowerPC, but -malign-natural and -malign-power look
 like good initial guesses.  That's probably going to cause more
 problems than it solves, but if you have a screwy BLAS/LAPACK or MPI,
 it might help. I only suggest it because if you've already tried #1,
 it will only take another half hour or hour of your time to recompile
 GROMACS again.  Other than that, tinkering with alignment flags is a
 really easy way to REALLY break code, so you might consider skipping
 this and moving straight on to #3.

 3) Snag GCC 4.4.4 or 4.4.5 and compile it, and use that to compile
 GROMACS, again with -O2.  GCC takes forever to compile, but beyond
 that, it's not as difficult as it could be.  Nothing preventing you
 from installing it in your home directory, either, assuming you set
 PATH and LD_LIBRARY_PATH (or DYLD_LIBRARY_PATH on OS X) properly.  You
 might need to snag a new copy of binutils as well, if gcc refuses to
 compile with the system ld.  This option would also probably get you
 threading, since you certainly have hardware support for it.

 4) Rebuild your entire GROMACS stack, including FFTW, BLAS/LAPACK,
 MPI, and GROMACS itself with the same compiler (preferably GCC from
 #3) and the same compiler options (which again should be -O2, and
 definitely NOT any sort of alignment flag).  Put them in their own
 tree (like /opt/sci), and definitely not in /usr (which is generally
 managed by the system) or /usr/local (which tends to accumulate
 cruft).  ATLAS is a good choice for BLAS, and there are directions on
 the ATLAS website for building a complete and optimized LAPACK based
 on BLAS.

 In practice, I've found I've had to do #4 for every piece of
 scientific software our group uses, because pretty much nothing works
 right with OS-installed versions of compilers, BLAS/LAPACK, and MPI.
 It takes 

Re: [gmx-users] free energy

2011-02-15 Thread Michael Shirts
One other thing I would point out is that the solvation free energy is
dependent on concentration.  you will get a different result with 4
polymer chains vs 3 vs 3, etc.  Make sure you understand the
dependence.  Also, the free energy will depend on the polymer chain
length.

Polymer and finite concentration calculations are harder to interpret
than monomer and infinite dilution calculations.  Make sure you
understand the differences.   I'm not sure understand all of them,
though, so you'll have to think about it yourself!  Basically, you
need to make sure the physical picture of the molecules in gromacs is
the physical picture of the realistic molecular system itself.


On Mon, Feb 14, 2011 at 6:28 PM, Moeed lecie...@googlemail.com wrote:

 Dear experts,

 I am going to do solvation FE of polymer (polyethylene) in a hydrocarbon
 solvent. I have prepared a system consisting of 4 polymer chains and 480
 hexane molecules with the actual density of polymer solution (~ 0.5 g/cm3).

 1- For such a study I dont know how many polymers I need to have in my
 system. If FE can be done with only one chain, am I making system bigger in
 vain? Does this matter affect the accuracy of results?

 2- I have switched off electrostatics so I am using

 free_energy  =   yes
 init_lambda  =   0
 delta_lambda =   0
 sc_alpha =   0.5
 sc-power =   1
 sc_sigma =   0.3
 couple-lambda0   =   vdw
 couple-lambda1   =   none
 couple-intramol  =   no

 In David Mobley's turorial the last three lines are not included. I wanted
 to know if I am to run say 10 simulations for different lambda, what purpose
 does the last three lines serve in 4.0.7  ? I got very close values in that
 tutorial without these settings. ( I know what these lines mean, just
 curious how these three lines affect the results in 4 X +).

 Please let me know your comments/point of view about the system and setting
 I am using.

 Thanks
 Moeed

 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Secondary structure loss in implicit solvent simulations

2011-01-25 Thread Michael Shirts
OK, that is what I was trying to figure out -- is the problem
reproducible on both GPU and CPU.  Now, you havent answered the direct
question, if the energies are the same for at least the first 5 steps
are so -- without that knowledge, then here might be different errors
occurring in GPU vs CPU.

At this point, the question is then whether this works with another
code with the same input parameters.  Sander (part of the AmberTools
system)is downloadable (I believe to anyone, not just academics), so
you can try running the same system on Sander, using the best fit to
the implicit solvent model in gromacs.

If THAT works, and Gromacs GPU fails, then it would appear to be a
problem with Gromacs implicit solvent implementation, and should be
posted to redmine as a bug, as well as described on the list with full
details (more than you have provided so far!) so that it can be
reproduced by the developers.

Best,

On Tue, Jan 25, 2011 at 5:05 AM, K. Singhal k.sing...@uva.nl wrote:
 Hi

 It's not necessarily GPU-specific, it's implicit solvent-specific. I don't 
 get these problems in explicit solvent simulations on CPU, only in implicit 
 solvent simulations both on GPU as well as CPU. One of the problems that I 
 can think of is unbalanced charges that I would have balanced out using NaCl 
 ions, but not any more.

 Regards
 Kush



 --
 Kushagra Singhal
 Promovendus, Computational Chemistry
 van 't Hoff Institute of Molecular Sciences
 Science Park 904, room C2.119
 1098 XH Amsterdam, The Netherlands
 +31 205256965
 Universiteit van Amsterdam
 k.sing...@uva.nl


 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at 
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Re: Secondary structure loss in implicit solvent simulations

2011-01-24 Thread Michael Shirts
 Do you get the same effects if you run a normal simulation on CPU and not
 GPU?  That information would be critical for properly diagnosing what's
 going on. If it's not GPU-specific, in all likelihood whatever you're doing
 is incorrect somewhere along the way.

Folllowing up on this -- if it's not the same trajectory as the CPU to
within single precision, then something is wrong.

So the energies should be the same to 5 digits for the first time
steps, and then gradually drift away. If it's not the same for the
first 5-10 steps, then there is something going wrong with the GPU
code or the setup for the GPU code.
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Secondary structure loss in implicit solvent simulations

2011-01-17 Thread Michael Shirts
A few questions:

1) What force field are you using?
2) do you get the same answers with and without GPU acceleration?
3) How long does it take for secondary structure to disappear?  100's
of ps?  10's of ns?


On Mon, Jan 17, 2011 at 8:45 AM, K. Singhal k.sing...@uva.nl wrote:
 Hi

 I have a question regarding the implicit solvent implementation in GROMACS 
 4.5.1 (with GPU acceleration). I have been trying to simulate molecules with 
 varying size (PYP with 125 AAs, BBA5 with ~20 AAs, and and Trigger Factor 
 with 432 AAs), but can't seem to be able to maintain the secondary structure 
 in any of them. Can any one suggest a particular set of parameters than need 
 to be taken care of or used with certain fixed values for better results?

 While I have used a varied range of values for almost all parameters, here is 
 the list of parameters recently and unsuccessfully used:
   integrator           = sd
   nsteps               = 1000
   init_step            = 0
   ns_type              = Grid
   nstlist              = 10
   ndelta               = 2
   nstcomm              = 10
   comm_mode            = Linear
   delta_t              = 0.002
   pme_order            = 4
   ewald_rtol           = 1e-05
   ewald_geometry       = 0
   epsilon_surface      = 0
   optimize_fft         = FALSE
   ePBC                 = xyz
   bPeriodicMols        = FALSE
   bContinuation        = FALSE
   bShakeSOR            = FALSE
   etc                  = No
   nsttcouple           = -1
   epc                  = No
   epctype              = Isotropic
   nstpcouple           = -1
   tau_p                = 1
   andersen_seed        = 815131
   rlist                = 1.1
   rlistlong            = 2.5
   rtpi                 = 0.05
   coulombtype          = Cut-off
   rcoulomb_switch      = 0
   rcoulomb             = 2.5
   vdwtype              = Cut-off
   rvdw_switch          = 0
   rvdw                 = 2.5
   epsilon_r            = 1
   epsilon_rf           = 1
   implicit_solvent     = GBSA
   gb_algorithm         = OBC
   gb_epsilon_solvent   = 80
   nstgbradii           = 1
   rgbradii             = 1.1
   gb_saltconc          = 0.02
   gb_obc_alpha         = 1
   gb_obc_beta          = 0.8
   gb_obc_gamma         = 4.85
   gb_dielectric_offset = 0.009
   sa_algorithm         = Still
   sa_surface_tension   = 2.092
   shake_tol            = 7.5e-06
   bd_fric              = 0.5
   ld_seed              = 2010

 Thanks and Regards
 Kush



 --
 Kushagra Singhal
 Promovendus, Computational Chemistry
 Universiteit van Amsterdam

 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at 
 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-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] 5 identical simulations -- 5 different results, why?

2010-06-16 Thread Michael Shirts
Things to check: were the random numbers initializing the velocities
the same?  Looking at the log files for each simulation should be able
to tell you (look at the initialized kinetic energy at the beginning).
 Some options choose a different random number automatically for each
simulation.

On Wed, Jun 16, 2010 at 12:09 PM, Ricardo O. S. Soares
ross_...@yahoo.com.br wrote:
 Hello dear users,

 in order to introduce Gromacs to a few students, I'm taking five computers
 with the same operational system (ubuntu), to run the same simulation
 independently. I prepared a box o water and ran all standard steps until the
 grompp. This is when I copied and pasted the folder to the other four
 computers. So I ran the same grompp command in each unit, and then the mdrun
 step. At the end we obtained 5 slightly different results for total energy,
 pressure, volume etc. Since MD is a deterministic process, one may expect to
 find the same results for the same simulations. So does anyone have a clue
 at which step is it possible that those little differences are being brought
 into the systems? I'm guessing the grompp step...

 Thanks a lot.


 ---



 Ricardo O. S. Soares , MsC.
 Group of Biological Physics - Department of Physics  Chemistry
 Faculty of Pharmaceutical Sciences at Ribeirão Preto - University of São
 Paulo.
 Av.do Café, S/N - ZIP:14040-903 - Ribeirão Preto, São Paulo,  Brazil.
 Phone: +55 16 36024840.

 ross_...@yahoo.com.br,rsoa...@fcfrp.usp.br





 --
 gmx-users mailing list    gmx-us...@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at http://www.gromacs.org/search before posting!
 Please don't post (un)subscribe requests to the list. Use the
 www interface or send it to gmx-users-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/mailing_lists/users.php

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


Re: [gmx-users] Gromacs commercial license

2010-01-27 Thread Michael Shirts
 Gromacs is under GPLv2 license (http://www.gnu.org/licenses/gpl-2.0.html) and
 you can use it for commercial purpose without any problem.

EXCEPT distributing new commercial binaries that use any of the
Gromacs code, without also distributing the source for the code.

(did I get all the GPL correct there? ;)


Best,

Michael Shirts
Assistant Professor
Department of Chemical Engineering
University of Virginia
michael.shi...@virginia.edu
(434)-243-1821
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


Re: [gmx-users] free energy calculations using MBAR and gromacs

2009-11-30 Thread Michael Shirts
This looks like more of a pymbar question!

  I am trying to calculate free energy of a system that involves
 disappearance of LJ particle at lambda=1 in explicit solvent. I ran the
 simulation at 20 different lambda points ranging from 0 to 1, using soft
 core potential. In order to use MBAR method (python implementation from
 Michael shirts ), I have to rerun the simulations, in order to get
 potential energies for each configuration at every lambda point.

That's correct.  Very soon (i.e., the code is in the CVS) it will be
possible to generate these energy differences during single runs, by
setting the keyword foreign-lambda.

 Now,
 here i get into trouble. If i evaluate the trajectories generated at
 lambda = 1 (particle completely decoupled) using the potential function
 at lambda = 0, i get very large energy values,

Yes; if you look at the example in pymbar, there are lots of values of
E0 - E1 of around 10^14.

 which quite likely are
 responsible for MBAR's failure to converge.

That shouldn't be the case.  These energy values turn into weights,
which are essentially zero, and thus don't affect MBAR.  This sounds
like is more of a pymbar question that a Gromacs question!

 The program works well if i
 calculate free energy to/from an intermediate. ie lambda 0 - 0.5 and
 0.5-1.0

Hmm. It should work well in either case!

 As the potential energies for the configurations generated at lambda=1,
 evaluated using potential energy function at lambda =0, is expected to
 be infinite (due to Van der Waals overlaps between solvent and LJ
 particle),

Very large but finite, but yes.

 Is it even possible to apply MBAR method to such cases,
 without splitting the analysis to an intermediate lambda value? Did
 anyone ran into similar problems, or am i missing something?

Yes, it is possible.  If you send a copy of 1) the data and 2) the
scripts you use to interface with MBAR, I can take a look to try to
see what the matter will be.

Best,
Michael
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


Re: [gmx-users] Free energy calculations

2009-11-09 Thread Michael Shirts
Can you be more specific in the question of what property you want to
compute for which molecule or molecules? Linear Interaction energy
approaches are not always very efficient for coulomb, and is
definitely not exact for the LJ term.  So LIE might not be what you
want to be doing.

Michael Shirts
Assistant Professor
Department of Chemical Engineering
University of Virginia
michael.shi...@virginia.edu
(434)-243-1821

On Mon, Nov 9, 2009 at 1:49 PM,  jorge_quint...@ciencias.uis.edu.co wrote:
 I'm interesting to calculate Gibbs energy using g_lie application.
 Unfourtunately, I don't know what kind of state equation is the
 sub-program using.  Can you help me with that.


 Best regards


 Jorge R. Quintero
 Grupo de Investigación en Fisicoquímica Teórica y Experimental
 Universidad Industrial de Santander
 Bucaramanga, Santander - Colombia

 --
 gmx-users mailing list    gmx-us...@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at http://www.gromacs.org/search before posting!
 Please don't post (un)subscribe requests to the list. Use the
 www interface or send it to gmx-users-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/mailing_lists/users.php

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


  1   2   >