Coverage for qml_essentials/entanglement.py: 90%
185 statements
« prev ^ index » next coverage.py v7.9.2, created at 2025-10-02 13:10 +0000
1from typing import Optional, Any, List
2import pennylane as qml
3import pennylane.numpy as np
4from copy import deepcopy
5from qml_essentials.utils import logm_v
6from qml_essentials.model import Model
7import logging
9log = logging.getLogger(__name__)
12class Entanglement:
14 @staticmethod
15 def meyer_wallach(
16 model: Model,
17 n_samples: Optional[int | None],
18 seed: Optional[int],
19 scale: bool = False,
20 **kwargs: Any,
21 ) -> float:
22 """
23 Calculates the entangling capacity of a given quantum circuit
24 using Meyer-Wallach measure.
26 Args:
27 model (Model): The quantum circuit model.
28 n_samples (Optional[int]): Number of samples per qubit.
29 If None or < 0, the current parameters of the model are used.
30 seed (Optional[int]): Seed for the random number generator.
31 scale (bool): Whether to scale the number of samples.
32 kwargs (Any): Additional keyword arguments for the model function.
34 Returns:
35 float: Entangling capacity of the given circuit, guaranteed
36 to be between 0.0 and 1.0.
37 """
38 if "noise_params" in kwargs:
39 log.warning(
40 "Meyer-Wallach measure not suitable for noisy circuits.\
41 Consider 'relative_entropy' instead."
42 )
44 if scale:
45 n_samples = np.power(2, model.n_qubits) * n_samples
47 rng = np.random.default_rng(seed)
48 if n_samples is not None and n_samples > 0:
49 assert seed is not None, "Seed must be provided when samples > 0"
50 # TODO: maybe switch to JAX rng
51 model.initialize_params(rng=rng, repeat=n_samples)
52 else:
53 if seed is not None:
54 log.warning("Seed is ignored when samples is 0")
56 # implicitly set input to none in case it's not needed
57 kwargs.setdefault("inputs", None)
58 # explicitly set execution type because everything else won't work
59 rhos = model(execution_type="density", **kwargs).reshape(
60 -1, 2**model.n_qubits, 2**model.n_qubits
61 )
63 measure = np.zeros(len(rhos))
65 for i, rho in enumerate(rhos):
66 measure[i] = Entanglement._compute_meyer_wallach_meas(rho, model.n_qubits)
68 # Average all iterated states
69 entangling_capability = min(max(measure.mean(), 0.0), 1.0)
70 log.debug(f"Variance of measure: {measure.var()}")
72 # catch floating point errors
73 return float(entangling_capability)
75 @staticmethod
76 def _compute_meyer_wallach_meas(rho: np.ndarray, n_qubits: int):
77 qb = list(range(n_qubits))
78 entropy = 0
79 for j in range(n_qubits):
80 # Formula 6 in https://doi.org/10.48550/arXiv.quant-ph/0305094
81 density = qml.math.partial_trace(rho, qb[:j] + qb[j + 1 :])
82 # only real values, because imaginary part will be separate
83 # in all following calculations anyway
84 # entropy should be 1/2 <= entropy <= 1
85 entropy += np.trace((density @ density).real)
87 # inverse averaged entropy and scale to [0, 1]
88 return 2 * (1 - entropy / n_qubits)
90 @staticmethod
91 def bell_measurements(
92 model: Model, n_samples: int, seed: int, scale: bool = False, **kwargs: Any
93 ) -> float:
94 """
95 Compute the Bell measurement for a given model.
97 Args:
98 model (Model): The quantum circuit model.
99 n_samples (int): The number of samples to compute the measure for.
100 seed (int): The seed for the random number generator.
101 scale (bool): Whether to scale the number of samples
102 according to the number of qubits.
103 **kwargs (Any): Additional keyword arguments for the model function.
105 Returns:
106 float: The Bell measurement value.
107 """
108 if "noise_params" in kwargs:
109 log.warning(
110 "Bell Measurements not suitable for noisy circuits.\
111 Consider 'relative_entropy' instead."
112 )
114 if scale:
115 n_samples = np.power(2, model.n_qubits) * n_samples
117 def _circuit(
118 params: np.ndarray,
119 inputs: np.ndarray,
120 pulse_params: Optional[np.ndarray] = None,
121 enc_params: Optional[np.ndarray] = None,
122 gate_mode: str = "unitary",
123 ) -> List[np.ndarray]:
124 """
125 Compute the Bell measurement circuit.
127 Args:
128 params (np.ndarray): The model parameters.
129 inputs (np.ndarray): The input to the model.
130 pulse_params (np.ndarray): The model pulse parameters.
131 enc_params (Optional[np.ndarray]): The frequency encoding parameters.
133 Returns:
134 List[np.ndarray]: The probabilities of the Bell measurement.
135 """
136 model._variational(params, inputs, pulse_params, enc_params, gate_mode)
138 qml.map_wires(
139 model._variational,
140 {i: i + model.n_qubits for i in range(model.n_qubits)},
141 )(params, inputs)
143 for q in range(model.n_qubits):
144 qml.CNOT(wires=[q, q + model.n_qubits])
145 qml.H(q)
147 obs_wires = [(q, q + model.n_qubits) for q in range(model.n_qubits)]
148 return [qml.probs(wires=w) for w in obs_wires]
150 model.circuit = qml.QNode(
151 _circuit,
152 qml.device(
153 "default.qubit",
154 shots=model.shots,
155 wires=model.n_qubits * 2,
156 ),
157 )
159 rng = np.random.default_rng(seed)
160 if n_samples is not None and n_samples > 0:
161 assert seed is not None, "Seed must be provided when samples > 0"
162 # TODO: maybe switch to JAX rng
163 model.initialize_params(rng=rng, repeat=n_samples)
164 params = model.params
165 else:
166 if seed is not None:
167 log.warning("Seed is ignored when samples is 0")
169 if len(model.params.shape) <= 2:
170 params = model.params.reshape(*model.params.shape, 1)
171 else:
172 log.info(f"Using sample size of model params: {model.params.shape[-1]}")
173 params = model.params
175 n_samples = params.shape[-1]
176 measure = np.zeros(n_samples)
178 # implicitly set input to none in case it's not needed
179 kwargs.setdefault("inputs", None)
180 exp = model(params=params, **kwargs)
181 exp = 1 - 2 * exp[..., -1]
182 measure = 2 * (1 - exp.mean(axis=0))
183 entangling_capability = min(max(measure.mean(), 0.0), 1.0)
184 log.debug(f"Variance of measure: {measure.var()}")
186 return float(entangling_capability)
188 @staticmethod
189 def relative_entropy(
190 model: Model,
191 n_samples: int,
192 n_sigmas: int,
193 seed: Optional[int],
194 scale: bool = False,
195 **kwargs: Any,
196 ) -> float:
197 """
198 Calculates the relative entropy of entanglement of a given quantum
199 circuit. This measure is also applicable to mixed state, albeit it
200 might me not fully accurate in this simplified case.
202 As the relative entropy is generally defined as the smallest relative
203 entropy from the state in question to the set of separable states.
204 However, as computing the nearest separable state is NP-hard, we select
205 n_sigmas of random separable states to compute the distance to, which
206 is not necessarily the nearest. Thus, this measure of entanglement
207 presents an upper limit of entanglement.
209 As the relative entropy is not necessarily between zero and one, this
210 function also normalises by the relative entroy to the GHZ state.
212 Args:
213 model (Model): The quantum circuit model.
214 n_samples (int): Number of samples per qubit.
215 If <= 0, the current parameters of the model are used.
216 n_sigmas (int): Number of random separable pure states to compare against.
217 seed (Optional[int]): Seed for the random number generator.
218 scale (bool): Whether to scale the number of samples.
219 kwargs (Any): Additional keyword arguments for the model function.
221 Returns:
222 float: Entangling capacity of the given circuit, guaranteed
223 to be between 0.0 and 1.0.
224 """
225 dim = np.power(2, model.n_qubits)
226 if scale:
227 n_samples = dim * n_samples
228 n_sigmas = dim * n_sigmas
230 rng = np.random.default_rng(seed)
232 # Random separable states
233 log_sigmas = sample_random_separable_states(
234 model.n_qubits, n_samples=n_sigmas, rng=rng, take_log=True
235 )
237 if n_samples is not None and n_samples > 0:
238 assert seed is not None, "Seed must be provided when samples > 0"
239 model.initialize_params(rng=rng, repeat=n_samples)
240 else:
241 if seed is not None:
242 log.warning("Seed is ignored when samples is 0")
244 if len(model.params.shape) <= 2:
245 model.params = model.params.reshape(*model.params.shape, 1)
246 else:
247 log.info(f"Using sample size of model params: {model.params.shape[-1]}")
249 ghz_model = Model(model.n_qubits, 1, "GHZ", data_reupload=False)
251 normalised_entropies = np.zeros((n_sigmas, model.params.shape[-1]))
252 for j, log_sigma in enumerate(log_sigmas):
254 # Entropy of GHZ states should be maximal
255 ghz_entropy = Entanglement._compute_rel_entropies(
256 ghz_model,
257 log_sigma,
258 )
260 rel_entropy = Entanglement._compute_rel_entropies(
261 model, log_sigma, **kwargs
262 )
264 normalised_entropies[j] = rel_entropy / ghz_entropy
266 # Average all iterated states
267 entangling_capability = normalised_entropies.min(axis=0).mean()
268 log.debug(f"Variance of measure: {normalised_entropies.var()}")
270 return entangling_capability
272 @staticmethod
273 def _compute_rel_entropies(
274 model: Model,
275 log_sigma: np.ndarray,
276 **kwargs,
277 ) -> np.ndarray:
278 """
279 Compute the relative entropy for a given model.
281 Args:
282 model (Model): The model for which to compute entanglement
283 log_sigma (np.ndarray): Density matrix of next separable state
285 Returns:
286 np.ndarray: Relative Entropy for each sample
287 """
288 # implicitly set input to none in case it's not needed
289 kwargs.setdefault("inputs", None)
290 # explicitly set execution type because everything else won't work
291 rho = model(execution_type="density", **kwargs)
292 rho = rho.reshape(-1, 2**model.n_qubits, 2**model.n_qubits)
293 log_rho = logm_v(rho) / np.log(2)
295 rel_entropies = np.abs(np.trace(rho @ (log_rho - log_sigma), axis1=1, axis2=2))
297 return rel_entropies
299 @staticmethod
300 def entanglement_of_formation(
301 model: Model,
302 n_samples: int,
303 seed: Optional[int],
304 scale: bool = False,
305 always_decompose: bool = False,
306 **kwargs: Any,
307 ) -> float:
308 """
309 This function implements the entanglement of formation for mixed
310 quantum systems.
311 In that a mixed state gets decomposed into pure states with respective
312 probabilities using the eigendecomposition of the density matrix.
313 Then, the Meyer-Wallach measure is computed for each pure state,
314 weighted by the eigenvalue.
315 See e.g. https://doi.org/10.48550/arXiv.quant-ph/0504163
317 Note that the decomposition is *not unique*! Therefore, this measure
318 presents the entanglement for *some* decomposition into pure states,
319 not necessarily the one that is anticipated when applying the Kraus
320 channels.
321 If a pure state is provided, this results in the same value as the
322 Entanglement.meyer_wallach function if `always_decompose` flag is not set.
324 Args:
325 model (Model): The quantum circuit model.
326 n_samples (int): Number of samples per qubit.
327 seed (Optional[int]): Seed for the random number generator.
328 scale (bool): Whether to scale the number of samples.
329 always_decompose (bool): Whether to explicitly compute the
330 entantlement of formation for the eigendecomposition of a pure
331 state.
332 kwargs (Any): Additional keyword arguments for the model function.
334 Returns:
335 float: Entangling capacity of the given circuit, guaranteed
336 to be between 0.0 and 1.0.
337 """
339 if scale:
340 n_samples = np.power(2, model.n_qubits) * n_samples
342 rng = np.random.default_rng(seed)
343 if n_samples is not None and n_samples > 0:
344 assert seed is not None, "Seed must be provided when samples > 0"
345 model.initialize_params(rng=rng, repeat=n_samples)
346 else:
347 if seed is not None:
348 log.warning("Seed is ignored when samples is 0")
350 if len(model.params.shape) <= 2:
351 model.params = model.params.reshape(*model.params.shape, 1)
352 else:
353 log.info(f"Using sample size of model params: {model.params.shape[-1]}")
355 # implicitly set input to none in case it's not needed
356 kwargs.setdefault("inputs", None)
357 rhos = model(execution_type="density", **kwargs)
358 rhos = rhos.reshape(-1, 2**model.n_qubits, 2**model.n_qubits)
359 entanglement = np.zeros(len(rhos))
360 for i, rho in enumerate(rhos):
361 entanglement[i] = Entanglement._compute_entanglement_of_formation(
362 rho, model.n_qubits, always_decompose
363 )
364 entangling_capability = min(max(entanglement.mean(), 0.0), 1.0)
365 return float(entangling_capability)
367 @staticmethod
368 def _compute_entanglement_of_formation(
369 rho: np.ndarray, n_qubits: int, always_decompose: bool
370 ) -> float:
371 """
372 Computes the entanglement of formation for a given density matrix rho.
374 Args:
375 rho (np.ndarray): The density matrix
376 n_qubits (int): Number of qubits
377 always_decompose (bool): Whether to explicitly compute the
378 entantlement of formation for the eigendecomposition of a pure
379 state.
381 Returns:
382 float: Entanglement for the provided state.
383 """
384 eigenvalues, eigenvectors = np.linalg.eigh(rho)
385 if any(np.isclose(eigenvalues, 1.0)) and not always_decompose: # Pure state
386 return Entanglement._compute_meyer_wallach_meas(rho, n_qubits)
387 ent = 0
388 for prob, ev in zip(eigenvalues, eigenvectors):
389 ev = ev.reshape(-1, 1)
390 rho = ev @ np.conjugate(ev).T
391 measure = Entanglement._compute_meyer_wallach_meas(rho, n_qubits)
392 ent += prob * measure
393 return ent
395 @staticmethod
396 def concentratable_entanglement(
397 model: Model, n_samples: int, seed: int, scale: bool = False, **kwargs: Any
398 ) -> float:
399 """
400 Computes the concentratable entanglement of a given model.
402 This method utilizes the Concentratable Entanglement measure from
403 https://arxiv.org/abs/2104.06923.
405 Args:
406 model (Model): The quantum circuit model.
407 n_samples (int): The number of samples to compute the measure for.
408 seed (int): The seed for the random number generator.
409 scale (bool): Whether to scale the number of samples according to
410 the number of qubits.
411 **kwargs (Any): Additional keyword arguments for the model function.
413 Returns:
414 float: Entangling capability of the given circuit, guaranteed
415 to be between 0.0 and 1.0.
416 """
417 if "noise_params" in kwargs:
418 log.warning(
419 "Concentratable entanglement is not suitable for noisy circuits.\
420 Consider 'relative_entropy' instead."
421 )
423 n = model.n_qubits
425 if scale:
426 n_samples = np.power(2, model.n) * n_samples
428 def _circuit(
429 params: np.ndarray,
430 inputs: np.ndarray,
431 pulse_params: Optional[np.ndarray] = None,
432 enc_params: Optional[np.ndarray] = None,
433 gate_mode: str = "unitary",
434 ) -> List[np.ndarray]:
435 """
436 Constructs a circuit to compute the concentratable entanglement using the
437 swap test by creating two copies of the models circuit and map the output
438 wires accordingly
440 Args:
441 params (np.ndarray): The model parameters.
442 inputs (np.ndarray): The input data for the model.
443 pulse_params (np.ndarray): The model pulse parameters.
444 enc_params (Optional[np.ndarray]): Optional encoding parameters.
446 Returns:
447 List[np.ndarray]: Probabilities obtained from the swap test circuit.
448 """
450 qml.map_wires(model._variational, {i: i + n for i in range(n)})(
451 params, inputs, pulse_params, enc_params, gate_mode
452 )
453 qml.map_wires(model._variational, {i: i + 2 * n for i in range(n)})(
454 params, inputs, pulse_params, enc_params, gate_mode
455 )
457 # Perform swap test
458 for i in range(n):
459 qml.H(i)
461 for i in range(n):
462 qml.CSWAP([i, i + n, i + 2 * n])
464 for i in range(n):
465 qml.H(i)
467 return qml.probs(wires=[i for i in range(n)])
469 model.circuit = qml.QNode(
470 _circuit,
471 qml.device(
472 "default.qubit",
473 shots=model.shots,
474 wires=n * 3,
475 ),
476 )
478 rng = np.random.default_rng(seed)
479 if n_samples is not None and n_samples > 0:
480 assert seed is not None, "Seed must be provided when samples > 0"
481 model.initialize_params(rng=rng, repeat=n_samples)
482 params = model.params
483 else:
484 if seed is not None:
485 log.warning("Seed is ignored when samples is 0")
487 if len(model.params.shape) <= 2:
488 params = model.params.reshape(*model.params.shape, 1)
489 else:
490 log.info(f"Using sample size of model params: {model.params.shape[-1]}")
491 params = model.params
493 n_samples = params.shape[-1]
495 # implicitly set input to none in case it's not needed
496 kwargs.setdefault("inputs", None)
498 samples_probs = model(params=params, execution_type="probs", **kwargs)
499 if n_samples == 1:
500 samples_probs = [samples_probs]
502 ce_measure = np.zeros(len(samples_probs))
504 for i, probs in enumerate(samples_probs):
505 ce_measure[i] = 1 - probs[0]
507 # Average all iterated states
508 entangling_capability = min(max(ce_measure.mean(), 0.0), 1.0)
509 log.debug(f"Variance of measure: {ce_measure.var()}")
511 # catch floating point errors
512 return float(entangling_capability)
515def sample_random_separable_states(
516 n_qubits: int, n_samples: int, rng: np.random.Generator, take_log: bool = False
517) -> np.ndarray:
518 """
519 Sample random separable states (density matrix).
521 Args:
522 n_qubits (int): number of qubits in the state
523 n_samples (int): number of states
524 rng (np.random.Generator): random number generator
525 take_log (bool): if the matrix logarithm of the density matrix should be taken.
527 Returns:
528 np.ndarray: Density matrices of shape (n_samples, 2**n_qubits, 2**n_qubits)
529 """
530 model = Model(n_qubits, 1, "No_Entangling", data_reupload=False)
531 model.initialize_params(rng=rng, repeat=n_samples)
532 # explicitly set execution type because everything else won't work
533 sigmas = model(execution_type="density", inputs=None)
534 if take_log:
535 sigmas = logm_v(sigmas) / np.log(2)
537 return sigmas