sisl.Geometry

class sisl.Geometry(xyz: numpy.typing.ArrayLike, atoms: sisl.typing.AtomsLike | None = None, lattice: sisl.typing.LatticeLike | None = None, names=None)

Bases: LatticeChild, _Dispatchs

Holds atomic information, coordinates, species, lattice vectors

The Geometry class holds information regarding atomic coordinates, the atomic species, the corresponding lattice-vectors.

It enables the interaction and conversion of atomic structures via simple routine methods.

All lengths are assumed to be in units of Angstrom, however, as long as units are kept same the exact units are irrespective.

>>> square = Geometry([[0.5, 0.5, 0.5]], Atom(1),
...                   lattice=Lattice([1, 1, 10], nsc=[3, 3, 1]))
>>> print(square)
Geometry{na: 1, no: 1,
 Atoms{species: 1,
  Atom{H, Z: 1, mass(au): 1.00794, maxR: -1.00000,
   Orbital{R: -1.00000, q0: 0.0}
  }: 1,
 },
 maxR: -1.00000,
 Lattice{volume: 1.0000e+01, nsc: [3 3 1]}
}
Parameters:
  • xyz – atomic coordinates xyz[i, :] is the atomic coordinate of the i’th atom.

  • atoms – atomic species retrieved from the PeriodicTable

  • lattice – the unit-cell describing the atoms in a periodic super-cell

Examples

An atomic cubic lattice of Hydrogen atoms

>>> xyz = [[0, 0, 0],
...        [1, 1, 1]]
>>> sc = Lattice([2,2,2])
>>> g = Geometry(xyz, Atom('H'), sc)

The following estimates the lattice vectors from the atomic coordinates, although possible, it is not recommended to be used.

>>> xyz = [[0, 0, 0],
...        [1, 1, 1]]
>>> g = Geometry(xyz, Atom('H'))

Conversion of geometries to other projects instances can be done via sisl’s dispatch functionality

>>> g.to.ase()
Atoms(...)

converts to an ASE ase.Atoms object.

See also

Atoms

contained atoms self.atoms

Atom

contained atoms are each an object of this

Conversion

new

new.Sile

new.ase

new.pymatgen

to

A dispatcher for classes, using __get__ it converts into ObjectDispatcher upon invocation from an object, or a TypeDispatcher when invoked from a class

to.Path(*args, **kwargs)

Writes the geometry to a sile with any optional arguments.

to.Sile(*args, **kwargs)

Writes the geometry to a sile with any optional arguments.

to.ase(**kwargs)

Conversion of Geometry to an ase.Atoms object

to.dataframe(*args, **kwargs)

Convert the geometry to a pandas.DataFrame with values stored in columns

to.pymatgen(**kwargs)

Conversion of Geometry to a pymatgen object.

to.str(*args, **kwargs)

Writes the geometry to a sile with any optional arguments.

Plotting

plot

Plotting functions for the Geometry class.

plot.geometry([axes, atoms, ...])

Builds a GeometryPlot by setting the value of "geometry" to the current object.

Methods

Rij(ia, ja)

Vector between atom ia and ja, atoms can be in super-cell indices

a2isc(atoms)

Super-cell indices for a specific/list atom

a2o(atoms[, all])

Returns an orbital index of the first orbital of said atom.

a2sc(atoms)

Returns the super-cell offset for a specific atom

a2transpose(atoms1[, atoms2])

Transposes connections from atoms1 to atoms2 such that supercell connections are transposed

add(other[, offset])

Merge two geometries (or a Geometry and Lattice) by adding the two atoms together

add_vacuum(vacuum, axis[, offset])

Add vacuum along the axis lattice vector

angle(atoms[, dir, ref, rad])

The angle between atom atoms and the direction dir, with possibility of a reference coordinate ref

append(other, axis[, offset])

Appends two structures along axis

area(ax0, ax1)

Calculate the area spanned by the two axis ax0 and ax1

as_primary(na_primary[, axes, ret_super])

Reduce the geometry to the primary unit-cell comprising na_primary atoms

as_supercell()

Create a new geometry equivalent to self * self.nsc, where the indices are ordered as the supercells

asc2uc(atoms[, unique])

Returns atoms from supercell indices to unit-cell indices, possibly removing dublicates

attach(atom, other, other_atom[, dist, axis])

Attaches another Geometry at the atom index with respect to other_atom using different methods.

auc2sc(atoms[, unique])

Returns atom from unit-cell indices to supercell indices, possibly removing dublicates

axyz([atoms, isc])

Return the atomic coordinates in the supercell of a given atom.

bond_correct(ia, atoms[, method])

Corrects the bond between ia and the atoms.

center([atoms, what])

Returns the center of the geometry

close(xyz_ia[, R, atoms, atoms_xyz, ...])

Indices of atoms in the entire supercell within a given radius from a given coordinate

close_sc(xyz_ia[, isc, R, atoms, atoms_xyz, ...])

Indices of atoms in a given supercell within a given radius from a given coordinate

copy()

Create a new object with the same content (a copy).

dihedral(atoms[, rad])

Calculate the dihedral angle defined by four atoms.

distance([atoms, R, atol, method])

Calculate the distances for all atoms in shells of radius tol within max_R

equal(other[, R, atol])

Whether two geometries are the same (optional not check of the orbital radius)

find_nsc([axes, R, method])

Find number of supercells for the geometry, depending on certain criteria

iR([na, iR, R])

Return an integer number of maximum radii (self.maxR()) which holds approximately na atoms

insert(atom, other)

Inserts other atoms right before index

iter()

An iterator over all atomic indices

iter_block([iR, R, atoms, method])

Iterator for performance critical loops

iter_block_rand([iR, R, atoms])

Perform the random block-iteration by randomly selecting the next center of block

iter_block_shape([shape, iR, atoms])

Perform the grid block-iteration by looping a grid

iter_orbitals([atoms, local])

Returns an iterator over all atoms and their associated orbitals

iter_species([atoms])

Iterator over all atoms (or a subset) and species as a tuple in this geometry

maxR([all])

Maximum orbital range of the atoms

mirror(method[, atoms, point])

Mirrors the atomic coordinates about a plane given by its normal vector

move(*args, **kwargs)

See translate for details

o2a(orbitals[, unique])

Atomic index corresponding to the orbital indicies.

o2isc(orbitals)

Returns the super-cell index for a specific orbital.

o2sc(orbitals)

Returns the super-cell offset for a specific orbital.

o2transpose(orb1[, orb2])

Transposes connections from orb1 to orb2 such that supercell connections are transposed

oRij(orbitals1, orbitals2)

Vector between orbital orbitals1 and orbitals2, orbitals can be in super-cell indices

optimize_nsc([axes, R])

Optimize the number of supercell connections based on self.maxR()

orij(orbitals1, orbitals2)

Distance between orbital orbitals1 and orbitals2, orbitals can be in super-cell indices

osc2uc(orbitals[, unique])

Orbitals from supercell indices to unit-cell indices, possibly removing dublicates

ouc2sc(orbitals[, unique])

Orbitals from unit-cell indices to supercell indices, possibly removing dublicates

overlap(other[, atol, offset, offset_other])

Calculate the overlapping indices between two geometries

prepend(other, axis[, offset])

Prepend two structures along axis

read(sile, *args, **kwargs)

Reads geometry from the Sile using Sile.read_geometry

reduce()

Remove all atoms not currently used in the self.atoms object

remove(atoms)

Remove atoms from the geometry.

remove_orbital(atoms, orbitals)

Remove a subset of orbitals on atoms according to orbitals

reorder()

Reorders atoms according to first occurence in the geometry

repeat(reps, axis)

Create a repeated geometry

replace(atoms, other[, offset])

Create a new geometry from self and replace atoms with other

reverse([atoms])

Returns a reversed geometry

rij(ia, ja)

Distance between atom ia and ja, atoms can be in super-cell indices

rotate(angle, v[, origin, atoms, rad, what])

Rotate geometry around vector and return a new geometry

rotate_miller(m, v)

Align Miller direction along v

sc2uc(atoms[, unique])

Returns atoms from supercell indices to unit-cell indices, possibly removing dublicates

sc_index(*args, **kwargs)

Call local Lattice.sc_index function

scale(scale[, what, scale_basis])

Scale coordinates and unit-cell to get a new geometry with proper scaling

set_lattice(lattice)

Overwrites the local lattice.

set_nsc(*args, **kwargs)

Set the number of super-cells in the Lattice object

set_supercell(lattice)

Overwrites the local lattice.

sort(**kwargs)

Sort atoms in a nested fashion according to various criteria

sparserij([dtype, na_iR, method])

Return the sparse matrix with all distances in the matrix The sparse matrix will only be defined for the elements which have orbitals overlapping with other atoms.

sub(atoms)

Create a new Geometry with a subset of this Geometry

sub_orbital(atoms, orbitals)

Retain only a subset of the orbitals on atoms according to orbitals

swap(atoms1, atoms2)

Swap a set of atoms in the geometry and return a new one

swapaxes(axes1, axes2[, what])

Swap the axes components by either lattice vectors (only cell), or Cartesian coordinates

tile(reps, axis)

Tile the geometry to create a bigger one

translate(v[, atoms])

Translates the geometry by v

translate2uc([atoms, axes])

Translates atoms in the geometry into the unit cell

uc2sc(atoms[, unique])

Returns atom from unit-cell indices to supercell indices, possibly removing dublicates

unrepeat(reps, axis, *args, **kwargs)

Unrepeats the geometry similarly as untile

untile(reps, axis[, segment, rtol, atol])

A subset of atoms from the geometry by cutting the geometry into reps parts along the direction axis.

within(shapes[, atoms, atoms_xyz, ret_xyz, ...])

Indices of atoms in the entire supercell within a given shape from a given coordinate

within_inf(lattice[, periodic, atol, origin])

Find all atoms within a provided supercell

within_sc(shapes[, isc, atoms, atoms_xyz, ...])

Indices of atoms in a given supercell within a given shape from a given coordinate

write(sile, *args, **kwargs)

Writes geometry to the sile using sile.write_geometry

Attributes

atoms

Atoms associated with the geometry

boundary_condition

cell

Returns the inherent Lattice.cell

firsto

The first orbital on the corresponding atom

fxyz

Returns geometry coordinates in fractional coordinates

icell

Returns the inherent Lattice.icell

isc_off

Returns the inherent Lattice.isc_off

lasto

The last orbital on the corresponding atom

length

Returns the inherent Lattice.length

mass

The mass of all atoms as an array

n_s

Returns the inherent Lattice.n_s

na

Number of atoms in geometry

na_s

Number of supercell atoms

names

The named index specifier

no

Number of orbitals in unit cell

no_s

Number of supercell orbitals

nsc

Returns the inherent Lattice.nsc

orbitals

List of orbitals per atom

origin

Returns the inherent Lattice.origin

pbc

q0

Total initial charge in this geometry (sum of q0 off all atoms)

rcell

Returns the inherent Lattice.rcell

sc

[deprecated] Return the lattice object associated with the Lattice.

sc_off

Returns the inherent Lattice.sc_off

volume

Returns the inherent Lattice.volume

Rij(ia: sisl.typing.AtomsIndex, ja: sisl.typing.AtomsIndex) ndarray[source]

Vector between atom ia and ja, atoms can be in super-cell indices

Returns the vector between two atoms:

\[\mathbf r^{IJ} = \mathbf r^J - \mathbf r^I\]
Parameters:
  • ia – atomic index of first atom

  • ja – atomic indices

a2isc(atoms: sisl.typing.AtomsIndex) ndarray[source]

Super-cell indices for a specific/list atom

Returns a vector of 3 numbers with integers. Any multi-dimensional input will be flattened before return.

The returned indices will thus always be a 2D matrix or a 1D vector.

Parameters:

atoms – atom indices to extract the supercell locations of

a2o(atoms: sisl.typing.AtomsIndex, all: bool = False) ndarray[source]

Returns an orbital index of the first orbital of said atom. This is particularly handy if you want to create TB models with more than one orbital per atom.

Note that this will preserve the super-cell offsets.

Parameters:
  • atoms – Atomic indices

  • allFalse, return only the first orbital corresponding to the atom, True, returns list of the full atom(s), will always return a 1D array.

a2sc(atoms: sisl.typing.AtomsIndex) ndarray[source]

Returns the super-cell offset for a specific atom

Parameters:

atoms – atom indices to extract the supercell offsets of

a2transpose(atoms1: sisl.typing.AtomsIndex, atoms2: sisl.typing.AtomsIndex = None) tuple[ndarray, ndarray][source]

Transposes connections from atoms1 to atoms2 such that supercell connections are transposed

When handling supercell indices it is useful to get the transposed connection. I.e. if you have a connection from site i (in unit cell indices) to site j (in supercell indices) it may be useful to get the equivalent supercell connection such for site j (in unit cell indices) to site i (in supercell indices) such that they correspond to the transposed coupling.

Note that since this transposes couplings the indices returned are always expanded to the full length if either of the inputs are a single index.

Examples

>>> gr = geom.graphene()
>>> atoms = gr.close(0, 1.5)
>>> atoms
array([0, 1, 5, 9], dtype=int32)
>>> gr.a2transpose(0, atoms)
(array([0, 1, 1, 1], dtype=int32), array([ 0,  0, 14, 10], dtype=int32))
Parameters:
  • atoms1 – atomic indices must have same length as atoms2 or length 1

  • atoms2 – atomic indices must have same length as atoms1 or length 1. If not present then only atoms1 will be returned in transposed indices.

Returns:

  • atoms2 (array_like) – transposed indices for atoms2 (only returned if atoms2 is not None)

  • atoms1 (array_like) – transposed indices for atoms1

add(other: sisl.typing.LatticeOrGeometryLike, offset: sisl.typing.Coord = (0, 0, 0)) Geometry

Merge two geometries (or a Geometry and Lattice) by adding the two atoms together

If other is a Geometry only the atoms gets added, to also add the supercell vectors simply do geom.add(other).add(other.lattice).

Parameters:
  • other – Other geometry class which is added

  • offset – offset in geometry of other when adding the atoms. Otherwise it is the offset of the geometry atoms.

See also

Geometry.append

appending geometries

Geometry.prepend

prepending geometries

Geometry.attach

attach a geometry

Geometry.insert

insert a geometry

Examples

>>> first = Geometry(...)
>>> second = Geometry(...)
>>> lattice = Lattice(...)
>>> added = first.add(second, offset=(0, 0, 2))
>>> assert np.allclose(added.xyz[:len(first)], first.xyz)
>>> assert np.allclose(added.xyz[len(first):] - [0, 0, 2], second.xyz)
add_vacuum(vacuum: float, axis: int, offset: Sequence[float] = (0, 0, 0)) Geometry[source]

Add vacuum along the axis lattice vector

When the vacuum is bigger than the maximum orbital ranges the number of supercells along that axis will be truncated to 1 (de-couple images).

Parameters:
  • vacuum – amount of vacuum added, in Ang

  • axis – the lattice vector to add vacuum along

  • offset – offset in geometry when adding the vacuum.

Returns:

Geometry (a new geometry with added vacuum)

angle(atoms: sisl.typing.AtomsIndex, dir: str | int | Sequence[int] = (1.0, 0, 0), ref: int | Sequence[float] | None = None, rad: bool = False) float | ndarray[source]

The angle between atom atoms and the direction dir, with possibility of a reference coordinate ref

The calculated angle can be written as this

\[\theta = \arccos \frac{(\mathbf r^I - \mathbf{r^{\mathrm{ref}}})\cdot \mathbf{d}} {|\mathbf r^I-\mathbf{r^{\mathrm{ref}}}||\mathbf{d}|}\]

and thus lies in the interval \([0 ; \pi]\) as one cannot distinguish orientation without additional vectors.

Parameters:
  • atoms – indices/boolean of all atoms where angles should be calculated on

  • dir – the direction from which the angle is calculated from, default to x. An integer specifies the corresponding lattice vector as the direction.

  • ref – the reference point from which the vectors are drawn, default to origin An integer specifies an atomic index.

  • rad – whether the returned value is in radians

append(other: sisl.typing.LatticeOrGeometryLike, axis: sisl.typing.CellAxis, offset: Literal['none', 'min'] | sisl.typing.Coord = 'none') Geometry

Appends two structures along axis

This will automatically add the geometry.cell[axis] to all atomic coordiates in the other structure before appending.

The basic algorithm is this:

>>> oxa = other.xyz + geometry.cell[axis][None,:]
>>> geometry.xyz = np.append(geometry.xyz,oxa)
>>> geometry.cell[axis] += other.cell[axis]

NOTE: The cell appended is only in the axis that is appended, which means that the other cell directions need not conform.

Parameters:
  • other – Other geometry class which needs to be appended If a Lattice only the super cell will be extended

  • axis – Cell direction to which the other geometry should be appended.

  • offset ({'none', 'min', (3,)}) – By default appending two structures will simply use the coordinates, as is. With ‘min’, the routine will shift both the structures along the cell axis of geometry such that they coincide at the first atom, lastly one may use a specified offset to manually select how other is displaced. NOTE: That geometry.cell[axis] will be added to offset if other is a geometry.

See also

Geometry.add

add geometries

Geometry.prepend

prepending geometries

Geometry.attach

attach a geometry

Geometry.insert

insert a geometry

apply(data, func, mapper: Callable[[int], int] | str, axis: int = 0, segments: Literal['atoms', 'orbitals', 'all'] | Iterator[int] = 'atoms') ndarray

Apply a function func to the data along axis axis using the method specified

This can be useful for applying conversions from orbital data to atomic data through sums or other functions.

The data may be of any shape but it is expected the function can handle arguments as func(data, axis=axis).

Parameters:
  • data (array_like) – the data to be converted

  • func (callable or str) – a callable function that transforms the data in some way. If a str, will use getattr(numpy, func).

  • mapper (func, optional) – a function transforming the segments into some other segments that is present in data. It can accept anything the segments returns. If a str, it will be equivalent to getattr(geometry, mapper)

  • axis – axis selector for data along which func will be applied

  • segments – which segments the mapper will recieve, if atoms, each atom index will be passed to the mapper(ia). If 'all', it will be range(data.shape[axis]).

Examples

Convert orbital data into summed atomic data

>>> g = sisl.geom.diamond(atoms=sisl.Atom(6, R=(1, 2)))
>>> orbital_data = np.random.rand(10, g.no, 3)
>>> atomic_data = g.apply(orbital_data, np.sum, mapper=partial(g.a2o, all=True), axis=1)

The same can be accomblished by passing an explicit segment iterator, note that iter(g) == range(g.na)

>>> atomic_data = g.apply(orbital_data, np.sum, mapper=partial(g.a2o, all=True), axis=1,
...                       segments=g)

To only take out every 2nd orbital:

>>> alternate_data = g.apply(orbital_data, np.sum, mapper=lambda idx: idx[::2], axis=1,
...                          segments="all")
area(ax0, ax1) float

Calculate the area spanned by the two axis ax0 and ax1

as_primary(na_primary: int, axes: Sequence[int] = (0, 1, 2), ret_super: bool = False) Geometry | tuple[Geometry, Lattice][source]

Reduce the geometry to the primary unit-cell comprising na_primary atoms

This will basically try and find the tiling/repetitions required for the geometry to only have na_primary atoms in the unit cell.

Parameters:
  • na_primary – number of atoms in the primary unit cell

  • axes – only search the given directions for supercells, default to all directions

  • ret_super – also return the number of supercells used in each direction

Returns:

  • Geometry – the primary unit cell

  • Lattice – the tiled supercell numbers used to find the primary unit cell (only if ret_super is true)

Raises:

SislError – If the algorithm fails.

as_supercell() Geometry[source]

Create a new geometry equivalent to self * self.nsc, where the indices are ordered as the supercells

Returns:

Geometry – the supercell expanded and reordered Geometry

asc2uc(atoms: sisl.typing.AtomsIndex, unique: bool = False) ndarray[source]

Returns atoms from supercell indices to unit-cell indices, possibly removing dublicates

Parameters:
  • atoms – the atomic supercell indices to be converted to unit-cell indices

  • unique – If True the returned indices are unique and sorted.

attach(atom: int, other: sisl.typing.GeometryLike, other_atom: int, dist='calc', axis: int | None = None) Geometry[source]

Attaches another Geometry at the atom index with respect to other_atom using different methods.

The attached geometry will be inserted at the end of the geometry via add.

Parameters:
  • atom (int) – atomic index which is the base position of the attachment. The distance between atom and other_atom is dist.

  • other (Geometry) – the other Geometry to attach at the given point. In this case dist from atom.

  • other_atom (int) – the index of the atom in other that is inserted at atom.

  • dist (array_like or float or str, optional) – the distance (in Ang) between the attached coordinates. If dist is array_like it should be the vector between the atoms; if dist is float the argument axis is required and the vector will be calculated along the corresponding latticevector; else if dist is str this will correspond to the method argument of the Atom.radius class of the two atoms. Here axis is also required.

  • axis (int) – specify the direction of the lattice vectors used. Not used if dist is an array-like argument.

auc2sc(atoms: sisl.typing.AtomsIndex, unique: bool = False) ndarray[source]

Returns atom from unit-cell indices to supercell indices, possibly removing dublicates

Parameters:
  • atoms – the atomic unit-cell indices to be converted to supercell indices

  • unique – If True the returned indices are unique and sorted.

axyz(atoms: sisl.typing.AtomsIndex = None, isc=None) ndarray[source]

Return the atomic coordinates in the supercell of a given atom.

The Geometry[...] slicing is calling this function with appropriate options.

Parameters:
  • atoms – atom(s) from which we should return the coordinates, the atomic indices may be in supercell format.

  • isc (array_like, optional) – Returns the atomic coordinates shifted according to the integer parts of the cell. Defaults to the unit-cell

Examples

>>> geom = Geometry([[0, 0, 0], [0.5, 0, 0]], lattice=1.)
>>> print(geom.axyz(isc=[1,0,0]))
[[1.   0.   0. ]
 [1.5  0.   0. ]]
>>> geom = Geometry([[0, 0, 0], [0.5, 0, 0]], lattice=1.)
>>> print(geom.axyz(0))
[0.  0.  0.]
bond_correct(ia: int, atoms: sisl.typing.AtomsIndex, method: str | float = 'calc') None[source]

Corrects the bond between ia and the atoms.

Corrects the bond-length between atom ia and atoms in such a way that the atomic radius is preserved. I.e. the sum of the bond-lengths minimizes the distance matrix.

Only atom ia is moved.

Parameters:
  • ia – The atom to be displaced according to the atomic radius

  • atoms – The atom(s) from which the radius should be reduced.

  • method – If str will use that as lookup in Atom.radius. Else it will be the new bond-length.

center(atoms: sisl.typing.AtomsIndex = None, what: Literal['COP|xyz|position', 'mm:xyz', 'mm:lattice|mm:cell', 'COM|mass', 'COMM:pbc|mass:pbc', 'COU|lattice|cell'] = 'xyz') ndarray

Returns the center of the geometry

By specifying what one can control whether it should be:

  • COP|xyz|position: Center of coordinates (default)

  • mm:xyz: Center of minimum+maximum of Cartesian coordinates

  • mm:lattice|mm:cell: Center of minimum+maximum of lattice vectors Cartesian coordinates

  • COM|mass: Center of mass

  • COM:pbc|mass:pbc: Center of mass using periodicity, if the point 0, 0, 0 is returned it

    may likely be because of a completely periodic system with no true center of mass

  • COU|lattice|cell: Center of lattice vectors

Parameters:
  • atoms – list of atomic indices to find center of

  • what – determine which center to calculate

close(xyz_ia, R=None, atoms: sisl.typing.AtomsIndex = None, atoms_xyz=None, ret_xyz: bool = False, ret_rij: bool = False, ret_isc: bool = False)[source]

Indices of atoms in the entire supercell within a given radius from a given coordinate

This heavily relies on the close_sc method.

Note that if a connection is made in a neighboring super-cell then the atomic index is shifted by the super-cell index times number of atoms. This allows one to decipher super-cell atoms from unit-cell atoms.

Parameters:
  • xyz_ia (coordinate/index) – Either a point in space or an index of an atom. If an index is passed it is the equivalent of passing the atomic coordinate close_sc(self.xyz[xyz_ia,:]).

  • R ((None), float/tuple of float) – The radii parameter to where the atomic connections are found. If R is an array it will return the indices: in the ranges:

    >>> ( x <= R[0] , R[0] < x <= R[1], R[1] < x <= R[2] )
    

    If a single float it will return:

    >>> x <= R
    
  • atoms – List of indices for atoms that are to be considered

  • atoms_xyz (array_like, optional) – The atomic coordinates of the equivalent atoms variable (atoms must also be passed)

  • ret_xyz – If true this method will return the coordinates for each of the couplings.

  • ret_rij – If true this method will return the distances from the xyz_ia for each of the couplings.

  • ret_isc – If true this method will return the lattice offset from xyz_ia for each of the couplings.

Returns:

  • index – indices of atoms (in supercell indices) within the shells of radius R

  • xyz – atomic coordinates of the indexed atoms (only for true ret_xyz)

  • rij – distance of the indexed atoms to the center coordinate (only for true ret_rij)

  • isc – integer lattice offsets for the couplings (related to rij without atomic coordinates)

close_sc(xyz_ia, isc=(0, 0, 0), R=None, atoms: sisl.typing.AtomsIndex = None, atoms_xyz=None, ret_xyz: bool = False, ret_rij: bool = False)[source]

Indices of atoms in a given supercell within a given radius from a given coordinate

This returns a set of atomic indices which are within a sphere of radius R.

If R is a tuple/list/array it will return the indices: in the ranges:

>>> ( x <= R[0] , R[0] < x <= R[1], R[1] < x <= R[2] )
Parameters:
  • xyz_ia (array_like of floats or int) – Either a point in space or an index of an atom. If an index is passed it is the equivalent of passing the atomic coordinate close_sc(self.xyz[xyz_ia,:]).

  • isc ((3,), optional) – Integer super-cell offsets in which the coordinates are checked in. I.e. isc=[0, 0, 0] is the primary cell (default).

  • R (float or array_like, optional) – The radii parameter to where the atomic connections are found. If R is an array it will return the indices: in the ranges ( x <= R[0] , R[0] < x <= R[1], R[1] < x <= R[2] ). If a single float it will return x <= R.

  • atoms – List of atoms that will be considered. This can be used to only take out a certain atom.

  • atoms_xyz (array_like of float, optional) – The atomic coordinates of the equivalent atoms variable (atoms must also be passed)

  • ret_xyz – If True this method will return the coordinates for each of the couplings.

  • ret_rij – If True this method will return the distance for each of the couplings.

Returns:

  • index – indices of atoms (in supercell indices) within the shells of radius R

  • xyz – atomic coordinates of the indexed atoms (only for true ret_xyz)

  • rij – distance of the indexed atoms to the center coordinate (only for true ret_rij)

copy() Geometry

Create a new object with the same content (a copy).

dihedral(atoms: sisl.typing.AtomsIndex, rad: bool = False) float | ndarray[source]

Calculate the dihedral angle defined by four atoms.

The dihehral angle is defined between 2 half-planes.

The first 3 atoms define the first plane The last 3 atoms define the second.

The dihedral angle is calculated using this formula:

\[\begin{split}\mathbf u_0 &= \mathbf r_1 - \mathbf r_0 \\ \mathbf u_1 &= \mathbf r_2 - \mathbf r_1 \\ \mathbf u_2 &= \mathbf r_3 - \mathbf r_2 \\ \phi &= \operatorname{atan2}\Big( \hat{\mathbf u}_0\cdot (\hat{\mathbf u}_1\times\hat{\mathbf u}_2), (\hat{\mathbf u}_0\times\hat{\mathbf u}_1)\cdot (\hat{\mathbf u}_1\times\hat{\mathbf u}_2) \Big)\end{split}\]

Where \(\hat{\cdot}\) means the unit-vector.

Parameters:
  • atoms – An array of shape (4,) or (*, 4) representing the indices of 4 atoms forming the dihedral angle

  • rad – whether the returned value is in radians

distance(atoms: AtomsIndex = None, R: float | None = None, atol: float | Sequence[float] = 0.1, method: Callable[[Sequence[float]], float] | Literal['average', 'mode', '<numpy.method>'] = 'average') float | ndarray[source]

Calculate the distances for all atoms in shells of radius tol within max_R

Parameters:
  • atoms – only create list of distances from the given atoms, default to all atoms

  • R – the maximum radius to consider, default to self.maxR(). To retrieve all distances for atoms within the supercell structure you can pass numpy.inf.

  • atol – the tolerance for grouping a set of atoms. This parameter sets the shell radius for each shell. I.e. the returned distances between two shells will be maximally 2*atol, but only if atoms are within two consecutive lists. If this is a list, the shells will be of unequal size.

    The first shell size will be atol * .5 or atol[0] * .5 if atol is a list.

  • method – How the distance in each shell is determined. A list of distances within each shell is gathered and the equivalent method will be used to extract a single quantity from the list of distances in the shell. If ‘mode’ is chosen it will use scipy.stats.mode. If another string is given it will correspond to getattr(numpy, method), while any callable function may be passed. The passed function will only be passed a list of unsorted distances that needs to be processed.

Notes

Using method='mode' requires scipy>=1.9.

Examples

>>> geom = Geometry([0]*3, Atom(1, R=1.), lattice=Lattice(1., nsc=[5, 5, 1]))
>>> geom.distance()
array([1.])
>>> geom.distance(atol=[0.5, 0.4, 0.3, 0.2])
array([1.])
>>> geom.distance(R=2, atol=[0.5, 0.4, 0.3, 0.2])
array([1.        ,  1.41421356,  2.        ])
>>> geom.distance(R=2, atol=[0.5, 0.7]) # the R = 1 and R = 2 ** .5 gets averaged
array([1.20710678,  2.        ])
Returns:

numpy.ndarray – an array of positive numbers yielding the distances from the atoms in reduced form

See also

sparserij

return a sparse matrix will all distances between atoms

equal(other: sisl.typing.GeometryLike, R: bool = True, atol: float = 1e-4) bool[source]

Whether two geometries are the same (optional not check of the orbital radius)

Parameters:
  • other – the other Geometry to check against

  • R – if True also check if the orbital radii are the same (see Atom.equal)

  • atol – tolerance for checking the atomic coordinates

find_nsc(axes: CellAxes | None = None, R: float | None = None, method: Literal['atoms', 'cell', 'overlap'] = 'atoms') ndarray[source]

Find number of supercells for the geometry, depending on certain criteria

This can find the optimal nsc values for a given method.

The important parameter, method determines how nsc is found. The method are shown here, from method that produces the smallest nsc, up to the largest nsc.

method=atoms

here only the atoms ranges are taken into account, and only whether atoms in the primary unit cell can connect to others in neigboring cells.

method=cell

only the atoms ranges are taken into account. For instance if a lattice vector is as long as the orbital range it will have 3 supercells (it can only connect to its neighboring cells).

method=overlap

determine nsc by examining at what range two orbitals overlaps.

Parameters:
  • axes – only discover new nsc the specified axes (defaults to all)

  • R – the maximum connection radius for each atom, defaults to self.maxR().

  • method – See discussion above.

Returns:

numpy.ndarray (the found nsc that obeys `method`)

See also

optimize_nsc

same as this, but equivalent to also doing self.set_nsc(self.find_nsc(...))

iR(na: int = 1000, iR: int = 20, R: float | None = None) int[source]

Return an integer number of maximum radii (self.maxR()) which holds approximately na atoms

Parameters:
  • na – number of atoms within the radius

  • iR – initial iR value, which the sphere is estitametd from

  • R – the value used for atomic range (defaults to self.maxR())

Returns:

int – number of radius needed to contain na atoms. Minimally 2 will be returned.

insert(atom: sisl.typing.AtomsIndex, other: sisl.typing.GeometryLike) Geometry

Inserts other atoms right before index

We insert the geometry Geometry before atom. Note that this will not change the unit cell.

Parameters:
  • atom – the atomic index at which the other geometry is inserted

  • other – the other geometry to be inserted

See also

Geometry.add

add geometries

Geometry.append

appending geometries

Geometry.prepend

prepending geometries

Geometry.attach

attach a geometry

iter() Iterator[int][source]

An iterator over all atomic indices

This iterator is the same as:

>>> for ia in range(len(self)):
...    <do something>

or equivalently

>>> for ia in self:
...    <do something>

See also

iter_species

iterate across indices and atomic species

iter_orbitals

iterate across atomic indices and orbital indices

iter_block(iR: int = 20, R: float | None = None, atoms: sisl.typing.AtomsIndex = None, method: str = 'rand') Iterator[tuple[ndarray, ndarray]][source]

Iterator for performance critical loops

NOTE: This requires that R has been set correctly as the maximum interaction range.

I.e. the loop would look like this:

>>> for ias, idxs in self.iter_block():
...    for ia in ias:
...        idx_a = self.close(ia, R = R, idx = idxs)

This iterator is intended for systems with more than 1000 atoms.

Remark that the iterator used is non-deterministic, i.e. any two iterators need not return the same atoms in any way.

Parameters:
  • iR – the number of R ranges taken into account when doing the iterator

  • R – enables overwriting the local R quantity. Defaults to self.maxR() + 0.001

  • atoms – enables only effectively looping a subset of the full geometry

  • method ({'rand', 'sphere', 'cube'}) – select the method by which the block iteration is performed. Possible values are:

    rand: a spherical object is constructed with a random center according to the internal atoms sphere: a spherical equispaced shape is constructed and looped cube: a cube shape is constructed and looped

Yields:
iter_block_rand(iR: int = 20, R: float | None = None, atoms: sisl.typing.AtomsIndex = None) Iterator[tuple[ndarray, ndarray]][source]

Perform the random block-iteration by randomly selecting the next center of block

iter_block_shape(shape=None, iR: int = 20, atoms: sisl.typing.AtomsIndex = None) Iterator[tuple[ndarray, ndarray]][source]

Perform the grid block-iteration by looping a grid

iter_orbitals(atoms: sisl.typing.AtomsIndex = None, local: bool = True) Iterator[int, int][source]

Returns an iterator over all atoms and their associated orbitals

>>> for ia, io in self.iter_orbitals():

with ia being the atomic index, io the associated orbital index on atom ia. Note that io will start from 0.

Parameters:
  • atoms – only loop on the given atoms, default to all atoms

  • local – whether the orbital index is the global index, or the local index relative to the atom it resides on.

Yields:
  • ia – atomic index

  • io – orbital index

See also

iter

iterate over atomic indices

iter_species

iterate across indices and atomic species

iter_species(atoms: sisl.typing.AtomsIndex = None) Iterator[int, Atom, int][source]

Iterator over all atoms (or a subset) and species as a tuple in this geometry

>>> for ia, a, idx_species in self.iter_species():
...     isinstance(ia, int) == True
...     isinstance(a, Atom) == True
...     isinstance(idx_species, int) == True

with ia being the atomic index, a the Atom object, idx_species is the index of the specie

Parameters:

atoms – only loop on the given atoms, default to all atoms

See also

iter

iterate over atomic indices

iter_orbitals

iterate across atomic indices and orbital indices

maxR(all: bool = False) float[source]

Maximum orbital range of the atoms

mirror(method, atoms: sisl.typing.AtomsIndex = None, point: Sequence[float] = (0, 0, 0)) Geometry[source]

Mirrors the atomic coordinates about a plane given by its normal vector

This will typically move the atomic coordinates outside of the unit-cell. This method should be used with care.

Parameters:
  • method ({'xy'/'z', ..., 'ab', ..., v}) – mirror the structure about a Cartesian direction (x, y, z), plane (xy, xz, yz) or about user defined vectors (v). A vector may also be specified by 'ab' which is the vector normal to the plane spanned by the first and second lattice vector. or user defined vector (v) which is defining a plane.

  • atoms – only mirror a subset of atoms

  • point – mirror coordinates around the plane that intersects the method vector and this point

Examples

>>> geom = geom.graphene()
>>> out = geom.mirror('x')
>>> out.xyz[:, 0]
[0.  -1.42]
>>> out = geom.mirror('x', point=(1.42/2, 0, 0))
>>> out.xyz[:, 0]
[1.42  0.]
move(*args, **kwargs) Geometry

See translate for details

new.Sile(**kwargs)
new.ase(**kwargs)
new.pymatgen(**kwargs)
o2a(orbitals: sisl.typing.OrbitalsIndex, unique: bool = False) ndarray[source]

Atomic index corresponding to the orbital indicies.

Note that this will preserve the super-cell offsets.

Parameters:
  • orbitals – List of orbital indices to return the atoms for

  • unique – If True only return the unique atoms.

o2isc(orbitals: sisl.typing.OrbitalsIndex) ndarray[source]

Returns the super-cell index for a specific orbital.

Returns a vector of 3 numbers with integers.

o2sc(orbitals: sisl.typing.OrbitalsIndex) ndarray[source]

Returns the super-cell offset for a specific orbital.

o2transpose(orb1: sisl.typing.OrbitalsIndex, orb2: sisl.typing.OrbitalsIndex | None = None) tuple[ndarray, ndarray][source]

Transposes connections from orb1 to orb2 such that supercell connections are transposed

When handling supercell indices it is useful to get the transposed connection. I.e. if you have a connection from site i (in unit cell indices) to site J (in supercell indices) it may be useful to get the equivalent supercell connection such for site j (in unit cell indices) to site I (in supercell indices) such that they correspond to the transposed coupling.

Note that since this transposes couplings the indices returned are always expanded to the full length if either of the inputs are a single index.

Examples

>>> gr = geom.graphene() # one orbital per site
>>> atoms = gr.close(0, 1.5)
>>> atoms
array([0, 1, 5, 9], dtype=int32)
>>> gr.o2transpose(0, atoms)
(array([0, 1, 1, 1], dtype=int32), array([ 0,  0, 14, 10], dtype=int32))
Parameters:
  • orb1 – orbital indices must have same length as orb2 or length 1

  • orb2 – orbital indices must have same length as orb1 or length 1. If not present then only orb1 will be returned in transposed indices.

Returns:

  • orb2 (array_like) – transposed indices for orb2 (only returned if orb2 is not None)

  • orb1 (array_like) – transposed indices for orb1

oRij(orbitals1: sisl.typing.OrbitalsIndex, orbitals2: sisl.typing.OrbitalsIndex) ndarray[source]

Vector between orbital orbitals1 and orbitals2, orbitals can be in super-cell indices

Returns the vector between two orbitals:

\[\mathbf r^{ij} = \mathbf r^j - \mathbf r^i\]
Parameters:
  • orbitals1 – orbital index of first orbital

  • orbitals2 – orbital indices

optimize_nsc(axes: sisl.typing.CellAxes | None = None, R: float | None = None) ndarray[source]

Optimize the number of supercell connections based on self.maxR()

After this routine the number of supercells may not necessarily be the same.

This is an in-place operation.

Deprecated method!

Parameters:
  • axes – only optimize the specified axes (default to all)

  • R – the maximum connection radius for each atom

orij(orbitals1: sisl.typing.OrbitalsIndex, orbitals2: sisl.typing.OrbitalsIndex) ndarray[source]

Distance between orbital orbitals1 and orbitals2, orbitals can be in super-cell indices

Returns the distance between two orbitals:

\[r^{ij} = |\mathbf r^j - \mathbf r^i|\]
Parameters:
  • orbitals1 – orbital index of first orbital

  • orbitals2 – orbital indices

osc2uc(orbitals: sisl.typing.OrbitalsIndex, unique: bool = False) ndarray[source]

Orbitals from supercell indices to unit-cell indices, possibly removing dublicates

Parameters:
  • orbitals – the orbital supercell indices to be converted to unit-cell indices

  • unique – If True the returned indices are unique and sorted.

ouc2sc(orbitals: sisl.typing.OrbitalsIndex, unique: bool = False) ndarray[source]

Orbitals from unit-cell indices to supercell indices, possibly removing dublicates

Parameters:
  • orbitals – the orbital unit-cell indices to be converted to supercell indices

  • unique – If True the returned indices are unique and sorted.

overlap(other: GeometryLikeType, atol: float = 0.1, offset: Sequence[float] = (0.0, 0.0, 0.0), offset_other: Sequence[float] = (0.0, 0.0, 0.0)) tuple[ndarray, ndarray][source]

Calculate the overlapping indices between two geometries

Find equivalent atoms (in the primary unit-cell only) in two geometries. This routine finds which atoms have the same atomic positions in self and other.

Note that this will return duplicate overlapping atoms if one atoms lies within eps of more than 1 atom in other.

Parameters:
  • other – Geometry to compare with self

  • atol – atoms within this distance will be considered equivalent

  • offset – offset for self.xyz before comparing

  • offset_other – offset for other.xyz before comparing

Examples

>>> gr22 = sisl.geom.graphene().tile(2, 0).tile(2, 1)
>>> gr44 = gr22.tile(2, 0).tile(2, 1)
>>> offset = np.array([0.2, 0.4, 0.4])
>>> gr22 = gr22.translate(offset)
>>> idx = np.arange(len(gr22))
>>> np.random.shuffle(idx)
>>> gr22 = gr22.sub(idx)
>>> idx22, idx44 = gr22.overlap(gr44, offset=-offset)
>>> assert idx22 == np.arange(len(gr22))
>>> assert idx44 == idx
Returns:

  • idx_self (numpy.ndarray of int) – indices in self that are equivalent with idx_other

  • idx_other (numpy.ndarray of int) – indices in other that are equivalent with idx_self

plot.geometry(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

Builds a GeometryPlot by setting the value of “geometry” to the current object.

Plots a geometry structure, with plentiful of customization options.

Parameters:
  • 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.

prepend(other: sisl.typing.LatticeOrGeometryLike, axis: sisl.typing.CellAxis, offset: Literal['none', 'min'] | sisl.typing.Coord = 'none') Geometry

Prepend two structures along axis

This will automatically add the geometry.cell[axis,:] to all atomic coordiates in the other structure before appending.

The basic algorithm is this:

>>> oxa = other.xyz
>>> geometry.xyz = np.append(oxa, geometry.xyz + other.cell[axis,:][None,:])
>>> geometry.cell[axis,:] += other.cell[axis,:]

NOTE: The cell prepended is only in the axis that is prependend, which means that the other cell directions need not conform.

Parameters:
  • other – Other geometry class which needs to be prepended If a Lattice only the super cell will be extended

  • axis – Cell direction to which the other geometry should be prepended

  • offset ({'none', 'min', (3,)}) – By default appending two structures will simply use the coordinates, as is. With ‘min’, the routine will shift both the structures along the cell axis of other such that they coincide at the first atom, lastly one may use a specified offset to manually select how geometry is displaced. NOTE: That other.cell[axis, :] will be added to offset if other is a geometry.

See also

Geometry.add

add geometries

Geometry.append

appending geometries

Geometry.attach

attach a geometry

Geometry.insert

insert a geometry

static read(sile: sisl.typing.SileLike, *args, **kwargs) Geometry[source]

Reads geometry from the Sile using Sile.read_geometry

Parameters:

sile – a Sile object which will be used to read the geometry if it is a string it will create a new sile using get_sile.

See also

write

writes a Geometry to a given Sile/file

reduce() None[source]

Remove all atoms not currently used in the self.atoms object

Notes

This is an in-place operation.

remove(atoms: sisl.typing.AtomsIndex) Geometry

Remove atoms from the geometry.

Indices passed MUST be unique.

Negative indices are wrapped and thus works.

Parameters:

atoms – indices/boolean of all atoms to be removed

See also

Geometry.sub

the negative of this routine, i.e. retain a subset of atoms

remove_orbital(atoms: sisl.typing.AtomsIndex, orbitals: sisl.typing.OrbitalsIndex) Geometry[source]

Remove a subset of orbitals on atoms according to orbitals

For more detailed examples, please see the equivalent (but opposite) method sub_orbital.

Parameters:
  • atoms (array_like of int or Atom) – indices of atoms or Atom that will be reduced in size according to orbitals

  • orbitals (array_like of int or Orbital) – indices of the orbitals on atoms that are removed from the geometry.

See also

sub_orbital

retaining a set of orbitals (see here for examples)

reorder() None[source]

Reorders atoms according to first occurence in the geometry

The atoms gets reordered according to their placement in the geometry. For instance, if the first atom is the 2nd species in the geometry. Then this routine will swap the 2nd and 1st species in the self.atoms object.

Notes

This is an in-place operation.

repeat(reps: int, axis: sisl.typing.CellAxis) Geometry

Create a repeated geometry

The atomic indices are NOT retained from the base structure.

The expansion of the atoms are basically performed using this algorithm:

>>> ja = 0
>>> for ia in range(geometry.na):
...     for id,r in args:
...        for i in range(r):
...           ja = ia + cell[id,:] * i

For geometries with a single atom this routine returns the same as tile.

Tiling and repeating a geometry will result in the same geometry. The only difference between the two is the final ordering of the atoms.

Parameters:
  • reps – number of repetitions

  • axis – direction of repetition, 0, 1, 2 according to the cell-direction

Examples

>>> geom = Geometry([[0, 0, 0], [0.5, 0, 0]], lattice=1)
>>> g = geom.repeat(2,axis=0)
>>> print(g.xyz) 
[[0.   0.   0. ]
 [1.   0.   0. ]
 [0.5  0.   0. ]
 [1.5  0.   0. ]]
>>> g = geom.repeat(2,0).repeat(2,1)
>>> print(g.xyz) 
[[0.   0.   0. ]
 [0.   1.   0. ]
 [1.   0.   0. ]
 [1.   1.   0. ]
 [0.5  0.   0. ]
 [0.5  1.   0. ]
 [1.5  0.   0. ]
 [1.5  1.   0. ]]

In functional form: >>> repeat(geom, 2, axis=0)

See also

Geometry.tile

equivalent geometry as repeat but different ordering of final structure

Geometry.unrepeat

opposite method of this

replace(atoms: sisl.typing.AtomsIndex, other: sisl.typing.GeometryLike, offset: Sequence[float] = (0.0, 0.0, 0.0)) Geometry[source]

Create a new geometry from self and replace atoms with other

Parameters:
  • atoms – atoms in self to be removed and replaced by other other will be placed in the geometry at the lowest index of atoms

  • other – the other Geometry to insert instead, the unit-cell will not be used.

  • offset – the offset for other when adding its coordinates, default to no offset

reverse(atoms: sisl.typing.AtomsIndex = None) Geometry[source]

Returns a reversed geometry

Also enables reversing a subset of the atoms.

Parameters:

atoms – only reverse the given atomic indices, if not specified, all atoms will be reversed

rij(ia: sisl.typing.AtomsIndex, ja: sisl.typing.AtomsIndex) ndarray[source]

Distance between atom ia and ja, atoms can be in super-cell indices

Returns the distance between two atoms:

\[r^{IJ} = |\mathbf r^J - \mathbf r^I|\]
Parameters:
  • ia – atomic index of first atom

  • ja – atomic indices

rotate(angle: float, v: str | int | sisl.typing.Coord, origin: int | sisl.typing.Coord | None = None, atoms: sisl.typing.AtomsIndex = None, rad: bool = False, what: Literal['xyz', 'abc', 'abc+xyz', 'x', 'a', Ellipsis] | None = None) Geometry

Rotate geometry around vector and return a new geometry

Per default will the entire geometry be rotated, such that everything is aligned as before rotation.

However, by supplying what = 'abc|xyz' one can designate which part of the geometry that will be rotated.

Parameters:
  • angle – the angle in degrees to rotate the geometry. Set the rad argument to use radians.

  • v – the normal vector to the rotated plane, i.e. [1, 0, 0] will rotate around the \(yz\) plane. If a str it refers to the Cartesian direction (xyz), or the lattice vectors (abc). Providing several is the combined direction.

  • origin – the origin of rotation. Anything but [0, 0, 0] is equivalent to a geometry.translate(-origin).rotate(…).translate(origin). If this is an int it corresponds to the atomic index.

  • atoms – only rotate the given atomic indices, if not specified, all atoms will be rotated.

  • rad – if True the angle is provided in radians (rather than degrees)

  • what ({'xyz', 'abc', 'abc+xyz', <or combinations of "xyzabc">}) – which coordinate subject should be rotated, if any of abc is in this string the corresponding cell vector will be rotated if any of xyz is in this string the corresponding coordinates will be rotated If atoms is None, this defaults to “abc+xyz”, otherwise it defaults to “xyz”. See Examples.

Examples

rotate coordinates around the \(x\)-axis >>> geom_x45 = geom.rotate(45, [1, 0, 0])

rotate around the (1, 1, 0) direction but project the rotation onto the \(x\) axis >>> geom_xy_x = geom.rotate(45, “xy”, what=’x’)

See also

Quaternion

class to rotate

Lattice.rotate

rotation for a Lattice object

rotate_miller(m, v) Geometry[source]

Align Miller direction along v

Rotate geometry and cell such that the Miller direction points along the Cartesian vector v.

sc2uc(atoms: sisl.typing.AtomsIndex, unique: bool = False) ndarray

Returns atoms from supercell indices to unit-cell indices, possibly removing dublicates

Parameters:
  • atoms – the atomic supercell indices to be converted to unit-cell indices

  • unique – If True the returned indices are unique and sorted.

sc_index(*args, **kwargs) int | Sequence[int]

Call local Lattice.sc_index function

scale(scale: sisl.typing.CoordOrScalar, what: Literal['abc', 'xyz'] = 'abc', scale_basis: bool = True) Geometry

Scale coordinates and unit-cell to get a new geometry with proper scaling

Parameters:
  • scale – the scale factor for the new geometry (lattice vectors, coordinates and the atomic radii are scaled).

  • what ({"abc", "xyz"}) –

    abc

    Is applied on the corresponding lattice vector and the fractional coordinates.

    xyz

    Is applied only to the atomic coordinates.

    If three different scale factors are provided, each will correspond to the Cartesian direction/lattice vector.

  • scale_basis – if true, the atoms basis-sets will be also be scaled. The scaling of the basis-sets will be done based on the largest scaling factor.

set_lattice(lattice: sisl.typing.LatticeLike)

Overwrites the local lattice.

set_nsc(*args, **kwargs)

Set the number of super-cells in the Lattice object

See set_nsc for allowed parameters.

See also

Lattice.set_nsc

the underlying called method

set_supercell(lattice: sisl.typing.LatticeLike)

Overwrites the local lattice.

sort(**kwargs) Geometry | tuple[Geometry, list[list[int]]]

Sort atoms in a nested fashion according to various criteria

There are many ways to sort a Geometry. - by Cartesian coordinates, axis - by lattice vectors, lattice - by user defined vectors, vector - by grouping atoms, group - by a user defined function, func - by a user defined function using internal sorting algorithm, func_sort

  • a combination of the above in arbitrary order

Additionally one may sort ascending or descending.

This method allows nested sorting based on keyword arguments.

Parameters:
  • atoms (AtomsIndex, optional) – only perform sorting algorithm for subset of atoms. This is NOT a positional dependent argument. All sorting algorithms will _only_ be performed on these atoms. Default, all atoms will be sorted.

  • ret_atoms (bool, optional) – return a list of list for the groups of atoms that have been sorted.

  • axis (int or tuple of int, optional) – sort coordinates according to Cartesian coordinates, if a tuple of ints is passed it will be equivalent to sort(axis0=axis[0], axis1=axis[1]). This behaves differently than numpy.lexsort!

  • lattice (int or tuple of int, optional) – sort coordinates according to lattice vectors, if a tuple of ints is passed it will be equivalent to sort(lattice0=lattice[0], lattice1=lattice[1]). Note that before sorting we multiply the fractional coordinates by the length of the lattice vector. This ensures that atol is meaningful for both axis and lattice since they will be on the same order of magnitude. This behaves differently than numpy.lexsort!

  • vector (Coord, optional) – sort along a user defined vector, similar to lattice but with a user defined direction. Note that lattice sorting and vector sorting are only equivalent when the lattice vector is orthogonal to the other lattice vectors.

  • group ({'Z', 'symbol', 'tag', 'species'} or (str, ), optional) – group together a set of atoms by various means. group may be one of the listed strings. For 'Z' atoms will be grouped in atomic number For 'symbol' atoms will be grouped by their atomic symbol. For 'tag' atoms will be grouped by their atomic tag. For 'species' atoms will be sorted according to their species index. If a tuple/list is passed the first item is described. All subsequent items are a list of groups, where each group comprises elements that should be sorted on an equal footing. If one of the groups is None, that group will be replaced with all non-mentioned elements. See examples.

  • func (callable or list-like of callable, optional) – pass a sorting function which should have an interface like func(geometry, atoms, **kwargs). The first argument is the geometry to sort. The 2nd argument is a list of indices in the current group of sorted atoms. And **kwargs are any optional arguments currently collected, i.e. ascend, atol etc. The function should return either a list of atoms, or a list of list of atoms (in which case the atomic indices have been split into several groups that will be sorted individually for subsequent sorting methods). In either case the returned indices must never hold any other indices but the ones passed as atoms. If a list/tuple of functions, they will be processed in that order.

  • func_sort (callable or list-like of callable, optional) – pass a function returning a 1D array corresponding to all atoms in the geometry. The interface should simply be: func(geometry). Those values will be passed down to the internal sorting algorithm. To be compatible with atol the returned values from func_sort should be on the scale of coordinates (in Ang).

  • ascend, descend (bool, optional) – control ascending or descending sorting for all subsequent sorting methods. Default ascend=True.

  • atol (float, optional) – absolute tolerance when sorting numerical arrays for subsequent sorting methods. When a selection of sorted coordinates are grouped via atol, we ensure such a group does not alter its indices. I.e. the group is always ascending indices. Note this may have unwanted side-effects if atol is very large compared to the difference between atomic coordinates. Default 1e-9.

Notes

The order of arguments is also the sorting order. sort(axis=0, lattice=0) is different from sort(lattice=0, axis=0)

All arguments may be suffixed with integers. This allows multiple keyword arguments to control sorting algorithms in different order. It also allows changing of sorting settings between different calls. Note that the integers have no relevance to the order of execution! See examples.

Returns:

  • geometry (Geometry) – sorted geometry

  • index (list of list of indices) – indices that would sort the original structure to the output, only returned if ret_atoms is True

Examples

>>> geom = sisl.geom.bilayer(top_atoms=sisl.Atom[5, 7], bottom_atoms=sisl.Atom(6))
>>> geom = geom.tile(5, 0).repeat(7, 1)

Sort according to \(x\) coordinate

>>> geom.sort(axis=0)

Sort according to \(z\), then \(x\) for each group created from first sort

>>> geom.sort(axis=(2, 0))

Sort according to \(z\), then first lattice vector

>>> geom.sort(axis=2, lattice=0)

Sort according to \(z\) (ascending), then first lattice vector (descending)

>>> geom.sort(axis=2, ascend=False, lattice=0)

Sort according to \(z\) (descending), then first lattice vector (ascending) Note how integer suffixes has no importance.

>>> geom.sort(ascend1=False, axis=2, ascend0=True, lattice=0)

Sort only atoms range(1, 5) first by \(z\), then by first lattice vector

>>> geom.sort(axis=2, lattice=0, atoms=np.arange(1, 5))

Sort two groups of atoms [range(1, 5), range(5, 10)] (individually) by \(z\) coordinate

>>> geom.sort(axis=2, atoms=[np.arange(1, 5), np.arange(5, 10)])

The returned sorting indices may be used for manual sorting. Note however, that this requires one to perform a sorting for all atoms. In such a case the following sortings are equal.

>>> geom0, atoms0 = geom.sort(axis=2, lattice=0, ret_atoms=True)
>>> _, atoms1 = geom.sort(axis=2, ret_atoms=True)
>>> geom1, atoms1 = geom.sort(lattice=0, atoms=atoms1, ret_atoms=True)
>>> geom2 = geom.sub(np.concatenate(atoms0))
>>> geom3 = geom.sub(np.concatenate(atoms1))
>>> assert geom0 == geom1
>>> assert geom0 == geom2
>>> assert geom0 == geom3

Default sorting is equivalent to axis=(0, 1, 2)

>>> assert geom.sort() == geom.sort(axis=(0, 1, 2))

Sort along a user defined vector [2.2, 1., 0.]

>>> geom.sort(vector=[2.2, 1., 0.])

Integer specification has no influence on the order of operations. It is _always_ the keyword argument order that determines the operation.

>>> assert geom.sort(axis2=1, axis0=0, axis1=2) == geom.sort(axis=(1, 0, 2))

Sort by atomic numbers

>>> geom.sort(group='Z') # 5, 6, 7

One may group several elements together on an equal footing (None means all non-mentioned elements) The order of the groups are important (the first two are _not_ equal, the last three _are_ equal)

>>> geom.sort(group=('symbol', 'C'), axis=2) # C will be sorted along z
>>> geom.sort(axis=1, atoms='C', axis1=2) # all along y, then C sorted along z
>>> geom.sort(group=('symbol', 'C', None)) # C, [B, N]
>>> geom.sort(group=('symbol', None, 'C')) # [B, N], C
>>> geom.sort(group=('symbol', ['N', 'B'], 'C')) # [B, N], C (B and N unaltered order)
>>> geom.sort(group=('symbol', ['B', 'N'], 'C')) # [B, N], C (B and N unaltered order)

A group based sorting can use anything that can be fetched from the Atom object, sort first according to mass, then for all with the same mass, sort according to atomic tag:

>>> geom.sort(group0='mass', group1='tag')

A too high atol may have unexpected side-effects. This is because of the way the sorting algorithm splits the sections for nested sorting. So for coordinates with a continuous displacement the sorting may break and group a large range into 1 group. Consider the following array to be split in groups while sorting.

An example would be a linear chain with a middle point with a much shorter distance.

>>> x = np.arange(5) * 0.1
>>> x[3:] -= 0.095
y = z = np.zeros(5)
geom = si.Geometry(np.stack((x, y, z), axis=1))
>>> geom.xyz[:, 0]
[0.    0.1   0.2   0.205 0.305]

In this case a high tolerance (atol>0.005) would group atoms 2 and 3 together

>>> geom.sort(atol=0.01, axis=0, ret_atoms=True)[1]
[[0], [1], [2, 3], [4]]

However, a very low tolerance will not find these two as atoms close to each other.

>>> geom.sort(atol=0.001, axis=0, ret_atoms=True)[1]
[[0], [1], [2], [3], [4]]
sparserij(dtype=np.float64, na_iR: int = 1000, method: str = 'rand')[source]

Return the sparse matrix with all distances in the matrix The sparse matrix will only be defined for the elements which have orbitals overlapping with other atoms.

Parameters:
  • dtype (numpy.dtype, numpy.float64) – the data-type of the sparse matrix

  • na_iR – number of atoms within the sphere for speeding up the iter_block loop.

  • method – see iter_block for details

Returns:

SparseAtom – sparse matrix with all rij elements

See also

iter_block

the method for looping the atoms

distance

create a list of distances

sub(atoms: sisl.typing.AtomsIndex) Geometry

Create a new Geometry with a subset of this Geometry

Indices passed MUST be unique.

Negative indices are wrapped and thus works.

Parameters:

atoms – indices/boolean of all atoms to be removed

See also

Lattice.fit

update the supercell according to a reference supercell

Geometry.remove

the negative of this routine, i.e. remove a subset of atoms

sub_orbital(atoms: sisl.typing.AtomsIndex, orbitals: sisl.typing.OrbitalsIndex) Geometry[source]

Retain only a subset of the orbitals on atoms according to orbitals

This allows one to retain only a given subset of geometry.

Parameters:
  • atoms (array_like of int or Atom) – indices of atoms or Atom that will be reduced in size according to orbitals

  • orbitals (array_like of int or Orbital) – indices of the orbitals on atoms that are retained in the geometry, the list of orbitals will be sorted.

Notes

Future implementations may allow one to re-arange orbitals using this method.

When using this method the internal species list will be populated by another specie that is named after the orbitals removed. This is to distinguish different atoms.

Examples

>>> # a Carbon atom with 2 orbitals
>>> C = sisl.Atom('C', [1., 2.])
>>> # an oxygen atom with 3 orbitals
>>> O = sisl.Atom('O', [1., 2., 2.4])
>>> geometry = sisl.Geometry([[0] * 3, [1] * 3]], 2, [C, O])

Now geometry is a geometry with 2 different species and 6 atoms (3 of each). They are ordered [C, O, C, O, C, O]. In the following we will note species that are different from the original by a ' in the list.

Retain 2nd orbital on the 2nd atom: [C, O', C, O, C, O]

>>> new_geom = geometry.sub_orbital(1, 1)

Retain 2nd orbital on 1st and 2nd atom: [C', O', C, O, C, O]

>>> new_geom = geometry.sub_orbital([0, 1], 1)

Retain 2nd orbital on the 1st atom and 3rd orbital on 4th atom: [C', O, C, O', C, O]

>>> new_geom = geometry.sub_orbital(0, 1).sub_orbital(3, 2)

Retain 2nd orbital on all atoms equivalent to the first atom: [C', O, C', O, C', O]

>>> new_geom = geometry.sub_orbital(obj.geometry.atoms[0], 1)

Retain 1st orbital on 1st atom, and 2nd orbital on 3rd and 5th atom: [C', O, C'', O, C'', O]

>>> new_geom = geometry.sub_orbital(0, 0).sub_orbital([2, 4], 1)

See also

remove_orbital

removing a set of orbitals (opposite of this)

swap(atoms1: sisl.typing.AtomsIndex, atoms2: sisl.typing.AtomsIndex) Geometry

Swap a set of atoms in the geometry and return a new one

This can be used to reorder elements of a geometry.

Parameters:
  • atoms1 – the first list of atomic coordinates

  • atoms2 – the second list of atomic coordinates

swapaxes(axes1: sisl.typing.AnyAxes, axes2: sisl.typing.AnyAxes, what: Literal['abc', 'xyz', 'abc+xyz'] = 'abc') Geometry

Swap the axes components by either lattice vectors (only cell), or Cartesian coordinates

See Lattice.swapaxes for details.

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 by xyz. If str, then what is not used.

  • axes2 – the new axis indices, same as axes1 old axis indices (or labels)

  • what – what to swap, lattice vectors (abc) or Cartesian components (xyz), or both. Neglected for integer axes arguments.

See also

Lattice.swapaxes

Examples

Only swap lattice vectors

>>> g_ba = g.swapaxes(0, 1)
>>> assert np.allclose(g.xyz, g_ba.xyz)

Only swap Cartesian coordinates

>>> g_ba = g.swapaxes(0, 1, "xyz")
>>> assert np.allclose(g.xyz[:, [1, 0, 2]], g_ba.xyz)

Consecutive swappings (what will be neglected if provided):

  1. abc, xyz -> bac, xyz

  2. bac, xyz -> bca, xyz

  3. bac, xyz -> bca, zyx

>>> g_s = g.swapaxes("abx", "bcz")
>>> assert np.allclose(g.xyz[:, [2, 1, 0]], g_s.xyz)
>>> assert np.allclose(g.cell[[1, 2, 0]][:, [2, 1, 0]], g_s.cell)
tile(reps: int, axis: sisl.typing.CellAxis) Geometry

Tile the geometry to create a bigger one

The atomic indices are retained for the base structure.

Tiling and repeating a geometry will result in the same geometry. The only difference between the two is the final ordering of the atoms.

Parameters:
  • reps – number of tiles (repetitions)

  • axis – direction of tiling, 0, 1, 2 according to the cell-direction

Examples

>>> geom = Geometry([[0, 0, 0], [0.5, 0, 0]], lattice=1.)
>>> g = geom.tile(2,axis=0)
>>> print(g.xyz) 
[[0.   0.   0. ]
 [0.5  0.   0. ]
 [1.   0.   0. ]
 [1.5  0.   0. ]]
>>> g = tile(geom, 2,0).tile(2,axis=1)
>>> print(g.xyz) 
[[0.   0.   0. ]
 [0.5  0.   0. ]
 [1.   0.   0. ]
 [1.5  0.   0. ]
 [0.   1.   0. ]
 [0.5  1.   0. ]
 [1.   1.   0. ]
 [1.5  1.   0. ]]

In functional form: >>> tile(geom, 2, axis=0)

See also

Geometry.repeat

equivalent but different ordering of final structure

Geometry.untile

opposite method of this

to.Path(*args, **kwargs) None

Writes the geometry to a sile with any optional arguments.

Examples

>>> geom = si.geom.graphene()
>>> geom.to("hello.xyz")
>>> geom.to(pathlib.Path("hello.xyz"))
to.Sile(*args, **kwargs) None

Writes the geometry to a sile with any optional arguments.

Examples

>>> geom = si.geom.graphene()
>>> geom.to("hello.xyz")
>>> geom.to(pathlib.Path("hello.xyz"))
to.ase(**kwargs) Atoms

Conversion of Geometry to an ase.Atoms object

to.dataframe(*args, **kwargs) DataFrame

Convert the geometry to a pandas.DataFrame with values stored in columns

to.pymatgen(**kwargs) pymatgen.core.Molecule | pymatgen.core.Structure

Conversion of Geometry to a pymatgen object.

Depending on the periodicity, it can be Molecule or Structure.

to.str(*args, **kwargs) None

Writes the geometry to a sile with any optional arguments.

Examples

>>> geom = si.geom.graphene()
>>> geom.to("hello.xyz")
>>> geom.to(pathlib.Path("hello.xyz"))
translate(v: sisl.typing.CoordOrScalar, atoms: sisl.typing.AtomsIndex = None) Geometry

Translates the geometry by v

move is a shorthand for this function.

One can translate a subset of the atoms by supplying atoms.

Returns a copy of the structure translated by v.

Parameters:
  • v – the value or vector to displace all atomic coordinates It should just be broad-castable with the geometry’s coordinates.

  • atoms – only displace the given atomic indices, if not specified, all atoms will be displaced

translate2uc(atoms: sisl.typing.AtomsIndex = None, axes: int | bool | Sequence[int] | None = None) Geometry[source]

Translates atoms in the geometry into the unit cell

One can translate a subset of the atoms or axes by appropriate arguments.

Warning

When coordinates are lying on one of the edges, they may move to the other side of the unit-cell due to small rounding errors. In such situations you are encouraged to shift all coordinates by a small amount to remove numerical errors, in the following case we have atomic coordinates lying close to the lower side of each lattice vector.

>>> geometry.translate(1e-8).translate2uc().translate(-1e-8)

Notes

By default only the periodic axes (self.pbc) will be translated to the UC. If translation is required for all axes, supply them directly.

Parameters:
  • atoms – only translate the given atomic indices, if not specified, all atoms will be translated

  • axes – only translate certain lattice directions, None specifies only the directions with supercells, True specifies all directions.

uc2sc(atoms: sisl.typing.AtomsIndex, unique: bool = False) ndarray

Returns atom from unit-cell indices to supercell indices, possibly removing dublicates

Parameters:
  • atoms – the atomic unit-cell indices to be converted to supercell indices

  • unique – If True the returned indices are unique and sorted.

unrepeat(reps: int, axis: sisl.typing.CellAxis, *args, **kwargs) Geometry

Unrepeats the geometry similarly as untile

This is the opposite of Geometry.repeat.

Please see Geometry.untile for argument details. This algorithm first re-arranges atoms as though they had been tiled, then subsequently calls untile.

See also

Geometry.repeat

opposite method of this

untile(reps: int, axis: sisl.typing.CellAxis, segment: int = 0, rtol: float = 1e-4, atol: float = 1e-4) Geometry

A subset of atoms from the geometry by cutting the geometry into reps parts along the direction axis.

This will effectively change the unit-cell in the axis as-well as removing geometry.na/reps atoms. It requires that geometry.na % reps == 0.

REMARK: You need to ensure that all atoms within the first cut out region are within the primary unit-cell.

Doing geom.untile(2, 1).tile(2, 1), could for symmetric setups, be equivalent to a no-op operation. A UserWarning will be issued if this is not the case.

This method is the reverse of tile.

Parameters:
  • reps – number of times the structure will be cut (untiled)

  • axis – the axis that will be cut

  • segment – returns the i’th segment of the untiled structure Currently the atomic coordinates are not translated, this may change in the future.

  • rtol – directly passed to numpy.allclose

  • atol – directly passed to numpy.allclose

Examples

>>> g = sisl.geom.graphene()
>>> gxyz = g.tile(4, 0).tile(3, 1).tile(2, 2)
>>> G = gxyz.untile(2, 2).untile(3, 1).untile(4, 0)
>>> np.allclose(g.xyz, G.xyz)
True

In functional form: >>> untile(geom, 2, axis=0)

See also

Geometry.tile

opposite method of this

Geometry.repeat

equivalent geometry as tile but different ordering of final structure

within(shapes, atoms: sisl.typing.AtomsIndex = None, atoms_xyz=None, ret_xyz: bool = False, ret_rij: bool = False, ret_isc: bool = False)[source]

Indices of atoms in the entire supercell within a given shape from a given coordinate

This heavily relies on the within_sc method.

Note that if a connection is made in a neighboring super-cell then the atomic index is shifted by the super-cell index times number of atoms. This allows one to decipher super-cell atoms from unit-cell atoms.

Parameters:
  • shapes (Shape, list of Shape)

  • atoms – List of indices for atoms that are to be considered

  • atoms_xyz (array_like, optional) – The atomic coordinates of the equivalent atoms variable (atoms must also be passed)

  • ret_xyz – If true this method will return the coordinates for each of the couplings.

  • ret_rij – If true this method will return the distances from the xyz_ia for each of the couplings.

  • ret_isc – If true this method will return the supercell offsets for each of the couplings.

Returns:

  • index – indices of atoms (in supercell indices) within the shape

  • xyz – atomic coordinates of the indexed atoms (only for true ret_xyz)

  • rij – distance of the indexed atoms to the center of the shape (only for true ret_rij)

  • isc – supercell indices of the couplings (only for true ret_isc)

within_inf(lattice: Lattice, periodic: Sequence[bool] | sisl.typing.CellAxes | None = None, atol: float = 1e-5, origin: Sequence[float] = (0.0, 0.0, 0.0)) tuple[ndarray, ndarray, ndarray][source]

Find all atoms within a provided supercell

Note this function is rather different from close and within. Specifically this routine is returning all indices for the infinite periodic system. The default periodic directions are self.pbc, unless periodic is provided.

Atomic coordinates lying on the boundary of the supercell will be duplicated on the neighboring supercell images. Thus performing geom.within_inf(geom.lattice) may result in more atoms than in the structure.

Notes

The name of this function may change. Currently it should only be used internally in sisl.

Parameters:
  • lattice (LatticeLike) – the supercell in which this geometry should be expanded into.

  • periodic – explicitly define the periodic directions, by default the periodic directions are only where self.pbc.

  • atol – length tolerance for the coordinates to be on a duplicate site (in Ang). This allows atoms within atol of the cell boundaries to be taken as inside the cell.

  • origin – origin that is the basis for comparison, default to 0.

Returns:

  • numpy.ndarray – unit-cell atomic indices which are inside the lattice cell

  • numpy.ndarray – atomic coordinates for the ia atoms (including supercell offsets)

  • numpy.ndarray – integer supercell offsets for ia atoms

within_sc(shapes, isc=None, atoms: sisl.typing.AtomsIndex = None, atoms_xyz=None, ret_xyz: bool = False, ret_rij: bool = False)[source]

Indices of atoms in a given supercell within a given shape from a given coordinate

This returns a set of atomic indices which are within a sphere of radius R.

If R is a tuple/list/array it will return the indices: in the ranges:

>>> ( x <= R[0] , R[0] < x <= R[1], R[1] < x <= R[2] )
Parameters:
  • shapes (Shape or list of Shape) – A list of increasing shapes that define the extend of the geometric volume that is searched. It is vital that:

    shapes[0] in shapes[1] in shapes[2] ...
    
  • isc (array_like, optional) – The super-cell which the coordinates are checked in. Defaults to [0, 0, 0]

  • atoms – List of atoms that will be considered. This can be used to only take out a certain atoms.

  • atoms_xyz (array_like, optional) – The atomic coordinates of the equivalent idx variable (idx must also be passed)

  • ret_xyz – If True this method will return the coordinates for each of the couplings.

  • ret_rij – If True this method will return the distance to the center of the shapes

Returns:

  • index – indices of atoms (in supercell indices) within the shape

  • xyz – atomic coordinates of the indexed atoms (only for true ret_xyz)

  • rij – distance of the indexed atoms to the center of the shape (only for true ret_rij)

write(sile: sisl.typing.SileLike, *args, **kwargs) None

Writes geometry to the sile using sile.write_geometry

Parameters:
  • sile – a Sile object which will be used to write the geometry if it is a string it will create a new sile using get_sile

  • *args, **kwargs – Any other args will be passed directly to the underlying routine

See also

Geometry.read

reads a Geometry from a given Sile/file

write

generic sisl function dispatcher

property atoms: Atoms

Atoms associated with the geometry

property boundary_condition: ndarray
property cell: ndarray

Returns the inherent Lattice.cell

property firsto: npt.NDArray[np.int32]

The first orbital on the corresponding atom

property fxyz: npt.NDArray[np.float64]

Returns geometry coordinates in fractional coordinates

property icell: ndarray

Returns the inherent Lattice.icell

property isc_off: ndarray

Returns the inherent Lattice.isc_off

property lasto: npt.NDArray[np.int32]

The last orbital on the corresponding atom

property length: float

Returns the inherent Lattice.length

property mass: ndarray

The mass of all atoms as an array

property n_s: int

Returns the inherent Lattice.n_s

property na: int

Number of atoms in geometry

property na_s: int

Number of supercell atoms

property names

The named index specifier

new = <TypeDispatcher{obj=<class 'sisl.Geometry'>}>
property no: int

Number of orbitals in unit cell

property no_s: int

Number of supercell orbitals

property nsc: ndarray

Returns the inherent Lattice.nsc

property orbitals: list[Orbital]

List of orbitals per atom

property origin: ndarray

Returns the inherent Lattice.origin

property pbc: ndarray
plot

Plotting functions for the Geometry class.

property q0: float

Total initial charge in this geometry (sum of q0 off all atoms)

property rcell: ndarray

Returns the inherent Lattice.rcell

property sc

[deprecated] Return the lattice object associated with the Lattice.

property sc_off: ndarray

Returns the inherent Lattice.sc_off

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 class

  • dispatchs (dict, optional) – dictionary of dispatch methods

  • obj_getattr (callable, optional) – method with 2 arguments, an obj and the attr 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.

property volume: float

Returns the inherent Lattice.volume