sisl.SparseAtom

class sisl.SparseAtom(geometry, dim=1, dtype=None, nnzpr=None, **kwargs)

Bases: _SparseGeometry

Sparse object with number of rows equal to the total number of atoms in the Geometry

Methods

Rij([dtype])

Create a sparse matrix with vectors between atoms

construct(func[, na_iR, method, eta])

Automatically construct the sparse model based on a function that does the setting up of the elements

copy([dtype])

A copy of this object

create_construct(R, params)

Create a simple function for passing to the construct function.

edges(atoms[, exclude])

Retrieve edges (connections) for all atoms

eliminate_zeros(*args, **kwargs)

Removes all zero elements from the sparse matrix

empty([keep_nnz])

See empty for details

finalize()

Finalizes the model

fromsp(geometry, P, **kwargs)

Create a sparse model from a preset Geometry and a list of sparse matrices

iter_nnz([atoms])

Iterations of the non-zero elements

nonzero([atoms, only_cols])

Indices row and column indices where non-zero elements exists

remove(atoms)

Create a subset of this sparse matrix by removing the atoms corresponding to atoms

repeat(reps, axis)

Create a repeated sparse atom object, equivalent to Geometry.repeat

reset([dim, dtype, nnzpr])

The sparsity pattern has all elements removed and everything is reset.

rij([dtype])

Create a sparse matrix with the distance between atoms

set_nsc(*args, **kwargs)

Reset the number of allowed supercells in the sparse atom

spalign(other)

See align for details

spsame(other)

Compare two sparse objects and check whether they have the same entries.

sub(atoms)

Create a subset of this sparse matrix by only retaining the elements corresponding to the atoms

swap(a, b)

Swaps atoms in the sparse geometry to obtain a new order of atoms

tile(reps, axis)

Create a tiled sparse atom object, equivalent to Geometry.tile

tocsr([dim, isc])

Return a csr_matrix for the specified dimension

translate2uc([atoms, axes])

Translates all primary atoms to the unit cell.

transpose([sort])

Create the transposed sparse geometry by interchanging supercell indices

unrepeat(reps, axis[, segment, sym])

Unrepeats the sparse model into different parts (retaining couplings)

untile(reps, axis[, segment, sym])

Untiles the sparse model into different parts (retaining couplings)

dim

Number of components per element

dkind

Data type of sparse elements (in str)

dtype

Data type of sparse elements

finalized

Whether the contained data is finalized and non-used elements have been removed

geometry

Associated geometry

nnz

Number of non-zero elements

shape

Shape of sparse matrix

Rij(dtype=<class 'numpy.float64'>)[source]

Create a sparse matrix with vectors between atoms

Parameters:

dtype (numpy.dtype, optional) – the data-type of the sparse matrix.

Notes

The returned sparse matrix with vectors are taken from the current sparse pattern. I.e. a subsequent addition of sparse elements will make them inequivalent. It is thus important to only create the sparse vector matrix when the sparse structure is completed.

__init__(geometry, dim=1, dtype=None, nnzpr=None, **kwargs)

Create sparse object with element between orbitals

construct(func, na_iR=1000, method='rand', eta=None)

Automatically construct the sparse model based on a function that does the setting up of the elements

This may be called in two variants.

  1. Pass a function (func), see e.g. create_construct which does the setting up.

  2. Pass a tuple/list in func which consists of two elements, one is R the radii parameters for the corresponding parameters. The second is the parameters corresponding to the R[i] elements. In this second case all atoms must only have one orbital.

Parameters:
  • func (callable or array_like) –

    this function must take 4 arguments. 1. Is this object (self) 2. Is the currently examined atom (ia) 3. Is the currently bounded indices (idxs) 4. Is the currently bounded indices atomic coordinates (idxs_xyz) An example func could be:

    >>> def func(self, ia, atoms, atoms_xyz=None):
    ...     idx = self.geometry.close(ia, R=[0.1, 1.44], atoms=atoms, atoms_xyz=atoms_xyz)
    ...     self[ia, idx[0]] = 0
    ...     self[ia, idx[1]] = -2.7
    

  • na_iR (int, optional) – number of atoms within the sphere for speeding up the iter_block loop.

  • method ({'rand', str}) – method used in Geometry.iter_block, see there for details

  • eta (bool, optional) – whether an ETA will be printed

See also

create_construct

a generic function used to create a generic function which this routine requires

tile

tiling after construct is much faster for very large systems

repeat

repeating after construct is much faster for very large systems

copy(dtype=None)

A copy of this object

Parameters:

dtype (numpy.dtype, optional) – it is possible to convert the data to a different data-type If not specified, it will use self.dtype

create_construct(R, params)

Create a simple function for passing to the construct function.

This is simply to leviate the creation of simplistic functions needed for setting up the sparse elements.

Basically this returns a function:

>>> def func(self, ia, atoms, atoms_xyz=None):
...     idx = self.geometry.close(ia, R=R, atoms=atoms, atoms_xyz=atoms_xyz)
...     for ix, p in zip(idx, params):
...         self[ia, ix] = p

Notes

This function only works for geometry sparse matrices (i.e. one element per atom). If you have more than one element per atom you have to implement the function your-self.

Parameters:
  • R (array_like) – radii parameters for different shells. Must have same length as params or one less. If one less it will be extended with R[0]/100

  • params (array_like) – coupling constants corresponding to the R ranges. params[0, :] are the elements for the all atoms within R[0] of each atom.

See also

construct

routine to create the sparse matrix from a generic function (as returned from create_construct)

property dim

Number of components per element

property dkind

Data type of sparse elements (in str)

property dtype

Data type of sparse elements

edges(atoms, exclude=None)

Retrieve edges (connections) for all atoms

The returned edges are unique and sorted (see numpy.unique) and are returned in supercell indices (i.e. 0 <= edge < self.geometry.na_s).

Parameters:
  • atoms (int or list of int) – the edges are returned only for the given atom

  • exclude (int or list of int or None, optional) – remove edges which are in the exclude list.

See also

SparseCSR.edges

the underlying routine used for extracting the edges

eliminate_zeros(*args, **kwargs)

Removes all zero elements from the sparse matrix

This is an in-place operation.

See also

SparseCSR.eliminate_zeros

method called, see there for parameters

empty(keep_nnz=False)

See empty for details

finalize()

Finalizes the model

Finalizes the model so that all non-used elements are removed. I.e. this simply reduces the memory requirement for the sparse matrix.

Note that adding more elements to the sparse matrix is more time-consuming than for a non-finalized sparse matrix due to the internal data-representation.

property finalized

Whether the contained data is finalized and non-used elements have been removed

classmethod fromsp(geometry, P, **kwargs)

Create a sparse model from a preset Geometry and a list of sparse matrices

The passed sparse matrices are in one of scipy.sparse formats.

Parameters:
  • geometry (Geometry) – geometry to describe the new sparse geometry

  • P (list of scipy.sparse or scipy.sparse) – the new sparse matrices that are to be populated in the sparse matrix

  • **kwargs (optional) – any arguments that are directly passed to the __init__ method of the class.

Returns:

a new sparse matrix that holds the passed geometry and the elements of P

Return type:

SparseGeometry

property geometry

Associated geometry

iter_nnz(atoms=None)[source]

Iterations of the non-zero elements

An iterator on the sparse matrix with, row and column

Examples

>>> for i, j in self.iter_nnz():
...    self[i, j] # is then the non-zero value
Parameters:

atoms (int or array_like) – only loop on the non-zero elements coinciding with the atoms

property nnz

Number of non-zero elements

nonzero(atoms=None, only_cols=False)[source]

Indices row and column indices where non-zero elements exists

Parameters:
  • atoms (int or array_like of int, optional) – only return the tuples for the requested atoms, default is all atoms

  • only_cols (bool, optional) – only return then non-zero columns

See also

SparseCSR.nonzero

the equivalent function call

remove(atoms)

Create a subset of this sparse matrix by removing the atoms corresponding to atoms

Negative indices are wrapped and thus works.

Parameters:

atoms (array_like of int) – indices of removed atoms

See also

Geometry.remove

equivalent to the resulting Geometry from this routine

Geometry.sub

the negative of Geometry.remove

sub

the opposite of remove, i.e. retain a subset of atoms

repeat(reps, axis)[source]

Create a repeated sparse atom object, equivalent to Geometry.repeat

The already existing sparse elements are extrapolated to the new supercell by repeating them in blocks like the coordinates.

Parameters:
  • reps (int) – number of repetitions along cell-vector axis

  • axis (int) – 0, 1, 2 according to the cell-direction

See also

Geometry.repeat

the same ordering as the final geometry

Geometry.tile

a different ordering of the final geometry

tile

a different ordering of the final geometry

reset(dim=None, dtype=<class 'numpy.float64'>, nnzpr=None)

The sparsity pattern has all elements removed and everything is reset.

The object will be the same as if it had been initialized with the same geometry as it were created with.

Parameters:
  • dim (int, optional) – number of dimensions per element, default to the current number of elements per matrix element.

  • dtype (numpy.dtype, optional) – the datatype of the sparse elements

  • nnzpr (int, optional) – number of non-zero elements per row

rij(dtype=<class 'numpy.float64'>)[source]

Create a sparse matrix with the distance between atoms

Parameters:

dtype (numpy.dtype, optional) – the data-type of the sparse matrix.

Notes

The returned sparse matrix with distances are taken from the current sparse pattern. I.e. a subsequent addition of sparse elements will make them inequivalent. It is thus important to only create the sparse distance when the sparse structure is completed.

set_nsc(*args, **kwargs)[source]

Reset the number of allowed supercells in the sparse atom

If one reduces the number of supercells any sparse element that references the supercell will be deleted.

See Lattice.set_nsc for allowed parameters.

See also

Lattice.set_nsc

the underlying called method

property shape

Shape of sparse matrix

spalign(other)

See align for details

spsame(other)

Compare two sparse objects and check whether they have the same entries.

This does not necessarily mean that the elements are the same

sub(atoms)[source]

Create a subset of this sparse matrix by only retaining the elements corresponding to the atoms

Indices passed MUST be unique.

Negative indices are wrapped and thus works.

Parameters:

atoms (array_like of int) – indices of retained atoms

See also

Geometry.remove

the negative of Geometry.sub

Geometry.sub

equivalent to the resulting Geometry from this routine

remove

the negative of sub, i.e. remove a subset of atoms

swap(a, b)

Swaps atoms in the sparse geometry to obtain a new order of atoms

This can be used to reorder elements of a geometry.

Parameters:
  • a (array_like) – the first list of atomic coordinates

  • b (array_like) – the second list of atomic coordinates

tile(reps, axis)[source]

Create a tiled sparse atom object, equivalent to Geometry.tile

The already existing sparse elements are extrapolated to the new supercell by repeating them in blocks like the coordinates.

Notes

Calling this routine will automatically finalize the SparseAtom. This is required to greatly increase performance.

Parameters:
  • reps (int) – number of repetitions along cell-vector axis

  • axis (int) – 0, 1, 2 according to the cell-direction

See also

repeat

a different ordering of the final geometry

untile

opposite of this method

Geometry.tile

the same ordering as the final geometry

Geometry.repeat

a different ordering of the final geometry

tocsr(dim=0, isc=None, **kwargs)

Return a csr_matrix for the specified dimension

Parameters:
  • dim (int, optional) – the dimension in the sparse matrix (for non-orthogonal cases the last dimension is the overlap matrix)

  • isc (int, optional) – the supercell index, or all (if isc=None)

translate2uc(atoms=None, axes=None)

Translates all primary atoms to the unit cell.

With this, the coordinates of the geometry are translated to the unit cell and the supercell connections in the matrix are updated accordingly.

Parameters:
  • atoms (AtomsArgument, optional) – only translate the specified atoms. If not specified, all atoms will be translated.

  • axes (int or array_like or None, optional) – only translate certain lattice directions, None species only the periodic directions

Returns:

A new sparse matrix with the updated connections and a new associated geometry.

Return type:

SparseOrbital or SparseAtom

transpose(sort=True)

Create the transposed sparse geometry by interchanging supercell indices

Sparse geometries are (typically) relying on symmetry in the supercell picture. Thus when one transposes a sparse geometry one should ideally get the same matrix. This is true for the Hamiltonian, density matrix, etc.

This routine transposes all rows and columns such that any interaction between row, r, and column c in a given supercell (i,j,k) will be transposed into row c, column r in the supercell (-i,-j,-k).

Parameters:

sort (bool, optional) – the returned columns for the transposed structure will be sorted if this is true, default

Notes

The components for each sparse element are not changed in this method.

Examples

Force a sparse geometry to be Hermitian:

>>> sp = SparseOrbital(...)
>>> sp = (sp + sp.transpose()) / 2
Returns:

an equivalent sparse geometry with transposed matrix elements

Return type:

object

unrepeat(reps, axis, segment=0, *args, sym=True, **kwargs)

Unrepeats the sparse model into different parts (retaining couplings)

Please see untile for details, the algorithm and arguments are the same however, this is the opposite of repeat.

untile(reps, axis, segment=0, *args, sym=True, **kwargs)[source]

Untiles the sparse model into different parts (retaining couplings)

Recreates a new sparse object with only the cutted atoms in the structure. This will preserve matrix elements in the supercell.

Parameters:
  • reps (int) – number of repetitions the tiling function created (opposite meaning as in untile)

  • axis (int) – which axis to untile along

  • segment (int, optional) – which segment to retain. Generally each segment should be equivalent, however requesting individiual segments can help uncover inconsistencies in the sparse matrix

  • *args – arguments passed directly to Geometry.untile

  • sym (bool, optional) – if True, the algorithm will ensure the returned matrix is symmetrized (i.e. return (M + M.transpose())/2, else return data as is. False should generally only be used for debugging precision of the matrix elements, or if one wishes to check the warnings.

  • **kwargs – keyword arguments passed directly to Geometry.untile

Notes

Untiling structures with nsc == 1 along axis are assumed to have periodic boundary conditions.

When untiling structures with nsc == 1 along axis it is important to untile as much as possible. This is because otherwise the algorithm cannot determine the correct couplings. Therefore to create a geometry of 3 times a unit-cell, one should untile to the unit-cell, and subsequently tile 3 times.

Consider for example a system of 4 atoms, each atom connects to its 2 neighbours. Due to the PBC atom 0 will connect to 1 and 3. Untiling this structure in 2 will group couplings of atoms 0 and 1. As it will only see one coupling to the right it will halve the coupling and use the same coupling to the left, which is clearly wrong.

In the following the latter is the correct way to do it.

>>> SPM.untile(2, 0) != SPM.untile(4, 0).tile(2, 0)
Raises:

ValueError : – in case the matrix elements are not conseuctive when determining the new supercell structure. This may often happen if untiling a matrix too few times, and then untiling it again.

See also

tile

opposite of this method

Geometry.untile

same as this method, see details about parameters here