Merge pull request #1345 from tnemoz/fix_1088

Added Rényi entropy of a quantum state

Author: vprusso

docs/refs.bib

@@ -987,6 +987,19 @@ @misc{Molina_2012_Optimal
   primaryClass={quant-ph}
 }
 
+@article{Muller_2013_Renyi_Generalization,
+   title={On quantum Rényi entropies: A new generalization and some properties},
+   volume={54},
+   ISSN={1089-7658},
+   url={http://dx.doi.org/10.1063/1.4838856},
+   DOI={10.1063/1.4838856},
+   number={12},
+   journal={Journal of Mathematical Physics},
+   publisher={AIP Publishing},
+   author={Müller-Lennert, Martin and Dupuis, Frédéric and Szehr, Oleg and Fehr, Serge and Tomamichel, Marco},
+   year={2013},
+   month=dec }
+
 % Last name begins with N
 @article{Navascues_2008_AConvergent,
   title={A convergent hierarchy of semidefinite programs characterizing the set of quantum correlations},

toqito/state_props/__init__.py

@@ -26,3 +26,4 @@
 from toqito.state_props.is_distinguishable import is_distinguishable
 from toqito.state_props.is_unextendible_product_basis import is_unextendible_product_basis
 from toqito.state_props.common_quantum_overlap import common_quantum_overlap
+from toqito.state_props.renyi_entropy import renyi_entropy

toqito/state_props/renyi_entropy.py

@@ -0,0 +1,107 @@
+"""Calculates the Rényi entropy metric of a quantum state."""
+
+import numpy as np
+
+from toqito.matrix_props import is_density
+from toqito.state_props import von_neumann_entropy
+
+
+def renyi_entropy(rho: np.ndarray, alpha: float) -> float:
+    r"""Compute the Rényi entropy of a density matrix :footcite:`Muller_2013_Renyi_Generalization`.
+
+    Let :math:`P \in \text{Pos}(\mathcal{X})` be a positive semidefinite operator, for a complex
+    Euclidean space :math:`\mathcal{X}`. Then one defines the *Rényi entropy of order*
+    :math:`\alpha\geqslant0` as
+
+    .. math::
+        H_{\alpha}(P) = H_{\alpha}(\lambda(P)),
+
+    where :math:`\lambda(P)` is the vector of eigenvalues of :math:`P` and where the function
+    :math:`H(\cdot)` is the classical Rényi entropy of order :math:`\alpha` defined as
+
+    .. math::
+        H_{\alpha}(u) = \frac{1}{1-\alpha}\log\left(\sum_{a \in \Sigma} u(a)^{\alpha}\right),
+
+    where the :math:`\log` function is assumed to be the base-2 logarithm, and where
+    :math:`\Sigma` is an alphabet where :math:`u \in [0, \infty)^{\Sigma}` is a vector of
+    nonnegative real numbers indexed by :math:`\Sigma`. It recovers the von Neumann entropy for
+    :math:`\alpha=1` and the min-entropy for :math:`\alpha=+\infty`.
+
+    Examples
+    ==========
+
+    Consider the following Bell state:
+
+    .. math::
+        u = \frac{1}{\sqrt{2}} \left(|00 \rangle + |11 \rangle \right) \in \mathcal{X}.
+
+    The corresponding density matrix of :math:`u` may be calculated by:
+
+    .. math::
+        \rho = u u^* = \frac{1}{2} \begin{pmatrix}
+                         1 & 0 & 0 & 1 \\
+                         0 & 0 & 0 & 0 \\
+                         0 & 0 & 0 & 0 \\
+                         1 & 0 & 0 & 1
+                       \end{pmatrix} \in \text{D}(\mathcal{X}).
+
+    Calculating the Rényi entropy of order :math:`2` of :math:`\rho` in :code:`|toqito⟩` can be
+    done as follows.
+
+    .. jupyter-execute::
+
+        from toqito.state_props import renyi_entropy
+        import numpy as np
+        test_input_mat = np.array(
+                [[1 / 2, 0, 0, 1 / 2], [0, 0, 0, 0],
+                [0, 0, 0, 0], [1 / 2, 0, 0, 1 / 2]]
+            )
+        renyi_entropy(test_input_mat, 2)
+
+    Consider the density operator corresponding to the maximally mixed state of dimension two
+
+    .. math::
+        \rho = \frac{1}{2}
+        \begin{pmatrix}
+            1 & 0 \\
+            0 & 1
+        \end{pmatrix}.
+
+    As this state is maximally mixed, the Rényi entropy of :math:`\rho` is
+    equal to one for all orders :math:`\alpha`. We can see this in :code:`|toqito⟩` as follows.
+
+    .. jupyter-execute::
+
+        from toqito.state_props import renyi_entropy
+        import numpy as np
+        rho = 1/2 * np.identity(2)
+        renyi_entropy(rho, 3/2)
+
+    References
+    ==========
+    .. footbibliography::
+
+
+
+    :param rho: Density operator.
+    :param alpha: Order for the Rényi entropy. Note that numerical instability may happen for small
+        positive values because of the computation of the spectral decomposition.
+    :return: The Rényi entropy of order :code:`alpha` of :code:`rho`.
+
+    """
+    if not is_density(rho):
+        raise ValueError("Rényi entropy is only defined for density operators.")
+    if alpha < 0:
+        raise ValueError("Rényi entropy is only defined for positive orders.")
+    if alpha == 0:
+        return np.log2(np.linalg.matrix_rank(rho))
+    if alpha == 1:
+        return von_neumann_entropy(rho)
+
+    eigs = np.linalg.eigvalsh(rho)
+    eigs = eigs[eigs > 0]
+
+    if alpha == float("inf"):
+        return -np.log2(eigs.max())
+
+    return np.log2(pow(eigs, alpha).sum()) / (1 - alpha)

toqito/state_props/tests/test_renyi_entropy.py

@@ -0,0 +1,46 @@
+"""Tests for renyi_entropy."""
+
+import numpy as np
+import pytest
+
+from toqito.state_props import purity, renyi_entropy, von_neumann_entropy
+from toqito.states import bell, max_mixed
+
+RHO_TEST = np.array([[0.8, 0.0], [0.0, 0.2]])
+
+
+@pytest.mark.parametrize(
+    "rho, alpha, expected_result",
+    [
+        # Entangled state Rényi entropy should be zero.
+        (bell(0) @ bell(0).conj().T, 3 / 2, 0),
+        # Rényi entropy of the maximally mixed state should be one.
+        (max_mixed(2, is_sparse=False), 5 / 2, 1),
+        # For alpha=0, log of the number of outcomes
+        (RHO_TEST, 0.0, 1.0),
+        # For alpha=1, recovers von Neumann entropy
+        (RHO_TEST, 1.0, von_neumann_entropy(RHO_TEST)),
+        # For alpha=2, collision entropy
+        (RHO_TEST, 2.0, -np.log2(purity(RHO_TEST))),
+        # For alpha=+inf, min-entropy
+        (RHO_TEST, float("inf"), -np.log2(0.8)),
+    ],
+)
+def test_renyi_entropy(rho, alpha, expected_result):
+    """Test function works as expected for a valid input."""
+    np.testing.assert_allclose(renyi_entropy(rho, alpha), expected_result, atol=1e-5)
+
+
+@pytest.mark.parametrize(
+    "rho, alpha",
+    [
+        # Test Rényi entropy on non-density matrix.
+        (np.array([[1, 2], [3, 4]]), 3 / 2),
+        # Test Rényi entropy on non-positive order.
+        (RHO_TEST, -1.0),
+    ],
+)
+def test_renyi_invalid_input(rho, alpha):
+    """Test function works as expected for an invalid input."""
+    with np.testing.assert_raises(ValueError):
+        renyi_entropy(rho, alpha)