Behaviour of `kwant.operator` is inconsistent
Thinking about this a bit more, I am not sure that the behaviour of the operators in kwant.operator
is consistent.
The following should example should illustrate what I mean.
Let us first consider an operator representing the "charge":
Q = kwant.operator.LocalOperator(sys)
This thing can be used with wavefunctions like so:
psi, phi = get_wavefunctions()
result = Q(phi, psi)
and the result is a vector of the length of the number of sites in the system.
This means that Q
actually represents some 3-tensor Q_iαβ
where
Greek indices run over the Hilbert space (orbitals) and Latin indices
run over the sites.
The code Q(phi, psi)
therefore calculates \sum_αβ Q_iαβ phi_α psi_β
.
However, if the act
method is used:
result_2 = Q.act(psi)
Then what is calculated is actually \sum_iβ Q_iαβ psi_β
.
This means that
np.dot(phi.conj(), Q.act(psi))
is actually equal to
sum(Q(phi, psi))
rather than simply
Q(phi, psi)
as one may have expected.
This feels inconsistent and is (I now realise) why I had trouble giving @anton-akhmerov a sparse matrix representation of the operator: we actually have 3-tensors.
Maybe the act
method should instead calculate \sum_β Q_iαβ psi_β
and the
result should be a sparse matrix...?