diff --git a/.pickle b/.pickle
index 431a18e..8b8a276 100644
Binary files a/.pickle and b/.pickle differ
diff --git a/src/grogupy/grogu.py b/src/grogupy/grogu.py
index 225f3f3..7c0aa85 100644
--- a/src/grogupy/grogu.py
+++ b/src/grogupy/grogu.py
@@ -454,8 +454,10 @@ def main(simulation_parameters, magnetic_entities, pairs):
                 site_j = magnetic_entities[pair["aj"]]
                 storage = []
                 # iterate over the first order local perturbations in all possible orientations for the two sites
-                for Vui in site_i["Vu1"][i]:
-                    for Vuj in site_j["Vu1"][i]:
+                # actually all possible orientations without the orientation for the off-diagonal anisotropy
+                # that is why we only take the first two of each Vu1
+                for Vui in site_i["Vu1"][i][:2]:
+                    for Vuj in site_j["Vu1"][i][:2]:
                         # The Szunyogh-Lichtenstein formula
                         traced = np.trace(
                             (Vui @ Gij @ Vuj @ Gji), axis1=1, axis2=2
diff --git a/tests/test_io.py b/tests/test_io.py
index 5caa549..be895fe 100644
--- a/tests/test_io.py
+++ b/tests/test_io.py
@@ -27,8 +27,8 @@ from pathlib import Path
 import numpy as np
 import pytest
 
+from grogupy.globals import DEFAULT_ARGUMENTS
 from grogupy.io import (
-    default_args,
     load_pickle,
     print_atoms_and_pairs,
     print_job_description,
@@ -42,7 +42,7 @@ from grogupy.io import (
 @pytest.fixture
 def simulation_parameters():
     """Create sample simulation parameters for testing"""
-    params = default_args.copy()
+    params = DEFAULT_ARGUMENTS.copy()
     params.update(
         {
             "infile": "test.fdf",
@@ -196,23 +196,23 @@ def test_print_job_description(simulation_parameters):
 
 
 # Test default arguments
-def test_default_args_structure():
+def test_DEFAULT_ARGUMENTS_structure():
     """Test the structure and values of default arguments"""
-    assert isinstance(default_args, dict)
-    assert "infile" in default_args
-    assert "outfile" in default_args
-    assert "kset" in default_args
-    assert "kdirs" in default_args
-    assert "eset" in default_args
-    assert "esetp" in default_args
+    assert isinstance(DEFAULT_ARGUMENTS, dict)
+    assert "infile" in DEFAULT_ARGUMENTS
+    assert "outfile" in DEFAULT_ARGUMENTS
+    assert "kset" in DEFAULT_ARGUMENTS
+    assert "kdirs" in DEFAULT_ARGUMENTS
+    assert "eset" in DEFAULT_ARGUMENTS
+    assert "esetp" in DEFAULT_ARGUMENTS
 
     # Test specific default values
-    assert default_args["kset"] == 2
-    assert default_args["kdirs"] == "xyz"
-    assert default_args["eset"] == 42
-    assert default_args["esetp"] == 1000
-    assert default_args["parallel_solver_for_Gk"] is False
-    assert default_args["padawan_mode"] is True
+    assert DEFAULT_ARGUMENTS["kset"] == 2
+    assert DEFAULT_ARGUMENTS["kdirs"] == "xyz"
+    assert DEFAULT_ARGUMENTS["eset"] == 42
+    assert DEFAULT_ARGUMENTS["esetp"] == 1000
+    assert DEFAULT_ARGUMENTS["parallel_solver_for_Gk"] is False
+    assert DEFAULT_ARGUMENTS["padawan_mode"] is True
 
 
 def test_simulation_parameters_validation(simulation_parameters):
diff --git a/tests/test_utilities.py b/tests/test_utilities.py
index 5ddd221..44f01db 100644
--- a/tests/test_utilities.py
+++ b/tests/test_utilities.py
@@ -24,6 +24,7 @@ from hypothesis import given
 from hypothesis import strategies as st
 from numpy.testing import assert_allclose, assert_array_almost_equal
 
+from grogupy.globals import TAU_0, TAU_X, TAU_Y, TAU_Z
 from grogupy.utilities import (
     RotM,
     RotMa2b,
@@ -33,11 +34,7 @@ from grogupy.utilities import (
     make_contour,
     make_kset,
     read_siesta_emin,
-    tau_0,
     tau_u,
-    tau_x,
-    tau_y,
-    tau_z,
 )
 
 
@@ -45,19 +42,19 @@ from grogupy.utilities import (
 def test_pauli_matrices_properties():
     """Test fundamental properties of Pauli matrices"""
     # Test anticommutation relations
-    assert_array_almost_equal(tau_x @ tau_y + tau_y @ tau_x, np.zeros((2, 2)))
-    assert_array_almost_equal(tau_y @ tau_z + tau_z @ tau_y, np.zeros((2, 2)))
-    assert_array_almost_equal(tau_z @ tau_x + tau_x @ tau_z, np.zeros((2, 2)))
+    assert_array_almost_equal(TAU_X @ TAU_Y + TAU_Y @ TAU_X, np.zeros((2, 2)))
+    assert_array_almost_equal(TAU_Y @ TAU_Z + TAU_Z @ TAU_Y, np.zeros((2, 2)))
+    assert_array_almost_equal(TAU_Z @ TAU_X + TAU_X @ TAU_Z, np.zeros((2, 2)))
 
     # Test square of Pauli matrices equals identity
-    assert_array_almost_equal(tau_x @ tau_x, tau_0)
-    assert_array_almost_equal(tau_y @ tau_y, tau_0)
-    assert_array_almost_equal(tau_z @ tau_z, tau_0)
+    assert_array_almost_equal(TAU_X @ TAU_X, TAU_0)
+    assert_array_almost_equal(TAU_Y @ TAU_Y, TAU_0)
+    assert_array_almost_equal(TAU_Z @ TAU_Z, TAU_0)
 
     # Test Hermiticity
-    assert_array_almost_equal(tau_x, tau_x.conj().T)
-    assert_array_almost_equal(tau_y, tau_y.conj().T)
-    assert_array_almost_equal(tau_z, tau_z.conj().T)
+    assert_array_almost_equal(TAU_X, TAU_X.conj().T)
+    assert_array_almost_equal(TAU_Y, TAU_Y.conj().T)
+    assert_array_almost_equal(TAU_Z, TAU_Z.conj().T)
 
 
 # Test hsk function