Linear Algebra Operations#

Collection of linear algebra operations that are commonly used in quantum information theory.

Anticommutator#

qibo.quantum_info.anticommutator(operator_1, operator_2)[source]#

Returns the anticommutator of operator_1 and operator_2.

The anticommutator of two matrices \(A\) and \(B\) is given by

\[\{A, B\} = A \, B + B \, A \,.\]
Parameters:

  • operator_1 (ArrayLike) – First operator.

  • operator_2 (ArrayLike) – Second operator.

Returns:

Anticommutator of operator_1 and operator_2.

Return type:

ArrayLike

Commutator#

qibo.quantum_info.commutator(operator_1, operator_2)[source]#

Returns the commutator of operator_1 and operator_2.

The commutator of two matrices \(A\) and \(B\) is given by

\[[A, B] = A \, B - B \, A \,.\]
Parameters:

  • operator_1 (ArrayLike) – First operator.

  • operator_2 (ArrayLike) – Second operator.

Returns:

Commutator of operator_1 and operator_2.

Return type:

ArrayLike

Double commutator#

qibo.quantum_info.double_commutator(operator_1, operator_2, operator_3=None, nesting='right')[source]#

Returns the double commutator of operator_1, operator_2, and operator_3.

The double commutator of three matrices \(A\), \(B\), and \(C\) can be nested to the right,

\[[A, [B, C]] = A \, [B, C] - [B, C] \, A \, ,\]

or to the left,

\[[[A, B], C] = [A, B] \, C - C \, [A, B] \, .\]

The latter is the generator of the Heisenberg equation of a double-bracket rotation [1]. If operator_3 is None, only two operators are involved and the result is \([A, [A, B]]\) for right nesting and \([[A, B], B]\) for left nesting.

Parameters:

  • operator_1 (ArrayLike) – First operator.

  • operator_2 (ArrayLike) – Second operator.

  • operator_3 (ArrayLike, optional) – Third operator. If None, the result is \([A, [A, B]]\) or \([[A, B], B]\) for nesting="right" or nesting="left", respectively. Defaults to None.

  • nesting (str, optional) – either "right", for \([A, [B, C]]\), or "left", for \([[A, B], C]\). Defaults to "right".

Returns:

Double commutator of the three operators.

Return type:

ArrayLike

References

1. M. Gluza, Double-bracket quantum algorithms for diagonalization, Quantum (2024), arXiv:2206.11772.

Group commutator#

qibo.quantum_info.group_commutator(operator_1, operator_2, backend=None)[source]#

Returns the group commutator of two hermitian operators.

The group commutator of two hermitian matrices \(A\) and \(B\) is the unitary

\[V(A, B) = e^{i \, A} \, e^{i \, B} \, e^{-i \, A} \, e^{-i \, B} \, .\]

It is a product formula that approximates the evolution generated by the commutator of \(A\) and \(B\) [1, 2],

\[e^{-[A, B]} = V(A, B) + E \, , \quad \|E\| \leq \|[A, [A, B]]\| + \|[B, [B, A]]\| \, ,\]

in any unitarily invariant norm. Hence, \(V(\sqrt{s} \, A, \sqrt{s} \, B)\) approximates \(e^{-s \, [A, B]}\) with an error of order \(s^{3/2}\).

Parameters:

  • operator_1 (ArrayLike) – First hermitian operator \(A\).

  • operator_2 (ArrayLike) – Second hermitian operator \(B\).

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

Group commutator \(V(A, B)\).

Return type:

ArrayLike

References

1. M. Gluza, Double-bracket quantum algorithms for diagonalization, Quantum (2024), arXiv:2206.11772.

2. Y.-A. Chen, A. M. Childs, M. Hafezi, Z. Jiang, H. Kim, and Y. Xu, Efficient product formulas for commutators and applications to quantum simulation, Phys. Rev. Research 4, 013191 (2022).

Partial trace#

qibo.quantum_info.partial_trace(state, traced_qubits, backend=None)[source]#

Return the density matrix resulting from tracing out traced_qubits from state.

Total number of qubits is inferred by the shape of state.

Parameters:

  • state (ArrayLike) – density matrix or statevector.

  • traced_qubits (list[int] or tuple[int, ...]) – indices of qubits to be traced out.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

Density matrix of the remaining qubit(s).

Return type:

ArrayLike

Partial transpose#

qibo.quantum_info.partial_transpose(operator, partition, backend=None)[source]#

Return matrix after the partial transposition of partition qubits in operator.

Given a \(n\)-qubit operator \(O \in \mathcal{H}_{A} \otimes \mathcal{H}_{B}\), the partial transpose with respect to partition \(B\) is given by

\[\begin{split}\begin{align} O^{T_{B}} &= \sum_{jklm} \, O_{lm}^{jk} \, \ketbra{j}{k} \otimes \left(\ketbra{l}{m}\right)^{T} \\ &= \sum_{jklm} \, O_{lm}^{jk} \, \ketbra{j}{k} \otimes \ketbra{m}{l} \\ &= \sum_{jklm} \, O_{ml}^{jk} \, \ketbra{j}{k} \otimes \ketbra{l}{m} \, , \end{align}\end{split}\]

where the superscript \(T\) indicates the transposition operation, and \(T_{B}\) indicates transposition on partition \(B\). The total number of qubits is inferred by the shape of operator.

Parameters:

  • operator (ArrayLike) – \(1\)- or \(2\)-dimensional operator, or an array of \(1\)- or \(2\)-dimensional operators,

  • partition (Union[List[int], Tuple[int, ...]]) – indices of qubits to be transposed.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

Partially transposed operator(s) \(O^{T_{B}}\).

Return type:

ArrayLike

Matrix exponentiation#

qibo.quantum_info.matrix_exponentiation(matrix, phase=None, eigenvectors=None, eigenvalues=None, backend=None)[source]#

Calculates the exponential of a matrix.

Given a matrix \(H\) and a phase \(\theta\), it returns the exponential of the form

\[\exp\left(\theta \, H \right) \, .\]

If the eigenvectors and eigenvalues are given, the matrix diagonalization is used for the exponentiation.

Parameters:

  • matrix (ArrayLike) – matrix to be exponentiated.

  • phase (float or int or complex, optional) – phase that multiplies the matrix. If None, defaults to \(1\). Defaults to None.

  • eigenvectors (ArrayLike, optional) – if not None, eigenvectors are used to calculate matrix exponentiation as part of diagonalization. Must be used together with eigenvalues. Defaults to None.

  • eigenvalues (ArrayLike, optional) – if not None, eigenvalues are used to calculate matrix exponentiation as part of diagonalization. Must be used together with eigenvectors. Defaults to None.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

matrix exponential \(\exp\left(\theta \, H \right)\).

Return type:

ArrayLike

Matrix logarithm#

qibo.quantum_info.matrix_logarithm(matrix, base=2, eigenvectors=None, eigenvalues=None, backend=None)[source]#

Calculates the logarithm of a matrix.

Given a matrix \(A\) and a log base \(b\), it returns the logarithm of the form

\[\log_{b}\left(A\right) \, .\]

If the eigenvectors and eigenvalues are given, the matrix diagonalization is used for the calculation.

Parameters:

  • matrix (ArrayLike) – matrix to be logarithmed.

  • eigenvectors (ArrayLike, optional) – if not None, eigenvectors are used to calculate the matrix logarithm as part of diagonalization. Must be used together with eigenvalues. Defaults to None.

  • eigenvalues (ArrayLike, optional) – if not None, eigenvalues are used to calculate the matrix logarithm as part of diagonalization. Must be used together with eigenvectors. Defaults to None.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

  • base (float) –

Returns:

Matrix logarithm \(\log_{b}(A)\).

Return type:

ArrayLike

Matrix power#

qibo.quantum_info.matrix_power(matrix, power, precision_singularity=1e-14, backend=None)[source]#

Given a matrix \(A\) and power \(\alpha\), calculate \(A^{\alpha}\).

Parameters:

  • matrix (ArrayLike) – matrix whose power to calculate.

  • power (float or int) – power to raise matrix to.

  • precision_singularity (float, optional) – If determinant of matrix is smaller than precision_singularity, then matrix is considered to be singular. Used when power is negative.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

matrix power \(A^{\alpha}\).

Return type:

ArrayLike

Matrix square root#

qibo.quantum_info.matrix_sqrt(matrix, backend=None)[source]#

Given a matrix \(A\), calculate \(A^{1/2}\).

Parameters:

  • matrix (ArrayLike) – matrix whose power to calculate.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

Matrix power \(A^{1/2}\).

Return type:

ArrayLike

Cosine-Sine decomposition#

qibo.quantum_info.cosine_sin_decomposition(matrix, block_nrows=None, block_ncols=None, backend=None)[source]#

Calculate the Cosine-Sine Decomposition (CSD) of a unitary matrix.

Given an \(N \times N\) unitary matrix \(U\), partitioned into blocks such that the upper left block has shape \(n_{r} \times n_{c}\), with \(n_{r}\) given by block_nrows and \(n_{c}\) by block_ncols, its CSD is given by

\[U = L \, M \, R^{\dagger} \, ,\]

where \(L = L_{1} \oplus L_{2}\) and \(R = R_{1} \oplus R_{2}\) are block-diagonal unitary matrices, with \(L_{1}\), \(L_{2}\), \(R_{1}\), and \(R_{2}\) of shapes \(n_{r} \times n_{r}\), \((N - n_{r}) \times (N - n_{r})\), \(n_{c} \times n_{c}\), and \((N - n_{c}) \times (N - n_{c})\), respectively. \(M\) is an \(N \times N\) orthogonal matrix built from the diagonal matrices \(C = \text{diag}(\cos\theta_{j})\) and \(S = \text{diag}(\sin\theta_{j})\), with \(0 \leq \theta_{j} \leq \pi / 2\). For \(n_{r} = n_{c} = N / 2\),

\[\begin{split}M = \begin{pmatrix} C & -S \\ S & C \end{pmatrix} \, .\end{split}\]

For other partitions, \(M\) also contains identity blocks, which are described in the documentation of scipy.linalg.cossin().

Parameters:

  • matrix (ArrayLike) – unitary matrix whose CSD to calculate.

  • block_nrows (int, optional) – number of rows of the upper left block of matrix. If None, defaults to \(N / 2\). Defaults to None.

  • block_ncols (int, optional) – number of columns of the upper left block of matrix. If None, defaults to \(N / 2\). Defaults to None.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

Cosine-Sine decomposition of \(U\), i.e. \(L\), \(M\), and \(R^{\dagger}\), in that order.

Return type:

Tuple[ArrayLike, ArrayLike, ArrayLike]

References

1. C. C. Paige and M. Wei, History and generality of the CS decomposition, Linear Algebra Appl. 208-209, 303 (1994).

Schur decomposition#

qibo.quantum_info.schur_decomposition(matrix, output='complex', backend=None)[source]#

Calculate the Schur decomposition of a square matrix.

Given an \(N \times N\) matrix \(A\), its Schur decomposition is given by

\[A = Z \, T \, Z^{\dagger} \, ,\]

where \(Z\) is an \(N \times N\) unitary matrix and \(T\) is an upper triangular matrix with the eigenvalues of \(A\) in its diagonal. If \(A\) is normal, which is the case of unitary and Hermitian matrices, then \(T\) is diagonal.

Parameters:

  • matrix (ArrayLike) – square matrix whose Schur decomposition to calculate.

  • output (str, optional) – if "complex", \(T\) is upper triangular. If "real" and matrix is real, \(Z\) is orthogonal and \(T\) is quasi upper triangular, i.e. with \(1 \times 1\) and \(2 \times 2\) blocks in its diagonal. Defaults to "complex".

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

Schur form \(T\) and unitary matrix \(Z\), in that order.

Return type:

Tuple[ArrayLike, ArrayLike]

Singular value decomposition#

qibo.quantum_info.singular_value_decomposition(matrix, backend=None)[source]#

Calculate the Singular Value Decomposition (SVD) of matrix.

Given an \(M \times N\) complex matrix \(A\), its SVD is given by

\[A = U \, S \, V^{\dagger} \, ,\]

where \(U\) and \(V\) are, respectively, an \(M \times M\) and an \(N \times N\) complex unitary matrices, and \(S\) is an \(M \times N\) diagonal matrix with the singular values of \(A\).

Parameters:

  • matrix (ArrayLike) – matrix whose SVD to calculate.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

Singular value decomposition of \(A\), i.e. \(U\), \(S\), and \(V^{\dagger}\), in that order.

Return type:

ArrayLike, ArrayLike, ArrayLike

Schmidt decomposition#

qibo.quantum_info.schmidt_decomposition(state, partition, backend=None)[source]#

Return the Schmidt decomposition of a \(n\)-qubit bipartite pure quantum state.

Given a bipartite pure state \(\ket{\psi}\in\mathcal{H}_{A}\otimes\mathcal{H}_{B}\), its Schmidt decomposition is given by

\[\ket{\psi} = \sum_{k = 1}^{\min\{a, \, b\}} \, c_{k} \, \ket{\phi_{k}} \otimes \ket{\nu_{k}} \, ,\]

with \(a\) and \(b\) being the respective cardinalities of \(\mathcal{H}_{A}\) and \(\mathcal{H}_{B}\), and \(\{\phi_{k}\}_{k\in[\min\{a, \, b\}]} \subset \mathcal{H}_{A}\) and \(\{\nu_{k}\}_{k\in[\min\{a, \, b\}]} \subset \mathcal{H}_{B}\) being orthonormal sets. The coefficients \(\{c_{k}\}_{k\in[\min\{a, \, b\}]}\) are real, non-negative, and unique up to re-ordering.

The decomposition is calculated using qibo.quantum_info.singular_value_decomposition() of the \(a \times b\) matrix of coefficients \(C\), with \(\ket{\psi} = \sum_{jk} \, C_{jk} \, \ket{j} \otimes \ket{k}\), resulting in

\[C = U \, S \, V^{\dagger} \, ,\]

where \(U\) is an \(a \times a\) unitary matrix, \(V\) is an \(b \times b\) unitary matrix, and \(S\) is the array of singular values of \(C\), i.e. the Schmidt coefficients \(c_{k}\).

Parameters:

  • state (ArrayLike) – statevector.

  • partition (Union[List[int], Tuple[int, ...]]) – indices of qubits in one of the two partitions. The other partition is inferred as the remaining qubits.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses qibo.backends.GlobalBackend. Defaults to None.

Returns:

Respectively, the unitary matrix \(U\), the array \(S\) of Schmidt coefficients, and the unitary matrix \(V^{\dagger}\).

Return type:

ArrayLike, ArrayLike, ArrayLike

Lanczos algorithm#

qibo.quantum_info.lanczos(matrix, steps=None, initial_vector=None, precision_tol=1e-08, seed=None, backend=None)[source]#

Lanczos iterative method to tridiagonalize a Hermitian matrix.

Given a \(N \times N\) Hermitian matrix \(H\) and a number of iterations \(m \leq N\), the Lanczos algorithm outputs a \(N \times m\) orthonormal matrix \(U\) and a \(m \times m\) tridiagonal real symmetric matrix \(T = U^{\dagger} \, H \, U\). If \(m = N\), then \(U\) is an unitary matrix. The eigenvalues of \(T\) and \(H\) coincide, while \(U \ket{\mathbf{x}}\) are the eigenvectors of \(H\), with \(\ket{\mathbf{x}}\) being the eigenvectors of \(T\).

This reduces the problem of diagonalization of \(H\) to constructing the matrix \(U\) and diagonalizing \(T\).

With \(\|\cdot\|_{2}\) being the Euclidean norm, the algorithm goes as follows:

  1. Generate random \(\ket{v_{1}} \in \mathbb{C}^{N}\) such that \(\|\ket{v_{1}}\|_{2} = 1\)

  2. \(\ket{\omega_{1}^{\prime}} = H \ket{v_{1}}\)

  3. \(\alpha_{1} = \braket{\omega_{1}^{\prime} | v_{1}}\)

  4. \(\ket{\omega_{1}} = H \ket{v_{1}} - \alpha_{1} \ket{v_{1}}\)

  5. For \(j = 2, \dots, m - 1\):

    1. \(\beta_{j} = \|\omega_{j-1}\|_{2}\)

    2. \(\ket{v_{j}} = \ket{\omega_{j-1}} \, / \, \beta_{j}\) If \(\beta_{j} \neq 0\) else generate random \(\ket{v_{j}}\) such that \(\ket{v_{j}} \perp \{\ket{v_{j^{\prime}}}\}_{j^{\prime} \in [1, j-1]}\)

    3. \(\ket{\omega_{j}^{\prime}} = H \ket{v_{j}}\)

    4. \(\alpha_{j} = \braket{\omega_{j}^{\prime} | v_{j}}\)

    5. \(\ket{\omega_{j}} = \ket{\omega_{j}^{\prime}} - \alpha_{j} \ket{v_{j}} - \beta_{j} \ket{v_{j-1}}\)

The columns of the orthogonal matrix \(U\) are the Lanczos vectors \(\{\ket{v_{j}}\}_{j\in[1, m]}\).

Parameters:

  • matrix (ArrayLike) – square Hermitian matrix to be tridiagonalized.

  • steps (int, optional) – number of iterations \(m\), between \(1\) and the size of matrix. If None, defaults to the size of matrix. Defaults to None.

  • initial_vector (ArrayLike, optional) – vector to be used as the initial Lanczos vector \(\ket{v_{1}}\), which is normalized. If None, array is uniformly sampled. Defaults to None.

  • precision_tol (float, optional) – precision threshold such that for \(\beta_{j}\) smaller than precision_tol, it is considered to be zero.

  • seed (int or numpy.random.Generator, optional) – Seed for the initial random vector \(\ket{v_{1}}\) Either a generator of random numbers or a fixed seed to initialize a generator. If None, initializes a generator with a random seed. Defaults to None.

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

Tridiagonal matrix and the orthogonal matrix of Lanczos vectors, respectively.

Return type:

(ArrayLike, ArrayLike)

References

1. Lanczos, C. An iteration method for the solution of the eigenvalue problem of linear differential and integral operators, Journal of Research of the National Bureau of Standards. 45 (4): 255–282 (1950).

Lie closure#

qibo.quantum_info.lie_closure(generators, max_iterations=10000, tol=1e-10, backend=None)[source]#

Compute the dynamical Lie algebra (DLA) generated by a set of operators.

The generators are either matrices, Pauli strings (e.g. "XXI"), or real linear combinations of Pauli strings given as dictionaries (e.g. {"XX": 1.0, "ZI": 0.5}). Given generators \(\{G_{j}\}_{j}\), each either Hermitian or skew-Hermitian, the DLA is the real span of all nested commutators \([G_{j}, [G_{k}, \dots, [G_{l}, G_{m}]]]\) of the generators. Following the physics convention, the algebra is \(\{i \, G_{\alpha}\}_{\alpha}\) and the returned operators \(G_{\alpha}\) are Hermitian. Commutators are added round by round: round \(r\) contains the commutators of the operators added in round \(r - 1\) with the original generators, i.e. nesting depth \(r\).

The commutator of two Pauli strings is either zero or proportional to a Pauli string. Hence, Pauli generators are handled on the phase-space (tableau) representation of Paulis, in which the product of two Paulis is the XOR of their bit vectors \((x | z)\), and they anticommute if and only if their symplectic inner product is \(1\). This avoids building \(2^{n} \times 2^{n}\) matrices. If all generators are Pauli strings, the DLA is itself spanned by Pauli strings. Otherwise, operators are stored as sparse coefficient vectors over the Pauli strings encountered, and commutators are computed term by term.

Parameters:

  • generators (list[ArrayLike] or list[str] or list[dict]) – either square matrices of identical shape, each Hermitian or skew-Hermitian, or Pauli strings of identical length composed of "I", "X", "Y" and "Z", or dictionaries mapping such Pauli strings to real coefficients. Pauli strings and dictionaries can be mixed, but not combined with matrices.

  • max_iterations (int, optional) – maximum nesting depth of commutators. If reached before the algebra is closed, a warning is logged. Defaults to \(10000\).

  • tol (float, optional) – threshold on the norm below which a candidate is considered linearly dependent on the current basis (and, for dictionaries, below which a coefficient is set to zero). Ignored if all generators are Pauli strings, which are handled exactly. If an operator is found to be independent with a residual norm close to tol, it may be numerical noise, and a warning is logged. This can happen for dense, ill-conditioned combinations of generators, for which a larger tol may be needed. Defaults to \(10^{-10}\).

  • backend (qibo.backends.abstract.Backend, optional) – backend to be used in the execution. If None, it uses the current backend. Defaults to None.

Returns:

For matrices, array of shape (dim, d, d) with a Hermitian basis of the DLA, orthonormal with respect to the Hilbert-Schmidt inner product \(\text{Tr}(A^{\dagger} \, B)\). For Pauli strings only, the list of Pauli strings (all with positive sign) that span the DLA. Otherwise, a list of dictionaries, each mapping Pauli strings to real coefficients with unit Euclidean norm (orthonormal basis).

Return type:

ArrayLike or list[str] or list[dict]

Example

Transverse-field Ising model on two qubits, whose DLA has dimension \(6\):

from qibo import matrices
from qibo.backends import NumpyBackend
from qibo.quantum_info import lie_closure

backend = NumpyBackend()
I, X, Z = matrices.I, matrices.X, matrices.Z
generators = [
    backend.kron(X, X),
    backend.kron(Z, I),
    backend.kron(I, Z),
]

dla = lie_closure(generators, backend=backend)
print(dla.shape)  # (6, 4, 4)

The same algebra from Pauli strings:

from qibo.quantum_info import lie_closure

dla = lie_closure(["XX", "ZI", "IZ"])
print(dla)  # ['XX', 'ZI', 'IZ', 'YX', 'XY', 'YY']

Heisenberg chain on four qubits, with generators that are sums of Pauli strings:

from qibo.quantum_info import lie_closure

generators = [
    {"XXII": 1.0, "YYII": 1.0, "ZZII": 1.0},
    {"IXXI": 1.0, "IYYI": 1.0, "IZZI": 1.0},
    {"IIXX": 1.0, "IIYY": 1.0, "IIZZ": 1.0},
]

dla = lie_closure(generators)
print(len(dla))  # 12

References

1. M. Larocca et al., Diagnosing barren plateaus with tools from quantum optimal control, Quantum 6, 824 (2022). 2. S. Aaronson and D. Gottesman, Improved Simulation of Stabilizer Circuits, Phys. Rev. A 70, 052328 (2004).