diff --git a/pyiron_lammps/compatibility/file.py b/pyiron_lammps/compatibility/file.py index fd3a5c1..709b873 100644 --- a/pyiron_lammps/compatibility/file.py +++ b/pyiron_lammps/compatibility/file.py @@ -88,18 +88,24 @@ def lammps_file_interface_function( calc_kwargs = {} os.makedirs(working_directory, exist_ok=True) - if isinstance(potential, str): - potential_dataframe = get_potential_by_name( - potential_name=potential, resource_path=resource_path - ) - elif isinstance(potential, pandas.DataFrame): - potential_dataframe = potential.iloc[0] - elif isinstance(potential, pandas.Series): - potential_dataframe = potential - else: - raise TypeError() - lmp_str_lst = lammps_file_initialization(structure=structure) - lmp_str_lst += potential_dataframe["Config"] + potential_lst, potential_replace, species = _get_potential( + potential=potential, resource_path=resource_path + ) + + lmp_str_lst = [] + atom_type = "atomic" + for l in lammps_file_initialization(structure=structure): + if l.startswith("units") and "units" in potential_replace: + lmp_str_lst.append(potential_replace["units"]) + elif l.startswith("atom_style") and "atom_style" in potential_replace: + lmp_str_lst.append(potential_replace["atom_style"]) + atom_type = potential_replace["atom_style"].split()[-1] + elif l.startswith("dimension") and "dimension" in potential_replace: + lmp_str_lst.append(potential_replace["dimension"]) + else: + lmp_str_lst.append(l) + + lmp_str_lst += potential_lst lmp_str_lst += ["variable dumptime equal {} ".format(calc_kwargs.get("n_print", 1))] lmp_str_lst += [ "dump 1 all custom ${dumptime} dump.out id type xsu ysu zsu fx fy fz vx vy vz", @@ -161,11 +167,12 @@ def lammps_file_interface_function( write_lammps_datafile( structure=structure, - potential_elements=potential_dataframe["Species"], + potential_elements=species, bond_dict=None, units=units, file_name="lammps.data", working_directory=working_directory, + atom_type=atom_type, ) shell = subprocess.check_output( @@ -178,7 +185,7 @@ def lammps_file_interface_function( output = parse_lammps_output( working_directory=working_directory, structure=structure, - potential_elements=potential_dataframe["Species"], + potential_elements=species, units=units, prism=None, dump_h5_file_name="dump.h5", @@ -222,3 +229,30 @@ def _modify_input_dict( return lmp_tmp_lst else: return lmp_str_lst + + +def _get_potential(potential, resource_path: Optional[str] = None): + if isinstance(potential, str): + potential_dataframe = get_potential_by_name( + potential_name=potential, resource_path=resource_path + ) + elif isinstance(potential, pandas.DataFrame): + potential_dataframe = potential.iloc[0] + elif isinstance(potential, pandas.Series): + potential_dataframe = potential + else: + raise TypeError() + + potential_replace = {} + potential_lst = [] + for l in potential_dataframe["Config"]: + if l.startswith("units"): + potential_replace["units"] = l + elif l.startswith("atom_style"): + potential_replace["atom_style"] = l + elif l.startswith("dimension"): + potential_replace["dimension"] = l + else: + potential_lst.append(l) + + return potential_lst, potential_replace, potential_dataframe["Species"] diff --git a/pyiron_lammps/potential.py b/pyiron_lammps/potential.py index e9df255..2f338f1 100644 --- a/pyiron_lammps/potential.py +++ b/pyiron_lammps/potential.py @@ -1,6 +1,6 @@ import os from pathlib import Path -from typing import Union +from typing import Optional, Union import pandas from ase.atoms import Atoms @@ -346,7 +346,7 @@ def get_potential_dataframe(structure: Atoms, resource_path=None): ) -def get_potential_by_name(potential_name: str, resource_path=None): +def get_potential_by_name(potential_name: str, resource_path: Optional[str] = None): if resource_path is None: resource_path = get_resource_path_from_conda() df = LammpsPotentialFile(resource_path=resource_path).list() diff --git a/tests/test_compatibility_file.py b/tests/test_compatibility_file.py index 57dde17..93cef2b 100644 --- a/tests/test_compatibility_file.py +++ b/tests/test_compatibility_file.py @@ -3,7 +3,10 @@ import shutil from ase.build import bulk import pandas -from pyiron_lammps.compatibility.file import lammps_file_interface_function +from pyiron_lammps.compatibility.file import ( + lammps_file_interface_function, + _get_potential, +) from pyiron_lammps.potential import get_potential_by_name @@ -271,10 +274,43 @@ def test_calc_md_nvt_langevin(self): def test_calc_md_nve(self): calc_kwargs = {"n_print": 100, "langevin": True} + potential = [ + "# Bouhadja et al., J. Chem. Phys. 138, 224510 (2013) \n", + "units metal\n", + "dimension 3\n", + "atom_style charge\n", + "\n", + "# create groups ###\n", + "group Al type 1\n", + "group Ca type 2\n", + "group O type 3\n", + "group Si type 4\n", + "\n### set charges ###\n", + "set type 1 charge 1.8\n", + "set type 2 charge 1.2\n", + "set type 3 charge -1.2\n", + "set type 4 charge 2.4\n", + "\n### Bouhadja Born-Mayer-Huggins + Coulomb Potential Parameters ###\n", + "pair_style born/coul/dsf 0.25 8.0\n", + "pair_coeff 1 1 0.002900 0.068000 1.570400 14.049800 0.000000\n", + "pair_coeff 1 2 0.003200 0.074000 1.957200 17.171000 0.000000\n", + "pair_coeff 1 3 0.007500 0.164000 2.606700 34.574700 0.000000\n", + "pair_coeff 1 4 0.002500 0.057000 1.505600 18.811600 0.000000\n", + "pair_coeff 2 2 0.003500 0.080000 2.344000 20.985600 0.000000\n", + "pair_coeff 2 3 0.007700 0.178000 2.993500 42.255600 0.000000\n", + "pair_coeff 2 4 0.002700 0.063000 1.892400 22.990700 0.000000\n", + "pair_coeff 3 3 0.012000 0.263000 3.643000 85.084000 0.000000\n", + "pair_coeff 3 4 0.007000 0.156000 2.541900 46.293000 0.000000\n", + "pair_coeff 4 4 0.001200 0.046000 1.440800 25.187300 0.000000\n", + "\npair_modify shift yes\n", + ] + element_lst = ["Al", "Ca", "O", "Si"] shell_output, parsed_output, job_crashed = lammps_file_interface_function( working_directory=self.working_dir, structure=self.structure, - potential=self.potential, + potential=pandas.DataFrame( + {"Config": [potential], "Species": [element_lst]} + ), calc_mode="md", calc_kwargs=calc_kwargs, units=self.units, @@ -292,15 +328,40 @@ def test_calc_md_nve(self): "units metal\n", "dimension 3\n", "boundary p p p\n", - "atom_style atomic\n", + "atom_style charge\n", + "\n", "read_data lammps.data\n", - "pair_style eam/alloy\n", + "# Bouhadja et al., J. Chem. Phys. 138, 224510 (2013) \n", + "# create groups ###\n", + "group Al type 1\n", + "group Ca type 2\n", + "group O type 3\n", + "group Si type 4\n", + "### set charges ###\n", + "set type 1 charge 1.8\n", + "set type 2 charge 1.2\n", + "set type 3 charge -1.2\n", + "set type 4 charge 2.4\n", + "### Bouhadja Born-Mayer-Huggins + Coulomb Potential Parameters ###\n", + "pair_style born/coul/dsf 0.25 8.0\n", + "pair_coeff 1 1 0.002900 0.068000 1.570400 14.049800 0.000000\n", + "pair_coeff 1 2 0.003200 0.074000 1.957200 17.171000 0.000000\n", + "pair_coeff 1 3 0.007500 0.164000 2.606700 34.574700 0.000000\n", + "pair_coeff 1 4 0.002500 0.057000 1.505600 18.811600 0.000000\n", + "pair_coeff 2 2 0.003500 0.080000 2.344000 20.985600 0.000000\n", + "pair_coeff 2 3 0.007700 0.178000 2.993500 42.255600 0.000000\n", + "pair_coeff 2 4 0.002700 0.063000 1.892400 22.990700 0.000000\n", + "pair_coeff 3 3 0.012000 0.263000 3.643000 85.084000 0.000000\n", + "pair_coeff 3 4 0.007000 0.156000 2.541900 46.293000 0.000000\n", + "pair_coeff 4 4 0.001200 0.046000 1.440800 25.187300 0.000000\n", + "pair_modify shift yes\n", "variable dumptime equal 100 \n", "dump 1 all custom ${dumptime} dump.out id type xsu ysu zsu fx fy fz vx vy vz\n", 'dump_modify 1 sort id format line "%d %d %20.15g %20.15g %20.15g %20.15g %20.15g %20.15g %20.15g %20.15g %20.15g"\n', "fix ensemble all nve\n", "variable thermotime equal 100 \n", "timestep 0.001\n", + "velocity all create None 80996 dist gaussian\n", "thermo_style custom step temp pe etotal pxx pxy pxz pyy pyz pzz vol\n", "thermo_modify format float %20.15g\n", "thermo ${thermotime}\n", @@ -460,3 +521,63 @@ def test_calc_minimize_pressure_3d(self): ] for line in content_expected: self.assertIn(line, content) + + +class TestGlassPotential(unittest.TestCase): + def setUp(self): + self.static_path = os.path.abspath( + os.path.join("..", os.path.dirname(__file__), "static") + ) + + def test_bouhadja(self): + potential = [ + "# Bouhadja et al., J. Chem. Phys. 138, 224510 (2013) \n", + "units metal\n", + "dimension 3\n", + "atom_style charge\n", + "\n", + "# create groups ###\n", + "group Al type 1\n", + "group Ca type 2\n", + "group O type 3\n", + "group Si type 4\n", + "\n### set charges ###\n", + "set type 1 charge 1.8\n", + "set type 2 charge 1.2\n", + "set type 3 charge -1.2\n", + "set type 4 charge 2.4\n", + "\n### Bouhadja Born-Mayer-Huggins + Coulomb Potential Parameters ###\n", + "pair_style born/coul/dsf 0.25 8.0\n", + "pair_coeff 1 1 0.002900 0.068000 1.570400 14.049800 0.000000\n", + "pair_coeff 1 2 0.003200 0.074000 1.957200 17.171000 0.000000\n", + "pair_coeff 1 3 0.007500 0.164000 2.606700 34.574700 0.000000\n", + "pair_coeff 1 4 0.002500 0.057000 1.505600 18.811600 0.000000\n", + "pair_coeff 2 2 0.003500 0.080000 2.344000 20.985600 0.000000\n", + "pair_coeff 2 3 0.007700 0.178000 2.993500 42.255600 0.000000\n", + "pair_coeff 2 4 0.002700 0.063000 1.892400 22.990700 0.000000\n", + "pair_coeff 3 3 0.012000 0.263000 3.643000 85.084000 0.000000\n", + "pair_coeff 3 4 0.007000 0.156000 2.541900 46.293000 0.000000\n", + "pair_coeff 4 4 0.001200 0.046000 1.440800 25.187300 0.000000\n", + "\npair_modify shift yes\n", + ] + element_lst = ["Al", "Ca", "O", "Si"] + potential_lst, potential_replace, species = _get_potential( + potential=pandas.DataFrame( + {"Config": [potential], "Species": [element_lst]} + ), + resource_path=os.path.join(self.static_path, "potential"), + ) + for i, l in enumerate(potential): + if i in [1, 2, 3]: + self.assertFalse(l in potential_lst) + else: + self.assertTrue(l in potential_lst) + + for k, v in { + "units": "units metal\n", + "dimension": "dimension 3\n", + "atom_style": "atom_style charge\n", + }.items(): + self.assertEqual(potential_replace[k], v) + + self.assertEqual(species, element_lst)