🐛 ProjwfcParser: Fix DOS parsing for spin-polarised case

Currently, the `ProjwfcParser` does not correctly parse the total DOS data in the case of spin-polarised calculations. Only the first column is parsed, which corresponds to spin-up density only. This issue passed the tests unnoticed since the regression tests only check the attributes of the `Dos` output, not the actual content. Here we adapt the `ProjwfcParser` to treat the spin-polarised case `nspin = 2` case correctly, constructing an `XyData` node that has two `y` arrays, one for each spin. The tests for the spin-polarised case are also adapted to add the `XyData` contents to the regression files.
aiidateam · Jan 30, 2025 · bcd5508 · bcd5508
1 parent 29f9c98
commit bcd5508
Show file tree

Hide file tree

Showing 3 changed files with 192 additions and 23 deletions.
diff --git a/src/aiida_quantumespresso/parsers/projwfc.py b/src/aiida_quantumespresso/parsers/projwfc.py
@@ -63,19 +63,16 @@ def parse(self, **kwargs):
 
         orbitals = self._parse_orbitals(header, structure, non_collinear, spinorbit)
         bands, projections = self._parse_bands_and_projections(kpoint_blocks, len(orbitals))
-        energy, dos, pdos_array = self._parse_pdos_files(retrieved_temporary_folder, nspin, spinorbit, logs)
+        energy, dos_node, pdos_array = self._parse_pdos_files(retrieved_temporary_folder, nspin, spinorbit, logs)
+
+        self.out('Dos', dos_node)
 
         output_node_dict = self._build_bands_and_projections(
             kpoints, bands, energy, orbitals, projections, pdos_array, nspin
         )
         for linkname, node in output_node_dict.items():
             self.out(linkname, node)
 
-        dos_out = XyData()
-        dos_out.set_x(energy, 'Energy', 'eV')
-        dos_out.set_y(dos, 'Dos', 'states/eV')
-        self.out('Dos', dos_out)
-
         for exit_code in [
             'ERROR_OUTPUT_STDOUT_MISSING',
             'ERROR_OUTPUT_STDOUT_READ',
@@ -261,7 +258,7 @@ def _parse_bands_and_projections(kpoint_blocks: list, num_orbitals: int) -> Tupl
         return bands, projections
 
     def _parse_pdos_files(self, retrieved_temporary_folder: Path, nspin: int,
-                          spinorbit: bool, logs: AttributeDict) -> Tuple[ArrayLike, ArrayLike, ArrayLike]:
+                          spinorbit: bool, logs: AttributeDict) -> Tuple[ArrayLike, XyData, ArrayLike]:
         """Parse the PDOS files and convert them into arrays.
 
         Reads in all of the ``*.pdos*`` files and converts the data into arrays. The PDOS arrays are then concatenated
@@ -276,7 +273,7 @@ def _parse_pdos_files(self, retrieved_temporary_folder: Path, nspin: int,
         :param nspin: nspin value of the parent calculation.
         :param spinorbit: True if the calculation used spin-orbit coupling.
 
-        :return: tuple of three arrays containing the energy grid, the total DOS and the PDOS
+        :return: tuple of three containing the energy grid, the total DOS as a node and the PDOS
         """
 
         def natural_sort_key(sort_key, _nsre=re.compile('([0-9]+)')):
@@ -300,11 +297,24 @@ def natural_sort_key(sort_key, _nsre=re.compile('([0-9]+)')):
             with pdostot_filepath.open('r') as pdostot_file:
                 # Columns: Energy(eV), Ldos, Pdos
                 pdostot_array = np.atleast_2d(np.genfromtxt(pdostot_file))
-                energy = pdostot_array[:, 0]
-                dos = pdostot_array[:, 1]
         except (OSError, KeyError):
             logs.error.append('ERROR_READING_PDOSTOT_FILE')
-            return np.array([]), np.array([]), np.array([])
+            return np.array([]), XyData(), np.array([])
+
+        energy = pdostot_array[:, 0]
+
+        dos_node = XyData()
+        dos_node.set_x(energy, 'Energy', 'eV')
+
+        # For spin-polarised calculations the total DOS is split up in spin up and down
+        if nspin == 2:
+            dos_node.set_y(
+                (pdostot_array[:, 1], pdostot_array[:, 2]),
+                ('dos_up', 'dos_down'),
+                ('states/eV', 'states/eV')
+            )
+        else:
+            dos_node.set_y(pdostot_array[:, 1], 'Dos', 'states/eV')
 
         # We're only interested in the PDOS, so we skip the first columns corresponding to the energy and LDOS
         if nspin == 1 or spinorbit:
@@ -340,7 +350,7 @@ def natural_sort_key(sort_key, _nsre=re.compile('([0-9]+)')):
                 pdos_list.append(np.concatenate([wfc_pdos_array[:, 0::2], wfc_pdos_array[:, 1::2]], axis=1))
             pdos_array = np.concatenate(pdos_list, axis=1)
 
-        return energy, dos, pdos_array
+        return energy, dos_node, pdos_array
 
     @classmethod
     def _build_bands_and_projections(

diff --git a/tests/parsers/test_projwfc.py b/tests/parsers/test_projwfc.py
@@ -65,7 +65,8 @@ def test_projwfc_spinpolarised(generate_projwfc_node, generate_parser, data_regr
         assert link_name in results, list(results.keys())
 
     data_regression.check({
-        'Dos': results['Dos'].base.attributes.all,
+        'Dos':
+        {array_name: results['Dos'].get_array(array_name).tolist() for array_name in results['Dos'].get_arraynames()},
         'bands_up': results['bands_up'].base.attributes.all,
         'bands_down': results['bands_down'].base.attributes.all,
         'projections_up': {

diff --git a/tests/parsers/test_projwfc/test_projwfc_spinpolarised.yml b/tests/parsers/test_projwfc/test_projwfc_spinpolarised.yml
@@ -1,14 +1,172 @@
 Dos:
-  array|x_array:
-  - 55
-  array|y_array_0:
-  - 55
-  x_name: Energy
-  x_units: eV
-  y_names:
-  - Dos
-  y_units:
-  - states/eV
+  x_array:
+  - -70.062
+  - -68.062
+  - -66.062
+  - -64.062
+  - -62.062
+  - -60.062
+  - -58.062
+  - -56.062
+  - -54.062
+  - -52.062
+  - -50.062
+  - -48.062
+  - -46.062
+  - -44.062
+  - -42.062
+  - -40.062
+  - -38.062
+  - -36.062
+  - -34.062
+  - -32.062
+  - -30.062
+  - -28.062
+  - -26.062
+  - -24.062
+  - -22.062
+  - -20.062
+  - -18.062
+  - -16.062
+  - -14.062
+  - -12.062
+  - -10.062
+  - -8.062
+  - -6.062
+  - -4.062
+  - -2.062
+  - -0.062
+  - 1.938
+  - 3.938
+  - 5.938
+  - 7.938
+  - 9.938
+  - 11.938
+  - 13.938
+  - 15.938
+  - 17.938
+  - 19.938
+  - 21.938
+  - 23.938
+  - 25.938
+  - 27.938
+  - 29.938
+  - 31.938
+  - 33.938
+  - 35.938
+  - 37.938
+  y_array_0:
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 6.48
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.256
+  - 0.256
+  - 0.518
+  - 4.55
+  - 1.58
+  - 0.108
+  - 0.108
+  - 0.108
+  - 0.108
+  - 0.278
+  - 0.278
+  - 0.278
+  - 1.8
+  - 0.0
+  - 0.0
+  y_array_1:
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.0
+  - 0.238
+  - 0.238
+  - 0.438
+  - 0.438
+  - 1.42
+  - 0.131
+  - 0.131
+  - 0.131
+  - 0.131
+  - 0.264
+  - 0.264
+  - 0.264
+  - 2.15
+  - 0.0
 bands_down:
   array|bands:
   - 2