diff --git a/cirq-core/cirq/ops/dense_pauli_string.py b/cirq-core/cirq/ops/dense_pauli_string.py index 0be1c66e350..dfeab3098d3 100644 --- a/cirq-core/cirq/ops/dense_pauli_string.py +++ b/cirq-core/cirq/ops/dense_pauli_string.py @@ -319,7 +319,7 @@ def sparse(self, qubits: Sequence[cirq.Qid] | None = None) -> cirq.PauliString: return pauli_string.PauliString( coefficient=self.coefficient, - qubit_pauli_map={q: PAULI_GATES[p] for q, p in zip(qubits, self.pauli_mask) if p}, + qubit_pauli_map={q: PAULI_GATES[p] for q, p in zip(qubits, self.pauli_mask)}, ) def __str__(self) -> str: diff --git a/cirq-core/cirq/ops/linear_combinations.py b/cirq-core/cirq/ops/linear_combinations.py index e7afeecd16d..5fcfb4da867 100644 --- a/cirq-core/cirq/ops/linear_combinations.py +++ b/cirq-core/cirq/ops/linear_combinations.py @@ -35,7 +35,7 @@ if TYPE_CHECKING: import cirq -UnitPauliStringT = frozenset[tuple[raw_types.Qid, pauli_gates.Pauli]] +UnitPauliStringT = frozenset[tuple[raw_types.Qid, pauli_gates.Pauli | identity.IdentityGate]] PauliSumLike = Union[ complex, PauliString, 'PauliSum', pauli_string.SingleQubitPauliStringGateOperation ] @@ -336,7 +336,7 @@ def _is_linear_dict_of_unit_pauli_string(linear_dict: value.LinearDict[UnitPauli for qid, pauli in k: if not isinstance(qid, raw_types.Qid): return False - if not isinstance(pauli, pauli_gates.Pauli): + if not isinstance(pauli, pauli_gates.Pauli) and pauli != identity.I: return False return True @@ -468,7 +468,10 @@ def from_pauli_strings(cls, terms: PauliString | list[PauliString]) -> PauliSum: terms = [terms] termdict: defaultdict[UnitPauliStringT, value.Scalar] = defaultdict(lambda: 0) for pstring in terms: - key = frozenset(pstring._qubit_pauli_map.items()) + key = frozenset( + list(pstring._qubit_pauli_map.items()) + + [(qubit, identity.I) for qubit in pstring.identity_qubits] + ) termdict[key] += pstring.coefficient return cls(linear_dict=value.LinearDict(termdict)) diff --git a/cirq-core/cirq/ops/pauli_string.py b/cirq-core/cirq/ops/pauli_string.py index 9bdcbcf931a..ac95da34a66 100644 --- a/cirq-core/cirq/ops/pauli_string.py +++ b/cirq-core/cirq/ops/pauli_string.py @@ -154,7 +154,7 @@ class PauliString(raw_types.Operation, Generic[TKey]): def __init__( self, *contents: cirq.PAULI_STRING_LIKE, - qubit_pauli_map: dict[TKey, cirq.Pauli] | None = None, + qubit_pauli_map: Mapping[TKey, cirq.Pauli | cirq.IdentityGate] | None = None, coefficient: cirq.TParamValComplex = 1, ): """Initializes a new `PauliString` operation. @@ -166,9 +166,7 @@ def __init__( values is given, they are each individually converted and then multiplied from left to right in order. qubit_pauli_map: Initial dictionary mapping qubits to pauli - operations. Defaults to the empty dictionary. Note that, unlike - dictionaries passed to contents, this dictionary must not - contain any identity gate values. Further note that this + operations. Defaults to the empty dictionary. Note that this argument specifies values that are logically *before* factors specified in `contents`; `contents` are *right* multiplied onto the values in this dictionary. @@ -177,12 +175,22 @@ def __init__( Raises: TypeError: If the `qubit_pauli_map` has values that are not Paulis. """ + self._qubit_pauli_map: dict[TKey, cirq.Pauli] = {} + # Maintain a list of identity qubits. + # Historically this class retains only qubits with Pauli operations applied + # to them. However, downstream calculations (e.g. applying exponentiation, etc.) + # relies on knowledge of where identity qubits are placed as well. To maintain + # existing behavior, track identity qubits under a separate field. + self._identity_qubits: Sequence[cirq.Qid] = [] if _compat.__cirq_debug__.get() and qubit_pauli_map is not None: - for v in qubit_pauli_map.values(): - if not isinstance(v, pauli_gates.Pauli): - raise TypeError(f'{v} is not a Pauli') + for k, v in qubit_pauli_map.items(): + if isinstance(v, pauli_gates.Pauli): + self._qubit_pauli_map[k] = v + elif v == identity.I: + self._identity_qubits.append(k) + else: + raise TypeError(f'{v} is not a Pauli or identity') - self._qubit_pauli_map: dict[TKey, cirq.Pauli] = qubit_pauli_map or {} self._coefficient: cirq.TParamValComplex | sympy.Expr = ( coefficient if isinstance(coefficient, sympy.Expr) else complex(coefficient) ) @@ -329,6 +337,11 @@ def qubits(self) -> tuple[TKey, ...]: """Returns a tuple of qubits on which this pauli string acts.""" return tuple(self.keys()) + @property + def identity_qubits(self) -> tuple[cirq.Qid, ...]: + """Returns a tuple of qubits that this pauli string is applying I on.""" + return tuple(self._identity_qubits) + def _circuit_diagram_info_(self, args: cirq.CircuitDiagramInfoArgs) -> list[str]: if not len(self._qubit_pauli_map): return NotImplemented diff --git a/cirq-core/cirq/ops/pauli_sum_exponential.py b/cirq-core/cirq/ops/pauli_sum_exponential.py index ae0a1885c70..dcde5e57b36 100644 --- a/cirq-core/cirq/ops/pauli_sum_exponential.py +++ b/cirq-core/cirq/ops/pauli_sum_exponential.py @@ -108,9 +108,35 @@ def matrix(self) -> np.ndarray: """ if protocols.is_parameterized(self._exponent): raise ValueError("Exponent should not parameterized.") - ret = np.ones(1) - for pauli_string_exp in self: - ret = np.kron(ret, protocols.unitary(pauli_string_exp)) + + # Retrieve all qubits that are involved in producing the PauliSum used by + # this class, including identity qubits. These qubits will be used to + # compute the appropriate total size of the final unitary. + qubits = set() + for pauli_string in self._pauli_sum: + qubits.update(pauli_string.qubits) + qubits.update(pauli_string.identity_qubits) + + # Calculate the identity e^{i \theta \sum P} = \prod e^{i \theta P_n} + ret = np.eye(2 ** len(qubits)) + for pauli_string in self._pauli_sum: + # For each e^{i \theta P_n}, compute its value using the identity + # e^{i \theta P_n} = I \cos{(\theta)} + i P_n \sin{(\theta)} + # in the context of the qubits involved in the PauliSum. + pauli_string_qubits = set(pauli_string.qubits) + pauli_string_unitary = np.ones(1, dtype=float) + for qubit in sorted(qubits): + if qubit in pauli_string_qubits: + pauli_string_unitary = np.kron( + pauli_string_unitary, protocols.unitary(pauli_string[qubit]) + ) + else: + pauli_string_unitary = np.kron(pauli_string_unitary, np.eye(2)) + theta = pauli_string.coefficient * self._multiplier + theta *= self._exponent + cos_term = np.cos(theta) * np.eye(2 ** len(qubits)) + sin_term = np.sin(theta) * 1j * pauli_string_unitary + ret = ret @ (cos_term + sin_term) return ret @_compat.cached_method diff --git a/cirq-core/cirq/ops/pauli_sum_exponential_test.py b/cirq-core/cirq/ops/pauli_sum_exponential_test.py index f685fb6cc9a..34bf344bccd 100644 --- a/cirq-core/cirq/ops/pauli_sum_exponential_test.py +++ b/cirq-core/cirq/ops/pauli_sum_exponential_test.py @@ -102,6 +102,14 @@ def test_pauli_sum_exponential_parameterized_matrix_raises() -> None: cirq.PauliSumExponential(2j * cirq.X(q0) + 3j * cirq.Z(q1), np.pi / 2), np.array([[1j, 0, 0, 0], [0, -1j, 0, 0], [0, 0, 1j, 0], [0, 0, 0, -1j]]), ), + ( + cirq.PauliSumExponential(cirq.DensePauliString('XI')(q0, q1), exponent=np.pi / 4), + np.array([[1, 0, 1j, 0], [0, 1, 0, 1j], [1j, 0, 1, 0], [0, 1j, 0, 1]]) / (2**0.5), + ), + ( + cirq.PauliSumExponential(cirq.DensePauliString('IY')(q0, q1), exponent=3 * np.pi / 4), + np.array([[-1, 1, 0, 0], [-1, -1, 0, 0], [0, 0, -1, 1], [0, 0, -1, -1]]) / (2**0.5), + ), ), ) def test_pauli_sum_exponential_has_correct_unitary(psum_exp, expected_unitary) -> None: