This brings a much needed clean-up of comments and
array indices in the TS approach.
Basically many of the comments were wrong based on
a transposed view of the matrix + phases.

Additionally we had the dm_update routine doing
the transposing of the matrix calculations.
This was very awkward. This is only for the EQ
with

 DM_eq \propto (Gf - Gf^\dagger)

Now we do it a little easier and do it in the
respective add_DM routines, this makes
dm_update much easier to comprehend and it
should also be faster since searching the sparse
elements should be slower than indexing the transpose
array element in add_DM.
This means we have a problem for ts_dq in the MUMPS@k
solver. The difference in minor but may result in
slower dEf convergence.

Cleaned up all phase calculations, now they are
all based on dot_product(k, sc_off) which is easier
to comprehend.

Cleaned up MUMPS scattering matrix calculation. It is
easier to understand now.

For the self-energy calculations we now rely on Siesta
to provide a Hermitian matrix. I.e. we do not calculate
H = (H + H^A) / 2. Generally Siesta behaves very well and
we should trust this. This may introduce slight noises
in the calculations. We could, if large noises occur make
this feature optional. Secondly, the transfer matrix
is not symmetrized, so one could argue that this needs to
be done in any case.
SPEEDUP.

For Hk we also remove this symmetrization and shift the
Fermi-level in one go. This removes a complete loop of
the sparse matrix elements in Hk.
SPEEDUP.

ts_dq calculations now have more is_nan checks and also
prints out the charge accummulated at Ef.
The interpolation from the file also tries to decide
whether the change in Ef co-incides with the sign of dq.
If not the file interpolation will be turned off for
that iteration. Generally a tight dq tolerance requires
a smaller TS.dQ.Factor.

All phase changes have been applied to the tbtrans code.
This also encouraged a small clean-up in the orb_current
calculation, merged some if-statements.

I ran all BTD, LAPACK and MUMPS calculations
for gamma and k-points. All give the same results.