From dd1dbe81bca99cf1d3e9d506a8f685be6bf7f1ba Mon Sep 17 00:00:00 2001 From: Atsushi Togo Date: Mon, 29 Apr 2024 11:52:46 +0900 Subject: [PATCH] Supercell energies in type-2 dataset --- phono3py/api_phono3py.py | 16 ++++-- phono3py/cui/create_force_sets.py | 83 +++++++++++++++++++++-------- phono3py/interface/phono3py_yaml.py | 60 ++++++++++----------- 3 files changed, 101 insertions(+), 58 deletions(-) diff --git a/phono3py/api_phono3py.py b/phono3py/api_phono3py.py index 8de413d4..1cb98b7f 100644 --- a/phono3py/api_phono3py.py +++ b/phono3py/api_phono3py.py @@ -582,6 +582,7 @@ def dataset(self): {'number': atom index of second displaced atom, 'displacement': displacement in Cartesian coordinates}, 'forces': forces on atoms in supercell, + 'supercell_energy': energy of supercell, 'pair_distance': distance between paired atoms, 'included': with cutoff pair distance in displacement pair generation, this indicates if this @@ -591,9 +592,11 @@ def dataset(self): ... ] }, ... ] } Type 2. All atomic displacements in each supercell: {'displacements': ndarray, dtype='double', order='C', - shape=(supercells, atoms in supercell, 3) + shape=(supercells, atoms in supercell, 3), 'forces': ndarray, dtype='double',, order='C', - shape=(supercells, atoms in supercell, 3)} + shape=(supercells, atoms in supercell, 3), + 'supercell_energies': ndarray, dtype='double', + shape=(supercells,)} In type 2, displacements and forces can be given by numpy array with different shape but that can be reshaped to (supercells, natom, 3). @@ -622,12 +625,15 @@ def phonon_dataset(self): 'first_atoms': [ {'number': atom index of first displaced atom, 'displacement': displacement in Cartesian coordinates, - 'forces': forces on atoms in supercell} ... ]} + 'forces': forces on atoms in supercell, + 'supercell_energy': energy of supercell}, ... ]} Type 2. All atomic displacements in each supercell: {'displacements': ndarray, dtype='double', order='C', - shape=(supercells, atoms in supercell, 3) + shape=(supercells, atoms in supercell, 3), 'forces': ndarray, dtype='double',, order='C', - shape=(supercells, atoms in supercell, 3)} + shape=(supercells, atoms in supercell, 3), + 'supercell_energies': ndarray, dtype='double', + shape=(supercells,)} In type 2, displacements and forces can be given by numpy array with different shape but that can be reshaped to (supercells, natom, 3). diff --git a/phono3py/cui/create_force_sets.py b/phono3py/cui/create_force_sets.py index 5173d4b2..8916aa17 100644 --- a/phono3py/cui/create_force_sets.py +++ b/phono3py/cui/create_force_sets.py @@ -40,6 +40,7 @@ import sys from typing import Optional +import numpy as np from phonopy.cui.create_force_sets import check_number_of_force_files from phonopy.cui.load_helper import get_nac_params from phonopy.cui.phonopy_script import file_exists, files_exist, print_error @@ -227,8 +228,16 @@ def _get_force_sets_fc2( interface_mode = settings.calculator if log_level: print(f'FC2 displacement dataset was read from "{disp_filename}".') - num_atoms = disp_dataset["natom"] - num_disps = len(disp_dataset["first_atoms"]) + + if "first_atoms" in disp_dataset: # type-1 + num_atoms = disp_dataset["natom"] + num_disps = len(disp_dataset["first_atoms"]) + elif "displacements" in disp_dataset: # type-2: + num_disps = len(disp_dataset["displacements"]) + num_atoms = len(disp_dataset["displacements"][0]) + else: + raise RuntimeError("FC2 displacement dataset is broken.") + force_filenames = settings.create_forces_fc2 for filename in force_filenames: file_exists(filename, log_level) @@ -276,12 +285,18 @@ def _get_force_sets_fc3( print("") print(f'FC3 Displacement dataset was read from "{disp_filename}".') - num_atoms = disp_dataset["natom"] - num_disps = len(disp_dataset["first_atoms"]) - for d1 in disp_dataset["first_atoms"]: - for d2 in d1["second_atoms"]: - if "included" not in d2 or d2["included"]: - num_disps += 1 + if "first_atoms" in disp_dataset: # type-1 + num_atoms = disp_dataset["natom"] + num_disps = len(disp_dataset["first_atoms"]) + for d1 in disp_dataset["first_atoms"]: + for d2 in d1["second_atoms"]: + if "included" not in d2 or d2["included"]: + num_disps += 1 + elif "displacements" in disp_dataset: # type-2: + num_disps = len(disp_dataset["displacements"]) + num_atoms = len(disp_dataset["displacements"][0]) + else: + raise RuntimeError("FC3 displacement dataset is broken.") if settings.create_forces_fc3_file: file_exists(settings.create_forces_fc3_file, log_level) @@ -336,23 +351,45 @@ def _get_force_sets_fc3( def _set_forces_and_nac_params( ph3py_yaml: Phono3pyYaml, settings, calc_dataset_fc3: dict, calc_dataset_fc2: dict ): - count = len(ph3py_yaml.dataset["first_atoms"]) - for i, d1 in enumerate(ph3py_yaml.dataset["first_atoms"]): - d1["forces"] = calc_dataset_fc3["forces"][i] - if "energies" in calc_dataset_fc3: - d1["energy"] = float(calc_dataset_fc3["energies"][i]) - for d2 in d1["second_atoms"]: - if "included" not in d2 or d2["included"]: - d2["forces"] = calc_dataset_fc3["forces"][count] - if "energies" in calc_dataset_fc3: - d2["energy"] = float(calc_dataset_fc3["energies"][count]) - count += 1 + if "first_atoms" in ph3py_yaml.dataset: + count = len(ph3py_yaml.dataset["first_atoms"]) + for i, d1 in enumerate(ph3py_yaml.dataset["first_atoms"]): + d1["forces"] = calc_dataset_fc3["forces"][i] + if "supercell_energies" in calc_dataset_fc3: + d1["supercell_energy"] = float( + calc_dataset_fc3["supercell_energies"][i] + ) + for d2 in d1["second_atoms"]: + if "included" not in d2 or d2["included"]: + d2["forces"] = calc_dataset_fc3["forces"][count] + if "supercell_energies" in calc_dataset_fc3: + d2["supercell_energy"] = float( + calc_dataset_fc3["supercell_energies"][count] + ) + count += 1 + else: + ph3py_yaml.dataset["forces"] = np.array( + calc_dataset_fc3["forces"], dtype="double", order="C" + ) + ph3py_yaml.dataset["supercell_energies"] = np.array( + calc_dataset_fc3["supercell_energies"], dtype="double" + ) if settings.create_forces_fc2: - for i, d in enumerate(ph3py_yaml.phonon_dataset["first_atoms"]): - d["forces"] = calc_dataset_fc2["forces"][i] - if "energies" in calc_dataset_fc2: - d["energy"] = float(calc_dataset_fc2["energies"][i]) + if "first_atoms" in ph3py_yaml.phonon_dataset: + for i, d in enumerate(ph3py_yaml.phonon_dataset["first_atoms"]): + d["forces"] = calc_dataset_fc2["forces"][i] + if "supercell_energies" in calc_dataset_fc2: + d["supercell_energy"] = float( + calc_dataset_fc2["supercell_energies"][i] + ) + else: + ph3py_yaml.phonon_dataset["forces"] = np.array( + calc_dataset_fc2["forces"], dtype="double", order="C" + ) + ph3py_yaml.phonon_dataset["supercell_energies"] = np.array( + calc_dataset_fc2["supercell_energies"], dtype="double" + ) nac_params = get_nac_params(primitive=ph3py_yaml.primitive) if nac_params: diff --git a/phono3py/interface/phono3py_yaml.py b/phono3py/interface/phono3py_yaml.py index 81d972f0..dc86c5d2 100644 --- a/phono3py/interface/phono3py_yaml.py +++ b/phono3py/interface/phono3py_yaml.py @@ -133,14 +133,11 @@ def _parse_dataset(self): This method override PhonopyYaml._parse_dataset. - key="phonon_displacements" at v2.2 or later. - key="displacements" at older version than v2.2. - """ self._data.phonon_dataset = self._get_dataset( - self._data.phonon_supercell, key="phonon_displacements" + self._data.phonon_supercell, key_prefix="phonon_" ) - if self._data.phonon_dataset is None: # key="displacements" + if self._data.phonon_dataset is None: self._data.phonon_dataset = self._get_dataset(self._data.phonon_supercell) def _parse_fc3_dataset(self): @@ -192,8 +189,8 @@ def _parse_fc3_dataset_type1(self, natom): data1["id"] = d1_id if "forces" in d1: data1["forces"] = np.array(d1["forces"], dtype="double", order="C") - if "energy" in d1: - data1["energy"] = d1["energy"] + if "supercell_energy" in d1: + data1["supercell_energy"] = d1["supercell_energy"] d2_list = d1.get("paired_with") if d2_list is None: # backward compatibility d2_list = d1.get("second_atoms") @@ -228,8 +225,8 @@ def _parse_fc3_dataset_type1_with_forces(self, data1, d2, disp2_id): second_atom_dict["id"] = d2_id if "pair_distance" in d2: second_atom_dict["pair_distance"] = d2["pair_distance"] - if "energy" in d2: - second_atom_dict["energy"] = d2["energy"] + if "supercell_energy" in d2: + second_atom_dict["supercell_energy"] = d2["supercell_energy"] disp2_id += 1 data1["second_atoms"].append(second_atom_dict) return disp2_id @@ -298,15 +295,6 @@ def _cell_info_yaml_lines(self): lines += self._phonon_supercell_yaml_lines() return lines - def _phonon_supercell_matrix_yaml_lines(self): - lines = [] - if self._data.phonon_supercell_matrix is not None: - lines.append("phonon_supercell_matrix:") - for v in self._data.supercell_matrix: - lines.append("- [ %3d, %3d, %3d ]" % tuple(v)) - lines.append("") - return lines - def _phonon_supercell_yaml_lines(self): lines = [] if self._data.phonon_supercell is not None: @@ -330,7 +318,7 @@ def _nac_yaml_lines(self): else: return self._nac_yaml_lines_given_symbols(self._data.primitive.symbols) - def _displacements_yaml_lines(self, with_forces=False, key="displacements"): + def _displacements_yaml_lines(self, with_forces: bool = False) -> list: """Get YAML lines for phonon_dataset and dataset. This method override PhonopyYaml._displacements_yaml_lines. @@ -341,16 +329,18 @@ def _displacements_yaml_lines(self, with_forces=False, key="displacements"): lines = [] if self._data.phonon_supercell_matrix is not None: lines += self._displacements_yaml_lines_2types( - self._data.phonon_dataset, with_forces=with_forces, key=f"phonon_{key}" + self._data.phonon_dataset, + with_forces=with_forces, + key_prefix="phonon_", ) lines += self._displacements_yaml_lines_2types( - self._data.dataset, with_forces=with_forces, key=key + self._data.dataset, with_forces=with_forces ) return lines def _displacements_yaml_lines_type1( - self, dataset, with_forces=False, key="displacements" - ): + self, dataset: dict, with_forces: bool = False, key_prefix: str = "" + ) -> list: """Get YAML lines for type1 phonon_dataset and dataset. This method override PhonopyYaml._displacements_yaml_lines_type1. @@ -358,7 +348,9 @@ def _displacements_yaml_lines_type1( Phono3pyYaml._displacements_yaml_lines_type1. """ - return displacements_yaml_lines_type1(dataset, with_forces=with_forces, key=key) + return displacements_yaml_lines_type1( + dataset, with_forces=with_forces, key_prefix=key_prefix + ) class Phono3pyYaml(PhonopyYaml): @@ -445,7 +437,9 @@ def set_phonon_info(self, phono3py: "Phono3py"): self._data.phonon_supercell = phono3py.phonon_supercell -def displacements_yaml_lines_type1(dataset, with_forces=False, key="displacements"): +def displacements_yaml_lines_type1( + dataset: dict, with_forces: bool = False, key_prefix: str = "" +) -> list: """Get YAML lines for type1 phonon_dataset and dataset. This is a function but not class method because used by other function. @@ -456,7 +450,7 @@ def displacements_yaml_lines_type1(dataset, with_forces=False, key="displacement if "second_atoms" in dataset["first_atoms"][0]: lines = ["displacement_pairs:"] else: - lines = [f"{key}:"] + lines = [f"{key_prefix}displacements:"] for disp1_id, d1 in enumerate(dataset["first_atoms"]): lines.append("- atom: %4d" % (d1["number"] + 1)) lines.append(" displacement:") @@ -468,8 +462,12 @@ def displacements_yaml_lines_type1(dataset, with_forces=False, key="displacement lines.append(" forces:") for v in d1["forces"]: lines.append(" - [ %19.16f, %19.16f, %19.16f ]" % tuple(v)) - if "energy" in d1: - lines.append(" energy: {energy:.8f}".format(energy=d1["energy"])) + if "supercell_energy" in d1: + lines.append( + " supercell_energy: {energy:.8f}".format( + energy=d1["supercell_energy"] + ) + ) if "second_atoms" in d1: ret_lines, id_offset = _second_displacements_yaml_lines( d1["second_atoms"], id_offset, with_forces=with_forces @@ -548,9 +546,11 @@ def _second_displacements_yaml_lines(dataset2, id_offset, with_forces=False): lines.append(" forces:") for v in dataset2[j]["forces"]: lines.append(" - [ %19.16f, %19.16f, %19.16f ]" % tuple(v)) - if "energy" in dataset2[j]: + if "supercell_energy" in dataset2[j]: lines.append( - " energy: {energy:.8f}".format(energy=dataset2[j]["energy"]) + " supercell_energy: {energy:.8f}".format( + energy=dataset2[j]["supercell_energy"] + ) ) else: lines.append(" - atom: %4d" % (i + 1))