• 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

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

4
"""
1✔
5
Classes for reading/manipulating/writing VASP input files. All major VASP input
6
files.
7
"""
8

9
from __future__ import annotations
1✔
10

11
import glob
1✔
12
import itertools
1✔
13
import json
1✔
14
import logging
1✔
15
import math
1✔
16
import os
1✔
17
import re
1✔
18
import subprocess
1✔
19
import warnings
1✔
20
from collections import namedtuple
1✔
21
from enum import Enum
1✔
22
from hashlib import md5, sha256
1✔
23
from typing import Any, Literal, Sequence
1✔
24

25
import numpy as np
1✔
26
import scipy.constants as const
1✔
27
from monty.io import zopen
1✔
28
from monty.json import MontyDecoder, MSONable
1✔
29
from monty.os import cd
1✔
30
from monty.os.path import zpath
1✔
31
from monty.serialization import loadfn
1✔
32
from tabulate import tabulate
1✔
33

34
from pymatgen.core import SETTINGS
1✔
35
from pymatgen.core.lattice import Lattice
1✔
36
from pymatgen.core.periodic_table import Element, get_el_sp
1✔
37
from pymatgen.core.structure import Structure
1✔
38
from pymatgen.electronic_structure.core import Magmom
1✔
39
from pymatgen.util.io_utils import clean_lines
1✔
40
from pymatgen.util.string import str_delimited
1✔
41
from pymatgen.util.typing import ArrayLike, PathLike
1✔
42

43
__author__ = "Shyue Ping Ong, Geoffroy Hautier, Rickard Armiento, Vincent L Chevrier, Stephen Dacek"
1✔
44
__copyright__ = "Copyright 2011, The Materials Project"
1✔
45

46
logger = logging.getLogger(__name__)
1✔
47

48

49
class Poscar(MSONable):
1✔
50
    """
51
    Object for representing the data in a POSCAR or CONTCAR file.
52
    Please note that this current implementation. Most attributes can be set
53
    directly.
54

55
    .. attribute:: structure
56

57
        Associated Structure.
58

59
    .. attribute:: comment
60

61
        Optional comment string.
62

63
    .. attribute:: true_names
64

65
        Boolean indication whether Poscar contains actual real names parsed
66
        from either a POTCAR or the POSCAR itself.
67

68
    .. attribute:: selective_dynamics
69

70
        Selective dynamics attribute for each site if available. A Nx3 array of
71
        booleans.
72

73
    .. attribute:: velocities
74

75
        Velocities for each site (typically read in from a CONTCAR). A Nx3
76
        array of floats.
77

78
    .. attribute:: predictor_corrector
79

80
        Predictor corrector coordinates and derivatives for each site; i.e.
81
        a list of three 1x3 arrays for each site (typically read in from a MD
82
        CONTCAR).
83

84
    .. attribute:: predictor_corrector_preamble
85

86
        Predictor corrector preamble contains the predictor-corrector key,
87
        POTIM, and thermostat parameters that precede the site-specic predictor
88
        corrector data in MD CONTCAR
89

90
    .. attribute:: temperature
91

92
        Temperature of velocity Maxwell-Boltzmann initialization. Initialized
93
        to -1 (MB hasn"t been performed).
94
    """
95

96
    def __init__(
1✔
97
        self,
98
        structure: Structure,
99
        comment: str | None = None,
100
        selective_dynamics=None,
101
        true_names: bool = True,
102
        velocities: ArrayLike | None = None,
103
        predictor_corrector: ArrayLike | None = None,
104
        predictor_corrector_preamble: str | None = None,
105
        sort_structure: bool = False,
106
    ):
107
        """
108
        :param structure: Structure object.
109
        :param comment: Optional comment line for POSCAR. Defaults to unit
110
            cell formula of structure. Defaults to None.
111
        :param selective_dynamics: bool values for selective dynamics,
112
            where N is number of sites. Defaults to None.
113
        :param true_names: Set to False if the names in the POSCAR are not
114
            well-defined and ambiguous. This situation arises commonly in
115
            vasp < 5 where the POSCAR sometimes does not contain element
116
            symbols. Defaults to True.
117
        :param velocities: Velocities for the POSCAR. Typically parsed
118
            in MD runs or can be used to initialize velocities.
119
        :param predictor_corrector: Predictor corrector for the POSCAR.
120
            Typically parsed in MD runs.
121
        :param predictor_corrector_preamble: Preamble to the predictor
122
            corrector.
123
        :param sort_structure: Whether to sort structure. Useful if species
124
            are not grouped properly together.
125
        """
126
        if structure.is_ordered:
1✔
127
            site_properties = {}
1✔
128
            if selective_dynamics:
1✔
129
                site_properties["selective_dynamics"] = selective_dynamics
1✔
130
            if velocities:
1✔
131
                site_properties["velocities"] = velocities
1✔
132
            if predictor_corrector:
1✔
133
                site_properties["predictor_corrector"] = predictor_corrector
1✔
134
            structure = Structure.from_sites(structure)
1✔
135
            self.structure = structure.copy(site_properties=site_properties)
1✔
136
            if sort_structure:
1✔
137
                self.structure = self.structure.get_sorted_structure()
×
138
            self.true_names = true_names
1✔
139
            self.comment = structure.formula if comment is None else comment
1✔
140
            self.predictor_corrector_preamble = predictor_corrector_preamble
1✔
141
        else:
142
            raise ValueError("Structure with partial occupancies cannot be converted into POSCAR!")
×
143

144
        self.temperature = -1.0
1✔
145

146
    @property
1✔
147
    def velocities(self):
1✔
148
        """Velocities in Poscar"""
149
        return self.structure.site_properties.get("velocities")
1✔
150

151
    @property
1✔
152
    def selective_dynamics(self):
1✔
153
        """Selective dynamics in Poscar"""
154
        return self.structure.site_properties.get("selective_dynamics")
1✔
155

156
    @property
1✔
157
    def predictor_corrector(self):
1✔
158
        """Predictor corrector in Poscar"""
159
        return self.structure.site_properties.get("predictor_corrector")
1✔
160

161
    @velocities.setter  # type: ignore
1✔
162
    def velocities(self, velocities):
1✔
163
        """Setter for Poscar.velocities"""
164
        self.structure.add_site_property("velocities", velocities)
×
165

166
    @selective_dynamics.setter  # type: ignore
1✔
167
    def selective_dynamics(self, selective_dynamics):
1✔
168
        """Setter for Poscar.selective_dynamics"""
169
        self.structure.add_site_property("selective_dynamics", selective_dynamics)
1✔
170

171
    @predictor_corrector.setter  # type: ignore
1✔
172
    def predictor_corrector(self, predictor_corrector):
1✔
173
        """Setter for Poscar.predictor_corrector"""
174
        self.structure.add_site_property("predictor_corrector", predictor_corrector)
×
175

176
    @property
1✔
177
    def site_symbols(self):
1✔
178
        """
179
        Sequence of symbols associated with the Poscar. Similar to 6th line in
180
        vasp 5+ POSCAR.
181
        """
182
        syms = [site.specie.symbol for site in self.structure]
1✔
183
        return [a[0] for a in itertools.groupby(syms)]
1✔
184

185
    @property
1✔
186
    def natoms(self):
1✔
187
        """
188
        Sequence of number of sites of each type associated with the Poscar.
189
        Similar to 7th line in vasp 5+ POSCAR or the 6th line in vasp 4 POSCAR.
190
        """
191
        syms = [site.specie.symbol for site in self.structure]
1✔
192
        return [len(tuple(a[1])) for a in itertools.groupby(syms)]
1✔
193

194
    def __setattr__(self, name, value):
1✔
195
        if name in ("selective_dynamics", "velocities"):
1✔
196
            if value is not None and len(value) > 0:
1✔
197
                value = np.array(value)
1✔
198
                dim = value.shape
1✔
199
                if dim[1] != 3 or dim[0] != len(self.structure):
1✔
200
                    raise ValueError(name + " array must be same length as the structure.")
1✔
201
                value = value.tolist()
1✔
202
        super().__setattr__(name, value)
1✔
203

204
    @staticmethod
1✔
205
    def from_file(filename, check_for_POTCAR=True, read_velocities=True) -> Poscar:
1✔
206
        """
207
        Reads a Poscar from a file.
208

209
        The code will try its best to determine the elements in the POSCAR in
210
        the following order:
211

212
        1. If check_for_POTCAR is True, the code will try to check if a POTCAR
213
        is in the same directory as the POSCAR and use elements from that by
214
        default. (This is the VASP default sequence of priority).
215
        2. If the input file is VASP5-like and contains element symbols in the
216
        6th line, the code will use that if check_for_POTCAR is False or there
217
        is no POTCAR found.
218
        3. Failing (2), the code will check if a symbol is provided at the end
219
        of each coordinate.
220

221
        If all else fails, the code will just assign the first n elements in
222
        increasing atomic number, where n is the number of species, to the
223
        Poscar. For example, H, He, Li, .... This will ensure at least a
224
        unique element is assigned to each site and any analysis that does not
225
        require specific elemental properties should work fine.
226

227
        Args:
228
            filename (str): File name containing Poscar data.
229
            check_for_POTCAR (bool): Whether to check if a POTCAR is present
230
                in the same directory as the POSCAR. Defaults to True.
231
            read_velocities (bool): Whether to read or not velocities if they
232
                are present in the POSCAR. Default is True.
233

234
        Returns:
235
            Poscar object.
236
        """
237
        dirname = os.path.dirname(os.path.abspath(filename))
1✔
238
        names = None
1✔
239
        if check_for_POTCAR:
1✔
240
            potcars = glob.glob(os.path.join(dirname, "*POTCAR*"))
1✔
241
            if potcars:
1✔
242
                try:
1✔
243
                    potcar = Potcar.from_file(sorted(potcars)[0])
1✔
244
                    names = [sym.split("_")[0] for sym in potcar.symbols]
1✔
245
                    [get_el_sp(n) for n in names]  # ensure valid names
1✔
246
                except Exception:
×
247
                    names = None
×
248
        with zopen(filename, "rt") as f:
1✔
249
            return Poscar.from_string(f.read(), names, read_velocities=read_velocities)
1✔
250

251
    @staticmethod
1✔
252
    def from_string(data, default_names=None, read_velocities=True):
1✔
253
        """
254
        Reads a Poscar from a string.
255

256
        The code will try its best to determine the elements in the POSCAR in
257
        the following order:
258

259
        1. If default_names are supplied and valid, it will use those. Usually,
260
        default names comes from an external source, such as a POTCAR in the
261
        same directory.
262

263
        2. If there are no valid default names but the input file is VASP5-like
264
        and contains element symbols in the 6th line, the code will use that.
265

266
        3. Failing (2), the code will check if a symbol is provided at the end
267
        of each coordinate.
268

269
        If all else fails, the code will just assign the first n elements in
270
        increasing atomic number, where n is the number of species, to the
271
        Poscar. For example, H, He, Li, .... This will ensure at least a
272
        unique element is assigned to each site and any analysis that does not
273
        require specific elemental properties should work fine.
274

275
        Args:
276
            data (str): String containing Poscar data.
277
            default_names ([str]): Default symbols for the POSCAR file,
278
                usually coming from a POTCAR in the same directory.
279
            read_velocities (bool): Whether to read or not velocities if they
280
                are present in the POSCAR. Default is True.
281

282
        Returns:
283
            Poscar object.
284
        """
285
        # "^\s*$" doesn't match lines with no whitespace
286
        chunks = re.split(r"\n\s*\n", data.rstrip(), flags=re.MULTILINE)
1✔
287
        try:
1✔
288
            if chunks[0] == "":
1✔
289
                chunks.pop(0)
1✔
290
                chunks[0] = "\n" + chunks[0]
1✔
291
        except IndexError:
1✔
292
            raise ValueError("Empty POSCAR")
1✔
293

294
        # Parse positions
295
        lines = tuple(clean_lines(chunks[0].split("\n"), False))
1✔
296
        comment = lines[0]
1✔
297
        scale = float(lines[1])
1✔
298
        lattice = np.array([[float(i) for i in line.split()] for line in lines[2:5]])
1✔
299
        if scale < 0:
1✔
300
            # In vasp, a negative scale factor is treated as a volume. We need
301
            # to translate this to a proper lattice vector scaling.
302
            vol = abs(np.linalg.det(lattice))
1✔
303
            lattice *= (-scale / vol) ** (1 / 3)
1✔
304
        else:
305
            lattice *= scale
1✔
306

307
        vasp5_symbols = False
1✔
308
        try:
1✔
309
            natoms = [int(i) for i in lines[5].split()]
1✔
310
            ipos = 6
1✔
311
        except ValueError:
1✔
312
            vasp5_symbols = True
1✔
313
            symbols = lines[5].split()
1✔
314

315
            """
316
            Atoms and number of atoms in POSCAR written with vasp appear on
317
            multiple lines when atoms of the same type are not grouped together
318
            and more than 20 groups are then defined ...
319
            Example :
320
            Cr16 Fe35 Ni2
321
               1.00000000000000
322
                 8.5415010000000002   -0.0077670000000000   -0.0007960000000000
323
                -0.0077730000000000    8.5224019999999996    0.0105580000000000
324
                -0.0007970000000000    0.0105720000000000    8.5356889999999996
325
               Fe   Cr   Fe   Cr   Fe   Cr   Fe   Cr   Fe   Cr   Fe   Cr   Fe   Cr   Fe   Ni   Fe   Cr   Fe   Cr
326
               Fe   Ni   Fe   Cr   Fe
327
                 1   1   2   4   2   1   1   1     2     1     1     1     4     1     1     1     5     3     6     1
328
                 2   1   3   2   5
329
            Direct
330
              ...
331
            """
332
            nlines_symbols = 1
1✔
333
            for nlines_symbols in range(1, 11):
1✔
334
                try:
1✔
335
                    int(lines[5 + nlines_symbols].split()[0])
1✔
336
                    break
1✔
337
                except ValueError:
1✔
338
                    pass
1✔
339
            for iline_symbols in range(6, 5 + nlines_symbols):
1✔
340
                symbols.extend(lines[iline_symbols].split())
1✔
341
            natoms = []
1✔
342
            iline_natoms_start = 5 + nlines_symbols
1✔
343
            for iline_natoms in range(iline_natoms_start, iline_natoms_start + nlines_symbols):
1✔
344
                natoms.extend([int(i) for i in lines[iline_natoms].split()])
1✔
345
            atomic_symbols = []
1✔
346
            for i, nat in enumerate(natoms):
1✔
347
                atomic_symbols.extend([symbols[i]] * nat)
1✔
348
            ipos = 5 + 2 * nlines_symbols
1✔
349

350
        pos_type = lines[ipos].split()[0]
1✔
351

352
        has_selective_dynamics = False
1✔
353
        # Selective dynamics
354
        if pos_type[0] in "sS":
1✔
355
            has_selective_dynamics = True
1✔
356
            ipos += 1
1✔
357
            pos_type = lines[ipos].split()[0]
1✔
358

359
        cart = pos_type[0] in "cCkK"
1✔
360
        n_sites = sum(natoms)
1✔
361

362
        # If default_names is specified (usually coming from a POTCAR), use
363
        # them. This is in line with VASP's parsing order that the POTCAR
364
        # specified is the default used.
365
        if default_names:
1✔
366
            try:
1✔
367
                atomic_symbols = []
1✔
368
                for i, nat in enumerate(natoms):
1✔
369
                    atomic_symbols.extend([default_names[i]] * nat)
1✔
370
                vasp5_symbols = True
1✔
371
            except IndexError:
×
372
                pass
×
373

374
        if not vasp5_symbols:
1✔
375
            ind = 3 if not has_selective_dynamics else 6
1✔
376
            try:
1✔
377
                # Check if names are appended at the end of the coordinates.
378
                atomic_symbols = [l.split()[ind] for l in lines[ipos + 1 : ipos + 1 + n_sites]]
1✔
379
                # Ensure symbols are valid elements
380
                if not all(Element.is_valid_symbol(sym) for sym in atomic_symbols):
1✔
381
                    raise ValueError("Non-valid symbols detected.")
×
382
                vasp5_symbols = True
1✔
383
            except (ValueError, IndexError):
1✔
384
                # Defaulting to false names.
385
                atomic_symbols = []
1✔
386
                for i, nat in enumerate(natoms):
1✔
387
                    sym = Element.from_Z(i + 1).symbol
1✔
388
                    atomic_symbols.extend([sym] * nat)
1✔
389
                warnings.warn(f"Elements in POSCAR cannot be determined. Defaulting to false names {atomic_symbols}.")
1✔
390

391
        # read the atomic coordinates
392
        coords = []
1✔
393
        selective_dynamics = [] if has_selective_dynamics else None
1✔
394
        for i in range(n_sites):
1✔
395
            toks = lines[ipos + 1 + i].split()
1✔
396
            crd_scale = scale if cart else 1
1✔
397
            coords.append([float(j) * crd_scale for j in toks[:3]])
1✔
398
            if has_selective_dynamics:
1✔
399
                selective_dynamics.append([tok.upper()[0] == "T" for tok in toks[3:6]])
1✔
400

401
        struct = Structure(
1✔
402
            lattice,
403
            atomic_symbols,
404
            coords,
405
            to_unit_cell=False,
406
            validate_proximity=False,
407
            coords_are_cartesian=cart,
408
        )
409

410
        if read_velocities:
1✔
411
            # Parse velocities if any
412
            velocities = []
1✔
413
            if len(chunks) > 1:
1✔
414
                for line in chunks[1].strip().split("\n"):
1✔
415
                    velocities.append([float(tok) for tok in line.split()])
1✔
416

417
            # Parse the predictor-corrector data
418
            predictor_corrector = []
1✔
419
            predictor_corrector_preamble = None
1✔
420

421
            if len(chunks) > 2:
1✔
422
                lines = chunks[2].strip().split("\n")
1✔
423
                # There are 3 sets of 3xN Predictor corrector parameters
424
                # So can't be stored as a single set of "site_property"
425

426
                # First line in chunk is a key in CONTCAR
427
                # Second line is POTIM
428
                # Third line is the thermostat parameters
429
                predictor_corrector_preamble = lines[0] + "\n" + lines[1] + "\n" + lines[2]
1✔
430
                # Rest is three sets of parameters, each set contains
431
                # x, y, z predictor-corrector parameters for every atom in order
432
                lines = lines[3:]
1✔
433
                for st in range(n_sites):
1✔
434
                    d1 = [float(tok) for tok in lines[st].split()]
1✔
435
                    d2 = [float(tok) for tok in lines[st + n_sites].split()]
1✔
436
                    d3 = [float(tok) for tok in lines[st + 2 * n_sites].split()]
1✔
437
                    predictor_corrector.append([d1, d2, d3])
1✔
438
        else:
439
            velocities = None
1✔
440
            predictor_corrector = None
1✔
441
            predictor_corrector_preamble = None
1✔
442

443
        return Poscar(
1✔
444
            struct,
445
            comment,
446
            selective_dynamics,
447
            vasp5_symbols,
448
            velocities=velocities,
449
            predictor_corrector=predictor_corrector,
450
            predictor_corrector_preamble=predictor_corrector_preamble,
451
        )
452

453
    def get_string(self, direct: bool = True, vasp4_compatible: bool = False, significant_figures: int = 16) -> str:
1✔
454
        """
455
        Returns a string to be written as a POSCAR file. By default, site
456
        symbols are written, which means compatibility is for vasp >= 5.
457

458
        Args:
459
            direct (bool): Whether coordinates are output in direct or
460
                Cartesian. Defaults to True.
461
            vasp4_compatible (bool): Set to True to omit site symbols on 6th
462
                line to maintain backward vasp 4.x compatibility. Defaults
463
                to False.
464
            significant_figures (int): No. of significant figures to
465
                output all quantities. Defaults to 16. Note that positions are
466
                output in fixed point, while velocities are output in
467
                scientific format.
468

469
        Returns:
470
            String representation of POSCAR.
471
        """
472
        # This corrects for VASP really annoying bug of crashing on lattices
473
        # which have triple product < 0. We will just invert the lattice
474
        # vectors.
475
        latt = self.structure.lattice
1✔
476
        if np.linalg.det(latt.matrix) < 0:
1✔
477
            latt = Lattice(-latt.matrix)
1✔
478

479
        format_str = f"{{:{significant_figures+5}.{significant_figures}f}}"
1✔
480
        lines = [self.comment, "1.0"]
1✔
481
        for v in latt.matrix:
1✔
482
            lines.append(" ".join(format_str.format(c) for c in v))
1✔
483

484
        if self.true_names and not vasp4_compatible:
1✔
485
            lines.append(" ".join(self.site_symbols))
1✔
486
        lines.append(" ".join(map(str, self.natoms)))
1✔
487
        if self.selective_dynamics:
1✔
488
            lines.append("Selective dynamics")
1✔
489
        lines.append("direct" if direct else "cartesian")
1✔
490

491
        selective_dynamics = self.selective_dynamics
1✔
492
        for i, site in enumerate(self.structure):
1✔
493
            coords = site.frac_coords if direct else site.coords
1✔
494
            line = " ".join(format_str.format(c) for c in coords)
1✔
495
            if selective_dynamics is not None:
1✔
496
                sd = ["T" if j else "F" for j in selective_dynamics[i]]
1✔
497
                line += f" {sd[0]} {sd[1]} {sd[2]}"
1✔
498
            line += " " + site.species_string
1✔
499
            lines.append(line)
1✔
500

501
        if self.velocities:
1✔
502
            try:
1✔
503
                lines.append("")
1✔
504
                for v in self.velocities:
1✔
505
                    lines.append(" ".join(format_str.format(i) for i in v))
1✔
506
            except Exception:
×
507
                warnings.warn("Velocities are missing or corrupted.")
×
508

509
        if self.predictor_corrector:
1✔
510
            lines.append("")
1✔
511
            if self.predictor_corrector_preamble:
1✔
512
                lines.append(self.predictor_corrector_preamble)
1✔
513
                pred = np.array(self.predictor_corrector)
1✔
514
                for col in range(3):
1✔
515
                    for z in pred[:, col]:
1✔
516
                        lines.append(" ".join(format_str.format(i) for i in z))
1✔
517
            else:
518
                warnings.warn(
×
519
                    "Preamble information missing or corrupt. Writing Poscar with no predictor corrector data."
520
                )
521

522
        return "\n".join(lines) + "\n"
1✔
523

524
    def __repr__(self):
1✔
525
        return self.get_string()
×
526

527
    def __str__(self):
1✔
528
        """
529
        String representation of Poscar file.
530
        """
531
        return self.get_string()
1✔
532

533
    def write_file(self, filename: PathLike, **kwargs):
1✔
534
        """
535
        Writes POSCAR to a file. The supported kwargs are the same as those for
536
        the Poscar.get_string method and are passed through directly.
537
        """
538
        with zopen(filename, "wt") as f:
1✔
539
            f.write(self.get_string(**kwargs))
1✔
540

541
    def as_dict(self) -> dict:
1✔
542
        """
543
        :return: MSONable dict.
544
        """
545
        return {
1✔
546
            "@module": type(self).__module__,
547
            "@class": type(self).__name__,
548
            "structure": self.structure.as_dict(),
549
            "true_names": self.true_names,
550
            "selective_dynamics": np.array(self.selective_dynamics).tolist(),
551
            "velocities": self.velocities,
552
            "predictor_corrector": self.predictor_corrector,
553
            "comment": self.comment,
554
        }
555

556
    @classmethod
1✔
557
    def from_dict(cls, d: dict) -> Poscar:
1✔
558
        """
559
        :param d: Dict representation.
560
        :return: Poscar
561
        """
562
        return Poscar(
1✔
563
            Structure.from_dict(d["structure"]),
564
            comment=d["comment"],
565
            selective_dynamics=d["selective_dynamics"],
566
            true_names=d["true_names"],
567
            velocities=d.get("velocities", None),
568
            predictor_corrector=d.get("predictor_corrector", None),
569
        )
570

571
    def set_temperature(self, temperature: float):
1✔
572
        """
573
        Initializes the velocities based on Maxwell-Boltzmann distribution.
574
        Removes linear, but not angular drift (same as VASP)
575

576
        Scales the energies to the exact temperature (microcanonical ensemble)
577
        Velocities are given in A/fs. This is the vasp default when
578
        direct/cartesian is not specified (even when positions are given in
579
        direct coordinates)
580

581
        Overwrites imported velocities, if any.
582

583
        Args:
584
            temperature (float): Temperature in Kelvin.
585
        """
586
        # mean 0 variance 1
587
        velocities = np.random.randn(len(self.structure), 3)
1✔
588

589
        # in AMU, (N,1) array
590
        atomic_masses = np.array([site.specie.atomic_mass.to("kg") for site in self.structure])
1✔
591
        dof = 3 * len(self.structure) - 3
1✔
592

593
        # remove linear drift (net momentum)
594
        velocities -= np.average(atomic_masses[:, np.newaxis] * velocities, axis=0) / np.average(atomic_masses)
1✔
595

596
        # scale velocities due to atomic masses
597
        # mean 0 std proportional to sqrt(1/m)
598
        velocities /= atomic_masses[:, np.newaxis] ** (1 / 2)
1✔
599

600
        # scale velocities to get correct temperature
601
        energy = np.sum(1 / 2 * atomic_masses * np.sum(velocities**2, axis=1))
1✔
602
        scale = (temperature * dof / (2 * energy / const.k)) ** (1 / 2)
1✔
603

604
        velocities *= scale * 1e-5  # these are in A/fs
1✔
605

606
        self.temperature = temperature
1✔
607
        try:
1✔
608
            del self.structure.site_properties["selective_dynamics"]
1✔
609
        except KeyError:
1✔
610
            pass
1✔
611

612
        try:
1✔
613
            del self.structure.site_properties["predictor_corrector"]
1✔
614
        except KeyError:
1✔
615
            pass
1✔
616
        # returns as a list of lists to be consistent with the other
617
        # initializations
618

619
        self.structure.add_site_property("velocities", velocities.tolist())
1✔
620

621

622
cwd = os.path.abspath(os.path.dirname(__file__))
1✔
623
with open(os.path.join(cwd, "incar_parameters.json")) as incar_params:
1✔
624
    incar_params = json.loads(incar_params.read())
1✔
625

626

627
class BadIncarWarning(UserWarning):
1✔
628
    """
629
    Warning class for bad Incar parameters.
630
    """
631

632

633
class Incar(dict, MSONable):
1✔
634
    """
635
    INCAR object for reading and writing INCAR files. Essentially consists of
636
    a dictionary with some helper functions
637
    """
638

639
    def __init__(self, params: dict[str, Any] | None = None):
1✔
640
        """
641
        Creates an Incar object.
642

643
        Args:
644
            params (dict): A set of input parameters as a dictionary.
645
        """
646
        super().__init__()
1✔
647
        if params:
1✔
648
            # if Incar contains vector-like magmoms given as a list
649
            # of floats, convert to a list of lists
650
            if (params.get("MAGMOM") and isinstance(params["MAGMOM"][0], (int, float))) and (
1✔
651
                params.get("LSORBIT") or params.get("LNONCOLLINEAR")
652
            ):
653
                val = []
1✔
654
                for i in range(len(params["MAGMOM"]) // 3):
1✔
655
                    val.append(params["MAGMOM"][i * 3 : (i + 1) * 3])
1✔
656
                params["MAGMOM"] = val
1✔
657

658
            self.update(params)
1✔
659

660
    def __setitem__(self, key: str, val: Any):
1✔
661
        """
662
        Add parameter-val pair to Incar. Warns if parameter is not in list of
663
        valid INCAR tags. Also cleans the parameter and val by stripping
664
        leading and trailing white spaces.
665
        """
666
        super().__setitem__(
1✔
667
            key.strip(),
668
            Incar.proc_val(key.strip(), val.strip()) if isinstance(val, str) else val,
669
        )
670

671
    def as_dict(self) -> dict:
1✔
672
        """
673
        :return: MSONable dict.
674
        """
675
        d = dict(self)
1✔
676
        d["@module"] = type(self).__module__
1✔
677
        d["@class"] = type(self).__name__
1✔
678
        return d
1✔
679

680
    @classmethod
1✔
681
    def from_dict(cls, d) -> Incar:
1✔
682
        """
683
        :param d: Dict representation.
684
        :return: Incar
685
        """
686
        if d.get("MAGMOM") and isinstance(d["MAGMOM"][0], dict):
1✔
687
            d["MAGMOM"] = [Magmom.from_dict(m) for m in d["MAGMOM"]]
1✔
688
        return Incar({k: v for k, v in d.items() if k not in ("@module", "@class")})
1✔
689

690
    def get_string(self, sort_keys: bool = False, pretty: bool = False) -> str:
1✔
691
        """
692
        Returns a string representation of the INCAR. The reason why this
693
        method is different from the __str__ method is to provide options for
694
        pretty printing.
695

696
        Args:
697
            sort_keys (bool): Set to True to sort the INCAR parameters
698
                alphabetically. Defaults to False.
699
            pretty (bool): Set to True for pretty aligned output. Defaults
700
                to False.
701
        """
702
        keys = list(self)
1✔
703
        if sort_keys:
1✔
704
            keys = sorted(keys)
1✔
705
        lines = []
1✔
706
        for k in keys:
1✔
707
            if k == "MAGMOM" and isinstance(self[k], list):
1✔
708
                value = []
1✔
709

710
                if isinstance(self[k][0], (list, Magmom)) and (self.get("LSORBIT") or self.get("LNONCOLLINEAR")):
1✔
711
                    value.append(" ".join(str(i) for j in self[k] for i in j))
1✔
712
                elif self.get("LSORBIT") or self.get("LNONCOLLINEAR"):
1✔
713
                    for m, g in itertools.groupby(self[k]):
1✔
714
                        value.append(f"3*{len(tuple(g))}*{m}")
1✔
715
                else:
716
                    # float() to ensure backwards compatibility between
717
                    # float magmoms and Magmom objects
718
                    for m, g in itertools.groupby(self[k], lambda x: float(x)):
1✔
719
                        value.append(f"{len(tuple(g))}*{m}")
1✔
720

721
                lines.append([k, " ".join(value)])
1✔
722
            elif isinstance(self[k], list):
1✔
723
                lines.append([k, " ".join(map(str, self[k]))])
1✔
724
            else:
725
                lines.append([k, self[k]])
1✔
726

727
        if pretty:
1✔
728
            return str(tabulate([[l[0], "=", l[1]] for l in lines], tablefmt="plain"))
1✔
729
        return str_delimited(lines, None, " = ") + "\n"
1✔
730

731
    def __str__(self):
1✔
732
        return self.get_string(sort_keys=True, pretty=False)
1✔
733

734
    def write_file(self, filename: PathLike):
1✔
735
        """
736
        Write Incar to a file.
737

738
        Args:
739
            filename (str): filename to write to.
740
        """
741
        with zopen(filename, "wt") as f:
1✔
742
            f.write(str(self))
1✔
743

744
    @staticmethod
1✔
745
    def from_file(filename: PathLike) -> Incar:
1✔
746
        """
747
        Reads an Incar object from a file.
748

749
        Args:
750
            filename (str): Filename for file
751

752
        Returns:
753
            Incar object
754
        """
755
        with zopen(filename, "rt") as f:
1✔
756
            return Incar.from_string(f.read())
1✔
757

758
    @staticmethod
1✔
759
    def from_string(string: str) -> Incar:
1✔
760
        """
761
        Reads an Incar object from a string.
762

763
        Args:
764
            string (str): Incar string
765

766
        Returns:
767
            Incar object
768
        """
769
        lines = list(clean_lines(string.splitlines()))
1✔
770
        params = {}
1✔
771
        for line in lines:
1✔
772
            for sline in line.split(";"):
1✔
773
                m = re.match(r"(\w+)\s*=\s*(.*)", sline.strip())
1✔
774
                if m:
1✔
775
                    key = m.group(1).strip()
1✔
776
                    val = m.group(2).strip()
1✔
777
                    val = Incar.proc_val(key, val)
1✔
778
                    params[key] = val
1✔
779
        return Incar(params)
1✔
780

781
    @staticmethod
1✔
782
    def proc_val(key: str, val: Any):
1✔
783
        """
784
        Static helper method to convert INCAR parameters to proper types, e.g.,
785
        integers, floats, lists, etc.
786

787
        Args:
788
            key: INCAR parameter key
789
            val: Actual value of INCAR parameter.
790
        """
791
        list_keys = (
1✔
792
            "LDAUU",
793
            "LDAUL",
794
            "LDAUJ",
795
            "MAGMOM",
796
            "DIPOL",
797
            "LANGEVIN_GAMMA",
798
            "QUAD_EFG",
799
            "EINT",
800
        )
801
        bool_keys = (
1✔
802
            "LDAU",
803
            "LWAVE",
804
            "LSCALU",
805
            "LCHARG",
806
            "LPLANE",
807
            "LUSE_VDW",
808
            "LHFCALC",
809
            "ADDGRID",
810
            "LSORBIT",
811
            "LNONCOLLINEAR",
812
        )
813
        float_keys = (
1✔
814
            "EDIFF",
815
            "SIGMA",
816
            "TIME",
817
            "ENCUTFOCK",
818
            "HFSCREEN",
819
            "POTIM",
820
            "EDIFFG",
821
            "AGGAC",
822
            "PARAM1",
823
            "PARAM2",
824
        )
825
        int_keys = (
1✔
826
            "NSW",
827
            "NBANDS",
828
            "NELMIN",
829
            "ISIF",
830
            "IBRION",
831
            "ISPIN",
832
            "ICHARG",
833
            "NELM",
834
            "ISMEAR",
835
            "NPAR",
836
            "LDAUPRINT",
837
            "LMAXMIX",
838
            "ENCUT",
839
            "NSIM",
840
            "NKRED",
841
            "NUPDOWN",
842
            "ISPIND",
843
            "LDAUTYPE",
844
            "IVDW",
845
        )
846

847
        def smart_int_or_float(numstr):
1✔
848
            if numstr.find(".") != -1 or numstr.lower().find("e") != -1:
1✔
849
                return float(numstr)
1✔
850
            return int(numstr)
1✔
851

852
        try:
1✔
853
            if key in list_keys:
1✔
854
                output = []
1✔
855
                toks = re.findall(r"(-?\d+\.?\d*)\*?(-?\d+\.?\d*)?\*?(-?\d+\.?\d*)?", val)
1✔
856
                for tok in toks:
1✔
857
                    if tok[2] and "3" in tok[0]:
1✔
858
                        output.extend([smart_int_or_float(tok[2])] * int(tok[0]) * int(tok[1]))
1✔
859
                    elif tok[1]:
1✔
860
                        output.extend([smart_int_or_float(tok[1])] * int(tok[0]))
1✔
861
                    else:
862
                        output.append(smart_int_or_float(tok[0]))
1✔
863
                return output
1✔
864
            if key in bool_keys:
1✔
865
                m = re.match(r"^\.?([T|F|t|f])[A-Za-z]*\.?", val)
1✔
866
                if m:
1✔
867
                    return m.group(1).lower() == "t"
1✔
868

869
                raise ValueError(key + " should be a boolean type!")
×
870

871
            if key in float_keys:
1✔
872
                return float(re.search(r"^-?\d*\.?\d*[e|E]?-?\d*", val).group(0))  # type: ignore
1✔
873

874
            if key in int_keys:
1✔
875
                return int(re.match(r"^-?[0-9]+", val).group(0))  # type: ignore
1✔
876

877
        except ValueError:
×
878
            pass
×
879

880
        # Not in standard keys. We will try a hierarchy of conversions.
881
        try:
1✔
882
            val = int(val)
1✔
883
            return val
1✔
884
        except ValueError:
1✔
885
            pass
1✔
886

887
        try:
1✔
888
            val = float(val)
1✔
889
            return val
×
890
        except ValueError:
1✔
891
            pass
1✔
892

893
        if "true" in val.lower():
1✔
894
            return True
1✔
895

896
        if "false" in val.lower():
1✔
897
            return False
×
898

899
        return val.strip().capitalize()
1✔
900

901
    def diff(self, other: Incar) -> dict[str, dict[str, Any]]:
1✔
902
        """
903
        Diff function for Incar. Compares two Incars and indicates which
904
        parameters are the same and which are not. Useful for checking whether
905
        two runs were done using the same parameters.
906

907
        Args:
908
            other (Incar): The other Incar object to compare to.
909

910
        Returns:
911
            Dict of the following format:
912
            {"Same" : parameters_that_are_the_same,
913
            "Different": parameters_that_are_different}
914
            Note that the parameters are return as full dictionaries of values.
915
            E.g. {"ISIF":3}
916
        """
917
        similar_param = {}
1✔
918
        different_param = {}
1✔
919
        for k1, v1 in self.items():
1✔
920
            if k1 not in other:
1✔
921
                different_param[k1] = {"INCAR1": v1, "INCAR2": None}
1✔
922
            elif v1 != other[k1]:
1✔
923
                different_param[k1] = {"INCAR1": v1, "INCAR2": other[k1]}
1✔
924
            else:
925
                similar_param[k1] = v1
1✔
926
        for k2, v2 in other.items():
1✔
927
            if k2 not in similar_param and k2 not in different_param:
1✔
928
                if k2 not in self:
1✔
929
                    different_param[k2] = {"INCAR1": None, "INCAR2": v2}
1✔
930
        return {"Same": similar_param, "Different": different_param}
1✔
931

932
    def __add__(self, other):
1✔
933
        """
934
        Add all the values of another INCAR object to this object.
935
        Facilitates the use of "standard" INCARs.
936
        """
937
        params = dict(self.items())
×
938
        for k, v in other.items():
×
939
            if k in self and v != self[k]:
×
940
                raise ValueError("Incars have conflicting values!")
×
941
            params[k] = v
×
942
        return Incar(params)
×
943

944
    def check_params(self):
1✔
945
        """
946
        Raises a warning for nonsensical or non-existent INCAR tags and
947
        parameters. If a keyword doesn't exist (e.g. there's a typo in a
948
        keyword), your calculation will still run, however VASP will ignore the
949
        parameter without letting you know, hence why we have this Incar method.
950
        """
951
        for k, v in self.items():
1✔
952
            # First check if this parameter even exists
953
            if k not in incar_params:
1✔
954
                warnings.warn(
1✔
955
                    f"Cannot find {k} in the list of INCAR flags",
956
                    BadIncarWarning,
957
                    stacklevel=2,
958
                )
959

960
            if k in incar_params:
1✔
961
                if type(incar_params[k]).__name__ == "str":
1✔
962
                    # Now we check if this is an appropriate parameter type
963
                    if incar_params[k] == "float":
1✔
964
                        if not type(v) not in ["float", "int"]:
1✔
965
                            warnings.warn(
×
966
                                f"{k}: {v} is not real",
967
                                BadIncarWarning,
968
                                stacklevel=2,
969
                            )
970
                    elif type(v).__name__ != incar_params[k]:
1✔
971
                        warnings.warn(
1✔
972
                            f"{k}: {v} is not a {incar_params[k]}",
973
                            BadIncarWarning,
974
                            stacklevel=2,
975
                        )
976

977
                # if we have a list of possible parameters, check
978
                # if the user given parameter is in this list
979
                elif type(incar_params[k]).__name__ == "list":
1✔
980
                    if v not in incar_params[k]:
1✔
981
                        warnings.warn(
1✔
982
                            f"{k}: Cannot find {v} in the list of parameters",
983
                            BadIncarWarning,
984
                            stacklevel=2,
985
                        )
986

987

988
class Kpoints_supported_modes(Enum):
1✔
989
    """
990
    Enum type of all supported modes for Kpoint generation.
991
    """
992

993
    Automatic = 0
1✔
994
    Gamma = 1
1✔
995
    Monkhorst = 2
1✔
996
    Line_mode = 3
1✔
997
    Cartesian = 4
1✔
998
    Reciprocal = 5
1✔
999

1000
    def __str__(self):
1✔
1001
        return str(self.name)
1✔
1002

1003
    @staticmethod
1✔
1004
    def from_string(s: str) -> Kpoints_supported_modes:
1✔
1005
        """
1006
        :param s: String
1007
        :return: Kpoints_supported_modes
1008
        """
1009
        c = s.lower()[0]
1✔
1010
        for m in Kpoints_supported_modes:
1✔
1011
            if m.name.lower()[0] == c:
1✔
1012
                return m
1✔
1013
        raise ValueError(f"Can't interpret Kpoint mode {s}")
×
1014

1015

1016
class Kpoints(MSONable):
1✔
1017
    """
1018
    KPOINT reader/writer.
1019
    """
1020

1021
    supported_modes = Kpoints_supported_modes
1✔
1022

1023
    def __init__(
1✔
1024
        self,
1025
        comment: str = "Default gamma",
1026
        num_kpts: int = 0,
1027
        style: Kpoints_supported_modes = supported_modes.Gamma,
1028
        kpts: Sequence[float | int | Sequence] = ((1, 1, 1),),
1029
        kpts_shift: tuple[float, float, float] = (0, 0, 0),
1030
        kpts_weights=None,
1031
        coord_type=None,
1032
        labels=None,
1033
        tet_number: int = 0,
1034
        tet_weight: float = 0,
1035
        tet_connections=None,
1036
    ):
1037
        """
1038
        Highly flexible constructor for Kpoints object. The flexibility comes
1039
        at the cost of usability and in general, it is recommended that you use
1040
        the default constructor only if you know exactly what you are doing and
1041
        requires the flexibility. For most usage cases, the three automatic
1042
        schemes can be constructed far more easily using the convenience static
1043
        constructors (automatic, gamma_automatic, monkhorst_automatic) and it
1044
        is recommended that you use those.
1045

1046
        Args:
1047
            comment (str): String comment for Kpoints. Defaults to "Default gamma".
1048
            num_kpts: Following VASP method of defining the KPOINTS file, this
1049
                parameter is the number of kpoints specified. If set to 0
1050
                (or negative), VASP automatically generates the KPOINTS.
1051
            style: Style for generating KPOINTS. Use one of the
1052
                Kpoints.supported_modes enum types.
1053
            kpts (2D array): 2D array of kpoints. Even when only a single
1054
                specification is required, e.g. in the automatic scheme,
1055
                the kpts should still be specified as a 2D array. e.g.,
1056
                [[20]] or [[2,2,2]].
1057
            kpts_shift (3x1 array): Shift for Kpoints.
1058
            kpts_weights: Optional weights for kpoints. Weights should be
1059
                integers. For explicit kpoints.
1060
            coord_type: In line-mode, this variable specifies whether the
1061
                Kpoints were given in Cartesian or Reciprocal coordinates.
1062
            labels: In line-mode, this should provide a list of labels for
1063
                each kpt. It is optional in explicit kpoint mode as comments for
1064
                k-points.
1065
            tet_number: For explicit kpoints, specifies the number of
1066
                tetrahedrons for the tetrahedron method.
1067
            tet_weight: For explicit kpoints, specifies the weight for each
1068
                tetrahedron for the tetrahedron method.
1069
            tet_connections: For explicit kpoints, specifies the connections
1070
                of the tetrahedrons for the tetrahedron method.
1071
                Format is a list of tuples, [ (sym_weight, [tet_vertices]),
1072
                ...]
1073

1074
        The default behavior of the constructor is for a Gamma centered,
1075
        1x1x1 KPOINTS with no shift.
1076
        """
1077
        if num_kpts > 0 and (not labels) and (not kpts_weights):
1✔
1078
            raise ValueError("For explicit or line-mode kpoints, either the labels or kpts_weights must be specified.")
×
1079

1080
        self.comment = comment
1✔
1081
        self.num_kpts = num_kpts
1✔
1082
        self.kpts = kpts
1✔
1083
        self.style = style
1✔
1084
        self.coord_type = coord_type
1✔
1085
        self.kpts_weights = kpts_weights
1✔
1086
        self.kpts_shift = kpts_shift
1✔
1087
        self.labels = labels
1✔
1088
        self.tet_number = tet_number
1✔
1089
        self.tet_weight = tet_weight
1✔
1090
        self.tet_connections = tet_connections
1✔
1091

1092
    @property
1✔
1093
    def style(self):
1✔
1094
        """
1095
        :return: Style for kpoint generation. One of Kpoints_supported_modes
1096
            enum.
1097
        """
1098
        return self._style
1✔
1099

1100
    @style.setter
1✔
1101
    def style(self, style):
1✔
1102
        """
1103
        :param style: Style
1104
        :return: Sets the style for the Kpoints. One of Kpoints_supported_modes
1105
            enum.
1106
        """
1107
        if isinstance(style, str):
1✔
1108
            style = Kpoints.supported_modes.from_string(style)
1✔
1109

1110
        if (
1✔
1111
            style
1112
            in (
1113
                Kpoints.supported_modes.Automatic,
1114
                Kpoints.supported_modes.Gamma,
1115
                Kpoints.supported_modes.Monkhorst,
1116
            )
1117
            and len(self.kpts) > 1
1118
        ):
1119
            raise ValueError(
×
1120
                "For fully automatic or automatic gamma or monk "
1121
                "kpoints, only a single line for the number of "
1122
                "divisions is allowed."
1123
            )
1124

1125
        self._style = style
1✔
1126

1127
    @staticmethod
1✔
1128
    def automatic(subdivisions):
1✔
1129
        """
1130
        Convenient static constructor for a fully automatic Kpoint grid, with
1131
        gamma centered Monkhorst-Pack grids and the number of subdivisions
1132
        along each reciprocal lattice vector determined by the scheme in the
1133
        VASP manual.
1134

1135
        Args:
1136
            subdivisions: Parameter determining number of subdivisions along
1137
                each reciprocal lattice vector.
1138

1139
        Returns:
1140
            Kpoints object
1141
        """
1142
        return Kpoints(
1✔
1143
            "Fully automatic kpoint scheme",
1144
            0,
1145
            style=Kpoints.supported_modes.Automatic,
1146
            kpts=[[subdivisions]],
1147
        )
1148

1149
    @staticmethod
1✔
1150
    def gamma_automatic(kpts: tuple[int, int, int] = (1, 1, 1), shift: tuple[float, float, float] = (0, 0, 0)):
1✔
1151
        """
1152
        Convenient static constructor for an automatic Gamma centered Kpoint
1153
        grid.
1154

1155
        Args:
1156
            kpts: Subdivisions N_1, N_2 and N_3 along reciprocal lattice
1157
                vectors. Defaults to (1,1,1)
1158
            shift: Shift to be applied to the kpoints. Defaults to (0,0,0).
1159

1160
        Returns:
1161
            Kpoints object
1162
        """
1163
        return Kpoints(
1✔
1164
            "Automatic kpoint scheme",
1165
            0,
1166
            Kpoints.supported_modes.Gamma,
1167
            kpts=[kpts],
1168
            kpts_shift=shift,
1169
        )
1170

1171
    @staticmethod
1✔
1172
    def monkhorst_automatic(kpts: tuple[int, int, int] = (2, 2, 2), shift: tuple[float, float, float] = (0, 0, 0)):
1✔
1173
        """
1174
        Convenient static constructor for an automatic Monkhorst pack Kpoint
1175
        grid.
1176

1177
        Args:
1178
            kpts: Subdivisions N_1, N_2 and N_3 along reciprocal lattice
1179
                vectors. Defaults to (2,2,2)
1180
            shift: Shift to be applied to the kpoints. Defaults to (0,0,0).
1181

1182
        Returns:
1183
            Kpoints object
1184
        """
1185
        return Kpoints(
1✔
1186
            "Automatic kpoint scheme",
1187
            0,
1188
            Kpoints.supported_modes.Monkhorst,
1189
            kpts=[kpts],
1190
            kpts_shift=shift,
1191
        )
1192

1193
    @staticmethod
1✔
1194
    def automatic_density(structure: Structure, kppa: float, force_gamma: bool = False):
1✔
1195
        """
1196
        Returns an automatic Kpoint object based on a structure and a kpoint
1197
        density. Uses Gamma centered meshes for hexagonal cells and
1198
        Monkhorst-Pack grids otherwise.
1199

1200
        Algorithm:
1201
            Uses a simple approach scaling the number of divisions along each
1202
            reciprocal lattice vector proportional to its length.
1203

1204
        Args:
1205
            structure (Structure): Input structure
1206
            kppa (float): Grid density
1207
            force_gamma (bool): Force a gamma centered mesh (default is to
1208
                use gamma only for hexagonal cells or odd meshes)
1209

1210
        Returns:
1211
            Kpoints
1212
        """
1213
        comment = f"pymatgen with grid density = {kppa:.0f} / number of atoms"
1✔
1214
        if math.fabs((math.floor(kppa ** (1 / 3) + 0.5)) ** 3 - kppa) < 1:
1✔
1215
            kppa += kppa * 0.01
1✔
1216
        latt = structure.lattice
1✔
1217
        lengths = latt.abc
1✔
1218
        ngrid = kppa / structure.num_sites
1✔
1219
        mult = (ngrid * lengths[0] * lengths[1] * lengths[2]) ** (1 / 3)
1✔
1220

1221
        num_div = [int(math.floor(max(mult / l, 1))) for l in lengths]
1✔
1222

1223
        is_hexagonal = latt.is_hexagonal()
1✔
1224

1225
        has_odd = any(i % 2 == 1 for i in num_div)
1✔
1226
        if has_odd or is_hexagonal or force_gamma:
1✔
1227
            style = Kpoints.supported_modes.Gamma
1✔
1228
        else:
1229
            style = Kpoints.supported_modes.Monkhorst
1✔
1230

1231
        return Kpoints(comment, 0, style, [num_div], (0, 0, 0))
1✔
1232

1233
    @staticmethod
1✔
1234
    def automatic_gamma_density(structure: Structure, kppa: float):
1✔
1235
        """
1236
        Returns an automatic Kpoint object based on a structure and a kpoint
1237
        density. Uses Gamma centered meshes always. For GW.
1238

1239
        Algorithm:
1240
            Uses a simple approach scaling the number of divisions along each
1241
            reciprocal lattice vector proportional to its length.
1242

1243
        Args:
1244
            structure: Input structure
1245
            kppa: Grid density
1246
        """
1247
        latt = structure.lattice
×
1248
        lengths = latt.abc
×
1249
        ngrid = kppa / structure.num_sites
×
1250

1251
        mult = (ngrid * lengths[0] * lengths[1] * lengths[2]) ** (1 / 3)
×
1252
        num_div = [int(round(mult / l)) for l in lengths]
×
1253

1254
        # ensure that numDiv[i] > 0
1255
        num_div = [i if i > 0 else 1 for i in num_div]
×
1256

1257
        # VASP documentation recommends to use even grids for n <= 8 and odd
1258
        # grids for n > 8.
1259
        num_div = [i + i % 2 if i <= 8 else i - i % 2 + 1 for i in num_div]
×
1260

1261
        style = Kpoints.supported_modes.Gamma
×
1262

1263
        comment = f"pymatgen with grid density = {kppa:.0f} / number of atoms"
×
1264

1265
        num_kpts = 0
×
1266
        return Kpoints(comment, num_kpts, style, [num_div], (0, 0, 0))
×
1267

1268
    @staticmethod
1✔
1269
    def automatic_density_by_vol(structure: Structure, kppvol: int, force_gamma: bool = False):
1✔
1270
        """
1271
        Returns an automatic Kpoint object based on a structure and a kpoint
1272
        density per inverse Angstrom^3 of reciprocal cell.
1273

1274
        Algorithm:
1275
            Same as automatic_density()
1276

1277
        Args:
1278
            structure (Structure): Input structure
1279
            kppvol (int): Grid density per Angstrom^(-3) of reciprocal cell
1280
            force_gamma (bool): Force a gamma centered mesh
1281

1282
        Returns:
1283
            Kpoints
1284
        """
1285
        vol = structure.lattice.reciprocal_lattice.volume
1✔
1286
        kppa = kppvol * vol * structure.num_sites
1✔
1287
        return Kpoints.automatic_density(structure, kppa, force_gamma=force_gamma)
1✔
1288

1289
    @staticmethod
1✔
1290
    def automatic_density_by_lengths(
1✔
1291
        structure: Structure, length_densities: Sequence[float], force_gamma: bool = False
1292
    ):
1293
        """
1294
        Returns an automatic Kpoint object based on a structure and a k-point
1295
        density normalized by lattice constants.
1296

1297
        Algorithm:
1298
            For a given dimension, the # of k-points is chosen as
1299
            length_density = # of kpoints * lattice constant, e.g. [50.0, 50.0, 1.0] would
1300
            have k-points of 50/a x 50/b x 1/c.
1301

1302
        Args:
1303
            structure (Structure): Input structure
1304
            length_densities (list[floats]): Defines the density of k-points in each
1305
            dimension, e.g. [50.0, 50.0, 1.0].
1306
            force_gamma (bool): Force a gamma centered mesh
1307

1308
        Returns:
1309
            Kpoints
1310
        """
1311
        comment = f"k-point density of {length_densities}/[a, b, c]"
1✔
1312
        lattice = structure.lattice
1✔
1313
        abc = lattice.abc
1✔
1314
        num_div = [
1✔
1315
            np.ceil(length_densities[0] / abc[0]),
1316
            np.ceil(length_densities[1] / abc[1]),
1317
            np.ceil(length_densities[2] / abc[2]),
1318
        ]
1319
        is_hexagonal = lattice.is_hexagonal()
1✔
1320
        has_odd = any(i % 2 == 1 for i in num_div)
1✔
1321
        if has_odd or is_hexagonal or force_gamma:
1✔
1322
            style = Kpoints.supported_modes.Gamma
1✔
1323
        else:
1324
            style = Kpoints.supported_modes.Monkhorst
×
1325
        style = Kpoints.supported_modes.Gamma
1✔
1326

1327
        return Kpoints(comment, 0, style, [num_div], (0, 0, 0))
1✔
1328

1329
    @staticmethod
1✔
1330
    def automatic_linemode(divisions, ibz):
1✔
1331
        """
1332
        Convenient static constructor for a KPOINTS in mode line_mode.
1333
        gamma centered Monkhorst-Pack grids and the number of subdivisions
1334
        along each reciprocal lattice vector determined by the scheme in the
1335
        VASP manual.
1336

1337
        Args:
1338
            divisions: Parameter determining the number of k-points along each high symmetry line.
1339
            ibz: HighSymmKpath object (pymatgen.symmetry.bandstructure)
1340

1341
        Returns:
1342
            Kpoints object
1343
        """
1344
        kpoints = []
×
1345
        labels = []
×
1346
        for path in ibz.kpath["path"]:
×
1347
            kpoints.append(ibz.kpath["kpoints"][path[0]])
×
1348
            labels.append(path[0])
×
1349
            for i in range(1, len(path) - 1):
×
1350
                kpoints.append(ibz.kpath["kpoints"][path[i]])
×
1351
                labels.append(path[i])
×
1352
                kpoints.append(ibz.kpath["kpoints"][path[i]])
×
1353
                labels.append(path[i])
×
1354

1355
            kpoints.append(ibz.kpath["kpoints"][path[-1]])
×
1356
            labels.append(path[-1])
×
1357

1358
        return Kpoints(
×
1359
            "Line_mode KPOINTS file",
1360
            style=Kpoints.supported_modes.Line_mode,
1361
            coord_type="Reciprocal",
1362
            kpts=kpoints,
1363
            labels=labels,
1364
            num_kpts=int(divisions),
1365
        )
1366

1367
    @staticmethod
1✔
1368
    def from_file(filename):
1✔
1369
        """
1370
        Reads a Kpoints object from a KPOINTS file.
1371

1372
        Args:
1373
            filename (str): filename to read from.
1374

1375
        Returns:
1376
            Kpoints object
1377
        """
1378
        with zopen(filename, "rt") as f:
1✔
1379
            return Kpoints.from_string(f.read())
1✔
1380

1381
    @staticmethod
1✔
1382
    def from_string(string):
1✔
1383
        """
1384
        Reads a Kpoints object from a KPOINTS string.
1385

1386
        Args:
1387
            string (str): KPOINTS string.
1388

1389
        Returns:
1390
            Kpoints object
1391
        """
1392
        lines = [line.strip() for line in string.splitlines()]
1✔
1393

1394
        comment = lines[0]
1✔
1395
        num_kpts = int(lines[1].split()[0].strip())
1✔
1396
        style = lines[2].lower()[0]
1✔
1397

1398
        # Fully automatic KPOINTS
1399
        if style == "a":
1✔
1400
            return Kpoints.automatic(int(lines[3].split()[0].strip()))
1✔
1401

1402
        coord_pattern = re.compile(r"^\s*([\d+.\-Ee]+)\s+([\d+.\-Ee]+)\s+" r"([\d+.\-Ee]+)")
1✔
1403

1404
        # Automatic gamma and Monk KPOINTS, with optional shift
1405
        if style in ["g", "m"]:
1✔
1406
            kpts = [int(i) for i in lines[3].split()]
1✔
1407
            kpts_shift = (0, 0, 0)
1✔
1408
            if len(lines) > 4 and coord_pattern.match(lines[4]):
1✔
1409
                try:
1✔
1410
                    kpts_shift = [float(i) for i in lines[4].split()]
1✔
1411
                except ValueError:
×
1412
                    pass
×
1413
            return (
1✔
1414
                Kpoints.gamma_automatic(kpts, kpts_shift)
1415
                if style == "g"
1416
                else Kpoints.monkhorst_automatic(kpts, kpts_shift)
1417
            )
1418

1419
        # Automatic kpoints with basis
1420
        if num_kpts <= 0:
1✔
1421
            style = Kpoints.supported_modes.Cartesian if style in "ck" else Kpoints.supported_modes.Reciprocal
1✔
1422
            kpts = [[float(j) for j in lines[i].split()] for i in range(3, 6)]
1✔
1423
            kpts_shift = [float(i) for i in lines[6].split()]
1✔
1424
            return Kpoints(
1✔
1425
                comment=comment,
1426
                num_kpts=num_kpts,
1427
                style=style,
1428
                kpts=kpts,
1429
                kpts_shift=kpts_shift,
1430
            )
1431

1432
        # Line-mode KPOINTS, usually used with band structures
1433
        if style == "l":
1✔
1434
            coord_type = "Cartesian" if lines[3].lower()[0] in "ck" else "Reciprocal"
1✔
1435
            style = Kpoints.supported_modes.Line_mode
1✔
1436
            kpts = []
1✔
1437
            labels = []
1✔
1438
            patt = re.compile(r"([e0-9.\-]+)\s+([e0-9.\-]+)\s+([e0-9.\-]+)" r"\s*!*\s*(.*)")
1✔
1439
            for i in range(4, len(lines)):
1✔
1440
                line = lines[i]
1✔
1441
                m = patt.match(line)
1✔
1442
                if m:
1✔
1443
                    kpts.append([float(m.group(1)), float(m.group(2)), float(m.group(3))])
1✔
1444
                    labels.append(m.group(4).strip())
1✔
1445
            return Kpoints(
1✔
1446
                comment=comment,
1447
                num_kpts=num_kpts,
1448
                style=style,
1449
                kpts=kpts,
1450
                coord_type=coord_type,
1451
                labels=labels,
1452
            )
1453

1454
        # Assume explicit KPOINTS if all else fails.
1455
        style = Kpoints.supported_modes.Cartesian if style in "ck" else Kpoints.supported_modes.Reciprocal
1✔
1456
        kpts = []
1✔
1457
        kpts_weights = []
1✔
1458
        labels = []
1✔
1459
        tet_number = 0
1✔
1460
        tet_weight = 0
1✔
1461
        tet_connections = None
1✔
1462

1463
        for i in range(3, 3 + num_kpts):
1✔
1464
            toks = lines[i].split()
1✔
1465
            kpts.append([float(j) for j in toks[0:3]])
1✔
1466
            kpts_weights.append(float(toks[3]))
1✔
1467
            if len(toks) > 4:
1✔
1468
                labels.append(toks[4])
1✔
1469
            else:
1470
                labels.append(None)
1✔
1471
        try:
1✔
1472
            # Deal with tetrahedron method
1473
            if lines[3 + num_kpts].strip().lower()[0] == "t":
1✔
1474
                toks = lines[4 + num_kpts].split()
1✔
1475
                tet_number = int(toks[0])
1✔
1476
                tet_weight = float(toks[1])
1✔
1477
                tet_connections = []
1✔
1478
                for i in range(5 + num_kpts, 5 + num_kpts + tet_number):
1✔
1479
                    toks = lines[i].split()
1✔
1480
                    tet_connections.append((int(toks[0]), [int(toks[j]) for j in range(1, 5)]))
1✔
1481
        except IndexError:
1✔
1482
            pass
1✔
1483

1484
        return Kpoints(
1✔
1485
            comment=comment,
1486
            num_kpts=num_kpts,
1487
            style=Kpoints.supported_modes[str(style)],
1488
            kpts=kpts,
1489
            kpts_weights=kpts_weights,
1490
            tet_number=tet_number,
1491
            tet_weight=tet_weight,
1492
            tet_connections=tet_connections,
1493
            labels=labels,
1494
        )
1495

1496
    def write_file(self, filename):
1✔
1497
        """
1498
        Write Kpoints to a file.
1499

1500
        Args:
1501
            filename (str): Filename to write to.
1502
        """
1503
        with zopen(filename, "wt") as f:
1✔
1504
            f.write(str(self))
1✔
1505

1506
    def __repr__(self):
1✔
1507
        return str(self)
×
1508

1509
    def __str__(self):
1✔
1510
        lines = [self.comment, str(self.num_kpts), self.style.name]
1✔
1511
        style = self.style.name.lower()[0]
1✔
1512
        if style == "l":
1✔
1513
            lines.append(self.coord_type)
1✔
1514
        for idx, kpt in enumerate(self.kpts):
1✔
1515
            lines.append(" ".join(map(str, kpt)))
1✔
1516
            if style == "l":
1✔
1517
                lines[-1] += " ! " + self.labels[idx]
1✔
1518
                if idx % 2 == 1:
1✔
1519
                    lines[-1] += "\n"
1✔
1520
            elif self.num_kpts > 0:
1✔
1521
                if self.labels is not None:
1✔
1522
                    lines[-1] += f" {int(self.kpts_weights[idx])} {self.labels[idx]}"
1✔
1523
                else:
1524
                    lines[-1] += f" {int(self.kpts_weights[idx])}"
×
1525

1526
        # Print tetrahedron parameters if the number of tetrahedrons > 0
1527
        if style not in "lagm" and self.tet_number > 0:
1✔
1528
            lines.append("Tetrahedron")
×
1529
            lines.append(f"{self.tet_number} {self.tet_weight:f}")
×
1530
            for sym_weight, vertices in self.tet_connections:
×
1531
                a, b, c, d = vertices
×
1532
                lines.append(f"{sym_weight} {a} {b} {c} {d}")
×
1533

1534
        # Print shifts for automatic kpoints types if not zero.
1535
        if self.num_kpts <= 0 and tuple(self.kpts_shift) != (0, 0, 0):
1✔
1536
            lines.append(" ".join(map(str, self.kpts_shift)))
×
1537
        return "\n".join(lines) + "\n"
1✔
1538

1539
    def as_dict(self):
1✔
1540
        """
1541
        :return: MSONable dict.
1542
        """
1543
        d = {
1✔
1544
            "comment": self.comment,
1545
            "nkpoints": self.num_kpts,
1546
            "generation_style": self.style.name,
1547
            "kpoints": self.kpts,
1548
            "usershift": self.kpts_shift,
1549
            "kpts_weights": self.kpts_weights,
1550
            "coord_type": self.coord_type,
1551
            "labels": self.labels,
1552
            "tet_number": self.tet_number,
1553
            "tet_weight": self.tet_weight,
1554
            "tet_connections": self.tet_connections,
1555
        }
1556

1557
        optional_paras = ["genvec1", "genvec2", "genvec3", "shift"]
1✔
1558
        for para in optional_paras:
1✔
1559
            if para in self.__dict__:
1✔
1560
                d[para] = self.__dict__[para]
1✔
1561
        d["@module"] = type(self).__module__
1✔
1562
        d["@class"] = type(self).__name__
1✔
1563
        return d
1✔
1564

1565
    @classmethod
1✔
1566
    def from_dict(cls, d):
1✔
1567
        """
1568
        :param d: Dict representation.
1569
        :return: Kpoints
1570
        """
1571
        comment = d.get("comment", "")
1✔
1572
        generation_style = d.get("generation_style")
1✔
1573
        kpts = d.get("kpoints", [[1, 1, 1]])
1✔
1574
        kpts_shift = d.get("usershift", [0, 0, 0])
1✔
1575
        num_kpts = d.get("nkpoints", 0)
1✔
1576
        return cls(
1✔
1577
            comment=comment,
1578
            kpts=kpts,
1579
            style=generation_style,
1580
            kpts_shift=kpts_shift,
1581
            num_kpts=num_kpts,
1582
            kpts_weights=d.get("kpts_weights"),
1583
            coord_type=d.get("coord_type"),
1584
            labels=d.get("labels"),
1585
            tet_number=d.get("tet_number", 0),
1586
            tet_weight=d.get("tet_weight", 0),
1587
            tet_connections=d.get("tet_connections"),
1588
        )
1589

1590

1591
def _parse_string(s):
1✔
1592
    return f"{s.strip()}"
1✔
1593

1594

1595
def _parse_bool(s):
1✔
1596
    m = re.match(r"^\.?([TFtf])[A-Za-z]*\.?", s)
1✔
1597
    if m:
1✔
1598
        return m.group(1) in ["T", "t"]
1✔
1599
    raise ValueError(s + " should be a boolean type!")
×
1600

1601

1602
def _parse_float(s):
1✔
1603
    return float(re.search(r"^-?\d*\.?\d*[eE]?-?\d*", s).group(0))
1✔
1604

1605

1606
def _parse_int(s):
1✔
1607
    return int(re.match(r"^-?[0-9]+", s).group(0))
1✔
1608

1609

1610
def _parse_list(s):
1✔
1611
    return [float(y) for y in re.split(r"\s+", s.strip()) if not y.isalpha()]
1✔
1612

1613

1614
Orbital = namedtuple("Orbital", ["n", "l", "j", "E", "occ"])
1✔
1615
OrbitalDescription = namedtuple("OrbitalDescription", ["l", "E", "Type", "Rcut", "Type2", "Rcut2"])
1✔
1616

1617

1618
class UnknownPotcarWarning(UserWarning):
1✔
1619
    """
1620
    Warning raised when POTCAR hashes do not pass validation
1621
    """
1622

1623

1624
class PotcarSingle:
1✔
1625
    """
1626
    Object for a **single** POTCAR. The builder assumes the POTCAR contains
1627
    the complete untouched data in "data" as a string and a dict of keywords.
1628

1629
    .. attribute:: data
1630

1631
        POTCAR data as a string.
1632

1633
    .. attribute:: keywords
1634

1635
        Keywords parsed from the POTCAR as a dict. All keywords are also
1636
        accessible as attributes in themselves. E.g., potcar.enmax,
1637
        potcar.encut, etc.
1638

1639
    md5 hashes of the entire POTCAR file and the actual data are validated
1640
    against a database of known good hashes. Appropriate warnings or errors
1641
    are raised if a POTCAR hash fails validation.
1642
    """
1643

1644
    functional_dir = {
1✔
1645
        "PBE": "POT_GGA_PAW_PBE",
1646
        "PBE_52": "POT_GGA_PAW_PBE_52",
1647
        "PBE_54": "POT_GGA_PAW_PBE_54",
1648
        "LDA": "POT_LDA_PAW",
1649
        "LDA_52": "POT_LDA_PAW_52",
1650
        "LDA_54": "POT_LDA_PAW_54",
1651
        "PW91": "POT_GGA_PAW_PW91",
1652
        "LDA_US": "POT_LDA_US",
1653
        "PW91_US": "POT_GGA_US_PW91",
1654
        "Perdew-Zunger81": "POT_LDA_PAW",
1655
    }
1656

1657
    functional_tags = {
1✔
1658
        "pe": {"name": "PBE", "class": "GGA"},
1659
        "91": {"name": "PW91", "class": "GGA"},
1660
        "rp": {"name": "revPBE", "class": "GGA"},
1661
        "am": {"name": "AM05", "class": "GGA"},
1662
        "ps": {"name": "PBEsol", "class": "GGA"},
1663
        "pw": {"name": "PW86", "class": "GGA"},
1664
        "lm": {"name": "Langreth-Mehl-Hu", "class": "GGA"},
1665
        "pb": {"name": "Perdew-Becke", "class": "GGA"},
1666
        "ca": {"name": "Perdew-Zunger81", "class": "LDA"},
1667
        "hl": {"name": "Hedin-Lundquist", "class": "LDA"},
1668
        "wi": {"name": "Wigner Interpoloation", "class": "LDA"},
1669
    }
1670

1671
    parse_functions = {
1✔
1672
        "LULTRA": _parse_bool,
1673
        "LUNSCR": _parse_bool,
1674
        "LCOR": _parse_bool,
1675
        "LPAW": _parse_bool,
1676
        "EATOM": _parse_float,
1677
        "RPACOR": _parse_float,
1678
        "POMASS": _parse_float,
1679
        "ZVAL": _parse_float,
1680
        "RCORE": _parse_float,
1681
        "RWIGS": _parse_float,
1682
        "ENMAX": _parse_float,
1683
        "ENMIN": _parse_float,
1684
        "EMMIN": _parse_float,
1685
        "EAUG": _parse_float,
1686
        "DEXC": _parse_float,
1687
        "RMAX": _parse_float,
1688
        "RAUG": _parse_float,
1689
        "RDEP": _parse_float,
1690
        "RDEPT": _parse_float,
1691
        "QCUT": _parse_float,
1692
        "QGAM": _parse_float,
1693
        "RCLOC": _parse_float,
1694
        "IUNSCR": _parse_int,
1695
        "ICORE": _parse_int,
1696
        "NDATA": _parse_int,
1697
        "VRHFIN": _parse_string,
1698
        "LEXCH": _parse_string,
1699
        "TITEL": _parse_string,
1700
        "STEP": _parse_list,
1701
        "RRKJ": _parse_list,
1702
        "GGA": _parse_list,
1703
        "SHA256": _parse_string,
1704
        "COPYR": _parse_string,
1705
    }
1706

1707
    def __init__(self, data, symbol=None):
1✔
1708
        """
1709
        Args:
1710
            data:
1711
                Complete and single potcar file as a string.
1712
            symbol:
1713
                POTCAR symbol corresponding to the filename suffix
1714
                e.g. "Tm_3" for POTCAR.TM_3". If not given, pymatgen
1715
                will attempt to extract the symbol from the file itself.
1716
                However, this is not always reliable!
1717
        """
1718
        self.data = data  # raw POTCAR as a string
1✔
1719

1720
        # VASP parses header in vasprun.xml and this differs from the titel
1721
        self.header = data.split("\n")[0].strip()
1✔
1722

1723
        search_lines = re.search(
1✔
1724
            r"(?s)(parameters from PSCTR are:" r".*?END of PSCTR-controll parameters)",
1725
            data,
1726
        ).group(1)
1727

1728
        self.keywords = {}
1✔
1729
        for key, val in re.findall(r"(\S+)\s*=\s*(.*?)(?=;|$)", search_lines, flags=re.MULTILINE):
1✔
1730
            try:
1✔
1731
                self.keywords[key] = self.parse_functions[key](val)
1✔
1732
            except KeyError:
×
1733
                warnings.warn(f"Ignoring unknown variable type {key}")
×
1734

1735
        PSCTR = {}
1✔
1736

1737
        array_search = re.compile(r"(-*[0-9.]+)")
1✔
1738
        orbitals = []
1✔
1739
        descriptions = []
1✔
1740
        atomic_configuration = re.search(
1✔
1741
            r"(?s)Atomic configuration(.*?)Description",
1742
            search_lines,
1743
        )
1744
        if atomic_configuration:
1✔
1745
            lines = atomic_configuration.group(1).splitlines()
1✔
1746
            num_entries = re.search(r"([0-9]+)", lines[1]).group(1)
1✔
1747
            num_entries = int(num_entries)
1✔
1748
            PSCTR["nentries"] = num_entries
1✔
1749
            for line in lines[3:]:
1✔
1750
                orbit = array_search.findall(line)
1✔
1751
                if orbit:
1✔
1752
                    orbitals.append(
1✔
1753
                        Orbital(
1754
                            int(orbit[0]),
1755
                            int(orbit[1]),
1756
                            float(orbit[2]),
1757
                            float(orbit[3]),
1758
                            float(orbit[4]),
1759
                        )
1760
                    )
1761
            PSCTR["Orbitals"] = tuple(orbitals)
1✔
1762

1763
        description_string = re.search(
1✔
1764
            r"(?s)Description\s*\n" r"(.*?)Error from kinetic" r" energy argument \(eV\)",
1765
            search_lines,
1766
        )
1767
        if description_string:
1✔
1768
            for line in description_string.group(1).splitlines():
1✔
1769
                description = array_search.findall(line)
1✔
1770
                if description:
1✔
1771
                    descriptions.append(
1✔
1772
                        OrbitalDescription(
1773
                            int(description[0]),
1774
                            float(description[1]),
1775
                            int(description[2]),
1776
                            float(description[3]),
1777
                            int(description[4]) if len(description) > 4 else None,
1778
                            float(description[5]) if len(description) > 4 else None,
1779
                        )
1780
                    )
1781

1782
        if descriptions:
1✔
1783
            PSCTR["OrbitalDescriptions"] = tuple(descriptions)
1✔
1784

1785
        rrkj_kinetic_energy_string = re.search(
1✔
1786
            r"(?s)Error from kinetic energy argument \(eV\)\s*\n" r"(.*?)END of PSCTR-controll parameters",
1787
            search_lines,
1788
        )
1789
        rrkj_array = []
1✔
1790
        if rrkj_kinetic_energy_string:
1✔
1791
            for line in rrkj_kinetic_energy_string.group(1).splitlines():
1✔
1792
                if "=" not in line:
1✔
1793
                    rrkj_array += _parse_list(line.strip("\n"))
1✔
1794
            if rrkj_array:
1✔
1795
                PSCTR["RRKJ"] = tuple(rrkj_array)
1✔
1796

1797
        PSCTR.update(self.keywords)
1✔
1798
        self.PSCTR = dict(sorted(PSCTR.items()))
1✔
1799

1800
        if symbol:
1✔
1801
            self._symbol = symbol
1✔
1802
        else:
1803
            try:
1✔
1804
                self._symbol = self.keywords["TITEL"].split(" ")[1].strip()
1✔
1805
            except IndexError:
×
1806
                self._symbol = self.keywords["TITEL"].strip()
×
1807

1808
        # Compute the POTCAR hashes to check them against the database of known
1809
        # VASP POTCARs and possibly SHA256 hashes contained in the file itself.
1810
        self.hash = self.get_potcar_hash()
1✔
1811
        self.file_hash = self.get_potcar_file_hash()
1✔
1812
        if hasattr(self, "SHA256"):
1✔
1813
            self.hash_sha256_from_file = self.SHA256.split()[0]
1✔
1814
            self.hash_sha256_computed = self.get_sha256_file_hash()
1✔
1815

1816
        if not self.identify_potcar(mode="data")[0]:
1✔
1817
            warnings.warn(
1✔
1818
                f"POTCAR with symbol {self.symbol} has metadata that does\n"
1819
                "not match any VASP POTCAR known to pymatgen. The data in this\n"
1820
                "POTCAR is known to match the following functionals:\n"
1821
                f"{self.identify_potcar(mode='data')[0]}",
1822
                UnknownPotcarWarning,
1823
            )
1824

1825
        has_sh256, hash_check_passed = self.verify_potcar()
1✔
1826
        if not has_sh256 and not hash_check_passed:
1✔
1827
            warnings.warn(
1✔
1828
                f"POTCAR data with symbol { self.symbol} does not match any VASP "
1829
                "POTCAR known to pymatgen. There is a possibility your "
1830
                "POTCAR is corrupted or that the pymatgen database is incomplete.",
1831
                UnknownPotcarWarning,
1832
            )
1833
        elif has_sh256 and not hash_check_passed:
1✔
1834
            raise ValueError(
1✔
1835
                f"POTCAR with symbol {self.symbol} and functional\n"
1836
                f"{self.functional} has a SHA256 hash defined,\n"
1837
                "but the computed hash differs.\n"
1838
                "YOUR POTCAR FILE HAS BEEN CORRUPTED AND SHOULD NOT BE USED!"
1839
            )
1840

1841
    def __str__(self):
1✔
1842
        return self.data + "\n"
1✔
1843

1844
    @property
1✔
1845
    def electron_configuration(self):
1✔
1846
        """
1847
        :return: Electronic configuration of the PotcarSingle.
1848
        """
1849
        if not self.nelectrons.is_integer():
1✔
1850
            warnings.warn("POTCAR has non-integer charge, electron configuration not well-defined.")
×
1851
            return None
×
1852
        el = Element.from_Z(self.atomic_no)
1✔
1853
        full_config = el.full_electronic_structure
1✔
1854
        nelect = self.nelectrons
1✔
1855
        config = []
1✔
1856
        while nelect > 0:
1✔
1857
            e = full_config.pop(-1)
1✔
1858
            config.append(e)
1✔
1859
            nelect -= e[-1]
1✔
1860
        return config
1✔
1861

1862
    def write_file(self, filename: str) -> None:
1✔
1863
        """
1864
        Writes PotcarSingle to a file.
1865
        :param filename: Filename
1866
        """
1867
        with zopen(filename, "wt") as f:
×
1868
            f.write(str(self))
×
1869

1870
    @staticmethod
1✔
1871
    def from_file(filename: str) -> PotcarSingle:
1✔
1872
        """
1873
        Reads PotcarSingle from file.
1874

1875
        :param filename: Filename.
1876
        :return: PotcarSingle.
1877
        """
1878
        match = re.search(r"(?<=POTCAR\.)(.*)(?=.gz)", str(filename))
1✔
1879
        if match:
1✔
1880
            symbol = match.group(0)
1✔
1881
        else:
1882
            symbol = ""
1✔
1883

1884
        try:
1✔
1885
            with zopen(filename, "rt") as f:
1✔
1886
                return PotcarSingle(f.read(), symbol=symbol or None)
1✔
1887
        except UnicodeDecodeError:
1✔
1888
            warnings.warn("POTCAR contains invalid unicode errors. We will attempt to read it by ignoring errors.")
×
1889
            import codecs
×
1890

1891
            with codecs.open(filename, "r", encoding="utf-8", errors="ignore") as f:
×
1892
                return PotcarSingle(f.read(), symbol=symbol or None)
×
1893

1894
    @staticmethod
1✔
1895
    def from_symbol_and_functional(symbol: str, functional: str | None = None):
1✔
1896
        """
1897
        Makes a PotcarSingle from a symbol and functional.
1898

1899
        :param symbol: Symbol, e.g., Li_sv
1900
        :param functional: E.g., PBE
1901
        :return: PotcarSingle
1902
        """
1903
        functional = functional or SETTINGS.get("PMG_DEFAULT_FUNCTIONAL", "PBE")
1✔
1904
        assert isinstance(functional, str)  # mypy type narrowing
1✔
1905
        funcdir = PotcarSingle.functional_dir[functional]
1✔
1906
        d = SETTINGS.get("PMG_VASP_PSP_DIR")
1✔
1907
        if d is None:
1✔
1908
            raise ValueError(
×
1909
                f"No POTCAR for {symbol} with functional {functional} found. Please set the PMG_VASP_PSP_DIR "
1910
                "environment in .pmgrc.yaml, or you may need to set PMG_DEFAULT_FUNCTIONAL to PBE_52 or "
1911
                "PBE_54 if you are using newer psps from VASP."
1912
            )
1913
        paths_to_try = [
1✔
1914
            os.path.join(d, funcdir, f"POTCAR.{symbol}"),
1915
            os.path.join(d, funcdir, symbol, "POTCAR"),
1916
        ]
1917
        for p in paths_to_try:
1✔
1918
            p = os.path.expanduser(p)
1✔
1919
            p = zpath(p)
1✔
1920
            if os.path.exists(p):
1✔
1921
                psingle = PotcarSingle.from_file(p)
1✔
1922
                return psingle
1✔
1923
        raise OSError(
×
1924
            f"You do not have the right POTCAR with functional {functional} and label {symbol} "
1925
            f"in your VASP_PSP_DIR. Paths tried: {paths_to_try}"
1926
        )
1927

1928
    @property
1✔
1929
    def element(self):
1✔
1930
        """
1931
        Attempt to return the atomic symbol based on the VRHFIN keyword.
1932
        """
1933
        element = self.keywords["VRHFIN"].split(":")[0].strip()
1✔
1934
        try:
1✔
1935
            return Element(element).symbol
1✔
1936
        except ValueError:
×
1937
            # VASP incorrectly gives the element symbol for Xe as "X"
1938
            # Some potentials, e.g., Zr_sv, gives the symbol as r.
1939
            if element == "X":
×
1940
                return "Xe"
×
1941
            return Element(self.symbol.split("_")[0]).symbol
×
1942

1943
    @property
1✔
1944
    def atomic_no(self) -> int:
1✔
1945
        """
1946
        Attempt to return the atomic number based on the VRHFIN keyword.
1947
        """
1948
        return Element(self.element).Z
1✔
1949

1950
    @property
1✔
1951
    def nelectrons(self):
1✔
1952
        """
1953
        :return: Number of electrons
1954
        """
1955
        return self.zval
1✔
1956

1957
    @property
1✔
1958
    def symbol(self):
1✔
1959
        """
1960
        :return: The POTCAR symbol, e.g. W_pv
1961
        """
1962
        return self._symbol
1✔
1963

1964
    @property
1✔
1965
    def potential_type(self) -> str:
1✔
1966
        """
1967
        :return: Type of PSP. E.g., US, PAW, etc.
1968
        """
1969
        if self.lultra:
1✔
1970
            return "US"
×
1971
        if self.lpaw:
1✔
1972
            return "PAW"
1✔
1973
        return "NC"
×
1974

1975
    @property
1✔
1976
    def functional(self):
1✔
1977
        """
1978
        :return: Functional associated with PotcarSingle.
1979
        """
1980
        return self.functional_tags.get(self.LEXCH.lower(), {}).get("name")
1✔
1981

1982
    @property
1✔
1983
    def functional_class(self):
1✔
1984
        """
1985
        :return: Functional class associated with PotcarSingle.
1986
        """
1987
        return self.functional_tags.get(self.LEXCH.lower(), {}).get("class")
1✔
1988

1989
    def verify_potcar(self) -> tuple[bool, bool]:
1✔
1990
        """
1991
        Attempts to verify the integrity of the POTCAR data.
1992

1993
        This method checks the whole file (removing only the SHA256
1994
        metadata) against the SHA256 hash in the header if this is found.
1995
        If no SHA256 hash is found in the file, the file hash (md5 hash of the
1996
        whole file) is checked against all POTCAR file hashes known to pymatgen.
1997

1998
        Returns
1999
        -------
2000
        (bool, bool)
2001
            has_sh256 and passed_hash_check are returned.
2002

2003
        """
2004
        if hasattr(self, "SHA256"):
1✔
2005
            has_sha256 = True
1✔
2006
            if self.hash_sha256_from_file == self.hash_sha256_computed:
1✔
2007
                passed_hash_check = True
1✔
2008
            else:
2009
                passed_hash_check = False
1✔
2010
        else:
2011
            has_sha256 = False
1✔
2012
            # if no sha256 hash is found in the POTCAR file, compare the whole
2013
            # file with known potcar file hashes.
2014
            md5_file_hash = self.file_hash
1✔
2015
            hash_db = loadfn(os.path.join(cwd, "vasp_potcar_file_hashes.json"))
1✔
2016
            if md5_file_hash in hash_db:
1✔
2017
                passed_hash_check = True
1✔
2018
            else:
2019
                passed_hash_check = False
1✔
2020
        return (has_sha256, passed_hash_check)
1✔
2021

2022
    def identify_potcar(self, mode: Literal["data", "file"] = "data"):
1✔
2023
        """
2024
        Identify the symbol and compatible functionals associated with this PotcarSingle.
2025

2026
        This method checks the md5 hash of either the POTCAR metadadata (PotcarSingle.hash)
2027
        or the entire POTCAR file (PotcarSingle.file_hash) against a database
2028
        of hashes for POTCARs distributed with VASP 5.4.4.
2029

2030
        Args:
2031
            mode ('data' | 'file'): 'data' mode checks the hash of the POTCAR metadata in self.PSCTR,
2032
                while 'file' mode checks the hash of the entire POTCAR file.
2033

2034
        Returns:
2035
            symbol (List): List of symbols associated with the PotcarSingle
2036
            potcar_functionals (List): List of potcar functionals associated with
2037
                                       the PotcarSingle
2038
        """
2039
        # Dict to translate the sets in the .json file to the keys used in
2040
        # DictSet
2041
        mapping_dict = {
1✔
2042
            "potUSPP_GGA": {
2043
                "pymatgen_key": "PW91_US",
2044
                "vasp_description": "Ultrasoft pseudo potentials\
2045
                                         for LDA and PW91 (dated 2002-08-20 and 2002-04-08,\
2046
                                         respectively). These files are outdated, not\
2047
                                         supported and only distributed as is.",
2048
            },
2049
            "potUSPP_LDA": {
2050
                "pymatgen_key": "LDA_US",
2051
                "vasp_description": "Ultrasoft pseudo potentials\
2052
                                         for LDA and PW91 (dated 2002-08-20 and 2002-04-08,\
2053
                                         respectively). These files are outdated, not\
2054
                                         supported and only distributed as is.",
2055
            },
2056
            "potpaw_GGA": {
2057
                "pymatgen_key": "PW91",
2058
                "vasp_description": "The LDA, PW91 and PBE PAW datasets\
2059
                                        (snapshot: 05-05-2010, 19-09-2006 and 06-05-2010,\
2060
                                        respectively). These files are outdated, not\
2061
                                        supported and only distributed as is.",
2062
            },
2063
            "potpaw_LDA": {
2064
                "pymatgen_key": "Perdew-Zunger81",
2065
                "vasp_description": "The LDA, PW91 and PBE PAW datasets\
2066
                                        (snapshot: 05-05-2010, 19-09-2006 and 06-05-2010,\
2067
                                        respectively). These files are outdated, not\
2068
                                        supported and only distributed as is.",
2069
            },
2070
            "potpaw_LDA.52": {
2071
                "pymatgen_key": "LDA_52",
2072
                "vasp_description": "LDA PAW datasets version 52,\
2073
                                           including the early GW variety (snapshot 19-04-2012).\
2074
                                           When read by VASP these files yield identical results\
2075
                                           as the files distributed in 2012 ('unvie' release).",
2076
            },
2077
            "potpaw_LDA.54": {
2078
                "pymatgen_key": "LDA_54",
2079
                "vasp_description": "LDA PAW datasets version 54,\
2080
                                           including the GW variety (original release 2015-09-04).\
2081
                                           When read by VASP these files yield identical results as\
2082
                                           the files distributed before.",
2083
            },
2084
            "potpaw_PBE": {
2085
                "pymatgen_key": "PBE",
2086
                "vasp_description": "The LDA, PW91 and PBE PAW datasets\
2087
                                        (snapshot: 05-05-2010, 19-09-2006 and 06-05-2010,\
2088
                                        respectively). These files are outdated, not\
2089
                                        supported and only distributed as is.",
2090
            },
2091
            "potpaw_PBE.52": {
2092
                "pymatgen_key": "PBE_52",
2093
                "vasp_description": "PBE PAW datasets version 52,\
2094
                                           including early GW variety (snapshot 19-04-2012).\
2095
                                           When read by VASP these files yield identical\
2096
                                           results as the files distributed in 2012.",
2097
            },
2098
            "potpaw_PBE.54": {
2099
                "pymatgen_key": "PBE_54",
2100
                "vasp_description": "PBE PAW datasets version 54,\
2101
                                           including the GW variety (original release 2015-09-04).\
2102
                                           When read by VASP these files yield identical results as\
2103
                                           the files distributed before.",
2104
            },
2105
            "unvie_potpaw.52": {
2106
                "pymatgen_key": "unvie_LDA_52",
2107
                "vasp_description": "files released previously\
2108
                                             for vasp.5.2 (2012-04) and vasp.5.4 (2015-09-04)\
2109
                                             by univie.",
2110
            },
2111
            "unvie_potpaw.54": {
2112
                "pymatgen_key": "unvie_LDA_54",
2113
                "vasp_description": "files released previously\
2114
                                             for vasp.5.2 (2012-04) and vasp.5.4 (2015-09-04)\
2115
                                             by univie.",
2116
            },
2117
            "unvie_potpaw_PBE.52": {
2118
                "pymatgen_key": "unvie_PBE_52",
2119
                "vasp_description": "files released previously\
2120
                                                for vasp.5.2 (2012-04) and vasp.5.4 (2015-09-04)\
2121
                                                by univie.",
2122
            },
2123
            "unvie_potpaw_PBE.54": {
2124
                "pymatgen_key": "unvie_PBE_52",
2125
                "vasp_description": "files released previously\
2126
                                                for vasp.5.2 (2012-04) and vasp.5.4 (2015-09-04)\
2127
                                                by univie.",
2128
            },
2129
        }
2130

2131
        cwd = os.path.abspath(os.path.dirname(__file__))
1✔
2132

2133
        if mode == "data":
1✔
2134
            hash_db = loadfn(os.path.join(cwd, "vasp_potcar_pymatgen_hashes.json"))
1✔
2135
            potcar_hash = self.hash
1✔
2136
        elif mode == "file":
×
2137
            hash_db = loadfn(os.path.join(cwd, "vasp_potcar_file_hashes.json"))
×
2138
            potcar_hash = self.file_hash
×
2139
        else:
2140
            raise ValueError("Bad 'mode' argument. Specify 'data' or 'file'.")
×
2141

2142
        identity = hash_db.get(potcar_hash)
1✔
2143

2144
        if identity:
1✔
2145
            # convert the potcar_functionals from the .json dict into the functional
2146
            # keys that pymatgen uses
2147
            potcar_functionals = []
1✔
2148
            for i in identity["potcar_functionals"]:
1✔
2149
                potcar_functionals.append(mapping_dict[i]["pymatgen_key"])
1✔
2150
            potcar_functionals = list(set(potcar_functionals))
1✔
2151

2152
            return potcar_functionals, identity["potcar_symbols"]
1✔
2153
        return [], []
1✔
2154

2155
    def get_sha256_file_hash(self):
1✔
2156
        """
2157
        Computes a SHA256 hash of the PotcarSingle EXCLUDING lines starting with 'SHA256' and 'CPRY'
2158

2159
        This hash corresponds to the sha256 hash printed in the header of modern POTCAR files.
2160

2161
        :return: Hash value.
2162
        """
2163
        # we have to remove lines with the hash itself and the copyright
2164
        # notice to get the correct hash.
2165
        potcar_list = self.data.split("\n")
1✔
2166
        potcar_to_hash = [l for l in potcar_list if not l.strip().startswith(("SHA256", "COPYR"))]
1✔
2167
        potcar_to_hash_str = "\n".join(potcar_to_hash)
1✔
2168
        return sha256(potcar_to_hash_str.encode("utf-8")).hexdigest()
1✔
2169

2170
    def get_potcar_file_hash(self):
1✔
2171
        """
2172
        Computes a md5 hash of the entire PotcarSingle.
2173

2174
        This hash corresponds to the md5 hash of the POTCAR file itself.
2175

2176
        :return: Hash value.
2177
        """
2178
        return md5(self.data.encode("utf-8")).hexdigest()
1✔
2179

2180
    def get_potcar_hash(self):
1✔
2181
        """
2182
        Computes a md5 hash of the metadata defining the PotcarSingle.
2183

2184
        :return: Hash value.
2185
        """
2186
        hash_str = ""
1✔
2187
        for k, v in self.PSCTR.items():
1✔
2188
            # for newer POTCARS we have to exclude 'SHA256' and 'COPYR lines
2189
            # since they were not used in the initial hashing
2190
            if k in ("nentries", "Orbitals", "SHA256", "COPYR"):
1✔
2191
                continue
1✔
2192
            hash_str += f"{k}"
1✔
2193
            if isinstance(v, bool):
1✔
2194
                hash_str += f"{v}"
1✔
2195
            elif isinstance(v, float):
1✔
2196
                hash_str += f"{v:.3f}"
1✔
2197
            elif isinstance(v, int):
1✔
2198
                hash_str += f"{v}"
1✔
2199
            elif isinstance(v, (tuple, list)):
1✔
2200
                for item in v:
1✔
2201
                    if isinstance(item, float):
1✔
2202
                        hash_str += f"{item:.3f}"
1✔
2203
                    elif isinstance(item, (Orbital, OrbitalDescription)):
1✔
2204
                        for item_v in item:
1✔
2205
                            if isinstance(item_v, (int, str)):
1✔
2206
                                hash_str += f"{item_v}"
1✔
2207
                            elif isinstance(item_v, float):
1✔
2208
                                hash_str += f"{item_v:.3f}"
1✔
2209
                            else:
2210
                                hash_str += f"{item_v}" if item_v else ""
1✔
2211
            else:
2212
                hash_str += v.replace(" ", "")
1✔
2213

2214
        self.hash_str = hash_str
1✔
2215
        return md5(hash_str.lower().encode("utf-8")).hexdigest()
1✔
2216

2217
    def __getattr__(self, a):
1✔
2218
        """
2219
        Delegates attributes to keywords. For example, you can use
2220
        potcarsingle.enmax to get the ENMAX of the POTCAR.
2221

2222
        For float type properties, they are converted to the correct float. By
2223
        default, all energies in eV and all length scales are in Angstroms.
2224
        """
2225
        try:
1✔
2226
            return self.keywords[a.upper()]
1✔
2227
        except Exception:
1✔
2228
            raise AttributeError(a)
1✔
2229

2230

2231
class Potcar(list, MSONable):
1✔
2232
    """
2233
    Object for reading and writing POTCAR files for calculations. Consists of a
2234
    list of PotcarSingle.
2235
    """
2236

2237
    FUNCTIONAL_CHOICES = list(PotcarSingle.functional_dir)
1✔
2238

2239
    def __init__(self, symbols=None, functional=None, sym_potcar_map=None):
1✔
2240
        """
2241
        Args:
2242
            symbols ([str]): Element symbols for POTCAR. This should correspond
2243
                to the symbols used by VASP. E.g., "Mg", "Fe_pv", etc.
2244
            functional (str): Functional used. To know what functional options
2245
                there are, use Potcar.FUNCTIONAL_CHOICES. Note that VASP has
2246
                different versions of the same functional. By default, the old
2247
                PBE functional is used. If you want the newer ones, use PBE_52 or
2248
                PBE_54. Note that if you intend to compare your results with the
2249
                Materials Project, you should use the default setting. You can also
2250
                override the default by setting PMG_DEFAULT_FUNCTIONAL in your
2251
                .pmgrc.yaml.
2252
            sym_potcar_map (dict): Allows a user to specify a specific element
2253
                symbol to raw POTCAR mapping.
2254
        """
2255
        if functional is None:
1✔
2256
            functional = SETTINGS.get("PMG_DEFAULT_FUNCTIONAL", "PBE")
1✔
2257
        super().__init__()
1✔
2258
        self.functional = functional
1✔
2259
        if symbols is not None:
1✔
2260
            self.set_symbols(symbols, functional, sym_potcar_map)
1✔
2261

2262
    def as_dict(self):
1✔
2263
        """
2264
        :return: MSONable dict representation
2265
        """
2266
        return {
1✔
2267
            "functional": self.functional,
2268
            "symbols": self.symbols,
2269
            "@module": type(self).__module__,
2270
            "@class": type(self).__name__,
2271
        }
2272

2273
    @classmethod
1✔
2274
    def from_dict(cls, d):
1✔
2275
        """
2276
        :param d: Dict representation
2277
        :return: Potcar
2278
        """
2279
        return Potcar(symbols=d["symbols"], functional=d["functional"])
1✔
2280

2281
    @staticmethod
1✔
2282
    def from_file(filename: str):
1✔
2283
        """
2284
        Reads Potcar from file.
2285

2286
        :param filename: Filename
2287
        :return: Potcar
2288
        """
2289
        try:
1✔
2290
            with zopen(filename, "rt") as f:
1✔
2291
                fdata = f.read()
1✔
2292
        except UnicodeDecodeError:
×
2293
            warnings.warn("POTCAR contains invalid unicode errors. We will attempt to read it by ignoring errors.")
×
2294
            import codecs
×
2295

2296
            with codecs.open(filename, "r", encoding="utf-8", errors="ignore") as f:
×
2297
                fdata = f.read()
×
2298

2299
        potcar = Potcar()
1✔
2300
        potcar_strings = re.compile(r"\n?(\s*.*?End of Dataset\n)", re.S).findall(fdata)
1✔
2301
        functionals = []
1✔
2302
        for p in potcar_strings:
1✔
2303
            single = PotcarSingle(p)
1✔
2304
            potcar.append(single)
1✔
2305
            functionals.append(single.functional)
1✔
2306
        if len(set(functionals)) != 1:
1✔
2307
            raise ValueError("File contains incompatible functionals!")
×
2308
        potcar.functional = functionals[0]
1✔
2309
        return potcar
1✔
2310

2311
    def __str__(self):
1✔
2312
        return "\n".join(str(potcar).strip("\n") for potcar in self) + "\n"
1✔
2313

2314
    def write_file(self, filename):
1✔
2315
        """
2316
        Write Potcar to a file.
2317

2318
        Args:
2319
            filename (str): filename to write to.
2320
        """
2321
        with zopen(filename, "wt") as f:
1✔
2322
            f.write(str(self))
1✔
2323

2324
    @property
1✔
2325
    def symbols(self):
1✔
2326
        """
2327
        Get the atomic symbols of all the atoms in the POTCAR file.
2328
        """
2329
        return [p.symbol for p in self]
1✔
2330

2331
    @symbols.setter
1✔
2332
    def symbols(self, symbols):
1✔
2333
        self.set_symbols(symbols, functional=self.functional)
1✔
2334

2335
    @property
1✔
2336
    def spec(self):
1✔
2337
        """
2338
        Get the atomic symbols and hash of all the atoms in the POTCAR file.
2339
        """
2340
        return [{"symbol": p.symbol, "hash": p.get_potcar_hash()} for p in self]
1✔
2341

2342
    def set_symbols(self, symbols, functional=None, sym_potcar_map=None):
1✔
2343
        """
2344
        Initialize the POTCAR from a set of symbols. Currently, the POTCARs can
2345
        be fetched from a location specified in .pmgrc.yaml. Use pmg config
2346
        to add this setting.
2347

2348
        Args:
2349
            symbols ([str]): A list of element symbols
2350
            functional (str): The functional to use. If None, the setting
2351
                PMG_DEFAULT_FUNCTIONAL in .pmgrc.yaml is used, or if this is
2352
                not set, it will default to PBE.
2353
            sym_potcar_map (dict): A map of symbol:raw POTCAR string. If
2354
                sym_potcar_map is specified, POTCARs will be generated from
2355
                the given map data rather than the config file location.
2356
        """
2357
        del self[:]
1✔
2358
        if sym_potcar_map:
1✔
2359
            for el in symbols:
1✔
2360
                self.append(PotcarSingle(sym_potcar_map[el]))
1✔
2361
        else:
2362
            for el in symbols:
1✔
2363
                p = PotcarSingle.from_symbol_and_functional(el, functional)
1✔
2364
                self.append(p)
1✔
2365

2366

2367
class VaspInput(dict, MSONable):
1✔
2368
    """
2369
    Class to contain a set of vasp input objects corresponding to a run.
2370
    """
2371

2372
    def __init__(self, incar, kpoints, poscar, potcar, optional_files=None, **kwargs):
1✔
2373
        """
2374
        Args:
2375
            incar: Incar object.
2376
            kpoints: Kpoints object.
2377
            poscar: Poscar object.
2378
            potcar: Potcar object.
2379
            optional_files: Other input files supplied as a dict of {
2380
                filename: object}. The object should follow standard pymatgen
2381
                conventions in implementing a as_dict() and from_dict method.
2382
        """
2383
        super().__init__(**kwargs)
1✔
2384
        self.update({"INCAR": incar, "KPOINTS": kpoints, "POSCAR": poscar, "POTCAR": potcar})
1✔
2385
        if optional_files is not None:
1✔
2386
            self.update(optional_files)
1✔
2387

2388
    def __str__(self):
1✔
2389
        output = []
×
2390
        for k, v in self.items():
×
2391
            output.append(k)
×
2392
            output.append(str(v))
×
2393
            output.append("")
×
2394
        return "\n".join(output)
×
2395

2396
    def as_dict(self):
1✔
2397
        """
2398
        :return: MSONable dict.
2399
        """
2400
        d = {k: v.as_dict() for k, v in self.items()}
1✔
2401
        d["@module"] = type(self).__module__
1✔
2402
        d["@class"] = type(self).__name__
1✔
2403
        return d
1✔
2404

2405
    @classmethod
1✔
2406
    def from_dict(cls, d):
1✔
2407
        """
2408
        :param d: Dict representation.
2409
        :return: VaspInput
2410
        """
2411
        dec = MontyDecoder()
1✔
2412
        sub_d = {"optional_files": {}}
1✔
2413
        for k, v in d.items():
1✔
2414
            if k in ["INCAR", "POSCAR", "POTCAR", "KPOINTS"]:
1✔
2415
                sub_d[k.lower()] = dec.process_decoded(v)
1✔
2416
            elif k not in ["@module", "@class"]:
1✔
2417
                sub_d["optional_files"][k] = dec.process_decoded(v)
1✔
2418
        return cls(**sub_d)
1✔
2419

2420
    def write_input(self, output_dir=".", make_dir_if_not_present=True):
1✔
2421
        """
2422
        Write VASP input to a directory.
2423

2424
        Args:
2425
            output_dir (str): Directory to write to. Defaults to current
2426
                directory (".").
2427
            make_dir_if_not_present (bool): Create the directory if not
2428
                present. Defaults to True.
2429
        """
2430
        if make_dir_if_not_present and not os.path.exists(output_dir):
1✔
2431
            os.makedirs(output_dir)
1✔
2432
        for k, v in self.items():
1✔
2433
            if v is not None:
1✔
2434
                with zopen(os.path.join(output_dir, k), "wt") as f:
1✔
2435
                    f.write(str(v))
1✔
2436

2437
    @staticmethod
1✔
2438
    def from_directory(input_dir, optional_files=None):
1✔
2439
        """
2440
        Read in a set of VASP input from a directory. Note that only the
2441
        standard INCAR, POSCAR, POTCAR and KPOINTS files are read unless
2442
        optional_filenames is specified.
2443

2444
        Args:
2445
            input_dir (str): Directory to read VASP input from.
2446
            optional_files (dict): Optional files to read in as well as a
2447
                dict of {filename: Object type}. Object type must have a
2448
                static method from_file.
2449
        """
2450
        sub_d = {}
1✔
2451
        for fname, ftype in [
1✔
2452
            ("INCAR", Incar),
2453
            ("KPOINTS", Kpoints),
2454
            ("POSCAR", Poscar),
2455
            ("POTCAR", Potcar),
2456
        ]:
2457
            try:
1✔
2458
                fullzpath = zpath(os.path.join(input_dir, fname))
1✔
2459
                sub_d[fname.lower()] = ftype.from_file(fullzpath)
1✔
2460
            except FileNotFoundError:  # handle the case where there is no KPOINTS file
×
2461
                sub_d[fname.lower()] = None
×
2462

2463
        sub_d["optional_files"] = {}
1✔
2464
        if optional_files is not None:
1✔
2465
            for fname, ftype in optional_files.items():
1✔
2466
                sub_d["optional_files"][fname] = ftype.from_file(os.path.join(input_dir, fname))
1✔
2467
        return VaspInput(**sub_d)
1✔
2468

2469
    def run_vasp(
1✔
2470
        self,
2471
        run_dir: PathLike = ".",
2472
        vasp_cmd: list | None = None,
2473
        output_file: PathLike = "vasp.out",
2474
        err_file: PathLike = "vasp.err",
2475
    ):
2476
        """
2477
        Write input files and run VASP.
2478

2479
        :param run_dir: Where to write input files and do the run.
2480
        :param vasp_cmd: Args to be supplied to run VASP. Otherwise, the
2481
            PMG_VASP_EXE in .pmgrc.yaml is used.
2482
        :param output_file: File to write output.
2483
        :param err_file: File to write err.
2484
        """
2485
        self.write_input(output_dir=run_dir)
1✔
2486
        vasp_cmd = vasp_cmd or SETTINGS.get("PMG_VASP_EXE")  # type: ignore[assignment]
1✔
2487
        if not vasp_cmd:
1✔
2488
            raise ValueError("No VASP executable specified!")
×
2489
        vasp_cmd = [os.path.expanduser(os.path.expandvars(t)) for t in vasp_cmd]
1✔
2490
        if not vasp_cmd:
1✔
2491
            raise RuntimeError("You need to supply vasp_cmd or set the PMG_VASP_EXE in .pmgrc.yaml to run VASP.")
×
2492
        with cd(run_dir):
1✔
2493
            with open(output_file, "w") as f_std, open(err_file, "w", buffering=1) as f_err:
1✔
2494
                subprocess.check_call(vasp_cmd, stdout=f_std, stderr=f_err)
1✔
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