sisl.Lattice
- class sisl.Lattice(cell: sisl.typing.CellLike, nsc: numpy.typing.ArrayLike = None, origin=None, boundary_condition: BoundaryCondition | int | str | bool | Sequence[BoundaryCondition | int | str | bool] = BoundaryCondition.PERIODIC)
Bases:
_Dispatchs
A cell class to retain lattice vectors and a supercell structure
The supercell structure is comprising the primary unit-cell and neighboring unit-cells. The number of supercells is given by the attribute
nsc
which is a vector with 3 elements, one per lattice vector. It describes how many times the primary unit-cell is extended along the i’th lattice vector. Fornsc[i] == 3
the supercell is made up of 3 unit-cells. One behind, the primary unit-cell and one after.- Parameters:
cell – the lattice parameters of the unit cell (the actual cell is returned from
tocell
.nsc – number of supercells along each lattice vector
origin (
(3,)
offloat
, optional) – the origin of the supercell.boundary_condition – the boundary conditions for each of the cell’s planes. Defaults to periodic boundary condition. See
BoundaryCondition
for valid enumerations.
Conversion
A dispatcher for classes, using __get__ it converts into ObjectDispatcher upon invocation from an object, or a TypeDispatcher when invoked from a class
to.Cuboid
([center, origin, orthogonal])Convert lattice parameters to a
Cuboid
to.Path
(*args, **kwargs)Lattice
writing to a sile.to.Sile
(*args, **kwargs)Lattice
writing to a sile.to.ase
(**kwargs)Lattice
conversion to an ase.Cell object.to.str
(*args, **kwargs)Lattice
writing to a sile.Methods
add
(other)Add two supercell lattice vectors to each other
add_vacuum
(vacuum, axis[, orthogonal_to_plane])Returns a new object with vacuum along the axis lattice vector
angle
(axis1, axis2[, rad])The angle between two of the cell vectors
append
(other, axis)Appends other
Lattice
to this grid along axisarea
(axis1, axis2)Calculate the area spanned by the two axis ax0 and ax1
cell2length
(length[, axes])Calculate cell vectors such that they each have length
length
center
([axes])Returns center of the
Lattice
, possibly with respect to axescopy
([cell])A deepcopy of the object
equal
(other[, atol])Check whether two lattices are equivalent
fit
(xyz[, axes, atol])Fit the supercell to xyz such that the unit-cell becomes periodic in the specified directions
is_cartesian
([atol])Checks if cell vectors a,b,c are multiples of the cartesian axis vectors (x, y, z)
is_orthogonal
([rtol])Returns true if the cell vectors are orthogonal.
lengthf
([axes])Length of specific lattice vectors (as a function) :Parameters: axes -- only calculate the volume based on a subset of axes
offset
([isc])Returns the supercell offset of the supercell index
parallel
(other[, axes])Returns true if the cell vectors are parallel to other
parameters
([rad])Cell parameters of this cell in 3 lengths and 3 angles
plane
(axis1, axis2[, origin])Query point and plane-normal for the plane spanning ax1 and ax2
prepend
(other, axis)Prepends other
Lattice
to this grid along axisread
(sile, *args, **kwargs)Reads the supercell from the
Sile
usingSile.read_lattice
repeat
(reps, axis)Extend the unit-cell reps times along the axis lattice vector
rotate
(angle, v[, rad, what])Rotates the supercell, in-place by the angle around the vector
sc_index
(sc_off)Returns the integer index in the sc_off list that corresponds to
sc_off
scale
(scale[, what])Scale lattice vectors
set_boundary_condition
([boundary, a, b, c])Set the boundary conditions on the grid
set_nsc
([nsc, a, b, c])Sets the number of supercells in the 3 different cell directions
swapaxes
(axes1, axes2[, what])Swaps axes axes1 and axes2
tile
(reps, axis)Extend the unit-cell reps times along the axis lattice vector
toCuboid
(*args, **kwargs)A cuboid with vectors as this unit-cell and center with respect to its origin
tocell
(*args)Returns a 3x3 unit-cell dependent on the input
unrepeat
(reps, axis)Reverses a
Lattice.tile
and returns the segmented versionuntile
(reps, axis)Reverses a
Lattice.tile
and returns the segmented versionvertices
()Vertices of the cell
volumef
([axes])Volume of cell (as a function)
write
(sile, *args, **kwargs)Writes latticey to the sile using sile.write_lattice
Attributes
Boundary conditions for each lattice vector (lower/upper) sides
(3, 2)
Returns the reciprocal (inverse) cell for the
Lattice
.Internal indexed supercell
[ia, ib, ic] == i
Length of each lattice vector
Origin for the cell
Boolean array to specify whether the boundary conditions are periodic`
Returns the reciprocal cell for the
Lattice
with2*np.pi
Integer supercell offsets
Volume of cell
- BC
alias of
BoundaryCondition
- add(other) Lattice
Add two supercell lattice vectors to each other
- Parameters:
other (
Lattice
,array_like
) – the lattice vectors of the other supercell to add
- add_vacuum(vacuum: float, axis: sisl.typing.CellAxis, orthogonal_to_plane: bool = False) Lattice [source]
Returns a new object with vacuum along the axis lattice vector
- angle(axis1: sisl.typing.CellAxis, axis2: sisl.typing.CellAxis, rad: bool = False) float [source]
The angle between two of the cell vectors
- Parameters:
axis1 – the first cell vector
axis2 – the second cell vector
rad (
bool
, optional) – whether the returned value is in radians
- area(axis1: sisl.typing.CellAxis, axis2: sisl.typing.CellAxis) float [source]
Calculate the area spanned by the two axis ax0 and ax1
- cell2length(length, axes: sisl.typing.CellAxes = (0, 1, 2)) ndarray [source]
Calculate cell vectors such that they each have length
length
- Parameters:
- Returns:
numpy.ndarray
– cell-vectors with prescribed length, same order as axes
- center(axes: sisl.typing.CellAxes = (0, 1, 2)) ndarray
Returns center of the
Lattice
, possibly with respect to axes
- copy(cell=None, **kwargs) Lattice
A deepcopy of the object
- Parameters:
cell (
array_like
) – the new cell parameters
- fit(xyz, axes: sisl.typing.CellAxes = (0, 1, 2), atol: float = 0.05) Lattice [source]
Fit the supercell to xyz such that the unit-cell becomes periodic in the specified directions
The fitted supercell tries to determine the unit-cell parameters by solving a set of linear equations corresponding to the current supercell vectors.
>>> numpy.linalg.solve(self.cell.T, xyz.T)
It is important to know that this routine will only work if at least some of the atoms are integer offsets of the lattice vectors. I.e. the resulting fit will depend on the translation of the coordinates.
- Parameters:
xyz (
array_like ``shape(*
,3)``
) – the coordinates that we will wish to encompass and analyze.axes – only the cell-vectors along the provided axes will be used
atol – tolerance (in Angstrom) of the positions. I.e. we neglect coordinates which are not within the radius of this magnitude
- Raises:
RuntimeError : – when the cell-parameters does not fit within the given tolerance (atol).
- is_cartesian(atol: float = 0.001) bool [source]
Checks if cell vectors a,b,c are multiples of the cartesian axis vectors (x, y, z)
- Parameters:
atol (
float
, optional) – the threshold above which an off diagonal term will be considered non-zero.
- is_orthogonal(rtol: float = 0.001) bool [source]
Returns true if the cell vectors are orthogonal.
Internally this will be done on the normalized lattice vectors to ensure no numerical instability.
- Parameters:
rtol (
float
, optional) – the threshold above which the scalar product of two normalized cell vectors will be considered non-zero.
- lengthf(axes: sisl.typing.CellAxes = (0, 1, 2)) ndarray [source]
Length of specific lattice vectors (as a function) :Parameters: axes – only calculate the volume based on a subset of axes
Examples
Only get lengths of two lattice vectors:
>>> lat = Lattice(1) >>> lat.lengthf([0, 1])
- new.Sile(**kwargs)
- new.ase(**kwargs)
- new.ndarray(**kwargs)
- offset(isc=None) tuple[float, float, float] [source]
Returns the supercell offset of the supercell index
- parallel(other, axes: sisl.typing.CellAxes = (0, 1, 2)) bool [source]
Returns true if the cell vectors are parallel to other
- Parameters:
other (
Lattice
) – the other object to check whether the axis are parallelaxes – only check the specified axes (default to all)
- parameters(rad: bool = False) tuple[float, float, float, float, float, float] [source]
Cell parameters of this cell in 3 lengths and 3 angles
Notes
Since we return the length and angles between vectors it may not be possible to recreate the same cell. Only in the case where the first lattice vector only has a Cartesian \(x\) component will this be the case.
- Parameters:
rad (
bool
, optional) – whether the angles are returned in radians (otherwise in degree)- Returns:
length (
numpy.ndarray
) – length of each lattice vectorangles (
numpy.ndarray
) – angles between the lattice vectors (in Voigt notation)[0]
is between 2nd and 3rd lattice vector, etc.
- plane(axis1: sisl.typing.CellAxis, axis2: sisl.typing.CellAxis, origin: bool = True) tuple[ndarray, ndarray] [source]
Query point and plane-normal for the plane spanning ax1 and ax2
- Parameters:
axis1 – the first axis vector
axis2 – the second axis vector
origin (
bool
, optional) – whether the plane intersects the origin or the opposite corner of the unit-cell.
- Returns:
normal_V (
numpy.ndarray
) – planes normal vector (pointing outwards with regards to the cell)p (
numpy.ndarray
) – a point on the plane
Examples
All 6 faces of the supercell can be retrieved like this:
>>> lattice = Lattice(4) >>> n1, p1 = lattice.plane(0, 1, True) >>> n2, p2 = lattice.plane(0, 1, False) >>> n3, p3 = lattice.plane(0, 2, True) >>> n4, p4 = lattice.plane(0, 2, False) >>> n5, p5 = lattice.plane(1, 2, True) >>> n6, p6 = lattice.plane(1, 2, False)
However, for performance critical calculations it may be advantageous to do this:
>>> lattice = Lattice(4) >>> uc = lattice.cell.sum(0) >>> n1, p1 = lattice.plane(0, 1) >>> n2 = -n1 >>> p2 = p1 + uc >>> n3, p3 = lattice.plane(0, 2) >>> n4 = -n3 >>> p4 = p3 + uc >>> n5, p5 = lattice.plane(1, 2) >>> n6 = -n5 >>> p6 = p5 + uc
Secondly, the variables
p1
,p3
andp5
are always[0, 0, 0]
andp2
,p4
andp6
are alwaysuc
. Hence this may be used to further reduce certain computations.
- static read(sile, *args, **kwargs) Lattice [source]
Reads the supercell from the
Sile
usingSile.read_lattice
- Parameters:
sile (
Sile
,str
orpathlib.Path
) – aSile
object which will be used to read the supercell if it is a string it will create a new sile usingsisl.io.get_sile
.
- repeat(reps: int, axis: sisl.typing.CellAxis) Lattice
Extend the unit-cell reps times along the axis lattice vector
Notes
This is exactly equivalent to the
tile
routine.- Parameters:
reps – number of times the unit-cell is repeated along the specified lattice vector
axis – the lattice vector along which the repetition is performed
- rotate(angle: float, v: str | int | Coord, rad: bool = False, what: Literal['abc', 'a', ...] = 'abc') Lattice
Rotates the supercell, in-place by the angle around the vector
One can control which cell vectors are rotated by designating them individually with
only='[abc]'
.- Parameters:
angle – the angle of which the geometry should be rotated
v – the vector around the rotation is going to happen
[1, 0, 0]
will rotate in theyz
planerad – Whether the angle is in radians (True) or in degrees (False)
what – only rotate the designated cell vectors.
- sc_index(sc_off) int | Sequence[int] [source]
Returns the integer index in the sc_off list that corresponds to
sc_off
Returns the index for the supercell in the global offset.
- Parameters:
sc_off (
(3,)
orlist
of(3,)
) – super cell specification. For each axis having valueNone
all supercells along that axis is returned.
- scale(scale: CoordOrScalar, what: Literal['abc', 'xyz'] = 'abc') Lattice
Scale lattice vectors
Does not scale
origin
.- Parameters:
scale – the scale factor for the new lattice vectors.
what – If three different scale factors are provided, whether each scaling factor is to be applied on the corresponding lattice vector (“abc”) or on the corresponding cartesian coordinate (“xyz”).
- set_boundary_condition(boundary: BoundaryCondition | int | str | bool | Sequence[BoundaryCondition | int | str | bool] | None = None, a: BoundaryCondition | int | str | bool | Sequence[BoundaryCondition | int | str | bool] | None = None, b: BoundaryCondition | int | str | bool | Sequence[BoundaryCondition | int | str | bool] | None = None, c: BoundaryCondition | int | str | bool | Sequence[BoundaryCondition | int | str | bool] | None = None)[source]
Set the boundary conditions on the grid
- Parameters:
boundary (
(3
,2)
or(3
,)
orint
, optional) – boundary condition for all boundaries (or the same for all)a (
int
orlist
ofint
, optional) – boundary condition for the first unit-cell vector directionb (
int
orlist
ofint
, optional) – boundary condition for the second unit-cell vector directionc (
int
orlist
ofint
, optional) – boundary condition for the third unit-cell vector direction
- Raises:
ValueError – if specifying periodic one one boundary, so must the opposite side.
- set_nsc(nsc=None, a: int | None = None, b: int | None = None, c: int | None = None) None [source]
Sets the number of supercells in the 3 different cell directions
- Parameters:
nsc (
list
ofint
, optional) – number of supercells in each directiona (
int
, optional) – number of supercells in the first unit-cell vector directionb (
int
, optional) – number of supercells in the second unit-cell vector directionc (
int
, optional) – number of supercells in the third unit-cell vector direction
- swapaxes(axes1: AnyAxes, axes2: AnyAxes, what: Literal['abc', 'xyz', 'abc+xyz'] = 'abc') Lattice
Swaps axes axes1 and axes2
Swapaxes is a versatile method for changing the order of axes elements, either lattice vector order, or Cartesian coordinate orders.
- Parameters:
axes1 – the old axis indices (or labels if str) A string will translate each character as a specific axis index. Lattice vectors are denoted by
abc
while the Cartesian coordinates are denote byxyz
. If str, then what is not used.axes2 – the new axis indices, same as axes1
what – which elements to swap, lattice vectors (
abc
), or Cartesian coordinates (xyz
), or both. This argument is only used if the axes arguments are ints.
Examples
Swap the first two axes
>>> sc_ba = sc.swapaxes(0, 1) >>> assert np.allclose(sc_ba.cell[(1, 0, 2)], sc.cell)
Swap the Cartesian coordinates of the lattice vectors
>>> sc_yx = sc.swapaxes(0, 1, what="xyz") >>> assert np.allclose(sc_ba.cell[:, (1, 0, 2)], sc.cell)
Consecutive swapping: 1. abc -> bac 2. bac -> bca
>>> sc_bca = sc.swapaxes("ab", "bc") >>> assert np.allclose(sc_ba.cell[:, (1, 0, 2)], sc.cell)
- tile(reps: int, axis: sisl.typing.CellAxis) Lattice
Extend the unit-cell reps times along the axis lattice vector
Notes
This is exactly equivalent to the
repeat
routine.- Parameters:
reps – number of times the unit-cell is repeated along the specified lattice vector
axis – the lattice vector along which the repetition is performed
- to.Cuboid(center=None, origin=None, orthogonal=False) Cuboid
Convert lattice parameters to a
Cuboid
- to.Path(*args, **kwargs) Any
Lattice
writing to a sile.Examples
>>> geom = si.geom.graphene() >>> geom.lattice.to("hello.xyz") >>> geom.lattice.to(pathlib.Path("hello.xyz"))
- to.Sile(*args, **kwargs) Any
Lattice
writing to a sile.Examples
>>> geom = si.geom.graphene() >>> geom.lattice.to("hello.xyz") >>> geom.lattice.to(pathlib.Path("hello.xyz"))
- to.str(*args, **kwargs) Any
Lattice
writing to a sile.Examples
>>> geom = si.geom.graphene() >>> geom.lattice.to("hello.xyz") >>> geom.lattice.to(pathlib.Path("hello.xyz"))
- toCuboid(*args, **kwargs)[source]
A cuboid with vectors as this unit-cell and center with respect to its origin
- Parameters:
orthogonal (
bool
, optional) – if true the cuboid has orthogonal sides such that the entire cell is contained
- classmethod tocell(*args) Lattice [source]
Returns a 3x3 unit-cell dependent on the input
- 1 argument
a unit-cell along Cartesian coordinates with side-length equal to the argument.
- 3 arguments
the diagonal components of a Cartesian unit-cell
- 6 arguments
the cell parameters given by \(a\), \(b\), \(c\), \(\alpha\), \(\beta\) and \(\gamma\) (angles in degrees).
- 9 arguments
a 3x3 unit-cell.
- Parameters:
*args (
float
) – May be either, 1, 3, 6 or 9 elements. Note that the arguments will be put into an array and flattened before checking the number of arguments.
Examples
>>> cell_1_1_1 = Lattice.tocell(1.) >>> cell_1_2_3 = Lattice.tocell(1., 2., 3.) >>> cell_1_2_3 = Lattice.tocell([1., 2., 3.]) # same as above
- unrepeat(reps: int, axis: sisl.typing.CellAxis) Lattice
Reverses a
Lattice.tile
and returns the segmented versionNotes
Untiling will not correctly re-calculate nsc since it has no knowledge of connections.
See also
Lattice.tile
opposite of this method
- untile(reps: int, axis: sisl.typing.CellAxis) Lattice
Reverses a
Lattice.tile
and returns the segmented versionNotes
Untiling will not correctly re-calculate nsc since it has no knowledge of connections.
See also
Lattice.tile
opposite of this method
- vertices() ndarray [source]
Vertices of the cell
- Returns:
array
ofshape (2
,2
,2
,3)
– The coordinates of the vertices of the cell. The first three dimensions correspond to each cell axis (off, on), and the last one contains the xyz coordinates.
- volumef(axes: sisl.typing.CellAxes = (0, 1, 2)) float [source]
Volume of cell (as a function)
Default to the 3D volume. For axes with only 2 elements, it corresponds to an area. For axes with only 1 element, it corresponds to a length.
- Parameters:
axes – only calculate the volume based on a subset of axes
Examples
Only get the volume of the periodic directions:
>>> lat = Lattice(1) >>> lat.pbc = (True, False, True) >>> lat.volumef(lat.pbc.nonzero()[0])
- write(sile: sisl.typing.SileLike, *args, **kwargs) None
Writes latticey to the sile using sile.write_lattice
- Parameters:
sile – a
Sile
object which will be used to write the lattice if it is a string it will create a new sile usingget_sile
*args, **kwargs – Any other args will be passed directly to the underlying routine
See also
Lattice.read
reads a
Lattice
from a givenSile
/file
- property boundary_condition: ndarray
Boundary conditions for each lattice vector (lower/upper) sides
(3, 2)
- cell
- property icell: ndarray
Returns the reciprocal (inverse) cell for the
Lattice
.Note: The returned vectors are still in
[0, :]
format and not as returned by an inverse LAPACK algorithm.
- n_s
- new = <TypeDispatcher{obj=<class 'sisl.Lattice'>}>
- nsc
- property rcell: ndarray
Returns the reciprocal cell for the
Lattice
with2*np.pi
Note: The returned vectors are still in [0, :] format and not as returned by an inverse LAPACK algorithm.
- to
A dispatcher for classes, using __get__ it converts into ObjectDispatcher upon invocation from an object, or a TypeDispatcher when invoked from a class
This is a class-placeholder allowing a dispatcher to be a class attribute and converted into an ObjectDispatcher when invoked from an object.
If it is called on the class, it will return a TypeDispatcher.
This class should be an attribute of a class. It heavily relies on the __get__ special method.
- Parameters:
name (
str
) – name of the attribute in the classdispatchs (
dict
, optional) – dictionary of dispatch methodsobj_getattr (
callable
, optional) – method with 2 arguments, anobj
and theattr
which may be used to control how the attribute is called.instance_dispatcher (
AbstractDispatcher
, optional) – control how instance dispatchers are handled through __get__ method. This controls the dispatcher used if called from an instance.type_dispatcher (
AbstractDispatcher
, optional) – control how class dispatchers are handled through __get__ method. This controls the dispatcher used if called from a class.
Examples
>>> class A: ... new = ClassDispatcher("new", obj_getattr=lambda obj, attr: getattr(obj.sub, attr))
The above defers any attributes to the contained A.sub attribute.