split vasp sets into package with smaller more navigable submodules

matches atomate2 code structure with base.py, mp.py, matpes.py, etc.

Author: janosh

dev_scripts/potcar_scrambler.py

@@ -12,7 +12,7 @@
 
 from pymatgen.core import SETTINGS
 from pymatgen.io.vasp import Potcar, PotcarSingle
-from pymatgen.io.vasp.sets import _load_yaml_config
+from pymatgen.io.vasp.sets.base import _load_yaml_config
 from pymatgen.util.testing import VASP_IN_DIR
 
 if TYPE_CHECKING:

pymatgen/io/vasp/sets.py

@@ -0,0 +1,72 @@
+"""
+This module defines the VaspInputSet abstract base class and a concrete implementation for the parameters developed
+and tested by the core team of pymatgen, including the Materials Virtual Lab, Materials Project and the MIT high
+throughput project. The basic concept behind an input set is to specify a scheme to generate a consistent set of VASP
+inputs from a structure without further user intervention. This ensures comparability across runs.
+
+Read the following carefully before implementing new input sets:
+
+1. 99% of what needs to be done can be done by specifying user_incar_settings to override some of the defaults of
+   various input sets. Unless there is an extremely good reason to add a new set, **do not** add one. e.g. if you want
+   to turn the Hubbard U off, just set "LDAU": False as a user_incar_setting.
+2. All derivative input sets should inherit appropriate configurations (e.g., from MPRelaxSet), and more often than
+   not, VaspInputSet should be the superclass. Superclass delegation should be used where possible. In particular,
+   you are not supposed to implement your own as_dict or from_dict for derivative sets unless you know what you are
+   doing. Improper overriding the as_dict and from_dict protocols is the major cause of implementation headaches. If
+   you need an example, look at how the MPStaticSet is initialized.
+
+The above are recommendations. The following are **UNBREAKABLE** rules:
+
+1. All input sets must take in a structure, list of structures or None as the first argument. If None, the input set
+   should perform a stateless initialization and before any output can be written, a structure must be set.
+2. user_incar_settings, user_kpoints_settings and user_<whatever>_settings are ABSOLUTE. Any new sets you implement
+   must obey this. If a user wants to override your settings, you assume he knows what he is doing. Do not
+   magically override user supplied settings. You can issue a warning if you think the user is wrong.
+3. All input sets must save all supplied args and kwargs as instance variables. e.g. self.arg = arg and
+   self.kwargs = kwargs in the __init__. This ensures the as_dict and from_dict work correctly.
+"""
+
+from __future__ import annotations
+
+from pymatgen.io.vasp.sets.base import DictSet, VaspInputGenerator, VaspInputSet
+from pymatgen.io.vasp.sets.lobster import LobsterSet
+from pymatgen.io.vasp.sets.matpes import MatPESStaticSet
+from pymatgen.io.vasp.sets.mit import MITMDSet, MITNEBSet, MITRelaxSet
+from pymatgen.io.vasp.sets.mp import (
+    MPAbsorptionSet,
+    MPHSEBSSet,
+    MPHSERelaxSet,
+    MPMDSet,
+    MPMetalRelaxSet,
+    MPNMRSet,
+    MPNonSCFSet,
+    MPRelaxSet,
+    MPScanRelaxSet,
+    MPScanStaticSet,
+    MPSOCSet,
+    MPStaticSet,
+)
+from pymatgen.io.vasp.sets.mvl import (
+    MVLElasticSet,
+    MVLGBSet,
+    MVLGWSet,
+    MVLNPTMDSet,
+    MVLRelax52Set,
+    MVLScanRelaxSet,
+    MVLSlabSet,
+)
+from pymatgen.io.vasp.sets.utils import (
+    MODULE_DIR,
+    BadInputSetWarning,
+    _apply_incar_updates,
+    _combine_kpoints,
+    _get_ispin,
+    _get_recommended_lreal,
+    _remove_unused_incar_params,
+    auto_kspacing,
+    batch_write_input,
+    get_structure_from_prev_run,
+    get_valid_magmom_struct,
+    get_vasprun_outcar,
+    standardize_structure,
+)

pymatgen/io/vasp/sets/__init__.py

@@ -0,0 +1,292 @@
+from __future__ import annotations
+
+from collections.abc import Sequence
+from glob import glob
+from pathlib import Path
+from typing import TYPE_CHECKING, cast
+
+import numpy as np
+
+from pymatgen.analysis.structure_matcher import StructureMatcher
+from pymatgen.io.vasp.inputs import Kpoints
+from pymatgen.io.vasp.outputs import Outcar, Vasprun
+from pymatgen.symmetry.analyzer import SpacegroupAnalyzer
+from pymatgen.util.typing import Kpoint
+
+if TYPE_CHECKING:
+    from typing import Literal, Union
+
+    from pymatgen.core import Structure
+
+    UserPotcarFunctional = Union[
+        Literal["PBE", "PBE_52", "PBE_54", "LDA", "LDA_52", "LDA_54", "PW91", "LDA_US", "PW91_US"], None
+    ]
+
+
+class BadInputSetWarning(UserWarning):
+    """Warning class for bad but legal VASP inputs."""
+
+
+def _get_ispin(vasprun: Vasprun | None, outcar: Outcar | None) -> int:
+    """Get value of ISPIN depending on the magnetisation in the OUTCAR and vasprun."""
+    if outcar is not None and outcar.magnetization is not None:
+        # Turn off spin when magmom for every site is smaller than 0.02.
+        site_magmom = np.array([i["tot"] for i in outcar.magnetization])
+        return 2 if np.any(np.abs(site_magmom) > 0.02) else 1
+    if vasprun is not None:
+        return 2 if vasprun.is_spin else 1
+    return 2
+
+
+def _get_recommended_lreal(structure: Structure) -> str | bool:
+    """Get recommended LREAL flag based on the structure."""
+    return "Auto" if structure.num_sites > 16 else False
+
+
+def _combine_kpoints(*kpoints_objects: Kpoints) -> Kpoints:
+    """Combine multiple Kpoints objects."""
+    _labels: list[list[str]] = []
+    _kpoints: list[Sequence[Kpoint]] = []
+    _weights = []
+
+    for kpoints_object in filter(None, kpoints_objects):
+        if kpoints_object.style != Kpoints.supported_modes.Reciprocal:
+            raise ValueError("Can only combine kpoints with style=Kpoints.supported_modes.Reciprocal")
+        if kpoints_object.labels is None:
+            _labels.append([""] * len(kpoints_object.kpts))
+        else:
+            _labels.append(kpoints_object.labels)
+
+        _kpoints.append(kpoints_object.kpts)
+        _weights.append(kpoints_object.kpts_weights)
+
+    labels = np.concatenate(_labels).tolist()
+    kpoints = np.concatenate(_kpoints).tolist()
+    weights = np.concatenate(_weights).tolist()
+
+    return Kpoints(
+        comment="Combined k-points",
+        style=Kpoints.supported_modes.Reciprocal,
+        num_kpts=len(kpoints),
+        kpts=cast(Sequence[Kpoint], kpoints),
+        labels=labels,
+        kpts_weights=weights,
+    )
+
+
+def _apply_incar_updates(incar, updates, skip: Sequence[str] = ()) -> None:
+    """
+    Apply updates to an INCAR file.
+
+    Args:
+        incar (dict): An incar.
+        updates (dict): Updates to apply.
+        skip (list of str): Keys to skip.
+    """
+    for k, v in updates.items():
+        if k in skip:
+            continue
+
+        if v is None:
+            incar.pop(k, None)
+        else:
+            incar[k] = v
+
+
+def _remove_unused_incar_params(incar, skip: Sequence[str] = ()) -> None:
+    """
+    Remove INCAR parameters that are not actively used by VASP.
+
+    Args:
+        incar (dict): An incar.
+        skip (list of str): Keys to skip.
+    """
+    # Turn off IBRION/ISIF/POTIM if NSW = 0
+    opt_flags = ["EDIFFG", "IBRION", "ISIF", "POTIM"]
+    if incar.get("NSW", 0) == 0:
+        for opt_flag in opt_flags:
+            if opt_flag not in skip:
+                incar.pop(opt_flag, None)
+
+    # Remove MAGMOMs if they aren't used
+    if incar.get("ISPIN", 1) == 1 and "MAGMOM" not in skip:
+        incar.pop("MAGMOM", None)
+
+    # Turn off +U flags if +U is not even used
+    ldau_flags = ["LDAUU", "LDAUJ", "LDAUL", "LDAUTYPE"]
+    if incar.get("LDAU", False) is False:
+        for ldau_flag in ldau_flags:
+            if ldau_flag not in skip:
+                incar.pop(ldau_flag, None)
+
+
+def _get_nedos(vasprun: Vasprun | None, dedos: float) -> int:
+    """Automatic setting of nedos using the energy range and the energy step."""
+    if vasprun is None:
+        return 2000
+
+    emax = max(eigs.max() for eigs in vasprun.eigenvalues.values())
+    emin = min(eigs.min() for eigs in vasprun.eigenvalues.values())
+    return int((emax - emin) / dedos)
+
+
+def auto_kspacing(bandgap, bandgap_tol):
+    """Automatically set kspacing based on the bandgap"""
+    if bandgap is None or bandgap <= bandgap_tol:  # metallic
+        return 0.22
+
+    rmin = max(1.5, 25.22 - 2.87 * bandgap)  # Eq. 25
+    kspacing = 2 * np.pi * 1.0265 / (rmin - 1.0183)  # Eq. 29
+
+    # cap kspacing at a max of 0.44, per internal benchmarking
+    return min(kspacing, 0.44)
+
+
+def get_valid_magmom_struct(structure: Structure, inplace: bool = True, spin_mode: str = "auto") -> Structure:
+    """
+    Make sure that the structure has valid magmoms based on the kind of calculation.
+
+    Fill in missing Magmom values.
+
+    Args:
+        structure: The input structure
+        inplace: True: edit magmoms of the input structure; False: return new structure
+        spin_mode: "scalar"/"vector"/"none"/"auto" only first letter (s/v/n) is needed.
+            dictates how the spin configuration will be determined.
+
+            - auto: read the existing magmom values and decide
+            - scalar: use a single scalar value (for spin up/down)
+            - vector: use a vector value for spin-orbit systems
+            - none: Remove all the magmom information
+
+    Returns:
+        New structure if inplace is False
+    """
+    default_values = {"s": 1.0, "v": [1.0, 1.0, 1.0], "n": None}
+    if spin_mode[0].lower() == "a":
+        mode = "n"
+        for site in structure:
+            if "magmom" not in site.properties or site.properties["magmom"] is None:
+                pass
+            elif isinstance(site.properties["magmom"], (float, int)):
+                if mode == "v":
+                    raise TypeError("Magmom type conflict")
+                mode = "s"
+                if isinstance(site.properties["magmom"], int):
+                    site.properties["magmom"] = float(site.properties["magmom"])
+            elif len(site.properties["magmom"]) == 3:
+                if mode == "s":
+                    raise TypeError("Magmom type conflict")
+                mode = "v"
+            else:
+                raise TypeError("Unrecognized Magmom Value")
+    else:
+        mode = spin_mode[0].lower()
+
+    ret_struct = structure if inplace else structure.copy()
+    for site in ret_struct:
+        if mode == "n":
+            if "magmom" in site.properties:
+                site.properties.pop("magmom")
+        elif "magmom" not in site.properties or site.properties["magmom"] is None:
+            site.properties["magmom"] = default_values[mode]
+
+    return ret_struct
+
+
+def get_structure_from_prev_run(vasprun, outcar=None) -> Structure:
+    """
+    Process structure from previous run.
+
+    Args:
+        vasprun (Vasprun): Vasprun that contains the final structure from previous run.
+        outcar (Outcar): Outcar that contains the magnetization info from previous run.
+
+    Returns:
+        Structure: The magmom-decorated structure that can be passed to get VASP input files, e.g.
+            get_kpoints().
+    """
+    structure = vasprun.final_structure
+
+    site_properties = {}
+    # magmom
+    if vasprun.is_spin:
+        if outcar and outcar.magnetization:
+            site_properties["magmom"] = [i["tot"] for i in outcar.magnetization]
+        else:
+            site_properties["magmom"] = vasprun.parameters["MAGMOM"]
+    # LDAU
+    if vasprun.parameters.get("LDAU", False):
+        for key in ("LDAUU", "LDAUJ", "LDAUL"):
+            vals = vasprun.incar[key]
+            m = {}
+            l_val = []
+            s = 0
+            for site in structure:
+                if site.specie.symbol not in m:
+                    m[site.specie.symbol] = vals[s]
+                    s += 1
+                l_val.append(m[site.specie.symbol])
+            if len(l_val) == len(structure):
+                site_properties.update({key.lower(): l_val})
+            else:
+                raise ValueError(f"length of list {l_val} not the same as structure")
+
+    return structure.copy(site_properties=site_properties)
+
+
+def standardize_structure(structure, sym_prec=0.1, international_monoclinic=True):
+    """Get the symmetrically standardized structure.
+
+    Args:
+        structure (Structure): The structure.
+        sym_prec (float): Tolerance for symmetry finding for standardization.
+        international_monoclinic (bool): Whether to use international
+            convention (vs Curtarolo) for monoclinic. Defaults True.
+
+    Returns:
+        The symmetrized structure.
+    """
+    sym_finder = SpacegroupAnalyzer(structure, symprec=sym_prec)
+    new_structure = sym_finder.get_primitive_standard_structure(international_monoclinic=international_monoclinic)
+
+    # the primitive structure finding has had several bugs in the past
+    # defend through validation
+    vpa_old = structure.volume / len(structure)
+    vpa_new = new_structure.volume / len(new_structure)
+
+    if abs(vpa_old - vpa_new) / vpa_old > 0.02:
+        raise ValueError(f"Standardizing cell failed! VPA old: {vpa_old}, VPA new: {vpa_new}")
+
+    matcher = StructureMatcher()
+    if not matcher.fit(structure, new_structure):
+        raise ValueError("Standardizing cell failed! Old structure doesn't match new.")
+
+    return new_structure
+
+
+def get_vasprun_outcar(path: str | Path, parse_dos: bool = True, parse_eigen: bool = True) -> tuple[Vasprun, Outcar]:
+    """Get a Vasprun and Outcar from a directory.
+
+    Args:
+        path: Path to get the vasprun.xml and OUTCAR.
+        parse_dos: Whether to parse dos. Defaults to True.
+        parse_eigen: Whether to parse eigenvalue. Defaults to True.
+
+    Returns:
+        Vasprun and Outcar files.
+    """
+    path = Path(path)
+    vruns = list(glob(str(path / "vasprun.xml*")))
+    outcars = list(glob(str(path / "OUTCAR*")))
+
+    if len(vruns) == 0 or len(outcars) == 0:
+        raise ValueError(f"Unable to get vasprun.xml/OUTCAR from prev calculation in {path}")
+    vs_file_fullpath = str(path / "vasprun.xml")
+    outcar_fullpath = str(path / "OUTCAR.gz")
+    vs_file = vs_file_fullpath if vs_file_fullpath in vruns else max(vruns)
+    outcar_file = outcar_fullpath if outcar_fullpath in outcars else max(outcars)
+    return (
+        Vasprun(vs_file, parse_dos=parse_dos, parse_eigen=parse_eigen),
+        Outcar(outcar_file),
+    )