Source code for pennylane.labs.dla.involutions
# Copyright 2024 Xanadu Quantum Technologies Inc.
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
# http://www.apache.org/licenses/LICENSE-2.0
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
"""Cartan involutions"""
from functools import singledispatch
# pylint: disable=missing-function-docstring
from typing import Optional, Union
import numpy as np
import pennylane as qml
from pennylane import Y, Z, matrix
from pennylane.operation import Operator
from pennylane.pauli import PauliSentence
[docs]def khaneja_glaser_involution(op: Union[np.ndarray, PauliSentence, Operator], wire: int = None):
r"""Khaneja-Glaser involution, which is a type-:func:`~AIII` Cartan involution with p=q.
Args:
op (Union[np.ndarray, PauliSentence, Operator]): Input operator
wire (int): Qubit wire on which to perform Khaneja-Glaser decomposition
Returns:
bool: Whether ``op`` should go to the even or odd subspace of the decomposition
.. seealso:: :func:`~cartan_decomp`
See `Khaneja and Glaser (2000) <https://arxiv.org/abs/quant-ph/0010100>`__ for reference.
**Example**
Let us first perform a single Khaneja-Glaser decomposition of
:math:`\mathfrak{g} = \mathfrak{su}(8)`, i.e. the Lie algebra of all Pauli words on 3 qubits.
>>> g = list(qml.pauli.pauli_group(3)) # u(8)
>>> g = g[1:] # remove identity
>>> g = [op.pauli_rep for op in g]
We perform the first iteration on the first qubit. We use :func:`~cartan_decomp`.
>>> from functools import partial
>>> from pennylane.labs.dla import cartan_decomp, khaneja_glaser_involution
>>> k0, m0 = cartan_decomp(g, partial(khaneja_glaser_involution, wire=0))
>>> print(f"First iteration, AIII: {len(k0)}, {len(m0)}")
First iteration, AIII: 31, 32
>>> assert qml.labs.dla.check_cartan_decomp(k0, m0) # check Cartan commutation relations
We see that we split the 63-dimensional algebra :math:`\mathfrak{su}(8)` into a 31-dimensional
and a 32-dimensional subspace, and :func:`~check_cartan_decomp` verified that they satisfy
the commutation relations of a Cartan decomposition.
To expand on this example, we continue process recursively on the :math:`\mathfrak{k}_0`
subalgebra. Before we can apply the Khaneja-Glaser (type AIII) decomposition a second time,
we need to a) remove the :math:`\mathfrak{u}(1)` center from
:math:`\mathfrak{k}_0=\mathfrak{su}(4)\oplus\mathfrak{su}(4)\oplus\mathfrak{u}(1)` to make it
semisimple and b) perform a different Cartan decomposition, which is sometimes termed class B.
>>> center_k0 = qml.center(k0, pauli=True) # Compute center of k0
>>> k0_semi = [op for op in k0 if op not in center_k0] # Remove center from k0
>>> print(f"Removed operators {center_k0}, new length: {len(k0_semi)}")
Removed operators [1.0 * Z(0)], new length: 30
>>> from pennylane.labs.dla import ClassB
>>> k1, m1 = cartan_decomp(k0_semi, partial(ClassB, wire=0))
>>> assert qml.labs.dla.check_cartan_decomp(k1, m1)
>>> print(f"First iteration, class B: {len(k1)}, {len(m1)}")
First iteration, class B: 15, 15
Now we arrived at the subalgebra :math:`\mathfrak{k}_1=\mathfrak{su}(4)` and can perform the
next iteration of the recursive decomposition.
>>> k2, m2 = cartan_decomp(k1, partial(khaneja_glaser_involution, wire=1))
>>> assert qml.labs.dla.check_cartan_decomp(k2, m2)
>>> print(f"Second iteration, AIII: {len(k2)}, {len(m2)}")
Second iteration, AIII: 7, 8
We can follow up on this by removing the center from :math:`\mathfrak{k}_1` again, and
applying a final class B decomposition:
>>> center_k2 = qml.center(k2, pauli=True)
>>> k2_semi = [op for op in k2 if op not in center_k2]
>>> print(f"Removed operators {center_k2}, new length: {len(k2_semi)}")
Removed operators [1.0 * Z(1)], new length: 6
>>> k3, m3 = cartan_decomp(k2_semi, partial(ClassB, wire=1))
>>> assert qml.labs.dla.check_cartan_decomp(k3, m3)
>>> print(f"Second iteration, class B: {len(k3)}, {len(m3)}")
Second iteration, class B: 3, 3
This concludes our example of the Khaneja-Glaser recursive decomposition, as we arrived at
easy to implement single-qubit operations.
"""
if wire is None:
raise ValueError(
"Please specify the wire for the Khaneja-Glaser involution via "
"functools.partial(khaneja_glaser_involution, wire=wire)"
)
if isinstance(op, np.ndarray):
p = q = op.shape[0] // 2
elif isinstance(op, (PauliSentence, Operator)):
p = q = 2 ** (len(op.wires) - 1)
else:
raise ValueError(f"Can't infer p and q from operator of type {type(op)}.")
return AIII(op, p, q, wire)
# Canonical involutions
# see https://arxiv.org/abs/2406.04418 appendix C
# matrices
def int_log2(x):
"""Compute the integer closest to log_2(x)."""
return int(np.round(np.log2(x)))
def is_qubit_case(p, q):
"""Return whether p and q are the same and a power of 2."""
return p == q and 2 ** int_log2(p) == p
def J(n, wire=None):
"""This is the standard choice for the symplectic transformation operator.
For an :math:`N`-qubit system (:math:`n=2^N`), it equals :math:`Y_0`."""
N = int_log2(n)
if 2**N == n:
if wire is None:
wire = 0
return Y(wire).matrix(wire_order=range(N + 1))
if wire is not None:
raise ValueError("The wire argument is only supported for n=2**N for some integer N.")
zeros = np.zeros((n, n))
eye = np.eye(n)
return np.block([[zeros, -1j * eye], [1j * eye, zeros]])
def Ipq(p, q, wire=None):
"""This is the canonical transformation operator for AIII and BDI Cartan
decompositions. For an :math:`N`-qubit system (:math:`n=2^N`) and
:math:`p=q=n/2`, it equals :math:`Z_0`."""
# If p = q and they are a power of two, use Pauli representation
if is_qubit_case(p, q):
if wire is None:
wire = 0
return Z(wire).matrix(wire_order=range(int_log2(p) + 1))
if wire is not None:
raise ValueError("The wire argument is only supported for p=q=2**N for some integer N.")
return np.block([[np.eye(p), np.zeros((p, q))], [np.zeros((q, p)), -np.eye(q)]])
def Kpq(p, q, wire=None):
"""This is the canonical transformation operator for CII Cartan
decompositions. For an :math:`N`-qubit system (:math:`n=2^N`) and
:math:`p=q=n/2`, it equals :math:`Z_1`."""
# If p = q and they are a power of two, use Pauli representation
if is_qubit_case(p, q):
if wire is None:
wire = 1
return Z(wire).matrix(wire_order=range(int_log2(p) + 1))
if wire is not None:
raise ValueError("The wire argument is only supported for p=q=2**N for some integer N.")
zeros = np.zeros((p + q, p + q))
d = np.diag(np.concatenate([np.ones(p), -np.ones(q)]))
KKm = np.block([[d, zeros], [zeros, d]])
return KKm
[docs]def AI(op: Union[np.ndarray, PauliSentence, Operator]) -> bool:
r"""Canonical Cartan decomposition of type AI, given by :math:`\theta: x \mapsto x^\ast`.
.. note:: Note that we work with Hermitian
operators internally, so that the input will be multiplied by :math:`i` before
evaluating the involution.
Args:
op (Union[np.ndarray, PauliSentence, Operator]): Operator on which the involution is
evaluated and for which the parity under the involution is returned.
Returns:
bool: Whether or not the input operator (times :math:`i`) is in the eigenspace of the
involution :math:`\theta` with eigenvalue :math:`+1`.
"""
return _AI(op)
@singledispatch
def _AI(op): # pylint:disable=unused-argument
r"""Default implementation of the canonical form of the AI involution
:math:`\theta: x \mapsto x^\ast`.
"""
raise NotImplementedError(f"Involution not implemented for operator {op} of type {type(op)}")
@_AI.register(np.ndarray)
def _AI_matrix(op: np.ndarray) -> bool:
r"""Matrix implementation of the canonical form of the AI involution
:math:`\theta: x \mapsto x^\ast`.
"""
op = op * 1j
return np.allclose(op, op.conj())
@_AI.register(PauliSentence)
def _AI_ps(op: PauliSentence) -> bool:
r"""PauliSentence implementation of the canonical form of the AI involution
:math:`\theta: x \mapsto x^\ast`.
"""
parity = []
for pw in op:
result = sum(el == "Y" for el in pw.values())
parity.append(bool(result % 2))
# only makes sense if parity is the same for all terms
assert all(parity) or not any(parity)
return parity[0]
@_AI.register(Operator)
def _AI_op(op: Operator) -> bool:
r"""Operator implementation of the canonical form of the AI involution
:math:`\theta: x \mapsto x^\ast`.
"""
return _AI_ps(op.pauli_rep)
[docs]def AII(op: Union[np.ndarray, PauliSentence, Operator], wire: Optional[int] = None) -> bool:
r"""Canonical Cartan decomposition of type AII, given by :math:`\theta: x \mapsto Y_0 x^\ast Y_0`.
.. note:: Note that we work with Hermitian
operators internally, so that the input will be multiplied by :math:`i` before
evaluating the involution.
Args:
op (Union[np.ndarray, PauliSentence, Operator]): Operator on which the involution is
evaluated and for which the parity under the involution is returned.
wire (int): The wire on which the Pauli-:math:`Y` operator acts to implement the
involution. Will default to ``0`` if ``None``.
Returns:
bool: Whether or not the input operator (times :math:`i`) is in the eigenspace of the
involution :math:`\theta` with eigenvalue :math:`+1`.
"""
return _AII(op, wire)
@singledispatch
def _AII(op, wire=None): # pylint:disable=unused-argument
r"""Default implementation of the canonical form of the AII involution
:math:`\theta: x \mapsto Y_0 x^\ast Y_0`.
"""
raise NotImplementedError(f"Involution not implemented for operator {op} of type {type(op)}")
@_AII.register(np.ndarray)
def _AII_matrix(op: np.ndarray, wire: Optional[int] = None) -> bool:
r"""Matrix implementation of the canonical form of the AII involution
:math:`\theta: x \mapsto Y_0 x^\ast Y_0`.
"""
op = op * 1j
y = J(op.shape[-1] // 2, wire=wire)
return np.allclose(op, y @ op.conj() @ y)
@_AII.register(PauliSentence)
def _AII_ps(op: PauliSentence, wire: Optional[int] = None) -> bool:
r"""PauliSentence implementation of the canonical form of the AII involution
:math:`\theta: x \mapsto Y_0 x^\ast Y_0`.
"""
if wire is None:
wire = 0
parity = []
for pw in op:
result = sum(el == "Y" for el in pw.values()) + (pw.get(wire, "I") in "XZ")
parity.append(bool(result % 2))
# only makes sense if parity is the same for all terms
assert all(parity) or not any(parity)
return parity[0]
@_AII.register(Operator)
def _AII_op(op: Operator, wire: Optional[int] = None) -> bool:
r"""Operator implementation of the canonical form of the AII involution
:math:`\theta: x \mapsto Y_0 x^\ast Y_0`.
"""
return _AII_ps(op.pauli_rep, wire)
[docs]def AIII(
op: Union[np.ndarray, PauliSentence, Operator],
p: int = None,
q: int = None,
wire: Optional[int] = None,
) -> bool:
r"""Canonical Cartan decomposition of type AIII, given by :math:`\theta: x \mapsto I_{p,q} x I_{p,q}`.
The matrix :math:`I_{p,q}` is given by
.. math::
I_{p,q}=\text{diag}(\underset{p \text{times}}{\underbrace{1, \dots 1}},
\underset{q \text{times}}{\underbrace{-1, \dots -1}}).
For :math:`p=q=2^N` for some integer :math:`N`, we have :math:`I_{p,q}=Z_0`.
.. note:: Note that we work with Hermitian operators internally, so that the input will be
multiplied by :math:`i` before evaluating the involution.
Args:
op (Union[np.ndarray, PauliSentence, Operator]): Operator on which the involution is
evaluated and for which the parity under the involution is returned.
p (int): Dimension of first subspace.
q (int): Dimension of second subspace.
wire (int): The wire on which the Pauli-:math:`Z` operator acts to implement the
involution. Will default to ``0`` if ``None``.
Returns:
bool: Whether or not the input operator (times :math:`i`) is in the eigenspace of the
involution :math:`\theta` with eigenvalue :math:`+1`.
"""
if p is None or q is None:
raise ValueError(
"please specify p and q for the involution via functools.partial(AIII, p=p, q=q)"
)
return _AIII(op, p, q, wire)
@singledispatch
def _AIII(op, p=None, q=None, wire=None): # pylint:disable=unused-argument
r"""Default implementation of the canonical form of the AIII involution
:math:`\theta: x \mapsto I_{p,q} x I_{p,q}`.
"""
raise NotImplementedError(f"Involution not implemented for operator {op} of type {type(op)}")
@_AIII.register(np.ndarray)
def _AIII_matrix(op: np.ndarray, p: int = None, q: int = None, wire: Optional[int] = None) -> bool:
r"""Matrix implementation of the canonical form of the AIII involution
:math:`\theta: x \mapsto I_{p,q} x I_{p,q}`.
"""
op = op * 1j
z = Ipq(p, q, wire=wire)
return np.allclose(op, z @ op @ z)
@_AIII.register(PauliSentence)
def _AIII_ps(op: PauliSentence, p: int = None, q: int = None, wire: Optional[int] = None) -> bool:
r"""PauliSentence implementation of the canonical form of the AIII involution
:math:`\theta: x \mapsto I_{p,q} x I_{p,q}`.
"""
if is_qubit_case(p, q):
if wire is None:
wire = 0
parity = [pw.get(wire, "I") in "IZ" for pw in op]
# only makes sense if parity is the same for all terms
assert all(parity) or not any(parity)
return parity[0]
# If it is not a qubit case, use the matrix representation
return _AIII_matrix(matrix(op, wire_order=sorted(op.wires)), p, q, wire)
@_AIII.register(Operator)
def _AIII_op(op: Operator, p: int = None, q: int = None, wire: Optional[int] = None) -> bool:
r"""Operator implementation of the canonical form of the AIII involution
:math:`\theta: x \mapsto I_{p,q} x I_{p,q}`.
"""
return _AIII_ps(op.pauli_rep, p, q, wire)
[docs]def BDI(
op: Union[np.ndarray, PauliSentence, Operator],
p: int = None,
q: int = None,
wire: Optional[int] = None,
) -> bool:
r"""Canonical Cartan decomposition of type BDI, given by :math:`\theta: x \mapsto I_{p,q} x I_{p,q}`.
The matrix :math:`I_{p,q}` is given by
.. math::
I_{p,q}=\text{diag}(\underset{p \text{times}}{\underbrace{1, \dots 1}},
\underset{q \text{times}}{\underbrace{-1, \dots -1}}).
For :math:`p=q=2^N` for some integer :math:`N`, we have :math:`I_{p,q}=Z_0`.
.. note:: Note that we work with Hermitian operators internally, so that the input will be
multiplied by :math:`i` before evaluating the involution.
Args:
op (Union[np.ndarray, PauliSentence, Operator]): Operator on which the involution is
evaluated and for which the parity under the involution is returned.
p (int): Dimension of first subspace.
q (int): Dimension of second subspace.
wire (int): The wire on which the Pauli-:math:`Z` operator acts to implement the
involution. Will default to ``0`` if ``None``.
Returns:
bool: Whether or not the input operator (times :math:`i`) is in the eigenspace of the
involution :math:`\theta` with eigenvalue :math:`+1`.
"""
return AIII(op, p, q, wire)
[docs]def CI(op: Union[np.ndarray, PauliSentence, Operator]) -> bool:
r"""Canonical Cartan decomposition of type CI, given by :math:`\theta: x \mapsto x^\ast`.
.. note:: Note that we work with Hermitian
operators internally, so that the input will be multiplied by :math:`i` before
evaluating the involution.
Args:
op (Union[np.ndarray, PauliSentence, Operator]): Operator on which the involution is
evaluated and for which the parity under the involution is returned.
Returns:
bool: Whether or not the input operator (times :math:`i`) is in the eigenspace of the
involution :math:`\theta` with eigenvalue :math:`+1`.
"""
return AI(op)
[docs]def CII(
op: Union[np.ndarray, PauliSentence, Operator],
p: int = None,
q: int = None,
wire: Optional[int] = None,
) -> bool:
r"""Canonical Cartan decomposition of type CII, given by :math:`\theta: x \mapsto K_{p,q} x K_{p,q}`.
The matrix :math:`K_{p,q}` is given by
.. math::
K_{p,q}=\text{diag}(
\underset{p \text{times}}{\underbrace{1, \dots 1}},
\underset{q \text{times}}{\underbrace{-1, \dots -1}},
\underset{p \text{times}}{\underbrace{1, \dots 1}},
\underset{q \text{times}}{\underbrace{-1, \dots -1}},
).
For :math:`p=q=2^N` for some integer :math:`N`, we have :math:`K_{p,q}=Z_1`.
.. note:: Note that we work with Hermitian operators internally, so that the input will be
multiplied by :math:`i` before evaluating the involution.
Args:
op (Union[np.ndarray, PauliSentence, Operator]): Operator on which the involution is
evaluated and for which the parity under the involution is returned.
p (int): Dimension of first subspace.
q (int): Dimension of second subspace.
wire (int): The wire on which the Pauli-:math:`Z` operator acts to implement the
involution. Will default to ``1`` if ``None``.
Returns:
bool: Whether or not the input operator (times :math:`i`) is in the eigenspace of the
involution :math:`\theta` with eigenvalue :math:`+1`.
"""
if p is None or q is None:
raise ValueError(
"please specify p and q for the involution via functools.partial(CII, p=p, q=q)"
)
return _CII(op, p, q, wire)
@singledispatch
def _CII(op, p=None, q=None, wire=None): # pylint:disable=unused-argument
r"""Default implementation of the canonical form of the CII involution
:math:`\theta: x \mapsto K_{p,q} x K_{p,q}`.
"""
raise NotImplementedError(f"Involution not implemented for operator {op} of type {type(op)}")
@_CII.register(np.ndarray)
def _CII_matrix(op: np.ndarray, p: int = None, q: int = None, wire: Optional[int] = None) -> bool:
r"""Matrix implementation of the canonical form of the CII involution
:math:`\theta: x \mapsto K_{p,q} x K_{p,q}`.
"""
op = op * 1j
z = Kpq(p, q, wire=wire)
return np.allclose(op, z @ op @ z)
@_CII.register(PauliSentence)
def _CII_ps(op: PauliSentence, p: int = None, q: int = None, wire: Optional[int] = None) -> bool:
r"""PauliSentence implementation of the canonical form of the CII involution
:math:`\theta: x \mapsto K_{p,q} x K_{p,q}`.
"""
if is_qubit_case(p, q):
if wire is None:
wire = 1
parity = [pw.get(wire, "I") in "IZ" for pw in op]
# only makes sense if parity is the same for all terms
assert all(parity) or not any(parity)
return parity[0]
# If it is not a qubit case, use the matrix representation
return _CII(matrix(op, wire_order=sorted(op.wires)), p, q, wire)
@_CII.register(Operator)
def _CII_op(op: Operator, p: int = None, q: int = None, wire: Optional[int] = None) -> bool:
r"""Operator implementation of the canonical form of the CII involution
:math:`\theta: x \mapsto K_{p,q} x K_{p,q}`.
"""
return _CII_ps(op.pauli_rep, p, q, wire)
[docs]def DIII(op: Union[np.ndarray, PauliSentence, Operator], wire: Optional[int] = None) -> bool:
r"""Canonical Cartan decomposition of type DIII, given by :math:`\theta: x \mapsto Y_0 x Y_0`.
Args:
op (Union[np.ndarray, PauliSentence, Operator]): Operator on which the involution is
evaluated and for which the parity under the involution is returned.
wire (int): The wire on which the Pauli-:math:`Y` operator acts to implement the
involution. Will default to ``0`` if ``None``.
Returns:
bool: Whether or not the input operator (times :math:`i`) is in the eigenspace of the
involution :math:`\theta` with eigenvalue :math:`+1`.
"""
return _DIII(op, wire)
@singledispatch
def _DIII(op, wire=None): # pylint:disable=unused-argument
r"""Default implementation of the canonical form of the DIII involution
:math:`\theta: x \mapsto Y_0 x Y_0`.
"""
raise NotImplementedError(f"Involution not implemented for operator {op} of type {type(op)}")
@_DIII.register(np.ndarray)
def _DIII_matrix(op: np.ndarray, wire: Optional[int] = None) -> bool:
r"""Matrix implementation of the canonical form of the DIII involution
:math:`\theta: x \mapsto Y_0 x Y_0`.
"""
y = J(op.shape[-1] // 2, wire=wire)
return np.allclose(op, y @ op @ y)
@_DIII.register(PauliSentence)
def _DIII_ps(op: PauliSentence, wire: Optional[int] = None) -> bool:
r"""PauliSentence implementation of the canonical form of the DIII involution
:math:`\theta: x \mapsto Y_0 x Y_0`.
"""
if wire is None:
wire = 0
parity = [pw.get(wire, "I") in "IY" for pw in op]
# only makes sense if parity is the same for all terms
assert all(parity) or not any(parity)
return parity[0]
@_DIII.register(Operator)
def _DIII_op(op: Operator, wire: Optional[int] = None) -> bool:
r"""Operator implementation of the canonical form of the DIII involution
:math:`\theta: x \mapsto Y_0 x Y_0`.
"""
return _DIII_ps(op.pauli_rep, wire)
[docs]def ClassB(op: Union[np.ndarray, PauliSentence, Operator], wire: Optional[int] = None) -> bool:
r"""Canonical Cartan decomposition of class B, given by :math:`\theta: x \mapsto Y_0 x Y_0`.
Args:
op (Union[np.ndarray, PauliSentence, Operator]): Operator on which the involution is
evaluated and for which the parity under the involution is returned.
wire (int): The wire on which the Pauli-:math:`Y` operator acts to implement the
involution. Will default to ``0`` if ``None``.
Returns:
bool: Whether or not the input operator (times :math:`i`) is in the eigenspace of the
involution :math:`\theta` with eigenvalue :math:`+1`.
"""
return DIII(op, wire)
# dispatch to different input types
[docs]def even_odd_involution(op: Union[PauliSentence, np.ndarray, Operator]) -> bool:
r"""The Even-Odd involution.
This is defined in `quant-ph/0701193 <https://arxiv.org/abs/quant-ph/0701193>`__.
For Pauli words and sentences, it comes down to counting non-trivial Paulis in Pauli words.
Args:
op ( Union[PauliSentence, np.ndarray, Operator]): Input operator
Returns:
bool: Boolean output ``True`` or ``False`` for odd (:math:`\mathfrak{k}`) and even parity subspace (:math:`\mathfrak{m}`), respectively
.. seealso:: :func:`~cartan_decomp`
**Example**
>>> from pennylane import X, Y, Z
>>> from pennylane.labs.dla import even_odd_involution
>>> ops = [X(0), X(0) @ Y(1), X(0) @ Y(1) @ Z(2)]
>>> [even_odd_involution(op) for op in ops]
[True, False, True]
Operators with an odd number of non-identity Paulis yield ``1``, whereas even ones yield ``0``.
The function also works with dense matrix representations.
>>> ops_m = [qml.matrix(op, wire_order=range(3)) for op in ops]
>>> [even_odd_involution(op_m) for op_m in ops_m]
[True, False, True]
"""
return _even_odd_involution(op)
@singledispatch
def _even_odd_involution(op): # pylint:disable=unused-argument, missing-function-docstring
return NotImplementedError(f"Involution not defined for operator {op} of type {type(op)}")
@_even_odd_involution.register(PauliSentence)
def _even_odd_involution_ps(op: PauliSentence):
# Generalization to sums of Paulis: check each term and assert they all have the same parity
parity = []
for pw in op.keys():
parity.append(len(pw) % 2)
# only makes sense if parity is the same for all terms, e.g. Heisenberg model
assert all(
parity[0] == p for p in parity
), f"The Even-Odd involution is not well-defined for operator {op} as individual terms have different parity"
return bool(parity[0])
@_even_odd_involution.register(np.ndarray)
def _even_odd_involution_matrix(op: np.ndarray):
"""see Table CI in https://arxiv.org/abs/2406.04418"""
n = int(np.round(np.log2(op.shape[-1])))
YYY = qml.prod(*[Y(i) for i in range(n)])
YYY = qml.matrix(YYY, range(n))
transformed = YYY @ op.conj() @ YYY
if np.allclose(transformed, op):
return False
if np.allclose(transformed, -op):
return True
raise ValueError(f"The Even-Odd involution is not well-defined for operator {op}.")
@_even_odd_involution.register(Operator)
def _even_odd_involution_op(op: Operator):
"""use pauli representation"""
return _even_odd_involution_ps(op.pauli_rep)
# dispatch to different input types
[docs]def concurrence_involution(op: Union[PauliSentence, np.ndarray, Operator]) -> bool:
r"""The Concurrence Canonical Decomposition :math:`\Theta(g) = -g^T` as a Cartan involution function.
This is defined in `quant-ph/0701193 <https://arxiv.org/abs/quant-ph/0701193>`__.
For Pauli words and sentences, it comes down to counting Pauli-Y operators.
Args:
op ( Union[PauliSentence, np.ndarray, Operator]): Input operator
Returns:
bool: Boolean output ``True`` or ``False`` for odd (:math:`\mathfrak{k}`) and even parity subspace (:math:`\mathfrak{m}`), respectively
.. seealso:: :func:`~cartan_decomp`
**Example**
>>> from pennylane import X, Y, Z
>>> from pennylane.labs.dla import concurrence_involution
>>> ops = [X(0), X(0) @ Y(1), X(0) @ Y(1) @ Z(2), Y(0) @ Y(2)]
>>> [concurrence_involution(op) for op in ops]
[False, True, True, False]
Operators with an odd number of ``Y`` operators yield ``1``, whereas even ones yield ``0``.
The function also works with dense matrix representations.
>>> ops_m = [qml.matrix(op, wire_order=range(3)) for op in ops]
>>> [even_odd_involution(op_m) for op_m in ops_m]
[False, True, True, False]
"""
return _concurrence_involution(op)
@singledispatch
def _concurrence_involution(op):
return NotImplementedError(f"Involution not defined for operator {op} of type {type(op)}")
@_concurrence_involution.register(PauliSentence)
def _concurrence_involution_pauli(op: PauliSentence):
# Generalization to sums of Paulis: check each term and assert they all have the same parity
parity = []
for pw in op.keys():
result = sum(1 if el == "Y" else 0 for el in pw.values())
parity.append(result % 2)
# only makes sense if parity is the same for all terms, e.g. Heisenberg model
assert all(
parity[0] == p for p in parity
), f"The concurrence canonical decomposition is not well-defined for operator {op} as individual terms have different parity"
return bool(parity[0])
@_concurrence_involution.register(Operator)
def _concurrence_involution_operator(op: Operator):
return _concurrence_involution_matrix(op.matrix())
@_concurrence_involution.register(np.ndarray)
def _concurrence_involution_matrix(op: np.ndarray):
if np.allclose(op, -op.T):
return True
if np.allclose(op, op.T):
return False
raise ValueError(
f"The concurrence canonical decomposition is not well-defined for operator {op}"
)
_modules/pennylane/labs/dla/involutions
Download Python script
Download Notebook
View on GitHub