Change nearest-neighbour analysis and parser

The self-consistent workflow now checks the maximum number of neighbours
and defines the hp inputs so to print all the possible intersites couples. The
parser is then instructed to only keep the interactions that were defined in
the initial HubbardStructureData. This is useful in case the actual first neighbours
are different then the activated intersite interactions, thus a screening
of the HUBBARD.dat file cannot be avoided. This methodology allows for
having a more consistent definition and parsing of the nearest neighbours
using Voronoi tessellation.

pyproject.toml

@@ -27,7 +27,7 @@ keywords = ['aiida', 'workflows']
 requires-python = '>=3.8'
 dependencies = [
     'aiida-core~=2.2',
-    'aiida-quantumespresso~=4.3',
+    'aiida-quantumespresso@git+https://github.com/bastonero/aiida-quantumespresso.git@feat/voronoi',
 ]
 
 [project.urls]

src/aiida_quantumespresso_hp/parsers/hp.py

@@ -5,7 +5,7 @@
 from aiida import orm
 from aiida.common import exceptions
 from aiida.parsers import Parser
-import numpy
+import numpy as np
 
 from aiida_quantumespresso_hp.calculations.hp import HpCalculation
 
@@ -114,8 +114,8 @@ def parse_stdout(self):
             parsed_data, logs = parse_raw_output(stdout)
         except Exception:  # pylint: disable=broad-except
             return self.exit_codes.ERROR_OUTPUT_STDOUT_PARSE
-        else:
-            self.out('parameters', orm.Dict(parsed_data))
+
+        self.out('parameters', orm.Dict(parsed_data))
 
         # If the stdout was incomplete, most likely the job was interrupted before it could cleanly finish, so the
         # output files are most likely corrupt and cannot be restarted from
@@ -202,17 +202,38 @@ def parse_hubbard_dat(self, folder_path):
 
         :return: optional exit code in case of an error
         """
+        from aiida_quantumespresso.common.hubbard import Hubbard
         from aiida_quantumespresso.utils.hubbard import HubbardUtils
         filename = HpCalculation.filename_output_hubbard_dat
 
         filepath = os.path.join(folder_path, filename)
 
         hubbard_structure = self.node.inputs.hubbard_structure.clone()
+
+        intersites = None
+        if 'settings' in self.node.inputs:
+            if 'radial_analysis' in self.node.inputs.settings.get_dict():
+                kwargs = self.node.inputs.settings.dict.radial_analysis
+                intersites = HubbardUtils(hubbard_structure).get_intersites_list(**kwargs)
+
         hubbard_structure.clear_hubbard_parameters()
         hubbard_utils = HubbardUtils(hubbard_structure)
         hubbard_utils.parse_hubbard_dat(filepath=filepath)
 
-        self.out('hubbard_structure', hubbard_utils.hubbard_structure)
+        if intersites is None:
+            self.out('hubbard_structure', hubbard_utils.hubbard_structure)
+        else:
+            hubbard_list = np.array(hubbard_utils.hubbard_structure.hubbard.to_list(), dtype='object')
+            parsed_intersites = hubbard_list[:, [0, 2, 5]].tolist()
+            selected_indices = []
+
+            for i, intersite in enumerate(parsed_intersites):
+                if intersite in intersites:
+                    selected_indices.append(i)
+
+            hubbard = Hubbard.from_list(hubbard_list[selected_indices])
+            hubbard_structure.hubbard = hubbard
+            self.out('hubbard_structure', hubbard_structure)
 
     def get_hubbard_structure(self):
         """Set in output an ``HubbardStructureData`` with standard Hubbard U formulation."""
@@ -264,7 +285,7 @@ def parse_chi_content(self, handle):
         for matrix_name in ('chi0', 'chi'):
             matrix_block = blocks[matrix_name]
             matrix_data = data[matrix_block[0]:matrix_block[1]]
-            matrix = numpy.array(self.parse_hubbard_matrix(matrix_data))
+            matrix = np.array(self.parse_hubbard_matrix(matrix_data))
             result[matrix_name] = matrix
 
         return result
@@ -377,4 +398,4 @@ def parse_hubbard_matrix(data):
         if row:
             matrix.append(row)
 
-        return numpy.array(matrix)
+        return np.array(matrix)

src/aiida_quantumespresso_hp/workflows/hubbard.py

@@ -594,17 +594,20 @@ def run_hp(self):
         inputs = AttributeDict(self.exposed_inputs(HpWorkChain, namespace='hubbard'))
 
         if 'radial_analysis' in self.inputs:
-            from qe_tools import CONSTANTS
-
             kwargs = self.inputs.radial_analysis.get_dict()
             hubbard_utils = HubbardUtils(self.ctx.current_hubbard_structure)
-            radius = hubbard_utils.get_intersites_radius(**kwargs)  # in Angstrom
+            num_neigh = hubbard_utils.get_max_number_of_neighbours(**kwargs)
 
             parameters = inputs.hp.parameters.get_dict()
-            parameters['INPUTHP'].pop('num_neigh', None)
-            parameters['INPUTHP']['rmax'] = radius / CONSTANTS.bohr_to_ang
+            parameters['INPUTHP']['num_neigh'] = num_neigh
+
+            settings = {'radial_analysis': self.inputs.radial_analysis.get_dict()}
+            if 'settings' in inputs.hp:
+                settings = inputs.hp.settings.get_dict()
+                settings['radial_analysis'] = self.inputs.radial_analysis.get_dict()
 
             inputs.hp.parameters = orm.Dict(parameters)
+            inputs.hp.settings = orm.Dict(settings)
 
         inputs.clean_workdir = self.inputs.clean_workdir
         inputs.hp.parent_scf = workchain.outputs.remote_folder

tests/fixtures/structures/MnCoS.cif

@@ -0,0 +1,28 @@
+data_image0
+_chemical_formula_structural       MnCoS
+_chemical_formula_sum              "Mn1 Co1 S1"
+_cell_length_a       4.02254
+_cell_length_b       4.02254
+_cell_length_c       4.02254
+_cell_angle_alpha    60
+_cell_angle_beta     60
+_cell_angle_gamma    60
+
+_space_group_name_H-M_alt    "P 1"
+_space_group_IT_number       1
+
+loop_
+  _space_group_symop_operation_xyz
+  'x, y, z'
+
+loop_
+  _atom_site_type_symbol
+  _atom_site_label
+  _atom_site_symmetry_multiplicity
+  _atom_site_fract_x
+  _atom_site_fract_y
+  _atom_site_fract_z
+  _atom_site_occupancy
+  Mn  Mn1       1.0  0.00000  0.00000  0.00000  1.0000
+  Co  Co1       1.0  0.75000  0.75000  0.75000  1.0000
+  S   S1        1.0  0.50000  0.50000  0.50000  1.0000

tests/parsers/fixtures/hp/voronoi_hubbard_structure/HUBBARD.dat

@@ -0,0 +1,24 @@
+# Copy this data in the pw.x input file for DFT+Hubbard calculations
+HUBBARD	{ortho-atomic}
+V    Mn-3d     Mn-3d    1   1    5.3443
+V    Mn-3d     Co-3d    1   14   0.4967
+V    Mn-3d     Co-3d    1   8    0.4967
+V    Mn-3d     Co-3d    1   5    0.4967
+V    Mn-3d     Co-3d    1   32   0.4967
+V    Mn-3d     S-3p     1   42   0.0256
+V    Mn-3d     S-3p     1   36   0.0256
+V    Mn-3d     S-3p     1   18   0.0256
+V    Mn-3d     S-3p     1   33   0.0256
+V    Mn-3d     S-3p     1   15   0.0256
+V    Mn-3d     S-3p     1   9    0.0256
+V    Co-3d     Co-3d    2   2    6.6772
+V    Co-3d     Mn-3d    2   79   0.4967
+V    Co-3d     S-3p     2   45   0.1165
+V    Co-3d     S-3p     2   69   0.1165
+V    Co-3d     S-3p     2   51   0.1165
+V    Co-3d     Mn-3d    2   76   0.4967
+V    Co-3d     Mn-3d    2   70   0.4967
+V    Co-3d     Mn-3d    2   52   0.4967
+V    Co-3d     S-3p     2   3    0.1165
+V    Co-3d     Co-3d    2   17   0.0948
+V    Co-3d     Co-3d    2   44   0.0948