[gmx-users] Re: Enthalpy Confusion

2013-06-14 Thread Jeffery Perkins
 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. 
 
 sorry I always have issues deciding how to talk about this stuff, so
 thanks 
 for putting up with my terrible notation =) 
 
 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. 
 
 I hadn't really thought about that... but the system is a monatomic 
 Lennard-Jones particle (uncharged sigma = 0.35 nm epsilon = 2 kJ/mol mass
 = 
 40 amu) which should be a liquid at the conditions I was looking at,
 P=1000 
 bar, T = 300 K using phase diagram in: Equation of state for the 
 Lennard-Jones fluid, J. J. Nicolas et al., MOLECULAR PHYSICS, 1979, VOL.
 37, 
 No. 5, 1429-1454 
 
 the U is -1600 and the pV is 880 when manually done, around 400 from 
 g_energy 
Here's the code: 
 /* This is pV (in kJ/mol).  The pressure is the reference
 pressure, 
not the instantaneous pressure */ 
 pv = vol*md-ref_p/PRESFAC; 

 add_ebin(md-ebin, md-ipv, 1, pv, bSum); 
 enthalpy = pv + enerd-term[F_ETOT]; 

What is your volume? 
What is Etot?

a bit messy but this is data from g_energy for the system with increasing
temperature (300 K-1000 K), and constant pressure of 1000 bar (so 100 000
000 Pascals):

Temperature VolumeEnthalpy  Potential   Total EnergypV  
Kinetic En.
K nm^3   kJ/molkJ/mol
kJ/mol   kJ/mol kJ/mol
299.914 14.7226  -1168.83   -2565.44-1611.63442.791 
953.817
399.901 17.462   -339.795   -2080.31-808.501468.706 
1271.81
499.898 21.0263   421.763   -1666.7 -76.8794498.643 
  
1589.83
599.899 25.1149   1088.79   -1348.14559.718 529.07  
1907.86
699.914 29.4219   1668.2-1115.471110.46 557.733 
2225.94
799.925 33.7642184.31   -943.6131600.39 583.922 
2544
899.923 38.0865  2660.3 -809.5722052.45 607.847 
2862.03
999.931 42.2798  3101.56-707.9082472.17 629.383 
3180.08

all uncorrected for the number of particles in the box (so Etot is in kJ/256
mol of particles)

when i use this data to get H myself i get:

H
kJ/mol
-725.32948
242.7114
1188.90386
2071.63498
2881.65838
3632.9828
4345.2573
5017.41396

again uncorrected for the number of particles in the box

Thoughts?



--
View this message in context: 
http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009164.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


Re: [gmx-users] Re: Enthalpy Confusion

2013-06-14 Thread David van der Spoel

On 2013-06-14 19:28, Jeffery Perkins wrote:

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.


sorry I always have issues deciding how to talk about this stuff, so
thanks
for putting up with my terrible notation =)


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.


I hadn't really thought about that... but the system is a monatomic
Lennard-Jones particle (uncharged sigma = 0.35 nm epsilon = 2 kJ/mol mass
=
40 amu) which should be a liquid at the conditions I was looking at,
P=1000
bar, T = 300 K using phase diagram in: Equation of state for the
Lennard-Jones fluid, J. J. Nicolas et al., MOLECULAR PHYSICS, 1979, VOL.
37,
No. 5, 1429-1454

the U is -1600 and the pV is 880 when manually done, around 400 from
g_energy

Here's the code:
 /* This is pV (in kJ/mol).  The pressure is the reference
pressure,
not the instantaneous pressure */
 pv = vol*md-ref_p/PRESFAC;

 add_ebin(md-ebin, md-ipv, 1, pv, bSum);
 enthalpy = pv + enerd-term[F_ETOT];

What is your volume?
What is Etot?


a bit messy but this is data from g_energy for the system with increasing
temperature (300 K-1000 K), and constant pressure of 1000 bar (so 100 000
000 Pascals):

Temperature VolumeEnthalpy  Potential   Total EnergypV  
Kinetic En.
K nm^3   kJ/molkJ/mol
kJ/mol   kJ/mol kJ/mol
299.914 14.7226  -1168.83   -2565.44-1611.63442.791 
953.817
399.901 17.462   -339.795   -2080.31-808.501468.706 
1271.81
499.898 21.0263   421.763   -1666.7 -76.8794498.643
1589.83
599.899 25.1149   1088.79   -1348.14559.718 529.07  
1907.86
699.914 29.4219   1668.2-1115.471110.46 557.733 
2225.94
799.925 33.7642184.31   -943.6131600.39 583.922 
2544
899.923 38.0865  2660.3 -809.5722052.45 607.847 
2862.03
999.931 42.2798  3101.56-707.9082472.17 629.383 
3180.08

all uncorrected for the number of particles in the box (so Etot is in kJ/256
mol of particles)

when i use this data to get H myself i get:

H
kJ/mol
-725.32948
242.7114
1188.90386
2071.63498
2881.65838
3632.9828
4345.2573
5017.41396

again uncorrected for the number of particles in the box

Thoughts?


Your calculation seems correct. Which gmx version did you use?
The correlation between the numbers is almost 100% so there must be a 
simple explanation.





--
View this message in context: 
http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009164.html
Sent from the GROMACS Users Forum mailing list archive at Nabble.com.




--
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
* 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: Enthalpy Confusion

2013-06-14 Thread Jeffery Perkins
Your calculation seems correct. Which gmx version did you use? 
The correlation between the numbers is almost 100% so there must be a 
simple explanation.

gmx version is 4.5.4, and yeah the correlation is odd, in the code you
listed:

  pv = vol*md-ref_p/PRESFAC; 
 
  add_ebin(md-ebin, md-ipv, 1, pv, bSum); 
  enthalpy = pv + enerd-term[F_ETOT]

so pv = vol * md-ref_p/PRESFAC; is the vol * ref_p/conversion factor,
right?
the only thing there that seems like it could be off is if PRESFAC isn't
correct for some reason

(correct me if I'm wrong of course!)

Jeffery



--
View this message in context: 
http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009168.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


Re: [gmx-users] Re: Enthalpy Confusion

2013-06-14 Thread David van der Spoel

On 2013-06-14 21:39, Jeffery Perkins wrote:

Your calculation seems correct. Which gmx version did you use?
The correlation between the numbers is almost 100% so there must be a
simple explanation.


gmx version is 4.5.4, and yeah the correlation is odd, in the code you
listed:


Then that is the reason. I forgot about it, but the calculation of pv is 
incorrect in that version. If you use 4.5.7 or 4.6.2 it should be fine.



  pv = vol*md-ref_p/PRESFAC;

  add_ebin(md-ebin, md-ipv, 1, pv, bSum);
  enthalpy = pv + enerd-term[F_ETOT]


so pv = vol * md-ref_p/PRESFAC; is the vol * ref_p/conversion factor,
right?
the only thing there that seems like it could be off is if PRESFAC isn't
correct for some reason

(correct me if I'm wrong of course!)

Jeffery



--
View this message in context: 
http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009168.html
Sent from the GROMACS Users Forum mailing list archive at Nabble.com.




--
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
* 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-12 Thread David van der Spoel

On 2013-06-11 23:31, Jeffery Perkins wrote:

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.


sorry I always have issues deciding how to talk about this stuff, so thanks
for putting up with my terrible notation =)


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.


I hadn't really thought about that... but the system is a monatomic
Lennard-Jones particle (uncharged sigma = 0.35 nm epsilon = 2 kJ/mol mass =
40 amu) which should be a liquid at the conditions I was looking at, P=1000
bar, T = 300 K using phase diagram in: Equation of state for the
Lennard-Jones fluid, J. J. Nicolas et al., MOLECULAR PHYSICS, 1979, VOL. 37,
No. 5, 1429-1454

the U is -1600 and the pV is 880 when manually done, around 400 from
g_energy

Here's the code:
/* This is pV (in kJ/mol).  The pressure is the reference 
pressure,

   not the instantaneous pressure */
pv = vol*md-ref_p/PRESFAC;

add_ebin(md-ebin, md-ipv, 1, pv, bSum);
enthalpy = pv + enerd-term[F_ETOT];

What is your volume?
What is Etot?

PRESFAC  ~ 16.6










--
View this message in context: 
http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009068.html
Sent from the GROMACS Users Forum mailing list archive at Nabble.com.




--
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
* 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: Enthalpy Confusion

2013-06-11 Thread Jeffery Perkins
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


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] Re: Enthalpy Confusion

2013-06-11 Thread Jeffery Perkins
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. 

Right, that's an easy one to miss but i don't think that's my problem here. 

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

I'm running in NPT and was calcaulating H from ave. P, ave. U and ave. V for
the run
while i understand that this doesn't exactly equal average H it should be
close enough for me at this point, going back and redoing it with P set to
the reference value (which the average hits with small fluctuations)

the resultant H still doesn't line up properly for some example data:

pV, g_energyV, m^3  pV, J pV, kJ/mol
442.42  1.46E-26  1.46E-18884.07

P = 1000 bar = 100,000,000 Pa
and i see that my pV is around 2x g_energies pV, no values are scaled for
the # of items in the box at this point, this is all for moles of the system

 i think that's alright though



--
View this message in context: 
http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009062.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


Re: [gmx-users] Re: Enthalpy Confusion

2013-06-11 Thread David van der Spoel

On 2013-06-11 21:57, Jeffery Perkins wrote:

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.


Right, that's an easy one to miss but i don't think that's my problem here.


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


I'm running in NPT and was calcaulating H from ave. P, ave. U and ave. V for
the run
while i understand that this doesn't exactly equal average H it should be
close enough for me at this point, going back and redoing it with P set to
the reference value (which the average hits with small fluctuations)

the resultant H still doesn't line up properly for some example data:

pV, g_energyV, m^3  pV, J pV, kJ/mol
442.42  1.46E-26  1.46E-18884.07

P = 1000 bar = 100,000,000 Pa
and i see that my pV is around 2x g_energies pV, no values are scaled for
the # of items in the box at this point, this is all for moles of the system

  i think that's alright though


You should not use pV from g_energy though, as Michael explained, rather 
you need ref_p times V. This precludes that your system is in 
equilibrium of course.






--
View this message in context: 
http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009062.html
Sent from the GROMACS Users Forum mailing list archive at Nabble.com.




--
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
* 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: Enthalpy Confusion

2013-06-11 Thread Jeffery Perkins
You should not use pV from g_energy though, as Michael explained, rather 
you need ref_p times V. This precludes that your system is in 
equilibrium of course.

That is what I had initially thought, then take that with the U to get
H?
or should i be doing  U+V*ref_p  = H?

But even so I still get different values for H then g_energy returns

example system gives H = -1168 kJ/mol and i find H = -725 kJ/mol either
way
I was mostly using pV as a stand in for H as it seems to be where the
difference between manual calculation and g_energy lies

thanks again for taking the time to assist,

Jeffery



--
View this message in context: 
http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009065.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


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


[gmx-users] Re: Enthalpy Confusion

2013-06-11 Thread Jeffery Perkins
 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. 

sorry I always have issues deciding how to talk about this stuff, so thanks
for putting up with my terrible notation =)

 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.

I hadn't really thought about that... but the system is a monatomic
Lennard-Jones particle (uncharged sigma = 0.35 nm epsilon = 2 kJ/mol mass =
40 amu) which should be a liquid at the conditions I was looking at, P=1000
bar, T = 300 K using phase diagram in: Equation of state for the
Lennard-Jones fluid, J. J. Nicolas et al., MOLECULAR PHYSICS, 1979, VOL. 37,
No. 5, 1429-1454

the U is -1600 and the pV is 880 when manually done, around 400 from
g_energy





--
View this message in context: 
http://gromacs.5086.x6.nabble.com/Enthalpy-Confusion-tp5009053p5009068.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