zerothi / sisl

Electronic structure Python package for post analysis and large scale tight-binding DFT/NEGF calculations
https://zerothi.github.io/sisl
Mozilla Public License 2.0
181 stars 58 forks source link

Hamiltonian to netCDF without all diagonal elements set #389

Closed tfrederiksen closed 2 years ago

tfrederiksen commented 2 years ago

I encountered some unexpected behaviour trying to write a simple Hamiltonian to netcdf, say with this simple script:

import sisl
g = sisl.geom.graphene()
H = sisl.Hamiltonian(g)
#H[0, 0] = 0
H[1, 1] = 0
H.write('test.nc')

which produces this error:

Traceback (most recent call last):
  File "test_write.py", line 6, in <module>
    H.write('test.nc')
  File "sisl-0.12.0.dev0-py3.8-linux-x86_64.egg/sisl/physics/hamiltonian.py", line 380, in write
    fh.write_hamiltonian(self, *args, **kwargs)
  File "sisl-0.12.0.dev0-py3.8-linux-x86_64.egg/sisl/io/siesta/siesta_nc.py", line 620, in write_hamiltonian
    self._write_overlap(spgroup, csr, H.orthogonal, H.S_idx)
  File "sisl-0.12.0.dev0-py3.8-linux-x86_64.egg/sisl/io/siesta/siesta_nc.py", line 553, in _write_overlap
    raise ValueError(f'{self}._write_overlap '
ValueError: ncSileSiesta(test.nc, base=.)._write_overlap is trying to write an Overlap in Siesta format with
missing diagonal terms on rows []. Please explicitly add *all* diagonal overlap terms.

As is clearly stated in this message, the behaviour is related to the fact that the whole overlap matrix is not explicitly given.

If I uncomment the line in my script above (ie. setting the matrix element also for the first atom) then things are OK and the file is correctly produced.

Now, wouldn't it not be better if sisl handled this situation internally? Is there a reason to force the user to set all diagonal matrix elements even when working with an orthogonal basis?

zerothi commented 2 years ago

True, I think this was a leftover from an old limitation when writing the file... So definitely should be looked into!

zerothi commented 2 years ago

As a side note. This is because tbtrans expects diagonal values. And if you change the sparsity and read it in again, they won't be the same! Quite unexpectedly for some users I guess?

tfrederiksen commented 2 years ago

Could the solution be a flag tbtrans=True (or fill_diag?) for write where eventual missing diagonal elements are added (at the expense of changing sparsity)?

In my case I only want to construct a Hamiltonian once to be able to quickly read it in many subsequent runs. It is no problem if sparsity is modified, but on the other hand it seems unnecessary.

zerothi commented 2 years ago

Could the solution be a flag tbtrans=True (or fill_diag?) for write where eventual missing diagonal elements are added (at the expense of changing sparsity)?

In my case I only want to construct a Hamiltonian once to be able to quickly read it in many subsequent runs. It is no problem if sparsity is modified, but on the other hand it seems unnecessary.

I am somewhat hesitant coming to think of it. The file formats TSHS and nc are generally Siesta/TBtrans specific. And there is no error handling from tbtrans side to check that the diagonal entries are set. Every calculation done with such a file would be wrong. So while the fileformats could be used by 3rd party libraries I think they should still uphold the standard of defining all diagonal entries. This is about safe-guarding users ;)

zerothi commented 2 years ago

If anything it should always change the sparsity pattern to contain the diagonals, without possibilities to not store them.