:hammer_and_wrench: pyQSP

A Python Package for Quantum Signal Processing

Recent news and updates

:round_pushpin: (10/24/2024) The codebase has recently been updated to work internally almost entirely with polynomials in the Chebyshev, rather than the monomial, basis. This improves the numerical stability of the laurent method substantially, and should present almost no visible differences for an end-user. Currently the poly2angles command-line function is the only major input for monomial basis targets, while the main QuantumSignalProcessingPhases method and various helper and plotting methods all accept only Chebyshev basis inputs with internal checks. For conversion between bases, see below.

Introduction

Quantum signal processing (QSP) is a flexible quantum algorithm subsuming Hamiltonian simulation, quantum linear systems solvers, amplitude amplification, and many other common quantum algorithms. Moreover, for each of these applications, QSP often exhibits state-of-the-art space and time complexity, with comparatively simple analysis. QSP and its related algorithms, e.g., Quantum singular value transformation (QSVT), using basic alternating circuit ansätze originally inspired by composite pulse techniques in NMR, permit one to transform the spectrum of linear operators by near arbitrary polynomial functions, with the aforementioned good numerical properties and simple analytic treatment.

In their most basic forms, QSP/QSVT give a recipe for a desired spectral transformation of a unitary[^1] $U$, requiring access to an auxiliary qubit, a controlled version of $U$, and single-qubit rotations on the auxiliary qubit. A standard application of QSP/QSVT might look like the following:

• Given a function one wants to apply to the spectrum of a unitary, classically generate a (Laurent) polynomial that suitably approximates this function over a desired spectral range.
• Having computed a good polynomial approximation, and checking that it obeys certain mild conditions, use one among many efficient classical algorithms to compute the sequence of single-qubit rotations (the QSP phases) interspersing applications of the controlled unitary (the QSP signal) corresponding to the polynomial approximation.
• Run the corresponding sequence of gates as a quantum circuit, interleaving QSP signals and phases, followed by a measurement in a chosen basis.

:round_pushpin: The theory of QSP is not only under active development, but comprises multiple subtly different conventions, each of which can use different terminology compared the barebones outline given here. These included conventions for how the signal is encoded, how the phases are applied, the basis to measure in, whether one desires to transform eigenvalues or singular values, whether the classical algorithm to find these phases is exact or iterative, and so on.

Regardless, the basic scheme of QSP and QSVT is relatively fixed: given a specific circuit ansatz and a theory for the polynomial transformations achievable for that ansatz, generate those conditions and algorithms relating the achieved function and circuit parameterization. Understanding the bidirectional map between phases and polynomial transforms, as well as the efficiency of loading linear systems into quantum processes, constitutes most of the theory of these algorithms.

This package provides such conditions and algorithms, and automatically treats a few common conventions, with options to modify the code in basic ways to encompass others. These conventions are enumerated in the pedagogical work A Grand Unification of Quantum Algorithms, and the QSP phase-finding algorithms we treat can be broken roughly into three types:

• :hammer: The laurent method employs techniques originated in Finding Angles for Quantum Signal Processing with Machine Precision, and extends code from its attached repository. This method exactly computes phases by studying the properties of the desired polynomials using a divide-and-conquer approach.
• :sparkles: The tf method employs TensorFlow + Keras, and finds QSP phase angles using straightforward (but sometimes slow) optimization techniques.
• :key: The sym_qsp method employs an iterative, quasi-Newton method to find QSP phases for a special, lightly-restricted sub-class of QSP protocols whose phases are symmetric. In comparison with the above two techniques, this method is almost invariably quick, numerically stable, and should suit most applications. It is based off work in Efficient phase-factor evaluation in quantum signal processing, and Matlab implementations in the QSPPACK repository. See the affiliated CITATION file for a full list of references.

:warning: As methods for numerically handling QSP protocols have been refined, we have tried to update this repository to reflect leading methods. Along the way, we have also had to lightly deprecate older methods (which may still be in use by others who use this repository). In the sections that follow, we try to give special attention to where a new user might begin to use this repository in the most frictionless way.

---

[^1]: In general, QSVT and related algorithms can compute matrix polynomials in wide classes of linear operators, from normal operators to non-square operators to infinite-dimensional operators, each with corresponding limitations and a related circuit form. For our purposes, treating normal operators will nearly always be sufficient, and captures most of the essence of QSP/QSVT.

An overview of the QSP ansatz and conventions

As organized in this introductory pedagogical overview, there are two major conventions used in QSP literature: $W_x$ convention, where the signal $W(a)$ is an $X$-rotation and the QSP phases correspond to $Z$-rotations, and the $W_z$ convention, where the signal $W(a)$ is a $Z$-rotation and the QSP phases correspond to $X$-rotations.

Concretely, in the $W_x$ convention, the overall unitary for some list of QSP phases $\Phi \in \mathbb{R}^{n + 1}$ is:

    U_x = e^{i\phi_0 Z} \prod_{k=1}^n W(a) e^{i\phi_k Z} \;\;\;\;\text{where}\;\;\;\; W(x)= \begin{bmatrix} a & i\sqrt{1-a^2} \\ i\sqrt{1-a^2} & a \end{bmatrix}

while in the $W_z$ convention, the overall unitary is:

    U_z = e^{i\phi_0 X} \prod_{k=1}^n W(a) e^{i\phi_k X} \;\;\;\;\text{where}\;\;\;\; W(a)=\begin{bmatrix} a + i\sqrt{1-a^2} & 0 \\ 0 & a - i\sqrt{1-a^2} \end{bmatrix}

As one might guess, these conventions are related by a Hadamard transform:

    U_x = H U_z H

The $W_z$ convention is convenient for and employed in Laurent polynomial formulations of QSP, while the $W_x$ convention is older and perhaps more widespread currently, e.g. as employed in quantum singular value transform.

In the first convention, the resulting QSP unitary is also given a standard form (which can be shown by induction), namely that

    U_x = e^{i\phi_0 Z} \prod_{k=1}^n W(a) e^{i\phi_k Z} = \begin{bmatrix} P_x(a) & iQ_x(a)\sqrt{1-a^2} \\ i Q_x^*(a)\sqrt{1-a^2} & P_x^*(a) \end{bmatrix},

where $P$ and $Q$ are polynomials of degree $n$ and $n-1$ respectively, with definite parity, and satisfying the condition $|P|^2 + |Q|^2(1 - a^2) = 1$ for all $a \in [-1,1]$. Evidently in a different basis the parity and degree of these polynomials may change, and moreover sometimes the substitution $a \rightarrow (b + 1/b)/2$ is sometimes made, moving from polynomials over $[-1,1]$ to Laurent polynomials over the unit circle in the complex plane. Each of these choices has benefits and limitations, but for most initial presentations, the $U_x$ convention with polynomials over $[-1,1]$ is used.

As stated, this package can generate QSP phases in both conventions. The impediment to immediately freely working between both conventions is that if one wants a certain polynomial $P_x(a) = \langle 0|U_x|0\rangle$ in the $W_x$ convention, one cannot just use the phases generated for this polynomial in the $W_z$ convention. Instead, first the $Q_x(a)$ corresponding to $P_x(a)$ is needed to complete the full $U_x$. From this one can compute $P_z(a) = \langle 0|U_z|0\rangle = P_x(a) + Q_x(a)$. Computing the QSP phases for this $P_z(a)$ in the $W_z$ convention yields the desired QSP phases for $P_x(a)$ in the $W_x$ convention.

:warning: In addition to specifying the signal unitary and the signal processing phase rotations, the implicit measurement basis must also be specified. In this codebase the default basis is $|\pm\rangle$ (i.e., the $X$ Pauli's eigenbasis), though the code also allows for the computational basis (e.g., when using tf optimization method to generate the phase angles). For information on this, see the --measurement option below.

The guiding principle to take away from the discussion above is the following: the choice of circuit convention can change the conditions required of the achieved polynomial transforms (e.g., their parity, degree, boundary conditions, etc.). That said, the classical subroutines used to find good polynomial approximations remain mostly unchanged, and for most applications the restrictions on achieved functions are not crucial for performance.

A few quick-and-dirty examples

In this section we aim to teach basic use-cases for this package through a series of related examples. Given the heterogeneous research efforts in QSP, there are multiple ways to do most