Implement AlternativeSpinHamiltonian for Issue #179

check_my_spin.py

@@ -0,0 +1,26 @@
+import numpy as np
+from moha.hamiltonians import AlternativeSpinHamiltonian
+
+# Initialize 2-spin system with homogeneous coupling
+n_spins = 2
+couplings = [1.0] 
+model = AlternativeSpinHamiltonian(n_spins, couplings)
+
+print(f"Testing AlternativeSpinHamiltonian (Spins: {n_spins})")
+H_matrix = model.generate_integrals()
+
+print("\nHamiltonian Matrix:")
+print(H_matrix)
+
+# Verify eigenvalues against analytical solutions for 2-spin Heisenberg model
+energies = np.linalg.eigvals(H_matrix)
+sorted_energies = np.sort(energies.real)
+
+print("\nCalculated Eigenvalues:")
+print(sorted_energies)
+
+expected = np.array([-0.75, 0.25, 0.25, 0.25])
+if np.allclose(sorted_energies, expected):
+    print("\nStatus: Passed ")
+else:
+    print("\nStatus: Failed - Eigenvalues deviate from analytical expectations.")

moha/hamiltonians.py

@@ -684,3 +684,70 @@ def __init__(self,
             J_ax=J_ax,
             connectivity=connectivity
         )
+
+
+class AlternativeSpinHamiltonian:
+    """
+    Implementation of the alternative approach to Spin Hamiltonians.
+    Reference: Issue #179
+    """
+    def __init__(self, num_spins, coupling_constants, connectivity=None):
+        self.num_spins = int(num_spins)
+        # Convert to numpy array immediately for safety
+        self.J = np.array(coupling_constants)
+
+        # Default to a 1D chain if no connectivity is provided.
+        if connectivity is None or len(connectivity) == 0:
+            self.connectivity = [(i, i + 1) for i in range(self.num_spins - 1)]
+        else:
+            self.connectivity = connectivity
+
+    def get_operator(self, site_index, op_type='z'):
+        """
+        Generates Sx, Sy, or Sz for a specific site using the Kronecker product.
+        op_type: 'x', 'y', or 'z'
+        """
+        if op_type == 'x':
+            base_op = 0.5 * np.array([[0, 1], [1, 0]])
+        elif op_type == 'y':
+            base_op = 0.5 * np.array([[0, -1j], [1j, 0]])
+        else:
+            base_op = 0.5 * np.array([[1, 0], [0, -1]])
+
+        identity = np.eye(2)
+        op = base_op if site_index == 0 else identity
+
+        # Chain the rest of the sites
+        for i in range(1, self.num_spins):
+            next_op = base_op if i == site_index else identity
+            op = np.kron(op, next_op)
+
+        return op
+
+    def generate_integrals(self):
+        """
+        Builds the full Heisenberg Hamiltonian matrix.
+        H = sum over <i,j> of J * (Sx_i*Sx_j + Sy_i*Sy_j + Sz_i*Sz_j)
+        """
+        dim = 2 ** self.num_spins
+        h_matrix = np.zeros((dim, dim), dtype=np.complex128)
+
+        # Safely extract the J value
+        try:
+            J_val = float(self.J[0] if self.J.size > 0 else 1.0)
+        except (TypeError, IndexError):
+            J_val = 1.0
+
+
+        for (i, j) in self.connectivity:
+            # Calculate Sx_i*Sx_j + Sy_i*Sy_j + Sz_i*Sz_j
+            for term in ['x', 'y', 'z']:
+                op_i = self.get_operator(i, term)
+                op_j = self.get_operator(j, term)
+
+
+                h_matrix += J_val * (op_i @ op_j)
+
+        # The Heisenberg Hamiltonian is Hermitian, so we return the real part.
+        return h_matrix.real
+

moha/test/test_alternative_spin.py

@@ -0,0 +1,23 @@
+import numpy as np
+import pytest
+from moha.hamiltonians import AlternativeSpinHamiltonian
+
+def test_heisenberg_2spin_energies():
+    """Verify the 2-spin Heisenberg model produces exact analytical eigenvalues."""
+    n_spins = 2
+    couplings = [1.0]
+    model = AlternativeSpinHamiltonian(n_spins, couplings)
+
+    H = model.generate_integrals()
+    energies = np.linalg.eigvals(H)
+    sorted_energies = np.sort(energies.real)
+
+    # Expected eigenvalues for J=1: [-0.75, 0.25, 0.25, 0.25]
+    expected = np.array([-0.75, 0.25, 0.25, 0.25])
+    assert np.allclose(sorted_energies, expected), f"Expected {expected}, got {sorted_energies}"
+
+def test_hermiticity():
+    """The Hamiltonian must be Hermitian (H = H.T)."""
+    model = AlternativeSpinHamiltonian(3, [1.0, 1.0])
+    H = model.generate_integrals()
+    assert np.allclose(H, H.T), "Hamiltonian is not Hermitian!"