You have two ways to do that: 1. Calculate the total energy as a function of the volume for points close to the optimized lattice parameter and then fit the results with the Birch-Murnaghan equation of state. As you know the equation of state, you can obtain the pressure by doing p=-dEt/dV. So, with the obtained pressure you can get the enthalpy H=U+PV. In order to get the Gibbs energy, you must also obtain the phonon's density of states (DOS) and, by using statistical mechanics, get the partition function Z and finally, the Gibbs energy, as well as the Helmholtz free energy. 2. Do the first calculation mentioned in 1., and, by using the quasi-harmonic approximation, within the Gibbs2 code (http://azufre.quimica.uniovi.es/software.html), you can get all the thermodynamic functions, as well as the Debye temperature.

Thank You Dr. Horacio for nice discussion .

I will try this!

Good answer! Thank you, too.

Dear Dr. Horacio W. Leite Alves 

In context my earlier post I want to ask here, you suggested me two methods to calculate Gibbs free energy

1. Any publication or reference regarding method first

2. Any publication or reference regarding thermodynamic properties using gibbs2 code

3. I tried gibbs2 code it is giving many thermodynamic parameter including Gibbs free energy, EOS and many others. Through gibbs2 code calculating all thermodynamic parameter is very easy compare to method1.

I want to ask both methods will give same Gibbs free energy and Helmholtz free energy. Are both methods comparable ?

Thank you for your kind help and support !

Shilendra

calculations are performed at temperature T = 0 K, the Gibbs free energy is equal to the enthalpy H.We can obtain pressure versus volume by Murnaghan equation, and the Gibbs free energy expressed as (per formula unit).  At a given pressure a stable structure is one for which enthalpy has its lowest value and the transition pressures are calculated at which the enthalpies for the phases are equal. 

Dear Shilendra Kumar Sharma 

May be I am late here, but I would like to know how can I get the Gibbs free energy of a system in VASP using the Gibbs2 code. As you have already mentioned that you use the Gibbs2 code so how can I also use it? I have compiled the Makefile.inc file of the 'src' directory. Could you please give me some hints? 

Thank you very much.

Hi Arindam,

Just calculate E-V data from VASP around equilibrium volume.Convert unit, for Gibbs2 code  VASP gives energy in eV and volume in Angstrom but  Gibbs2 code  supports atomic unit for example energy unit Hartree.

Then execute by command and it will generate 3-4 output files in which one file will have thermodynamic data.

Shilendra

Dear Shilendra Kumar Sharma

Thank you very much for your concern for a very old post like this. As you can assume that I am beginner for VASP I suspect that the E-V data from VASP you mentioned is related to volume relaxation along with the ionic relaxation. Am I right? Here is my following two inquiry:  

1. If E-V data is related to the volume relaxation then I am afraid that for my system I only performed the ionic relaxation. So I would like to know in that case how can I generate the 'ing' file of the E-V data for my system.

2. Is there any other way to calculate the Gibbs free energy for ionic relaxation only within the Gibbs2 code?

Thank you very much again.

Hi Arindam,

After relaxation system which will be your optimized volume and E0 equilibrium energy. Now you have to vary your volume around your optimized volume for example your optimized volume is 300 Angs3 and now you will vary volume below optimized volume as 295 280 275 and above optimized volume 305 310 315. You can take any interval between volume like 300 302 304 check your system energy is varying significantly according that decide volume interval for energy generally 5 angs3 is OK.

How does write POSCAR for example

1.00000000000000

20.0746845202259685 0.0000000000000000 0.0000000000000000

0.0000000000000000 19.5534067032932128 0.0000000000000000

0.0000000000000000 0.0000000000000000 13.1181589915946599

These lattice vector for your initial POSCAR after optimization you find volume is 300 Angs3 for that your system energy is minimum.

Now CONTCAR to POSCAR and change POSCAR as for example 305 volume

-305

20.0746845202259685 0.0000000000000000 0.0000000000000000

0.0000000000000000 19.5534067032932128 0.0000000000000000

0.0000000000000000 0.0000000000000000 13.1181589915946599

So scaling factor multiply by required volume with -ve sign that after run you will find volume by command

grep volume OUTCAR it will give 305 volume so corresponding energy will be energy for volume 305.

Repeat for other volume and plot E-V data. It will be in parabola shape, you can check in literature such type of graph.

All the best

Shilendra

how do I calculate the molecule capture or uptake capacity using the software. I only know how to calculate the adsorption energy from the energy calculation but dont know how to calculate the quantity of the molecule adsorbed. Please any info on this ?

Entropy: A concept that is not a physical quantity

https://www.researchgate.net/publication/230554936_Entropy_A_concept_that_is_not_a_physical_quantity

Shilendra Kumar Sharma , Arindam Sannyal can you guys tell me in detail procedure to calculate the gibbs free energy for adsorbed H in the system? Since you guys have done this type of work , i am expecting a good reply with details.

We have already done the structural relaxation and DFPT frequency calculations. Gibbs Energy for an adsorbed hydrogen is calculated by: ΔGH = ΔEH + ΔEZPE - TΔS Where, * T is temperature in Kelvin * ΔS is entropy, * ΔEH is the binding enegy, * ΔEZPE (Zero-point-energy).

Can you any one tell me in details how to get all the parameters from the calculation means T, ΔS , ΔEH , ΔEZPE? Can anyone tell me the step by step calculations by sharing the INCAR file for this as I am new for such type of calculations?

Your help is highly appreciated.

I don't know whether it is required the DFPT calculation. If not, can you tell me the step by step process?

（http://azufre.quimica.uniovi.es/software.html）is no longer valid，Do you still have this Gibbss2 ，and upload a copy？Thanks！@Shilendra Kumar Sharma

Hi Sourav Bhakta

It seems that you're on the right track. Intending to know the right procedure to calculate the Gibbs free energy for adsorbed species, you must know the reaction (the elementary step) so you can find out the chemical potentials associated. Thus, the Gibbs free energy is deduced from it.

For example, suppose the adsorption of OH (minus) onto a surface. The reaction goes like * + OH(minus) -> *OH + e (Where * is the surface available site). Having said that, the Gibbs free energy is DeltaG = μ_*OH - μ_* - μ_H2O(l) + 1/2μ_H2 -eU_RHE (μ = chemical potential). The appearance of μ_H2O(l), μ_H2 and eU_RHE goes beyond the scope of this discussion. Essentially, μ_H2O(l), μ_H2 come into play to avoid the calculation of μ_OH(minus) and μ_e (electron chemical potential). This equation is used H2O(l) H(plus) + OH(minus) for the preceding purpose. On the other hand, eU_RHE shows up naturally when the H2(g) dissociation is used to determine the chemical potentials of H(plus) and μ_e.

Anyhow, it is simpler than look like and I had the intention to show you how important the reaction equation really is.

Back to the Gibbs free energy expression,

DeltaG = μ_*OH - μ_* - μ_H2O(l) + 1/2μ_H2 -eU_RHE

μ_*OH = E_*OH + ZPE_*OH - TS_*OH

μ_* = E_* (Equilibrium energy of the surface without the adsorbed specie; run an individual calculation for the surface only; take E0 from OSZICAR file)

μ_H2 = E_H2 (Energy formation of H2; run an individual calculation for two atoms of hydrogen; take the energy of equilibrium E0 in OSZICAR file once the relaxation is done)

μ_H2O(l) = E_H2O + ZPE_H2O - TS_H2O(l)

When it comes to VASP, the following quantities can be obtained

E_*OH, E_*, E_H2, E_H2O, ZPE_H2O, ZPE_*OH

For adsorption energies, E_*OH and E_*, after full relaxation, go to OSZICAR file and look for E0 in the last cycle.

For the formation energies, E_H2 and E_H2O, run a calculation individually for each one containing only the atoms of interest and repeat the same thing said above (go to OSZICAR [...]).

For the ZPE's, ZPE_H2O, ZPE_*OH, you run a frequency job (IBRION = 5 or IBRION =6, both yield the same results. Just use IBRION =5 if your system has no symmetry or IBRION =6, otherwise). After you reach accuracy (full relaxation), go to OUTCAR file, scroll down to the very end and look for the frequencies. Add them up and divide by 2 (Harmonic approximation is adopted). For IBRION usage, see INCAR file attached.

The OUTCAR file attached brings these frequencies (ZPE is roughly 0.33 eV, check it yourself).

1 f = 101.940914 THz 640.513654 2PiTHz 3400.382877 cm-1 421.593555 meV

X Y Z dx dy dz

0.000000 0.000000 0.000000 0 0 0

1.380811 1.952762 0.000000 0 0 0

2.755834 1.947668 2.122186 0 0 0

2.758525 3.900519 1.792325 0 0 0

0.056896 1.946426 4.120966 -0.227477 -0.000513 0.086096

1.883641 1.945120 4.449031 0.910349 0.002398 -0.334815

2 f = 33.059548 THz 207.719263 2PiTHz 1102.747807 cm-1 136.723241 meV

X Y Z dx dy dz

0.000000 0.000000 0.000000 0 0 0

1.380811 1.952762 0.000000 0 0 0

2.755834 1.947668 2.122186 0 0 0

2.758525 3.900519 1.792325 0 0 0

0.056896 1.946426 4.120966 0.133124 0.001101 0.169405

1.883641 1.945120 4.449031 -0.321840 0.000546 -0.921953

3 f = 15.653855 THz 98.356071 2PiTHz 522.156393 cm-1 64.739113 meV

X Y Z dx dy dz

0.000000 0.000000 0.000000 0 0 0

1.380811 1.952762 0.000000 0 0 0

2.755834 1.947668 2.122186 0 0 0

2.758525 3.900519 1.792325 0 0 0

0.056896 1.946426 4.120966 0.012187 -0.001006 0.981653

1.883641 1.945120 4.449031 -0.020218 0.002072 0.189193

4 f = 4.114407 THz 25.851581 2PiTHz 137.241842 cm-1 17.015812 meV

X Y Z dx dy dz

0.000000 0.000000 0.000000 0 0 0

1.380811 1.952762 0.000000 0 0 0

2.755834 1.947668 2.122186 0 0 0

2.758525 3.900519 1.792325 0 0 0

0.056896 1.946426 4.120966 0.963908 -0.028957 -0.015439

1.883641 1.945120 4.449031 0.259212 -0.022694 0.045810

5 f = 3.003346 THz 18.870577 2PiTHz 100.180824 cm-1 12.420834 meV

X Y Z dx dy dz

0.000000 0.000000 0.000000 0 0 0

1.380811 1.952762 0.000000 0 0 0

2.755834 1.947668 2.122186 0 0 0

2.758525 3.900519 1.792325 0 0 0

0.056896 1.946426 4.120966 -0.034327 -0.915027 0.000716

1.883641 1.945120 4.449031 -0.009187 -0.401814 -0.002948

6 f/i= 5.893296 THz 37.028671 2PiTHz 196.579197 cm-1 24.372703 meV

X Y Z dx dy dz

0.000000 0.000000 0.000000 0 0 0

1.380811 1.952762 0.000000 0 0 0

2.755834 1.947668 2.122186 0 0 0

2.758525 3.900519 1.792325 0 0 0

0.056896 1.946426 4.120966 -0.009317 0.402350 0.002617

1.883641 1.945120 4.449031 -0.000246 -0.915434 -0.000840

Obs.: Frequency job is performed on fully relaxed structures. Use selective dynamics to let only atoms of interest to move (See POSCAR file attached for reference).

Finally, TS_*OH and TS_H2O(l) are found in thermodynamics books or related articles.

I hope it helps,

Best,

Gabriel Vinicius Thank you very much, Sir, for your help. It will help me a lot. You did a great help. I will take care all of the points you mentioned.

Thanks a lot, Sir.

Linshuang Zhang It is available here https://aoterodelaroza.github.io/gibbs2/