The value of chemical potential for pure water calculated by fix-nvt and fix-gcmc combination

Dear Lammps users and developers,

I have been learning fix-gcmc for a few days. Now I am testing a target system: pure water.

I use fix-nvt + fix-gcmc to calculate the SPC/E water at 300 K. There are totally about 800 waters in the simulation box. When I set mu=-4.8 kcal/mol, I got a density of water, 0.990. Thus it is estimated that for density of 1.00, the chemical potential should be less than -4.8 kcal/mol.

However, my result of mu for water at 298 K is different with some publications. See below.

ref.1 DYNAMIC SIMULATIONS OF WATER AT CONSTANT CHEMICAL-POTENTIAL
http://scitation.aip.org/content/aip/journal/jcp/96/2/10.1063/1.462169

The value of -8.0 kcal/mol is given.

ref.2 Molecular simulation of aqueous electrolytes: Water chemical potential results and Gibbs-Duhem equation consistency tests
http://scitation.aip.org/content/aip/journal/jcp/139/12/10.1063/1.4821153
The value of -240.0 kJ/mol (-56.7 kcal/mol) is given.

My question is:

What is the difference between ref. 1 and ref. 2?

Why fix gcmc give a even different value?

What is the correct mu value by gix-gcmc?

Somebody must have calculated the mu for pure water. Please help!

Best wishes,

Yongbiao

Dear Lammps users and developers,

     I have been learning fix-gcmc for a few days. Now I am testing a target
system: pure water.
     I use fix-nvt + fix-gcmc to calculate the SPC/E water at 300 K. There
are totally about 800 waters in the simulation box. When I set mu=-4.8
kcal/mol, I got a density of water, 0.990. Thus it is estimated that for
density of 1.00, the chemical potential should be less than -4.8 kcal/mol.

     However, my result of mu for water at 298 K is different with some
publications. See below.

   ref.1 DYNAMIC SIMULATIONS OF WATER AT CONSTANT CHEMICAL-POTENTIAL

http://scitation.aip.org/content/aip/journal/jcp/96/2/10.1063/1.462169
         The value of -8.0 kcal/mol is given.

  ref.2 Molecular simulation of aqueous electrolytes: Water chemical
potential results and Gibbs-Duhem equation consistency tests

http://scitation.aip.org/content/aip/journal/jcp/139/12/10.1063/1.4821153
        The value of -240.0 kJ/mol (-56.7 kcal/mol) is given.

My question is:

   What is the difference between ref. 1 and ref. 2?

you have the two papers, you can make the comparison of the two papers yourself.

   Why fix gcmc give a even different value?

why not? are your simulation conditions and force field settings
*exactly* like in any of the two papers? do you have any estimates
about convergence and degree of equilibration? are you certain, that
there are no typos or rounding issues anywhere?

   What is the correct mu value by gix-gcmc?

   Somebody must have calculated the mu for pure water. Please help!

on the level of classical MD simulations, there is not going to be
*the* correct value for pure water, as there are a gazillion of
different water potentials, each being different from the others in
some more-or-less subtle ways. before even thinking about reproducing
numbers for grand canonical ensemble simulations, you have to confirm
that you can reproduce data for other, simpler ensembles accurately.

axel.

Dear Axel,

Thanks for your help.

I have tested the fix-npt simulations and a density of ~1.0 is perfectly found for SPC/E water at 300K.

As for the chemical potential for pure water at 298 K, it is found around -4.6 kcal/mol. As the manual says, the input is the total chemical potential for water. I have read some papers based on SPCE or SPC or TIP4P models, their results are also around -5.7.

I feel the results should be consistent to some extent… It can not be an arbitrary value…

I am very confused now.

Best,

Yongbiao

Dear Axel,

     Thanks for your help.

      I have tested the fix-npt simulations and a density of ~1.0 is
perfectly found for SPC/E water at 300K.

density is only one of many properties and not a very sensitive one.

      As for the chemical potential for pure water at 298 K, it is found
around -4.6 kcal/mol. As the manual says, the input is the total chemical
potential for water. I have read some papers based on SPCE or SPC or TIP4P
models, their results are also around -5.7.

       I feel the results should be consistent to some extent....... It can
not be an arbitrary value...

no arbitrary, but there are many factors that can impact it. outside
of the issues i mentioned before there are:
- system size
- the exact variant of the water potential (e.g. whether it was
"tuned" for a specific kind of long-range electrostatics)
- what kind of cutoffs (particularly LJ) and if there is an energy
correction for the tail of the LJ potential
- choice of simulation parameters (e.g. for fix gcmc)
- how well a simulation is equilibrated and converged.
- whether you have a COM drift.
and many more.

there also is the fact, that a significant chunk of bad simulation
data gets published (some people claim it is more these days than what
it used to be), and that the threshold for what is considered a "good"
result was different when there was much less computer power available
to people.

       I am very confused now.

welcome to the wonderful world of computer simulations, where things
get worse, the closer you look. i think at this point, you should
really discuss with your adviser and your senior colleagues. you seem
to have a very optimistic approach to what kind of information a
simulation can provide and what is required to make sure what you
compute is consistent with what you are expecting. sometimes subtle
differences can have a large impact, sometimes gross mistakes have
minimal consequences. none of that is specific to LAMMPS, but just the
reality of the "business". in general you need to learn that basic
assumption should not be "i did everything correct, so the results
must be correct", but rather "what else could i have done wrong or
what could i have overlooked?"

axel.