Skip to content

random_povm

Generates a random POVM.

random_povm

random_povm(
    dim: int,
    num_inputs: int,
    num_outputs: int,
    seed: int | None = None,
) -> ndarray

Generate random positive operator valued measurements (POVMs) 1.

Randomness model

For each input we draw \(n_{\text{out}}\) matrices from the real Ginibre ensemble, i.e., each entry is sampled independently from the standard normal distribution using numpy's default_rng. We interpret these matrices as Kraus operators \(A_{x,a}\) and normalize them so that the measurement is complete. Concretely, for each input \(x\) we form

\[ G_x = \sum_a A_{x,a}^\dagger A_{x,a}, \qquad B_{x,a} = G_x^{-1/2} A_{x,a}, \qquad M_{x,a} = B_{x,a}^\dagger B_{x,a}. \]

The matrices \(M_{x,a}\) constitute a POVM satisfying \(\sum_a M_{x,a} = \mathbb{I}\). This procedure induces the (Hilbert–Schmidt) normalized Wishart measure on the POVM effects. Supplying seed reproduces the same sample sequence.

Parameters:

  • dim (int) –

    The dimensions of the measurements.

  • num_inputs (int) –

    The number of inputs for the measurement.

  • num_outputs (int) –

    The number of outputs for the measurement.

  • seed (int | None, default: None ) –

    A seed used to instantiate numpy's random number generator (Ginibre sampling).

Returns:

  • ndarray

    A set of dim-by-dim POVMs of shape (dim, dim, num_inputs, num_outputs).

Examples:

We can generate a set of dim-by-dim POVMs consisting of a specific dimension along with a given number of measurement inputs and measurement outputs. As an example, we can construct a random set of \(2\)-by-\(2\) POVMs of dimension with \(2\) inputs and \(2\) outputs.

import numpy as np
from toqito.rand import random_povm

dim, num_inputs, num_outputs = 2, 2, 2

povms = random_povm(dim, num_inputs, num_outputs)

print(povms)
[[[[ 0.86242361+0.j  0.13757639+0.j]
   [ 0.77505005+0.j  0.22494995+0.j]]

  [[-0.16453336+0.j  0.16453336+0.j]
   [-0.3192279 +0.j  0.3192279 +0.j]]]


 [[[-0.16453336+0.j  0.16453336+0.j]
   [-0.3192279 +0.j  0.3192279 +0.j]]

  [[ 0.10289044+0.j  0.89710956+0.j]
   [ 0.40225942+0.j  0.59774058+0.j]]]]

We can verify that this constitutes a valid set of POVM elements as checking that these operators all sum to the identity operator.

print(np.round(povms[:, :, 0, 0] + povms[:, :, 0, 1]))
[[1.+0.j 0.+0.j]
 [0.+0.j 1.+0.j]]

It is also possible to add a seed for reproducibility.

import numpy as np
from toqito.rand import random_povm

dim, num_inputs, num_outputs = 2, 2, 2

povms = random_povm(dim, num_inputs, num_outputs, seed=42)

print(povms)
[[[[ 0.22988028+0.j  0.77011972+0.j]
   [ 0.45021752+0.j  0.54978248+0.j]]

  [[-0.23938341+0.j  0.23938341+0.j]
   [ 0.32542956+0.j -0.32542956+0.j]]]


 [[[-0.23938341+0.j  0.23938341+0.j]
   [ 0.32542956+0.j -0.32542956+0.j]]

  [[ 0.83184406+0.j  0.16815594+0.j]
   [ 0.61323275+0.j  0.38676725+0.j]]]]

We can once again verify that this constitutes a valid set of POVM elements as checking that these operators all sum to the identity operator.

print(np.round(povms[:, :, 0, 0] + povms[:, :, 0, 1]))
[[ 1.+0.j -0.+0.j]
 [-0.+0.j  1.+0.j]]

References

1 Wikipedia. {POVM}. link.

Source code in toqito/rand/random_povm.py
def random_povm(dim: int, num_inputs: int, num_outputs: int, seed: int | None = None) -> np.ndarray:
    r"""Generate random positive operator valued measurements (POVMs) [@wikipediapovm].

    Randomness model
    ----------------

    For each input we draw \(n_{\text{out}}\) matrices from the real Ginibre ensemble, i.e., each
    entry is sampled independently from the standard normal distribution using ``numpy``'s
    ``default_rng``.  We interpret these matrices as Kraus operators \(A_{x,a}\) and normalize
    them so that the measurement is complete.  Concretely, for each input \(x\) we form

    \[
        G_x = \sum_a A_{x,a}^\dagger A_{x,a}, \qquad
        B_{x,a} = G_x^{-1/2} A_{x,a}, \qquad
        M_{x,a} = B_{x,a}^\dagger B_{x,a}.
    \]

    The matrices \(M_{x,a}\) constitute a POVM satisfying
    \(\sum_a M_{x,a} = \mathbb{I}\).  This procedure induces the (Hilbert–Schmidt) normalized
    Wishart measure on the POVM effects.  Supplying ``seed`` reproduces the same sample sequence.

    Args:
        dim: The dimensions of the measurements.
        num_inputs: The number of inputs for the measurement.
        num_outputs: The number of outputs for the measurement.
        seed: A seed used to instantiate numpy's random number generator (Ginibre sampling).

    Returns:
        A set of `dim`-by-`dim` POVMs of shape `(dim, dim, num_inputs, num_outputs)`.

    Examples:
        We can generate a set of `dim`-by-`dim` POVMs consisting of a specific dimension along with a given number of
        measurement inputs and measurement outputs. As an example, we can construct a random set of \(2\)-by-\(2\)
        POVMs of dimension with \(2\) inputs and \(2\) outputs.

        ```python exec="1" source="above" result="text" session="povm_example"
        import numpy as np
        from toqito.rand import random_povm

        dim, num_inputs, num_outputs = 2, 2, 2

        povms = random_povm(dim, num_inputs, num_outputs)

        print(povms)
        ```


        We can verify that this constitutes a valid set of POVM elements as checking that these operators all sum to the
        identity operator.

        ```python exec="1" source="above" result="text" session="povm_example"
        print(np.round(povms[:, :, 0, 0] + povms[:, :, 0, 1]))
        ```

        It is also possible to add a seed for reproducibility.

        ```python exec="1" source="above" result="text" session="povm_example"
        import numpy as np
        from toqito.rand import random_povm

        dim, num_inputs, num_outputs = 2, 2, 2

        povms = random_povm(dim, num_inputs, num_outputs, seed=42)

        print(povms)
        ```

        We can once again verify that this constitutes a valid set of POVM elements as checking that
        these operators all sum to the identity operator.

        ```python exec="1" source="above" result="text" session="povm_example"
        print(np.round(povms[:, :, 0, 0] + povms[:, :, 0, 1]))
        ```

    """
    povms = []
    gen = np.random.default_rng(seed=seed)
    gram_vectors = gen.normal(size=(num_inputs, num_outputs, dim, dim))
    for input_block in gram_vectors:
        normalizer = sum(np.array(output_block).T.conj() @ output_block for output_block in input_block)
        u_mat, d_mat, _ = np.linalg.svd(normalizer)

        output_povms = []
        for output_block in input_block:
            partial = np.array(output_block, dtype=complex).dot(u_mat).dot(np.diag(d_mat ** (-1 / 2.0)))
            internal = partial.dot(np.diag(np.ones(dim)) ** (1 / 2.0))
            output_povms.append(internal.T.conj() @ internal)
        povms.append(output_povms)

    # This allows us to index the POVMs as [dim, dim, num_inputs, num_outputs].
    povms = np.swapaxes(np.array(povms), 0, 2)
    povms = np.swapaxes(povms, 1, 3)

    return povms