Fix the basis when reading from XV

[pfebrer]

When reading the geometry from an XV file using fdf.read_geometry(output=True), the basis was not correctly set. It was completely broken, I don't know why I (or anyone) hadn't noticed until now :sweat_smile:

[codecov[bot]]

Codecov Report

Attention: Patch coverage is 92.30769% with 1 lines in your changes are missing coverage. Please review.

Project coverage is 87.32%. Comparing base (087767d) to head (bf0f877).

Files	Patch %	Lines
src/sisl/io/siesta/fdf.py	90.90%	1 Missing :warning:

Additional details and impacted files

```diff @@ Coverage Diff @@ ## main #778 +/- ## ======================================= Coverage 87.32% 87.32% ======================================= Files 396 396 Lines 50579 50575 -4 ======================================= - Hits 44168 44166 -2 + Misses 6411 6409 -2 ```

:umbrella: View full report in Codecov by Sentry.
:loudspeaker: Have feedback on the report? Share it here.

[pfebrer]

Minimal example of how it failed:

XV_fail.zip

From fdf:

python -c 'import sisl; print(sisl.get_sile("siesta.fdf").read_geometry().atoms.atom)'

[<sisl.Atom Au, Z=79, M=196.97, maxR=3.3210000000000024, no=15>, <sisl.Atom O, Z=8, M=16.0, maxR=3.1449000000000105, no=13>]

From XV:

python -c 'import sisl; print(sisl.get_sile("siesta.fdf").read_geometry(output=True).atoms.atom)'

[<sisl.Atom Au, Z=79, M=196.97, maxR=3.3210000000000024, no=15>]

[tfrederiksen]

I can confirm this issue. For one of my systems I find:

# read from fdf
g = sisl.get_sile("./RUN.fdf").read_geometry()
print(g.atoms.atom)
[<sisl.Atom H.opt88, Z=1, M=1.01, maxR=2.8520000000000043, no=5>, 
 <sisl.Atom Au.blyp, Z=79, M=196.97, maxR=3.328400000000018, no=15>,
 <sisl.Atom Aus.blyp, Z=79, M=196.97, maxR=4.7660999999999865, no=16>,
 <sisl.Atom C.blyp, Z=6, M=12.01, maxR=2.954100000000009, no=13>]

# read from XV
g = sisl.get_sile("./RUN.fdf").read_geometry(True)
print(g.atoms.atom)
[<sisl.Atom H.opt88, Z=1, M=1.01, maxR=2.8520000000000043, no=5>,
 <sisl.Atom Li, Z=3, M=6.941, maxR=-1.0, no=1>]

The first case has the correct atoms (and input coordinates), the latter wrong atoms (but correct, relaxed coordinates).

[pfebrer]

Actually the fix here isn't fully correct either.

It seems that using species_as_Z is not a good idea. When reading the basis from the fdf file the atoms list is not sorted by species index but by appearance.

So something like:

In the fdf file

Species:
1 Au
2 O

Coordinates:
X X X 2
X X X 1

In the XV file

Coordinates:
2 8   X X X
1 79 X X X

Will lead to reversed species if you read from fdf.read_geometry(output=True)

So either:

• Species must be sorted when reading the basis from fdf.
• Atoms from the XV file must be retreived from the basis based on their atomic number, not their index.

I think the second one makes more sense, but I don't know if other files that use _replace_with_species also contain the atomic numbers.

[pfebrer]

Hmm getting the basis by Z won't work either if you have multiple atoms with the same Z, as in Thomas' example :sweat_smile:

So the only option is that fdf.read_basis() returns the atoms sorted by species.

[zerothi]

nice catches. Indeed these things should be resolved. Let me have a look

[pfebrer]

Ok, now everything should work fine. You might want to change the name of the argument that I introduced on read_basis (sort_by_fdf_appearance) :sweat_smile:

[zerothi]

Ok, now everything should work fine. You might want to change the name of the argument that I introduced on read_basis (sort_by_fdf_appearance) 😅

I'll probably revert that one, it should be handled differently, the whole idea of using species indices like that would be wrong... Other parts of sisl uses atoms argument to pass a reference basis, here I would do the same (wait for it ;))

Ideally this issues should be resolved with some generic methods that resolves merges of geometries, atoms, orbitals, in a very strict manner, see e.g. #775. This kind of procedure is used all-over, so it could be streamlined quite a bit I would assume.

[zerothi]

Thanks for discovering this!

Could you both try out the #779 branch? It should solve this problem (and more!). Also added some more tests.

[zerothi]

Closed via #779