From 7d467befc2010a38adfc7420fb6992ddf8025b02 Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Sun, 11 Aug 2024 19:45:04 -0700 Subject: [PATCH 01/14] Add quasi-harmonic phonon calculation to matcalc MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit The harmonic Hamiltonian of an ideal system is independent of volume, leading to a zero volume expansion coefficient. In order to accurately describe lattice expansion, it is necessary to incorporate anharmonic effects. The quasi-harmonic approximation (QHA) extends beyond the harmonic model by accounting for the variation of phonon frequencies with changes in the equilibrium volume induced by temperature. QHA provides a reasonable approximation over a broad temperature range below the melting point, balancing computational efficiency with accuracy when compared to fully anharmonic methods. In quasi-harmonic phonon calculation, it is of the most importance to compute the phonon properties of specific material for slightly expanded and reduced volumes and find the equilibrium volume at the given temperature by minimizing the Gibbs free energy equation with respect to volume at each temperature. However, first-principles harmonic phonon calculations have been already very expensive, not to mention that for QHA, these calculations are required to perform for multiple volumes. Thanks to the advent and development of universal machine interatomic potentials (uMLIPs), it becomes possible to efficiently calculate materials properties with near-ab initio accuracy from potential energy surface. The code developed here firstly achieves the open-source application of uMLIPs in quasi-harmonic phonon calculations, allowing users to efficiently compute temperature-dependent thermodynamic properties, including but not limited to thermal expansion coefficient, Gibbs free energy, Grüneisen parameters, bulk modulus and heat capacity at constant pressure with one line. Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- qha.py | 183 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 183 insertions(+) create mode 100644 qha.py diff --git a/qha.py b/qha.py new file mode 100644 index 0000000..8d7492f --- /dev/null +++ b/qha.py @@ -0,0 +1,183 @@ +"""Calculator for phonon properties under quasi-harmonic approxiamtion.""" + +from __future__ import annotations + +import numpy as np + +from dataclasses import dataclass, field +from typing import TYPE_CHECKING + +import numpy as np +import copy +from phonopy import PhonopyQHA + +from .base import PropCalc +from .relaxation import RelaxCalc +from .phonon import PhononCalc + +if TYPE_CHECKING: + from ase.calculators.calculator import Calculator + from pymatgen.core import Structure + +@dataclass +class QhaPhononCalc(PropCalc): + """Calculator for phonon properties under quasi-harmonic approxiamtion. + + Args: + calculator (Calculator): ASE Calculator to use. + t_step (float): Temperature step. Defaults to 10 (in Kelvin) + t_max (float): Max temperature (in Kelvin). Defaults to 1000 (in Kelvin) + t_min (float): Min temperature (in Kelvin). Defaults to 0 (in Kelvin) + fmax (float): Max forces. This criterion is more stringent than for simple relaxation. + Defaults to 0.1 (in eV/Angstrom) + optimizer (str): Optimizer used for RelaxCalc. Default to "FIRE" + eos (str): Equation of state used to fit F vs V, including "vinet", "murnaghan" or + "birch_murnaghan". Default to "vinet". + relax_structure (bool): Whether to first relax the structure. Set to False if structures + provided are pre-relaxed with the same calculator. + relax_calc_kwargs (dict): Arguments to be passed to the RelaxCalc, if relax_structure is True. + phonon_calc_kwargs (dict): Arguments to be passed to the PhononCalc. + scale_factors (list[float]): Factors to scale the lattice constants of the structure. + write_helmholtz_volume (bool | str | Path): Whether to save Helmholtz free energy vs volume in file. + Pass string or Path for custom filename. Defaults to False. + write_volume_temperature (bool | str | Path): Whether to save equilibrium volume vs temperature in file. + Pass string or Path for custom filename. Defaults to False. + write_thermal_expansion (bool | str | Path): Whether to save thermal expansion vs temperature in file. + Pass string or Path for custom filename. Defaults to False. + write_gibbs_temperature (bool | str | Path): Whether to save Gibbs free energy vs temperature in file. + Pass string or Path for custom filename. Defaults to False. + write_bulk_modulus_temperature (bool | str | Path): Whether to save bulk modulus vs temperature in file. + Pass string or Path for custom filename. Defaults to False. + write_heat_capacity_P_numerical (bool | str | Path): Whether to save heat capacity at constant pressure + by numerical difference vs temperature in file. Pass string or Path for custom filename. + Defaults to False. + write_heat_capacity_P_polyfit (bool | str | Path): Whether to save heat capacity at constant pressure + by fitting vs temperature in file. Pass string or Path for custom filename. Defaults to False. + write_gruneisen_temperature (bool | str | Path): Whether to save Grueneisen parameter vs temperature in + file. Pass string or Path for custom filename. Defaults to False. + """ + + calculator: Calculator + t_step: float = 10 + t_max: float = 1000 + t_min: float = 0 + fmax: float = 0.1 + optimizer: str = "FIRE" + eos: str = "vinet" + relax_structure: bool = True + relax_calc_kwargs: dict | None = None + phonon_calc_kwargs: dict | None = None + scale_factors: list[float] = field(default_factory=lambda: np.arange(0.95, 1.05, 0.01).tolist()) + write_helmholtz_volume: bool | str | Path = False + write_volume_temperature: bool | str | Path = False + write_thermal_expansion: bool | str | Path = False + write_gibbs_temperature: bool | str | Path = False + write_bulk_modulus_temperature: bool | str | Path = False + write_heat_capacity_P_numerical: bool | str | Path = False + write_heat_capacity_P_polyfit: bool | str | Path = False + write_gruneisen_temperature: bool | str | Path = False + + def __post_init__(self) -> None: + """Set default paths for where to save output files.""" + # map True to canonical default path, False to "" and Path to str + for key, val, default_path in ( + ("write_helmholtz_volume", self.write_helmholtz_volume, "helmholtz_volume.dat"), + ("write_volume_temperature", self.write_volume_temperature, "volume_temperature.dat"), + ("write_thermal_expansion", self.write_thermal_expansion, "thermal_expansion.dat"), + ("write_write_gibbs_temperature", self.write_gibbs_temperature, "gibbs_temperature.dat"), + ("write_bulk_modulus_temperature", self.write_bulk_modulus_temperature, "bulk_modulus_temperature.dat"), + ("write_heat_capacity_P_numerical", self.write_heat_capacity_P_numerical, "Cp_temperature.dat"), + ("write_heat_capacity_P_polyfit", self.write_heat_capacity_P_polyfit, "Cp_temperature_polyfit.dat"), + ("write_gruneisen_temperature", self.write_gruneisen_temperature, "gruneisen_temperature.dat"), + ): + setattr(self, key, str({True: default_path, False: ""}.get(val, val))) # type: ignore[arg-type] + + def calc(self, structure: Structure) -> dict: + """Calculates thermal properties of Pymatgen structure with phonopy under quasi-harmonic approxiamtion. + + Args: + structure: Pymatgen structure. + + Returns: + { + "scale_factors": list of scale factors of lattice constants, + "volumes": list of unit cell volumes at corresponding scale factors (in Angstrom^3), + "electronic_energies": list of electronic energies at corresponding volumes (in eV), + "temperatures": list of temperatures in ascending order (in Kelvin), + "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding + temperatures (in Kelvin^-1), + "gibbs_free_energies": list of Gibbs free energies at corresponding temperatures (in eV), + "bulk_modulus_P": list of bulk modulus at constant presure at corresponding temperatures (in GPa), + "heat_capacity_P": list of heat capacities at constant pressure at corresponding temperatures (in J/K/mol), + "gruneisen_parameters": list of Gruneisen parameters at corresponding temperatures, + } + """ + volumes = [] + electronic_energies = [] + temperatures = np.arange(self.t_min, self.t_max+self.t_step, self.t_step) + + free_energies = [] + entropies = [] + heat_capacities = [] + + if self.relax_structure: + relaxer = RelaxCalc( + self.calculator, fmax=self.fmax, optimizer=self.optimizer, **(self.relax_calc_kwargs or {}) + ) + structure = relaxer.calc(structure)["final_structure"] + + for scale_factor in self.scale_factors: + struct = copy.deepcopy(structure) + struct.scale_lattice(struct.volume*scale_factor**3) + + static_calc = RelaxCalc( + self.calculator, relax_atoms=False, relax_cell=False) + volumes.append(struct.volume) + electronic_energies.append(static_calc.calc(struct)["energy"]) + + phonon_calc = PhononCalc( + self.calculator, t_step=self.t_step, t_max=self.t_max, t_min=self.t_min, relax_structure=False, write_phonon=False, **(self.phonon_calc_kwargs or {}) + ) + thermal_properties = phonon_calc.calc(struct)["thermal_properties"] + free_energies.append(thermal_properties["free_energy"]) + entropies.append(thermal_properties["entropy"]) + heat_capacities.append(thermal_properties["heat_capacity"]) + + qha = PhonopyQHA( + volumes=volumes, + electronic_energies=electronic_energies, + temperatures=temperatures, + free_energy=np.transpose(free_energies), + entropy=np.transpose(entropies), + cv=np.transpose(heat_capacities), + eos=self.eos, + t_max=self.t_max + ) + + if self.write_helmholtz_volume: + qha.write_helmholtz_volume(filename=self.write_helmholtz_volume) + if self.write_volume_temperature: + qha.write_volume_temperature(filename=self.write_volume_temperature) + if self.write_thermal_expansion: + qha.write_thermal_expansion(filename=self.write_thermal_expansion) + if self.write_gibbs_temperature: + qha.write_gibbs_temperature(filename=self.write_gibbs_temperature) + if self.write_bulk_modulus_temperature: + qha.write_bulk_modulus_temperature(filename=self.write_bulk_modulus_temperature) + if self.write_heat_capacity_P_numerical: + qha.write_heat_capacity_P_numerical(filename=self.write_heat_capacity_P_numerical) + if self.write_heat_capacity_P_polyfit: + qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_P_polyfit) + if self.write_gruneisen_temperature: + qha.write_gruneisen_temperature(filename=self.write_gruneisen_temperature) + + return {"scale_factors": self.scale_factors, + "volumes": volumes, + "electronic_energies": electronic_energies, + "temperatures": temperatures, + "thermal_expansion_coefficients": qha.thermal_expansion, + "gibbs_free_energies": qha.gibbs_temperature, + "bulk_modulus_P": qha.bulk_modulus_temperature, + "heat_capacity_P": qha.heat_capacity_P_polyfit, + "gruneisen_parameters": qha.gruneisen_temperature, + } From 4ac7f5862266975b3c6af2638f03c0acd8c4dccf Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Sun, 11 Aug 2024 20:24:12 -0700 Subject: [PATCH 02/14] Add quasi-harmonic phonon calculation to matcalc Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- matcalc/qha.py | 183 +++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 183 insertions(+) create mode 100644 matcalc/qha.py diff --git a/matcalc/qha.py b/matcalc/qha.py new file mode 100644 index 0000000..8d7492f --- /dev/null +++ b/matcalc/qha.py @@ -0,0 +1,183 @@ +"""Calculator for phonon properties under quasi-harmonic approxiamtion.""" + +from __future__ import annotations + +import numpy as np + +from dataclasses import dataclass, field +from typing import TYPE_CHECKING + +import numpy as np +import copy +from phonopy import PhonopyQHA + +from .base import PropCalc +from .relaxation import RelaxCalc +from .phonon import PhononCalc + +if TYPE_CHECKING: + from ase.calculators.calculator import Calculator + from pymatgen.core import Structure + +@dataclass +class QhaPhononCalc(PropCalc): + """Calculator for phonon properties under quasi-harmonic approxiamtion. + + Args: + calculator (Calculator): ASE Calculator to use. + t_step (float): Temperature step. Defaults to 10 (in Kelvin) + t_max (float): Max temperature (in Kelvin). Defaults to 1000 (in Kelvin) + t_min (float): Min temperature (in Kelvin). Defaults to 0 (in Kelvin) + fmax (float): Max forces. This criterion is more stringent than for simple relaxation. + Defaults to 0.1 (in eV/Angstrom) + optimizer (str): Optimizer used for RelaxCalc. Default to "FIRE" + eos (str): Equation of state used to fit F vs V, including "vinet", "murnaghan" or + "birch_murnaghan". Default to "vinet". + relax_structure (bool): Whether to first relax the structure. Set to False if structures + provided are pre-relaxed with the same calculator. + relax_calc_kwargs (dict): Arguments to be passed to the RelaxCalc, if relax_structure is True. + phonon_calc_kwargs (dict): Arguments to be passed to the PhononCalc. + scale_factors (list[float]): Factors to scale the lattice constants of the structure. + write_helmholtz_volume (bool | str | Path): Whether to save Helmholtz free energy vs volume in file. + Pass string or Path for custom filename. Defaults to False. + write_volume_temperature (bool | str | Path): Whether to save equilibrium volume vs temperature in file. + Pass string or Path for custom filename. Defaults to False. + write_thermal_expansion (bool | str | Path): Whether to save thermal expansion vs temperature in file. + Pass string or Path for custom filename. Defaults to False. + write_gibbs_temperature (bool | str | Path): Whether to save Gibbs free energy vs temperature in file. + Pass string or Path for custom filename. Defaults to False. + write_bulk_modulus_temperature (bool | str | Path): Whether to save bulk modulus vs temperature in file. + Pass string or Path for custom filename. Defaults to False. + write_heat_capacity_P_numerical (bool | str | Path): Whether to save heat capacity at constant pressure + by numerical difference vs temperature in file. Pass string or Path for custom filename. + Defaults to False. + write_heat_capacity_P_polyfit (bool | str | Path): Whether to save heat capacity at constant pressure + by fitting vs temperature in file. Pass string or Path for custom filename. Defaults to False. + write_gruneisen_temperature (bool | str | Path): Whether to save Grueneisen parameter vs temperature in + file. Pass string or Path for custom filename. Defaults to False. + """ + + calculator: Calculator + t_step: float = 10 + t_max: float = 1000 + t_min: float = 0 + fmax: float = 0.1 + optimizer: str = "FIRE" + eos: str = "vinet" + relax_structure: bool = True + relax_calc_kwargs: dict | None = None + phonon_calc_kwargs: dict | None = None + scale_factors: list[float] = field(default_factory=lambda: np.arange(0.95, 1.05, 0.01).tolist()) + write_helmholtz_volume: bool | str | Path = False + write_volume_temperature: bool | str | Path = False + write_thermal_expansion: bool | str | Path = False + write_gibbs_temperature: bool | str | Path = False + write_bulk_modulus_temperature: bool | str | Path = False + write_heat_capacity_P_numerical: bool | str | Path = False + write_heat_capacity_P_polyfit: bool | str | Path = False + write_gruneisen_temperature: bool | str | Path = False + + def __post_init__(self) -> None: + """Set default paths for where to save output files.""" + # map True to canonical default path, False to "" and Path to str + for key, val, default_path in ( + ("write_helmholtz_volume", self.write_helmholtz_volume, "helmholtz_volume.dat"), + ("write_volume_temperature", self.write_volume_temperature, "volume_temperature.dat"), + ("write_thermal_expansion", self.write_thermal_expansion, "thermal_expansion.dat"), + ("write_write_gibbs_temperature", self.write_gibbs_temperature, "gibbs_temperature.dat"), + ("write_bulk_modulus_temperature", self.write_bulk_modulus_temperature, "bulk_modulus_temperature.dat"), + ("write_heat_capacity_P_numerical", self.write_heat_capacity_P_numerical, "Cp_temperature.dat"), + ("write_heat_capacity_P_polyfit", self.write_heat_capacity_P_polyfit, "Cp_temperature_polyfit.dat"), + ("write_gruneisen_temperature", self.write_gruneisen_temperature, "gruneisen_temperature.dat"), + ): + setattr(self, key, str({True: default_path, False: ""}.get(val, val))) # type: ignore[arg-type] + + def calc(self, structure: Structure) -> dict: + """Calculates thermal properties of Pymatgen structure with phonopy under quasi-harmonic approxiamtion. + + Args: + structure: Pymatgen structure. + + Returns: + { + "scale_factors": list of scale factors of lattice constants, + "volumes": list of unit cell volumes at corresponding scale factors (in Angstrom^3), + "electronic_energies": list of electronic energies at corresponding volumes (in eV), + "temperatures": list of temperatures in ascending order (in Kelvin), + "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding + temperatures (in Kelvin^-1), + "gibbs_free_energies": list of Gibbs free energies at corresponding temperatures (in eV), + "bulk_modulus_P": list of bulk modulus at constant presure at corresponding temperatures (in GPa), + "heat_capacity_P": list of heat capacities at constant pressure at corresponding temperatures (in J/K/mol), + "gruneisen_parameters": list of Gruneisen parameters at corresponding temperatures, + } + """ + volumes = [] + electronic_energies = [] + temperatures = np.arange(self.t_min, self.t_max+self.t_step, self.t_step) + + free_energies = [] + entropies = [] + heat_capacities = [] + + if self.relax_structure: + relaxer = RelaxCalc( + self.calculator, fmax=self.fmax, optimizer=self.optimizer, **(self.relax_calc_kwargs or {}) + ) + structure = relaxer.calc(structure)["final_structure"] + + for scale_factor in self.scale_factors: + struct = copy.deepcopy(structure) + struct.scale_lattice(struct.volume*scale_factor**3) + + static_calc = RelaxCalc( + self.calculator, relax_atoms=False, relax_cell=False) + volumes.append(struct.volume) + electronic_energies.append(static_calc.calc(struct)["energy"]) + + phonon_calc = PhononCalc( + self.calculator, t_step=self.t_step, t_max=self.t_max, t_min=self.t_min, relax_structure=False, write_phonon=False, **(self.phonon_calc_kwargs or {}) + ) + thermal_properties = phonon_calc.calc(struct)["thermal_properties"] + free_energies.append(thermal_properties["free_energy"]) + entropies.append(thermal_properties["entropy"]) + heat_capacities.append(thermal_properties["heat_capacity"]) + + qha = PhonopyQHA( + volumes=volumes, + electronic_energies=electronic_energies, + temperatures=temperatures, + free_energy=np.transpose(free_energies), + entropy=np.transpose(entropies), + cv=np.transpose(heat_capacities), + eos=self.eos, + t_max=self.t_max + ) + + if self.write_helmholtz_volume: + qha.write_helmholtz_volume(filename=self.write_helmholtz_volume) + if self.write_volume_temperature: + qha.write_volume_temperature(filename=self.write_volume_temperature) + if self.write_thermal_expansion: + qha.write_thermal_expansion(filename=self.write_thermal_expansion) + if self.write_gibbs_temperature: + qha.write_gibbs_temperature(filename=self.write_gibbs_temperature) + if self.write_bulk_modulus_temperature: + qha.write_bulk_modulus_temperature(filename=self.write_bulk_modulus_temperature) + if self.write_heat_capacity_P_numerical: + qha.write_heat_capacity_P_numerical(filename=self.write_heat_capacity_P_numerical) + if self.write_heat_capacity_P_polyfit: + qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_P_polyfit) + if self.write_gruneisen_temperature: + qha.write_gruneisen_temperature(filename=self.write_gruneisen_temperature) + + return {"scale_factors": self.scale_factors, + "volumes": volumes, + "electronic_energies": electronic_energies, + "temperatures": temperatures, + "thermal_expansion_coefficients": qha.thermal_expansion, + "gibbs_free_energies": qha.gibbs_temperature, + "bulk_modulus_P": qha.bulk_modulus_temperature, + "heat_capacity_P": qha.heat_capacity_P_polyfit, + "gruneisen_parameters": qha.gruneisen_temperature, + } From 5379de40b6aa11aeba7527f6634b24427acecf76 Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Sun, 11 Aug 2024 20:26:24 -0700 Subject: [PATCH 03/14] Delete the misplaced file Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- qha.py | 183 --------------------------------------------------------- 1 file changed, 183 deletions(-) delete mode 100644 qha.py diff --git a/qha.py b/qha.py deleted file mode 100644 index 8d7492f..0000000 --- a/qha.py +++ /dev/null @@ -1,183 +0,0 @@ -"""Calculator for phonon properties under quasi-harmonic approxiamtion.""" - -from __future__ import annotations - -import numpy as np - -from dataclasses import dataclass, field -from typing import TYPE_CHECKING - -import numpy as np -import copy -from phonopy import PhonopyQHA - -from .base import PropCalc -from .relaxation import RelaxCalc -from .phonon import PhononCalc - -if TYPE_CHECKING: - from ase.calculators.calculator import Calculator - from pymatgen.core import Structure - -@dataclass -class QhaPhononCalc(PropCalc): - """Calculator for phonon properties under quasi-harmonic approxiamtion. - - Args: - calculator (Calculator): ASE Calculator to use. - t_step (float): Temperature step. Defaults to 10 (in Kelvin) - t_max (float): Max temperature (in Kelvin). Defaults to 1000 (in Kelvin) - t_min (float): Min temperature (in Kelvin). Defaults to 0 (in Kelvin) - fmax (float): Max forces. This criterion is more stringent than for simple relaxation. - Defaults to 0.1 (in eV/Angstrom) - optimizer (str): Optimizer used for RelaxCalc. Default to "FIRE" - eos (str): Equation of state used to fit F vs V, including "vinet", "murnaghan" or - "birch_murnaghan". Default to "vinet". - relax_structure (bool): Whether to first relax the structure. Set to False if structures - provided are pre-relaxed with the same calculator. - relax_calc_kwargs (dict): Arguments to be passed to the RelaxCalc, if relax_structure is True. - phonon_calc_kwargs (dict): Arguments to be passed to the PhononCalc. - scale_factors (list[float]): Factors to scale the lattice constants of the structure. - write_helmholtz_volume (bool | str | Path): Whether to save Helmholtz free energy vs volume in file. - Pass string or Path for custom filename. Defaults to False. - write_volume_temperature (bool | str | Path): Whether to save equilibrium volume vs temperature in file. - Pass string or Path for custom filename. Defaults to False. - write_thermal_expansion (bool | str | Path): Whether to save thermal expansion vs temperature in file. - Pass string or Path for custom filename. Defaults to False. - write_gibbs_temperature (bool | str | Path): Whether to save Gibbs free energy vs temperature in file. - Pass string or Path for custom filename. Defaults to False. - write_bulk_modulus_temperature (bool | str | Path): Whether to save bulk modulus vs temperature in file. - Pass string or Path for custom filename. Defaults to False. - write_heat_capacity_P_numerical (bool | str | Path): Whether to save heat capacity at constant pressure - by numerical difference vs temperature in file. Pass string or Path for custom filename. - Defaults to False. - write_heat_capacity_P_polyfit (bool | str | Path): Whether to save heat capacity at constant pressure - by fitting vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_gruneisen_temperature (bool | str | Path): Whether to save Grueneisen parameter vs temperature in - file. Pass string or Path for custom filename. Defaults to False. - """ - - calculator: Calculator - t_step: float = 10 - t_max: float = 1000 - t_min: float = 0 - fmax: float = 0.1 - optimizer: str = "FIRE" - eos: str = "vinet" - relax_structure: bool = True - relax_calc_kwargs: dict | None = None - phonon_calc_kwargs: dict | None = None - scale_factors: list[float] = field(default_factory=lambda: np.arange(0.95, 1.05, 0.01).tolist()) - write_helmholtz_volume: bool | str | Path = False - write_volume_temperature: bool | str | Path = False - write_thermal_expansion: bool | str | Path = False - write_gibbs_temperature: bool | str | Path = False - write_bulk_modulus_temperature: bool | str | Path = False - write_heat_capacity_P_numerical: bool | str | Path = False - write_heat_capacity_P_polyfit: bool | str | Path = False - write_gruneisen_temperature: bool | str | Path = False - - def __post_init__(self) -> None: - """Set default paths for where to save output files.""" - # map True to canonical default path, False to "" and Path to str - for key, val, default_path in ( - ("write_helmholtz_volume", self.write_helmholtz_volume, "helmholtz_volume.dat"), - ("write_volume_temperature", self.write_volume_temperature, "volume_temperature.dat"), - ("write_thermal_expansion", self.write_thermal_expansion, "thermal_expansion.dat"), - ("write_write_gibbs_temperature", self.write_gibbs_temperature, "gibbs_temperature.dat"), - ("write_bulk_modulus_temperature", self.write_bulk_modulus_temperature, "bulk_modulus_temperature.dat"), - ("write_heat_capacity_P_numerical", self.write_heat_capacity_P_numerical, "Cp_temperature.dat"), - ("write_heat_capacity_P_polyfit", self.write_heat_capacity_P_polyfit, "Cp_temperature_polyfit.dat"), - ("write_gruneisen_temperature", self.write_gruneisen_temperature, "gruneisen_temperature.dat"), - ): - setattr(self, key, str({True: default_path, False: ""}.get(val, val))) # type: ignore[arg-type] - - def calc(self, structure: Structure) -> dict: - """Calculates thermal properties of Pymatgen structure with phonopy under quasi-harmonic approxiamtion. - - Args: - structure: Pymatgen structure. - - Returns: - { - "scale_factors": list of scale factors of lattice constants, - "volumes": list of unit cell volumes at corresponding scale factors (in Angstrom^3), - "electronic_energies": list of electronic energies at corresponding volumes (in eV), - "temperatures": list of temperatures in ascending order (in Kelvin), - "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding - temperatures (in Kelvin^-1), - "gibbs_free_energies": list of Gibbs free energies at corresponding temperatures (in eV), - "bulk_modulus_P": list of bulk modulus at constant presure at corresponding temperatures (in GPa), - "heat_capacity_P": list of heat capacities at constant pressure at corresponding temperatures (in J/K/mol), - "gruneisen_parameters": list of Gruneisen parameters at corresponding temperatures, - } - """ - volumes = [] - electronic_energies = [] - temperatures = np.arange(self.t_min, self.t_max+self.t_step, self.t_step) - - free_energies = [] - entropies = [] - heat_capacities = [] - - if self.relax_structure: - relaxer = RelaxCalc( - self.calculator, fmax=self.fmax, optimizer=self.optimizer, **(self.relax_calc_kwargs or {}) - ) - structure = relaxer.calc(structure)["final_structure"] - - for scale_factor in self.scale_factors: - struct = copy.deepcopy(structure) - struct.scale_lattice(struct.volume*scale_factor**3) - - static_calc = RelaxCalc( - self.calculator, relax_atoms=False, relax_cell=False) - volumes.append(struct.volume) - electronic_energies.append(static_calc.calc(struct)["energy"]) - - phonon_calc = PhononCalc( - self.calculator, t_step=self.t_step, t_max=self.t_max, t_min=self.t_min, relax_structure=False, write_phonon=False, **(self.phonon_calc_kwargs or {}) - ) - thermal_properties = phonon_calc.calc(struct)["thermal_properties"] - free_energies.append(thermal_properties["free_energy"]) - entropies.append(thermal_properties["entropy"]) - heat_capacities.append(thermal_properties["heat_capacity"]) - - qha = PhonopyQHA( - volumes=volumes, - electronic_energies=electronic_energies, - temperatures=temperatures, - free_energy=np.transpose(free_energies), - entropy=np.transpose(entropies), - cv=np.transpose(heat_capacities), - eos=self.eos, - t_max=self.t_max - ) - - if self.write_helmholtz_volume: - qha.write_helmholtz_volume(filename=self.write_helmholtz_volume) - if self.write_volume_temperature: - qha.write_volume_temperature(filename=self.write_volume_temperature) - if self.write_thermal_expansion: - qha.write_thermal_expansion(filename=self.write_thermal_expansion) - if self.write_gibbs_temperature: - qha.write_gibbs_temperature(filename=self.write_gibbs_temperature) - if self.write_bulk_modulus_temperature: - qha.write_bulk_modulus_temperature(filename=self.write_bulk_modulus_temperature) - if self.write_heat_capacity_P_numerical: - qha.write_heat_capacity_P_numerical(filename=self.write_heat_capacity_P_numerical) - if self.write_heat_capacity_P_polyfit: - qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_P_polyfit) - if self.write_gruneisen_temperature: - qha.write_gruneisen_temperature(filename=self.write_gruneisen_temperature) - - return {"scale_factors": self.scale_factors, - "volumes": volumes, - "electronic_energies": electronic_energies, - "temperatures": temperatures, - "thermal_expansion_coefficients": qha.thermal_expansion, - "gibbs_free_energies": qha.gibbs_temperature, - "bulk_modulus_P": qha.bulk_modulus_temperature, - "heat_capacity_P": qha.heat_capacity_P_polyfit, - "gruneisen_parameters": qha.gruneisen_temperature, - } From ac6015992bb74c43be11c1855739c651877ec56c Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Sun, 11 Aug 2024 20:35:40 -0700 Subject: [PATCH 04/14] fix the ruff Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- matcalc/qha.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/matcalc/qha.py b/matcalc/qha.py index 8d7492f..1963367 100644 --- a/matcalc/qha.py +++ b/matcalc/qha.py @@ -76,7 +76,7 @@ class QhaPhononCalc(PropCalc): write_heat_capacity_P_numerical: bool | str | Path = False write_heat_capacity_P_polyfit: bool | str | Path = False write_gruneisen_temperature: bool | str | Path = False - + def __post_init__(self) -> None: """Set default paths for where to save output files.""" # map True to canonical default path, False to "" and Path to str @@ -129,7 +129,7 @@ def calc(self, structure: Structure) -> dict: for scale_factor in self.scale_factors: struct = copy.deepcopy(structure) struct.scale_lattice(struct.volume*scale_factor**3) - + static_calc = RelaxCalc( self.calculator, relax_atoms=False, relax_cell=False) volumes.append(struct.volume) @@ -142,7 +142,7 @@ def calc(self, structure: Structure) -> dict: free_energies.append(thermal_properties["free_energy"]) entropies.append(thermal_properties["entropy"]) heat_capacities.append(thermal_properties["heat_capacity"]) - + qha = PhonopyQHA( volumes=volumes, electronic_energies=electronic_energies, @@ -170,7 +170,7 @@ def calc(self, structure: Structure) -> dict: qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_P_polyfit) if self.write_gruneisen_temperature: qha.write_gruneisen_temperature(filename=self.write_gruneisen_temperature) - + return {"scale_factors": self.scale_factors, "volumes": volumes, "electronic_energies": electronic_energies, From c064e667ce963de948a90444f0e0ecd439909bdc Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Sun, 11 Aug 2024 20:57:35 -0700 Subject: [PATCH 05/14] fix the ruff Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- matcalc/qha.py | 40 ++++++++++++++++++++-------------------- 1 file changed, 20 insertions(+), 20 deletions(-) diff --git a/matcalc/qha.py b/matcalc/qha.py index 1963367..6d54d96 100644 --- a/matcalc/qha.py +++ b/matcalc/qha.py @@ -2,8 +2,6 @@ from __future__ import annotations -import numpy as np - from dataclasses import dataclass, field from typing import TYPE_CHECKING @@ -16,6 +14,8 @@ from .phonon import PhononCalc if TYPE_CHECKING: + from pathlib import Path + from ase.calculators.calculator import Calculator from pymatgen.core import Structure @@ -31,7 +31,7 @@ class QhaPhononCalc(PropCalc): fmax (float): Max forces. This criterion is more stringent than for simple relaxation. Defaults to 0.1 (in eV/Angstrom) optimizer (str): Optimizer used for RelaxCalc. Default to "FIRE" - eos (str): Equation of state used to fit F vs V, including "vinet", "murnaghan" or + eos (str): Equation of state used to fit F vs V, including "vinet", "murnaghan" or "birch_murnaghan". Default to "vinet". relax_structure (bool): Whether to first relax the structure. Set to False if structures provided are pre-relaxed with the same calculator. @@ -42,18 +42,18 @@ class QhaPhononCalc(PropCalc): Pass string or Path for custom filename. Defaults to False. write_volume_temperature (bool | str | Path): Whether to save equilibrium volume vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_thermal_expansion (bool | str | Path): Whether to save thermal expansion vs temperature in file. + write_thermal_expansion (bool | str | Path): Whether to save thermal expansion vs temperature in file. Pass string or Path for custom filename. Defaults to False. write_gibbs_temperature (bool | str | Path): Whether to save Gibbs free energy vs temperature in file. Pass string or Path for custom filename. Defaults to False. write_bulk_modulus_temperature (bool | str | Path): Whether to save bulk modulus vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_heat_capacity_P_numerical (bool | str | Path): Whether to save heat capacity at constant pressure - by numerical difference vs temperature in file. Pass string or Path for custom filename. + write_heat_capacity_p_numerical (bool | str | Path): Whether to save heat capacity at constant pressure + by numerical difference vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_heat_capacity_P_polyfit (bool | str | Path): Whether to save heat capacity at constant pressure + write_heat_capacity_p_polyfit (bool | str | Path): Whether to save heat capacity at constant pressure by fitting vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_gruneisen_temperature (bool | str | Path): Whether to save Grueneisen parameter vs temperature in + write_gruneisen_temperature (bool | str | Path): Whether to save Grueneisen parameter vs temperature in file. Pass string or Path for custom filename. Defaults to False. """ @@ -73,8 +73,8 @@ class QhaPhononCalc(PropCalc): write_thermal_expansion: bool | str | Path = False write_gibbs_temperature: bool | str | Path = False write_bulk_modulus_temperature: bool | str | Path = False - write_heat_capacity_P_numerical: bool | str | Path = False - write_heat_capacity_P_polyfit: bool | str | Path = False + write_heat_capacity_p_numerical: bool | str | Path = False + write_heat_capacity_p_polyfit: bool | str | Path = False write_gruneisen_temperature: bool | str | Path = False def __post_init__(self) -> None: @@ -86,8 +86,8 @@ def __post_init__(self) -> None: ("write_thermal_expansion", self.write_thermal_expansion, "thermal_expansion.dat"), ("write_write_gibbs_temperature", self.write_gibbs_temperature, "gibbs_temperature.dat"), ("write_bulk_modulus_temperature", self.write_bulk_modulus_temperature, "bulk_modulus_temperature.dat"), - ("write_heat_capacity_P_numerical", self.write_heat_capacity_P_numerical, "Cp_temperature.dat"), - ("write_heat_capacity_P_polyfit", self.write_heat_capacity_P_polyfit, "Cp_temperature_polyfit.dat"), + ("write_heat_capacity_p_numerical", self.write_heat_capacity_p_numerical, "Cp_temperature.dat"), + ("write_heat_capacity_p_polyfit", self.write_heat_capacity_p_polyfit, "Cp_temperature_polyfit.dat"), ("write_gruneisen_temperature", self.write_gruneisen_temperature, "gruneisen_temperature.dat"), ): setattr(self, key, str({True: default_path, False: ""}.get(val, val))) # type: ignore[arg-type] @@ -104,7 +104,7 @@ def calc(self, structure: Structure) -> dict: "volumes": list of unit cell volumes at corresponding scale factors (in Angstrom^3), "electronic_energies": list of electronic energies at corresponding volumes (in eV), "temperatures": list of temperatures in ascending order (in Kelvin), - "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding + "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding temperatures (in Kelvin^-1), "gibbs_free_energies": list of Gibbs free energies at corresponding temperatures (in eV), "bulk_modulus_P": list of bulk modulus at constant presure at corresponding temperatures (in GPa), @@ -131,7 +131,7 @@ def calc(self, structure: Structure) -> dict: struct.scale_lattice(struct.volume*scale_factor**3) static_calc = RelaxCalc( - self.calculator, relax_atoms=False, relax_cell=False) + self.calculator, relax_atoms=False, relax_cell=False) volumes.append(struct.volume) electronic_energies.append(static_calc.calc(struct)["energy"]) @@ -164,10 +164,10 @@ def calc(self, structure: Structure) -> dict: qha.write_gibbs_temperature(filename=self.write_gibbs_temperature) if self.write_bulk_modulus_temperature: qha.write_bulk_modulus_temperature(filename=self.write_bulk_modulus_temperature) - if self.write_heat_capacity_P_numerical: - qha.write_heat_capacity_P_numerical(filename=self.write_heat_capacity_P_numerical) - if self.write_heat_capacity_P_polyfit: - qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_P_polyfit) + if self.write_heat_capacity_p_numerical: + qha.write_heat_capacity_P_numerical(filename=self.write_heat_capacity_p_numerical) + if self.write_heat_capacity_p_polyfit: + qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_p_polyfit) if self.write_gruneisen_temperature: qha.write_gruneisen_temperature(filename=self.write_gruneisen_temperature) @@ -177,7 +177,7 @@ def calc(self, structure: Structure) -> dict: "temperatures": temperatures, "thermal_expansion_coefficients": qha.thermal_expansion, "gibbs_free_energies": qha.gibbs_temperature, - "bulk_modulus_P": qha.bulk_modulus_temperature, - "heat_capacity_P": qha.heat_capacity_P_polyfit, + "bulk_modulus": qha.bulk_modulus_temperature, + "heat_capacity": qha.heat_capacity_p_polyfit, "gruneisen_parameters": qha.gruneisen_temperature, } From 20534f27e2cae2e5e4faf39d272b5644ff1530ca Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Thu, 29 Aug 2024 21:44:44 -0700 Subject: [PATCH 06/14] fix the ruff Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- matcalc/qha.py | 132 ++++++++++++++++++++++++++++--------------------- 1 file changed, 77 insertions(+), 55 deletions(-) diff --git a/matcalc/qha.py b/matcalc/qha.py index 6d54d96..6b5c58d 100644 --- a/matcalc/qha.py +++ b/matcalc/qha.py @@ -2,58 +2,58 @@ from __future__ import annotations +import copy from dataclasses import dataclass, field from typing import TYPE_CHECKING import numpy as np -import copy from phonopy import PhonopyQHA from .base import PropCalc -from .relaxation import RelaxCalc from .phonon import PhononCalc +from .relaxation import RelaxCalc if TYPE_CHECKING: from pathlib import Path - from ase.calculators.calculator import Calculator + from numpy.typing import ArrayLike from pymatgen.core import Structure @dataclass -class QhaPhononCalc(PropCalc): +class QHACalc(PropCalc): """Calculator for phonon properties under quasi-harmonic approxiamtion. Args: calculator (Calculator): ASE Calculator to use. - t_step (float): Temperature step. Defaults to 10 (in Kelvin) - t_max (float): Max temperature (in Kelvin). Defaults to 1000 (in Kelvin) - t_min (float): Min temperature (in Kelvin). Defaults to 0 (in Kelvin) + t_step (float): Temperature step. Defaults to 10 (in Kelvin). + t_max (float): Max temperature (in Kelvin). Defaults to 1000 (in Kelvin). + t_min (float): Min temperature (in Kelvin). Defaults to 0 (in Kelvin). fmax (float): Max forces. This criterion is more stringent than for simple relaxation. - Defaults to 0.1 (in eV/Angstrom) - optimizer (str): Optimizer used for RelaxCalc. Default to "FIRE" - eos (str): Equation of state used to fit F vs V, including "vinet", "murnaghan" or + Defaults to 0.1 (in eV/Angstrom). + optimizer (str): Optimizer used for RelaxCalc. Default to "FIRE". + eos (str): Equation of state used to fit F vs V, including "vinet", "murnaghan" or "birch_murnaghan". Default to "vinet". relax_structure (bool): Whether to first relax the structure. Set to False if structures provided are pre-relaxed with the same calculator. relax_calc_kwargs (dict): Arguments to be passed to the RelaxCalc, if relax_structure is True. phonon_calc_kwargs (dict): Arguments to be passed to the PhononCalc. - scale_factors (list[float]): Factors to scale the lattice constants of the structure. + scale_factors (ArrayLike): Factors to scale the lattice constants of the structure. write_helmholtz_volume (bool | str | Path): Whether to save Helmholtz free energy vs volume in file. Pass string or Path for custom filename. Defaults to False. write_volume_temperature (bool | str | Path): Whether to save equilibrium volume vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_thermal_expansion (bool | str | Path): Whether to save thermal expansion vs temperature in file. + write_thermal_expansion (bool | str | Path): Whether to save thermal expansion vs temperature in file. Pass string or Path for custom filename. Defaults to False. write_gibbs_temperature (bool | str | Path): Whether to save Gibbs free energy vs temperature in file. Pass string or Path for custom filename. Defaults to False. write_bulk_modulus_temperature (bool | str | Path): Whether to save bulk modulus vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_heat_capacity_p_numerical (bool | str | Path): Whether to save heat capacity at constant pressure - by numerical difference vs temperature in file. Pass string or Path for custom filename. + write_heat_capacity_P_numerical (bool | str | Path): Whether to save heat capacity at constant pressure + by numerical difference vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_heat_capacity_p_polyfit (bool | str | Path): Whether to save heat capacity at constant pressure + write_heat_capacity_P_polyfit (bool | str | Path): Whether to save heat capacity at constant pressure by fitting vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_gruneisen_temperature (bool | str | Path): Whether to save Grueneisen parameter vs temperature in + write_gruneisen_temperature (bool | str | Path): Whether to save Grueneisen parameter vs temperature in file. Pass string or Path for custom filename. Defaults to False. """ @@ -67,14 +67,14 @@ class QhaPhononCalc(PropCalc): relax_structure: bool = True relax_calc_kwargs: dict | None = None phonon_calc_kwargs: dict | None = None - scale_factors: list[float] = field(default_factory=lambda: np.arange(0.95, 1.05, 0.01).tolist()) + scale_factors: ArrayLike = np.arange(0.95, 1.06, 0.01) write_helmholtz_volume: bool | str | Path = False write_volume_temperature: bool | str | Path = False write_thermal_expansion: bool | str | Path = False write_gibbs_temperature: bool | str | Path = False write_bulk_modulus_temperature: bool | str | Path = False - write_heat_capacity_p_numerical: bool | str | Path = False - write_heat_capacity_p_polyfit: bool | str | Path = False + write_heat_capacity_P_numerical: bool | str | Path = False + write_heat_capacity_P_polyfit: bool | str | Path = False write_gruneisen_temperature: bool | str | Path = False def __post_init__(self) -> None: @@ -86,8 +86,8 @@ def __post_init__(self) -> None: ("write_thermal_expansion", self.write_thermal_expansion, "thermal_expansion.dat"), ("write_write_gibbs_temperature", self.write_gibbs_temperature, "gibbs_temperature.dat"), ("write_bulk_modulus_temperature", self.write_bulk_modulus_temperature, "bulk_modulus_temperature.dat"), - ("write_heat_capacity_p_numerical", self.write_heat_capacity_p_numerical, "Cp_temperature.dat"), - ("write_heat_capacity_p_polyfit", self.write_heat_capacity_p_polyfit, "Cp_temperature_polyfit.dat"), + ("write_heat_capacity_P_numerical", self.write_heat_capacity_P_numerical, "Cp_temperature.dat"), + ("write_heat_capacity_P_polyfit", self.write_heat_capacity_P_polyfit, "Cp_temperature_polyfit.dat"), ("write_gruneisen_temperature", self.write_gruneisen_temperature, "gruneisen_temperature.dat"), ): setattr(self, key, str({True: default_path, False: ""}.get(val, val))) # type: ignore[arg-type] @@ -100,11 +100,12 @@ def calc(self, structure: Structure) -> dict: Returns: { + "qha": Phonopy.qha object, "scale_factors": list of scale factors of lattice constants, "volumes": list of unit cell volumes at corresponding scale factors (in Angstrom^3), "electronic_energies": list of electronic energies at corresponding volumes (in eV), "temperatures": list of temperatures in ascending order (in Kelvin), - "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding + "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding temperatures (in Kelvin^-1), "gibbs_free_energies": list of Gibbs free energies at corresponding temperatures (in eV), "bulk_modulus_P": list of bulk modulus at constant presure at corresponding temperatures (in GPa), @@ -112,13 +113,6 @@ def calc(self, structure: Structure) -> dict: "gruneisen_parameters": list of Gruneisen parameters at corresponding temperatures, } """ - volumes = [] - electronic_energies = [] - temperatures = np.arange(self.t_min, self.t_max+self.t_step, self.t_step) - - free_energies = [] - entropies = [] - heat_capacities = [] if self.relax_structure: relaxer = RelaxCalc( @@ -126,27 +120,51 @@ def calc(self, structure: Structure) -> dict: ) structure = relaxer.calc(structure)["final_structure"] - for scale_factor in self.scale_factors: - struct = copy.deepcopy(structure) - struct.scale_lattice(struct.volume*scale_factor**3) + volumes, electronic_energies, free_energies, entropies, heat_capacities = self._collect_properties(structure) - static_calc = RelaxCalc( - self.calculator, relax_atoms=False, relax_cell=False) - volumes.append(struct.volume) - electronic_energies.append(static_calc.calc(struct)["energy"]) + qha = self._create_qha(volumes, electronic_energies, free_energies, entropies, heat_capacities) - phonon_calc = PhononCalc( - self.calculator, t_step=self.t_step, t_max=self.t_max, t_min=self.t_min, relax_structure=False, write_phonon=False, **(self.phonon_calc_kwargs or {}) - ) - thermal_properties = phonon_calc.calc(struct)["thermal_properties"] + self._write_output_files(qha) + + return self._generate_output_dict(qha, volumes, electronic_energies) + + def _collect_properties(self, structure: Structure): + volumes = [] + electronic_energies = [] + free_energies = [] + entropies = [] + heat_capacities = [] + for scale_factor in self.scale_factors: + struct = self._scale_structure(structure, scale_factor) + volumes.append(struct.volume) + electronic_energies.append(self._calculate_energy(struct)) + thermal_properties = self._calculate_thermal_properties(struct) free_energies.append(thermal_properties["free_energy"]) entropies.append(thermal_properties["entropy"]) heat_capacities.append(thermal_properties["heat_capacity"]) + return volumes, electronic_energies, free_energies, entropies, heat_capacities + + def _scale_structure(self, structure: Structure, scale_factor: float) -> Structure: + struct = copy.deepcopy(structure) + struct.scale_lattice(struct.volume * scale_factor ** 3) + return struct + + def _calculate_energy(self, structure: Structure) -> float: + static_calc = RelaxCalc(self.calculator, relax_atoms=False, relax_cell=False) + return static_calc.calc(structure)["energy"] + + def _calculate_thermal_properties(self, structure: Structure) -> dict: + phonon_calc = PhononCalc( + self.calculator, t_step=self.t_step, t_max=self.t_max, t_min=self.t_min, + relax_structure=False, write_phonon=False, **(self.phonon_calc_kwargs or {}) + ) + return phonon_calc.calc(structure)["thermal_properties"] - qha = PhonopyQHA( + def _create_qha(self, volumes, electronic_energies, free_energies, entropies, heat_capacities) -> PhonopyQHA: + return PhonopyQHA( volumes=volumes, electronic_energies=electronic_energies, - temperatures=temperatures, + temperatures=np.arange(self.t_min, self.t_max + self.t_step, self.t_step), free_energy=np.transpose(free_energies), entropy=np.transpose(entropies), cv=np.transpose(heat_capacities), @@ -154,6 +172,7 @@ def calc(self, structure: Structure) -> dict: t_max=self.t_max ) + def _write_output_files(self, qha: PhonopyQHA) -> None: if self.write_helmholtz_volume: qha.write_helmholtz_volume(filename=self.write_helmholtz_volume) if self.write_volume_temperature: @@ -164,20 +183,23 @@ def calc(self, structure: Structure) -> dict: qha.write_gibbs_temperature(filename=self.write_gibbs_temperature) if self.write_bulk_modulus_temperature: qha.write_bulk_modulus_temperature(filename=self.write_bulk_modulus_temperature) - if self.write_heat_capacity_p_numerical: - qha.write_heat_capacity_P_numerical(filename=self.write_heat_capacity_p_numerical) - if self.write_heat_capacity_p_polyfit: - qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_p_polyfit) + if self.write_heat_capacity_P_numerical: + qha.write_heat_capacity_P_numerical(filename=self.write_heat_capacity_P_numerical) + if self.write_heat_capacity_P_polyfit: + qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_P_polyfit) if self.write_gruneisen_temperature: qha.write_gruneisen_temperature(filename=self.write_gruneisen_temperature) - return {"scale_factors": self.scale_factors, - "volumes": volumes, - "electronic_energies": electronic_energies, - "temperatures": temperatures, - "thermal_expansion_coefficients": qha.thermal_expansion, - "gibbs_free_energies": qha.gibbs_temperature, - "bulk_modulus": qha.bulk_modulus_temperature, - "heat_capacity": qha.heat_capacity_p_polyfit, - "gruneisen_parameters": qha.gruneisen_temperature, + def _generate_output_dict(self, qha: PhonopyQHA, volumes, electronic_energies) -> dict: + return { + "qha": qha, + "scale_factors": self.scale_factors, + "volumes": volumes, + "electronic_energies": electronic_energies, + "temperatures": qha.temperatures, + "thermal_expansion_coefficients": qha.thermal_expansion, + "gibbs_free_energies": qha.gibbs_temperature, + "bulk_modulus_P": qha.bulk_modulus_temperature, + "heat_capacity_P": qha.heat_capacity_P_polyfit, + "gruneisen_parameters": qha.gruneisen_temperature, } From ca2da298210787ce7deded01c0096e9a507277dd Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Thu, 29 Aug 2024 23:19:02 -0700 Subject: [PATCH 07/14] fix the ruff Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- matcalc/qha.py | 130 +++++++++++++++++++++++++++++++++++++++---------- 1 file changed, 103 insertions(+), 27 deletions(-) diff --git a/matcalc/qha.py b/matcalc/qha.py index 6b5c58d..0cb3e41 100644 --- a/matcalc/qha.py +++ b/matcalc/qha.py @@ -1,9 +1,9 @@ -"""Calculator for phonon properties under quasi-harmonic approxiamtion.""" +"""Calculator for phonon properties under quasi-harmonic approximation.""" from __future__ import annotations import copy -from dataclasses import dataclass, field +from dataclasses import dataclass from typing import TYPE_CHECKING import numpy as np @@ -14,9 +14,10 @@ from .relaxation import RelaxCalc if TYPE_CHECKING: + from collections.abc import Sequence from pathlib import Path + from ase.calculators.calculator import Calculator - from numpy.typing import ArrayLike from pymatgen.core import Structure @dataclass @@ -31,29 +32,29 @@ class QHACalc(PropCalc): fmax (float): Max forces. This criterion is more stringent than for simple relaxation. Defaults to 0.1 (in eV/Angstrom). optimizer (str): Optimizer used for RelaxCalc. Default to "FIRE". - eos (str): Equation of state used to fit F vs V, including "vinet", "murnaghan" or + eos (str): Equation of state used to fit F vs V, including "vinet", "murnaghan" or "birch_murnaghan". Default to "vinet". relax_structure (bool): Whether to first relax the structure. Set to False if structures provided are pre-relaxed with the same calculator. relax_calc_kwargs (dict): Arguments to be passed to the RelaxCalc, if relax_structure is True. phonon_calc_kwargs (dict): Arguments to be passed to the PhononCalc. - scale_factors (ArrayLike): Factors to scale the lattice constants of the structure. + scale_factors (Sequence[float]): Factors to scale the lattice constants of the structure. write_helmholtz_volume (bool | str | Path): Whether to save Helmholtz free energy vs volume in file. Pass string or Path for custom filename. Defaults to False. write_volume_temperature (bool | str | Path): Whether to save equilibrium volume vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_thermal_expansion (bool | str | Path): Whether to save thermal expansion vs temperature in file. + write_thermal_expansion (bool | str | Path): Whether to save thermal expansion vs temperature in file. Pass string or Path for custom filename. Defaults to False. write_gibbs_temperature (bool | str | Path): Whether to save Gibbs free energy vs temperature in file. Pass string or Path for custom filename. Defaults to False. write_bulk_modulus_temperature (bool | str | Path): Whether to save bulk modulus vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_heat_capacity_P_numerical (bool | str | Path): Whether to save heat capacity at constant pressure - by numerical difference vs temperature in file. Pass string or Path for custom filename. + write_heat_capacity_p_numerical (bool | str | Path): Whether to save heat capacity at constant pressure + by numerical difference vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_heat_capacity_P_polyfit (bool | str | Path): Whether to save heat capacity at constant pressure + write_heat_capacity_p_polyfit (bool | str | Path): Whether to save heat capacity at constant pressure by fitting vs temperature in file. Pass string or Path for custom filename. Defaults to False. - write_gruneisen_temperature (bool | str | Path): Whether to save Grueneisen parameter vs temperature in + write_gruneisen_temperature (bool | str | Path): Whether to save Grueneisen parameter vs temperature in file. Pass string or Path for custom filename. Defaults to False. """ @@ -67,14 +68,14 @@ class QHACalc(PropCalc): relax_structure: bool = True relax_calc_kwargs: dict | None = None phonon_calc_kwargs: dict | None = None - scale_factors: ArrayLike = np.arange(0.95, 1.06, 0.01) + scale_factors: Sequence[float] = (0.95, 0.96, 0.97, 0.98, 0.99, 1.00, 1.01, 1.02, 1.03, 1.04, 1.05) write_helmholtz_volume: bool | str | Path = False write_volume_temperature: bool | str | Path = False write_thermal_expansion: bool | str | Path = False write_gibbs_temperature: bool | str | Path = False write_bulk_modulus_temperature: bool | str | Path = False - write_heat_capacity_P_numerical: bool | str | Path = False - write_heat_capacity_P_polyfit: bool | str | Path = False + write_heat_capacity_p_numerical: bool | str | Path = False + write_heat_capacity_p_polyfit: bool | str | Path = False write_gruneisen_temperature: bool | str | Path = False def __post_init__(self) -> None: @@ -84,10 +85,10 @@ def __post_init__(self) -> None: ("write_helmholtz_volume", self.write_helmholtz_volume, "helmholtz_volume.dat"), ("write_volume_temperature", self.write_volume_temperature, "volume_temperature.dat"), ("write_thermal_expansion", self.write_thermal_expansion, "thermal_expansion.dat"), - ("write_write_gibbs_temperature", self.write_gibbs_temperature, "gibbs_temperature.dat"), + ("write_gibbs_temperature", self.write_gibbs_temperature, "gibbs_temperature.dat"), ("write_bulk_modulus_temperature", self.write_bulk_modulus_temperature, "bulk_modulus_temperature.dat"), - ("write_heat_capacity_P_numerical", self.write_heat_capacity_P_numerical, "Cp_temperature.dat"), - ("write_heat_capacity_P_polyfit", self.write_heat_capacity_P_polyfit, "Cp_temperature_polyfit.dat"), + ("write_heat_capacity_p_numerical", self.write_heat_capacity_p_numerical, "Cp_temperature.dat"), + ("write_heat_capacity_p_polyfit", self.write_heat_capacity_p_polyfit, "Cp_temperature_polyfit.dat"), ("write_gruneisen_temperature", self.write_gruneisen_temperature, "gruneisen_temperature.dat"), ): setattr(self, key, str({True: default_path, False: ""}.get(val, val))) # type: ignore[arg-type] @@ -115,10 +116,10 @@ def calc(self, structure: Structure) -> dict: """ if self.relax_structure: - relaxer = RelaxCalc( - self.calculator, fmax=self.fmax, optimizer=self.optimizer, **(self.relax_calc_kwargs or {}) - ) - structure = relaxer.calc(structure)["final_structure"] + relaxer = RelaxCalc( + self.calculator, fmax=self.fmax, optimizer=self.optimizer, **(self.relax_calc_kwargs or {}) + ) + structure = relaxer.calc(structure)["final_structure"] volumes, electronic_energies, free_energies, entropies, heat_capacities = self._collect_properties(structure) @@ -128,7 +129,17 @@ def calc(self, structure: Structure) -> dict: return self._generate_output_dict(qha, volumes, electronic_energies) - def _collect_properties(self, structure: Structure): + def _collect_properties(self, structure: Structure) -> tuple[list, list, list, list, list]: + """Helper to collect properties like volumes, electronic energies, and thermal properties. + + Args: + structure: Pymatgen structure for which the properties need to be calculated. + + Returns: + Tuple containing lists of volumes, electronic energies, free energies, entropies, + and heat capacities for different scale factors. + """ + volumes = [] electronic_energies = [] free_energies = [] @@ -145,22 +156,70 @@ def _collect_properties(self, structure: Structure): return volumes, electronic_energies, free_energies, entropies, heat_capacities def _scale_structure(self, structure: Structure, scale_factor: float) -> Structure: + """Helper to scale the lattice of a structure. + + Args: + structure: Pymatgen structure to be scaled. + scale_factor: Factor by which the lattice constants are scaled. + + Returns: + Pymatgen structure with scaled lattice constants. + """ + struct = copy.deepcopy(structure) struct.scale_lattice(struct.volume * scale_factor ** 3) return struct def _calculate_energy(self, structure: Structure) -> float: + """Helper to calculate the electronic energy of a structure. + + Args: + structure: Pymatgen structure for which the energy is calculated. + + Returns: + Electronic energy of the structure. + """ + static_calc = RelaxCalc(self.calculator, relax_atoms=False, relax_cell=False) return static_calc.calc(structure)["energy"] def _calculate_thermal_properties(self, structure: Structure) -> dict: + """Helper to calculate the thermal properties of a structure. + + Args: + structure: Pymatgen structure for which the thermal properties are calculated. + + Returns: + Dictionary of thermal properties containing free energies, entropies and heat capacities. + """ + phonon_calc = PhononCalc( self.calculator, t_step=self.t_step, t_max=self.t_max, t_min=self.t_min, relax_structure=False, write_phonon=False, **(self.phonon_calc_kwargs or {}) ) return phonon_calc.calc(structure)["thermal_properties"] - def _create_qha(self, volumes, electronic_energies, free_energies, entropies, heat_capacities) -> PhonopyQHA: + def _create_qha( + self, + volumes: list, + electronic_energies: list, + free_energies: list, + entropies: list, + heat_capacities: list + ) -> PhonopyQHA: + """Helper to create a PhonopyQHA object for quasi-harmonic approximation. + + Args: + volumes: List of volumes corresponding to different scale factors. + electronic_energies: List of electronic energies corresponding to different volumes. + free_energies: List of free energies corresponding to different volumes and temperatures. + entropies: List of entropies corresponding to different volumes and temperatures. + heat_capacities: List of heat capacities corresponding to different volumes and temperatures. + + Returns: + Phonopy.qha object. + """ + return PhonopyQHA( volumes=volumes, electronic_energies=electronic_energies, @@ -173,6 +232,12 @@ def _create_qha(self, volumes, electronic_energies, free_energies, entropies, he ) def _write_output_files(self, qha: PhonopyQHA) -> None: + """Helper to write various output files based on the QHA calculation. + + Args: + qha: Phonopy.qha object + """ + if self.write_helmholtz_volume: qha.write_helmholtz_volume(filename=self.write_helmholtz_volume) if self.write_volume_temperature: @@ -183,14 +248,25 @@ def _write_output_files(self, qha: PhonopyQHA) -> None: qha.write_gibbs_temperature(filename=self.write_gibbs_temperature) if self.write_bulk_modulus_temperature: qha.write_bulk_modulus_temperature(filename=self.write_bulk_modulus_temperature) - if self.write_heat_capacity_P_numerical: - qha.write_heat_capacity_P_numerical(filename=self.write_heat_capacity_P_numerical) - if self.write_heat_capacity_P_polyfit: - qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_P_polyfit) + if self.write_heat_capacity_p_numerical: + qha.write_heat_capacity_P_numerical(filename=self.write_heat_capacity_p_numerical) + if self.write_heat_capacity_p_polyfit: + qha.write_heat_capacity_P_polyfit(filename=self.write_heat_capacity_p_polyfit) if self.write_gruneisen_temperature: qha.write_gruneisen_temperature(filename=self.write_gruneisen_temperature) - def _generate_output_dict(self, qha: PhonopyQHA, volumes, electronic_energies) -> dict: + def _generate_output_dict(self, qha: PhonopyQHA, volumes: list, electronic_energies: list) -> dict: + """Helper to generate the output dictionary after QHA calculation. + + Args: + qha: Phonopy.qha object. + volumes: List of volumes corresponding to different scale factors. + electronic_energies: List of electronic energies corresponding to different volumes. + + Returns: + Dictionary containing the results of QHA calculation. + """ + return { "qha": qha, "scale_factors": self.scale_factors, From 9195132c107f9cfd94395bfba129c0f74f6db700 Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Thu, 29 Aug 2024 23:32:43 -0700 Subject: [PATCH 08/14] fix the ruff Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- matcalc/qha.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/matcalc/qha.py b/matcalc/qha.py index 0cb3e41..b636e48 100644 --- a/matcalc/qha.py +++ b/matcalc/qha.py @@ -22,7 +22,7 @@ @dataclass class QHACalc(PropCalc): - """Calculator for phonon properties under quasi-harmonic approxiamtion. + """Calculator for phonon properties under quasi-harmonic approximation. Args: calculator (Calculator): ASE Calculator to use. @@ -94,7 +94,7 @@ def __post_init__(self) -> None: setattr(self, key, str({True: default_path, False: ""}.get(val, val))) # type: ignore[arg-type] def calc(self, structure: Structure) -> dict: - """Calculates thermal properties of Pymatgen structure with phonopy under quasi-harmonic approxiamtion. + """Calculates thermal properties of Pymatgen structure with phonopy under quasi-harmonic approximation. Args: structure: Pymatgen structure. @@ -106,7 +106,7 @@ def calc(self, structure: Structure) -> dict: "volumes": list of unit cell volumes at corresponding scale factors (in Angstrom^3), "electronic_energies": list of electronic energies at corresponding volumes (in eV), "temperatures": list of temperatures in ascending order (in Kelvin), - "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding + "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding temperatures (in Kelvin^-1), "gibbs_free_energies": list of Gibbs free energies at corresponding temperatures (in eV), "bulk_modulus_P": list of bulk modulus at constant presure at corresponding temperatures (in GPa), @@ -136,8 +136,8 @@ def _collect_properties(self, structure: Structure) -> tuple[list, list, list, l structure: Pymatgen structure for which the properties need to be calculated. Returns: - Tuple containing lists of volumes, electronic energies, free energies, entropies, - and heat capacities for different scale factors. + Tuple containing lists of volumes, electronic energies, free energies, entropies, + and heat capacities for different scale factors. """ volumes = [] @@ -215,7 +215,7 @@ def _create_qha( free_energies: List of free energies corresponding to different volumes and temperatures. entropies: List of entropies corresponding to different volumes and temperatures. heat_capacities: List of heat capacities corresponding to different volumes and temperatures. - + Returns: Phonopy.qha object. """ @@ -232,7 +232,7 @@ def _create_qha( ) def _write_output_files(self, qha: PhonopyQHA) -> None: - """Helper to write various output files based on the QHA calculation. + """Helper to write various output files based on the QHA calculation. Args: qha: Phonopy.qha object @@ -256,7 +256,7 @@ def _write_output_files(self, qha: PhonopyQHA) -> None: qha.write_gruneisen_temperature(filename=self.write_gruneisen_temperature) def _generate_output_dict(self, qha: PhonopyQHA, volumes: list, electronic_energies: list) -> dict: - """Helper to generate the output dictionary after QHA calculation. + """Helper to generate the output dictionary after QHA calculation. Args: qha: Phonopy.qha object. From 9127d7ff77bde00c9fbf9f7a1c2f0d10e6fe2910 Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Thu, 29 Aug 2024 23:45:16 -0700 Subject: [PATCH 09/14] fix the ruff Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- matcalc/qha.py | 8 -------- 1 file changed, 8 deletions(-) diff --git a/matcalc/qha.py b/matcalc/qha.py index b636e48..15d2f64 100644 --- a/matcalc/qha.py +++ b/matcalc/qha.py @@ -114,7 +114,6 @@ def calc(self, structure: Structure) -> dict: "gruneisen_parameters": list of Gruneisen parameters at corresponding temperatures, } """ - if self.relax_structure: relaxer = RelaxCalc( self.calculator, fmax=self.fmax, optimizer=self.optimizer, **(self.relax_calc_kwargs or {}) @@ -139,7 +138,6 @@ def _collect_properties(self, structure: Structure) -> tuple[list, list, list, l Tuple containing lists of volumes, electronic energies, free energies, entropies, and heat capacities for different scale factors. """ - volumes = [] electronic_energies = [] free_energies = [] @@ -165,7 +163,6 @@ def _scale_structure(self, structure: Structure, scale_factor: float) -> Structu Returns: Pymatgen structure with scaled lattice constants. """ - struct = copy.deepcopy(structure) struct.scale_lattice(struct.volume * scale_factor ** 3) return struct @@ -179,7 +176,6 @@ def _calculate_energy(self, structure: Structure) -> float: Returns: Electronic energy of the structure. """ - static_calc = RelaxCalc(self.calculator, relax_atoms=False, relax_cell=False) return static_calc.calc(structure)["energy"] @@ -192,7 +188,6 @@ def _calculate_thermal_properties(self, structure: Structure) -> dict: Returns: Dictionary of thermal properties containing free energies, entropies and heat capacities. """ - phonon_calc = PhononCalc( self.calculator, t_step=self.t_step, t_max=self.t_max, t_min=self.t_min, relax_structure=False, write_phonon=False, **(self.phonon_calc_kwargs or {}) @@ -219,7 +214,6 @@ def _create_qha( Returns: Phonopy.qha object. """ - return PhonopyQHA( volumes=volumes, electronic_energies=electronic_energies, @@ -237,7 +231,6 @@ def _write_output_files(self, qha: PhonopyQHA) -> None: Args: qha: Phonopy.qha object """ - if self.write_helmholtz_volume: qha.write_helmholtz_volume(filename=self.write_helmholtz_volume) if self.write_volume_temperature: @@ -266,7 +259,6 @@ def _generate_output_dict(self, qha: PhonopyQHA, volumes: list, electronic_energ Returns: Dictionary containing the results of QHA calculation. """ - return { "qha": qha, "scale_factors": self.scale_factors, From 58cf5de7cd7dfdbba4313c542468ab27bc208306 Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Fri, 30 Aug 2024 15:23:14 -0700 Subject: [PATCH 10/14] fix the ruff Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- matcalc/qha.py | 117 +++++++++++++++++++++++++++++++++++++------------ 1 file changed, 88 insertions(+), 29 deletions(-) diff --git a/matcalc/qha.py b/matcalc/qha.py index 15d2f64..f315256 100644 --- a/matcalc/qha.py +++ b/matcalc/qha.py @@ -20,6 +20,7 @@ from ase.calculators.calculator import Calculator from pymatgen.core import Structure + @dataclass class QHACalc(PropCalc): """Calculator for phonon properties under quasi-harmonic approximation. @@ -68,7 +69,19 @@ class QHACalc(PropCalc): relax_structure: bool = True relax_calc_kwargs: dict | None = None phonon_calc_kwargs: dict | None = None - scale_factors: Sequence[float] = (0.95, 0.96, 0.97, 0.98, 0.99, 1.00, 1.01, 1.02, 1.03, 1.04, 1.05) + scale_factors: Sequence[float] = ( + 0.95, + 0.96, + 0.97, + 0.98, + 0.99, + 1.00, + 1.01, + 1.02, + 1.03, + 1.04, + 1.05, + ) write_helmholtz_volume: bool | str | Path = False write_volume_temperature: bool | str | Path = False write_thermal_expansion: bool | str | Path = False @@ -82,14 +95,46 @@ def __post_init__(self) -> None: """Set default paths for where to save output files.""" # map True to canonical default path, False to "" and Path to str for key, val, default_path in ( - ("write_helmholtz_volume", self.write_helmholtz_volume, "helmholtz_volume.dat"), - ("write_volume_temperature", self.write_volume_temperature, "volume_temperature.dat"), - ("write_thermal_expansion", self.write_thermal_expansion, "thermal_expansion.dat"), - ("write_gibbs_temperature", self.write_gibbs_temperature, "gibbs_temperature.dat"), - ("write_bulk_modulus_temperature", self.write_bulk_modulus_temperature, "bulk_modulus_temperature.dat"), - ("write_heat_capacity_p_numerical", self.write_heat_capacity_p_numerical, "Cp_temperature.dat"), - ("write_heat_capacity_p_polyfit", self.write_heat_capacity_p_polyfit, "Cp_temperature_polyfit.dat"), - ("write_gruneisen_temperature", self.write_gruneisen_temperature, "gruneisen_temperature.dat"), + ( + "write_helmholtz_volume", + self.write_helmholtz_volume, + "helmholtz_volume.dat", + ), + ( + "write_volume_temperature", + self.write_volume_temperature, + "volume_temperature.dat", + ), + ( + "write_thermal_expansion", + self.write_thermal_expansion, + "thermal_expansion.dat", + ), + ( + "write_gibbs_temperature", + self.write_gibbs_temperature, + "gibbs_temperature.dat", + ), + ( + "write_bulk_modulus_temperature", + self.write_bulk_modulus_temperature, + "bulk_modulus_temperature.dat", + ), + ( + "write_heat_capacity_p_numerical", + self.write_heat_capacity_p_numerical, + "Cp_temperature.dat", + ), + ( + "write_heat_capacity_p_polyfit", + self.write_heat_capacity_p_polyfit, + "Cp_temperature_polyfit.dat", + ), + ( + "write_gruneisen_temperature", + self.write_gruneisen_temperature, + "gruneisen_temperature.dat", + ), ): setattr(self, key, str({True: default_path, False: ""}.get(val, val))) # type: ignore[arg-type] @@ -102,31 +147,35 @@ def calc(self, structure: Structure) -> dict: Returns: { "qha": Phonopy.qha object, - "scale_factors": list of scale factors of lattice constants, - "volumes": list of unit cell volumes at corresponding scale factors (in Angstrom^3), - "electronic_energies": list of electronic energies at corresponding volumes (in eV), - "temperatures": list of temperatures in ascending order (in Kelvin), - "thermal_expansion_coefficients": list of volumetric thermal expansion coefficients at corresponding + "scale_factors": List of scale factors of lattice constants, + "volumes": List of unit cell volumes at corresponding scale factors (in Angstrom^3), + "electronic_energies": List of electronic energies at corresponding volumes (in eV), + "temperatures": List of temperatures in ascending order (in Kelvin), + "thermal_expansion_coefficients": List of volumetric thermal expansion coefficients at corresponding temperatures (in Kelvin^-1), - "gibbs_free_energies": list of Gibbs free energies at corresponding temperatures (in eV), - "bulk_modulus_P": list of bulk modulus at constant presure at corresponding temperatures (in GPa), - "heat_capacity_P": list of heat capacities at constant pressure at corresponding temperatures (in J/K/mol), - "gruneisen_parameters": list of Gruneisen parameters at corresponding temperatures, + "gibbs_free_energies": List of Gibbs free energies at corresponding temperatures (in eV), + "bulk_modulus_P": List of bulk modulus at constant presure at corresponding temperatures (in GPa), + "heat_capacity_P": List of heat capacities at constant pressure at corresponding temperatures (in J/K/mol), + "gruneisen_parameters": List of Gruneisen parameters at corresponding temperatures, } """ if self.relax_structure: relaxer = RelaxCalc( - self.calculator, fmax=self.fmax, optimizer=self.optimizer, **(self.relax_calc_kwargs or {}) + self.calculator, + fmax=self.fmax, + optimizer=self.optimizer, + **(self.relax_calc_kwargs or {}), ) structure = relaxer.calc(structure)["final_structure"] + temperatures = np.arange(self.t_min, self.t_max + self.t_step, self.t_step) volumes, electronic_energies, free_energies, entropies, heat_capacities = self._collect_properties(structure) - qha = self._create_qha(volumes, electronic_energies, free_energies, entropies, heat_capacities) + qha = self._create_qha(volumes, electronic_energies, temperatures, free_energies, entropies, heat_capacities) self._write_output_files(qha) - return self._generate_output_dict(qha, volumes, electronic_energies) + return self._generate_output_dict(qha, volumes, electronic_energies, temperatures) def _collect_properties(self, structure: Structure) -> tuple[list, list, list, list, list]: """Helper to collect properties like volumes, electronic energies, and thermal properties. @@ -164,7 +213,7 @@ def _scale_structure(self, structure: Structure, scale_factor: float) -> Structu Pymatgen structure with scaled lattice constants. """ struct = copy.deepcopy(structure) - struct.scale_lattice(struct.volume * scale_factor ** 3) + struct.scale_lattice(struct.volume * scale_factor**3) return struct def _calculate_energy(self, structure: Structure) -> float: @@ -189,8 +238,13 @@ def _calculate_thermal_properties(self, structure: Structure) -> dict: Dictionary of thermal properties containing free energies, entropies and heat capacities. """ phonon_calc = PhononCalc( - self.calculator, t_step=self.t_step, t_max=self.t_max, t_min=self.t_min, - relax_structure=False, write_phonon=False, **(self.phonon_calc_kwargs or {}) + self.calculator, + t_step=self.t_step, + t_max=self.t_max, + t_min=self.t_min, + relax_structure=False, + write_phonon=False, + **(self.phonon_calc_kwargs or {}), ) return phonon_calc.calc(structure)["thermal_properties"] @@ -198,15 +252,17 @@ def _create_qha( self, volumes: list, electronic_energies: list, + temperatures: list, free_energies: list, entropies: list, - heat_capacities: list + heat_capacities: list, ) -> PhonopyQHA: """Helper to create a PhonopyQHA object for quasi-harmonic approximation. Args: volumes: List of volumes corresponding to different scale factors. electronic_energies: List of electronic energies corresponding to different volumes. + temperatures: List of temperatures in ascending order (in Kelvin). free_energies: List of free energies corresponding to different volumes and temperatures. entropies: List of entropies corresponding to different volumes and temperatures. heat_capacities: List of heat capacities corresponding to different volumes and temperatures. @@ -217,12 +273,12 @@ def _create_qha( return PhonopyQHA( volumes=volumes, electronic_energies=electronic_energies, - temperatures=np.arange(self.t_min, self.t_max + self.t_step, self.t_step), + temperatures=temperatures, free_energy=np.transpose(free_energies), entropy=np.transpose(entropies), cv=np.transpose(heat_capacities), eos=self.eos, - t_max=self.t_max + t_max=self.t_max, ) def _write_output_files(self, qha: PhonopyQHA) -> None: @@ -248,13 +304,16 @@ def _write_output_files(self, qha: PhonopyQHA) -> None: if self.write_gruneisen_temperature: qha.write_gruneisen_temperature(filename=self.write_gruneisen_temperature) - def _generate_output_dict(self, qha: PhonopyQHA, volumes: list, electronic_energies: list) -> dict: + def _generate_output_dict( + self, qha: PhonopyQHA, volumes: list, electronic_energies: list, temperatures: list + ) -> dict: """Helper to generate the output dictionary after QHA calculation. Args: qha: Phonopy.qha object. volumes: List of volumes corresponding to different scale factors. electronic_energies: List of electronic energies corresponding to different volumes. + temperatures: List of temperatures in ascending order (in Kelvin). Returns: Dictionary containing the results of QHA calculation. @@ -264,7 +323,7 @@ def _generate_output_dict(self, qha: PhonopyQHA, volumes: list, electronic_energ "scale_factors": self.scale_factors, "volumes": volumes, "electronic_energies": electronic_energies, - "temperatures": qha.temperatures, + "temperatures": temperatures, "thermal_expansion_coefficients": qha.thermal_expansion, "gibbs_free_energies": qha.gibbs_temperature, "bulk_modulus_P": qha.bulk_modulus_temperature, From 4afa6864d4866556ba75f7aba9c600658597fcda Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Sat, 7 Sep 2024 20:31:19 -0700 Subject: [PATCH 11/14] Add unittests for qha.py Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- tests/test_qha.py | 120 ++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 120 insertions(+) create mode 100644 tests/test_qha.py diff --git a/tests/test_qha.py b/tests/test_qha.py new file mode 100644 index 0000000..5102652 --- /dev/null +++ b/tests/test_qha.py @@ -0,0 +1,120 @@ +"""Tests for QHACalc class""" + +from __future__ import annotations + +import inspect +import os +from typing import TYPE_CHECKING + +import pytest + +from matcalc.qha import QHACalc + +if TYPE_CHECKING: + from pathlib import Path + + from matgl.ext.ase import M3GNetCalculator + from pymatgen.core import Structure + + +@pytest.mark.parametrize( + ( + "helmholtz_file", + "volume_temp_file", + "thermal_exp_file", + "gibbs_file", + "bulk_mod_file", + "cp_numerical_file", + "cp_polyfit_file", + "gruneisen_file", + ), + [ + ("", "", "", "", "", "", "", ""), + ( + "helmholtz.dat", + "volume_temp.dat", + "thermal_expansion.dat", + "gibbs.dat", + "bulk_mod.dat", + "cp_numerical.dat", + "cp_polyfit.dat", + "gruneisen.dat", + ), + ], +) +def test_qha_calc( + Li2O: Structure, + M3GNetCalc: M3GNetCalculator, + tmp_path: Path, + helmholtz_file: str, + volume_temp_file: str, + thermal_exp_file: str, + gibbs_file: str, + bulk_mod_file: str, + cp_numerical_file: str, + cp_polyfit_file: str, + gruneisen_file: str, +) -> None: + """Tests for QHACalc class.""" + # Note that the fmax is probably too high. This is for testing purposes only. + + # change dir to tmp_path + helmholtz_volume = tmp_path / helmholtz_file if helmholtz_file else False + volume_temperature = tmp_path / volume_temp_file if volume_temp_file else False + thermal_expansion = tmp_path / thermal_exp_file if thermal_exp_file else False + gibbs_temperature = tmp_path / gibbs_file if gibbs_file else False + bulk_modulus_temperature = tmp_path / bulk_mod_file if bulk_mod_file else False + cp_numerical = tmp_path / cp_numerical_file if cp_numerical_file else False + cp_polyfit = tmp_path / cp_polyfit_file if cp_polyfit_file else False + gruneisen_temperature = tmp_path / gruneisen_file if gruneisen_file else False + + write_kwargs = { + "write_helmholtz_volume": helmholtz_volume, + "write_volume_temperature": volume_temperature, + "write_thermal_expansion": thermal_expansion, + "write_gibbs_temperature": gibbs_temperature, + "write_bulk_modulus_temperature": bulk_modulus_temperature, + "write_heat_capacity_P_numerical": cp_numerical, + "write_heat_capacity_P_polyfit": cp_polyfit, + "write_gruneisen_temperature": gruneisen_temperature, + } + + # Initialize QHACalc + qha_calc = QHACalc( + calculator=M3GNetCalc, + t_step=50, + t_max=1000, + scale_factors=[0.97, 0.98, 0.99, 1.00, 1.01, 1.02, 1.03], + **write_kwargs, # type: ignore[arg-type] + ) + + result = qha_calc.calc(Li2O) + + # Test values corresponding to different scale factors + assert result["volumes"] == pytest.approx( + [23.07207, 23.79302, 24.52884, 25.27967, 26.04567, 26.82699, 27.62378], + rel=1e-3, + ) + assert result["electronic_energies"] == pytest.approx( + [-14.12540, -14.15696, -14.17359, -14.17671, -14.17039, -14.15557, -14.13327], + rel=1e-3, + ) + + # Test values at 300 K + ind = result["temperatures"].tolist().index(300) + assert result["thermal_expansion_coefficients"][ind] == pytest.approx(1.03973e-04, rel=1e-2) + assert result["gibbs_free_energies"][ind] == pytest.approx(-14.04472, rel=1e-2) + assert result["bulk_modulus_P"][ind] == pytest.approx(54.25954, rel=1e-2) + assert result["heat_capacity_P"][ind] == pytest.approx(62.27455, rel=1e-2) + assert result["gruneisen_parameters"][ind] == pytest.approx(1.49330, rel=1e-2) + + qha_calc_params = inspect.signature(QHACalc).parameters + # get all keywords starting with write_ and their default values + file_write_defaults = {key: val.default for key, val in qha_calc_params.items() if key.startswith("write_")} + assert len(file_write_defaults) == 8 + + for keyword, default_path in file_write_defaults.items(): + if instance_val := write_kwargs[keyword]: + assert os.path.isfile(str(instance_val)) + elif not default_path and not instance_val: + assert not os.path.isfile(default_path) From cee4998022de139528f60fdd593d17b4e2082503 Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Sat, 7 Sep 2024 20:39:22 -0700 Subject: [PATCH 12/14] fix the pytest Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- tests/test_qha.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/test_qha.py b/tests/test_qha.py index 5102652..f347f0d 100644 --- a/tests/test_qha.py +++ b/tests/test_qha.py @@ -74,8 +74,8 @@ def test_qha_calc( "write_thermal_expansion": thermal_expansion, "write_gibbs_temperature": gibbs_temperature, "write_bulk_modulus_temperature": bulk_modulus_temperature, - "write_heat_capacity_P_numerical": cp_numerical, - "write_heat_capacity_P_polyfit": cp_polyfit, + "write_heat_capacity_p_numerical": cp_numerical, + "write_heat_capacity_p_polyfit": cp_polyfit, "write_gruneisen_temperature": gruneisen_temperature, } From b8988c52cb37fe0df108458feca62d8710bdaf35 Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Thu, 12 Sep 2024 16:32:29 -0700 Subject: [PATCH 13/14] Update requirements.txt Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- requirements.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/requirements.txt b/requirements.txt index 439c065..a703ec7 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,7 +1,7 @@ pymatgen==2024.2.20 ase==3.23.0 joblib==1.3.1 -phonopy==2.20.0 +phonopy==2.27.0 scikit-learn>=1.3.0 seekpath==2.1.0 numpy==1.26.4 From 916ebaedc4d891428d3c1ff1640bdac1cf1de0b6 Mon Sep 17 00:00:00 2001 From: Runze Liu <146490083+rul048@users.noreply.github.com> Date: Thu, 12 Sep 2024 17:42:15 -0700 Subject: [PATCH 14/14] Update qha.py Signed-off-by: Runze Liu <146490083+rul048@users.noreply.github.com> --- matcalc/qha.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/matcalc/qha.py b/matcalc/qha.py index f315256..a3984a4 100644 --- a/matcalc/qha.py +++ b/matcalc/qha.py @@ -2,7 +2,6 @@ from __future__ import annotations -import copy from dataclasses import dataclass from typing import TYPE_CHECKING @@ -212,8 +211,8 @@ def _scale_structure(self, structure: Structure, scale_factor: float) -> Structu Returns: Pymatgen structure with scaled lattice constants. """ - struct = copy.deepcopy(structure) - struct.scale_lattice(struct.volume * scale_factor**3) + struct = structure.copy() + struct.apply_strain(scale_factor - 1) return struct def _calculate_energy(self, structure: Structure) -> float: