Kpoints line mode for VASP

pyiron/atomistics/structure/atoms.py

@@ -94,6 +94,7 @@ def __init__(
         elements=None,
         dimension=None,
         species=None,
+        high_symmetry_points=None,
         **qwargs
     ):
         if symbols is not None:
@@ -227,6 +228,7 @@ def __init__(
             else:
                 self.pbc = pbc
         self.set_initial_magnetic_moments(magmoms)
+        self._high_symmetry_points = high_symmetry_points
 
     @property
     def cell(self):
@@ -306,6 +308,24 @@ def elements(self):
         """
         return np.array([self.species[el] for el in self.indices])
 
+    def get_high_symmetry_points(self):
+        """
+        Returns:
+            dict: high_symmetry_points
+        """
+        return self._high_symmetry_points
+
+    def set_high_symmetry_points(self, new_high_symmetry_points):
+        """
+        Sets new high symmetry points dictionary.
+
+        Args:
+            new_high_symmetry_points (dict): new high symmetry points
+        """
+        if not isinstance(new_high_symmetry_points, dict):
+            raise ValueError("has to be dict!")
+        self._high_symmetry_points = new_high_symmetry_points
+
     def new_array(self, name, a, dtype=None, shape=None):
         """
         Adding a new array to the instance. This function is for the purpose of compatibility with the ASE package
@@ -467,6 +487,9 @@ def to_hdf(self, hdf, group_name="structure"):
 
             # print ('time in atoms.to_hdf: ', time.time() - time_start)
 
+            if self._high_symmetry_points is not None:
+                hdf_structure["high_symmetry_points"] = self._high_symmetry_points
+
     def from_hdf(self, hdf, group_name="structure"):
         """
         Retrieve the object from a HDF5 file
@@ -544,6 +567,10 @@ def from_hdf(self, hdf, group_name="structure"):
 
                 if "bonds" in hdf_atoms.list_nodes():
                     self.bonds = hdf_atoms["explicit_bonds"]
+
+                self._high_symmetry_points = None
+                if "high_symmetry_points" in hdf_atoms.list_nodes():
+                    self._high_symmetry_points = hdf_atoms["high_symmetry_points"]
                 return self
 
         else:
@@ -603,6 +630,10 @@ def _from_hdf_old(self, hdf, group_name="structure"):
 
             if "bonds" in hdf_atoms.list_nodes():
                 self.bonds = hdf_atoms["explicit_bonds"]
+
+            self._high_symmetry_points = None
+            if "high_symmetry_points" in hdf_atoms.list_nodes():
+                self._high_symmetry_points = hdf_atoms["high_symmetry_points"]
             return self
 
     def center(self, vacuum=None, axis=(0, 1, 2)):

pyiron/base/generic/parameters.py

@@ -471,9 +471,9 @@ def set_value(self, line, val):
             new_array.append(val)
             new_comments.append("")
             new_params.append("")
-            new_array = np.array(new_array)
-            new_comments = np.array(new_comments)
-            new_params = np.array(new_params)
+            new_array = np.array(new_array).tolist()
+            new_comments = np.array(new_comments).tolist()
+            new_params = np.array(new_params).tolist()
             new_dict = OrderedDict()
             new_dict["Value"] = new_array
             new_dict["Comment"] = new_comments

pyiron/dft/job/generic.py

@@ -152,6 +152,8 @@ def _set_kpoints(
         manual_kpoints=None,
         weights=None,
         reciprocal=True,
+        n_trace=None,
+        trace=None,
     ):
         raise NotImplementedError(
             "The set_kpoints function is not implemented for this code."
@@ -167,6 +169,8 @@ def set_kpoints(
         weights=None,
         reciprocal=True,
         kpoints_per_angstrom=None,
+        n_trace=None,
+        trace=None,
     ):
         """
         Function to setup the k-points
@@ -181,6 +185,8 @@ def set_kpoints(
             reciprocal (bool): Tells if the supplied values are in reciprocal (direct) or cartesian coordinates (in
             reciprocal space)
             kpoints_per_angstrom (float): Number of kpoint per angstrom in each direction
+            n_trace (int): Number of points per trace part for line mode
+            trace (list): ordered list of high symmetry points for line mode
         """
         if kpoints_per_angstrom is not None:
             if mesh is not None:
@@ -200,6 +206,8 @@ def set_kpoints(
             manual_kpoints=manual_kpoints,
             weights=weights,
             reciprocal=reciprocal,
+            n_trace=n_trace,
+            trace=trace,
         )
 
     def calc_static(

pyiron/vasp/base.py

@@ -1015,6 +1015,8 @@ def _set_kpoints(
         manual_kpoints=None,
         weights=None,
         reciprocal=True,
+        n_trace=None,
+        trace=None,
     ):
         """
         Function to setup the k-points for the VASP job
@@ -1029,8 +1031,10 @@ def _set_kpoints(
             reciprocal (bool): Tells if the supplied values are in reciprocal (direct) or cartesian coordinates (in
             reciprocal space)
             kmesh_density (float): Value of the required density
+            n_trace (int): Number of points per trace part for line mode
+            trace (list): ordered list of high symmetry points for line mode
         """
-
+        self.structure.get_high_symmetry_points()
         if not symmetry_reduction:
             self.input.incar["ISYM"] = -1
         scheme_list = ["MP", "GP", "Line", "Manual"]
@@ -1043,7 +1047,19 @@ def _set_kpoints(
         if scheme == "GP":
             self.input.kpoints.set(size_of_mesh=[1, 1, 1], method="Gamma Point")
         if scheme == "Line":
-            raise NotImplementedError("The line mode is not implemented as yet")
+            if n_trace is None:
+                raise ValueError("n_trace has to be defined")
+            high_symmetry_points = self.structure.get_high_symmetry_points()
+            if high_symmetry_points is None:
+                raise ValueError("high_symmetry_points has to be defined")
+            if trace is None:
+                raise ValueError("trace_points has to be defined")
+            self.input.kpoints._set_trace(trace)
+            self.input.kpoints.set(
+                method="Line",
+                n_trace=n_trace,
+                trace_coord=self._get_trace_for_kpoints(trace)
+            )
         if scheme == "Manual":
             if manual_kpoints is None:
                 raise ValueError(
@@ -1068,6 +1084,28 @@ def _set_kpoints(
                         val=" ".join([str(kpt[0]), str(kpt[1]), str(kpt[2]), str(wt)]),
                     )
 
+    def _get_trace_for_kpoints(self, trace):
+        """
+        gets the trace for k-points line mode in a VASP readable form.
+
+        Args:
+            trace (list): ordered list of names for k-points trace
+
+        Returns:
+            list: trace points coordinates for VASP
+        """
+        for t in trace:
+            if not t in self.structure.get_high_symmetry_points().keys():
+                raise ValueError("trace point '{}' is not in high symmetry points".format(t))
+
+        trace_roll = np.roll(trace, -1)
+        k_trace = []
+        for i, t in enumerate(trace):
+            k_trace.append(self.structure.get_high_symmetry_points()[t])
+            k_trace.append(self.structure.get_high_symmetry_points()[trace_roll[i]])
+
+        return k_trace[:-2]
+
     def set_for_band_structure_calc(
         self, num_points, structure=None, read_charge_density=True
     ):
@@ -1648,8 +1686,13 @@ def to_hdf(self, hdf):
             self.kpoints.to_hdf(hdf5_input)
             self.potcar.to_hdf(hdf5_input)
 
-            vasp_dict = {"eddrmm_handling": self._eddrmm}
-            hdf5_input["vasp_dict"] = vasp_dict
+            if "vasp_dict" in hdf5_input.list_nodes():
+                vasp_dict = hdf5_input["vasp_dict"]
+                vasp_dict.update({"eddrmm_handling": self._eddrmm})
+                hdf5_input["vasp_dict"] = vasp_dict
+            else:
+                vasp_dict = {"eddrmm_handling": self._eddrmm}
+                hdf5_input["vasp_dict"] = vasp_dict
 
     def from_hdf(self, hdf):
         """
@@ -1662,11 +1705,12 @@ def from_hdf(self, hdf):
             self.incar.from_hdf(hdf5_input)
             self.kpoints.from_hdf(hdf5_input)
             self.potcar.from_hdf(hdf5_input)
+
+            self._eddrmm = "ignore"
             if "vasp_dict" in hdf5_input.list_nodes():
                 vasp_dict = hdf5_input["vasp_dict"]
-                self._eddrmm = vasp_dict["eddrmm_handling"]
-            else:
-                self._eddrmm = "ignore"
+                if "eddrmm_handling" in vasp_dict.keys():
+                    self._eddrmm = vasp_dict["eddrmm_handling"]
 
 
 class Output:
@@ -2135,15 +2179,52 @@ def __init__(self, input_file_name=None, table_name="kpoints"):
             val_only=True,
             comment_char="!",
         )
+        self._trace = None
+
+    def _set_trace(self, trace):
+        """
+        Sets high symmetry points names of k-points trace (line mode only)
+
+        Args:
+            trace (list): new trace
+        """
+        self._trace = trace
 
-    def set(self, method=None, size_of_mesh=None, shift=None):
+    def _get_trace(self):
+        """
+        Returns high symmetry points name of k-points trace (Line mode only)
+
+        Returns:
+            list: trace values
+        """
+        return self._trace
+
+    def set(self, method=None, size_of_mesh=None, shift=None, n_trace=None, trace_coord=None):
         """
         Sets appropriate tags and values in the KPOINTS file
         Args:
             method (str): Type of meshing scheme (Gamma Point, MP, Manual or Line)
             size_of_mesh (list/numpy.ndarray): List of size 1x3 specifying the required mesh size
             shift (list): List of size 1x3 specifying the user defined shift from the Gamma point
+            n_trace (int): Number of points per trace for line mode
+            trace_coord (list): coordinates of k-points trace in VASP KPOINTS format
         """
+        if n_trace is not None:
+            if self._trace is None or trace_coord is None:
+                raise ValueError("trace have to be defined")
+
+            self.set_value(line=1, val=n_trace)
+            self.set_value(line=3, val="rec")
+
+            trace_names = []
+            for i in range(len(self._trace) - 1):
+                trace_names.append(self._trace[i])
+                trace_names.append(self._trace[i + 1])
+
+            for i, t in enumerate(trace_coord):
+                val = " ".join([str(ii) for ii in t])
+                val = val + " !" + trace_names[i]
+                self.set_value(line=i + 4, val=val)
         if method is not None:
             self.set_value(line=2, val=method)
         if size_of_mesh is not None:
@@ -2178,6 +2259,45 @@ def set_kmesh_by_density(self, structure):
                 )
                 self.set(size_of_mesh=k_mesh)
 
+    def to_hdf(self, hdf, group_name=None):
+        """
+        Store the GenericParameters in an HDF5 file
+
+        Args:
+            hdf (ProjectHDFio): HDF5 group object
+            group_name (str): HDF5 subgroup name - optional
+        """
+        super(Kpoints, self).to_hdf(
+            hdf=hdf,
+            group_name=group_name
+        )
+        if self._trace is not None:
+            if "vasp_dict" in hdf.list_nodes():
+                vasp_dict = hdf["vasp_dict"]
+                vasp_dict.update({"trace": self._trace})
+                hdf["vasp_dict"] = vasp_dict
+            else:
+                vasp_dict = {"trace": self._trace}
+                hdf["vasp_dict"] = vasp_dict
+
+    def from_hdf(self, hdf, group_name=None):
+        """
+        Restore the GenericParameters from an HDF5 file
+
+        Args:
+            hdf (ProjectHDFio): HDF5 group object
+            group_name (str): HDF5 subgroup name - optional
+        """
+        super(Kpoints, self).from_hdf(
+            hdf=hdf,
+            group_name=group_name
+        )
+        self._trace = None
+        if "vasp_dict" in hdf.list_nodes():
+            vasp_dict = hdf["vasp_dict"]
+            if "trace" in vasp_dict.keys():
+                self._trace = vasp_dict["trace"]
+
 
 class Potcar(GenericParameters):
     pot_path_dict = {"GGA": "paw-gga-pbe", "PBE": "paw-gga-pbe", "LDA": "paw-lda"}