R: RE: R: RE: R: RE: R: RE: R: RE: R: Re: R:[gmx-users] Tabulated potential - Problem
Mark Abraham
Mark.Abraham at anu.edu.au
Thu Oct 1 22:10:37 CEST 2009
Berk Hess wrote:
> Ah, I thought the sentence in one of my previous mails was clear enough:
> "The points in the table_b.xvg file should be equally spaced.
> grompp assumes equal spacing, but does not check this."
>
> I think "should be" and "must be" are equivalent in this context.
> But we can clarify the manual, and add a check in grompp.
It's close, but in general there is a shade of certainty conveyed by
"must" that is additional to that of "should". The dictionary I have to
hand contrasts "to be obliged or compelled to, as by some constraining
force or necessity" with "to denote duty, propriety, expediency"
respectively. Wiktionary agrees - http://en.wiktionary.org/wiki/must. I
am aware of usage in a body of law where failure to follow a direction
qualified by "must" carries a more severe consequence than failure to
follow a direction qualified by "should".
Berk is correct in this context, however. The only relevant qualities of
the table are that it be accurate and reasonably efficient. A user
failing to follow "should" would be making an unwarranted assumption
that accuracy and/or efficiency was compromised only to a negligible
degree. "must" is slightly superior in such a context where no purpose
is served by the distinction with "should," and the reader might err
more often when "should" is used. Manual section 6.7.1 certainly
expresses the idea that equal spacing of table points is normal and
expected, but there is no statement which conveys a suitable degree of
definiteness. That could be improved.
Even the use of language can be as exacting as the methods described by
it :-) I tied myself in knots last weekend with an incorrect use of
"abjure" instead of "adjure"...
Mark
> But at the connecting point the potentials should be equal by
> definition, right?
> Anyhow, for numerical work you should smooth you potentials in some way.
> You can do this and generate a lot of extra (intermediate, smoothed) points
> or use a small amount of points. But in the end they should be equally
> spaced
> for computational efficiency.
> mdrun could, in principle, also do this job, but it is better if the
> user thinks about this.
>
> Berk
>
> ------------------------------------------------------------------------
> Date: Thu, 1 Oct 2009 17:13:30 +0200
> From: albitauro at virgilio.it
> To: gmx3 at hotmail.com
> Subject: R: RE: R: RE: R: RE: R: RE: R: RE: R: Re: R:[gmx-users]
> Tabulated potential - Problem
>
> Yes I know about it.
> This is the reason for my previous questions if "should be uniformly
> spaced" means indeed "must be equally spaced".
> When generating a numerical potential from boltzmann inversion I need to
> "complete" on a wider range of distance with the left and right
> repulsive branches. At the connecting points of the different branches
> there is usually the need for smoothing to avoid numerical "dirt"...the
> program I use causes a uniform spacing to be lost.
> Here's the problem: if grompp works only if points are equally spaced I
> have to find a new way to prepare my table but if it can work also
> without uniform spacing then there would not have been the need for this
> effort. On the manual I interpreted that "should be" as a "it's better
> but not required".
> Probably it's better to add a message to warn that it is required.
> Thank you very much for your help
>
> AM
>
> ----Messaggio originale----
> Da: gmx3 at hotmail.com
> Data: 1-ott-2009 2.56 PM
> A: <albitauro at virgilio.it>, "Discussion list for GROMACS
> users"<gmx-users at gromacs.org>
> Ogg: RE: R: RE: R: RE: R: RE: R: RE: R: Re: R:[gmx-users] Tabulated
> potential - Problem
>
> Have you ever looked at the table file?
> It looks like this:
>
> 0.000000 112.891000 531.46
> 0.002000 111.831000 530.234883
> 0.002001 111.830000 529.000000
> 0.003000 111.302000 527.736132
> 0.003001 111.302000 526.000000
> 0.004000 110.776000 526.736632
> 0.004002 110.775000 526.000000
> 0.005000 110.250000 525.737131
>
> This is certainly not equally spaced.
>
> I guess it might be useful to add a check for this in grompp.
> But the user also has a responsibility for checking the input.
>
>
> Berk
>
> ------------------------------------------------------------------------
> Date: Thu, 1 Oct 2009 16:48:01 +0200
> From: albitauro at virgilio.it
> To: gmx3 at hotmail.com
> Subject: R: RE: R: RE: R: RE: R: RE: R: Re: R:[gmx-users] Tabulated
> potential - Problem
>
> Sure!
> Thank you!
>
>
> ----Messaggio originale----
> Da: gmx3 at hotmail.com
> Data: 1-ott-2009 2.36 PM
> A: <albitauro at virgilio.it>, "Discussion list for GROMACS
> users"<gmx-users at gromacs.org>
> Ogg: RE: R: RE: R: RE: R: RE: R: Re: R:[gmx-users] Tabulated potential -
> Problem
>
> Could you mail me your table file?
>
> Berk
>
> ------------------------------------------------------------------------
> Date: Thu, 1 Oct 2009 16:29:43 +0200
> From: albitauro at virgilio.it
> To: gmx3 at hotmail.com
> Subject: R: RE: R: RE: R: RE: R: Re: R:[gmx-users] Tabulated potential -
> Problem
>
> Hi Berk,
> I tried such a test on only one bond (already at its equilibrium value)
> and I got:
>
> -For analytic potential:
>
> Polak-Ribiere Conjugate Gradients converged to Fmax < 0.1 in -1 steps
> Potential Energy = 5.5511151e-13
> Maximum force = 3.7252903e-05 on atom 1
> Norm of force = 3.7252903e-05
>
> -for tabulated potential:
> Polak-Ribiere Conjugate Gradients converged to machine precision in -1
> steps,
> but did not reach the requested Fmax < 0.1.
> Potential Energy = 7.0767384e-04
> Maximum force = 5.3229499e-01 on atom 1
> Norm of force = 5.3229499e-01
>
> I am sure that there are no problems in the code...it's something
> related to my input but I cannot understand what!
> Could be due to the forces specified in table.xvg? Since I would like to
> use a general numerical potential generated by Boltzmann inversion I
> calculated them in any case (also for the harmonic potential)
> numerically as F(xi)= - [V(i+1)-V(i-1)] / [x(i+1) - x(i-1)]
> Thanks!
>
> Alberto
>
>
> ----Messaggio originale----
> Da: gmx3 at hotmail.com
>
> Data: 1-ott-2009 8.04 AM
> A: <albitauro at virgilio.it>
> Ogg: RE: R: RE: R: RE: R: Re: R:[gmx-users] Tabulated potential - Problem
>
> Hi,
>
> I would start checking a simple system.
> Take a system with only one of these bonds.
> Run a simulation of 0 steps with the analytic and tabulated bond
> and compare the energy and the force.
>
> I am 99.9% sure there are no problems in the code,
> so it should be something in your input.
>
> Berk
>
> ------------------------------------------------------------------------
> Date: Thu, 1 Oct 2009 09:57:10 +0200
> From: albitauro at virgilio.it
> To: gmx-users at gromacs.org
> Subject: R: RE: R: RE: R: Re: R:[gmx-users] Tabulated potential - Problem
>
> Ok! Thanks.
> My table should be correct...yes in the third column I have the forces.
> In any case, my simulations with the same potential (analytical or
> numerical) still give results which are largely different.
> Do you have any further suggestions?
> Thank you very very much!
>
> AM
>
>
>
> ----Messaggio originale----
> Da: gmx3 at hotmail.com
> Data: 1-ott-2009 7.00 AM
> A: <albitauro at virgilio.it>, "Discussion list for GROMACS
> users"<gmx-users at gromacs.org>
> Ogg: RE: R: RE: R: Re: R:[gmx-users] Tabulated potential - Problem
>
> Hi,
>
> You can put the force constant in the table and use a force constant of
> 1 in the topology,
> or only put the 1/2 in the table and put the force constant in the
> topolgy, the effect is the same.
>
> The points in the table_b.xvg file should be equally spaced.
> grompp assumes equal spacing, but does not check this.
>
> Also note that you have to give the force, the negative of the derivative.
> But I assume you have done this, as I have built a check in grompp to give
> a warning when the potential and the force are inconsistent.
>
> Berk
>
> ------------------------------------------------------------------------
> Date: Thu, 1 Oct 2009 08:53:51 +0200
> From: albitauro at virgilio.it
> To: gmx-users at gromacs.org
> Subject: R: RE: R: Re: R:[gmx-users] Tabulated potential - Problem
>
> Hi Berk,
> the tabulated potential I used (the file table_b.xvg) is exactly the
> function
> V(R)=1/2 *k*(R-R0)^2 in numerical form. I printed in the first column
> the R variable from 0 to 2 nm in step of 0.001 nm,
> on the second column the associated V(R) value (KJ/mol) and in the
> third column its derivatives. With this choice I chose the force constant
> k (that must be defined for tabulated function) to 1.0.
> About this point (please forgive me...it's probably a silly question!)
> I am a bit confused about this force constant since I don't understand
> its meaning...
> my numerical potential should be already correct in dimensions and units
> and so I guessed that k should have been equal to 1.0.
> Further, I use the occasion to ask you also if the x points in the table
> "should be uniformly spaced" as reported in the manual so that there are
> no problems even if it's not so or they "must be uniformly spaced" to
> obtain the correct results.
> Thank you,
>
> AM
>
>
> ----Messaggio originale----
> Da: gmx3 at hotmail.com
> Data: 30-set-2009 1.36 PM
> A: <albitauro at virgilio.it>, "Discussion list for GROMACS
> users"<gmx-users at gromacs.org>
> Ogg: RE: R: Re: R: RE: R: RE: R: RE: R: RE: [gmx-users] Tabulated
> potential - Problem
>
> Hi,
>
> Are you really sure you entered the tabulated potential correctly?
> There is for instance the pre-factor 1/2 in front of the harmonic potential,
> which you will have to add explicitly, either in the table or in the
> force constant.
>
> Berk
>
> ------------------------------------------------------------------------
> Date: Wed, 30 Sep 2009 14:33:31 +0200
> From: albitauro at virgilio.it
> To: gmx-users at gromacs.org
> Subject: R: Re: R: RE: R: RE: R: RE: R: RE: [gmx-users] Tabulated
> potential - Problem
>
> Dear Mark,
> unfortunately the differences I was speaking about are very consistent.
> I am simulating a short chain with MARTINI CG force field:
> with an analytical harmonic potential the chain becomes a random coil,
> when using the same potential but in tabulated numerical form the chain
> remains extended and only slightly distorted from linearity. All the
> parameters of the simulation (box sizes, mdp option...) are the same.
> I will try the options for mdrun.
> Thank you,
>
> Alberto
>
> ----Messaggio originale----
> Da: Mark.Abraham at anu.edu.au
> Data: 30-set-2009 12.15 PM
> A: "Discussion list for GROMACS users"<gmx-users at gromacs.org>
> Ogg: Re: R: RE: R: RE: R: RE: R: RE: [gmx-users] Tabulated potential -
> Problem
>
> albitauro at virgilio.it wrote:
> > I changed bond type to 8.
> > In any case the two MD simulations (same harmonic potential but in
> > analytical form vs tabulated form) gives different results.
> > All the pararameters of the simulation are the same.
>
> Slight differences will be normal unless you're running with "mdrun
> -reprod." Even then, there will be small differences between tabulated
> and non-tabulated. If you really want to be sure, you should consider
> doing an "mdrun -rerun -reprod" so that you are calculating the same
> quantities on the same inputs.
>
> Mark
>
> > ----Messaggio originale----
> > Da: gmx3 at hotmail.com
> > Data: 29-set-2009 3.54 PM
> > A: <albitauro at virgilio.it>, "Discussion list for GROMACS
> > users"<gmx-users at gromacs.org>
> > Ogg: RE: R: RE: R: RE: R: RE: [gmx-users] Tabulated potential - Problem
> >
> > No.
> >
> > The ONLY difference between bonds type 8 and type 9 is that type 8
> > generates exclusions
> > while type 9 does not (see table 5.4 in the manual).
> > Simply changing from type 9 to 8 will generate the exclusions.
> >
> > Berk
> >
> > ------------------------------------------------------------------------
> > Date: Tue, 29 Sep 2009 17:52:03 +0200
> > From: albitauro at virgilio.it
> > To: gmx-users at gromacs.org
> > Subject: R: RE: R: RE: R: RE: [gmx-users] Tabulated potential - Problem
> >
> > In my simulation I used bond type 9 because I was convinced that 1-2 LJ
> > interaction would have been excluded in any case as default for bond
> > stretching interaction. After my simulations it seems not and so I
> > suppose that I must use bond type 8 and list the exclusion in the
> > itp...Right?
> >
> > Alberto
> >
> > ----Messaggio originale----
> > Da: gmx3 at hotmail.com
> > Data: 29-set-2009 2.15 PM
> > A: "Discussion list for GROMACS users"<gmx-users at gromacs.org>
> > Ogg: RE: R: RE: R: RE: [gmx-users] Tabulated potential - Problem
> >
> > Are you really sure about this and that this is with bond type 8?
> >
> > The whole point of having a tabulated bond type 8 and 9
> > is that 8 does generate exclusions and 9 does not.
> >
> > Berk
> >
> > ------------------------------------------------------------------------
> > Date: Tue, 29 Sep 2009 14:17:37 +0200
> > From: albitauro at virgilio.it
> > To: gmx-users at gromacs.org
> > Subject: R: RE: R: RE: [gmx-users] Tabulated potential - Problem
> >
> > Hi,
> >
> > in order to check further which kind of problems are present when using
> > tabulated potential, I carried out two simulations on the same system by
> > using the same harmonic
> > potential for bond stretching but either in analytical form or in
> > tabulated form respectively. The
> > results of the two simulations are different!
> > I repeated the calculationis on just two bonded particles and I verified
> > that when using
> > analytical stretching potentials 1-2 LJ interactions are excluded
> while they
> > are not excluded when using tabulated potential.
> > Am I right?
> >
> > This should mean that I have to use a function type 8 and list
> > explicitly the elements of the [
> > exclusions ] field or is there another method?
> >
> > Thanks
> >
> > Alberto
> >
> > ----Messaggio originale----
> > Da: gmx3 at hotmail.com
> > Data: 25-set-2009 12.27 PM
> > A: "Discussion list for GROMACS users"<gmx-users at gromacs.org>
> > Ogg: RE: R: RE: [gmx-users] Tabulated potential - Problem
> >
> > Your system could be unstable.
> > You can check for large forces with mdrun -pforce
> > I don't know what a reasonable range of forces is, you can try 5000.
> > If you have instabilities, you should get large forces printed
> > before you get the fatal error.
> >
> > Berk
> >
> >
> > ------------------------------------------------------------------------
> > Date: Fri, 25 Sep 2009 14:10:08 +0200
> > From: albitauro at virgilio.it
> > To: gmx-users at gromacs.org
> > Subject: R: RE: [gmx-users] Tabulated potential - Problem
> >
> > Unfortunately, my box sizes are not close to 23. I also carried out
> > calculations switching off PBC or on much smaller systems.
> > I received always the same error.
> > I tried also a geometry optimization. It finished without warnings nor
> > errors: anyway the potential energy changed only very slightly during
> > the simulation with too large values.
> >
> > Thanks
> >
> > AM
> >
> > ----Messaggio originale----
> > Da: gmx3 at hotmail.com
> > Data: 24-set-2009 11.29 AM
> > A: "Discussion list for GROMACS users"<gmx-users at gromacs.org>
> > Ogg: RE: [gmx-users] Tabulated potential - Problem
> >
> > This is not nonsense, it is exactly what is says.
> > The distance between two atoms is more than 10 times as large as your
> > table length.
> >
> > Maybe you are somehow having issues with periodic boundary conditions.
> > Is you box size close to 23?
> >
> > Berk
> >
> >
> > ------------------------------------------------------------------------
> > Date: Thu, 24 Sep 2009 12:32:36 +0200
> > From: albitauro at virgilio.it
> > To: gmx-users at gromacs.org
> > Subject: [gmx-users] Tabulated potential - Problem
> >
> > Hi,
> >
> > I'm trying to carry out a CG simulation and I'm using
> > a tabulated potential for a bond stretching term.
> > My MD simulations stops immediately with the error message:
> >
> > -------------------------------------------------------
> > Program mdrun_mpi, VERSION 4.0.5
> > Source code file: bondfree.c, line: 1772
> >
> > Fatal error:
> > A tabulated bond interaction table number 0 is out of the table range: r
> > 23.678833, between table indices 12069 and 12070, table length 1020
> > -------------------------------------------------------
> >
> > This should mean that some distances are beyond table length (as
> > reported in the manual) but this is
> > nonsense considering my input files and topology.
> >
> > Do you have any suggestion?
> > Thanks!
> >
> > AM
> >
> >
> > ------------------------------------------------------------------------
> > Express yourself instantly with MSN Messenger! MSN Messenger
> > <http://clk.atdmt.com/AVE/go/onm00200471ave/direct/01/>
> >
> >
> > ------------------------------------------------------------------------
> > What can you do with the new Windows Live? Find out
> > <http://www.microsoft.com/windows/windowslive/default.aspx>
> >
> >
> > ------------------------------------------------------------------------
> > Express yourself instantly with MSN Messenger! MSN Messenger
> > <http://clk.atdmt.com/AVE/go/onm00200471ave/direct/01/>
> >
> >
> > ------------------------------------------------------------------------
> > Express yourself instantly with MSN Messenger! MSN Messenger
> > <http://clk.atdmt.com/AVE/go/onm00200471ave/direct/01/>
> >
> >
> > ------------------------------------------------------------------------
> >
> > _______________________________________________
> > gmx-users mailing list gmx-users at 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-request at gromacs.org.
> > Can't post? Read http://www.gromacs.org/mailing_lists/users.php
> _______________________________________________
> gmx-users mailing list gmx-users at 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-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>
>
>
> ------------------------------------------------------------------------
> See all the ways you can stay connected to friends and family
> <http://www.microsoft.com/windows/windowslive/default.aspx>
>
>
> ------------------------------------------------------------------------
> Express yourself instantly with MSN Messenger! MSN Messenger
> <http://clk.atdmt.com/AVE/go/onm00200471ave/direct/01/>
>
>
> ------------------------------------------------------------------------
> See all the ways you can stay connected to friends and family
> <http://www.microsoft.com/windows/windowslive/default.aspx>
>
>
> ------------------------------------------------------------------------
> Express yourself instantly with MSN Messenger! MSN Messenger
> <http://clk.atdmt.com/AVE/go/onm00200471ave/direct/01/>
>
>
> ------------------------------------------------------------------------
> Express yourself instantly with MSN Messenger! MSN Messenger
> <http://clk.atdmt.com/AVE/go/onm00200471ave/direct/01/>
>
>
> ------------------------------------------------------------------------
> See all the ways you can stay connected to friends and family
> <http://www.microsoft.com/windows/windowslive/default.aspx>
>
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> gmx-users mailing list gmx-users at 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-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
More information about the gromacs.org_gmx-users
mailing list