# [gmx-users] Test Particle Insertion

Javier Cerezo jcb1 at um.es
Wed May 16 11:28:28 CEST 2012

```About the red curve, I guess fluctuations might be directly related to
volume fluctuations, you can extract the volume over time from g_energy
(boxXX*boxYY**boxZZ) and compare. (just another comment, now I am not
very sure about the "f." that precedes the red line legend..)

About the interpretation of the quantities, the Widom technique does not
provide you with an absolute value of the chemical potential but
directly with the excess chemical potential. So, mu=-kTlog(Ve ^
(U*B)/(V))n+1 is the excess chemical potential, where (if I recall
correctly) U_{n+1} is the the interaction energy  between the inserted
particle and the rest of the system. You don't need (and should not do)
such post-processing operations that you proposed to get the excess
chemical potential.

Javier

El 16/05/12 11:06, Steven Neumann escribió:
> Thank you very much! I just saw your response.
>
> As I run it in NPT ensemble the plot with volume is important for me.
>
> http://speedy.sh/CJn5b/tpiN.jpg
>
> So does the fluctuating red curve make any sesnse then if it does not
> consider volume?
>
> Another thing: this is chemical potential of the system with extra
> water molecule (N+1), right (u=-kTlog(Ve ^ (U*B)/(V))n+1? So if I want
> to obtain the excess chemical potential:
> u=-kTlog(Ve ^ (-deltaU*B)/(V)) I should calculate it for the system
> with N molecules and then substract it.
> Is it calculated somewhere or I should use g_energy of my previosu
> system and calculate the total potential energy then -kTlog... of this
> values and then substract it? Please correct me if I am wrong.
>
> Steven
>
>
> On Tue, May 15, 2012 at 11:59 AM, Javier Cerezo <jcb1 at um.es
> <mailto:jcb1 at um.es>> wrote:
>
>     Hi Steven.
>
>
>>     1. Why this value is divided by nm3? Shall I multiply it by the
>>     simulation box?
>     It is not not divided by nm3. The legend for "y" axis is not
>     appropriate for your plot. Keep in mind that the same graph is
>     used to represent lots of quantities (you can plot all of them
>     with xmgrace -nxy tpi.xvg). The "y" axis is not the same for all,
>     but only one label is possible, so developers have to chose which
>     label to place on the axis. But this is just a label, don't give
>     much importance to it and analyse you results (including units)
>     according to the equations and the standard units in gromacs.
>
>>     2. Why e^(-BU) is multiplied by V? I just want to have the excess
>>     chemical potential: u=-kTlog(e ^ (-deltaU*B) - so how can I get
>>     deltaU?
>     The volume appears in the expression of the excess chemical
>     potential if you are running a NpT ensemble. The second plot (if
>     you use xmgrace -nxy tpi.xvg) does not contain the volume.
>
>>     3. The value corresponds to the plateau so I should run it for
>>     longer time?
>     You are getting a time&ensemble average and for large sampling
>     (and large simulation times), this average should converge. So,
>     the final value you will get is the last point of the graph, it up
>     to you to say if it is converged. So you can try to enlarge the
>     number of points sampled, if the shape does not change you are
>     sampling correctly every snapshot, then take longer simulation
>     times if you want to converge your results.
>
>     Javier
>
>
>     El 15/05/12 09:57, Steven Neumann escribió:
>>
>>
>>     On Mon, May 14, 2012 at 5:05 PM, Justin A. Lemkul
>>     <jalemkul at vt.edu <mailto:jalemkul at vt.edu>> wrote:
>>
>>
>>
>>         On 5/14/12 11:53 AM, Steven Neumann wrote:
>>
>>             Dear Gmx Users,
>>
>>             Did anyone use TPI method for the calculation of chemical
>>             potential? The tpi.xvg
>>             files consists of:
>>
>>             @ s0 legend "-kT log(<Ve\S-\xb\f{}U\N>/<V>)"
>>             @ s1 legend "f. -kT log<e\S-\xb\f{}U\N>"
>>             @ s2 legend "f. <e\S-\xb\f{}U\N>"
>>             @ s3 legend "f. V"
>>             @ s4 legend "f. <Ue\S-\xb\f{}U\N>"
>>             @ s5 legend "f. <U\sVdW System\Ne\S-\xb\f{}U\N>"
>>             @ s6 legend "f. <U\sdisp c\Ne\S-\xb\f{}U\N>"
>>             @ s7 legend "f. <U\sCoul System\Ne\S-\xb\f{}U\N>"
>>             @ s8 legend "f. <U\sCoul recip\Ne\S-\xb\f{}U\N>"
>>
>>             @    xaxis  label "Time (ps)"
>>             @    yaxis  label "(kJ mol\S-1\N) / (nm\S3\N)"
>>
>>             Can anyone explain me these legends? I just want obtain a
>>             value of the excess
>>             chemical potential according to the equation:
>>             u=-kT log (-deltaV/kT), Which legend is responsible for
>>             this and what are the
>>             units? kJ/mol? Please, explain as the above letters does
>>             not mean to me anything?
>>
>>
>>         These strings are formatted for XmGrace.  Have you tried
>>         plotting the file to see what it contains?  The legends will
>>         be far more obvious if you do.
>>
>>         -Justin
>>
>>
>>     Thank you Justin.
>>     Can anyone explain me from the plot:
>>
>>     http://speedy.sh/Xpnws/tpi.JPG
>>
>>     1. Why this value is divided by nm3? Shall I multiply it by the
>>     simulation box?
>>     2. Why e^(-BU) is multiplied by V? I just want to have the excess
>>     chemical potential: u=-kTlog(e ^ (-deltaU*B) - so how can I get
>>     deltaU?
>>     3. The value corresponds to the plateau so I should run it for
>>     longer time?
>>
>>
>>     Thank you,
>>
>>     Steven
>>
>>
>>         --
>>         ========================================
>>
>>         Justin A. Lemkul, Ph.D.
>>         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
>>
>>         ========================================
>>         --
>>         gmx-users mailing list gmx-users at gromacs.org
>>         <mailto:gmx-users at 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-request at gromacs.org
>>         <mailto:gmx-users-request at gromacs.org>.
>>
>>
>>
>>
>
>     --
>     Javier CEREZO BASTIDA
>     PhD Student
>     Physical Chemistry
>     Murcia (Spain)
>     Tel: (+34)868887434 <tel:%28%2B34%29868887434>
>
>     --
>     gmx-users mailing list gmx-users at gromacs.org
>     <mailto:gmx-users at 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-request at gromacs.org
>     <mailto:gmx-users-request at gromacs.org>.
>
>
>
>

--
Javier CEREZO BASTIDA
Ph.D. Student
Physical Chemistry