Coverage for qml_essentials/entanglement.py: 79%
145 statements
« prev ^ index » next coverage.py v7.9.2, created at 2025-07-04 21:13 +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 mw_measure = np.zeros(len(rhos))
65 for i, rho in enumerate(rhos):
66 mw_measure[i] = Entanglement._compute_meyer_wallach_meas(
67 rho, model.n_qubits
68 )
70 # Average all iterated states
71 entangling_capability = min(max(mw_measure.mean(), 0.0), 1.0)
72 log.debug(f"Variance of measure: {mw_measure.var()}")
74 # catch floating point errors
75 return float(entangling_capability)
77 @staticmethod
78 def _compute_meyer_wallach_meas(rho: np.ndarray, n_qubits: int):
79 qb = list(range(n_qubits))
80 entropy = 0
81 for j in range(n_qubits):
82 # Formula 6 in https://doi.org/10.48550/arXiv.quant-ph/0305094
83 density = qml.math.partial_trace(rho, qb[:j] + qb[j + 1 :])
84 # only real values, because imaginary part will be separate
85 # in all following calculations anyway
86 # entropy should be 1/2 <= entropy <= 1
87 entropy += np.trace((density @ density).real)
89 # inverse averaged entropy and scale to [0, 1]
90 return 2 * (1 - entropy / n_qubits)
92 @staticmethod
93 def bell_measurements(
94 model: Model, n_samples: int, seed: int, scale: bool = False, **kwargs: Any
95 ) -> float:
96 """
97 Compute the Bell measurement for a given model.
99 Args:
100 model (Model): The quantum circuit model.
101 n_samples (int): The number of samples to compute the measure for.
102 seed (int): The seed for the random number generator.
103 scale (bool): Whether to scale the number of samples
104 according to the number of qubits.
105 **kwargs (Any): Additional keyword arguments for the model function.
107 Returns:
108 float: The Bell measurement value.
109 """
110 if "noise_params" in kwargs:
111 log.warning(
112 "Bell Measurements not suitable for noisy circuits.\
113 Consider 'relative_entropy' instead."
114 )
116 if scale:
117 n_samples = np.power(2, model.n_qubits) * n_samples
119 def _circuit(
120 params: np.ndarray, inputs: np.ndarray,
121 enc_params: Optional[np.ndarray] = None
122 ) -> List[np.ndarray]:
123 """
124 Compute the Bell measurement circuit.
126 Args:
127 params (np.ndarray): The model parameters.
128 inputs (np.ndarray): The input to the model.
129 enc_params (Optional[np.ndarray]): The frequency encoding parameters.
131 Returns:
132 List[np.ndarray]: The probabilities of the Bell measurement.
133 """
134 model._variational(params, inputs, enc_params)
136 qml.map_wires(
137 model._variational,
138 {i: i + model.n_qubits for i in range(model.n_qubits)},
139 )(params, inputs)
141 for q in range(model.n_qubits):
142 qml.CNOT(wires=[q, q + model.n_qubits])
143 qml.H(q)
145 obs_wires = [(q, q + model.n_qubits) for q in range(model.n_qubits)]
146 return [qml.probs(wires=w) for w in obs_wires]
148 model.circuit = qml.QNode(
149 _circuit,
150 qml.device(
151 "default.qubit",
152 shots=model.shots,
153 wires=model.n_qubits * 2,
154 ),
155 )
157 rng = np.random.default_rng(seed)
158 if n_samples is not None and n_samples > 0:
159 assert seed is not None, "Seed must be provided when samples > 0"
160 # TODO: maybe switch to JAX rng
161 model.initialize_params(rng=rng, repeat=n_samples)
162 params = model.params
163 else:
164 if seed is not None:
165 log.warning("Seed is ignored when samples is 0")
167 if len(model.params.shape) <= 2:
168 params = model.params.reshape(*model.params.shape, 1)
169 else:
170 log.info(f"Using sample size of model params: {model.params.shape[-1]}")
171 params = model.params
173 n_samples = params.shape[-1]
174 mw_measure = np.zeros(n_samples)
176 # implicitly set input to none in case it's not needed
177 kwargs.setdefault("inputs", None)
178 exp = model(params=params, **kwargs)
179 exp = 1 - 2 * exp[:, :, -1]
180 mw_measure = 2 * (1 - exp.mean(axis=0))
181 entangling_capability = min(max(mw_measure.mean(), 0.0), 1.0)
182 log.debug(f"Variance of measure: {mw_measure.var()}")
184 return float(entangling_capability)
186 @staticmethod
187 def relative_entropy(
188 model: Model,
189 n_samples: int,
190 n_sigmas: int,
191 seed: Optional[int],
192 scale: bool = False,
193 **kwargs: Any,
194 ) -> float:
195 """
196 Calculates the relative entropy of entanglement of a given quantum
197 circuit. This measure is also applicable to mixed state, albeit it
198 might me not fully accurate in this simplified case.
200 As the relative entropy is generally defined as the smallest relative
201 entropy from the state in question to the set of separable states.
202 However, as computing the nearest separable state is NP-hard, we select
203 n_sigmas of random separable states to compute the distance to, which
204 is not necessarily the nearest. Thus, this measure of entanglement
205 presents an upper limit of entanglement.
207 As the relative entropy is not necessarily between zero and one, this
208 function also normalises by the relative entroy to the GHZ state.
210 Args:
211 model (Model): The quantum circuit model.
212 n_samples (int): Number of samples per qubit.
213 If <= 0, the current parameters of the model are used.
214 n_sigmas (int): Number of random separable pure states to compare against.
215 seed (Optional[int]): Seed for the random number generator.
216 scale (bool): Whether to scale the number of samples.
217 kwargs (Any): Additional keyword arguments for the model function.
219 Returns:
220 float: Entangling capacity of the given circuit, guaranteed
221 to be between 0.0 and 1.0.
222 """
223 dim = np.power(2, model.n_qubits)
224 if scale:
225 n_samples = dim * n_samples
226 n_sigmas = dim * n_sigmas
228 rng = np.random.default_rng(seed)
230 # Random separable states
231 log_sigmas = sample_random_separable_states(
232 model.n_qubits, n_samples=n_sigmas, rng=rng, take_log=True
233 )
235 if n_samples > 0:
236 assert seed is not None, "Seed must be provided when samples > 0"
237 model.initialize_params(rng=rng, repeat=n_samples)
238 else:
239 if seed is not None:
240 log.warning("Seed is ignored when samples is 0")
242 if len(model.params.shape) <= 2:
243 model.params = model.params.reshape(*model.params.shape, 1)
244 else:
245 log.info(f"Using sample size of model params: {model.params.shape[-1]}")
247 ghz_model = Model(model.n_qubits, 1, "GHZ", data_reupload=False)
249 normalised_entropies = np.zeros((n_sigmas, n_samples))
250 for j, log_sigma in enumerate(log_sigmas):
252 # Entropy of GHZ states should be maximal
253 ghz_entropy = Entanglement._compute_rel_entropies(
254 ghz_model,
255 log_sigma,
256 )
258 rel_entropy = Entanglement._compute_rel_entropies(
259 model, log_sigma, **kwargs
260 )
262 normalised_entropies[j] = rel_entropy / ghz_entropy
264 # Average all iterated states
265 entangling_capability = normalised_entropies.min(axis=0).mean()
266 log.debug(f"Variance of measure: {normalised_entropies.var()}")
268 return entangling_capability
270 @staticmethod
271 def _compute_rel_entropies(
272 model: Model,
273 log_sigma: np.ndarray,
274 **kwargs,
275 ) -> np.ndarray:
276 """
277 Compute the relative entropy for a given model.
279 Args:
280 model (Model): The model for which to compute entanglement
281 log_sigma (np.ndarray): Density matrix of next separable state
283 Returns:
284 np.ndarray: Relative Entropy for each sample
285 """
286 # implicitly set input to none in case it's not needed
287 kwargs.setdefault("inputs", None)
288 # explicitly set execution type because everything else won't work
289 rho = model(execution_type="density", **kwargs)
290 rho = rho.reshape(-1, 2**model.n_qubits, 2**model.n_qubits)
291 log_rho = logm_v(rho) / np.log(2)
293 rel_entropies = np.abs(np.trace(rho @ (log_rho - log_sigma), axis1=1, axis2=2))
295 return rel_entropies
297 @staticmethod
298 def entanglement_of_formation(
299 model: Model,
300 n_samples: int,
301 seed: Optional[int],
302 scale: bool = False,
303 always_decompose: bool = False,
304 **kwargs: Any,
305 ) -> float:
306 """
307 This function implements the entanglement of formation for mixed
308 quantum systems.
309 In that a mixed state gets decomposed into pure states with respective
310 probabilities using the eigendecomposition of the density matrix.
311 Then, the Meyer-Wallach measure is computed for each pure state,
312 weighted by the eigenvalue.
313 See e.g. https://doi.org/10.48550/arXiv.quant-ph/0504163
315 Note that the decomposition is *not unique*! Therefore, this measure
316 presents the entanglement for *some* decomposition into pure states,
317 not necessarily the one that is anticipated when applying the Kraus
318 channels.
319 If a pure state is provided, this results in the same value as the
320 Entanglement.meyer_wallach function if `always_decompose` flag is not set.
322 Args:
323 model (Model): The quantum circuit model.
324 n_samples (int): Number of samples per qubit.
325 seed (Optional[int]): Seed for the random number generator.
326 scale (bool): Whether to scale the number of samples.
327 always_decompose (bool): Whether to explicitly compute the
328 entantlement of formation for the eigendecomposition of a pure
329 state.
330 kwargs (Any): Additional keyword arguments for the model function.
332 Returns:
333 float: Entangling capacity of the given circuit, guaranteed
334 to be between 0.0 and 1.0.
335 """
337 if scale:
338 n_samples = np.power(2, model.n_qubits) * n_samples
340 rng = np.random.default_rng(seed)
341 if n_samples > 0:
342 assert seed is not None, "Seed must be provided when samples > 0"
343 model.initialize_params(rng=rng, repeat=n_samples)
344 else:
345 if seed is not None:
346 log.warning("Seed is ignored when samples is 0")
348 if len(model.params.shape) <= 2:
349 model.params = model.params.reshape(*model.params.shape, 1)
350 else:
351 log.info(f"Using sample size of model params: {model.params.shape[-1]}")
353 # implicitly set input to none in case it's not needed
354 kwargs.setdefault("inputs", None)
355 rhos = model(execution_type="density", **kwargs)
356 rhos = rhos.reshape(-1, 2**model.n_qubits, 2**model.n_qubits)
357 entanglement = np.zeros(len(rhos))
358 for i, rho in enumerate(rhos):
359 entanglement[i] = Entanglement._compute_entanglement_of_formation(
360 rho, model.n_qubits, always_decompose
361 )
362 entangling_capability = min(max(entanglement.mean(), 0.0), 1.0)
363 return float(entangling_capability)
365 @staticmethod
366 def _compute_entanglement_of_formation(
367 rho: np.ndarray, n_qubits: int, always_decompose: bool
368 ) -> float:
369 """
370 Computes the entanglement of formation for a given density matrix rho.
372 Args:
373 rho (np.ndarray): The density matrix
374 n_qubits (int): Number of qubits
375 always_decompose (bool): Whether to explicitly compute the
376 entantlement of formation for the eigendecomposition of a pure
377 state.
379 Returns:
380 float: Entanglement for the provided state.
381 """
382 eigenvalues, eigenvectors = np.linalg.eigh(rho)
383 if any(np.isclose(eigenvalues, 1.0)) and not always_decompose: # Pure state
384 return Entanglement._compute_meyer_wallach_meas(rho, n_qubits)
385 ent = 0
386 for prob, ev in zip(eigenvalues, eigenvectors):
387 ev = ev.reshape(-1, 1)
388 rho = ev @ np.conjugate(ev).T
389 mw_measure = Entanglement._compute_meyer_wallach_meas(rho, n_qubits)
390 ent += prob * mw_measure
391 return ent
394def sample_random_separable_states(
395 n_qubits: int, n_samples: int, rng: np.random.Generator, take_log: bool = False
396) -> np.ndarray:
397 """
398 Sample random separable states (density matrix).
400 Args:
401 n_qubits (int): number of qubits in the state
402 n_samples (int): number of states
403 rng (np.random.Generator): random number generator
404 take_log (bool): if the matrix logarithm of the density matrix should be taken.
406 Returns:
407 np.ndarray: Density matrices of shape (n_samples, 2**n_qubits, 2**n_qubits)
408 """
409 model = Model(n_qubits, 1, "No_Entangling", data_reupload=False)
410 model.initialize_params(rng=rng, repeat=n_samples)
411 # explicitly set execution type because everything else won't work
412 sigmas = model(execution_type="density", inputs=None)
413 if take_log:
414 sigmas = logm_v(sigmas) / np.log(2)
416 return sigmas