# An Efficient Quantum Algorithm for Linear System Problem in Tensor Format

Zeguan Wu  Lehigh University  Los Alamos National Laboratory  Sidhant Misra
Los Alamos National Laboratory  Tamás Terlaky  Lehigh University  Xiu Yang
Lehigh University  Marc Vuffray  Los Alamos National Laboratory

###### Abstract

Solving linear systems is at the foundation of many algorithms. Recently,
quantum linear system algorithms (QLSAs) have attracted great attention since
they converge to a solution exponentially faster than classical algorithms in
terms of the problem dimension. However, low-complexity circuit
implementations of the oracles assumed in these QLSAs constitute the major
bottleneck for practical quantum speed-up in solving linear systems. In this
work, we focus on the application of QLSAs for linear systems that are
expressed as a low rank tensor sums, which arise in solving discretized PDEs.
Previous works uses modified Krylov subspace methods to solve such linear
systems with a per-iteration complexity being polylogarithmic of the dimension
but with no guarantees on the total convergence cost. We propose a quantum
algorithm based on the recent advances on adiabatic-inspired QLSA and perform
a detailed analysis of the circuit depth of its implementation. We rigorously
show that the _total complexity_ of our implementation is polylogarithmic in
the dimension, which is comparable to the _per-iteration complexity_ of the
classical heuristic methods.

Keywords— quantum computing; linear system; tensor format

##  1 Introduction

Quantum computing has obtained great attention for its potential in solving
certain problems faster than classical algorithms. Starting from Deutsch’s
theory, see [15], a large number of quantum algorithms have been proposed,
including Shor’s algorithm for integer factorization, [29], Grover’s algorithm
for Database Search, [19], Quantum Approximate Optimization Algorithm (QAOA)
for combinatorial optimization problems, [17], QLSAs for solving quantum
linear system problem (QLSP), [20, 11, 9, 30, 13, 2, 22], and many other
algorithms. We refer readers to the survey paper [14] for more details. Among
these quantum algorithms, some provide polynomial speed-up while some provide
exponential speed-up. QLSAs are a family of algorithms that show exponential
convergence speed-up in terms of problem dimension when compared with
classical methods including Cholesky factorization, [20, 11, 9, 30, 13, 2,
22]. Due to the significance of solving linear systems in numerical
computation, a significant amount of research have been invested in using
QLSAs to speed up classical algorithms. These efforts include using QLSAs to
speed up classical PDE algorithms [24, 23, 21] and classical optimization
algorithms [8, 31, 26, 4, 3]. We refer the reader to the survey paper [1] for
more details on quantum optimization algorithms. The main caveat in speeding
up classical algorithms with QLSAs is the efficient quantum circuit
implementation of the oracles assumed for preprocessing the input data and
building quantum circuit from the classical system description. Finding such
efficient implementations is an area of active research and is crucial for
building quantum-classical hybrid algorithms that can provide quantum speed-
up.

In this work, we study the use of QLSAs for solving a class of linear systems,
whose coefficient matrix and right-hand-side vector can be represented as a
linear combination of a few tensor product of 2-by-2 matrices and
2-dimensional vectors, respectively. The problem is a special case of the so-
called linear system problem in tensor format (LSP-TF), where matrices/vectors
in tensor product are not necessarily 2-dimensional. LSP-TF is frequently
encountered in discretized PDE problems [18]. The problem size grows
exponentially as the length of the tensor product chain grows, which makes it
difficult to solve using general classic linear system solvers, including
Cholesky factorization and Krylov subspaces methods. Some classical algorithms
have been proposed to solve such LSP-TF [5]. The main idea of these algorithms
is to modify the Krylov subspaces methods by taking into account the tensor
format of the problem. The resulted Krylov subspaces are constructed
approximately, which sacrifices accuracy while keeps the complexity in each
iteration polylogarithmic of the problem size. Despite the polylogarithmic
per-iteration complexity of these algorithms, their total complexity is
unknown and is unlikely to be better than that of the original Krylov
subspaces for the Krylov subspaces in these algorithms are constructed
approximately. Our main contribution is to show that such linear systems can
be efficiently solved using a quantum computer and we provide a full and
explicit circuit implementation of the algorithm.

##  2 Problem Definition

In this section, we start with notations and then introduce LSP-TF, QLSP, and
the Trotterization method. Finally, we summarize our contributions.

###  2.1 Notation

Vectors are denoted by lower case letters and matrices are denoted by upper
case letters. We use $e_{i}$ to denote the unit basis vector with the $i$th
entry being $1$. We use $I_{n}$ to denote the identity matrix with dimension
$n\times n$, or simply $I$ if the dimension is obvious from the context.
Single-qubit Pauli matrices are $\\{I,\ X,\ Y,Z\\}$, where

$\displaystyle X=\begin{bmatrix}0&1\\\ 1&0\end{bmatrix},\
Y=\begin{bmatrix}0&-i\\\ i&0\end{bmatrix},\ Z=\begin{bmatrix}1&0\\\
0&-1\end{bmatrix}.$

We use

  * •

$\|\cdot\|_{1}$ to denote the trace norm of matrices;

  * •

$\|\cdot\|_{2}$ to denote the $\ell_{2}$ norm for vectors and spectral norm
for matrices;

  * •

$\|\cdot\|_{F}$ to denote the Frobenius norm for matrices.

The condition number of a general matrix $M$ is denoted by $\kappa_{M}$. Let
$p$ be a positive integer. We use $[p]$ to represent the set
$\\{1,2,\dots,p\\}$.

We use $\ket{\cdot}$ to represent quantum state, which can be taken as a
column vector in this work. Let $\psi=(\psi_{1},\dots,\psi_{N})$ be a column
vector. We use
$\ket{\psi}=\sum_{i=1}^{N}\psi_{i}\ket{i}/\sqrt{\sum_{i=1}^{N}|\psi_{i}|^{2}}$
to denote the quantum state representation of $\psi$. We also take the
convention that

$\ket{0}=\begin{bmatrix}1\\\ 0\end{bmatrix}{\rm
and}\ket{1}=\begin{bmatrix}0\\\ 1\end{bmatrix}.$

###  2.2 LSP-TF

The general linear system problem (LSP) is to find $x\in\mathbb{R}^{N}$ such
that

$Ax=b,$ (LSP)

where $A\in\mathbb{R}^{N\times N}$ and $b\in\mathbb{R}^{N}$. For matrix $A$,
its condition number is denoted by $\kappa_{A}$. For simplicity, we set
$\kappa=\kappa_{A}$. The LSP-TF problem is to find $x\in\mathbb{R}^{N}$ such
that

$Ax=b,\ A=\sum_{i=1}^{m}\otimes_{k=1}^{n}A_{ik},\
b=\sum_{j=1}^{d}\otimes_{k=1}^{n}b_{jk},$ (LSP-TF)

where $A\in\mathbb{R}^{N\times N}$ and $b\in\mathbb{R}^{N}$. This type of
problems arises frequently in high-dimensional discretized PDEs [5]. In [5], a
combination of classical projection method and low-rank tensor format
approximation was proposed to solve LSP-TF. Their main idea of their method is
to use low rank tensors to construct subspaces and maintain the low rank
tensor structure using hierarchical Tucker format. Their method obtains per-
iteration complexity $\mathcal{O}(mn)$ while assuming $d$ is small. We refer
readers to the survey paper [18] on these low tensor rank iterative methods.

###  2.3 Contributions

In this work, we apply Algorithm 1 to solve problem LSP-TF. We focus on the
special case when $A=\sum_{i=1}^{m}\otimes_{k=1}^{n}A_{ik}$ and
$b=\sum_{j=1}^{d}\otimes_{k=1}^{n}b_{jk}$ with $A_{ik}\in\mathbb{R}^{2\times
2}$ being Hermitian and $b_{jk}\in\mathbb{R}^{2}$. Without loss of generality,
we make the following assumption.

######  Assumption 2.1.

$\|A\|_{2}\leq 1$, $\|\otimes_{k=1}^{n}A_{ik}\|_{2}\leq 1$ for all $i\in[m]$,
and $\|\otimes_{k=1}^{n}b_{jk}\|_{2}\leq 1$ for all $j\in[d]$.

We also assume $A$ and $b$ are sparse in tensor product strings.

######  Assumption 2.2.

$m=\mathcal{O}(1)$ and $d=\mathcal{O}(1)$.

We propose a Hamiltonian decomposition for the Hamiltonian used in Algorithm 1
when solving LSP-TF and a detailed circuit implementation for the simulation
of the decomposed Hamiltonians. We combine the circuits and the Trotterization
method to implement Algorithm 1. We show the circuit depth in our
implementation of the QLSA is polylogarithmic of the problem dimension. Here
we provide a simplified version of our main result. The full statement is
provided in Theorem 4.1.

######  Theorem 2.1.

Let $0<\epsilon\leq 1/(3n)$,
$\epsilon_{0}=\epsilon^{2}/(\kappa\log^{2}\kappa)$ and $p=\lceil
1-\log_{2}\epsilon_{0}\rceil$. Our implementation of Algorithm 1 prepares an
$\mathcal{O}(\epsilon)$-approximate solution using $\mathcal{T}$ classical
arithmetic operations, $\mathcal{T}$ single qubit unitary circuits  and their
controlled version with gate depth $\mathcal{O}(1)$, and $\mathcal{T}$ calls
to $p$-qubit quantum multiplier  and their controlled version with circuit
depth $\mathcal{O}(\mathcal{T}\mathcal{T}_{p})$, where
$\mathcal{T}=\mathcal{O}({\rm poly}(\kappa\log(N)/\epsilon))$ and
$\mathcal{T}_{p}=\mathcal{O}(p^{3})$ is the gate depth of $p$-qubit quantum
multiplier.

######  Remark 2.1.

The final result obtained from our implementation is a quantum state. One
needs to apply quantum tomography algorithm to read out any entries.

The results suggest that the time complexity of our implementation is
polylogarithmic in the problem dimension if the time complexity of each gate
is $\mathcal{O}(1)$. This is an exponential speed-up compared with classical
general linear system algorithms, including Gaussian elimination, Cholesky
factorization, and Krylov subspaces methods. When compared with classical
algorithms designed for LSP-TF, for example the algorithm introduced in [5],
our total time complexity is comparable to their per-iteration complexity.

In the remaining of this section, we introduce QLSA and Trotterization method.

###  2.4 QLSA

In this section, we start with the definition of QLSP and give a brief
introduction of the QLSA proposed in [30].

######  Definition 2.1 (QLSP).

Given a Hermitian matrix $A\in\mathbb{R}^{N\times N}$ and
$b\in\mathbb{R}^{N}$, the QLSP is to find an approximation of the quantum
state

$\ket{x}=\frac{\sum_{i=1}^{N}x_{i}\ket{i}}{\sqrt{\sum_{i=1}^{N}|x_{i}|^{2}}},$

where $(x_{1},\dots,x_{N})^{\top}=A^{-1}b$. Quantum state $\ket{\tilde{x}}$ is
called an $\epsilon$-approximate solution to the QLSP if

$\left\|\ket{\tilde{x}}-\ket{x}\right\|_{2}\leq\epsilon.$

In [20], the authors propose the first QLSA to solve QLSP with polylogarithmic
dependence on dimension when the problem is sparse. Later, researchers
proposed different QLSAs with the sparsity condition relaxed and better
dependence on other parameters including condition number of linear system and
solution accuracy [11, 9, 30, 2, 13, 22]. In [30], the authors propose a QLSA
inspired by adiabatic quantum computing to solve QLSP. They study two types of
Hamiltonians and use the randomization method introduced by [7] to design two
discretized adiabatic-like quantum algorithms. Algorithm 1 is the one of the
two QLSAs with better complexity. Later in [22], an improved version of
Algorithm 1 is proposed with better dependence on the solution accuracy. In
this work, we opt to study Algorithm 1 for its simplicity. For the remainder
of this section,we give a brief explanation of Algorithm 1 and we refer the
reader to [30] for further details.

Given a linear system problem $Ax=b$, let

$\displaystyle A(s)=(1-s)Z\otimes I+sX\otimes A$

and

$\displaystyle\ket{\bar{b}}$
$\displaystyle=\frac{\sqrt{2}}{2}(\ket{0}+\ket{1})\otimes\ket{b},\
P_{\bar{b}}^{\perp}$ $\displaystyle=I-\ket{\bar{b}}\bra{\bar{b}},$

where $s\in[0,1]$ is the evolution parameter. With these, they introduce the
following family of quantum states

$\ket{x(s)}=\frac{\sum_{i=1}^{2N}x_{i}(s)\ket{i}}{\sqrt{\sum_{i=1}^{2N}|x_{i}(s)|^{2}}},$

where $(x_{1}(s),\dots,x_{2N}(s))^{\top}=A(s)^{-1}\bar{b}$. They then study
the following Hamiltonian

$\displaystyle H(s)=\left(\frac{X+iY}{2}\right)\otimes
A(s)P_{\bar{b}}^{\perp}+\left(\frac{X-iY}{2}\right)\otimes
P_{\bar{b}}^{\perp}A(s),\ s\in[0,1].$ (1)

They prove that, if one starts with quantum state $\ket{0}\otimes\ket{x(0)}$
and lets the quantum system evolve adiabatically, then the final state will be
$\ket{0}\otimes\ket{x(1)}$ with sufficiently high probability. Notice that, in
our framework, the initial state is
$\frac{\sqrt{2}}{2}(\ket{0}-\ket{1})\otimes\ket{b}$, where $\ket{b}$ has low
tensor rank. This allows us to prepare the initial state efficiently, see the
proof of Theorem 4.1 for details. They also prove that the relevant spectral
gap of the Hamiltonian is bounded from below by the following quantity

$\sqrt{\Delta^{\ast}(s)}=\sqrt{(1-s)^{2}+(s/\kappa_{A})^{2}}.$

They propose to use the phase randomization method proposed in [7] to turn the
continuous adiabatic evolution into piece-wise time-independent evolution. The
continuous evolution parameter $s\in[0,1]$, which is time dependent, is
discretized into $q$ discretization points $s_{j}$ for $j\in[q]$. To determine
the discretization of $s$, they parametrize $s$ using

$s(v)=\frac{e^{v(\sqrt{1+\kappa^{2}}/\sqrt{2\kappa^{2}})}+2\kappa^{2}-\kappa^{2}e^{-v(\sqrt{1+\kappa^{2}}/\sqrt{2\kappa^{2}})}}{2(1+\kappa^{2})}.$
(2)

Then they discretize $s$ by uniformly discretizing the value of $v$ between
$v_{a}$ and $v_{b}$, where

$\displaystyle v_{a}$
$\displaystyle=\frac{\sqrt{2}\kappa}{\sqrt{1+\kappa^{2}}}\log\left(\kappa\sqrt{1+\kappa^{2}}-\kappa^{2}\right)$
$\displaystyle v_{b}$
$\displaystyle=\frac{\sqrt{2}\kappa}{\sqrt{1+\kappa^{2}}}\log\left(\sqrt{1+\kappa^{2}}+1\right).$

At each discretization point $s_{j}$, the quantum system is evolved using
Hamiltonian $H(s_{j})$ for time $t_{j}$, which is sampled from the uniform
distribution $t_{j}\sim[0,2\pi/\sqrt{\Delta^{\ast}(s_{j})}]$. The pseudocode
of their algorithm is provided in Algorithm 1.

Algorithm 1 QLSA for QLSP by [30]

1:Given $A,\ b,\ \kappa,$ and $\epsilon$

2:Compute $v_{a}$ and $v_{b}$ by Eq. (2.4)

3:Set $q=\Theta\left(\log^{2}(\kappa)/\epsilon\right)$ and
$\delta=(v_{b}-v_{a})/q$

4:For $j=1,\dots,q$, let $v_{j}=v_{a}+j\delta,\ s_{j}=s(v_{j})$, and $t_{j}$
be sampled uniformly from $[0,2\pi/\sqrt{\Delta^{\ast}(s_{j})}]$

5:Apply $e^{-it_{q}H(s_{q})}\cdots e^{-it_{1}H(s_{1})}$ to
$\ket{0}\otimes\ket{x(0)}$, discard the ancilla

As mentioned in [30], the algorithm is inspired by adiabatic quantum computing
and thus its run time can be measured by the evolution time needed. The time
complexity is $T=\mathcal{O}(\kappa\log(\kappa)/\epsilon)$. However, as
pointed out by the authors, the Hamiltonian $H(s)$ easily encodes the
information of QLSP but not necessarily corresponds to any physical model,
which makes the implementation of the algorithm on analog quantum computers
difficult. The authors then introduce a gate-model implementation of Algorithm
1 and analyze its complexity. With the oracle access to problem information
$A$ and $b$, the authors obtain the query complexity
$\tilde{\mathcal{O}}(d\kappa/\epsilon)$, with polylogarithmic factors of
$d\kappa/\epsilon$ hidden.

Despite the achievement on complexities, Algorithm 1 assumes the existence of
oracle access to $A$, $b$, and several other non-trivial unitaries.
Constructing these oracles from classical input data $(A,b)$ could be more
difficult than solving the problem using the QLSA. In this work, we propose a
simpler implementation of Algorithm 1 and estimate the implementation cost
when it is applied to LSP-TF.

###  2.5 Trotterization Method

In this section, we give a brief introduction to the Trotterization method and
its impact on Hamiltonian evolution. The detailed impact of Trotterization
method on Algorithm 1 is analyzed in Section 4.2.

Consider the time evolution of a time-independent Hamiltonian $M$ , i.e.,
$e^{-iMt}.$ Since $M$ is a general Hermitian matrix, it is difficult to design
the circuit for $e^{-iMt}$. The main idea of Trotterization methods is that
one can consider a decomposition of Hamiltonian $M$, e.g.,

$M=\sum_{\gamma=1}^{\Gamma}M_{\gamma},$

where the time evolution of each $M_{j}$ are are assumed to be easy to
compile. Then, for short enough evolution time, one can approximate the
Hamiltonian evolution using

$\displaystyle
e^{-it\sum_{\gamma=1}^{\Gamma}M_{\gamma}}\approx\prod_{\gamma=1}^{\Gamma}e^{-itM_{\gamma}},$

which is the first-order Lie-Trotter formula. The accuracy of the
Trotterization methods is extensively studied in the literature, see e.g.,
[25, 6, 10, 12]. Here we cite Theorem 6 from [12] on the accuracy of general
Trotterization methods. We restate it by considering $e^{-iMt}$ instead of
$e^{Mt}$ and only for a Hermitian matrix $M$.

######  Lemma 2.1 (Restated Theorem 6 of [12]).

Let $M=\sum_{\gamma=1}^{\Gamma}M_{\gamma}$ be an Hermitian operator consisting
of $\Gamma$ summands with each $M_{\gamma}$ Hermitian and $t\geq 0$. Let
$\mathcal{S}(t)$ be a $p$-th order $\Upsilon$-stage product formula as

$\displaystyle\mathcal{S}(t)=\prod_{\upsilon=1}^{\Upsilon}\prod_{\gamma=1}^{\Gamma}e^{-ita_{\upsilon,\gamma}M_{\upsilon,\gamma}},$

where $a_{\upsilon,\gamma}$ and $M_{\upsilon,\gamma}$ are to be specified by
certain product formula. Define $\tilde{\alpha}_{\rm
comm}=\sum_{\gamma_{1}=1}^{\Gamma}\cdots\sum_{\gamma_{p+1}=1}^{\Gamma}\left\|\left[M_{\gamma_{p+1}},\cdots\left[M_{\gamma_{2}},M_{\gamma_{1}}\right]\cdots\right]\right\|_{2}$.
Then, the additive error $\mathcal{E}_{\mathcal{A}}(t)$ and the multiplicative
error $\mathcal{E}_{\mathcal{M}}(t)$, defined, respectively, by
$\mathcal{S}(t)=e^{-iMt}+\mathcal{E}_{\mathcal{A}}(t)$ and
$\mathcal{S}(t)=e^{-iMt}\left(I+\mathcal{E}_{\mathcal{M}}(t)\right)$, can be
asymptotically bounded as

$\displaystyle\left\|\mathcal{E}_{\mathcal{A}}(t)\right\|_{2},\left\|\mathcal{E}_{\mathcal{M}}(t)\right\|_{2}=\mathcal{O}\left(\tilde{\alpha}_{\rm
comm}t^{p+1}\right).$

######  Remark 2.2.

In Lemma 2.1, $a_{\upsilon,\gamma}$ and $M_{\upsilon,\gamma}$ are determined
in the corresponding Trotterization formula. In this paper, we work with the
aforementioned first-order Lie-Trotter formula and thus

$\Upsilon=1,\ a_{\upsilon,\gamma}=1,\ M_{\upsilon,\gamma}=M_{\gamma},$

and

$\displaystyle\mathcal{S}(t)=\prod_{\gamma=1}^{\Gamma}e^{-itM_{\gamma}}.$

The commutator scaling factor becomes

$\displaystyle\tilde{\alpha}_{\rm comm}$
$\displaystyle=\sum_{\gamma_{1}=1}^{\Gamma}\sum_{\gamma_{2}=1}^{\Gamma}\left\|[M_{\gamma_{1}},M_{\gamma_{2}}]\right\|_{2}$
$\displaystyle=\sum_{\gamma_{1}=1}^{\Gamma}\sum_{\gamma_{2}=1}^{\Gamma}\left\|M_{\gamma_{1}}M_{\gamma_{2}}-M_{\gamma_{2}}M_{\gamma_{1}}\right\|_{2}$
$\displaystyle\leq\sum_{\gamma_{1}=1}^{\Gamma}\sum_{\gamma_{2}=1}^{\Gamma}\left\|M_{\gamma_{1}}M_{\gamma_{2}}\right\|_{2}+\left\|M_{\gamma_{2}}M_{\gamma_{1}}\right\|_{2}$
$\displaystyle\leq
2\sum_{\gamma_{1}=1}^{\Gamma}\sum_{\gamma_{2}=1}^{\Gamma}\left\|M_{\gamma_{1}}\right\|_{2}\left\|M_{\gamma_{2}}\right\|_{2}.$

Lemma 2.1 applies for general Hamiltonians and thus tighter bounds are
possible for well-structured Hamiltonians. In Section 4.2, we discuss the
bound of the commutator scaling factor $\tilde{\alpha}_{\rm comm}$ when the
first-order Lie-Trotter formula is applied in the QLSA for LSP-TF. With
$\tilde{\alpha}_{\rm comm}$ quantified, one can find the Trotter number. We
conclude this section by restating Corollary 7 on Trotter number in [12].

######  Lemma 2.2 (Restated Corollary 7 in [12]).

Let $M=\sum_{\gamma=1}^{\Gamma}M_{\gamma}$ be a Hermitian operator consisting
of $\Gamma$ summands with each $M_{\gamma}$ Hermitian and $t\geq 0$. Let
$\mathcal{S}(t)$ be a $p$-th order $\Upsilon$-stage product formula as

$\displaystyle\mathcal{S}(t)=\prod_{\upsilon=1}^{\Upsilon}\prod_{\gamma=1}^{\Gamma}e^{-ita_{\upsilon,\gamma}M_{\upsilon,\gamma}}.$

Define $\tilde{\alpha}_{\rm
comm}=\sum_{\gamma_{1}=1}^{\Gamma}\cdots\sum_{\gamma_{p+1}=1}^{\Gamma}\left\|\left[M_{\gamma_{p+1}},\cdots\left[M_{\gamma_{2}},M_{\gamma_{1}}\right]\cdots\right]\right\|_{2}$.
Then, we have $\|\mathcal{S}^{r}(t/r)-e^{-itM}\|_{2}=\mathcal{O}(\epsilon)$,
provided that

$r=\mathcal{O}\left(\frac{\tilde{\alpha}_{\rm
comm}^{1/p}t^{1+1/p}}{\epsilon^{1/p}}\right).$

The remaining of the paper is organized as follows. In Section 3, we describe
our circuit design in the implementation. In Section 4, we analyze the cost to
implement the proposed circuits, and finally give the total cost of the
implementation.

##  3 QLSA Circuit Design for LSP-TF

In this section, we show how to decompose Hamiltonian (1) into two types of
structured Hamiltonians when Algorithm 1 is applied on LSP-TF. Then we show
how to implement the simulation of these two types of Hamiltonians and give
their cost estimations.

Following the notation introduced in Section 2.2 and Section 2.4, we have

$\displaystyle A(s)$
$\displaystyle=(1-s)Z\otimes(\otimes_{k=1}^{n}I)+sX\otimes\sum_{i=1}^{m}\otimes_{k=1}^{n}A_{ik}$
$\displaystyle\ket{\bar{b}}$
$\displaystyle=\frac{1}{\sqrt{2}}\begin{bmatrix}b\\\
b\end{bmatrix}=\frac{1}{\sqrt{2}}\begin{bmatrix}\sum_{j=1}^{d}\otimes_{k=1}^{n}b_{jk}\\\
\sum_{j=1}^{d}\otimes_{k=1}^{n}b_{jk}\end{bmatrix}$ $\displaystyle
P_{\bar{b}}^{\perp}$
$\displaystyle=I_{2N}-\ket{\bar{b}}\bra{\bar{b}}=I_{2N}-\frac{1}{2}\begin{bmatrix}bb^{\dagger}&bb^{\dagger}\\\
bb^{\dagger}&bb^{\dagger}\end{bmatrix},$

where

$bb^{\dagger}=\sum_{j_{1},j_{2}=1}^{d}b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}=\sum_{j_{1},j_{2}=1}^{d}\otimes_{k=1}^{n}b_{j_{1}k}b_{j_{2}k}^{\dagger}.$

Here $b_{j_{1}\cdot}$ and $b_{j_{2}\cdot}$ represent the corresponding tensor
product strings. We also use the similar notation $A_{i\cdot}$. With these
notations, the Hamiltonian $H(s)$ can be written as

$\displaystyle H(s)$ $\displaystyle=\left(\frac{X+iY}{2}\right)\otimes
A(s)P_{\bar{b}}^{\perp}+\left(\frac{X-iY}{2}\right)\otimes
P_{\bar{b}}^{\perp}A(s)$ (3)
$\displaystyle=\begin{bmatrix}0&A(s)P_{\bar{b}}^{\perp}\\\
P_{\bar{b}}^{\perp}A(s)&0\end{bmatrix}$
$\displaystyle=H_{1}(s)+H_{2}(s)+H_{3}(s)+H_{4}(s),$

where

$\displaystyle H_{1}(s)$ $\displaystyle=(1-s)X\otimes
Z\otimes(\otimes_{l=1}^{n}I)$ $\displaystyle H_{2}(s)$
$\displaystyle=sX\otimes X\otimes A=s\sum_{i=1}^{m}X\otimes
X\otimes(\otimes_{l=1}^{n}A_{il})$ $\displaystyle H_{3}(s)$
$\displaystyle=-\frac{1}{2}(1-s)\begin{bmatrix}0&(Z+iY)\otimes bb^{\dagger}\\\
(Z+iY)^{\dagger}\otimes bb^{\dagger}&0\end{bmatrix}$
$\displaystyle=\frac{1}{2}(1-s)\sum_{j_{1},j_{2}}\left(\begin{bmatrix}0&Z\otimes
b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}\\\ Z\otimes
b_{j_{2}\cdot}b_{j_{1}\cdot}^{\dagger}&0\end{bmatrix}+\begin{bmatrix}0&iY\otimes
b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}\\\ (iY)^{\dagger}\otimes
b_{j_{2}\cdot}b_{j_{1}\cdot}^{\dagger}&0\end{bmatrix}\right)$ $\displaystyle
H_{4}(s)$ $\displaystyle=-\frac{1}{4}s(X+iY)\otimes(I+X)\otimes
Abb^{\dagger}-\frac{1}{4}s(X-iY)\otimes(I+X)\otimes bb^{\dagger}A$
$\displaystyle=\frac{1}{2}s\sum_{ij_{1}j_{2}}\left(\begin{bmatrix}0&X\otimes\left(A_{i\cdot}b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}\right)\\\
X\otimes\left(b_{j_{2}\cdot}b_{j_{1}\cdot}^{\dagger}A_{i\cdot}\right)&0\end{bmatrix}+\begin{bmatrix}0&I\otimes\left(A_{i\cdot}b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}\right)\\\
I\otimes\left(b_{j_{2}\cdot}b_{j_{1}\cdot}^{\dagger}A_{i\cdot}\right)&0\end{bmatrix}\right).$

The summands in the four Hamiltonians can be categorized into the two types of
Hamiltonian defined in the following. Let us define the two types of
Hamiltonian

$\displaystyle\text{Type-1: }H^{1}$ $\displaystyle=P_{1}\otimes P_{2}\otimes
C_{1}\otimes\cdots\otimes C_{n}$ $\displaystyle\text{Type-2: }H^{2}$
$\displaystyle=\begin{bmatrix}0&P_{0}\otimes D_{1}\otimes\cdots\otimes
D_{n}\\\ P_{0}^{\dagger}\otimes D_{1}^{\dagger}\otimes\cdots\otimes
D_{n}^{\dagger}&0\end{bmatrix},$

where $P_{j}\in\\{I,X,Y,Z\\}$, $C_{i}$ and $D_{i}$ are both
$\mathbb{C}^{2\times 2}$, and $C_{i}$ is Hermitian. It is obvious that the
summands of $H_{1}(s)$ and $H_{2}(s)$ are Type 1 Hamiltonians; the summands of
$H_{3}(s)$ and $H_{4}(s)$ are Type 2 Hamiltonians. One can also verify that
the total amount of such summands in $H(s)$ is a polynomial of $m$ and $d$.

######  Lemma 3.1.

$H(s)$ can be represented as the summation of $(m+1)$ Type-1 Hamiltonians and
$(2d^{2}+2md^{2})$ Type-2 Hamiltonians.

###### Proof.

We see that $H_{1}(s)$ is Type-1; $H_{2}(s)$ is summation of $m$ Type-1
Hamiltonians; $H_{3}(s)$ is summation of $2d^{2}$ Type-2 Hamiltonians; and
$H_{4}(s)$ is summation of $2md^{2}$ Type-2 Hamiltonians. ∎

######  Remark 3.1.

The fact that $H(s)$ can be decomposed into $\mathcal{O}(md^{2})$ Type-1 and
Type-2 Hamiltonians allows us to apply Trotterization Methods to implement
Algorithm 1 because it will only contribute polylogarithmic overhead to
complexities.

We apply the Trotterization method to the evolution of Hamiltonian $H(s)$ and
approximate it by the product of time evolution of all the summands, which are
Type 1 and Type 2 Hamiltonians. In the next section, we introduce our circuit
design for the two types of circuits and discuss their properties. In Section
4, we analyze the circuit approximation accuracy needed for the Troterization
Method and the total circuit depth of the whole implementation.

###  3.1 Circuit for Type-1 Hamiltonian

A Type-1 Hamiltonian is a tensor product of Pauli and Hermitian matrices.
$H_{1}(s)$ is a special case of Type-1 Hamiltonian since it is a scaled Pauli
matrix. For Pauli matrices, their time evolution is discussed in [27] Section
4.7.3. The idea is to apply single qubit operations on each qubit to turn any
$X$ or $Y$ to $Z$; then apply the matrix exponential for the new Hamiltonian
with only $I$ and $Z$; finally undo the operations to turn certain $Z$ back to
$X$ or $Y$. We refer the readers to [27] for details about this case.

For general Type-1 Hamiltonians, the following lemma is helpful for our
circuit design.

######  Lemma 3.2.

For any $t\geq 0$ and any Hermitian matrix $M$ with decomposition
$M=U_{P}\Sigma_{M}U_{P}^{\dagger}$, where $U_{P}$ is a unitary matrix, the
following identity holds

$e^{-itM}=U_{P}e^{-it\Sigma_{M}}U_{P}^{\dagger}.$

###### Proof.

Using Taylor expansion, it is obvious that

$\displaystyle e^{-itM}$
$\displaystyle=I+(-it)M+\frac{(-it)^{2}}{2!}M^{2}+\cdots$
$\displaystyle=U_{P}\left(I+(-it)\Sigma_{M}+\frac{(-it)^{2}}{2!}\Sigma_{M}^{2}+\cdots\right)U_{P}^{\dagger}$
$\displaystyle=U_{P}e^{-it\Sigma_{M}}U_{P}^{\dagger}.$

∎

######  Remark 3.2.

For Type-1 Hamiltonian $H^{1}$, it can be decomposed as

$\displaystyle H^{1}$
$\displaystyle=U_{H^{1}}\Sigma_{H^{1}}U_{H^{1}}^{\dagger}$

with

$\displaystyle U_{H^{1}}$ $\displaystyle=U_{P_{1}}\otimes U_{P_{2}}\otimes
U_{C_{1}}\otimes\cdots\otimes U_{C_{n}}$ $\displaystyle\Sigma_{H^{1}}$
$\displaystyle=\Sigma_{P_{1}}\otimes\Sigma_{P_{2}}\otimes\Sigma_{C_{1}}\otimes\cdots\otimes\Sigma_{C_{n}},$

where $U_{P_{i}}\Sigma_{P_{i}}U_{P_{i}}^{\dagger}$ is the eigenvalue
decomposition of $P_{i}$ and $U_{C_{j}}\Sigma_{C_{j}}U_{C_{j}}^{\dagger}$ is
the singular value decomposition of $C_{j}$. We choose the order of the
singular values such that $\Sigma_{H^{1}}$ can be represented as

$\displaystyle\Sigma_{H^{1}}$
$\displaystyle=\left\|\Sigma_{H^{1}}\right\|_{2}\Sigma_{P_{1}}\otimes\Sigma_{P_{2}}\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{1},1}\end{bmatrix}\otimes\cdots\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{1},n}\end{bmatrix},$ (4)

where $\sigma_{H^{1},i}\in[0,1]$ for all $i\in\\{1,\dots,n\\}.$ Then the time
evolution can be written as

$e^{-itH^{1}}=U_{H^{1}}e^{-it\Sigma_{H^{1}}}U_{H^{1}}^{\dagger}.$

For $U_{H^{1}}$, each matrix in the tensor product is a unitary matrix in
$\mathbb{C}^{2\times 2}$, thus we can implement the circuit for them
efficiently. For $e^{-it\Sigma_{H^{1}}}$, we introduce an algorithm to
implement the circuit. For the simplicity and the ease to read, we start with
a generic version, see Algorithm 2, and give explanations afterwards.

Algorithm 2 Generic Implementation of $e^{-it\Sigma}$

1:Given $t\geq 0$ and $\Sigma$ in the format of Eq. (4)

2:Represent diagonal element $\sigma_{i}$ in binary format,

3:Encode the binary representation of $\sigma_{i}$ as quantum state
$\ket{\sigma_{i}}$,

4:Apply quantum multiplier to $\ket{\sigma_{i}}$ controlled by the input
quantum states,

5:Apply a series of Phase Gate controlled by the result of quantum multiplier.

Step 2-4 in Algorithm 2 are inspired by [28]. We encourage the readers to read
[28] to learn their quantum adder and quantum multiplier. Their quantum adder
and quantum multiplier are designed for integer numbers. In the followings, we
give a thorough explanation of each step. For Step 4, we start with a brief
introduction of the quantum adder and multiplier. Then, we show that the
quantum adder and multiplier also work for rational numbers in the specific
form described in Lemma 3.3 and explain how they can be used in the
implementation.

####  3.1.1 Algorithm 2 Step 2 & 3

By the definition of $\sigma_{i}$ in Eq. (4), we have $\sigma_{i}\in[0,1]$.
Here we temporarily assume $\sigma_{i}$ can be represented using $p$ binary
digits

$\displaystyle\sigma_{i}=\sum_{j=1}^{p}\nu_{j}^{i}2^{-j},$

where $\nu_{j}^{i}$ are binary numbers. This assumption is relaxed in Section
4.1, where inexact binary representation is discussed. Then, we use
$\ket{\nu_{j}^{i}}=\ket{0}$ to represent $\nu_{j}^{i}=0$ and
$\ket{\nu_{j}^{i}}=\ket{1}$ to represent $\nu_{j}^{i}=1$, which gives us a
quantum representation of $\sigma_{i}$ as

$\ket{\sigma_{i}}=\ket{\nu_{j}^{1}}\ket{\nu_{j}^{2}}\cdots\ket{\nu_{j}^{p}}.$

Such binary presentation of numbers is used to construct quantum adder and
multiplier in [28].

####  3.1.2 Algorithm 2 Step 4 & 5

In [28], the authors introduce a non-modular quantum Fourier transform (QFT)
based quantum adder following the QFT adder introduced in [16]. Then they
introduce a quantum multiplier based on their quantum adder. The quantum adder
and the quantum multiplier are designed for positive integers originally and
easily accommodate signed integers by adding an extra qubit for the sign.
Staring from Lemma 3.3, we show that the quantum adder and the quantum
multiplier also work for rational numbers between $-1$ and $1$ with mild
modification of the original version. With the QFT-based quantum adder and
quantum multiplier, we are ready to explain Step 4 of Algorithm 2.

Notice that $\Sigma$ is a diagonal matrix, it is known that $e^{-it\Sigma}$ is
also a diagonal matrix with the $j$th diagonal entry being
$e^{-it\Sigma_{j}}$, where $\Sigma_{j}$ is the $j$th entry of $\Sigma$. When
we apply the circuit $e^{-it\Sigma}$ on quantum state $\ket{j}$, we expect
$e^{-it\Sigma_{j}}\ket{j}$ being the resulting state, i.e.,
$e^{-it\Sigma}\ket{j}=e^{-it\Sigma_{j}}\ket{j}$. Once we have a binary
representation of $\Sigma_{j}$ encoded as a quantum state, we apply a gate
phase circuit to turn $\ket{j}$ into $e^{-it\Sigma_{j}}\ket{j}$. Specifically,
for $\ket{\Sigma_{j}}=\ket{\nu_{j}^{1}}\cdots\ket{\nu_{j}^{p}}$, if
$\nu_{j}^{k}=1$, then we apply the phase shift $e^{-it2^{-k}}$ to $\ket{j}$.
This is the controlled Phase gate in Step 5.

What is left is the computation of the binary representation of $\Sigma_{j}$.
Notice that $\Sigma$ is a tensor product of $(n+2)$ $2$-by-$2$ diagonal
matrices with all the diagonal entries between $-1$ and $1$, we can get
$\Sigma_{j}$ by multiplying $n+2$ numbers. Specifically, for
$\ket{j}=\ket{j_{1}}\ket{j_{2}}\cdots\ket{j_{n+2}}$ with $j_{k}$ being binary
numbers, if $j_{k}=0$, then we take the first diagonal entry of the $k$th
$2$-by-$2$ diagonal matrix into the multiplication; otherwise, we take the
second diagonal entry of the matrix. This is the controlled quantum multiplier
in Step 4. The quantum circuit we describe maps $\ket{j}$ to
$e^{-it\Sigma_{j}}\ket{j}$ and thus implements the circuit $e^{-it\Sigma}$
because of the linearity of quantum circuit.

In the rest of this section, we prove our aforementioned claim that the
quantum adder and quantum multiplier can be generalized to certain non-integer
numbers in Lemma 3.3 to 3.6. Finally, we summarize our circuit for Type-1
Hamiltonian afterwards in Lemma 3.8.

######  Lemma 3.3.

Let $K$ be an integer. Let
$\nu^{1}=\nu^{1}_{1}2^{K-1}+\cdots+\nu^{1}_{p}2^{K-p}$ and
$\nu^{2}=\nu^{2}_{1}2^{K-1}+\cdots+\nu^{2}_{p}2^{K-p}$, where $\nu^{1}_{i}$
and $\nu^{2}_{j}$ are binary numbers for all $i\in[p]$ and $j\in[p]$. Let
$\ket{\nu^{1}}=\ket{\nu^{1}_{1}}\otimes\cdots\otimes\ket{\nu^{1}_{p}}$ and
$\ket{\nu^{2}}=\ket{\nu^{2}_{1}}\otimes\cdots\otimes\ket{\nu^{2}_{p}}$. If $C$
is a circuit such that, when $K=p$,

$C\ket{0}\otimes\ket{\nu^{1}}\otimes\ket{\nu^{2}}=\ket{\nu^{3}}\otimes\ket{\nu^{2}},$
(5)

where
$\ket{\nu^{3}}=\ket{\nu^{3}_{1}}\otimes\cdots\otimes\ket{\nu^{3}_{p+1}}$,
$\nu^{3}_{i}\in\\{0,1\\}$ for all $i\in[p+1]$, and

$\nu^{1}+\nu^{2}=\nu^{3}_{1}2^{K}+\cdots+\nu^{3}_{p+1}2^{K-p},$ (6)

then, Eq. (5) holds for any integer $K$.

###### Proof.

Proof. When $K=p$, the circuit $C$ is the quantum adder introduced in [28],
which applies to integer numbers. Since Eq. (6) holds for $K=p$, we have

$2^{K}\left(\sum_{i=1}^{p}\nu^{1}_{i}2^{-i}+\sum_{i=1}^{p}\nu^{2}_{i}2^{-i}-\sum_{i=1}^{p+1}\nu^{3}_{i}2^{1-i}\right)=0,$

which holds for any integer $K$. So the circuit $C$ is an adder for any number
pair $\nu^{1}$ and $\nu^{2}$. ∎

We have shown the quantum adder works for positive rational numbers in the
specific form described in Lemma 3.3. In fact, it also works for corresponding
signed rational numbers after mild modification following the argument in
[28]. This modification only adds one extra qubit for the sign and does not
change the asymptotic number of the gates needed to implement the adder. In
this work, we only care about the addition of numbers between $-1$ and $1$,
which is the case when $K=0$. For the adder to accomplish this job, we
summarize its properties in the following lemma. Readers are referred to
Section 4 in [28] for the proof.

######  Lemma 3.4.

Let $p$ be an integer. Let $\nu^{1}$ and $\nu^{2}$ be signed numbers with
$|\nu^{1}|=\nu^{1}_{1}2^{-1}+\cdots+\nu^{1}_{p}2^{-p}$ and
$|\nu^{2}|=\nu^{2}_{1}2^{-1}+\cdots+\nu^{2}_{p}2^{-p}$, where $\nu^{1}_{i}$
and $\nu^{2}_{j}$ are binary numbers for all $i\in[p]$ and $j\in[p]$. The
quantum adder introduced in [28] computes $\nu^{1}+\nu^{2}$ using
$\mathcal{O}(p^{2})$ single qubit gates and their controlled versions with
gate depth $\mathcal{O}(p^{2})$. Specifically, the adder uses one
$(p+2)$-qubit QFT, one $(p+2)$-qubit IQFT, and $(p+1)(p+2)/2$ controlled
single-qubit gates.

Similar argument can be made for the quantum multiplier. We summarize them in
the following two lemmas.

######  Lemma 3.5.

Let $K$ be an integer. Let
$\nu^{1}=\nu^{1}_{1}2^{K-1}+\cdots+\nu^{1}_{p}2^{K-p}$ and
$\nu^{2}=\nu^{2}_{1}2^{K-1}+\cdots+\nu^{2}_{p}2^{K-p}$, where $\nu^{1}_{i}$
and $\nu^{2}_{j}$ are binary numbers for all $i\in[p]$ and $j\in[p]$. Let
$\ket{\nu^{1}}=\ket{\nu^{1}_{1}}\otimes\cdots\otimes\ket{\nu^{1}_{p}}$ and
$\ket{\nu^{2}}=\ket{\nu^{2}_{1}}\otimes\cdots\otimes\ket{\nu^{2}_{p}}$. If $C$
is a circuit such that, when $K=p$,

$C\ket{0}_{2p}\otimes\ket{\nu^{1}}\otimes\ket{\nu^{2}}=\ket{\nu^{3}}\otimes\ket{\nu^{1}}\otimes\ket{\nu^{2}},$
(7)

where $\ket{\nu^{3}}=\ket{\nu^{3}_{1}}\otimes\cdots\otimes\ket{\nu^{3}_{2p}}$,
$\nu^{3}_{i}\in\\{0,1\\}$ for all $i\in[2p]$, and

$\nu^{1}\nu^{2}=\nu^{3}_{1}2^{2K-1}+\cdots+\nu^{3}_{2p}2^{2K-2p},$ (8)

then, Eq. (7) holds for any integer $K$.

###### Proof.

Proof. When $K=p$, the circuit $C$ is the original quantum multiplier
introduced in [28], which applies to integer numbers. Since Eq. (8) holds for
$K=p$, we have

$2^{2K}\left(\sum_{i=1}^{p}\nu^{1}_{i}2^{-i}\sum_{i=1}^{p}\nu_{i}^{2}2^{-i}-\sum_{i=1}^{2p}\nu_{i}^{3}2^{-i}\right)=0,$

which holds for any integer $K$. So the circuit $C$ is an multiplier for any
number pair $\nu^{1}$ and $\nu^{2}$. ∎

In [28], the authors also discuss how to multiply signed integer numbers by
adding qubits for signs. We refer the readers to Section 7 of [28] and
summarize the main result in the following lemma.

######  Lemma 3.6.

Let $\nu^{1}$ and $\nu^{2}$ be two signed numbers with
$|\nu^{1}|=\nu^{1}_{1}2^{-1}+\cdots+\nu^{1}_{p}2^{-p}$ and
$|\nu^{2}|=\nu^{2}_{1}2^{-1}+\cdots+\nu^{2}_{p}2^{-p}$, where $\nu^{1}_{i}$
and $\nu^{2}_{j}$ are binary numbers for all $i\in[p]$ and $j\in[p]$. The
quantum multiplier introduced in [28] computes $\nu^{1}\nu^{2}$ using
$\mathcal{O}(p^{3})$ single-qubit gates and their controlled version with gate
depth $\mathcal{O}(p^{3})$. Specifically, the multiplier uses one
$(2p+2)$-qubit QFT, one $(2p+2)$-qubit IQFT, and $p$ controlled adders.

For both the quantum adder and quantum multiplier, when their input numbers
are represented by $p$-qubit quantum states, we call the adder and multiplier
by $p$-qubit adder and $p$-qubit multiplier, respectively. The following lemma
discusses the complexity of Algorithm 2.

######  Lemma 3.7.

If $\sigma_{H^{1},i}=\sum_{j=1}^{p}\nu_{j}^{i}2^{-j}$ with $\nu_{j}^{i}$’s
being binary for all $i\in[n]$, there exists a quantum circuit that prepares
$e^{-it\Sigma_{H^{1}}}$ as described in Algorithm 2, for any $t\geq 0$, using
at most $\mathcal{O}(n)$ calls to the $p$-qubit quantum multiplier with gate
depth $\mathcal{O}(np^{3})$.

###### Proof.

Proof. Step 4 in Algorithm 2 computes the multiplication of at most $n+2$
$p$-qubit numbers and thus needs at most $n+2$ calls to the $p$-qubit
multiplier. ∎

The following lemma summarizes the cost to perform time evolution of Type-1
Hamiltonian $H^{1}$.

######  Lemma 3.8.

For Type-1 Hamiltonian $H^{1}$ with decomposition

$\displaystyle H^{1}$
$\displaystyle=U_{H^{1}}\Sigma_{H^{1}}U_{H^{1}}^{\dagger}$

with

$\displaystyle U_{H^{1}}$ $\displaystyle=U_{P_{1}}\otimes U_{P_{2}}\otimes
U_{C_{1}}\otimes\cdots\otimes U_{C_{n}}$ $\displaystyle\Sigma_{H^{1}}$
$\displaystyle=\Sigma_{P_{1}}\otimes\Sigma_{P_{2}}\otimes\Sigma_{C_{1}}\otimes\cdots\otimes\Sigma_{C_{n}},$

where $U_{P_{i}}\Sigma_{P_{i}}U_{P_{i}}^{\dagger}$ is the eigenvalue
decomposition of $P_{i}$ and $U_{C_{j}}\Sigma_{C_{j}}U_{C_{j}}^{\dagger}$ is
the singular value decomposition of $C_{j}$. We can rewrite $\Sigma_{H^{1}}$
as

$\displaystyle\Sigma_{H^{1}}$
$\displaystyle=\left\|H^{1}\right\|_{2}\Sigma_{P_{1}}\otimes\Sigma_{P_{2}}\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{1},1}\end{bmatrix}\otimes\cdots\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{1},n}\end{bmatrix},$

where $\sigma_{H^{1},i}\in[0,1]$ for all $i\in\\{1,\dots,n\\}.$ If
$\sigma_{H^{1},i}=\sum_{j=1}^{p}\nu_{j}^{i}2^{-j}$ with $\nu_{j}^{i}$’s being
binary for all $i\in[n]$, then there exists a quantum circuit that prepares
$e^{-itH^{1}}$ for any $t\geq 0$, using $\mathcal{O}(n)$ classical arithmetic
operations, $\mathcal{O}(n)$ single qubit unitary circuits with gate depth
$\mathcal{O}(1)$, and $\mathcal{O}(n)$ calls to the $p$-qubit quantum
multiplier with gate depth $\mathcal{O}(np^{3})$.

###### Proof.

Proof. We start by computing the singular value decomposition of $H^{1}$.
Considering $H^{1}$ is a tensor product of $n+2$ two-by-two matrices, the
singular value decomposition of $H^{1}$ is the tensor product of the singular
value decomposition of the $2$-by-$2$ matrices. For each $2$-by-$2$ matrix,
the number of classical arithmetic operations needed for computing singular
value decomposition is a constant. So in total $\mathcal{O}(n)$ classical
arithmetic operations are needed.

Then, we need to construct the circuit for unitary operation $U_{H^{1}}$,
which is the tensor product of $n+2$ two-by-two matrices. So we need
$\mathcal{O}(n)$ single qubit unitary circuits.

Finally, we need to construct the circuit for $e^{-it\Sigma_{H^{1}}}$.
According to Lemma 3.7, we need $\mathcal{O}(n)$ calls to the $p$-qubit
quantum multiplier. ∎

###  3.2 Circuit for Type-2 Hamiltonian

Unlike Type-1 Hamiltonians, in general, Type-2 Hamiltonians are not in the
tensor format, which brings extra challenges for implementation. In this
section, we introduce a method to avoid the issue and show that circuit for
Type-2 Hamiltonian can be implemented in a similar way as Type-1 Hamiltonian.

First, we classically compute the singular value decomposition of
$D_{1}\otimes\cdots\otimes D_{n}$. This step takes at most
$\mathcal{O}\left(n\right)$ classical arithmetic operations. Denote the
singular value decomposition of $P_{0}\otimes D_{1}\otimes\cdots\otimes D_{n}$
by

$\displaystyle P_{0}\otimes D_{1}\otimes\cdots\otimes D_{n}$
$\displaystyle=U_{H^{2}}\Sigma_{H^{2}}V_{H^{2}}^{\dagger},$

where

$\displaystyle U_{H^{2}}$ $\displaystyle=P_{0}\otimes
U_{D_{1}}\otimes\cdots\otimes U_{D_{n}}$ $\displaystyle V_{H^{2}}$
$\displaystyle=I\otimes V_{D_{1}}\otimes\cdots\otimes V_{D_{n}}$
$\displaystyle\Sigma_{H^{2}}$
$\displaystyle=\left\|H^{2}\right\|_{2}I\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{2},1}\end{bmatrix}\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{2},2}\end{bmatrix}\otimes\cdots\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{2},n}\end{bmatrix}$

with $\sigma_{H^{2},i}\in[0,1]$ for all $i\in[n]$. Then, we have

$H^{2}=\begin{bmatrix}U_{H^{2}}&0\\\
0&V_{H^{2}}\end{bmatrix}\times\begin{bmatrix}0&\Sigma_{H^{2}}\\\
\Sigma_{H^{2}}&0\end{bmatrix}\times\begin{bmatrix}U_{H^{2}}^{\dagger}&0\\\
0&V_{H^{2}}^{\dagger}\end{bmatrix}.$

Note that the decomposition is not a singular value decomposition but it is
useful for the implementation because the outer matrices are unitaries and the
inner matrix is in tensor format, where Lemma 3.2 applies. The lemma indicates
that we need to implement the circuit of the outer unitary matrix and the
circuit for the time evolution of the inner matrix.

For the outer matrix, it is unitary but it is not in the tensor format. We can
further decompose it as follows

$\begin{bmatrix}U_{H^{2}}&0\\\ 0&V_{H^{2}}\end{bmatrix}=L_{H^{2}}R_{H^{2}},$

where

$\displaystyle L_{H^{2}}$ $\displaystyle=\begin{bmatrix}U_{H^{2}}&0\\\
0&I_{2^{n+1}}\end{bmatrix},\ R_{H^{2}}=\begin{bmatrix}I_{2^{n+1}}&0\\\
0&V_{H^{2}}\end{bmatrix}.$

$L_{H^{2}}$ and $R_{H^{2}}$ are still unitary and represent the controlled
version of $U_{H^{2}}$ and $V_{H^{2}}$, respectively. Since $U_{H^{2}}$ and
$V_{H^{2}}$ are tensor product of $2$-by-$2$ unitaries, we can further
decompose $L_{H^{2}}$ and $R_{H^{2}}$ into

$\displaystyle L_{H^{2}}$ $\displaystyle=\begin{bmatrix}P_{0}\otimes
I\otimes\cdots\otimes I&0\\\ 0&I_{2^{n+1}}\end{bmatrix}\begin{bmatrix}I\otimes
U_{D_{1}}\otimes\cdots\otimes I&0\\\
0&I_{2^{n+1}}\end{bmatrix}\cdots\begin{bmatrix}I\otimes I\otimes\cdots\otimes
U_{D_{n}}&0\\\ 0&I_{2^{n+1}}\end{bmatrix}$ $\displaystyle R_{H^{2}}$
$\displaystyle=\begin{bmatrix}I_{2^{n+1}}&0\\\ 0&I\otimes
I\otimes\cdots\otimes I\end{bmatrix}\begin{bmatrix}I_{2^{n+1}}&0\\\ 0&I\otimes
V_{D_{1}}\otimes\cdots\otimes
I\end{bmatrix}\cdots\begin{bmatrix}I_{2^{n+1}}&0\\\ 0&I\otimes
I\otimes\cdots\otimes V_{D_{n}}\end{bmatrix}.$

Each matrix in the decomposition of $L_{H^{2}}$ and $R_{H^{2}}$ is a single
qubit unitary controlled by $n+1$ qubits and can be implemented efficiently.

According to Lemma 3.2, we still need to implement
$e^{-itX\otimes\Sigma_{H^{2}}}$. However, $X\otimes\Sigma_{H^{2}}$ is not a
diagonal matrix and thus cannot be implemented directly using Algorithm 2. To
void this, we consider the eigenvalue decomposition of $X$ and have

$\displaystyle X\otimes\Sigma_{H^{2}}$
$\displaystyle=(U_{X}\Sigma_{X}U_{X}^{\dagger})\otimes\Sigma_{H^{2}}$ (9)
$\displaystyle=(U_{X}\otimes
I_{2^{n+1}})(\Sigma_{X}\otimes\Sigma_{H^{2}})(U_{X}\otimes
I_{2^{n+1}})^{\dagger}.$

Now we apply Lemma 3.2 again since $U_{X}\otimes I_{2^{n+1}}$ is a unitary.
Also, since $\Sigma_{X}=\begin{bmatrix}1&0\\\ 0&-1\end{bmatrix}$, the matrix
$\Sigma_{X}\otimes\Sigma_{H^{2}}$ is a tensor product of $2$-by-$2$ diagonal
matrices, where Algorithm 2 applies. We summarize the analysis above in the
following lemma.

######  Lemma 3.9.

For Type-2 Hamiltonian $H^{2}$, it can be decomposed into

$H^{2}=\begin{bmatrix}U_{H^{2}}&0\\\
0&V_{H^{2}}\end{bmatrix}\times\begin{bmatrix}0&\Sigma_{H^{2}}\\\
\Sigma_{H^{2}}&0\end{bmatrix}\times\begin{bmatrix}U_{H^{2}}^{\dagger}&0\\\
0&V_{H^{2}}^{\dagger}\end{bmatrix}$

with

$\displaystyle U_{H^{2}}$ $\displaystyle=P_{0}\otimes
U_{D_{1}}\otimes\cdots\otimes U_{D_{n}}$ $\displaystyle V_{H^{2}}$
$\displaystyle=I\otimes V_{D_{1}}\otimes\cdots\otimes V_{D_{n}}$
$\displaystyle\Sigma_{H^{2}}$
$\displaystyle=\left\|H^{2}\right\|_{2}I\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{2},1}\end{bmatrix}\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{2},2}\end{bmatrix}\otimes\cdots\otimes\begin{bmatrix}1&0\\\
0&\sigma_{H^{2},n}\end{bmatrix},$ (10)

where $U_{H^{2}}\Sigma_{H^{2}}V_{H^{2}}^{\dagger}$ is the singular value
decomposition of $P_{0}\otimes D_{1}\otimes\cdots\otimes D_{n}$ and
$\sigma_{H^{2},i}\in[0,1]$ for all $i\in[n]$. If
$\sigma_{H^{2},i}=\sum_{j=1}^{p}\nu_{j}^{i}2^{-j}$ with $\nu_{j}^{i}$’s being
binary for all $i\in[n]$, then there exists a quantum circuit that prepares
$e^{-itH^{2}}$ for any $t\geq 0$, using $\mathcal{O}(n)$ classical arithmetic
operations, $\mathcal{O}(n)$ single qubit unitary circuits, and
$\mathcal{O}(n)$ calls to $p$-qubit quantum multiplier.

###### Proof.

Proof. According to the analysis in Section 3, the decomposition exists.
According to Lemma 3.2, we need to implement unitary circuit
$\begin{bmatrix}U_{H^{2}}&0\\\ 0&V_{H^{2}}\end{bmatrix}$ and unitary circuit
$exp\left\\{-it\begin{bmatrix}0&\Sigma_{H^{2}}\\\
\Sigma_{H^{2}}&0\end{bmatrix}\right\\}$. The first one can be implemented
using $\mathcal{O}(n)$ controlled single qubit unitary circuits. The second
one can be decomposed using Eq. (9). After decomposition, one can implement
$U_{X}\otimes I_{2^{n+1}}$ using one controlled single qubit unitary circuit.
As for the time evolution circuit
$exp\left\\{-it\Sigma_{X}\otimes\Sigma_{H^{2}}\right\\}$, the Hamiltonian
$\Sigma_{X}\otimes\Sigma_{H^{2}}$ is in the same form of $\Sigma_{H^{1}}$, so
we can use Algorithm 2 to implement the time evolution circuit, which comes
with at most $\mathcal{O}(n)$ calls to $p$-qubit quantum multiplier according
to Lemma 3.7. ∎

In this section, we showed how to implement the time evolution for Type-1 and
Type-2 Hamiltonian and give cost estimation for both of them. However, we
assume that all the numbers $\sigma_{H^{j},i}$ can be represented using $p$
qubits, which does not hold for the general case. We discuss this issue in
Section 4.

##  4 Circuit Cost Analysis

In this section, we analyze the total cost needed to use the circuits proposed
in the previous section to implement Algorithm 1. Our analysis are twofold: we
first estimate the cost needed to implement the circuits when we do not have
exact $p$-qubit binary representation of data, and then we estimate the
Trotterization steps needed for Algorithm 1 to converge.

###  4.1 Cost Analysis for Approximate Data

For the analysis in Section 3, we assume $\sigma_{H^{1},i}$ in Eq. (4) and
$\sigma_{H^{2},j}$ in Eq. (10) can be represented by $p$-qubit quantum states.
This assumption does not hold for the general case. In fact, we do not need an
exact binary representation of these numbers. Here we introduce the definition
of $\epsilon$-approximate binary representation of a number $\nu\in[0,1]$ and
then discuss the behavior of quantum adder and multiplier on approximate
binary representations. Finally, we discuss the relationship between the
inexactness of data and the inexactness of the circuit constructed by
Algorithm 2.

######  Definition 4.1.

Let $\nu\in[0,1]$, $\epsilon>0$, and $p\in\mathbb{Z}_{+}$. Let
$\nu_{i}\in\\{0,1\\}$ for all $i\in\\{1,2,\dots,p\\}$. We call

$\tilde{\nu}=\nu_{1}2^{-1}+\nu_{2}2^{-2}+\dots+\nu_{2}^{-p}$

a $p^{\rm th}$-order $\epsilon$-approximation of $\nu$, if

$|\nu-\tilde{\nu}|\leq\epsilon.$

The corresponding binary representation is called a $p^{\rm th}$-order
$\epsilon$-approximate binary representation of $\nu$.

In this work, we use $p$ qubits to store a $p^{\rm th}$-order binary
representation, so we also call the approximation a $p$-qubit approximation
and the binary representation a $p$-qubit binary representation.

######  Lemma 4.1.

Let $\nu\in[0,1]$, $\epsilon>0$, and $p_{0}=\lceil
1-\log_{2}(\epsilon)\rceil$. There exists $\nu_{i}\in\\{0,1\\}$ for all
$i\in\\{1,2,\dots,p_{0}\\}$ such that

$\tilde{\nu}=\nu_{1}2^{-1}+\nu_{2}2^{-2}+\dots+\nu_{p_{0}}2^{-p_{0}}$

is a $p_{0}$-qubit $\epsilon$-approximation $\nu$.

###### Proof.

Proof. Consider the exact binary representation of $\nu$:

$\displaystyle\nu=\sum_{i=1}^{+\infty}\nu_{i}2^{-i}=\tilde{\nu}+\sum_{i=p_{0}+1}^{+\infty}\nu_{i}2^{-i}.$

It is obvious that

$|\nu-\tilde{\nu}|=\sum_{i=K+1}^{+\infty}\nu_{i}2^{-i}\leq\sum_{i=p_{0}+1}^{+\infty}2^{-i}\leq
2^{-p_{0}}\leq\epsilon.$

∎

Lemma 4.1 indicates that we can obtain an $\epsilon$-approximate binary
representation of a number using at most $\lceil 1-\log_{2}(\epsilon)\rceil$
qubits. In the next lemma, we discuss the quantum multiplication of two
$p$-qubit approximate numbers.

######  Lemma 4.2.

Let $0<\epsilon<\epsilon^{j}<1$ for $j\in[2]$, $p=\lceil
1-\log_{2}(\epsilon)\rceil$, and $0\leq\nu^{1},\nu^{2}\leq 1$. Let
$\tilde{\nu}^{j}=\sum_{i=1}^{p}\nu_{i}^{j}2^{-i}$ be the $p$-qubit
$\epsilon^{j}$-approximation of $\nu^{j}$ for $j\in[2]$. Let
$\nu^{3}=\sum_{i=1}^{2p}\nu^{3}_{i}2^{-i}$ be $\tilde{\nu}^{1}\tilde{\nu}^{2}$
computed by the quantum multiplier introduced in [28]. Then, $\nu^{3}$ is a
$2p$-qubit
$(\epsilon^{1}+\epsilon^{2}+\epsilon^{1}\epsilon^{2})$-approximation of
$\nu^{1}\nu^{2}$. Let $\tilde{\nu}^{3}=\sum_{i=1}^{p}\nu^{3}_{i}2^{-i}$ be an
approximation of $\nu^{3}$. Then, $\tilde{\nu}^{3}$ is a $p$-qubit
$(\epsilon+\epsilon^{1}+\epsilon^{2}+\epsilon^{1}\epsilon^{2})$-approximation
of $\nu^{1}\nu^{2}$.

###### Proof.

Proof. According to Lemma 4.1, both $\tilde{\nu}^{1}$ and $\tilde{\nu}^{2}$
exist. Since $\tilde{\nu}^{1}$ is a $p$-qubit $\epsilon$-approximation of
$\nu^{1}$, we have $\tilde{\nu}^{1}=\nu^{1}+\delta^{1}$ with
$|\delta^{1}|\leq\epsilon^{1}$. Similarly for $\nu^{2}$, we have
$\tilde{\nu}^{2}=\nu^{2}+\delta^{2}$ with $|\delta^{2}|\leq\epsilon^{2}$. The
quantum multiplier gives an exact $2p$-qubit binary representation of
$\nu^{3}=\tilde{\nu}^{1}\tilde{\nu}^{2}$. By their definitions, we have

$\displaystyle\left|\nu^{3}-\nu^{1}\nu^{2}\right|$
$\displaystyle=\left|(\nu^{1}+\delta^{1})(\nu^{2}+\delta^{2})-\nu^{1}\nu^{2}\right|$
$\displaystyle=\left|\nu^{1}\delta^{2}+\nu^{2}\delta^{1}+\delta^{1}\delta^{2}\right|$
$\displaystyle\leq\epsilon^{1}+\epsilon^{2}+\epsilon^{1}\epsilon^{2}.$

Since $\nu^{3}$ is represented by $2p$ qubits, one can discard the last $p$
less-important qubits, leaving a $p$-qubit approximation of $\nu^{3}$, i.e.,
$\tilde{\nu}^{3}$. As a result, we have

$\displaystyle\left|\tilde{\nu}^{3}-\nu^{1}\nu^{2}\right|$
$\displaystyle=\left|\tilde{\nu}^{3}-\nu^{3}+\nu^{3}-\nu^{1}\nu^{2}\right|$
$\displaystyle\leq\epsilon+\left|\nu^{3}-\nu^{1}\nu^{2}\right|,$

where the last inequality holds because

$\left|\nu^{3}-\tilde{\nu}^{3}\right|=\sum_{i=p+1}^{2p}\nu^{3}_{i}2^{-i}\leq\sum_{i=p+1}^{2p}2^{-i}\leq
2^{-p}\leq\epsilon.$

∎

As indicated in the previous lemma, when multiplying two $p$-qubit numbers,
one can either keep the $2p$-qubit result or approximate it by a $p$-qubit
approximation. However, when considering the multiplication of multiple
$p$-qubit numbers, since we stick with the quantum multiplier introduced in
[28], we have to consider the following issue. For example, if we add an extra
$p$-qubit number $\nu^{4}$ and we want the result of $\nu^{1}\nu^{2}\nu^{4}$.
Then, after we get $\nu^{3}$ from the quantum multiplier, we have to decide
whether to use $\tilde{\nu}^{3}\nu^{4}$ as an approximate result, or to turn
$\nu^{4}$ into a $2p$ qubit number and obtain a $4p$-qubit number as the
result. If we apply the second strategy, then number of qubits needed would
grow exponentially as the amount of numbers grow. Instead, if we apply the
first strategy, for every extra number, we need to use $p$ more qubits, which
grows linearly. We call the first strategy by truncation strategy. In this
work, we stick with the quantum multiplier introduced in [28] and choose the
truncation strategy described above.

In the next lemma, we discuss the quantum multiplication of multiple $p$-qubit
approximate numbers.

######  Lemma 4.3.

Let $K$ be a positive integer larger than 1. Let $0<\epsilon<1$ such that
$K\epsilon<1/3$. Let $p=\lceil 1-\log_{2}(\epsilon)\rceil$. Let
$0\leq\nu^{j}\leq 1$ for $j\in[K]$ and
$\tilde{\nu}^{j}=\sum_{i=1}^{p}\nu_{i}^{j}2^{-i}$ be the $p$-qubit
$\epsilon$-approximation of $\nu^{j}$ for $j\in[K]$. Let $\nu^{\ast}$ be the
result computed by the quantum multiplier introduced in [28] following the
truncation strategy. Then, $\nu^{\ast}$ is a $p$-qubit
$\epsilon_{K}$-approximation of $\prod_{j=1}^{K}\nu^{j}$, where

$\displaystyle\epsilon_{K}\leq 3K\epsilon.$

###### Proof.

Proof. We prove the lemma using induction. Define function
$f(k)=2\sum_{i=1}^{k}(k\epsilon)^{i}$, for $k=2,\dots.$ It is obvious that
$f(k)\leq 2k\epsilon/(1-k\epsilon)$ when $0<k\epsilon<1$. When $k=2$,
according to Lemma 4.2, we have

$\epsilon_{2}=3\epsilon+\epsilon^{2}\leq 2(2\epsilon+4\epsilon^{2})=f(2).$

Let us assume that when $k=K$, $\nu^{+}$ is a $p$-qubit
$\epsilon_{K}$-approximation of $\prod_{j=1}^{K}\nu^{j}$ with
$\epsilon_{K}\leq f(K)$. Then, when $k=K+1$, according to Lemma 4.2, it
follows that

$\displaystyle\epsilon_{K+1}$
$\displaystyle=\epsilon+\epsilon_{K}+\epsilon+\epsilon_{K}\epsilon$
$\displaystyle\leq 2\epsilon+f(K)+f(K)\epsilon$
$\displaystyle=2\epsilon+\left(2K\epsilon+2\sum_{i=2}^{K}(K\epsilon)^{i}\right)+\left(2\sum_{i=1}^{K-1}(K\epsilon)^{i}\epsilon+2K^{K}\epsilon^{K+1}\right)$
$\displaystyle=2(K+1)\epsilon+2\sum_{i=2}^{K}(K^{i}+K^{i-1})\epsilon^{i}+2K^{K}\epsilon^{K+1}$
$\displaystyle\leq
2(K+1)\epsilon+2\sum_{i=2}^{K}(K+1)^{i}\epsilon^{i}+2(K+1)^{K+1}\epsilon^{K+1}$
$\displaystyle=f(K+1).$

So we conclude $\epsilon_{K}\leq f(K)\leq 2K\epsilon/(1-K\epsilon)\leq
3K\epsilon$. ∎

######  Lemma 4.4.

Let $0<\epsilon<\frac{1}{3n}$ and $p=\lceil 1-\log_{2}(\epsilon)\rceil$. Let
$\tilde{\sigma}_{H^{1},i}$ be a $p$-qubit $\epsilon$-approximation of
$\sigma_{H^{1},i}$ for all $i\in[n]$. When we use $\tilde{\sigma}_{H^{1},i}$
in Algorithm 1 and obtain an inexact circuit $e^{-it\tilde{\Sigma}_{H^{1}}}$,
the following bound holds

$\left\|e^{-it\tilde{\Sigma}_{H^{1}}}-e^{-it\Sigma_{H^{1}}}\right\|_{2}\leq
t\epsilon.$

###### Proof.

Proof. Both $e^{-it\tilde{\Sigma}_{H^{1}}}$ and $e^{-it\Sigma_{H^{1}}}$ are
diagonal matrices, it follows that

$\displaystyle\left\|e^{-it\tilde{\Sigma}_{H^{1}}}-e^{-it\Sigma_{H^{1}}}\right\|_{2}$
$\displaystyle\leq\max_{i\in[n]}\left|e^{-it\tilde{\sigma}_{H^{1},i}}-e^{-it\sigma_{H^{1},i}}\right|$
$\displaystyle=\max_{i\in[n]}\left|(e^{-it(\tilde{\sigma}_{H^{1},i}-\sigma_{H^{1},i})}-1)e^{-it\sigma_{H^{1},i}}\right|$
$\displaystyle=\max_{i\in[n]}\left|e^{-it(\tilde{\sigma}_{H^{1},i}-\sigma_{H^{1},i})}-1\right|$
$\displaystyle=\max_{i\in[n]}\sqrt{\left(e^{-it(\tilde{\sigma}_{H^{1},i}-\sigma_{H^{1},i})}-1\right)\left(e^{-it(\tilde{\sigma}_{H^{1},i}-\sigma_{H^{1},i})}-1\right)^{\dagger}}$
$\displaystyle=\max_{i\in[n]}\sqrt{2-2\cos\left(t\tilde{\sigma}_{H^{1},i}-t\sigma_{H^{1},i}\right)}$
$\displaystyle\leq\max_{i\in[n]}t\left|\tilde{\sigma}_{H^{1},i}-\sigma_{H^{1},i}\right|$
$\displaystyle\leq t\epsilon.$

∎

###  4.2 Trotterization Cost Analysis

In Section 3, we introduce how to implement the time evolution circuits of two
types of Hamiltonian. When we use those circuits in the implementation of
Algorithm 1, we decompose the Hamiltonian (1) into some structured
Hamiltonians as shown in Eq.(3). Then we use the time evolution circuits of
these structured Hamiltonians to approximate the time evolution circuits used
in Algorithm 1, which is the Trotterization method. As mentioned in Section
2.5, Trotterization method introduces error into circuits and makes circuits
inexact. Lemma 2.1 provides a bound for the Trotterization error in terms of
communicative factor $\tilde{\alpha}_{\rm comm}$. In the next lemma, we
provide a bound for $\tilde{\alpha}_{\rm comm}$ when the Hamiltonian
decomposition in Eq. (3) is used by the specific Trotterization method – the
first-order Lie-Trotter formula.

######  Lemma 4.5.

$\tilde{\alpha}_{\rm comm}\leq 2(1+m+2d^{2}+2md^{2})^{2}=\mathcal{O}(1)$.

###### Proof.

Proof. As discussed in Remark 2.2 and Lemma 3.1, we need to find the spectral
norm of the $(m+1)$ Type-1 Hamiltonians and the $(2d^{2}+2md)$ Type-2
Hamiltonians. According to Eq. (3), $H(s)$ is decomposed into $4$
Hamiltonians, each Hamiltonian is further decomposed into Type-1 and Type-2
Hamiltonians. For the single Type-1 Hamiltonian in $H_{1}(s)$, we have

$\|X\otimes Z\otimes(\otimes_{l=1}^{n}I)\|_{2}=1.$

For the $m$ Type-1 Hamiltonians in $H_{2}(s)$, we have

$\|X\otimes X\otimes(\otimes_{l=1}^{n}A_{il})\|_{2}\leq 1$

according to Assumption 2.1.

For the $2d^{2}$ Type-2 Hamiltonians in $H_{3}(s)$, denote the singular value
decomposition of $Z\otimes b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}$ by
$U_{Zb}\Sigma_{Zb}V_{Zb}^{\dagger}$, it follows that

$\displaystyle\left\|\begin{bmatrix}0&Z\otimes
b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}\\\ Z\otimes
b_{j_{2}\cdot}b_{j_{1}\cdot}^{\dagger}&0\end{bmatrix}\right\|_{2}$
$\displaystyle=\left\|\begin{bmatrix}U_{Zb}&0\\\
0&V_{Zb}\end{bmatrix}\begin{bmatrix}0&\Sigma_{Zb}\\\
\Sigma_{Zb}^{\dagger}&0\end{bmatrix}\begin{bmatrix}U_{Zb}^{\dagger}&0\\\
0&V_{Zb}^{\dagger}\end{bmatrix}\right\|_{2}$
$\displaystyle=\left\|\begin{bmatrix}0&\Sigma_{Zb}\\\
\Sigma_{Zb}^{\dagger}&0\end{bmatrix}\right\|_{2}$
$\displaystyle=\left\|\Sigma_{Zb}\right\|_{2}$
$\displaystyle\leq\left\|b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}\right\|_{2}$
$\displaystyle=\left\|b_{j_{1}\cdot}\right\|_{2}\left\|b_{j_{2}\cdot}^{\dagger}\right\|_{2}$
$\displaystyle\leq 1,$

where the last inequality holds because of Assumption 2.1. Similar for the
remaining Type-2 Hamiltonians in $H_{3}(s)$.

For the $2md^{2}$ Type-2 Hamiltonians in $H_{4}(s)$, notice that

$\left\|A_{i\cdot}b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}\right\|_{2}\leq\left\|A_{i\cdot}\right\|_{2}\left\|b_{j_{1}\cdot}b_{j_{2}\cdot}^{\dagger}\right\|_{2}\leq
1,$

it follows that all the Type-2 Hamiltonians in $H_{4}(s)$ have spectral norm
bounded by $1$. Combine these bounds with Lemma 3.1, we have

$\tilde{\alpha}_{\rm comm}\leq 2(1+m+2d^{2}+2md^{2})^{2}.$

∎

###  4.3 Total Cost Analysis

Both inexact data and Trotterization method introduce errors into the unitary
circuits. In this section, we discuss the relationship between the inexactness
of circuits and the accuracy of final state and then estimate the total cost
needed to implement Algorithm 1.

######  Lemma 4.6.

Let $K$ be a positive integer and $\epsilon>0$. Let $U_{i}$ for $i\in[K]$ be
unitary circuits and $\tilde{U}_{i}$ for $i\in[K]$ be the corresponding
inexact circuits such that $\tilde{U}_{i}=U_{i}+\mathcal{E}_{i}$ for all
$i=[K]$. Then, the following inequality holds

$\left\|\prod_{i=1}^{K}\tilde{U}_{i}-\prod_{i=1}^{K}U_{i}\right\|_{2}\leq\sum_{i=1}^{K}\|\mathcal{E}_{i}\|_{2}.$

###### Proof.

Proof. Although the circuits $\tilde{U}_{i}$ are inexact, they are still
unitary operators. So the following identity holds

$\prod_{i=1}^{k}\left(U_{i}+\mathcal{E}_{i}\right)-\prod_{i=1}^{k}U_{i}=\left(\prod_{i=1}^{k-1}\left(U_{i}+\mathcal{E}_{i}\right)-\prod_{i=1}^{k}U_{i}(U_{k}+\mathcal{E}_{k})^{\dagger}\right)(U_{k}+\mathcal{E}_{k})$

for all $k\in[K]$ because $U_{k}+\mathcal{E}_{k}$ is unitary. From the above
identity, we have

$\displaystyle\left\|\prod_{i=1}^{K}\tilde{U}_{i}-\prod_{i=1}^{K}U_{i}\right\|_{2}$
$\displaystyle=\left\|\prod_{i=1}^{K}\left(U_{i}+\mathcal{E}_{i}\right)-\prod_{i=1}^{K}U_{i}\right\|_{2}$
$\displaystyle=\left\|\left(\prod_{i=1}^{K-1}\left(U_{i}+\mathcal{E}_{i}\right)-\prod_{i=1}^{K-1}U_{i}(U_{K}+\mathcal{E}_{K})^{\dagger}\right)(U_{K}+\mathcal{E}_{K})\right\|_{2}$
$\displaystyle=\left\|\prod_{i=1}^{K-1}\left(U_{i}+\mathcal{E}_{i}\right)-\prod_{i=1}^{K}U_{i}(U_{K}+\mathcal{E}_{K})^{\dagger}\right\|_{2}$
$\displaystyle=\left\|\prod_{i=1}^{K-1}\left(U_{i}+\mathcal{E}_{i}\right)-\prod_{i=1}^{K-1}U_{i}-\prod_{i=1}^{K}U_{i}\mathcal{E}_{K}^{\dagger}\right\|_{2}$
$\displaystyle\leq\left\|\prod_{i=1}^{K-1}\left(U_{i}+\mathcal{E}_{i}\right)-\prod_{i=1}^{K-1}U_{i}\right\|_{2}+\left\|\prod_{i=1}^{K}U_{i}\mathcal{E}_{K}^{\dagger}\right\|_{2}$
$\displaystyle\leq\left\|\prod_{i=1}^{K-1}\left(U_{i}+\mathcal{E}_{i}\right)-\prod_{i=1}^{K-1}U_{i}\right\|_{2}+\left\|\mathcal{E}_{K}\right\|_{2}$
$\displaystyle\quad\vdots$
$\displaystyle\leq\sum_{i=1}^{K}\|\mathcal{E}_{i}\|_{2},$

where the third equality and the second inequality hold because unitary
operators do not change spectral norm and the last inequality holds because of
recursion. ∎

The following lemma comes from [22] discussing the trace distance between
exact and inexact matrices when applied to the same density operator.

######  Lemma 4.7 (Lemma 9 in [22]).

Let $U_{i}$ for $i\in[2]$ be square matrices and $\tilde{U}_{i}$ for $i\in[2]$
be the corresponding approximate matrices such that
$\tilde{U}_{i}=U_{i}+\mathcal{E}_{i}$ with
$\|\mathcal{E}_{i}\|_{2}\leq\epsilon$ for $i\in[2]$. Then for every density
operator $\rho$, the following inequality holds

$\left\|\tilde{U}_{1}\rho\tilde{U}_{2}-U_{1}\rho
U_{2}\right\|_{1}\leq\epsilon(\|U_{1}\|_{2}+\|U_{2}||_{2})+\epsilon^{2}.$

Now we are ready to provide the main theorem of this work.

######  Theorem 4.1.

Let $0<\epsilon\leq 1/(3n)$,
$\epsilon_{0}=\epsilon^{2}/(\kappa\log^{2}\kappa)$ and $p=\lceil
1-\log_{2}\epsilon_{0}\rceil$. When using Trotter number
$r=\mathcal{O}(m^{2}d^{4}\kappa^{2}/\epsilon)$ for Trotterization method and
$p$-qubit $\epsilon_{0}$-approximation of $\sigma_{i}$ in Algorithm 2, our
implementation of Algorithm 1 prepares an $\mathcal{O}(\epsilon)$-approximate
solution using $\mathcal{T}$ classical arithmetic operations, $\mathcal{T}$
single qubit unitary circuits  and their controlled version with gate depth
$\mathcal{O}(1)$, and $\mathcal{T}$ calls to $p$-qubit quantum multiplier  and
their controlled version with gate depth $\mathcal{O}(\mathcal{T}p^{3})$,
where
$\mathcal{T}=\mathcal{O}({\kappa^{2}log(N)\log^{2}(\kappa)/\epsilon^{2}})$.

###### Proof.

Proof. In Algorithm 1, denote the initial state by $\ket{\psi_{0}}$ and the
ideal final state by $\ket{\psi^{\ast}}$. We first discuss the cost to
implement Algorithm 1 and finally demonstrate that the cost to prepare the
initial is negligible. Notice that evolution time $t=(t_{1},\dots,t_{q})$ is a
$q$-dimensional random variable, we denote the circuits by
$U_{j}(t_{j})=e^{-it_{j}H(s^{j})}$ and $U(t)=\prod_{j=1}^{q}U_{j}(t_{j})$.
According to [30], when $U(t)$ is implemented exactly, then the expected trace
distance between the actual and ideal final states is bounded by

$\mathbb{E}_{t}\left\|U(t)\ket{\psi_{0}}\bra{\psi_{0}}U(t)^{\dagger}-\ket{\psi^{\ast}}\bra{\psi^{\ast}}\right\|_{1}\leq\epsilon.$

In our implementation, each unitary circuit $U_{j}(t_{j})$ is approximated by
a unitary circuit $\tilde{U}_{j}(t_{j})$ satisfying
$\tilde{U}(t)=U(t)+\mathcal{E}(t)$. It follows that

$\displaystyle\mathbb{E}_{t}\left\|\tilde{U}(t)\ket{\psi_{0}}\bra{\psi_{0}}\tilde{U}(t)^{\dagger}-\ket{\psi^{\ast}}\bra{\psi^{\ast}}\right\|_{1}$
$\displaystyle=$
$\displaystyle\mathbb{E}_{t}\left\|\tilde{U}(t)\ket{\psi_{0}}\bra{\psi_{0}}\tilde{U}(t)^{\dagger}-U(t)\ket{\psi_{0}}\bra{\psi_{0}}U(t)^{\dagger}+U(t)\ket{\psi_{0}}\bra{\psi_{0}}U(t)^{\dagger}-\ket{\psi^{\ast}}\bra{\psi^{\ast}}\right\|_{1}$
$\displaystyle\leq$
$\displaystyle\mathbb{E}_{t}\left(\left\|\tilde{U}(t)\ket{\psi_{0}}\bra{\psi_{0}}\tilde{U}(t)^{\dagger}-U(t)\ket{\psi_{0}}\bra{\psi_{0}}U(t)^{\dagger}\right\|_{1}+\left\|U(t)\ket{\psi_{0}}\bra{\psi_{0}}U(t)^{\dagger}-\ket{\psi^{\ast}}\bra{\psi^{\ast}}\right\|_{1}\right)$
$\displaystyle\leq$
$\displaystyle\mathbb{E}_{t}\left\|\tilde{U}(t)\ket{\psi_{0}}\bra{\psi_{0}}\tilde{U}(t)^{\dagger}-U(t)\ket{\psi_{0}}\bra{\psi_{0}}U(t)^{\dagger}\right\|_{1}+\epsilon$
$\displaystyle\leq$
$\displaystyle\mathbb{E}_{t}\left(2\|\mathcal{E}(t)\|_{2}+\|\mathcal{E}(t)\|_{2}^{2}\right)+\epsilon,$

where the last inequality holds because of Lemma 4.7.

As we know, $U(t)$ is the product of $q$ circuits $U_{j}(t_{j})$. Each
$U_{j}(t_{j})$ is approximated by $\tilde{U}_{j}(t_{j})$, whose inexactness
comes from Trotterization method and inexact data. We use $\bar{U}_{j}(t_{j})$
to denote the circuit that approximates $U_{j}(t_{j})$ using the same
Trotterization method with accurate data. Let $r_{j}$ be the Trotter number
for $U_{j}(t_{j})$. Then

$\displaystyle\left\|\tilde{U}_{j}(t_{j})-U_{j}(t_{j})\right\|_{2}$
$\displaystyle=\left\|\tilde{U}_{j}(t_{j})-\bar{U}_{j}(t_{j})+\bar{U}_{j}(t_{j})-U_{j}(t_{j})\right\|_{2}$
$\displaystyle\leq\left\|\tilde{U}_{j}(t_{j})-\bar{U}_{j}(t_{j})\right\|_{2}+\left\|\bar{U}_{j}(t_{j})-U_{j}(t_{j})\right\|_{2}$

The first term measures the error introduced by inexact data. There are
$r_{j}$ circuits in both $\tilde{U}_{j}(t_{j})$ and $\bar{U}_{j}(t_{j})$ since
they are both Trotterized with Trotter number $r_{j}$, i.e.,

$\tilde{U}_{j}(t_{j})=\prod_{l=1}^{r_{j}}\tilde{U}_{j}(t_{j}/r_{j}),\
\bar{U}_{j}(t_{j})=\prod_{l=1}^{r_{j}}\bar{U}_{j}(t_{j}/{r_{j}}).$

According to Lemma 4.6 and Lemma 4.4, we have

$\displaystyle\left\|\tilde{U}_{j}(t_{j})-\bar{U}_{j}(t_{j})\right\|_{2}$
$\displaystyle\leq
r_{j}\left\|\tilde{U}_{j}(t_{j}/r_{j})-\bar{U}_{j}(t_{j}/r_{j})\right\|_{2}$
$\displaystyle\leq t_{j}\epsilon_{0}.$

By the definition of $t_{j}$ in Algorithm 1, $t_{j}\leq 2\pi\kappa$. It
follows that

$\left\|\tilde{U}_{j}(t_{j})-\bar{U}_{j}(t_{j})\right\|_{2}\leq
2\pi\epsilon^{2}/\log^{2}\kappa.$

As for the second term, it measures the error introduced by Trotterization
method. According to Lemma 2.2 and Lemma 4.5, we have
$\|\bar{U}_{j}(t_{j})-U_{j}(t_{j})\|_{2}=\mathcal{O}(\epsilon^{2}/\log^{2}\kappa)$,
provided that Trotter number is bounded by
$\mathcal{O}\left({m^{2}d^{4}t_{j}^{2}}/{\epsilon^{2}}\right)$. In this work,
we choose Trotter number be the upper bound of this bound, which guarantees
$\|\bar{U}_{j}(t_{j})-U_{j}(t_{j})\|_{2}=\mathcal{O}(\epsilon^{2}/\log^{2}\kappa)$.
Put these together, we have

$\left\|\tilde{U}_{j}(t_{j})-U_{j}(t_{j})\right\|_{2}=\mathcal{O}(\epsilon^{2}/\log^{2}\kappa).$

Following Lemma 4.6, we have

$\displaystyle\mathbb{E}_{t}\left\|\mathcal{E}(t)\right\|_{2}$
$\displaystyle=\mathbb{E}_{t}\left\|\tilde{U}(t)-U(t)\right\|_{2}$
$\displaystyle\leq\mathbb{E}_{t}\sum_{j=1}^{q}\left\|\tilde{U}_{j}(t_{j})-U_{j}(t_{j})\right\|_{2}$
$\displaystyle=\mathcal{O}(\epsilon),$

where the last inequality holds by the definition of $q$ in Algorithm 1.
Finally, we have

$\mathbb{E}_{t}\left\|\tilde{U}(t)\ket{\psi_{0}}\bra{\psi_{0}}\tilde{U}(t)^{\dagger}-\ket{\psi^{\ast}}\bra{\psi^{\ast}}\right\|_{1}\leq\mathbb{E}_{t}\left(2\|\mathcal{E}(t)\|_{2}+\|\mathcal{E}(t)\|_{2}^{2}\right)+\epsilon=\mathcal{O}(\epsilon).$

In this implementation of Algorithm 1, there are $q$ circuits
$e^{-it_{j}H(s_{j})}$ and each of them is approximated using $r$ Trottered
circuits. Each Trottered circuits consists of $\mathcal{O}(md^{2})$ time
evolution circuits of Type-1 or Type-2 Hamiltonians. So the total cost of our
implementation is $\mathcal{T}$ classical arithmetic operations, $\mathcal{T}$
single qubit unitary circuits that can be performed in parallel, and
$\mathcal{T}$ calls to $p$-qubit quantum multiplier, where
$\mathcal{T}=\mathcal{O}(qrmnd^{2})=\mathcal{O}(nm^{3}\kappa^{2}d^{6}\log^{2}\kappa/\epsilon^{2})$.

Finally, we discuss the cost to prepare the initial state. We first discuss
the cost to prepare the initial state. The initial state
$\ket{\psi_{0}}=\frac{\sqrt{2}}{2}(\ket{0}-\ket{1})\otimes\ket{b}$ has
$\mathcal{O}(1)$ low tensor rank and thus can be implemented efficiently. We
provide two methods to prepare the initial state and show that in both cases
we can ignore the cost to prepare the initial state. The first method starts
from the fact that, for each $b_{j}=\otimes_{k=1}^{n}b_{jk}$, we have

$\ket{b_{j}}=\otimes_{k=1}^{n}\left(\frac{1}{\|b_{jk}\|_{2}}\begin{bmatrix}b_{jk}&\tilde{b}_{jk}\end{bmatrix}\ket{0}\right),$

where $\tilde{b}_{jk}$ is the complex conjugate of $-iYb_{jk}$ and
$\frac{1}{\|b_{jk}\|_{2}}\begin{bmatrix}b_{jk}&\tilde{b}_{jk}\end{bmatrix}\in\mathbb{C}^{2\times
2}$ is unitary. When $d=1$, we only need $\mathcal{O}(n)$ classical arithmetic
operations and $\mathcal{O}(n)$ single qubit gates to prepare the initial
state, which is negligible compared with the cost to implement Algorithm 2.
When $d=\mathcal{O}(1)$, instead of solving $Ax=b$ directly, we can solve $d$
linear systems $Ax=\otimes_{k=1}^{n}b_{jk}$ for $j\in[d]$. The total
complexity would be $d$ times the complexity to implement Algorithm 2 when
$d=1$. This does not change the asymptotic total complexity since
$d=\mathcal{O}(1)$. One can also prepare the initial state directly without
solving $d$ linear systems. Similar as earlier, it takes $\mathcal{O}(dn)$
classical arithmetic operations and $\mathcal{O}(dn)$ single qubit gates with
gate depth $\mathcal{O}(d)$ to prepare $d$ unitaries $U_{b_{j}}$ such that
$U_{b_{j}}\ket{0}=\ket{b_{j}}$ for $j\in[d]$. We claim the following circuit
can be used to prepare the initial state.

$\Qcircuit<EMAIL_ADDRESS>\lstick{\ket{0}_{n}}\qw\gate{U_{1}}\qw\cdots\gate{U_{d}}\qw\qw\qw\lstick{\ket{\psi}}}$

The first register has $d_{0}=\lceil\log_{2}(d)\rceil$ qubits and represents
states $\ket{1}$ to $\ket{d_{0}}$. For the controlled unitaries, if the first
register is in state $\ket{j}$ with $j\in[d]$, then $U_{b_{j}}$ is applied on
the second register. Finally, when state $\ket{1}$ is measured in the first
register, the second register is in state $\ket{b}$. It takes $\mathcal{O}(d)$
tries to get $\ket{b}$, with which, the initial state can be easily prepared.
In this case, the total cost to prepare the initial state is still negligible.

∎

Notice that, when $d>1$, one can solve $d$ linear systems
$Ax=\otimes_{k=1}^{n}b_{jk}$ for $j\in[d]$ and use the $d$ solutions to
construct the solution of the original linear systems. So we have the
following corollary.

######  Corollary 4.1.

Let $0<\epsilon\leq 1/(3n)$,
$\epsilon_{0}=\epsilon^{2}/(\kappa\log^{2}\kappa)$ and $p=\lceil
1-\log_{2}\epsilon_{0}\rceil$. When using Trotter number
$r=\mathcal{O}(m^{2}d^{4}\kappa^{2}/\epsilon)$ for Trotterization method and
$p$-qubit $\epsilon_{0}$-approximation of $\sigma_{i}$ in Algorithm 2, one can
use our implementation of Algorithm 1 to prepare an
$\mathcal{O}(\epsilon)$-approximate solution using $\mathcal{T}$ classical
arithmetic operations, $\mathcal{T}$ single qubit unitary circuits  and their
controlled version with gate depth $\mathcal{O}(1)$, and $\mathcal{T}$ calls
to $p$-qubit quantum multiplier  and their controlled version with gate depth
$\mathcal{O}(Tp^{3})$, where
$\mathcal{T}=\mathcal{O}(nm^{3}\kappa^{2}d\log^{2}\kappa/\epsilon^{2})$.

##  5 Conclusion

In this work, we study using an adiabatic inspired quantum linear system
algorithm to solve linear system problem in tensor format. We focus on a class
of linear systems whose matrices consist of a linear combination of tensor
products of Hermitian $2$-by-$2$ matrices and linear system vector consists of
a linear combinations of tensor products of $2$-dimensional vectors. We
explicitly describe all the quantum circuits components used in the
implementation of the QLSA. The implementation only uses single qubit unitary
circuits, $p$-qubit quantum multipliers, and their controlled versions. The
number of classical arithmetic operations and the number of these gates are
polynomial of $n$, $m$, and $d$, and thus polylogarithmic in the dimension of
the linear system. Considering the execution time for each of these gates is
$\mathcal{O}(1)$, the total time complexity of the implementation is
polylogarithmic in the problem dimension, which is better than any classical
algorithm for general linear system problems. When compared with the classical
algorithm designed for linear system in tensor format, our total complexity is
comparable to the single step complexity of the classical algorithm proposed
by [5] in terms of the problem dimension.

Possible interesting generalizations of our method would be to extend these
results to higher dimensional $A_{ik}$ and $b_{jk}$ and to cases where
$A_{ik}$ and/or $b_{jk}$ are in different sub-spaces. It is also worth
continuing investigating how to obtain better dependence on condition number
and accuracy for the structured problems.

## Acknowledgement

The authors thank Faisal Alam and Muqing Zheng for helpful discussions. This
work was supported by Defense Advanced Research Projects Agency as part of the
project W911NF2010022: The Quantum Computing Revolution and Optimization:
Challenges and Opportunities. This work was also supported by National Science
Foundation CAREER DMS-2143915 and by U.S. Department of Energy/Office of
Electricity Advanced Grid Modeling program. This research used resources of
the Oak Ridge Leadership Computing Facility, which is a DOE Office of Science
User Facility supported under Contract DE-AC05-00OR22725.

## References

  * Abbas et al. [2023] A. Abbas, A. Ambainis, B. Augustino, A. Bärtschi, H. Buhrman, C. Coffrin, G. Cortiana, V. Dunjko, D. J. Egger, B. G. Elmegreen, et al.  Quantum optimization: Potential, challenges, and the path forward.  _arXiv preprint arXiv:2312.02279_ , 2023.
  * An and Lin [2022] D. An and L. Lin.  Quantum linear system solver based on time-optimal adiabatic quantum computing and quantum approximate optimization algorithm.  _ACM Transactions on Quantum Computing_ , 3(2):1–28, 2022.
  * Apers and Gribling [2023] S. Apers and S. Gribling.  Quantum speedups for linear programming via interior point methods.  _arXiv preprint arXiv:2311.03215_ , 2023.
  * Augustino et al. [2023] B. Augustino, G. Nannicini, T. Terlaky, and L. F. Zuluaga.  Quantum interior point methods for semidefinite optimization.  _Quantum_ , 7:1110, 2023.
  * Ballani and Grasedyck [2013] J. Ballani and L. Grasedyck.  A projection method to solve linear systems in tensor format.  _Numerical Linear Algebra with Applications_ , 20(1):27–43, 2013.
  * Berry et al. [2007] D. W. Berry, G. Ahokas, R. Cleve, and B. C. Sanders.  Efficient quantum algorithms for simulating sparse hamiltonians.  _Communications in Mathematical Physics_ , 270:359–371, 2007.
  * Boixo et al. [2009] S. Boixo, E. Knill, and R. D. Somma.  Eigenpath traversal by phase randomization.  _Quantum Info. Comput._ , 9(9&10):833–855, 2009.
  * Casares and Martin-Delgado [2020] P. A. Casares and M. A. Martin-Delgado.  A quantum interior-point predictor–corrector algorithm for linear programming.  _Journal of Physics A: Mathematical and Theoretical_ , 53(44):445305, 2020.
  * Chakraborty et al. [2018] S. Chakraborty, A. Gilyén, and S. Jeffery.  The power of block-encoded matrix powers: improved regression techniques via faster Hamiltonian simulation.  _arXiv preprint arXiv:1804.01973_ , 2018.
  * Childs and Su [2019] A. M. Childs and Y. Su.  Nearly optimal lattice simulation by product formulas.  _Physical review letters_ , 123(5):050503, 2019.
  * Childs et al. [2017] A. M. Childs, R. Kothari, and R. D. Somma.  Quantum algorithm for systems of linear equations with exponentially improved dependence on precision.  _SIAM Journal on Computing_ , 46(6):1920–1950, 2017.
  * Childs et al. [2021] A. M. Childs, Y. Su, M. C. Tran, N. Wiebe, and S. Zhu.  Theory of Trotter error with commutator scaling.  _Physical Review X_ , 11(1):011020, 2021.
  * Costa et al. [2022] P. C. Costa, D. An, Y. R. Sanders, Y. Su, R. Babbush, and D. W. Berry.  Optimal scaling quantum linear-systems solver via discrete adiabatic theorem.  _PRX Quantum_ , 3(4):040303, 2022.
  * Dalzell et al. [2023] A. M. Dalzell, S. McArdle, M. Berta, P. Bienias, C.-F. Chen, A. Gilyén, C. T. Hann, M. J. Kastoryano, E. T. Khabiboulline, A. Kubica, et al.  Quantum algorithms: A survey of applications and end-to-end complexities.  _arXiv preprint arXiv:2310.03011_ , 2023.
  * Deutsch and Jozsa [1992] D. Deutsch and R. Jozsa.  Rapid solution of problems by quantum computation.  _Proceedings of the Royal Society of London. Series A: Mathematical and Physical Sciences_ , 439(1907):553–558, 1992.
  * Draper [2000] T. G. Draper.  Addition on a quantum computer.  _arXiv preprint quant-ph/0008033_ , 2000.
  * Farhi et al. [2014] E. Farhi, J. Goldstone, and S. Gutmann.  A quantum approximate optimization algorithm.  _arXiv preprint arXiv:1411.4028_ , 2014.
  * Grasedyck et al. [2013] L. Grasedyck, D. Kressner, and C. Tobler.  A literature survey of low-rank tensor approximation techniques.  _GAMM-Mitteilungen_ , 36(1):53–78, 2013.
  * Grover [1996] L. K. Grover.  A fast quantum mechanical algorithm for database search.  In _Proceedings of the twenty-eighth annual ACM symposium on Theory of computing_ , pages 212–219, 1996.
  * Harrow et al. [2009] A. W. Harrow, A. Hassidim, and S. Lloyd.  Quantum algorithm for linear systems of equations.  _Physical Review Letters_ , 103(15):150502, 2009.
  * Jennings et al. [2023a] D. Jennings, M. Lostaglio, R. B. Lowrie, S. Pallister, and A. T. Sornborger.  The cost of solving linear differential equations on a quantum computer: fast-forwarding to explicit resource counts.  _arXiv preprint arXiv:2309.07881_ , 2023a.
  * Jennings et al. [2023b] D. Jennings, M. Lostaglio, S. Pallister, A. T. Sornborger, and Y. Subaşı.  Efficient quantum linear solver algorithm with detailed running costs.  _arXiv preprint arXiv:2305.11352_ , 2023b.
  * Krovi [2023] H. Krovi.  Improved quantum algorithms for linear and nonlinear differential equations.  _Quantum_ , 7:913, 2023.
  * Liu et al. [2021] J.-P. Liu, H. Ø. Kolden, H. K. Krovi, N. F. Loureiro, K. Trivisa, and A. M. Childs.  Efficient quantum algorithm for dissipative nonlinear differential equations.  _Proceedings of the National Academy of Sciences_ , 118(35):e2026805118, 2021.
  * Lloyd [1996] S. Lloyd.  Universal quantum simulators.  _Science_ , 273(5278):1073–1078, 1996.
  * Mohammadisiahroudi et al. [2023] M. Mohammadisiahroudi, R. Fakhimi, Z. Wu, and T. Terlaky.  An inexact feasible interior point method for linear optimization with high adaptability to quantum computers.  _arXiv preprint arXiv:2307.14445_ , 2023.
  * Nielsen and Chuang [2010] M. A. Nielsen and I. L. Chuang.  _Quantum Computation and Quantum Information_.  Cambridge University Press, 2010.
  * Ruiz-Perez and Garcia-Escartin [2017] L. Ruiz-Perez and J. C. Garcia-Escartin.  Quantum arithmetic with the quantum Fourier transform.  _Quantum Information Processing_ , 16:1–14, 2017.
  * Shor [1994] P. W. Shor.  Algorithms for quantum computation: discrete logarithms and factoring.  In _Proceedings 35th Annual Symposium on Foundations of Computer Science_ , pages 124–134. IEEE, 1994.
  * Subaşı et al. [2019] Y. Subaşı, R. D. Somma, and D. Orsucci.  Quantum algorithms for systems of linear equations inspired by adiabatic quantum computing.  _Physical Review Letters_ , 122(6):060504, 2019.
  * Wu et al. [2023] Z. Wu, M. Mohammadisiahroudi, B. Augustino, X. Yang, and T. Terlaky.  An inexact feasible quantum interior point method for linearly constrained quadratic optimization.  _Entropy_ , 25(2), 2023.