diff --git a/phono3py/api_phono3py.py b/phono3py/api_phono3py.py index e758b027..2b14b669 100644 --- a/phono3py/api_phono3py.py +++ b/phono3py/api_phono3py.py @@ -40,7 +40,6 @@ from collections.abc import Sequence from typing import Literal, Optional, Union import numpy as np -from phonopy import Phonopy from phonopy.exception import ForceCalculatorRequiredError from phonopy.harmonic.displacement import ( directions_to_displacement_dataset, @@ -56,7 +55,13 @@ from phonopy.harmonic.force_constants import ( symmetrize_force_constants, ) from phonopy.interface.fc_calculator import get_fc2 -from phonopy.interface.pypolymlp import PypolymlpParams +from phonopy.interface.pypolymlp import ( + PypolymlpData, + PypolymlpParams, + develop_polymlp, + evalulate_polymlp, + parse_mlp_params, +) from phonopy.structure.atoms import PhonopyAtoms from phonopy.structure.cells import ( Primitive, @@ -295,8 +300,9 @@ class Phono3py: self._fc2 = None self._fc3 = None - # MLP container - self._mlp_phonopy = None + # MLP + self._mlp = None + self._mlp_dataset = None # Setup interaction self._interaction = None @@ -681,6 +687,35 @@ class Phono3py: self._phonon_supercells_with_displacements = None + @property + def mlp_dataset(self) -> Optional[dict]: + """Return displacement-force dataset. + + The supercell matrix is equal to that of usual displacement-force + dataset. Only type 2 format is supported. "displacements", + "forces", and "supercell_energies" should be contained. + + """ + return self._mlp_dataset + + @mlp_dataset.setter + def mlp_dataset(self, mlp_dataset: dict): + if not isinstance(mlp_dataset, dict): + raise TypeError("mlp_dataset has to be a dictionary.") + if "displacements" not in mlp_dataset: + raise RuntimeError("Displacements have to be given.") + if "forces" not in mlp_dataset: + raise RuntimeError("Forces have to be given.") + if "supercell_energy" in mlp_dataset: + raise RuntimeError("Supercell energies have to be given.") + if len(mlp_dataset["displacements"]) != len(mlp_dataset["forces"]): + raise RuntimeError("Length of displacements and forces are different.") + if len(mlp_dataset["displacements"]) != len(mlp_dataset["supercell_energies"]): + raise RuntimeError( + "Length of displacements and supercell_energies are different." + ) + self._mlp_dataset = mlp_dataset + @property def band_indices(self) -> list[np.ndarray]: """Setter and getter of band indices. @@ -811,7 +846,7 @@ class Phono3py: for disp1 in dataset["first_atoms"]: num_scells += len(disp1["second_atoms"]) displacements = np.zeros( - (num_scells, self._supercell.get_number_of_atoms(), 3), + (num_scells, len(self._supercell), 3), dtype="double", order="C", ) @@ -990,30 +1025,6 @@ class Phono3py: def phonon_supercell_energies(self, values): self._set_phonon_forces_energies(values, target="supercell_energies") - @property - def mlp_dataset(self) -> Optional[dict]: - """Return displacement-force dataset. - - The supercell matrix is equal to that of usual displacement-force - dataset. Only type 2 format is supported. "displacements", - "forces", and "supercell_energies" should be contained. - - """ - if self._mlp_phonopy is None: - return None - return self._mlp_phonopy.mlp_dataset - - @mlp_dataset.setter - def mlp_dataset(self, mlp_dataset: dict): - if self._mlp_phonopy is None: - self._mlp_phonopy = Phonopy( - self._unitcell, - supercell_matrix=self._supercell_matrix, - symprec=self._symprec, - log_level=self._log_level, - ) - self._mlp_phonopy.mlp_dataset = mlp_dataset - @property def phph_interaction(self): """Return Interaction instance.""" @@ -1385,53 +1396,31 @@ class Phono3py: symmetrize_fc3r: bool = False, is_compact_fc: bool = False, fc_calculator: Optional[str] = None, - fc_calculator_options: Optional[str] = None, - use_pypolymlp: bool = False, - mlp_params: Optional[Union[PypolymlpParams, dict]] = None, + fc_calculator_options: Optional[Union[str, dict]] = None, ): """Calculate fc3 from displacements and forces. Parameters ---------- - symmetrize_fc3r : bool + symmetrize_fc3r : bool, optional Only for type 1 displacement_dataset, translational and permutation symmetries are applied after creating fc3. This symmetrization is not very sophisticated and can break space group symmetry, but often useful. If better symmetrization is expected, it is recommended to use external force constants calculator such as ALM. Default is False. - is_compact_fc : bool + is_compact_fc : bool, optional fc3 shape is False: (supercell, supercell, supecell, 3, 3, 3) True: (primitive, supercell, supecell, 3, 3, 3) where 'supercell' and 'primitive' indicate number of atoms in these cells. Default is False. - fc_calculator : str or None + fc_calculator : str, optional Force constants calculator given by str. - fc_calculator_options : dict + fc_calculator_options : dict or str, optional Options for external force constants calculator. - use_pypolymlp : bool - Use MLP of pypolymlp to calculate fc3. Default is False. mlp_dataset - has to be set before calling this method. - mlp_params : PypolymlpParams or dict, optional - Parameters for developing MLP. Default is None. When dict is given, - PypolymlpParams instance is created from the dict. """ - if use_pypolymlp: - if self._mlp_phonopy is None: - msg = "mlp_dataset has to be set before calling this method." - raise RuntimeError(msg) - if self._dataset is None or "displacements" not in self._dataset: - raise RuntimeError( - "Type 2 displacements are not set. Run generate_displacements." - ) - self._mlp_phonopy.develop_mlp(params=mlp_params) - self._mlp_phonopy.displacements = self._dataset["displacements"] - self._mlp_phonopy.evaluate_mlp() - self._dataset["forces"] = self._mlp_phonopy.forces - self._dataset["supercell_energies"] = self._mlp_phonopy.supercell_energies - fc3_calculator, fc3_calculator_options = self._extract_fc2_fc3_calculators( fc_calculator, fc_calculator_options, 3 ) @@ -2172,6 +2161,65 @@ class Phono3py: with open(filename, "w") as w: w.write(str(ph3py_yaml)) + def develop_mlp(self, params: Optional[Union[PypolymlpParams, dict, str]] = None): + """Develop MLP of pypolymlp. + + Parameters + ---------- + params : PypolymlpParams or dict, optional + Parameters for developing MLP. Default is None. When dict is given, + PypolymlpParams instance is created from the dict. + + """ + if self._mlp_dataset is None: + raise RuntimeError("MLP dataset is not set.") + + if params is not None: + _params = parse_mlp_params(params) + + disps = self._mlp_dataset["displacements"] + forces = self._mlp_dataset["forces"] + energies = self._mlp_dataset["supercell_energies"] + n = int(len(disps) * 0.9) + train_data = PypolymlpData( + displacements=disps[:n], forces=forces[:n], supercell_energies=energies[:n] + ) + test_data = PypolymlpData( + displacements=disps[n:], forces=forces[n:], supercell_energies=energies[n:] + ) + self._mlp = develop_polymlp( + self._supercell, + train_data, + test_data, + params=_params, + verbose=self._log_level - 1 > 0, + ) + + def evaluate_mlp(self): + """Evaluate the machine learning potential of pypolymlp. + + This method calculates the supercell energies and forces from the MLP + for the displacements in self._dataset of type 2. The results are stored + in self._dataset. + + The displacements may be generated by the produce_force_constants method + with number_of_snapshots > 0. With MLP, a small distance parameter, such + as 0.001, can be numerically stable for the computation of force + constants. + + """ + if self._mlp is None: + raise RuntimeError("MLP is not developed yet.") + + if self.supercells_with_displacements is None: + raise RuntimeError("Displacements are not set. Run generate_displacements.") + + energies, forces, _ = evalulate_polymlp( + self._mlp, self.supercells_with_displacements + ) + self.supercell_energies = energies + self.forces = forces + ################### # private methods # ################### @@ -2252,32 +2300,44 @@ class Phono3py: raise RuntimeError def _build_phonon_supercells_with_displacements( - self, supercell: PhonopyAtoms, displacement_dataset + self, supercell: PhonopyAtoms, dataset ): supercells = [] + positions = supercell.positions magmoms = supercell.magnetic_moments masses = supercell.masses numbers = supercell.numbers lattice = supercell.cell - for disp1 in displacement_dataset["first_atoms"]: - disp_cart1 = disp1["displacement"] - positions = supercell.positions - positions[disp1["number"]] += disp_cart1 - supercells.append( - PhonopyAtoms( - numbers=numbers, - masses=masses, - magnetic_moments=magmoms, - positions=positions, - cell=lattice, + if "displacements" in dataset: + for disp in dataset["displacements"]: + supercells.append( + PhonopyAtoms( + numbers=numbers, + masses=masses, + magnetic_moments=magmoms, + positions=positions + disp, + cell=lattice, + ) + ) + else: + for disp1 in dataset["first_atoms"]: + disp_cart1 = disp1["displacement"] + positions = supercell.positions + positions[disp1["number"]] += disp_cart1 + supercells.append( + PhonopyAtoms( + numbers=numbers, + masses=masses, + magnetic_moments=magmoms, + positions=positions, + cell=lattice, + ) ) - ) return supercells def _build_supercells_with_displacements(self): - supercells = [] magmoms = self._supercell.magnetic_moments masses = self._supercell.masses numbers = self._supercell.numbers @@ -2287,28 +2347,29 @@ class Phono3py: self._supercell, self._dataset ) - for disp1 in self._dataset["first_atoms"]: - disp_cart1 = disp1["displacement"] - for disp2 in disp1["second_atoms"]: - if "included" in disp2: - included = disp2["included"] - else: - included = True - if included: - positions = self._supercell.positions - positions[disp1["number"]] += disp_cart1 - positions[disp2["number"]] += disp2["displacement"] - supercells.append( - PhonopyAtoms( - numbers=numbers, - masses=masses, - magnetic_moments=magmoms, - positions=positions, - cell=lattice, + if "first_atoms" in self._dataset: + for disp1 in self._dataset["first_atoms"]: + disp_cart1 = disp1["displacement"] + for disp2 in disp1["second_atoms"]: + if "included" in disp2: + included = disp2["included"] + else: + included = True + if included: + positions = self._supercell.positions + positions[disp1["number"]] += disp_cart1 + positions[disp2["number"]] += disp2["displacement"] + supercells.append( + PhonopyAtoms( + numbers=numbers, + masses=masses, + magnetic_moments=magmoms, + positions=positions, + cell=lattice, + ) ) - ) - else: - supercells.append(None) + else: + supercells.append(None) self._supercells_with_displacements = supercells diff --git a/phono3py/cui/create_force_constants.py b/phono3py/cui/create_force_constants.py index d3bae975..64e75bd8 100644 --- a/phono3py/cui/create_force_constants.py +++ b/phono3py/cui/create_force_constants.py @@ -40,7 +40,8 @@ import copy import os import pathlib import sys -from typing import Optional +from dataclasses import asdict +from typing import Optional, Union import numpy as np from phonopy.cui.phonopy_script import file_exists, print_error @@ -52,6 +53,7 @@ from phonopy.harmonic.force_constants import ( ) from phonopy.interface.calculator import get_default_physical_units from phonopy.interface.fc_calculator import fc_calculator_names +from phonopy.interface.pypolymlp import PypolymlpParams, parse_mlp_params from phono3py import Phono3py from phono3py.cui.show_log import show_phono3py_force_constants_settings @@ -122,6 +124,11 @@ def create_phono3py_force_constants( settings.cutoff_pair_distance, fc_calculator, fc_calculator_options, + settings.use_pypolymlp, + settings.mlp_params, + settings.displacement_distance, + settings.random_displacements, + settings.random_seed, log_level, ) @@ -249,7 +256,7 @@ def parse_forces( # Try to read FORCES_FC* if type-2 and return dataset. # None is returned unless type-2. # can emit FileNotFoundError. - if dataset is None or dataset is not None and not forces_in_dataset(dataset): + if dataset is None or (dataset is not None and not forces_in_dataset(dataset)): _dataset = read_type2_dataset( natom, filename=force_filename, log_level=log_level ) @@ -321,8 +328,10 @@ def forces_in_dataset(dataset: dict) -> bool: ) -def displacements_in_dataset(dataset: dict) -> bool: +def displacements_in_dataset(dataset: Optional[dict]) -> bool: """Return whether displacements in dataset or not.""" + if dataset is None: + return False return "displacements" in dataset or "first_atoms" in dataset @@ -443,6 +452,11 @@ def _create_phono3py_fc3( cutoff_pair_distance: Optional[float], fc_calculator: Optional[str], fc_calculator_options: Optional[str], + use_pypoymplp: bool, + mlp_params: Union[str, dict, PypolymlpParams], + displacement_distance: Optional[float], + number_of_snapshots: Optional[int], + random_seed: Optional[int], log_level: int, ): """Read or calculate fc3. @@ -481,7 +495,18 @@ def _create_phono3py_fc3( # from _get_type2_dataset file_exists(e.filename, log_level) - phono3py.dataset = dataset + if use_pypoymplp: + phono3py.mlp_dataset = dataset + run_pypolymlp_to_compute_forces( + phono3py, + mlp_params=mlp_params, + displacement_distance=displacement_distance, + number_of_snapshots=number_of_snapshots, + random_seed=random_seed, + log_level=log_level, + ) + else: + phono3py.dataset = dataset phono3py.produce_fc3( symmetrize_fc3r=symmetrize_fc3r, is_compact_fc=is_compact_fc, @@ -490,6 +515,68 @@ def _create_phono3py_fc3( ) +def run_pypolymlp_to_compute_forces( + ph3py: Phono3py, + mlp_params: Union[str, dict, PypolymlpParams], + displacement_distance: float, + number_of_snapshots: int, + random_seed: int, + log_level: int, +): + """Run pypolymlp to compute forces.""" + if displacement_distance is None: + _displacement_distance = 0.001 + else: + _displacement_distance = displacement_distance + + if log_level: + print("-" * 29 + " pypolymlp start " + "-" * 30) + print("MLP parameters:") + if mlp_params: + for k, v in asdict(parse_mlp_params(mlp_params)).items(): + if v is not None: + print(f" {k}: {v}") + + if log_level: + if number_of_snapshots: + print("Generate random displacements") + print( + " Twice of number of snapshots will be generated " + "for plus-minus displacements." + ) + else: + print("Generate displacements") + print( + f" displacement distance: {_displacement_distance:.5f}".rstrip("0").rstrip( + "." + ) + ) + ph3py.generate_displacements( + distance=_displacement_distance, + is_plusminus=True, + number_of_snapshots=number_of_snapshots, + random_seed=random_seed, + ) + + if log_level: + print( + " Number of supercells for computing forces: " + f"{ph3py.displacements.shape[0]}", + flush=True, + ) + + if ph3py.mlp_dataset is None: + msg = "mlp_dataset has to be set before calling this method." + raise RuntimeError(msg) + if ph3py.supercells_with_displacements is None: + raise RuntimeError("Displacements are not set. Run generate_displacements.") + ph3py.develop_mlp(params=mlp_params) + ph3py.evaluate_mlp() + + if log_level: + print("-" * 30 + " pypolymlp end " + "-" * 31, flush=True) + + def _create_phono3py_fc2( phono3py: Phono3py, ph3py_yaml: Optional[Phono3pyYaml], diff --git a/phono3py/cui/load.py b/phono3py/cui/load.py index e51ec844..a768a6c6 100644 --- a/phono3py/cui/load.py +++ b/phono3py/cui/load.py @@ -49,10 +49,9 @@ from phonopy.structure.cells import determinant from phono3py import Phono3py from phono3py.cui.create_force_constants import ( - displacements_in_dataset, forces_in_dataset, parse_forces, - read_type2_dataset, + run_pypolymlp_to_compute_forces, ) from phono3py.file_IO import read_fc2_from_hdf5, read_fc3_from_hdf5 from phono3py.interface.phono3py_yaml import Phono3pyYaml @@ -86,6 +85,8 @@ def load( symmetrize_fc: bool = True, is_mesh_symmetry: bool = True, is_compact_fc: bool = False, + use_pypolymlp: bool = False, + mlp_params: Optional[dict] = None, use_grg: bool = False, make_r0_average: bool = True, symprec: float = 1e-5, @@ -236,6 +237,10 @@ def load( True: (primitive, supecell, 3, 3) False: (supercell, supecell, 3, 3) where 'supercell' and 'primitive' indicate number of atoms in these cells. Default is False. + use_pypolymlp : bool, optional + Use pypolymlp for generating force constants. Default is False. + mlp_params : dict, optional + A set of parameters used by machine learning potentials. use_grg : bool, optional Use generalized regular grid when True. Default is False. make_r0_average : bool, optional @@ -337,6 +342,7 @@ def load( forces_fc3_filename=forces_fc3_filename, forces_fc2_filename=forces_fc2_filename, phono3py_yaml_filename=phono3py_yaml, + use_pypolymlp=use_pypolymlp, log_level=log_level, ) @@ -348,6 +354,8 @@ def load( fc_calculator_options=fc_calculator_options, symmetrize_fc=symmetrize_fc, is_compact_fc=is_compact_fc, + use_pypolymlp=use_pypolymlp, + mlp_params=mlp_params, log_level=log_level, ) @@ -370,6 +378,7 @@ def set_dataset_and_force_constants( forces_fc2_filename: Optional[Union[os.PathLike, Sequence]] = None, phono3py_yaml_filename: Optional[os.PathLike] = None, cutoff_pair_distance: Optional[float] = None, + use_pypolymlp: bool = False, log_level: int = 0, ) -> dict: """Set displacements, forces, and create force constants. @@ -387,7 +396,7 @@ def set_dataset_and_force_constants( """ read_fc = {"fc2": False, "fc3": False} - read_fc["fc3"] = _set_dataset_or_fc3( + read_fc["fc3"], dataset = _get_dataset_or_fc3( ph3py, ph3py_yaml=ph3py_yaml, fc3_filename=fc3_filename, @@ -396,20 +405,20 @@ def set_dataset_and_force_constants( cutoff_pair_distance=cutoff_pair_distance, log_level=log_level, ) - read_fc["fc2"] = _set_dataset_phonon_dataset_or_fc2( + if not read_fc["fc3"]: + if use_pypolymlp: + ph3py.mlp_dataset = dataset + else: + ph3py.dataset = dataset + read_fc["fc2"], phonon_dataset = _get_dataset_phonon_dataset_or_fc2( ph3py, ph3py_yaml=ph3py_yaml, fc2_filename=fc2_filename, forces_fc2_filename=forces_fc2_filename, log_level=log_level, ) - - # Cases that dataset is in phono3py.yaml but not forces. - if ph3py.dataset is None: - if ph3py_yaml is not None and ph3py_yaml.dataset is not None: - ph3py.dataset = ph3py_yaml.dataset - if ph3py_yaml is not None and ph3py_yaml.phonon_dataset is not None: - ph3py.phonon_dataset = ph3py_yaml.phonon_dataset + if not read_fc["fc2"]: + ph3py.phonon_dataset = phonon_dataset return read_fc @@ -418,9 +427,14 @@ def compute_force_constants_from_datasets( ph3py: Phono3py, read_fc: dict, fc_calculator: Optional[str] = None, - fc_calculator_options: Optional[dict] = None, + fc_calculator_options: Optional[Union[dict, str]] = None, symmetrize_fc: bool = True, is_compact_fc: bool = True, + use_pypolymlp: bool = False, + mlp_params: Optional[Union[dict, str]] = None, + displacement_distance: Optional[float] = None, + number_of_snapshots: Optional[int] = None, + random_seed: Optional[int] = None, log_level: int = 0, ): """Compute force constants from datasets. @@ -435,15 +449,27 @@ def compute_force_constants_from_datasets( fc2 : bool """ - if not read_fc["fc3"] and ph3py.dataset: - ph3py.produce_fc3( - symmetrize_fc3r=symmetrize_fc, - is_compact_fc=is_compact_fc, - fc_calculator=fc_calculator, - fc_calculator_options=fc_calculator_options, - ) - if log_level and symmetrize_fc and fc_calculator is None: - print("fc3 was symmetrized.") + if not read_fc["fc3"]: + if ph3py.dataset or ph3py.mlp_dataset: + if ph3py.mlp_dataset and use_pypolymlp: + run_pypolymlp_to_compute_forces( + ph3py, + mlp_params=mlp_params, + displacement_distance=displacement_distance, + number_of_snapshots=number_of_snapshots, + random_seed=random_seed, + log_level=log_level, + ) + + ph3py.produce_fc3( + symmetrize_fc3r=symmetrize_fc, + is_compact_fc=is_compact_fc, + fc_calculator=fc_calculator, + fc_calculator_options=fc_calculator_options, + ) + + if log_level and symmetrize_fc and fc_calculator is None: + print("fc3 was symmetrized.") if not read_fc["fc2"] and (ph3py.dataset or ph3py.phonon_dataset): if ( @@ -467,7 +493,7 @@ def compute_force_constants_from_datasets( print("fc2 was symmetrized.") -def _set_dataset_or_fc3( +def _get_dataset_or_fc3( ph3py: Phono3py, ph3py_yaml: Optional[Phono3pyYaml] = None, fc3_filename: Optional[os.PathLike] = None, @@ -475,19 +501,40 @@ def _set_dataset_or_fc3( phono3py_yaml_filename: Optional[os.PathLike] = None, cutoff_pair_distance: Optional[float] = None, log_level: int = 0, -) -> bool: +) -> tuple[bool, dict]: p2s_map = ph3py.primitive.p2s_map read_fc3 = False - if fc3_filename is not None: - fc3 = read_fc3_from_hdf5(filename=fc3_filename, p2s_map=p2s_map) - _check_fc3_shape(ph3py, fc3, filename=fc3_filename) + dataset = None + if fc3_filename is not None or pathlib.Path("fc3.hdf5").exists(): + if fc3_filename is None: + _fc3_filename = "fc3.hdf5" + else: + _fc3_filename = fc3_filename + fc3 = read_fc3_from_hdf5(filename=_fc3_filename, p2s_map=p2s_map) + _check_fc3_shape(ph3py, fc3, filename=_fc3_filename) ph3py.fc3 = fc3 read_fc3 = True if log_level: - print('fc3 was read from "%s".' % fc3_filename) - elif forces_fc3_filename is not None: - force_filename = forces_fc3_filename - _set_dataset_for_fc3( + print(f'fc3 was read from "{_fc3_filename}".') + elif ( + ph3py_yaml is not None + and ph3py_yaml.dataset is not None + and forces_in_dataset(ph3py_yaml.dataset) + ): + dataset = _get_dataset_for_fc3( + ph3py, + ph3py_yaml, + None, + phono3py_yaml_filename, + cutoff_pair_distance, + log_level, + ) + elif forces_fc3_filename is not None or pathlib.Path("FORCES_FC3").exists(): + if forces_fc3_filename is None: + force_filename = "FORCES_FC3" + else: + force_filename = forces_fc3_filename + dataset = _get_dataset_for_fc3( ph3py, ph3py_yaml, force_filename, @@ -495,132 +542,74 @@ def _set_dataset_or_fc3( cutoff_pair_distance, log_level, ) - elif pathlib.Path("fc3.hdf5").exists(): - fc3 = read_fc3_from_hdf5(filename="fc3.hdf5", p2s_map=p2s_map) - _check_fc3_shape(ph3py, fc3) - ph3py.fc3 = fc3 - read_fc3 = True - if log_level: - print('fc3 was read from "fc3.hdf5".') - elif dataset := read_type2_dataset( - len(ph3py.supercell), "FORCES_FC3", log_level=log_level - ): # := Assignment Expressions (Python>=3.8) - ph3py.dataset = dataset - elif ph3py_yaml is not None and ph3py_yaml.dataset is not None: - if pathlib.Path("FORCES_FC3").exists() and displacements_in_dataset( - ph3py_yaml.dataset - ): - _set_dataset_for_fc3( - ph3py, - ph3py_yaml, - "FORCES_FC3", - phono3py_yaml_filename, - cutoff_pair_distance, - log_level, - ) - if forces_in_dataset(ph3py_yaml.dataset): - _set_dataset_for_fc3( - ph3py, - ph3py_yaml, - None, - phono3py_yaml_filename, - cutoff_pair_distance, - log_level, - ) - return read_fc3 + if not forces_in_dataset(dataset): + dataset = None + + return read_fc3, dataset -def _set_dataset_phonon_dataset_or_fc2( +def _get_dataset_phonon_dataset_or_fc2( ph3py: Phono3py, ph3py_yaml: Optional[Phono3pyYaml] = None, fc2_filename: Optional[os.PathLike] = None, forces_fc2_filename: Optional[Union[os.PathLike, Sequence]] = None, log_level: int = 0, -) -> bool: +) -> tuple[bool, dict, dict]: phonon_p2s_map = ph3py.phonon_primitive.p2s_map read_fc2 = False - if fc2_filename is not None: - fc2 = read_fc2_from_hdf5(filename=fc2_filename, p2s_map=phonon_p2s_map) - _check_fc2_shape(ph3py, fc2, filename=fc2_filename) + phonon_dataset = None + if fc2_filename is not None or pathlib.Path("fc2.hdf5").exists(): + if fc2_filename is None: + _fc2_filename = "fc2.hdf5" + else: + _fc2_filename = fc2_filename + fc2 = read_fc2_from_hdf5(filename=_fc2_filename, p2s_map=phonon_p2s_map) + _check_fc2_shape(ph3py, fc2, filename=_fc2_filename) ph3py.fc2 = fc2 read_fc2 = True if log_level: - print('fc2 was read from "%s".' % fc2_filename) - elif forces_fc2_filename is not None: - force_filename = forces_fc2_filename - _set_dataset_for_fc2( + print(f'fc2 was read from "{fc2_filename}".') + elif ( + ph3py_yaml is not None + and ph3py_yaml.phonon_dataset is not None + and forces_in_dataset(ph3py_yaml.phonon_dataset) + ): + phonon_dataset = _get_dataset_for_fc2( + ph3py, + ph3py_yaml, + None, + "phonon_fc2", + log_level, + ) + elif ( + forces_fc2_filename is not None or pathlib.Path("FORCES_FC2").exists() + ) and ph3py.phonon_supercell_matrix: + if forces_fc2_filename is None: + force_filename = forces_fc2_filename + else: + force_filename = forces_fc2_filename + phonon_dataset = _get_dataset_for_fc2( ph3py, ph3py_yaml, force_filename, "phonon_fc2", log_level, ) - elif pathlib.Path("fc2.hdf5").exists(): - fc2 = read_fc2_from_hdf5(filename="fc2.hdf5", p2s_map=phonon_p2s_map) - _check_fc2_shape(ph3py, fc2) - ph3py.fc2 = fc2 - read_fc2 = True - if log_level: - print('fc2 was read from "fc2.hdf5".') - elif ( - pathlib.Path("FORCES_FC2").exists() - and ph3py.phonon_supercell_matrix is not None - ): - _set_dataset_for_fc2( - ph3py, - ph3py_yaml, - "FORCES_FC2", - "phonon_fc2", - log_level, - ) - elif ( - ph3py_yaml is not None - and ph3py_yaml.phonon_dataset is not None - and forces_in_dataset(ph3py_yaml.phonon_dataset) - ): - _set_dataset_for_fc2( - ph3py, - ph3py_yaml, - None, - "phonon_fc2", - log_level, - ) - elif ( - ph3py_yaml is not None - and ph3py_yaml.dataset is not None - and forces_in_dataset(ph3py_yaml.dataset) - ): - _set_dataset_for_fc2( - ph3py, - ph3py_yaml, - None, - "fc2", - log_level, - ) - elif ( - pathlib.Path("FORCES_FC3").exists() - and ph3py.phonon_supercell_matrix is not None - ): - # suppose fc3.hdf5 is read but fc2.hdf5 doesn't exist. - _set_dataset_for_fc2( - ph3py, - ph3py_yaml, - "FORCES_FC3", - "fc2", - log_level, - ) - return read_fc2 + if not forces_in_dataset(phonon_dataset): + phonon_dataset = None + + return read_fc2, phonon_dataset -def _set_dataset_for_fc3( +def _get_dataset_for_fc3( ph3py: Phono3py, ph3py_yaml: Optional[Phono3pyYaml], force_filename, phono3py_yaml_filename, cutoff_pair_distance, log_level, -): - ph3py.dataset = parse_forces( +) -> dict: + dataset = parse_forces( ph3py, ph3py_yaml=ph3py_yaml, cutoff_pair_distance=cutoff_pair_distance, @@ -629,9 +618,10 @@ def _set_dataset_for_fc3( fc_type="fc3", log_level=log_level, ) + return dataset -def _set_dataset_for_fc2( +def _get_dataset_for_fc2( ph3py: Phono3py, ph3py_yaml: Optional[Phono3pyYaml], force_filename, @@ -645,10 +635,7 @@ def _set_dataset_for_fc2( fc_type=fc_type, log_level=log_level, ) - if fc_type == "phonon_fc2": - ph3py.phonon_dataset = dataset - else: - ph3py.dataset = dataset + return dataset def _check_fc2_shape(ph3py: Phono3py, fc2, filename="fc2.hdf5"): diff --git a/phono3py/cui/phono3py_argparse.py b/phono3py/cui/phono3py_argparse.py index 48bf9820..6c51cb8b 100644 --- a/phono3py/cui/phono3py_argparse.py +++ b/phono3py/cui/phono3py_argparse.py @@ -471,6 +471,15 @@ def get_parser(fc_symmetry=False, is_nac=False, load_phono3py_yaml=False): parser.add_argument( "--mesh", nargs="+", dest="mesh_numbers", default=None, help="Mesh numbers" ) + parser.add_argument( + "--mlp-params", + dest="mlp_params", + default=None, + help=( + "Parameters for machine learning potentials as comma separated " + "string with the style of key = values" + ), + ) parser.add_argument( "--mv", "--mass-variances", @@ -622,7 +631,7 @@ def get_parser(fc_symmetry=False, is_nac=False, load_phono3py_yaml=False): "--pypolymlp", dest="use_pypolymlp", action="store_true", - default=None, + default=False, help="Use pypolymlp and symfc for generating force constants", ) parser.add_argument( @@ -647,6 +656,21 @@ def get_parser(fc_symmetry=False, is_nac=False, load_phono3py_yaml=False): default=False, help="Print out smallest information", ) + parser.add_argument( + "--random-seed", + dest="random_seed", + type=int, + default=None, + help="Random seed by a 32 bit unsigned integer", + ) + parser.add_argument( + "--rd", + "--random-displacements", + dest="random_displacements", + type=int, + default=None, + help="Number of supercells with random displacements", + ) parser.add_argument( "--read-collision", dest="read_collision", diff --git a/phono3py/cui/phono3py_script.py b/phono3py/cui/phono3py_script.py index 4cbacee1..2fda8b0c 100644 --- a/phono3py/cui/phono3py_script.py +++ b/phono3py/cui/phono3py_script.py @@ -528,6 +528,7 @@ def store_force_constants( ph3py_yaml=ph3py_yaml, phono3py_yaml_filename=phono3py_yaml_filename, cutoff_pair_distance=settings.cutoff_pair_distance, + use_pypolymlp=settings.use_pypolymlp, log_level=log_level, ) try: @@ -538,6 +539,11 @@ def store_force_constants( fc_calculator_options=fc_calculator_options, symmetrize_fc=settings.fc_symmetry, is_compact_fc=settings.is_compact_fc, + use_pypolymlp=settings.use_pypolymlp, + mlp_params=settings.mlp_params, + displacement_distance=settings.displacement_distance, + number_of_snapshots=settings.random_displacements, + random_seed=settings.random_seed, log_level=log_level, ) except ForceCalculatorRequiredError: @@ -605,8 +611,8 @@ def _show_fc_calculator_not_found(log_level): print( "Built-in force constants calculator doesn't support the " "dispalcements-forces dataset. " - "An external force calculator, e.g., ALM (--alm), has to be used " - "to compute force constants." + "An external force calculator, e.g., symfc (--symfc_ or ALM (--alm), " + "has to be used to compute force constants." ) print_error() sys.exit(1) diff --git a/phono3py/cui/show_log.py b/phono3py/cui/show_log.py index 2ef4f1ba..bc731832 100644 --- a/phono3py/cui/show_log.py +++ b/phono3py/cui/show_log.py @@ -107,7 +107,7 @@ def show_phono3py_cells(phono3py: Phono3py): print_cell(phonon_primitive) print("-" * 21 + " supercell for harmonic phonon " + "-" * 22) print_cell(phonon_supercell, mapping=phonon_primitive.s2p_map) - print("-" * 76) + print("-" * 76, flush=True) def show_phono3py_force_constants_settings(settings: Phono3pySettings): diff --git a/pyproject.toml b/pyproject.toml index 1827bc6e..52340c86 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -3,14 +3,14 @@ requires = ["setuptools", "wheel", "numpy"] [tool.ruff] line-length = 88 -select = [ +lint.select = [ "F", # Flake8 "B", # Black "I", # isort "E", # pycodestyle-errors "D", # pydocstyle ] -extend-ignore = [ +lint.extend-ignore = [ "D417", "D100", ] diff --git a/test/api/test_api_phono3py.py b/test/api/test_api_phono3py.py index 35a1b2fb..23ad9471 100644 --- a/test/api/test_api_phono3py.py +++ b/test/api/test_api_phono3py.py @@ -170,7 +170,7 @@ def test_type2_forces_energies_setter_Si(si_111_222_rd: Phono3py): np.testing.assert_allclose(ph3_in.phonon_forces + 1, ph3.phonon_forces) -def test_use_pypolymlp_mgs(mgo_222rd_444rd: Phono3py): +def test_use_pypolymlp_mgo(mgo_222rd_444rd: Phono3py): """Test use_pypolymlp in produce_fc3.""" pytest.importorskip("pypolymlp") @@ -193,12 +193,16 @@ def test_use_pypolymlp_mgs(mgo_222rd_444rd: Phono3py): atom_energies = {"Mg": -0.00896717, "O": -0.95743902} params = PypolymlpParams(gtinv_maxl=(4, 4), atom_energies=atom_energies) - ph3.produce_fc3(use_pypolymlp=True, mlp_params=params, fc_calculator="symfc") + + ph3.generate_displacements(distance=0.001, is_plusminus=True) + ph3.develop_mlp(params=params) + ph3.evaluate_mlp() + ph3.produce_fc3(fc_calculator="symfc") ph3.produce_fc2(fc_calculator="symfc") ph3.mesh_numbers = 30 ph3.init_phph_interaction() ph3.run_thermal_conductivity(temperatures=[300]) assert ( - pytest.approx(63.6001137, abs=1e-3) == ph3.thermal_conductivity.kappa[0, 0, 0] + pytest.approx(63.0018546, abs=1e-3) == ph3.thermal_conductivity.kappa[0, 0, 0] ) diff --git a/test/conftest.py b/test/conftest.py index 91fda478..61c2408d 100644 --- a/test/conftest.py +++ b/test/conftest.py @@ -415,7 +415,6 @@ def nacl_pbe_cutoff_fc3(request) -> Phono3py: count += 1 ph3.dataset = dataset ph3.produce_fc3() - # ph3.produce_fc3(symmetrize_fc3r=True) return ph3 @@ -504,7 +503,7 @@ def mgo_222rd_444rd() -> Phono3py: 4 and 400 supercells for fc2 and fc3, respectively. """ - yaml_filename = cwd / "phono3py_params_MgO-222rd-444fd.yaml.xz" + yaml_filename = cwd / "phono3py_params_MgO-222rd-444rd.yaml.xz" return phono3py.load(yaml_filename, produce_fc=False, log_level=1) diff --git a/test/phono3py_params_MgO-222rd-444fd.yaml.xz b/test/phono3py_params_MgO-222rd-444rd.yaml.xz similarity index 100% rename from test/phono3py_params_MgO-222rd-444fd.yaml.xz rename to test/phono3py_params_MgO-222rd-444rd.yaml.xz