VASP does not write charge component for orbital moments

vasp-dev · Jun 17, 2024 · f0c8a10 · f0c8a10
1 parent ac9f063
commit f0c8a10
Show file tree

Hide file tree

Showing 3 changed files with 7 additions and 4 deletions.
diff --git a/src/py4vasp/_calculation/magnetism.py b/src/py4vasp/_calculation/magnetism.py
@@ -230,7 +230,7 @@ def _collinear_moments(self):
     def _noncollinear_moments(self, selection):
         spin_moments = self._raw_data.spin_moments[self._steps, 1:]
         if self._has_orbital_moments:
-            orbital_moments = self._raw_data.orbital_moments[self._steps, 1:]
+            orbital_moments = self._raw_data.orbital_moments[self._steps]
         else:
             orbital_moments = np.zeros_like(spin_moments)
         if selection == "orbital":
@@ -246,7 +246,7 @@ def _add_spin_and_orbital_moments(self):
         if not self._has_orbital_moments:
             return {}
         spin_moments = self._raw_data.spin_moments[self._steps, 1:]
-        orbital_moments = self._raw_data.orbital_moments[self._steps, 1:]
+        orbital_moments = self._raw_data.orbital_moments[self._steps]
         direction_axis = 1 if spin_moments.ndim == 4 else 0
         return {
             "spin_moments": np.moveaxis(spin_moments, direction_axis, -1),

diff --git a/tests/calculation/test_magnetism.py b/tests/calculation/test_magnetism.py
@@ -63,7 +63,7 @@ def __getitem__(self, step):
         reference.moments = np.moveaxis(raw_magnetism.spin_moments[:, 1:4], 1, 3)
     else:
         spin_moments = np.moveaxis(raw_magnetism.spin_moments[:, 1:4], 1, 3)
-        orbital_moments = np.moveaxis(raw_magnetism.orbital_moments[:, 1:4], 1, 3)
+        orbital_moments = np.moveaxis(raw_magnetism.orbital_moments, 1, 3)
         reference.moments = spin_moments + orbital_moments
         reference.spin_moments = spin_moments
         reference.orbital_moments = orbital_moments

diff --git a/tests/conftest.py b/tests/conftest.py
@@ -543,7 +543,10 @@ def _magnetism(selection):
         spin_moments=_make_data(np.arange(np.prod(shape)).reshape(shape)),
     )
     if selection == "orbital_moments":
-        magnetism.orbital_moments = _make_data(np.sqrt(magnetism.spin_moments))
+        remove_charge_component = magnetism.spin_moments[:, 1:]
+        magnetism.orbital_moments = _make_data(np.sqrt(remove_charge_component))
+        print("spin_moments", magnetism.spin_moments.shape)
+        print("orb_moments", magnetism.orbital_moments.shape)
     return magnetism