# Question about performing thermodynamical integration using LAMMPS to calculate free energy of an solid system with EAM potential

Hi, everyone!
I want to calculate the free energy of a solid with EAM potential using lammps peforming thermodynamical integration method, however the EAM potential is not coded on fix adapt & compute ti combo, and I am incapable of programing.

I have noticed that someone has already met the same problem and Rodrigo Freitas may provide a fix command to fix it(http://lammps.sandia.gov/threads/msg28675.html),
So, dear Steve and Rodrigo,could you add this fix command in the new version?

Yang LI2013年01月29日 星期二

What attribute of the EAM potential are you wanting to adapt?
The EAM potential is input and stored in tabulated form, so
there isn't a simple parameter like "epsilon" that could be time
varied.

Steve

I have been working on a thermodynamical integration framework using
the Hamiltonian-Switching method. My package can deal with solid and
liquid phases. There are a bunch of additions I made to my LAMMPS
repo, such as a new pair/hybrid method that can switch gradually
between two different potentials and minor stuff like a restore_state
command and a new function that returns the mass of a certain species
in the simulation (rather than the mass of a certain group). The
package also includes a python analysis script which process the data
and output free energy surfaces for alloys (as function of temperature
and concentration).
Since I work in a National Lab I cannot simply publish the whole thing
online right now (I'll have to go through a release process). But I'm
sure something can be arranged in terms of a collaboration.
Let me know if you are interested.
Daniel

Steve,

sometimes it is necessary to perform an interpolation between two
different potentials, this is slightly different from simply changing
one parameter of the potential.

What I've done is to write a fix that perform an linear interpolation
between the interatomic potential and an potential where each atom is
connected to its equilibrium position by a spring, so my effective
potential is

U(x) = (1-x) * U_int + x * U_harm

where U_int is the interatomic potential, U_harm is the harmonic
potential and x is the parameter that I change during my simulation.

In this sense the thermodynamic integration that I've programmed (and
probably the same for Daniel Schwen's code) is more general than the
current fix adapt + compute ti, it allows the complete interpolation
between any potential and a 'reference system', as reference system I
used the harmonic crystal, in other cases it is possible to use a
Gaussian potential or a ideal gas. In the end of the day all we want is
to compute the free energy difference between something that we know the
free energy and the system with the interatomic potential.

There are other details about the method (eg: it is very useful to be
able to perform this interpolation in both directions) that require some
attention.

I have been eager to see such kind of methods implemented in the
official version of LAMMPS, just like Daniel I worked out my own LAMMPS
modifications and have been successfully using them over the last year
(I have some analytical results to gain confidence on the code). I don't
have any issues with the publication of the code, despite my email being
@llnl.gov all these codes where developed during my Master's in Brazil,
but I recognize that it may need some changes and discussion to be sure
that it works according to the LAMMPS code standards. I've been
interest on it we can talk about implementing this on LAMMPS.

Rodrigo Freitas

What I've done is to write a fix that perform an linear interpolation
between the interatomic potential and an potential where each atom is
connected to its equilibrium position by a spring, so my effective
potential is

U(x) = (1-x) * U_int + x * U_harm

where U_int is the interatomic potential, U_harm is the harmonic
potential and x is the parameter that I change during my simulation.

In this sense the thermodynamic integration that I've programmed (and
probably the same for Daniel Schwen's code) is more general than the

I'm using fix spring/self and a new fix force/mul to do the switching
to the einstein crystal. Pretty much everything else can be done in
the input script (calculating the spring constants from Debye
temperatures, obtaining the Hamiltonians of the pure systems etc.)
Of course this leaves the switching to an ideal gas to be able to
obtain the free energy of a liquid. As I mentioned this I do with a
new pair style hybrid/weighted. Coincidentally I've been doing this
for about a year now myself and had some good success (good results
for melting points, free energy calculation is independent of
switching temperature etc.).
It might be interesting to collaborate on making this fit for
integration into LAMMPS. If the whole budget mess gets sorted out and
we finally receive DoE approval, I'll be attending TMS in a month,
giving a talk on my thermo package.

I'll take this off-line with the two of you ...

Steve