git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@13881 f3b2605a-c512-4ea7-a41b-209d697bcdaa

This commit is contained in:
sjplimp 2015-08-14 14:56:57 +00:00
parent 36315f3a96
commit 24bf97204b
3 changed files with 262 additions and 164 deletions

View File

@ -17,8 +17,8 @@ d = dump("dump.*") wildcard expands to multiple files
d = dump("dump.*",0) two args = store filenames, but don't read d = dump("dump.*",0) two args = store filenames, but don't read
incomplete and duplicate snapshots are deleted incomplete and duplicate snapshots are deleted
if atoms have 5 or 8 columns, assign id,type,x,y,z (ix,iy,iz)
atoms will be unscaled if stored in files as scaled atoms will be unscaled if stored in files as scaled
self-describing column names assigned
time = d.next() read next snapshot from dump files time = d.next() read next snapshot from dump files
@ -28,7 +28,7 @@ time = d.next() read next snapshot from dump files
return -1 if no snapshots left or last snapshot is incomplete return -1 if no snapshots left or last snapshot is incomplete
no column name assignment or unscaling is performed no column name assignment or unscaling is performed
d.map(1,"id",3,"x") assign names to atom columns (1-N) d.map(1,"id",3,"x") assign names to columns (1-N)
not needed if dump file is self-describing not needed if dump file is self-describing
@ -61,7 +61,7 @@ d.write("file",head,app) write selected steps/atoms to dump file
d.scatter("tmp") write selected steps/atoms to multiple files d.scatter("tmp") write selected steps/atoms to multiple files
write() can be specified with 2 additional flags write() can be specified with 2 additional flags
headd = 0/1 for no/yes snapshot header, app = 0/1 for write vs append head = 0/1 for no/yes snapshot header, app = 0/1 for write vs append
scatter() files are given timestep suffix: e.g. tmp.0, tmp.100, etc scatter() files are given timestep suffix: e.g. tmp.0, tmp.100, etc
d.scale() scale x,y,z to 0-1 for all timesteps d.scale() scale x,y,z to 0-1 for all timesteps
@ -115,10 +115,9 @@ fx,fy,... = d.vecs(1000,"fx","fy",...) return vector(s) for timestep N
vecs() returns vectors with one value for each selected atom in the timestep vecs() returns vectors with one value for each selected atom in the timestep
index,time,flag = d.iterator(0/1) loop over dump snapshots index,time,flag = d.iterator(0/1) loop over dump snapshots
time,box,atoms,bonds,tris = d.viz(index) return list of viz objects time,box,atoms,bonds,tris,lines = d.viz(index) return list of viz objects
d.atype = "color" set column returned as "type" by viz d.atype = "color" set column returned as "type" by viz
d.extra("dump.bond") read bond list from dump file d.extra(obj) extract bond/tri/line info from obj
d.extra(data) extract bond/tri/line list from data
iterator() loops over selected timesteps iterator() loops over selected timesteps
iterator() called with arg = 0 first time, with arg = 1 on subsequent calls iterator() called with arg = 0 first time, with arg = 1 on subsequent calls
@ -126,17 +125,21 @@ d.extra(data) extract bond/tri/line list from data
time = timestep value time = timestep value
flag = -1 when iteration is done, 1 otherwise flag = -1 when iteration is done, 1 otherwise
viz() returns info for selected atoms for specified timestep index viz() returns info for selected atoms for specified timestep index
can also call as viz(time,1) and will find index of preceding snapshot
time = timestep value time = timestep value
box = [xlo,ylo,zlo,xhi,yhi,zhi] box = \[xlo,ylo,zlo,xhi,yhi,zhi\]
atoms = id,type,x,y,z for each atom as 2d array atoms = id,type,x,y,z for each atom as 2d array
bonds = id,type,x1,y1,z1,x2,y2,z2,t1,t2 for each bond as 2d array bonds = id,type,x1,y1,z1,x2,y2,z2,t1,t2 for each bond as 2d array
if bonds() was used to define bonds, else empty list if extra() used to define bonds, else NULL
tris = id,type,x1,y1,z1,x2,y2,z2,x3,y3,z3,nx,ny,nz for each tri as 2d array tris = id,type,x1,y1,z1,x2,y2,z2,x3,y3,z3,nx,ny,nz for each tri as 2d array
if extra() was used to define tris, else empty list if extra() used to define tris, else NULL
lines = id,type,x1,y1,z1,x2,y2,z2 for each line as 2d array lines = id,type,x1,y1,z1,x2,y2,z2 for each line as 2d array
if extra() was used to define lines, else empty list if extra() used to define lines, else NULL
atype is column name viz() will return as atom type (def = "type") atype is column name viz() will return as atom type (def = "type")
extra() stores list of bonds/tris/lines to return each time viz() is called extra() extracts bonds/tris/lines from obj each time viz() is called
obj can be data object for bonds, cdata object for tris and lines,
bdump object for bonds, tdump object for tris, ldump object for lines.
mdump object for tris
""" """
# History # History
@ -153,6 +156,7 @@ d.extra(data) extract bond/tri/line list from data
# increment = 1 if reading snapshots one-at-a-time # increment = 1 if reading snapshots one-at-a-time
# nextfile = which file to read from via next() # nextfile = which file to read from via next()
# eof = ptr into current file for where to read via next() # eof = ptr into current file for where to read via next()
# scale_original = 0/1/-1 if coords were read in as unscaled/scaled/unknown
# nsnaps = # of snapshots # nsnaps = # of snapshots
# nselect = # of selected snapshots # nselect = # of selected snapshots
# snaps = list of snapshots # snaps = list of snapshots
@ -161,20 +165,22 @@ d.extra(data) extract bond/tri/line list from data
# tselect = class for time selection # tselect = class for time selection
# aselect = class for atom selection # aselect = class for atom selection
# atype = name of vector used as atom type by viz extract # atype = name of vector used as atom type by viz extract
# bondflag = 0 if no bonds, 1 if they are defined statically # bondflag = 0 if no bonds, 1 if they are defined statically, 2 if dynamic
# bondlist = static list of bonds to viz() return for all snapshots # bondlist = static list of bonds to return w/ viz() for all snapshots
# only a list of atom pairs, coords have to be created for each snapshot
# triflag = 0 if no tris, 1 if they are defined statically, 2 if dynamic # triflag = 0 if no tris, 1 if they are defined statically, 2 if dynamic
# trilist = static list of tris to return via viz() for all snapshots # trilist = static list of tris to return w/ viz() for all snapshots
# lineflag = 0 if no lines, 1 if they are defined statically # lineflag = 0 if no lines, 1 if they are defined statically, 2 if dynamic
# linelist = static list of lines to return via viz() for all snapshots # linelist = static list of lines to return w/ viz() for all snapshots
# objextra = object to get bonds,tris,lines from dynamically
# Snap = one snapshot # Snap = one snapshot
# time = time stamp # time = time stamp
# tselect = 0/1 if this snapshot selected # tselect = 0/1 if this snapshot selected
# natoms = # of atoms # natoms = # of atoms
# boxstr = format string after BOX BOUNDS, if it exists
# triclinic = 0/1 for orthogonal/triclinic based on BOX BOUNDS fields
# nselect = # of selected atoms in this snapshot # nselect = # of selected atoms in this snapshot
# aselect[i] = 0/1 for each atom # aselect[i] = 0/1 for each atom
# xlo,xhi,ylo,yhi,zlo,zhi = box bounds (float) # xlo,xhi,ylo,yhi,zlo,zhi,xy,xz,yz = box bounds (float)
# atoms[i][j] = 2d array of floats, i = 0 to natoms-1, j = 0 to ncols-1 # atoms[i][j] = 2d array of floats, i = 0 to natoms-1, j = 0 to ncols-1
# Imports and external programs # Imports and external programs
@ -184,11 +190,11 @@ from os import popen
from math import * # any function could be used by set() from math import * # any function could be used by set()
try: try:
import numpy as np import numpy as np
oldnumeric = False oldnumeric = False
except: except:
import Numeric as np import Numeric as np
oldnumeric = True oldnumeric = True
try: from DEFAULTS import PIZZA_GUNZIP try: from DEFAULTS import PIZZA_GUNZIP
except: PIZZA_GUNZIP = "gunzip" except: PIZZA_GUNZIP = "gunzip"
@ -210,9 +216,9 @@ class dump:
self.bondlist = [] self.bondlist = []
self.triflag = 0 self.triflag = 0
self.trilist = [] self.trilist = []
self.triobj = 0
self.lineflag = 0 self.lineflag = 0
self.linelist = [] self.linelist = []
self.objextra = None
# flist = list of all dump file names # flist = list of all dump file names
@ -263,19 +269,9 @@ class dump:
self.tselect.all() self.tselect.all()
# set default names for atom columns if file wasn't self-describing # print column assignments
if len(self.snaps) == 0: if len(self.names):
print "no column assignments made"
elif len(self.names):
print "assigned columns:",self.names2str()
elif self.snaps[0].atoms == None:
print "no column assignments made"
elif len(self.snaps[0].atoms[0]) == 5:
self.map(1,"id",2,"type",3,"x",4,"y",5,"z")
print "assigned columns:",self.names2str()
elif len(self.snaps[0].atoms[0]) == 8:
self.map(1,"id",2,"type",3,"x",4,"y",5,"z",6,"ix",7,"iy",8,"iz")
print "assigned columns:",self.names2str() print "assigned columns:",self.names2str()
else: else:
print "no column assignments made" print "no column assignments made"
@ -285,10 +281,11 @@ class dump:
if (not self.names.has_key("x")) or \ if (not self.names.has_key("x")) or \
(not self.names.has_key("y")) or \ (not self.names.has_key("y")) or \
(not self.names.has_key("z")): (not self.names.has_key("z")):
print "no unscaling could be performed" print "dump scaling status is unknown"
elif self.nsnaps > 0: elif self.nsnaps > 0:
if self.scaled(self.nsnaps-1): self.unscale() if self.scale_original == 1: self.unscale()
else: print "dump is already unscaled" elif self.scale_original == 0: print "dump is already unscaled"
else: print "dump scaling status is unknown"
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# read next snapshot from list of files # read next snapshot from list of files
@ -333,8 +330,10 @@ class dump:
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# read a single snapshot from file f # read a single snapshot from file f
# return snapshot or 0 if failed # return snapshot or 0 if failed
# assign column names if not already done and file is self-describing # for first snapshot only:
# convert xs,xu to x # assign column names (file must be self-describing)
# set scale_original to 0/1/-1 for unscaled/scaled/unknown
# convert xs,xu to x in names
def read_snapshot(self,f): def read_snapshot(self,f):
try: try:
@ -347,26 +346,62 @@ class dump:
snap.aselect = np.zeros(snap.natoms) snap.aselect = np.zeros(snap.natoms)
item = f.readline() item = f.readline()
words = item.split("BOUNDS ")
if len(words) == 1: snap.boxstr = ""
else: snap.boxstr = words[1].strip()
if "xy" in snap.boxstr: snap.triclinic = 1
else: snap.triclinic = 0
words = f.readline().split() words = f.readline().split()
snap.xlo,snap.xhi = float(words[0]),float(words[1]) if len(words) == 2:
words = f.readline().split() snap.xlo,snap.xhi,snap.xy = float(words[0]),float(words[1]),0.0
snap.ylo,snap.yhi = float(words[0]),float(words[1]) else:
words = f.readline().split() snap.xlo,snap.xhi,snap.xy = \
snap.zlo,snap.zhi = float(words[0]),float(words[1]) float(words[0]),float(words[1]),float(words[2])
words = f.readline().split()
if len(words) == 2:
snap.ylo,snap.yhi,snap.xz = float(words[0]),float(words[1]),0.0
else:
snap.ylo,snap.yhi,snap.xz = \
float(words[0]),float(words[1]),float(words[2])
words = f.readline().split()
if len(words) == 2:
snap.zlo,snap.zhi,snap.yz = float(words[0]),float(words[1]),0.0
else:
snap.zlo,snap.zhi,snap.yz = \
float(words[0]),float(words[1]),float(words[2])
item = f.readline() item = f.readline()
if len(self.names) == 0: if len(self.names) == 0:
self.scale_original = -1
xflag = yflag = zflag = -1
words = item.split()[2:] words = item.split()[2:]
if len(words): if len(words):
for i in range(len(words)): for i in range(len(words)):
if words[i] == "xs" or words[i] == "xu": if words[i] == "x" or words[i] == "xu":
xflag = 0
self.names["x"] = i self.names["x"] = i
elif words[i] == "ys" or words[i] == "yu": elif words[i] == "xs" or words[i] == "xsu":
xflag = 1
self.names["x"] = i
elif words[i] == "y" or words[i] == "yu":
yflag = 0
self.names["y"] = i self.names["y"] = i
elif words[i] == "zs" or words[i] == "zu": elif words[i] == "ys" or words[i] == "ysu":
yflag = 1
self.names["y"] = i
elif words[i] == "z" or words[i] == "zu":
zflag = 0
self.names["z"] = i
elif words[i] == "zs" or words[i] == "zsu":
zflag = 1
self.names["z"] = i self.names["z"] = i
else: self.names[words[i]] = i else: self.names[words[i]] = i
if xflag == 0 and yflag == 0 and zflag == 0: self.scale_original = 0
if xflag == 1 and yflag == 1 and zflag == 1: self.scale_original = 1
if snap.natoms: if snap.natoms:
words = f.readline().split() words = f.readline().split()
ncol = len(words) ncol = len(words)
@ -387,27 +422,6 @@ class dump:
except: except:
return 0 return 0
# --------------------------------------------------------------------
# decide if snapshot i is scaled/unscaled from coords of first and last atom
def scaled(self,i):
ix = self.names["x"]
iy = self.names["y"]
iz = self.names["z"]
natoms = self.snaps[i].natoms
if natoms == 0: return 0
x1 = self.snaps[i].atoms[0][ix]
y1 = self.snaps[i].atoms[0][iy]
z1 = self.snaps[i].atoms[0][iz]
x2 = self.snaps[i].atoms[natoms-1][ix]
y2 = self.snaps[i].atoms[natoms-1][iy]
z2 = self.snaps[i].atoms[natoms-1][iz]
if x1 >= -0.1 and x1 <= 1.1 and y1 >= -0.1 and y1 <= 1.1 and \
z1 >= -0.1 and z1 <= 1.1 and x2 >= -0.1 and x2 <= 1.1 and \
y2 >= -0.1 and y2 <= 1.1 and z2 >= -0.1 and z2 <= 1.1:
return 1
else: return 0
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# map atom column names # map atom column names
@ -418,9 +432,8 @@ class dump:
j = i + 1 j = i + 1
self.names[pairs[j]] = pairs[i]-1 self.names[pairs[j]] = pairs[i]-1
# delete unselected snapshots
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# delete unselected snapshots
def delete(self): def delete(self):
ndel = i = 0 ndel = i = 0
@ -435,6 +448,7 @@ class dump:
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# scale coords to 0-1 for all snapshots or just one # scale coords to 0-1 for all snapshots or just one
# use 6 params as h-matrix to treat orthongonal or triclinic boxes
def scale(self,*list): def scale(self,*list):
if len(list) == 0: if len(list) == 0:
@ -453,17 +467,52 @@ class dump:
# -------------------------------------------------------------------- # --------------------------------------------------------------------
def scale_one(self,snap,x,y,z): def scale_one(self,snap,x,y,z):
if snap.atoms == None: return if snap.xy == 0.0 and snap.xz == 0.0 and snap.yz == 0.0:
xprdinv = 1.0 / (snap.xhi - snap.xlo) xprdinv = 1.0 / (snap.xhi - snap.xlo)
yprdinv = 1.0 / (snap.yhi - snap.ylo) yprdinv = 1.0 / (snap.yhi - snap.ylo)
zprdinv = 1.0 / (snap.zhi - snap.zlo) zprdinv = 1.0 / (snap.zhi - snap.zlo)
atoms = snap.atoms atoms = snap.atoms
atoms[:,x] = (atoms[:,x] - snap.xlo) * xprdinv if atoms != None:
atoms[:,y] = (atoms[:,y] - snap.ylo) * yprdinv atoms[:,x] = (atoms[:,x] - snap.xlo) * xprdinv
atoms[:,z] = (atoms[:,z] - snap.zlo) * zprdinv atoms[:,y] = (atoms[:,y] - snap.ylo) * yprdinv
atoms[:,z] = (atoms[:,z] - snap.zlo) * zprdinv
else:
xlo_bound = snap.xlo; xhi_bound = snap.xhi
ylo_bound = snap.ylo; yhi_bound = snap.yhi
zlo_bound = snap.zlo; zhi_bound = snap.zhi
xy = snap.xy
xz = snap.xz
yz = snap.yz
xlo = xlo_bound - min((0.0,xy,xz,xy+xz))
xhi = xhi_bound - max((0.0,xy,xz,xy+xz))
ylo = ylo_bound - min((0.0,yz))
yhi = yhi_bound - max((0.0,yz))
zlo = zlo_bound
zhi = zhi_bound
h0 = xhi - xlo
h1 = yhi - ylo
h2 = zhi - zlo
h3 = yz
h4 = xz
h5 = xy
h0inv = 1.0 / h0
h1inv = 1.0 / h1
h2inv = 1.0 / h2
h3inv = yz / (h1*h2)
h4inv = (h3*h5 - h1*h4) / (h0*h1*h2)
h5inv = xy / (h0*h1)
atoms = snap.atoms
if atoms != None:
atoms[:,x] = (atoms[:,x] - snap.xlo)*h0inv + \
(atoms[:,y] - snap.ylo)*h5inv + \
(atoms[:,z] - snap.zlo)*h4inv
atoms[:,y] = (atoms[:,y] - snap.ylo)*h1inv + \
(atoms[:,z] - snap.zlo)*h3inv
atoms[:,z] = (atoms[:,z] - snap.zlo)*h2inv
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# unscale coords from 0-1 to box size for all snapshots or just one # unscale coords from 0-1 to box size for all snapshots or just one
# use 6 params as h-matrix to treat orthongonal or triclinic boxes
def unscale(self,*list): def unscale(self,*list):
if len(list) == 0: if len(list) == 0:
@ -482,15 +531,40 @@ class dump:
# -------------------------------------------------------------------- # --------------------------------------------------------------------
def unscale_one(self,snap,x,y,z): def unscale_one(self,snap,x,y,z):
if snap.atoms == None: return if snap.xy == 0.0 and snap.xz == 0.0 and snap.yz == 0.0:
xprd = snap.xhi - snap.xlo xprd = snap.xhi - snap.xlo
yprd = snap.yhi - snap.ylo yprd = snap.yhi - snap.ylo
zprd = snap.zhi - snap.zlo zprd = snap.zhi - snap.zlo
atoms = snap.atoms atoms = snap.atoms
atoms[:,x] = snap.xlo + atoms[:,x]*xprd if atoms != None:
atoms[:,y] = snap.ylo + atoms[:,y]*yprd atoms[:,x] = snap.xlo + atoms[:,x]*xprd
atoms[:,z] = snap.zlo + atoms[:,z]*zprd atoms[:,y] = snap.ylo + atoms[:,y]*yprd
atoms[:,z] = snap.zlo + atoms[:,z]*zprd
else:
xlo_bound = snap.xlo; xhi_bound = snap.xhi
ylo_bound = snap.ylo; yhi_bound = snap.yhi
zlo_bound = snap.zlo; zhi_bound = snap.zhi
xy = snap.xy
xz = snap.xz
yz = snap.yz
xlo = xlo_bound - min((0.0,xy,xz,xy+xz))
xhi = xhi_bound - max((0.0,xy,xz,xy+xz))
ylo = ylo_bound - min((0.0,yz))
yhi = yhi_bound - max((0.0,yz))
zlo = zlo_bound
zhi = zhi_bound
h0 = xhi - xlo
h1 = yhi - ylo
h2 = zhi - zlo
h3 = yz
h4 = xz
h5 = xy
atoms = snap.atoms
if atoms != None:
atoms[:,x] = snap.xlo + atoms[:,x]*h0 + atoms[:,y]*h5 + atoms[:,z]*h4
atoms[:,y] = snap.ylo + atoms[:,y]*h1 + atoms[:,z]*h3
atoms[:,z] = snap.zlo + atoms[:,z]*h2
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# wrap coords from outside box to inside # wrap coords from outside box to inside
@ -537,7 +611,8 @@ class dump:
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# wrap coords to same image as atom ID stored in "other" column # wrap coords to same image as atom ID stored in "other" column
# if dynamic extra lines or triangles defined, owrap them as well
def owrap(self,other): def owrap(self,other):
print "Wrapping to other ..." print "Wrapping to other ..."
@ -549,7 +624,7 @@ class dump:
iy = self.names["iy"] iy = self.names["iy"]
iz = self.names["iz"] iz = self.names["iz"]
iother = self.names[other] iother = self.names[other]
for snap in self.snaps: for snap in self.snaps:
xprd = snap.xhi - snap.xlo xprd = snap.xhi - snap.xlo
yprd = snap.yhi - snap.ylo yprd = snap.yhi - snap.ylo
@ -563,14 +638,17 @@ class dump:
atoms[i][x] += (atoms[i][ix]-atoms[j][ix])*xprd atoms[i][x] += (atoms[i][ix]-atoms[j][ix])*xprd
atoms[i][y] += (atoms[i][iy]-atoms[j][iy])*yprd atoms[i][y] += (atoms[i][iy]-atoms[j][iy])*yprd
atoms[i][z] += (atoms[i][iz]-atoms[j][iz])*zprd atoms[i][z] += (atoms[i][iz]-atoms[j][iz])*zprd
# should bonds also be owrapped ?
if self.lineflag == 2 or self.triflag == 2:
self.objextra.owrap(snap.time,xprd,yprd,zprd,ids,atoms,iother,ix,iy,iz)
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# convert column names assignment to a string, in column order # convert column names assignment to a string, in column order
def names2str(self): def names2str(self):
ncol = len(self.snaps[0].atoms[0])
pairs = self.names.items() pairs = self.names.items()
values = self.names.values() values = self.names.values()
ncol = len(pairs)
str = "" str = ""
for i in xrange(ncol): for i in xrange(ncol):
if i in values: str += pairs[values.index(i)][0] + ' ' if i in values: str += pairs[values.index(i)][0] + ' '
@ -614,6 +692,12 @@ class dump:
if len(self.snaps): namestr = self.names2str() if len(self.snaps): namestr = self.names2str()
if not append: f = open(file,"w") if not append: f = open(file,"w")
else: f = open(file,"a") else: f = open(file,"a")
if "id" in self.names: id = self.names["id"]
else: id = -1
if "type" in self.names: type = self.names["type"]
else: type = -1
for snap in self.snaps: for snap in self.snaps:
if not snap.tselect: continue if not snap.tselect: continue
print snap.time, print snap.time,
@ -624,10 +708,16 @@ class dump:
print >>f,snap.time print >>f,snap.time
print >>f,"ITEM: NUMBER OF ATOMS" print >>f,"ITEM: NUMBER OF ATOMS"
print >>f,snap.nselect print >>f,snap.nselect
print >>f,"ITEM: BOX BOUNDS" if snap.boxstr: print >>f,"ITEM: BOX BOUNDS",snap.boxstr
print >>f,snap.xlo,snap.xhi else: print >>f,"ITEM: BOX BOUNDS"
print >>f,snap.ylo,snap.yhi if snap.triclinic:
print >>f,snap.zlo,snap.zhi print >>f,snap.xlo,snap.xhi,snap.xy
print >>f,snap.ylo,snap.yhi,snap.xz
print >>f,snap.zlo,snap.zhi,snap.yz
else:
print >>f,snap.xlo,snap.xhi
print >>f,snap.ylo,snap.yhi
print >>f,snap.zlo,snap.zhi
print >>f,"ITEM: ATOMS",namestr print >>f,"ITEM: ATOMS",namestr
atoms = snap.atoms atoms = snap.atoms
@ -636,7 +726,7 @@ class dump:
if not snap.aselect[i]: continue if not snap.aselect[i]: continue
line = "" line = ""
for j in xrange(nvalues): for j in xrange(nvalues):
if (j < 2): if j == id or j == type:
line += str(int(atoms[i][j])) + " " line += str(int(atoms[i][j])) + " "
else: else:
line += str(atoms[i][j]) + " " line += str(atoms[i][j]) + " "
@ -660,10 +750,16 @@ class dump:
print >>f,snap.time print >>f,snap.time
print >>f,"ITEM: NUMBER OF ATOMS" print >>f,"ITEM: NUMBER OF ATOMS"
print >>f,snap.nselect print >>f,snap.nselect
print >>f,"ITEM: BOX BOUNDS" if snap.boxstr: print >>f,"ITEM: BOX BOUNDS",snap.boxstr
print >>f,snap.xlo,snap.xhi else: print >>f,"ITEM: BOX BOUNDS"
print >>f,snap.ylo,snap.yhi if snap.triclinic:
print >>f,snap.zlo,snap.zhi print >>f,snap.xlo,snap.xhi,snap.xy
print >>f,snap.ylo,snap.yhi,snap.xz
print >>f,snap.zlo,snap.zhi,snap.yz
else:
print >>f,snap.xlo,snap.xhi
print >>f,snap.ylo,snap.yhi
print >>f,snap.zlo,snap.zhi
print >>f,"ITEM: ATOMS",namestr print >>f,"ITEM: ATOMS",namestr
atoms = snap.atoms atoms = snap.atoms
@ -897,9 +993,20 @@ class dump:
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# return list of atoms to viz for snapshot isnap # return list of atoms to viz for snapshot isnap
# if called with flag, then index is timestep, so convert to snapshot index
# augment with bonds, tris, lines if extra() was invoked # augment with bonds, tris, lines if extra() was invoked
def viz(self,isnap): def viz(self,index,flag=0):
if not flag: isnap = index
else:
times = self.time()
n = len(times)
i = 0
while i < n:
if times[i] > index: break
i += 1
isnap = i - 1
snap = self.snaps[isnap] snap = self.snaps[isnap]
time = snap.time time = snap.time
@ -919,18 +1026,22 @@ class dump:
atom = snap.atoms[i] atom = snap.atoms[i]
atoms.append([atom[id],atom[type],atom[x],atom[y],atom[z]]) atoms.append([atom[id],atom[type],atom[x],atom[y],atom[z]])
# create list of current bond coords from static bondlist # create list of bonds from static or dynamic bond list
# then generate bond coords from bondlist
# alist = dictionary of atom IDs for atoms list # alist = dictionary of atom IDs for atoms list
# lookup bond atom IDs in alist and grab their coords # lookup bond atom IDs in alist and grab their coords
# try is used since some atoms may be unselected # try is used since some atoms may be unselected
# any bond with unselected atom is not returned to viz caller # any bond with unselected atom is not added to bonds
# need Numeric/Numpy mode here # need Numeric/Numpy mode here
bonds = [] bonds = []
if self.bondflag: if self.bondflag:
if self.bondflag == 1: bondlist = self.bondlist
elif self.bondflag == 2:
tmp1,tmp2,tmp3,bondlist,tmp4,tmp5 = self.objextra.viz(time,1)
alist = {} alist = {}
for i in xrange(len(atoms)): alist[int(atoms[i][0])] = i for i in xrange(len(atoms)): alist[int(atoms[i][0])] = i
for bond in self.bondlist: for bond in bondlist:
try: try:
i = alist[bond[2]] i = alist[bond[2]]
j = alist[bond[3]] j = alist[bond[3]]
@ -940,15 +1051,23 @@ class dump:
atom2[2],atom2[3],atom2[4],atom1[1],atom2[1]]) atom2[2],atom2[3],atom2[4],atom1[1],atom2[1]])
except: continue except: continue
# create list of tris from static or dynamic tri list
# if dynamic, could eliminate tris for unselected atoms
tris = [] tris = []
if self.triflag: if self.triflag:
if self.triflag == 1: tris = self.trilist if self.triflag == 1: tris = self.trilist
elif self.triflag == 2: elif self.triflag == 2:
timetmp,boxtmp,atomstmp,bondstmp, \ tmp1,tmp2,tmp3,tmp4,tris,tmp5 = self.objextra.viz(time,1)
tris,linestmp = self.triobj.viz(time,1)
# create list of lines from static or dynamic tri list
# if dynamic, could eliminate lines for unselected atoms
lines = [] lines = []
if self.lineflag: lines = self.linelist if self.lineflag:
if self.lineflag == 1: lines = self.linelist
elif self.lineflag == 2:
tmp1,tmp2,tmp3,tmp4,tmp5,lines = self.objextra.viz(time,1)
return time,box,atoms,bonds,tris,lines return time,box,atoms,bonds,tris,lines
@ -991,49 +1110,14 @@ class dump:
# -------------------------------------------------------------------- # --------------------------------------------------------------------
# grab bonds/tris/lines from another object # grab bonds/tris/lines from another object
# if static, grab once, else store obj to grab dynamically
def extra(self,arg): def extra(self,arg):
# read bonds from bond dump file # data object, grab bonds statically
if type(arg) is types.StringType: if type(arg) is types.InstanceType and ".data" in str(arg.__class__):
try: self.bondflag = 0
f = open(arg,'r')
item = f.readline()
time = int(f.readline())
item = f.readline()
nbonds = int(f.readline())
item = f.readline()
if not re.search("BONDS",item):
raise StandardError, "could not read bonds from dump file"
words = f.readline().split()
ncol = len(words)
for i in xrange(1,nbonds):
words += f.readline().split()
f.close()
# convert values to int and absolute value since can be negative types
if oldnumeric: bondlist = np.zeros((nbonds,4),np.Int)
else: bondlist = np.zeros((nbonds,4),np.int)
ints = [abs(int(value)) for value in words]
start = 0
stop = 4
for i in xrange(nbonds):
bondlist[i] = ints[start:stop]
start += ncol
stop += ncol
if bondlist:
self.bondflag = 1
self.bondlist = bondlist
except:
raise StandardError,"could not read from bond dump file"
# request bonds from data object
elif type(arg) is types.InstanceType and ".data" in str(arg.__class__):
try: try:
bondlist = [] bondlist = []
bondlines = arg.sections["Bonds"] bondlines = arg.sections["Bonds"]
@ -1047,9 +1131,10 @@ class dump:
except: except:
raise StandardError,"could not extract bonds from data object" raise StandardError,"could not extract bonds from data object"
# request tris/lines from cdata object # cdata object, grab tris and lines statically
elif type(arg) is types.InstanceType and ".cdata" in str(arg.__class__): elif type(arg) is types.InstanceType and ".cdata" in str(arg.__class__):
self.triflag = self.lineflag = 0
try: try:
tmp,tmp,tmp,tmp,tris,lines = arg.viz(0) tmp,tmp,tmp,tmp,tris,lines = arg.viz(0)
if tris: if tris:
@ -1061,14 +1146,29 @@ class dump:
except: except:
raise StandardError,"could not extract tris/lines from cdata object" raise StandardError,"could not extract tris/lines from cdata object"
# request tris from mdump object # mdump object, grab tris dynamically
elif type(arg) is types.InstanceType and ".mdump" in str(arg.__class__): elif type(arg) is types.InstanceType and ".mdump" in str(arg.__class__):
try: self.triflag = 2
self.triflag = 2 self.objextra = arg
self.triobj = arg
except: # bdump object, grab bonds dynamically
raise StandardError,"could not extract tris from mdump object"
elif type(arg) is types.InstanceType and ".bdump" in str(arg.__class__):
self.bondflag = 2
self.objextra = arg
# ldump object, grab lines dynamically
elif type(arg) is types.InstanceType and ".ldump" in str(arg.__class__):
self.lineflag = 2
self.objextra = arg
# tdump object, grab tris dynamically
elif type(arg) is types.InstanceType and ".tdump" in str(arg.__class__):
self.triflag = 2
self.objextra = arg
else: else:
raise StandardError,"unrecognized argument to dump.extra()" raise StandardError,"unrecognized argument to dump.extra()"

View File

@ -19,7 +19,7 @@ g.plot(a,b) plot B against A
g.plot(a,b,c,d,...) plot B against A, D against C, etc g.plot(a,b,c,d,...) plot B against A, D against C, etc
g.mplot(M,N,S,"file",a,b,...) multiple plots saved to file0000.eps, etc g.mplot(M,N,S,"file",a,b,...) multiple plots saved to file0000.eps, etc
each plot argument can be a tuple, list, or Numeric vector each plot argument can be a tuple, list, or Numeric/NumPy vector
mplot loops over range(M,N,S) and create one plot per iteration mplot loops over range(M,N,S) and create one plot per iteration
last args are same as list of vectors for plot(), e.g. 1, 2, 4 vectors last args are same as list of vectors for plot(), e.g. 1, 2, 4 vectors
each plot is made from a portion of the vectors, depending on loop index i each plot is made from a portion of the vectors, depending on loop index i

View File

@ -219,6 +219,8 @@ class log:
if self.style == 1: if self.style == 1:
s1 = txt.find(self.firststr) s1 = txt.find(self.firststr)
s2 = txt.find("\n--",s1) s2 = txt.find("\n--",s1)
if (s2 == -1):
s2 = txt.find("\nLoop time of",s1)
pattern = "\s(\S*)\s*=" pattern = "\s(\S*)\s*="
keywords = re.findall(pattern,txt[s1:s2]) keywords = re.findall(pattern,txt[s1:s2])
keywords.insert(0,"Step") keywords.insert(0,"Step")
@ -270,9 +272,7 @@ class log:
s1 = txt.find(self.firststr,start) s1 = txt.find(self.firststr,start)
s2 = txt.find("Loop time of",start+1) s2 = txt.find("Loop time of",start+1)
if s2 == -1:
s2 = txt.find("ERROR",start+1)
if s1 >= 0 and s2 >= 0 and s1 < s2: # found s1,s2 with s1 before s2 if s1 >= 0 and s2 >= 0 and s1 < s2: # found s1,s2 with s1 before s2
if self.style == 2: if self.style == 2:
s1 = txt.find("\n",s1) + 1 s1 = txt.find("\n",s1) + 1
@ -296,9 +296,6 @@ class log:
if txt.find("Loop time of",start) == start: # end of file, so exit if txt.find("Loop time of",start) == start: # end of file, so exit
eof -= len(txt) - start # reset eof to "Loop" eof -= len(txt) - start # reset eof to "Loop"
break break
if txt.find("ERROR",start) == start: # end of file, so exit
eof -= len(txt) - start # reset eof to "ERROR"
break
last = 1 # entire read is a chunk last = 1 # entire read is a chunk
s1 = 0 s1 = 0
@ -324,6 +321,7 @@ class log:
word2 = re.findall(pat2,section) word2 = re.findall(pat2,section)
words = word1 + word2 words = word1 + word2
self.data.append(map(float,words)) self.data.append(map(float,words))
else: else:
lines = chunk.split("\n") lines = chunk.split("\n")
for line in lines: for line in lines: