sisl.physics.phonon.ModeCPhonon
- class sisl.physics.phonon.ModeCPhonon(state, c, parent=None, **info)[source]
Bases:
StateC
A mode describing a physical quantity related to phonons, with associated coefficients of the mode
Methods
align_norm
(other[, ret_index, inplace])Align self with the site-norms of other, a copy may optionally be returned
align_phase
(other[, ret_index, inplace])Align self with the phases for other, a copy may be returned
asState
()change_gauge
(gauge[, offset])In-place change of the gauge of the state coefficients
copy
()Return a copy (only the coefficients and states are copied),
parent
andinfo
are passed by referencedegenerate
(atol)Find degenerate coefficients with a specified precision
derivative
([order, matrix, axes, operator])Calculate the derivative with respect to \(\mathbf k\) for a set of states up to a given order
inner
([ket, matrix, projection])Calculate the inner product as \(\mathbf A_{ij} = \langle\psi_i|\mathbf M|\psi'_j\rangle\)
ipr
([q])Calculate the inverse participation ratio (IPR) for arbitrary q values
iter
([asarray])An iterator looping over the states in this system
norm
()Return a vector with the Euclidean norm of each state \(\sqrt{\langle\psi|\psi\rangle}\)
norm2
([projection])Return a vector with the norm of each state \(\langle\psi|\psi\rangle\)
Return a normalized state where each state has \(|\psi|^2=1\)
outer
([ket, matrix])Return the outer product by \(\sum_\alpha|\psi_\alpha\rangle\langle\psi'_\alpha|\)
phase
([method, ret_index])Calculate the Euler angle (phase) for the elements of the state, in the range \(]-\pi;\pi]\)
remove
(index[, inplace])Return a new state without the specified indices
rotate
([phi, individual, inplace])Rotate all states to rotate the largest component to be along the angle phi
sort
([ascending])Sort and return a new
StateC
by sorting the coefficients (default to ascending)sub
(index[, inplace])Return a new state with only the specified states
tile
(reps, axis[, normalize, offset])Tile the state vectors for a new supercell
translate
(isc)Translate the vectors to a new unit-cell position
velocity
(*args, **kwargs)Calculate velocity of the modes
Attributes
The data-type of the state (in str)
Data-type for the state
Eigenmodes (states)
Returns the shape of the state
- align_norm(other: State, ret_index: bool = False, inplace: bool = False)
Align self with the site-norms of other, a copy may optionally be returned
To determine the new ordering of self first calculate the residual norm of the site-norms.
\[\delta N_{\alpha\beta} = \sum_i \big(\langle \psi^\alpha_i | \psi^\alpha_i\rangle - \langle \psi^\beta_i | \psi^\beta_i\rangle\big)^2\]where \(\alpha\) and \(\beta\) correspond to state indices in self and other, respectively. The new states (from self) returned is then ordered such that the index \(\alpha \equiv \beta'\) where \(\delta N_{\alpha\beta}\) is smallest.
- Parameters:
other (
State
) – the other state to align againstret_index – also return indices for the swapped indices
inplace – swap states in-place
- Returns:
self_swap (
State
) – A swapped instance of self, only if inplace is Falseindex (
array
ofint
) – the indices that swaps self to beself_swap
, i.e.self_swap = self.sub(index)
Only if inplace is False and ret_index is True
Notes
The input state and output state have the same number of states, but their ordering is not necessarily the same.
See also
align_phase
rotate states such that their phases align
- align_phase(other: State, ret_index: bool = False, inplace: bool = False)
Align self with the phases for other, a copy may be returned
States will be rotated by \(\pi\) provided the phase difference between the states are above \(|\Delta\theta| > \pi/2\).
- Parameters:
other (
State
) – the other state to align onto this stateret_index – return which indices got swapped
inplace – rotate the states in-place
See also
align_norm
re-order states such that site-norms have a smaller residual
- asCoefficient()
- asState()
- change_gauge(gauge: sisl.typing.GaugeType, offset=(0, 0, 0))
In-place change of the gauge of the state coefficients
The two gauges are related through:
\[\tilde C_\alpha = e^{i\mathbf k\mathbf r_\alpha} C_\alpha\]where \(C_\alpha\) and \(\tilde C_\alpha\) belongs to the
atom
andcell
gauge, respectively.- Parameters:
gauge – specify the new gauge for the mode coefficients
offset (
array_like
, optional) – whether the coordinates should be offset by another phase-factor
- copy() StateC
Return a copy (only the coefficients and states are copied),
parent
andinfo
are passed by reference
- degenerate(atol: float)
Find degenerate coefficients with a specified precision
- Parameters:
atol – the precision above which coefficients are not considered degenerate
- Returns:
list
ofnumpy.ndarray
– a list of indices
- derivative(order: Literal[1, 2] = 1, matrix: bool = False, axes: sisl.typing.CartesianAxes = 'xyz', operator: Callable[[Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], Literal['x', 'y', 'z', 'xx', 'yy', 'zz', 'yz', 'xz', 'xy'] | None], Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes]] = lambda M, d=None: ...)
Calculate the derivative with respect to \(\mathbf k\) for a set of states up to a given order
These are calculated using the analytic expression (\(\alpha\) corresponding to the Cartesian directions), here only shown for the 1st order derivative:
\[\mathbf{d}_{\alpha ij} = \langle \psi_i | \frac{\partial}{\partial\mathbf k_\alpha} \mathbf H(\mathbf k) | \psi_j \rangle\]In case of non-orthogonal basis the equations substitutes \(\mathbf H(\mathbf k)\) by \(\mathbf H(\mathbf k) - \epsilon_i\mathbf S(\mathbf k)\).
The 2nd order derivatives are calculated with the Berry curvature correction:
\[\mathbf d^2_{\alpha \beta ij} = \langle\psi_i| \frac{\partial^2}{\partial\mathbf k_\alpha\partial\mathbf k_\beta} \mathbf H(\mathbf k) | \psi_j\rangle - \frac12\frac{\mathbf{d}_{\alpha ij}\mathbf{d}_{\beta ij}} {\epsilon_i - \epsilon_j}\]Notes
When requesting 2nd derivatives it will not be advisable to use a
sub
before calculating the derivatives since the 1st order perturbation uses the energy differences (Berry contribution) and the 1st derivative matrix for correcting the curvature.For states at the \(\Gamma\) point you may get warnings about casting complex numbers to reals. In these cases you should force the state at the \(\Gamma\) point to be calculated in complex numbers to enable the correct decoupling.
- Parameters:
order – an integer specifying which order of the derivative is being calculated.
matrix – whether the full matrix or only the diagonal components are returned
axes – NOTE: this argument may change in future versions. only calculate the derivative(s) along specified Cartesian directions. The axes argument will be sorted internally, so the order will always be xyz. For the higher order derivatives all those involving only the provided axes will be used.
operator – an operator that translates the \(\delta\) matrices to another operator. The same operator will be applied to both
P
andS
matrices.
See also
SparseOrbitalBZ.dPk
function for generating the matrix derivatives
SparseOrbitalBZ.dSk
function for generating the matrix derivatives in non-orthogonal basis
- Returns:
dv
– the 1st derivative, has shape(3, state.shape[0])
formatrix=False
, else has shape(3, state.shape[0], state.shape[0])
Also returned fororder >= 2
since it is used in the higher order derivativesddv
– the 2nd derivative, has shape(6, state.shape[0])
formatrix=False
, else has shape(6, state.shape[0], state.shape[0])
, the first dimension is in the Voigt representation Only returned fororder >= 2
- inner(ket=None, matrix=None, projection: Literal['diag', 'atoms', 'basis', 'matrix'] = 'diag')
Calculate the inner product as \(\mathbf A_{ij} = \langle\psi_i|\mathbf M|\psi'_j\rangle\)
Inner product calculation allows for a variety of things.
for
matrix
it will compute off-diagonal elements as well
\[\mathbf A_{\alpha\beta} = \langle\psi_\alpha|\mathbf M|\psi'_\beta\rangle\]for
diag
only the diagonal components will be returned
\[\mathbf a_\alpha = \langle\psi_\alpha|\mathbf M|\psi_\alpha\rangle\]for
basis
, only do inner products for individual states, but return them basis-resolved
\[\mathbf A_{\alpha\beta} = \psi^*_{\alpha,\beta} \mathbf M|\psi_\alpha\rangle_\beta\]for
atoms
, only do inner products for individual states, but return them atom-resolved
- Parameters:
ket (
State
, optional) – the ket object to calculate the inner product with, if not passed it will do the inner product with itself. The object itself will always be the bra \(\langle\psi_i|\)matrix (
array_like
, optional) – whether a matrix is sandwiched between the bra and ket, defaults to the identity matrix. 1D arrays will be treated as a diagonal matrix.projection – how to perform the final projection. This can be used to sum specific sub-elements, return the diagonal, or the full matrix.
diag
only return the diagonal of the inner productmatrix
a matrix with diagonals and the off-diagonalsbasis
only do inner products for individual states, but return them basis-resolvedatoms
only do inner products for individual states, but return them atom-resolved
Notes
This does not take into account a possible overlap matrix when non-orthogonal basis sets are used. One have to add the overlap matrix in the matrix argument, if needed.
- Raises:
ValueError – if the number of state coefficients are different for the bra and ket
RuntimeError – if the matrix shapes are incompatible with an atomic resolution conversion
- Returns:
numpy.ndarray
– a matrix with the sum of inner state products
- ipr(q: int = 2)
Calculate the inverse participation ratio (IPR) for arbitrary q values
The inverse participation ratio is defined as
\[I_{q,\alpha} = \frac{\sum_i |\psi_{\alpha,i}|^{2q}}{ \big[\sum_i |\psi_{\alpha,i}|^2\big]^q}\]where \(\alpha\) is the band index and \(i\) is the orbital. The order of the IPR is defaulted to \(q=2\), see (1) for details. The IPR may be used to distinguish Anderson localization and extended states:
\begin{align} \lim_{L\to\infty} I_{2,\alpha} = \left\{ \begin{aligned} 1/L^d &\quad \text{extended state} \\ \text{const.} &\quad \text{localized state} \end{aligned}\right. \end{align}For further details see [7]. Note that for eigen states the IPR reduces to:
\[I_{q,\alpha} = \sum_i |\psi_{\alpha,i}|^{2q}\]since the denominator is \(1^{q} = 1\).
- Parameters:
q – order parameter for the IPR
- iter(asarray: bool = False)
An iterator looping over the states in this system
- Parameters:
asarray (
bool
, optional) – if true the yielded values are the state vectors, i.e. a numpy array. Otherwise an equivalent object is yielded.- Yields:
state (
State
) – a state only containing individual elements, if asarray is falsestate (
numpy.ndarray
) – a state only containing individual elements, if asarray is true
- norm()
Return a vector with the Euclidean norm of each state \(\sqrt{\langle\psi|\psi\rangle}\)
- Returns:
numpy.ndarray
– the Euclidean norm for each state
- norm2(projection: Literal['sum', 'atoms', 'basis'] = 'sum')
Return a vector with the norm of each state \(\langle\psi|\psi\rangle\)
- Parameters:
projection – whether to compute the norm per state as a single number, atom-resolved or per basis dimension.
See also
inner
used method for calculating the squared norm.
- Returns:
numpy.ndarray
– the squared norm for each state
- normalize()
Return a normalized state where each state has \(|\psi|^2=1\)
This is roughly equivalent to:
>>> state = StateC(np.arange(10), 1) >>> n = state.norm() >>> norm_state = StateC(state.state / n.reshape(-1, 1), state.c.copy()) >>> norm_state.c[0] == 1
- Returns:
State
– a new state with all states normalized, otherwise equal to this
- outer(ket=None, matrix=None)
Return the outer product by \(\sum_\alpha|\psi_\alpha\rangle\langle\psi'_\alpha|\)
- Parameters:
ket (
State
, optional) – the ket object to calculate the outer product of, if not passed it will do the outer product with itself. The object itself will always be the bra \(|\psi_\alpha\rangle\)matrix (
array_like
, optional) – whether a matrix is sandwiched between the ket and bra, defaults to the identity matrix. 1D arrays will be treated as a diagonal matrix.
Notes
This does not take into account a possible overlap matrix when non-orthogonal basis sets are used.
- Returns:
numpy.ndarray
– a matrix with the sum of outer state products
- phase(method: Literal['max', 'all'] = 'max', ret_index: bool = False)
Calculate the Euler angle (phase) for the elements of the state, in the range \(]-\pi;\pi]\)
- Parameters:
method (
{'max', 'all'}
) – for max, the phase for the element which has the largest absolute magnitude is returned, for all, all phases are calculatedret_index – return indices for the elements used when
method=='max'
- remove(index: sisl.typing.SimpleIndex, inplace: bool = False) StateC | None
Return a new state without the specified indices
- Parameters:
index – indices that are removed in the returned object
inplace – whether the values will be removed inplace
- Returns:
StateC
– a new state without containing the requested elements, only if inplace is false
- rotate(phi: float = 0.0, individual: bool = False, inplace: bool = False) State | None
Rotate all states to rotate the largest component to be along the angle phi
The states will be rotated according to:
\[\mathbf S' = \mathbf S / \mathbf S^\dagger_{\phi-\mathrm{max}} \exp (i \phi),\]where \(\mathbf S^\dagger_{\phi-\mathrm{max}}\) is the phase of the component with the largest amplitude and \(\phi\) is the angle to align on.
- Parameters:
phi (
float
, optional) – angle to align the state at (in radians), 0 is the positive real axisindividual (
bool
, optional) – whether the rotation is per state, or a single maximum component is chosen.inplace – whether to do the rotation on the object it-self (True), or return a copy with the rotated states (False).
- sort(ascending: bool = True)
Sort and return a new
StateC
by sorting the coefficients (default to ascending)- Parameters:
ascending – sort the contained elements ascending, else they will be sorted descending
- sub(index: sisl.typing.SimpleIndex, inplace: bool = False) StateC | None
Return a new state with only the specified states
- Parameters:
index – indices that are retained in the returned object
inplace – whether the values will be retained inplace
- Returns:
StateC
– a new object with a subset of the states, only if inplace is false
- tile(reps: int, axis: int, normalize: bool = False, offset: float = 0) State
Tile the state vectors for a new supercell
Tiling a state vector makes use of the Bloch factors for a state by utilizing
\[\psi_{\mathbf k}(\mathbf r + \mathbf T) \propto e^{i\mathbf k\cdot \mathbf T}\]where \(\mathbf T = i\mathbf a_0 + j\mathbf a_1 + l\mathbf a_2\). Note that axis selects which of the \(\mathbf a_i\) vectors that are translated and reps corresponds to the \(i\), \(j\) and \(l\) variables. The offset moves the individual states by said amount, i.e. \(i\to i+\mathrm{offset}\).
- Parameters:
reps – number of repetitions along a specific lattice vector
axis – lattice vector to tile along
normalize – whether the states are normalized upon return, may be useful for eigenstates, equivalent to
state.tile().normalize()
offset – the offset for the phase factors
See also
Geometry.tile
,Grid.tile
,Lattice.tile
- translate(isc)
Translate the vectors to a new unit-cell position
The method is thoroughly explained in
tile
while this one only selects the corresponding state vector- Parameters:
isc (
(3,)
) – number of offsets for the statevector
See also
tile
equivalent method for generating more cells simultaneously
- velocity(*args, **kwargs)[source]
Calculate velocity of the modes
This routine calls
derivative
with appropriate arguments (1st order derivative) and returns the velocity for the modes.Note that the coefficients associated with the
ModeCPhonon
must correspond to the energies of the modes.See
derivative
for details and possible arguments. One cannot pass theorder
argument as that is fixed to1
in this call.Notes
The states and energies for the modes may have changed after calling this routine. This is because of the velocity un-folding for degenerate modes. I.e. calling displacement and/or
PDOS
after this method may change the result.See also
derivative
for details of the implementation
- c
- property dkind
The data-type of the state (in str)
- property dtype
Data-type for the state
- info
- property mode
Eigenmodes (states)
- parent
- property shape
Returns the shape of the state
- state