Commit 82b152b6 by Francois Gygi

reorganized utility scripts


git-svn-id: http://qboxcode.org/svn/qb/trunk@1415 cba15fb0-1239-40c8-b417-11db7ca47a34
parent 755db8fc
#!/usr/bin/python
# qbox_dos.py: extract dos from Qbox output
# generate dos plot in gnuplot format
# qbox_dos.py: extract electronic DOS from Qbox output
# generate DOS plot in gnuplot format
# use: qbox_dos.py emin emax width file.r
# emin, emax: bounds of plot in [eV]
# width: gaussian broadening in [eV]
# the DOS is accumulated separately for each spin
import xml.sax
......
#!/usr/bin/python
# qbox_msd.py: compute mean-square displacement in an MD simulation
# generate plot of <r^2>(t) in gnuplot format
# use: qbox_msd.py species file.r [file.r ...]
import xml.sax
import sys
import math
if len(sys.argv) < 3:
print "use: ",sys.argv[0]," species file [file ...]"
sys.exit()
species = sys.argv[1]
# Qbox output handler to extract and process data
class QboxOutputHandler(xml.sax.handler.ContentHandler):
def __init__(self):
self.step = 0
self.readData = 0
self.readAtom = 0
self.tau0 = []
def startElement(self, name, attributes):
if name == "atomset":
self.tau=[]
elif name == "atom":
self.readAtom = (attributes["species"] == species)
elif (name == "position") & self.readAtom:
self.readData = 1
self.buffer = ""
def characters(self, data):
if self.readData:
self.buffer += data
def endElement(self, name):
if (name == "position") & self.readAtom:
self.readData = 0
pos = self.buffer.split()
x = float(pos[0])
y = float(pos[1])
z = float(pos[2])
self.tau.append([x,y,z])
elif name == "atomset":
if self.step == 0:
# copy initial positions to tau0
self.tau0 = self.tau
# compute square displacement
disp2sum = 0.0
for i in range(len(self.tau)):
dx = self.tau[i][0]-self.tau0[i][0]
dy = self.tau[i][1]-self.tau0[i][1]
dz = self.tau[i][2]-self.tau0[i][2]
disp2sum += dx*dx + dy*dy + dz*dz
print self.step,'%12.6f'%(disp2sum/len(self.tau))
self.step += 1
print "# ",species
parser = xml.sax.make_parser()
handler = QboxOutputHandler()
parser.setContentHandler(handler)
for i in range(len(sys.argv)-2):
infile = sys.argv[i+2]
parser.parse(infile)
print
print
#!/usr/bin/python
# qbox_repair_h2o.py: repair broken h2o molecules in a Qbox sys file
# move hydrogen atoms across periodic boundaries to repair molecules
# use: qbox_repair_h2o.py file.sys
import sys
import math
olist = []
hlist = []
def distance(a,b,sx,sy,sz):
return math.sqrt((a[3]-b[3]-sx)**2+(a[4]-b[4]-sy)**2+(a[5]-b[5]-sz)**2)
f = open(sys.argv[1])
for line in f:
l = line.split()
if ( l[0] == "set" ) & ( l[1] == "cell" ):
a_cell = float(l[2])
b_cell = float(l[6])
c_cell = float(l[10])
print line
elif ( l[0] == "species" ):
print line
elif ( l[0] == "atom" ) & ( l[2] == "oxygen" ):
olist.append([l[0],l[1],l[2],float(l[3]),float(l[4]),float(l[5])])
elif ( l[0] == "atom" ) & ( l[2] == "hydrogen" ):
hlist.append([l[0],l[1],l[2],float(l[3]),float(l[4]),float(l[5])])
for h in hlist:
# find nearest oxygen atom in olist
mindist = 1.e10
for o in olist:
# compute minimal distance o-h
for sx in [-a_cell,0,a_cell]:
for sy in [-b_cell,0,b_cell]:
for sz in [-c_cell,0,c_cell]:
d = distance(o,h,sx,sy,sz)
# print "dist(",sx,sy,sz,") = ",d
if d < mindist:
mindist = d
sx_min = sx
sy_min = sy
sz_min = sz
nearest_o = o;
# print "min shift is: ",sx_min,sy_min,sz_min
if ( sx_min != 0 ) | ( sy_min != 0 ) | ( sz_min != 0 ):
print "# current ",h[1]," at ", h[3],h[4],h[5]
print "# nearest O is at ", nearest_o[3],nearest_o[4],nearest_o[5]
print "# move ",h[1]," by ", sx_min, sy_min, sz_min
h[3] += sx_min
h[4] += sy_min
h[5] += sz_min
for o in olist:
print o[0],o[1],o[2],'%10.5f'%o[3],'%10.5f'%o[4],'%10.5f'%o[5]
for h in hlist:
print h[0],h[1],h[2],'%10.5f'%h[3],'%10.5f'%h[4],'%10.5f'%h[5]
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment