Quick Guide

To introduce the main functionality of Haarpy consider the following problem: imagine that you can generate unitary matrices at random (from the Haar measure); you would like to estimate the average of \(|U_{i,j}|^2\) which we write mathematically as \(\int dU |U_{i,j}|^2 = \int dU U_{i,j} U_{i,j}^*\). We could obtain this average by using the random matrix functionality of SciPy as follows:

import numpy as np
from scipy.stats import unitary_group

np.random.seed(137)
shots = 1000

# For unitary matrices of size between 2 and 4 produce 1000 shots and calculate the average
# of the absolute value squared of the 0,1 entry
np.array(
    [
        np.mean([np.abs(unitary_group.rvs(dim=dim)[0, 1]) ** 2 for _ in range(shots)])
        for dim in range(2, 5)
    ]
)

array([0.4964599 , 0.32742463, 0.25429793])

Haarpy allows you to obtain this average (and many others!) analytically. We first recall that the expression we are trying to calculate is \(\int dU \left|U_{i,j}\right|^2 = \int dU U_{i,j} U_{i,j}^*\). With this expression in mind we can use Sympy to create a symbolic variable \(d\) for the dimension of the unitary and write

from sympy import Symbol
from haarpy import haar_integral_unitary

d = Symbol("d")
haar_integral_unitary(("i", "j", "i", "j"), d)

$\displaystyle \frac{1}{d}$

We can also put integers

haar_integral_unitary(("i", "j", "i", "j"), 3)

Fraction(1, 3)

Notice the order of the indices! The first "i" and "j" are the indices of \(U\) while the second pair of "i" and "j" are the indices of \(U^*\).

Imagine that now we want to calculate something like \(\int dU U_{i,m} U_{j,n} U_{k,o} U_{i,m}^* U_{j,n}^* U_{k,p}^*\) \(= \int dU \left|U_{i,m} U_{j,n} U_{k,o}\right|^2\) we simply do

haar_integral_unitary(("ijk", "mno", "ijk", "mno"), d)

$\displaystyle \frac{d^{2} - 2}{d \left(d - 2\right) \left(d - 1\right) \left(d + 1\right) \left(d + 2\right)}$

The averages we are calculating are obtained by using so-called Weingarten calculus.