sisl: Kernel sometimes dies when reading the SIESTA output hamiltonian.

Describe the issue

For electronic structure calculations, I run SIESTA with the following lines in the input:

CDF.Save true CDF.Compress 9 SaveHS true SaveRho true

Afterwards, if the calculations finish correctly, I use the following commands to make the post-processing in SISL:

import sisl
import numpy as np

H = sisl.get_sile('input.fdf').read_hamiltonian()

It always worked correctly so far, but I am working with a molecule with relatively few atoms (56, carbons and hydrogens), SIESTA finishes the job but SISL fails to read the hamiltonian. In jupyter notebook, it says: “The kernel appears to have died. It will restart automatically.” In python 3.9: “STOP Error in reading data, not allocated”.

About this issue

Original URL
State: closed
Created 2 years ago
Comments: 17 (11 by maintainers)

Commits related to this issue

bug: fixed #512 due to compiler bug A compiler bug has been found when doing: Gamma .neqv. no_u == no_s I have tried to create a small report: program t integer :: n, m logical :: Gam call... — committed to nils-wittemeier/sisl by zerothi 2 years ago

Most upvoted comments

Ok, this is really weird.

The problem stems from https://github.com/zerothi/sisl/blob/730155dd3fd13970937fa7f116a5f6eb6a880e22/sisl/io/siesta/_src/hsx_read.f90#L372

However, that line is coded correct. The problem seems to be a compiler bug that I can also reproduce with GCC 11.2.0.

For instance, when adding this line:

  print *, Gamma, no_u == no_s, Gamma .neqv. (no_u == no_s)

I get the output

T T T

which should never occur (neqv means logical non-equivalence). I have tried to create a small minimal reproducing code, but I can’t…

program t

  integer :: n, m
  logical :: Gam

  call wrt()

  open(238, form="unformatted", status="old")
  call rd(Gam, n, m)
  print *, Gam, n == m, Gam .neqv. n == m
  call rd(Gam, n, m)
  print *, Gam, n == m, Gam .neqv. n == m
  call rd(Gam, n, m)
  print *, Gam, n == m, Gam .neqv. n == m
  call rd(Gam, n, m)
  print *, Gam, n == m, Gam .neqv. n == m

  close(238)
  
contains

  subroutine wrt()

    open(238, form="unformatted")
    write(238) .true.
    write(238) 14, 352
    write(238) .true.
    write(238) 14, 14
    write(238) .false.
    write(238) 14, 352
    write(238) .false.
    write(238) 14, 14
    close(238)
    
  end subroutine wrt

  subroutine rd(Gam, n, m)
    logical, intent(out) :: Gam
    integer, intent(out) :: n, m

    read(238) Gam
    read(238) n, m

  end subroutine rd
    
end program t

I think I’ll have to report this upstream… But this is just weird…

I have now tried to circumvent this by spelling out the checks, I could easily read your output file.

In any case, I would still recommend a different format than the HSX files. 😃

Thanks for reporting this!

zerothi on Nov 15, 2022

Did you try executing the above lines as a stand-alone script (instead of in a jupyter-notebook)?

tfrederiksen on Nov 14, 2022