sumo.cli package

Module contents

Package containing command line tools.

Submodules

sumo.cli.bandplot module

A script to plot electronic band structure diagrams.

Todo

  • Replace the elements and project formats with the dream syntax

sumo.cli.bandplot.bandplot(filenames=None, code='vasp', prefix=None, directory=None, vbm_cbm_marker=False, projection_selection=None, mode='rgb', normalise='all', interpolate_factor=4, color1='#FF0000', color2='#0000FF', color3='#00FF00', colorspace='lab', circle_size=150, dos_file=None, cart_coords=False, scissor=None, ylabel='Energy (eV)', dos_label=None, zero_line=False, zero_energy=None, elements=None, lm_orbitals=None, atoms=None, spin=None, total_only=False, plot_total=True, legend_cutoff=3, gaussian=None, height=None, width=None, ymin=-6.0, ymax=6.0, colours=None, yscale=1, style=None, no_base_style=False, image_format='pdf', dpi=400, plt=None, fonts=None, title=None)[source]

Plot electronic band structure diagrams from vasprun.xml files.

Parameters:
  • filenames (str or list, optional) –

    Path to input files:

    Vasp:

    Use vasprun.xml or vasprun.xml.gz file.

    Questaal:

    Path to a bnds.ext file. The extension will also be used to find site.ext and syml.ext files in the same directory.

    Castep:

    Path to a seedname.bands file. The prefix (“seedname”) is used to locate a seedname.cell file in the same directory and read in the positions of high-symmetry points.

    If no filenames are provided, sumo will search for vasprun.xml or vasprun.xml.gz files in folders named ‘split-0*’. Failing that, the code will look for a vasprun in the current directory. If a list of vasprun files is provided, these will be combined into a single band structure.

  • code (str, optional) – Calculation type. Default is ‘vasp’; ‘questaal’ and ‘castep’ also supported (with a reduced feature-set).

  • prefix (str, optional) – Prefix for file names.

  • directory (str, optional) – The directory in which to save files.

  • vbm_cbm_marker (bool, optional) – Plot markers to indicate the VBM and CBM locations.

  • projection_selection (list) –

    A list of tuple or string identifying which elements and orbitals to project on to the band structure. These can be specified by both element and orbital, for example, the following will project the Bi s, p and S p orbitals:

    [('Bi', 's'), ('Bi', 'p'), ('S', 'p')]
    

    If just the element is specified then all the orbitals of that element are combined. For example, to sum all the S orbitals:

    [('Bi', 's'), ('Bi', 'p'), 'S']
    

    You can also choose to sum particular orbitals by supplying a tuple of orbitals. For example, to sum the S s, p, and d orbitals into a single projection:

    [('Bi', 's'), ('Bi', 'p'), ('S', ('s', 'p', 'd'))]
    

    If mode = 'rgb', a maximum of 3 orbital/element combinations can be plotted simultaneously (one for red, green and blue), otherwise an unlimited number of elements/orbitals can be selected.

  • mode (str, optional) –

    Type of projected band structure to plot. Options are:

    ”rgb”

    The band structure line color depends on the character of the band. Each element/orbital contributes either red, green or blue with the corresponding line colour a mixture of all three colours. This mode only supports up to 3 elements/orbitals combinations. The order of the selection tuple determines which colour is used for each selection.

    ”stacked”

    The element/orbital contributions are drawn as a series of stacked circles, with the colour depending on the composition of the band. The size of the circles can be scaled using the circle_size option.

  • normalise (str, optional) –

    Normalisation the projections. Options are:

    • 'all': Projections normalised against the sum of all

      other projections.

    • 'select': Projections normalised against the sum of the

      selected projections.

    • None: No normalisation performed.

  • color1 (str) – A color specified in any way supported by matplotlib. Used when mode = 'rgb'.

  • color2 (str) – A color specified in any way supported by matplotlib. Used when mode = 'rgb'.

  • color3 (str) – A color specified in any way supported by matplotlib. Used when mode = 'rgb'.

  • colorspace (str) – The colorspace in which to perform the interpolation. The allowed values are rgb, hsv, lab, luvlc, lablch, and xyz. Used when mode = 'rgb'.

  • circle_size (float, optional) – The area of the circles used when mode = 'stacked'.

  • cart_coords (bool, optional) – Whether the k-points are read as cartesian or reciprocal coordinates. This is only required for Questaal output; Vasp output is less ambiguous. Defaults to False (fractional coordinates).

  • scissor (float, optional) – Apply a scissor operator (rigid shift of the CBM), use with caution if applying to metals.

  • dos_file (str, optional) – Path to vasprun.xml file from which to read the density of states information. If set, the density of states will be plotted alongside the bandstructure.

  • zero_line (bool, optional) – If true, draw a horizontal line at zero energy. (To adjust where zero energy sits, use zero_energy.)

  • zero_energy (float, optional) – Energy offset determining position of zero energy. By default, this is the VBM. It may be useful to set this to e.g. a calculated self-consistent Fermi energy.

  • elements (dict, optional) –

    The elements and orbitals to extract from the projected density of states. Should be provided as a dict with the keys as the element names and corresponding values as a tuple of orbitals. For example, the following would extract the Bi s, px, py and d orbitals:

    {'Bi': ('s', 'px', 'py', 'd')}
    

    If an element is included with an empty tuple, all orbitals for that species will be extracted. If elements is not set or set to None, all elements for all species will be extracted.

  • lm_orbitals (dict, optional) –

    The orbitals to decompose into their lm contributions (e.g. p -> px, py, pz). Should be provided as a dict, with the elements names as keys and a tuple of orbitals as the corresponding values. For example, the following would be used to decompose the oxygen p and d orbitals:

    {'O': ('p', 'd')}
    

  • atoms (dict, optional) –

    Which atomic sites to use when calculating the projected density of states. Should be provided as a dict, with the element names as keys and a tuple of int specifying the atomic indices as the corresponding values. The elemental projected density of states will be summed only over the atom indices specified. If an element is included with an empty tuple, then all sites for that element will be included. The indices are 0 based for each element specified in the POSCAR. For example, the following will calculate the density of states for the first 4 Sn atoms and all O atoms in the structure:

    {'Sn': (1, 2, 3, 4), 'O': (, )}
    

    If atoms is not set or set to None then all atomic sites for all elements will be considered.

  • spin (Spin, optional) – Plot only one spin channel from a spin-polarised calculation; “up” or “1” for spin up only, “down” or “-1” for spin down only. Defaults to None.

  • total_only (bool, optional) – Only extract the total density of states. Defaults to False.

  • plot_total (bool, optional) – Plot the total density of states. Defaults to True.

  • legend_cutoff (float, optional) – The cut-off (in % of the maximum density of states within the plotting range) for an elemental orbital to be labelled in the legend. This prevents the legend from containing labels for orbitals that have very little contribution in the plotting range.

  • gaussian (float, optional) – Broaden the density of states using convolution with a gaussian function. This parameter controls the sigma or standard deviation of the gaussian distribution.

  • height (float, optional) – The height of the plot.

  • width (float, optional) – The width of the plot.

  • ymin (float, optional) – The minimum energy on the y-axis.

  • ymax (float, optional) – The maximum energy on the y-axis.

  • style (list or str, optional) – (List of) matplotlib style specifications, to be composed on top of Sumo base style.

  • no_base_style (bool, optional) – Prevent use of sumo base style. This can make alternative styles behave more predictably.

  • colours (dict, optional) –

    Use custom colours for specific element and orbital combinations. Specified as a dict of dict of the colours. For example:

    {
        'Sn': {'s': 'r', 'p': 'b'},
        'O': {'s': '#000000'}
    }
    

    The colour can be a hex code, series of rgb value, or any other format supported by matplotlib.

  • yscale (float, optional) – Scaling factor for the y-axis.

  • image_format (str, optional) – The image file format. Can be any format supported by matplotlib, including: png, jpg, pdf, and svg. Defaults to pdf.

  • dpi (int, optional) – The dots-per-inch (pixel density) for the image.

  • plt (matplotlib.pyplot, optional) – A matplotlib.pyplot object to use for plotting.

  • fonts (list, optional) – Fonts to use in the plot. Can be a a single font, specified as a str, or several fonts, specified as a list of str.

Returns:

If plt set then the plt object will be returned. Otherwise, the method will return a list of filenames written to disk.

sumo.cli.bandplot.find_vasprun_files()[source]

Search for vasprun files from the current directory.

The precedence order for file locations is:

  1. First search for folders named: ‘split-0*’

  2. Else, look in the current directory.

The split folder names should always be zero based, therefore easily sortable.

sumo.cli.bandplot.save_data_files(bs, prefix=None, directory=None)[source]

Write the band structure data files to disk.

Parameters:
  • bs (BandStructureSymmLine) – Calculated band structure.

  • prefix (str, optional) – Prefix for data file.

  • directory (str, optional) – Directory in which to save the data.

Returns:

The filename of the written data file.

sumo.cli.bandstats module

A script to calculate effective masses and band degeneracies.

Todo

  • Would be good to reimplement get_vbm and get_cbm methods, with the ability to give degeneracies + other band edges within a certain tolerance.

  • Sample custom k-point, band, spin combinations.

sumo.cli.bandstats.bandstats(filenames=None, num_sample_points=3, temperature=None, degeneracy_tol=0.0001, parabolic=True)[source]

Calculate the effective masses of the bands of a semiconductor.

Parameters:
  • filenames (str or list, optional) – Path to vasprun.xml or vasprun.xml.gz file. If no filenames are provided, the code will search for vasprun.xml or vasprun.xml.gz files in folders named ‘split-0*’. Failing that, the code will look for a vasprun in the current directory. If a list of vasprun files is provided, these will be combined into a single band structure.

  • num_sample_points (int, optional) – Number of k-points to sample when fitting the effective masses.

  • temperature (int, optional) – Find band edges within kB * T of the valence band maximum and conduction band minimum. Not currently implemented.

  • degeneracy_tol (float, optional) – Tolerance for determining the degeneracy of the valence band maximum and conduction band minimum.

  • parabolic (bool, optional) – Use a parabolic fit of the band edges. If False then nonparabolic fitting will be attempted. Defaults to True.

Returns:

The hole and electron effective masses. Formatted as a dict with keys: 'hole_data' and 'electron_data'. The data is a list of dict with the keys:

’effective_mass’ (float)

The effective mass in units of electron rest mass, \(m_0\).

’energies’ (numpy.ndarray)

Band eigenvalues in eV.

’distances’ (numpy.ndarray)

Distances of the k-points in reciprocal space.

’band_id’ (int)

The index of the band,

’spin’ (Spin)

The spin channel

’start_kpoint’ (int)

The index of the k-point at which the band extrema occurs

’end_kpoint’ (int)

Return type:

dict

sumo.cli.dosplot module

A script to plot density of states diagrams.

Todo

  • Add ability to scale an orbitals density of states

sumo.cli.dosplot.dosplot(filename=None, code='vasp', prefix=None, directory=None, elements=None, lm_orbitals=None, atoms=None, spin=None, subplot=False, shift=True, total_only=False, plot_total=True, legend_on=True, legend_frame_on=False, legend_cutoff=3.0, gaussian=None, colours=None, height=6.0, width=8.0, xmin=-6.0, xmax=6.0, num_columns=2, xlabel='Energy (eV)', ylabel='DOS', yscale=1, zero_energy=None, zero_line=False, style=None, no_base_style=False, image_format='pdf', dpi=400, plt=None, fonts=None)[source]

A script to plot the density of states from a vasprun.xml file.

Parameters:
  • filename (str, optional) – Path to a DOS data file (can be gzipped). The preferred file type depends on the electronic structure code: vasprun.xml (VASP); .bands (CASTEP); dos. (Questaal).

  • code (str, optional) – Electronic structure code used (‘vasp’, ‘castep’ or ‘questaal’). Note that for Castep only a rough TDOS is available, assembled by sampling the eigenvalues.

  • prefix (str, optional) – Prefix for file names.

  • directory (str, optional) – The directory in which to save files.

  • elements (dict, optional) –

    The elements and orbitals to extract from the projected density of states. Should be provided as a dict with the keys as the element names and corresponding values as a tuple of orbitals. For example, the following would extract the Bi s, px, py and d orbitals:

    {'Bi': ('s', 'px', 'py', 'd')}
    

    If an element is included with an empty tuple, all orbitals for that species will be extracted. If elements is not set or set to None, all elements for all species will be extracted.

  • lm_orbitals (dict, optional) –

    The orbitals to decompose into their lm contributions (e.g. p -> px, py, pz). Should be provided as a dict, with the elements names as keys and a tuple of orbitals as the corresponding values. For example, the following would be used to decompose the oxygen p and d orbitals:

    {'O': ('p', 'd')}
    

  • atoms (dict, optional) –

    Which atomic sites to use when calculating the projected density of states. Should be provided as a dict, with the element names as keys and a tuple of int specifying the atomic indices as the corresponding values. The elemental projected density of states will be summed only over the atom indices specified. If an element is included with an empty tuple, then all sites for that element will be included. The indices are 0 based for each element specified in the POSCAR. For example, the following will calculate the density of states for the first 4 Sn atoms and all O atoms in the structure:

    {'Sn': (1, 2, 3, 4), 'O': (, )}
    

    If atoms is not set or set to None then all atomic sites for all elements will be considered.

  • spin (Spin, optional) – Plot only one spin channel from a spin-polarised calculation; “up” or “1” for spin up only, “down” or “-1” for spin down only. Defaults to None.

  • subplot (bool, optional) – Plot the density of states for each element on separate subplots. Defaults to False.

  • shift (bool, optional) – Shift the energies such that the valence band maximum (or Fermi level for metals) is at 0 eV. Defaults to True.

  • total_only (bool, optional) – Only extract the total density of states. Defaults to False.

  • plot_total (bool, optional) – Plot the total density of states. Defaults to True.

  • legend_on (bool, optional) – Plot the graph legend. Defaults to True.

  • legend_frame_on (bool, optional) – Plot a frame around the graph legend. Defaults to False.

  • legend_cutoff (float, optional) – The cut-off (in % of the maximum density of states within the plotting range) for an elemental orbital to be labelled in the legend. This prevents the legend from containing labels for orbitals that have very little contribution in the plotting range.

  • gaussian (float, optional) – Broaden the density of states using convolution with a gaussian function. This parameter controls the sigma or standard deviation of the gaussian distribution.

  • height (float, optional) – The height of the plot.

  • width (float, optional) – The width of the plot.

  • xmin (float, optional) – The minimum energy on the x-axis.

  • xmax (float, optional) – The maximum energy on the x-axis.

  • num_columns (int, optional) – The number of columns in the legend.

  • colours (dict, optional) –

    Use custom colours for specific element and orbital combinations. Specified as a dict of dict of the colours. For example:

    {
        'Sn': {'s': 'r', 'p': 'b'},
        'O': {'s': '#000000'}
    }
    

    The colour can be a hex code, series of rgb value, or any other format supported by matplotlib.

  • xlabel (str, optional) – Label/units for x-axis (i.e. energy)

  • ylabel (str, optional) – Label/units for y-axis (i.e. DOS)

  • yscale (float, optional) – Scaling factor for the y-axis.

  • zero_line (bool, optional) – Plot vertical line at energy zero.

  • zero_energy (float, optional) – Zero energy reference (e.g. Fermi energy from sc-fermi.) If not given, behaviour is determined by boolean shift.

  • style (list or str, optional) – (List of) matplotlib style specifications, to be composed on top of Sumo base style.

  • no_base_style (bool, optional) – Prevent use of sumo base style. This can make alternative styles behave more predictably.

  • image_format (str, optional) – The image file format. Can be any format supported by matplotlib, including: png, jpg, pdf, and svg. Defaults to pdf.

  • dpi (int, optional) – The dots-per-inch (pixel density) for the image.

  • plt (matplotlib.pyplot, optional) – A matplotlib.pyplot object to use for plotting.

  • fonts (list, optional) – Fonts to use in the plot. Can be a a single font, specified as a str, or several fonts, specified as a list of str.

Returns:

A matplotlib pyplot object.

sumo.cli.kgen module

A script to generate input files for band structure calculations

Todo

  • Add support for CDML labels

  • Save Brillouin zone diagram

  • Return a list of filenames/folders

sumo.cli.kgen.kgen(filename='POSCAR', code='vasp', directory=None, make_folders=False, symprec=0.01, kpts_per_split=None, ibzkpt=None, spg=None, density=60, phonon=False, mode='bradcrack', cart_coords=False, kpt_list=None, labels=None)[source]

Generate KPOINTS files for VASP band structure calculations.

This script provides a wrapper around several frameworks used to generate k-points along a high-symmetry path. The paths found in Bradley and Cracknell, SeeK-path, and pymatgen are all supported.

It is important to note that the standard primitive cell symmetry is different between SeeK-path and pymatgen. If the correct the structure is not used, the high-symmetry points (and band path) may be invalid.

Parameters:
  • filename (str, optional) – Path to VASP structure file. Default is POSCAR.

  • code (str, optional) – Calculation type. Default is ‘vasp’; ‘questaal’ also supported.

  • directory (str, optional) – The output file directory.

  • make_folders (bool, optional) – Generate folders and copy in required files (INCAR, POTCAR, POSCAR, and possibly CHGCAR) from the current directory.

  • symprec (float, optional) – The precision used for determining the cell symmetry.

  • kpts_per_split (int, optional) – If set, the k-points are split into separate k-point files (or folders) each containing the number of k-points specified. This is useful for hybrid band structure calculations where it is often intractable to calculate all k-points in the same calculation.

  • ibzkpt (str, optional) – Path to IBZKPT file. If set, the generated k-points will be appended to the k-points in this file and given a weight of 0. This is necessary for hybrid band structure calculations.

  • spg (str or int, optional) – The space group international number or symbol to override the symmetry determined by spglib. This is not recommended and only provided for testing purposes. This option will only take effect when mode = 'bradcrack'.

  • line_density (int, optional) – Density of k-points along the path.

  • phonon (bool, optional) – Write phonon q-point path instead of k-points. (Not appropriate for most codes.)

  • mode (str, optional) –

    Method used for calculating the high-symmetry path. The options are:

    bradcrack

    Use the paths from Bradley and Cracknell. See [brad].

    pymatgen

    Use the paths from pymatgen. See [curt].

    latimer-munro

    Use the paths from Latimer & Munro. See [lm].

    seekpath

    Use the paths from SeeK-path. See [seek].

  • cart_coords (bool, optional) – Whether the k-points are returned in cartesian or reciprocal coordinates. Defaults to False (fractional coordinates).

  • kpt_list (list, optional) –

    List of k-points to use, formatted as a list of subpaths, each containing a list of fractional k-points. For example:

    [ [[0., 0., 0.], [0., 0., 0.5]],
      [[0.5, 0., 0.], [0.5, 0.5, 0.]] ]
    

    Will return points along 0 0 0 -> 0 0 1/2 | 1/2 0 0 -> 1/2 1/2 0

  • path_labels (list, optional) –

    The k-point labels. These should be provided as a list of str for each subpath of the overall path. For example:

    [ ['Gamma', 'Z'], ['X', 'M'] ]
    

    combined with the above example for kpt_list would indicate the path: Gamma -> Z | X -> M. If no labels are provided, letters from A -> Z will be used instead. If a label begins with ‘@’ it will be concealed when plotting with sumo-bandplot.

sumo.cli.optplot module

A script to calculate and plot optical spectra from ab initio calculations.

sumo.cli.optplot.optplot(modes=('absorption',), filenames=None, codes='vasp', prefix=None, directory=None, gaussian=None, band_gaps=None, labels=None, average=True, height=6, width=6, xmin=0, xmax=None, ymin=0, ymax=100000.0, colours=None, style=None, no_base_style=None, image_format='pdf', dpi=400, plt=None, fonts=None, units='eV')[source]

A script to plot optical absorption spectra from VASP calculations.

Parameters:
  • modes (list or tuple) – Ordered list of str determining properties to plot. Accepted options are ‘absorption’ (default), ‘eps_real’, ‘eps_imag’, ‘n_real’, ‘n_imag’, ‘loss’ (equivalent to n_imag).

  • filenames (str or list, optional) – Path to data file. For VASP this is a vasprun.xml file (can be gzipped); for Questaal the opt.ext file from lmf or eps_BSE.out from bethesalpeter may be used. Alternatively, a list of paths can be provided, in which case the absorption spectra for each will be plotted concurrently.

  • codes (str or list, optional) – Original calculator. Accepted values are ‘vasp’ and ‘questaal’. Items should correspond to filenames.

  • prefix (str, optional) – Prefix for file names.

  • directory (str, optional) – The directory in which to save files.

  • gaussian (float) – Standard deviation for gaussian broadening.

  • band_gaps (float, str or list, optional) – The band gap as a float, in eV, plotted as a dashed line. If plotting multiple spectra then a list of band gaps can be provided. Band gaps can be provided as a floating-point number or as a path to a vasprun.xml file. To skip over a line, set its bandgap to zero or a negative number to place it outside the visible range.

  • labels (str or list) – A label to identify the spectra. If plotting multiple spectra then a list of labels can be provided.

  • average (bool, optional) – Average the dielectric response across all lattice directions. Defaults to True.

  • height (float, optional) – The height of the plot.

  • width (float, optional) – The width of the plot.

  • xmin (float, optional) – The minimum energy on the x-axis.

  • xmax (float, optional) – The maximum energy on the x-axis.

  • ymin (float, optional) – The minimum absorption intensity on the y-axis.

  • ymax (float, optional) – The maximum absorption intensity on the y-axis.

  • colours (list, optional) – A list of colours to use in the plot. The colours can be specified as a hex code, set of rgb values, or any other format supported by matplotlib.

  • style (list or str, optional) – (List of) matplotlib style specifications, to be composed on top of Sumo base style.

  • no_base_style (bool, optional) – Prevent use of sumo base style. This can make alternative styles behave more predictably.

  • image_format (str, optional) – The image file format. Can be any format supported by matplotlib, including: png, jpg, pdf, and svg. Defaults to pdf.

  • dpi (int, optional) – The dots-per-inch (pixel density) for the image.

  • plt (matplotlib.pyplot, optional) – A matplotlib.pyplot object to use for plotting.

  • fonts (list, optional) – Fonts to use in the plot. Can be a a single font, specified as a str, or several fonts, specified as a list of str.

  • units (str, optional) – X-axis units for the plot. ‘eV’ for energy in electronvolts or ‘nm’ for wavelength in nanometers. Defaults to ‘eV’.

Returns:

A matplotlib pyplot object.

sumo.cli.phonon_bandplot module

A script to plot phonon band structure diagrams.

Todo

  • automatically plot dos if present in band.yaml

  • make band structure from vasprun displacement/DFPT files

  • deal with magnetic moments

  • Read force_constants.hdf5

  • read settings from phonopy config file

  • prefix file names

sumo.cli.phonon_bandplot.phonon_bandplot(filename, poscar=None, prefix=None, directory=None, dim=None, born=None, qmesh=None, spg=None, primitive_axis=None, line_density=60, units='THz', symprec=0.01, mode='bradcrack', kpt_list=None, eigenvectors=None, labels=None, height=6.0, width=6.0, style=None, no_base_style=False, ymin=None, ymax=None, image_format='pdf', dpi=400, plt=None, fonts=None, dos=None, to_json=None, from_json=None, to_web=None, legend=None)[source]

A script to plot phonon band structure diagrams.

Parameters:
  • filename (str) – Path to phonopy output. Can be a band structure yaml file, FORCE_SETS, FORCE_CONSTANTS, or force_constants.hdf5.

  • poscar (str, optional) – Path to POSCAR file of unitcell. Not required if plotting the phonon band structure from a yaml file. If not specified, the script will search for a POSCAR file in the current directory.

  • prefix (str, optional) – Prefix for file names.

  • directory (str, optional) – The directory in which to save files.

  • dim (list, optional) – supercell matrix

  • born (str, optional) – Path to file containing Born effective charges. Should be in the same format as the file produced by the phonopy-vasp-born script provided by phonopy.

  • qmesh (list of int, optional) – Q-point mesh to use for calculating the density of state. Formatted as a 3x1 list of int.

  • spg (str or int, optional) – The space group international number or symbol to override the symmetry determined by spglib. This is not recommended and only provided for testing purposes. This option will only take effect when mode = 'bradcrack'.

  • primitive_axis (list or str, optional) – The transformation matrix from the conventional to primitive cell. Only required when the conventional cell was used as the starting structure. Should be provided as a 3x3 list of float. Alternatively the string ‘auto’ may be used to request that a primitive matrix is determined automatically.

  • line_density (int, optional) – Density of k-points along the path.

  • units (str, optional) – Units of phonon frequency. Accepted (case-insensitive) values are Thz, cm-1, eV, meV.

  • symprec (float, optional) – Tolerance for space-group-finding operations

  • mode (str, optional) –

    Method used for calculating the high-symmetry path. The options are:

    bradcrack

    Use the paths from Bradley and Cracknell. See [brad].

    pymatgen

    Use the paths from pymatgen. See [curt].

    seekpath

    Use the paths from SeeK-path. See [seek].

  • kpt_list (list, optional) –

    List of k-points to use, formatted as a list of subpaths, each containing a list of fractional k-points. For example:

    [ [[0., 0., 0.], [0., 0., 0.5]],
      [[0.5, 0., 0.], [0.5, 0.5, 0.]] ]
    

    Will return points along 0 0 0 -> 0 0 1/2 | 1/2 0 0 -> 1/2 1/2 0

  • labels (list, optional) –

    The k-point labels. These should be provided as a list of str for each subpath of the overall path. For example:

    [ ['Gamma', 'Z'], ['X', 'M'] ]
    

    combined with the above example for kpt_list would indicate the path: Gamma -> Z | X -> M. If no labels are provided, letters from A -> Z will be used instead.

  • eigenvectors (bool, optional) – Write the eigenvectors to the yaml file. (Always True if to_web is set.)

  • dos (str) – Path to Phonopy total dos .dat file

  • height (float, optional) – The height of the plot.

  • width (float, optional) – The width of the plot.

  • ymin (float, optional) – The minimum energy on the y-axis.

  • ymax (float, optional) – The maximum energy on the y-axis.

  • style (list or str, optional) – (List of) matplotlib style specifications, to be composed on top of Sumo base style.

  • no_base_style (bool, optional) – Prevent use of sumo base style. This can make alternative styles behave more predictably.

  • image_format (str, optional) – The image file format. Can be any format supported by matplotlib, including: png, jpg, pdf, and svg. Defaults to pdf.

  • dpi (int, optional) – The dots-per-inch (pixel density) for the image.

  • plt (matplotlib.pyplot, optional) – A matplotlib.pyplot object to use for plotting.

  • fonts (list, optional) – Fonts to use in the plot. Can be a a single font, specified as a str, or several fonts, specified as a list of str.

  • to_json (str, optional) – JSON file to dump the Pymatgen band path object.

  • from_json (list or str, optional) – (List of) JSON bandpath data filename(s) to import and overlay.

  • to_web (str, optional) – JSON file to write data for visualisation with http://henriquemiranda.github.io/phononwebsite

  • legend (list or None, optional) – Legend labels. If None, don’t add a legend. With a list length equal to from_json, label json inputs only. With one extra list entry, label all lines beginning with new plot.

Returns:

A matplotlib pyplot object.

sumo.cli.phonon_bandplot.save_data_files(bs, prefix=None, directory=None)[source]

Write the phonon band structure data files to disk.

Parameters:
  • bs (PhononBandStructureSymmLine) – The phonon band structure.

  • prefix (str, optional) – Prefix for data file.

  • directory (str, optional) – Directory in which to save the data.

Returns:

The filename of the written data file.

Return type:

str