sisl.physics.DensityMatrix¶
- class sisl.physics.DensityMatrix(geometry, dim=1, dtype=None, nnzpr=None, **kwargs)¶
Bases:
sisl.physics.densitymatrix._densitymatrix
Sparse density matrix object
Assigning or changing elements is as easy as with standard
numpy
assignments:>>> DM = DensityMatrix(...) >>> DM.D[1,2] = 0.1
which assigns 0.1 as the density element between orbital 2 and 3. (remember that Python is 0-based elements).
- Parameters
geometry (Geometry) – parent geometry to create a density matrix from. The density matrix will have size equivalent to the number of orbitals in the geometry
dim (int or Spin, optional) – number of components per element, may be a
Spin
objectdtype (np.dtype, optional) – data type contained in the density matrix. See details of
Spin
for default values.nnzpr (int, optional) – number of initially allocated memory per orbital in the density matrix. For increased performance this should be larger than the actual number of entries per orbital.
spin (Spin, optional) – equivalent to
dim
argument. This keyword-only argument has precedence overdim
.orthogonal (bool, optional) – whether the density matrix corresponds to a non-orthogonal basis. In this case the dimensionality of the density matrix is one more than
dim
. This is a keyword-only argument.
Methods
Dk
([k, dtype, gauge, format])Setup the density matrix for a given k-point
Rij
([what, dtype])Create a sparse matrix with the vectors between atoms/orbitals
Sk
([k, dtype, gauge, format])Setup the overlap matrix for a given k-point
add
(other[, axis, offset])Add two sparse matrices by adding the parameters to one set.
append
(other, axis[, eps, scale])Append other along axis to construct a new connected sparse matrix
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, param)Create a simple function for passing to the
construct
function.cut
(seps, axis, *args, **kwargs)Cuts the sparse orbital model into different parts.
dDk
([k, dtype, gauge, format])Setup the density matrix derivative for a given k-point
dSk
([k, dtype, gauge, format])Setup the \(k\)-derivatie of the overlap matrix for a given k-point
ddDk
([k, dtype, gauge, format])Setup the density matrix double derivative for a given k-point
ddSk
([k, dtype, gauge, format])Setup the double \(k\)-derivatie of the overlap matrix for a given k-point
density
(grid[, spinor, tol, eta])Expand the density matrix to the charge density on a grid
edges
([atoms, exclude, orbitals])Retrieve edges (connections) of a given atom or list of atom’s
eig
([k, gauge, eigvals_only])Returns the eigenvalues of the physical quantity (using the non-Hermitian solver)
eigh
([k, gauge, eigvals_only])Returns the eigenvalues of the physical quantity
eigsh
([k, n, gauge, eigvals_only])Calculates a subset of eigenvalues of the physical quantity (default 10)
eliminate_zeros
(*args, **kwargs)Removes all zero elements from the sparse matrix
empty
([keep_nnz])See
empty
for detailsfinalize
()Finalizes the model
fromsp
(geometry, P[, S])Create a sparse model from a preset Geometry and a list of sparse matrices
iter_nnz
([atoms, orbitals])Iterations of the non-zero elements
iter_orbitals
([atoms, local])Iterations of the orbital space in the geometry, two indices from loop
mulliken
([projection])Calculate Mulliken charges from the density matrix
nonzero
([atoms, only_col])Indices row and column indices where non-zero elements exists
orbital_momentum
([projection, method])Calculate orbital angular momentum on either atoms or orbitals
prepend
(other, axis[, eps])See
append
for detailsread
(sile, *args, **kwargs)Reads density matrix from Sile using read_density_matrix.
remove
(atoms)Remove a subset of this sparse matrix by only retaining the atoms corresponding to atom
remove_orbital
(atoms, orbitals)Remove a subset of orbitals on atom according to orbital
repeat
(reps, axis)Create a repeated sparse orbital object, equivalent to Geometry.repeat
reset
([dim, dtype, nnzpr])The sparsity pattern has all elements removed and everything is reset.
rij
([what, dtype])Create a sparse matrix with the distance between atoms/orbitals
set_nsc
(*args, **kwargs)Reset the number of allowed supercells in the sparse orbital
spalign
(other)See
align
for detailsspin_align
(vec)Aligns all spin along the vector vec
spin_rotate
(angles[, rad])Rotates spin-boxes by fixed angles around the \(x\), \(y\) and \(z\) axis, respectively.
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 atoms corresponding to atom
sub_orbital
(atoms, orbitals)Retain only a subset of the orbitals on atom according to orbital
swap
(a, b)Swaps atoms in the sparse geometry to obtain a new order of atoms
tile
(reps, axis)Create a tiled sparse orbital object, equivalent to Geometry.tile
toSparseAtom
([dim, dtype])Convert the sparse object (without data) to a new sparse object with equivalent but reduced sparse pattern
tocsr
([dim, isc])Return a
csr_matrix
for the specified dimensiontranspose
([hermitian, spin, sort])A transpose copy of this object, possibly apply the Hermitian conjugate as well
trs
()Create a new matrix with applied time-reversal-symmetry
write
(sile, *args, **kwargs)Writes a density matrix to the Sile as implemented in the
Sile.write_density_matrix
methodAccess the density matrix elements
Access the overlap elements associated with the sparse matrix
Number of components per element
Data type of sparse elements (in str)
Data type of sparse elements
Whether the contained data is finalized and non-used elements have been removed
deprecated geometry
Associated geometry
Number of non-zero elements
True if the object is using a non-orthogonal basis
True if the object is using an orthogonal basis
Shape of sparse matrix
Associated spin class
- property D¶
Access the density matrix elements
- Dk(k=(0, 0, 0), dtype=None, gauge='R', format='csr', *args, **kwargs)[source]¶
Setup the density matrix for a given k-point
Creation and return of the density matrix for a given k-point (default to Gamma).
Notes
Currently the implemented gauge for the k-point is the cell vector gauge:
\[\mathbf D(k) = \mathbf D_{\nu\mu} e^{i k R}\]where \(R\) is an integer times the cell vector and \(\nu\), \(\mu\) are orbital indices.
Another possible gauge is the orbital distance which can be written as
\[\mathbf D(k) = \mathbf D_{\nu\mu} e^{i k r}\]where \(r\) is the distance between the orbitals.
- Parameters
k (array_like) – the k-point to setup the density matrix at
dtype (numpy.dtype , optional) – the data type of the returned matrix. Do NOT request non-complex data-type for non-Gamma k. The default data-type is
numpy.complex128
gauge ({'R', 'r'}) – the chosen gauge, R for cell vector gauge, and r for orbital distance gauge.
format ({'csr', 'array', 'dense', 'coo', ...}) – the returned format of the matrix, defaulting to the
scipy.sparse.csr_matrix
, however if one always requires operations on dense matrices, one can always return innumpy.ndarray
(‘array’/’dense’/’matrix’). Prefixing with ‘sc:’, or simply ‘sc’ returns the matrix in supercell format with phases.spin (int, optional) – if the density matrix is a spin polarized one can extract the specific spin direction matrix by passing an integer (0 or 1). If the density matrix is not
Spin.POLARIZED
this keyword is ignored.
See also
- Returns
matrix – the density matrix at \(k\). The returned object depends on format.
- Return type
numpy.ndarray or scipy.sparse.*_matrix
- Rij(what='orbital', dtype=<class 'numpy.float64'>)¶
Create a sparse matrix with the vectors between atoms/orbitals
- Parameters
what ({'orbital', 'atom'}) – which kind of sparse vector matrix to return, either an atomic vector matrix or an orbital vector matrix. The orbital matrix is equivalent to the atomic one with the same vectors repeated for the same atomic orbitals. The default is the same type as the parent class.
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.
- property S¶
Access the overlap elements associated with the sparse matrix
- Sk(k=(0, 0, 0), dtype=None, gauge='R', format='csr', *args, **kwargs)¶
Setup the overlap matrix for a given k-point
Creation and return of the overlap matrix for a given k-point (default to Gamma).
Notes
Currently the implemented gauge for the k-point is the cell vector gauge:
\[\mathbf S(k) = \mathbf S_{\nu\mu} e^{i k R}\]where \(R\) is an integer times the cell vector and \(\nu\), \(\mu\) are orbital indices.
Another possible gauge is the orbital distance which can be written as
\[\mathbf S(k) = \mathbf S_{\nu\mu} e^{i k r}\]where \(r\) is the distance between the orbitals.
- Parameters
k (array_like, optional) – the k-point to setup the overlap at (default Gamma point)
dtype (numpy.dtype, optional) – the data type of the returned matrix. Do NOT request non-complex data-type for non-Gamma k. The default data-type is
numpy.complex128
gauge ({'R', 'r'}) – the chosen gauge, R for cell vector gauge, and r for orbital distance gauge.
format ({'csr', 'array', 'matrix', 'coo', ...}) – the returned format of the matrix, defaulting to the
scipy.sparse.csr_matrix
, however if one always requires operations on dense matrices, one can always return innumpy.ndarray
(‘array’/’dense’/’matrix’). Prefixing with ‘sc:’, or simply ‘sc’ returns the matrix in supercell format with phases. This is useful for e.g. bond-current calculations where individual hopping + phases are required.
See also
- Returns
matrix – the overlap matrix at \(k\). The returned object depends on format.
- Return type
numpy.ndarray or scipy.sparse.*_matrix
- add(other, axis=None, offset=(0, 0, 0))¶
Add two sparse matrices by adding the parameters to one set. The final matrix will have no couplings between self and other
The final sparse matrix will not have any couplings between self and other. Not even if they have commensurate overlapping regions. If you want to create couplings you have to use
append
but that requires the structures are commensurate in the coupling region.- Parameters
other (SparseGeometry) – the other sparse matrix to be added, all atoms will be appended
axis (int or None, optional) – whether a specific axis of the cell will be added to the final geometry. For
None
the final cell will be that of self, otherwise the lattice vector corresponding to axis will be appended.offset ((3,), optional) – offset in geometry of other when adding the atoms.
- append(other, axis, eps=0.01, scale=1)¶
Append other along axis to construct a new connected sparse matrix
This method tries to append two sparse geometry objects together by the following these steps:
Create the new extended geometry
Use neighbor cell couplings from self as the couplings to other This may cause problems if the coupling atoms are not exactly equi-positioned. If the coupling coordinates and the coordinates in other differ by more than 0.01 Ang, a warning will be issued. If this difference is above eps the couplings will be removed.
When appending sparse matrices made up of atoms, this method assumes that the orbitals on the overlapping atoms have the same orbitals, as well as the same orbital ordering.
Examples
>>> sporb = SparseOrbital(....) >>> sporb2 = sporb.append(sporb, 0) >>> sporbt = sporb.tile(2, 0) >>> sporb2.spsame(sporbt) True
To retain couplings only from the left sparse matrix, do:
>>> sporb = left.append(right, 0, scale=(2, 0)) >>> sporb = (sporb + sporb.transpose()) * 0.5
To retain couplings only from the right sparse matrix, do:
>>> sporb = left.append(right, 0, scale=(0, 2.)) >>> sporb = (sporb + sporb.transpose()) * 0.5
Notes
The current implentation does not preserve the hermiticity of the matrix. If you want to preserve hermiticity of the matrix you have to do the following:
>>> h = (h + h.transpose()) / 2
- Parameters
other (object) – must be an object of the same type as self
axis (int) – axis to append the two sparse geometries along
eps (float, optional) – tolerance that all coordinates must be within to allow an append. It is important that this value is smaller than half the distance between the two closests atoms such that there is no ambiguity in selecting equivalent atoms. An internal stricter eps is used as a baseline, see above.
scale (float or array_like, optional) – the scale used for the overlapping region. For scalar values it corresponds to passing:
(scale, scale)
. For array-like inputscale[0]
refers to the scale of the matrix elements coupling from self, whilescale[1]
is the scale of the matrix elements in other.
See also
prepend
equivalent scheme as this method
add
merge two matrices without considering overlap or commensurability
transpose
ensure hermiticity by using this routine
Geometry.append
,Geometry.prepend
SparseCSR.scale_columns
method used to scale the two matrix elements values
- Raises
ValueError if the two geometries are not compatible for either coordinate, orbital or supercell errors –
- Returns
a new instance with two sparse matrices joined and appended together
- Return type
- construct(func, na_iR=1000, method='rand', eta=False)¶
Automatically construct the sparse model based on a function that does the setting up of the elements
This may be called in two variants.
Pass a function (func), see e.g.
create_construct
which does the setting up.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 theR[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, param)¶
Create a simple function for passing to the
construct
function.This is to relieve the creation of simplistic functions needed for setting up sparse elements.
For simple matrices 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, param): ... self[ia, ix] = p
In the non-colinear case the matrix element \(M_{ij}\) will be set to input values param if \(i \le j\) and the Hermitian conjugated values for \(j < i\).
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.
This method issues warnings if the on-site terms are not Hermitian for spin-orbit systems. Do note that it still creates the matrices based on the input.
- Parameters
R (array_like) – radii parameters for different shells. Must have same length as param or one less. If one less it will be extended with
R[0]/100
param (array_like) – coupling constants corresponding to the R ranges.
param[0,:]
are the elements for the all atoms withinR[0]
of each atom.
See also
construct
routine to create the sparse matrix from a generic function (as returned from
create_construct
)
- cut(seps, axis, *args, **kwargs)¶
Cuts the sparse orbital model into different parts.
Recreates a new sparse orbital object with only the cutted atoms in the structure.
Cutting is the opposite of tiling.
- dDk(k=(0, 0, 0), dtype=None, gauge='R', format='csr', *args, **kwargs)[source]¶
Setup the density matrix derivative for a given k-point
Creation and return of the density matrix derivative for a given k-point (default to Gamma).
Notes
Currently the implemented gauge for the k-point is the cell vector gauge:
\[\nabla_k \mathbf D_\alpha(k) = i R_\alpha \mathbf D_{\nu\mu} e^{i k R}\]where \(R\) is an integer times the cell vector and \(\nu\), \(\mu\) are orbital indices. And \(\alpha\) is one of the Cartesian directions.
Another possible gauge is the orbital distance which can be written as
\[\nabla_k \mathbf D_\alpha(k) = i r_\alpha \mathbf D_{\nu\mu} e^{i k r}\]where \(r\) is the distance between the orbitals.
- Parameters
k (array_like) – the k-point to setup the density matrix at
dtype (numpy.dtype , optional) – the data type of the returned matrix. Do NOT request non-complex data-type for non-Gamma k. The default data-type is
numpy.complex128
gauge ({'R', 'r'}) – the chosen gauge, R for cell vector gauge, and r for orbital distance gauge.
format ({'csr', 'array', 'dense', 'coo', ...}) – the returned format of the matrix, defaulting to the
scipy.sparse.csr_matrix
, however if one always requires operations on dense matrices, one can always return innumpy.ndarray
(‘array’/’dense’/’matrix’).spin (int, optional) – if the density matrix is a spin polarized one can extract the specific spin direction matrix by passing an integer (0 or 1). If the density matrix is not
Spin.POLARIZED
this keyword is ignored.
- Returns
for each of the Cartesian directions a \(\partial \mathbf D(k)/\partial k\) is returned.
- Return type
- dSk(k=(0, 0, 0), dtype=None, gauge='R', format='csr', *args, **kwargs)¶
Setup the \(k\)-derivatie of the overlap matrix for a given k-point
Creation and return of the derivative of the overlap matrix for a given k-point (default to Gamma).
Notes
Currently the implemented gauge for the k-point is the cell vector gauge:
\[\nabla_k \mathbf S_\alpha(k) = i R_\alpha \mathbf S_{\nu\mu} e^{i k R}\]where \(R\) is an integer times the cell vector and \(\nu\), \(\mu\) are orbital indices. And \(\alpha\) is one of the Cartesian directions.
Another possible gauge is the orbital distance which can be written as
\[\nabla_k \mathbf S_\alpha(k) = i r_\alpha \mathbf S_{ij} e^{i k r}\]where \(r\) is the distance between the orbitals.
- Parameters
k (array_like, optional) – the k-point to setup the overlap at (default Gamma point)
dtype (numpy.dtype, optional) – the data type of the returned matrix. Do NOT request non-complex data-type for non-Gamma k. The default data-type is
numpy.complex128
gauge ({'R', 'r'}) – the chosen gauge, R for cell vector gauge, and r for orbital distance gauge.
format ({'csr', 'array', 'matrix', 'coo', ...}) – the returned format of the matrix, defaulting to the
scipy.sparse.csr_matrix
, however if one always requires operations on dense matrices, one can always return innumpy.ndarray
(‘array’/’dense’/’matrix’).
- Returns
for each of the Cartesian directions a \(\partial \mathbf S(k)/\partial k\) is returned.
- Return type
- ddDk(k=(0, 0, 0), dtype=None, gauge='R', format='csr', *args, **kwargs)[source]¶
Setup the density matrix double derivative for a given k-point
Creation and return of the density matrix double derivative for a given k-point (default to Gamma).
Notes
Currently the implemented gauge for the k-point is the cell vector gauge:
\[\nabla_k^2 \mathbf D_{\alpha\beta}(k) = - R_\alpha R_\beta \mathbf D_{\nu\mu} e^{i k R}\]where \(R\) is an integer times the cell vector and \(\nu\), \(\mu\) are orbital indices. And \(\alpha\) and \(\beta\) are one of the Cartesian directions.
Another possible gauge is the orbital distance which can be written as
\[\nabla_k^2 \mathbf D_{\alpha\beta}(k) = - r_\alpha r_\beta \mathbf D_{\nu\mu} e^{i k r}\]where \(r\) is the distance between the orbitals.
- Parameters
k (array_like) – the k-point to setup the density matrix at
dtype (numpy.dtype , optional) – the data type of the returned matrix. Do NOT request non-complex data-type for non-Gamma k. The default data-type is
numpy.complex128
gauge ({'R', 'r'}) – the chosen gauge, R for cell vector gauge, and r for orbital distance gauge.
format ({'csr', 'array', 'dense', 'coo', ...}) – the returned format of the matrix, defaulting to the
scipy.sparse.csr_matrix
, however if one always requires operations on dense matrices, one can always return innumpy.ndarray
(‘array’/’dense’/’matrix’).spin (int, optional) – if the density matrix is a spin polarized one can extract the specific spin direction matrix by passing an integer (0 or 1). If the density matrix is not
Spin.POLARIZED
this keyword is ignored.
- Returns
for each of the Cartesian directions (in Voigt representation); xx, yy, zz, zy, xz, xy
- Return type
list of matrices
- ddSk(k=(0, 0, 0), dtype=None, gauge='R', format='csr', *args, **kwargs)¶
Setup the double \(k\)-derivatie of the overlap matrix for a given k-point
Creation and return of the double derivative of the overlap matrix for a given k-point (default to Gamma).
Notes
Currently the implemented gauge for the k-point is the cell vector gauge:
\[\nabla_k^2 \mathbf S_{\alpha\beta}(k) = - R_\alpha R_\beta \mathbf S_{\nu\mu} e^{i k R}\]where \(R\) is an integer times the cell vector and \(\nu\), \(\mu\) are orbital indices. And \(\alpha\) and \(\beta\) are one of the Cartesian directions.
Another possible gauge is the orbital distance which can be written as
\[\nabla_k^2 \mathbf S_{\alpha\beta}(k) = - r_\alpha r_\beta \mathbf S_{ij} e^{i k r}\]where \(r\) is the distance between the orbitals.
- Parameters
k (array_like, optional) – the k-point to setup the overlap at (default Gamma point)
dtype (numpy.dtype, optional) – the data type of the returned matrix. Do NOT request non-complex data-type for non-Gamma k. The default data-type is
numpy.complex128
gauge ({'R', 'r'}) – the chosen gauge, R for cell vector gauge, and r for orbital distance gauge.
format ({'csr', 'array', 'matrix', 'coo', ...}) – the returned format of the matrix, defaulting to the
scipy.sparse.csr_matrix
, however if one always requires operations on dense matrices, one can always return innumpy.ndarray
(‘array’/’dense’/’matrix’).
- Returns
for each of the Cartesian directions (in Voigt representation); xx, yy, zz, zy, xz, xy
- Return type
list of matrices
- density(grid, spinor=None, tol=1e-07, eta=False)¶
Expand the density matrix to the charge density on a grid
This routine calculates the real-space density components on a specified grid.
This is an in-place operation that adds to the current values in the grid.
Note: To calculate \(\rho(\mathbf r)\) in a unit-cell different from the originating geometry, simply pass a grid with a unit-cell different than the originating supercell.
The real-space density is calculated as:
\[\rho(\mathbf r) = \sum_{\nu\mu}\phi_\nu(\mathbf r)\phi_\mu(\mathbf r) D_{\nu\mu}\]While for non-collinear/spin-orbit calculations the density is determined from the spinor component (spinor) by
\[\rho_{\boldsymbol\sigma}(\mathbf r) = \sum_{\nu\mu}\phi_\nu(\mathbf r)\phi_\mu(\mathbf r) \sum_\alpha [\boldsymbol\sigma \mathbf \rho_{\nu\mu}]_{\alpha\alpha}\]Here \(\boldsymbol\sigma\) corresponds to a spinor operator to extract relevant quantities. By passing the identity matrix the total charge is added. By using the Pauli matrix \(\boldsymbol\sigma_x\) only the \(x\) component of the density is added to the grid (see
Spin.X
).- Parameters
grid (Grid) – the grid on which to add the density (the density is in
e/Ang^3
)spinor ((2,) or (2, 2), optional) – the spinor matrix to obtain the diagonal components of the density. For un-polarized density matrices this keyword has no influence. For spin-polarized it has to be either 1 integer or a vector of length 2 (defaults to total density). For non-collinear/spin-orbit density matrices it has to be a 2x2 matrix (defaults to total density).
tol (float, optional) – DM tolerance for accepted values. For all density matrix elements with absolute values below the tolerance, they will be treated as strictly zeros.
eta (bool, optional) – show a progressbar on stdout
- 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=None, exclude=None, orbitals=None)¶
Retrieve edges (connections) of a given atom or list of atom’s
The returned edges are unique and sorted (see
numpy.unique
) and are returned in supercell indices (i.e.0 <= edge < self.geometry.no_s
).- Parameters
atoms (int or list of int) – the edges are returned only for the given atom (but by using all orbitals of the requested atom). The returned edges are also atoms.
exclude (int or list of int or None, optional) – remove edges which are in the exclude list, this list refers to orbitals.
orbital (int or list of int) – the edges are returned only for the given orbital. The returned edges are orbitals.
See also
SparseCSR.edges
the underlying routine used for extracting the edges
- eig(k=(0, 0, 0), gauge='R', eigvals_only=True, **kwargs)¶
Returns the eigenvalues of the physical quantity (using the non-Hermitian solver)
Setup the system and overlap matrix with respect to the given k-point and calculate the eigenvalues.
All subsequent arguments gets passed directly to
scipy.linalg.eig
- Parameters
spin (int, optional) – the spin-component to calculate the eigenvalue spectrum of, note that this parameter is only valid for
Spin.POLARIZED
matrices.
- eigh(k=(0, 0, 0), gauge='R', eigvals_only=True, **kwargs)¶
Returns the eigenvalues of the physical quantity
Setup the system and overlap matrix with respect to the given k-point and calculate the eigenvalues.
All subsequent arguments gets passed directly to
scipy.linalg.eigh
- Parameters
spin (int, optional) – the spin-component to calculate the eigenvalue spectrum of, note that this parameter is only valid for
Spin.POLARIZED
matrices.
- eigsh(k=(0, 0, 0), n=10, gauge='R', eigvals_only=True, **kwargs)¶
Calculates a subset of eigenvalues of the physical quantity (default 10)
Setup the quantity and overlap matrix with respect to the given k-point and calculate a subset of the eigenvalues using the sparse algorithms.
All subsequent arguments gets passed directly to
scipy.linalg.eigsh
- Parameters
spin (int, optional) – the spin-component to calculate the eigenvalue spectrum of, note that this parameter is only valid for
Spin.POLARIZED
matrices.
- 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, S=None, **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
S (scipy.sparse, optional) – if provided this refers to the overlap matrix and will force the returned sparse matrix to be non-orthogonal
**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 and optionally being non-orthogonal if
S
is not none- Return type
SparseGeometry
- property geom¶
deprecated geometry
- property geometry¶
Associated geometry
- iter_nnz(atoms=None, orbitals=None)¶
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
- iter_orbitals(atoms=None, local=False)¶
Iterations of the orbital space in the geometry, two indices from loop
An iterator returning the current atomic index and the corresponding orbital index.
>>> for ia, io in self.iter_orbitals():
In the above case
io
always belongs to atom ia and ia may be repeated according to the number of orbitals associated with the atom ia.- Parameters
- Yields
ia – atomic index
io – orbital index
See also
Geometry.iter_orbitals
method used to iterate orbitals
- mulliken(projection='orbital')¶
Calculate Mulliken charges from the density matrix
In the following \(\nu\) and \(\mu\) are orbital indices. Atomic indices noted by \(\alpha\), \(\beta\). Matrices \(\boldsymbol\rho\) and \(\mathbf S\) are density and overlap matrices, respectively.
For polarized calculations the Mulliken charges are calculated as (for each spin-channel)
\[\begin{split}M_{\nu} &= \sum_mu [\boldsymbol\rho \mathbf S]_{\nu\mu} \\ M_{\alpha} &= \sum_{\nu\in\alpha} M_{\nu}\end{split}\]For non-colinear calculations (including spin-orbit) they are calculated as above but using the spin-box per orbital (\(\sigma\) is spin)
\[\begin{split}M_{\nu} &= \sum_\sigma\sum_mu [\boldsymbol\rho \mathbf S]_{\nu\mu,\sigma\sigma} \\ S_{\nu}^x &= \sum_mu \Re [\boldsymbol\rho \mathbf S]_{\nu\mu,\uparrow\downarrow} + \Re [\boldsymbol\rho \mathbf S]_{\nu\mu,\downarrow\uparrow} \\ S_{\nu}^y &= \sum_mu \Im [\boldsymbol\rho \mathbf S]_{\nu\mu,\uparrow\downarrow} - \Im [\boldsymbol\rho \mathbf S]_{\nu\mu,\downarrow\uparrow} \\ S_{\nu}^z &= \sum_mu \Re [\boldsymbol\rho \mathbf S]_{\nu\mu,\uparrow\uparrow} - \Re [\boldsymbol\rho \mathbf S]_{\nu\mu,\downarrow\downarrow}\end{split}\]- Parameters
projection ({'orbital', 'atom'}) – how the Mulliken charges are returned. Can be atom-resolved, orbital-resolved or the charge matrix (off-diagonal elements)
- Returns
if projection does not contain matrix, the first dimension contains the orbitals, and the 2nd the spin information
- Return type
- property nnz¶
Number of non-zero elements
- property non_orthogonal¶
True if the object is using a non-orthogonal basis
- nonzero(atoms=None, only_col=False)¶
Indices row and column indices where non-zero elements exists
- Parameters
See also
SparseCSR.nonzero
the equivalent function call
- orbital_momentum(projection='orbital', method='onsite')[source]¶
Calculate orbital angular momentum on either atoms or orbitals
Currently this implementation equals the Siesta implementation in that the on-site approximation is enforced thus limiting the calculated quantities to obey the following conditions:
Same atom
\(l>0\)
\(l_\nu \equiv l_\mu\)
\(m_\nu \neq m_\mu\)
\(\zeta_\nu \equiv \zeta_\mu\)
This allows one to sum the orbital angular moments on a per atom site.
- Parameters
projection ({'orbital', 'atom'}) – whether the angular momentum is resolved per atom, or per orbital
method ({'onsite'}) – method used to calculate the angular momentum
- Returns
orbital angular momentum with the last dimension equalling the \(L_x\), \(L_y\) and \(L_z\) components
- Return type
- property orthogonal¶
True if the object is using an orthogonal basis
- prepend(other, axis, eps=0.01)¶
See
append
for detailsThis is currently equivalent to:
>>> other.append(self, axis, eps)
- static read(sile, *args, **kwargs)[source]¶
Reads density matrix from Sile using read_density_matrix.
- Parameters
sile (Sile, str or pathlib.Path) – a Sile object which will be used to read the density matrix and the overlap matrix (if any) if it is a string it will create a new sile using get_sile.
* (args passed directly to
read_density_matrix(,**)
) –
- remove(atoms)¶
Remove a subset of this sparse matrix by only retaining the atoms corresponding to atom
- Parameters
atoms (array_like of int or Atom) – indices of removed atoms or Atom for direct removal of all atoms
- remove_orbital(atoms, orbitals)¶
Remove a subset of orbitals on atom according to orbital
- Parameters
Examples
>>> obj = SparseOrbital(...) >>> # remove the second orbital on the 2nd atom >>> # all other orbitals are retained >>> obj.remove_orbital(1, 1)
- repeat(reps, axis)¶
Create a repeated sparse orbital 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
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(what='orbital', dtype=<class 'numpy.float64'>)¶
Create a sparse matrix with the distance between atoms/orbitals
- Parameters
what ({'orbital', 'atom'}) – which kind of sparse distance matrix to return, either an atomic distance matrix or an orbital distance matrix. The orbital matrix is equivalent to the atomic one with the same distance repeated for the same atomic orbitals. The default is the same type as the parent class.
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)¶
Reset the number of allowed supercells in the sparse orbital
If one reduces the number of supercells any sparse element that references the supercell will be deleted.
See SuperCell.set_nsc for allowed parameters.
See also
SuperCell.set_nsc
the underlying called method
- property shape¶
Shape of sparse matrix
- spalign(other)¶
See
align
for details
- property spin¶
Associated spin class
- spin_align(vec)¶
Aligns all spin along the vector vec
In case the matrix is polarized and vec is not aligned at the z-axis, the returned matrix will be a non-collinear spin configuration.
- Parameters
vec ((3,)) – vector to align the spin boxes against
See also
spin_rotate
rotate spin-boxes by a fixed amount (does not align spins)
- Returns
a new object with aligned spins
- Return type
- spin_rotate(angles, rad=False)¶
Rotates spin-boxes by fixed angles around the \(x\), \(y\) and \(z\) axis, respectively.
The angles are with respect to each spin-boxes initial angle. One should use
spin_align
to fix all angles along a specific direction.Notes
For a polarized matrix: The returned matrix will be in non-collinear spin-configuration in case the angles does not reflect a pure flip of spin in the \(z\)-axis.
- Parameters
angles ((3,)) – angle to rotate spin boxes, \(x\), \(y\) and :math`z`, respectively
rad (bool, optional) – Determines the unit of angles, for true it is in radians
See also
spin_align
align all spin-boxes along a specific direction
- Returns
a new object with rotated spins
- Return type
- 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)¶
Create a subset of this sparse matrix by only retaining the atoms corresponding to atom
Negative indices are wrapped and thus works, supercell atoms are also wrapped to the unit-cell.
- Parameters
atoms (array_like of int or Atom) – indices of retained atoms or Atom for retaining only that atom
Examples
>>> obj = SparseOrbital(...) >>> obj.sub(1) # only retain the second atom in the SparseGeometry >>> obj.sub(obj.atoms.atom[0]) # retain all atoms which is equivalent to >>> # the first atomic specie
- sub_orbital(atoms, orbitals)¶
Retain only a subset of the orbitals on atom according to orbital
This allows one to retain only a given subset of the sparse matrix elements.
- Parameters
atoms (array_like of int or Atom) – indices of atoms or Atom that will be reduced in size according to orbital
orbitals (array_like of int or Orbital) – indices of the orbitals on atom that are retained in the sparse matrix, the list of orbitals will be sorted. One cannot re-arrange matrix elements currently.
Notes
Future implementations may allow one to re-arange orbitals using this method.
Examples
>>> obj = SparseOrbital(...) >>> # only retain the second orbital on the 2nd atom >>> # all other orbitals are retained >>> obj.sub_orbital(1, 1)
- 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)¶
Create a tiled sparse orbital object, equivalent to Geometry.tile
The already existing sparse elements are extrapolated to the new supercell by repeating them in blocks like the coordinates.
- Parameters
See also
Geometry.tile
the same ordering as the final geometry
Geometry.repeat
a different ordering of the final geometry
repeat
a different ordering of the final geometry
- toSparseAtom(dim=None, dtype=None)¶
Convert the sparse object (without data) to a new sparse object with equivalent but reduced sparse pattern
This converts the orbital sparse pattern to an atomic sparse pattern.
- Parameters
dim (int, optional) – number of dimensions allocated in the SparseAtom object, default to the same
dtype (numpy.dtype, optional) – used data-type for the sparse object. Defaults to the same.
- tocsr(dim=0, isc=None, **kwargs)¶
Return a
csr_matrix
for the specified dimension
- transpose(hermitian=False, spin=True, sort=True)¶
A transpose copy of this object, possibly apply the Hermitian conjugate as well
- Parameters
hermitian (bool, optional) – if true, also emply a spin-box Hermitian operator to ensure TRS, otherwise only return the transpose values.
spin (bool, optional) – whether the spin-box is also transposed if this is false, and hermitian is true, then only imaginary values will change sign.
sort (bool, optional) – the returned columns for the transposed structure will be sorted if this is true, default
- trs()¶
Create a new matrix with applied time-reversal-symmetry
Time reversal symmetry is applied using the following equality:
\[2\mathbf M^{\mathrm{TRS}} = \mathbf M + \boldsymbol\sigma_y \mathbf M^* \boldsymbol\sigma_y\]where \(*\) is the conjugation operator.