ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Trying to remove symlink

Browse Source
This commit is contained in:
rohskopf
2022-07-07 17:12:40 -06:00
parent 8f431c5904
commit 139ecad13c
10 changed files with 0 additions and 4384 deletions
Download Patch File Download Diff File Expand all files Collapse all files

49
python/lammps/__init__.py
View File

@ -1,49 +0,0 @@
"""
LAMMPS module global members:
.. data:: __version__
Numerical representation of the LAMMPS version this
module was taken from. Has the same format as the
result of :py:func:`lammps.version`.
"""
from .constants import * # lgtm [py/polluting-import]
from .core import * # lgtm [py/polluting-import]
from .data import * # lgtm [py/polluting-import]
from .pylammps import * # lgtm [py/polluting-import]
# convert installed module string version to numeric version
def get_version_number():
import time
from os.path import join
from sys import version_info
# must report 0 when inside LAMMPS source tree
if __file__.find(join('python', 'lammps', '__init__.py')) > 0:
return 0
vstring = None
if version_info.major == 3 and version_info.minor >= 8:
from importlib.metadata import version, PackageNotFoundError
try:
vstring = version('lammps')
except PackageNotFoundError:
# nothing to do, ignore
pass
else:
from pkg_resources import get_distribution, DistributionNotFound
try:
vstring = get_distribution('lammps').version
except DistributionNotFound:
# nothing to do, ignore
pass
if not vstring:
return 0
t = time.strptime(vstring, "%Y.%m.%d")
return t.tm_year*10000 + t.tm_mon*100 + t.tm_mday
__version__ = get_version_number()

48
python/lammps/constants.py
View File

@ -1,48 +0,0 @@
# ----------------------------------------------------------------------
# LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
# https://www.lammps.org/ Sandia National Laboratories
# Steve Plimpton, sjplimp@sandia.gov
#
# Copyright (2003) 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.
#
# See the README file in the top-level LAMMPS directory.
# -------------------------------------------------------------------------
# various symbolic constants to be used
# in certain calls to select data formats
LAMMPS_AUTODETECT = None
LAMMPS_INT = 0
LAMMPS_INT_2D = 1
LAMMPS_DOUBLE = 2
LAMMPS_DOUBLE_2D = 3
LAMMPS_INT64 = 4
LAMMPS_INT64_2D = 5
LAMMPS_STRING = 6
# these must be kept in sync with the enums in library.h
LMP_STYLE_GLOBAL = 0
LMP_STYLE_ATOM = 1
LMP_STYLE_LOCAL = 2
LMP_TYPE_SCALAR = 0
LMP_TYPE_VECTOR = 1
LMP_TYPE_ARRAY = 2
LMP_SIZE_VECTOR = 3
LMP_SIZE_ROWS = 4
LMP_SIZE_COLS = 5
LMP_VAR_EQUAL = 0
LMP_VAR_ATOM = 1
# -------------------------------------------------------------------------
def get_ctypes_int(size):
from ctypes import c_int, c_int32, c_int64
if size == 4:
return c_int32
elif size == 8:
return c_int64
return c_int

2097
python/lammps/core.py
View File

File diff suppressed because it is too large Load Diff

92
python/lammps/data.py
View File

@ -1,92 +0,0 @@
# ----------------------------------------------------------------------
# LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
# https://www.lammps.org/ Sandia National Laboratories
# Steve Plimpton, sjplimp@sandia.gov
#
# Copyright (2003) 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.
#
# See the README file in the top-level LAMMPS directory.
# -------------------------------------------------------------------------
################################################################################
# LAMMPS data structures
# Written by Richard Berger <richard.berger@temple.edu>
################################################################################
class NeighList(object):
"""This is a wrapper class that exposes the contents of a neighbor list.
It can be used like a regular Python list. Each element is a tuple of:
* the atom local index
* its number of neighbors
* and a pointer to an c_int array containing local atom indices of its
neighbors
Internally it uses the lower-level LAMMPS C-library interface.
:param lmp: reference to instance of :py:class:`lammps`
:type lmp: lammps
:param idx: neighbor list index
:type idx: int
"""
def __init__(self, lmp, idx):
self.lmp = lmp
self.idx = idx
def __str__(self):
return "Neighbor List ({} atoms)".format(self.size)
def __repr__(self):
return self.__str__()
@property
def size(self):
"""
:return: number of elements in neighbor list
"""
return self.lmp.get_neighlist_size(self.idx)
def get(self, element):
"""
Access a specific neighbor list entry. "element" must be a number from 0 to the size-1 of the list
:return: tuple with atom local index, number of neighbors and ctypes pointer to neighbor's local atom indices
:rtype: (int, int, ctypes.POINTER(c_int))
"""
iatom, numneigh, neighbors = self.lmp.get_neighlist_element_neighbors(self.idx, element)
return iatom, numneigh, neighbors
# the methods below implement the iterator interface, so NeighList can be used like a regular Python list
def __getitem__(self, element):
return self.get(element)
def __len__(self):
return self.size
def __iter__(self):
inum = self.size
for ii in range(inum):
yield self.get(ii)
def find(self, iatom):
"""
Find the neighbor list for a specific (local) atom iatom.
If there is no list for iatom, (-1, None) is returned.
:return: tuple with number of neighbors and ctypes pointer to neighbor's local atom indices
:rtype: (int, ctypes.POINTER(c_int))
"""
inum = self.size
for ii in range(inum):
idx, numneigh, neighbors = self.get(ii)
if idx == iatom:
return numneigh, neighbors
return -1, None

227
python/lammps/formats.py
View File

@ -1,227 +0,0 @@
# ----------------------------------------------------------------------
# LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
# https://www.lammps.org/ Sandia National Laboratories
# Steve Plimpton, sjplimp@sandia.gov
#
# Copyright (2003) 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.
#
# See the README file in the top-level LAMMPS directory.
# -------------------------------------------------------------------------
################################################################################
# LAMMPS output formats
# Written by Richard Berger <richard.berger@temple.edu>
# and Axel Kohlmeyer <akohlmey@gmail.com>
################################################################################
import re
has_yaml = False
try:
import yaml
has_yaml = True
try:
from yaml import CSafeLoader as Loader
except ImportError:
from yaml import SafeLoader as Loader
except ImportError:
# ignore here, raise an exception when trying to parse yaml instead
pass
class LogFile:
"""Reads LAMMPS log files and extracts the thermo information
It supports the line, multi, and yaml thermo output styles.
:param filename: path to log file
:type filename: str
:ivar runs: List of LAMMPS runs in log file. Each run is a dictionary with
thermo fields as keys, storing the values over time
:ivar errors: List of error lines in log file
"""
STYLE_DEFAULT = 0
STYLE_MULTI = 1
STYLE_YAML = 2
def __init__(self, filename):
alpha = re.compile(r'[a-df-zA-DF-Z]') # except e or E for floating-point numbers
kvpairs = re.compile(r'([a-zA-Z_0-9]+)\s+=\s*([0-9\.eE\-]+)')
style = LogFile.STYLE_DEFAULT
yamllog = ""
self.runs = []
self.errors = []
with open(filename, 'rt') as f:
in_thermo = False
in_data_section = False
for line in f:
if "ERROR" in line or "exited on signal" in line:
self.errors.append(line)
elif re.match(r'^ *Step ', line):
in_thermo = True
in_data_section = True
keys = line.split()
current_run = {}
for k in keys:
current_run[k] = []
elif re.match(r'^(keywords:.*$|data:$|---$| - \[.*\]$)', line):
if not has_yaml:
raise Exception('Cannot process YAML format logs without the PyYAML Python module')
style = LogFile.STYLE_YAML
yamllog += line;
current_run = {}
elif re.match(r'^\.\.\.$', line):
thermo = yaml.load(yamllog, Loader=Loader)
for k in thermo['keywords']:
current_run[k] = []
for step in thermo['data']:
icol = 0
for k in thermo['keywords']:
current_run[k].append(step[icol])
icol += 1
self.runs.append(current_run)
yamllog = ""
elif re.match(r'^------* Step ', line):
if not in_thermo:
current_run = {'Step': [], 'CPU': []}
in_thermo = True
in_data_section = True
style = LogFile.STYLE_MULTI
str_step, str_cpu = line.strip('-\n').split('-----')
step = float(str_step.split()[1])
cpu = float(str_cpu.split('=')[1].split()[0])
current_run["Step"].append(step)
current_run["CPU"].append(cpu)
elif line.startswith('Loop time of'):
in_thermo = False
if style != LogFile.STYLE_YAML:
self.runs.append(current_run)
elif in_thermo and in_data_section:
if style == LogFile.STYLE_DEFAULT:
if alpha.search(line):
continue
for k, v in zip(keys, map(float, line.split())):
current_run[k].append(v)
elif style == LogFile.STYLE_MULTI:
if '=' not in line:
in_data_section = False
continue
for k,v in kvpairs.findall(line):
if k not in current_run:
current_run[k] = [float(v)]
else:
current_run[k].append(float(v))
class AvgChunkFile:
"""Reads files generated by fix ave/chunk
:param filename: path to ave/chunk file
:type filename: str
:ivar timesteps: List of timesteps stored in file
:ivar total_count: total count over time
:ivar chunks: List of chunks. Each chunk is a dictionary containing its ID, the coordinates, and the averaged quantities
"""
def __init__(self, filename):
with open(filename, 'rt') as f:
timestep = None
chunks_read = 0
self.timesteps = []
self.total_count = []
self.chunks = []
for lineno, line in enumerate(f):
if lineno == 0:
if not line.startswith("# Chunk-averaged data for fix"):
raise Exception("Chunk data reader only supports default avg/chunk headers!")
parts = line.split()
self.fix_name = parts[5]
self.group_name = parts[8]
continue
elif lineno == 1:
if not line.startswith("# Timestep Number-of-chunks Total-count"):
raise Exception("Chunk data reader only supports default avg/chunk headers!")
continue
elif lineno == 2:
if not line.startswith("#"):
raise Exception("Chunk data reader only supports default avg/chunk headers!")
columns = line.split()[1:]
ndim = line.count("Coord")
compress = 'OrigID' in line
if ndim > 0:
coord_start = columns.index("Coord1")
coord_end = columns.index("Coord%d" % ndim)
ncount_start = coord_end + 1
data_start = ncount_start + 1
else:
coord_start = None
coord_end = None
ncount_start = 2
data_start = 3
continue
parts = line.split()
if timestep is None:
timestep = int(parts[0])
num_chunks = int(parts[1])
total_count = float(parts[2])
self.timesteps.append(timestep)
self.total_count.append(total_count)
for i in range(num_chunks):
self.chunks.append({
'coord' : [],
'ncount' : []
})
elif chunks_read < num_chunks:
chunk = int(parts[0])
ncount = float(parts[ncount_start])
if compress:
chunk_id = int(parts[1])
else:
chunk_id = chunk
current = self.chunks[chunk_id - 1]
current['id'] = chunk_id
current['ncount'].append(ncount)
if ndim > 0:
coord = tuple(map(float, parts[coord_start:coord_end+1]))
current['coord'].append(coord)
for i, data_column in list(enumerate(columns))[data_start:]:
value = float(parts[i])
if data_column in current:
current[data_column].append(value)
else:
current[data_column] = [value]
chunks_read += 1
assert chunk == chunks_read
else:
# do not support changing number of chunks
if not (num_chunks == int(parts[1])):
raise Exception("Currently, changing numbers of chunks are not supported.")
timestep = int(parts[0])
total_count = float(parts[2])
chunks_read = 0
self.timesteps.append(timestep)
self.total_count.append(total_count)

20
python/lammps/mliap/__init__.py
View File

@ -1,20 +0,0 @@
# Check compatiblity of this build with the python shared library.
# If this fails, lammps will segfault because its library will
# try to improperly start up a new interpreter.
import sysconfig
import ctypes
library = sysconfig.get_config_vars('INSTSONAME')[0]
try:
pylib = ctypes.CDLL(library)
except OSError as e:
if pylib.endswith(".a"):
pylib.strip(".a") + ".so"
pylib = ctypes.CDLL(library)
else:
raise e
if not pylib.Py_IsInitialized():
raise RuntimeError("This interpreter is not compatible with python-based mliap for LAMMPS.")
del sysconfig, ctypes, library, pylib
from .loader import load_model, activate_mliappy

52
python/lammps/mliap/loader.py
View File

@ -1,52 +0,0 @@
# ----------------------------------------------------------------------
# LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
# https://www.lammps.org/ Sandia National Laboratories
# Steve Plimpton, sjplimp@sandia.gov
#
# Copyright (2003) 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.
#
# See the README file in the top-level LAMMPS directory.
# -------------------------------------------------------------------------
# ----------------------------------------------------------------------
# Contributing author: Nicholas Lubbers (LANL)
# -------------------------------------------------------------------------
import sys
import importlib.util
import importlib.machinery
def activate_mliappy(lmp):
try:
# Begin Importlib magic to find the embedded python module
# This is needed because the filename for liblammps does not
# match the spec for normal python modules, wherein
# file names match with PyInit function names.
# Also, python normally doesn't look for extensions besides '.so'
# We fix both of these problems by providing an explict
# path to the extension module 'mliap_model_python_couple' in
path = lmp.lib._name
loader = importlib.machinery.ExtensionFileLoader('mliap_model_python_couple', path)
spec = importlib.util.spec_from_loader('mliap_model_python_couple', loader)
module = importlib.util.module_from_spec(spec)
sys.modules['mliap_model_python_couple'] = module
spec.loader.exec_module(module)
# End Importlib magic to find the embedded python module
except Exception as ee:
raise ImportError("Could not load ML-IAP python coupling module.") from ee
def load_model(model):
try:
import mliap_model_python_couple
except ImportError as ie:
raise ImportError("ML-IAP python module must be activated before loading\n"
"the pair style. Call lammps.mliap.activate_mliappy(lmp)."
) from ie
mliap_model_python_couple.load_from_python(model)

326
python/lammps/mliap/pytorch.py
View File

@ -1,326 +0,0 @@
# ----------------------------------------------------------------------
# LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
# https://www.lammps.org/ Sandia National Laboratories
# Steve Plimpton, sjplimp@sandia.gov
#
# Copyright (2003) 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.
#
# See the README file in the top-level LAMMPS directory.
# -------------------------------------------------------------------------
# ----------------------------------------------------------------------
# Contributing author: Nicholas Lubbers (LANL)
# -------------------------------------------------------------------------
import numpy as np
import torch
def calc_n_params(model):
"""
Returns the sum of two decimal numbers in binary digits.
Parameters:
model (torch.nn.Module): Network model that maps descriptors to a per atom attribute
Returns:
n_params (int): Number of NN model parameters
"""
return sum(p.nelement() for p in model.parameters())
class TorchWrapper(torch.nn.Module):
"""
A class to wrap Modules to ensure lammps mliap compatability.
...
Attributes
----------
model : torch.nn.Module
Network model that maps descriptors to a per atom attribute
device : torch.nn.Module (None)
Accelerator device
dtype : torch.dtype (torch.float64)
Dtype to use on device
n_params : torch.nn.Module (None)
Number of NN model parameters
n_descriptors : int
Max number of per atom descriptors
n_elements : int
Max number of elements
Methods
-------
forward(descriptors, elems):
Feeds descriptors to network model to produce per atom energies and forces.
"""
def __init__(self, model, n_descriptors, n_elements, n_params=None, device=None, dtype=torch.float64):
"""
Constructs all the necessary attributes for the network module.
Parameters
----------
model : torch.nn.Module
Network model that maps descriptors to a per atom attribute
n_descriptors : int
Max number of per atom descriptors
n_elements : int
Max number of elements
n_params : torch.nn.Module (None)
Number of NN model parameters
device : torch.nn.Module (None)
Accelerator device
dtype : torch.dtype (torch.float64)
Dtype to use on device
"""
super().__init__()
self.model = model
self.device = device
self.dtype = dtype
# Put model on device and convert to dtype
self.to(self.dtype)
self.to(self.device)
if n_params is None:
n_params = calc_n_params(model)
self.n_params = n_params
self.n_descriptors = n_descriptors
self.n_elements = n_elements
def forward(self, elems, descriptors, beta, energy):
"""
Takes element types and descriptors calculated via lammps and
calculates the per atom energies and forces.
Parameters
----------
elems : numpy.array
Per atom element types
descriptors : numpy.array
Per atom descriptors
beta : numpy.array
Expired beta array to be filled with new betas
energy : numpy.array
Expired per atom energy array to be filled with new per atom energy
(Note: This is a pointer to the lammps per atom energies)
Returns
-------
None
"""
descriptors = torch.from_numpy(descriptors).to(dtype=self.dtype, device=self.device).requires_grad_(True)
elems = torch.from_numpy(elems).to(dtype=torch.long, device=self.device) - 1
with torch.autograd.enable_grad():
energy_nn = self.model(descriptors, elems)
if energy_nn.ndim > 1:
energy_nn = energy_nn.flatten()
beta_nn = torch.autograd.grad(energy_nn.sum(), descriptors)[0]
beta[:] = beta_nn.detach().cpu().numpy().astype(np.float64)
energy[:] = energy_nn.detach().cpu().numpy().astype(np.float64)
class IgnoreElems(torch.nn.Module):
"""
A class to represent a NN model agnostic of element typing.
...
Attributes
----------
subnet : torch.nn.Module
Network model that maps descriptors to a per atom attribute
Methods
-------
forward(descriptors, elems):
Feeds descriptors to network model
"""
def __init__(self, subnet):
"""
Constructs all the necessary attributes for the network module.
Parameters
----------
subnet : torch.nn.Module
Network model that maps descriptors to a per atom attribute
"""
super().__init__()
self.subnet = subnet
def forward(self, descriptors, elems):
"""
Feeds descriptors to network model
Parameters
----------
descriptors : torch.tensor
Per atom descriptors
elems : torch.tensor
Per atom element types
Returns
-------
self.subnet(descriptors) : torch.tensor
Per atom attribute computed by the network model
"""
return self.subnet(descriptors)
class UnpackElems(torch.nn.Module):
"""
A class to represent a NN model pseudo-agnostic of element typing for
systems with multiple element typings.
...
Attributes
----------
subnet : torch.nn.Module
Network model that maps descriptors to a per atom attribute
n_types : int
Number of atom types used in training the NN model.
Methods
-------
forward(descriptors, elems):
Feeds descriptors to network model after adding zeros into
descriptor columns relating to different atom types
"""
def __init__(self, subnet, n_types):
"""
Constructs all the necessary attributes for the network module.
Parameters
----------
subnet : torch.nn.Module
Network model that maps descriptors to a per atom attribute.
n_types : int
Number of atom types used in training the NN model.
"""
super().__init__()
self.subnet = subnet
self.n_types = n_types
def forward(self, descriptors, elems):
"""
Feeds descriptors to network model after adding zeros into
descriptor columns relating to different atom types
Parameters
----------
descriptors : torch.tensor
Per atom descriptors
elems : torch.tensor
Per atom element types
Returns
-------
self.subnet(descriptors) : torch.tensor
Per atom attribute computed by the network model
"""
unpacked_descriptors = torch.zeros(elems.shape[0], self.n_types, descriptors.shape[1], dtype=torch.float64)
for i, ind in enumerate(elems):
unpacked_descriptors[i, ind, :] = descriptors[i]
return self.subnet(torch.reshape(unpacked_descriptors, (elems.shape[0], -1)), elems)
class ElemwiseModels(torch.nn.Module):
"""
A class to represent a NN model dependent on element typing.
...
Attributes
----------
subnets : list of torch.nn.Modules
Per element type network models that maps per element type
descriptors to a per atom attribute.
n_types : int
Number of atom types used in training the NN model.
Methods
-------
forward(descriptors, elems):
Feeds descriptors to network model after adding zeros into
descriptor columns relating to different atom types
"""
def __init__(self, subnets, n_types):
"""
Constructs all the necessary attributes for the network module.
Parameters
----------
subnets : list of torch.nn.Modules
Per element type network models that maps per element
type descriptors to a per atom attribute.
n_types : int
Number of atom types used in training the NN model.
"""
super().__init__()
self.subnets = subnets
self.n_types = n_types
def forward(self, descriptors, elems):
"""
Feeds descriptors to network model after adding zeros into
descriptor columns relating to different atom types
Parameters
----------
descriptors : torch.tensor
Per atom descriptors
elems : torch.tensor
Per atom element types
Returns
-------
self.subnets(descriptors) : torch.tensor
Per atom attribute computed by the network model
"""
per_atom_attributes = torch.zeros(elems.size[0])
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])
return per_atom_attributes

483
python/lammps/numpy_wrapper.py
View File

@ -1,483 +0,0 @@
# ----------------------------------------------------------------------
# LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
# https://www.lammps.org/ Sandia National Laboratories
# Steve Plimpton, sjplimp@sandia.gov
#
# Copyright (2003) 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.
#
# See the README file in the top-level LAMMPS directory.
# -------------------------------------------------------------------------
################################################################################
# NumPy additions
# Written by Richard Berger <richard.berger@temple.edu>
################################################################################
import warnings
from ctypes import POINTER, c_void_p, c_char_p, c_double, c_int, c_int32, c_int64, cast
from .constants import * # lgtm [py/polluting-import]
from .data import NeighList
class numpy_wrapper:
"""lammps API NumPy Wrapper
This is a wrapper class that provides additional methods on top of an
existing :py:class:`lammps` instance. The methods transform raw ctypes
pointers into NumPy arrays, which give direct access to the
original data while protecting against out-of-bounds accesses.
There is no need to explicitly instantiate this class. Each instance
of :py:class:`lammps` has a :py:attr:`numpy <lammps.numpy>` property
that returns an instance.
:param lmp: instance of the :py:class:`lammps` class
:type lmp: lammps
"""
def __init__(self, lmp):
self.lmp = lmp
# -------------------------------------------------------------------------
def _ctype_to_numpy_int(self, ctype_int):
import numpy as np
if ctype_int == c_int32:
return np.int32
elif ctype_int == c_int64:
return np.int64
return np.intc
# -------------------------------------------------------------------------
def extract_atom(self, name, dtype=LAMMPS_AUTODETECT, nelem=LAMMPS_AUTODETECT, dim=LAMMPS_AUTODETECT):
"""Retrieve per-atom properties from LAMMPS as NumPy arrays
This is a wrapper around the :py:meth:`lammps.extract_atom()` method.
It behaves the same as the original method, but returns NumPy arrays
instead of ``ctypes`` pointers.
.. note::
While the returned arrays of per-atom data are dimensioned
for the range [0:nmax] - as is the underlying storage -
the data is usually only valid for the range of [0:nlocal],
unless the property of interest is also updated for ghost
atoms. In some cases, this depends on a LAMMPS setting, see
for example :doc:`comm_modify vel yes <comm_modify>`.
:param name: name of the property
:type name: string
:param dtype: type of the returned data (see :ref:`py_datatype_constants`)
:type dtype: int, optional
:param nelem: number of elements in array
:type nelem: int, optional
:param dim: dimension of each element
:type dim: int, optional
:return: requested data as NumPy array with direct access to C data or None
:rtype: numpy.array or NoneType
"""
if dtype == LAMMPS_AUTODETECT:
dtype = self.lmp.extract_atom_datatype(name)
if nelem == LAMMPS_AUTODETECT:
if name == "mass":
nelem = self.lmp.extract_global("ntypes") + 1
else:
nelem = self.lmp.extract_global("nlocal")
if dim == LAMMPS_AUTODETECT:
if dtype in (LAMMPS_INT_2D, LAMMPS_DOUBLE_2D, LAMMPS_INT64_2D):
# TODO add other fields
if name in ("x", "v", "f", "x0","omega", "angmom", "torque", "csforce", "vforce", "vest"):
dim = 3
elif name == "smd_data_9":
dim = 9
elif name == "smd_stress":
dim = 6
else:
dim = 2
else:
dim = 1
raw_ptr = self.lmp.extract_atom(name, dtype)
if dtype in (LAMMPS_DOUBLE, LAMMPS_DOUBLE_2D):
return self.darray(raw_ptr, nelem, dim)
elif dtype in (LAMMPS_INT, LAMMPS_INT_2D):
return self.iarray(c_int32, raw_ptr, nelem, dim)
elif dtype in (LAMMPS_INT64, LAMMPS_INT64_2D):
return self.iarray(c_int64, raw_ptr, nelem, dim)
return raw_ptr
# -------------------------------------------------------------------------
def extract_atom_iarray(self, name, nelem, dim=1):
warnings.warn("deprecated, use extract_atom instead", DeprecationWarning)
if name in ['id', 'molecule']:
c_int_type = self.lmp.c_tagint
elif name in ['image']:
c_int_type = self.lmp.c_imageint
else:
c_int_type = c_int
if dim == 1:
raw_ptr = self.lmp.extract_atom(name, LAMMPS_INT)
else:
raw_ptr = self.lmp.extract_atom(name, LAMMPS_INT_2D)
return self.iarray(c_int_type, raw_ptr, nelem, dim)
# -------------------------------------------------------------------------
def extract_atom_darray(self, name, nelem, dim=1):
warnings.warn("deprecated, use extract_atom instead", DeprecationWarning)
if dim == 1:
raw_ptr = self.lmp.extract_atom(name, LAMMPS_DOUBLE)
else:
raw_ptr = self.lmp.extract_atom(name, LAMMPS_DOUBLE_2D)
return self.darray(raw_ptr, nelem, dim)
# -------------------------------------------------------------------------
def extract_compute(self, cid, cstyle, ctype):
"""Retrieve data from a LAMMPS compute
This is a wrapper around the
:py:meth:`lammps.extract_compute() <lammps.lammps.extract_compute()>` method.
It behaves the same as the original method, but returns NumPy arrays
instead of ``ctypes`` pointers.
:param cid: compute ID
:type cid: string
:param cstyle: style of the data retrieve (global, atom, or local), see :ref:`py_style_constants`
:type cstyle: int
:param ctype: type of the returned data (scalar, vector, or array), see :ref:`py_type_constants`
:type ctype: int
:return: requested data either as float, as NumPy array with direct access to C data, or None
:rtype: float, numpy.array, or NoneType
"""
value = self.lmp.extract_compute(cid, cstyle, ctype)
if cstyle == LMP_STYLE_GLOBAL:
if ctype == LMP_TYPE_VECTOR:
nrows = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_VECTOR)
return self.darray(value, nrows)
elif ctype == LMP_TYPE_ARRAY:
nrows = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_ROWS)
ncols = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_COLS)
return self.darray(value, nrows, ncols)
elif cstyle == LMP_STYLE_LOCAL:
nrows = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_ROWS)
ncols = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_COLS)
if ncols == 0:
return self.darray(value, nrows)
else:
return self.darray(value, nrows, ncols)
elif cstyle == LMP_STYLE_ATOM:
if ctype == LMP_TYPE_VECTOR:
nlocal = self.lmp.extract_global("nlocal")
return self.darray(value, nlocal)
elif ctype == LMP_TYPE_ARRAY:
nlocal = self.lmp.extract_global("nlocal")
ncols = self.lmp.extract_compute(cid, cstyle, LMP_SIZE_COLS)
return self.darray(value, nlocal, ncols)
return value
# -------------------------------------------------------------------------
def extract_fix(self, fid, fstyle, ftype, nrow=0, ncol=0):
"""Retrieve data from a LAMMPS fix
This is a wrapper around the :py:meth:`lammps.extract_fix() <lammps.lammps.extract_fix()>` method.
It behaves the same as the original method, but returns NumPy arrays
instead of ``ctypes`` pointers.
:param fid: fix ID
:type fid: string
:param fstyle: style of the data retrieve (global, atom, or local), see :ref:`py_style_constants`
:type fstyle: int
:param ftype: type or size of the returned data (scalar, vector, or array), see :ref:`py_type_constants`
:type ftype: int
:param nrow: index of global vector element or row index of global array element
:type nrow: int
:param ncol: column index of global array element
:type ncol: int
:return: requested data
:rtype: integer or double value, pointer to 1d or 2d double array or None
"""
value = self.lmp.extract_fix(fid, fstyle, ftype, nrow, ncol)
if fstyle == LMP_STYLE_ATOM:
if ftype == LMP_TYPE_VECTOR:
nlocal = self.lmp.extract_global("nlocal")
return self.darray(value, nlocal)
elif ftype == LMP_TYPE_ARRAY:
nlocal = self.lmp.extract_global("nlocal")
ncols = self.lmp.extract_fix(fid, fstyle, LMP_SIZE_COLS, 0, 0)
return self.darray(value, nlocal, ncols)
elif fstyle == LMP_STYLE_LOCAL:
if ftype == LMP_TYPE_VECTOR:
nrows = self.lmp.extract_fix(fid, fstyle, LMP_SIZE_ROWS, 0, 0)
return self.darray(value, nrows)
elif ftype == LMP_TYPE_ARRAY:
nrows = self.lmp.extract_fix(fid, fstyle, LMP_SIZE_ROWS, 0, 0)
ncols = self.lmp.extract_fix(fid, fstyle, LMP_SIZE_COLS, 0, 0)
return self.darray(value, nrows, ncols)
return value
# -------------------------------------------------------------------------
def extract_variable(self, name, group=None, vartype=LMP_VAR_EQUAL):
""" Evaluate a LAMMPS variable and return its data
This function is a wrapper around the function
:py:meth:`lammps.extract_variable() <lammps.lammps.extract_variable()>`
method. It behaves the same as the original method, but returns NumPy arrays
instead of ``ctypes`` pointers.
:param name: name of the variable to execute
:type name: string
:param group: name of group for atom-style variable (ignored for equal-style variables)
:type group: string
:param vartype: type of variable, see :ref:`py_vartype_constants`
:type vartype: int
:return: the requested data or None
:rtype: c_double, numpy.array, or NoneType
"""
import numpy as np
value = self.lmp.extract_variable(name, group, vartype)
if vartype == LMP_VAR_ATOM:
return np.ctypeslib.as_array(value)
return value
# -------------------------------------------------------------------------
def gather_bonds(self):
"""Retrieve global list of bonds as NumPy array
This is a wrapper around :py:meth:`lammps.gather_bonds() <lammps.lammps.gather_bonds()>`
It behaves the same as the original method, but returns a NumPy array instead
of a ``ctypes`` list.
.. versionadded:: 28Jul2021
:return: the requested data as a 2d-integer numpy array
:rtype: numpy.array(nbonds,3)
"""
import numpy as np
nbonds, value = self.lmp.gather_bonds()
return np.ctypeslib.as_array(value).reshape(nbonds,3)
# -------------------------------------------------------------------------
def fix_external_get_force(self, fix_id):
"""Get access to the array with per-atom forces of a fix external instance with a given fix ID.
This function is a wrapper around the
:py:meth:`lammps.fix_external_get_force() <lammps.lammps.fix_external_get_force()>`
method. It behaves the same as the original method, but returns a NumPy array instead
of a ``ctypes`` pointer.
.. versionchanged:: 28Jul2021
:param fix_id: Fix-ID of a fix external instance
:type: string
:return: requested data
:rtype: numpy.array
"""
import numpy as np
nlocal = self.lmp.extract_setting('nlocal')
value = self.lmp.fix_external_get_force(fix_id)
return self.darray(value,nlocal,3)
# -------------------------------------------------------------------------
def fix_external_set_energy_peratom(self, fix_id, eatom):
"""Set the per-atom energy contribution for a fix external instance with the given ID.
This function is an alternative to
:py:meth:`lammps.fix_external_set_energy_peratom() <lammps.lammps.fix_external_set_energy_peratom()>`
method. It behaves the same as the original method, but accepts a NumPy array
instead of a list as argument.
.. versionadded:: 28Jul2021
:param fix_id: Fix-ID of a fix external instance
:type: string
:param eatom: per-atom potential energy
:type: numpy.array
"""
import numpy as np
nlocal = self.lmp.extract_setting('nlocal')
if len(eatom) < nlocal:
raise Exception('per-atom energy dimension must be at least nlocal')
c_double_p = POINTER(c_double)
value = eatom.astype(np.double)
return self.lmp.lib.lammps_fix_external_set_energy_peratom(self.lmp.lmp, fix_id.encode(),
value.ctypes.data_as(c_double_p))
# -------------------------------------------------------------------------
def fix_external_set_virial_peratom(self, fix_id, vatom):
"""Set the per-atom virial contribution for a fix external instance with the given ID.
This function is an alternative to
:py:meth:`lammps.fix_external_set_virial_peratom() <lammps.lammps.fix_external_set_virial_peratom()>`
method. It behaves the same as the original method, but accepts a NumPy array
instead of a list as argument.
.. versionadded:: 28Jul2021
:param fix_id: Fix-ID of a fix external instance
:type: string
:param eatom: per-atom potential energy
:type: numpy.array
"""
import numpy as np
nlocal = self.lmp.extract_setting('nlocal')
if len(vatom) < nlocal:
raise Exception('per-atom virial first dimension must be at least nlocal')
if len(vatom[0]) != 6:
raise Exception('per-atom virial second dimension must be 6')
c_double_pp = np.ctypeslib.ndpointer(dtype=np.uintp, ndim=1, flags='C')
# recast numpy array to be compatible with library interface
value = (vatom.__array_interface__['data'][0]
+ np.arange(vatom.shape[0])*vatom.strides[0]).astype(np.uintp)
# change prototype to our custom type
self.lmp.lib.lammps_fix_external_set_virial_peratom.argtypes = [ c_void_p, c_char_p, c_double_pp ]
self.lmp.lib.lammps_fix_external_set_virial_peratom(self.lmp.lmp, fix_id.encode(), value)
# -------------------------------------------------------------------------
def get_neighlist(self, idx):
"""Returns an instance of :class:`NumPyNeighList` which wraps access to the neighbor list with the given index
:param idx: index of neighbor list
:type idx: int
:return: an instance of :class:`NumPyNeighList` wrapping access to neighbor list data
:rtype: NumPyNeighList
"""
if idx < 0:
return None
return NumPyNeighList(self.lmp, idx)
# -------------------------------------------------------------------------
def get_neighlist_element_neighbors(self, idx, element):
"""Return data of neighbor list entry
This function is a wrapper around the function
:py:meth:`lammps.get_neighlist_element_neighbors() <lammps.lammps.get_neighlist_element_neighbors()>`
method. It behaves the same as the original method, but returns a NumPy array containing the neighbors
instead of a ``ctypes`` pointer.
:param element: neighbor list index
:type element: int
:param element: neighbor list element index
:type element: int
:return: tuple with atom local index and numpy array of neighbor local atom indices
:rtype: (int, numpy.array)
"""
iatom, numneigh, c_neighbors = self.lmp.get_neighlist_element_neighbors(idx, element)
neighbors = self.iarray(c_int, c_neighbors, numneigh, 1)
return iatom, neighbors
# -------------------------------------------------------------------------
def iarray(self, c_int_type, raw_ptr, nelem, dim=1):
if raw_ptr is None:
return None
import numpy as np
np_int_type = self._ctype_to_numpy_int(c_int_type)
if dim == 1:
ptr = cast(raw_ptr, POINTER(c_int_type * nelem))
else:
ptr = cast(raw_ptr[0], POINTER(c_int_type * nelem * dim))
a = np.frombuffer(ptr.contents, dtype=np_int_type)
if dim > 1:
a.shape = (nelem, dim)
else:
a.shape = (nelem)
return a
# -------------------------------------------------------------------------
def darray(self, raw_ptr, nelem, dim=1):
if raw_ptr is None:
return None
import numpy as np
if dim == 1:
ptr = cast(raw_ptr, POINTER(c_double * nelem))
else:
ptr = cast(raw_ptr[0], POINTER(c_double * nelem * dim))
a = np.frombuffer(ptr.contents)
if dim > 1:
a.shape = (nelem, dim)
else:
a.shape = (nelem)
return a
# -------------------------------------------------------------------------
class NumPyNeighList(NeighList):
"""This is a wrapper class that exposes the contents of a neighbor list.
It can be used like a regular Python list. Each element is a tuple of:
* the atom local index
* a NumPy array containing the local atom indices of its neighbors
Internally it uses the lower-level LAMMPS C-library interface.
:param lmp: reference to instance of :py:class:`lammps`
:type lmp: lammps
:param idx: neighbor list index
:type idx: int
"""
def __init__(self, lmp, idx):
super(NumPyNeighList, self).__init__(lmp, idx)
def get(self, element):
"""
Access a specific neighbor list entry. "element" must be a number from 0 to the size-1 of the list
:return: tuple with atom local index, numpy array of neighbor local atom indices
:rtype: (int, numpy.array)
"""
iatom, neighbors = self.lmp.numpy.get_neighlist_element_neighbors(self.idx, element)
return iatom, neighbors
def find(self, iatom):
"""
Find the neighbor list for a specific (local) atom iatom.
If there is no list for iatom, None is returned.
:return: numpy array of neighbor local atom indices
:rtype: numpy.array or None
"""
inum = self.size
for ii in range(inum):
idx, neighbors = self.get(ii)
if idx == iatom:
return neighbors
return None

990
python/lammps/pylammps.py
View File

@ -1,990 +0,0 @@
# ----------------------------------------------------------------------
# LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
# https://www.lammps.org/ Sandia National Laboratories
# Steve Plimpton, sjplimp@sandia.gov
#
# Copyright (2003) 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.
#
# See the README file in the top-level LAMMPS directory.
# -------------------------------------------------------------------------
################################################################################
# Alternative Python Wrapper
# Written by Richard Berger <richard.berger@temple.edu>
################################################################################
# for python2/3 compatibility
from __future__ import print_function
import io
import os
import re
import sys
import tempfile
from collections import namedtuple
from .core import lammps
# -------------------------------------------------------------------------
class OutputCapture(object):
""" Utility class to capture LAMMPS library output """
def __init__(self):
self.stdout_fd = 1
self.captured_output = ""
def __enter__(self):
self.tmpfile = tempfile.TemporaryFile(mode='w+b')
sys.stdout.flush()
# make copy of original stdout
self.stdout_orig = os.dup(self.stdout_fd)
# replace stdout and redirect to temp file
os.dup2(self.tmpfile.fileno(), self.stdout_fd)
return self
def __exit__(self, exc_type, exc_value, traceback):
os.dup2(self.stdout_orig, self.stdout_fd)
os.close(self.stdout_orig)
self.tmpfile.close()
@property
def output(self):
sys.stdout.flush()
self.tmpfile.flush()
self.tmpfile.seek(0, io.SEEK_SET)
self.captured_output = self.tmpfile.read().decode('utf-8')
return self.captured_output
# -------------------------------------------------------------------------
class Variable(object):
def __init__(self, pylammps_instance, name, style, definition):
self._pylmp = pylammps_instance
self.name = name
self.style = style
self.definition = definition.split()
@property
def value(self):
if self.style == 'atom':
return list(self._pylmp.lmp.extract_variable(self.name, "all", 1))
else:
value = self._pylmp.lmp_print('"${%s}"' % self.name).strip()
try:
return float(value)
except ValueError:
return value
# -------------------------------------------------------------------------
class AtomList(object):
"""
A dynamic list of atoms that returns either an :py:class:`Atom` or
:py:class:`Atom2D` instance for each atom. Instances are only allocated
when accessed.
:ivar natoms: total number of atoms
:ivar dimensions: number of dimensions in system
"""
def __init__(self, pylammps_instance):
self._pylmp = pylammps_instance
self.natoms = self._pylmp.system.natoms
self.dimensions = self._pylmp.system.dimensions
self._loaded = {}
def __getitem__(self, index):
"""
Return Atom with given local index
:param index: Local index of atom
:type index: int
:rtype: Atom or Atom2D
"""
if index not in self._loaded:
if self.dimensions == 2:
atom = Atom2D(self._pylmp, index)
else:
atom = Atom(self._pylmp, index)
self._loaded[index] = atom
return self._loaded[index]
def __len__(self):
return self.natoms
# -------------------------------------------------------------------------
class Atom(object):
"""
A wrapper class then represents a single atom inside of LAMMPS
It provides access to properties of the atom and allows you to change some of them.
"""
def __init__(self, pylammps_instance, index):
self._pylmp = pylammps_instance
self.index = index
def __dir__(self):
return [k for k in super().__dir__() if not k.startswith('_')]
def get(self, name, index):
prop = self._pylmp.lmp.numpy.extract_atom(name)
if prop is not None:
return prop[index]
return None
@property
def id(self):
"""
Return the atom ID
:type: int
"""
return self.get("id", self.index)
@property
def type(self):
"""
Return the atom type
:type: int
"""
return self.get("type", self.index)
@property
def mol(self):
"""
Return the atom molecule index
:type: int
"""
return self.get("mol", self.index)
@property
def mass(self):
"""
Return the atom mass
:type: float
"""
return self.get("mass", self.index)
@property
def radius(self):
"""
Return the particle radius
:type: float
"""
return self.get("radius", self.index)
@property
def position(self):
"""
:getter: Return position of atom
:setter: Set position of atom
:type: numpy.array (float, float, float)
"""
return self.get("x", self.index)
@position.setter
def position(self, value):
current = self.position
current[:] = value
@property
def velocity(self):
"""
:getter: Return velocity of atom
:setter: Set velocity of atom
:type: numpy.array (float, float, float)
"""
return self.get("v", self.index)
@velocity.setter
def velocity(self, value):
current = self.velocity
current[:] = value
@property
def force(self):
"""
Return the total force acting on the atom
:type: numpy.array (float, float, float)
"""
return self.get("f", self.index)
@force.setter
def force(self, value):
current = self.force
current[:] = value
@property
def torque(self):
"""
Return the total torque acting on the atom
:type: numpy.array (float, float, float)
"""
return self.get("torque", self.index)
@force.setter
def torque(self, value):
current = self.torque
current[:] = value
@property
def omega(self):
"""
Return the rotational velocity of the particle
:type: numpy.array (float, float, float)
"""
return self.get("torque", self.index)
@omega.setter
def omega(self, value):
current = self.torque
current[:] = value
@property
def torque(self):
"""
Return the total torque acting on the particle
:type: numpy.array (float, float, float)
"""
return self.get("torque", self.index)
@torque.setter
def torque(self, value):
current = self.torque
current[:] = value
@property
def angular_momentum(self):
"""
Return the angular momentum of the particle
:type: numpy.array (float, float, float)
"""
return self.get("angmom", self.index)
@angular_momentum.setter
def angular_momentum(self, value):
current = self.angular_momentum
current[:] = value
@property
def charge(self):
"""
Return the atom charge
:type: float
"""
return self.get("q", self.index)
# -------------------------------------------------------------------------
class Atom2D(Atom):
"""
A wrapper class then represents a single 2D atom inside of LAMMPS
Inherits all properties from the :py:class:`Atom` class, but returns 2D versions
of position, velocity, and force.
It provides access to properties of the atom and allows you to change some of them.
"""
def __init__(self, pylammps_instance, index):
super(Atom2D, self).__init__(pylammps_instance, index)
@property
def position(self):
"""Access to coordinates of an atom
:getter: Return position of atom
:setter: Set position of atom
:type: numpy.array (float, float)
"""
return super(Atom2D, self).position[0:2]
@position.setter
def position(self, value):
current = self.position
current[:] = value
@property
def velocity(self):
"""Access to velocity of an atom
:getter: Return velocity of atom
:setter: Set velocity of atom
:type: numpy.array (float, float)
"""
return super(Atom2D, self).velocity[0:2]
@velocity.setter
def velocity(self, value):
current = self.velocity
current[:] = value
@property
def force(self):
"""Access to force of an atom
:getter: Return force of atom
:setter: Set force of atom
:type: numpy.array (float, float)
"""
return super(Atom2D, self).force[0:2]
@force.setter
def force(self, value):
current = self.force
current[:] = value
# -------------------------------------------------------------------------
class variable_set:
def __init__(self, name, variable_dict):
self._name = name
array_pattern = re.compile(r"(?P<arr>.+)\[(?P<index>[0-9]+)\]")
for key, value in variable_dict.items():
m = array_pattern.match(key)
if m:
g = m.groupdict()
varname = g['arr']
idx = int(g['index'])
if varname not in self.__dict__:
self.__dict__[varname] = {}
self.__dict__[varname][idx] = value
else:
self.__dict__[key] = value
def __str__(self):
return "{}({})".format(self._name, ','.join(["{}={}".format(k, self.__dict__[k]) for k in self.__dict__.keys() if not k.startswith('_')]))
def __dir__(self):
return [k for k in self.__dict__.keys() if not k.startswith('_')]
def __repr__(self):
return self.__str__()
# -------------------------------------------------------------------------
def get_thermo_data(output):
""" traverse output of runs and extract thermo data columns """
if isinstance(output, str):
lines = output.splitlines()
else:
lines = output
runs = []
columns = []
in_run = False
current_run = {}
for line in lines:
if line.startswith("Per MPI rank memory allocation"):
in_run = True
elif in_run and len(columns) == 0:
# first line after memory usage are column names
columns = line.split()
current_run = {}
for col in columns:
current_run[col] = []
elif line.startswith("Loop time of "):
in_run = False
columns = []
thermo_data = variable_set('ThermoData', current_run)
r = {'thermo' : thermo_data }
runs.append(namedtuple('Run', list(r.keys()))(*list(r.values())))
elif in_run and len(columns) > 0:
items = line.split()
# Convert thermo output and store it.
# It must have the same number of columns and
# all of them must be convertible to floats.
# Otherwise we ignore the line
if len(items) == len(columns):
try:
values = [float(x) for x in items]
for i, col in enumerate(columns):
current_run[col].append(values[i])
except ValueError:
# cannot convert. must be a non-thermo output. ignore.
pass
return runs
# -------------------------------------------------------------------------
# -------------------------------------------------------------------------
class PyLammps(object):
"""
This is a Python wrapper class around the lower-level
:py:class:`lammps` class, exposing a more Python-like,
object-oriented interface for prototyping system inside of IPython and
Jupyter notebooks.
It either creates its own instance of :py:class:`lammps` or can be
initialized with an existing instance. The arguments are the same of the
lower-level interface. The original interface can still be accessed via
:py:attr:`PyLammps.lmp`.
:param name: "machine" name of the shared LAMMPS library ("mpi" loads ``liblammps_mpi.so``, "" loads ``liblammps.so``)
:type name: string
:param cmdargs: list of command line arguments to be passed to the :cpp:func:`lammps_open` function. The executable name is automatically added.
:type cmdargs: list
:param ptr: pointer to a LAMMPS C++ class instance when called from an embedded Python interpreter. None means load symbols from shared library.
:type ptr: pointer
:param comm: MPI communicator (as provided by `mpi4py <mpi4py_docs_>`_). ``None`` means use ``MPI_COMM_WORLD`` implicitly.
:type comm: MPI_Comm
:param verbose: print all LAMMPS output to stdout
:type verbose: bool
:ivar lmp: instance of original LAMMPS Python interface
:vartype lmp: :py:class:`lammps`
:ivar runs: list of completed runs, each storing the thermo output
:vartype run: list
"""
def __init__(self, name="", cmdargs=None, ptr=None, comm=None, verbose=False):
self.has_echo = False
self.verbose = verbose
if cmdargs:
if '-echo' in cmdargs:
idx = cmdargs.index('-echo')
# ensures that echo line is ignored during output capture
self.has_echo = idx+1 < len(cmdargs) and cmdargs[idx+1] in ('screen', 'both')
if ptr:
if isinstance(ptr,PyLammps):
self.lmp = ptr.lmp
elif isinstance(ptr,lammps):
self.lmp = ptr
else:
self.lmp = lammps(name=name,cmdargs=cmdargs,ptr=ptr,comm=comm)
else:
self.lmp = lammps(name=name,cmdargs=cmdargs,ptr=None,comm=comm)
print("LAMMPS output is captured by PyLammps wrapper")
self._cmd_history = []
self._enable_cmd_history = False
self.runs = []
def __enter__(self):
return self
def __exit__(self, ex_type, ex_value, ex_traceback):
self.close()
def __del__(self):
if self.lmp: self.lmp.close()
self.lmp = None
def close(self):
"""Explicitly delete a LAMMPS instance
This is a wrapper around the :py:meth:`lammps.close` of the Python interface.
"""
if self.lmp: self.lmp.close()
self.lmp = None
def version(self):
"""Return a numerical representation of the LAMMPS version in use.
This is a wrapper around the :py:meth:`lammps.version` function of the Python interface.
:return: version number
:rtype: int
"""
return self.lmp.version()
def file(self, file):
"""Read LAMMPS commands from a file.
This is a wrapper around the :py:meth:`lammps.file` function of the Python interface.
:param path: Name of the file/path with LAMMPS commands
:type path: string
"""
self.lmp.file(file)
@property
def enable_cmd_history(self):
"""
:getter: Return whether command history is saved
:setter: Set if command history should be saved
:type: bool
"""
return self._enable_cmd_history
@enable_cmd_history.setter
def enable_cmd_history(self, value):
"""
:getter: Return whether command history is saved
:setter: Set if command history should be saved
:type: bool
"""
self._enable_cmd_history = (value == True)
def write_script(self, filepath):
"""
Write LAMMPS script file containing all commands executed up until now
:param filepath: path to script file that should be written
:type filepath: string
"""
with open(filepath, "w") as f:
for cmd in self._cmd_history:
print(cmd, file=f)
def clear_cmd_history(self):
"""
Clear LAMMPS command history up to this point
"""
self._cmd_history = []
def command(self, cmd):
"""
Execute LAMMPS command
If :py:attr:`PyLammps.enable_cmd_history` is set to ``True``, commands executed
will be recorded. The entire command history can be written to a file using
:py:meth:`PyLammps.write_script()`. To clear the command history, use
:py:meth:`PyLammps.clear_cmd_history()`.
:param cmd: command string that should be executed
:type: cmd: string
"""
self.lmp.command(cmd)
if self.enable_cmd_history:
self._cmd_history.append(cmd)
def run(self, *args, **kwargs):
"""
Execute LAMMPS run command with given arguments
All thermo output during the run is captured and saved as new entry in
:py:attr:`PyLammps.runs`. The latest run can be retrieved by
:py:attr:`PyLammps.last_run`.
"""
output = self.__getattr__('run')(*args, **kwargs)
self.runs += get_thermo_data(output)
return output
@property
def last_run(self):
"""
Return data produced of last completed run command
:getter: Returns an object containing information about the last run command
:type: dict
"""
if len(self.runs) > 0:
return self.runs[-1]
return None
@property
def atoms(self):
"""
All atoms of this LAMMPS instance
:getter: Returns a list of atoms currently in the system
:type: AtomList
"""
return AtomList(self)
@property
def system(self):
"""
The system state of this LAMMPS instance
:getter: Returns an object with properties storing the current system state
:type: namedtuple
"""
output = self.lmp_info("system")
output = output[output.index("System information:")+1:]
d = self._parse_info_system(output)
return namedtuple('System', d.keys())(*d.values())
@property
def communication(self):
"""
The communication state of this LAMMPS instance
:getter: Returns an object with properties storing the current communication state
:type: namedtuple
"""
output = self.lmp_info("communication")
output = output[output.index("Communication information:")+1:]
d = self._parse_info_communication(output)
return namedtuple('Communication', d.keys())(*d.values())
@property
def computes(self):
"""
The list of active computes of this LAMMPS instance
:getter: Returns a list of computes that are currently active in this LAMMPS instance
:type: list
"""
output = self.lmp_info("computes")
output = output[output.index("Compute information:")+1:]
return self._parse_element_list(output)
@property
def dumps(self):
"""
The list of active dumps of this LAMMPS instance
:getter: Returns a list of dumps that are currently active in this LAMMPS instance
:type: list
"""
output = self.lmp_info("dumps")
output = output[output.index("Dump information:")+1:]
return self._parse_element_list(output)
@property
def fixes(self):
"""
The list of active fixes of this LAMMPS instance
:getter: Returns a list of fixes that are currently active in this LAMMPS instance
:type: list
"""
output = self.lmp_info("fixes")
output = output[output.index("Fix information:")+1:]
return self._parse_element_list(output)
@property
def groups(self):
"""
The list of active atom groups of this LAMMPS instance
:getter: Returns a list of atom groups that are currently active in this LAMMPS instance
:type: list
"""
output = self.lmp_info("groups")
output = output[output.index("Group information:")+1:]
return self._parse_groups(output)
@property
def variables(self):
"""
Returns a dictionary of all variables defined in the current LAMMPS instance
:getter: Returns a dictionary of all variables that are defined in this LAMMPS instance
:type: dict
"""
output = self.lmp_info("variables")
output = output[output.index("Variable information:")+1:]
variables = {}
for v in self._parse_element_list(output):
variables[v['name']] = Variable(self, v['name'], v['style'], v['def'])
return variables
def eval(self, expr):
"""
Evaluate expression
:param expr: the expression string that should be evaluated inside of LAMMPS
:type expr: string
:return: the value of the evaluated expression
:rtype: float if numeric, string otherwise
"""
value = self.lmp_print('"$(%s)"' % expr).strip()
try:
return float(value)
except ValueError:
return value
def _split_values(self, line):
return [x.strip() for x in line.split(',')]
def _get_pair(self, value):
return [x.strip() for x in value.split('=')]
def _parse_info_system(self, output):
system = {}
for line in output:
if line.startswith("Units"):
system['units'] = self._get_pair(line)[1]
elif line.startswith("Atom style"):
system['atom_style'] = self._get_pair(line)[1]
elif line.startswith("Atom map"):
system['atom_map'] = self._get_pair(line)[1]
elif line.startswith("Atoms"):
parts = self._split_values(line)
system['natoms'] = int(self._get_pair(parts[0])[1])
system['ntypes'] = int(self._get_pair(parts[1])[1])
system['style'] = self._get_pair(parts[2])[1]
elif line.startswith("Kspace style"):
system['kspace_style'] = self._get_pair(line)[1]
elif line.startswith("Dimensions"):
system['dimensions'] = int(self._get_pair(line)[1])
elif line.startswith("Orthogonal box"):
system['orthogonal_box'] = [float(x) for x in self._get_pair(line)[1].split('x')]
elif line.startswith("Boundaries"):
system['boundaries'] = self._get_pair(line)[1]
elif line.startswith("xlo"):
keys, values = [self._split_values(x) for x in self._get_pair(line)]
for key, value in zip(keys, values):
system[key] = float(value)
elif line.startswith("ylo"):
keys, values = [self._split_values(x) for x in self._get_pair(line)]
for key, value in zip(keys, values):
system[key] = float(value)
elif line.startswith("zlo"):
keys, values = [self._split_values(x) for x in self._get_pair(line)]
for key, value in zip(keys, values):
system[key] = float(value)
elif line.startswith("Molecule type"):
system['molecule_type'] = self._get_pair(line)[1]
elif line.startswith("Bonds"):
parts = self._split_values(line)
system['nbonds'] = int(self._get_pair(parts[0])[1])
system['nbondtypes'] = int(self._get_pair(parts[1])[1])
system['bond_style'] = self._get_pair(parts[2])[1]
elif line.startswith("Angles"):
parts = self._split_values(line)
system['nangles'] = int(self._get_pair(parts[0])[1])
system['nangletypes'] = int(self._get_pair(parts[1])[1])
system['angle_style'] = self._get_pair(parts[2])[1]
elif line.startswith("Dihedrals"):
parts = self._split_values(line)
system['ndihedrals'] = int(self._get_pair(parts[0])[1])
system['ndihedraltypes'] = int(self._get_pair(parts[1])[1])
system['dihedral_style'] = self._get_pair(parts[2])[1]
elif line.startswith("Impropers"):
parts = self._split_values(line)
system['nimpropers'] = int(self._get_pair(parts[0])[1])
system['nimpropertypes'] = int(self._get_pair(parts[1])[1])
system['improper_style'] = self._get_pair(parts[2])[1]
return system
def _parse_info_communication(self, output):
comm = {}
for line in output:
if line.startswith("MPI library"):
comm['mpi_version'] = line.split(':')[1].strip()
elif line.startswith("Comm style"):
parts = self._split_values(line)
comm['comm_style'] = self._get_pair(parts[0])[1]
comm['comm_layout'] = self._get_pair(parts[1])[1]
elif line.startswith("Processor grid"):
comm['proc_grid'] = [int(x) for x in self._get_pair(line)[1].split('x')]
elif line.startswith("Communicate velocities for ghost atoms"):
comm['ghost_velocity'] = (self._get_pair(line)[1] == "yes")
elif line.startswith("Nprocs"):
parts = self._split_values(line)
comm['nprocs'] = int(self._get_pair(parts[0])[1])
comm['nthreads'] = int(self._get_pair(parts[1])[1])
return comm
def _parse_element_list(self, output):
elements = []
for line in output:
if not line or (":" not in line): continue
element_info = self._split_values(line.split(':')[1].strip())
element = {'name': element_info[0]}
for key, value in [self._get_pair(x) for x in element_info[1:]]:
element[key] = value
elements.append(element)
return elements
def _parse_groups(self, output):
groups = []
group_pattern = re.compile(r"(?P<name>.+) \((?P<type>.+)\)")
for line in output:
m = group_pattern.match(line.split(':')[1].strip())
group = {'name': m.group('name'), 'type': m.group('type')}
groups.append(group)
return groups
def lmp_print(self, s):
""" needed for Python2 compatibility, since print is a reserved keyword """
return self.__getattr__("print")(s)
def __dir__(self):
return sorted(set(['angle_coeff', 'angle_style', 'atom_modify', 'atom_style', 'atom_style',
'bond_coeff', 'bond_style', 'boundary', 'change_box', 'communicate', 'compute',
'create_atoms', 'create_box', 'delete_atoms', 'delete_bonds', 'dielectric',
'dihedral_coeff', 'dihedral_style', 'dimension', 'dump', 'fix', 'fix_modify',
'group', 'improper_coeff', 'improper_style', 'include', 'kspace_modify',
'kspace_style', 'lattice', 'mass', 'minimize', 'min_style', 'neighbor',
'neigh_modify', 'newton', 'nthreads', 'pair_coeff', 'pair_modify',
'pair_style', 'processors', 'read', 'read_data', 'read_restart', 'region',
'replicate', 'reset_timestep', 'restart', 'run', 'run_style', 'thermo',
'thermo_modify', 'thermo_style', 'timestep', 'undump', 'unfix', 'units',
'variable', 'velocity', 'write_restart'] + self.lmp.available_styles("command")))
def lmp_info(self, s):
# skip anything before and after Info-Info-Info
# also skip timestamp line
output = self.__getattr__("info")(s)
indices = [index for index, line in enumerate(output) if line.startswith("Info-Info-Info-Info")]
start = indices[0]
end = indices[1]
return [line for line in output[start+2:end] if line]
def __getattr__(self, name):
"""
This method is where the Python 'magic' happens. If a method is not
defined by the class PyLammps, it assumes it is a LAMMPS command. It takes
all the arguments, concatinates them to a single string, and executes it using
:py:meth:`lammps.PyLammps.command()`.
:param verbose: Print output of command
:type verbose: bool
:return: line or list of lines of output, None if no output
:rtype: list or string
"""
def handler(*args, **kwargs):
cmd_args = [name] + [str(x) for x in args]
self.lmp.flush_buffers()
with OutputCapture() as capture:
cmd = ' '.join(cmd_args)
self.command(cmd)
self.lmp.flush_buffers()
output = capture.output
comm = self.lmp.get_mpi_comm()
if comm:
output = self.lmp.comm.bcast(output, root=0)
if self.verbose or ('verbose' in kwargs and kwargs['verbose']):
print(output, end = '')
lines = output.splitlines()
if self.has_echo:
lines = lines[1:]
if len(lines) > 1:
return lines
elif len(lines) == 1:
return lines[0]
return None
return handler
class IPyLammps(PyLammps):
"""
IPython wrapper for LAMMPS which adds embedded graphics capabilities to PyLammmps interface
It either creates its own instance of :py:class:`lammps` or can be
initialized with an existing instance. The arguments are the same of the
lower-level interface. The original interface can still be accessed via
:py:attr:`PyLammps.lmp`.
:param name: "machine" name of the shared LAMMPS library ("mpi" loads ``liblammps_mpi.so``, "" loads ``liblammps.so``)
:type name: string
:param cmdargs: list of command line arguments to be passed to the :cpp:func:`lammps_open` function. The executable name is automatically added.
:type cmdargs: list
:param ptr: pointer to a LAMMPS C++ class instance when called from an embedded Python interpreter. None means load symbols from shared library.
:type ptr: pointer
:param comm: MPI communicator (as provided by `mpi4py <mpi4py_docs_>`_). ``None`` means use ``MPI_COMM_WORLD`` implicitly.
:type comm: MPI_Comm
"""
def __init__(self,name="",cmdargs=None,ptr=None,comm=None):
super(IPyLammps, self).__init__(name=name,cmdargs=cmdargs,ptr=ptr,comm=comm)
def image(self, filename="snapshot.png", group="all", color="type", diameter="type",
size=None, view=None, center=None, up=None, zoom=1.0, background_color="white"):
""" Generate image using write_dump command and display it
See :doc:`dump image <dump_image>` for more information.
:param filename: Name of the image file that should be generated. The extension determines whether it is PNG or JPEG
:type filename: string
:param group: the group of atoms write_image should use
:type group: string
:param color: name of property used to determine color
:type color: string
:param diameter: name of property used to determine atom diameter
:type diameter: string
:param size: dimensions of image
:type size: tuple (width, height)
:param view: view parameters
:type view: tuple (theta, phi)
:param center: center parameters
:type center: tuple (flag, center_x, center_y, center_z)
:param up: vector pointing to up direction
:type up: tuple (up_x, up_y, up_z)
:param zoom: zoom factor
:type zoom: float
:param background_color: background color of scene
:type background_color: string
:return: Image instance used to display image in notebook
:rtype: :py:class:`IPython.core.display.Image`
"""
cmd_args = [group, "image", filename, color, diameter]
if size is not None:
width = size[0]
height = size[1]
cmd_args += ["size", width, height]
if view is not None:
theta = view[0]
phi = view[1]
cmd_args += ["view", theta, phi]
if center is not None:
flag = center[0]
Cx = center[1]
Cy = center[2]
Cz = center[3]
cmd_args += ["center", flag, Cx, Cy, Cz]
if up is not None:
Ux = up[0]
Uy = up[1]
Uz = up[2]
cmd_args += ["up", Ux, Uy, Uz]
if zoom is not None:
cmd_args += ["zoom", zoom]
cmd_args.append("modify backcolor " + background_color)
self.write_dump(*cmd_args)
from IPython.core.display import Image
return Image(filename)
def video(self, filename):
"""
Load video from file
Can be used to visualize videos from :doc:`dump movie <dump_image>`.
:param filename: Path to video file
:type filename: string
:return: HTML Video Tag used by notebook to embed a video
:rtype: :py:class:`IPython.display.HTML`
"""
from IPython.display import HTML
return HTML("<video controls><source src=\"" + filename + "\"></video>")