sisl.physics.BandStructure
- class sisl.physics.BandStructure(parent, *args, **kwargs)
Bases:
BrillouinZone
Create a path in the Brillouin zone for plotting band-structures etc.
- Parameters:
parent (
object
orarray_like
) – An object with associated parent.cell and parent.rcell or an array of floats which may be turned into a Latticepoints (
array_like
offloat
) – a list of points that are the corners of the path. Define a discontinuity in the points by adding a None in the list.divisions (
int
orarray_like
ofint
) – number of divisions in each segment. If a single integer is passed it is the total number of points on the path (equally separated). If it is an array_like input it must have length one less than point, in this case the total number of points will besum(divisions) + 1
due to the end-point constraint.names (
array_like
ofstr
) – the associated names of the points on the Brillouin Zone pathjump_dk (
float
orarray_like
, optional) – Percentage ofself.lineark()[-1]
that is used as separation between discontinued jumps in the band-structure. For band-structures with disconnected jumps thelineark
andlineartick
methods returns a separation between the disconnected points according to this percentage. Default value is 5% of the total distance. Alternatively an array equal to the number of discontinuity jumps may be passed for individual percentages. Keyword only, argument.
Examples
>>> lattice = Lattice(10) >>> bs = BandStructure(lattice, [[0] * 3, [0.5] * 3], 200) >>> bs = BandStructure(lattice, [[0] * 3, [0.5] * 3, [1.] * 3], 200) >>> bs = BandStructure(lattice, [[0] * 3, [0.5] * 3, [1.] * 3], 200, ['Gamma', 'M', 'Gamma'])
A disconnected band structure may be created by having None as the element. Note that the number of names does not contain the empty points (they are simply removed). Such a band-structure may be useful when one is interested in a discontinuous band structure.
>>> bs = BandStructure(lattice, [[0, 0, 0], [0, 0.5, 0], None, [0.5, 0, 0], [0.5, 0.5, 0]], 200)
Plotting
Plotting functions for the
BrillouinZone
class.plot.bands
([H, extra_vars, ...])Creates a
BandsData
object and then plots aBandsPlot
from it.plot.fatbands
([H, extra_vars, ...])Creates a
BandsData
object and then plots aFatbandsPlot
from it.plot.sites
([axes, sites, ...])Builds a
SitesPlot
by setting the value of "sites_obj" to the current object.\(k\)-point calculations
Loop over all k-points by applying parent methods for all k.
apply.array
(method[, eta_key])Dispatch the method by one array
apply.average
(method)Dispatch the method by averaging
apply.dataarray
(method)Dispatch the method by returning a DataArray or data-set
apply.grid
(method[, eta_key])Dispatch the method by putting values on the grid
apply.iter
(method[, eta_key])Dispatch the method by iterating values
apply.list
(method)Dispatch the method by returning list of values
apply.ndarray
(method[, eta_key])Dispatch the method by one array
apply.none
(method)Dispatch the method by doing nothing (mostly useful if wrapped)
apply.oplist
(method)Dispatch the method by returning oplist of values
apply.sum
(method)Dispatch the method by summing
apply.xarray
(method)Dispatch the method by returning a DataArray or data-set
Methods
copy
([parent])Create a copy of this object, optionally changing the parent
in_primitive
(k)Move the k-point into the primitive point(s) ]-0.5 ; 0.5]
insert_jump
(*arrays[, value])Return a copy of arrays filled with value at indices of discontinuity jumps
iter
([ret_weight])An iterator for the k-points and (possibly) the weights
lineark
([ticks])A 1D array which corresponds to the delta-k values of the path
The tick-marks corresponding to the linear-k values
merge
(bzs[, weight_scale, parent])Merge several BrillouinZone objects into one
param_circle
(parent, N_or_dk, kR, normal, origin)Create a parameterized k-point list where the k-points are generated on a circle around an origin
parametrize
(parent, func, N, *args, **kwargs)Generate a new
BrillouinZone
object with k-points parameterized via the function func in N separationsset_parent
(parent)Update the parent associated to this object
tocartesian
([k])Transfer a k-point in reduced coordinates to the Cartesian coordinates
tolinear
(k[, ret_index, atol])Convert a k-point into the equivalent linear k-point via the distance
toreduced
(k)Transfer a k-point in Cartesian coordinates to the reduced coordinates
volume
([ret_dim, axes])Calculate the volume of the BrillouinZone, optionally only on some axes axes
write
(sile, *args, **kwargs)Writes k-points to a
tableSile
.Attributes
A list of all k-points (if available)
Weight of the k-points in the
BrillouinZone
object- apply.array(method, eta_key='ndarray')
Dispatch the method by one array
- apply.average(method)
Dispatch the method by averaging
- apply.dataarray(method)
Dispatch the method by returning a DataArray or data-set
- apply.grid(method, eta_key='grid')
Dispatch the method by putting values on the grid
- apply.iter(method, eta_key='iter')
Dispatch the method by iterating values
- apply.list(method)
Dispatch the method by returning list of values
- apply.ndarray(method, eta_key='ndarray')
Dispatch the method by one array
- apply.none(method)
Dispatch the method by doing nothing (mostly useful if wrapped)
- apply.oplist(method)
Dispatch the method by returning oplist of values
- apply.sum(method)
Dispatch the method by summing
- apply.xarray(method)
Dispatch the method by returning a DataArray or data-set
- copy(parent=None) BandStructure
Create a copy of this object, optionally changing the parent
- Parameters:
parent (optional) – change the parent
- static in_primitive(k: numpy.typing.ArrayLike) ndarray
Move the k-point into the primitive point(s) ]-0.5 ; 0.5]
- Parameters:
k (
array_like
) – k-point(s) to move into the primitive cell- Returns:
numpy.ndarray
– all k-points moved into the primitive cell
- insert_jump(*arrays, value=np.nan)[source]
Return a copy of arrays filled with value at indices of discontinuity jumps
Arrays with value in jumps is easier to plot since those lines will be naturally discontinued. For band structures without discontinuity jumps in the Brillouin zone the arrays will be return as is.
It will insert value along the first dimension matching the length of self. For each discontinuity jump an element will be inserted.
This may be useful for plotting since np.nan gets interpreted as a discontinuity in the graph thus removing connections between the segments.
- Parameters:
*arrays (
array_like
) – arrays will get value inserted where there are jumps in the band structurevalue (optional) – the value to be inserted at the jump points in the data array
Examples
Create a bandrstructure with a discontinuity.
>>> gr = geom.graphene() >>> bs = BandStructure(gr, [[0, 0, 0], [0.5, 0, 0], None, [0, 0, 0], [0, 0.5, 0]], 4) >>> data = np.zeros([len(bs), 10]) >>> data_with_jump = bs.insert_jump(data) >>> assert data_with_jump.shape == (len(bs)+1, 10) >>> np.all(data_with_jump[2] == np.nan) True
- lineark(ticks: bool = False)[source]
A 1D array which corresponds to the delta-k values of the path
This is mainly meant for plotting but may be useful for finding out distances in the reciprocal lattice.
Examples
>>> p = BandStructure(...) >>> eigs = Hamiltonian.eigh(p) >>> for i in range(len(Hamiltonian)): ... plt.plot(p.lineark(), eigs[:, i])
>>> p = BandStructure(...) >>> eigs = Hamiltonian.eigh(p) >>> lk, kt, kl = p.lineark(True) >>> plt.xticks(kt, kl) >>> for i in range(len(Hamiltonian)): ... plt.plot(lk, eigs[:, i])
- Parameters:
ticks – if True the ticks for the points are also returned
See also
linspace_bz
converts k-points into a linear distance parameterization
- Returns:
linear_k (
numpy.ndarray
) – the positions in reciprocal space determined by the distance between pointsticks (
numpy.ndarray
) – linear k-positions of the points, only returned if ticks isTrue
ticklabels (
list
ofstr
) – labels at ticks, only returned if ticks isTrue
- lineartick()[source]
The tick-marks corresponding to the linear-k values
- Returns:
numpy.ndarray
– the positions in reciprocal space determined by the distance between points
See also
lineark
Routine used to calculate the tick-marks.
- static merge(bzs, weight_scale: Sequence[float] | float = 1.0, parent=None)
Merge several BrillouinZone objects into one
The merging strategy only stores the new list of k-points and weights. Information retained in the merged objects will not be stored.
- Parameters:
bzs (
list-like
ofBrillouinZone objects
) – each element is a BrillouinZone object withbzs[i].k
andbzs[i].weight
fields.weight_scale (
list-like
orfloat
) – these are matched item-wise with bzs and applied to. Internallyitertools.zip_longest(fillvalue=weight_scale[-1])
will be used to extend for all bzs.parent (
object
, optional) – Associated parent in the returned object, will default tobzs[0].parent
- Returns:
BrillouinZone
– even if all objects are not BrillouinZone objects the returned object will be.
- classmethod param_circle(parent, N_or_dk: int | float, kR: float, normal, origin, loop: bool = False)
Create a parameterized k-point list where the k-points are generated on a circle around an origin
The generated circle is a perfect circle in the reciprocal space (Cartesian coordinates). To generate a perfect circle in units of the reciprocal lattice vectors one can generate the circle for a diagonal supercell with side-length \(2\pi\), see example below.
- Parameters:
parent (
Lattice
, orLatticeChild
) – the parent objectN_or_dk (
int
) – number of k-points generated using the parameterization (if an integer), otherwise it specifies the discretization length on the circle (in 1/Ang), If the latter case will use less than 2 points a warning will be raised and the number of points increased to 2.kR (
float
) – radius of the k-point. In 1/Angnormal (
array_like
offloat
) – normal vector to determine the circle planeorigin (
array_like
offloat
) – origin of the circle used to generate the circular parameterizationloop – whether the first and last point are equal
Examples
>>> lattice = Lattice([1, 1, 10, 90, 90, 60]) >>> bz = BrillouinZone.param_circle(lattice, 10, 0.05, [0, 0, 1], [1./3, 2./3, 0])
To generate a circular set of k-points in reduced coordinates (reciprocal
>>> lattice = Lattice([1, 1, 10, 90, 90, 60]) >>> bz = BrillouinZone.param_circle(lattice, 10, 0.05, [0, 0, 1], [1./3, 2./3, 0]) >>> bz_rec = BrillouinZone.param_circle(2*np.pi, 10, 0.05, [0, 0, 1], [1./3, 2./3, 0]) >>> bz.k[:, :] = bz_rec.k[:, :]
- Returns:
BrillouinZone
– with the parameterized k-points.
- static parametrize(parent, func, N: Sequence[int] | int, *args, **kwargs) BrillouinZone
Generate a new
BrillouinZone
object with k-points parameterized via the function func in N separationsGenerator of a parameterized Brillouin zone object that contains a parameterized k-point list.
- Parameters:
parent (
Lattice
, orLatticeChild
) – the object that the returned object will contain as parentfunc (
callable
) – method that parameterizes the k-points, must at least accept three arguments, 1.parent
: object 2.N
: total number of k-points 3.i
: current index of the k-point (starting from 0)the function must return a k-point in 3 dimensions.
N (
int
orlist
ofint
) – number of k-points generated using the parameterization, or a list of integers that will be looped over. In this case argumentsN
andi
in func will be lists accordingly.*args – additional arguments passed directly to func
**kwargs – additional keyword arguments passed directly to func
Examples
Simple linear k-points
>>> def func(sc, N, i): ... return [i/N, 0, 0] >>> bz = BrillouinZone.parametrize(1, func, 10) >>> assert len(bz) == 10 >>> assert np.allclose(bz.k[-1, :], [9./10, 0, 0])
For double looping, say to create your own grid
>>> def func(sc, N, i): ... return [i[0]/N[0], i[1]/N[1], 0] >>> bz = BrillouinZone.parametrize(1, func, [10, 5]) >>> assert len(bz) == 50
- plot.bands(H: Hamiltonian | None = None, extra_vars: Sequence[dict | str] = (), *, Erange: tuple[float, float] | None = None, E0: float = 0.0, E_axis: Literal['x', 'y'] = 'y', bands_range: tuple[int, int] | None = None, spin: Literal[0, 1] | None = None, bands_style: StyleSpec = {'color': 'black', 'dash': 'solid', 'opacity': 1, 'width': 1}, spindown_style: StyleSpec = {'color': 'blue', 'width': 1}, colorscale: Colorscale | None = None, gap: bool = False, gap_tol: float = 0.01, gap_color: str = 'red', gap_marker: dict = {'size': 7}, direct_gaps_only: bool = False, custom_gaps: Sequence[dict] = [], line_mode: Literal['line', 'scatter', 'area_line'] = 'line', group_legend: bool = True, backend: str = 'plotly') BandsPlot
Creates a
BandsData
object and then plots aBandsPlot
from it.- Parameters:
H – The Hamiltonian to use for the calculations. If None, the parent of the Brillouin zone will be used, which is typically what you want!
extra_vars – Additional variables to calculate for each eigenstate, apart from their energy. Each item of the list should be a dictionary with the following keys:
‘name’, str: The name of the variable.
‘getter’, callable: A function that gets 3 arguments: eigenstate, plot and
spin index, and returns the values of the variable in a numpy array. This function will be called for each eigenstate object separately. That is, once for each (k-point, spin) combination. - ‘coords’, tuple of str: The names of the dimensions of the returned array. The number of coordinates should match the number of dimensions. - ‘coords_values’, dict: If this variable introduces a new coordinate, you should pass the values for that coordinate here. If the coordinates were already defined by another variable, they will already have values. If you are unsure that the coordinates are new, just pass the values for them, they will get overwritten.
Each item can also be a string indicating the name of a known variable: ‘norm2’, ‘spin_moment’, ‘ipr’.
Erange – The energy range to plot. If None, the range is determined by
bands_range
.E0 – The energy reference.
E_axis – Axis to plot the energies.
bands_range – The bands to plot. Only used if
Erange
is None. If None, the 15 bands above and below the Fermi level are plotted.spin – Which spin channel to display. Only meaningful for spin-polarized calculations. If None and the calculation is spin polarized, both are plotted.
bands_style – Styling attributes for bands.
spindown_style – Styling attributes for the spin down bands (if present). Any missing attribute will be taken from
bands_style
.colorscale – Colorscale to use for the bands in case the color attribute is an array of values. If None, the default colorscale is used for each backend.
gap – Whether to display the gap.
gap_tol – Tolerance in k for determining whether two gaps are the same.
gap_color – Color of the gap.
gap_marker – Marker styles for the gap (as plotly marker’s styles).
direct_gaps_only – Whether to only display direct gaps.
custom_gaps – List of custom gaps to display. See the showcase notebooks for examples.
line_mode – The method used to draw the band lines.
group_legend – Whether to group all bands in the legend to show a single legend item.
If the bands are spin polarized, bands are grouped by spin channel.
backend – The backend to use to generate the figure.
See also
BandsPlot
The plot class used to generate the plot.
BandsData
The class to which data is converted.
- plot.fatbands(H: Hamiltonian | None = None, extra_vars: Sequence[dict | str] = (), *, Erange: tuple[float, float] | None = None, E0: float = 0.0, E_axis: Literal['x', 'y'] = 'y', bands_range: tuple[int, int] | None = None, spin: Literal[0, 1] | None = None, bands_style: StyleSpec = {'color': 'black', 'opacity': 1, 'width': 1}, spindown_style: StyleSpec = {'color': 'blue', 'width': 1}, gap: bool = False, gap_tol: float = 0.01, gap_color: str = 'red', gap_marker: dict = {'size': 7}, direct_gaps_only: bool = False, custom_gaps: Sequence[dict] = [], bands_mode: Literal['line', 'scatter', 'area_line'] = 'line', bands_group_legend: bool = True, groups: OrbitalQueries = [], fatbands_var: str = 'norm2', fatbands_mode: Literal['line', 'scatter', 'area_line'] = 'area_line', fatbands_scale: float = 1.0, backend: str = 'plotly') FatbandsPlot
Creates a
BandsData
object and then plots aFatbandsPlot
from it.- Parameters:
H – The Hamiltonian to use for the calculations. If None, the parent of the Brillouin zone will be used, which is typically what you want!
extra_vars – Additional variables to calculate for each eigenstate, apart from their energy. Each item of the list should be a dictionary with the following keys:
‘name’, str: The name of the variable.
‘getter’, callable: A function that gets 3 arguments: eigenstate, plot and
spin index, and returns the values of the variable in a numpy array. This function will be called for each eigenstate object separately. That is, once for each (k-point, spin) combination. - ‘coords’, tuple of str: The names of the dimensions of the returned array. The number of coordinates should match the number of dimensions. - ‘coords_values’, dict: If this variable introduces a new coordinate, you should pass the values for that coordinate here. If the coordinates were already defined by another variable, they will already have values. If you are unsure that the coordinates are new, just pass the values for them, they will get overwritten.
Each item can also be a string indicating the name of a known variable: ‘norm2’, ‘spin_moment’, ‘ipr’.
Erange – The energy range to plot. If None, the range is determined by
bands_range
.E0 – The energy reference.
E_axis – Axis to plot the energies.
bands_range – The bands to plot. Only used if
Erange
is None. If None, the 15 bands above and below the Fermi level are plotted.spin – Which spin channel to display. Only meaningful for spin-polarized calculations. If None and the calculation is spin polarized, both are plotted.
bands_style – Styling attributes for bands.
spindown_style – Styling attributes for the spin down bands (if present). Any missing attribute will be taken from
bands_style
.gap – Whether to display the gap.
gap_tol – Tolerance in k for determining whether two gaps are the same.
gap_color – Color of the gap.
gap_marker – Marker styles for the gap (as plotly marker’s styles).
direct_gaps_only – Whether to only display direct gaps.
custom_gaps – List of custom gaps to display. See the showcase notebooks for examples.
bands_mode – The method used to draw the band lines.
bands_group_legend – Whether to group all bands in the legend to show a single legend item.
If the bands are spin polarized, bands are grouped by spin channel.
groups – Orbital groups to plots. See showcase notebook for examples.
fatbands_var – The variable to use from bands_data to determine the width of the fatbands. This variable must have as coordinates (k, band, orb, [spin]).
fatbands_mode – The method used to draw the fatbands.
fatbands_scale – Factor that scales the size of all fatbands.
backend – The backend to use to generate the figure.
See also
FatbandsPlot
The plot class used to generate the plot.
BandsData
The class to which data is converted.
- plot.sites(axes: Axes = ['x', 'y', 'z'], sites: AtomsIndex = None, sites_style: Sequence[AtomsStyleSpec] = [], sites_scale: float = 1.0, sites_name: str = 'Sites', sites_colorscale: Colorscale | None = None, drawing_mode: Literal['scatter', 'balls', 'line', None] = None, show_cell: Literal['box', 'axes', False] = False, cell_style: StyleSpec = {}, nsc: tuple[int, int, int] = (1, 1, 1), sites_ndim_scale: tuple[float, float, float] = (1, 1, 1), dataaxis_1d: np.ndarray | Callable | None = None, arrows: Sequence[AtomArrowSpec] = (), backend='plotly') SitesPlot
Builds a
SitesPlot
by setting the value of “sites_obj” to the current object.Plots sites from an object that can be parsed into a geometry.
The only differences between this plot and a geometry plot is the naming of the inputs and the fact that there are no options to plot bonds.
- Parameters:
axes – The axes to project the sites to.
sites – The sites to plot. If None, all sites are plotted.
sites_style – List of style specifications for the sites. See the showcase notebooks for examples.
sites_scale – Scaling factor for the size of all sites.
sites_name – Name to give to the trace that draws the sites.
sites_colorscale – Colorscale to use for the sites in case the color attribute is an array of values. If None, the default colorscale is used for each backend.
drawing_mode – The method used to draw the sites.
show_cell – Mode to display the reciprocal cell. If False, the cell is not displayed.
cell_style – Style specification for the reciprocal cell. See the showcase notebooks for examples.
nsc – Number of unit cells to display in each direction.
sites_ndim_scale – Scaling factor for the size of the sites for different dimensionalities (1D, 2D, 3D).
dataaxis_1d – Only meaningful for 1D plots. The data to plot on the Y axis.
arrows – List of arrow specifications to display. See the showcase notebooks for examples.
backend – The backend to use to generate the figure.
- set_parent(parent) None
Update the parent associated to this object
- Parameters:
parent (
object
orarray_like
) – an object containing cell vectors
- tocartesian(k: npt.ArrayLike | None = None) np.ndarray
Transfer a k-point in reduced coordinates to the Cartesian coordinates
- Parameters:
k – k-point in reduced coordinates, defaults to this objects k-points.
- Returns:
numpy.ndarray
– in units of 1/Ang
- tolinear(k, ret_index: bool = False, atol: float = 1e-4)[source]
Convert a k-point into the equivalent linear k-point via the distance
Finds the index of the k-point in self.k that is closests to
k
. The returned value is then the equivalent index inlineark
.This is very useful for extracting certain points along the band structure.
- Parameters:
k (
array_like
) – the k-point(s) to locate in the linear valuesret_index – whether the indices are also returned
atol – when the found k-point has a distance (in Cartesian coordinates) is differing by more than tol a warning will be issued. The tolerance is in units 1/Ang.
- toreduced(k: numpy.typing.ArrayLike) ndarray
Transfer a k-point in Cartesian coordinates to the reduced coordinates
- Parameters:
- Returns:
numpy.ndarray
– in units of reciprocal lattice vectors ]-0.5 ; 0.5] (if k is in the primitive cell)
- volume(ret_dim: bool = False, axes: CellAxes | None = None) float | tuple[float, int]
Calculate the volume of the BrillouinZone, optionally only on some axes axes
This will return the volume of the Brillouin zone, depending on the dimensions of the system. Here the dimensions of the system is determined by how many dimensions have auxilliary supercells that can contribute to Brillouin zone integrals. Therefore the returned value will have differing units depending on dimensionality.
- Parameters:
ret_dim – also return the dimensionality of the system
axes – estimate the volume using only the directions indexed by this array. The default axes are only the periodic ones (
self.parent.pbc.nonzero()[0]
). Hence the units might not necessarily be 1/Ang^3.
- Returns:
vol
– the volume of the Brillouin zone. Units are 1/Ang^D with D being the dimensionality. For 0D it will return 0.dimensionality (
int
) – the dimensionality of the volume
- write(sile: sisl.typing.SileLike, *args, **kwargs) None
Writes k-points to a
tableSile
.This allows one to pass a tableSile or a file-name.
- apply
Loop over all k-points by applying parent methods for all k.
This allows potential for running and collecting various computationally heavy methods from a single point on all k-points.
The
apply
method will dispatch the parent methods through all k-points and passingk
as arguments to the parent methods in a straight-forward manner.For instance to iterate over all eigenvalues of a Hamiltonian
>>> H = Hamiltonian(...) >>> bz = BrillouinZone(H) >>> for ik, eigh in enumerate(bz.apply.eigh()): ... # do something with eigh which corresponds to bz.k[ik]
By default the
apply
method exposes a set of dispatch methods:apply.iter
, the default iterator moduleapply.average
reduced result by averaging (usingBrillouinZone.weight
as the weight per k-point.apply.sum
reduced result without weighingapply.array
return a single array with all values; has len equal to number of k-pointsapply.none
, specialized method that is mainly useful when wrapping methodsapply.list
same asapply.array
but using Python list as return valueapply.oplist
usingsisl.oplist
allows greater flexibility for mathematical operations element wiseapply.datarray
ifxarray
is available one can retrieve anxarray.DataArray
instance
Please see Brillouin zone for further examples.
- plot
Plotting functions for the
BrillouinZone
class.
- property weight: ndarray
Weight of the k-points in the
BrillouinZone
object