VASP 2.0

Introduction

This module introduces an updated version of the ASE VASP calculator, which adds the functionality of the Calculator. This allows a more general usage of the other ASE methods, such as BandStructure.

For a general introduction please refer to the vasp calculator documentation, as this is just a list of things which have changed.

Warning

This calculator is currently in BETA testing. If you are not comfortable testing new software, please use the old calculator object, see vasp.

Note

If you encounter any bugs using this calculator, please report it as an issue on the ASE github or on the ASE IRC.

Environment variables

The calculator needs to know how to execute VASP. One way of doing this, is by using the command() method, with instructions on how to execute vasp, e.g.:

Vasp2(command='mpiexec vasp_std')

which requires that the executable vasp_std is in your PATH. Alternatively, similar to the original implementation, one of the following environment variables can be set: ASE_VASP_COMMAND, VASP_COMMAND or VASP_SCRIPT - note, that the environment variables are prioritized in that order, so if ASE_VASP_COMMAND is set, the two others are ignored. The variables ASE_VASP_COMMAND or VASP_COMMAND should be commands which executes vasp. Additionally, remember to set the VASP_PP_PATH. An example shell configuration could contain

$ export ASE_VASP_COMMAND="mpiexec vasp_std"
$ export VASP_PP_PATH=$HOME/vasp/mypps

Alternatively, the VASP_SCRIPT could be used, as described in the original VASP calculator documentation.

Vasp 2.0 Calculator

The VASP specific keywords are unchanged, and should be included as described in VASP calculator. See the official VASP manual for more details on the VASP specific keywords.

class ase.calculators.vasp.vasp2.Vasp2(atoms=None, restart=None, directory='', label='vasp', ignore_bad_restart_file=False, command=None, txt=None, **kwargs)[source]

ASE interface for the Vienna Ab initio Simulation Package (VASP), with the Calculator interface.

Parameters:

atoms: object

Attach an atoms object to the calculator.

label: str

Prefix for the output file, and sets the working directory. Default is ‘vasp’.

directory: str

Set the working directory. Is prepended to label.

restart: str or bool

Sets a label for the directory to load files from. if restart=True, the working directory from label is used.

txt: bool, None, str or writable object
  • If txt is None, default ouput stream will be to PREFIX.out, where PREFIX is determined by label, i.e. the default would be vasp.out.

  • If txt is False or ‘-‘ the output will be sent through stdout

  • If txt is a string a file will be opened, and the output will be sent to that file.

  • Finally, txt can also be a an output stream, which has a ‘write’ attribute.

  • Example:

    >>> Vasp2(label='mylabel', txt=None) # Redirect stdout to :file:`mylabel.out`
    >>> Vasp2(txt='myfile.txt') # Redirect stdout to :file:`myfile.txt`
    >>> Vasp2(txt='-') # Print vasp output to stdout
    
command: str

Custom instructions on how to execute VASP. Has priority over environment variables.

Note

Parameters can be changed after the calculator has been constructed by using the set() method:

>>> calc.set(prec='Accurate', ediff=1E-5)

This would set the precision to Accurate and the break condition for the electronic SC-loop to 1E-5 eV.

Storing the calculator state

The results from the Vasp2 calculator can exported as a dictionary, which can then be saved in a JSON format, which enables easy and compressed sharing and storing of the input & outputs of a VASP calculation. The following methods of Vasp2 can be used for this purpose:

Vasp2.asdict()[source]

Return a dictionary representation of the calculator state. Does NOT contain information on the command, txt or directory keywords. Contains the following keys:

  • ase_version

  • vasp_version

  • inputs

  • results

  • atoms (Only if the calculator has an Atoms object)

Vasp2.fromdict(dct)[source]

Restore calculator from a asdicti() dictionary.

Parameters:

dct: Dictionary

The dictionary which is used to restore the calculator state.

Vasp2.write_json(filename)[source]

Dump calculator state to JSON file.

Parameters:

filename: string

The filename which the JSON file will be stored to. Prepends the directory path to the filename.

Vasp2.read_json(filename)[source]

Load Calculator state from an exported JSON Vasp2 file.

First we can dump the state of the calculation using the write_json() method:

# After a calculation
calc.write_json('mystate.json')

# This is equivalent to
from ase.io import jsonio
dct = calc.asdict()  # Get the calculator in a dictionary format
jsonio.write_json('mystate.json', dct)

At a later stage, that file can be used to restore a the input and (simple) output parameters of a calculation, without the need to copy around all the VASP specific files, using either the ase.io.jsonio.read_json() function or the Vasp2 fromdict() method.

calc = Vasp2()
calc.read_json('mystate.json')
atoms = calc.get_atoms()  # Get the atoms object

# This is equivalent to
from ase.calculators.vasp import Vasp2
from ase.io import jsonio
dct = jsonio.read_json('mystate.json')  # Load exported dict object from the JSON file
calc = Vasp2()
calc.fromdict(dct)
atoms = calc.get_atoms()  # Get the atoms object

The dictionary object, which is created from the todict() method, also contains information about the ASE and VASP version which was used at the time of the calculation, through the ase_version and vasp_version keys.

import json
with open('mystate.json', 'r') as f:
    dct = json.load(f)
print('ASE version: {}, VASP version: {}'.format(dct['ase_version'], dct['vasp_version']))

Note

The ASE calculator contains no information about the wavefunctions or charge densities, so these are NOT stored in the dictionary or JSON file, and therefore results may vary on a restarted calculation.

Examples

The Vasp 2 calculator now integrates with existing ASE functions, such as BandStructure or bandgap.

Band structure with VASP

The VASP manual has an example of creating a Si band structure - we can easily reproduce a similar result, by using the ASE Vasp2 calculator.

We can use the directory keyword to control the folder in which the calculations take place, and keep a more structured folder structure. The following script does the initial calculations, in order to construct the band structure for silicon

from ase.build import bulk
from ase.calculators.vasp import Vasp2

si = bulk('Si')

mydir = 'bandstructure'    # Directory where we will do the calculations

# Make self-consistent ground state
calc = Vasp2(kpts=(4, 4, 4), directory=mydir)

si.set_calculator(calc)
si.get_potential_energy()  # Run the calculation

# Non-SC calculation along band path
kpts = {'path': 'WGX',     # The BS path
        'npoints': 30}     # Number of points along the path

calc.set(isym=0,           # Turn off kpoint symmetry reduction
         icharg=11,        # Non-SC calculation
         kpts=kpts)

# Run the calculation
si.get_potential_energy()

As this calculation might be longer, depending on your system, it may be more convenient to split the plotting into a separate file, as all of the VASP data is written to files. The plotting can then be achieved by using the restart keyword, in a second script

from ase.calculators.vasp import Vasp2

mydir = 'bandstructure'    # Directory where we did the calculations

# Load the calculator from the VASP output files
calc_load = Vasp2(restart=True, directory=mydir)

bs = calc_load.band_structure() # ASE Band structure object
bs.plot(emin=-13)               # Plot the band structure

Which results in the following image

../../_images/vasp_si_bandstructure.png

We could also find the band gap in the same calculation,

>>> from ase.dft.bandgap import bandgap
>>> bandgap(calc_load)
Gap: 0.474 eV
Transition (v -> c):
  (s=0, k=15, n=3, [0.000, 0.000, 0.000]) -> (s=0, k=27, n=4, [0.429, 0.000, 0.429])

Note

When using hybrids, due to the exact-exchange calculations, one needs to treat the k-point sampling more carefully, see VASP HSE band structure wiki.

Currently, we have no functions to easily handle this issue, but may be added in the future.

Density of States

Vasp2 also allows for quick access to the Density of States (DOS), through the ASE DOS module, see DOS. Quick access to this function, however, can be found by using the get_dos() function:

>>> energies, dos = calc.get_dos()