Fix qubit order for 2-qubit LieTrotter evolution synthesis (backport #7551)

qiskit/circuit/library/pauli_evolution.py

@@ -13,6 +13,7 @@
 """A gate to implement time-evolution of operators."""
 
 from typing import Union, Optional
+import numpy as np
 
 from qiskit.circuit.gate import Gate
 from qiskit.circuit.parameterexpression import ParameterExpression
@@ -127,14 +128,17 @@ def _to_sparse_pauli_op(operator):
     if isinstance(operator, PauliSumOp):
         sparse_pauli = operator.primitive
         sparse_pauli._coeffs *= operator.coeff
-        return sparse_pauli
-    if isinstance(operator, PauliOp):
+    elif isinstance(operator, PauliOp):
         sparse_pauli = SparsePauliOp(operator.primitive)
         sparse_pauli._coeffs *= operator.coeff
-        return sparse_pauli
-    if isinstance(operator, Pauli):
-        return SparsePauliOp(operator)
-    if isinstance(operator, SparsePauliOp):
-        return operator
+    elif isinstance(operator, Pauli):
+        sparse_pauli = SparsePauliOp(operator)
+    elif isinstance(operator, SparsePauliOp):
+        sparse_pauli = operator
+    else:
+        raise ValueError(f"Unsupported operator type for evolution: {type(operator)}.")
 
-    raise ValueError(f"Unsupported operator type for evolution: {type(operator)}.")
+    if any(np.iscomplex(sparse_pauli.coeffs)):
+        raise ValueError("Operator contains complex coefficients, which are not supported.")
+
+    return sparse_pauli

qiskit/synthesis/evolution/product_formula.py

@@ -138,6 +138,8 @@ def evolve_pauli(
 
 def _single_qubit_evolution(pauli, time):
     definition = QuantumCircuit(pauli.num_qubits)
+    # Note that all phases are removed from the pauli label and are only in the coefficients.
+    # That's because the operators we evolved have all been translated to a SparsePauliOp.
     for i, pauli_i in enumerate(reversed(pauli.to_label())):
         if pauli_i == "X":
             definition.rx(2 * time, i)
@@ -150,8 +152,10 @@ def _single_qubit_evolution(pauli, time):
 
 
 def _two_qubit_evolution(pauli, time, cx_structure):
-    # get the Paulis and the qubits they act on
-    labels_as_array = np.array(list(pauli.to_label()))
+    # Get the Paulis and the qubits they act on.
+    # Note that all phases are removed from the pauli label and are only in the coefficients.
+    # That's because the operators we evolved have all been translated to a SparsePauliOp.
+    labels_as_array = np.array(list(reversed(pauli.to_label())))
     qubits = np.where(labels_as_array != "I")[0]
     labels = np.array([labels_as_array[idx] for idx in qubits])
 
@@ -165,9 +169,9 @@ def _two_qubit_evolution(pauli, time, cx_structure):
     elif all(labels == "Z"):  # RZZ
         definition.rzz(2 * time, qubits[0], qubits[1])
     elif labels[0] == "Z" and labels[1] == "X":  # RZX
-        definition.rzx(2 * time, qubits[1], qubits[0])
-    elif labels[0] == "X" and labels[1] == "Z":  # RXZ
         definition.rzx(2 * time, qubits[0], qubits[1])
+    elif labels[0] == "X" and labels[1] == "Z":  # RXZ
+        definition.rzx(2 * time, qubits[1], qubits[0])
     else:  # all the others are not native in Qiskit, so use default the decomposition
         definition = _multi_qubit_evolution(pauli, time, cx_structure)
 
@@ -186,6 +190,8 @@ def _multi_qubit_evolution(pauli, time, cx_structure):
 
     # determine qubit to do the rotation on
     target = None
+    # Note that all phases are removed from the pauli label and are only in the coefficients.
+    # That's because the operators we evolved have all been translated to a SparsePauliOp.
     for i, pauli_i in enumerate(reversed(pauli.to_label())):
         if pauli_i != "I":
             target = i

releasenotes/notes/fix-lietrotter-2q-61d5cd66e0bf7359.yaml

@@ -0,0 +1,36 @@
+---
+fixes:
+  - |
+    Fix the qubit order of 2-qubit evolutions in the
+    :class:`qiskit.circuit.library.PauliEvolutionGate`, if used with a product formula synthesis.
+    For instance, before, the evolution of ``IIZ + IZI + IZZ``
+
+    .. code-block:: python
+
+        from qiskit.circuit.library import PauliEvolutionGate
+        from qiskit.opflow import I, Z
+        operator = (I ^ I ^ Z) + (I ^ Z ^ I) + (I ^ Z ^ Z)
+        print(PauliEvolutionGate(operator).definition.decompose())
+
+    produced
+
+    .. code-block::
+
+             ┌───────┐
+        q_0: ┤ Rz(2) ├────────
+             ├───────┤
+        q_1: ┤ Rz(2) ├─■──────
+             └───────┘ │ZZ(2)
+        q_2: ──────────■──────
+
+
+    whereas now it correctly yields
+
+    .. code-block::
+
+             ┌───────┐
+        q_0: ┤ Rz(2) ├─■──────
+             ├───────┤ │ZZ(2)
+        q_1: ┤ Rz(2) ├─■──────
+             └───────┘
+        q_2: ─────────────────

test/python/circuit/library/test_evolution_gate.py

@@ -57,19 +57,21 @@ def test_lie_trotter(self):
         self.assertEqual(decomposed.count_ops()["cx"], reps * 3 * 4)
 
     def test_rzx_order(self):
-        """Test ZX is mapped onto the correct qubits."""
-        evo_gate = PauliEvolutionGate(X ^ Z)
-        decomposed = evo_gate.definition.decompose()
+        """Test ZX and XZ is mapped onto the correct qubits."""
+        for op, indices in zip([X ^ Z, Z ^ X], [(0, 1), (1, 0)]):
+            with self.subTest(op=op, indices=indices):
+                evo_gate = PauliEvolutionGate(op)
+                decomposed = evo_gate.definition.decompose()
 
-        ref = QuantumCircuit(2)
-        ref.h(1)
-        ref.cx(1, 0)
-        ref.rz(2.0, 0)
-        ref.cx(1, 0)
-        ref.h(1)
+                ref = QuantumCircuit(2)
+                ref.h(indices[1])
+                ref.cx(indices[1], indices[0])
+                ref.rz(2.0, indices[0])
+                ref.cx(indices[1], indices[0])
+                ref.h(indices[1])
 
-        # don't use circuit equality since RZX here decomposes with RZ on the bottom
-        self.assertTrue(Operator(decomposed).equiv(ref))
+                # don't use circuit equality since RZX here decomposes with RZ on the bottom
+                self.assertTrue(Operator(decomposed).equiv(ref))
 
     def test_suzuki_trotter(self):
         """Test constructing the circuit with Lie Trotter decomposition."""
@@ -262,6 +264,23 @@ def test_paulisumop_coefficients_respected(self):
         ]
         self.assertListEqual(rz_angles, [20, 30, -10])
 
+    def test_lie_trotter_two_qubit_correct_order(self):
+        """Test that evolutions on two qubit operators are in the right order.
+
+        Regression test of Qiskit/qiskit-terra#7544.
+        """
+        operator = I ^ Z ^ Z
+        time = 0.5
+        exact = scipy.linalg.expm(-1j * time * operator.to_matrix())
+        lie_trotter = PauliEvolutionGate(operator, time, synthesis=LieTrotter())
+
+        self.assertTrue(Operator(lie_trotter).equiv(exact))
+
+    def test_complex_op_raises(self):
+        """Test an operator with complex coefficient raises an error."""
+        with self.assertRaises(ValueError):
+            _ = PauliEvolutionGate(Pauli("iZ"))
+
     @data(LieTrotter, MatrixExponential)
     def test_inverse(self, synth_cls):
         """Test calculating the inverse is correct."""