[labs.dla 4] Cartan decomposition

pennylane/labs/dla/__init__.py

@@ -22,6 +22,7 @@
 
     ~lie_closure_dense
     ~structure_constants_dense
+    ~cartan_decomposition
 
 
 Utility functions
@@ -39,11 +40,23 @@
     ~check_orthonormal
     ~trace_inner_product
 
+Involutions
+~~~~~~~~~~~
+
+.. currentmodule:: pennylane.labs.dla
+
+.. autosummary::
+    :toctree: api
+
+    ~even_odd_involution
+    ~concurrence_involution
+
 
 """
 
 from .lie_closure_dense import lie_closure_dense
 from .structure_constants_dense import structure_constants_dense
+from .cartan import cartan_decomposition, even_odd_involution, concurrence_involution
 from .dense_util import (
     pauli_coefficients,
     pauli_decompose,

pennylane/labs/dla/cartan.py

@@ -0,0 +1,165 @@
+# 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.
+"""Functionality for Cartan decomposition"""
+from functools import singledispatch
+from typing import Union
+
+import numpy as np
+
+import pennylane as qml
+from pennylane import Y
+from pennylane.operation import Operator
+from pennylane.pauli import PauliSentence
+
+
+def cartan_decomposition(g, involution):
+    r"""Cartan Decomposition :math:`\mathfrak{g} = \mathfrak{k} \plus \mathfrak{m}`.
+
+    Given a Lie algebra :math:`\mathfrak{g}`, the Cartan decomposition is a decomposition
+    :math:`\mathfrak{g} = \mathfrak{k} \oplus \mathfrak{m}` into orthogonal complements.
+    This is realized by an involution :math:`\Theta(g)` that maps each operator :math:`g \in \mathfrak{g}`
+    back to itself after two consecutive applications, i.e., :math:`\Theta(\Theta(g)) = g \forall g \in \mathfrak{g}`.
+
+    The ``involution`` argument can be any function that maps the operators in the provided ``g`` to a boolean output.
+    ``True`` for operators that go into :math:`\mathfrak{k}` and ``False`` for operators in :math:`\mathfrak{m}`.
+
+    The resulting subspaces fulfill the Cartan commutation relations
+
+    .. math:: [\mathfrak{k}, \mathfrak{k}] \subseteq \mathfrak{k} \text{ ; } [\mathfrak{k}, \mathfrak{m}] \subseteq \mathfrak{m} \text{ ; } [\mathfrak{m}, \mathfrak{m}] \subseteq \mathfrak{k}
+
+    Args:
+        g (List[Union[PauliSentence, Operator]]): the (dynamical) Lie algebra to decompose
+        involution (callable): Involution function :math:`\Theta(\cdot)` to act on the input operator, should return ``0/1`` or ``True/False``.
+            E.g., :func:`~even_odd_involution` or :func:`~concurrence_involution`.
+
+    Returns:
+        k (List[Union[PauliSentence, Operator]]): the even parity subspace :math:`\Theta(\mathfrak{k}) = \mathfrak{k}`
+        m (List[Union[PauliSentence, Operator]]): the odd parity subspace :math:`\Theta(\mathfrak{m}) = \mathfrak{m}`
+
+    .. seealso:: :func:`~even_odd_involution`, :func:`~concurrence_involution`
+    """
+    # simple implementation assuming all elements in g are already either in k and m
+    # TODO: Figure out more general way to do this when the above is not the case
+    m = []
+    k = []
+
+    for op in g:
+        if involution(op):  # odd parity
+            k.append(op)
+        else:  # even parity
+            m.append(op)
+
+    return k, m
+
+
+# dispatch to different input types
+def even_odd_involution(op: Union[PauliSentence, np.ndarray, Operator]):
+    r"""The Even-Odd involution
+
+    This is defined in `quant-ph/0701193 <https://arxiv.org/pdf/quant-ph/0701193>`__, and for Pauli words and sentences 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_decomposition`
+    """
+    return _even_odd_involution(op)
+
+
+@singledispatch
+def _even_odd_involution(op):  # pylint:disable=unused-argument
+    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 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
+    return not np.allclose(transformed, 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
+def concurrence_involution(op: Union[PauliSentence, np.ndarray, Operator]):
+    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/pdf/quant-ph/0701193>`__, and for Pauli words and sentences 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_decomposition`
+
+    """
+    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_operation(op: Operator):
+    op = op.matrix()
+    return np.allclose(op, -op.T)
+
+
+@_concurrence_involution.register(np.ndarray)
+def _concurrence_involution_matrix(op: np.ndarray):
+    return np.allclose(op, -op.T)

pennylane/labs/tests/dla/test_cartan.py

@@ -0,0 +1,104 @@
+# 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.
+"""Tests for pennylane/dla/lie_closure_dense.py functionality"""
+# pylint: disable=no-self-use,too-few-public-methods,missing-class-docstring
+import pytest
+
+import pennylane as qml
+from pennylane import X, Y, Z
+from pennylane.labs.dla import cartan_decomposition, concurrence_involution, even_odd_involution
+
+
+def check_commutation(ops1, ops2, vspace):
+    """Helper function to check things like [k, m] subspace m; expensive"""
+    for o1 in ops1:
+        for o2 in ops2:
+            com = o1.commutator(o2)
+            assert not vspace.is_independent(com)
+
+    return True
+
+
+Ising2 = qml.lie_closure([X(0), X(1), Z(0) @ Z(1)])
+Ising3 = qml.lie_closure([X(0), X(1), X(2), Z(0) @ Z(1), Z(1) @ Z(2)])
+Heisenberg3 = qml.lie_closure(
+    [X(0) @ X(1), X(1) @ X(2), Y(0) @ Y(1), Y(1) @ Y(2), Z(0) @ Z(1), Z(1) @ Z(2)]
+)
+
+
+class TestCartanDecomposition:
+    @pytest.mark.parametrize("involution", [even_odd_involution, concurrence_involution])
+    @pytest.mark.parametrize("g", [Ising2, Ising3, Heisenberg3])
+    def test_cartan_decomposition(self, g, involution):
+        """Test basic properties and Cartan decomposition definitions"""
+
+        g = [op.pauli_rep for op in g]
+        k, m = cartan_decomposition(g, involution)
+
+        assert all(involution(op) == 1 for op in k)
+        assert all(involution(op) == 0 for op in m)
+
+        k_space = qml.pauli.PauliVSpace(k)
+        m_space = qml.pauli.PauliVSpace(m)
+
+        # Commutation relations for Cartan pair
+        assert check_commutation(k, k, k_space)
+        assert check_commutation(k, m, m_space)
+        assert check_commutation(m, m, k_space)
+
+    @pytest.mark.parametrize("involution", [even_odd_involution, concurrence_involution])
+    @pytest.mark.parametrize("g", [Ising2, Ising3, Heisenberg3])
+    def test_cartan_decomposition_dense(self, g, involution):
+        """Test basic properties and Cartan decomposition definitions using dense representations"""
+
+        g = [op.pauli_rep for op in g]
+        k, m = cartan_decomposition(g, involution)
+
+        assert all(involution(op) == 1 for op in k)
+        assert all(involution(op) == 0 for op in m)
+
+        k_space = qml.pauli.PauliVSpace(k)
+        m_space = qml.pauli.PauliVSpace(m)
+
+        # Commutation relations for Cartan pair
+        assert check_commutation(k, k, k_space)
+        assert check_commutation(k, m, m_space)
+        assert check_commutation(m, m, k_space)
+
+
+involution_ops = [X(0) @ X(1), X(0) @ X(1) + Y(0) @ Y(1)]
+
+
+class TestInvolutions:
+    """Test involutions"""
+
+    @pytest.mark.parametrize("op", involution_ops)
+    def test_concurrence_involution_inputs(self, op):
+        """Test different input types yield consistent results"""
+        res_op = concurrence_involution(op)
+        res_ps = concurrence_involution(op.pauli_rep)
+        res_m = concurrence_involution(op.matrix())
+
+        assert res_op == res_ps
+        assert res_op == res_m
+
+    @pytest.mark.parametrize("op", involution_ops)
+    def test_even_odd_involution_inputs(self, op):
+        """Test different input types yield consistent results"""
+        res_op = even_odd_involution(op)
+        res_ps = even_odd_involution(op.pauli_rep)
+        res_m = even_odd_involution(op.matrix())
+
+        assert res_op == res_ps
+        assert res_op == res_m