Coverage for /home/runner/work/tket/tket/pytket/pytket/utils/operators.py: 91%
135 statements
« prev ^ index » next coverage.py v7.6.12, created at 2025-03-14 11:30 +0000
« prev ^ index » next coverage.py v7.6.12, created at 2025-03-14 11:30 +0000
1# Copyright Quantinuum
2#
3# Licensed under the Apache License, Version 2.0 (the "License");
4# you may not use this file except in compliance with the License.
5# You may obtain a copy of the License at
6#
7# http://www.apache.org/licenses/LICENSE-2.0
8#
9# Unless required by applicable law or agreed to in writing, software
10# distributed under the License is distributed on an "AS IS" BASIS,
11# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12# See the License for the specific language governing permissions and
13# limitations under the License.
15import copy
16from typing import TYPE_CHECKING, Any, Union
18import numpy
19import numpy as np
20from sympy import Expr, Symbol, im, re, sympify
22from pytket.circuit import Qubit
23from pytket.pauli import QubitPauliString, pauli_string_mult
24from pytket.utils.serialization import complex_to_list, list_to_complex
26CoeffTypeAccepted = Union[int, float, complex, Expr]
28if TYPE_CHECKING: 28 ↛ 29line 28 didn't jump to line 29 because the condition on line 28 was never true
29 from scipy.sparse import csc_matrix
32def _coeff_convert(coeff: CoeffTypeAccepted | str) -> Expr:
33 sympy_val = sympify(coeff)
34 if not isinstance(sympy_val, Expr): 34 ↛ 35line 34 didn't jump to line 35 because the condition on line 34 was never true
35 raise ValueError("Unsupported value for QubitPauliString coefficient")
36 return sympy_val
39class QubitPauliOperator:
40 """
41 Generic data structure for generation of circuits and expectation
42 value calculation. Contains a dictionary from QubitPauliString to
43 sympy Expr. Capacity for symbolic expressions allows the operator
44 to be used to generate ansätze for variational algorithms.
46 Represents a mathematical object :math:`\\sum_j \\alpha_j P_j`,
47 where each :math:`\\alpha_j` is a complex symbolic expression and
48 :math:`P_j` is a Pauli string, i.e. :math:`P_j \\in \\{ I, X, Y,
49 Z\\}^{\\otimes n}`.
51 A prototypical example is a molecular Hamiltonian, for which one
52 may wish to calculate the expectation value :math:`\\langle \\Psi
53 | H | \\Psi \\rangle` by decomposing :math:`H` into individual
54 Pauli measurements. Alternatively, one may wish to evolve a state
55 by the operator :math:`e^{-iHt}` for digital quantum simulation.
56 In this case, the whole operator must be decomposed into native
57 operations.
59 In both cases, :math:`H` may be represented by a
60 QubitPauliOperator.
61 """
63 def __init__(
64 self,
65 dictionary: dict[QubitPauliString, CoeffTypeAccepted] | None = None,
66 ) -> None:
67 self._dict: dict[QubitPauliString, Expr] = dict()
68 if dictionary:
69 for key, value in dictionary.items():
70 self._dict[key] = _coeff_convert(value)
71 self._collect_qubits()
73 def __repr__(self) -> str:
74 return self._dict.__repr__()
76 def __getitem__(self, key: QubitPauliString) -> Expr:
77 return self._dict[key]
79 def get(self, key: QubitPauliString, default: CoeffTypeAccepted) -> Expr:
80 return self._dict.get(key, _coeff_convert(default))
82 def __setitem__(self, key: QubitPauliString, value: CoeffTypeAccepted) -> None:
83 """Update value in dictionary ([]). Automatically converts value into sympy
84 Expr.
86 :param key: String to use as key
87 :type key: QubitPauliString
88 :param value: Associated coefficient
89 :type value: Union[int, float, complex, Expr]
90 """
91 self._dict[key] = _coeff_convert(value)
92 self._all_qubits.update(key.map.keys())
94 def __getstate__(self) -> dict[QubitPauliString, Expr]:
95 return self._dict
97 def __setstate__(self, _dict: dict[QubitPauliString, Expr]) -> None:
98 # values assumed to be already sympified
99 self._dict = _dict
100 self._collect_qubits()
102 def __eq__(self, other: object) -> bool:
103 if isinstance(other, QubitPauliOperator): 103 ↛ 105line 103 didn't jump to line 105 because the condition on line 103 was always true
104 return self._dict == other._dict
105 return False
107 def __iadd__(self, addend: "QubitPauliOperator") -> "QubitPauliOperator":
108 """In-place addition (+=) of QubitPauliOperators.
110 :param addend: The operator to add
111 :type addend: QubitPauliOperator
112 :return: Updated operator (self)
113 :rtype: QubitPauliOperator
114 """
115 if isinstance(addend, QubitPauliOperator): 115 ↛ 120line 115 didn't jump to line 120 because the condition on line 115 was always true
116 for key, value in addend._dict.items():
117 self[key] = self.get(key, 0.0) + value
118 self._all_qubits.update(addend._all_qubits)
119 else:
120 raise TypeError(f"Cannot add {type(addend)} to QubitPauliOperator.")
122 return self
124 def __add__(self, addend: "QubitPauliOperator") -> "QubitPauliOperator":
125 """Addition (+) of QubitPauliOperators.
127 :param addend: The operator to add
128 :type addend: QubitPauliOperator
129 :return: Sum operator
130 :rtype: QubitPauliOperator
131 """
132 summand = copy.deepcopy(self)
133 summand += addend
134 return summand
136 def __imul__(
137 self, multiplier: Union[float, Expr, "QubitPauliOperator"]
138 ) -> "QubitPauliOperator":
139 """In-place multiplication (*=) with QubitPauliOperator or scalar.
140 Multiply coefficients and terms.
142 :param multiplier: The operator or scalar to multiply
143 :type multiplier: Union[QubitPauliOperator, int, float, complex, Expr]
144 :return: Updated operator (self)
145 :rtype: QubitPauliOperator
146 """
148 # Handle operator of the same type
149 if isinstance(multiplier, QubitPauliOperator):
150 result_terms: dict = dict()
151 for left_key, left_value in self._dict.items():
152 for right_key, right_value in multiplier._dict.items():
153 new_term, bonus_coeff = pauli_string_mult(left_key, right_key)
154 new_coefficient = bonus_coeff * left_value * right_value
156 # Update result dict.
157 if new_term in result_terms: 157 ↛ 158line 157 didn't jump to line 158 because the condition on line 157 was never true
158 result_terms[new_term] += new_coefficient
159 else:
160 result_terms[new_term] = new_coefficient
161 self._dict = result_terms
162 self._all_qubits.update(multiplier._all_qubits)
163 return self
165 # Handle scalars.
166 if isinstance(multiplier, (float, Expr)): 166 ↛ 172line 166 didn't jump to line 172 because the condition on line 166 was always true
167 for key in self._dict:
168 self[key] *= multiplier
169 return self
171 # Invalid multiplier type
172 raise TypeError(f"Cannot multiply QubitPauliOperator with {type(multiplier)}")
174 def __mul__(
175 self, multiplier: Union[float, Expr, "QubitPauliOperator"]
176 ) -> "QubitPauliOperator":
177 """Multiplication (*) by QubitPauliOperator or scalar.
179 :param multiplier: The scalar to multiply by
180 :type multiplier: Union[int, float, complex, Expr, QubitPauliOperator]
181 :return: Product operator
182 :rtype: QubitPauliOperator
183 """
184 product = copy.deepcopy(self)
185 product *= multiplier
186 return product
188 def __rmul__(self, multiplier: CoeffTypeAccepted) -> "QubitPauliOperator":
189 """Multiplication (*) by a scalar.
190 We only define __rmul__ for scalars because left multiply is
191 queried as default behaviour, and is used for
192 QubitPauliOperator*QubitPauliOperator.
194 :param multiplier: The scalar to multiply by
195 :type multiplier: Union[int, float, complex, Expr]
196 :return: Product operator
197 :rtype: QubitPauliOperator
198 """
199 return self.__mul__(_coeff_convert(multiplier))
201 @property
202 def all_qubits(self) -> set[Qubit]:
203 """
204 :return: The set of all qubits the operator ranges over (including qubits
205 that were provided explicitly as identities)
207 :rtype: Set[Qubit]
208 """
209 return self._all_qubits
211 def subs(self, symbol_dict: dict[Symbol, complex]) -> None:
212 """Substitutes any matching symbols in the QubitPauliOperator.
214 :param symbol_dict: A dictionary of symbols to fixed values.
215 :type symbol_dict: Dict[Symbol, complex]
216 """
217 for key, value in self._dict.items():
218 self._dict[key] = value.subs(symbol_dict)
220 def get_dict(self) -> dict[QubitPauliString, Expr]:
221 """Generate a dict representation of QubitPauliOperator,
222 mapping each :py:class:`QubitPauliString` in the support
223 to its corresponding value.
225 :return: A dict of Pauli strings and their coefficients
226 as key-value pairs
227 """
228 return self._dict
230 def to_list(self) -> list[dict[str, Any]]:
231 """Generate a list serialized representation of QubitPauliOperator,
232 suitable for writing to JSON.
234 :return: JSON serializable list of dictionaries.
235 :rtype: List[Dict[str, Any]]
236 """
237 ret: list[dict[str, Any]] = []
238 for k, v in self._dict.items():
239 try:
240 coeff = complex_to_list(complex(v))
241 except TypeError:
242 assert isinstance(Expr(v), Expr)
243 coeff = str(v)
244 ret.append(
245 {
246 "string": k.to_list(),
247 "coefficient": coeff,
248 }
249 )
250 return ret
252 @classmethod
253 def from_list(cls, pauli_list: list[dict[str, Any]]) -> "QubitPauliOperator":
254 """Construct a QubitPauliOperator from a serializable JSON list format,
255 as returned by QubitPauliOperator.to_list()
257 :return: New QubitPauliOperator instance.
258 :rtype: QubitPauliOperator
259 """
261 def get_qps(obj: dict[str, Any]) -> QubitPauliString:
262 return QubitPauliString.from_list(obj["string"])
264 def get_coeff(obj: dict[str, Any]) -> Expr:
265 coeff = obj["coefficient"]
266 if type(coeff) is str:
267 return _coeff_convert(coeff)
268 return _coeff_convert(list_to_complex(coeff))
270 return QubitPauliOperator({get_qps(obj): get_coeff(obj) for obj in pauli_list})
272 def to_sparse_matrix(self, qubits: list[Qubit] | int | None = None) -> "csc_matrix":
273 """Represents the sparse operator as a dense operator under the ordering
274 scheme specified by ``qubits``, and generates the corresponding matrix.
276 - When ``qubits`` is an explicit list, the qubits are ordered with
277 ``qubits[0]`` as the most significant qubit for indexing into the matrix.
278 - If ``None``, then no padding qubits are introduced and we use the ILO-BE
279 convention, e.g. ``Qubit("a", 0)`` is more significant than
280 ``Qubit("a", 1)`` or ``Qubit("b")``.
281 - Giving a number specifies the number of qubits to use in the final
282 operator, treated as sequentially indexed from 0 in the default register
283 (padding with identities as necessary) and ordered by ILO-BE so
284 ``Qubit(0)`` is the most significant.
286 :param qubits: Sequencing of qubits in the matrix, either as an explicit
287 list, number of qubits to pad to, or infer from the operator.
288 Defaults to None
289 :type qubits: Union[List[Qubit], int, None], optional
290 :return: A sparse matrix representation of the operator.
291 :rtype: csc_matrix
292 """
293 if qubits is None:
294 qubits_ = sorted(list(self._all_qubits))
295 return sum(
296 complex(coeff) * pauli.to_sparse_matrix(qubits_)
297 for pauli, coeff in self._dict.items()
298 )
299 return sum(
300 complex(coeff) * pauli.to_sparse_matrix(qubits)
301 for pauli, coeff in self._dict.items()
302 )
304 def dot_state(
305 self, state: np.ndarray, qubits: list[Qubit] | None = None
306 ) -> np.ndarray:
307 """Applies the operator to the given state, mapping qubits to indexes
308 according to ``qubits``.
310 - When ``qubits`` is an explicit list, the qubits are ordered with
311 ``qubits[0]`` as the most significant qubit for indexing into ``state``.
312 - If ``None``, qubits sequentially indexed from 0 in the default register
313 and ordered by ILO-BE so ``Qubit(0)`` is the most significant.
315 :param state: The initial statevector
316 :type state: numpy.ndarray
317 :param qubits: Sequencing of qubits in ``state``, if not mapped to the
318 default register. Defaults to None
319 :type qubits: Union[List[Qubit], None], optional
320 :return: The dot product of the operator with the statevector
321 :rtype: numpy.ndarray
322 """
323 if qubits:
324 product_sum = sum(
325 complex(coeff) * pauli.dot_state(state, qubits)
326 for pauli, coeff in self._dict.items()
327 )
328 else:
329 product_sum = sum(
330 complex(coeff) * pauli.dot_state(state)
331 for pauli, coeff in self._dict.items()
332 )
333 return product_sum if isinstance(product_sum, numpy.ndarray) else state
335 def state_expectation(
336 self, state: np.ndarray, qubits: list[Qubit] | None = None
337 ) -> complex:
338 """Calculates the expectation value of the given statevector with respect
339 to the operator, mapping qubits to indexes according to ``qubits``.
341 - When ``qubits`` is an explicit list, the qubits are ordered with
342 ``qubits[0]`` as the most significant qubit for indexing into ``state``.
343 - If ``None``, qubits sequentially indexed from 0 in the default register
344 and ordered by ILO-BE so ``Qubit(0)`` is the most significant.
346 :param state: The initial statevector
347 :type state: numpy.ndarray
348 :param qubits: Sequencing of qubits in ``state``, if not mapped to the
349 default register. Defaults to None
350 :type qubits: Union[List[Qubit], None], optional
351 :return: The expectation value of the statevector and operator
352 :rtype: complex
353 """
354 if qubits:
355 return sum(
356 complex(coeff) * pauli.state_expectation(state, qubits)
357 for pauli, coeff in self._dict.items()
358 )
359 return sum(
360 complex(coeff) * pauli.state_expectation(state)
361 for pauli, coeff in self._dict.items()
362 )
364 def compress(self, abs_tol: float = 1e-10) -> None:
365 """Substitutes all free symbols in the QubitPauliOperator with
366 1, and then removes imaginary and real components which have
367 magnitudes below the tolerance. If the resulting expression is
368 0, the term is removed entirely.
370 Warning: This methods assumes significant expression structure
371 is known a priori, and is best suited to operators which have
372 simple product expressions, such as excitation operators for
373 VQE ansätze and digital quantum simulation. Otherwise, it may
374 remove terms relevant to computation. Each expression is of
375 the form :math:`f(a_1,a_2,\\ldots,a_n)` for some symbols
376 :math:`a_i`. :math:`|f(a_1,a_2,\\ldots,a_n)|` is assumed to
377 monotonically increase in both real and imaginary components
378 for all :math:`a_i \\in [0, 1]`.
380 :param abs_tol: The threshold below which to remove values.
381 :type abs_tol: float
382 """
384 to_delete = []
385 for key, value in self._dict.items():
386 placeholder = value.subs({s: 1 for s in value.free_symbols})
387 if abs(re(placeholder)) <= abs_tol:
388 if abs(im(placeholder)) <= abs_tol:
389 to_delete.append(key)
390 else:
391 self._dict[key] = im(value) * 1j
392 elif abs(im(placeholder)) <= abs_tol: 392 ↛ 385line 392 didn't jump to line 385 because the condition on line 392 was always true
393 self._dict[key] = re(value)
395 for key in to_delete:
396 del self._dict[key]
398 def _collect_qubits(self) -> None:
399 self._all_qubits: set[Qubit] = set()
400 for key in self._dict.keys():
401 for q in key.map.keys():
402 self._all_qubits.add(q)