1 Introduction↩︎

1.1 Quantum linear system solvers↩︎

Solving systems of linear equations constitutes a foundational task across scientific and engineering disciplines. Classical solvers quickly become computationally prohibitive as the number of constraints and unknown variables increases. This intractability led Harrow, Hassidim, and Lloyd to develop the first quantum linear system algorithm in the seminal work [1], where the BQP-completeness of matrix inversion was established. Remarkably, the query complexity of quantum linear system solvers does not depend explicitly on the problem dimensionality, enabling a potential exponential speedup over classical solvers. It remains challenging to identify end-to-end applications that fully realize this exponential quantum speedup; however, recent studies point to physics-inspired problems such as estimating many-body Green’s functions [2] as promising avenues.

Beyond direct applications, quantum linear system solvers also serve as a key primitive for constructing more advanced quantum algorithms. This was demonstrated in previous work on solving differential equations [3]–[5], where a system of linear equations was introduced whose inverse encodes the matrix exponential via a truncated Taylor series. More generally, quantum linear system solvers have found applications in transforming eigenvalues of non-normal operators, achieved through rational generating functions for Faber polynomials [6], contour integration of resolvents [7]–[9], and rational approximation [10].

A critical parameter determining the runtime of quantum linear system solvers is the spectral condition number. Specifically, consider the linear system \(Ax=|b\rangle\), where \(A\) is the coefficient matrix with spectral norm \(\left\lVert A\right\rVert\leq1\) accessed through a block encoding oracle \(O_A\), and \(|b\rangle\) is the initial vector/state with unit Euclidean norm \(\left\lVert|b\rangle\right\rVert=1\) prepared by a unitary \(O_b\) from a standard reference state. If \(A\) is invertible, then the solution to the linear system is given by \[x=A^{-1}|b\rangle,\] and we let the condition number be the largest norm of the solution vector maximized over all possible initial states \[\kappa=\max_{\left\lVert|b\rangle\right\rVert=1}\left\lVert A^{-1}|b\rangle\right\rVert=\left\lVert A^{-1}\right\rVert.\] In practice, the condition number is rarely known exactly; we instead work with an upper bound \(\kappa \geq \left\lVert A^{-1}\right\rVert\). To produce the normalized state \(|x\rangle=\frac{x}{\left\lVert x\right\rVert}\) with accuracy \(\epsilon\), the Harrow–Hassidim–Lloyd solver has query complexity \(\operatorname{\mathbf{O}}\left(\mathop{\mathrm{poly}}\left(\kappa,\frac{1}{\epsilon}\right)\right)\). This was improved through a series of works, ultimately yielding solvers with query complexity \(\operatorname{\mathbf{O}}\left(\kappa\log\left(\frac{1}{\epsilon}\right)\right)\) [11]–[14]. This tradeoff is optimal, which is established by an unpublished lower bound of Harrow and Kothari and confirmed recently in [15] and [16].

All known solvers attaining the optimal scaling treat the coefficient matrix and the initial vector on equal footing, and thus make the same asymptotic number of queries to \(O_A\) and \(O_b\): \[\operatorname{\mathbf{O}}\left(\kappa\log\left(\frac{1}{\epsilon}\right)\operatorname{\mathbf{Cost}}\left(O_A\right) +\kappa\log\left(\frac{1}{\epsilon}\right)\operatorname{\mathbf{Cost}}\left(O_b\right)\right).\] Alternatively, one may solve the quantum linear system problem using Variable Time Amplitude Amplification (VTAA) [17] incurring a cost of [14] \[\operatorname{\mathbf{O}}\left(\kappa\log\left(\frac{\kappa}{\left\lVert x\right\rVert}\right) \left(\log\log\left(\frac{\kappa}{\left\lVert x\right\rVert}\right)+\log\left(\frac{1}{\epsilon}\right)\right)\operatorname{\mathbf{Cost}}\left(O_A\right) +\frac{\kappa}{\left\lVert x\right\rVert}\operatorname{\mathbf{Cost}}\left(O_b\right)\right).\] Here, the solution norm satisfies \(1\leq\left\lVert x\right\rVert\leq\kappa\), so the scaling for the state preparation cost \(\frac{\kappa}{\left\lVert x\right\rVert}\) can be much better than the worst-case scaling \(\kappa\log\left(\frac{1}{\epsilon}\right)\) in terms of the condition number. Moreover, this cost of \(O_b\) can be shown to be optimal through a reduction from the quantum search problem. That said, the number of queries to \(O_A\) is still constrained by the condition number \(\kappa\)—the block encoding cost is actually logarithmic factors away from the optimal tradeoff and, even with the deterministic amplification schedule [14], VTAA-based solvers are still considerably more complex to implement than those from [11]–[14].

It is therefore natural to ask whether the condition number barrier can be lifted, especially with respect to the block encoding cost. This question was raised by Harrow, Hassidim, and Lloyd in [1] and answered affirmatively in the recent work of Li [18]. Li’s original solver handles only sparse linear systems, though it is fairly straightforward to generalize to the block encoding model—we sketch this generalization in Appendix [[append:li]](#append:li){reference-type=“ref” reference=“append:li”}. The generalized solver then has the query cost \[\operatorname{\mathbf{O}}\left(\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert\epsilon^2}\operatorname{\mathbf{Cost}}\left(O_A\right) +\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert\epsilon^2}\operatorname{\mathbf{Cost}}\left(O_b\right)\right).\] Achieving the above scaling requires the nontrivial assumption that the solution norm \(\left\lVert x\right\rVert\) is known up to a constant multiplicative factor. If \(\left\lVert x\right\rVert\) is unknown, then the complexity worsens to \[\operatorname{\mathbf{O}}\left(\frac{\left\lVert A^{-1\dagger}x\right\rVert\left\lVert x\right\rVert^2}{\epsilon^2}\operatorname{\mathbf{Cost}}\left(O_A\right) +\frac{\left\lVert A^{-1\dagger}x\right\rVert\left\lVert x\right\rVert^2}{\epsilon^2}\operatorname{\mathbf{Cost}}\left(O_b\right)\right).\] Regardless, the cost now depends on vector norms as opposed to the condition number \(\kappa\), satisfying the inequality \(\left\lVert x\right\rVert^2\leq\left\lVert A^{-1\dagger}x\right\rVert\leq\kappa\left\lVert x\right\rVert\). We shall henceforth call these the beyond-\(\kappa\) solvers.

Li demonstrated that for applications such as solving multivariate polynomial systems, these vector norm inequalities are far from saturated, and beyond-\(\kappa\) solvers remain efficient even in the presence of exponentially large condition numbers. This opens the door to solving ill-conditioned linear systems on quantum computers, a task out of reach for all prior approaches. However, at the opposite extreme where \(\left\lVert A^{-1\dagger}x\right\rVert\approx\kappa\left\lVert x\right\rVert\approx\kappa^2\), the query complexity of Li’s solver is \(\operatorname{\mathbf{O}}\left(\frac{\kappa}{\epsilon^2}\right)\) with a known \(\left\lVert x\right\rVert\) and degrades to \(\operatorname{\mathbf{O}}\left(\frac{\kappa^4}{\epsilon^2}\right)\) otherwise, falling far short of optimality. Even in the intermediate regime, these unfavorable exponents could offset any practical benefit of the beyond-\(\kappa\) solvers, undermining their intended advantage over conventional methods.

1.2 Main result↩︎

We present two quantum linear system solvers whose complexity goes beyond the condition number barrier: one based on an effective truncation of the linear system, and one based on eigenstate filtering over the unit circle with an effective gap. We summarize our result and compare it against previous results in Table 1.

Our truncation-based solver has query complexity \[\label{eq:vtaa95trunc95intro} \operatorname{\mathbf{O}}\left(\kappa_{\text{eff}}\log\left(\frac{\kappa_{\text{eff}}}{\left\lVert x\right\rVert}\right) \left(\log\log\left(\frac{\kappa_{\text{eff}}}{\left\lVert x\right\rVert}\right)+\log\left(\frac{1}{\epsilon}\right)\right)\operatorname{\mathbf{Cost}}\left(O_A\right) +\frac{\kappa_{\text{eff}}}{\left\lVert x\right\rVert}\operatorname{\mathbf{Cost}}\left(O_b\right)\right).\tag{1}\] Here, the effective condition number \(1\leq\kappa_{\text{eff}}\leq\kappa\) is chosen such that the truncated solvers still produce the normalized solution state with the target accuracy \(\epsilon\). We show that the following family of bounds hold for the effective condition number \[\label{eq:keff95t95intro} \kappa_{\text{eff}}\leq\left(\frac{\left\lVert(A)_{\mathbf{sv}}^{-t}x\right\rVert}{\left\lVert x\right\rVert\epsilon}\right)^{\frac{1}{t}},\qquad 0<t<\infty,\tag{2}\] where \((A)_{\mathbf{sv}}^{-t}\) denotes the singular value transformation [25] with scalar function \(\sigma\mapsto\sigma^{-t}\). For a positive integer \(t\), the singular value transformation is given by the standard matrix inversion as \[\label{eq:keff95integer95intro} \kappa_{\text{eff}}\leq\begin{cases} \left(\frac{\left\lVert A^{-1\dagger}(A^\dagger A)^{-(t-1)/2}x\right\rVert}{\left\lVert x\right\rVert\epsilon}\right)^{1/t},\quad&t\text{ is odd},\\ \left(\frac{\left\lVert(A^\dagger A)^{-t/2}x\right\rVert}{\left\lVert x\right\rVert\epsilon}\right)^{1/t},\quad&t\text{ is even}. \end{cases}\tag{3}\] This reduces further to the succinct bound \(\kappa_{\text{eff}}\leq\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert\epsilon}\) when \(t=1\), but it also gives alternative bounds with logarithmic scaling in the inverse accuracy when \(t=\operatorname{\mathbf{\Theta}}\left(\frac{\log\left(\frac{1}{\epsilon}\right)}{\log\log\left(\frac{1}{\epsilon}\right)}\right)\). As \(t\) increases, the vector norm \((\left\lVert(A)_{\mathbf{sv}}^{-t}x\right\rVert/\left\lVert x\right\rVert)^{1/t}\nearrow\kappa\) increases monotonically to the condition number while the inverse error \((1/\epsilon)^{1/t}\searrow1\) decreases to unity, recovering the scaling of conventional VTAA.

Our filtering-based solver is remarkably simple and has a favorable runtime prefactor. In particular, we show that the solver has query complexity \[\label{eq:filter95intro} 6\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert\epsilon}\ln\left(\frac{1}{\epsilon}\right)\operatorname{\mathbf{Cost}}\left(O_A\right) +6\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert\epsilon}\ln\left(\frac{1}{\epsilon}\right)\operatorname{\mathbf{Cost}}\left(O_b\right)\tag{4}\] to leading order. This is achieved assuming that relative approximations of the solution norm \(\left\lVert x\right\rVert\) are available with an approximation ratio approaching \(1\). If \(\left\lVert x\right\rVert\) is unknown, then we present a comparably simple solution norm estimator achieving a constant relative precision with query cost \[\operatorname{\mathbf{O}}\left(\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert}\log^2\left(\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert}\right)\operatorname{\mathbf{Cost}}\left(O_A\right) +\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert}\log^2\left(\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert}\right)\operatorname{\mathbf{Cost}}\left(O_b\right)\right).\] The approximation ratio can then be refined using standard methods such as amplitude estimation.

Our beyond-\(\kappa\) solvers thus substantially improve upon the prior art, while our vector norm scalings degrade to almost linear in \(\kappa\) in the worst case, reproducing the performance of conventional solvers. This constitutes a concrete step toward efficiently solving ill-conditioned linear systems on a quantum computer.

1.3 Block encoding and affine dilation↩︎

Quantum linear system solvers can be formulated in an abstract setting where the coefficient matrix \(A\) has norm \(\left\lVert A\right\rVert\leq1\) and is accessed through a block encoding as \[A=G_1^\dagger O_AG_0,\] with unitary \(O_A\) and isometries \(G_0,G_1\), while the initial state \(|b\rangle\) is prepared by a unitary \(O_b\) as \[|b\rangle=O_b|0\rangle.\] Assuming the implementation cost of \(G_0\) and \(G_1\) is negligible without loss of generality, the complexity of these solvers is then measured in terms of the number of queries to \(O_A\) and \(O_b\) respectively. This applies to a range of concrete settings, including sparse linear systems, and extends to cases where \(A\) and \(|b\rangle\) are constructed indirectly via common linear algebraic operations. We term this the standard input model for the remainder of the paper.

The beyond-\(\kappa\) solver of [18] operates on an augmented matrix of the form \(\begin{bmatrix} A & -|b\rangle \end{bmatrix}\). Using oracles \(O_A\) and \(O_b\), one can construct a block encoding of this augmented matrix with normalization factor \(\sqrt{2}\) and thereby reformulate the solver in the standard input model; see Appendix [[append:li]](#append:li){reference-type=“ref” reference=“append:li”} for details. However, as shown in [18], the query complexity can be further refined through a direct block encoding of \(\begin{bmatrix} A & -|b\rangle \end{bmatrix}\), bypassing the standard input model.

We introduce the affine dilation model to facilitate comparison of Li’s solver [18] with other quantum linear system solvers while also exploiting this refined block encoding. Specifically, we consider \(\left[\begin{smallmatrix} A & -|b\rangle\\ 0 & c \end{smallmatrix}\right]\), where the augmented matrix is recovered by setting the scalar \(c=0\). However, if \(c\neq0\), then the affine dilated matrix is invertible and \[\begin{bmatrix} A & -|b\rangle\\ 0 & c \end{bmatrix}^{-1} =\begin{bmatrix} A^{-1} & \frac{1}{c}A^{-1}|b\rangle\\ 0 & \frac{1}{c} \end{bmatrix}.\] Consequently, the solution \(x = A^{-1}|b\rangle\) can be extracted from the top-right entry of the inverse matrix, while \(A\) and \(|b\rangle\) are handled jointly in a single block encoding analogous to [18]. We formally introduce the quantum linear system problem and its input models in Section [[sec:input]](#sec:input){reference-type=“ref” reference=“sec:input”}. While we primarily focus on the standard input model for simplicity, we also analyze the query complexity of beyond-\(\kappa\) solvers in the affine dilation model in Appendix [[append:affine]](#append:affine){reference-type=“ref” reference=“append:affine”}.

1.4 Solver based on effective truncation of linear system↩︎

We now present our first beyond-\(\kappa\) solver based on an effective truncation of the linear system. To describe this truncation, we define the projections \(\Pi_{\text{left},\mathcal{S}}=\sum_{\sigma_j\in\mathcal{S}}|u_j\rangle\!\langle u_j|\) and \(\Pi_{\text{right},\mathcal{S}}=\sum_{\sigma_j\in\mathcal{S}}|v_j\rangle\!\langle v_j|\), for \(A=\sum_j\sigma_j|u_j\rangle\!\langle v_j|\) the singular value decomposition with singular values \(\kappa^{-1}\leq\sigma_j\leq1\) and \(\mathcal{S}\subseteq\mathbb{R}\) a subset of real numbers. Then the goal of our truncation-based solver is to produce a quantum state proportional to the truncated solution \(\Pi_{\text{right},\left[\kappa_{\text{eff}}^{-1},1\right]}A^{-1}|b\rangle\). Here, we aim to choose an effective condition number \(1\leq\kappa_{\text{eff}}=\kappa_{\text{eff}}(\epsilon)\leq\kappa\) as the smallest value satisfying \[\left\lVert\frac{\Pi_{\text{right},\left[\kappa_{\text{eff}}^{-1},1\right]}A^{-1}|b\rangle}{\left\lVert\Pi_{\text{right},\left[\kappa_{\text{eff}}^{-1},1\right]}A^{-1}|b\rangle\right\rVert} -\frac{A^{-1}|b\rangle}{\left\lVert A^{-1}|b\rangle\right\rVert}\right\rVert =\operatorname{\mathbf{O}}(\epsilon),\] so that the solver still outputs the desired solution state up to error \(\operatorname{\mathbf{O}}(\epsilon)\).

We show how to fulfill this truncation requirement by characterizing \(\kappa_{\text{eff}}\) through the generalized inverse (quantile function) of the weighted singular value distribution, where the weights are determined by the squared overlaps of \(A^{-1}|b\rangle\) with the singular vector subspaces: \[\frac{\left\lVert\Pi_{\text{right},\left[0,\kappa_{\text{eff}}^{-1}\right)}A^{-1}|b\rangle\right\rVert^2}{\left\lVert A^{-1}|b\rangle\right\rVert^2}\leq\epsilon^2,\qquad \frac{\left\lVert\Pi_{\text{right},\left[0,\kappa_{\text{eff}}^{-1}\right]}A^{-1}|b\rangle\right\rVert^2}{\left\lVert A^{-1}|b\rangle\right\rVert^2}>\epsilon^2.\] This definition is illustrated in Figure [[fig:quantile_keff]](#fig:quantile_keff){reference-type=“ref” reference=“fig:quantile_keff”}. In practice, determining the exact value of the effective condition number \(\kappa_{\text{eff}}\) can be computationally prohibitive. We address this by developing a family of upper bounds on \(\kappa_{\text{eff}}\) such as [eq:keff_integer_intro], which can be analyzed using standard linear algebraic techniques. See Theorem [1](#thm:keff){reference-type=“ref” reference=“thm:keff”} for details.

Note that although the singular vector projection \(\Pi_{\text{right},\left[\kappa_{\text{eff}}^{-1},1\right]}\) may be realized to precision \(\epsilon\) using Quantum Singular Value Transformation (QSVT) [25] with cost \(\operatorname{\mathbf{O}}\left(\kappa_{\text{eff}}\log\left(\frac{1}{\epsilon}\right)\right)\), we do not pursue this approach, as this would incur query overhead and significantly complicate the algorithm. Instead, we invoke a quantum linear system solver with the condition number set to \(\kappa_{\text{eff}}\), without modifying the input oracles \(O_A\) and \(O_b\). Since our \(\kappa_{\text{eff}}\) is no longer an upper bound on the condition number, this constitutes a somewhat “abnormal” use of the quantum linear system solver, and the resulting behavior is generally unpredictable. To address this, we demand the solver to effectively truncate the coefficient matrix in the following sense.

Under this definition, any quantum linear system solver with the strong truncation property can produce a truncated solution state, with a query complexity determined by the effective condition number \(\kappa_{\text{eff}}\). Recall from [eq:keff_t_intro] and [eq:keff_integer_intro] that \(\kappa_{\text{eff}}\) can be bounded by vector norms that are independent of the condition number; consequently, we have obtained a beyond-\(\kappa\) solver. We formally state this result in Theorem [[thm:beyondk_strong_trunc]](#thm:beyondk_strong_trunc){reference-type=“ref” reference=“thm:beyondk_strong_trunc”}.

It remains to identify solvers that fulfill the strong truncation requirement. One simple solver can be constructed by block encoding \(A^{-1}\) using QSVT with condition number set to \(\kappa_{\text{eff}}\), followed by a (fixed-point) amplitude amplification. This QSVT-based solver makes \(\operatorname{\mathbf{O}}\left(\frac{\kappa_{\text{eff}}^2}{\left\lVert x\right\rVert}\log\left(\frac{\kappa_{\text{eff}}}{\left\lVert x\right\rVert\epsilon}\right)\right)\) queries to \(O_A\) and \(\operatorname{\mathbf{O}}\left(\frac{\kappa_{\text{eff}}}{\left\lVert x\right\rVert}\right)\) queries to \(O_b\). Running VTAA with \(\kappa_{\text{eff}}\) yields an improved query complexity, as shown in [eq:vtaa_trunc_intro]. However, rigorously establishing the strong truncation property for VTAA is highly nontrivial. The difficulty arises because VTAA recursively amplifies all components of the state that potentially lead to success. When configuring with \(\kappa_{\mathrm{eff}}\), we intend to zero out components corresponding to singular values \(\sigma_j < \kappa_{\mathrm{eff}}^{-1}\). However, VTAA still treats these components as potentially successful and amplifies them, requiring an additional postprocessing step. In effect, VTAA amplifies first and then truncates, rather than truncating first and then amplifying. This nuance leads to a nonstandard amplification schedule that is difficult to analyze directly. To circumvent this technical difficulty, we introduce the weak truncation property. Combined with an effective truncation of the initial state, this admits a significantly cleaner proof.

Specifically, with the same effective condition number \(\kappa_{\text{eff}}=\kappa_{\text{eff}}(\epsilon)\) as above, we consider the truncated initial state \(|b_{\text{eff}}\rangle\) proportional to \(\Pi_{\text{left},\left[\kappa_{\text{eff}}^{-1},1\right]}|b\rangle\). We then construct an oracle \(O_{b_{\text{eff}}}\) that prepares \(|b_{\text{eff}}\rangle\) with one query to \(O_b\) and approximates \(O_b\), i.e., \[O_{b_{\text{eff}}}|0\rangle=\frac{\Pi_{\text{left},\left[\kappa_{\text{eff}}^{-1},1\right]}|b\rangle}{\left\lVert\Pi_{\text{left},\left[\kappa_{\text{eff}}^{-1},1\right]}|b\rangle\right\rVert}=|b_{\text{eff}}\rangle,\qquad \left\lVert O_{b_{\text{eff}}}-O_b\right\rVert=\operatorname{\mathbf{O}}\left(\frac{\left\lVert x\right\rVert\epsilon}{\kappa_{\text{eff}}}\right).\] Once again, the singular vector projection is not performed in the actual implementation of the algorithm. Instead, we run a quantum linear system solver configured with condition number \(\kappa_{\text{eff}}\), while retaining the original input oracles \(O_A\) and \(O_b\). Suppose that the solver we use has an optimal query complexity of initial state preparation, i.e., it makes \(\operatorname{\mathbf{O}}\left(\frac{\kappa_{\text{eff}}}{\left\lVert x\right\rVert}\right)\) queries to \(O_b\). Then replacing every occurrence of \(O_b\) by \(O_{b_{\text{eff}}}\) introduces an error at most \(\operatorname{\mathbf{O}}(\epsilon)\). Up to this error tolerance, our solver is conceptually an algorithm with condition number \(\kappa_{\text{eff}}\) and oracles \(O_A\) and \(O_{b_{\text{eff}}}\); in this sense the initial state is effectively truncated.

As before, this is a non-standard use of the quantum linear system solver, and its correct execution requires the following additional assumption on the underlying solver.

Note that the strong truncation property naturally implies the weak truncation property. This follows from the fact that \(A^{-1}|b\rangle = A^{-1}\Pi_{\text{left},[\kappa_{\text{eff}}^{-1},1]}|b\rangle = \Pi_{\text{right},[\kappa_{\text{eff}}^{-1},1]}A^{-1}|b\rangle\) whenever \(|b\rangle \in \operatorname{\mathbf{Im}}\left(\Pi_{\text{left},[\kappa_{\text{eff}}^{-1},1]}\right)\) holds. Conversely, we show in Theorem [[thm:beyondk_weak_trunc]](#thm:beyondk_weak_trunc){reference-type=“ref” reference=“thm:beyondk_weak_trunc”} that a solver satisfying the weak truncation property with optimal query complexity of initial state preparation necessarily fulfills the strong truncation requirement. Given the additional promise that \(|b\rangle \in \operatorname{\mathbf{Im}}(\Pi_{\text{left},[\kappa_{\mathrm{eff}}^{-1},1]})\), it is routine to check that VTAA satisfies the weak truncation property. As VTAA also has an optimal query complexity of initial state preparation, it can be used to construct beyond-\(\kappa\) solvers based on effective state truncation (Corollary [[cor:beyondk_vtaa]](#cor:beyondk_vtaa){reference-type=“ref” reference=“cor:beyondk_vtaa”}). We illustrate this algorithmic template in Figure 2 and present a detailed analysis of the truncation-based solver in Section [[sec:trunc]](#sec:trunc){reference-type=“ref” reference=“sec:trunc”}.

To obtain the beyond-\(\kappa\) solvers, we require that VTAA produces the truncated solution state when configured with the effective condition number \(\kappa_{\text{eff}}\). Interestingly, VTAA can be modified to produce the truncated solution state for any \(\alpha\) supplied as the condition number parameter. This is a more general truncation property of VTAA and is discussed further in Appendix [[append:trunc_general]](#append:trunc_general){reference-type=“ref” reference=“append:trunc_general”}.

1.5 Solver based on filtering with effective gap↩︎

Our second beyond-\(\kappa\) solver is based on an improved eigenstate filtering on the unit circle with an effective gap. Specifically, suppose that the augmented matrix \(\begin{bmatrix} \beta A & -|b\rangle \end{bmatrix}\) is block encoded by the overlap of isometries \(G_1^\dagger G_0\) with some scalar \(\beta>0\). We then start with the initial state \(G_0\begin{bmatrix} 0\\ 1 \end{bmatrix}\) partitioned conformally to the augmented matrix, which decomposes as \[\begin{align} G_0\begin{bmatrix} 0\\1 \end{bmatrix} &=\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2} G_0\begin{bmatrix} x\\\beta \end{bmatrix} -\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2}G_0\begin{bmatrix} x\\ -\frac{\left\lVert x\right\rVert^2}{\beta} \end{bmatrix}\\ &=\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2}G_0\begin{bmatrix} x\\ \beta \end{bmatrix} -\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2} G_0G_0^\dagger G_1A^{-1\dagger}x. \end{align}\] Our goal is to retain the first term while discarding the second, from which a state proportional to the solution \(x=A^{-1}|b\rangle\) can be extracted.

A direct calculation shows that the first term belongs to \(\operatorname{\mathbf{Ker}}\left(\Pi_1 \Pi_0\right)\) while the second term is supported in \(\operatorname{\mathbf{Im}}\left(\Pi_0 \Pi_1\right)\), where \(\Pi_0=G_0G_0^\dagger\) and \(\Pi_1=G_1G_1^\dagger\). As the smallest nonzero singular value of \(\Pi_1\Pi_0\) is at least \(\operatorname{\mathbf{\Omega}}(\kappa^{-1})\) [18], one can perform a right singular vector filtering targeting at the singular value \(0\) with \(\Pi_1\Pi_0\) as the underlying operator. This preserves the first term while arbitrarily suppressing the second to magnitude below \(\operatorname{\mathbf{O}}(\delta)\), although the query cost would be \(\operatorname{\mathbf{O}}\left(\kappa\log\left(\frac{1}{\delta}\right)\right)\) scaling linearly with the condition number.

In contrast, the beyond-\(\kappa\) solver uses the operator \(W=\left(2\Pi_0-I\right)\left(I-2\Pi_1\right)\). A key observation is that the first term from the decomposition is an eigenvector of \(W\) with eigenvalue \(e^{i0}=+1\), while the second term corresponds to eigenvalues \(e^{i\theta}\) with nonzero phase angles \(\theta\neq0\). The effective gap lemma [26] then allows us to reduce the magnitude of the second term below \(\operatorname{\mathbf{O}}(\delta)\) with a query cost \(\operatorname{\mathbf{O}}\left(\frac{1}{\delta}\right)\) independent of the condition number. This suggests performing an eigenstate filtering over the unit circle that preserves the \(\theta=0\) component while suppressing the remaining ones with the help of effective gap lemma.

This eigenstate filtering can be realized in a variety of ways. For instance, one can use generalized quantum signal processing [27], [28] with the target polynomial chosen as the Dolph-Chebyshev polynomial over the unit circle [29]. Prior work [29] also proposed a filtering method based on the conventional quantum signal processing [30], [31] and QSVT. However, their method relies on \(\sin(\theta)\) and appears incapable of distinguishing between eigenphases \(\theta=0\) and \(\theta=\pm\pi\), a critical issue left unresolved. We describe a corrected QSVT-based filtering algorithm with a matching asymptotic complexity and a reduced constant prefactor (Appendix [[append:lambertw]](#append:lambertw){reference-type=“ref” reference=“append:lambertw”}). Regardless, this improves upon the prior result of [18] based on quantum phase estimation and repeated sampling.

If \(\left\lVert x\right\rVert\) is known a priori, we can run the beyond-\(\kappa\) solver with \(\beta=\operatorname{\mathbf{\Theta}}(\left\lVert x\right\rVert)\), successfully producing the solution state with a constant probability. This filter-then-measure strategy is exceptionally simple. In fact, we show that its query complexity has a constant prefactor \(6\) to leading order; see [eq:filter_intro]. If \(\left\lVert x\right\rVert\) is unknown, then we develop a similarly simple solution norm estimation algorithm by applying an exponential search starting from \(\beta=1\) and terminating at \(\beta=\operatorname{\mathbf{\Theta}}(\left\lVert x\right\rVert)\) with a high probability. This yields a constant multiplicative approximation of the solution norm, and the approximation ratio can be further refined using amplitude estimation as in [12]. We present the filtering-based solver in detail in Section [[sec:gap]](#sec:gap){reference-type=“ref” reference=“sec:gap”} with the result formally stated in Theorem [[thm:beyondk_filtering]](#thm:beyondk_filtering){reference-type=“ref” reference=“thm:beyondk_filtering”}.

We remark that although both solvers approximate the solution state \(|x\rangle\) to an accuracy of \(\operatorname{\mathbf{O}}(\epsilon)\), their error profiles are structurally different. The truncation-based solver produces a state approximately supported on the right singular vectors of \(A\) with singular values of at least \(\kappa_{\mathrm{eff}}^{-1}\), discarding all components corresponding to smaller singular values. The filtering-based solver preserves the full solution but incurs errors from incomplete suppression of an auxiliary component in a constrained orthogonal decomposition induced by the operator \(W\). Despite this structural difference, the two solvers provide equivalent approximation guarantees: both produce a quantum state within Euclidean distance \(\operatorname{\mathbf{O}}(\epsilon)\) of the normalized \(|x\rangle\). Which solver to use then depends on the regime of interest: the truncation-based solver (via VTAA) achieves nearly linear scaling in the effective condition number with an optimal number of queries to the initial state preparation, while the filtering-based solver (with the effective spectral gap) provides a much simpler approach with a more favorable runtime prefactor.

We conclude in Section [[sec:discuss]](#sec:discuss){reference-type=“ref” reference=“sec:discuss”} with a summary of our main contributions and a number of open questions for future investigation.

2 Quantum linear system problem and input models↩︎

In this section, we formally introduce the quantum linear system problem and discuss its input models. We first define the standard input model in Section [[sec:input_standard]](#sec:input_standard){reference-type=“ref” reference=“sec:input_standard”} based on block encoding, and then introduce the affine dilation model in Section [[sec:input_affine]](#sec:input_affine){reference-type=“ref” reference=“sec:input_affine”} which allows further fine-tuning of the query complexity. Finally, we formulate the quantum linear system problem in Section [[sec:input_problem]](#sec:input_problem){reference-type=“ref” reference=“sec:input_problem”} and discuss its solvers.

2.1 Standard input model↩︎

Block encoding provides a natural input model for quantum linear algebra algorithms [21], [32]. It can be defined in full generality using unitaries and isometries. We say an operator \(G:\mathcal{G}\rightarrow\mathcal{H}\) is an isometry between spaces \(\mathcal{G}\) and \(\mathcal{H}\) if \(G^\dagger G=I\). By definition, \(G\) is injective with image \(\operatorname{\mathbf{Im}}(G)=\operatorname{\mathbf{Im}}(GG^\dagger)\), whereas \(G^\dagger\) is surjective with kernel \(\operatorname{\mathbf{Ker}}(G^\dagger)=\operatorname{\mathbf{Ker}}(GG^\dagger)\). We thus obtain the isometric embedding \(\mathcal{G} \xrightleftharpoons[G^\dagger]{G} \operatorname{\mathbf{Im}}(GG^\dagger)\subseteq\mathcal{H}\). Fixing any bases with respect to this embedding, the operator \(G\) can be represented in matrix form as \[G=\begin{bmatrix} I\\ 0 \end{bmatrix}.\] Examples of isometries include: (i) quantum state \(|\psi\rangle\), where \(\left\lVert|\psi\rangle\right\rVert=1\); (ii) unitary operator \(U\), where \(UU^\dagger=U^\dagger U=I\); (iii) tensor product \(G_0\otimes G_1\), where \(G_0\) and \(G_1\) are both isometries; (iv) operator composition \(G_1G_0\), where \(G_0\) and \(G_1\) are both isometries and the composition is well defined; and (v) Stinespring dilation \(\sum_j|j\rangle\otimes E_j\), where \(\sum_jE_j^\dagger E_j=I\) and \(\{|j\rangle\}\) are orthonormal.

We say that \(A\) is block encoded by isometries \(G_0\), \(G_1\) and unitary \(O_A\), if \[A=G_1^\dagger O_AG_0.\] In bases compatible with the isometric embeddings, we can write \(O_A\) in the matrix form \[O_A=\begin{bmatrix} A & \cdot\\ \cdot & \cdot \end{bmatrix},\] where \(A\) appears as the top-left block; hence the name block encoding. Note that such a block encoding is feasible only when \(\left\lVert A\right\rVert\leq\left\lVert G_1^\dagger\right\rVert\left\lVert O_A\right\rVert\left\lVert G_0\right\rVert=1\). Conversely, for any \(\left\lVert A\right\rVert\leq1\), one can mathematically define a dilation \(O_A=\left[\begin{smallmatrix} A & -\sqrt{I-AA^\dagger}\\ \sqrt{I-A^\dagger A} & A^\dagger \end{smallmatrix}\right]\) as in [33] and [34], although \(A\) may need to be further rescaled by some normalization factor for \(O_A\) to be efficiently realizable as a quantum circuit.

For a given linear system \(Ax=|b\rangle\), we assume without loss of generality that the initial vector \(|b\rangle\) is a normalized quantum state, prepared by oracle \(O_b\) from a standard reference state as \[|b\rangle=O_b|0\rangle.\] We further assume that the coefficient matrix \(A\) is accessed through a block encoding, with any normalization factor absorbed into the matrix so that \(\left\lVert A\right\rVert \leq 1\). Depending on the quantum linear system solver under consideration, we adopt one of two simplified forms of block encoding. In the first form, we neglect the implementation cost of the isometries \(G_0\), \(G_1\) and simply say that \(A\) is block encoded by the unitary \(O_A\). This can happen for a block encoding like \(A=\left(\langle 0|\otimes I\right)O_A\left(|0\rangle\otimes I\right)\) where isometries are given by trivial state preparations. In the second form, we absorb \(O_A\) into the definition of either \(G_0\) or \(G_1\) and say that \(A\) is block encoded by the overlap of isometries \(G_1^\dagger G_0\). The query complexity is then measured by the number of calls to either \(G_0\) or \(G_1\).

2.2 Affine dilation model↩︎

In the standard input model, we assumed that the coefficient matrix \(A\) has spectral norm \(\left\lVert A\right\rVert \leq 1\) and is accessed through a block encoding \(A = G_1^\dagger O_A G_0\). However, this often cannot be efficiently realized by a quantum circuit unless \(A\) is rescaled by an additional factor \(\alpha > 1\), giving \(\frac{A}{\alpha} = G_1^\dagger O_A G_0\). The normalization factor \(\alpha>1\) pushes the singular values of \(A\) toward \(0\), making the linear system harder to solve quantumly; minimizing \(\alpha\) is thus critical for reducing the overall query complexity.

To this end, we introduce the affine dilation model, where the input is given by the block encoding of \[\widetilde{A}=\frac{1}{\widetilde{\alpha}} \begin{bmatrix} A & -|b\rangle\\ 0 & c \end{bmatrix}\] for some scalar \(c\) and normalization factor \(\widetilde{\alpha}\geq1\). If \(c=0\), then the affine dilation reduces to the augmented matrix \(\frac{1}{\alpha}\begin{bmatrix} A & -|b\rangle \end{bmatrix}\) used as the input model by the beyond-\(\kappa\) solver of [18]. However, we can also choose \(c>0\). Then the affine dilated matrix becomes invertible and the solution \(x = A^{-1}|b\rangle\) can be extracted from the top-right entry of the inverse matrix \[\widetilde{A}^{-1} =\left(\frac{1}{\widetilde{\alpha}}\begin{bmatrix} A & -|b\rangle\\ 0 & c \end{bmatrix}\right)^{-1} =\widetilde{\alpha}\begin{bmatrix} A^{-1} & \frac{1}{c}A^{-1}|b\rangle\\ 0 & \frac{1}{c} \end{bmatrix}.\] This means the affine dilation model can be paired with any quantum linear system solver, allowing the query complexity to be optimized over different choices of \(A\), \(|b\rangle\), \(c\), and \(\widetilde{\alpha}\).

Let us show that, even with the worst-case estimate, it never incurs significant overhead converting from the standard input model to the affine dilation model. To this end, we set \(c=1\) and decompose the affine dilation into \[\begin{bmatrix} A & -|b\rangle\\ 0 & 1 \end{bmatrix} =\begin{bmatrix} A & 0\\ 0 & 1 \end{bmatrix} +\frac{1}{2}\begin{bmatrix} 0 & -|b\rangle\\ -\langle b| & 0 \end{bmatrix} +\frac{1}{2}\begin{bmatrix} 0 & -|b\rangle\\ \langle b| & 0 \end{bmatrix},\] where each summand can be block encoded respectively as \[\begin{align} \begin{bmatrix} A & 0\\ 0 & 1 \end{bmatrix} &=\begin{bmatrix} G_1^\dagger & 0\\ 0 & \langle 0| \end{bmatrix} \begin{bmatrix} O_A & 0\\ 0 & I \end{bmatrix} \begin{bmatrix} G_0 & 0\\ 0 & |0\rangle \end{bmatrix},\\ \begin{bmatrix} 0 & -|b\rangle\\ -\langle b| & 0 \end{bmatrix} &=\begin{bmatrix} I & 0\\ 0 & \langle 0| \end{bmatrix} \begin{bmatrix} 0 & -O_b\\ -O_b^\dagger & 0 \end{bmatrix} \begin{bmatrix} I & 0\\ 0 & |0\rangle \end{bmatrix},\\ \begin{bmatrix} 0 & -|b\rangle\\ \langle b| & 0 \end{bmatrix} &=\begin{bmatrix} I & 0\\ 0 & \langle 0| \end{bmatrix} \begin{bmatrix} 0 & -O_b\\ O_b^\dagger & 0 \end{bmatrix} \begin{bmatrix} I & 0\\ 0 & |0\rangle \end{bmatrix}. \end{align}\] Taking their linear combination [25] then yields a block encoding of the affine dilation matrix \(\widetilde{A}=\frac{1}{2}\left[\begin{smallmatrix} A & -|b\rangle\\ 0 & 1 \end{smallmatrix}\right]\) with normalization factor \(\widetilde{\alpha} = 2\). This implies that the condition number increases to at most \[\left\lVert\widetilde{A}^{-1}\right\rVert = \left\lVert 2\begin{bmatrix} A^{-1} & A^{-1}|b\rangle\\ 0 & 1 \end{bmatrix}\right\rVert \leq\left\lVert 2\begin{bmatrix} A^{-1} & 0\\ 0 & 1 \end{bmatrix}\right\rVert +\left\lVert 2\begin{bmatrix} 0 & A^{-1}|b\rangle\\ 0 & 0 \end{bmatrix}\right\rVert \leq4\kappa,\] while the norm of the solution encoded into the upper-right block rescales to \(\left\lVert\widetilde{\alpha}A^{-1}|b\rangle\right\rVert=2\left\lVert x\right\rVert\). Since the cost of conventional quantum linear system solvers is governed by the condition number and the solution norm, switching from the standard input model to the affine dilation model thus increases the runtime by at most a constant factor. An exception to this conclusion may occur if \(\operatorname{\mathbf{Cost}}\left(O_b\right)\) is asymptotically larger than \(\operatorname{\mathbf{Cost}}\left(O_A\right)\). In this case, quantum linear system solvers that make asymptotically fewer calls to the initial vector than to the coefficient matrix may have asymptotically better cost in the standard input model than in the affine dilation model, since block-encoding the dilated coefficient matrix \(\widetilde{A}\) via linear combination requires querying not only \(O_A\) (cheap) but also \(O_b\) (expensive). As shown in Table 1, this category of quantum linear system solver includes the VTAA method, but excludes all methods with strictly optimal query complexity to \(O_A\).

Li [18] showed that for sparse linear systems, one can adjust the row norms in the affine dilation model, substantially reducing the query complexity. In general, searching over \(A\), \(|b\rangle\), \(c\), and \(\widetilde{\alpha}\) to minimize the query complexity is a nontrivial task. Accordingly, we primarily state our results in the standard input model, though we also analyze the runtime in the affine dilation model in Appendix [[append:affine]](#append:affine){reference-type=“ref” reference=“append:affine”}. A complete exploration of the affine dilation model—and the potential improvements it may yield—is an interesting direction for future work.

2.3 Quantum linear system problem↩︎

We now introduce notation and terminology for the quantum linear system problem. Specifically, consider the linear system \(Ax = |b\rangle\) which has the solution \[x=A^{-1}|b\rangle\] for an invertible square matrix \(A\). Then the goal of the quantum linear system problem is to produce the normalized solution vector \[|x\rangle=\frac{x}{\left\lVert x\right\rVert}=\frac{A^{-1}|b\rangle}{\left\lVert A^{-1}|b\rangle\right\rVert}\] as a quantum state. Since both \(A^{-1}\) and \(|b\rangle\) appear linearly in the numerator and denominator of the solution state, we assume without loss of generality that the coefficient matrix satisfies \(\left\lVert A\right\rVert \leq 1\) and the initial vector has unit Euclidean norm \(\left\lVert|b\rangle\right\rVert=1\). We shall use the Dirac notation only for unit vectors.

We consider the abstract setting where a quantum linear system solver is described by a unitary \(U\), acting jointly on the ancilla and system register, with the behavior \[\left\lVert\left(\langle 0|\otimes I\right)U|00\rangle\right\rVert^2>\frac{1}{2},\qquad \left\lVert\frac{\left(\langle 0|\otimes I\right)U|00\rangle}{\left\lVert\left(\langle 0|\otimes I\right)U|00\rangle\right\rVert}-\frac{A^{-1}|b\rangle}{\left\lVert A^{-1}|b\rangle\right\rVert}\right\rVert= \operatorname{\mathbf{O}}(\epsilon). \label{eq:QLSP95formal95specs}\tag{5}\] Here, the first condition guarantees that the solver succeeds with probability strictly greater than \(\frac{1}{2}\), while the second ensures that the normalized output state approximates the desired solution state to within a fixed constant factor of the target accuracy \(\epsilon>0\). We focus primarily on the standard input model where \(A\) is block encoded by an oracle \(O_A\) and \(|b\rangle\) is prepared by a unitary \(O_b\). The complexity of the solver is then measured in terms of the number of queries to \(O_A\) and \(O_b\) made by \(U\). When making non-asymptotic comparisons between two solvers that both produce solutions to \(\operatorname{\mathbf{O}}(\epsilon)\) error as in [eq:QLSP_formal_specs]—but with different prefactors in front of the \(\epsilon\)—care must be taken to rescale \(\epsilon\) accordingly.

As discussed in Section [[sec:intro_solver]](#sec:intro_solver){reference-type=“ref” reference=“sec:intro_solver”}, the spectral condition number is a key parameter determining the runtime of quantum linear system solvers, which is mathematically defined by \(\left\lVert A\right\rVert\left\lVert A^{-1}\right\rVert\). However, the precise value of the condition number is often inaccessible in practice; instead, we work with a known upper bound \[\kappa\geq\left\lVert A^{-1}\right\rVert =\max_{\left\lVert|b\rangle\right\rVert=1}\left\lVert A^{-1}|b\rangle\right\rVert.\] As this involves a maximization over all initial states \(|b\rangle\), the condition number \(\kappa\) serves as a worst-case cost measure for quantum linear system solvers. The main goal of this paper is to develop beyond-\(\kappa\) solvers whose query complexity does not scale with the condition number.

Another parameter of importance is the solution norm \(\left\lVert x\right\rVert = \left\lVert A^{-1}|b\rangle\right\rVert\). Since \(1 \leq \left\lVert x\right\rVert \leq \left\lVert A^{-1}\right\rVert\), this quantity reduces to the condition number in the worst case, but typically provides a more refined cost measure for quantum linear system solvers. Like the condition number, the exact solution norm is often unavailable in practice. Instead, we assume knowledge of a value \(\alpha_x\) that approximates \(\left\lVert x\right\rVert\) to within a multiplicative factor \(\mu \geq 1\), i.e., \[\frac{1}{\mu} \leq \frac{\alpha_x}{\left\lVert x\right\rVert} \leq \mu.\] We will often assume that the solution norm is known to constant multiplicative accuracy, i.e., \(\mu = \operatorname{\mathbf{O}}(1)\). When such an estimate is not available a priori, we provide solution norm estimation algorithms to obtain one—this estimation only needs to be performed once as a preprocessing step. Furthermore, we may assume without loss of generality that \(\mu \to 1\), as this refinement can be achieved using amplitude estimation incurring multiplicative complexity overhead \(1/(\mu-1)\) [12].

To summarize, a quantum linear system solver is described by a unitary \(U(O_A, O_b, \kappa, \alpha_x, \epsilon)\) depending on the input oracles \(O_A\) and \(O_b\), an upper bound on the condition number \(\kappa\), a constant multiplicative estimate of the solution norm \(\alpha_x\), and the target accuracy \(\epsilon\). With this notation, the query complexity of the solver takes the form \[\operatorname{\mathbf{Cost}}(U(O_A, O_b, \kappa, \alpha_x, \epsilon)) =h_A(\kappa, \alpha_x, \epsilon)\operatorname{\mathbf{Cost}}(O_A) +h_b(\kappa, \alpha_x, \epsilon)\operatorname{\mathbf{Cost}}(O_b),\] where \(h_A\) and \(h_b\) are positive-integer-valued functions denoting the number of calls to the respective oracles. To represent them in closed form, we introduce \(\operatorname{\mathbf{Floor}}(\cdot)\) to denote the largest integer not exceeding a real number and \(\operatorname{\mathbf{Ceil}}(\cdot)\) for the smallest integer no less than a number.

The setting above is somewhat restricted for clarity of exposition, but the framework readily accommodates several generalizations. First, the solver need not be a unitary—more generally, it can be described by a quantum channel, with the error measured in trace distance rather than Euclidean distance. This generalization is especially relevant when the solver uses classical randomness, in which case it naturally implements a mixed unitary channel. Second, we have specified the solver’s performance through its success probability and output accuracy. It is also possible to characterize the solver instead by its expected behavior and expected query complexity, assuming that the solver is repeated until a successful run is obtained. Third, we have focused on the standard input model, where \(A\) and \(|b\rangle\) are accessed through separate oracles. One can also consider the affine dilation model, where \(A\) and \(|b\rangle\) are block encoded jointly by a single oracle, allowing further refinements of the query complexity. Fourth, the solver’s performance may depend on additional parameters beyond \(\kappa\) and \(\left\lVert x\right\rVert\). In particular, we will introduce a family of parameters \(\left\lVert(A)_{\mathbf{sv}}^{-t}|x\rangle\right\rVert\) with \(0<t<\infty\) based on the singular value transformation of \(A\), which upper bounds the query complexity of beyond-\(\kappa\) solvers. Fifth, we have focused on the case where \(A\) is an invertible square matrix, but the framework generalizes to non-square and non-invertible matrices. In particular, when the linear system has no solution—i.e., when \(|b\rangle\notin\operatorname{\mathbf{Im}}(A)\)—a natural output is the least-squares solution, obtained by applying the Moore–Penrose pseudoinverse of \(A\). Our truncation-based solvers can be adapted to achieve this. Finally, we have focused primarily on query complexity. Other resource measures—such as space complexity, time complexity, and query depth—are also of interest and merit further investigation.

3 Beyond-\(\kappa\) solver based on effective truncation of linear system↩︎

In this section, we present an algorithmic template to construct beyond-\(\kappa\) solvers based on effective truncation of the linear system. We begin with a formal definition of the effective condition number \(\kappa_{\text{eff}}\) in Section [[sec:trunc_keff_def]](#sec:trunc_keff_def){reference-type=“ref” reference=“sec:trunc_keff_def”} and present a family of upper bounds on \(\kappa_{\text{eff}}\) in Section [[sec:trunc_keff_bound]](#sec:trunc_keff_bound){reference-type=“ref” reference=“sec:trunc_keff_bound”}, proving Theorem [1](#thm:keff){reference-type=“ref” reference=“thm:keff”}. We then describe an algorithm based on effective truncation of the coefficient matrix and analyze its runtime in Theorem [[thm:beyondk_strong_trunc]](#thm:beyondk_strong_trunc){reference-type=“ref” reference=“thm:beyondk_strong_trunc”} in Section [[sec:trunc_matrix]](#sec:trunc_matrix){reference-type=“ref” reference=“sec:trunc_matrix”}. Next, we construct an oracle in Section [[sec:trunc_oracle]](#sec:trunc_oracle){reference-type=“ref” reference=“sec:trunc_oracle”} to realize effective truncation of the initial state. Finally, we propose and analyze a solver based on effective state truncation in Theorem [[thm:beyondk_weak_trunc]](#thm:beyondk_weak_trunc){reference-type=“ref” reference=“thm:beyondk_weak_trunc”} in Section [[sec:trunc_state]](#sec:trunc_state){reference-type=“ref” reference=“sec:trunc_state”}.

3.1 Effective condition number↩︎

Let \(A=\sum_j\sigma_j|u_j\rangle\!\langle v_j|\) be the singular value decomposition of the coefficient matrix with singular values \(\kappa^{-1}\leq\sigma_j\leq1\). Define the singular vector projections \(\Pi_{\text{left},\mathcal{S}}=\sum_{\sigma_j\in\mathcal{S}}|u_j\rangle\!\langle u_j|\) and \(\Pi_{\text{right},\mathcal{S}}=\sum_{\sigma_j\in\mathcal{S}}|v_j\rangle\!\langle v_j|\) for \(\mathcal{S}\subseteq\mathbb{R}\) a subset of real numbers. Our aim is to choose an effective condition number \(\kappa_{\text{eff}}\) such that \[\frac{\left\lVert\Pi_{\text{right},\left[0,\kappa_{\text{eff}}^{-1}\right)}A^{-1}|b\rangle\right\rVert}{\left\lVert A^{-1}|b\rangle\right\rVert}\leq\epsilon,\qquad \frac{\left\lVert\Pi_{\text{right},\left[0,\kappa_{\text{eff}}^{-1}\right]}A^{-1}|b\rangle\right\rVert}{\left\lVert A^{-1}|b\rangle\right\rVert}>\epsilon.\] Intuitively, we use \(\kappa_{\text{eff}}^{-1}\) to truncate an \(\epsilon\)-fraction of weight from the solution vector \(x=A^{-1}|b\rangle\), starting from the smallest singular values. In cases where multiple \(\kappa_{\mathrm{eff}}\) qualify, we shall select the smallest to minimize the cost of our beyond-\(\kappa\) solvers (corresponding to the largest possible \(\kappa_{\text{eff}}^{-1}\)). This defines the effective condition number as a function \(\kappa_{\mathrm{eff}}(A, |b\rangle, \epsilon)\); we suppress the arguments when clear from context.

Although the effective condition number can be defined and analyzed directly, we find it more natural to work with cumulative distribution functions and their generalized inverses. Specifically, for any non-decreasing function \(f:\mathbb{R}\rightarrow\mathbb{R}\), we denote its one-sided limits by \[f(x-) = \lim_{\varepsilon \searrow 0} f(x - \varepsilon), \qquad f(x+) = \lim_{\varepsilon \searrow 0} f(x + \varepsilon).\] Since \(f\) is non-decreasing, both limits exist and satisfy \(f(x-) \leq f(x) \leq f(x+)\) for all \(x\). We then say that \(x\) is a \(y\)-quantile of \(f\) if \(f(x-) \leq y \leq f(x+)\). In general, the \(y\)-quantiles are not uniquely determined when \(f\) is non-invertible. We will use the largest possible value to define the generalized inverse \(f^+\).

Similarly, one can choose the smallest \(y\)-quantile and write \(f^-(y) = \inf\{x : f(x-) \leq y \leq f(x+)\}\), which is more commonly used in probability. These two definitions are conceptually related though not exactly identical. For a detailed treatment of this topic, we refer the reader to the excellent note of [35]; see also Figure [[fig:quantile_all]](#fig:quantile_all){reference-type=“ref” reference=“fig:quantile_all”}. For our analysis, we need the following properties of the generalized inverse.

The quantile relation \(f(f^+(y)-) \leq y \leq f(f^+(y))\) cannot be refined in general as both inequalities are saturated when \(f\) is an invertible continuous function. However, it can be further strengthened when \(f\) is a step function. Specifically, given \(n\) discrete points \(x_1<x_2<\cdots<x_n\) with associated weights \(w_k>0\), define the step function \(f(x)=\sum_{x_k\leq x}w_k\). Then it holds \(f(f^+(y)-) \leq y < f(f^+(y))\).

Suppose that the initial vector expands under the left singular vector basis as \(|b\rangle=\sum_jw_j|u_j\rangle\). We now consider the function \[f(\sigma)=\frac{\left\lVert\Pi_{\text{right},\left[0,\sigma\right]}A^{-1}|b\rangle\right\rVert^2}{\left\lVert A^{-1}|b\rangle\right\rVert^2} =\frac{\sum_{\sigma_j\leq \sigma}\sigma_j^{-2}\left\lvert w_j\right\rvert^2}{\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2}.\] This function can be seen as the cumulative distribution function of the discrete weights \(\frac{\sigma_j^{-2}\left\lvert w_j\right\rvert^2}{\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2}\). Applying the above analysis, we obtain \[\begin{align} \frac{\sum_{\sigma_j<f^{+}(\epsilon^2)}\sigma_j^{-2}\left\lvert w_j\right\rvert^2}{\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2} =f\left(f^{+}(\epsilon^2)-\right) \leq\epsilon^2,\qquad \frac{\sum_{\sigma_j\leq f^{+}(\epsilon^2)}\sigma_j^{-2}\left\lvert w_j\right\rvert^2}{\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2} =f\left(f^{+}(\epsilon^2)\right) >\epsilon^2. \end{align}\] Finally, let us define \[\begin{align} \kappa_{\text{eff}}=\frac{1}{f^{+}(\epsilon^2)} =\sup\left\{\sigma^{-1}\;\bigg|\;\frac{\sum_{\sigma_j\leq \sigma}\sigma_j^{-2}\left\lvert w_j\right\rvert^2}{\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2}>\epsilon^2\right\} =\inf\left\{\sigma^{-1}\;\bigg|\;\frac{\sum_{\sigma_j\leq \sigma}\sigma_j^{-2}\left\lvert w_j\right\rvert^2}{\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2}\leq\epsilon^2\right\}. \end{align}\] This then gives the desired inequalities \[\begin{align} \frac{\sum_{\sigma_j<\kappa_{\text{eff}}^{-1}}\sigma_j^{-2}\left\lvert w_j\right\rvert^2}{\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2} \leq\epsilon^2,\qquad \frac{\sum_{\sigma_j\leq \kappa_{\text{eff}}^{-1}}\sigma_j^{-2}\left\lvert w_j\right\rvert^2}{\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2} >\epsilon^2. \end{align}\]

3.2 Bounds on the effective condition number↩︎

The effective condition number \(\kappa_{\mathrm{eff}}\) from the previous subsection serves as a threshold for truncating the linear system, enabling quantum linear system solving beyond the condition number. However, evaluating \(\kappa_{\mathrm{eff}}\) directly can be computationally difficult. We now present a family of upper bounds on \(\kappa_{\mathrm{eff}}\) that can be used to bound the runtime of truncation-based beyond-\(\kappa\) solvers.

Specifically, for any tunable parameter \(t>0\), we have \[\sum_j\sigma_j^{-2-2t}\left\lvert w_j\right\rvert^2 \geq\sum_{\sigma_j\leq\kappa_{\text{eff}}^{-1}}\sigma_j^{-2-2t}\left\lvert w_j\right\rvert^2 \geq\kappa_{\text{eff}}^{2t}\sum_{\sigma_j\leq\kappa_{\text{eff}}^{-1}}\sigma_j^{-2}\left\lvert w_j\right\rvert^2 \geq\kappa_{\text{eff}}^{2t}\epsilon^2\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2,\] which gives a family of upper bounds indexed by \(t\): \[\kappa_{\text{eff}}\leq\left(\frac{\sum_j\sigma_j^{-2-2t}\left\lvert w_j\right\rvert^2}{\sum_{k}\sigma_k^{-2}\left\lvert w_k\right\rvert^2\epsilon^2}\right)^{\frac{1}{2t}}.\] The hierarchy of these bounds can be understood through the following weighted power mean inequality.

To apply this inequality, we rewrite our bounds on \(\kappa_{\text{eff}}\) as \[\kappa_{\text{eff}} \leq\frac{\left(\sum_jp_ja_j^{2t}\right)^{\frac{1}{2t}}}{\epsilon^{\frac{1}{t}}},\qquad a_j=\sigma_j^{-1},\qquad p_j=\frac{\sigma_j^{-2}\left\lvert w_j\right\rvert^2}{\sum_k\sigma_k^{-2}\left\lvert w_k\right\rvert^2}.\] Here, we have \(a_j>0\) and we can assume without loss of generality that \(p_j>0\) for all \(j\) as well, for otherwise we simply drop the zero terms from the summand. By the power mean inequality, we see that the numerator monotonically increases to \[\left(\sum_jp_ja_j^{2t}\right)^{\frac{1}{2t}}\nearrow\max_j\sigma_j^{-1}.\] This is at most the condition number \(\kappa\), with equality when the initial state \(|b\rangle\) has nonzero support on the left singular vector of \(A\) with the smallest singular value. On the other hand, the denominator monotonically increases to \[\epsilon^{\frac{1}{t}}\nearrow 1.\]

Recall that \(A=\sum_j\sigma_j|u_j\rangle\!\langle v_j|\), \(|b\rangle=\sum_jw_j|u_j\rangle\), and \(x=A^{-1}|b\rangle=\sum_j\sigma_j^{-1}w_j|v_j\rangle\). Using the singular value transformation with an odd function extension [25], \[\left(\sum_j\sigma_j|u_j\rangle\!\langle v_j|\right)_{\mathbf{sv}}^{-t} =\sum_j\sigma_j^{-t}|u_j\rangle\!\langle v_j|,\] we can then represent the bounds compactly as \[\kappa_{\text{eff}}\leq\left(\frac{\left\lVert(A)_{\mathbf{sv}}^{-t}x\right\rVert}{\left\lVert x\right\rVert\epsilon}\right)^{\frac{1}{t}}.\] For positive integers \(t\), the singular value transformation reduces to the standard matrix inversion as: \[\kappa_{\text{eff}}\leq\begin{cases} \left(\frac{\left\lVert A^{-1\dagger}(A^\dagger A)^{-(t-1)/2}x\right\rVert}{\left\lVert x\right\rVert\epsilon}\right)^{1/t},\quad&t\text{ is odd},\\ \left(\frac{\left\lVert(A^\dagger A)^{-t/2}x\right\rVert}{\left\lVert x\right\rVert\epsilon}\right)^{1/t},\quad&t\text{ is even}. \end{cases}\] These matrix inversions may then be analyzed using standard linear algebraic techniques.

In particular, when \(t=1\), we can bound the effective condition number by \[\kappa_{\text{eff}}\leq\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert\epsilon}.\] This yields beyond-\(\kappa\) solvers whose runtime is comparable to that of the solver from Section [[sec:gap]](#sec:gap){reference-type=“ref” reference=“sec:gap”} based on filtering with effective gap. However, other choices of \(t\) may yield solvers with even lower query complexity. For instance, we can choose \[t=\frac{\ln\left(\frac{1}{\epsilon}\right)}{\ln\ln\left(\frac{1}{\epsilon}\right)}\] such that \[\left(\frac{1}{\epsilon}\right)^{\frac{1}{t}} =\ln\left(\frac{1}{\epsilon}\right).\] The solver thus achieves a logarithmic dependence on the inverse accuracy without increasing the vector norm factor beyond the condition number. We collect our findings on the effective condition number in the following theorem.

3.3 Analysis of beyond-\(\kappa\) solver based on effective matrix truncation↩︎

We are now ready to present and analyze our beyond-\(\kappa\) quantum linear system solver based on an effective matrix truncation. The output of the solver will approximate a quantum state proportional to the truncated solution \(\Pi_{\text{right},\left[\kappa_{\text{eff}}^{-1},1\right]}x\). The following corollary confirms that this truncated solution, after normalization, is close to the true solution state.

As noted in Section [[sec:intro_trunc]](#sec:intro_trunc){reference-type=“ref” reference=“sec:intro_trunc”}, we do not perform the singular vector projection, but instead invoke a quantum linear system solver with effective condition number \(\kappa_{\mathrm{eff}}\). Specifically, we denote the solver by the unitary \(U\left(O_A,O_b,\kappa,\epsilon\right)\), where \(\kappa\) is an upper bound on the condition number and \(\epsilon\) is the target accuracy. The solver then acts on the ancilla and system register jointly and has the behavior \[\begin{align} \left\lVert\left(\langle 0|\otimes I\right)U\left(O_A,O_b,\kappa,\epsilon\right)|00\rangle\right\rVert^2 &>\frac{1}{2},\\ \left\lVert\frac{\left(\langle 0|\otimes I\right)U\left(O_A,O_b,\kappa,\epsilon\right)|00\rangle}{\left\lVert\left(\langle 0|\otimes I\right)U\left(O_A,O_b,\kappa,\epsilon\right)|00\rangle\right\rVert}-|x\rangle\right\rVert &=\operatorname{\mathbf{O}}(\epsilon), \end{align}\] guaranteed to succeed with a constant probability and produce a quantum state close to the normalized solution state. Since \(\kappa_{\mathrm{eff}}\) is not an upper bound on the true condition number, the solver’s output is no longer well defined in the usual sense. Instead, we require it to correctly output the projected solution state supported on terms whose singular values exceed \(\kappa_{\mathrm{eff}}^{-1}\).

We can now construct a beyond-\(\kappa\) solver by invoking a quantum linear system algorithm with the strong truncation property, using a suitable effective condition number.

Theorem [[thm:beyondk_strong_trunc]](#thm:beyondk_strong_trunc){reference-type=“ref” reference=“thm:beyondk_strong_trunc”} applies broadly, converting any conventional quantum linear system solver with the strong truncation property into a beyond-\(\kappa\) solver via truncation. In particular, this applies to VTAA-based solvers, although a direct proof of this claim is somewhat challenging. The difficulty is that, when configured with the effective condition number \(\kappa_{\mathrm{eff}}\) and the objective of preparing a state proportional to \(\Pi_{\text{right}\left[\kappa_{\text{eff}}^{-1},1\right]}|x\rangle\), the solver aims to apply matrix inversion on components of the initial vector corresponding to singular values above \(\kappa_{\mathrm{eff}}^{-1}\) while suppressing those below \(\kappa_{\mathrm{eff}}^{-1}\). However, both components are amplified in VTAA, transforming \(|b\rangle=\sum_jw_j|u_j\rangle\) into a state proportional to \[\frac{1}{\operatorname{\mathbf{\Theta}}(\alpha_{\kappa_{\text{eff}}})}|0\rangle\sum_{\sigma_j\geq\alpha_{\kappa_{\text{eff}}}^{-1}}\sigma_j^{-1}w_j|v_j\rangle +|1\rangle\sum_{\sigma_j<\alpha_{\kappa_{\text{eff}}}^{-1}}w_j|v_j\rangle,\] where we have omitted all irrelevant ancilla registers used by VTAA, and used singular vectors as opposed to eigenvectors of the walk operator [14] for presentational purpose. It may be possible to show directly that a sufficiently good approximation of \(|x\rangle\) can still be extracted from the output of truncated VTAA; however, this nuance may lead to an unconventional amplification schedule that is cumbersome to analyze. To circumvent this, we take a different route: we construct a beyond-\(\kappa\) solver based on an effective truncation of the initial state, which admits a cleaner, more general analysis that applies to VTAA and implies the above matrix truncation result as a corollary.

3.4 Oracle for effective state truncation↩︎

We now develop a beyond-\(\kappa\) solver based on an effective truncation of the initial state. Effectively, the solver replaces every instance of oracle \(O_b\) with \(O_{b_{\text{eff}}}\) that prepares a truncated initial state. The construction of this truncated state preparation oracle relies on the following rotation within the 2D subspace \(\operatorname{\mathbf{Span}}(|b\rangle,|b_{\text{eff}}\rangle)\) spanned by \(|b\rangle\) and \(|b_{\text{eff}}\rangle\).

When \(|b_0\rangle\approx|b_1\rangle\), we expect that the Palais matrix rotating \(|b_0\rangle\) to \(|b_1\rangle\) is close to identity. Explicitly, \[\begin{align} \left\lVert V-I\right\rVert &=\left\lVert-2|w\rangle\!\langle w|+2|b_1\rangle\!\langle b_0|\right\rVert =\left\lVert 2\begin{bmatrix} -|w\rangle & |b_1\rangle \end{bmatrix} \begin{bmatrix} \langle w|\\ \langle b_0| \end{bmatrix}\right\rVert\\ &=\left\lVert 2\begin{bmatrix} \langle w|\\ \langle b_0| \end{bmatrix} \begin{bmatrix} -|w\rangle & |b_1\rangle \end{bmatrix}\right\rVert =\left\lVert 2\begin{bmatrix} -1 & \langle w|b_1\rangle\\ -\langle b_0|w\rangle & \langle b_0|b_1\rangle \end{bmatrix}\right\rVert\\ &=\left\lvert\langle b_0|b_1\rangle-1\pm i\sqrt{1-\langle b_0|b_1\rangle^2}\right\rvert =\sqrt{2-2\langle b_0|b_1\rangle}, \end{align}\] where the third equality uses the fact that \(V-I\) is unitarily diagonalizable with rank \(2\) [33], and the fifth equality holds as \(2\left[\begin{smallmatrix} -1 & \langle w|b_1\rangle\\ -\langle b_0|w\rangle & \langle b_0|b_1\rangle \end{smallmatrix}\right]\) has eigenvalues \(\langle b_0|b_1\rangle-1\pm i\sqrt{1-\langle b_0|b_1\rangle^2}\). When applied to the oracles for state preparation, this gives the following distance bound.

3.5 Analysis of beyond-\(\kappa\) solver based on effective state truncation↩︎

We are now ready to present and analyze our beyond-\(\kappa\) quantum linear system solver based on an effective state truncation. As mentioned in Section [[sec:intro_trunc]](#sec:intro_trunc){reference-type=“ref” reference=“sec:intro_trunc”}, we do not perform the state truncation in the actual implementation of the algorithm. Instead, we run a quantum linear system solver configured with condition number \(\kappa_{\text{eff}}\), while retaining the original input oracles \(O_A\) and \(O_b\) as well as the constant multiplicative approximation of the solution norm \(\left\lVert x\right\rVert\). Suppose that the solver has an optimal query complexity of initial state preparation, i.e., it makes \(\operatorname{\mathbf{O}}\left(\frac{\kappa_{\text{eff}}}{\left\lVert x\right\rVert}\right)\) queries to \(O_b\). Then replacing every occurrence of \(O_b\) by \(O_{b_{\text{eff}}}\) introduces an error at most \(\operatorname{\mathbf{O}}(\epsilon)\). Up to this error, our solver is thus effectively configured with condition number \(\kappa_{\mathrm{eff}}\) and oracles \(O_A\) and \(O_{b_{\text{eff}}}\).

As \(\kappa_{\text{eff}}\) is not an upper bound on the true condition number, this is an unconventional application of quantum linear system solvers. We require that the solver outputs the correct solution state, given the promise that the initial vector has been truncated a priori.

It follows directly from the definition that the strong truncation property implies the weak truncation property. This is because, when given the promise \(|b\rangle\in\operatorname{\mathbf{Im}}\left(\Pi_{\text{left},\left[\alpha_{\kappa_{\text{eff}}}^{-1},1\right]}\right)\), we can say that \[\frac{\Pi_{\text{right},\left[\alpha_{\kappa_{\text{eff}}}^{-1},1\right]}x}{\left\lVert\Pi_{\text{right},\left[\alpha_{\kappa_{\text{eff}}}^{-1},1\right]}x\right\rVert} =\frac{\Pi_{\text{right},\left[\alpha_{\kappa_{\text{eff}}}^{-1},1\right]}A^{-1}|b\rangle}{\left\lVert\Pi_{\text{right},\left[\alpha_{\kappa_{\text{eff}}}^{-1},1\right]}A^{-1}|b\rangle\right\rVert} =\frac{A^{-1}\Pi_{\text{left},\left[\alpha_{\kappa_{\text{eff}}}^{-1},1\right]}|b\rangle}{\left\lVert A^{-1}\Pi_{\text{left},\left[\alpha_{\kappa_{\text{eff}}}^{-1},1\right]}|b\rangle\right\rVert} =\frac{A^{-1}|b\rangle}{\left\lVert A^{-1}|b\rangle\right\rVert} =\frac{x}{\left\lVert x\right\rVert} =|x\rangle.\] Conversely, if the solver has optimal initial state preparation complexity, then the weak truncation property implies the strong truncation property. This gives rise to the following algorithmic template for constructing beyond-\(\kappa\) solvers.

Like Theorem [[thm:beyondk_strong_trunc]](#thm:beyondk_strong_trunc){reference-type=“ref” reference=“thm:beyondk_strong_trunc”}, Theorem [[thm:beyondk_weak_trunc]](#thm:beyondk_weak_trunc){reference-type=“ref” reference=“thm:beyondk_weak_trunc”} is broadly applicable: it converts any solver with the weak truncation property and optimal query cost for state preparation into a beyond-\(\kappa\) solver. However, the weak truncation property has an additional promise on the support of the initial state. Thanks to this promise, verifying the weak truncation property for VTAA is straightforward: the proof from [14] carries over line by line within \(\operatorname{\mathbf{Im}}\left(\Pi_{\text{left},\left[\alpha_{\kappa_{\mathrm{eff}}}^{-1},1\right]}\right)\). We now obtain:

We re-emphasize that the complexity of preparing the solution state \(|x\rangle\) has an implicit \(\epsilon\) dependence through \(\kappa_{\text{eff}}\). On the other hand, for estimating the solution norm \(\left\lVert x\right\rVert\) to constant multiplicative error, one may take \(\epsilon = \operatorname{\pmb{\Theta}}(1)\) when computing \(\kappa_{\text{eff}}\), and there is no implicit \(\epsilon\) dependence. Additionally, our estimator has a cost that scales with any known lower bound \(\alpha_x\) on the solution norm. This may be improved to the exact value of \(\left\lVert x\right\rVert\) using techniques from [39] and [40]. However, we leave the detailed analysis for future work.

4 Beyond-\(\kappa\) solver based on filtering with effective gap↩︎

In this section, we present a beyond-\(\kappa\) solver based on eigenstate filtering with effective gap. We begin by discussing the input model in Section [[sec:gap_decomposition]](#sec:gap_decomposition){reference-type=“ref” reference=“sec:gap_decomposition”} and a constrained orthogonal decomposition induced by the given linear system. We then introduce the effective gap lemma in Section [[sec:gap_effective]](#sec:gap_effective){reference-type=“ref” reference=“sec:gap_effective”}. Next, we present a QSVT-based algorithm in Section [[sec:gap_filtering]](#sec:gap_filtering){reference-type=“ref” reference=“sec:gap_filtering”} for eigenstate filtering over the unit circle. In Section [[sec:gap_solver]](#sec:gap_solver){reference-type=“ref” reference=“sec:gap_solver”}, we propose the filtering-based beyond-\(\kappa\) solver and analyze the constant prefactor of its query complexity. Finally, we develop a solution norm estimation algorithm in Section [[sec:gap_est]](#sec:gap_est){reference-type=“ref” reference=“sec:gap_est”} with a query complexity independent of the condition number, establishing Theorem [[thm:beyondk_filtering]](#thm:beyondk_filtering){reference-type=“ref” reference=“thm:beyondk_filtering”}.

4.1 Block encoding and constrained orthogonal decomposition↩︎

Let \(Ax=|b\rangle\) be a system of linear equations. Recall that the filtering-based solver operates on an augmented matrix of the form \(\begin{bmatrix} A & -|b\rangle \end{bmatrix}\).

Starting from the standard input model—where \(O_A\) block encodes \(A\) with \(\left\lVert A\right\rVert \leq 1\) and \(O_b\) prepares \(|b\rangle\)—one can set \(c = 0\) in the analysis of Section [[sec:input_affine]](#sec:input_affine){reference-type=“ref” reference=“sec:input_affine”} to obtain a block encoding of \(\frac{1}{2}\begin{bmatrix} A & -|b\rangle \end{bmatrix}\) with normalization factor \(2\). However, the normalization factor can be improved to \(\sqrt{2}\) through the following construction \[\frac{1}{\sqrt{2}}\begin{bmatrix} A & -|b\rangle \end{bmatrix} =\frac{1}{\sqrt{2}}\begin{bmatrix} \langle 0|\otimes I & -\langle 0|\otimes I \end{bmatrix} \begin{bmatrix} O_A & 0\\ 0 & I\otimes O_b \end{bmatrix} \begin{bmatrix} |0\rangle\otimes I & 0\\ 0 & |0\rangle\otimes|0\rangle \end{bmatrix}.\] We then denote this block encoding in terms of an overlap of isometries as \[\frac{1}{\sqrt{2}}\begin{bmatrix} A & -|b\rangle \end{bmatrix} =G_1^\dagger G_0,\] where \[G_1^\dagger=\frac{1}{\sqrt{2}}\begin{bmatrix} \left(\langle 0|\otimes I\right)O_A & -\langle 0|\otimes O_b \end{bmatrix},\qquad G_0=\begin{bmatrix} |0\rangle\otimes I & 0\\ 0 & |0\rangle\otimes |0\rangle \end{bmatrix}\] are isometries satisfying \(G_0^\dagger G_0=\left[\begin{smallmatrix} I & 0\\ 0 & 1 \end{smallmatrix}\right]\) and \(G_1^\dagger G_1=I\). Within this block encoding of the augmented matrix, the isometry \(G_1\) makes \(1\) controlled query to \(O_A\) and \(O_b\) in the standard input model, while \(G_0\) requires no query to implement.

There are two projections \(\Pi_0=G_0G_0^\dagger\) and \(\Pi_1=G_1G_1^\dagger\) associated naturally with the block encoding. As \(\Pi_0\) and \(\Pi_1\) do not commute in general, these operators are not simultaneously diagonalizable and the orthogonal decomposition \(\operatorname{\mathbf{Im}}(\Pi_0) = \left(\operatorname{\mathbf{Im}}(\Pi_0) \cap \operatorname{\mathbf{Ker}}(\Pi_1)\right) \obot \left(\operatorname{\mathbf{Im}}(\Pi_0) \cap \operatorname{\mathbf{Im}}(\Pi_1)\right)\) need not hold. Nevertheless, we can utilize the following pair of decompositions.

Taking the smaller subspace from each containment then yields the following constrained orthogonal decomposition \[\left(\operatorname{\mathbf{Im}}(\Pi_0)\cap\operatorname{\mathbf{Ker}}(\Pi_1)\right) \obot\operatorname{\mathbf{Im}}\left(\Pi_0\Pi_1\right).\] Here, a generic vector from \(\operatorname{\mathbf{Im}}(\Pi_0)\) has the form \(G_0\begin{bmatrix} y\\ a \end{bmatrix}\) for some \(y\) and \(a\), partitioned conformally to the augmented matrix. For it to also lie in \(\operatorname{\mathbf{Ker}}(\Pi_1)\), we require that \(G_1^\dagger G_0\begin{bmatrix} y\\ a \end{bmatrix}=\frac{1}{\sqrt{2}}\begin{bmatrix} A & -|b\rangle \end{bmatrix}\begin{bmatrix} y\\ a \end{bmatrix}=0\), which is satisfied by \(G_0\begin{bmatrix} x\\ 1 \end{bmatrix}\).

On the other hand, a generic vector from \(\operatorname{\mathbf{Im}}\left(\Pi_0\Pi_1\right)\) has the form \(G_0G_0^\dagger G_1z =G_0\frac{1}{\sqrt{2}}\begin{bmatrix} A^\dagger z\\ -\langle b|z \end{bmatrix}\) for some vector \(z\). Suppose that we initialize the algorithm with \(G_0\begin{bmatrix} 0\\ 1 \end{bmatrix}\). To find its constrained orthogonal decomposition, we then choose \(A^\dagger z=x\) to match the first component which yields \(z = A^{-1\dagger} x\). This analysis can be slightly generalized to give:

4.2 Effective gap lemma↩︎

In the above decomposition, the first term encodes the solution \(x=A^{-1}|b\rangle\) and is contained in \(\operatorname{\mathbf{Ker}}\left(G_1G_1^\dagger G_0G_0^\dagger\right)\), whereas the second term is enclosed by \(\operatorname{\mathbf{Im}}\left(G_0G_0^\dagger G_1G_1^\dagger\right)\). Using the fact that the smallest nonzero singular value of \(G_1^\dagger G_0\) is gapped from zero by \(\operatorname{\mathbf{\Omega}}\left(\frac{1}{\kappa}\right)\), one can apply eigenstate filtering with the underlying operator \(G_1^\dagger G_0\) and extract the first term, although doing so incurs a query cost scaling linearly in the coniditon number \(\kappa\) [18].

In contrast, we now present the following adaption of the effective gap lemma, which can be applied with the constrained orthogonal decomposition to give a filtering with query cost independent of the condition number.

4.3 Eigenstate filtering over unit circle↩︎

In light of the constrained orthogonal decomposition, our goal is to choose a filter function to preserve the eigenvalue \(+1\) of the underlying unitary operator \(W\) and suppress eigenvalues with gap \(\delta\) on the unit circle. That is, we want to perform \(\operatorname{\mathbf{F}}(W)\) such that \[\operatorname{\mathbf{F}}(1)=1,\qquad \max_{\delta\leq\left\lvert\theta\right\rvert\leq\pi}\left\lvert\operatorname{\mathbf{F}}(e^{i\theta})\right\rvert\leq\xi,\qquad \max_{\theta\in[-\pi,\pi]}\left\lvert\operatorname{\mathbf{F}}(e^{i\theta})\right\rvert=1,\] with circular gap \(\delta>0\) and precision \(\xi>0\). All the remaining eigenvalues with \(\left\lvert\theta\right\rvert<\delta\) are then suppressed by the effective gap lemma.

In this subsection, we present a QSVT-based algorithm to realize this eigenstate filtering over the unit circle, correcting the phase ambiguity in the method of [29]. The core idea is to shift the input unitary operator, so that the eigenvalue of interest is mapped to the origin. See also [41]. Specifically, starting from the unitary \(W\), we construct the block encoding \[\frac{W-I}{2} =\left(\frac{\langle 0|+\langle 1|}{\sqrt{2}}\otimes I\right) \left(|0\rangle\!\langle 0|\otimes W-|1\rangle\!\langle 1|\otimes I\right) \left(\frac{|0\rangle+|1\rangle}{\sqrt{2}}\otimes I\right).\] Their eigenvalues correspond as follows.

1. Eigenvalues of \(W\) are on the unit circle. The eigenvalue of interest is \(+1\), while eigenvalues to be suppressed are \(\delta\)-gapped away in circular distance.

2. Eigenvalues of \(\frac{W-I}{2}\) are on the circle of radius \(\frac{1}{2}\) centered at \(-\frac{1}{2}\). The eigenvalue of interest is \(0\), while eigenvalues to be suppressed are \(\sin\left(\frac{\delta}{2}\right)\)-gapped away in Euclidean distance.

See Figure 4 for an illustration of this shifting.

After the shifting, \(\frac{W-I}{2}\) is still unitarily diagonalizable, and its eigenvalues coincide with singular values in magnitude. The singular value of interest is \(0\), while the singular values to be suppressed are \(\sin\left(\frac{\delta}{2}\right)\)-gapped away. We achieve this filtering using QSVT to perform the Dolph-Chebyshev polynomials, whose relevant properties are summarized in the following lemma.

4.4 Constant-prefactor analysis of beyond-\(\kappa\) solver based on filtering↩︎

We now present our beyond-\(\kappa\) solver based on eigenstate filtering over the unit circle with effective gap. We also estimate the constant prefactor of its query complexity to leading order.

Let \(A\) be a matrix with \(\left\lVert A\right\rVert\leq1\) block encoded by \(O_A\), and \(|b\rangle\) be a normalized quantum state prepared by \(O_b\). Assume that \(A\) is invertible and denote \(x=A^{-1}|b\rangle\). Applying Corollary [[cor:gap_block]](#cor:gap_block){reference-type=“ref” reference=“cor:gap_block”}, we obtain a block encoding of the augmented matrix as an overlap of isometries \[\frac{1}{\sqrt{\beta^2+1}}\begin{bmatrix} \beta A & -|b\rangle \end{bmatrix} =G_1^\dagger G_0\] for any scalar \(\beta>0\), using \(1\) query to \(O_A\) and \(O_b\). Moreover, this induces the constrained orthogonal decomposition \[\begin{align} G_0\begin{bmatrix} 0\\ 1 \end{bmatrix} &=\frac{\beta}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}} \frac{1}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}} G_0\begin{bmatrix} x\\\beta \end{bmatrix} -\frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}\frac{\beta}{\left\lVert x\right\rVert\sqrt{\left\lVert x\right\rVert^2+\beta^2}}G_0\begin{bmatrix} x\\ -\frac{\left\lVert x\right\rVert^2}{\beta} \end{bmatrix}\\ &=\frac{\beta}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}} \frac{1}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}G_0\begin{bmatrix} x\\ \beta \end{bmatrix} -\frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}} \Pi_0 \Pi_1 \frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert\sqrt{\left\lVert x\right\rVert^2+\beta^2}}G_1A^{-1\dagger}x,\\ \end{align}\] where \(\Pi_0=G_0G_0^\dagger\) and \(\Pi_1=G_1G_1^\dagger\).

Define the unitary \(W=\left(2\Pi_0-I\right)\left(I-2\Pi_1\right)=\sum_ue^{i\theta_u}|\phi_u\rangle\!\langle\phi_u|\) and perform the filtering \(\operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right)\) with circular gap \(0<\delta<\pi\) and precision \(0<\xi<1\), where we have omitted the subscripts for notational convenience. By the effective gap lemma (Lemma [[lem:eff_gap]](#lem:eff_gap){reference-type=“ref” reference=“lem:eff_gap”}), \[\begin{align} \operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right)G_0\begin{bmatrix} 0\\ 1 \end{bmatrix} &=\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2} G_0\begin{bmatrix} x\\\beta \end{bmatrix} -\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2}\operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right)G_0\begin{bmatrix} x\\ -\frac{\left\lVert x\right\rVert^2}{\beta} \end{bmatrix}\\ &=\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2}G_0\begin{bmatrix} x\\ \beta \end{bmatrix} -\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2}\operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right) \Pi_0 \Pi_1G_1A^{-1\dagger}x, \end{align}\] where Corollary [[cor:dolph_cheby]](#cor:dolph_cheby){reference-type=“ref” reference=“cor:dolph_cheby”} shows \[\operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right) =\sum_u\operatorname{\mathbf{F}}\left(\sin\left(\frac{\theta_u}{2}\right)\right)|\phi_u\rangle\!\langle\phi_u|.\] We now split the second term further into two sub terms. To this end, we introduce the circular indicator function \[\operatorname{\mathbf{Ind}}_{\left[0,\sin\left(\frac{\delta}{2}\right)\right]}\left(\sin\left(\frac{\theta}{2}\right)\right) =\begin{cases} 1,\quad&\left\lvert\theta\right\rvert<\delta,\\ 0,&\delta\leq\left\lvert\theta\right\rvert\leq\pi. \end{cases}\] This yields the final decomposition \[\begin{align} \operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right)G_0\begin{bmatrix} 0\\ 1 \end{bmatrix} &=\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2} G_0\begin{bmatrix} x\\\beta \end{bmatrix}\\ &\quad-\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2} \sum_u\left(\operatorname{\mathbf{F}}\operatorname{\mathbf{Ind}}\right)\left(\sin\left(\frac{\theta_u}{2}\right)\right)|\phi_u\rangle\!\langle\phi_u| \Pi_0 \Pi_1G_1A^{-1\dagger}x\\ &\quad-\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2} \sum_u\left(\operatorname{\mathbf{F}}(1-\operatorname{\mathbf{Ind}})\right)\left(\sin\left(\frac{\theta_u}{2}\right)\right)|\phi_u\rangle\!\langle\phi_u| G_0\begin{bmatrix} x\\ -\frac{\left\lVert x\right\rVert^2}{\beta} \end{bmatrix}.\\ \end{align}\]

In this new decomposition, the third term can be bounded directly using Corollary [[cor:dolph_cheby]](#cor:dolph_cheby){reference-type=“ref” reference=“cor:dolph_cheby”} as \[\begin{align} &\left\lVert\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2} \sum_u\left(\operatorname{\mathbf{F}}(1-\operatorname{\mathbf{Ind}})\right)\left(\sin\left(\frac{\theta_u}{2}\right)\right)|\phi_u\rangle\!\langle\phi_u| G_0\begin{bmatrix} x\\ -\frac{\left\lVert x\right\rVert^2}{\beta} \end{bmatrix}\right\rVert\\ &\leq\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2} \sqrt{\left\lVert x\right\rVert^2+\frac{\left\lVert x\right\rVert^4}{\beta^2}} \max_{\theta\in[-\pi,\pi]}\left\lvert\operatorname{\mathbf{F}}\left(\sin\left(\frac{\theta}{2}\right)\right) (1-\operatorname{\mathbf{Ind}})\left(\sin\left(\frac{\theta}{2}\right)\right)\right\rvert\\ &=\frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}} \max_{\delta\leq\left\lvert\theta\right\rvert\leq\pi}\left\lvert\operatorname{\mathbf{F}}\left(\sin\left(\frac{\theta}{2}\right)\right)\right\rvert \leq\frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}\xi, \end{align}\] whereas for the second term we apply the effective gap lemma (Lemma [[lem:eff_gap]](#lem:eff_gap){reference-type=“ref” reference=“lem:eff_gap”}) \[\begin{align} &\left\lVert\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2} \sum_u\left(\operatorname{\mathbf{F}}\operatorname{\mathbf{Ind}}\right)\left(\sin\left(\frac{\theta_u}{2}\right)\right)|\phi_u\rangle\!\langle\phi_u| \Pi_0 \Pi_1G_1A^{-1\dagger}x\right\rVert\\ &\leq\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2} \max_{\theta\in[-\pi,\pi]}\left\lvert\operatorname{\mathbf{F}}\left(\sin\left(\frac{\theta}{2}\right)\right) \operatorname{\mathbf{Ind}}\left(\sin\left(\frac{\theta}{2}\right)\right) \frac{\theta}{2}\right\rvert \left\lVert G_1A^{-1\dagger}x\right\rVert\\ &\leq\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2} \max_{\left\lvert\theta\right\rvert<\delta}\left\lvert\operatorname{\mathbf{F}}\left(\sin\left(\frac{\theta}{2}\right)\right) \frac{\theta}{2}\right\rvert \left\lVert A^{-1\dagger}x\right\rVert \leq\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2}\left\lVert A^{-1\dagger}x\right\rVert\frac{\delta}{2}. \end{align}\] We thus obtain the following error bound for the eigenstate filtering \[\left\lVert\operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right)G_0\begin{bmatrix} 0\\ 1 \end{bmatrix} -\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2} G_0\begin{bmatrix} x\\\beta \end{bmatrix}\right\rVert \leq\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2}\left\lVert A^{-1\dagger}x\right\rVert\frac{\delta}{2} +\frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}\xi.\] This filtering uses \(2\operatorname{\mathbf{Ceil}}\left(\frac{1}{\delta}\ln\left(\frac{2}{\xi}\right)\right)\) queries to controlled \(W\) and \(W^\dagger\). With \(W=(2G_0G_0^\dagger-I)(I-2G_1G_1^\dagger)\), each \(W\) or \(W^\dagger\) can be realized using \(2\) (uncontrolled) queries to the block encoding isometries \(G_0\) and \(G_1\). Each \(G_1\) then requires \(1\) controlled application of \(O_A\) and \(O_b\) whereas \(G_0\) can be realized with zero query ([eq:def_g]), giving the cost \[4\operatorname{\mathbf{Ceil}}\left(\frac{1}{\delta}\ln\left(\frac{2}{\xi}\right)\right)\] in the standard input model.

Let us choose the effective gap \(\delta\) and filtering accuracy \(\xi\) so that the solution vector \(x\) has error at most \(\left\lVert x\right\rVert\epsilon\). Equivalently, we require that \[\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2}\left\lVert A^{-1\dagger}x\right\rVert\frac{\delta}{2} +\frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}\xi \leq\frac{\beta\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+\beta^2}\epsilon.\] A standard strategy is to divide the total error budget into two parts of comparable size as in [29]. Here, we present an alternative approach that matches their asymptotic scaling but improves the constant prefactor. Specifically, we let \[\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2}\left\lVert A^{-1\dagger}x\right\rVert\frac{\delta}{2} =\gamma\frac{\beta\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+\beta^2}\epsilon,\qquad \frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}\xi =(1-\gamma)\frac{\beta\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+\beta^2}\epsilon,\] for some tunable parameter \(0<\gamma<1\) to be chosen momentarily. Then the query cost becomes \[\begin{align} \operatorname{\mathbf{Ceil}}\left(\frac{4}{\delta}\ln\left(\frac{2}{\xi}\right)\right) &=\operatorname{\mathbf{Ceil}}\left(\frac{2\sqrt{\beta^2+1}\left\lVert A^{-1\dagger}x\right\rVert}{\beta\left\lVert x\right\rVert} \frac{1}{\gamma\epsilon} \ln\left(\frac{2\sqrt{\left\lVert x\right\rVert^2+\beta^2}}{\beta}\frac{1}{(1-\gamma)\epsilon}\right)\right)\\ &=\operatorname{\mathbf{Ceil}}\left(\frac{\sqrt{\beta^2+1}\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert\sqrt{\left\lVert x\right\rVert^2+\beta^2}}\frac{1}{\gamma\widetilde{\epsilon}}\ln\left(\frac{1}{(1-\gamma)\widetilde{\epsilon}}\right)\right), \end{align}\] where \[\widetilde{\epsilon} =\frac{\beta\epsilon}{2\sqrt{\left\lVert x\right\rVert^2+\beta^2}}.\]

In Appendix [[append:lambertw]](#append:lambertw){reference-type=“ref” reference=“append:lambertw”}, we prove that the function \[g(\gamma)=\frac{1}{\gamma\widetilde{\epsilon}}\ln\left(\frac{1}{(1-\gamma)\widetilde{\epsilon}}\right),\qquad 0<\gamma<1.\] is minimized at \(\gamma_0=1+\frac{1}{\operatorname{\mathbf{W}_{-1}}\left(-\frac{\widetilde{\epsilon}}{e}\right)}\), where \(\operatorname{\mathbf{W}_{-1}}(\cdot)\) is the (\(-1\))-branch of the Lambert-W function. Correspondingly, \[\begin{align} \min g(\gamma) =g(\gamma_0) &=\frac{\operatorname{\mathbf{W}_{-1}}\left(-\frac{\widetilde{\epsilon}}{e}\right)}{1+\operatorname{\mathbf{W}_{-1}}\left(-\frac{\widetilde{\epsilon}}{e}\right)} \frac{1}{\widetilde{\epsilon}} \ln\left(\frac{-\operatorname{\mathbf{W}_{-1}}\left(-\frac{\widetilde{\epsilon}}{e}\right)}{\widetilde{\epsilon}}\right)\\ \end{align}\] which reduces to \[\frac{1}{\widetilde{\epsilon}}\ln\left(\frac{1}{\widetilde{\epsilon}}\right)\] when higher order terms are dropped. However, the exact values of \(\left\lVert x\right\rVert\) and \(\left\lVert A^{-1\dagger}x\right\rVert\) are typically unavailable in practice. Instead, we assume that \(\alpha_{A^{-1\dagger}|x\rangle}\geq\left\lVert A^{-1\dagger}|x\rangle\right\rVert\) is a norm upper bound and \(\frac{1}{\mu}\leq\frac{\alpha_x}{\left\lVert x\right\rVert}\leq\mu\) is a \(\mu\)-multiplicative approximation of the solution norm \(\left\lVert x\right\rVert\) for some \(\mu\geq1\). Setting \(\beta=\alpha_x\), we can say that the query complexity is at most

since shrinking the error budget from \(\widetilde{\epsilon}=\frac{\alpha_x\epsilon}{2\sqrt{\left\lVert x\right\rVert^2+\alpha_x^2}}\) to \(\widetilde{\epsilon}=\frac{\epsilon}{2\sqrt{\mu^2+1}}\) never decreases the cost. To leading order, this simplifies to \[\frac{\mu\sqrt{\alpha_x^2+1}\alpha_{A^{-1\dagger}|x\rangle}}{\sqrt{\mu^2+1}\alpha_x} \frac{2\sqrt{\mu^2+1}}{\epsilon} \ln\left(\frac{2\sqrt{\mu^2+1}}{\epsilon}\right) \simeq\frac{2\alpha_{A^{-1\dagger}|x\rangle}}{\epsilon} \ln\left(\frac{1}{\epsilon}\right).\]

With the above choice of effective gap \(\delta\), filtering accuracy \(\xi\) and scalar \(\beta\), the success amplitude is lower bounded by \[\frac{\beta\left\lVert x\right\rVert(1-\epsilon)}{\left\lVert x\right\rVert^2+\beta^2} =\frac{\alpha_x\left\lVert x\right\rVert(1-\epsilon)}{\left\lVert x\right\rVert^2+\alpha_x^2} =\frac{(1-\epsilon)}{\frac{\left\lVert x\right\rVert}{\alpha_x}+\frac{\alpha_x}{\left\lVert x\right\rVert}} \geq\frac{1-\epsilon}{\mu+\frac{1}{\mu}}.\] Applying the fixed-point amplitude amplification [12], we obtain the solution state \(|x\rangle=\frac{x}{\left\lVert x\right\rVert}\) with accuracy \(\epsilon\) and failure probability \(p_{\text{fail}}\) using \[2\operatorname{\mathbf{Floor}}\left(\frac{1}{2}\operatorname{\mathbf{Ceil}}\left(\frac{\mu+\frac{1}{\mu}}{1-\epsilon}\ln\left(\frac{2}{\sqrt{p_{\text{fail}}}}\right)\right)\right)+1\] repetitions of filtering. Here, we have used \(2\operatorname{\mathbf{Floor}}\left(\frac{1}{2}\operatorname{\mathbf{Ceil}}\left(\cdot\right)\right)+1\) to denote the smallest odd integer greater than or equal to its argument.

Altogether, to succeed with a probability at least \(\frac{1}{2}\), the solver makes

queries to the oracles \(O_A\) and \(O_b\) (recall from [eq:def_g] that undoing the isometry \(G_0^\dagger\) has zero query cost). In the limit \(\mu\rightarrow1\), \(\epsilon\rightarrow0\), \(\alpha_x\rightarrow\infty\), \(\alpha_{A^{-1\dagger}|x\rangle}\rightarrow\infty\), the leading-order query complexity is \[\left(2\operatorname{\mathbf{Floor}}\left(\frac{1}{2}\operatorname{\mathbf{Ceil}}\left(3\ln\left(2\right)\right)\right)+1\right) \frac{2\alpha_{A^{-1\dagger}|x\rangle}}{\epsilon} \ln\left(\frac{1}{\epsilon}\right) =\frac{6\alpha_{A^{-1\dagger}|x\rangle}}{\epsilon} \ln\left(\frac{1}{\epsilon}\right).\]

4.5 Beyond-\(\kappa\) solution norm estimation↩︎

The simple beyond-\(\kappa\) solver from the previous subsection relies on the assumption that the solution norm \(\left\lVert x\right\rVert\) is known to a constant relative accuracy. If we do not have this prior knowledge, then we can turn the solver into a solution norm estimation algorithm with essentially the same query complexity up to logarithmic factors.

To start, we recall that the eigenstate filtering introduces an error at most \[\left\lVert\operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right)G_0\begin{bmatrix} 0\\ 1 \end{bmatrix} -\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2} G_0\begin{bmatrix} x\\\beta \end{bmatrix}\right\rVert \leq\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2}\left\lVert A^{-1\dagger}x\right\rVert\frac{\delta}{2} +\frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}\xi.\] Unlike before, we now aim to choose \(\delta\) and \(\xi\), so that when projected onto the first component, the error in the success amplitude is sufficiently small. Additionally, we perform an exponential search over the rescaling factor \(\beta=1,2,4,\ldots\), so \(\beta\geq1\) is guaranteed. This means we can choose \[\delta=\operatorname{\mathbf{\Theta}}\left(\frac{1}{\alpha_{A^{-1\dagger}|x\rangle}}\right),\qquad \xi=\operatorname{\mathbf{\Theta}}(1),\] where \(\alpha_{A^{-1\dagger}x}\geq\left\lVert A^{-1\dagger}|x\rangle\right\rVert\) is a known upper bound, so that \[\begin{align} \frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2}\left\lVert A^{-1\dagger}x\right\rVert\frac{\delta}{2} &\leq\frac{\sqrt{2}\beta\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+\beta^2} \left\lVert A^{-1\dagger}|x\rangle\right\rVert\frac{\delta}{2} \leq\frac{1}{8}\frac{\beta\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+\beta^2},\\ \frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}\xi &\leq\xi \leq\frac{1}{16}. \end{align}\] Then, the filtering has asymptotic query complexity \[\operatorname{\mathbf{O}}\left(\alpha_{A^{-1\dagger}|x\rangle}\right).\]

Note that the success amplitude is bounded within the interval \[\left[\frac{7}{8}\frac{\beta\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+\beta^2}-\frac{1}{16},\; \frac{9}{8}\frac{\beta\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+\beta^2}+\frac{1}{16}\right] \subseteq\left[\frac{7}{8}\frac{\beta\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+\beta^2}-\frac{1}{16},\; \frac{9}{8}\frac{\beta}{\left\lVert x\right\rVert}+\frac{1}{16}\right].\] Let us choose the rescaling factor \(\beta_j=2^j\) for \(j=0,1,\ldots,\operatorname{\mathbf{Ceil}}\left(\log_2\left(\alpha_{A^{-1\dagger}|x\rangle}\right)\right)\) and consider two cases. If \(2^j\leq\frac{\left\lVert x\right\rVert}{16}\), then the success amplitude is at most \[\frac{9}{8}\frac{2^j}{\left\lVert x\right\rVert}+\frac{1}{16} \leq\frac{9}{8}\frac{\frac{\left\lVert x\right\rVert}{16}}{\left\lVert x\right\rVert}+\frac{1}{16} =\frac{17}{128}.\] On the other hand, if \(\frac{\left\lVert x\right\rVert}{2}\leq2^j\leq\left\lVert x\right\rVert\), then the success amplitude is at least \[\frac{7}{8}\frac{2^j\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+4^{j}}-\frac{1}{16} \geq\frac{7}{8}\frac{\frac{\left\lVert x\right\rVert^2}{2}}{\left\lVert x\right\rVert^2+\left\lVert x\right\rVert^2}-\frac{1}{16} =\frac{5}{32} =\frac{20}{128}.\]

To distinguish these two cases, we loop through \(j=0,1,\ldots,\operatorname{\mathbf{Ceil}}\left(\log_2\left(\alpha_{A^{-1\dagger}|x\rangle}\right)\right)\), and perform amplitude estimation with a sufficiently small constant precision and failure probability \(p_j=\frac{1}{6\cdot 2^j}\). Then with probability at least \(\frac{2}{3}\), the algorithm produces estimates of amplitudes at most \(\frac{18}{128}\) for all \(j\) such that \(2^j\leq\frac{\left\lVert x\right\rVert}{16}\). On the other hand, if \(\frac{\left\lVert x\right\rVert}{2}\leq2^j\leq\left\lVert x\right\rVert\), the algorithm outputs an estimate at least \(\frac{19}{128}\) with probability at least \(1-\frac{1}{6\cdot 2^j}\). Suppose the first time we find the amplitude above \(\frac{19}{128}\) is at \(j_0\). This means \[2^{j_0-1}<\frac{\left\lVert x\right\rVert}{2},\qquad 2^{j_0}>\frac{\left\lVert x\right\rVert}{16},\] which gives a constant multiplicative approximation of solution norm: \[\frac{1}{4}<\frac{\left\lVert x\right\rVert}{2^{j_0}\cdot4}<4.\] The total query complexity is then at most \[\operatorname{\mathbf{O}}\left(\sum_{j=0}^{\operatorname{\mathbf{Ceil}}\left(\log_2\left(\alpha_{A^{-1\dagger}|x\rangle}\right)\right)} \alpha_{A^{-1\dagger}|x\rangle} \log\left(6\cdot 2^j\right) \right) =\operatorname{\mathbf{O}}\left( \alpha_{A^{-1\dagger}|x\rangle}\log^2\left(\alpha_{A^{-1\dagger}|x\rangle}\right) \right).\] The approximation ratio can be refined to \(\mu\) arbitrarily close to \(1\) using amplitude estimation incurring multiplicative complexity overhead \(1/(\mu-1)\) [12]. We now summarize the filtering-based beyond-\(\kappa\) solver and the corresponding solution norm estimator in the following theorem.

5 Discussion↩︎

In this work, we have developed algorithms that solve the quantum linear system problem with query complexity independent of the spectral condition number. We have proposed an algorithmic template based on effective truncation of the linear system, which converts conventional quantum linear system solvers into beyond-\(\kappa\) solvers. Our template applies to any solver satisfying the strong truncation property, or the weak truncation property with optimal query complexity for initial state preparation. Additionally, we have presented a beyond-\(\kappa\) solver using eigenstate filtering with effective gap. We have shown that the filtering-based solver is particularly simple and achieves a favorable runtime prefactor when the solution norm is known. When it is not, we have designed a similarly straightforward solution norm estimator—making our result fully self-contained.

We have introduced the effective condition number \(\kappa_{\mathrm{eff}}\) as a measure of complexity for the truncation-based solvers, and established a family of upper bounds \(\kappa_{\text{eff}}\leq\left(\frac{\left\lVert(A)_{\mathbf{sv}}^{-t}|x\rangle\right\rVert}{\epsilon}\right)^{\frac{1}{t}}\) in terms of the singular transformation of the coefficient matrix \(A\), normalized solution vector \(|x\rangle\) and adjustable parameter \(0<t<\infty\). The case where \(t\) is a positive integer is especially straightforward, as the transformation reduces to \((A)_{\mathbf{sv}}^{-t}=(A^\dagger A)^{-t/2}\) for even \(t\) and \((A)_{\mathbf{sv}}^{-t}=A^{-1\dagger}(A^\dagger A)^{-(t-1)/2}\) for odd \(t\) and can be analyzed via standard linear algebraic techniques. In particular, choosing \(t=1\) gives \((A)_{\mathbf{sv}}^{-t}=A^{-1\dagger}\). This choice introduces a polynomial dependence on \(1/\epsilon\), which can be improved to a polylogarithmic scaling by setting \(t=\operatorname{\mathbf{\Theta}}\left(\frac{\log\left(\frac{1}{\epsilon}\right)}{\log\log\left(\frac{1}{\epsilon}\right)}\right)\). Finding efficiently computable bounds on such transformations or alternative tight estimates of \(\kappa_{\mathrm{eff}}\) remains an interesting open problem. Our effective condition number \(\kappa_{\mathrm{eff}}\) captures both the strength and the limitation of truncation-based solvers. Its reciprocal is the smallest cutoff for which the truncated linear system can still be solved to accuracy \(\epsilon\). However, it remains plausible that solvers based on different techniques could achieve scaling below \(\kappa_{\text{eff}}\). An interesting direction is to either construct such solvers or prove a lower bound ruling out their existence. Along similar lines, the filtering-based solver is much simpler but only achieves a scaling with \(\frac{\left\lVert A^{-1\dagger}|x\rangle\right\rVert}{\epsilon}\geq\kappa_{\mathrm{eff}}\). Whether this gap can be closed by a simple algorithm is an interesting question for future work.

In analyzing our solvers, we have focused primarily on the standard input model, where the coefficient matrix \(A\) and initial state \(|b\rangle\) are accessed through separate oracles. However, we have also defined the affine dilation model, which block encodes \(A\) and \(|b\rangle\) jointly, allowing further improvements in query complexity. To demonstrate that this improvement can be substantial, we now give a toy example where the separation between the two models diverges as \(\kappa \to \infty\). Specifically, consider \(A=\left[\begin{smallmatrix} \frac{1}{\kappa} & 0\\ 0 & 1 \end{smallmatrix}\right]\), \(|b\rangle=\left[\begin{smallmatrix} \frac{1}{\kappa}\\ \sqrt{1-\frac{1}{\kappa^2}} \end{smallmatrix}\right]\), \(x=A^{-1}|b\rangle=\left[\begin{smallmatrix} 1\\ \sqrt{1-\frac{1}{\kappa^2}} \end{smallmatrix}\right]\), \(A^{-1\dagger}x=\left[\begin{smallmatrix} \kappa\\ \sqrt{1-\frac{1}{\kappa^2}} \end{smallmatrix}\right]\). In the standard model, we can block encode \(A\) with normalization \(1\) and prepare the normalized state \(|b\rangle\). Then the quantity \(\left\lVert A^{-1\dagger}|x\rangle\right\rVert\) which determines the runtime of beyond-\(\kappa\) solvers scales like \[\left\lVert A^{-1\dagger}|x\rangle\right\rVert =\operatorname{\mathbf{\Theta}}\left(\frac{\sqrt{\kappa^2+1-\frac{1}{\kappa^2}}}{\sqrt{1+1-\frac{1}{\kappa^2}}}\right) =\operatorname{\mathbf{\Theta}}(\kappa).\] In contrast, we can define \(\widetilde{A}=\frac{1}{\widetilde{\alpha}} \left[\begin{smallmatrix} MA & -M|b\rangle\\ 0 & 1 \end{smallmatrix}\right]\), \(|\widetilde{b}\rangle=\left[\begin{smallmatrix} 0\\ 1 \end{smallmatrix}\right]\), \(\widetilde{x}=\widetilde{A}^{-1}|\widetilde{b}\rangle =\widetilde{\alpha}\left[\begin{smallmatrix} x\\ 1 \end{smallmatrix}\right]\), and \(\widetilde{A}^{-1\dagger}\widetilde{x} =\widetilde{\alpha}^2\left[\begin{smallmatrix} M^{-1\dagger}A^{-1\dagger}x\\ \left\lVert x\right\rVert^2+1 \end{smallmatrix}\right]\) in the affine dilation model, where \(M=\left[\begin{smallmatrix} \kappa & 0\\ 0 & 1 \end{smallmatrix}\right]\) is a diagonal matrix and \(\widetilde{\alpha}=\operatorname{\mathbf{\Theta}}(1)\) is the normalization factor for block encoding, giving \[\left\lVert\widetilde{A}^{-1\dagger}|\widetilde{x}\rangle\right\rVert =\widetilde{\alpha}\frac{\sqrt{\left\lVert M^{-1\dagger}A^{-1\dagger}x\right\rVert^2+\left(\left\lVert x\right\rVert^2+1\right)^2}}{\sqrt{\left\lVert x\right\rVert^2+1}} =\operatorname{\mathbf{\Theta}}(1).\] This shows that the affine dilation model can offer an arbitrarily large advantage as the condition number \(\kappa\) grows.

The preceding example is in fact a special case of a broader strategy known as preconditioning, which solves the modified linear system \(MAx=M|b\rangle\) within the affine dilation model. More generally, our framework opens a new regime for designing preconditioners for quantum linear system solvers: rather than targeting the condition number \(\kappa(A)\) of matrix \(A\), the purpose of this new preconditioning strategy is to redistribute how \(|b\rangle\) aligns with the left singular vectors of \(A\), thereby reducing the effective condition number. This can yield substantial complexity improvements even when \(\kappa(MA)\) itself remains large, since the gains are captured entirely by the effective condition number.

For sparse linear systems, [18] constructs a preconditioner using diagonal matrices to balance the row norms. Interestingly, the instance-aware parameter \(\operatorname{\mathbf{ET}} = \sum_i p_i^2\, d_i\) introduced in [18] can be recovered within our framework through diagonal preconditioning. Here \(p_i = \left\lvert\langle i|(AA^\dagger)^{-1}|b\rangle\right\rvert\) is the \(i\)-th component of vector \((AA^\dagger)^{-1}|b\rangle\) in absolute value and \(d_i = \|a_i\|^2 + |b_i|^2\) is the squared norm of the \(i\)-th row of the augmented matrix \(\begin{bmatrix} A & -|b\rangle \end{bmatrix}\). The sparse access model normalizes each row by its norm, implicitly block encoding the reweighted system \(\begin{bmatrix} MA & -M|b\rangle \end{bmatrix}\) with diagonal matrix \(M = \operatorname{\mathbf{Diag}}([\ldots,1/\sqrt{d_i},\ldots])\). A direct calculation then gives \((MA)^{-1\dagger}x = M^{-1}(AA^\dagger)^{-1}|b\rangle\), so that \(\left\|(MA)^{-1\dagger}x\right\|^2 = \sum_i p_i^2\, d_i = \operatorname{\mathbf{ET}}\) recovers the result of [18]. Moreover, by designing new input models, one can unlock preconditioning techniques beyond those available in the existing framework. To illustrate, consider a setting where \(\left\lVert|b\rangle\!\langle b|A\right\rVert\) and \(\left\lVert(I - |b\rangle\!\langle b|)A\right\rVert\) are substantially different. In this case, one can introduce an input model that block encodes \(|b\rangle\!\langle b|A\) and \((I - |b\rangle\!\langle b|)A\) separately, each with its own normalization factor. This effectively preconditions the linear system with respect to the orthogonal decomposition \(\operatorname{\mathbf{Im}}(|b\rangle\!\langle b|)\obot\operatorname{\mathbf{Ker}}(|b\rangle\!\langle b|)\). We expect that a systematic exploration of this preconditioning regime, identifying which classical preconditioners can be lifted into suitable input models and what complexity gains they unlock, is a fruitful direction for future work. As another example, [42] showed that for positive-definite matrices, the worst-case runtime lower bound of a quantum linear system solver remains linear in \(\kappa\), matching the indefinite case. Nevertheless, they identified broad classes of positive-definite systems where this lower bound can be circumvented, achieving a quadratic improvement to \(\sqrt{\kappa}\) via more efficient block-encodings or Cholesky-type preconditioners. It would be interesting to investigate whether our beyond-\(\kappa\) techniques can be combined with such structural assumptions to yield further speedups.

Setting aside these technical refinements, we believe the most impactful direction for future work is to identify practical applications where the condition number is prohibitively large for conventional solvers, yet beyond-\(\kappa\) solvers remain efficient. Our filtering-based solver is particularly appealing in this regard—its simplicity makes it a prime candidate for compilation into explicit quantum circuits. A concrete resource estimate for running this algorithm on a fault-tolerant quantum computer, even for a small-scale instance, would be an important step toward demonstrating practical quantum advantage for solving linear systems.

Acknowledgements↩︎

A.M.D. and Y.S. thank Fernando Brandão, Oskar Painter, James Hamilton, Nafea Bshara, Peter DeSantis, and Andy Jassy for their involvement and support of the research activities at the AWS Center for Quantum Computing. J.L. was supported by NSF Award FET-2243659 and NSF Career Award FET-2339116, Welch Foundation Grant no. A22-0307, a Microsoft Research Award, an Amazon Research Award, and a Ken Kennedy Research Cluster award (Rice K2I), in part through seed funding from the Ken Kennedy Institute at Rice University, the Rice CS Department, and the George R. Brown School of Engineering and Computing, by Rice University and the Department of Computer Science at Rice University, by the Ken Kennedy Institute and Rice Quantum Initiative, which is part of the Smalley-Curl Institute.

6 Analysis of Li’s beyond-\(\kappa\) solver in the standard input model↩︎

The beyond-\(\kappa\) solver of Li [18] was originally proposed for solving sparse linear systems. In this appendix, we show that it can be extended to work in the standard input model, where the coefficient matrix is accessed through a block encoding oracle and the initial vector is prepared by a unitary. To simplify the discussion, we use the precise values of vector norms such as \(\left\lVert x\right\rVert\) and \(\left\lVert A^{-1\dagger}x\right\rVert\), rather than their known upper and lower bounds.

Specifically, suppose that \(A\) is a matrix with \(\left\lVert A\right\rVert\leq1\) block encoded by \(O_A\), and \(|b\rangle\) be a normalized quantum state prepared by \(O_b\). Assume that \(A\) is invertible and define the solution vector by \(x=A^{-1}|b\rangle\). Then we know from Corollary [[cor:gap_block]](#cor:gap_block){reference-type=“ref” reference=“cor:gap_block”} that, for any scalar \(\beta>0\), the augmented matrix \[\frac{1}{\sqrt{\beta^2+1}}\begin{bmatrix} \beta A & -|b\rangle \end{bmatrix} =G_1^\dagger G_0\] can be block encoded in terms of an overlap of isometries using \(1\) query to \(O_A\) and \(O_b\). Moreover, this induces a constrained orthogonal decomposition \[\begin{align} G_0\begin{bmatrix} 0\\ 1 \end{bmatrix}&=\frac{\beta}{\left\lVert x\right\rVert^2+\beta^2} G_0\begin{bmatrix} x\\ \beta \end{bmatrix} -\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2} \Pi_0 \Pi_1 G_1A^{-1\dagger}x\\ &\in\left(\operatorname{\mathbf{Im}}(\Pi_0)\cap\operatorname{\mathbf{Ker}}(\Pi_1)\right) \obot\operatorname{\mathbf{Im}}\left(\Pi_0 \Pi_1\right), \end{align}\] where \(\Pi_0=G_0G_0^\dagger\) and \(\Pi_1=G_1G_1^\dagger\).

Here, the goal is to preserve the first term that encodes the solution vector \(x\) while suppressing the second term using the effective gap lemma. In the main text, we show that this can be achieved with QSVT-based eigenstate filtering over the unit circle. In [18], this was instead realized using quantum phase estimation. More concretely, [18] shows that applying quantum phase estimation produces the normalized state \(\frac{1}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}G_0\begin{bmatrix} x\\ \beta \end{bmatrix}\) with accuracy \(\xi\) in trace distance, using \[\operatorname{\mathbf{O}}\left(\frac{\left\lVert x\right\rVert^2+\beta^2}{\beta^2\xi^2} \frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2}\left\lVert A^{-1\dagger}x\right\rVert\right) =\operatorname{\mathbf{O}}\left(\frac{1}{\beta\xi^2} \left\lVert A^{-1\dagger}x\right\rVert\right)\] queries to the oracles \(O_A\) and \(O_b\) and succeeding with probability \[\operatorname{\mathbf{\Theta}}\left(\frac{\beta^2}{\left\lVert x\right\rVert^2+\beta^2}\right),\] assuming \(\beta=\operatorname{\mathbf{\Omega}}(1)\).

To ensure that the normalized solution \(|x\rangle=\frac{x}{\left\lVert x\right\rVert}\) is produced with accuracy \(\epsilon\), Li chose \(\xi=\operatorname{\mathbf{\Theta}}\left(\frac{\left\lVert x\right\rVert^2\epsilon}{\left\lVert x\right\rVert^2+\beta^2}\right)\). Then, he boosted the success probability close to unity by repeating the quantum phase estimation \[\operatorname{\mathbf{\Theta}}\left(\frac{\left\lVert x\right\rVert^2+\beta^2}{\beta^2}\frac{\left\lVert x\right\rVert^2+\beta^2}{\left\lVert x\right\rVert^2}\right) =\operatorname{\mathbf{\Theta}}\left(\frac{\left\lVert x\right\rVert^2}{\beta^2}+\frac{\beta^2}{\left\lVert x\right\rVert^2}\right)\] times. This leads to a total query complexity of \[\operatorname{\mathbf{O}}\left(\left(\frac{\left\lVert x\right\rVert^2}{\beta^2}+\frac{\beta^2}{\left\lVert x\right\rVert^2}\right) \frac{\left(\left\lVert x\right\rVert^2+\beta^2\right)^2}{\left\lVert x\right\rVert^4\epsilon^2} \frac{\left\lVert A^{-1\dagger}x\right\rVert}{\beta} \right) =\operatorname{\mathbf{O}}\left(\left(\frac{\left\lVert x\right\rVert^2}{\beta^2}+\frac{\beta^6}{\left\lVert x\right\rVert^6}\right) \frac{\left\lVert A^{-1\dagger}x\right\rVert}{\beta\epsilon^2} \right).\]

When the solution norm \(\left\lVert x\right\rVert\) is known, Li set \(\beta=\operatorname{\mathbf{\Theta}}(\left\lVert x\right\rVert)\), simplifying the query cost to \[\operatorname{\mathbf{O}}\left( \frac{\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert\epsilon^2} \right).\] Otherwise, he chose \(\beta=1\). Under the convention that \(\left\lVert x\right\rVert\geq1\), the query complexity then becomes \[\operatorname{\mathbf{O}}\left( \frac{\left\lVert A^{-1\dagger}x\right\rVert\left\lVert x\right\rVert^2}{\epsilon^2} \right).\]

7 Analysis of beyond-\(\kappa\) solvers in the affine dilation model↩︎

To facilitate a direct comparison with prior quantum linear system solvers, we have focused on the standard input model when presenting our beyond-\(\kappa\) solvers in the main text. In this appendix, we analyze beyond-\(\kappa\) solvers in the affine dilation model, which allows further refinement of the query complexity.

7.1 Analysis of truncation-based solvers in the affine dilation model↩︎

In the affine dilation model, the input to a quantum linear system is accessed through the block encoding \[\widetilde{A}=\frac{1}{\widetilde{\alpha}}\begin{bmatrix} A & -|b\rangle\\ 0 & c \end{bmatrix}.\] whereas the dilated initial vector is \(|\widetilde{b}\rangle=\begin{bmatrix} 0\\ 1 \end{bmatrix}\). Here, we assume that \(|b\rangle\) is normalized after a potential rescaling of the block encoding. For truncation-based solvers, we take \(c>0\) so that the dilated matrix is invertible. The normalization factor \(\widetilde{\alpha}\) then satisfies \(\widetilde{\alpha}\geq\max\{\left\lVert A\right\rVert,c,1\}\).

After matrix inversion, we have \[\widetilde{A}^{-1}=\widetilde{\alpha}\begin{bmatrix} A^{-1} & \frac{1}{c}x\\ 0 & \frac{1}{c} \end{bmatrix},\] which encodes the dilated solution vector \[\widetilde{x}=\widetilde{A}^{-1}|\widetilde{b}\rangle =\frac{\widetilde{\alpha}}{c}\begin{bmatrix} x\\ 1 \end{bmatrix}.\] The performance of the truncation-based solvers is still determined by the effective condition number from Theorem [1](#thm:keff){reference-type=“ref” reference=“thm:keff”}. The only difference is that \(\kappa_{\mathrm{eff}}\) is now defined in terms of the dilated \(\widetilde{A}\) and \(|\widetilde{b}\rangle\), and we denote it \(\widetilde{\kappa}_{\mathrm{eff}}\).

It is generally challenging to give a closed-form formula for \(\widetilde{\kappa}_{\mathrm{eff}}\) analogous to the one in Theorem [1](#thm:keff){reference-type=“ref” reference=“thm:keff”}, since the singular value decomposition of the dilated matrix \(\widetilde{A}\) does not relate simply to that of \(A\). However, the effective condition number admits the succinct upper bound \(\widetilde{\kappa}_{\mathrm{eff}}\leq\frac{\left\lVert\widetilde{A}^{-1\dagger}\widetilde{x}\right\rVert}{\left\lVert\widetilde{x}\right\rVert\epsilon}\) same as in the standard input model, where \(\left\lVert\widetilde{A}^{-1\dagger}\widetilde{x}\right\rVert\) can be evaluated as \[\left\lVert\widetilde{A}^{-1\dagger}\widetilde{x}\right\rVert =\widetilde{\alpha}^2\left\lVert\begin{bmatrix} \frac{A^{-1\dagger}x}{c}\\[.25em] \frac{\left\lVert x\right\rVert^2+1}{c^2} \end{bmatrix}\right\rVert =\widetilde{\alpha}^2\sqrt{\frac{\left\lVert A^{-1\dagger}x\right\rVert^2}{c^2}+\frac{\big(\left\lVert x\right\rVert^2+1\big)^2}{c^4}}.\] This implies that the effective condition number for the dilated system is at most \[\widetilde{\kappa}_{\mathrm{eff}}\leq\frac{\widetilde{\alpha}}{\epsilon}\frac{1}{\sqrt{\left\lVert x\right\rVert^2+1}} \sqrt{\left\lVert A^{-1\dagger}x\right\rVert^2+\frac{\big(\left\lVert x\right\rVert^2+1\big)^2}{c^2}}.\]

7.2 Analysis of filtering-based solver in the affine dilation model↩︎

For the filtering-based solver, we instead consider the augmented matrix \[\widetilde{A}=\frac{1}{\alpha} \begin{bmatrix} A & -|b\rangle \end{bmatrix}.\] This corresponds to setting \(c = 0\) in the generic affine dilation model and removing the redundant zero blocks at the bottom.

Suppose that this is block encoded by an overlap of isometries as \(\widetilde{A}=\frac{1}{\alpha} \begin{bmatrix} A & -|b\rangle \end{bmatrix}=G_1^\dagger G_0\). Then by an analysis similar to that in the main text, we have the constrained orthogonal decomposition \[\begin{align} G_0\begin{bmatrix} 0\\ 1 \end{bmatrix}&=\frac{1}{\left\lVert x\right\rVert^2+1} G_0\begin{bmatrix} x\\ 1 \end{bmatrix} -\frac{\alpha}{\left\lVert x\right\rVert^2+1} \Pi_0 \Pi_1 G_1A^{-1\dagger}x\\ &\in\left(\operatorname{\mathbf{Im}}(\Pi_0)\cap\operatorname{\mathbf{Ker}}(\Pi_1)\right) \obot\operatorname{\mathbf{Im}}\left(\Pi_0 \Pi_1\right), \end{align}\] where \(\Pi_0=G_0G_0^\dagger\) and \(\Pi_1=G_1G_1^\dagger\). Performing filtering \(\operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right)\) with \(W=\left(2\Pi_0-I\right)\left(I-2\Pi_1\right)\) then gives \[\operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right)G_0\begin{bmatrix} 0\\ 1 \end{bmatrix}=\frac{1}{\left\lVert x\right\rVert^2+1} G_0\begin{bmatrix} x\\ 1 \end{bmatrix} -\frac{\alpha}{\left\lVert x\right\rVert^2+1} \operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right) \Pi_0 \Pi_1 G_1A^{-1\dagger}x,\] where \[\left\lVert\operatorname{\mathbf{F}}_{\mathbf{sv}}\left(\frac{W-I}{2}\right) \Pi_0 \Pi_1G_1A^{-1\dagger}x\right\rVert \leq\frac{\delta}{2} \left\lVert A^{-1\dagger}x\right\rVert.\]

Assuming that the solution norm \(\left\lVert x\right\rVert=\left\lVert A^{-1}|b\rangle\right\rVert\) is estimated to a constant multiplicative accuracy, we set \[\delta=\operatorname{\mathbf{\Theta}}\left(\frac{\left\lVert x\right\rVert\epsilon}{\alpha\left\lVert A^{-1\dagger}x\right\rVert}\right).\] Note that \(\delta=\operatorname{\mathbf{O}}(1)\), since \(\epsilon=\operatorname{\mathbf{O}}(1)\) and \[\frac{\alpha\left\lVert A^{-1\dagger}x\right\rVert}{\left\lVert x\right\rVert} =\alpha\left\lVert A^{-1\dagger}|x\rangle\right\rVert \geq\left\lVert A\right\rVert\left\lVert A^{-1\dagger}|x\rangle\right\rVert =\left\lVert A^\dagger\right\rVert\left\lVert A^{-1\dagger}|x\rangle\right\rVert \geq\left\lVert A^\dagger A^{-1\dagger}|x\rangle\right\rVert =1.\] Finally, we perform \(\operatorname{\mathbf{O}}\left(\frac{\left\lVert x\right\rVert^2+1}{\left\lVert x\right\rVert}\right)\) steps of amplitude amplification, so the total query complexity is determined by \[\frac{\alpha}{\epsilon}\frac{\left\lVert x\right\rVert^2+1}{\left\lVert x\right\rVert^2}\left\lVert A^{-1\dagger}x\right\rVert\] up to logarithmic factors.

7.3 Complexity comparison↩︎

We now compare the runtimes of the truncation-based and filtering-based solvers in the affine dilation model. In particular, we show that the query complexity of the former is no worse than that of the latter up to logarithmic factors, as in the standard input model.

To this end, suppose we start with a block encoding of the augmented matrix \(\frac{1}{\alpha}\begin{bmatrix} A & -|b\rangle \end{bmatrix}=G_1^\dagger G_0\). Then, we can construct a block encoding of the affine dilation \[\widetilde{A}=\frac{1}{\sqrt{\alpha^2+c^2}}\begin{bmatrix} A & -|b\rangle\\ 0 & c \end{bmatrix}\] with normalization factor \(\sqrt{\alpha^2+c^2}\). Indeed, this can be achieved through the following factorization of the dilated matrix as an overlap of isometries \[\frac{1}{\sqrt{\alpha^2+c^2}}\begin{bmatrix} A & -|b\rangle\\ 0 & c \end{bmatrix} =\begin{bmatrix} G_1^\dagger & 0\\ 0 & \begin{bmatrix} 1 & 0 & 0 \end{bmatrix} \end{bmatrix} \frac{1}{\sqrt{\alpha^2+c^2}} \begin{bmatrix} \alpha G_0\\ c\begin{bmatrix} 0 & 1\\ I & 0\\ 0 & 0 \end{bmatrix} \end{bmatrix}.\] Then the effective condition number is bounded by \[\widetilde{\kappa}_{\mathrm{eff}}\leq\frac{\sqrt{\alpha^2+c^2}}{\epsilon\sqrt{\left\lVert x\right\rVert^2+1}} \sqrt{\left\lVert A^{-1\dagger}x\right\rVert^2+\frac{\big(\left\lVert x\right\rVert^2+1\big)^2}{c^2}} \leq\frac{\sqrt{\alpha^2+c^2}}{\epsilon}\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+1}} +\frac{\sqrt{\alpha^2+c^2}}{\epsilon}\frac{\sqrt{\left\lVert x\right\rVert^2+1}}{c}.\]

Let us choose \[c=\alpha\frac{\left\lVert x\right\rVert^2+1}{\left\lVert x\right\rVert} \geq2\alpha,\] so that \[\sqrt{\alpha^2+c^2} =\operatorname{\mathbf{O}}\left(\alpha\frac{\left\lVert x\right\rVert^2+1}{\left\lVert x\right\rVert}\right).\] Meanwhile, by our block encoding convention, \[1=\left\lVert|b\rangle\right\rVert=\left\lVert Ax\right\rVert \leq\alpha\left\lVert x\right\rVert \leq\alpha^2\left\lVert A^{-1\dagger}x\right\rVert,\] which gives \[c=\frac{\left\lVert x\right\rVert^2+1}{\frac{\left\lVert x\right\rVert}{\alpha}} \geq\frac{\left\lVert x\right\rVert^2+1}{\left\lVert A^{-1\dagger}x\right\rVert}.\]

Altogether, \[\widetilde{\kappa}_{\mathrm{eff}} =\operatorname{\mathbf{O}}\left(\frac{\alpha}{\epsilon}\frac{\left\lVert x\right\rVert^2+1}{\left\lVert x\right\rVert}\frac{\left\lVert A^{-1\dagger}x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+1}}\right) =\operatorname{\mathbf{O}}\left(\frac{\alpha}{\epsilon}\frac{\left\lVert x\right\rVert^2+1}{\left\lVert x\right\rVert^2}\left\lVert A^{-1\dagger}x\right\rVert\right).\] Hence, with this specific choice of \(c\) in the affine dilation model, the runtime of the truncation-based algorithm matches that of the filtering-based method up to logarithmic factors.

8 Generalized truncation property↩︎

In the main text, we have defined the strong and weak truncation property (Definition [[defn:strong_trunc]](#defn:strong_trunc){reference-type=“ref” reference=“defn:strong_trunc”} and Definition [[defn:weak_trunc]](#defn:weak_trunc){reference-type=“ref” reference=“defn:weak_trunc”}) for a quantum linear system solver. Informally, these capture the requirement that if one supplies \(\alpha\) as the condition number parameter, the solver correctly inverts singular value components above \(1/\alpha\) and discards those below. To obtain a beyond-\(\kappa\) solver, we have required that the truncation properties hold when \(\alpha\) is at least the effective condition number from Theorem [1](#thm:keff){reference-type=“ref” reference=“thm:keff”}. More generally, we can define the strong truncation property for an arbitrary \(\alpha\) as follows.

When configured with \(\alpha\), VTAA correctly inverts all singular values above \(1/\alpha\). Combining this with a standard amplitude amplification yields the truncated solution state supported on singular values above \(1/\alpha\). We thus conclude that VTAA satisfies the generalized strong truncation property with a query overhead from amplification. Similarly, one can show that VTAA satisfies the following generalized version of weak truncation property.

9 Constant-prefactor optimization for filtering with effective gap↩︎

In this appendix, we optimize the constant prefactor of the query complexity for eigenstate filtering with effective gap.

Recall from the main text (Section [[sec:gap_solver]](#sec:gap_solver){reference-type=“ref” reference=“sec:gap_solver”}) that, in order to prepare the normalized solution state with accuracy \(\epsilon\), we choose the effective gap \(\delta\) and filtering accuracy \(\xi\) so that \[\frac{\sqrt{\beta^2+1}}{\left\lVert x\right\rVert^2+\beta^2}\left\lVert A^{-1\dagger}x\right\rVert\frac{\delta}{2} +\frac{\left\lVert x\right\rVert}{\sqrt{\left\lVert x\right\rVert^2+\beta^2}}\xi \leq\frac{\beta\left\lVert x\right\rVert}{\left\lVert x\right\rVert^2+\beta^2}\epsilon.\] Correspondingly, the query complexity of filtering is \[\operatorname{\mathbf{Ceil}}\left(\frac{4}{\delta}\ln\left(\frac{2}{\xi}\right)\right).\] A naive strategy is to divide the error budget into equal halves. Here, we propose an alternative approach whose leading-order query complexity has a smaller constant prefactor.

To this end, we consider the function \[g(\gamma)=\frac{1}{\gamma\epsilon}\ln\left(\frac{1}{(1-\gamma)\epsilon}\right),\] where \(\epsilon\) is a fixed positive number and \(0<\gamma<1\) is a tunable parameter. See Figure 5 for an illustration of its qualitative behavior. This function corresponds to the query cost when we allocate a \(\gamma\)-fraction of the error budget to the effective gap and the remaining fraction to the filtering accuracy parameter. Intuitively, since \(\ln(\cdot)\) grows much slower than linearly, we expect the minimum of \(g(\gamma)\) to occur at some \(\gamma=\gamma_0\) close to \(1\).

To rigorize this intuition, we differentiate \(g(\gamma)\) to get \[\frac{\operatorname{d}}{\operatorname{d}\gamma}g(\gamma) =-\frac{(\gamma-1)\ln\left(\frac{1}{\epsilon(1-\gamma)}\right)+\gamma}{\epsilon(\gamma-1)\gamma^2}.\] Thus \(g(\gamma)\) is minimized at \(\gamma=\gamma_0\) where \(\frac{\operatorname{d}}{\operatorname{d}\gamma}g(\gamma_0)=0 \Leftrightarrow (\gamma_0-1)\ln\left(\frac{1}{\epsilon(1-\gamma_0)}\right)+\gamma_0=0\). Defining \(\upsilon_0=\gamma_0-1\) with \(-1<\upsilon_0<0\), we have \(\upsilon_0\ln\left(\frac{1}{-\epsilon\upsilon_0}\right)+\upsilon_0+1=0 \Leftrightarrow \frac{1}{\upsilon_0}e^{\frac{1}{\upsilon_0}} =-\frac{\epsilon}{e}\). Finally, letting \(\zeta_0=\frac{1}{\upsilon_0}\) with \(-\infty<\zeta_0<-1\), we get \(\zeta_0 e^{\zeta_0}=-\frac{\epsilon}{e} \Leftrightarrow \zeta_0=\operatorname{\mathbf{W}_{-1}}\left(-\frac{\epsilon}{e}\right)\), for \(\operatorname{\mathbf{W}_{-1}}(\cdot)\) the (\(-1\))-branch of the Lambert-W function. This implies \[\gamma_0=1+\frac{1}{\operatorname{\mathbf{W}_{-1}}\left(-\frac{\epsilon}{e}\right)}.\]

Now let us bound \(g(\gamma)\) at \(\gamma=\gamma_0\). To this end, we rewrite \(\operatorname{\mathbf{W}_{-1}}\left(-\frac{\epsilon}{e}\right) =\operatorname{\mathbf{W}_{-1}}\left(-e^{\ln\left(\frac{\epsilon}{e}\right)}\right) =\operatorname{\mathbf{W}_{-1}}\left(-e^{-\ln\left(\frac{1}{\epsilon}\right)-1}\right)\). Applying bounds on the (\(-1\))-branch of the Lambert-W function [43], we have \[-1-\sqrt{2\ln\left(\frac{1}{\epsilon}\right)}-\ln\left(\frac{1}{\epsilon}\right) <\operatorname{\mathbf{W}_{-1}}\left(-e^{-\ln\left(\frac{1}{\epsilon}\right)-1}\right) <-1-\sqrt{2\ln\left(\frac{1}{\epsilon}\right)}-\frac{2}{3}\ln\left(\frac{1}{\epsilon}\right).\] This implies \[g(\gamma_0) <\left(\frac{1+\sqrt{2\ln\left(\frac{1}{\epsilon}\right)}+\frac{2}{3}\ln\left(\frac{1}{\epsilon}\right)}{\sqrt{2\ln\left(\frac{1}{\epsilon}\right)}+\frac{2}{3}\ln\left(\frac{1}{\epsilon}\right)}\right) \frac{1}{\epsilon} \ln\left(\frac{1+\sqrt{2\ln\left(\frac{1}{\epsilon}\right)}+\ln\left(\frac{1}{\epsilon}\right)}{\epsilon}\right),\] which reduces to \[\frac{1}{\epsilon}\ln\left(\frac{1}{\epsilon}\right)\] when higher order terms are dropped.

References↩︎