From 1e347c42c01a4e926e15b910cca8964c1a0cc826 Mon Sep 17 00:00:00 2001 From: Jason Xie <48645456+jinlhr542@users.noreply.github.com> Date: Thu, 8 Feb 2024 22:47:02 +0800 Subject: [PATCH] Breaking: fix SubstrateAnalyzer film + substrate vectors not using original crystal coordinates (#3572) * Fix the problem that the film vectors and substrate vectors were not expressed in the original crystal coordinates when making the intreface This change is to help the users to be capable to track the real film vectors and substrate vectors expressed in their own original crystal coordinates. This is important when they make interface structures because they need to figure out the exact crystal vectors in along the a, b, c vectors in the supercell Signed-off-by: Jason Xie <48645456+jinlhr542@users.noreply.github.com> * pre-commit auto-fixes * breaking: make film + substrate args of generate_surface_vectors() before, method failed if called before calling calculate() and passing in film + substrate * add test_generate_surface_vectors * fix missing doc str white space --------- Signed-off-by: Jason Xie <48645456+jinlhr542@users.noreply.github.com> Co-authored-by: Janosh Riebesell --- .../analysis/interfaces/substrate_analyzer.py | 70 ++++++++++--------- pymatgen/analysis/interfaces/zsl.py | 4 +- pymatgen/analysis/piezo_sensitivity.py | 2 +- pymatgen/analysis/xas/spectrum.py | 2 +- pymatgen/io/aims/inputs.py | 2 +- pymatgen/io/lobster/inputs.py | 6 +- pymatgen/io/vasp/sets.py | 2 +- .../interfaces/test_substrate_analyzer.py | 32 +++++++-- 8 files changed, 72 insertions(+), 48 deletions(-) diff --git a/pymatgen/analysis/interfaces/substrate_analyzer.py b/pymatgen/analysis/interfaces/substrate_analyzer.py index 48e51d0d069..7e069aaea09 100644 --- a/pymatgen/analysis/interfaces/substrate_analyzer.py +++ b/pymatgen/analysis/interfaces/substrate_analyzer.py @@ -10,6 +10,8 @@ from pymatgen.core.surface import SlabGenerator, get_symmetrically_distinct_miller_indices if TYPE_CHECKING: + from numpy.typing import ArrayLike + from pymatgen.core import Structure @@ -95,11 +97,11 @@ def __init__(self, film_max_miller=1, substrate_max_miller=1, **kwargs): Initializes the substrate analyzer Args: - zslgen(ZSLGenerator): Defaults to a ZSLGenerator with standard + zslgen (ZSLGenerator): Defaults to a ZSLGenerator with standard tolerances, but can be fed one with custom tolerances - film_max_miller(int): maximum miller index to generate for film + film_max_miller (int): maximum miller index to generate for film surfaces - substrate_max_miller(int): maximum miller index to generate for + substrate_max_miller (int): maximum miller index to generate for substrate surfaces. """ self.film_max_miller = film_max_miller @@ -107,38 +109,47 @@ def __init__(self, film_max_miller=1, substrate_max_miller=1, **kwargs): self.kwargs = kwargs super().__init__(**kwargs) - def generate_surface_vectors(self, film_millers, substrate_millers): + def generate_surface_vectors( + self, film: Structure, substrate: Structure, film_millers: ArrayLike, substrate_millers: ArrayLike + ): """ Generates the film/substrate slab combinations for a set of given miller indices. Args: - film_millers(array): all miller indices to generate slabs for + film (Structure): film structure + substrate (Structure): substrate structure + film_millers (array): all miller indices to generate slabs for film - substrate_millers(array): all miller indices to generate slabs + substrate_millers (array): all miller indices to generate slabs for substrate """ vector_sets = [] - for f in film_millers: - film_slab = SlabGenerator(self.film, f, 20, 15, primitive=False).get_slab() - film_vectors = reduce_vectors(film_slab.lattice.matrix[0], film_slab.lattice.matrix[1]) + for f_miller in film_millers: + film_slab = SlabGenerator(film, f_miller, 20, 15, primitive=False).get_slab() + film_vectors = reduce_vectors( + film_slab.oriented_unit_cell.lattice.matrix[0], film_slab.oriented_unit_cell.lattice.matrix[1] + ) - for s in substrate_millers: - substrate_slab = SlabGenerator(self.substrate, s, 20, 15, primitive=False).get_slab() - substrate_vectors = reduce_vectors(substrate_slab.lattice.matrix[0], substrate_slab.lattice.matrix[1]) + for s_miller in substrate_millers: + substrate_slab = SlabGenerator(substrate, s_miller, 20, 15, primitive=False).get_slab() + substrate_vectors = reduce_vectors( + substrate_slab.oriented_unit_cell.lattice.matrix[0], + substrate_slab.oriented_unit_cell.lattice.matrix[1], + ) - vector_sets.append((film_vectors, substrate_vectors, f, s)) + vector_sets.append((film_vectors, substrate_vectors, f_miller, s_miller)) return vector_sets def calculate( self, - film, - substrate, + film: Structure, + substrate: Structure, elasticity_tensor=None, - film_millers=None, - substrate_millers=None, + film_millers: ArrayLike = None, + substrate_millers: ArrayLike = None, ground_state_energy=0, lowest=False, ): @@ -148,33 +159,28 @@ def calculate( ground state energy are provided: Args: - film(Structure): conventional standard structure for the film - substrate(Structure): conventional standard structure for the + film (Structure): conventional standard structure for the film + substrate (Structure): conventional standard structure for the substrate - elasticity_tensor(ElasticTensor): elasticity tensor for the film + elasticity_tensor (ElasticTensor): elasticity tensor for the film in the IEEE orientation - film_millers(array): film facets to consider in search as defined by + film_millers (array): film facets to consider in search as defined by miller indices - substrate_millers(array): substrate facets to consider in search as + substrate_millers (array): substrate facets to consider in search as defined by miller indices - ground_state_energy(float): ground state energy for the film - lowest(bool): only consider lowest matching area for each surface + ground_state_energy (float): ground state energy for the film + lowest (bool): only consider lowest matching area for each surface """ - self.film = film - self.substrate = substrate - # Generate miller indices if none specified for film if film_millers is None: - film_millers = sorted(get_symmetrically_distinct_miller_indices(self.film, self.film_max_miller)) + film_millers = sorted(get_symmetrically_distinct_miller_indices(film, self.film_max_miller)) # Generate miller indices if none specified for substrate if substrate_millers is None: - substrate_millers = sorted( - get_symmetrically_distinct_miller_indices(self.substrate, self.substrate_max_miller) - ) + substrate_millers = sorted(get_symmetrically_distinct_miller_indices(substrate, self.substrate_max_miller)) # Check each miller index combination - surface_vector_sets = self.generate_surface_vectors(film_millers, substrate_millers) + surface_vector_sets = self.generate_surface_vectors(film, substrate, film_millers, substrate_millers) for [ film_vectors, substrate_vectors, diff --git a/pymatgen/analysis/interfaces/zsl.py b/pymatgen/analysis/interfaces/zsl.py index 56d5e3b15df..cee9ab74194 100644 --- a/pymatgen/analysis/interfaces/zsl.py +++ b/pymatgen/analysis/interfaces/zsl.py @@ -116,8 +116,8 @@ def generate_sl_transformation_sets(self, film_area, substrate_area): lattices with a maximum area Args: - film_area(int): the unit cell area for the film - substrate_area(int): the unit cell area for the substrate + film_area (int): the unit cell area for the film + substrate_area (int): the unit cell area for the substrate Returns: transformation_sets: a set of transformation_sets defined as: diff --git a/pymatgen/analysis/piezo_sensitivity.py b/pymatgen/analysis/piezo_sensitivity.py index e7a86b3006a..1e50503df1c 100644 --- a/pymatgen/analysis/piezo_sensitivity.py +++ b/pymatgen/analysis/piezo_sensitivity.py @@ -229,7 +229,7 @@ def get_rand_IST(self, max_force=1): symmetry and the acoustic sum rule. Args: - max_force(float): maximum born effective charge value + max_force (float): maximum born effective charge value Return: InternalStrainTensor object diff --git a/pymatgen/analysis/xas/spectrum.py b/pymatgen/analysis/xas/spectrum.py index 8bab6ee6aef..93c7e670859 100644 --- a/pymatgen/analysis/xas/spectrum.py +++ b/pymatgen/analysis/xas/spectrum.py @@ -96,7 +96,7 @@ def stitch(self, other: XAS, num_samples: int = 500, mode: Literal["XAFS", "L23" Args: other: Another XAS object. - num_samples(int): Number of samples for interpolation. + num_samples (int): Number of samples for interpolation. mode("XAFS" | "L23"): Either XAFS mode for stitching XANES and EXAFS or L23 mode for stitching L2 and L3. diff --git a/pymatgen/io/aims/inputs.py b/pymatgen/io/aims/inputs.py index 68ecab11b0e..0052f1843c4 100644 --- a/pymatgen/io/aims/inputs.py +++ b/pymatgen/io/aims/inputs.py @@ -476,7 +476,7 @@ def get_aims_control_parameter_str(self, key: str, value: Any, fmt: str) -> str: fmt (str): The format string to apply to the value Returns: - The line to add to the control.in file + str: The line to add to the control.in file """ return f"{key:35s}{fmt % value}\n" diff --git a/pymatgen/io/lobster/inputs.py b/pymatgen/io/lobster/inputs.py index fb78d0bc470..95c520c3c5e 100644 --- a/pymatgen/io/lobster/inputs.py +++ b/pymatgen/io/lobster/inputs.py @@ -619,7 +619,7 @@ def _get_potcar_symbols(POTCAR_input: str) -> list: Will return the name of the species in the POTCAR. Args: - POTCAR_input(str): string to potcar file + POTCAR_input (str): string to potcar file Returns: list of the names of the species in string format @@ -662,8 +662,8 @@ def standard_calculations_from_vasp_files( Will generate Lobsterin with standard settings. Args: - POSCAR_input(str): path to POSCAR - INCAR_input(str): path to INCAR + POSCAR_input (str): path to POSCAR + INCAR_input (str): path to INCAR POTCAR_input (str): path to POTCAR dict_for_basis (dict): can be provided: it should look the following: dict_for_basis={"Fe":'3p 3d 4s 4f', "C": '2s 2p'} and will overwrite all settings from POTCAR_input diff --git a/pymatgen/io/vasp/sets.py b/pymatgen/io/vasp/sets.py index 3766db5fdfd..0b841500cd2 100644 --- a/pymatgen/io/vasp/sets.py +++ b/pymatgen/io/vasp/sets.py @@ -2162,7 +2162,7 @@ class MVLGBSet(DictSet): Class for writing a vasp input files for grain boundary calculations, slab or bulk. Args: - structure(Structure): provide the structure + structure (Structure): provide the structure k_product: Kpoint number * length for a & b directions, also for c direction in bulk calculations. Default to 40. slab_mode (bool): Defaults to False. Use default (False) for a bulk supercell. diff --git a/tests/analysis/interfaces/test_substrate_analyzer.py b/tests/analysis/interfaces/test_substrate_analyzer.py index 82650d2c60a..a842cd4dad7 100644 --- a/tests/analysis/interfaces/test_substrate_analyzer.py +++ b/tests/analysis/interfaces/test_substrate_analyzer.py @@ -1,20 +1,22 @@ from __future__ import annotations +from numpy.testing import assert_allclose + from pymatgen.analysis.elasticity.elastic import ElasticTensor from pymatgen.analysis.interfaces.substrate_analyzer import SubstrateAnalyzer from pymatgen.symmetry.analyzer import SpacegroupAnalyzer from pymatgen.util.testing import PymatgenTest +VO2 = PymatgenTest.get_structure("VO2") +TiO2 = PymatgenTest.get_structure("TiO2") -def test_substrate_analyzer_init(): - # Film VO2 - film = SpacegroupAnalyzer(PymatgenTest.get_structure("VO2"), symprec=0.1).get_conventional_standard_structure() +# Film VO2 +film = SpacegroupAnalyzer(VO2, symprec=0.1).get_conventional_standard_structure() +# Substrate TiO2 +substrate = SpacegroupAnalyzer(TiO2, symprec=0.1).get_conventional_standard_structure() - # Substrate TiO2 - substrate = SpacegroupAnalyzer( - PymatgenTest.get_structure("TiO2"), symprec=0.1 - ).get_conventional_standard_structure() +def test_substrate_analyzer_init(): film_elastic_tensor = ElasticTensor.from_voigt( [ [324.32, 187.3, 170.92, 0, 0, 0], @@ -32,3 +34,19 @@ def test_substrate_analyzer_init(): assert len(matches) == 296 for match in matches: assert isinstance(match.match_area, float) + + +def test_generate_surface_vectors(): + film_miller_indices = [(1, 0, 0)] + substrate_miller_indices = [(1, 1, 1)] + + vector_sets = SubstrateAnalyzer().generate_surface_vectors( + film, substrate, film_miller_indices, substrate_miller_indices + ) + assert len(vector_sets) == 1 + film_vectors, substrate_vectors, film_millers, substrate_millers = vector_sets[0] + + assert [film_millers] == film_miller_indices + assert [substrate_millers] == substrate_miller_indices + assert_allclose(film_vectors, [[0, 0, 3.035429], [-2.764654e-16, 4.515023, 2.764654e-16]], atol=1e-6) + assert_allclose(substrate_vectors, [[-3.766937, -1.928326, -6.328967], [3.766937, -12.307154, 0.0]], atol=1e-6)