1 #!/usr/bin/env python
   2 
   3 """Bulk Al(fcc) equation of state
   4 
   5 This script calculates the total energy of an fcc Al unit cell
   6 at several different volumes for fitting to an equation of state
   7 
   8 john kitchin May 20, 2005
   9 """
  10 import os
  11 from Dacapo import Dacapo
  12 from ASE import Atom,ListOfAtoms
  13 from ASE.Utilities.GeometricTransforms import SetUnitCellVolume
  14 
  15 bulk = ListOfAtoms([Atom('Al', (0, 0, 0))] )
  16 
  17 b = 4.05
  18 bulk.SetUnitCell([(0,b/2,b/2),
  19                   (b/2,0,b/2),
  20                   (b/2,b/2,0)])
  21 
  22 calc = Dacapo(kpts=(4,4,4),               
  23 		  planewavecutoff=350,
  24 		  nbands=6,                   
  25 		  usesymm=False)
  26     
  27 bulk.SetCalculator(calc)
  28 
  29 initial_volume = bulk.GetUnitCellVolume()
  30 
  31 energies = []
  32 volumes = []
  33 
  34 for fraction in [0.9, 0.95, 1.0, 1.05, 1.1]:
  35 
  36     ncfile = 'bulkAl_%1.2f.nc' % fraction
  37     #only run calculations if the files don't exist
  38     if not os.path.exists(ncfile):
  39         new_volume = fraction * initial_volume
  40 
  41         SetUnitCellVolume(bulk,new_volume)
  42 
  43         calc.SetNetCDFFile(ncfile)
  44 
  45         bulk.GetPotentialEnergy()
  46 
  47     atoms = Dacapo.ReadAtoms(ncfile)
  48     energies.append(atoms.GetPotentialEnergy())
  49     volumes.append(atoms.GetUnitCellVolume())
  50 
  51 from ASE.Utilities.EquationOfState import EquationOfState
  52 eos = EquationOfState('Murnaghan',volumes,energies)
  53 print eos
  54 eos.GetPlot()
  55 eos.SavePlot('Al_murn.png')
  56 
  57 V0 = eos.GetV0()
  58 
  59 '''
  60 the volume of a primitive fcc cell is equal to
  61 V=1/4*a^3 where a is the lattice constant.
  62 '''
  63 a = (V0*4)**(1./3.)
  64 print 'lattice constant = %1.2f angstroms' % a
  65 
  66 '''
  67 Typical output:
  68 
  69 V0(A^3)       B(Gpa)   E0 (eV)    pressure(GPa) Murnaghan
  70  16.78        85.91     -56.4668       -0.00
  71 chisq = 0.0000
  72 lattice constant = 4.05 angstroms
  73 '''
  74 
  75     
  76