Skip to content

Commit

Permalink
Merge pull request #295 from pyiron/interactive_units
Browse files Browse the repository at this point in the history
Interactive units
  • Loading branch information
sudarsan-surendralal authored Aug 6, 2021
2 parents 2957fc7 + 7751f51 commit b5b92da
Show file tree
Hide file tree
Showing 5 changed files with 221 additions and 14 deletions.
162 changes: 162 additions & 0 deletions notebooks/integration/interactive_units.ipynb
Original file line number Diff line number Diff line change
@@ -0,0 +1,162 @@
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This notebooks compares the outputs from two static LAMMPS jobs one of which is an interactive job. This is to ensure that the outputs from interactive and non-interactive jobs are consistent."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"from pyiron import Project\n",
"import numpy as np"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"pr = Project(\"water_interactive\")\n",
"pr.remove_jobs_silently()"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"dx = 0.7\n",
"cell = np.eye(3) * 10\n",
"r_O = [0, 0, 0]\n",
"r_H1 = [dx, dx, 0]\n",
"r_H2 = [-dx, dx, 0]\n",
"water = pr.create_atoms(elements=['H', 'H', 'O'],\n",
" positions=[r_H1, r_H2, r_O],\n",
" cell=cell, pbc=True)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/cmmc/u/chandu/programs/pyiron_mpie/pyiron_atomistics/pyiron_atomistics/lammps/base.py:210: UserWarning: WARNING: Non-'metal' units are not fully supported. Your calculation should run OK, but results may not be saved in pyiron units.\n",
" \"WARNING: Non-'metal' units are not fully supported. Your calculation should run OK, but \"\n",
"/cmmc/u/chandu/programs/pyiron_mpie/pyiron_atomistics/pyiron_atomistics/lammps/base.py:262: UserWarning: No potential set via job.potential - use default potential, H2O_tip3p\n",
" warnings.warn(\"No potential set via job.potential - use default potential, \" + lst_of_potentials[0])\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"The job test was saved and received the ID: 15458456\n"
]
}
],
"source": [
"# Interactive job\n",
"job_int = pr.create.job.Lammps(\"test\", delete_existing_job=True)\n",
"job_int.structure = water\n",
"job_int.interactive_open()\n",
"job_int.calc_static()\n",
"job_int.run()\n",
"job_int.interactive_close()"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"The job test_ni was saved and received the ID: 15458457\n"
]
}
],
"source": [
"# Non-interactive job\n",
"job = pr.create.job.Lammps(\"test_ni\", delete_existing_job=True)\n",
"job.structure = water\n",
"job.calc_static()\n",
"job.run()"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"cells\n",
"energy_pot\n",
"energy_tot\n",
"forces\n",
"indices\n",
"positions\n",
"pressures\n",
"steps\n",
"temperature\n",
"volume\n"
]
}
],
"source": [
"# Assert that the unit converstions work even in the interactive mode\n",
"\n",
"int_nodes = job_int[\"output/generic\"].list_nodes()\n",
"usual_nodes = job[\"output/generic\"].list_nodes()\n",
"for node in int_nodes:\n",
" if node in usual_nodes:\n",
" print(node)\n",
" assert np.allclose(job_int[\"output/generic/\" + node], job[\"output/generic/\" + node])"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.7.9"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
30 changes: 26 additions & 4 deletions pyiron_atomistics/lammps/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,7 @@
from pyiron_atomistics.lammps.control import LammpsControl
from pyiron_atomistics.lammps.potential import LammpsPotential
from pyiron_atomistics.lammps.structure import LammpsStructure, UnfoldingPrism
from pyiron_atomistics.lammps.units import UnitConverter
from pyiron_atomistics.lammps.units import UnitConverter, LAMMPS_UNIT_CONVERSIONS


__author__ = "Joerg Neugebauer, Sudarsan Surendralal, Jan Janssen"
Expand Down Expand Up @@ -61,6 +61,28 @@ def __init__(self, project, job_name):
self._prism = None
s.publication_add(self.publication)

@property
def units(self):
"""
Type of LAMMPS units used in the calculations. Can be either of 'metal', 'real', 'si', 'cgs', and 'lj'
Returns:
str: Type of LAMMPS unit
"""
if self.input.control["units"] is not None:
return self.input.control["units"]
else:
# Default to metal units
return "metal"

@units.setter
def units(self, val):
allowed_types = LAMMPS_UNIT_CONVERSIONS.keys()
if val in allowed_types:
self.input.control["units"] = val
else:
raise ValueError("'{}' is not a valid LAMMPS unit")

@property
def bond_dict(self):
"""
Expand Down Expand Up @@ -449,7 +471,7 @@ def collect_h5md_file(self, file_name="dump.h5", cwd=None):
Returns:
"""
uc = UnitConverter(self.input.control["units"])
uc = UnitConverter(self.units)
prism = UnfoldingPrism(self.structure.cell, digits=15)
if np.matrix.trace(prism.R) != 3:
raise RuntimeError("The Lammps output will not be mapped back to pyiron correctly.")
Expand Down Expand Up @@ -537,7 +559,7 @@ def collect_output_log(self, file_name="log.lammps", cwd=None):
Returns:
"""
uc = UnitConverter(self.input.control["units"])
uc = UnitConverter(self.units)
self.collect_errors(file_name=file_name, cwd=cwd)
file_name = self.job_file_name(file_name=file_name, cwd=cwd)
if os.path.exists(file_name):
Expand Down Expand Up @@ -905,7 +927,7 @@ def collect_dump_file(self, file_name="dump.out", cwd=None):
Returns:
"""
uc = UnitConverter(self.input.control["units"])
uc = UnitConverter(self.units)
file_name = self.job_file_name(file_name=file_name, cwd=cwd)
if os.path.exists(file_name):
output = {}
Expand Down
29 changes: 21 additions & 8 deletions pyiron_atomistics/lammps/interactive.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@
from pylammpsmpi import LammpsLibrary
except ImportError:
pass
from pyiron_atomistics.lammps.units import UnitConverter

__author__ = "Osamu Waseda, Jan Janssen"
__copyright__ = (
Expand Down Expand Up @@ -56,12 +57,14 @@ def _interactive_lib_command(self, command):
self._interactive_library.command(command)

def interactive_positions_getter(self):
uc = UnitConverter(units=self.units)
positions = np.reshape(
np.array(self._interactive_library.gather_atoms("x", 1, 3)),
(len(self.structure), 3),
)
if np.matrix.trace(self._prism.R) != 3:
positions = np.matmul(positions, self._prism.R.T)
positions = uc.convert_array_to_pyiron_units(positions, label="positions")
return positions.tolist()

def interactive_positions_setter(self, positions):
Expand All @@ -78,6 +81,7 @@ def interactive_positions_setter(self, positions):
self._interactive_lib_command("change_box all remap")

def interactive_cells_getter(self):
uc = UnitConverter(units=self.units)
cc = np.array(
[
[self._interactive_library.get_thermo("lx"), 0, 0],
Expand All @@ -93,7 +97,7 @@ def interactive_cells_getter(self):
],
]
)
return self._prism.unfold_cell(cc)
return uc.convert_array_to_pyiron_units(self._prism.unfold_cell(cc), label="cells")

def interactive_cells_setter(self, cell):
self._prism = UnfoldingPrism(cell)
Expand Down Expand Up @@ -132,15 +136,18 @@ def interactive_cells_setter(self, cell):
)

def interactive_volume_getter(self):
return self._interactive_library.get_thermo("vol")
uc = UnitConverter(units=self.units)
return uc.convert_array_to_pyiron_units(self._interactive_library.get_thermo("vol"), label="volume")

def interactive_forces_getter(self):
uc = UnitConverter(units=self.units)
ff = np.reshape(
np.array(self._interactive_library.gather_atoms("f", 1, 3)),
(len(self.structure), 3),
)
if np.matrix.trace(self._prism.R) != 3:
ff = np.matmul(ff, self._prism.R.T)
ff = uc.convert_array_to_pyiron_units(ff, label="forces")
return ff.tolist()

def interactive_execute(self):
Expand Down Expand Up @@ -534,8 +541,9 @@ def update_potential(self):
self._interactive_lib_command(self.potential.Config[0][1])

def interactive_indices_getter(self):
uc = UnitConverter(units=self.units)
lammps_indices = np.array(self._interactive_library.gather_atoms("type", 0, 1))
indices = self.remap_indices(lammps_indices)
indices = uc.convert_array_to_pyiron_units(self.remap_indices(lammps_indices), label="indices")
return indices.tolist()

def interactive_indices_setter(self, indices):
Expand All @@ -559,16 +567,20 @@ def interactive_indices_setter(self, indices):
self._interactive_library.scatter_atoms("type", elem_all)

def interactive_energy_pot_getter(self):
return self._interactive_library.get_thermo("pe")
uc = UnitConverter(units=self.units)
return uc.convert_array_to_pyiron_units(self._interactive_library.get_thermo("pe"), label="energy_pot")

def interactive_energy_tot_getter(self):
return self._interactive_library.get_thermo("etotal")
uc = UnitConverter(units=self.units)
return uc.convert_array_to_pyiron_units(self._interactive_library.get_thermo("etotal"), label="energy_tot")

def interactive_steps_getter(self):
return self._interactive_library.get_thermo("step")
uc = UnitConverter(units=self.units)
return uc.convert_array_to_pyiron_units(self._interactive_library.get_thermo("step"), label="steps")

def interactive_temperatures_getter(self):
return self._interactive_library.get_thermo("temp")
uc = UnitConverter(units=self.units)
return uc.convert_array_to_pyiron_units(self._interactive_library.get_thermo("temp"), label="temperature")

def interactive_stress_getter(self):
"""
Expand All @@ -595,6 +607,7 @@ def interactive_stress_getter(self):
return ss

def interactive_pressures_getter(self):
uc = UnitConverter(units=self.units)
pp = np.array(
[
[
Expand All @@ -617,7 +630,7 @@ def interactive_pressures_getter(self):
rotation_matrix = self._prism.R.T
if np.matrix.trace(rotation_matrix) != 3:
pp = rotation_matrix.T @ pp @ rotation_matrix
return pp / 10000 # bar -> GPa
return uc.convert_array_to_pyiron_units(pp, label="pressure")

def interactive_close(self):
if self.interactive_is_activated():
Expand Down
2 changes: 1 addition & 1 deletion pyiron_atomistics/lammps/units.py
Original file line number Diff line number Diff line change
Expand Up @@ -204,7 +204,7 @@ def convert_array_to_pyiron_units(self, array, label):
label (str): The label of the quantity (must be a key in the dictionary `quantity_dict`)
Returns:
numpy.ndarray: The array after conversion
ndarray: The array after conversion
"""
if label in quantity_dict.keys():
Expand Down
12 changes: 11 additions & 1 deletion tests/lammps/test_base.py
Original file line number Diff line number Diff line change
Expand Up @@ -332,7 +332,6 @@ def test_dump_parser_water(self):

self.assertAlmostEqual(self.job_water_dump["output/generic/pressures"][-2][0, 0], 515832.570508186 /
uc.pyiron_to_lammps("pressure"), 2)

self.job_water_dump.write_traj(filename="test.xyz",
file_format="xyz")
atom_indices = self.job_water_dump.structure.select_index("H")
Expand Down Expand Up @@ -427,6 +426,7 @@ def test_dump_parser_water(self):
self.assertTrue(np.allclose(neigh_traj_obj.indices, neigh_traj_obj_loaded.indices))
self.assertTrue(np.allclose(neigh_traj_obj.distances, neigh_traj_obj_loaded.distances))
self.assertTrue(np.allclose(neigh_traj_obj.vecs, neigh_traj_obj_loaded.vecs))
self.assertTrue(self.job_water_dump.units, "real")

def test_dump_parser(self):
structure = Atoms(
Expand Down Expand Up @@ -723,5 +723,15 @@ def test_potential_check(self):
potential['Species'][0][0] = 'Al'
self.job.potential = potential # shouldn't raise ValueError

def test_units(self):
self.assertTrue(self.job.units, "metal")
self.job.units = "real"
self.assertTrue(self.job.units, "real")

def setter(x):
self.job.units = x
self.assertRaises(ValueError, setter, "nonsense")


if __name__ == "__main__":
unittest.main()

0 comments on commit b5b92da

Please sign in to comment.