Some scripts for file format conversion
Jump to navigation
Jump to search
QE use the xsf file format produced by pwi2xsf.sh and pwo2xsf.sh. They can be opened only by Xcrysden. If you need to open your structure with other tools you can convert it to a cif file using the BABEL converter (http://openbabel.org/wiki/Main_Page) the only problem is that they do it wrong... the converter is unable to write the atom elements of the xsf files to the cif file.
Here is a small script that correct the errors of the wrong cif file obtained with BABEL.
#!/usr/bin/python2
#This script use the wrong output obtained from Babel when it tries to convert a xsf file to a cif file and convert it to a correct cif file
#tha can be opened by Material studio
import os
import sys
import math
#function to count number of lines in file
def file_len(fname):
with open(fname,"r") as f:
for i, l in enumerate(f):
pass
return i + 1
f.close
## program start here ##
#EDIT ONLY HERE (below)
name_prefix="ZnFeO2_spinel_001.out2"
#EDIT ONLY HERE (above)
filename=name_prefix+".xsf" #correct xsf
filename1=name_prefix+".xsf.cif1" #wrong cif from Babel
output=name_prefix+".xsf.cif" #right cif at the output
#Step 1 read from xsf the Element names
nmax=file_len(filename)
f0=open(filename,"r")
i=0 #used as line counter
T=0 #used as logical switch
k=0 #used for PRIMVEC
A=[] #store name
while i<nmax:
line=f0.readline().split()
#print line
if 'PRIMVEC' in line:
i=i+1
T=1
elif 'PRIMCOORD' in line:
i=i+1
T=2
elif T==1:
if k==0:
a=float(line[0])
b=a
k=k+1
elif k==1:
k=k+1
elif k==2:
c=float(line[2])
k=k+1
T=0
i=i+1
elif T==2:
numatoms=int(line[0])
T=3
i=i+1
elif T==3:
A.append([line[0]])
i=i+1
else:
i=i+1
f0.close
#Step 2 read the wrong cif file and write the right cif file
f1=open(filename1,"r")
out=open(output,"w")
i=0 #used as line counter
while i<1:
line=f1.readline()
print >> out, line,
if '_atom_site_occupancy' in line:
i=1 #coordinates start here, exit while loop
for i in range(numatoms):
line=f1.readline().split()
print >> out, " "+str(A[i][0])+str(i), A[i][0], line[2], line[3], line[4], line[5]
f1.close
out.close
Lazy script to convert cif file to QE input (must be updated)
#!/usr/bin/python2
#This script prepare the input file for quantum espresso ATOMIC_POSITIONS from a .cif file for a geometry optimisation
#Not for cell optimisation. Replace relax by vc-relax for that and add the CELL section
import os
import sys
import math
#function to count number of lines in file
def file_len(fname):
with open(fname,"r") as f:
for i, l in enumerate(f):
pass
return i + 1
f.close
## program start here ##
name_prefix="ZnFeO2_spinel_super222_001"
filename="fake.list" #changing the whole code is more bothersome than rewriting the whole code...
filename1=name_prefix+".cif"
#xyzname=name_prefix+".xyz" #not needed
output=name_prefix+".in"
#xyz=open(xyzname,"r")
#numatoms=int(xyz.readline())
#xyz.close
#numatoms=128 #Is now extracted from xyz file
#count number of atoms in cif file
f1=open(filename1,"r")
nmax1=file_len(filename1)
i=0
j=0
shift=1 #by default input list start at 1
for line in f1:
if '_atom_site_occupancy' in line:
i=1 #coordinates start here
elif i == 1:
j=j+1
if 'Fe0' in line:
print 'Fe0 found'
shift=0 #sometimes the list start at 0. I do not know why
#print i,j,numatoms
f1.close
numatoms=j #taken from the cif file, xyz file not needed
print numatoms,shift
nmax=file_len(filename)
#print nmax
#create the list
uplist=[]
downlist=[]
freezelist=[]
T=0 #used for list type
i=0 #used for line index
f=open(filename,"r") #should be read-only
while i<nmax:
line=f.readline().split()
#print T
#print line
if 'up' in line:
i=i+1
T=1
elif 'down' in line:
i=i+1
T=2
elif 'freeze' in line:
i=i+1
T=3
elif T==1:
b=[int(l) for l in line]
uplist = uplist + b
i=i+1
elif T==2:
b=[int(l) for l in line]
downlist = downlist + b
i=i+1
elif T==3:
b=[int(l) for l in line]
freezelist = freezelist + b
i=i+1
else:
i=i+1
#print i
f.close()
#Check if number is in both up & down list
for i in range(len(uplist)):
if uplist[i] in downlist:
print "BUG in list"
#print uplist
#print downlist
#print freezelist
f1=open(filename1,"r")
nmax1=file_len(filename1)
i=0
while i<1:
line=f1.readline()
if '_cell_length_a' in line:
a=float(line.split()[1])
if '_cell_length_b' in line:
b=float(line.split()[1])
if '_cell_length_c' in line:
c=float(line.split()[1])
if '_cell_angle_alpha' in line:
cosBC=math.cos(math.radians(float(line.split()[1])))
if -1.e-15 <cosBC < 1.e-15:
cosBC=0
if '_cell_angle_beta' in line:
cosAC=math.cos(math.radians(float(line.split()[1])))
if -1.e-15 < cosAC < 1.e-15:
cosAC=0
if '_cell_angle_gamma' in line:
cosAB=math.cos(math.radians(float(line.split()[1])))
if -1.e-15 < cosAB < 1.e-15:
cosAB=0
if '_atom_site_occupancy' in line:
i=1 #coordinates start here, exit while loop
out=open(output,"w")
print >> out, " &CONTROL"
print >> out, " calculation = 'relax' ,"
print >> out, " restart_mode = 'from_scratch' ,"
print >> out, " wf_collect = .true. ,"
print >> out, " prefix =","'"+name_prefix+"'", ","
print >> out, " verbosity = 'high' ,"
print >> out, " etot_conv_thr = 1.0D-8 ,"
print >> out, " forc_conv_thr = 1.0D-4 ,"
print >> out, " tstress = .true. ,"
print >> out, " tprnfor = .true. ,"
print >> out, " dipfield = .true. ,"
print >> out, " nstep = 400 ,"
print >> out, " max_seconds = 171000 , !47h30"
print >> out, " /"
print >> out, "&SYSTEM"
print >> out, " ibrav = 14,"
print >> out, " A =",a, ","
print >> out, " B =",b, ","
print >> out, " C =",c, ","
print >> out, " cosAB =",cosAB, ","
print >> out, " cosAC =",cosAC, ","
print >> out, " cosBC =",cosBC, ","
print >> out, " nat =",numatoms,","
print >> out, " ntyp = 3,"
print >> out, " ecutwfc = 40.0 , !J. Chem. Phys. 138, 194709" #34
print >> out, " ecutrho = 320.0 ,!J. Chem. Phys. 138, 194709" #323
print >> out, " input_dft = 'PBE' ,"
print >> out, " occupations = 'smearing' ,"
print >> out, " degauss = 0.02 ,"
print >> out, " smearing = 'gaussian' ,"
print >> out, " nspin = 2 ,"
print >> out, " starting_magnetization(1) = 0.0,"
print >> out, " starting_magnetization(2) = -0.0,"
print >> out, " starting_magnetization(3) = 0.0,"
print >> out, " lda_plus_u = .true. ,"
print >> out, " lda_plus_u_kind = 0 , !type 1 not implemented, ignores J"
print >> out, " Hubbard_U(1) = 0.0,"
print >> out, " Hubbard_U(2) = 4.2,"
print >> out, " Hubbard_J0(1) = 0.0, !Is ignored!"
print >> out, " Hubbard_J0(2) = 0.0, !Is ignored!"
print >> out, " !vdw_corr = Grimme-D2, ! DFT-D2 dispersion correction"
print >> out, "/"
print >> out, "&ELECTRONS"
print >> out, " electron_maxstep = 800,"
print >> out, " conv_thr_init = 1e-4 ,"
print >> out, " conv_thr = 1e-9 ,"
print >> out, " startingpot = 'atomic' ,"
print >> out, " startingwfc = 'random' ,"
print >> out, " adaptive_thr = .true. ,"
print >> out, " mixing_beta = 0.514,"
print >> out, " diagonalization = 'david' ,"
print >> out, "/"
print >> out, "&IONS"
print >> out, " ion_dynamics = 'bfgs' ,"
print >> out, " trust_radius_min = 1.D-5 ,"
print >> out, "/"
print >> out, "ATOMIC_SPECIES"
print >> out, " Zn 65.39000 Zn_pbe_v1.uspp.F.UPF"
print >> out, " Fe 55.84500 Fe.pbe-sp-van_ak.UPF"
print >> out, " O 15.99990 O.pbe-van_ak.UPF"
print >> out, "ATOMIC_POSITIONS crystal"
for i in range(numatoms):
line=f1.readline().split()
if (i+shift) in uplist and (i+shift) not in freezelist:
print >> out, line[1]+"1", line[2], line[3], line[4]
elif (i+shift) in uplist and (i+shift) in freezelist:
print >> out, line[1]+"1", line[2], line[3], line[4], 0 ,0, 0
elif (i+shift) in downlist and (i+shift) not in freezelist:
print >> out, line[1]+"2", line[2], line[3], line[4], 1 ,1, 1
elif (i+shift) in downlist and (i+shift) in freezelist:
print >> out, line[1]+"2", line[2], line[3], line[4], 0 ,0, 0
else:
if (i+shift) in freezelist:
print >> out, line[1], line[2], line[3], line[4], 0 ,0, 0
else:
print >> out, line[1], line[2], line[3], line[4]
print >> out, "K_POINTS Gamma"
#closing file
f1.close
f.close
out.close