diff --git a/pyiron_atomistics/sphinx/base.py b/pyiron_atomistics/sphinx/base.py index 706f9ecff..77e771def 100644 --- a/pyiron_atomistics/sphinx/base.py +++ b/pyiron_atomistics/sphinx/base.py @@ -29,7 +29,15 @@ ) from pyiron_atomistics.sphinx.structure import read_atoms from pyiron_atomistics.sphinx.potential import SphinxJTHPotentialFile -from pyiron_atomistics.sphinx.output_parser import SphinxLogParser, splitter +from pyiron_atomistics.sphinx.output_parser import ( + SphinxLogParser, + collect_energy_dat, + collect_residue_dat, + collect_spins_dat, + collect_relaxed_hist, + collect_energy_struct, + collect_eps_dat, +) from pyiron_atomistics.sphinx.potential import ( find_potential_file as find_potential_file_jth, ) @@ -2133,7 +2141,7 @@ def __init__(self, job): self.generic.create_group("dft") self.old_version = False - def collect_spins_dat(self, file_name="spins.dat", cwd=None): + def collect_spins_dat(self, file_name="spins.dat", cwd="."): """ Args: @@ -2143,14 +2151,16 @@ def collect_spins_dat(self, file_name="spins.dat", cwd=None): Returns: """ - if not os.path.isfile(posixpath.join(cwd, file_name)): - return None - spins = np.loadtxt(posixpath.join(cwd, file_name)) - self.generic.dft.atom_scf_spins = splitter( - np.array([ss[self._job.id_spx_to_pyi] for ss in spins[:, 1:]]), spins[:, 0] - ) + try: + results = collect_spins_dat( + file_name=file_name, cwd=cwd, index_permutation=self._job.id_spx_to_pyi + ) + except FileNotFoundError: + return + for k, v in results.items(): + self.generic.dft[k] = v - def collect_energy_dat(self, file_name="energy.dat", cwd=None): + def collect_energy_dat(self, file_name="energy.dat", cwd="."): """ Args: @@ -2160,26 +2170,14 @@ def collect_energy_dat(self, file_name="energy.dat", cwd=None): Returns: """ - if not os.path.isfile(posixpath.join(cwd, file_name)): - return None - energies = np.loadtxt(posixpath.join(cwd, file_name)) - self.generic.dft.scf_computation_time = splitter(energies[:, 1], energies[:, 0]) - self.generic.dft.scf_energy_int = splitter( - energies[:, 2] * HARTREE_TO_EV, energies[:, 0] - ) - - def en_split(e, counter=energies[:, 0]): - return splitter(e * HARTREE_TO_EV, counter) - - if len(energies[0]) == 7: - self.generic.dft.scf_energy_free = en_split(energies[:, 3]) - self.generic.dft.scf_energy_zero = en_split(energies[:, 4]) - self.generic.dft.scf_energy_band = en_split(energies[:, 5]) - self.generic.dft.scf_electronic_entropy = en_split(energies[:, 6]) - else: - self.generic.dft.scf_energy_band = en_split(energies[:, 3]) + try: + results = collect_energy_dat(file_name=file_name, cwd=cwd) + except FileNotFoundError: + return + for k, v in results.items(): + self.generic.dft[k] = v - def collect_residue_dat(self, file_name="residue.dat", cwd=None): + def collect_residue_dat(self, file_name="residue.dat", cwd="."): """ Args: @@ -2189,17 +2187,14 @@ def collect_residue_dat(self, file_name="residue.dat", cwd=None): Returns: """ - if not os.path.isfile(posixpath.join(cwd, file_name)): - return None - residue = np.loadtxt(posixpath.join(cwd, file_name)) - if len(residue) == 0: - return None - if len(residue[0]) == 2: - self.generic.dft.scf_residue = splitter(residue[:, 1], residue[:, 0]) - else: - self.generic.dft.scf_residue = splitter(residue[:, 1:], residue[:, 0]) + try: + results = collect_residue_dat(file_name=file_name, cwd=cwd) + except FileNotFoundError: + return + for k, v in results.items(): + self.generic.dft[k] = v - def collect_eps_dat(self, file_name="eps.dat", cwd=None): + def collect_eps_dat(self, file_name=None, cwd="."): """ Args: @@ -2209,17 +2204,15 @@ def collect_eps_dat(self, file_name="eps.dat", cwd=None): Returns: """ - if isinstance(file_name, str): - file_name = [file_name] - values = [] - for f in file_name: - file_tmp = posixpath.join(cwd, f) - if not os.path.isfile(file_tmp): - return - values.append(np.loadtxt(file_tmp)[..., 1:]) - values = np.stack(values, axis=0) - if "bands_eigen_values" not in self.generic.dft.list_nodes(): - self.generic.dft.bands_eigen_values = values.reshape((-1,) + values.shape) + try: + results = collect_eps_dat( + file_name=file_name, cwd=cwd, spins=self._job._spin_enabled + ) + for k, v in results.items(): + if k not in self.generic.dft: + self.generic.dft[k] = v + except FileNotFoundError: + return def collect_energy_struct(self, file_name="energy-structOpt.dat", cwd=None): """ @@ -2231,11 +2224,12 @@ def collect_energy_struct(self, file_name="energy-structOpt.dat", cwd=None): Returns: """ - file_name = posixpath.join(cwd, file_name) - if os.path.isfile(file_name): - self.generic.dft.energy_free = ( - np.loadtxt(file_name).reshape(-1, 2)[:, 1] * HARTREE_TO_EV - ) + try: + results = collect_energy_struct(file_name=file_name, cwd=cwd) + except FileNotFoundError: + return + for k, v in results.items(): + self.generic.dft[k] = v def collect_sphinx_log(self, file_name="sphinx.log", cwd="."): """ @@ -2271,7 +2265,7 @@ def collect_sphinx_log(self, file_name="sphinx.log", cwd="."): self.generic.dft.scf_energy_int = self._spx_log_parser.get_energy_int() if "scf_energy_free" not in self.generic.dft.list_nodes(): self.generic.dft.scf_energy_free = self._spx_log_parser.get_energy_free() - if "forces" in self.generic.list_nodes(): + if "forces" not in self.generic.list_nodes(): self.generic.forces = self._spx_log_parser.get_forces( index_permutation=self._job.id_spx_to_pyi ) @@ -2292,41 +2286,14 @@ def collect_relaxed_hist(self, file_name="relaxHist.sx", cwd=None): Returns: """ - file_name = posixpath.join(cwd, file_name) - if not os.path.isfile(file_name): - return None - with open(file_name, "r") as f: - file_content = "".join(f.readlines()) - natoms = len(self._job.id_spx_to_pyi) - - def get_value(term, file_content=file_content, natoms=natoms): - value = ( - np.array( - [ - re.split("\[|\]", line)[1].split(",") - for line in re.findall(term, file_content, re.MULTILINE) - ] - ) - .astype(float) - .reshape(-1, natoms, 3) - ) - return np.array([ff[self._job.id_spx_to_pyi] for ff in value]) - - self.generic.positions = get_value("coords.*$") * BOHR_TO_ANGSTROM - self.generic.forces = ( - get_value("force.*$") * HARTREE_OVER_BOHR_TO_EV_OVER_ANGSTROM - ) - self.generic.cells = ( - np.array( - [ - json.loads(line) - for line in re.findall( - "cell =(.*?);", file_content.replace("\n", ""), re.MULTILINE - ) - ] + try: + results = collect_relaxed_hist( + file_name=file_name, cwd=cwd, index_permutation=self._job.id_spx_to_pyi ) - * BOHR_TO_ANGSTROM - ) + except FileNotFoundError: + return + for k, v in results.items(): + self.generic[k] = v def collect_charge_density(self, file_name, cwd): if ( @@ -2364,23 +2331,20 @@ def _get_electronic_structure_object(self): es.generate_from_matrices() return es - def collect(self, directory=os.getcwd()): + def collect(self, directory=None): """ The collect function, collects all the output from a SPHInX simulation. Args: directory (str): the directory to collect the output from. """ + if directory is None: + directory = self._job.working_directory self.collect_energy_struct(file_name="energy-structOpt.dat", cwd=directory) self.collect_sphinx_log(file_name="sphinx.log", cwd=directory) self.collect_energy_dat(file_name="energy.dat", cwd=directory) self.collect_residue_dat(file_name="residue.dat", cwd=directory) - if self._job._spin_enabled: - self.collect_eps_dat(file_name="eps.dat", cwd=directory) - else: - self.collect_eps_dat( - file_name=[f"eps.{i}.dat" for i in [0, 1]], cwd=directory - ) + self.collect_eps_dat(file_name=None, cwd=directory) self.collect_spins_dat(file_name="spins.dat", cwd=directory) self.collect_relaxed_hist(file_name="relaxHist.sx", cwd=directory) self.collect_electrostatic_potential(file_name="vElStat-eV.sxb", cwd=directory) diff --git a/pyiron_atomistics/sphinx/output_parser.py b/pyiron_atomistics/sphinx/output_parser.py index 646641abc..27d1f967a 100644 --- a/pyiron_atomistics/sphinx/output_parser.py +++ b/pyiron_atomistics/sphinx/output_parser.py @@ -23,6 +23,100 @@ def splitter(arr, counter): return arr_new +def collect_energy_dat(file_name="energy.dat", cwd="."): + """ + + Args: + file_name (str): file name + cwd (str): directory path + + Returns: + (dict): results + + """ + if cwd is None: + cwd = "." + energies = np.loadtxt(str(Path(cwd) / Path(file_name))) + results = {"scf_computation_time": splitter(energies[:, 1], energies[:, 0])} + results["scf_energy_int"] = splitter(energies[:, 2] * HARTREE_TO_EV, energies[:, 0]) + + def en_split(e, counter=energies[:, 0]): + return splitter(e * HARTREE_TO_EV, counter) + + if len(energies[0]) == 7: + results["scf_energy_free"] = en_split(energies[:, 3]) + results["scf_energy_zero"] = en_split(energies[:, 4]) + results["scf_energy_band"] = en_split(energies[:, 5]) + results["scf_electronic_entropy"] = en_split(energies[:, 6]) + else: + results["scf_energy_band"] = en_split(energies[:, 3]) + return results + + +def collect_residue_dat(file_name="residue.dat", cwd="."): + """ + + Args: + file_name (str): file name + cwd (str): directory path + + Returns: + (dict): results + + """ + if cwd is None: + cwd = "." + residue = np.loadtxt(str(Path(cwd) / Path(file_name))) + if len(residue) == 0: + return {} + return {"scf_residue": splitter(residue[:, 1:].squeeze(), residue[:, 0])} + + +def _collect_eps_dat(file_name="eps.dat", cwd="."): + """ + + Args: + file_name: + cwd: + + Returns: + + """ + if cwd is None: + cwd = "." + return np.loadtxt(str(Path(cwd) / Path(file_name)))[..., 1:] + + +def collect_eps_dat(file_name=None, cwd=".", spins=True): + if file_name is not None: + values = [_collect_eps_dat(file_name=file_name, cwd=cwd)] + elif spins: + values = [_collect_eps_dat(file_name=f"eps.{i}.dat", cwd=cwd) for i in [0, 1]] + else: + values = [_collect_eps_dat(file_name="eps.dat", cwd=cwd)] + values = np.stack(values, axis=0) + return {"bands_eigen_values": values.reshape((-1,) + values.shape)} + + +def collect_energy_struct(file_name="energy-structOpt.dat", cwd="."): + """ + + Args: + file_name (str): file name + cwd (str): directory path + + Returns: + (dict): results + + """ + if cwd is None: + cwd = "." + return { + "energy_free": np.loadtxt(str(Path(cwd) / Path(file_name))).reshape(-1, 2)[:, 1] + * HARTREE_TO_EV + } + + def check_permutation(index_permutation): if index_permutation is None: return @@ -33,6 +127,69 @@ def check_permutation(index_permutation): raise ValueError("missing entries in the index_permutation") +def collect_spins_dat(file_name="spins.dat", cwd=".", index_permutation=None): + """ + + Args: + file_name (str): file name + cwd (str): directory path + index_permutation (numpy.ndarray): Indices for the permutation + + Returns: + (dict): results + + """ + check_permutation(index_permutation) + if cwd is None: + cwd = "." + spins = np.loadtxt(str(Path(cwd) / Path(file_name))) + if index_permutation is not None: + s = np.array([ss[index_permutation] for ss in spins[:, 1:]]) + else: + s = spins[:, 1:] + return {"atom_scf_spins": splitter(s, spins[:, 0])} + + +def collect_relaxed_hist(file_name="relaxHist.sx", cwd=None, index_permutation=None): + """ + + Args: + file_name (str): file name + cwd (str): directory path + index_permutation (numpy.ndarray): Indices for the permutation + + Returns: + (dict): results + + # TODO: parse movable, elements, species etc. + """ + check_permutation(index_permutation) + if cwd is None: + cwd = "." + with open(file_name, "r") as f: + file_content = "".join(f.readlines()) + n_steps = len(re.findall("// --- step \d", file_content, re.MULTILINE)) + f_v = ",".join(3 * [r"\s*([\d.-]+)"]) + + def get_value(term, f=file_content, n=n_steps, p=index_permutation): + value = np.array(re.findall(p, f, re.MULTILINE)).astype(float).reshape(n, -1, 3) + if p is not None: + value = np.array([ff[p] for ff in value]) + return value + + cell = re.findall( + r"cell = \[\[" + r"\],\n\s*\[".join(3 * [f_v]) + r"\]\];", + file_content, + re.MULTILINE, + ) + cell = np.array(cell).astype(float).reshape(n_steps, 3, 3) * BOHR_TO_ANGSTROM + return { + "positions": get_value(r"atom {coords = \[" + f_v + r"\];") * BOHR_TO_ANGSTROM, + "forces": get_value(r"force = \[" + f_v + r"\]; }"), + "cell": cell, + } + + class SphinxLogParser: def __init__(self, file_name="sphinx.log", cwd="."): """ @@ -165,9 +322,7 @@ def get_forces(self, index_permutation=None): arr = np.array(re.findall(pattern, self.log_file)) if len(arr) == 0: return [] - indices = arr[:, 0].astype(int) - indices = indices.reshape(-1, max(indices) + 1) - forces = arr[:, 1:].astype(float).reshape(indices.shape + (3,)) + forces = arr[:, 1:].astype(float).reshape(-1, self.n_atoms, 3) forces *= HARTREE_OVER_BOHR_TO_EV_OVER_ANGSTROM if index_permutation is not None: for ii, ff in enumerate(forces):