fixed hamiltonian simulation of projector operators

src/qrisp/operators/pauli/pauli_hamiltonian.py

@@ -19,14 +19,10 @@
 from qrisp.operators.pauli.pauli_term import PauliTerm
 from qrisp.operators.pauli.pauli_measurement import PauliMeasurement
 from qrisp.operators.pauli.measurement import get_measurement
-from qrisp import h, sx, IterationEnvironment, conjugate, merge
+from qrisp import h, s, x, IterationEnvironment, conjugate, merge
 
 import sympy as sp
 
-from sympy import init_printing
-# Initialize automatic LaTeX rendering
-init_printing()
-
 threshold = 1e-9
 
 #
@@ -405,8 +401,6 @@ def to_sparse_matrix(self, factor_amount = None):
         import scipy.sparse as sp
         from scipy.sparse import kron as TP, csr_matrix
 
-        I = csr_matrix([[1,0],[0,1]])
-
         def get_matrix(P):
             if P=="I":
                 return csr_matrix([[1,0],[0,1]])
@@ -440,9 +434,6 @@ def recursive_TP(keys,pauli_dict):
             coeffs.append(coeff)
             participating_indices = participating_indices.union(pauli.non_trivial_indices())
 
-
-
-
         if factor_amount is None:
             if len(participating_indices):
                 factor_amount = max(participating_indices) + 1
@@ -465,7 +456,7 @@ def recursive_TP(keys,pauli_dict):
         for k in range(m):
             M += complex(coeffs[k])*recursive_TP(keys.copy(),pauli_dicts[k])
 
-        res = ((M + M.transpose())/2)
+        res = ((M + M.transpose().conjugate())/2)
         res.sum_duplicates()
         return res
 
@@ -716,7 +707,10 @@ def change_of_basis(qarg, pauli_dict):
                 if axis=="X":
                     h(qarg[index])
                 if axis=="Y":
-                    sx(qarg[index])
+                    s(qarg[index])
+                    h(qarg[index])
+                    x(qarg[index])
+
 
         groups, bases = self.commuting_qw_groups(show_bases=True)
 

src/qrisp/operators/pauli/pauli_term.py

@@ -193,21 +193,22 @@ def inv_ghz_state(qb_list):
 
             projector_ctrl_state = ""
             projector_qubits = []
-            for i in range(len(projector_indices)):
-                if projector_state[i]:
-                    projector_ctrl_state += "1"
-                else:
-                    projector_ctrl_state += "0"
 
+
+            for i in range(len(projector_indices)):
+                projector_ctrl_state += str(int(projector_state[i]))
                 projector_qubits.append(qv[projector_indices[i]])
-
-            if len(ladder_indices) == 0:
+
+
+            control_qubit_available = False
+
+            if len(ladder_indices + projector_qubits) == 0:
                 env = QuantumEnvironment()
-                coeff *= 2
-            elif len(ladder_indices) == 1:
+            elif len(ladder_indices) == 1 and len(projector_indices) == 0:
                 env = QuantumEnvironment()
-            elif len(ladder_indices) == 2:
+            elif len(ladder_indices) == 2 and len(projector_indices) == 0:
                 hs_ancilla = qv[ladder_indices[0]]
+                control_qubit_available = True
                 if ladder_ctrl_state[0] == "0":
                     env = conjugate(x)(hs_ancilla)
                 else:
@@ -216,6 +217,7 @@ def inv_ghz_state(qb_list):
                 # We furthermore allocate an ancillae to perform an efficient
                 # multi controlled rz.
                 hs_ancilla = QuantumBool()
+                control_qubit_available = True
 
                 env = conjugate(mcx)(ladder_qubits[:-1] + projector_qubits, 
                                     hs_ancilla, 
@@ -246,15 +248,17 @@ def flip_anchor_qubit(qv, anchor_index, Z_indices):
                 # Perform the conjugation
                 with conjugate(flip_anchor_qubit)(qv, anchor_index, Z_indices):
 
-
-                    if len(ladder_indices) > 1:
+                    if control_qubit_available:
                         env = control(hs_ancilla)
                     else:
                         env = QuantumEnvironment()
-
+
+                    if len(ladder_indices) == 0:
+                        coeff *= 2
                     # Perform the controlled RZ
                     with env:
                             rz(-coeff, qv[anchor_index])
+
 
         if len(ladder_indices) > 2:
             # Delete ancilla

tests/hamiltonian_tests/test_pauli_hamiltonian_simulation.py

@@ -17,7 +17,7 @@
 """
 
 from qrisp import QuantumVariable, x, QPE
-from qrisp.operators.pauli import X, Z, A, C, P0, P1
+from qrisp.operators.pauli import X, Y, Z, A, C, P0, P1
 import numpy as np
 
 def test_pauli_hamiltonian_simulation():
@@ -103,4 +103,31 @@ def verify_trotterization(H):
     verify_trotterization(H)
 
     H = 0.1809312*X(0)*X(1)*A(2)
+    verify_trotterization(H)
+
+    H = 0.1809312*X(0)*X(1)*A(2)*Y(3)
+    verify_trotterization(H)
+
+    H = Z(0)*Z(1)*P0(2)
+    verify_trotterization(H)
+
+    H = Z(0)*Z(1)*A(2)*P1(1)
+    verify_trotterization(H)
+
+    H = Z(0)*Z(1)*A(2)*C(3)*P0(4)*P1(3)
+    verify_trotterization(H)
+
+    H = A(2)*P0(1)
+    verify_trotterization(H)
+
+    H = A(2)*C(1)*P0(0)*P1(2)
+    verify_trotterization(H)
+
+    H = 0.1809312*X(0)*X(1)*P1(3)
+    verify_trotterization(H)
+
+    H = 0.1809312*X(0)*X(1)*A(2)*P1(0)
+    verify_trotterization(H)
+
+    H = 0.1809312*X(0)*X(1)*A(2)*Y(3)*P1(5)
     verify_trotterization(H)