pymatgen.io.lobster module

Module for reading Lobster output files. For more information on LOBSTER see www.cohp.de.

class Bandoverlaps(filename: str = 'bandOverlaps.lobster')[source]

Bases: object

Class to read in bandOverlaps.lobster files. These files are not created during every Lobster run. .. attribute: bandoverlapsdict is a dict of the following form:

{spin:{“kpoint as string”: {“maxDeviation”: float that describes the max deviation, “matrix”: 2D array of the size number of bands times number of bands including the overlap matrices with } }}

Parameters

filename – filename of the “bandOverlaps.lobster” file

has_good_quality_check_occupied_bands(number_occ_bands_spin_up: int, number_occ_bands_spin_down: Optional[int] = None, spin_polarized: bool = False, limit_deviation: float = 0.1) → bool[source]

will check if the deviation from the ideal bandoverlap of all occupied bands is smaller or equal to limit_deviation

Args: number_occ_bands_spin_up (int): number of occupied bands of spin up number_occ_bands_spin_down (int): number of occupied bands of spin down spin_polarized (bool): If True, then it was a spin polarized calculation limit_deviation (float): limit of the maxDeviation :returns: Boolean that will give you information about the quality of the projection

has_good_quality_maxDeviation(limit_maxDeviation: float = 0.1) → bool[source]

will check if the maxDeviation from the ideal bandoverlap is smaller or equal to limit_maxDeviation :param limit_maxDeviation: limit of the maxDeviation

Returns

Boolean that will give you information about the quality of the projection

class Charge(filename: str = 'CHARGE.lobster')[source]

Bases: object

Class to read CHARGE files generated by LOBSTER

Parameters

filename – filename for the CHARGE file, typically “CHARGE.lobster”

get_structure_with_charges(structure_filename)[source]

get a Structure with Mulliken and Loewdin charges as site properties :param structure_filename: filename of POSCAR

Returns

Structure Object with Mulliken and Loewdin charges as site properties

class Cohpcar(are_coops: bool = False, filename: str = None)[source]

Bases: object

Class to read COHPCAR/COOPCAR files generated by LOBSTER.

Parameters
  • are_coops – Determines if the file is a list of COHPs or COOPs. Default is False for COHPs.

  • filename – Name of the COHPCAR file. If it is None, the default file name will be chosen, depending on the value of are_coops.

class Doscar(doscar: str = 'DOSCAR.lobster', structure_file: str = 'POSCAR', dftprogram: str = 'Vasp')[source]

Bases: object

Class to deal with Lobster’s projected DOS and local projected DOS. The beforehand quantum-chemical calculation was performed with VASP

completedos

LobsterCompleteDos Object

pdos
List of Dict including numpy arrays with pdos. Access as pdos[atomindex]['orbitalstring']['Spin.up/Spin.down']
tdos
Dos Object of the total density of states
energies
numpy array of the energies at which the DOS was calculated (in eV, relative to Efermi)
tdensities
tdensities[Spin.up]: numpy array of the total density of states for the Spin.up contribution at each of the
energies
tdensities[Spin.down]: numpy array of the total density of states for the Spin.down contribution at each of the
energies

if is_spin_polarized=False: tdensities[Spin.up]: numpy array of the total density of states

is_spin_polarized
Boolean. Tells if the system is spin polarized
Parameters
  • doscar – DOSCAR filename, typically “DOSCAR.lobster”

  • structure_file – for vasp, this is typically “POSCAR”

  • dftprogram – so far only “vasp” is implemented

property completedos

CompleteDos

Type

return

property energies

Energies

Type

return

property is_spin_polarized

Whether run is spin polarized.

Type

return

property pdos

Projected DOS

Type

return

property tdensities

total densities as a list

Type

return

property tdos

Total DOS

Type

return

class Fatband(filenames='.', vasprun='vasprun.xml', Kpointsfile='KPOINTS')[source]

Bases: object

Reads in FATBAND_x_y.lobster files

Parameters
  • filenames (list or string) – can be a list of file names or a path to a folder folder from which all “FATBAND_*” files will be read

  • vasprun – corresponding vasprun file

  • Kpointsfile – KPOINTS file for bandstructure calculation, typically “KPOINTS”

get_bandstructure()[source]

returns a LobsterBandStructureSymmLine object which can be plotted with a normal BSPlotter

class Grosspop(filename: str = 'GROSSPOP.lobster')[source]

Bases: object

Class to read in GROSSPOP.lobster files.

Parameters

filename – filename of the “GROSSPOP.lobster” file

get_structure_with_total_grosspop(structure_filename: str) → pymatgen.core.structure.Structure[source]

get a Structure with Mulliken and Loewdin total grosspopulations as site properties :param structure_filename: filename of POSCAR :type structure_filename: str

Returns

Structure Object with Mulliken and Loewdin total grosspopulations as site properties

class Icohplist(are_coops: bool = False, filename: str = None)[source]

Bases: object

Class to read ICOHPLIST/ICOOPLIST files generated by LOBSTER.

Parameters
  • are_coops – Determines if the file is a list of ICOHPs or ICOOPs. Defaults to False for ICOHPs.

  • filename – Name of the ICOHPLIST file. If it is None, the default file name will be chosen, depending on the value of are_coops.

property icohpcollection

IcohpCollection object

Type

Returns

property icohplist

icohplist compatible with older version of this class

Type

Returns

class Lobsterin(settingsdict: dict)[source]

Bases: dict, monty.json.MSONable

This class can handle and generate lobsterin files Furthermore, it can also modify INCAR files for lobster, generate KPOINT files for fatband calculations in Lobster, and generate the standard primitive cells in a POSCAR file that are needed for the fatband calculations. There are also several standard lobsterin files that can be easily generated.

Parameters

settingsdict – dict to initialize Lobsterin

AVAILABLEKEYWORDS = ['COHPstartEnergy', 'COHPendEnergy', 'basisSet', 'cohpGenerator', 'gaussianSmearingWidth', 'saveProjectionToFile', 'basisfunctions', 'skipdos', 'skipcohp', 'skipcoop', 'skipPopulationAnalysis', 'skipGrossPopulation', 'userecommendedbasisfunctions', 'loadProjectionFromFile', 'forceEnergyRange', 'DensityOfEnergy', 'BWDF', 'BWDFCOHP', 'skipProjection', 'createFatband', 'writeBasisFunctions', 'writeMatricesToFile', 'realspaceHamiltonian', 'realspaceOverlap', 'printPAWRealSpaceWavefunction', 'printLCAORealSpaceWavefunction', 'noFFTforVisualization', 'RMSp', 'onlyReadVasprun.xml', 'noMemoryMappedFiles', 'skipPAWOrthonormalityTest', 'doNotIgnoreExcessiveBands', 'doNotUseAbsoluteSpilling', 'skipReOrthonormalization', 'forceV1HMatrix', 'useOriginalTetrahedronMethod', 'useDecimalPlaces', 'kSpaceCOHP']
BOOLEANKEYWORDS = ['saveProjectionToFile', 'skipdos', 'skipcohp', 'skipcoop', 'loadProjectionFromFile', 'forceEnergyRange', 'DensityOfEnergy', 'BWDF', 'BWDFCOHP', 'skipPopulationAnalysis', 'skipGrossPopulation', 'userecommendedbasisfunctions', 'skipProjection', 'writeBasisFunctions', 'writeMatricesToFile', 'noFFTforVisualization', 'RMSp', 'onlyReadVasprun.xml', 'noMemoryMappedFiles', 'skipPAWOrthonormalityTest', 'doNotIgnoreExcessiveBands', 'doNotUseAbsoluteSpilling', 'skipReOrthonormalization', 'forceV1HMatrix', 'useOriginalTetrahedronMethod', 'forceEnergyRange', 'bandwiseSpilling', 'kpointwiseSpilling']
FLOATKEYWORDS = ['COHPstartEnergy', 'COHPendEnergy', 'gaussianSmearingWidth', 'useDecimalPlaces', 'COHPSteps']
LISTKEYWORDS = ['basisfunctions', 'cohpbetween', 'createFatband']
STRINGKEYWORDS = ['basisSet', 'cohpGenerator', 'realspaceHamiltonian', 'realspaceOverlap', 'printPAWRealSpaceWavefunction', 'printLCAORealSpaceWavefunction', 'kSpaceCOHP']
as_dict()[source]
Returns

MSONable dict

diff(other)[source]

Diff function for lobsterin. Compares two lobsterin and indicates which parameters are the same. Similar to the diff in INCAR. :param other: Lobsterin object to compare to :type other: Lobsterin

Returns

dict with differences and similarities

classmethod from_dict(d)[source]
Parameters

d – Dict representation

Returns

Lobsterin

classmethod from_file(lobsterin: str)[source]
Parameters

lobsterin (str) – path to lobsterin

Returns

Lobsterin object

classmethod standard_calculations_from_vasp_files(POSCAR_input: str = 'POSCAR', INCAR_input: str = 'INCAR', POTCAR_input: Optional[str] = None, dict_for_basis: Optional[dict] = None, option: str = 'standard')[source]

will generate Lobsterin with standard settings

Parameters
  • POSCAR_input (str) – path to POSCAR

  • INCAR_input (str) – path to INCAR

  • POTCAR_input (str) – path to POTCAR

  • dict_for_basis (dict) – can be provided: it should look the following: dict_for_basis={“Fe”:’3p 3d 4s 4f’, “C”: ‘2s 2p’} and will overwrite all settings from POTCAR_input

  • option (str) – ‘standard’ will start a normal lobster run where COHPs, COOPs, DOS, CHARGE etc. will be calculated ‘standard_from_projection’ will start a normal lobster run from a projection ‘standard_with_fatband’ will do a fatband calculation, run over all orbitals ‘onlyprojection’ will only do a projection ‘onlydos’ will only calculate a projected dos ‘onlycohp’ will only calculate cohp ‘onlycoop’ will only calculate coop ‘onlycohpcoop’ will only calculate cohp and coop

Returns

Lobsterin Object with standard settings

write_INCAR(incar_input: str = 'INCAR', incar_output: str = 'INCAR.lobster', poscar_input: str = 'POSCAR', isym: int = -1, further_settings: dict = None)[source]

Will only make the run static, insert nbands, make ISYM=-1, set LWAVE=True and write a new INCAR. You have to check for the rest. :param incar_input: path to input INCAR :type incar_input: str :param incar_output: path to output INCAR :type incar_output: str :param poscar_input: path to input POSCAR :type poscar_input: str :param isym: isym equal to -1 or 0 are possible. Current Lobster version only allow -1. :type isym: int :param further_settings: A dict can be used to include further settings, e.g. {“ISMEAR”:-5} :type further_settings: dict

static write_KPOINTS(POSCAR_input: str = 'POSCAR', KPOINTS_output='KPOINTS.lobster', reciprocal_density: int = 100, isym: int = -1, from_grid: bool = False, input_grid: list = [5, 5, 5], line_mode: bool = True, kpoints_line_density: int = 20, symprec: float = 0.01)[source]

writes a KPOINT file for lobster (only ISYM=-1 and ISYM=0 are possible), grids are gamma centered :param POSCAR_input: path to POSCAR :type POSCAR_input: str :param KPOINTS_output: path to output KPOINTS :type KPOINTS_output: str :param reciprocal_density: Grid density :type reciprocal_density: int :param isym: either -1 or 0. Current Lobster versions only allow -1. :type isym: int :param from_grid: If True KPOINTS will be generated with the help of a grid given in input_grid. Otherwise,

they will be generated from the reciprocal_density

Parameters
  • input_grid (list) – grid to generate the KPOINTS file

  • line_mode (bool) – If True, band structure will be generated

  • kpoints_line_density (int) – density of the lines in the band structure

  • symprec (float) – precision to determine symmetry

static write_POSCAR_with_standard_primitive(POSCAR_input='POSCAR', POSCAR_output='POSCAR.lobster', symprec=0.01)[source]

writes a POSCAR with the standard primitive cell. This is needed to arrive at the correct kpath :param POSCAR_input: filename of input POSCAR :type POSCAR_input: str :param POSCAR_output: filename of output POSCAR :type POSCAR_output: str :param symprec: precision to find symmetry :type symprec: float

write_lobsterin(path='lobsterin', overwritedict=None)[source]

writes a lobsterin file :param path: filename of the lobsterin file that will be written :type path: str :param overwritedict: dict that can be used to overwrite lobsterin, e.g. {“skipdos”: True} :type overwritedict: dict

class Lobsterout(filename='lobsterout')[source]

Bases: object

Class to read in the lobsterout and evaluate the spilling, save the basis, save warnings, save infos

Parameters

filename – filename of lobsterout

get_doc()[source]

Returns: LobsterDict with all the information stored in lobsterout