• Home
  • Features
  • Pricing
  • Docs
  • Announcements
  • Sign In

materialsproject / pymatgen / 4075885785

pending completion
4075885785

push

github

Shyue Ping Ong
Merge branch 'master' of github.com:materialsproject/pymatgen

96 of 96 new or added lines in 27 files covered. (100.0%)

81013 of 102710 relevant lines covered (78.88%)

0.79 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

26.64
/pymatgen/io/vasp/sets.py
1
# Copyright (c) Pymatgen Development Team.
2
# Distributed under the terms of the MIT License.
3

4
"""
1✔
5
This module defines the VaspInputSet abstract base class and a concrete
6
implementation for the parameters developed and tested by the core team
7
of pymatgen, including the Materials Virtual Lab, Materials Project and the MIT
8
high throughput project. The basic concept behind an input set is to specify
9
a scheme to generate a consistent set of VASP inputs from a structure
10
without further user intervention. This ensures comparability across
11
runs.
12

13
Read the following carefully before implementing new input sets:
14

15
1. 99% of what needs to be done can be done by specifying user_incar_settings
16
   to override some of the defaults of various input sets. Unless there is an
17
   extremely good reason to add a new set, DO NOT add one. E.g., if you want
18
   to turn the hubbard U off, just set "LDAU": False as a user_incar_setting.
19
2. All derivative input sets should inherit from one of the usual MPRelaxSet or
20
   MITRelaxSet, and proper superclass delegation should be used where possible.
21
   In particular, you are not supposed to implement your own as_dict or
22
   from_dict for derivative sets unless you know what you are doing.
23
   Improper overriding the as_dict and from_dict protocols is the major
24
   cause of implementation headaches. If you need an example, look at how the
25
   MPStaticSet or MPNonSCFSets are constructed.
26

27
The above are recommendations. The following are UNBREAKABLE rules:
28

29
1. All input sets must take in a structure or list of structures as the first
30
   argument.
31
2. user_incar_settings, user_kpoints_settings and user_<whatever>_settings are
32
   ABSOLUTE. Any new sets you implement must obey this. If a user wants to
33
   override your settings, you assume he knows what he is doing. Do not
34
   magically override user supplied settings. You can issue a warning if you
35
   think the user is wrong.
36
3. All input sets must save all supplied args and kwargs as instance variables.
37
   E.g., self.my_arg = my_arg and self.kwargs = kwargs in the __init__. This
38
   ensures the as_dict and from_dict work correctly.
39
"""
40

41
from __future__ import annotations
1✔
42

43
import abc
1✔
44
import glob
1✔
45
import itertools
1✔
46
import os
1✔
47
import re
1✔
48
import shutil
1✔
49
import warnings
1✔
50
from copy import deepcopy
1✔
51
from itertools import chain
1✔
52
from pathlib import Path
1✔
53
from zipfile import ZipFile
1✔
54

55
import numpy as np
1✔
56
from monty.io import zopen
1✔
57
from monty.json import MSONable
1✔
58
from monty.serialization import loadfn
1✔
59

60
from pymatgen.analysis.structure_matcher import StructureMatcher
1✔
61
from pymatgen.core.periodic_table import Element, Species
1✔
62
from pymatgen.core.sites import PeriodicSite
1✔
63
from pymatgen.core.structure import Structure
1✔
64
from pymatgen.io.vasp.inputs import Incar, Kpoints, Poscar, Potcar, VaspInput
1✔
65
from pymatgen.io.vasp.outputs import Outcar, Vasprun
1✔
66
from pymatgen.symmetry.analyzer import SpacegroupAnalyzer
1✔
67
from pymatgen.symmetry.bandstructure import HighSymmKpath
1✔
68

69
MODULE_DIR = Path(__file__).resolve().parent
1✔
70

71

72
class VaspInputSet(MSONable, metaclass=abc.ABCMeta):
1✔
73
    """
74
    Base class representing a set of VASP input parameters with a structure
75
    supplied as init parameters. Typically, you should not inherit from this
76
    class. Start from DictSet or MPRelaxSet or MITRelaxSet.
77
    """
78

79
    @property
1✔
80
    @abc.abstractmethod
1✔
81
    def incar(self):
1✔
82
        """Incar object"""
83

84
    @property
1✔
85
    @abc.abstractmethod
1✔
86
    def kpoints(self):
1✔
87
        """Kpoints object"""
88

89
    @property
1✔
90
    @abc.abstractmethod
1✔
91
    def poscar(self):
1✔
92
        """Poscar object"""
93

94
    @property
1✔
95
    def potcar_symbols(self):
1✔
96
        """
97
        List of POTCAR symbols.
98
        """
99
        # pylint: disable=E1101
100
        elements = self.poscar.site_symbols
1✔
101
        potcar_symbols = []
1✔
102
        settings = self._config_dict["POTCAR"]
1✔
103

104
        if isinstance(settings[elements[-1]], dict):
1✔
105
            for el in elements:
×
106
                potcar_symbols.append(settings[el]["symbol"] if el in settings else el)
×
107
        else:
108
            for el in elements:
1✔
109
                potcar_symbols.append(settings.get(el, el))
1✔
110

111
        return potcar_symbols
1✔
112

113
    @property
1✔
114
    def potcar(self):
1✔
115
        """
116
        Potcar object.
117
        """
118
        # pylint: disable=E1101
119
        potcar = Potcar(self.potcar_symbols, functional=self.potcar_functional)
1✔
120

121
        # warn if the selected POTCARs do not correspond to the chosen
122
        # potcar_functional
123
        for psingle in potcar:
1✔
124
            if self.potcar_functional not in psingle.identify_potcar()[0]:
1✔
125
                warnings.warn(
1✔
126
                    f"POTCAR data with symbol {psingle.symbol} is not known by pymatgen to\
127
                    correspond with the selected potcar_functional {self.potcar_functional}. This POTCAR\
128
                    is known to correspond with functionals {psingle.identify_potcar(mode='data')[0]}. "
129
                    f"Please verify that you are using the right POTCARs!",
130
                    BadInputSetWarning,
131
                )
132

133
        return potcar
1✔
134

135
    def get_vasp_input(self) -> VaspInput:
1✔
136
        """
137
        Returns:
138
            VaspInput
139
        """
140
        return VaspInput(
1✔
141
            incar=self.incar,
142
            kpoints=self.kpoints,
143
            poscar=self.poscar,
144
            potcar=self.potcar,
145
        )
146

147
    def write_input(
1✔
148
        self,
149
        output_dir,
150
        make_dir_if_not_present=True,
151
        include_cif=False,
152
        potcar_spec=False,
153
        zip_output=False,
154
    ):
155
        """
156
        Writes a set of VASP input to a directory.
157

158
        Args:
159
            output_dir (str): Directory to output the VASP input files
160
            make_dir_if_not_present (bool): Set to True if you want the
161
                directory (and the whole path) to be created if it is not
162
                present.
163
            include_cif (bool): Whether to write a CIF file in the output
164
                directory for easier opening by VESTA.
165
            potcar_spec (bool): Instead of writing the POTCAR, write a "POTCAR.spec".
166
                This is intended to help sharing an input set with people who might
167
                not have a license to specific Potcar files. Given a "POTCAR.spec",
168
                the specific POTCAR file can be re-generated using pymatgen with the
169
                "generate_potcar" function in the pymatgen CLI.
170
            zip_output (bool): If True, output will be zipped into a file with the
171
                same name as the InputSet (e.g., MPStaticSet.zip)
172
        """
173
        if potcar_spec:
×
174
            if make_dir_if_not_present and not os.path.exists(output_dir):
×
175
                os.makedirs(output_dir)
×
176

177
            with zopen(os.path.join(output_dir, "POTCAR.spec"), "wt") as f:
×
178
                f.write("\n".join(self.potcar_symbols))
×
179

180
            for k, v in {
×
181
                "INCAR": self.incar,
182
                "POSCAR": self.poscar,
183
                "KPOINTS": self.kpoints,
184
            }.items():
185
                if v is not None:
×
186
                    with zopen(os.path.join(output_dir, k), "wt") as f:
×
187
                        f.write(str(v))
×
188
        else:
189
            vinput = self.get_vasp_input()
×
190
            vinput.write_input(output_dir, make_dir_if_not_present=make_dir_if_not_present)
×
191

192
        cif_name = ""
×
193
        if include_cif:
×
194
            s = vinput["POSCAR"].structure
×
195
            cif_name = f"{output_dir}/{s.formula.replace(' ', '')}.cif"
×
196
            s.to(filename=cif_name)
×
197

198
        if zip_output:
×
199
            filename = type(self).__name__ + ".zip"
×
200
            with ZipFile(os.path.join(output_dir, filename), "w") as zip:
×
201
                for file in [
×
202
                    "INCAR",
203
                    "POSCAR",
204
                    "KPOINTS",
205
                    "POTCAR",
206
                    "POTCAR.spec",
207
                    cif_name,
208
                ]:
209
                    try:
×
210
                        zip.write(os.path.join(output_dir, file), arcname=file)
×
211
                    except FileNotFoundError:
×
212
                        pass
×
213
                    try:
×
214
                        os.remove(os.path.join(output_dir, file))
×
215
                    except (FileNotFoundError, PermissionError, IsADirectoryError):
×
216
                        pass
×
217

218
    def as_dict(self, verbosity=2):
1✔
219
        """
220
        Args:
221
            verbosity: Verbosity for generated dict. If 1, structure is
222
            excluded.
223

224
        Returns:
225
            MSONable dict
226
        """
227
        d = MSONable.as_dict(self)
×
228
        if verbosity == 1:
×
229
            d.pop("structure", None)
×
230
        return d
×
231

232

233
def _load_yaml_config(fname):
1✔
234
    config = loadfn(str(MODULE_DIR / (f"{fname}.yaml")))
1✔
235
    if "PARENT" in config:
1✔
236
        parent_config = _load_yaml_config(config["PARENT"])
1✔
237
        for k, v in parent_config.items():
1✔
238
            if k not in config:
1✔
239
                config[k] = v
×
240
            elif isinstance(v, dict):
1✔
241
                v_new = config.get(k, {})
1✔
242
                v_new.update(v)
1✔
243
                config[k] = v_new
1✔
244
    return config
1✔
245

246

247
class DictSet(VaspInputSet):
1✔
248
    """
249
    Concrete implementation of VaspInputSet that is initialized from a dict
250
    settings. This allows arbitrary settings to be input. In general,
251
    this is rarely used directly unless there is a source of settings in yaml
252
    format (e.g., from a REST interface). It is typically used by other
253
    VaspInputSets for initialization.
254

255
    Special consideration should be paid to the way the MAGMOM initialization
256
    for the INCAR is done. The initialization differs depending on the type of
257
    structure and the configuration settings. The order in which the magmom is
258
    determined is as follows:
259

260
    1. If the site itself has a magmom setting (i.e. site.properties["magmom"] = float),
261
        that is used. This can be set with structure.add_site_property().
262
    2. If the species of the site has a spin setting, that is used. This can be set
263
        with structure.add_spin_by_element().
264
    3. If the species itself has a particular setting in the config file, that
265
       is used, e.g., Mn3+ may have a different magmom than Mn4+.
266
    4. Lastly, the element symbol itself is checked in the config file. If
267
       there are no settings, a default value of 0.6 is used.
268
    """
269

270
    def __init__(
1✔
271
        self,
272
        structure,
273
        config_dict,
274
        files_to_transfer=None,
275
        user_incar_settings=None,
276
        user_kpoints_settings=None,
277
        user_potcar_settings=None,
278
        constrain_total_magmom=False,
279
        sort_structure=True,
280
        potcar_functional=None,
281
        user_potcar_functional=None,
282
        force_gamma=False,
283
        reduce_structure=None,
284
        vdw=None,
285
        use_structure_charge=False,
286
        standardize=False,
287
        sym_prec=0.1,
288
        international_monoclinic=True,
289
        validate_magmom=True,
290
    ):
291
        """
292
        Args:
293
            structure (Structure): The Structure to create inputs for.
294
            config_dict (dict): The config dictionary to use.
295
            files_to_transfer (dict): A dictionary of {filename: filepath}. This
296
                allows the transfer of files from a previous calculation.
297
            user_incar_settings (dict): User INCAR settings. This allows a user
298
                to override INCAR settings, e.g., setting a different MAGMOM for
299
                various elements or species. Note that in the new scheme,
300
                ediff_per_atom and hubbard_u are no longer args. Instead, the
301
                config_dict supports EDIFF_PER_ATOM and EDIFF keys. The former
302
                scales with # of atoms, the latter does not. If both are
303
                present, EDIFF is preferred. To force such settings, just supply
304
                user_incar_settings={"EDIFF": 1e-5, "LDAU": False} for example.
305
                The keys 'LDAUU', 'LDAUJ', 'LDAUL' are special cases since
306
                pymatgen defines different values depending on what anions are
307
                present in the structure, so these keys can be defined in one
308
                of two ways, e.g. either {"LDAUU":{"O":{"Fe":5}}} to set LDAUU
309
                for Fe to 5 in an oxide, or {"LDAUU":{"Fe":5}} to set LDAUU to
310
                5 regardless of the input structure.
311

312
                If a None value is given, that key is unset. For example,
313
                {"ENCUT": None} will remove ENCUT from the incar settings.
314
            user_kpoints_settings (dict or Kpoints): Allow user to override kpoints
315
                setting by supplying a dict E.g., {"reciprocal_density": 1000}.
316
                User can also supply Kpoints object. Default is None.
317
            user_potcar_settings (dict: Allow user to override POTCARs. E.g.,
318
                {"Gd": "Gd_3"}. This is generally not recommended. Default is None.
319
            constrain_total_magmom (bool): Whether to constrain the total magmom
320
                (NUPDOWN in INCAR) to be the sum of the expected MAGMOM for all
321
                species. Defaults to False.
322
            sort_structure (bool): Whether to sort the structure (using the
323
                default sort order of electronegativity) before generating input
324
                files. Defaults to True, the behavior you would want most of the
325
                time. This ensures that similar atomic species are grouped
326
                together.
327
            user_potcar_functional (str): Functional to use. Default (None) is to use
328
                the functional in the config dictionary. Valid values:
329
                "PBE", "PBE_52", "PBE_54", "LDA", "LDA_52", "LDA_54", "PW91",
330
                "LDA_US", "PW91_US".
331
            force_gamma (bool): Force gamma centered kpoint generation. Default
332
                (False) is to use the Automatic Density kpoint scheme, which
333
                will use the Gamma centered generation scheme for hexagonal
334
                cells, and Monkhorst-Pack otherwise.
335
            reduce_structure (None/str): Before generating the input files,
336
                generate the reduced structure. Default (None), does not
337
                alter the structure. Valid values: None, "niggli", "LLL".
338
            vdw: Adds default parameters for van-der-Waals functionals supported
339
                by VASP to INCAR. Supported functionals are: DFT-D2, undamped
340
                DFT-D3, DFT-D3 with Becke-Jonson damping, Tkatchenko-Scheffler,
341
                Tkatchenko-Scheffler with iterative Hirshfeld partitioning,
342
                MBD@rSC, dDsC, Dion's vdW-DF, DF2, optPBE, optB88, optB86b and
343
                rVV10.
344
            use_structure_charge (bool): If set to True, then the public
345
                variable used for setting the overall charge of the
346
                structure (structure.charge) is used to set the NELECT
347
                variable in the INCAR
348
                Default is False (structure's overall charge is not used)
349
            standardize (float): Whether to standardize to a primitive standard
350
                cell. Defaults to False.
351
            sym_prec (float): Tolerance for symmetry finding.
352
            international_monoclinic (bool): Whether to use international convention
353
                (vs Curtarolo) for monoclinic. Defaults True.
354
            validate_magmom (bool): Ensure that the missing magmom values are filled
355
                in with the VASP default value of 1.0
356
        """
357
        if reduce_structure:
1✔
358
            structure = structure.get_reduced_structure(reduce_structure)
×
359
        if sort_structure:
1✔
360
            structure = structure.get_sorted_structure()
1✔
361
        if validate_magmom:
1✔
362
            get_valid_magmom_struct(structure, spin_mode="auto", inplace=True)
1✔
363

364
        self._structure = structure
1✔
365
        self._config_dict = deepcopy(config_dict)
1✔
366
        self.files_to_transfer = files_to_transfer or {}
1✔
367
        self.constrain_total_magmom = constrain_total_magmom
1✔
368
        self.sort_structure = sort_structure
1✔
369
        self.force_gamma = force_gamma
1✔
370
        self.reduce_structure = reduce_structure
1✔
371
        self.user_incar_settings = user_incar_settings or {}
1✔
372
        self.user_kpoints_settings = user_kpoints_settings or {}
1✔
373
        self.user_potcar_settings = user_potcar_settings
1✔
374
        self.vdw = vdw.lower() if vdw is not None else None
1✔
375
        self.use_structure_charge = use_structure_charge
1✔
376
        self.standardize = standardize
1✔
377
        self.sym_prec = sym_prec
1✔
378
        self.international_monoclinic = international_monoclinic
1✔
379

380
        if self.user_incar_settings.get("KSPACING") and user_kpoints_settings is not None:
1✔
381
            warnings.warn(
×
382
                "You have specified KSPACING and also supplied kpoints "
383
                "settings. KSPACING only has effect when there is no "
384
                "KPOINTS file. Since both settings were given, pymatgen"
385
                "will generate a KPOINTS file and ignore KSPACING."
386
                "Remove the `user_kpoints_settings` argument to enable KSPACING.",
387
                BadInputSetWarning,
388
            )
389

390
        if self.vdw:
1✔
391
            vdw_par = loadfn(str(MODULE_DIR / "vdW_parameters.yaml"))
×
392
            try:
×
393
                self._config_dict["INCAR"].update(vdw_par[self.vdw])
×
394
            except KeyError:
×
395
                raise KeyError(
×
396
                    f"Invalid or unsupported van-der-Waals functional. Supported functionals are {', '.join(vdw_par)}."
397
                )
398
        # read the POTCAR_FUNCTIONAL from the .yaml
399
        self.potcar_functional = self._config_dict.get("POTCAR_FUNCTIONAL", "PBE")
1✔
400

401
        if potcar_functional is not None and user_potcar_functional is not None:
1✔
402
            raise ValueError(
×
403
                "Received both 'potcar_functional' and "
404
                "'user_potcar_functional arguments. 'potcar_functional "
405
                "is deprecated."
406
            )
407
        if potcar_functional:
1✔
408
            warnings.warn(
×
409
                "'potcar_functional' argument is deprecated. Use 'user_potcar_functional' instead.",
410
                FutureWarning,
411
            )
412
            self.potcar_functional = potcar_functional
×
413
        elif user_potcar_functional:
1✔
414
            self.potcar_functional = user_potcar_functional
×
415

416
        # warn if a user is overriding POTCAR_FUNCTIONAL
417
        if self.potcar_functional != self._config_dict.get("POTCAR_FUNCTIONAL"):
1✔
418
            warnings.warn(
×
419
                "Overriding the POTCAR functional is generally not recommended "
420
                " as it significantly affect the results of calculations and "
421
                "compatibility with other calculations done with the same "
422
                "input set. Note that some POTCAR symbols specified in "
423
                "the configuration file may not be available in the selected "
424
                "functional.",
425
                BadInputSetWarning,
426
            )
427

428
        if self.user_potcar_settings:
1✔
429
            warnings.warn(
×
430
                "Overriding POTCARs is generally not recommended as it "
431
                "significantly affect the results of calculations and "
432
                "compatibility with other calculations done with the same "
433
                "input set. In many instances, it is better to write a "
434
                "subclass of a desired input set and override the POTCAR in "
435
                "the subclass to be explicit on the differences.",
436
                BadInputSetWarning,
437
            )
438
            for k, v in self.user_potcar_settings.items():
×
439
                self._config_dict["POTCAR"][k] = v
×
440

441
    @property
1✔
442
    def structure(self) -> Structure:
1✔
443
        """
444
        :return: Structure
445
        """
446
        if self.standardize and self.sym_prec:
1✔
447
            return standardize_structure(
×
448
                self._structure,
449
                sym_prec=self.sym_prec,
450
                international_monoclinic=self.international_monoclinic,
451
            )
452
        return self._structure
1✔
453

454
    @property
1✔
455
    def incar(self) -> Incar:
1✔
456
        """
457
        :return: Incar
458
        """
459
        settings = dict(self._config_dict["INCAR"])
1✔
460
        for k, v in self.user_incar_settings.items():
1✔
461
            if v is None:
×
462
                try:
×
463
                    del settings[k]
×
464
                except KeyError:
×
465
                    settings[k] = v
×
466
            elif k == "KSPACING" and self.user_kpoints_settings != {}:
×
467
                pass  # Ignore KSPACING if user_kpoints_settings are given
×
468
            else:
469
                settings[k] = v
×
470
        structure = self.structure
1✔
471
        incar = Incar()
1✔
472
        comp = structure.composition
1✔
473
        elements = sorted((el for el in comp.elements if comp[el] > 0), key=lambda e: e.X)
1✔
474
        most_electroneg = elements[-1].symbol
1✔
475
        poscar = Poscar(structure)
1✔
476
        hubbard_u = settings.get("LDAU", False)
1✔
477

478
        for k, v in settings.items():
1✔
479
            if k == "MAGMOM":
1✔
480
                mag = []
1✔
481
                for site in structure:
1✔
482
                    if hasattr(site, "magmom"):
1✔
483
                        mag.append(site.magmom)
×
484
                    elif hasattr(site.specie, "spin"):
1✔
485
                        mag.append(site.specie.spin)
×
486
                    elif str(site.specie) in v:
1✔
487
                        if site.specie.symbol == "Co" and v[str(site.specie)] <= 1.0:
1✔
488
                            warnings.warn(
×
489
                                "Co without an oxidation state is initialized as low spin by default in Pymatgen. "
490
                                "If this default behavior is not desired, please set the spin on the magmom on the "
491
                                "site directly to ensure correct initialization."
492
                            )
493
                        mag.append(v.get(str(site.specie)))
1✔
494
                    else:
495
                        if site.specie.symbol == "Co":
1✔
496
                            warnings.warn(
×
497
                                "Co without an oxidation state is initialized as low spin by default in Pymatgen. "
498
                                "If this default behavior is not desired, please set the spin on the magmom on the "
499
                                "site directly to ensure correct initialization."
500
                            )
501
                        mag.append(v.get(site.specie.symbol, 0.6))
1✔
502
                incar[k] = mag
1✔
503
            elif k in ("LDAUU", "LDAUJ", "LDAUL"):
1✔
504
                if hubbard_u:
1✔
505
                    if hasattr(structure[0], k.lower()):
1✔
506
                        m = {site.specie.symbol: getattr(site, k.lower()) for site in structure}
×
507
                        incar[k] = [m[sym] for sym in poscar.site_symbols]
×
508
                        # lookup specific LDAU if specified for most_electroneg atom
509
                    elif most_electroneg in v and isinstance(v[most_electroneg], dict):
1✔
510
                        incar[k] = [v[most_electroneg].get(sym, 0) for sym in poscar.site_symbols]
1✔
511
                        # else, use fallback LDAU value if it exists
512
                    else:
513
                        incar[k] = [
×
514
                            v.get(sym, 0) if isinstance(v.get(sym, 0), (float, int)) else 0
515
                            for sym in poscar.site_symbols
516
                        ]
517
            elif k.startswith("EDIFF") and k != "EDIFFG":
1✔
518
                if "EDIFF" not in settings and k == "EDIFF_PER_ATOM":
1✔
519
                    incar["EDIFF"] = float(v) * structure.num_sites
1✔
520
                else:
521
                    incar["EDIFF"] = float(settings["EDIFF"])
×
522
            else:
523
                incar[k] = v
1✔
524
        has_u = hubbard_u and sum(incar["LDAUU"]) > 0
1✔
525
        if not has_u:
1✔
526
            for key in list(incar):
×
527
                if key.startswith("LDAU"):
×
528
                    del incar[key]
×
529

530
        # Modify LMAXMIX if you have d or f electrons present.
531
        # Note that if the user explicitly sets LMAXMIX in settings it will
532
        # override this logic.
533
        # Previously, this was only set if Hubbard U was enabled as per the
534
        # VASP manual but following an investigation it was determined that
535
        # this would lead to a significant difference between SCF -> NonSCF
536
        # even without Hubbard U enabled. Thanks to Andrew Rosen for
537
        # investigating and reporting.
538
        if "LMAXMIX" not in settings:
1✔
539
            # contains f-electrons
540
            if any(el.Z > 56 for el in structure.composition):
1✔
541
                incar["LMAXMIX"] = 6
×
542
            # contains d-electrons
543
            elif any(el.Z > 20 for el in structure.composition):
1✔
544
                incar["LMAXMIX"] = 4
1✔
545

546
        # Warn user about LASPH for +U, meta-GGAs, hybrids, and vdW-DF
547
        if not incar.get("LASPH", False) and (
1✔
548
            incar.get("METAGGA")
549
            or incar.get("LHFCALC", False)
550
            or incar.get("LDAU", False)
551
            or incar.get("LUSE_VDW", False)
552
        ):
553
            warnings.warn(
×
554
                "LASPH = True should be set for +U, meta-GGAs, hybrids, and vdW-DFT",
555
                BadInputSetWarning,
556
            )
557

558
        if self.constrain_total_magmom:
1✔
559
            nupdown = sum(mag if abs(mag) > 0.6 else 0 for mag in incar["MAGMOM"])
×
560
            if abs(nupdown - round(nupdown)) > 1e-5:
×
561
                warnings.warn(
×
562
                    "constrain_total_magmom was set to True, but the sum of MAGMOM "
563
                    "values is not an integer. NUPDOWN is meant to set the spin "
564
                    "multiplet and should typically be an integer. You are likely "
565
                    "better off changing the values of MAGMOM or simply setting "
566
                    "NUPDOWN directly in your INCAR settings.",
567
                    UserWarning,
568
                )
569
            incar["NUPDOWN"] = nupdown
×
570

571
        if self.use_structure_charge:
1✔
572
            incar["NELECT"] = self.nelect
×
573

574
        # Check that ALGO is appropriate
575
        if incar.get("LHFCALC", False) is True and incar.get("ALGO", "Normal") not in ["Normal", "All", "Damped"]:
1✔
576
            warnings.warn(
×
577
                "Hybrid functionals only support Algo = All, Damped, or Normal.",
578
                BadInputSetWarning,
579
            )
580

581
        # Ensure adequate number of KPOINTS are present for the tetrahedron
582
        # method (ISMEAR=-5). If KSPACING is in the INCAR file the number
583
        # of kpoints is not known before calling VASP, but a warning is raised
584
        # when the KSPACING value is > 0.5 (2 reciprocal Angstrom).
585
        # An error handler in Custodian is available to
586
        # correct overly large KSPACING values (small number of kpoints)
587
        # if necessary.
588
        # if "KSPACING" not in self.user_incar_settings:
589
        if self.kpoints is not None:
1✔
590
            if np.product(self.kpoints.kpts) < 4 and incar.get("ISMEAR", 0) == -5:
1✔
591
                incar["ISMEAR"] = 0
×
592

593
        if self.user_incar_settings.get("KSPACING", 0) > 0.5 and incar.get("ISMEAR", 0) == -5:
1✔
594
            warnings.warn(
×
595
                "Large KSPACING value detected with ISMEAR = -5. Ensure that VASP "
596
                "generates an adequate number of KPOINTS, lower KSPACING, or "
597
                "set ISMEAR = 0",
598
                BadInputSetWarning,
599
            )
600

601
        if all(k.is_metal for k in structure.composition):
1✔
602
            if incar.get("NSW", 0) > 0 and incar.get("ISMEAR", 1) < 1:
×
603
                warnings.warn(
×
604
                    "Relaxation of likely metal with ISMEAR < 1 "
605
                    "detected. Please see VASP recommendations on "
606
                    "ISMEAR for metals.",
607
                    BadInputSetWarning,
608
                )
609

610
        return incar
1✔
611

612
    @property
1✔
613
    def poscar(self) -> Poscar:
1✔
614
        """
615
        :return: Poscar
616
        """
617
        return Poscar(self.structure)
1✔
618

619
    @property
1✔
620
    def nelect(self) -> float:
1✔
621
        """
622
        Gets the default number of electrons for a given structure.
623
        """
624
        nelectrons_by_element = {p.element: p.nelectrons for p in self.potcar}
×
625
        nelect = sum(
×
626
            num_atoms * nelectrons_by_element[str(el)]
627
            for el, num_atoms in self.structure.composition.element_composition.items()
628
        )
629

630
        if self.use_structure_charge:
×
631
            return nelect - self.structure.charge
×
632
        return nelect
×
633

634
    @property
1✔
635
    def kpoints(self) -> Kpoints | None:
1✔
636
        """
637
        Returns a KPOINTS file using the fully automated grid method. Uses
638
        Gamma centered meshes for hexagonal cells and Monk grids otherwise.
639

640
        If KSPACING is set in user_incar_settings (or the INCAR file), no
641
        file is created because VASP will automatically generate the kpoints.
642

643
        Algorithm:
644
            Uses a simple approach scaling the number of divisions along each
645
            reciprocal lattice vector proportional to its length.
646
        """
647
        # Return None if KSPACING is present in the INCAR, because this will
648
        # cause VASP to generate the kpoints automatically
649
        if self.user_incar_settings.get("KSPACING") or self._config_dict["INCAR"].get("KSPACING"):
1✔
650
            if self.user_kpoints_settings == {}:
×
651
                return None
×
652

653
        settings = self.user_kpoints_settings or self._config_dict.get("KPOINTS")
1✔
654

655
        if isinstance(settings, Kpoints):
1✔
656
            return settings
×
657

658
        # Return None if KSPACING is present in the INCAR, because this will
659
        # cause VASP to generate the kpoints automatically
660
        if self.user_incar_settings.get("KSPACING") and self.user_kpoints_settings == {}:
1✔
661
            return None
×
662

663
        # If grid_density is in the kpoints_settings use
664
        # Kpoints.automatic_density
665
        if settings.get("grid_density"):
1✔
666
            return Kpoints.automatic_density(self.structure, int(settings["grid_density"]), self.force_gamma)
×
667

668
        # If reciprocal_density is in the kpoints_settings use
669
        # Kpoints.automatic_density_by_vol
670
        if settings.get("reciprocal_density"):
1✔
671
            return Kpoints.automatic_density_by_vol(
1✔
672
                self.structure, int(settings["reciprocal_density"]), self.force_gamma
673
            )
674

675
        # If length is in the kpoints_settings use Kpoints.automatic
676
        if settings.get("length"):
×
677
            return Kpoints.automatic(settings["length"])
×
678

679
        # Raise error. Unsure of which kpoint generation to use
680
        raise ValueError(
×
681
            "Invalid KPoint Generation algo : Supported Keys are "
682
            "grid_density: for Kpoints.automatic_density generation, "
683
            "reciprocal_density: for KPoints.automatic_density_by_vol "
684
            "generation, and length  : for Kpoints.automatic generation"
685
        )
686

687
    def estimate_nbands(self) -> int:
1✔
688
        """
689
        Estimate the number of bands that VASP will initialize a
690
        calculation with by default. Note that in practice this
691
        can depend on # of cores (if not set explicitly)
692
        """
693
        nions = len(self.structure)
×
694

695
        # from VASP's point of view, the number of magnetic atoms are
696
        # the number of atoms with non-zero magmoms, so use Incar as
697
        # source of truth
698
        nmag = len([m for m in self.incar["MAGMOM"] if not np.allclose(m, 0)])
×
699

700
        # by definition, if non-spin polarized ignore nmag
701
        if (not nmag) or (self.incar["ISPIN"] == 1):
×
702
            nbands = np.ceil(self.nelect / 2 + nions / 2)
×
703
        else:
704
            nbands = np.ceil(0.6 * self.nelect + nmag)
×
705

706
        return int(nbands)
×
707

708
    def __str__(self):
1✔
709
        return type(self).__name__
×
710

711
    def __repr__(self):
1✔
712
        return type(self).__name__
×
713

714
    def write_input(
1✔
715
        self,
716
        output_dir: str,
717
        make_dir_if_not_present: bool = True,
718
        include_cif: bool = False,
719
        potcar_spec: bool = False,
720
        zip_output: bool = False,
721
    ):
722
        """
723
        Writes out all input to a directory.
724

725
        Args:
726
            output_dir (str): Directory to output the VASP input files
727
            make_dir_if_not_present (bool): Set to True if you want the
728
                directory (and the whole path) to be created if it is not
729
                present.
730
            include_cif (bool): Whether to write a CIF file in the output
731
                directory for easier opening by VESTA.
732
            potcar_spec (bool): Instead of writing the POTCAR, write a "POTCAR.spec".
733
                This is intended to help sharing an input set with people who might
734
                not have a license to specific Potcar files. Given a "POTCAR.spec",
735
                the specific POTCAR file can be re-generated using pymatgen with the
736
                "generate_potcar" function in the pymatgen CLI.
737
            zip_output (bool): Whether to zip each VASP input file written to the output directory.
738
        """
739
        super().write_input(
×
740
            output_dir=output_dir,
741
            make_dir_if_not_present=make_dir_if_not_present,
742
            include_cif=include_cif,
743
            potcar_spec=potcar_spec,
744
            zip_output=zip_output,
745
        )
746
        for k, v in self.files_to_transfer.items():
×
747
            with zopen(v, "rb") as fin, zopen(str(Path(output_dir) / k), "wb") as fout:
×
748
                shutil.copyfileobj(fin, fout)
×
749

750
    def calculate_ng(self, max_prime_factor: int = 7, must_inc_2: bool = True) -> tuple:
1✔
751
        """
752
        Calculates the NGX, NGY, and NGZ values using the information available in the INCAR and POTCAR
753
        This is meant to help with making initial guess for the FFT grid so we can interact with the Charge density API
754

755
        Args:
756
            max_prime_factor (int): the valid prime factors of the grid size in each direction
757
                VASP has many different setting for this to handle many compiling options.
758
                For typical MPI options all prime factors up to 7 are allowed
759
            must_inc_2 (bool): Whether 2 must be a prime factor of the result. Defaults to True.
760
        """
761
        # TODO throw error for Ultrasoft potentials
762

763
        _RYTOEV = 13.605826
×
764
        _AUTOA = 0.529177249
×
765
        _PI = 3.141592653589793238
×
766

767
        # TODO Only do this for VASP 6 for now. Older version require more advanced logic
768

769
        # get the ENCUT val
770
        if "ENCUT" in self.incar and self.incar["ENCUT"] > 0:
×
771
            encut = self.incar["ENCUT"]
×
772
        else:
773
            encut = max(i_species.enmax for i_species in self.get_vasp_input()["POTCAR"])
×
774

775
        _CUTOF = [
×
776
            np.sqrt(encut / _RYTOEV) / (2 * _PI / (anorm / _AUTOA)) for anorm in self.poscar.structure.lattice.abc
777
        ]
778

779
        _PREC = "Normal"  # VASP default
×
780
        if "PREC" in self.incar:
×
781
            _PREC = self.incar["PREC"]
×
782

783
        if _PREC[0].lower() in {"l", "m", "h"}:
×
784
            raise NotImplementedError(
×
785
                "PREC = LOW/MEDIUM/HIGH from VASP 4.x and not supported, Please use NORMA/SINGLE/ACCURATE"
786
            )
787

788
        if _PREC[0].lower() in {"a", "s"}:  # TODO This only works in VASP 6.x
×
789
            _WFACT = 4
×
790
        else:
791
            _WFACT = 3
×
792

793
        def next_g_size(cur_g_size):
×
794
            g_size = int(_WFACT * cur_g_size + 0.5)
×
795
            return next_num_with_prime_factors(g_size, max_prime_factor, must_inc_2)
×
796

797
        ng_vec = [*map(next_g_size, _CUTOF)]
×
798

799
        if _PREC[0].lower() in {"a", "n"}:  # TODO This works for VASP 5.x and 6.x
×
800
            finer_g_scale = 2
×
801
        else:
802
            finer_g_scale = 1
×
803

804
        return ng_vec, [ng_ * finer_g_scale for ng_ in ng_vec]
×
805

806

807
# Helper functions to determine valid FFT grids for VASP
808
def next_num_with_prime_factors(n: int, max_prime_factor: int, must_inc_2: bool = True) -> int:
1✔
809
    """
810
    Return the next number greater than or equal to n that only has the desired prime factors
811

812
    Args:
813
        n (int): Initial guess at the grid density
814
        max_prime_factor (int): the maximum prime factor
815
        must_inc_2 (bool): 2 must be a prime factor of the result
816

817
    Returns:
818
        int: first product of the prime_factors that is >= n
819
    """
820
    if max_prime_factor < 2:
×
821
        raise ValueError("Must choose a maximum prime factor greater than 2")
×
822
    prime_factors = primes_less_than(max_prime_factor)
×
823
    for new_val in itertools.count(start=n):
×
824
        if must_inc_2 and new_val % 2 != 0:
×
825
            continue
×
826
        cur_val_ = new_val
×
827
        for j in prime_factors:
×
828
            while cur_val_ % j == 0:
×
829
                cur_val_ //= j
×
830
        if cur_val_ == 1:
×
831
            return new_val
×
832
    raise ValueError("No factorable number found, not possible.")
×
833

834

835
def primes_less_than(max_val: int) -> list[int]:
1✔
836
    """
837
    Get the primes less than or equal to the max value
838
    """
839
    res = []
×
840
    for i in range(2, max_val + 1):
×
841
        for j in range(2, i):
×
842
            if i % j == 0:
×
843
                break
×
844
        else:
845
            res.append(i)
×
846
    return res
×
847

848

849
class MITRelaxSet(DictSet):
1✔
850
    """
851
    Standard implementation of VaspInputSet utilizing parameters in the MIT
852
    High-throughput project.
853
    The parameters are chosen specifically for a high-throughput project,
854
    which means in general pseudopotentials with fewer electrons were chosen.
855

856
    Please refer::
857

858
        A Jain, G. Hautier, C. Moore, S. P. Ong, C. Fischer, T. Mueller,
859
        K. A. Persson, G. Ceder. A high-throughput infrastructure for density
860
        functional theory calculations. Computational Materials Science,
861
        2011, 50(8), 2295-2310. doi:10.1016/j.commatsci.2011.02.023
862
    """
863

864
    CONFIG = _load_yaml_config("MITRelaxSet")
1✔
865

866
    def __init__(self, structure: Structure, **kwargs):
1✔
867
        """
868
        :param structure: Structure
869
        :param kwargs: Same as those supported by DictSet.
870
        """
871
        super().__init__(structure, MITRelaxSet.CONFIG, **kwargs)
×
872
        self.kwargs = kwargs
×
873

874

875
class MPRelaxSet(DictSet):
1✔
876
    """
877
    Implementation of VaspInputSet utilizing parameters in the public
878
    Materials Project. Typically, the pseudopotentials chosen contain more
879
    electrons than the MIT parameters, and the k-point grid is ~50% more dense.
880
    The LDAUU parameters are also different due to the different psps used,
881
    which result in different fitted values.
882
    """
883

884
    CONFIG = _load_yaml_config("MPRelaxSet")
1✔
885

886
    def __init__(self, structure: Structure, **kwargs):
1✔
887
        """
888
        :param structure: Structure
889
        :param kwargs: Same as those supported by DictSet.
890
        """
891
        super().__init__(structure, MPRelaxSet.CONFIG, **kwargs)
1✔
892
        self.kwargs = kwargs
1✔
893

894

895
class MPScanRelaxSet(DictSet):
1✔
896
    """
897
    Class for writing a relaxation input set using the accurate and numerically
898
    efficient r2SCAN variant of the Strongly Constrained and Appropriately Normed
899
    (SCAN) metaGGA density functional.
900

901
    Notes:
902
        1. This functional is officially supported in VASP 6.0.0 and above. On older version,
903
        source code may be obtained by contacting the authors of the referenced manuscript.
904
        The original SCAN functional, available from VASP 5.4.3 onwards, maybe used instead
905
        by passing `user_incar_settings={"METAGGA": "SCAN"}` when instantiating this InputSet.
906
        r2SCAN and SCAN are expected to yield very similar results.
907

908
        2. Meta-GGA calculations require POTCAR files that include
909
        information on the kinetic energy density of the core-electrons,
910
        i.e. "PBE_52" or "PBE_54". Make sure the POTCARs include the
911
        following lines (see VASP wiki for more details):
912

913
            $ grep kinetic POTCAR
914
            kinetic energy-density
915
            mkinetic energy-density pseudized
916
            kinetic energy density (partial)
917

918
    References:
919
         James W. Furness, Aaron D. Kaplan, Jinliang Ning, John P. Perdew, and Jianwei Sun.
920
         Accurate and Numerically Efficient r2SCAN Meta-Generalized Gradient Approximation.
921
         The Journal of Physical Chemistry Letters 0, 11 DOI: 10.1021/acs.jpclett.0c02405
922
    """
923

924
    CONFIG = _load_yaml_config("MPSCANRelaxSet")
1✔
925

926
    def __init__(self, structure: Structure, bandgap=0, **kwargs):
1✔
927
        """
928
        Args:
929
            structure (Structure): Input structure.
930
            bandgap (int): Bandgap of the structure in eV. The bandgap is used to
931
                    compute the appropriate k-point density and determine the
932
                    smearing settings.
933

934
                    Metallic systems (default, bandgap = 0) use a KSPACING value of 0.22
935
                    and Methfessel-Paxton order 2 smearing (ISMEAR=2, SIGMA=0.2).
936

937
                    Non-metallic systems (bandgap > 0) use the tetrahedron smearing
938
                    method (ISMEAR=-5, SIGMA=0.05). The KSPACING value is
939
                    calculated from the bandgap via Eqs. 25 and 29 of Wisesa, McGill,
940
                    and Mueller [1] (see References). Note that if 'user_incar_settings'
941
                    or 'user_kpoints_settings' override KSPACING, the calculation from
942
                    bandgap is not performed.
943

944
            vdw (str): set "rVV10" to enable SCAN+rVV10, which is a versatile
945
                    van der Waals density functional by combing the SCAN functional
946
                    with the rVV10 non-local correlation functional. rvv10 is the only
947
                    dispersion correction available for SCAN at this time.
948
            **kwargs: Same as those supported by DictSet.
949

950
        References:
951
            [1] P. Wisesa, K.A. McGill, T. Mueller, Efficient generation of
952
            generalized Monkhorst-Pack grids through the use of informatics,
953
            Phys. Rev. B. 93 (2016) 1-10. doi:10.1103/PhysRevB.93.155109.
954
        """
955
        super().__init__(structure, MPScanRelaxSet.CONFIG, **kwargs)
×
956
        self.bandgap = bandgap
×
957
        self.kwargs = kwargs
×
958

959
        if self.potcar_functional not in ["PBE_52", "PBE_54"]:
×
960
            raise ValueError("SCAN calculations require PBE_52 or PBE_54!")
×
961

962
        # self.kwargs.get("user_incar_settings", {
963
        updates = {}
×
964
        # select the KSPACING and smearing parameters based on the bandgap
965
        if self.bandgap == 0:
×
966
            updates["KSPACING"] = 0.22
×
967
            updates["SIGMA"] = 0.2
×
968
            updates["ISMEAR"] = 2
×
969
        else:
970
            rmin = 25.22 - 2.87 * bandgap  # Eq. 25
×
971
            kspacing = 2 * np.pi * 1.0265 / (rmin - 1.0183)  # Eq. 29
×
972
            # cap the KSPACING at a max of 0.44, per internal benchmarking
973
            if 0.22 < kspacing < 0.44:
×
974
                updates["KSPACING"] = kspacing
×
975
            else:
976
                updates["KSPACING"] = 0.44
×
977
            updates["ISMEAR"] = -5
×
978
            updates["SIGMA"] = 0.05
×
979

980
        # Don't overwrite things the user has supplied
981
        if self.user_incar_settings.get("KSPACING"):
×
982
            del updates["KSPACING"]
×
983

984
        if self.user_incar_settings.get("ISMEAR"):
×
985
            del updates["ISMEAR"]
×
986

987
        if self.user_incar_settings.get("SIGMA"):
×
988
            del updates["SIGMA"]
×
989

990
        if self.vdw:
×
991
            if self.vdw != "rvv10":
×
992
                warnings.warn(
×
993
                    "Use of van der waals functionals other than rVV10 with SCAN is not supported at this time. "
994
                )
995
                # delete any vdw parameters that may have been added to the INCAR
996
                vdw_par = loadfn(str(MODULE_DIR / "vdW_parameters.yaml"))
×
997
                for k in vdw_par[self.vdw]:
×
998
                    try:
×
999
                        del self._config_dict["INCAR"][k]
×
1000
                    except KeyError:
×
1001
                        pass
×
1002

1003
        self._config_dict["INCAR"].update(updates)
×
1004

1005

1006
class MPMetalRelaxSet(MPRelaxSet):
1✔
1007
    """
1008
    Implementation of VaspInputSet utilizing parameters in the public
1009
    Materials Project, but with tuning for metals. Key things are a denser
1010
    k point density, and a
1011
    """
1012

1013
    CONFIG = _load_yaml_config("MPRelaxSet")
1✔
1014

1015
    def __init__(self, structure: Structure, **kwargs):
1✔
1016
        """
1017
        :param structure: Structure
1018
        :param kwargs: Same as those supported by DictSet.
1019
        """
1020
        super().__init__(structure, **kwargs)
×
1021
        self._config_dict["INCAR"].update({"ISMEAR": 1, "SIGMA": 0.2})
×
1022
        self._config_dict["KPOINTS"].update({"reciprocal_density": 200})
×
1023
        self.kwargs = kwargs
×
1024

1025

1026
class MPHSERelaxSet(DictSet):
1✔
1027
    """
1028
    Same as the MPRelaxSet, but with HSE parameters.
1029
    """
1030

1031
    CONFIG = _load_yaml_config("MPHSERelaxSet")
1✔
1032

1033
    def __init__(self, structure: Structure, **kwargs):
1✔
1034
        """
1035
        :param structure: Structure
1036
        :param kwargs: Same as those supported by DictSet.
1037
        """
1038
        super().__init__(structure, MPHSERelaxSet.CONFIG, **kwargs)
×
1039
        self.kwargs = kwargs
×
1040

1041

1042
class MPStaticSet(MPRelaxSet):
1✔
1043
    """
1044
    Creates input files for a static calculation.
1045
    """
1046

1047
    def __init__(
1✔
1048
        self,
1049
        structure,
1050
        prev_incar=None,
1051
        prev_kpoints=None,
1052
        lepsilon=False,
1053
        lcalcpol=False,
1054
        reciprocal_density=100,
1055
        small_gap_multiply=None,
1056
        **kwargs,
1057
    ):
1058
        """
1059
        Args:
1060
            structure (Structure): Structure from previous run.
1061
            prev_incar (Incar): Incar file from previous run.
1062
            prev_kpoints (Kpoints): Kpoints from previous run.
1063
            lepsilon (bool): Whether to add static dielectric calculation
1064
            lcalcpol (bool): Whether to turn on evaluation of the Berry phase approximations
1065
                for electronic polarization
1066
            reciprocal_density (int): For static calculations, we usually set the
1067
                reciprocal density by volume. This is a convenience arg to change
1068
                that, rather than using user_kpoints_settings. Defaults to 100,
1069
                which is ~50% more than that of standard relaxation calculations.
1070
            small_gap_multiply ([float, float]): If the gap is less than
1071
                1st index, multiply the default reciprocal_density by the 2nd
1072
                index.
1073
            **kwargs: kwargs supported by MPRelaxSet.
1074
        """
1075
        super().__init__(structure, **kwargs)
×
1076
        if isinstance(prev_incar, str):
×
1077
            prev_incar = Incar.from_file(prev_incar)
×
1078
        if isinstance(prev_kpoints, str):
×
1079
            prev_kpoints = Kpoints.from_file(prev_kpoints)
×
1080

1081
        self.prev_incar = prev_incar
×
1082
        self.prev_kpoints = prev_kpoints
×
1083
        self.reciprocal_density = reciprocal_density
×
1084
        self.kwargs = kwargs
×
1085
        self.lepsilon = lepsilon
×
1086
        self.lcalcpol = lcalcpol
×
1087
        self.small_gap_multiply = small_gap_multiply
×
1088

1089
    @property
1✔
1090
    def incar(self):
1✔
1091
        """
1092
        :return: Incar
1093
        """
1094
        parent_incar = super().incar
×
1095
        incar = Incar(self.prev_incar) if self.prev_incar is not None else Incar(parent_incar)
×
1096

1097
        incar.update(
×
1098
            {
1099
                "IBRION": -1,
1100
                "ISMEAR": -5,
1101
                "LAECHG": True,
1102
                "LCHARG": True,
1103
                "LORBIT": 11,
1104
                "LVHAR": True,
1105
                "LWAVE": False,
1106
                "NSW": 0,
1107
                "ALGO": "Normal",
1108
            }
1109
        )
1110

1111
        if self.lepsilon:
×
1112
            incar["IBRION"] = 8
×
1113
            incar["LEPSILON"] = True
×
1114

1115
            # LPEAD=T: numerical evaluation of overlap integral prevents
1116
            # LRF_COMMUTATOR errors and can lead to better expt. agreement
1117
            # but produces slightly different results
1118
            incar["LPEAD"] = True
×
1119

1120
            # Note that DFPT calculations MUST unset NSW. NSW = 0 will fail
1121
            # to output ionic.
1122
            incar.pop("NSW", None)
×
1123
            incar.pop("NPAR", None)
×
1124

1125
            # tighter ediff for DFPT
1126
            incar["EDIFF"] = 1e-5
×
1127

1128
        if self.lcalcpol:
×
1129
            incar["LCALCPOL"] = True
×
1130

1131
        for k in ["MAGMOM", "NUPDOWN"] + list(self.user_incar_settings):
×
1132
            # For these parameters as well as user specified settings, override
1133
            # the incar settings.
1134
            if parent_incar.get(k, None) is not None:
×
1135
                incar[k] = parent_incar[k]
×
1136
            else:
1137
                incar.pop(k, None)
×
1138

1139
        # use new LDAUU when possible b/c the Poscar might have changed
1140
        # representation
1141
        if incar.get("LDAU"):
×
1142
            u = incar.get("LDAUU", [])
×
1143
            j = incar.get("LDAUJ", [])
×
1144
            if sum(u[x] - j[x] for x, y in enumerate(u)) > 0:
×
1145
                for tag in ("LDAUU", "LDAUL", "LDAUJ"):
×
1146
                    incar.update({tag: parent_incar[tag]})
×
1147
            # ensure to have LMAXMIX for GGA+U static run
1148
            if "LMAXMIX" not in incar:
×
1149
                incar.update({"LMAXMIX": parent_incar["LMAXMIX"]})
×
1150

1151
        # Compare ediff between previous and staticinputset values,
1152
        # choose the tighter ediff
1153
        incar["EDIFF"] = min(incar.get("EDIFF", 1), parent_incar["EDIFF"])
×
1154
        return incar
×
1155

1156
    @property
1✔
1157
    def kpoints(self) -> Kpoints | None:
1✔
1158
        """
1159
        :return: Kpoints
1160
        """
1161
        self._config_dict["KPOINTS"]["reciprocal_density"] = self.reciprocal_density
×
1162
        kpoints = super().kpoints
×
1163

1164
        # Prefer to use k-point scheme from previous run
1165
        # except for when lepsilon = True is specified
1166
        if kpoints is not None:
×
1167
            if self.prev_kpoints and self.prev_kpoints.style != kpoints.style:
×
1168
                if (self.prev_kpoints.style == Kpoints.supported_modes.Monkhorst) and (not self.lepsilon):
×
1169
                    k_div = [kp + 1 if kp % 2 == 1 else kp for kp in kpoints.kpts[0]]  # type: ignore
×
1170
                    kpoints = Kpoints.monkhorst_automatic(k_div)  # type: ignore
×
1171
                else:
1172
                    kpoints = Kpoints.gamma_automatic(kpoints.kpts[0])  # type: ignore
×
1173
        return kpoints
×
1174

1175
    def override_from_prev_calc(self, prev_calc_dir="."):
1✔
1176
        """
1177
        Update the input set to include settings from a previous calculation.
1178

1179
        Args:
1180
            prev_calc_dir (str): The path to the previous calculation directory.
1181

1182
        Returns:
1183
            The input set with the settings (structure, k-points, incar, etc)
1184
            updated using the previous VASP run.
1185
        """
1186
        vasprun, outcar = get_vasprun_outcar(prev_calc_dir)
×
1187

1188
        self.prev_incar = vasprun.incar
×
1189
        self.prev_kpoints = vasprun.kpoints
×
1190

1191
        if self.standardize:
×
1192
            warnings.warn(
×
1193
                "Use of standardize=True with from_prev_run is not "
1194
                "recommended as there is no guarantee the copied "
1195
                "files will be appropriate for the standardized "
1196
                "structure."
1197
            )
1198

1199
        self._structure = get_structure_from_prev_run(vasprun, outcar)
×
1200

1201
        # multiply the reciprocal density if needed
1202
        if self.small_gap_multiply:
×
1203
            gap = vasprun.eigenvalue_band_properties[0]
×
1204
            if gap <= self.small_gap_multiply[0]:
×
1205
                self.reciprocal_density = self.reciprocal_density * self.small_gap_multiply[1]
×
1206

1207
        return self
×
1208

1209
    @classmethod
1✔
1210
    def from_prev_calc(cls, prev_calc_dir, **kwargs):
1✔
1211
        """
1212
        Generate a set of VASP input files for static calculations from a
1213
        directory of previous VASP run.
1214

1215
        Args:
1216
            prev_calc_dir (str): Directory containing the outputs(
1217
                vasprun.xml and OUTCAR) of previous vasp run.
1218
            **kwargs: All kwargs supported by MPStaticSet, other than prev_incar
1219
                and prev_structure and prev_kpoints which are determined from
1220
                the prev_calc_dir.
1221
        """
1222
        input_set = cls(_dummy_structure, **kwargs)
×
1223
        return input_set.override_from_prev_calc(prev_calc_dir=prev_calc_dir)
×
1224

1225

1226
class MPScanStaticSet(MPScanRelaxSet):
1✔
1227
    """
1228
    Creates input files for a static calculation using the accurate and numerically
1229
    efficient r2SCAN variant of the Strongly Constrained and Appropriately Normed
1230
    (SCAN) metaGGA functional.
1231
    """
1232

1233
    def __init__(self, structure: Structure, bandgap=0, prev_incar=None, lepsilon=False, lcalcpol=False, **kwargs):
1✔
1234
        """
1235
        Args:
1236
            structure (Structure): Structure from previous run.
1237
            bandgap (float): Bandgap of the structure in eV. The bandgap is used to
1238
                    compute the appropriate k-point density and determine the
1239
                    smearing settings.
1240
            prev_incar (Incar): Incar file from previous run.
1241
            lepsilon (bool): Whether to add static dielectric calculation
1242
            lcalcpol (bool): Whether to turn on evaluation of the Berry phase approximations
1243
                for electronic polarization.
1244
            **kwargs: kwargs supported by MPScanRelaxSet.
1245
        """
1246
        super().__init__(structure, bandgap, **kwargs)
×
1247
        if isinstance(prev_incar, str):
×
1248
            prev_incar = Incar.from_file(prev_incar)
×
1249

1250
        self.prev_incar = prev_incar
×
1251
        self.kwargs = kwargs
×
1252
        self.lepsilon = lepsilon
×
1253
        self.lcalcpol = lcalcpol
×
1254

1255
    @property
1✔
1256
    def incar(self):
1✔
1257
        """
1258
        :return: Incar
1259
        """
1260
        parent_incar = super().incar
×
1261
        incar = Incar(self.prev_incar) if self.prev_incar is not None else Incar(parent_incar)
×
1262

1263
        incar.update({"LREAL": False, "NSW": 0, "LORBIT": 11, "LVHAR": True, "ISMEAR": -5})
×
1264

1265
        if self.lepsilon:
×
1266
            incar["IBRION"] = 8
×
1267
            incar["LEPSILON"] = True
×
1268

1269
            # LPEAD=T: numerical evaluation of overlap integral prevents
1270
            # LRF_COMMUTATOR errors and can lead to better expt. agreement
1271
            # but produces slightly different results
1272
            incar["LPEAD"] = True
×
1273

1274
            # Note that DFPT calculations MUST unset NSW. NSW = 0 will fail
1275
            # to output ionic.
1276
            incar.pop("NSW", None)
×
1277
            incar.pop("NPAR", None)
×
1278

1279
        if self.lcalcpol:
×
1280
            incar["LCALCPOL"] = True
×
1281

1282
        for k in list(self.user_incar_settings):
×
1283
            # For user specified settings, override
1284
            # the incar settings.
1285
            if parent_incar.get(k, None) is not None:
×
1286
                incar[k] = parent_incar[k]
×
1287
            else:
1288
                incar.pop(k, None)
×
1289

1290
        return incar
×
1291

1292
    def override_from_prev_calc(self, prev_calc_dir="."):
1✔
1293
        """
1294
        Update the input set to include settings from a previous calculation.
1295

1296
        Args:
1297
            prev_calc_dir (str): The path to the previous calculation directory.
1298

1299
        Returns:
1300
            The input set with the settings (structure, k-points, incar, etc)
1301
            updated using the previous VASP run.
1302
        """
1303
        vasprun, outcar = get_vasprun_outcar(prev_calc_dir)
×
1304

1305
        self.prev_incar = vasprun.incar
×
1306

1307
        self._structure = get_structure_from_prev_run(vasprun, outcar)
×
1308

1309
        return self
×
1310

1311
    @classmethod
1✔
1312
    def from_prev_calc(cls, prev_calc_dir, **kwargs):
1✔
1313
        """
1314
        Generate a set of VASP input files for static calculations from a
1315
        directory of previous VASP run.
1316

1317
        Args:
1318
            prev_calc_dir (str): Directory containing the outputs(
1319
                vasprun.xml and OUTCAR) of previous vasp run.
1320
            **kwargs: All kwargs supported by MPScanStaticSet, other than prev_incar
1321
                which is determined from the prev_calc_dir.
1322
        """
1323
        input_set = cls(_dummy_structure, **kwargs)
×
1324
        return input_set.override_from_prev_calc(prev_calc_dir=prev_calc_dir)
×
1325

1326

1327
class MPHSEBSSet(MPHSERelaxSet):
1✔
1328
    """
1329
    Implementation of a VaspInputSet for HSE band structure computations.
1330
    Remember that HSE band structures must be self-consistent in VASP. A
1331
    band structure along symmetry lines for instance needs BOTH a uniform
1332
    grid with appropriate weights AND a path along the lines with weight 0.
1333

1334
    Thus, the "Uniform" mode is just like regular static SCF but allows
1335
    adding custom kpoints (e.g., corresponding to known VBM/CBM) to the
1336
    uniform grid that have zero weight (e.g., for better gap estimate).
1337

1338
    The "Gap" mode behaves just like the "Uniform" mode, however, if starting
1339
    from a previous calculation, the VBM and CBM k-points will automatically
1340
    be added to ``added_kpoints``.
1341

1342
    The "Line" mode is just like Uniform mode, but additionally adds
1343
    k-points along symmetry lines with zero weight.
1344
    """
1345

1346
    def __init__(
1✔
1347
        self,
1348
        structure,
1349
        user_incar_settings=None,
1350
        added_kpoints=None,
1351
        mode="Gap",
1352
        reciprocal_density=None,
1353
        copy_chgcar=True,
1354
        kpoints_line_density=20,
1355
        **kwargs,
1356
    ):
1357
        """
1358
        Args:
1359
            structure (Structure): Structure to compute
1360
            user_incar_settings (dict): A dict specifying additional incar
1361
                settings
1362
            added_kpoints (list): a list of kpoints (list of 3 number list)
1363
                added to the run. The k-points are in fractional coordinates
1364
            mode (str): "Line" - generate k-points along symmetry lines for
1365
                bandstructure. "Uniform" - generate uniform k-points grid.
1366
            reciprocal_density (int): k-point density to use for uniform mesh.
1367
            copy_chgcar (bool): Whether to copy the CHGCAR of a previous run.
1368
            kpoints_line_density (int): k-point density for high symmetry lines
1369
            **kwargs (dict): Any other parameters to pass into DictSet.
1370
        """
1371
        super().__init__(structure, **kwargs)
×
1372
        self.user_incar_settings = user_incar_settings or {}
×
1373
        self._config_dict["INCAR"].update(
×
1374
            {
1375
                "NSW": 0,
1376
                "ISMEAR": 0,
1377
                "SIGMA": 0.05,
1378
                "ISYM": 3,
1379
                "LCHARG": False,
1380
                "NELMIN": 5,
1381
            }
1382
        )
1383
        self.added_kpoints = added_kpoints if added_kpoints is not None else []
×
1384
        self.mode = mode
×
1385

1386
        if not reciprocal_density or "reciprocal_density" not in self.user_kpoints_settings:
×
1387
            self.reciprocal_density = 50
×
1388
        else:
1389
            self.reciprocal_density = reciprocal_density or self.user_kpoints_settings["reciprocal_density"]
×
1390

1391
        self.kpoints_line_density = kpoints_line_density
×
1392
        self.copy_chgcar = copy_chgcar
×
1393

1394
    @property
1✔
1395
    def kpoints(self) -> Kpoints:
1✔
1396
        """
1397
        :return: Kpoints
1398
        """
1399
        kpts: list[int | float | None] = []
×
1400
        weights: list[float | None] = []
×
1401
        all_labels: list[str | None] = []
×
1402
        structure = self.structure
×
1403

1404
        # for both modes, include the Uniform mesh w/standard weights
1405
        grid = Kpoints.automatic_density_by_vol(structure, self.reciprocal_density).kpts
×
1406
        ir_kpts = SpacegroupAnalyzer(structure, symprec=0.1).get_ir_reciprocal_mesh(grid[0])
×
1407
        for k in ir_kpts:
×
1408
            kpts.append(k[0])
×
1409
            weights.append(int(k[1]))
×
1410
            all_labels.append(None)
×
1411

1412
        # for both modes, include any user-added kpoints w/zero weight
1413
        for k in self.added_kpoints:
×
1414
            kpts.append(k)
×
1415
            weights.append(0.0)
×
1416
            all_labels.append("user-defined")
×
1417

1418
        # for line mode only, add the symmetry lines w/zero weight
1419
        if self.mode.lower() == "line":
×
1420
            kpath = HighSymmKpath(structure)
×
1421
            frac_k_points, labels = kpath.get_kpoints(
×
1422
                line_density=self.kpoints_line_density, coords_are_cartesian=False
1423
            )
1424

1425
            for k, f in enumerate(frac_k_points):
×
1426
                kpts.append(f)
×
1427
                weights.append(0.0)
×
1428
                all_labels.append(labels[k])
×
1429

1430
        comment = "HSE run along symmetry lines" if self.mode.lower() == "line" else "HSE run on uniform grid"
×
1431

1432
        return Kpoints(
×
1433
            comment=comment,
1434
            style=Kpoints.supported_modes.Reciprocal,
1435
            num_kpts=len(kpts),
1436
            kpts=kpts,  # type: ignore
1437
            kpts_weights=weights,
1438
            labels=all_labels,
1439
        )
1440

1441
    def override_from_prev_calc(self, prev_calc_dir="."):
1✔
1442
        """
1443
        Update the input set to include settings from a previous calculation.
1444

1445
        Args:
1446
            prev_calc_dir (str): The path to the previous calculation directory.
1447

1448
        Returns:
1449
            The input set with the settings (structure, k-points, incar, etc)
1450
            updated using the previous VASP run.
1451
        """
1452
        vasprun, outcar = get_vasprun_outcar(prev_calc_dir)
×
1453

1454
        self._structure = get_structure_from_prev_run(vasprun, outcar)
×
1455

1456
        # note: recommend not standardizing the cell because we want to retain
1457
        # k-points
1458
        if self.standardize:
×
1459
            warnings.warn(
×
1460
                "Use of standardize=True with from_prev_calc is not "
1461
                "recommended as there is no guarantee the copied "
1462
                "files will be appropriate for the standardized "
1463
                "structure."
1464
            )
1465

1466
        if self.mode.lower() == "gap":
×
1467
            added_kpoints = []
×
1468

1469
            bs = vasprun.get_band_structure()
×
1470
            vbm, cbm = bs.get_vbm()["kpoint"], bs.get_cbm()["kpoint"]
×
1471
            if vbm:
×
1472
                added_kpoints.append(vbm.frac_coords)
×
1473
            if cbm:
×
1474
                added_kpoints.append(cbm.frac_coords)
×
1475

1476
            self.added_kpoints.extend(added_kpoints)
×
1477

1478
        files_to_transfer = {}
×
1479
        if self.copy_chgcar:
×
1480
            chgcars = sorted(glob.glob(str(Path(prev_calc_dir) / "CHGCAR*")))
×
1481
            if chgcars:
×
1482
                files_to_transfer["CHGCAR"] = str(chgcars[-1])
×
1483

1484
        self.files_to_transfer.update(files_to_transfer)
×
1485

1486
        return self
×
1487

1488
    @classmethod
1✔
1489
    def from_prev_calc(cls, prev_calc_dir, **kwargs):
1✔
1490
        """
1491
        Generate a set of VASP input files for HSE calculations from a
1492
        directory of previous VASP run.
1493

1494
        Args:
1495
            prev_calc_dir (str): Directory containing the outputs
1496
                (vasprun.xml and OUTCAR) of previous vasp run.
1497
            **kwargs: All kwargs supported by MPHSEBSStaticSet, other than
1498
                prev_structure which is determined from the previous calc dir.
1499
        """
1500
        input_set = cls(_dummy_structure, **kwargs)
×
1501
        return input_set.override_from_prev_calc(prev_calc_dir=prev_calc_dir)
×
1502

1503

1504
class MPNonSCFSet(MPRelaxSet):
1✔
1505
    """
1506
    Init a MPNonSCFSet. Typically, you would use the classmethod
1507
    from_prev_calc to initialize from a previous SCF run.
1508
    """
1509

1510
    def __init__(
1✔
1511
        self,
1512
        structure,
1513
        prev_incar=None,
1514
        mode="line",
1515
        nedos=2001,
1516
        dedos=0.005,
1517
        reciprocal_density=100,
1518
        sym_prec=0.1,
1519
        kpoints_line_density=20,
1520
        optics=False,
1521
        copy_chgcar=True,
1522
        nbands_factor=1.2,
1523
        small_gap_multiply=None,
1524
        **kwargs,
1525
    ):
1526
        """
1527
        Args:
1528
            structure (Structure): Structure to compute
1529
            prev_incar (Incar/string): Incar file from previous run.
1530
            mode (str): Line, Uniform or Boltztrap mode supported.
1531
            nedos (int): nedos parameter. Default to 2001.
1532
            dedos (float): setting nedos=0 and uniform mode in from_prev_calc,
1533
                an automatic nedos will be calculated using the total energy range
1534
                divided by the energy step dedos
1535
            reciprocal_density (int): density of k-mesh by reciprocal
1536
                volume (defaults to 100)
1537
            sym_prec (float): Symmetry precision (for Uniform mode).
1538
            kpoints_line_density (int): Line density for Line mode.
1539
            optics (bool): whether to add dielectric function
1540
            copy_chgcar: Whether to copy the old CHGCAR when starting from a
1541
                previous calculation.
1542
            nbands_factor (float): Multiplicative factor for NBANDS when starting
1543
                from a previous calculation. Choose a higher number if you are
1544
                doing an LOPTICS calculation.
1545
            small_gap_multiply ([float, float]): When starting from a previous
1546
                calculation, if the gap is less than 1st index, multiply the default
1547
                reciprocal_density by the 2nd index.
1548
            **kwargs: kwargs supported by MPRelaxSet.
1549
        """
1550
        super().__init__(structure, **kwargs)
×
1551
        if isinstance(prev_incar, str):
×
1552
            prev_incar = Incar.from_file(prev_incar)
×
1553
        self.prev_incar = prev_incar
×
1554
        self.kwargs = kwargs
×
1555
        self.nedos = nedos
×
1556
        self.dedos = dedos
×
1557
        self.reciprocal_density = reciprocal_density
×
1558
        self.sym_prec = sym_prec
×
1559
        self.kpoints_line_density = kpoints_line_density
×
1560
        self.optics = optics
×
1561
        self.mode = mode.lower()
×
1562
        self.copy_chgcar = copy_chgcar
×
1563
        self.nbands_factor = nbands_factor
×
1564
        self.small_gap_multiply = small_gap_multiply
×
1565

1566
        if self.mode.lower() not in ["line", "uniform", "boltztrap"]:
×
1567
            raise ValueError("Supported modes for NonSCF runs are 'Line', 'Uniform' and 'Boltztrap!")
×
1568

1569
        if (self.mode.lower() != "uniform" or nedos < 2000) and optics:
×
1570
            warnings.warn("It is recommended to use Uniform mode with a high NEDOS for optics calculations.")
×
1571

1572
    @property
1✔
1573
    def incar(self) -> Incar:
1✔
1574
        """
1575
        :return: Incar
1576
        """
1577
        incar = super().incar
×
1578
        if self.prev_incar is not None:
×
1579
            incar.update(self.prev_incar.items())
×
1580

1581
        # Overwrite necessary INCAR parameters from previous runs
1582
        incar.update(
×
1583
            {
1584
                "IBRION": -1,
1585
                "LCHARG": False,
1586
                "LORBIT": 11,
1587
                "LWAVE": False,
1588
                "NSW": 0,
1589
                "ISYM": 0,
1590
                "ICHARG": 11,
1591
            }
1592
        )
1593

1594
        if self.mode.lower() == "uniform":
×
1595
            # use tetrahedron method for DOS and optics calculations
1596
            incar.update({"ISMEAR": -5, "ISYM": 2})
×
1597
        else:
1598
            # if line mode, can't use ISMEAR=-5; also use small sigma to avoid
1599
            # partial occupancies for small band gap materials.
1600
            # finally, explicit k-point generation (needed for bolztrap mode)
1601
            # is incompatible with ISMEAR = -5.
1602
            incar.update({"ISMEAR": 0, "SIGMA": 0.01})
×
1603

1604
        incar.update(self.user_incar_settings)
×
1605

1606
        if self.mode.lower() in "uniform":
×
1607
            # Set smaller steps for DOS and optics output
1608
            incar["NEDOS"] = self.nedos
×
1609

1610
        if self.optics:
×
1611
            incar["LOPTICS"] = True
×
1612

1613
        incar.pop("MAGMOM", None)
×
1614

1615
        return incar
×
1616

1617
    @property
1✔
1618
    def kpoints(self) -> Kpoints | None:
1✔
1619
        """
1620
        :return: Kpoints
1621
        """
1622
        # override pymatgen kpoints if provided
1623
        user_kpoints = self.user_kpoints_settings
×
1624
        if isinstance(user_kpoints, Kpoints):
×
1625
            return user_kpoints
×
1626

1627
        if self.mode.lower() == "line":
×
1628
            kpath = HighSymmKpath(self.structure)
×
1629
            frac_k_points, k_points_labels = kpath.get_kpoints(
×
1630
                line_density=self.kpoints_line_density, coords_are_cartesian=False
1631
            )
1632
            kpoints = Kpoints(
×
1633
                comment="Non SCF run along symmetry lines",
1634
                style=Kpoints.supported_modes.Reciprocal,
1635
                num_kpts=len(frac_k_points),
1636
                kpts=frac_k_points,
1637
                labels=k_points_labels,
1638
                kpts_weights=[1] * len(frac_k_points),
1639
            )
1640
        elif self.mode.lower() == "boltztrap":
×
1641
            kpoints = Kpoints.automatic_density_by_vol(self.structure, self.reciprocal_density)
×
1642
            mesh = kpoints.kpts[0]
×
1643
            ir_kpts = SpacegroupAnalyzer(self.structure, symprec=self.sym_prec).get_ir_reciprocal_mesh(mesh)
×
1644
            kpts = []
×
1645
            weights = []
×
1646
            for k in ir_kpts:
×
1647
                kpts.append(k[0])
×
1648
                weights.append(int(k[1]))
×
1649
            kpoints = Kpoints(
×
1650
                comment="Non SCF run on uniform grid",
1651
                style=Kpoints.supported_modes.Reciprocal,
1652
                num_kpts=len(ir_kpts),
1653
                kpts=kpts,
1654
                kpts_weights=weights,
1655
            )
1656
        else:
1657
            self._config_dict["KPOINTS"]["reciprocal_density"] = self.reciprocal_density
×
1658
            return super().kpoints
×
1659

1660
        return kpoints
×
1661

1662
    def override_from_prev_calc(self, prev_calc_dir="."):
1✔
1663
        """
1664
        Update the input set to include settings from a previous calculation.
1665

1666
        Args:
1667
            prev_calc_dir (str): The path to the previous calculation directory.
1668

1669
        Returns:
1670
            The input set with the settings (structure, k-points, incar, etc)
1671
            updated using the previous VASP run.
1672
        """
1673
        vasprun, outcar = get_vasprun_outcar(prev_calc_dir)
×
1674

1675
        self.prev_incar = vasprun.incar
×
1676

1677
        # Get a Magmom-decorated structure
1678
        self._structure = get_structure_from_prev_run(vasprun, outcar)
×
1679

1680
        if self.standardize:
×
1681
            warnings.warn(
×
1682
                "Use of standardize=True with from_prev_run is not "
1683
                "recommended as there is no guarantee the copied "
1684
                "files will be appropriate for the standardized"
1685
                " structure. copy_chgcar is enforced to be false."
1686
            )
1687
            self.copy_chgcar = False
×
1688

1689
        # Turn off spin when magmom for every site is smaller than 0.02.
1690
        if outcar and outcar.magnetization:
×
1691
            site_magmom = np.array([i["tot"] for i in outcar.magnetization])
×
1692
            ispin = 2 if np.any(site_magmom[np.abs(site_magmom) > 0.02]) else 1
×
1693

1694
        elif vasprun.is_spin:
×
1695
            ispin = 2
×
1696

1697
        else:
1698
            ispin = 1
×
1699

1700
        nbands = int(np.ceil(vasprun.parameters["NBANDS"] * self.nbands_factor))
×
1701
        self.prev_incar.update({"ISPIN": ispin, "NBANDS": nbands})
×
1702

1703
        files_to_transfer = {}
×
1704

1705
        if self.copy_chgcar:
×
1706
            chgcars = sorted(glob.glob(str(Path(prev_calc_dir) / "CHGCAR*")))
×
1707
            if chgcars:
×
1708
                files_to_transfer["CHGCAR"] = str(chgcars[-1])
×
1709

1710
        self.files_to_transfer.update(files_to_transfer)
×
1711

1712
        # multiply the reciprocal density if needed:
1713
        if self.small_gap_multiply:
×
1714
            gap = vasprun.eigenvalue_band_properties[0]
×
1715
            if gap <= self.small_gap_multiply[0]:
×
1716
                self.reciprocal_density = self.reciprocal_density * self.small_gap_multiply[1]
×
1717
                self.kpoints_line_density = self.kpoints_line_density * self.small_gap_multiply[1]
×
1718

1719
        # automatic setting of nedos using the energy range and the energy step dedos
1720
        if self.nedos == 0:
×
1721
            emax = max(eigs.max() for eigs in vasprun.eigenvalues.values())
×
1722
            emin = min(eigs.min() for eigs in vasprun.eigenvalues.values())
×
1723
            self.nedos = int((emax - emin) / self.dedos)
×
1724

1725
        return self
×
1726

1727
    @classmethod
1✔
1728
    def from_prev_calc(cls, prev_calc_dir, **kwargs):
1✔
1729
        """
1730
        Generate a set of VASP input files for NonSCF calculations from a
1731
        directory of previous static VASP run.
1732

1733
        Args:
1734
            prev_calc_dir (str): The directory contains the outputs(
1735
                vasprun.xml and OUTCAR) of previous vasp run.
1736
            **kwargs: All kwargs supported by MPNonSCFSet, other than structure,
1737
                prev_incar and prev_chgcar which are determined from the
1738
                prev_calc_dir.
1739
        """
1740
        input_set = cls(_dummy_structure, **kwargs)
×
1741
        return input_set.override_from_prev_calc(prev_calc_dir=prev_calc_dir)
×
1742

1743

1744
class MPSOCSet(MPStaticSet):
1✔
1745
    """
1746
    An input set for running spin-orbit coupling (SOC) calculations.
1747
    """
1748

1749
    def __init__(
1✔
1750
        self,
1751
        structure,
1752
        saxis=(0, 0, 1),
1753
        copy_chgcar=True,
1754
        nbands_factor=1.2,
1755
        reciprocal_density=100,
1756
        small_gap_multiply=None,
1757
        magmom=None,
1758
        **kwargs,
1759
    ):
1760
        """
1761
        Args:
1762
            structure (Structure): the structure must have the 'magmom' site
1763
                property and each magnetic moment value must have 3
1764
                components. eg: ``magmom = [[0,0,2], ...]``
1765
            saxis (tuple): magnetic moment orientation
1766
            copy_chgcar: Whether to copy the old CHGCAR. Defaults to True.
1767
            nbands_factor (float): Multiplicative factor for NBANDS. Choose a
1768
                higher number if you are doing an LOPTICS calculation.
1769
            reciprocal_density (int): density of k-mesh by reciprocal volume.
1770
            small_gap_multiply ([float, float]): If the gap is less than
1771
                1st index, multiply the default reciprocal_density by the 2nd
1772
                index.
1773
            magmom (list[list[float]]): Override for the structure magmoms.
1774
            **kwargs: kwargs supported by MPStaticSet.
1775
        """
1776
        if not hasattr(structure[0], "magmom") and not isinstance(structure[0].magmom, list):
×
1777
            raise ValueError(
×
1778
                "The structure must have the 'magmom' site "
1779
                "property and each magnetic moment value must have 3 "
1780
                "components. eg:- magmom = [0,0,2]"
1781
            )
1782

1783
        super().__init__(structure, reciprocal_density=reciprocal_density, **kwargs)
×
1784
        self.saxis = saxis
×
1785
        self.copy_chgcar = copy_chgcar
×
1786
        self.nbands_factor = nbands_factor
×
1787
        self.small_gap_multiply = small_gap_multiply
×
1788
        self.magmom = magmom
×
1789

1790
    @property
1✔
1791
    def incar(self) -> Incar:
1✔
1792
        """
1793
        :return: Incar
1794
        """
1795
        incar = super().incar
×
1796
        if self.prev_incar is not None:
×
1797
            incar.update(self.prev_incar.items())
×
1798

1799
        # Overwrite necessary INCAR parameters from previous runs
1800
        incar.update({"ISYM": -1, "LSORBIT": "T", "ICHARG": 11, "SAXIS": list(self.saxis)})
×
1801
        incar.update(self.user_incar_settings)
×
1802

1803
        return incar
×
1804

1805
    def override_from_prev_calc(self, prev_calc_dir="."):
1✔
1806
        """
1807
        Update the input set to include settings from a previous calculation.
1808

1809
        Args:
1810
            prev_calc_dir (str): The path to the previous calculation directory.
1811

1812
        Returns:
1813
            The input set with the settings (structure, k-points, incar, etc)
1814
            updated using the previous VASP run.
1815
        """
1816
        vasprun, outcar = get_vasprun_outcar(prev_calc_dir)
×
1817

1818
        self.prev_incar = vasprun.incar
×
1819

1820
        # Remove magmoms from previous INCAR, since we will prefer
1821
        # the final calculated magmoms
1822
        # TODO: revisit in context of MPStaticSet incar logic
1823
        if "MAGMOM" in self.prev_incar:
×
1824
            del self.prev_incar["magmom"]
×
1825

1826
        # Get a magmom-decorated structure
1827
        self._structure = get_structure_from_prev_run(vasprun, outcar)
×
1828
        if self.standardize:
×
1829
            warnings.warn(
×
1830
                "Use of standardize=True with from_prev_run is not "
1831
                "recommended as there is no guarantee the copied "
1832
                "files will be appropriate for the standardized"
1833
                " structure. copy_chgcar is enforced to be false."
1834
            )
1835
            self.copy_chgcar = False
×
1836

1837
        # override magmom if provided
1838
        if self.magmom:
×
1839
            self._structure = self._structure.copy(site_properties={"magmom": self.magmom})
×
1840

1841
        # magmom has to be 3D for SOC calculation.
1842
        if hasattr(self._structure[0], "magmom"):
×
1843
            if not isinstance(self._structure[0].magmom, list):
×
1844
                self._structure = self._structure.copy(
×
1845
                    site_properties={"magmom": [[0, 0, site.magmom] for site in self._structure]}
1846
                )
1847
        else:
1848
            raise ValueError("Neither the previous structure has magmom property nor magmom provided")
×
1849

1850
        nbands = int(np.ceil(vasprun.parameters["NBANDS"] * self.nbands_factor))
×
1851
        self.prev_incar.update({"NBANDS": nbands})
×
1852

1853
        files_to_transfer = {}
×
1854
        if self.copy_chgcar:
×
1855
            chgcars = sorted(glob.glob(str(Path(prev_calc_dir) / "CHGCAR*")))
×
1856
            if chgcars:
×
1857
                files_to_transfer["CHGCAR"] = str(chgcars[-1])
×
1858

1859
        self.files_to_transfer.update(files_to_transfer)
×
1860

1861
        # multiply the reciprocal density if needed:
1862
        if self.small_gap_multiply:
×
1863
            gap = vasprun.eigenvalue_band_properties[0]
×
1864
            if gap <= self.small_gap_multiply[0]:
×
1865
                self.reciprocal_density = self.reciprocal_density * self.small_gap_multiply[1]
×
1866

1867
        return self
×
1868

1869
    @classmethod
1✔
1870
    def from_prev_calc(cls, prev_calc_dir, **kwargs):
1✔
1871
        """
1872
        Generate a set of VASP input files for SOC calculations from a
1873
        directory of previous static VASP run. SOC calc requires all 3
1874
        components for MAGMOM for each atom in the structure.
1875

1876
        Args:
1877
            prev_calc_dir (str): The directory contains the outputs(
1878
                vasprun.xml and OUTCAR) of previous vasp run.
1879
            **kwargs: All kwargs supported by MPSOCSet, other than structure,
1880
                prev_incar and prev_chgcar which are determined from the
1881
                prev_calc_dir.
1882
        """
1883
        input_set = cls(_dummy_structure, **kwargs)
×
1884
        return input_set.override_from_prev_calc(prev_calc_dir=prev_calc_dir)
×
1885

1886

1887
class MPNMRSet(MPStaticSet):
1✔
1888
    """
1889
    Init a MPNMRSet.
1890
    """
1891

1892
    def __init__(
1✔
1893
        self, structure: Structure, mode="cs", isotopes=None, prev_incar=None, reciprocal_density=100, **kwargs
1894
    ):
1895
        """
1896
        Args:
1897
            structure (Structure): Structure to compute
1898
            mode (str): The NMR calculation to run
1899
                            "cs": for Chemical Shift
1900
                            "efg" for Electric Field Gradient
1901
            isotopes (list): list of Isotopes for quadrupole moments
1902
            prev_incar (Incar): Incar file from previous run.
1903
            reciprocal_density (int): density of k-mesh by reciprocal
1904
                                    volume (defaults to 100)
1905
            **kwargs: kwargs supported by MPStaticSet.
1906
        """
1907
        self.mode = mode
×
1908
        self.isotopes = isotopes or []
×
1909
        super().__init__(structure, prev_incar=prev_incar, reciprocal_density=reciprocal_density, **kwargs)
×
1910

1911
    @property
1✔
1912
    def incar(self):
1✔
1913
        """
1914
        :return: Incar
1915
        """
1916
        incar = super().incar
×
1917

1918
        if self.mode.lower() == "cs":
×
1919
            incar.update(
×
1920
                {
1921
                    "LCHIMAG": True,
1922
                    "EDIFF": -1.0e-10,
1923
                    "ISYM": 0,
1924
                    "LCHARG": False,
1925
                    "LNMR_SYM_RED": True,
1926
                    "NELMIN": 10,
1927
                    "NSLPLINE": True,
1928
                    "PREC": "ACCURATE",
1929
                    "SIGMA": 0.01,
1930
                }
1931
            )
1932
        elif self.mode.lower() == "efg":
×
1933
            isotopes = {ist.split("-")[0]: ist for ist in self.isotopes}
×
1934

1935
            quad_efg = [
×
1936
                float(Species(p).get_nmr_quadrupole_moment(isotopes.get(p, None))) for p in self.poscar.site_symbols
1937
            ]
1938

1939
            incar.update(
×
1940
                {
1941
                    "ALGO": "FAST",
1942
                    "EDIFF": -1.0e-10,
1943
                    "ISYM": 0,
1944
                    "LCHARG": False,
1945
                    "LEFG": True,
1946
                    "QUAD_EFG": quad_efg,
1947
                    "NELMIN": 10,
1948
                    "PREC": "ACCURATE",
1949
                    "SIGMA": 0.01,
1950
                }
1951
            )
1952
        incar.update(self.user_incar_settings)
×
1953

1954
        return incar
×
1955

1956

1957
class MVLElasticSet(MPRelaxSet):
1✔
1958
    """
1959
    MVL denotes VASP input sets that are implemented by the Materials Virtual
1960
    Lab (http://materialsvirtuallab.org) for various research.
1961

1962
    This input set is used to calculate elastic constants in VASP. It is used
1963
    in the following work::
1964

1965
        Z. Deng, Z. Wang, I.-H. Chu, J. Luo, S. P. Ong.
1966
        “Elastic Properties of Alkali Superionic Conductor Electrolytes
1967
        from First Principles Calculations”, J. Electrochem. Soc.
1968
        2016, 163(2), A67-A74. doi: 10.1149/2.0061602jes
1969

1970
    To read the elastic constants, you may use the Outcar class which parses the
1971
    elastic constants.
1972
    """
1973

1974
    def __init__(self, structure: Structure, potim=0.015, **kwargs):
1✔
1975
        """
1976
        Args:
1977
            scale (float): POTIM parameter. The default of 0.015 is usually fine,
1978
                but some structures may require a smaller step.
1979
            user_incar_settings (dict): A dict specifying additional incar
1980
                settings.
1981
            kwargs:
1982
                Parameters supported by MPRelaxSet.
1983
        """
1984
        super().__init__(structure, **kwargs)
×
1985
        self._config_dict["INCAR"].update({"IBRION": 6, "NFREE": 2, "POTIM": potim})
×
1986
        self._config_dict["INCAR"].pop("NPAR", None)
×
1987

1988

1989
class MVLGWSet(DictSet):
1✔
1990
    """
1991
    MVL denotes VASP input sets that are implemented by the Materials Virtual
1992
    Lab (http://materialsvirtuallab.org) for various research. This is a
1993
    flexible input set for GW calculations.
1994

1995
    Note that unlike all other input sets in this module, the PBE_54 series of
1996
    functional is set as the default. These have much improved performance for
1997
    GW calculations.
1998

1999
    A typical sequence is mode="STATIC" -> mode="DIAG" -> mode="GW" ->
2000
    mode="BSE". For all steps other than the first one (static), the
2001
    recommendation is to use from_prev_calculation on the preceding run in
2002
    the series.
2003
    """
2004

2005
    CONFIG = _load_yaml_config("MVLGWSet")
1✔
2006

2007
    SUPPORTED_MODES = ("DIAG", "GW", "STATIC", "BSE")
1✔
2008

2009
    def __init__(
1✔
2010
        self,
2011
        structure,
2012
        prev_incar=None,
2013
        nbands=None,
2014
        reciprocal_density=100,
2015
        mode="STATIC",
2016
        copy_wavecar=True,
2017
        nbands_factor=5,
2018
        ncores=16,
2019
        **kwargs,
2020
    ):
2021
        """
2022
        Args:
2023
            structure (Structure): Input structure.
2024
            prev_incar (Incar/string): Incar file from previous run.
2025
            mode (str): Supported modes are "STATIC" (default), "DIAG", "GW",
2026
                and "BSE".
2027
            nbands (int): For subsequent calculations, it is generally
2028
                recommended to perform NBANDS convergence starting from the
2029
                NBANDS of the previous run for DIAG, and to use the exact same
2030
                NBANDS for GW and BSE. This parameter is used by
2031
                from_previous_calculation to set nband.
2032
            copy_wavecar: Whether to copy the old WAVECAR, WAVEDER and associated
2033
                files when starting from a previous calculation.
2034
            nbands_factor (int): Multiplicative factor for NBANDS when starting
2035
                from a previous calculation. Only applies if mode=="DIAG".
2036
                Need to be tested for convergence.
2037
            ncores (int): Numbers of cores used for the calculation. VASP will alter
2038
                NBANDS if it was not dividable by ncores. Only applies if
2039
                mode=="DIAG".
2040
            **kwargs: All kwargs supported by DictSet. Typically,
2041
                user_incar_settings is a commonly used option.
2042
        """
2043
        super().__init__(structure, MVLGWSet.CONFIG, **kwargs)
×
2044
        self.prev_incar = prev_incar
×
2045
        self.nbands = nbands
×
2046
        self.reciprocal_density = reciprocal_density
×
2047
        self.mode = mode.upper()
×
2048
        if self.mode not in MVLGWSet.SUPPORTED_MODES:
×
2049
            raise ValueError(f"{self.mode} not one of the support modes : {MVLGWSet.SUPPORTED_MODES}")
×
2050
        self.kwargs = kwargs
×
2051
        self.copy_wavecar = copy_wavecar
×
2052
        self.nbands_factor = nbands_factor
×
2053
        self.ncores = ncores
×
2054

2055
    @property
1✔
2056
    def kpoints(self):
1✔
2057
        """
2058
        Generate gamma center k-points mesh grid for GW calc,
2059
        which is requested by GW calculation.
2060
        """
2061
        return Kpoints.automatic_density_by_vol(self.structure, self.reciprocal_density, force_gamma=True)
×
2062

2063
    @property
1✔
2064
    def incar(self):
1✔
2065
        """
2066
        :return: Incar
2067
        """
2068
        parent_incar = super().incar
×
2069
        incar = Incar(self.prev_incar) if self.prev_incar is not None else Incar(parent_incar)
×
2070

2071
        if self.mode == "DIAG":
×
2072
            # Default parameters for diagonalization calculation.
2073
            incar.update({"ALGO": "Exact", "NELM": 1, "LOPTICS": True, "LPEAD": True})
×
2074
        elif self.mode == "GW":
×
2075
            # Default parameters for GW calculation.
2076
            incar.update({"ALGO": "GW0", "NELM": 1, "NOMEGA": 80, "ENCUTGW": 250})
×
2077
            incar.pop("EDIFF", None)
×
2078
            incar.pop("LOPTICS", None)
×
2079
            incar.pop("LPEAD", None)
×
2080
        elif self.mode == "BSE":
×
2081
            # Default parameters for BSE calculation.
2082
            incar.update({"ALGO": "BSE", "ANTIRES": 0, "NBANDSO": 20, "NBANDSV": 20})
×
2083

2084
        if self.nbands:
×
2085
            incar["NBANDS"] = self.nbands
×
2086

2087
        # Respect user set INCAR.
2088
        incar.update(self.kwargs.get("user_incar_settings", {}))
×
2089

2090
        return incar
×
2091

2092
    def override_from_prev_calc(self, prev_calc_dir="."):
1✔
2093
        """
2094
        Update the input set to include settings from a previous calculation.
2095

2096
        Args:
2097
            prev_calc_dir (str): The path to the previous calculation directory.
2098

2099
        Returns:
2100
            The input set with the settings (structure, k-points, incar, etc)
2101
            updated using the previous VASP run.
2102
        """
2103
        vasprun, outcar = get_vasprun_outcar(prev_calc_dir)
×
2104
        self.prev_incar = vasprun.incar
×
2105
        self._structure = vasprun.final_structure
×
2106

2107
        if self.standardize:
×
2108
            warnings.warn(
×
2109
                "Use of standardize=True with from_prev_run is not "
2110
                "recommended as there is no guarantee the copied "
2111
                "files will be appropriate for the standardized "
2112
                "structure."
2113
            )
2114

2115
        self.nbands = int(vasprun.parameters["NBANDS"])
×
2116
        if self.mode.upper() == "DIAG":
×
2117
            self.nbands = int(np.ceil(self.nbands * self.nbands_factor / self.ncores) * self.ncores)
×
2118

2119
        # copy WAVECAR, WAVEDER (derivatives)
2120
        files_to_transfer = {}
×
2121
        if self.copy_wavecar:
×
2122
            for fname in ("WAVECAR", "WAVEDER", "WFULL"):
×
2123
                w = sorted(glob.glob(str(Path(prev_calc_dir) / (fname + "*"))))
×
2124
                if w:
×
2125
                    if fname == "WFULL":
×
2126
                        for f in w:
×
2127
                            fname = Path(f).name
×
2128
                            fname = fname.split(".")[0]
×
2129
                            files_to_transfer[fname] = f
×
2130
                    else:
2131
                        files_to_transfer[fname] = str(w[-1])
×
2132

2133
        self.files_to_transfer.update(files_to_transfer)
×
2134

2135
        return self
×
2136

2137
    @classmethod
1✔
2138
    def from_prev_calc(cls, prev_calc_dir, mode="DIAG", **kwargs):
1✔
2139
        """
2140
        Generate a set of VASP input files for GW or BSE calculations from a
2141
        directory of previous Exact Diag VASP run.
2142

2143
        Args:
2144
            prev_calc_dir (str): The directory contains the outputs(
2145
                vasprun.xml of previous vasp run.
2146
            mode (str): Supported modes are "STATIC", "DIAG" (default), "GW",
2147
                and "BSE".
2148
            **kwargs: All kwargs supported by MVLGWSet, other than structure,
2149
                prev_incar and mode, which are determined from the
2150
                prev_calc_dir.
2151
        """
2152
        input_set = cls(_dummy_structure, mode=mode, **kwargs)
×
2153
        return input_set.override_from_prev_calc(prev_calc_dir=prev_calc_dir)
×
2154

2155

2156
class MVLSlabSet(MPRelaxSet):
1✔
2157
    """
2158
    Class for writing a set of slab vasp runs,
2159
    including both slabs (along the c direction) and orient unit cells (bulk),
2160
    to ensure the same KPOINTS, POTCAR and INCAR criterion.
2161
    """
2162

2163
    def __init__(
1✔
2164
        self,
2165
        structure: Structure,
2166
        k_product=50,
2167
        bulk=False,
2168
        auto_dipole=False,
2169
        set_mix=True,
2170
        sort_structure=True,
2171
        **kwargs,
2172
    ):
2173
        """
2174
        :param structure: Structure
2175
        :param k_product: default to 50, kpoint number * length for a & b
2176
            directions, also for c direction in bulk calculations
2177
        :param bulk:
2178
        :param auto_dipole:
2179
        :param set_mix:
2180
        :param sort_structure:
2181
        :param kwargs: Other kwargs supported by :class:`DictSet`.
2182
        """
2183
        super().__init__(structure, **kwargs)
×
2184

2185
        if sort_structure:
×
2186
            structure = structure.get_sorted_structure()
×
2187

2188
        self.k_product = k_product
×
2189
        self.bulk = bulk
×
2190
        self.auto_dipole = auto_dipole
×
2191
        self.kwargs = kwargs
×
2192
        self.set_mix = set_mix
×
2193
        self.kpt_calc = None
×
2194

2195
        slab_incar = {
×
2196
            "EDIFF": 1e-4,
2197
            "EDIFFG": -0.02,
2198
            "ENCUT": 400,
2199
            "ISMEAR": 0,
2200
            "SIGMA": 0.05,
2201
            "ISIF": 3,
2202
        }
2203
        if not self.bulk:
×
2204
            slab_incar["ISIF"] = 2
×
2205
            slab_incar["LVTOT"] = True
×
2206
            if self.set_mix:
×
2207
                slab_incar["AMIN"] = 0.01
×
2208
                slab_incar["AMIX"] = 0.2
×
2209
                slab_incar["BMIX"] = 0.001
×
2210
            slab_incar["NELMIN"] = 8
×
2211
            if self.auto_dipole:
×
2212
                weights = [s.species.weight for s in structure]
×
2213
                center_of_mass = np.average(structure.frac_coords, weights=weights, axis=0)
×
2214

2215
                slab_incar["IDIPOL"] = 3
×
2216
                slab_incar["LDIPOL"] = True
×
2217
                slab_incar["DIPOL"] = center_of_mass
×
2218

2219
        self._config_dict["INCAR"].update(slab_incar)
×
2220

2221
    @property
1✔
2222
    def kpoints(self):
1✔
2223
        """
2224
        k_product, default to 50, is kpoint number * length for a & b
2225
            directions, also for c direction in bulk calculations
2226
        Automatic mesh & Gamma is the default setting.
2227
        """
2228
        # To get input sets, the input structure has to has the same number
2229
        # of required parameters as a Structure object (ie. 4). Slab
2230
        # attributes aren't going to affect the VASP inputs anyways so
2231
        # converting the slab into a structure should not matter
2232

2233
        kpt = super().kpoints
×
2234
        kpt.comment = "Automatic mesh"
×
2235
        kpt.style = "Gamma"
×
2236

2237
        # use k_product to calculate kpoints, k_product = kpts[0][0] * a
2238
        lattice_abc = self.structure.lattice.abc
×
2239
        kpt_calc = [
×
2240
            int(self.k_product / lattice_abc[0] + 0.5),
2241
            int(self.k_product / lattice_abc[1] + 0.5),
2242
            1,
2243
        ]
2244

2245
        self.kpt_calc = kpt_calc
×
2246
        # calculate kpts (c direction) for bulk. (for slab, set to 1)
2247
        if self.bulk:
×
2248
            kpt_calc[2] = int(self.k_product / lattice_abc[2] + 0.5)
×
2249

2250
        kpt.kpts[0] = kpt_calc
×
2251

2252
        return kpt
×
2253

2254
    def as_dict(self, verbosity=2):
1✔
2255
        """
2256
        :param verbosity: Verbosity of dict. E.g., whether to include Structure.
2257
        :return: MSONAble dict
2258
        """
2259
        d = MSONable.as_dict(self)
×
2260
        if verbosity == 1:
×
2261
            d.pop("structure", None)
×
2262
        return d
×
2263

2264

2265
class MVLGBSet(MPRelaxSet):
1✔
2266
    """
2267
    Class for writing a vasp input files for grain boundary calculations, slab
2268
    or bulk.
2269
    """
2270

2271
    def __init__(self, structure: Structure, k_product=40, slab_mode=False, is_metal=True, **kwargs):
1✔
2272
        """
2273
        Args:
2274
            structure(Structure): provide the structure
2275
            k_product: Kpoint number * length for a & b directions, also for c
2276
                direction in bulk calculations. Default to 40.
2277
            slab_mode (bool): Defaults to False. Use default (False) for a
2278
                bulk supercell. Use True if you are performing calculations on a
2279
                slab-like (i.e., surface) of the GB, for example, when you are
2280
                calculating the work of separation.
2281
            is_metal (bool): Defaults to True. This determines whether an ISMEAR of
2282
                1 is used (for metals) or not (for insulators and semiconductors)
2283
                by default. Note that it does *not* override user_incar_settings,
2284
                which can be set by the user to be anything desired.
2285
            **kwargs:
2286
                Other kwargs supported by :class:`MPRelaxSet`.
2287
        """
2288
        super().__init__(structure, **kwargs)
×
2289
        self.k_product = k_product
×
2290
        self.slab_mode = slab_mode
×
2291
        self.is_metal = is_metal
×
2292

2293
    @property
1✔
2294
    def kpoints(self):
1✔
2295
        """
2296
        k_product, default to 40, is kpoint number * length for a & b
2297
        directions, also for c direction in bulk calculations
2298
        Automatic mesh & Gamma is the default setting.
2299
        """
2300
        # To get input sets, the input structure has to has the same number
2301
        # of required parameters as a Structure object.
2302

2303
        kpt = super().kpoints
×
2304
        kpt.comment = "Generated by pymatgen's MVLGBSet"
×
2305
        kpt.style = "Gamma"
×
2306

2307
        # use k_product to calculate kpoints, k_product = kpts[0][0] * a
2308
        lengths = self.structure.lattice.abc
×
2309
        kpt_calc = [
×
2310
            int(self.k_product / lengths[0] + 0.5),
2311
            int(self.k_product / lengths[1] + 0.5),
2312
            int(self.k_product / lengths[2] + 0.5),
2313
        ]
2314

2315
        if self.slab_mode:
×
2316
            kpt_calc[2] = 1
×
2317

2318
        kpt.kpts[0] = kpt_calc
×
2319

2320
        return kpt
×
2321

2322
    @property
1✔
2323
    def incar(self):
1✔
2324
        """
2325
        :return: Incar
2326
        """
2327
        incar = super().incar
×
2328

2329
        # The default incar setting is used for metallic system, for
2330
        # insulator or semiconductor, ISMEAR need to be changed.
2331
        incar.update(
×
2332
            {
2333
                "LCHARG": False,
2334
                "NELM": 60,
2335
                "PREC": "Normal",
2336
                "EDIFFG": -0.02,
2337
                "ICHARG": 0,
2338
                "NSW": 200,
2339
                "EDIFF": 0.0001,
2340
            }
2341
        )
2342

2343
        if self.is_metal:
×
2344
            incar["ISMEAR"] = 1
×
2345
            incar["LDAU"] = False
×
2346

2347
        if self.slab_mode:
×
2348
            # for clean grain boundary and bulk relaxation, full optimization
2349
            # relaxation (ISIF=3) is used. For slab relaxation (ISIF=2) is used.
2350
            incar["ISIF"] = 2
×
2351
            incar["NELMIN"] = 8
×
2352

2353
        incar.update(self.user_incar_settings)
×
2354

2355
        return incar
×
2356

2357

2358
class MVLRelax52Set(DictSet):
1✔
2359
    """
2360
    Implementation of VaspInputSet utilizing the public Materials Project
2361
    parameters for INCAR & KPOINTS and VASP's recommended PAW potentials for
2362
    POTCAR.
2363

2364
    Keynotes from VASP manual:
2365
        1. Recommended potentials for calculations using vasp.5.2+
2366
        2. If dimers with short bonds are present in the compound (O2, CO,
2367
            N2, F2, P2, S2, Cl2), it is recommended to use the h potentials.
2368
            Specifically, C_h, O_h, N_h, F_h, P_h, S_h, Cl_h
2369
        3. Released on Oct 28, 2018 by VASP. Please refer to VASP
2370
            Manual 1.2, 1.3 & 10.2.1 for more details.
2371
    """
2372

2373
    CONFIG = _load_yaml_config("MVLRelax52Set")
1✔
2374

2375
    def __init__(self, structure: Structure, **kwargs):
1✔
2376
        """
2377
        Args:
2378
            structure (Structure): input structure.
2379
            potcar_functional (str): choose from "PBE_52" and "PBE_54".
2380
            **kwargs: Other kwargs supported by :class:`DictSet`.
2381
        """
2382
        if kwargs.get("potcar_functional") or kwargs.get("user_potcar_functional"):
×
2383
            super().__init__(structure, MVLRelax52Set.CONFIG, **kwargs)
×
2384
        else:
2385
            super().__init__(structure, MVLRelax52Set.CONFIG, user_potcar_functional="PBE_52", **kwargs)
×
2386
        if self.potcar_functional not in ["PBE_52", "PBE_54"]:
×
2387
            raise ValueError("Please select from PBE_52 and PBE_54!")
×
2388

2389
        self.kwargs = kwargs
×
2390

2391

2392
class MITNEBSet(MITRelaxSet):
1✔
2393
    """
2394
    Class for writing NEB inputs. Note that EDIFF is not on a per atom
2395
    basis for this input set.
2396
    """
2397

2398
    def __init__(self, structures, unset_encut=False, **kwargs):
1✔
2399
        """
2400
        Args:
2401
            structures: List of Structure objects.
2402
            unset_encut (bool): Whether to unset ENCUT.
2403
            **kwargs: Other kwargs supported by :class:`DictSet`.
2404
        """
2405
        if len(structures) < 3:
×
2406
            raise ValueError("You need at least 3 structures for an NEB.")
×
2407
        kwargs["sort_structure"] = False
×
2408
        super().__init__(structures[0], **kwargs)
×
2409
        self.structures = self._process_structures(structures)
×
2410

2411
        self.unset_encut = False
×
2412
        if unset_encut:
×
2413
            self._config_dict["INCAR"].pop("ENCUT", None)
×
2414

2415
        if "EDIFF" not in self._config_dict["INCAR"]:
×
2416
            self._config_dict["INCAR"]["EDIFF"] = self._config_dict["INCAR"].pop("EDIFF_PER_ATOM")
×
2417

2418
        # NEB specific defaults
2419
        defaults = {
×
2420
            "IMAGES": len(structures) - 2,
2421
            "IBRION": 1,
2422
            "ISYM": 0,
2423
            "LCHARG": False,
2424
            "LDAU": False,
2425
        }
2426
        self._config_dict["INCAR"].update(defaults)
×
2427

2428
    @property
1✔
2429
    def poscar(self):
1✔
2430
        """
2431
        :return: Poscar for structure of first end point.
2432
        """
2433
        return Poscar(self.structures[0])
×
2434

2435
    @property
1✔
2436
    def poscars(self):
1✔
2437
        """
2438
        :return: List of Poscars.
2439
        """
2440
        return [Poscar(s) for s in self.structures]
×
2441

2442
    @staticmethod
1✔
2443
    def _process_structures(structures):
1✔
2444
        """
2445
        Remove any atom jumps across the cell
2446
        """
2447
        input_structures = structures
×
2448
        structures = [input_structures[0]]
×
2449
        for s in input_structures[1:]:
×
2450
            prev = structures[-1]
×
2451
            for i, site in enumerate(s):
×
2452
                t = np.round(prev[i].frac_coords - site.frac_coords)
×
2453
                if np.any(np.abs(t) > 0.5):
×
2454
                    s.translate_sites([i], t, to_unit_cell=False)
×
2455
            structures.append(s)
×
2456
        return structures
×
2457

2458
    def write_input(
1✔
2459
        self,
2460
        output_dir,
2461
        make_dir_if_not_present=True,
2462
        write_cif=False,
2463
        write_path_cif=False,
2464
        write_endpoint_inputs=False,
2465
    ):
2466
        """
2467
        NEB inputs has a special directory structure where inputs are in 00,
2468
        01, 02, ....
2469

2470
        Args:
2471
            output_dir (str): Directory to output the VASP input files
2472
            make_dir_if_not_present (bool): Set to True if you want the
2473
                directory (and the whole path) to be created if it is not
2474
                present.
2475
            write_cif (bool): If true, writes a cif along with each POSCAR.
2476
            write_path_cif (bool): If true, writes a cif for each image.
2477
            write_endpoint_inputs (bool): If true, writes input files for
2478
                running endpoint calculations.
2479
        """
2480
        output_dir = Path(output_dir)
×
2481
        if make_dir_if_not_present and not output_dir.exists():
×
2482
            output_dir.mkdir(parents=True)
×
2483
        self.incar.write_file(str(output_dir / "INCAR"))
×
2484
        self.kpoints.write_file(str(output_dir / "KPOINTS"))
×
2485
        self.potcar.write_file(str(output_dir / "POTCAR"))
×
2486

2487
        for i, p in enumerate(self.poscars):
×
2488
            d = output_dir / str(i).zfill(2)
×
2489
            if not d.exists():
×
2490
                d.mkdir(parents=True)
×
2491
            p.write_file(str(d / "POSCAR"))
×
2492
            if write_cif:
×
2493
                p.structure.to(filename=str(d / f"{i}.cif"))
×
2494
        if write_endpoint_inputs:
×
2495
            end_point_param = MITRelaxSet(self.structures[0], user_incar_settings=self.user_incar_settings)
×
2496

2497
            for image in ["00", str(len(self.structures) - 1).zfill(2)]:
×
2498
                end_point_param.incar.write_file(str(output_dir / image / "INCAR"))
×
2499
                end_point_param.kpoints.write_file(str(output_dir / image / "KPOINTS"))
×
2500
                end_point_param.potcar.write_file(str(output_dir / image / "POTCAR"))
×
2501
        if write_path_cif:
×
2502
            sites = set()
×
2503
            lat = self.structures[0].lattice
×
2504
            for site in chain(*(s.sites for s in self.structures)):
×
2505
                sites.add(PeriodicSite(site.species, site.frac_coords, lat))
×
2506
            nebpath = Structure.from_sites(sorted(sites))
×
2507
            nebpath.to(filename=str(output_dir / "path.cif"))
×
2508

2509

2510
class MITMDSet(MITRelaxSet):
1✔
2511
    """
2512
    Class for writing a vasp md run. This DOES NOT do multiple stage
2513
    runs.
2514
    """
2515

2516
    def __init__(self, structure: Structure, start_temp, end_temp, nsteps, time_step=2, spin_polarized=False, **kwargs):
1✔
2517
        """
2518
        Args:
2519
            structure (Structure): Input structure.
2520
            start_temp (int): Starting temperature.
2521
            end_temp (int): Final temperature.
2522
            nsteps (int): Number of time steps for simulations. NSW parameter.
2523
            time_step (int): The time step for the simulation. The POTIM
2524
                parameter. Defaults to 2fs.
2525
            spin_polarized (bool): Whether to do spin polarized calculations.
2526
                The ISPIN parameter. Defaults to False.
2527
            **kwargs: Other kwargs supported by :class:`DictSet`.
2528
        """
2529
        # MD default settings
2530
        defaults = {
×
2531
            "TEBEG": start_temp,
2532
            "TEEND": end_temp,
2533
            "NSW": nsteps,
2534
            "EDIFF_PER_ATOM": 0.000001,
2535
            "LSCALU": False,
2536
            "LCHARG": False,
2537
            "LPLANE": False,
2538
            "LWAVE": True,
2539
            "ISMEAR": 0,
2540
            "NELMIN": 4,
2541
            "LREAL": True,
2542
            "BMIX": 1,
2543
            "MAXMIX": 20,
2544
            "NELM": 500,
2545
            "NSIM": 4,
2546
            "ISYM": 0,
2547
            "ISIF": 0,
2548
            "IBRION": 0,
2549
            "NBLOCK": 1,
2550
            "KBLOCK": 100,
2551
            "SMASS": 0,
2552
            "POTIM": time_step,
2553
            "PREC": "Low",
2554
            "ISPIN": 2 if spin_polarized else 1,
2555
            "LDAU": False,
2556
        }
2557

2558
        super().__init__(structure, **kwargs)
×
2559

2560
        self.start_temp = start_temp
×
2561
        self.end_temp = end_temp
×
2562
        self.nsteps = nsteps
×
2563
        self.time_step = time_step
×
2564
        self.spin_polarized = spin_polarized
×
2565
        self.kwargs = kwargs
×
2566

2567
        # use VASP default ENCUT
2568
        self._config_dict["INCAR"].pop("ENCUT", None)
×
2569

2570
        if defaults["ISPIN"] == 1:
×
2571
            self._config_dict["INCAR"].pop("MAGMOM", None)
×
2572
        self._config_dict["INCAR"].update(defaults)
×
2573

2574
    @property
1✔
2575
    def kpoints(self):
1✔
2576
        """
2577
        :return: Kpoints
2578
        """
2579
        return Kpoints.gamma_automatic()
×
2580

2581

2582
class MPMDSet(MPRelaxSet):
1✔
2583
    """
2584
    This a modified version of the old MITMDSet pre 2018/03/12.
2585

2586
    This set serves as the basis for the amorphous skyline paper.
2587

2588
    (1) Aykol, M.; Dwaraknath, S. S.; Sun, W.; Persson, K. A. Thermodynamic
2589
        Limit for Synthesis of Metastable Inorganic Materials. Sci. Adv. 2018,
2590
        4 (4).
2591

2592
    Class for writing a vasp md run. This DOES NOT do multiple stage runs.
2593
    Precision remains normal, to increase accuracy of stress tensor.
2594
    """
2595

2596
    def __init__(self, structure: Structure, start_temp, end_temp, nsteps, spin_polarized=False, **kwargs):
1✔
2597
        """
2598
        Args:
2599
            structure (Structure): Input structure.
2600
            start_temp (int): Starting temperature.
2601
            end_temp (int): Final temperature.
2602
            nsteps (int): Number of time steps for simulations. NSW parameter.
2603
            time_step (int): The time step for the simulation. The POTIM
2604
                parameter. Defaults to 2fs.
2605
            spin_polarized (bool): Whether to do spin polarized calculations.
2606
                The ISPIN parameter. Defaults to False.
2607
            **kwargs: Other kwargs supported by :class:`DictSet`.
2608
        """
2609
        # MD default settings
2610
        defaults = {
×
2611
            "TEBEG": start_temp,
2612
            "TEEND": end_temp,
2613
            "NSW": nsteps,
2614
            "EDIFF_PER_ATOM": 0.00001,
2615
            "LSCALU": False,
2616
            "LCHARG": False,
2617
            "LPLANE": False,
2618
            "LWAVE": True,
2619
            "ISMEAR": 0,
2620
            "NELMIN": 4,
2621
            "LREAL": True,
2622
            "BMIX": 1,
2623
            "MAXMIX": 20,
2624
            "NELM": 500,
2625
            "NSIM": 4,
2626
            "ISYM": 0,
2627
            "ISIF": 0,
2628
            "IBRION": 0,
2629
            "NBLOCK": 1,
2630
            "KBLOCK": 100,
2631
            "SMASS": 0,
2632
            "POTIM": 2,
2633
            "PREC": "Normal",
2634
            "ISPIN": 2 if spin_polarized else 1,
2635
            "LDAU": False,
2636
            "ADDGRID": True,
2637
        }
2638

2639
        if Element("H") in structure.species:
×
2640
            defaults["POTIM"] = 0.5
×
2641
            defaults["NSW"] = defaults["NSW"] * 4
×
2642

2643
        super().__init__(structure, **kwargs)
×
2644

2645
        self.start_temp = start_temp
×
2646
        self.end_temp = end_temp
×
2647
        self.nsteps = nsteps
×
2648
        self.spin_polarized = spin_polarized
×
2649
        self.kwargs = kwargs
×
2650

2651
        # use VASP default ENCUT
2652
        self._config_dict["INCAR"].pop("ENCUT", None)
×
2653

2654
        if defaults["ISPIN"] == 1:
×
2655
            self._config_dict["INCAR"].pop("MAGMOM", None)
×
2656
        self._config_dict["INCAR"].update(defaults)
×
2657

2658
    @property
1✔
2659
    def kpoints(self):
1✔
2660
        """
2661
        :return: Kpoints
2662
        """
2663
        return Kpoints.gamma_automatic()
×
2664

2665

2666
class MVLNPTMDSet(MITMDSet):
1✔
2667
    """
2668
    Class for writing a vasp md run in NPT ensemble.
2669

2670
    Notes:
2671
        To eliminate Pulay stress, the default ENCUT is set to a rather large
2672
        value of ENCUT, which is 1.5 * ENMAX.
2673
    """
2674

2675
    def __init__(self, structure: Structure, start_temp, end_temp, nsteps, time_step=2, spin_polarized=False, **kwargs):
1✔
2676
        """
2677
        Args:
2678
            structure (Structure): input structure.
2679
            start_temp (int): Starting temperature.
2680
            end_temp (int): Final temperature.
2681
            nsteps(int): Number of time steps for simulations. NSW parameter.
2682
            time_step (int): The time step for the simulation. The POTIM
2683
                parameter. Defaults to 2fs.
2684
            spin_polarized (bool): Whether to do spin polarized calculations.
2685
                The ISPIN parameter. Defaults to False.
2686
            **kwargs: Other kwargs supported by :class:`DictSet`.
2687
        """
2688
        user_incar_settings = kwargs.get("user_incar_settings", {})
×
2689

2690
        # NPT-AIMD default settings
2691
        defaults = {
×
2692
            "ALGO": "Fast",
2693
            "ISIF": 3,
2694
            "LANGEVIN_GAMMA": [10] * structure.ntypesp,
2695
            "LANGEVIN_GAMMA_L": 1,
2696
            "MDALGO": 3,
2697
            "PMASS": 10,
2698
            "PSTRESS": 0,
2699
            "SMASS": 0,
2700
        }
2701

2702
        defaults.update(user_incar_settings)
×
2703
        kwargs["user_incar_settings"] = defaults
×
2704

2705
        super().__init__(structure, start_temp, end_temp, nsteps, time_step, spin_polarized, **kwargs)
×
2706

2707
        # Set NPT-AIMD ENCUT = 1.5 * VASP_default
2708
        enmax = [self.potcar[i].keywords["ENMAX"] for i in range(structure.ntypesp)]
×
2709
        encut = max(enmax) * 1.5
×
2710
        self._config_dict["INCAR"]["ENCUT"] = encut
×
2711

2712

2713
class MVLScanRelaxSet(MPRelaxSet):
1✔
2714
    """
2715
    Class for writing a relax input set using Strongly Constrained and
2716
    Appropriately Normed (SCAN) semilocal density functional.
2717

2718
    Notes:
2719
        1. This functional is only available from VASP.5.4.3 upwards.
2720

2721
        2. Meta-GGA calculations require POTCAR files that include
2722
        information on the kinetic energy density of the core-electrons,
2723
        i.e. "PBE_52" or "PBE_54". Make sure the POTCAR including the
2724
        following lines (see VASP wiki for more details):
2725

2726
            $ grep kinetic POTCAR
2727
            kinetic energy-density
2728
            mkinetic energy-density pseudized
2729
            kinetic energy density (partial)
2730
    """
2731

2732
    def __init__(self, structure: Structure, **kwargs):
1✔
2733
        """
2734
        Args:
2735
            structure (Structure): input structure.
2736
            vdw (str): set "rVV10" to enable SCAN+rVV10, which is a versatile
2737
                van der Waals density functional by combing the SCAN functional
2738
                with the rVV10 non-local correlation functional.
2739
            **kwargs: Other kwargs supported by :class:`DictSet`.
2740
        """
2741
        # choose PBE_52 unless the user specifies something else
2742
        if kwargs.get("potcar_functional") or kwargs.get("user_potcar_functional"):
×
2743
            super().__init__(structure, **kwargs)
×
2744
        else:
2745
            super().__init__(structure, user_potcar_functional="PBE_52", **kwargs)
×
2746

2747
        if self.potcar_functional not in ["PBE_52", "PBE_54"]:
×
2748
            raise ValueError("SCAN calculations required PBE_52 or PBE_54!")
×
2749

2750
        updates = {
×
2751
            "ADDGRID": True,
2752
            "EDIFF": 1e-05,
2753
            "EDIFFG": -0.05,
2754
            "LASPH": True,
2755
            "LDAU": False,
2756
            "METAGGA": "SCAN",
2757
            "NELM": 200,
2758
        }
2759

2760
        if kwargs.get("vdw", "").lower() == "rvv10":
×
2761
            updates["BPARAM"] = 15.7  # This is the correct BPARAM for SCAN+rVV10
×
2762

2763
        self._config_dict["INCAR"].update(updates)
×
2764

2765

2766
class LobsterSet(MPRelaxSet):
1✔
2767
    """
2768
    Input set to prepare VASP runs that can be digested by Lobster (See cohp.de)
2769
    """
2770

2771
    CONFIG = _load_yaml_config("MPRelaxSet")
1✔
2772

2773
    def __init__(
1✔
2774
        self,
2775
        structure: Structure,
2776
        isym: int = 0,
2777
        ismear: int = -5,
2778
        reciprocal_density: int | None = None,
2779
        address_basis_file: str | None = None,
2780
        user_supplied_basis: dict | None = None,
2781
        user_potcar_settings: dict | None = None,
2782
        **kwargs,
2783
    ):
2784
        """
2785
        Args:
2786
            structure (Structure): input structure.
2787
            isym (int): ISYM entry for INCAR, only isym=-1 and isym=0 are allowed
2788
            ismear (int): ISMEAR entry for INCAR, only ismear=-5 and ismear=0 are allowed
2789
            reciprocal_density (int): density of k-mesh by reciprocal volume
2790
            user_supplied_basis (dict): dict including basis functions for all elements in structure,
2791
                e.g. {"Fe": "3d 3p 4s", "O": "2s 2p"}; if not supplied, a standard basis is used
2792
            address_basis_file (str): address to a file similar to "BASIS_PBE_54_standaard.yaml"
2793
                in pymatgen.io.lobster.lobster_basis
2794
            **kwargs: Other kwargs supported by :class:`DictSet`.
2795
        """
2796
        from pymatgen.io.lobster import Lobsterin
×
2797

2798
        warnings.warn("Make sure that all parameters are okay! This is a brand new implementation.")
×
2799

2800
        if isym not in (-1, 0):
×
2801
            raise ValueError("Lobster cannot digest WAVEFUNCTIONS with symmetry")
×
2802
        if ismear not in (-5, 0):
×
2803
            raise ValueError("Lobster usually works with ismear=-5 or ismear=0")
×
2804

2805
        # newest potcars are preferred
2806
        # Choose PBE_54 unless the user specifies a different potcar_functional
2807
        if kwargs.get("potcar_functional") or kwargs.get("user_potcar_functional"):
×
2808
            super().__init__(structure, **kwargs)
×
2809
        else:
2810
            super().__init__(structure, user_potcar_functional="PBE_54", user_potcar_settings={"W": "W_sv"}, **kwargs)
×
2811

2812
        # reciprocal density
2813
        if self.user_kpoints_settings is not None:
×
2814
            if not reciprocal_density or "reciprocal_density" not in self.user_kpoints_settings:
×
2815
                # test, if this is okay
2816
                self.reciprocal_density = 310
×
2817
            else:
2818
                self.reciprocal_density = reciprocal_density or self.user_kpoints_settings["reciprocal_density"]
×
2819
        else:
2820
            if not reciprocal_density:
×
2821
                # test, if this is okay
2822
                self.reciprocal_density = 310
×
2823
            else:
2824
                self.reciprocal_density = reciprocal_density
×
2825

2826
        self._config_dict["POTCAR"].update({"W": "W_sv"})
×
2827
        self.isym = isym
×
2828
        self.ismear = ismear
×
2829
        self.user_supplied_basis = user_supplied_basis
×
2830
        self.address_basis_file = address_basis_file
×
2831
        # predefined basis! Check if the basis is okay! (charge spilling and bandoverlaps!)
2832
        if user_supplied_basis is None and address_basis_file is None:
×
2833
            basis = Lobsterin.get_basis(structure=structure, potcar_symbols=self.potcar_symbols)
×
2834
        elif address_basis_file is not None:
×
2835
            basis = Lobsterin.get_basis(
×
2836
                structure=structure,
2837
                potcar_symbols=self.potcar_symbols,
2838
                address_basis_file=address_basis_file,
2839
            )
2840
        elif user_supplied_basis is not None:
×
2841
            # test if all elements from structure are in user_supplied_basis
2842
            for atomtype in structure.symbol_set:
×
2843
                if atomtype not in user_supplied_basis:
×
2844
                    raise ValueError("There are no basis functions for the atom type " + str(atomtype))
×
2845
            basis = [f"{key} {value}" for key, value in user_supplied_basis.items()]
×
2846

2847
        lobsterin = Lobsterin(settingsdict={"basisfunctions": basis})
×
2848
        nbands = lobsterin._get_nbands(structure=structure)
×
2849

2850
        update_dict = {
×
2851
            "EDIFF": 1e-6,
2852
            "NSW": 0,
2853
            "LWAVE": True,
2854
            "ISYM": isym,
2855
            "NBANDS": nbands,
2856
            "IBRION": -1,
2857
            "ISMEAR": ismear,
2858
            "LORBIT": 11,
2859
            "ICHARG": 0,
2860
            "ALGO": "Normal",
2861
        }
2862

2863
        self._config_dict["INCAR"].update(update_dict)
×
2864
        self._config_dict["KPOINTS"].update({"reciprocal_density": self.reciprocal_density})
×
2865

2866

2867
def get_vasprun_outcar(path, parse_dos=True, parse_eigen=True):
1✔
2868
    """
2869
    :param path: Path to get the vasprun.xml and OUTCAR.
2870
    :param parse_dos: Whether to parse dos. Defaults to True.
2871
    :param parse_eigen: Whether to parse eigenvalue. Defaults to True.
2872
    :return:
2873
    """
2874
    path = Path(path)
×
2875
    vruns = list(glob.glob(str(path / "vasprun.xml*")))
×
2876
    outcars = list(glob.glob(str(path / "OUTCAR*")))
×
2877

2878
    if len(vruns) == 0 or len(outcars) == 0:
×
2879
        raise ValueError(f"Unable to get vasprun.xml/OUTCAR from prev calculation in {path}")
×
2880
    vsfile_fullpath = str(path / "vasprun.xml")
×
2881
    outcarfile_fullpath = str(path / "OUTCAR")
×
2882
    vsfile = vsfile_fullpath if vsfile_fullpath in vruns else sorted(vruns)[-1]
×
2883
    outcarfile = outcarfile_fullpath if outcarfile_fullpath in outcars else sorted(outcars)[-1]
×
2884
    return (
×
2885
        Vasprun(vsfile, parse_dos=parse_dos, parse_eigen=parse_eigen),
2886
        Outcar(outcarfile),
2887
    )
2888

2889

2890
def get_structure_from_prev_run(vasprun, outcar=None):
1✔
2891
    """
2892
    Process structure from previous run.
2893

2894
    Args:
2895
        vasprun (Vasprun): Vasprun that contains the final structure
2896
            from previous run.
2897
        outcar (Outcar): Outcar that contains the magnetization info from
2898
            previous run.
2899

2900
    Returns:
2901
        Returns the magmom-decorated structure that can be passed to get
2902
        VASP input files, e.g. get_kpoints.
2903
    """
2904
    structure = vasprun.final_structure
×
2905

2906
    site_properties = {}
×
2907
    # magmom
2908
    if vasprun.is_spin:
×
2909
        if outcar and outcar.magnetization:
×
2910
            site_properties.update({"magmom": [i["tot"] for i in outcar.magnetization]})
×
2911
        else:
2912
            site_properties.update({"magmom": vasprun.parameters["MAGMOM"]})
×
2913
    # ldau
2914
    if vasprun.parameters.get("LDAU", False):
×
2915
        for k in ("LDAUU", "LDAUJ", "LDAUL"):
×
2916
            vals = vasprun.incar[k]
×
2917
            m = {}
×
2918
            l_val = []
×
2919
            s = 0
×
2920
            for site in structure:
×
2921
                if site.specie.symbol not in m:
×
2922
                    m[site.specie.symbol] = vals[s]
×
2923
                    s += 1
×
2924
                l_val.append(m[site.specie.symbol])
×
2925
            if len(l_val) == len(structure):
×
2926
                site_properties.update({k.lower(): l_val})
×
2927
            else:
2928
                raise ValueError(f"length of list {l_val} not the same as structure")
×
2929

2930
    return structure.copy(site_properties=site_properties)
×
2931

2932

2933
def standardize_structure(structure, sym_prec=0.1, international_monoclinic=True):
1✔
2934
    """
2935
    Get the symmetrically standardized structure.
2936

2937
    Args:
2938
        structure (Structure): The structure.
2939
        sym_prec (float): Tolerance for symmetry finding for standardization.
2940
        international_monoclinic (bool): Whether to use international
2941
            convention (vs Curtarolo) for monoclinic. Defaults True.
2942

2943
    Returns:
2944
        The symmetrized structure.
2945
    """
2946
    sym_finder = SpacegroupAnalyzer(structure, symprec=sym_prec)
×
2947
    new_structure = sym_finder.get_primitive_standard_structure(international_monoclinic=international_monoclinic)
×
2948

2949
    # the primitive structure finding has had several bugs in the past
2950
    # defend through validation
2951
    vpa_old = structure.volume / structure.num_sites
×
2952
    vpa_new = new_structure.volume / new_structure.num_sites
×
2953

2954
    if abs(vpa_old - vpa_new) / vpa_old > 0.02:
×
2955
        raise ValueError(f"Standardizing cell failed! VPA old: {vpa_old}, VPA new: {vpa_new}")
×
2956

2957
    sm = StructureMatcher()
×
2958
    if not sm.fit(structure, new_structure):
×
2959
        raise ValueError("Standardizing cell failed! Old structure doesn't match new.")
×
2960

2961
    return new_structure
×
2962

2963

2964
class BadInputSetWarning(UserWarning):
1✔
2965
    """
2966
    Warning class for bad but legal inputs.
2967
    """
2968

2969

2970
def batch_write_input(
1✔
2971
    structures,
2972
    vasp_input_set=MPRelaxSet,
2973
    output_dir=".",
2974
    make_dir_if_not_present=True,
2975
    subfolder=None,
2976
    sanitize=False,
2977
    include_cif=False,
2978
    potcar_spec=False,
2979
    zip_output=False,
2980
    **kwargs,
2981
):
2982
    """
2983
    Batch write vasp input for a sequence of structures to
2984
    output_dir, following the format output_dir/{group}/{formula}_{number}.
2985

2986
    Args:
2987
        structures ([Structure]): Sequence of Structures.
2988
        vasp_input_set (VaspInputSet): VaspInputSet class that creates
2989
            vasp input files from structures. Note that a class should be
2990
            supplied. Defaults to MPRelaxSet.
2991
        output_dir (str): Directory to output files. Defaults to current
2992
            directory ".".
2993
        make_dir_if_not_present (bool): Create the directory if not present.
2994
            Defaults to True.
2995
        subfolder (callable): Function to create subdirectory name from
2996
            structure. Defaults to simply "formula_count".
2997
        sanitize (bool): Boolean indicating whether to sanitize the
2998
            structure before writing the VASP input files. Sanitized output
2999
            are generally easier for viewing and certain forms of analysis.
3000
            Defaults to False.
3001
        include_cif (bool): Whether to output a CIF as well. CIF files are
3002
            generally better supported in visualization programs.
3003
        potcar_spec (bool): Instead of writing the POTCAR, write a "POTCAR.spec".
3004
                This is intended to help sharing an input set with people who might
3005
                not have a license to specific Potcar files. Given a "POTCAR.spec",
3006
                the specific POTCAR file can be re-generated using pymatgen with the
3007
                "generate_potcar" function in the pymatgen CLI.
3008
        zip_output (bool): If True, output will be zipped into a file with the
3009
            same name as the InputSet (e.g., MPStaticSet.zip)
3010
        **kwargs: Additional kwargs are passed to the vasp_input_set class
3011
            in addition to structure.
3012
    """
3013
    output_dir = Path(output_dir)
×
3014
    for i, s in enumerate(structures):
×
3015
        formula = re.sub(r"\s+", "", s.formula)
×
3016
        if subfolder is not None:
×
3017
            subdir = subfolder(s)
×
3018
            d = output_dir / subdir
×
3019
        else:
3020
            d = output_dir / f"{formula}_{i}"
×
3021
        if sanitize:
×
3022
            s = s.copy(sanitize=True)
×
3023
        v = vasp_input_set(s, **kwargs)
×
3024
        v.write_input(
×
3025
            str(d),
3026
            make_dir_if_not_present=make_dir_if_not_present,
3027
            include_cif=include_cif,
3028
            potcar_spec=potcar_spec,
3029
            zip_output=zip_output,
3030
        )
3031

3032

3033
_dummy_structure = Structure(
1✔
3034
    [1, 0, 0, 0, 1, 0, 0, 0, 1],
3035
    ["I"],
3036
    [[0, 0, 0]],
3037
    site_properties={"magmom": [[0, 0, 1]]},
3038
)
3039

3040

3041
def get_valid_magmom_struct(structure, inplace=True, spin_mode="auto"):
1✔
3042
    """
3043
    Make sure that the structure has valid magmoms based on the kind of calculation
3044
    Fill in missing Magmom values
3045

3046
    Args:
3047
        structure: The input structure
3048
        inplace: True - edit the magmom of the input structurel; False - return new structure
3049
        spin_mode: "scalar"/"vector"/"none"/"auto" only first letter (s/v/n) is needed.
3050
            dictates how the spin configuration will be determined.
3051

3052
            - auto: read the existing magmom values and decide
3053
            - scalar: use a single scalar value (for spin up/down)
3054
            - vector: use a vector value for spin-orbit systems
3055
            - none: Remove all the magmom information
3056

3057
    Returns:
3058
        New structure if inplace is False
3059
    """
3060
    default_values = {"s": 1.0, "v": [1.0, 1.0, 1.0], "n": None}
1✔
3061
    if spin_mode[0].lower() == "a":
1✔
3062
        mode = "n"
1✔
3063
        for isite in structure.sites:
1✔
3064
            if "magmom" not in isite.properties or isite.properties["magmom"] is None:
1✔
3065
                pass
×
3066
            elif isinstance(isite.properties["magmom"], (float, int)):
×
3067
                if mode == "v":
×
3068
                    raise TypeError("Magmom type conflict")
×
3069
                mode = "s"
×
3070
                if isinstance(isite.properties["magmom"], int):
×
3071
                    isite.properties["magmom"] = float(isite.properties["magmom"])
×
3072
            elif len(isite.properties["magmom"]) == 3:
×
3073
                if mode == "s":
×
3074
                    raise TypeError("Magmom type conflict")
×
3075
                mode = "v"
×
3076
            else:
3077
                raise TypeError("Unrecognized Magmom Value")
×
3078
    else:
3079
        mode = spin_mode[0].lower()
×
3080

3081
    if not inplace:
1✔
3082
        new_struct = structure.copy()
×
3083
    else:
3084
        new_struct = structure
1✔
3085
    for isite in new_struct.sites:
1✔
3086
        if mode == "n":
1✔
3087
            if "magmom" in isite.properties:
1✔
3088
                isite.properties.pop("magmom")
×
3089
        elif "magmom" not in isite.properties or isite.properties["magmom"] is None:
×
3090
            isite.properties["magmom"] = default_values[mode]
×
3091

3092
    if not inplace:
1✔
3093
        return new_struct
×
3094
    return None
1✔
3095

3096

3097
class MPAbsorptionSet(MPRelaxSet):
1✔
3098
    """
3099
    MP input set for generating frequency dependent dielectrics.
3100
    Two modes are supported: "IPA" or "RPA".
3101
    A typical sequence is mode="STATIC" -> mode="IPA" -> mode="RPA"(optional)
3102
    For all steps other than the first one (static), the
3103
    recommendation is to use from_prev_calculation on the preceding run in
3104
    the series. It is important to ensure Gamma centred kpoints for the RPA step.
3105

3106
    """
3107

3108
    # CONFIG = _load_yaml_config("MPAbsorptionSet")
3109

3110
    SUPPORTED_MODES = ("IPA", "RPA")
1✔
3111

3112
    def __init__(
1✔
3113
        self,
3114
        structure,
3115
        mode="IPA",
3116
        copy_wavecar=True,
3117
        nbands=None,
3118
        nbands_factor=2,
3119
        reciprocal_density=400,
3120
        nkred=None,
3121
        nedos=2001,
3122
        prev_incar=None,
3123
        **kwargs,
3124
    ):
3125
        """
3126
        Args:
3127
            structure (Structure): Input structure.
3128
            prev_incar (Incar/string): Incar file from previous run.
3129
            mode (str): Supported modes are "IPA", "RPA"
3130
            nbands (int): For subsequent calculations, it is generally
3131
                recommended to perform NBANDS convergence starting from the
3132
                NBANDS of the previous run for DIAG, and to use the exact same
3133
                NBANDS for RPA. This parameter is used by
3134
                from_previous_calculation to set nband.
3135
            nbands_factor (int): Multiplicative factor for NBANDS when starting
3136
                from a previous calculation. Only applies if mode=="IPA".
3137
                Need to be tested for convergence.
3138
            reciprocal_density: the k-points density
3139
            nkred: the reduced number of kpoints to calculate, equal to the k-mesh. Only applies in "RPA" mode
3140
                  because of the q->0 limit.
3141
            nedos: the density of DOS, default: 2001.
3142
            **kwargs: All kwargs supported by DictSet. Typically, user_incar_settings is a commonly used option.
3143
        """
3144
        # Initialize the input set (default: IPA absorption)
3145
        super().__init__(structure, **kwargs)
×
3146

3147
        self.prev_incar = prev_incar
×
3148
        self.nbands = nbands
×
3149
        self.reciprocal_density = reciprocal_density
×
3150
        self.mode = mode.upper()
×
3151
        if self.mode not in MPAbsorptionSet.SUPPORTED_MODES:
×
3152
            raise ValueError(f"{self.mode} not one of the support modes : {MPAbsorptionSet.SUPPORTED_MODES}")
×
3153
        self.copy_wavecar = copy_wavecar
×
3154
        self.nbands_factor = nbands_factor
×
3155
        self.nedos = nedos
×
3156
        self.nkred = nkred
×
3157
        self.kwargs = kwargs
×
3158

3159
    @property
1✔
3160
    def kpoints(self):
1✔
3161
        """
3162
        Generate gamma center k-points mesh grid for optical calculation. It is not mandatory for 'ALGO = Exact',
3163
        but is requested by 'ALGO = CHI' calculation.
3164
        """
3165
        return Kpoints.automatic_density_by_vol(self.structure, self.reciprocal_density, force_gamma=True)
×
3166

3167
    @property
1✔
3168
    def incar(self):
1✔
3169
        """
3170
        :return: Incar
3171
        """
3172
        parent_incar = super().incar
×
3173
        absorption_incar = {
×
3174
            "ALGO": "Exact",
3175
            "EDIFF": 1.0e-8,
3176
            "IBRION": -1,
3177
            "ICHARG": 1,
3178
            "ISMEAR": 0,
3179
            "SIGMA": 0.01,
3180
            "LWAVE": True,
3181
            "LREAL": False,  # for small cell it's more efficient to use reciprocal
3182
            "NELM": 100,
3183
            "NSW": 0,
3184
            "LOPTICS": True,
3185
            "CSHIFT": 0.1,
3186
            "NEDOS": 2001,
3187
        }
3188
        self._config_dict["INCAR"].update(absorption_incar)
×
3189

3190
        if self.mode == "IPA":
×
3191
            # use the incar from previous static calculation
3192
            if self.prev_incar is not None:
×
3193
                incar = Incar(self.prev_incar)
×
3194
                # Default parameters for diagonalization calculation.
3195
                incar.update({"ALGO": "Exact", "LOPTICS": True, "CSHIFT": 0.1, "LWAVE": "True"})
×
3196
                incar.update({"NEDOS": self.nedos}) if self.nedos is not None else incar.update({"NEDOS": 2001})
×
3197
            else:
3198
                incar = Incar(parent_incar)
×
3199

3200
        elif self.mode == "RPA":
×
3201
            # Default parameters for the response function calculation. NELM has to be set to 1.
3202
            # NOMEGA is set to 1000 in order to get smooth spectrum
3203
            incar = Incar(self.prev_incar) if self.prev_incar is not None else Incar(parent_incar)
×
3204
            incar.update({"ALGO": "CHI", "NELM": 1, "NOMEGA": 1000})
×
3205

3206
            if self.nkred is not None:
×
3207
                incar["NKREDX"] = self.nkred[0]
×
3208
                incar["NKREDY"] = self.nkred[1]
×
3209
                incar["NKREDZ"] = self.nkred[2]
×
3210

3211
            incar.pop("EDIFF", None)
×
3212
            incar.pop("LOPTICS", None)
×
3213
            incar.pop("LWAVE", False)
×
3214

3215
        else:
3216
            raise Exception("mode has to be from 'IPA' or 'RPA'")
×
3217

3218
        if self.nbands:
×
3219
            incar["NBANDS"] = self.nbands
×
3220

3221
        # Respect user set INCAR.
3222
        incar.update(self.kwargs.get("user_incar_settings", {}))
×
3223

3224
        return incar
×
3225

3226
    def override_from_prev_calc(self, prev_calc_dir=".", **kwargs):
1✔
3227
        """
3228
        Update the input set to include settings from a previous calculation.
3229

3230
        Args:
3231
            prev_calc_dir (str): The path to the previous calculation directory.
3232

3233
        Returns:
3234
            The input set with the settings (structure, k-points, incar, etc)
3235
            updated using the previous VASP run.
3236
        """
3237
        vasprun, outcar = get_vasprun_outcar(prev_calc_dir)
×
3238
        self.prev_incar = vasprun.incar
×
3239
        self._structure = vasprun.final_structure
×
3240

3241
        # The number of bands is multiplied by the factor
3242
        prev_nbands = int(vasprun.parameters["NBANDS"])
×
3243

3244
        if self.mode.upper() == "IPA":
×
3245
            self.nbands = int(np.ceil(prev_nbands * self.nbands_factor))
×
3246

3247
        # Since in the optical calculation, only the q->0 transition is of interests, we can reduce the number of q by
3248
        # the factor of the number of kpoints in each corresponding x, y, z directions. This will reduce the
3249
        # computational work by factor of 1/nkredx*nkredy*nkredz. An isotropic NKRED can be used for cubic
3250
        # lattice, but using NKREDX, NKREDY, NKREDZ is more sensible for other lattice.
3251
        if self.mode.upper() == "RPA":
×
3252
            # self.nbands = int(np.ceil(self.nbands * self.nbands_factor / self.ncores) * self.ncores)
3253
            self.nkred = vasprun.kpoints.kpts[0]
×
3254

3255
        files_to_transfer = {}
×
3256
        if self.copy_wavecar:
×
3257
            for fname in ("WAVECAR", "WAVEDER"):
×
3258
                w = sorted(glob.glob(str(Path(prev_calc_dir) / (fname + "*"))))
×
3259
                if w:
×
3260
                    files_to_transfer[fname] = str(w[-1])
×
3261

3262
        self.files_to_transfer.update(files_to_transfer)
×
3263

3264
        return self
×
3265

3266
    @classmethod
1✔
3267
    def from_prev_calc(cls, prev_calc_dir, mode, **kwargs):
1✔
3268
        """
3269
        Generate a set of VASP input files for absorption calculation
3270
        Args:
3271
            prev_calc_dir (str): The directory contains the outputs(
3272
                vasprun.xml of previous vasp run.
3273
            mode (str): Supported modes are "IPA", "RPA" (default)
3274
            **kwargs: All kwargs supported by MPAbsorptionsSet, other than structure
3275
        """
3276
        input_set = cls(_dummy_structure, mode, **kwargs)
×
3277
        return input_set.override_from_prev_calc(prev_calc_dir=prev_calc_dir)
×
STATUS · Troubleshooting · Open an Issue · Sales · Support · CAREERS · ENTERPRISE · START FREE · SCHEDULE DEMO
ANNOUNCEMENTS · TWITTER · TOS & SLA · Supported CI Services · What's a CI service? · Automated Testing

© 2025 Coveralls, Inc