Merge branch 'lammps:develop' into mliappy_unified

This commit is contained in:
Steven Anaya
2022-09-09 14:03:43 -06:00
committed by GitHub
407 changed files with 297686 additions and 8972 deletions

File diff suppressed because it is too large Load Diff

View File

@ -0,0 +1 @@
../../../tools/python/pizza/dump.py

View File

@ -1,390 +0,0 @@
# Pizza.py toolkit, https://lammps.github.io/pizza
# Steve Plimpton, sjplimp@sandia.gov, Sandia National Laboratories
#
# Copyright (2005) Sandia Corporation. Under the terms of Contract
# DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
# certain rights in this software. This software is distributed under
# the GNU General Public License.
# for python3 compatibility
from __future__ import print_function
# gnu tool
oneline = "Create plots via GnuPlot plotting program"
docstr = """
g = gnu() start up GnuPlot
g.stop() shut down GnuPlot process
g.plot(a) plot vector A against linear index
g.plot(a,b) plot B against A
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
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
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
Ith plot is of b[0:i] vs a[0:i], etc
series of plots saved as file0000.eps, file0001.eps, etc
if use xrange(),yrange() then plot axes will be same for all plots
g("plot 'file.dat' using 2:3 with lines") execute string in GnuPlot
g.enter() enter GnuPlot shell
gnuplot> plot sin(x) with lines type commands directly to GnuPlot
gnuplot> exit, quit exit GnuPlot shell
g.export("data",range(100),a,...) create file with columns of numbers
all vectors must be of equal length
could plot from file with GnuPlot command: plot 'data' using 1:2 with lines
g.select(N) figure N becomes the current plot
subsequent commands apply to this plot
g.hide(N) delete window for figure N
g.save("file") save current plot as file.eps
Set attributes for current plot:
g.erase() reset all attributes to default values
g.aspect(1.3) aspect ratio
g.xtitle("Time") x axis text
g.ytitle("Energy") y axis text
g.title("My Plot") title text
g.title("title","x","y") title, x axis, y axis text
g.xrange(xmin,xmax) x axis range
g.xrange() default x axis range
g.yrange(ymin,ymax) y axis range
g.yrange() default y axis range
g.xlog() toggle x axis between linear and log
g.ylog() toggle y axis between linear and log
g.label(x,y,"text") place label at x,y coords
g.curve(N,'r') set color of curve N
colors: 'k' = black, 'r' = red, 'g' = green, 'b' = blue
'm' = magenta, 'c' = cyan, 'y' = yellow
"""
# History
# 8/05, Matt Jones (BYU): original version
# 9/05, Steve Plimpton: added mplot() method
# ToDo list
# allow choice of JPG or PNG or GIF when saving ?
# can this be done from GnuPlot or have to do via ImageMagick convert ?
# way to trim EPS plot that is created ?
# hide does not work on Mac aqua
# select does not pop window to front on Mac aqua
# Variables
# current = index of current figure (1-N)
# figures = list of figure objects with each plot's attributes
# so they aren't lost between replots
# Imports and external programs
import os
import sys
try: from DEFAULTS import PIZZA_GNUPLOT
except ImportError: PIZZA_GNUPLOT = "gnuplot -p"
try: from DEFAULTS import PIZZA_GNUTERM
except ImportError: PIZZA_GNUTERM = "x11"
# Class definition
class gnu:
# --------------------------------------------------------------------
def __init__(self):
self.GNUPLOT = os.popen(PIZZA_GNUPLOT,'w')
self.file = "tmp.gnu"
self.figures = []
self.select(1)
# --------------------------------------------------------------------
def stop(self):
self.__call__("quit")
del self.GNUPLOT
# --------------------------------------------------------------------
def __call__(self,command):
self.GNUPLOT.write(command + '\n')
self.GNUPLOT.flush()
# --------------------------------------------------------------------
def enter(self):
while 1:
if sys.version_info[0] == 3:
command = input("gnuplot> ")
else:
command = raw_input("gnuplot> ")
if command == "quit" or command == "exit": return
self.__call__(command)
# --------------------------------------------------------------------
# write plot vectors to files and plot them
def plot(self,*vectors):
if len(vectors) == 1:
file = self.file + ".%d.1" % self.current
linear = range(len(vectors[0]))
self.export(file,linear,vectors[0])
self.figures[self.current-1].ncurves = 1
else:
if len(vectors) % 2: raise Exception("vectors must come in pairs")
for i in range(0,len(vectors),2):
file = self.file + ".%d.%d" % (self.current,i/2+1)
self.export(file,vectors[i],vectors[i+1])
self.figures[self.current-1].ncurves = len(vectors)/2
self.draw()
# --------------------------------------------------------------------
# create multiple plots from growing vectors, save to numbered files
# don't plot empty vector, create a [0] instead
def mplot(self,start,stop,skip,file,*vectors):
n = 0
for i in range(start,stop,skip):
partial_vecs = []
for vec in vectors:
if i: partial_vecs.append(vec[:i])
else: partial_vecs.append([0])
self.plot(*partial_vecs)
if n < 10: newfile = file + "000" + str(n)
elif n < 100: newfile = file + "00" + str(n)
elif n < 1000: newfile = file + "0" + str(n)
else: newfile = file + str(n)
self.save(newfile)
n += 1
# --------------------------------------------------------------------
# write list of equal-length vectors to filename
def export(self,filename,*vectors):
n = len(vectors[0])
for vector in vectors:
if len(vector) != n: raise Exception("vectors must be same length")
f = open(filename,'w')
nvec = len(vectors)
for i in range(n):
for j in range(nvec):
print(str(vectors[j][i])+" ",file=f,end='')
print ("",file=f)
f.close()
# --------------------------------------------------------------------
# select plot N as current plot
def select(self,n):
self.current = n
if len(self.figures) < n:
for i in range(n - len(self.figures)):
self.figures.append(figure())
cmd = "set term " + PIZZA_GNUTERM + ' ' + str(n)
self.__call__(cmd)
if self.figures[n-1].ncurves: self.draw()
# --------------------------------------------------------------------
# delete window for plot N
def hide(self,n):
cmd = "set term %s close %d" % (PIZZA_GNUTERM,n)
self.__call__(cmd)
# --------------------------------------------------------------------
# save plot to file.eps
# final re-select will reset terminal
# do not continue until plot file is written out
# else script could go forward and change data file
# use tmp.done as semaphore to indicate plot is finished
def save(self,file):
self.__call__("set terminal postscript enhanced solid lw 2 color portrait")
cmd = "set output '%s.eps'" % file
self.__call__(cmd)
if os.path.exists("tmp.done"): os.remove("tmp.done")
self.draw()
self.__call__("!touch tmp.done")
while not os.path.exists("tmp.done"): continue
self.__call__("set output")
self.select(self.current)
# --------------------------------------------------------------------
# restore default attributes by creating a new fig object
def erase(self):
fig = figure()
fig.ncurves = self.figures[self.current-1].ncurves
self.figures[self.current-1] = fig
self.draw()
# --------------------------------------------------------------------
def aspect(self,value):
self.figures[self.current-1].aspect = value
self.draw()
# --------------------------------------------------------------------
def xrange(self,*values):
if len(values) == 0:
self.figures[self.current-1].xlimit = 0
else:
self.figures[self.current-1].xlimit = (values[0],values[1])
self.draw()
# --------------------------------------------------------------------
def yrange(self,*values):
if len(values) == 0:
self.figures[self.current-1].ylimit = 0
else:
self.figures[self.current-1].ylimit = (values[0],values[1])
self.draw()
# --------------------------------------------------------------------
def label(self,x,y,text):
self.figures[self.current-1].labels.append((x,y,text))
self.figures[self.current-1].nlabels += 1
self.draw()
# --------------------------------------------------------------------
def nolabels(self):
self.figures[self.current-1].nlabel = 0
self.figures[self.current-1].labels = []
self.draw()
# --------------------------------------------------------------------
def title(self,*strings):
if len(strings) == 1:
self.figures[self.current-1].title = strings[0]
else:
self.figures[self.current-1].title = strings[0]
self.figures[self.current-1].xtitle = strings[1]
self.figures[self.current-1].ytitle = strings[2]
self.draw()
# --------------------------------------------------------------------
def xtitle(self,label):
self.figures[self.current-1].xtitle = label
self.draw()
# --------------------------------------------------------------------
def ytitle(self,label):
self.figures[self.current-1].ytitle = label
self.draw()
# --------------------------------------------------------------------
def xlog(self):
if self.figures[self.current-1].xlog:
self.figures[self.current-1].xlog = 0
else:
self.figures[self.current-1].xlog = 1
self.draw()
# --------------------------------------------------------------------
def ylog(self):
if self.figures[self.current-1].ylog:
self.figures[self.current-1].ylog = 0
else:
self.figures[self.current-1].ylog = 1
self.draw()
# --------------------------------------------------------------------
def curve(self,num,color):
fig = self.figures[self.current-1]
while len(fig.colors) < num: fig.colors.append(0)
fig.colors[num-1] = colormap[color]
self.draw()
# --------------------------------------------------------------------
# draw a plot with all its settings
# just return if no files of vectors defined yet
def draw(self):
fig = self.figures[self.current-1]
if not fig.ncurves: return
cmd = 'set size ratio ' + str(1.0/float(fig.aspect))
self.__call__(cmd)
cmd = 'set title ' + '"' + fig.title + '"'
self.__call__(cmd)
cmd = 'set xlabel ' + '"' + fig.xtitle + '"'
self.__call__(cmd)
cmd = 'set ylabel ' + '"' + fig.ytitle + '"'
self.__call__(cmd)
if fig.xlog: self.__call__("set logscale x")
else: self.__call__("unset logscale x")
if fig.ylog: self.__call__("set logscale y")
else: self.__call__("unset logscale y")
if fig.xlimit:
cmd = 'set xr [' + str(fig.xlimit[0]) + ':' + str(fig.xlimit[1]) + ']'
self.__call__(cmd)
else: self.__call__("set xr [*:*]")
if fig.ylimit:
cmd = 'set yr [' + str(fig.ylimit[0]) + ':' + str(fig.ylimit[1]) + ']'
self.__call__(cmd)
else: self.__call__("set yr [*:*]")
self.__call__("set nolabel")
for i in range(fig.nlabels):
x = fig.labels[i][0]
y = fig.labels[i][1]
text = fig.labels[i][2]
cmd = 'set label ' + '\"' + text + '\" at ' + str(x) + ',' + str(y)
self.__call__(cmd)
self.__call__("set key off")
cmd = 'plot '
for i in range(int(fig.ncurves)):
file = self.file + ".%d.%d" % (self.current,i+1)
if len(fig.colors) > i and fig.colors[i]:
cmd += "'" + file + "' using 1:2 with line %d, " % fig.colors[i]
else:
cmd += "'" + file + "' using 1:2 with lines, "
self.__call__(cmd[:-2])
# --------------------------------------------------------------------
# class to store settings for a single plot
class figure:
def __init__(self):
self.ncurves = 0
self.colors = []
self.title = ""
self.xtitle = ""
self.ytitle = ""
self.aspect = 1.3
self.xlimit = 0
self.ylimit = 0
self.xlog = 0
self.ylog = 0
self.nlabels = 0
self.labels = []
# --------------------------------------------------------------------
# line color settings
colormap = {'k':-1, 'r':1, 'g':2, 'b':3, 'm':4, 'c':5, 'y':7}

View File

@ -0,0 +1 @@
../../../tools/python/pizza/gnu.py

View File

@ -1,299 +0,0 @@
# Pizza.py toolkit, https://lammps.github.io/pizza
# Steve Plimpton, sjplimp@sandia.gov, Sandia National Laboratories
#
# Copyright (2005) Sandia Corporation. Under the terms of Contract
# DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
# certain rights in this software. This software is distributed under
# the GNU General Public License.
# for python3 compatibility
from __future__ import print_function
# pdb tool
oneline = "Read, write PDB files in combo with LAMMPS snapshots"
docstr = """
p = pdbfile("3CRO") create pdb object from PDB file or WWW
p = pdbfile("pep1 pep2") read in multiple PDB files
p = pdbfile("pep*") can use wildcards
p = pdbfile(d) read in snapshot data with no PDB file
p = pdbfile("3CRO",d) read in single PDB file with snapshot data
string arg contains one or more PDB files
don't need .pdb suffix except wildcard must expand to file.pdb
if only one 4-char file specified and it is not found,
it will be downloaded from http://www.rcsb.org as 3CRO.pdb
d arg is object with atom coordinates (dump, data)
p.one() write all output as one big PDB file to tmp.pdb
p.one("mine") write to mine.pdb
p.many() write one PDB file per snapshot: tmp0000.pdb, ...
p.many("mine") write as mine0000.pdb, mine0001.pdb, ...
p.single(N) write timestamp N as tmp.pdb
p.single(N,"new") write as new.pdb
how new PDB files are created depends on constructor inputs:
if no d: one new PDB file for each file in string arg (just a copy)
if only d specified: one new PDB file per snapshot in generic format
if one file in str arg and d: one new PDB file per snapshot
using input PDB file as template
multiple input PDB files with a d is not allowed
index,time,flag = p.iterator(0)
index,time,flag = p.iterator(1)
iterator = loop over number of PDB files
call first time with arg = 0, thereafter with arg = 1
N = length = # of snapshots or # of input PDB files
index = index of snapshot or input PDB file (0 to N-1)
time = timestep value (time stamp for snapshot, index for multiple PDB)
flag = -1 when iteration is done, 1 otherwise
typically call p.single(time) in iterated loop to write out one PDB file
"""
# History
# 8/05, Steve Plimpton (SNL): original version
# 3/17, Richard Berger (Temple U): improve Python 3 compatibility
# ToDo list
# for generic PDB file (no template) from a LJ unit system,
# the atoms in PDB file are too close together
# Variables
# files = list of input PDB files
# data = data object (ccell,data,dump) to read snapshots from
# atomlines = dict of ATOM lines in original PDB file
# key = atom id, value = tuple of (beginning,end) of line
# Imports and external programs
import sys, glob, urllib
PY3 = sys.version_info[0] == 3
if PY3:
string_types = str,
else:
string_types = basestring
# Class definition
class pdbfile:
# --------------------------------------------------------------------
def __init__(self,*args):
if len(args) == 1:
if type(args[0]) is string_types:
filestr = args[0]
self.data = None
else:
filestr = None
self.data = args[0]
elif len(args) == 2:
filestr = args[0]
self.data = args[1]
else: raise Exception("invalid args for pdb()")
# flist = full list of all PDB input file names
# append .pdb if needed
if filestr:
list = filestr.split()
flist = []
for file in list:
if '*' in file: flist += glob.glob(file)
else: flist.append(file)
for i in range(len(flist)):
if flist[i][-4:] != ".pdb": flist[i] += ".pdb"
if len(flist) == 0:
raise Exception("no PDB file specified")
self.files = flist
else: self.files = []
if len(self.files) > 1 and self.data:
raise Exception("cannot use multiple PDB files with data object")
if len(self.files) == 0 and not self.data:
raise Exception("no input PDB file(s)")
# grab PDB file from http://rcsb.org if not a local file
if len(self.files) == 1 and len(self.files[0]) == 8:
try:
open(self.files[0],'r').close()
except:
print("downloading %s from http://rcsb.org" % self.files[0])
fetchstr = "http://www.rcsb.org/pdb/cgi/export.cgi/%s?format=PDB&pdbId=2cpk&compression=None" % self.files[0]
urllib.urlretrieve(fetchstr,self.files[0])
if self.data and len(self.files): self.read_template(self.files[0])
# --------------------------------------------------------------------
# write a single large PDB file for concatenating all input data or files
# if data exists:
# only selected atoms returned by extract
# atoms written in order they appear in snapshot
# atom only written if its tag is in PDB template file
# if no data:
# concatenate all input files to one output file
def one(self,*args):
if len(args) == 0: file = "tmp.pdb"
elif args[0][-4:] == ".pdb": file = args[0]
else: file = args[0] + ".pdb"
f = open(file,'w')
# use template PDB file with each snapshot
if self.data:
n = flag = 0
while 1:
which,time,flag = self.data.iterator(flag)
if flag == -1: break
self.convert(f,which)
print("END",file=f)
print(time,end='')
sys.stdout.flush()
n += 1
else:
for file in self.files:
f.write(open(file,'r').read())
print("END",file=f)
print(file,end='')
sys.stdout.flush()
f.close()
print("\nwrote %d datasets to %s in PDB format" % (n,file))
# --------------------------------------------------------------------
# write series of numbered PDB files
# if data exists:
# only selected atoms returned by extract
# atoms written in order they appear in snapshot
# atom only written if its tag is in PDB template file
# if no data:
# just copy all input files to output files
def many(self,*args):
if len(args) == 0: root = "tmp"
else: root = args[0]
if self.data:
n = flag = 0
while 1:
which,time,flag = self.data.iterator(flag)
if flag == -1: break
if n < 10:
file = root + "000" + str(n)
elif n < 100:
file = root + "00" + str(n)
elif n < 1000:
file = root + "0" + str(n)
else:
file = root + str(n)
file += ".pdb"
f = open(file,'w')
self.convert(f,which)
f.close()
print(time,end='')
sys.stdout.flush()
n += 1
else:
n = 0
for infile in self.files:
if n < 10:
file = root + "000" + str(n)
elif n < 100:
file = root + "00" + str(n)
elif n < 1000:
file = root + "0" + str(n)
else:
file = root + str(n)
file += ".pdb"
f = open(file,'w')
f.write(open(infile,'r').read())
f.close()
print(file,end='')
sys.stdout.flush()
n += 1
print("\nwrote %d datasets to %s*.pdb in PDB format" % (n,root))
# --------------------------------------------------------------------
# write a single PDB file
# if data exists:
# time is timestamp in snapshot
# only selected atoms returned by extract
# atoms written in order they appear in snapshot
# atom only written if its tag is in PDB template file
# if no data:
# time is index into list of input PDB files
# just copy one input file to output file
def single(self,time,*args):
if len(args) == 0: file = "tmp.pdb"
elif args[0][-4:] == ".pdb": file = args[0]
else: file = args[0] + ".pdb"
f = open(file,'w')
if self.data:
which = self.data.findtime(time)
self.convert(f,which)
else:
f.write(open(self.files[time],'r').read())
f.close()
# --------------------------------------------------------------------
# iterate over list of input files or selected snapshots
# latter is done via data objects iterator
def iterator(self,flag):
if not self.data:
if not flag: self.iterate = 0
else:
self.iterate += 1
if self.iterate > len(self.files): return 0,0,-1
return self.iterate,self.iterate,1
return self.data.iterator(flag)
# --------------------------------------------------------------------
# read a PDB file and store ATOM lines
def read_template(self,file):
lines = open(file,'r').readlines()
self.atomlines = {}
for line in lines:
if line.find("ATOM") == 0:
tag = int(line[4:11])
begin = line[:30]
end = line[54:]
self.atomlines[tag] = (begin,end)
# --------------------------------------------------------------------
# convert one set of atoms to PDB format and write to f
def convert(self,f,which):
time,box,atoms,bonds,tris,lines = self.data.viz(which)
if len(self.files):
for atom in atoms:
id = atom[0]
if self.atomlines.has_key(id):
(begin,end) = self.atomlines[id]
line = "%s%8.3f%8.3f%8.3f%s" % (begin,atom[2],atom[3],atom[4],end)
print(line,file=f,end='')
else:
for atom in atoms:
begin = "ATOM %6d %2d R00 1 " % (atom[0],atom[1])
middle = "%8.3f%8.3f%8.3f" % (atom[2],atom[3],atom[4])
end = " 1.00 0.00 NONE"
print(begin+middle+end,file=f)

View File

@ -0,0 +1 @@
../../../tools/python/pizza/pdbfile.py

View File

@ -300,7 +300,7 @@ class ElemwiseModels(torch.nn.Module):
self.subnets = subnets
self.n_types = n_types
def forward(self, descriptors, elems):
def forward(self, descriptors, elems, dtype=torch.float64):
"""
Feeds descriptors to network model after adding zeros into
descriptor columns relating to different atom types
@ -319,8 +319,12 @@ class ElemwiseModels(torch.nn.Module):
Per atom attribute computed by the network model
"""
per_atom_attributes = torch.zeros(elems.size[0])
self.dtype=dtype
self.to(self.dtype)
per_atom_attributes = torch.zeros(elems.size(dim=0), dtype=self.dtype)
given_elems, elem_indices = torch.unique(elems, return_inverse=True)
for i, elem in enumerate(given_elems):
per_atom_attribute[elem_indices == i] = self.subnets[elem](descriptors[elem_indices == i])
self.subnets[elem].to(self.dtype)
per_atom_attributes[elem_indices == i] = self.subnets[elem](descriptors[elem_indices == i]).flatten()
return per_atom_attributes