✨ BandsOutput: parse .dat.gnu for cumulative k-path distances

src/qe_tools/outputs/bands.py

@@ -5,13 +5,14 @@
 from typing import Annotated, TextIO
 
 import numpy as np
-from glom import Spec
+from glom import Coalesce, Spec
 
 from dough import Unit
 from dough.outputs import BaseOutput, output_mapping
 
 from .parsers.bands import (
     BandsDatParser,
+    BandsGnuParser,
     BandsRapParser,
     BandsStdoutParser,
 )
@@ -21,10 +22,10 @@
 class _BandsMapping:
     """Typed outputs of a bands.x calculation."""
 
-    number_of_kpoints: Annotated[int, Spec("dat.nks")]
+    number_of_kpoints: Annotated[int, Spec(Coalesce("dat.nks", "gnu.nks"))]
     """Number of k-points along the band-structure path."""
 
-    number_of_bands: Annotated[int, Spec("dat.nbnd")]
+    number_of_bands: Annotated[int, Spec(Coalesce("dat.nbnd", "gnu.nbnd"))]
     """Number of bands written by bands.x."""
 
     k_points: Annotated[np.ndarray, Spec("dat.k_points")]
@@ -35,7 +36,11 @@ class _BandsMapping:
     coordinates for `crystal_b`). bands.x does not transform them.
     """
 
-    eigenvalues: Annotated[np.ndarray, Spec("dat.eigenvalues"), Unit("eV")]
+    eigenvalues: Annotated[
+        np.ndarray,
+        Spec(Coalesce("dat.eigenvalues", "gnu.eigenvalues")),
+        Unit("eV"),
+    ]
     """Kohn-Sham eigenvalues along the band path, in eV.
 
     Numpy array of shape `(n_kpoints, n_bands)`:
@@ -47,6 +52,14 @@ class _BandsMapping:
     array therefore covers a single spin channel.
     """
 
+    k_path_distances: Annotated[np.ndarray, Spec("gnu.k_path_distances")]
+    """Cumulative k-path distance per k-point, in `2π/alat`.
+
+    Numpy array of shape `(n_kpoints,)`. Suitable as the x-axis for a band-structure
+    plot. Parsed from the `*.dat.gnu` file written by bands.x; QE inserts zero-length
+    jumps for path discontinuities, which are preserved here.
+    """
+
     high_symmetry_points: Annotated[np.ndarray, Spec("stdout.high_symmetry_points")]
     """Crystal-momentum coordinates of the high-symmetry points along the path.
 
@@ -85,17 +98,18 @@ class BandsOutput(BaseOutput[_BandsMapping]):
 
     @classmethod
     def from_dir(cls, directory: str | Path):
-        """Locate filband (`*.dat`, `*.dat.rap`) and bands.x stdout in `directory`."""
+        """Locate filband (`*.dat`, `*.dat.gnu`, `*.dat.rap`) and bands.x stdout in `directory`."""
         directory = Path(directory)
 
         if not directory.is_dir():
             raise ValueError(f"Path `{directory}` is not a valid directory.")
 
         rap_file = next(directory.glob("*.dat.rap"), None)
+        gnu_file = next(directory.glob("*.dat.gnu"), None)
 
         dat_file = None
         for candidate in directory.glob("*.dat"):
-            if candidate.name.endswith(".dat.rap"):
+            if candidate.name.endswith((".dat.rap", ".dat.gnu")):
                 continue
             with candidate.open("r") as handle:
                 if "&plot" in handle.readline():
@@ -112,13 +126,16 @@ def from_dir(cls, directory: str | Path):
                 stdout_file = file
                 break
 
-        return cls.from_files(dat=dat_file, rap=rap_file, stdout=stdout_file)
+        return cls.from_files(
+            dat=dat_file, gnu=gnu_file, rap=rap_file, stdout=stdout_file
+        )
 
     @classmethod
     def from_files(
         cls,
         *,
         dat: None | str | Path | TextIO = None,
+        gnu: None | str | Path | TextIO = None,
         rap: None | str | Path | TextIO = None,
         stdout: None | str | Path | TextIO = None,
     ):
@@ -127,6 +144,8 @@ def from_files(
 
         if dat is not None:
             raw_outputs["dat"] = BandsDatParser.parse_from_file(dat)
+        if gnu is not None:
+            raw_outputs["gnu"] = BandsGnuParser.parse_from_file(gnu)
         if rap is not None:
             raw_outputs["rap"] = BandsRapParser.parse_from_file(rap)
         if stdout is not None:

src/qe_tools/outputs/parsers/bands.py

@@ -97,6 +97,46 @@ def parse(content: str) -> dict:
         }
 
 
+class BandsGnuParser(BaseOutputFileParser):
+    """Parse the gnuplot-friendly `.dat.gnu` output of bands.x.
+
+    The `.gnu` file lists `(k_path_distance, eigenvalue)` pairs grouped per band,
+    with bands separated by blank lines. The first column is the cumulative k-path
+    coordinate (in `2π/alat`); QE inserts zero-length jumps for path discontinuities.
+    """
+
+    @staticmethod
+    def parse(content: str) -> dict:
+        blocks = [b for b in (block.strip() for block in content.split("\n\n")) if b]
+        if not blocks:
+            raise ValueError("filband.gnu file is empty.")
+
+        per_band = []
+        for block in blocks:
+            data = np.fromstring(block, sep=" ")
+            if data.size % 2 != 0:
+                raise ValueError(
+                    f"filband.gnu band block has {data.size} numbers; expected an even count."
+                )
+            per_band.append(data.reshape(-1, 2))
+
+        nks = per_band[0].shape[0]
+        if any(block.shape[0] != nks for block in per_band):
+            raise ValueError(
+                "filband.gnu band blocks have inconsistent k-point counts."
+            )
+
+        k_path_distances = per_band[0][:, 0]
+        eigenvalues = np.column_stack([block[:, 1] for block in per_band])
+
+        return {
+            "nbnd": len(per_band),
+            "nks": nks,
+            "k_path_distances": k_path_distances,
+            "eigenvalues": eigenvalues,
+        }
+
+
 class BandsStdoutParser(BaseOutputFileParser):
     """Parse the stdout of bands.x for high-symmetry point markers."""
 

tests/outputs/fixtures/bands/mgo/MgO-bands.dat.gnu

@@ -0,0 +1,14 @@
+    0.0000  -10.0000
+    0.5000  -10.5000
+    1.3660  -12.0000
+    1.8660  -10.8000
+
+    0.0000    2.0000
+    0.5000    1.8000
+    1.3660    1.0000
+    1.8660    1.5000
+
+    0.0000    5.0000
+    0.5000    4.5000
+    1.3660    4.0000
+    1.8660    4.7000

tests/outputs/test_bands_output.py

@@ -14,9 +14,33 @@ def test_bands_mgo(robust_data_regression_check):
         "number_of_bands": out["number_of_bands"],
         "k_points": out["k_points"],
         "eigenvalues": out["eigenvalues"],
+        "k_path_distances": out["k_path_distances"],
         "high_symmetry_points": out["high_symmetry_points"],
         "high_symmetry_distances": out["high_symmetry_distances"],
         "is_high_symmetry": out["is_high_symmetry"].tolist(),
         "representations": out["representations"],
     }
     robust_data_regression_check(snapshot)
+
+
+def test_bands_eigenvalues_gnu_fallback(robust_data_regression_check):
+    """`eigenvalues` should resolve from the `.gnu` file when no `.dat` is provided.
+
+    The `.dat` file is the canonical source for `eigenvalues`, but `BandsOutput` also
+    exposes a `Coalesce` fallback to `.gnu` so that users who only have the gnuplot
+    file can still read eigenvalues.
+    """
+    gnu_file = (
+        Path(__file__).parent / "fixtures" / "bands" / "mgo" / "MgO-bands.dat.gnu"
+    )
+
+    bands = BandsOutput.from_files(gnu=gnu_file)
+
+    out = bands.get_output_dict()
+    snapshot = {
+        "number_of_kpoints": out["number_of_kpoints"],
+        "number_of_bands": out["number_of_bands"],
+        "eigenvalues": out["eigenvalues"],
+        "k_path_distances": out["k_path_distances"],
+    }
+    robust_data_regression_check(snapshot)

tests/outputs/test_bands_output/test_bands_eigenvalues_gnu_fallback.yml

@@ -0,0 +1,20 @@
+eigenvalues:
+- - -10.0
+  - 2.0
+  - 5.0
+- - -10.5
+  - 1.8
+  - 4.5
+- - -12.0
+  - 1.0
+  - 4.0
+- - -10.8
+  - 1.5
+  - 4.7
+k_path_distances:
+- 0.0
+- 0.5
+- 1.366
+- 1.866
+number_of_bands: 3.0
+number_of_kpoints: 4.0

tests/outputs/test_bands_output/test_bands_mgo.yml

@@ -30,6 +30,11 @@ is_high_symmetry:
 - false
 - true
 - true
+k_path_distances:
+- 0.0
+- 0.5
+- 1.366
+- 1.866
 k_points:
 - - 1.0
   - 0.5