Table Of Contents

Previous topic

Quantum Information Toolkit

Next topic

Linear maps (qit.lmap)

This Page

Utilities (qit.utils)

This module contains utility functions which do not logically fit anywhere else.

Mathematical utilities

gcd(a, b) Greatest common divisor.
lcm(a, b) Least common multiple.
majorize(x, y) Majorization partial order of real vectors.

Matrix functions

comm(A, B) Array commutator.
acomm(A, B) Array anticommutator.
mkron(*arg) This is how kron should work, dammit.
rank(A[, tol]) Matrix rank.
orth(A[, tol, kernel, spectrum]) Construct an orthonormal basis for the range of A using SVD
nullspace(A[, tol, spectrum]) Kernel of a matrix.
nullspace_hermitian(A[, tol]) Kernel of a superoperator matrix restricted to the Hermitian subspace.
projector(v) Projector corresponding to vector v.
eighsort(A) Returns eigenvalues and eigenvectors of a Hermitian matrix, sorted in nonincreasing order.
expv(t, A, v[, tol, m, iteration]) Multiply a vector by an exponentiated matrix.

Random matrices

rand_hermitian(n) Random Hermitian n*n matrix.
rand_U(n) Random U(n) matrix.
rand_SU(n) Random SU(n) matrix.
rand_U1(n) Random diagonal unitary matrix.
rand_pu(n) Random n-partition of unity.
rand_positive(n) Random n*n positive semidefinite matrix.
rand_SL(n) Random SL(n, C) matrix.

Superoperators

vec(rho) Flattens a matrix into a vector.
inv_vec(v[, dim]) Reshapes a vector into a matrix.
lmul(L[, q]) Superoperator equivalent for multiplying from the left.
rmul(R[, p]) Superoperator equivalent for multiplying from the right.
lrmul(L, R) Superoperator equivalent for multiplying both from left and right.
superop_lindblad(A[, H]) Liouvillian superoperator for a set of Lindblad operators.
superop_fp(L[, tol]) Fixed point states of a Liouvillian superoperator.

Physics

angular_momentum(*args) Angular momentum matrices.
boson_ladder(*args) Bosonic ladder operators.
fermion_ladder(*args) Fermionic ladder operators.

Bases, decompositions

spectral_decomposition(A) Spectral decomposition of a Hermitian matrix.
gellmann(*args) Gell-Mann matrices.
tensorbasis(n[, d, get_locality]) Hermitian tensor-product basis for End(H).

Miscellaneous

assert_o(actual, desired, tolerance) Octave-style assert.
copy_memoize(func) Memoization decorator for functions with immutable args, returns deep copies.
op_list(G, dim) Operator consisting of k-local terms, given as a list.
qubits(n) Dimension vector for an all-qubit system.
R_nmr(theta, phi) SU(2) rotation \theta_\phi (NMR notation).
R_x(theta) SU(2) x-rotation.
R_y(theta) SU(2) y-rotation.
R_z(theta) SU(2) z-rotation.
qit.utils.assert_o(actual, desired, tolerance)

Octave-style assert.

qit.utils.copy_memoize(func)

Memoization decorator for functions with immutable args, returns deep copies.

qit.utils.gcd(a, b)

Greatest common divisor.

Uses the Euclidean algorithm.

qit.utils.lcm(a, b)

Least common multiple.

qit.utils.majorize(x, y)

Majorization partial order of real vectors.

Parameters:x, y (vector) – real vectors, dimension d
Returns:True iff x \preceq y

Returns True iff x is majorized by y, denoted by x \preceq y. This is equivalent to

\sum_{k=1}^n x^{\downarrow}_k \le \sum_{k=1}^n y^{\downarrow}_k \quad \text{for all} \quad n \in \{1, 2, \ldots, d\},

where x^{\downarrow} is the vector x with the elements sorted in nonincreasing order.

x \preceq y if and only if x is in the convex hull of all the coordinate permutations of y.

qit.utils.comm(A, B)

Array commutator.

Returns [A, B] := A*B - B*A

qit.utils.acomm(A, B)

Array anticommutator.

Returns {A, B} := A*B + B*A

qit.utils.mkron(*arg)

This is how kron should work, dammit.

Returns the tensor (Kronecker) product X = A \otimes B \otimes \ldots

qit.utils.projector(v)

Projector corresponding to vector v.

qit.utils.eighsort(A)

Returns eigenvalues and eigenvectors of a Hermitian matrix, sorted in nonincreasing order.

qit.utils.expv(t, A, v, tol=1e-07, m=30, iteration=u'arnoldi')

Multiply a vector by an exponentiated matrix.

Approximates exp(t A) v using a Krylov subspace technique. Efficient for large sparse matrices. The basis for the Krylov subspace is constructed using either Arnoldi or Lanczos iteration.

Input: t vector of nondecreasing time instances >= 0 A n*n matrix (usually sparse) (as an (n,n)-shaped ndarray) v n-dimensional vector (as an (n,)-shaped ndarray) tol tolerance m Krylov subspace dimension, <= n iteration ‘arnoldi’ or ‘lanczos’. Lanczos is faster but requires a Hermitian A.

Output: W result matrix, W[i,:] \approx \exp(t[i] A) v error total truncation error estimate hump \max_{s \in [0, t]}  \| \exp(s A) \|

Uses the sparse algorithm from [21].

qit.utils.rank(A, tol=None)

Matrix rank.

qit.utils.orth(A, tol=None, kernel=False, spectrum=False)

Construct an orthonormal basis for the range of A using SVD

A : array, shape (M, N)

Q : array, shape (M, K)
Orthonormal basis for the range of A. K = effective rank of A, as determined by tolerance

Adaptation of scipy.linalg.orth with optional absolute tolerance.

qit.utils.nullspace(A, tol=None, spectrum=False)

Kernel of a matrix.

Given a matrix A (and optionally a tolerance tol), returns a basis for the kernel (null space) of A in the columns of Z.

qit.utils.nullspace_hermitian(A, tol=None)

Kernel of a superoperator matrix restricted to the Hermitian subspace.

Solves the intersection of the kernel of the superop A and the Hermitian subspace. A maps d*d matrices to whatever.

Singular values <= tol are considered zero.

qit.utils.rand_hermitian(n)

Random Hermitian n*n matrix.

Returns a random Hermitian matrix of size n*n. NOTE: The randomness is not defined in any deeply meaningful sense.

qit.utils.rand_U(n)

Random U(n) matrix.

Returns a random unitary n*n matrix. The matrix is random with respect to the Haar measure.

Uses the algorithm in [16].

qit.utils.rand_SU(n)

Random SU(n) matrix.

Returns a random special unitary n*n matrix. The matrix is random with respect to the Haar measure.

qit.utils.rand_U1(n)

Random diagonal unitary matrix.

Returns a random diagonal unitary n*n matrix. The matrix is random with respect to the Haar measure.

qit.utils.rand_pu(n)

Random n-partition of unity.

Returns the n-partition p, which is random with respect to the order-n Dirichlet distribution Dir(\alpha) with \alpha = (1, 1, ..., 1).

qit.utils.rand_positive(n)

Random n*n positive semidefinite matrix.

Normalized to Tr(A) = 1. Since A has all-real eigenvalues, it is Hermitian by construction.

qit.utils.rand_SL(n)

Random SL(n, C) matrix.

Returns a random special linear n*n matrix. NOTE: The randomness is not defined in any deeply meaningful sense.

qit.utils.vec(rho)

Flattens a matrix into a vector.

Matrix rho is flattened columnwise into a column vector v.

Used e.g. to convert state operators to superoperator representation.

qit.utils.inv_vec(v, dim=None)

Reshapes a vector into a matrix.

Given dim == (n, m), reshapes vector v (length n*m) into a matrix \rho (shape == dim), using column-major ordering. If dim is not given, \rho is assumed to be square.

Used e.g. to convert state operators from superoperator representation to standard matrix representation.

qit.utils.lmul(L, q=None)

Superoperator equivalent for multiplying from the left.

Parameters:
  • L (array) – matrix, shape (m, p)
  • q (int) – The shape of \rho is (p, q). If q is not given, \rho is assumed square.

Returns the superoperator that implements left multiplication of a vectorized matrix \rho by the matrix L.

L \rho = \text{inv\_vec}(\text{lmul}(L) \text{vec}(\rho))

qit.utils.rmul(R, p=None)

Superoperator equivalent for multiplying from the right.

Parameters:
  • R (array) – matrix, shape (q, r)
  • p (int) – The shape of \rho is (p, q). If p is not given, \rho is assumed square.

Returns the superoperator that implements right multiplication of a vectorized matrix \rho by the matrix R.

\rho R = \text{inv\_vec}(\text{rmul}(R) \text{vec}(\rho))

qit.utils.lrmul(L, R)

Superoperator equivalent for multiplying both from left and right.

Parameters:
  • L (array) – matrix, shape (m, p)
  • R (array) – matrix, shape (q, r)

Returns the superoperator that implements left multiplication by the matrix L and right multiplication by the matrix R of a vectorized matrix \rho.

L \rho R = \text{inv\_vec}(\text{lrmul}(L, R) \text{vec}(\rho))

qit.utils.superop_lindblad(A, H=None)

Liouvillian superoperator for a set of Lindblad operators.

A is a vector of traceless, orthogonal Lindblad operators. H is an optional Hamiltonian operator.

Returns the Liouvillian superoperator L corresponding to the diagonal-form Lindblad equation

\dot{\rho} = \text{inv\_vec}(L \text{vec}(\rho)) = -i [H, \rho] +\sum_k \left(A_k \rho A_k^\dagger -\frac{1}{2} \{A_k^\dagger A_k, \rho\}\right)

qit.utils.superop_fp(L, tol=None)

Fixed point states of a Liouvillian superoperator.

Finds the intersection of the kernel of the Liouvillian L with the set of valid state operators, thus giving the set of fixed point states for the quantum channel represented by the master equation

\dot{\rho} = \text{inv\_vec}(L \text{vec}(\rho)).

Let size(L) == [D, D] (and d = sqrt(D) be the dimension of the Hilbert space).

Returns the D*n array A, which contains as its columns a set of n vectorized orthogonal Hermitian matrices (with respect to the Hilbert-Schmidt inner product) which “span” the set of FP states in the following sense:

vec(\rho) = A c,  \quad \text{where} \quad c \in \text{R}^n \quad \text{and} \quad c_1 = 1.

A[:,0] is the shortest vector in the Hermitian kernel of L that has trace 1, the other columns of A are traceless and normalized. Hence, A defines an (n-1)-dimensional hyperplane.

A valid state operator also has to fulfill \rho \ge 0. These operators form a convex set in the Hermitian trace-1 hyperplane defined by A. Currently this function does nothing to enforce positivity, it is up to the user to choose the coefficients a_k such that this condition is satisfied.

Singular values of L less than or equal to the tolerance tol are treated as zero.

qit.utils.angular_momentum(*args)

Angular momentum matrices.

Parameters:d (int) – dimension, the corresponding angular momentum quantum number is j = (d-1)/2
Returns:3-tuple of angular momentum matrices (Jx, Jy, Jz)

Returns the dimensionless angular momentum matrices \vec{J} / \hbar for the d-dimensional subspace.

The angular momentum matrices fulfill the commutation relation [J_x, J_y] = i J_z and all its cyclic permutations.

qit.utils.boson_ladder(*args)

Bosonic ladder operators.

Parameters:d (int) – truncation dimension
Returns:bosonic annihilation operator b

Returns the d-dimensional approximation of the bosonic annihilation operator b for a single bosonic mode in the number basis \{|0\rangle, |1\rangle, ..., |d-1\rangle\}.

The corresponding creation operator is b^\dagger. The (infinite-dimensional) bosonic annihilation and creation operators fulfill the commutation relation [b, b^\dagger] = I. Due to the d-dimensional basis truncation, this does not hold for the last dimension of the b matrices returned by this function.

qit.utils.fermion_ladder(*args)

Fermionic ladder operators.

Parameters:n (int) – number of fermionic modes
Returns:array of fermionic annihilation operators f_k

Returns the fermionic annihilation operators for a system of n fermionic modes in the second quantization.

The annihilation operators are built using the Jordan-Wigner transformation for a chain of n qubits, where the state of each qubit denotes the occupation number of the corresponding mode. First define annihilation and number operators for a lone fermion mode:

\sigma^- &:= (\sigma_x + i \sigma_y)/2,\\
n &:= \sigma^+ \sigma^- = (I-\sigma_z)/2,\\
&\sigma^-|0\rangle = 0, \quad \sigma^-|1\rangle = |0\rangle, \quad n|k\rangle = k|k\rangle.

Then define a phase operator to keep track of sign changes when permuting the order of the operators: \phi_k := \sum_{j=1}^{k-1} n_j. Now, the fermionic annihilation operators for the n-mode system are given by

f_k := (-1)^{\phi_k} \sigma^-_k = \left(\bigotimes_{j=1}^{k-1} {\sigma_z}_j \right) \sigma^-_k.

These operators fulfill the required anticommutation relations:

\{f_j, f_k\}  &= 0,\\
\{f_j, f_k^\dagger\} &= I \delta_{jk},\\
f_k^\dagger f_k &= n_k.

qit.utils.R_nmr(theta, phi)

SU(2) rotation \theta_\phi (NMR notation).

Returns the one-qubit rotation by angle theta about the unit vector [\cos(\phi), \sin(\phi), 0], or \theta_\phi in the NMR notation.

qit.utils.R_x(theta)

SU(2) x-rotation.

Returns the one-qubit rotation about the x axis by the angle theta, e^{-i \sigma_x \theta/2}.

qit.utils.R_y(theta)

SU(2) y-rotation.

Returns the one-qubit rotation about the y axis by the angle theta, e^{-i \sigma_y \theta/2}.

qit.utils.R_z(theta)

SU(2) z-rotation.

Returns the one-qubit rotation about the z axis by the angle theta, e^{-i \sigma_z \theta/2}.

qit.utils.spectral_decomposition(A)

Spectral decomposition of a Hermitian matrix.

Returns the unique eigenvalues a and the corresponding projectors P for the Hermitian matrix A, such that A = \sum_k  a_k P_k.

qit.utils.gellmann(*args)

Gell-Mann matrices.

Returns the d**2 - 1 (traceless, Hermitian) Gell-Mann matrices of dimension d, normalized such that \mathrm{Tr}(G_i^\dagger G_j) = \delta_{ij}. They form an orthonormal basis for the real vector space of d*d traceless Hermitian matrices.

The matrices are returned in the array A, arranged such that the n:th Gell-Mann matrix is A[n,:,:].

qit.utils.tensorbasis(n, d=None, get_locality=False)

Hermitian tensor-product basis for End(H).

Returns a Hermitian basis for linear operators on the Hilbert space H which shares H’s tensor product structure. The basis elements are tensor products of Gell-Mann matrices (which in the case of qubits are equal to Pauli matrices). The basis elements are normalized such that \mathrm{Tr}(b_i^\dagger b_j) = \delta_{ij}.

Input is either two scalars, n and d, in which case the system consists of n qu(d)its, H = C_d^{\otimes n}, or the vector dim, which contains the dimensions of the individual subsystems: H = C_{dim[0]} \otimes ... \otimes C_{dim[n-1]}.

In addition to expanding Hermitian operators on H, this basis can be multiplied by the imaginary unit to obtain the antihermitian generators of U(prod(dim)).

qit.utils.op_list(G, dim)

Operator consisting of k-local terms, given as a list.

Parameters:
  • G (list) – list of k-local operator terms
  • dim (vector) – vector of subsystem dimensions
Returns:

operator O defined by the connection list

G is a list of arrays, G = [c_1, c_2, ..., c_n], where each array c_i corresponds to a term in O.

An array that has 2 columns and k rows, c_i = [(A1, s1), (A2, s2), ... , (Ak, sk)], where Aj are arrays and sj subsystem indices, corresponds to the k-local term given by the tensor product

A_1^{(s_1)} A_2^{(s_2)} \cdots A_k^{(s_k)}.

The dimensions of all operators acting on subsystem sj must match dim[sj].

Alternatively one can think of G as defining a hypergraph, where each subsystem corresponds to a vertex and each array c_i in the list describes a hyperedge connecting the vertices {s1, s2, ..., sk}.

Example: The connection list G = [[(sz,1)], [(sx,1), (sx,3)], [(sy,1), (sy,3)], [(sz,1), (sz,3)], [(sz,2)], [(A,2), (B+C,3)], [(2*sz,3)]] corresponds to the operator

\sigma_{z1} +\sigma_{z2} +2 \sigma_{z3} +\sigma_{x1} \sigma_{x3} +\sigma_{y1} \sigma_{y3} +\sigma_{z1} \sigma_{z3} +A_2 (B+C)_3.

qit.utils.qubits(n)

Dimension vector for an all-qubit system.

For the extemely lazy, returns (2,) * n