From 49b72ac36e14dbb8baab462df03c8621b999adbf Mon Sep 17 00:00:00 2001
From: sjplimp <sjplimp@f3b2605a-c512-4ea7-a41b-209d697bcdaa>
Date: Tue, 22 Feb 2011 22:32:28 +0000
Subject: [PATCH] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@5712
f3b2605a-c512-4ea7-a41b-209d697bcdaa
---
tools/eff/README | 5 +++-
tools/eff/VMD-input.py | 33 ++++++++++++++--------
tools/eff/cfg2lammps.py | 57 ++++++++++++++++++++++++++++++++------
tools/eff/lmp2radii.py | 10 +++----
tools/eff/lmp2radii.pyx | 8 +++---
tools/eff/lmp2radii_col.py | 6 ++--
tools/eff/lmp2xyz.py | 57 ++++++++++++++++++++++++++++++--------
tools/eff/radii.vmd | 35 +++++++++++++++++++----
8 files changed, 162 insertions(+), 49 deletions(-)
diff --git a/tools/eff/README b/tools/eff/README
index edc52edaa0..9c68472c5d 100644
--- a/tools/eff/README
+++ b/tools/eff/README
@@ -5,9 +5,12 @@ generating LAMMPS structure files:
- h2.pl: rectangular lattice of hydrogen atoms
- Diamond.pl: diamond A4 box
- Li-hydride: Lithium hydride solid box
-- Li-solid: Lithium solid box
+- Li-solid: Lithium solid box (Li-solid-angstromg produces data file in Angstroms)
- Uniform-electron-gas.pl: uniform electron gas on an NaCl lattice
+- bohr2ang.py: converts an eFF structure file described in bohr to
+ its corresponding version in Angstroms
+
And other useful scripts for processing pEFF related information are
included, such as:
diff --git a/tools/eff/VMD-input.py b/tools/eff/VMD-input.py
index 2d6c5f7233..bb1fa905aa 100644
--- a/tools/eff/VMD-input.py
+++ b/tools/eff/VMD-input.py
@@ -12,8 +12,8 @@ ajaramil@wag.caltech.edu
Project: pEFF
Version: August 2009
-Usage: python VMD-input.py lammps_dump_filename radii_column_number
-Example: python VMD-input.py dump.lammpstrj 6
+Usage: python VMD-input.py lammps_dump_filename radii_column xpos_column
+Example: python VMD-input.py dump.lammpstrj 5 6
1. Extracts the electron radii from a lammps trajectory dump into %s.out
2. Creates %s.xyz file
@@ -28,20 +28,31 @@ def printHelp(input):
Info%(input,input,input)
if __name__ == '__main__':
-
+
# if no input, print help and exit
if len(sys.argv) < 2:
- print "Usage: python VMD-input.py lammps_dump_filename radii_column_number\n"
+ print "Usage: python VMD-input.py lammps_dump_filename radii_column xpos_column\n"
sys.exit(1)
else:
infile=sys.argv[1]
+ workdir=os.getcwd()
+ tools=sys.argv[0].split('VMD-input.py')[0]
# set defaults
outfile = infile.split('.')[0]
- if len(sys.argv) == 2:
+ print sys.argv
+ if len(sys.argv) == 4:
column = int(sys.argv[2])
- else:
- column=6 # default = radius for dump -> id type x y z spin radius
+ xpos = int(sys.argv[3])
+ print "Assuming xpos=eradius+1"
+ elif len(sys.argv) == 3:
+ column = int(sys.argv[2])
+ xpos = column+1
+ print "Assuming xpos=eradius+1"
+ elif len(sys.argv) == 2:
+ column=5 # default = radius for dump -> id type q spin eradius x y z
+ xpos=6
+ else: print "Incorrect number of arguments"
# check for input:
opts, argv = getopt(sys.argv[1:], 'c:o:ha')
@@ -54,12 +65,12 @@ if __name__ == '__main__':
outfile=arg
if opt == '-c': # select column from lammpstrj file to tabulate
column=int(arg)
-
+ print column,xpos
makeradii(infile,outfile+".out",column,True)
- lmp2xyz(infile,outfile+".xyz")
+ lmp2xyz(infile,outfile+".xyz",xpos)
print "Creating %s file ..."%(outfile+".vmd")
- os.system("cat %s | sed 's/xyzfile/%s/' > %s"%("radii.vmd",outfile+".xyz","temp"))
- os.system("cat %s | sed 's/radiifile/%s/' > %s; rm temp"%("temp",outfile+".out",outfile+".vmd"))
+ os.system("cat %s | sed 's/xyzfile/%s/' > %s"%(tools+"radii.vmd",outfile+".xyz","temp"))
+ os.system("cat %s | sed 's/radiifile/%s/' > %s; rm temp"%("temp",outfile+".out",workdir+'/'+outfile+".vmd"))
print "Done !! (you can now source %s using VMD's console) \n"%(outfile+".vmd")
print "NOTE: In VMD, set graphics representation for electrons to transparency,"
diff --git a/tools/eff/cfg2lammps.py b/tools/eff/cfg2lammps.py
index 4e200c9723..1bf364aa93 100644
--- a/tools/eff/cfg2lammps.py
+++ b/tools/eff/cfg2lammps.py
@@ -12,6 +12,8 @@ Version: August 2009
Reads in an eff .cfg file and produces the corresponding lammps data and input files
NOTE: Unsupported functions will be reported in the output log
+
+12/2010: Added support for fixed-core and pseudo-core structures
"""
# import essentials:
@@ -38,15 +40,31 @@ units electron
newton on
boundary %s
-atom_style hybrid charge electron
+atom_style electron
read_data data.${sname}
pair_style eff/cut %s
pair_coeff * *
+compute energies all pair eff/cut
+variable eke equal c_energies[1]
+variable epauli equal c_energies[2]
+variable estatics equal c_energies[3]
+variable errestrain equal c_energies[4]
+
+communicate single vel yes
+
+compute peratom all stress/atom
+compute p all reduce sum c_peratom[1] c_peratom[2] c_peratom[3]
+variable press equal -(c_p[1]+c_p[2]+c_p[3])/(3*vol)
+
+compute effTemp all temp/eff
+compute effPress all pressure effTemp
+
thermo %s
-thermo_style multi
+thermo_style custom step etotal pe ke v_eke v_epauli v_estatics v_errestrain temp press v_press
+thermo_modify temp effTemp press effPress
"""
#%(date,name,boundary,cutoff,period)
@@ -55,7 +73,7 @@ minimize="""
min_style cg
dump 1 %s xyz %s ${sname}.min.xyz
-dump 2 %s custom %s ${sname}.min.lammpstrj id type x y z radius fx fy fz rf
+dump 2 %s custom %s ${sname}.min.lammpstrj id type q spin eradius x y z fx fy fz erforce
min_modify line quadratic
minimize 0 1.0e-5 %s %s
@@ -77,7 +95,7 @@ timestep %s
fix %s
-dump 1 %s custom %s ${sname}.%s.lammpstrj id type x y z spin radius
+dump 1 %s custom %s ${sname}.%s.lammpstrj id type q spin eradius x y z
dump 2 %s custom %s ${sname}.%s.xyz
run %s
@@ -114,8 +132,10 @@ def generate_lammps_input(infile):
lines = fin.xreadlines()
print 7*"\b"+"[DONE]"
+ numcores=0
numnuclei=0
numelec=0
+ cores={}
nuclei={}
electrons={}
masses=[]
@@ -132,6 +152,9 @@ def generate_lammps_input(infile):
if line.find("@params")==0:
flag='params'
continue
+ elif line.find("@cores")==0:
+ flag='cores'
+ continue
elif line.find("@nuclei")==0:
flag='nuclei'
continue
@@ -198,7 +221,7 @@ def generate_lammps_input(infile):
if line.find("set_limit_stiffness")>=0:
continue
if line.find("output_position")>=0:
- dump_pos="dump\t1 all custom %s ${sname}.lammpstrj id type x y z spin radius "%(period)
+ dump_pos="dump\t1 all custom %s ${sname}.lammpstrj id type q spin eradius x y z "%(period)
if line.find("output_velocities")>=0:
dump_pos+="vx vy vz "
if line.find("output_energy_forces")>=0:
@@ -232,6 +255,22 @@ def generate_lammps_input(infile):
cutoff=line.split()[2]
continue
+ if flag=='cores' and len(line)>1:
+ numcores+=1
+ ln=line.split()
+ nc=' '.join(ln[0:3])
+ q=ln[3]
+ spin='3'
+ radius=ln[4]
+ m=q2m[int(float(q))]
+ if m not in masses:
+ masses.append(m)
+ massstr+="%d %s\n"%(types,m)
+ q2type[q]=types
+ types+=1
+ cores[numcores]=[nc,q,spin,radius]
+ continue
+
if flag=='nuclei' and len(line)>1:
numnuclei+=1
ln=line.split()
@@ -282,15 +321,17 @@ def generate_lammps_input(infile):
print "Writing datafile to %s ... "%('data.'+infile),
sys.stdout.flush()
print "\b"*19+"General section ",
- datafile.writelines("Created using cfg2lammps (c) AJB-2009\n\n%d atoms\n%d atom types\n\n%s%s%s\n"%(numnuclei+numelec,types,xbound,ybound,zbound))
+ datafile.writelines("Created using cfg2lammps (c) AJB-2009\n\n%d atoms\n%d atom types\n\n%s%s%s\n"%(numcores+numnuclei+numelec,types,xbound,ybound,zbound))
print "\b"*19+"Masses section ",
datafile.writelines(massstr)
print "\b"*19+"Atoms section ",
datafile.writelines("Atoms\n\n")
+ for n in range(numcores):
+ datafile.writelines("%d %d %2.2f %s %s %s\n"%(n+1,q2type[cores[n+1][1]],float(cores[n+1][1]),cores[n+1][2],cores[n+1][3],cores[n+1][0]))
for n in range(numnuclei):
- datafile.writelines("%d %d %s %s 0 0.0\n"%(n+1,q2type[nuclei[n+1][1]],nuclei[n+1][0],nuclei[n+1][1]))
+ datafile.writelines("%d %d %2.2f 0 0.0 %s\n"%(n+numcores+1,q2type[nuclei[n+1][1]],float(nuclei[n+1][1]),nuclei[n+1][0]))
for e in range(numelec):
- datafile.write("%d %d %s 0.0 %s %s\n"%(e+numnuclei+1,types,electrons[e+1][0],electrons[e+1][1],electrons[e+1][2]))
+ datafile.write("%d %d 0.0 %s %s %s\n"%(e+numnuclei+numcores+1,types,electrons[e+1][1],electrons[e+1][2],electrons[e+1][0]))
print "\b"*19+"Velocities section\n",
datafile.writelines(vels)
datafile.writelines("\n")
diff --git a/tools/eff/lmp2radii.py b/tools/eff/lmp2radii.py
index 28c23d6289..488d19932c 100644
--- a/tools/eff/lmp2radii.py
+++ b/tools/eff/lmp2radii.py
@@ -11,9 +11,9 @@ Version: August 2009
Extracts the electron radii from a lammps trajectory dump of style custom:
-dump 1 all custom period dump_file id type x y z spin radius ...
+dump 1 all custom period dump_file id type q spin eradius x y z...
-NOTE: The radius must be the 6th column per trajectory entry in the dump file
+NOTE: The radius must be the 5th column per trajectory entry in the dump file
"""
@@ -45,7 +45,7 @@ def makeradii(infile):
tmp=open("atoms",'r')
atoms=int(tmp.readlines()[1].split()[0])
tmp.close()
- os.system("rm -rf frames atoms")
+ os.system("rm -rf frames atoms lines")
arry=numpy.zeros((atoms,frames),dtype=float)
framecnt=0
header=9
@@ -59,8 +59,8 @@ def makeradii(infile):
elif (i >= lo) and (i < hi):
lparse=line.split()
id=int(lparse[0])
- r=float(lparse[6])
- if (r!=0):
+ r=float(lparse[4])
+ if (r!=0.0):
arry[id-1][framecnt]=r
if (framecnt==0): ecount+=1
if (i==lo+1):
diff --git a/tools/eff/lmp2radii.pyx b/tools/eff/lmp2radii.pyx
index f114f243ca..5e109fea15 100644
--- a/tools/eff/lmp2radii.pyx
+++ b/tools/eff/lmp2radii.pyx
@@ -11,9 +11,9 @@ Version: August 2009
Extracts the electron radii from a lammps trajectory dump of style custom:
-dump 1 all custom period dump_file id type x y z spin radius ...
+dump 1 all custom period dump_file id type q spin eradius x y z...
-NOTE: The radius must be the 6th column per trajectory entry in the dump file
+NOTE: The radius must be the 5th column per trajectory entry in the dump file
"""
@@ -59,7 +59,7 @@ def makeradii(infile):
elif (i >= lo) and (i < hi):
lparse=line.split()
id=int(lparse[0])
- r=float(lparse[6])
+ r=float(lparse[4])
if (r!=0):
arry[id-1][framecnt]=r
if (framecnt==0): ecount+=1
@@ -87,7 +87,7 @@ def makeradii(infile):
fout.writelines("%f\t"%(arry[a][f]))
fout.writelines("\n")
print
- print "DONE .... GOODBYE !!"
+ print "Done !! (generated radii/frame table) \n"
fout.close()
fin.close()
diff --git a/tools/eff/lmp2radii_col.py b/tools/eff/lmp2radii_col.py
index d5505116e9..653fccbb66 100644
--- a/tools/eff/lmp2radii_col.py
+++ b/tools/eff/lmp2radii_col.py
@@ -11,7 +11,7 @@ Version: August 2009
Extracts the electron radii from a lammps trajectory dump of style custom:
-dump 1 all custom period dump_file id type x y z spin radius ...
+dump 1 all custom period dump_file id type spin eradius x y z ...
NOTE: The radius must be the "column" per trajectory entry in the dump file
@@ -59,7 +59,7 @@ def makeradii(infile,outfile,column,True):
elif (i >= lo) and (i < hi):
lparse=line.split()
id=int(lparse[0])
- r=float(lparse[column])
+ r=float(lparse[column-1])
if (r!=0):
arry[id-1][framecnt]=r
if (framecnt==0): ecount+=1
@@ -96,7 +96,7 @@ if __name__ == '__main__':
# set defaults
outfile = ""
flag_all = False
- column=6 # default = radius
+ column=5 # default = radius
# check for input:
opts, argv = getopt(sys.argv[1:], 'c:o:ha')
diff --git a/tools/eff/lmp2xyz.py b/tools/eff/lmp2xyz.py
index 704188118c..8fa9660ec1 100644
--- a/tools/eff/lmp2xyz.py
+++ b/tools/eff/lmp2xyz.py
@@ -15,15 +15,35 @@ Usage: python lmp2xyz.py lammps_dump_filename xyz_filename
"""
import os, sys
-from math import log10
+from math import log10,floor
from numpy import zeros
-mass={"1.00794":"H","4.002602":"He","6.941":"Li","9.012182":"Be","10.811":"B","12.0107":"C","1.00":"Au","0.0005486":"Au"}
+masses={"1.00794":"H","4.002602":"He","6.941":"Li","9.012182":"Be","10.811":"B","12.0107":"C","1.00":"Au","0.0005486":"Au"}
+mass_floor={1:"H",4:"He",6:"Li",9:"Be",10:"B",12:"C",0:"Au",28:"Si"}
-def lmp2xyz(lammps,xyz):
+def lmp2xyz(lammps,xyz,xpos):
print "\nGenerating %s file"%(xyz)
fin=open(lammps,'r')
fout=open(xyz,'w')
+ data=raw_input("Do you have a corresponding data file? please enter filename or 'n': ")
+ count=1
+ if data!='n':
+ dataf=open(data,'r')
+ datafile=dataf.readlines()
+ dataf.close()
+ for line in datafile:
+ if line.find("atom types")>=0:
+ numtypes=int(line.split()[0])
+ mass=zeros(numtypes,dtype=float)
+ elif line.find("Masses")>=0:
+ count+=1+datafile.index(line)
+ elif line.find("Atoms")>=0:
+ break
+ for i in range(numtypes):
+ mass[i]=float(datafile[count].split()[1])
+ count+=1
+ else:
+ print "\nWill continue without a data file specification"
header=9
lines=fin.readlines()
numatoms=lines[3].split()[0]
@@ -31,13 +51,15 @@ def lmp2xyz(lammps,xyz):
tmp=open('lines','r')
tlines=tmp.readline()
tmp.close()
+ os.system("rm lines")
flines=int(tlines.split()[0])
snaps=flines/(int(numatoms)+header)
countsnap=1
- coords=zeros((int(numatoms),4),dtype=float)
+ if data!='n': coords={}
+ else: coords=zeros((int(numatoms),4),dtype=float)
sys.stdout.write("Writing [%d]: "%(snaps))
sys.stdout.flush()
- read_atoms=1
+ read_atoms=0
for line in lines:
if line.find('ITEM: TIMESTEP')==0:
read_atom_flag=False
@@ -53,14 +75,23 @@ def lmp2xyz(lammps,xyz):
read_atoms+=1
parse=line.split()
if parse[0]!="":
- coords[int(parse[0])-1][0]=int(parse[1])
- coords[int(parse[0])-1][1]=float(parse[2])
- coords[int(parse[0])-1][2]=float(parse[3])
- coords[int(parse[0])-1][3]=float(parse[4])
+# print [mass_floor[int(floor(mass[int(parse[1])-1]))],float(parse[xpos-1]),float(parse[xpos]),float(parse[xpos+1])]
+ if data!='n':
+ if mass[int(parse[1])-1]==1.0:
+ type='Au'
+ else:
+ type=mass_floor[int(floor(mass[int(parse[1])-1]))]
+ coords[int(parse[0])-1]=[type,float(parse[xpos-1]),float(parse[xpos]),float(parse[xpos+1])]
+ else:
+ coords[int(parse[0])-1][0]=int(parse[1])
+ coords[int(parse[0])-1][1]=float(parse[xpos-1])
+ coords[int(parse[0])-1][2]=float(parse[xpos])
+ coords[int(parse[0])-1][3]=float(parse[xpos+1])
if read_atoms==int(numatoms):
read_atoms=0
for i in range(int(numatoms)):
- fout.writelines("%d %2.4f %2.4f %2.4f\n"%(coords[i][0],coords[i][1],coords[i][2],coords[i][3]))
+ if data!='n': fout.writelines("%s %2.4f %2.4f %2.4f\n"%(coords[i][0],coords[i][1],coords[i][2],coords[i][3]))
+ else: fout.writelines("%d %2.4f %2.4f %2.4f\n"%(coords[i][0],coords[i][1],coords[i][2],coords[i][3]))
print "\nDone converting to xyz!!\n"
fin.close()
@@ -76,6 +107,8 @@ if __name__ == '__main__':
inputfile=sys.argv[1]
outfile=sys.argv[2]
-
- lmp2xyz(inputfile,outfile.split()[0])
+ if len(sys.argv)==4:
+ xpos=sys.arv[3]-1
+ else: xpos=5
+ lmp2xyz(inputfile,outfile.split()[0],xpos)
diff --git a/tools/eff/radii.vmd b/tools/eff/radii.vmd
index d6d211233b..0dfc6e3dc0 100644
--- a/tools/eff/radii.vmd
+++ b/tools/eff/radii.vmd
@@ -171,7 +171,7 @@ proc do_radii {args} {
foreach elec $data_array(0) r $data_array($fr) {
set s [atomselect $molid "index [expr {$elec -1}]"]
#set nr [expr {exp($r)}]
- set nr $r
+ set nr [expr $r*1.5]
$s set radius $nr
$s delete
#puts stderr "$fr $elec $nr"
@@ -192,9 +192,6 @@ global data_array molid
set xyz xyzfile
set data radiifile
-# switch default rep to VDW.
-mol default style VDW
-
# load nuclear and electron xyz trajectory
#set xyz [lindex $::argv 0]
@@ -202,7 +199,35 @@ mol default style VDW
#set datafile [lindex $::argv 1]
set molid [mol new $xyz waitfor all]
-mol modstyle 0 [molinfo top] VDW 1.0 32.0
+#mol modstyle 0 [molinfo top] VDW 0.2 32.0
+
+set s [atomselect $molid "all"]
+$s set radius 0.1
+$s delete
+
+color Display Background white
+display projection orthographic
+
+# switch default rep to VDW.
+mol default style VDW
+mol delrep 0 $molid
+mol selection "name H He Li Be B C Si"
+mol addrep $molid
+mol modstyle 0 $molid VDW 2.0 32.0
+# spin +1
+mol material Transparent
+mol selection "name Au"
+mol addrep $molid
+mol modstyle 1 $molid VDW 0.4 32.0
+# spin -1
+#mol material Transparent
+mol selection "name Tl"
+mol addrep $molid
+mol modstyle 2 $molid VDW 0.4 32.0
+mol material Opaque
+mol selection "name H He Li Be B C Si"
+mol addrep $molid
+mol modstyle 3 $molid DynamicBonds 1.6 0.1 32
puts "Starting ..."
readFile $data
--
GitLab