sisl.physics.hamiltonian module¶
Tight-binding class to create tight-binding models.
-
class
sisl.physics.hamiltonian.
Hamiltonian
(geom, nnzpr=None, orthogonal=True, spin=1, dtype=None, *args, **kwargs)[source]¶ Bases:
object
Hamiltonian object containing the coupling constants between orbitals.
- The Hamiltonian object contains information regarding the
- geometry
- coupling constants between orbitals
It contains an intrinsic sparse matrix of the Hamiltonian elements.
Assigning or changing Hamiltonian elements is as easy as with standard
numpy
assignments:>>> ham = Hamiltonian(...) >>> ham.H[1,2] = 0.1
which assigns 0.1 as the coupling constant between orbital 2 and 3. (remember that Python is 0-based elements).
Attributes
H
S
dtype
Return data type of Hamiltonian (and overlap matrix) finalized
Whether the contained data is finalized and non-used elements have been removed geom
Return the attached geometry geometry
Return the attached geometry nnz
Returns number of non-zero elements in the tight-binding model no
Returns number of orbitals as used when the object was created orthogonal
Return whether the Hamiltonian is orthogonal spin
Return number of spin-components in Hamiltonian Methods
construct
(func[, na_iR, method, eta])Automatically construct the Hamiltonian model based on a function that does the setting up of the Hamiltonian copy
([dtype])Return a copy of the Hamiltonian
objectcreate_construct
(dR, param)Returns a simple function for passing to the construct
function.cut
(seps, axis, *args, **kwargs)Cuts the tight-binding model into different parts. eigh
([k, atoms, eigvals_only, overwrite_a, ...])Returns the eigenvalues of the Hamiltonian eigsh
([k, n, atoms, eigvals_only])Returns a subset of eigenvalues of the Hamiltonian (default 10) eliminate_zeros
()Removes all zero elememts from the sparse matrix empty
([keep])See SparseCSR.empty
for detailsfinalize
()Finalizes the tight-binding model fromsp
(geom, H[, S])Returns a tight-binding model from a preset H, S and Geometry iter
([local])Iterations of the orbital space in the geometry, two indices from loop iter_nnz
([atom, orbital])Iterations of the non-zero elements, returns a tuple of orbital and coupling orbital read
(sile, *args, **kwargs)Reads Hamiltonian from Sile
usingread_H
.remove
(atom)Remove atom from the Hamiltonian. repeat
(reps, axis)Returns a repeated tight-binding model for this, much like the Geometry
reset
([nnzpr, orthogonal, spin, dtype])The sparsity pattern is cleaned and every thing is reset. spsame
(other)Compare two Hamiltonians and check whether they have the same entries. sub
(atom)Returns a subset of atoms from the geometry. sub_supercell
(sc)Creates a new object with only the given supercells tile
(reps, axis)Returns a tiled tight-binding model for this, much like the Geometry
tocsr
(index[, isc])Return a scipy.sparse.csr_matrix
from the specified indexwrite
(sile, *args, **kwargs)Writes a tight-binding model to the Sile
as implemented in theSile.write_hamiltonian
methodCreate tight-binding model from geometry
Initializes a tight-binding model using the
geom
object as the underlying geometry for the tight-binding parameters.Attributes
H
S
dtype
Return data type of Hamiltonian (and overlap matrix) finalized
Whether the contained data is finalized and non-used elements have been removed geom
Return the attached geometry geometry
Return the attached geometry nnz
Returns number of non-zero elements in the tight-binding model no
Returns number of orbitals as used when the object was created orthogonal
Return whether the Hamiltonian is orthogonal spin
Return number of spin-components in Hamiltonian Methods
construct
(func[, na_iR, method, eta])Automatically construct the Hamiltonian model based on a function that does the setting up of the Hamiltonian copy
([dtype])Return a copy of the Hamiltonian
objectcreate_construct
(dR, param)Returns a simple function for passing to the construct
function.cut
(seps, axis, *args, **kwargs)Cuts the tight-binding model into different parts. eigh
([k, atoms, eigvals_only, overwrite_a, ...])Returns the eigenvalues of the Hamiltonian eigsh
([k, n, atoms, eigvals_only])Returns a subset of eigenvalues of the Hamiltonian (default 10) eliminate_zeros
()Removes all zero elememts from the sparse matrix empty
([keep])See SparseCSR.empty
for detailsfinalize
()Finalizes the tight-binding model fromsp
(geom, H[, S])Returns a tight-binding model from a preset H, S and Geometry iter
([local])Iterations of the orbital space in the geometry, two indices from loop iter_nnz
([atom, orbital])Iterations of the non-zero elements, returns a tuple of orbital and coupling orbital read
(sile, *args, **kwargs)Reads Hamiltonian from Sile
usingread_H
.remove
(atom)Remove atom from the Hamiltonian. repeat
(reps, axis)Returns a repeated tight-binding model for this, much like the Geometry
reset
([nnzpr, orthogonal, spin, dtype])The sparsity pattern is cleaned and every thing is reset. spsame
(other)Compare two Hamiltonians and check whether they have the same entries. sub
(atom)Returns a subset of atoms from the geometry. sub_supercell
(sc)Creates a new object with only the given supercells tile
(reps, axis)Returns a tiled tight-binding model for this, much like the Geometry
tocsr
(index[, isc])Return a scipy.sparse.csr_matrix
from the specified indexwrite
(sile, *args, **kwargs)Writes a tight-binding model to the Sile
as implemented in theSile.write_hamiltonian
method-
H
¶
-
S
¶
-
construct
(func, na_iR=1000, method='rand', eta=False)[source]¶ Automatically construct the Hamiltonian model based on a function that does the setting up of the Hamiltonian
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 isdR
the radii parameters for the corresponding tight-binding parameters. The second is the tight-binding parameters corresponding to thedR[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 the Hamiltonian object it-self (
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 examplefunc
could be:>>> def func(self, ia, idxs, idxs_xyz=None): >>> idx = self.geom.close(ia, dR=[0.1, 1.44], idx=idxs, idx_xyz=idxs_xyz) >>> self.H[ia, idx[0]] = 0. # on-site >>> self.H[ia, idx[1]] = -2.7 # nearest-neighbour
na_iR : int, 1000
number of atoms within the sphere for speeding up the
iter_block
loop.method : str, ‘rand’
method used in
Geometry.iter_block
, see there for detailseta: bool, False
whether an ETA will be printed
- Pass a function (
-
create_construct
(dR, param)[source]¶ Returns a simple function for passing to the
construct
function.This is simply to leviate the creation of simplistic functions needed for setting up the Hamiltonian.
Basically this returns a function: >>> def func(self, ia, idxs, idxs_xyz=None): >>> idx = self.geom.close(ia, dR=dR, idx=idxs) >>> for ix, p in zip(idx, param): >>> self[ia, ix] = p
Parameters: dR : array_like
radii parameters for tight-binding parameters. Must have same length as
param
or one less. If one less it will be extended withdR[0]/100
param : array_like
coupling constants corresponding to the
dR
ranges.param[0,:]
are the tight-binding parameter for the all atoms withindR[0]
of each atom.
-
cut
(seps, axis, *args, **kwargs)[source]¶ Cuts the tight-binding model into different parts.
Creates a tight-binding model by retaining the parameters for the cut-out region, possibly creating a super-cell.
Parameters: seps : integer, optional
number of times the structure will be cut.
axis : integer
the axis that will be cut
-
dtype
¶ Return data type of Hamiltonian (and overlap matrix)
-
eigh
(k=(0, 0, 0), atoms=None, eigvals_only=True, overwrite_a=True, overwrite_b=True, *args, **kwargs)[source]¶ Returns the eigenvalues of the Hamiltonian
Setup the Hamiltonian and overlap matrix with respect to the given k-point, then reduce the space to the specified atoms and calculate the eigenvalues.
All subsequent arguments gets passed directly to
scipy.linalg.eigh
-
eigsh
(k=(0, 0, 0), n=10, atoms=None, eigvals_only=True, *args, **kwargs)[source]¶ Returns a subset of eigenvalues of the Hamiltonian (default 10)
Setup the Hamiltonian and overlap matrix with respect to the given k-point, then reduce the space to the specified atoms and calculate a subset of the eigenvalues using the sparse algorithms.
All subsequent arguments gets passed directly to
scipy.linalg.eigsh
-
eliminate_zeros
()[source]¶ Removes all zero elememts from the sparse matrix
This is an in-place operation
-
empty
(keep=False)[source]¶ See
SparseCSR.empty
for details
-
finalize
()[source]¶ Finalizes the tight-binding model
Finalizes the tight-binding model so that no new sparse elements can be added.
Sparse elements can still be changed.
-
finalized
¶ Whether the contained data is finalized and non-used elements have been removed
-
classmethod
fromsp
(geom, H, S=None)[source]¶ Returns a tight-binding model from a preset H, S and Geometry
-
geom
¶ Return the attached geometry
-
geometry
¶ Return the attached geometry
-
iter
(local=False)[source]¶ 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:
In the above case
io
always belongs to atomia
andia
may be repeated according to the number of orbitals associated with the atomia
.Parameters: local :
bool=False
whether the orbital index is the global index, or the local index relative to the atom it resides on.
-
iter_nnz
(atom=None, orbital=None)[source]¶ Iterations of the non-zero elements, returns a tuple of orbital and coupling orbital
An iterator returning the current orbital index and the corresponding connected orbital where a non-zero is defined
>>> for io, jo in self.iter_nnz():
In the above case
io
andjo
are orbitals such that:>>> self.H[io,jo]
returns the non-zero element of the Hamiltonian.
One may reduce the iterated space by either requesting a specific set of atoms, or orbitals, _not_ both simultaneously.
Parameters: atom :
int
/array_like
iterate on couplings to the set of atoms (not compatible with
orbital
)orbital :
int
/array_like
iterate on couplings to the set of orbitals (not compatible with
atom
)Examples
Looping only on one or more atoms:
>>> for io, jo in self.iter_nnz(atom=[2, 3]): >>> # loop on all orbitals on atom 3 and 4 (0 indexing)
>>> for io, jo in self.iter_nnz(orbital=[2, 3]): >>> # loop on orbitals 3 and 4 (0 indexing)
-
nnz
¶ Returns number of non-zero elements in the tight-binding model
-
no
¶ Returns number of orbitals as used when the object was created
-
orthogonal
¶ Return whether the Hamiltonian is orthogonal
-
static
read
(sile, *args, **kwargs)[source]¶ Reads Hamiltonian from
Sile
usingread_H
.Parameters: sile :
Sile
, str* : args passed directly to
read_hamiltonian(,**)
-
remove
(atom)[source]¶ Remove atom from the Hamiltonian.
Indices passed MUST be unique.
Negative indices are wrapped and thus works.
Parameters: atom : array_like
indices of all atoms to be removed.
-
repeat
(reps, axis)[source]¶ Returns a repeated tight-binding model for this, much like the
Geometry
The already existing tight-binding parameters are extrapolated to the new supercell by repeating them in blocks like the coordinates.
Parameters: reps : number of repetitions
axis : direction of repetition
0, 1, 2 according to the cell-direction
-
reset
(nnzpr=None, orthogonal=True, spin=1, dtype=None)[source]¶ The sparsity pattern is cleaned and every thing is reset.
The object will be the same as if it had been initialized with the same geometry as it were created with.
Parameters: nnzpr: int
number of non-zero elements per row
orthogonal: boolean, True
if there is an overlap matrix associated with the Hamiltonian
spin: int, 1
number of spin-components
dtype: ``numpy.dtype``, `numpy.float64`
the datatype of the Hamiltonian
-
spin
¶ Return number of spin-components in Hamiltonian
-
spsame
(other)[source]¶ Compare two Hamiltonians and check whether they have the same entries.
This does not necessarily mean that the Hamiltonian values are the same
-
sub
(atom)[source]¶ Returns a subset of atoms from the geometry.
Indices passed MUST be unique.
Negative indices are wrapped and thus works.
Parameters: atom :
array_like
indices of all atoms to be removed.
-
tile
(reps, axis)[source]¶ Returns a tiled tight-binding model for this, much like the
Geometry
The already existing tight-binding parameters are extrapolated to the new supercell by repeating them in blocks like the coordinates.
Parameters: reps : number of tiles (repetitions)
axis : direction of tiling
0, 1, 2 according to the cell-direction
-
sisl.physics.hamiltonian.
TightBinding
¶ alias of
Hamiltonian