sisl.io.siesta.fdfSileSiesta
- class sisl.io.siesta.fdfSileSiesta(filename, *args, **kwargs)
Bases:
SileSiesta
FDF-input file
By supplying base you can reference files in other directories. By default the
base
is the directory given in the file name.- Parameters:
Examples
>>> fdf = fdfSileSiesta("tmp/RUN.fdf") # reads output files in "tmp/" folder >>> fdf = fdfSileSiesta("tmp/RUN.fdf", base=".") # reads output files in "./" folder
Plotting
Plotting functions for the
fdfSileSiesta
class.plot.bands
([bands_file, ...])Creates a
BandsData
object and then plots aBandsPlot
from it.plot.geometry
([output, ...])Calls
read_geometry
and creates aGeometryPlot
from its output.plot.grid
(name, *args[, ...])Calls
read_grid
and creates aGridPlot
from its output.plot.pdos
([source, ...])Creates a
PDOSData
object and then plots aPdosPlot
from it.plot.wavefunction
([source, k, ...])Creates a
EigenstateData
object and then plots aWavefunctionPlot
from it.Methods
base_directory
([relative_to])Retrieve the base directory of the file, relative to the path relative_to
close
()dir_file
([filename, filename_base])File of the current Sile
get
(label[, default, unit, with_unit])Retrieve fdf-keyword from the file
includes
()Return a list of all files that are included or otherwise necessary for reading the fdf file
print
(key, value)Return a string which is pretty-printing the key+value
read
(*args, **kwargs)Generic read method which should be overloaded in child-classes
read_basis
(*args, **kwargs)Read the atomic species and figure out the number of atomic orbitals in their basis
read_density_matrix
(*args, **kwargs)Try and read density matrix by reading the <>.nc, <>.TSDE files, <>.DM (in that order)
read_dynamical_matrix
(*args, **kwargs)Read dynamical matrix from output of the calculation
read_energy_density_matrix
(*args, **kwargs)Try and read energy density matrix by reading the <>.nc or <>.TSDE files (in that order)
read_fermi_level
(*args, **kwargs)Read Fermi-level from output of the calculation
read_force
(*args, **kwargs)Read forces from the output of the calculation (forces are not defined in the input)
read_force_constant
(*args, **kwargs)Read Hessian/force constant from the output of the calculation
read_geometry
([output])Returns Geometry object by reading fdf or Siesta output related files.
read_grid
(name, *args, **kwargs)Read grid related information from any of the output files
read_hamiltonian
(*args, **kwargs)Try and read the Hamiltonian by reading the <>.nc, <>.TSHS files, <>.HSX (in that order)
read_hessian
(*args, **kwargs)Read Hessian/force constant from the output of the calculation
read_lattice
([output])Returns Lattice object by reading fdf or Siesta output related files.
read_lattice_nsc
(*args, **kwargs)Read supercell size using any method available
read_overlap
(*args, **kwargs)Try and read the overlap matrix by reading the <>.nc, <>.TSHS files, <>.HSX, <>.onlyS (in that order)
set
(key, value[, keep])Add the key and value to the FDF file
type
(label)Return the type of the fdf-keyword
write
(*args, **kwargs)Generic write method which should be overloaded in child-classes
Writes Brillouinzone information to the fdf input file
write_geometry
(geometry[, fmt])Writes the geometry
write_lattice
(lattice[, fmt])Writes the supercell
Attributes
File of the current Sile
File of the current Sile
- base_directory(relative_to='.')
Retrieve the base directory of the file, relative to the path relative_to
- close()
- dir_file(filename=None, filename_base='')
File of the current Sile
- get(label: str, default: Any | None = None, unit: str | None = None, with_unit: bool = False)[source]
Retrieve fdf-keyword from the file
- Parameters:
label (
str
) – the fdf-label to search fordefault (optional) – if the label is not found, this will be the returned value (default to
None
)unit (
str
, optional) – unit of the physical quantity to returnwith_unit (
bool
, optional) – whether the physical quantity gets returned with the found unit in the fdf file.
- Returns:
value (
the value
ofthe fdf-label. If the label is a block
,a `list
is returned`,for
) – a real value a float (or if the default is of float), for an integer, an int is returned.unit (
if `with_unit
is true this will contain the associated unit if it is specified`)
Examples
>>> print(open(...).readlines()) LabeleV 1. eV LabelRy 1. Ry Label name FakeInt 1 %block Hello line 1 line2 %endblock >>> fdf.get("LabeleV") == 1. # default unit is eV >>> fdf.get("LabelRy") == unit.siesta.unit_convert("Ry", "eV") >>> fdf.get("LabelRy", unit="Ry") == 1. >>> fdf.get("LabelRy", with_unit=True) == (1., "Ry") >>> fdf.get("FakeInt", "0") == "1" >>> fdf.get("LabeleV", with_unit=True) == (1., "eV") >>> fdf.get("Label", with_unit=True) == "name" # no unit present on line >>> fdf.get("Hello") == ["line 1", "line2"]
- includes()[source]
Return a list of all files that are included or otherwise necessary for reading the fdf file
- plot.bands(bands_file: str | bandsSileSiesta | None = None, *, 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:
bands_file – Path to the bands file. If None, it will be assumed that the bands file is in the same directory as the fdf file and has the name <SystemLabel>.bands, with SystemLabel being retrieved from the fdf file.
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.geometry(output: bool = False, *args, data_kwargs={}, axes: Axes = ['x', 'y', 'z'], atoms: AtomsIndex = None, atoms_style: Sequence[AtomsStyleSpec] = [], atoms_scale: float = 1.0, atoms_colorscale: Colorscale | None = None, drawing_mode: Literal['scatter', 'balls', None] = None, bind_bonds_to_ats: bool = True, points_per_bond: int = 20, bonds_style: StyleSpec = {}, bonds_scale: float = 1.0, bonds_colorscale: Colorscale | None = None, show_atoms: bool = True, show_bonds: bool = True, show_cell: Literal['box', 'axes', False] = 'box', cell_style: StyleSpec = {}, nsc: tuple[int, int, int] = (1, 1, 1), atoms_ndim_scale: tuple[float, float, float] = (16, 16, 1), bonds_ndim_scale: tuple[float, float, float] = (1, 1, 10), dataaxis_1d: np.ndarray | Callable | None = None, arrows: Sequence[AtomArrowSpec] = (), backend='plotly') GeometryPlot
Calls
read_geometry
and creates aGeometryPlot
from its output.- Parameters:
output – whether to read geometry from output files (default to read from the fdf file).
order (
list
ofstr
, optional) – the order of which to try and read the geometry. By default this is[XV, nc, TSHS, STRUCT, HSX, fdf]
if output is true If order is present output is disregarded.axes – The axes to project the geometry to.
atoms – The atoms to plot. If None, all atoms are plotted.
atoms_style – List of style specifications for the atoms. See the showcase notebooks for examples.
atoms_scale – Scaling factor for the size of all atoms.
atoms_colorscale – Colorscale to use for the atoms 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 atoms.
bind_bonds_to_ats – Whether to display only bonds between atoms that are being displayed.
points_per_bond – When the points are drawn using points instead of lines (e.g. in some frameworks to draw multicolor bonds), the number of points used per bond.
bonds_style – Style specification for the bonds. See the showcase notebooks for examples.
bonds_scale – Scaling factor for the width of all bonds.
bonds_colorscale – Colorscale to use for the bonds in case the color attribute is an array of values. If None, the default colorscale is used for each backend.
show_atoms – Whether to display the atoms.
show_bonds – Whether to display the bonds.
show_cell – Mode to display the cell. If False, the cell is not displayed.
cell_style – Style specification for the cell. See the showcase notebooks for examples.
nsc – Number of unit cells to display in each direction.
atoms_ndim_scale – Scaling factor for the size of the atoms for different dimensionalities (1D, 2D, 3D).
bonds_ndim_scale – Scaling factor for the width of the bonds 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.
See also
GeometryPlot
The plot class used to generate the plot.
read_geometry
The method called to get the data.
- plot.grid(name: str, *args, data_kwargs={}, axes: Axes = ['z'], represent: Literal['real', 'imag', 'mod', 'phase', 'deg_phase', 'rad_phase'] = 'real', transforms: Sequence[str | Callable] = (), reduce_method: Literal['average', 'sum'] = 'average', boundary_mode: str = 'grid-wrap', nsc: tuple[int, int, int] = (1, 1, 1), interp: tuple[int, int, int] = (1, 1, 1), isos: Sequence[dict] = [], smooth: bool = False, colorscale: Colorscale | None = None, crange: tuple[float, float] | None = None, cmid: float | None = None, show_cell: Literal['box', 'axes', False] = 'box', cell_style: dict = {}, x_range: Sequence[float] | None = None, y_range: Sequence[float] | None = None, z_range: Sequence[float] | None = None, plot_geom: bool = False, geom_kwargs: dict = {}, backend: str = 'plotly') GridPlot
Calls
read_grid
and creates aGridPlot
from its output.- Parameters:
name – name of data to read. The list of names correspond to the Siesta output manual (Rho, TotalPotential, etc.), the strings are case insensitive.
order (
list
ofstr
, optional) – the order of which to try and read the geometry. By default this is["nc", "grid.nc", "bin"]
(bin refers to the binary files)axes – The axes to project the grid to.
represent – The representation of the grid to plot.
transforms – List of transforms to apply to the grid before plotting.
reduce_method – The method used to reduce the grid axes that are not displayed.
boundary_mode – The method used to deal with the boundary conditions. Only used if the grid is to be orthogonalized. See scipy docs for more info on the possible values.
nsc – The number of unit cells to display in each direction.
interp – The interpolation factor to use for each axis to make the grid smoother.
isos – List of isosurfaces or isocontours to plot. See the showcase notebooks for examples.
smooth – Whether to ask the plotting backend to make an attempt at smoothing the grid display.
colorscale – Colorscale to use for the grid display in the 2D representation. If None, the default colorscale is used for each backend.
crange – Min and max values for the colorscale.
cmid – The value at which the colorscale is centered.
show_cell – Method used to display the unit cell. If False, the cell is not displayed.
cell_style – Style specification for the cell. See the showcase notebooks for examples.
x_range – The range of the x axis to take into account. Even if the X axis is not displayed! This is important because the reducing operation will only be applied on this range.
y_range – The range of the y axis to take into account. Even if the Y axis is not displayed! This is important because the reducing operation will only be applied on this range.
z_range – The range of the z axis to take into account. Even if the Z axis is not displayed! This is important because the reducing operation will only be applied on this range.
plot_geom – Whether to plot the associated geometry (if any).
geom_kwargs – Keyword arguments to pass to the geometry plot of the associated geometry.
backend – The backend to use to generate the figure.
See also
GridPlot
The plot class used to generate the plot.
read_grid
The method called to get the data.
- plot.pdos(source: Literal['pdos', 'tbtnc', 'wfsx', 'hamiltonian'] = 'pdos', *, data_kwargs={}, groups: Sequence[OrbitalStyleQuery] = [{'name': 'DOS'}], Erange: tuple[float, float] = (-2, 2), E_axis: Literal['x', 'y'] = 'x', line_mode: Literal['line', 'scatter', 'area_line'] = 'line', line_scale: float = 1.0, backend: str = 'plotly') PdosPlot
Creates a
PDOSData
object and then plots aPdosPlot
from it.- Parameters:
source – The source to read the PDOS data from.
**kwargs – Extra arguments to be passed to the PDOSData constructor, which depends on the source requested.
Except for the hamiltonian source, no extra arguments are needed (and they won’t be used). See PDOSData.from_hamiltonian for the extra arguments accepted by the hamiltonian data constructor.
groups – List of orbital specifications to filter and accumulate the PDOS. The contribution of each group will be displayed in a different line. See showcase notebook for examples.
Erange – The energy range to plot.
E_axis – Axis to project the energies.
line_mode – Mode used to draw the PDOS lines.
line_scale – Scaling factor for the width of all lines.
backend – The backend to generate the figure.
See also
PdosPlot
The plot class used to generate the plot.
PDOSData
The class to which data is converted.
- plot.wavefunction(source: Literal['wfsx', 'hamiltonian'] = 'wfsx', k: tuple[float, float, float] = (0, 0, 0), spin: int = 0, *, i: int = 0, geometry: Geometry | None = None, grid_prec: float = 0.2, grid: Grid | None = None, axes: Axes = ['z'], represent: Literal['real', 'imag', 'mod', 'phase', 'deg_phase', 'rad_phase'] = 'real', transforms: Sequence[str | Callable] = (), reduce_method: Literal['average', 'sum'] = 'average', boundary_mode: str = 'grid-wrap', nsc: tuple[int, int, int] = (1, 1, 1), interp: tuple[int, int, int] = (1, 1, 1), isos: Sequence[dict] = [], smooth: bool = False, colorscale: Colorscale | None = None, crange: tuple[float, float] | None = None, cmid: float | None = None, show_cell: Literal['box', 'axes', False] = 'box', cell_style: dict = {}, x_range: Sequence[float] | None = None, y_range: Sequence[float] | None = None, z_range: Sequence[float] | None = None, plot_geom: bool = False, geom_kwargs: dict = {}, backend: str = 'plotly') WavefunctionPlot
Creates a
EigenstateData
object and then plots aWavefunctionPlot
from it.- Parameters:
source – The source of the eigenstates. It can be either from the WFSX file or from the Hamiltonian, in which case the eigenstates need to be calculated.
k – The k-point for which the eigenstate coefficients are desired.
spin – The spin for which the eigenstate coefficients are desired.
i – The index of the eigenstate to plot.
geometry – Geometry to use to project the eigenstate to real space. If None, the geometry associated with the eigenstate is used.
grid_prec – The precision of the grid where the wavefunction is projected.
grid – The grid to plot.
axes – The axes to project the grid to.
represent – The representation of the grid to plot.
transforms – List of transforms to apply to the grid before plotting.
reduce_method – The method used to reduce the grid axes that are not displayed.
boundary_mode – The method used to deal with the boundary conditions. Only used if the grid is to be orthogonalized. See scipy docs for more info on the possible values.
nsc – The number of unit cells to display in each direction.
interp – The interpolation factor to use for each axis to make the grid smoother.
isos – List of isosurfaces or isocontours to plot. See the showcase notebooks for examples.
smooth – Whether to ask the plotting backend to make an attempt at smoothing the grid display.
colorscale – Colorscale to use for the grid display in the 2D representation. If None, the default colorscale is used for each backend.
crange – Min and max values for the colorscale.
cmid – The value at which the colorscale is centered.
show_cell – Method used to display the unit cell. If False, the cell is not displayed.
cell_style – Style specification for the cell. See the showcase notebooks for examples.
x_range – The range of the x axis to take into account. Even if the X axis is not displayed! This is important because the reducing operation will only be applied on this range.
y_range – The range of the y axis to take into account. Even if the Y axis is not displayed! This is important because the reducing operation will only be applied on this range.
z_range – The range of the z axis to take into account. Even if the Z axis is not displayed! This is important because the reducing operation will only be applied on this range.
plot_geom – Whether to plot the associated geometry (if any).
geom_kwargs – Keyword arguments to pass to the geometry plot of the associated geometry.
backend – The backend to use to generate the figure.
See also
WavefunctionPlot
The plot class used to generate the plot.
EigenstateData
The class to which data is converted.
- read(*args, **kwargs)
Generic read method which should be overloaded in child-classes
- Parameters:
kwargs – keyword arguments will try and search for the attribute
read_<>
and call it with the remaining**kwargs
as arguments.
- read_basis(*args, **kwargs)[source]
Read the atomic species and figure out the number of atomic orbitals in their basis
The order of the read is shown below.
One can limit the tried files to only one file by passing only a single file ending.
- read_density_matrix(*args, **kwargs) DensityMatrix [source]
Try and read density matrix by reading the <>.nc, <>.TSDE files, <>.DM (in that order)
One can limit the tried files to only one file by passing only a single file ending.
- read_dynamical_matrix(*args, **kwargs)[source]
Read dynamical matrix from output of the calculation
Generally the mass is stored in the basis information output, but for dynamical matrices it makes sense to let the user control this, e.g. through the fdf file. By default the mass will be read from the AtomicMass key in the fdf file and not from the basis set information.
- Parameters:
order (
list
ofstr
, optional) – the order of which to try and read the dynamical matrix. By default this is["nc", "FC"]
.cutoff_dist (
float
, optional) – cutoff value for the distance of the force-constants (everything farther than cutoff_dist will be set to 0 Ang). Default, no cutoff.cutoff (
float
, optional) – absolute values below the cutoff are considered 0. Defaults to 0. eV/Ang**2.trans_inv (
bool
, optional) – if true (default), the force-constant matrix will be fixed so that translational invariance will be enforcedsum0 (
bool
, optional) – if true (default), the sum of forces on atoms for each displacement will be forced to 0.hermitian (
bool
, optional) – if true (default), the returned dynamical matrix will be hermitian
Notes
This is highly untested and may result in errorneous matrices. Please report back to developers about problems and suggestions.
- Returns:
dynamic_matrix (
DynamicalMatrix
) – the dynamical matrix
- read_energy_density_matrix(*args, **kwargs) EnergyDensityMatrix [source]
Try and read energy density matrix by reading the <>.nc or <>.TSDE files (in that order)
One can limit the tried files to only one file by passing only a single file ending.
- read_force(*args, **kwargs) ndarray [source]
Read forces from the output of the calculation (forces are not defined in the input)
- read_force_constant(*args, **kwargs)
Read Hessian/force constant from the output of the calculation
- Returns:
M (
numpy.ndarray
) – vector[*, 3, 2, *, 3]
with Hessian/force constant element for each of the atomic displacements
- read_geometry(output: bool = False, *args, **kwargs) Geometry [source]
Returns Geometry object by reading fdf or Siesta output related files.
One can limit the tried files to only one file by passing only a single file ending.
- Parameters:
Examples
>>> fdf = get_sile("RUN.fdf") >>> fdf.read_geometry() # read from fdf >>> fdf.read_geometry(True) # read from [XV, nc, fdf] >>> fdf.read_geometry(order=["nc"]) # read from [nc] >>> fdf.read_geometry(True, order=["nc"]) # read from [nc]
- read_grid(name: str, *args, **kwargs) Grid [source]
Read grid related information from any of the output files
The order of the readed data is shown below.
One can limit the tried files to only one file by passing only a single file ending.
- Parameters:
name – name of data to read. The list of names correspond to the Siesta output manual (Rho, TotalPotential, etc.), the strings are case insensitive.
order (
list
ofstr
, optional) – the order of which to try and read the geometry. By default this is["nc", "grid.nc", "bin"]
(bin refers to the binary files)
- read_hamiltonian(*args, **kwargs) Hamiltonian [source]
Try and read the Hamiltonian by reading the <>.nc, <>.TSHS files, <>.HSX (in that order)
One can limit the tried files to only one file by passing only a single file ending.
- read_hessian(*args, **kwargs)[source]
Read Hessian/force constant from the output of the calculation
- Returns:
M (
numpy.ndarray
) – vector[*, 3, 2, *, 3]
with Hessian/force constant element for each of the atomic displacements
- read_lattice(output: bool = False, *args, **kwargs) Lattice [source]
Returns Lattice object by reading fdf or Siesta output related files.
One can limit the tried files to only one file by passing only a single file ending.
- Parameters:
Examples
>>> fdf = get_sile("RUN.fdf") >>> fdf.read_lattice() # read from fdf >>> fdf.read_lattice(True) # read from [XV, nc, fdf] >>> fdf.read_lattice(order=["nc"]) # read from [nc] >>> fdf.read_lattice(True, order=["nc"]) # read from [nc]
- read_lattice_nsc(*args, **kwargs)[source]
Read supercell size using any method available
- Raises:
SislWarning – if none of the files can be read
- read_overlap(*args, **kwargs) Overlap [source]
Try and read the overlap matrix by reading the <>.nc, <>.TSHS files, <>.HSX, <>.onlyS (in that order)
One can limit the tried files to only one file by passing only a single file ending.
- set(key: str, value: Any, keep: bool = True)[source]
Add the key and value to the FDF file
- Parameters:
key (
str
) – the fdf-key value to be set in the fdf filevalue (
str
orlist
ofstr
) – the value of the string. If a str is passed a regular fdf-key is used, if a list it will be a %block.keep (
bool
, optional) – whether old flags will be kept in the fdf file. In this case a time-stamp will be written to show when the key was overwritten.
- type(label: str)[source]
Return the type of the fdf-keyword
- Parameters:
label (
str
) – the label to look-up
- write(*args, **kwargs)
Generic write method which should be overloaded in child-classes
- Parameters:
**kwargs – keyword arguments will try and search for the attribute write_ and call it with the remaining
**kwargs
as arguments.
- write_brillouinzone(bz: BrillouinZone)[source]
Writes Brillouinzone information to the fdf input file
The bz object will be written as options in the input file. The class of bz decides which options gets written. For instance a BandStructure class will write the corresponding
BandLines
block.- Parameters:
bz – which setting to write to the file
Notes
Currently, only the BandStructure class may be accepted as bz.
- write_geometry(geometry: Geometry, fmt: str = '.8f', *args, **kwargs)[source]
Writes the geometry
- Parameters:
geometry – geometry object to write
fmt – precision used to store the atomic coordinates
unit (
{"Ang", "Bohr", "fractional", "frac"}
) – the unit used when writing the data. fractional and frac are the same.
- write_lattice(lattice: Lattice, fmt: str = '.8f', *args, **kwargs)[source]
Writes the supercell
- Parameters:
lattice – supercell object to write
fmt – precision used to store the lattice vectors
unit (
{"Ang", "Bohr"}
) – the unit used when writing the data.
- property base_file
File of the current Sile
- property file
File of the current Sile
- plot
Plotting functions for the
fdfSileSiesta
class.