scnc
/
phono3py
mirror of https://github.com/phonopy/phono3py.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #402 from phonopy/symfc-projector 

Use symfc-projector for FD when phono3py-load
This commit is contained in:
Atsushi Togo 2025-06-24 15:38:18 +09:00 committed by GitHub
parent 3ee42f3b81 22fa0e1862
commit ac2653147e
No known key found for this signature in database
GPG Key ID: B5690EEEBB952194
6 changed files with 189 additions and 126 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

91
phono3py/api_phono3py.py
View File

@ -59,7 +59,7 @@ from phonopy.interface.mlp import PhonopyMLP
from phonopy.interface.pypolymlp import (
PypolymlpParams,
)
from phonopy.interface.symfc import SymfcFCSolver
from phonopy.interface.symfc import SymfcFCSolver, symmetrize_by_projector
from phonopy.physical_units import get_physical_units
from phonopy.structure.atoms import PhonopyAtoms
from phonopy.structure.cells import (
@ -744,16 +744,21 @@ class Phono3py:
self._phonon_mlp_dataset = mlp_dataset
@property
def mlp(self) -> Optional[PhonopyMLP]:
def mlp(self) -> PhonopyMLP | None:
"""Setter and getter of PhonopyMLP dataclass."""
return self._mlp
@mlp.setter
def mlp(self, mlp) -> Optional[PhonopyMLP]:
def mlp(self, mlp: PhonopyMLP):
self._mlp = mlp
@property
def band_indices(self) -> list[NDArray]:
def phonon_mlp(self):
"""Return MLP instance for fc2."""
return self._phonon_mlp
@property
def band_indices(self) -> list[NDArray] | None:
"""Setter and getter of band indices.
list[NDArray]
@ -811,7 +816,7 @@ class Phono3py:
return self._supercells_with_displacements
@property
def phonon_supercells_with_displacements(self):
def phonon_supercells_with_displacements(self) -> list[PhonopyAtoms] | None:
"""Return supercells with displacements for fc2.
list of PhonopyAtoms
@ -1084,16 +1089,6 @@ class Phono3py:
"""
return self._bz_grid
@property
def mlp(self):
"""Return MLP instance."""
return self._mlp
@property
def phonon_mlp(self):
"""Return MLP instance for fc2."""
return self._phonon_mlp
def init_phph_interaction(
self,
nac_q_direction=None,
@ -1496,8 +1491,9 @@ class Phono3py:
self,
symmetrize_fc3r: bool = False,
is_compact_fc: bool = False,
fc_calculator: Optional[Union[str, dict]] = None,
fc_calculator_options: Optional[Union[str, dict]] = None,
fc_calculator: Literal["traditional", "symfc", "alm"] | None = None,
fc_calculator_options: str | None = None,
use_symfc_projector: bool = False,
):
"""Calculate fc3 from displacements and forces.
@ -1518,8 +1514,11 @@ class Phono3py:
cells. Default is False.
fc_calculator : str, optional
Force constants calculator given by str.
fc_calculator_options : dict or str, optional
fc_calculator_options : str, optional
Options for external force constants calculator.
use_symfc_projector : bool, optional
If True, the force constants are symmetrized by symfc projector
instead of traditional approach.
"""
fc_solver_name = fc_calculator if fc_calculator is not None else "traditional"
@ -1537,8 +1536,30 @@ class Phono3py:
fc2 = fc_solver.force_constants[2]
fc3 = fc_solver.force_constants[3]
if fc_calculator is None or fc_calculator == "traditional":
if symmetrize_fc3r:
if symmetrize_fc3r and (
fc_calculator is None or fc_calculator == "traditional"
):
if use_symfc_projector:
if self._log_level:
print("Symmetrizing fc3 by symfc projector.", flush=True)
fc3 = symmetrize_by_projector(
self._supercell,
fc3,
3,
primitive=self._primitive,
log_level=self._log_level,
)
if self._fc2 is None:
if self._log_level:
print("Symmetrizing fc2 by symfc projector.", flush=True)
fc2 = symmetrize_by_projector(
self._supercell,
fc2,
2,
primitive=self._primitive,
log_level=self._log_level,
)
else:
if is_compact_fc:
set_translational_invariance_compact_fc3(fc3, self._primitive)
set_permutation_symmetry_compact_fc3(fc3, self._primitive)
@ -1585,10 +1606,11 @@ class Phono3py:
def produce_fc2(
self,
symmetrize_fc2=False,
is_compact_fc=False,
fc_calculator=None,
fc_calculator_options=None,
symmetrize_fc2: bool = False,
is_compact_fc: bool = False,
fc_calculator: Literal["traditional", "symfc", "alm"] | None = None,
fc_calculator_options: str | None = None,
use_symfc_projector: bool = False,
):
"""Calculate fc2 from displacements and forces.
@ -1609,8 +1631,11 @@ class Phono3py:
cells. Default is False.
fc_calculator : str or None
Force constants calculator given by str.
fc_calculator_options : dict
fc_calculator_options : str or None
Options for external force constants calculator.
use_symfc_projector : bool, optional
If True, the force constants are symmetrized by symfc projector
instead of traditional approach.
"""
if self._phonon_dataset is None:
@ -1618,6 +1643,8 @@ class Phono3py:
else:
disp_dataset = self._phonon_dataset
if disp_dataset is None:
raise RuntimeError("Displacement dataset is not set.")
if not forces_in_dataset(disp_dataset):
raise RuntimeError("Forces are not set in the dataset.")
@ -1632,8 +1659,16 @@ class Phono3py:
log_level=self._log_level,
)
if fc_calculator is None or fc_calculator == "traditional":
if symmetrize_fc2:
if symmetrize_fc2 and (fc_calculator is None or fc_calculator == "traditional"):
if use_symfc_projector:
self._fc2 = symmetrize_by_projector(
self._phonon_supercell,
self._fc2,
2,
primitive=self._primitive,
log_level=self._log_level,
)
else:
if is_compact_fc:
symmetrize_compact_force_constants(
self._fc2, self._phonon_primitive
@ -2485,7 +2520,7 @@ class Phono3py:
def _build_phonon_supercells_with_displacements(
self, supercell: PhonopyAtoms, dataset
):
) -> list[PhonopyAtoms]:
supercells = []
positions = supercell.positions
magmoms = supercell.magnetic_moments

5
phono3py/cui/create_force_constants.py
View File

@ -104,6 +104,9 @@ def parse_forces(
if dataset:
filename_read_from = force_filename
if dataset is None:
raise RuntimeError("Dataset is not found.")
# Units of displacements and forces are converted. If forces don't
# exist, the conversion will not be performed for forces.
if calculator is not None:
@ -112,8 +115,6 @@ def parse_forces(
distance_to_A=physical_units["distance_to_A"],
force_to_eVperA=physical_units["force_to_eVperA"],
)
if dataset is None:
raise RuntimeError("Dataset is not found.")
if "natom" in dataset and dataset["natom"] != natom:
raise RuntimeError(

40
phono3py/cui/create_supercells.py
View File

@ -37,10 +37,14 @@
from __future__ import annotations
import dataclasses
import os
from typing import cast
import numpy as np
from numpy.typing import ArrayLike
from phonopy.cui.collect_cell_info import CellInfoResult
from phonopy.cui.collect_cell_info import get_cell_info as phonopy_get_cell_info
from phonopy.cui.settings import PhonopySettings
from phonopy.interface.calculator import write_supercells_with_displacements
from phonopy.structure.cells import print_cell
@ -53,6 +57,7 @@ from phono3py.interface.calculator import (
get_default_displacement_distance,
)
from phono3py.interface.fc_calculator import determine_cutoff_pair_distance
from phono3py.interface.phono3py_yaml import Phono3pyYaml
@dataclasses.dataclass
@ -63,9 +68,39 @@ class Phono3pyCellInfoResult(CellInfoResult):
"""
phono3py_yaml: Phono3pyYaml | None = None
phonon_supercell_matrix: ArrayLike | None = None
def get_cell_info(
settings: Phono3pySettings,
cell_filename: str | os.PathLike | None,
log_level: int,
load_phonopy_yaml: bool = True,
) -> Phono3pyCellInfoResult:
"""Return calculator interface and crystal structure information."""
cell_info = phonopy_get_cell_info(
cast(PhonopySettings, settings),
cell_filename,
log_level=log_level,
load_phonopy_yaml=load_phonopy_yaml,
phonopy_yaml_cls=Phono3pyYaml,
)
cell_info_dict = dataclasses.asdict(cell_info)
cell_info_dict["phono3py_yaml"] = cell_info_dict.pop("phonopy_yaml")
cell_info = Phono3pyCellInfoResult(
**cell_info_dict,
phonon_supercell_matrix=settings.phonon_supercell_matrix,
)
ph3py_yaml = cell_info.phono3py_yaml
if cell_info.phonon_supercell_matrix is None and ph3py_yaml:
cell_info.phonon_supercell_matrix = ph3py_yaml.phonon_supercell_matrix
return cell_info
def create_phono3py_supercells(
cell_info: Phono3pyCellInfoResult,
settings: Phono3pySettings,
@ -171,7 +206,10 @@ def create_phono3py_supercells(
print(f"Cutoff distance for displacements: {cutoff_pair_distance}")
print(f"Number of displacement supercell files created: {num_disp_files}")
if ph3.phonon_supercell_matrix is not None:
if (
ph3.phonon_supercell_matrix is not None
and ph3.phonon_supercells_with_displacements is not None
):
num_disps = len(ph3.phonon_supercells_with_displacements)
additional_info = get_additional_info_to_write_fc2_supercells(
interface_mode, ph3.phonon_supercell_matrix

8
phono3py/cui/load.py
View File

@ -38,6 +38,7 @@ from __future__ import annotations
import os
import pathlib
from collections.abc import Sequence
from typing import Literal
import numpy as np
import phonopy.cui.load_helper as load_helper
@ -411,12 +412,13 @@ def load_fc2_and_fc3(
def compute_force_constants_from_datasets(
ph3py: Phono3py,
fc_calculator: str | None = None,
fc_calculator_options: dict | str | None = None,
fc_calculator: Literal["traditional", "symfc", "alm"] | None = None,
fc_calculator_options: str | None = None,
cutoff_pair_distance: float | None = None,
symmetrize_fc: bool = True,
is_compact_fc: bool = True,
log_level: int = 0,
load_phono3py_yaml: bool = False,
):
"""Compute force constants from datasets.
@ -444,6 +446,7 @@ def compute_force_constants_from_datasets(
is_compact_fc=is_compact_fc,
fc_calculator=fc3_calculator,
fc_calculator_options=fc3_calc_opts,
use_symfc_projector=load_phono3py_yaml,
)
if log_level and symmetrize_fc and fc_calculator is None:
@ -461,6 +464,7 @@ def compute_force_constants_from_datasets(
is_compact_fc=is_compact_fc,
fc_calculator=fc2_calculator,
fc_calculator_options=fc2_calc_opts,
use_symfc_projector=load_phono3py_yaml,
)
if log_level and symmetrize_fc and fc_calculator is None:
print("fc2 was symmetrized.")

101
phono3py/cui/phono3py_script.py
View File

@ -37,7 +37,6 @@
from __future__ import annotations
import argparse
import dataclasses
import os
import pathlib
import sys
@ -48,7 +47,6 @@ from typing import Optional, cast
import numpy as np
from numpy.typing import NDArray
from phonopy.api_phonopy import Phonopy
from phonopy.cui.collect_cell_info import collect_cell_info
from phonopy.cui.phonopy_argparse import show_deprecated_option_warnings
from phonopy.cui.phonopy_script import (
file_exists,
@ -57,7 +55,6 @@ from phonopy.cui.phonopy_script import (
print_error_message,
print_time,
print_version,
set_magnetic_moments,
store_nac_params,
)
from phonopy.cui.settings import PhonopySettings
@ -80,6 +77,7 @@ from phono3py.cui.create_force_sets import (
from phono3py.cui.create_supercells import (
Phono3pyCellInfoResult,
create_phono3py_supercells,
get_cell_info,
)
from phono3py.cui.load import (
compute_force_constants_from_datasets,
@ -342,40 +340,6 @@ def _get_input_output_filenames_from_args(args: argparse.Namespace):
return input_filename, output_filename
def _get_cell_info(
settings: Phono3pySettings, cell_filename: str | os.PathLike | None, log_level: int
) -> Phono3pyCellInfoResult:
"""Return calculator interface and crystal structure information."""
try:
cell_info = collect_cell_info(
supercell_matrix=settings.supercell_matrix,
primitive_matrix=settings.primitive_matrix,
interface_mode=settings.calculator,
cell_filename=cell_filename,
chemical_symbols=settings.chemical_symbols,
phonopy_yaml_cls=Phono3pyYaml,
)
except CellNotFoundError as e:
print_error_message(str(e))
if log_level:
print_error()
sys.exit(1)
cell_info = Phono3pyCellInfoResult(
**dataclasses.asdict(cell_info),
phonon_supercell_matrix=settings.phonon_supercell_matrix,
)
set_magnetic_moments(cell_info.unitcell, settings.magnetic_moments, log_level)
ph3py_yaml = cast(Phono3pyYaml, cell_info.phonopy_yaml)
if cell_info.phonon_supercell_matrix is None and ph3py_yaml:
ph_smat = ph3py_yaml.phonon_supercell_matrix
cell_info.phonon_supercell_matrix = ph_smat
return cell_info
def _get_default_values(settings: Phono3pySettings):
"""Set default values."""
# Brillouin zone integration: Tetrahedron (default) or smearing method
@ -457,12 +421,12 @@ def _check_supercell_in_yaml(
log_level: int,
):
"""Check consistency between generated cells and cells in yaml."""
if cell_info.phonopy_yaml is not None:
if cell_info.phono3py_yaml is not None:
if distance_to_A is None:
d2A = 1.0
else:
d2A = distance_to_A
phono3py_yaml = cast(Phono3pyYaml, cell_info.phonopy_yaml)
phono3py_yaml = cell_info.phono3py_yaml
if phono3py_yaml.supercell is not None and ph3.supercell is not None: # noqa E129
yaml_cell = phono3py_yaml.supercell.copy()
yaml_cell.cell = yaml_cell.cell * d2A
@ -585,7 +549,7 @@ def _create_supercells_with_displacements(
store_nac_params(
cast(Phonopy, phono3py),
cast(PhonopySettings, settings),
cell_info.phonopy_yaml,
cell_info.phono3py_yaml,
unitcell_filename,
log_level,
nac_factor=get_physical_units().Hartree * get_physical_units().Bohr,
@ -610,7 +574,10 @@ def _create_supercells_with_displacements(
def _produce_force_constants(
ph3py: Phono3py, settings: Phono3pySettings, log_level: int
ph3py: Phono3py,
settings: Phono3pySettings,
log_level: int,
load_phono3py_yaml: bool,
):
"""Calculate, read, and write force constants."""
if log_level:
@ -620,7 +587,7 @@ def _produce_force_constants(
read_fc2 = ph3py.fc2 is not None
cutoff_pair_distance = None
if settings.use_pypolymlp:
if settings.use_pypolymlp and ph3py.dataset is not None:
cutoff_pair_distance = ph3py.dataset.get("cutoff_distance")
if cutoff_pair_distance is None:
cutoff_pair_distance = determine_cutoff_pair_distance(
@ -651,10 +618,17 @@ def _produce_force_constants(
symmetrize_fc=settings.fc_symmetry,
is_compact_fc=settings.is_compact_fc,
log_level=log_level,
load_phono3py_yaml=load_phono3py_yaml,
)
except ForceCalculatorRequiredError:
if log_level:
print("Symfc will be used to handle general (or random) displacements.")
except ForceCalculatorRequiredError as e:
if load_phono3py_yaml:
if log_level:
print("Symfc will be used to handle general (or random) displacements.")
else:
print_error_message(str(e))
if log_level:
print_error()
sys.exit(1)
compute_force_constants_from_datasets(
ph3py,
@ -671,8 +645,6 @@ def _produce_force_constants(
print("fc3 could not be obtained.")
if not forces_in_dataset(ph3py.dataset):
print("Forces were not found.")
else:
show_drift_fc3(ph3py.fc3, primitive=ph3py.primitive)
if ph3py.fc2 is None:
print("fc2 could not be obtained.")
if ph3py.phonon_supercell_matrix is None:
@ -681,15 +653,21 @@ def _produce_force_constants(
else:
if not forces_in_dataset(ph3py.phonon_dataset):
print("Forces for dim-fc2 were not found.")
else:
show_drift_force_constants(
ph3py.fc2, primitive=ph3py.phonon_primitive, name="fc2"
)
if ph3py.fc3 is None or ph3py.fc2 is None:
print_error()
sys.exit(1)
if log_level:
if load_phono3py_yaml:
print("Max drift after symmetrization by symfc projector: ")
else:
print("Max drift after symmetrization by translation: ")
show_drift_fc3(ph3py.fc3, primitive=ph3py.primitive)
show_drift_force_constants(
ph3py.fc2, primitive=ph3py.phonon_primitive, name="fc2"
)
cutoff_distance = settings.cutoff_fc3_distance
if cutoff_distance is not None and cutoff_distance > 0:
if log_level:
@ -1055,10 +1033,21 @@ def main(**argparse_control):
else:
symprec = settings.symmetry_tolerance
cell_info = _get_cell_info(settings, cell_filename, log_level)
try:
cell_info = get_cell_info(
settings,
cell_filename,
log_level=log_level,
load_phonopy_yaml=load_phono3py_yaml,
)
except CellNotFoundError as e:
print_error_message(str(e))
if log_level:
print_error()
sys.exit(1)
unitcell_filename = cell_info.optional_structure_info[0]
interface_mode = cell_info.interface_mode
# ph3py_yaml = cell_info.phonopy_yaml
if run_mode is None:
run_mode = _get_run_mode(settings)
@ -1213,7 +1202,7 @@ def main(**argparse_control):
store_nac_params(
cast(Phonopy, ph3py),
cast(PhonopySettings, settings),
cell_info.phonopy_yaml,
cell_info.phono3py_yaml,
unitcell_filename,
log_level,
nac_factor=get_physical_units().Hartree * get_physical_units().Bohr,
@ -1234,7 +1223,7 @@ def main(**argparse_control):
############
_load_dataset_and_phonon_dataset(
ph3py,
ph3py_yaml=cast(Phono3pyYaml, cell_info.phonopy_yaml),
ph3py_yaml=cell_info.phono3py_yaml,
phono3py_yaml_filename=unitcell_filename,
cutoff_pair_distance=settings.cutoff_pair_distance,
calculator=interface_mode,
@ -1290,7 +1279,7 @@ def main(**argparse_control):
###########################
# Produce force constants #
###########################
_produce_force_constants(ph3py, settings, log_level)
_produce_force_constants(ph3py, settings, log_level, load_phono3py_yaml)
############################################
# Phonon Gruneisen parameter and then exit #

70
phono3py/interface/fc_calculator.py
View File

@ -36,7 +36,7 @@
from __future__ import annotations
from typing import Optional, Union
from typing import Literal
import numpy as np
from phonopy.interface.fc_calculator import FCSolver, fc_calculator_names
@ -108,7 +108,7 @@ class FDFC3Solver:
class FC3Solver(FCSolver):
"""Force constants solver for fc3."""
def _set_traditional_solver(self, solver_class: Optional[type] = FDFC3Solver):
def _set_traditional_solver(self, solver_class: type | None = FDFC3Solver):
return super()._set_traditional_solver(solver_class=solver_class)
def _set_symfc_solver(self):
@ -120,9 +120,9 @@ class FC3Solver(FCSolver):
def extract_fc2_fc3_calculators(
fc_calculator: Optional[Union[str, dict]],
fc_calculator: str | None,
order: int,
) -> Optional[Union[str, dict]]:
) -> str | None:
"""Extract fc_calculator and fc_calculator_options for fc2 and fc3.
fc_calculator : str
@ -132,43 +132,39 @@ def extract_fc2_fc3_calculators(
2 and 3 indicate fc2 and fc3, respectively.
"""
if isinstance(fc_calculator, dict) or fc_calculator is None:
if fc_calculator is None:
return fc_calculator
elif isinstance(fc_calculator, str):
if "|" in fc_calculator:
_fc_calculator = fc_calculator.split("|")[order - 2]
if _fc_calculator == "":
return None
return _fc_calculator
else:
if fc_calculator.strip() == "":
return None
return fc_calculator
else:
_fc_calculator = fc_calculator
return _fc_calculator
else:
raise RuntimeError("fc_calculator should be str, dict, or None.")
raise RuntimeError("fc_calculator should be str or None.")
def update_cutoff_fc_calculator_options(
fc_calc_opts: Optional[Union[str, dict]],
cutoff_pair_distance: Optional[float],
) -> Optional[Union[str, dict]]:
fc_calc_opts: str | None,
cutoff_pair_distance: float | None,
) -> str | None:
"""Update fc_calculator_options with cutoff distances.
Parameters
----------
fc_calc_opts : str or dict
fc_calc_opts : str or None
FC calculator options.
cutoff_pair_distance : float, optional
Cutoff distance for pair interaction.
"""
if cutoff_pair_distance is not None:
if not isinstance(fc_calc_opts, (str, dict)) and fc_calc_opts is not None:
raise RuntimeError("fc_calculator_options should be str, dict, or None.")
if isinstance(fc_calc_opts, dict) and "cutoff" not in fc_calc_opts:
fc_calc_opts["cutoff"] = float(cutoff_pair_distance)
elif isinstance(fc_calc_opts, str) and "cutoff" not in fc_calc_opts:
if isinstance(fc_calc_opts, str) and "cutoff" not in fc_calc_opts:
fc_calc_opts = f"{fc_calc_opts}, cutoff = {cutoff_pair_distance}"
elif fc_calc_opts is None:
fc_calc_opts = f"cutoff = {cutoff_pair_distance}"
@ -177,11 +173,11 @@ def update_cutoff_fc_calculator_options(
def get_fc_calculator_params(
fc_calculator: Optional[str],
fc_calculator_options: Optional[str],
cutoff_pair_distance: Optional[float],
fc_calculator: str | None,
fc_calculator_options: str | None,
cutoff_pair_distance: float | None,
log_level: int = 0,
) -> tuple[Optional[str], Optional[str]]:
) -> tuple[Literal["traditional", "symfc", "alm"] | None, str | None]:
"""Compile fc_calculator and fc_calculator_options from input settings."""
_fc_calculator = None
fc_calculator_list = []
@ -211,16 +207,16 @@ def get_fc_calculator_params(
def determine_cutoff_pair_distance(
fc_calculator: Optional[str] = None,
fc_calculator_options: Optional[str] = None,
cutoff_pair_distance: Optional[float] = None,
symfc_memory_size: Optional[float] = None,
random_displacements: Optional[Union[int, str]] = None,
supercell: Optional[PhonopyAtoms] = None,
primitive: Optional[Primitive] = None,
symmetry: Optional[Symmetry] = None,
fc_calculator: str | None = None,
fc_calculator_options: str | None = None,
cutoff_pair_distance: float | None = None,
symfc_memory_size: float | None = None,
random_displacements: int | str | None = None,
supercell: PhonopyAtoms | None = None,
primitive: Primitive | None = None,
symmetry: Symmetry | None = None,
log_level: int = 0,
) -> float:
) -> float | None:
"""Determine cutoff pair distance for displacements."""
_cutoff_pair_distance, _symfc_memory_size = _get_cutoff_pair_distance(
fc_calculator,
@ -248,7 +244,7 @@ def determine_cutoff_pair_distance(
def _set_cutoff_in_fc_calculator_options(
fc_calculator_options: Optional[str],
fc_calculator_options: str | None,
cutoff_str: str,
log_level: int,
):
@ -275,11 +271,11 @@ def _set_cutoff_in_fc_calculator_options(
def _get_cutoff_pair_distance(
fc_calculator: Optional[str],
fc_calculator_options: Optional[str],
cutoff_pair_distance: Optional[float],
symfc_memory_size: Optional[float] = None,
) -> Optional[float]:
fc_calculator: str | None,
fc_calculator_options: str | None,
cutoff_pair_distance: float | None,
symfc_memory_size: float | None,
) -> tuple[float | None, float | None]:
"""Return cutoff_pair_distance from settings."""
_, _fc_calculator_options = get_fc_calculator_params(
fc_calculator,