Fix SparsePauliOp initialization with dense Y labels

qiskit/quantum_info/operators/symplectic/sparse_pauli_op.py

@@ -167,7 +167,16 @@ def __init__(
             # move the phase of `pauli_list` to `self._coeffs`
             phase = pauli_list._phase
             count_y = pauli_list._count_y()
-            self._coeffs = np.asarray((-1j) ** (phase - count_y) * coeffs, dtype=coeffs.dtype)
+
+            # Compute exponentiation via integer arithmetic and lookup table to avoid floating point errors
+            exp_val = (phase - count_y) % 4
+            p = (3 * exp_val) % 4  # (−i)^(phase - count_y) = i^(3*(phase - count_y))
+            lookup = np.array([1 + 0j, 1j, -1 + 0j, -1j], dtype=coeffs.dtype)
+
+            vals = lookup[p]
+            self._coeffs = vals * coeffs
+
+            # Update pauli_list phase
             pauli_list._phase = np.mod(count_y, 4)
             self._pauli_list = pauli_list
 

test/python/quantum_info/operators/symplectic/test_sparse_pauli_op.py

@@ -140,6 +140,17 @@ def test_sparse_pauli_op_init(self):
             coeffs[:] = 0
             self.assertEqual(spp_op, ref_op)
 
+    def test_sparse_pauli_op_init_long_ys(self):
+        """Test SparsePauliOp initialization."""
+        with self.subTest(msg="make SparsePauliOp from SparsePauliOp"):
+            # produces non-zero imaginary component
+            yyy = SparsePauliOp("Y" * 100)
+            print(yyy.coeffs)
+            # correctly contains no imaginary component
+            yy = SparsePauliOp("Y" * 99)
+            print(yy.coeffs)
+            self.assertNotEqual(yyy, yy)
+
 
 @ddt.ddt
 class TestSparsePauliOpConversions(QiskitTestCase):