Generalized quantum singular value transformation with application in quantum conjugate gradient least squares algorithm

TL;DR

This paper extends generalized quantum signal processing to general matrices, creating GQSVT, and applies it to develop a hybrid quantum-classical conjugate gradient least squares algorithm.

Contribution

It introduces GQSVT, relaxing polynomial parity restrictions, and demonstrates its application in a hybrid quantum-classical algorithm for least squares problems.

Findings

GQSVT relaxes parity restrictions compared to QSVT.

The proposed hybrid algorithm applies GQSVT to solve least squares.

Extension of GQSP to general matrices broadens quantum signal processing capabilities.

Abstract

Quantum signal processing (QSP) and generalized quantum signal processing (GQSP) are essential tools for implementing the block encoding of matrix functions. The achievable polynomials of QSP have restrictions on parity, while GQSP eliminates these restrictions. But GQSP only constructs functions of unitary matrices. In this paper, we further investigate GQSP and extend it to general matrices. Compared with the quantum singular value transformation (QSVT), our proposed method relaxes the requirements on the parity of polynomials. We refer to this extension as generalized quantum singular value transformation (GQSVT). Subsequently, by utilizing the relationship between generalized matrix functions and standard matrix functions, we propose a classical-quantum hybrid quantum conjugate gradient least squares (CGLS) algorithm using GQSVT.

Equations

$$f(A)=Uf(\Lambda )U^{†},$$

$$f^{\diamond }(A)=Wf(\Sigma )V^{†},$$

$$f^L(A)=Wf(\Sigma )W^{†}\text{and}{f}^R(A)=Vf(\Sigma )V^{†},$$

$$q(A^{†}A)A^{†}=q(V{\Sigma }^2{V}^{†})V\Sigma W^{†}=Vq({\Sigma }^2)\Sigma W^{†}=Vp(\Sigma )W^{†}=p^{\diamond }(A^{†}).$$

$$T_j(x):=\cos (j\eta ),$$

$${\int }_{-1}^1{T}_j(x)T_k(x)\frac{dx}{\sqrt{1-x^2}}$$

$$f(x)=\sum _{j=0}^d{a}_j{T}_j(x),$$

$$a_j$$

$$f(\cos (\eta ))=\sum _{j=0}^d{a}_j\cos (j\eta )=\sum _{j=0}^d{a}_j\frac{e^{ij\eta }+e^{-ij\eta }}{2}=\sum _{j=-d}^d{c}_j{e}^{ij\eta }:={\tilde{f}}_d(e^{i\eta }),$$

$$M=\left(\begin{array}[]{cc}U&0\\ 0&I\\ \end{array}\right).$$

$$R(\theta,\phi,\lambda)=\left(\begin{array}[]{cc}e^{i(\lambda+\phi)}\cos\theta&e^{i\phi}\sin\theta\\ e^{i\lambda}\sin\theta&-\cos\theta\\ \end{array}\right)\otimes I.$$

$$R_0:=R({\theta }_0,{\varphi }_0,\lambda )={{\rm Y}}_0\otimes I,R_j:=R({\theta }_j,{\varphi }_j,0)={{\rm Y}}_j\otimes I(j\ge 1),$$

$$\Upsilon_{0}:=\left(\begin{array}[]{cc}e^{i(\lambda+\phi_{0})}\cos\theta_{0}&e^{i\phi_{0}}\sin\theta_{0}\\ e^{i\lambda}\sin\theta_{0}&-\cos\theta_{0}\\ \end{array}\right),~~\Upsilon_{j}:=\left(\begin{array}[]{cc}e^{i\phi_{j}}\cos\theta_{j}&e^{i\phi_{j}}\sin\theta_{j}\\ \sin\theta_{j}&-\cos\theta_{j}\\ \end{array}\right).$$

$$\left(\begin{array}[]{cc}P(U)&\cdot\\ \cdot&\cdot\\ \end{array}\right)=(\prod_{j=1}^{d}R(\theta_{j},\phi_{j},0)M)R(\theta_{0},\phi_{0},\lambda).$$

Full text

Key words. quantum signal processing, quantum singular value transformation, least squares, quantum algorithm.

1 Introduction

QSP is one of the most efficient and versatile methods for approximate matrix functions in recent years, requiring only a few ancillary qubits. QSP and its extensions have become important tools in quantum computing, with applications including ground-state preparation [YLY2022Ground, DL2024MultilevelQSP, MGB2025Coherent], solving linear systems [LT20filtering, TWY2024QCG], Hamiltonian simulation [LC2017Uniform, LC2017Optimal, LC19Hamiltonian], and block encoding [CLV2024Explicit, LNY2023Efficient, SCC2024Block-encoding]. It was originally proposed by Low et al. in 2016 to design composite quantum gates that can implement quantum response functions [LYC16Methodology]. These functions are polynomials with definite parity and satisfy the condition that magnitudes do not exceed $1$. In 2019, Low and Chuang presented a technique known as qubitization [LC19Hamiltonian]. This method embeds a Hamiltonian $H$ in an invariant subspace and can compute operator functions of $H$ with optimal query complexity. Utilizing qubitization and singular value decomposition (SVD), QSVT [GSL19QSVT] further extends QSP to embed the polynomial of general operators into a unitary. Subsequently, GQSP was proposed in 2024 by Motlagh and Wiebe [MW24GQSP]. It removes the restrictions on the parity of implementable polynomials in QSP, retaining only the condition that the polynomial $p$ satisfies $|p|\leq 1$. This relaxation further broadens the applicability of QSP.

However, GQSP only constructs the function of unitary $U=e^{iH}$, where $H$ is a Hermitian matrix. In this paper, we generalize GQSP from unitary matrices to general matrices by combining singular value decomposition and qubitization techniques. We adopt the framework of GQSP, using arbitrary $SU(2)$ rotations as the signal processing operator and controlled unitary matrices as the signal operator. We use two classes of [math]-controlled and $1$-controlled unitary operators to implement transformations between different spaces. Using the Chebyshev expansion of arbitrary polynomials, our method further extends the range of realizable block encoding matrix functions, that is, the block encoding of polynomials without definite parity of general operators.

QSP and its series of extensions have provided new and powerful approaches for solving the quantum linear system problem (QLSP). Solving linear systems is among the most frequent problems in scientific computing. The classical methods for solving linear systems are classified as either direct or iterative [Vorst2003, Golub2013Matrix, Meurant2020Nonsymmetric, Saad2003, Barrett1994Templates]. However, as the scale of the problem gradually increases, the computational resources and time required by classical computers to solve such linear systems increase greatly. Quantum computers hold the promise of enabling new algorithms that can solve problems requiring excessive resources on classical computers [Xiang2022, Nielsen2010, Shor97Polynomial, HHL09, Grover96fast].

The QLSP refers to the preparation of a quantum state $\ket{\bm{x}}=A^{-1}\ket{\bm{b}}$, which is proportional to the solution of a given linear system $A\bm{x}=\bm{b}$. Harrow, Hassidim, and Lloyd presented the HHL algorithm [HHL09] in 2009, which was the first quantum algorithm to solve linear systems. The computational complexity of the HHL algorithm is $O(\kappa^{2}\log(n)/\varepsilon)$ and scales logarithmically with the system size $n$, where $\kappa=\|A\|\|A^{-1}\|$ is the condition number of $A$, $\|\cdot\|$ is the usual $\ell^{2}$ matrix norm, and $\varepsilon$ is the desired precision. In subsequent studies, the HHL algorithm has been further improved [AA10VTAA, CKS17LCU, WZP18Dense]. The development of QSP has introduced new methods for solving QLSP. Gilyén et al. found the polynomial that approximates $1/x$ and employed QSVT to construct an operator that approximates $A^{-1}$ with maximum circuit depth $O(\kappa\log(\kappa/\varepsilon))$[GSL19QSVT, MRT21Grand]. Lin and Tong presented a quantum eigenstate filtering algorithm and used it to solve quantum linear system problems, achieving near-optimal $O(\kappa\log(1/\varepsilon))$ complexity for both $\kappa$ and $\varepsilon$ [LT20filtering]. In 2024, Toyoizumi et al. introduced a quantum conjugate gradient method that achieves a square-root improvement for $\kappa$ in the maximum circuit depth[TWY2024QCG].

As an application of GQSVT, we use it to solve the least squares (LS) problem, a fundamental and widely used approach for data fitting[Trefethen97Numerical]. For a matrix $A\in\mathbb{R}^{m\times n}$ and a vector $\bm{b}\in\mathbb{R}^{n}$ $(m\geq n)$, the LS procedure find the vector $\bm{x}$ that minimizes $\|A\bm{x}-\bm{b}\|$. The vector $\bm{x}$ is the LS solution if and only if $A^{\top}A\bm{x}=A^{\top}\bm{b}$ when $A$ has full column rank. This equation can be solved by the conjugate gradient (CG) method, which leads to the CGLS algorithm. The CG method works for symmetric positive (SPD) definite problems and it is an iterative Krylov subspace method, where the Krylov subspace with dimension $j$ is defined as $\mathcal{K}_{j}(A,\bm{v})=\operatorname{span}\{\bm{v},A\bm{v},A^{2}\bm{v},\ldots,A^{j-1}\bm{v}\}$.

This paper is organized as follows. In Sect.2, we present several fundamental definitions and a brief review of the GQSP method. In Sect.3, we provide a detailed process for extending GQSP to GQSVT. We also discuss its circuit implementation. In Sect.4, we first introduce the LS problem and the CGLS algorithm, then we present the quantum CGLS algorithm by employing GQSVT to compute vectors and using the swap test to calculate inner products. By using the convergence properties of CG, we derive the maximum circuit depth for our proposed algorithm. Conclusions are given in Sect.5.

2 Preliminaries

2.1 Functions of matrices

We first review the definitions of standard matrix functions and generalized matrix functions.

Standard matrix functions. For simplicity, let $A$ be an $n\times n$ Hermitian matrix with the eigendecomposition $A=U\Lambda U^{\dagger}$, where $U$ is unitary and $\Lambda=\operatorname{diag}(\lambda_{1},\lambda_{2},...,\lambda_{n})$ is a diagonal matrix containing the eigenvalues of $A$. For a function $f$, the standard matrix function $f(A)$ is defined as

[TABLE]

where $f(\Lambda)=\mathrm{diag}(f(\lambda_{1}),f(\lambda_{2}),...,f(\lambda_{n}))$.

Generalized matrix functions. Let $A$ be an $m\times n$ matrix with rank $r=\operatorname{rank}(A)$. The singular value decomposition of $A$ is given by $A=W\Sigma V^{{\dagger}}$, where $W$ and $V$ are unitary. Here, $\Sigma=\mathrm{diag}(\sigma_{1},...,\sigma_{r},0,...,0)$ is an $m\times n$ diagonal matrix whose diagonal entries $\sigma_{k}$ are the non-negative singular values of $A$. For a function $f$, the generalized matrix functions $f^{\diamond}(A)$ is defined as

[TABLE]

where $f(\Sigma)=\operatorname{diag}(f(\sigma_{1}),...,f(\sigma_{r}),0,...,0)$. The left and right generalized matrix polynomials are defined as

[TABLE]

respectively[HB73generalizedMF, BG2003Generalized].

In general, $f^{\diamond}(A)\neq f(A)$ for square indefinite matrices $A$, yet certain relations exist between $f^{\diamond}(A)$ and $f(A)$ [Aurentz19GMF]. Let $q(x)$ be a polynomial, and define $p(x)=q(x^{2})x$. We have

[TABLE]

Recall that when solving the normal equations using CG method, the approximate solution vector $\bm{x}_{j}$ lies in the Krylov subspace $\bm{x}_{0}+\mathcal{K}_{j}(A^{\dagger}A,A^{\dagger}\bm{b})$, which means that $\bm{x}_{j}$ can be represented as the generalized matrix polynomial $q(A^{\dagger}A)A^{\dagger}=p^{\diamond}(A^{\dagger})$ applied to vector $\bm{b}$ when choosing $\bm{x}_{0}=0$.

2.2 Chebyshev expansions

The Chebyshev polynomial $T_{j}(x)$ of the first kind is a polynomial in $x$ of degree $j$, defined by the relation

[TABLE]

where $\eta=\arccos(x)$, $x\in[-1,1]$. The Chebyshev polynomials are orthogonal in the sense that

[TABLE]

Then, a polynomial $f(x)$ of degree $d$ can be expanded in terms of Chebyshev polynomials as

[TABLE]

where the expansion coefficients are determined by

[TABLE]

We apply the substitution $\eta=\arccos(x)$ and re-express the Chebyshev expansion as follows:

[TABLE]

where $c_{j}=a_{j}/2$ for $j\geq 0$ and $c_{j}=a_{-j}/2$ for $j<0$. We define $\hat{f}_{2d}(e^{i\eta})=\sum\limits_{j=-d}^{d}c_{j}e^{i(d+j)\eta}$. Then we have $\hat{f}_{2d}(e^{i\eta})=e^{id\eta}\tilde{f}_{d}(e^{i\eta})$.

2.3 Generalized quantum signal processing

Let $H$ be a Hamiltonian. The signal operator used in generalized quantum signal processing is a [math]-controlled application of $U=e^{iH}$:

[TABLE]

The signal processing operations are defined as arbitrary $SU(2)$ rotations of the ancillary qubit:

[TABLE]

For the sake of brevity, we denote

[TABLE]

where

[TABLE]

The conclusion in Ref.[MW24GQSP] shows that the polynomial transformations of the unitary matrix $U$ can be block encoded by alternately performing the signal operator and signal processing operator. The quantum circuit is shown in Fig. LABEL:fig:GQSP. Importantly, these polynomials require no additional assumptions after being appropriately scaled, as shown in the following lemma.

Lemma 1**.**

(Generalized quantum signal processing, GQSP[MW24GQSP]) $\forall P\in\mathbb{C}[x]$, with $deg(P)=d$, if $\forall x\in\mathbb{R}$, $|P(e^{ix})|^{2}\leq 1$. Then $\exists\bm{\theta},\bm{\phi}\in\mathbb{R}^{d+1}$, $\lambda\in\mathbb{R}$ such that

[TABLE]