From 384b443b1415a1fe4d09d2499708700929eedd5f Mon Sep 17 00:00:00 2001 From: Harrison LaBollita Date: Fri, 4 Sep 2026 10:14:16 -0400 Subject: [PATCH 1/2] [wien2k] Store high-symmetry k-path labels in dft_bands_input Read the high-symmetry point block that dmftproj appends to case.outband and store it as kpts_labels and kpts_labels_idx, the same keys and format as the VASP KPOINTS_OPT conversion. The keys are only written when the block is present, so older files convert unchanged. --- doc/ChangeLog.md | 3 ++ doc/h5structure.rst | 9 +++++ python/triqs_dftkit/wien2k/converter.py | 51 +++++++++++++++++++++++++ 3 files changed, 63 insertions(+) diff --git a/doc/ChangeLog.md b/doc/ChangeLog.md index f70657d..0d15cc2 100644 --- a/doc/ChangeLog.md +++ b/doc/ChangeLog.md @@ -32,6 +32,9 @@ Find below an itemized list of changes in this release. * Add ABINIT support to `Wannier90Converter` for charge self-consistent calculations * Read VASP `ICHARG=5` miscellaneous input from `vaspout.h5` +### Wien2k +* Read the high-symmetry k-path labels from the end of `case.outband` and store them as `kpts_labels` / `kpts_labels_idx` in `dft_bands_input`, matching the VASP band conversion + ### Fix * Fix a bug in the `deltaN` write for the Quantum Espresso and Abinit interfaces diff --git a/doc/h5structure.rst b/doc/h5structure.rst index 0a9349c..1a429e1 100644 --- a/doc/h5structure.rst +++ b/doc/h5structure.rst @@ -603,6 +603,15 @@ counting the k-points along the path. ``max(n_parproj)``, ``max(shells['dim'])``, ``max(n_orbitals)``] - As in ``dft_parproj_input``, along the path. Elk writes a dummy ``array([0])``. + * - ``kpts_labels`` + - list of string + - Names of the high-symmetry points of the path (e.g. ``'GAMMA'``, ``'X'``), + for labelling the ticks of a band plot. Only written by Wien2k and VASP, + and only when the underlying DFT output provides the labels. + * - ``kpts_labels_idx`` + - numpy.array.int, dim [``len(kpts_labels)``] + - Position of each entry of ``kpts_labels`` along the path, as a 0-based + index into the ``n_k`` k-points. .. note:: diff --git a/python/triqs_dftkit/wien2k/converter.py b/python/triqs_dftkit/wien2k/converter.py index 2a6ac9a..4b0625f 100644 --- a/python/triqs_dftkit/wien2k/converter.py +++ b/python/triqs_dftkit/wien2k/converter.py @@ -28,6 +28,7 @@ from h5 import * from ..converter_tools import * import os.path +import re class Converter(ConverterTools): @@ -396,6 +397,51 @@ def convert_bands_input(self): if not (mpi.is_master_node()): return + def _read_kpath_labels(band_file, n_k): + """ + Read the high-symmetry k-path labels appended to the end of + case.outband and map them onto the flattened band k-point index. + + dmftproj writes one line per high-symmetry point after the projector + data, in fixed Fortran format (2i6,a): a running counter, the 1-based + position of the point along the band path, and the label, e.g. + + 1 1GAMMA + 2 122X + + Returns (labels, idx) where labels is a list of label strings and idx + is a 0-based numpy int array giving, for each label, the position of + that high-symmetry point in the n_k band path. Returns (None, None) + if no label block is present. + """ + with open(band_file, 'r') as R: + lines = R.readlines() + while lines and not lines[-1].strip(): + lines.pop() + + # Walk backwards from the end of the file: the label block is the + # trailing run of lines matching 'counter index label'. + labels = [] + idx = [] + for line in reversed(lines): + match = re.match(r'\s*(\d+)\s+(\d+)\s*([A-Za-z]\S*)\s*$', line) + if match is None: + break + labels.append(match.group(3)) + idx.append(int(match.group(2)) - 1) + labels.reverse() + idx.reverse() + + if not labels: + return None, None + + idx = numpy.array(idx, dtype=int) + if idx[0] < 0 or idx[-1] >= n_k or numpy.any(numpy.diff(idx) <= 0): + mpi.report("convert_bands_input : WARNING : inconsistent high-symmetry point indices in %s; skipping k-path labels." % band_file) + return None, None + + return labels, idx + try: # get needed data from hdf file with HDFArchive(self.hdf_file, 'a') as ar: @@ -483,6 +529,8 @@ def convert_bands_input(self): # Reading done! + kpts_labels, kpts_labels_idx = _read_kpath_labels(self.band_file, n_k) + # Save it to the HDF: with HDFArchive(self.hdf_file, 'a') as ar: if not (self.bands_subgrp in ar): @@ -491,6 +539,9 @@ def convert_bands_input(self): # created. If it exists, the data is overwritten! things_to_save = ['n_k', 'n_orbitals', 'proj_mat', 'hopping', 'n_parproj', 'proj_mat_all'] + if kpts_labels is not None: + things_to_save += ['kpts_labels', 'kpts_labels_idx'] + mpi.report(" Stored %i high-symmetry k-path labels: %s" % (len(kpts_labels), ', '.join(kpts_labels))) for it in things_to_save: ar[self.bands_subgrp][it] = locals()[it] From 37cbfc0940479f1b539cde0f5ffb033066496aba Mon Sep 17 00:00:00 2001 From: Harrison LaBollita Date: Mon, 7 Sep 2026 10:06:02 -0400 Subject: [PATCH 2/2] [wien2k] Harden the case.outband k-path label parser and test it --- doc/ChangeLog.md | 2 +- doc/h5structure.rst | 4 +- python/triqs_dftkit/wien2k/converter.py | 131 +++++++++++++++--------- test/python/wien2k/CMakeLists.txt | 6 ++ test/python/wien2k/test_kpath_labels.py | 74 +++++++++++++ 5 files changed, 169 insertions(+), 48 deletions(-) create mode 100644 test/python/wien2k/test_kpath_labels.py diff --git a/doc/ChangeLog.md b/doc/ChangeLog.md index 0d15cc2..bc151f3 100644 --- a/doc/ChangeLog.md +++ b/doc/ChangeLog.md @@ -27,10 +27,10 @@ Find below an itemized list of changes in this release. * Add a VASP driver for charge self-consistent DFT+DMFT calculations * Add `KPOINTS_OPT` band conversion from `vaspout.h5`: when `LOCPROJ_OPT` data are available, the converter writes `dft_bands_input` for band/spectral workflows, applying the same PLO config settings (`EWINDOW`, `TRANSFORM`, `NORMALIZE`, and optional `EFERMI`) as the regular VASP conversion, and stores the high-symmetry k-path labels * Warn on misplaced or unknown tags in the PLOVASP configuration +* Read VASP `ICHARG=5` miscellaneous input from `vaspout.h5` ### Wannier90 * Add ABINIT support to `Wannier90Converter` for charge self-consistent calculations -* Read VASP `ICHARG=5` miscellaneous input from `vaspout.h5` ### Wien2k * Read the high-symmetry k-path labels from the end of `case.outband` and store them as `kpts_labels` / `kpts_labels_idx` in `dft_bands_input`, matching the VASP band conversion diff --git a/doc/h5structure.rst b/doc/h5structure.rst index 1a429e1..b12206b 100644 --- a/doc/h5structure.rst +++ b/doc/h5structure.rst @@ -611,7 +611,9 @@ counting the k-points along the path. * - ``kpts_labels_idx`` - numpy.array.int, dim [``len(kpts_labels)``] - Position of each entry of ``kpts_labels`` along the path, as a 0-based - index into the ``n_k`` k-points. + index into the ``n_k`` k-points. Each converter validates these indices + as far as its own DFT output allows, so consumers should not rely on + them being sorted or in range without checking. .. note:: diff --git a/python/triqs_dftkit/wien2k/converter.py b/python/triqs_dftkit/wien2k/converter.py index 4b0625f..37e25f3 100644 --- a/python/triqs_dftkit/wien2k/converter.py +++ b/python/triqs_dftkit/wien2k/converter.py @@ -28,7 +28,91 @@ from h5 import * from ..converter_tools import * import os.path -import re + + +def _parse_kpath_label_line(line): + """ + Parse one line of the high-symmetry point block that dmftproj appends to + case.outband. The format is the fixed Fortran (2i6,a): columns 1-6 hold a + running counter, columns 7-12 the 1-based position of the point along the + band path, and columns 13 onwards the label. + + Returns (counter, pos, label), or None if the line does not have that form. + """ + line = line.rstrip('\n') + if len(line) < 13: + return None + try: + counter = int(line[0:6]) + pos = int(line[6:12]) + except ValueError: + return None + label = line[12:].strip() + if not label: + return None + return counter, pos, label + + +def _read_kpath_labels(band_file, n_k): + """ + Read the high-symmetry k-path labels appended to the end of case.outband + and map them onto the flattened band k-point index. + + dmftproj writes one line per high-symmetry point after the projector data, + in fixed Fortran format (2i6,a), e.g. + + 1 1GAMMA + 2 122X + + Returns (labels, idx) where labels is a list of label strings and idx is a + 0-based numpy int array giving, for each label, the position of that + high-symmetry point in the n_k band path. Returns (None, None) if no label + block is present or if the block is incomplete. + """ + # The block is a handful of lines at the very end of a file that holds all + # the projectors, so read the tail rather than the whole file. + n_bytes = 8192 + with open(band_file, 'rb') as R: + R.seek(0, os.SEEK_END) + size = R.tell() + R.seek(max(0, size - n_bytes)) + lines = R.read().decode('utf-8', 'replace').splitlines() + if size > n_bytes: + # the first line of the chunk is in general cut in the middle + lines = lines[1:] + while lines and not lines[-1].strip(): + lines.pop() + + # Walk backwards from the end of the file: the label block is the trailing + # run of (2i6,a) lines, and the counter of its first line is 1. + labels = [] + idx = [] + first_counter = None + for line in reversed(lines): + parsed = _parse_kpath_label_line(line) + if parsed is None: + break + first_counter, pos, label = parsed + labels.append(label) + idx.append(pos - 1) + if first_counter == 1: + break + labels.reverse() + idx.reverse() + + if not labels: + return None, None + + if first_counter != 1: + mpi.report("convert_bands_input : WARNING : the high-symmetry point block in %s does not start at counter 1, so it is truncated or corrupted; skipping k-path labels." % band_file) + return None, None + + idx = numpy.array(idx, dtype=int) + if idx[0] < 0 or idx[-1] >= n_k or numpy.any(numpy.diff(idx) < 0): + mpi.report("convert_bands_input : WARNING : inconsistent high-symmetry point indices in %s; skipping k-path labels." % band_file) + return None, None + + return labels, idx class Converter(ConverterTools): @@ -397,51 +481,6 @@ def convert_bands_input(self): if not (mpi.is_master_node()): return - def _read_kpath_labels(band_file, n_k): - """ - Read the high-symmetry k-path labels appended to the end of - case.outband and map them onto the flattened band k-point index. - - dmftproj writes one line per high-symmetry point after the projector - data, in fixed Fortran format (2i6,a): a running counter, the 1-based - position of the point along the band path, and the label, e.g. - - 1 1GAMMA - 2 122X - - Returns (labels, idx) where labels is a list of label strings and idx - is a 0-based numpy int array giving, for each label, the position of - that high-symmetry point in the n_k band path. Returns (None, None) - if no label block is present. - """ - with open(band_file, 'r') as R: - lines = R.readlines() - while lines and not lines[-1].strip(): - lines.pop() - - # Walk backwards from the end of the file: the label block is the - # trailing run of lines matching 'counter index label'. - labels = [] - idx = [] - for line in reversed(lines): - match = re.match(r'\s*(\d+)\s+(\d+)\s*([A-Za-z]\S*)\s*$', line) - if match is None: - break - labels.append(match.group(3)) - idx.append(int(match.group(2)) - 1) - labels.reverse() - idx.reverse() - - if not labels: - return None, None - - idx = numpy.array(idx, dtype=int) - if idx[0] < 0 or idx[-1] >= n_k or numpy.any(numpy.diff(idx) <= 0): - mpi.report("convert_bands_input : WARNING : inconsistent high-symmetry point indices in %s; skipping k-path labels." % band_file) - return None, None - - return labels, idx - try: # get needed data from hdf file with HDFArchive(self.hdf_file, 'a') as ar: diff --git a/test/python/wien2k/CMakeLists.txt b/test/python/wien2k/CMakeLists.txt index d9968e3..05cd5ad 100644 --- a/test/python/wien2k/CMakeLists.txt +++ b/test/python/wien2k/CMakeLists.txt @@ -11,3 +11,9 @@ add_test(NAME Py_wien2k_convert WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}) set_property(TEST Py_wien2k_convert APPEND PROPERTY ENVIRONMENT PYTHONPATH=${PROJECT_BINARY_DIR}/python:$ENV{PYTHONPATH} ${SANITIZER_RT_PRELOAD}) + +add_test(NAME Py_wien2k_kpath_labels + COMMAND ${TRIQS_PYTHON_EXECUTABLE} ${CMAKE_CURRENT_SOURCE_DIR}/test_kpath_labels.py + WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}) +set_property(TEST Py_wien2k_kpath_labels APPEND PROPERTY ENVIRONMENT + PYTHONPATH=${PROJECT_BINARY_DIR}/python:$ENV{PYTHONPATH} ${SANITIZER_RT_PRELOAD}) diff --git a/test/python/wien2k/test_kpath_labels.py b/test/python/wien2k/test_kpath_labels.py new file mode 100644 index 0000000..b2f1f45 --- /dev/null +++ b/test/python/wien2k/test_kpath_labels.py @@ -0,0 +1,74 @@ +################################################################################ +# +# TRIQS: a Toolbox for Research in Interacting Quantum Systems +# +# Copyright (C) 2011 by M. Aichhorn, L. Pourovskii, V. Vildosola +# +# TRIQS is free software: you can redistribute it and/or modify it under the +# terms of the GNU General Public License as published by the Free Software +# Foundation, either version 3 of the License, or (at your option) any later +# version. +# +# TRIQS is distributed in the hope that it will be useful, but WITHOUT ANY +# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS +# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more +# details. +# +# You should have received a copy of the GNU General Public License along with +# TRIQS. If not, see . +# +################################################################################ + +""" +Test the high-symmetry k-path label block that dmftproj appends to case.outband. +A real case.outband is far too large to ship as a fixture, so this writes a +synthetic tail in the exact Fortran (2i6,a) format of outband.f:277, including +the trailing blank the `a` descriptor emits for the declared length of kname. +""" + +import os +import tempfile + +import numpy as np + +from triqs_dftkit.wien2k.converter import _read_kpath_labels + +N_K = 501 + + +def line(counter, pos, label): + return '%6d%6d%s' % (counter, pos, label.ljust(10)) + + +def read(lines, n_k=N_K): + fd, path = tempfile.mkstemp(suffix='.outband') + try: + with os.fdopen(fd, 'w') as f: + f.write('\n'.join(lines) + '\n') + return _read_kpath_labels(path, n_k) + finally: + os.remove(path) + + +# One block exercising the awkward cases at once: '\xG' does not start with a +# letter (XCrySDen writes Gamma that way), X|Y puts two labels on the same +# k-point, the projector line in front of the block itself parses as (2i6,a), +# the padding pushes the block past the tail that is read back, and the file +# ends on blank lines. +padding = ['%20.14f%20.14f' % (0.5, 0.25)] * 2000 +block = [line(1, 1, '\\xG'), line(2, 122, 'X'), line(3, 122, 'Y'), + line(4, 268, 'L'), line(5, 501, 'K')] +labels, idx = read(padding + [line(3, 7, '0.52')] + block + ['', ' ']) + +assert labels == ['\\xG', 'X', 'Y', 'L', 'K'], labels +np.testing.assert_array_equal(idx, [0, 121, 121, 267, 500]) + +# A malformed line in the middle of the block truncates the backward walk, so +# the block no longer starts at counter 1 and is rejected rather than silently +# returning only its tail. +broken = list(block) +broken[2] = 'this line is not a label line' +assert read(padding + broken) == (None, None) + +# A position past the end of the band path is rejected as well. +assert read(padding + block, n_k=400) == (None, None)