Finding the volume of a unit cell at a fixed pressure

| categories: uncategorized | tags: | View Comments

A typical unit cell optimization in DFT is performed by minimizing the total energy with respect to variations in the unit cell parameters and atomic positions. In this approach, a pressure of 0 GPa is implied, as well as a temperature of 0K. For non-zero pressures, the volume that minimizes the total energy is not the same as the volume at P=0.

Let \(x\) be the unit cell parameters that can be varied. For P ≠ 0, and T = 0, we have the following

\(G(x; p) = E(x) + p V(x)\)

and we need to minimize this function to find the groundstate volume. We will do this for fcc Cu at 5 GPa of pressure. We will assume there is only one degree of freedom in the unit cell, the lattice constant. First we get the \(E(x)\) function, and then add the analytical correction.

from jasp import *
from ase import Atom, Atoms
from ase.utils.eos import EquationOfState

LC = [3.5, 3.55, 3.6, 3.65, 3.7, 3.75]
volumes, energies = [], []
ready = True

P = 5.0 / 160.2176487  # pressure in eV/ang**3

for a in LC:
    atoms = Atoms([Atom('Cu',(0, 0, 0))],
              cell=0.5 * a*np.array([[1.0, 1.0, 0.0],
                                     [0.0, 1.0, 1.0],
                                     [1.0, 0.0, 1.0]]))

    with jasp('../bulk/Cu-{0}'.format(a),
              xc='PBE',
              encut=350,
              kpts=(8,8,8),
              atoms=atoms) as calc:

        try:
            e = atoms.get_potential_energy()
            energies.append(e)
            volumes.append(atoms.get_volume())
        except (VaspSubmitted, VaspQueued):
            ready = False

if not ready:
    import sys; sys.exit()

import numpy as np
energies = np.array(energies)
volumes = np.array(volumes)

eos = EquationOfState(volumes, energies)
v0, e0, B = eos.fit()
print 'V0 at 0 GPa = {0:1.2f} ang^3'.format(v0)

eos5 = EquationOfState(volumes, energies + P * volumes)
v0_5, e0, B = eos5.fit()
print 'V0 at 5 GPa = {0:1.2f} ang^3'.format(v0_5)
V0 at 0 GPa = 12.02 ang^3
V0 at 5 GPa = 11.62 ang^3

You can see here that apply pressure decreases the equilibrium volume, and increases the total energy.

Copyright (C) 2013 by John Kitchin. See the License for information about copying.

org-mode source

blog comments powered by Disqus