# Quantum algorithm for linear systems of equations (HHL09): Step 1 - Confusion regarding the usage of phase estimation algorithm

I have been trying to get my head around the famous(?) paper Quantum algorithm for linear systems of equations (Harrow, Hassidim & Lloyd, 2009) (more popularly known as the HHL09 algorithm paper) for some time, now.

On the very first page, they say:

We sketch here the basic idea of our algorithm and then discuss it in more detail in the next section. Given a Hermitian $N\times N$ matrix $A$, and a unit vector $\vec{b}$, suppose we would like to find $\vec{x}$ satisfying $A\vec{x} = \vec{b}$. (We discuss later questions of efficiency as well as how the assumptions we have made about $A$ and $\vec{b}$ can be relaxed.) First, the algorithm represents $\vec{b}$ as a quantum state $|b\rangle = \sum_{i=1}^{N}b_i|i\rangle$. Next, we use techniques of Hamiltonian simulation [3, 4] to apply $e^{iAt}$ to $|b_i\rangle$ for a superposition of different times $t$. This ability to exponentiate $A$ translates, via the well-known technique of phase-estimation [5–7] into the ability to decompose $|b\rangle$ in the eigenbasis of $A$ and to find the corresponding eigenvalues $\lambda_j$ Informally, the state of the system after this stage is close to $\sum_{j=1}^{j=N} \beta_j |u_j\rangle|\lambda_j\rangle$, where $u_j$ is the eigenvector basis of $A$ and $|b\rangle = \sum_{j=1}^{j=N} \beta_j|u_j\rangle$.

So far so good. As described in Nielsen & Chuang in the chapter "The quantum Fourier transform and its applications", the phase estimation algorithm is used to estimate $\varphi$ in $e^{i2\pi \varphi}$ which is the eigenvalue corresponding to an eigenvector $|u\rangle$ of the unitary operator $U$.

Here's the relevant portion from Nielsen & Chuang:

The phase estimation algorithm uses two registers. The first register contains $t$ qubits initially in the state $|0\rangle$. How we choose $t$ depends on two things: the number of digits of accuracy we wish to have in our estimate for $\varphi$, and with what probability we wish the phase estimation procedure to be successful. The dependence of $t$ on these quantities emerges naturally from the following analysis.

The second register begins in the state $|u\rangle$ and contains as many qubits as is necessary to store $|u\rangle$. Phase estimation is performed in two stages. First, we apply the circuit shown in Figure 5.2. The circuit begins by applying a Hadamard transform to the first register, followed by application of controlled - $U$ operations on the second register, with $U$ raised to successive powers of two. The final state of the first register is easily seen to be:

$$\frac{1}{2^{t/2}}\left(|0\rangle+\text{exp}(2\pi i 2^{t-1}\varphi)|1\rangle)(|0\rangle+\text{exp}(2\pi i 2^{t-2}\varphi)|1\rangle)...(|0\rangle+\text{exp}(2\pi i 2^{0}\varphi)|1\rangle\right)= \frac{1}{2^{t/2}}\sum_{k=0}^{2^{t}-1}\text{exp}(2\pi i \varphi k)|k\rangle$$

The second stage of phase estimation is to apply the inverse quantum Fourier transform on the first register. This is obtained by reversing the circuit for the quantum Fourier transform in the previous section (Exercise 5.5) and can be done in $\Theta (t^2)$ steps. The third and final stage of phase estimation is to read out the state of the first register by doing a measurement in the computational basis. We will show that this provides a pretty good estimate of $\varphi$. An overall schematic of the algorithm is shown in Figure 5.3.

To sharpen our intuition as to why phase estimation works, suppose $\varphi$ may be expressed exactly int bits, as $\varphi = 0.\varphi_1 ... \varphi_t$. Then the state (5.20) resulting from the first stage of phase estimation may be rewritten

$$\frac{1}{2^{t/2}}(|0\rangle + \exp(2\pi i 0.\varphi_t|1\rangle)(|0\rangle + \exp(2\pi i 0.\varphi_{t-1}\varphi_t|1\rangle)...(|0\rangle + \exp(2\pi i 0.\varphi_1...\varphi_t|1\rangle)$$

The second stage of phase estimation is to apply the inverse quantum Fourier transform. But comparing the previous equation with the product form for the Fourier transform, Equation (5.4), we see that the output state from the second stage is the product state $|\varphi_1 ...\varphi_t\rangle$. A measurement in the computational basis, therefore, gives us $\varphi$ exactly!

Summarizing, the phase estimation algorithm allows one to estimate the phase $\varphi$ of an eigenvalue of a unitary operator $U$, given the corresponding eigenvector $|u\rangle$. An essential feature at the heart of this procedure is the ability of the inverse Fourier transform to perform the transformation

$$\frac{1}{2^{t/2}}\sum_{j = 0}^{2^t-1}\exp(2\pi i \varphi j)|j\rangle |u\rangle \to |\tilde \varphi \rangle |u\rangle$$

Let's proceed from here. I found a nice circuit diagram for the HHL09 algorithm here[$\dagger$]:

## Step 1 (Phase Estimation):

In the first step of the HHL09 algorithm the same concept (of the standard Quantum Phase Estimation algorithm as described in Nielsen and Chuang) is used. However, we must keep in mind that $A$ by itself isn't a unitary operator. However, if we assume that $A$ is Hermitian then the exponential $e^{iAt}$ is unitary (no worries, there's exists a workaround in case $A$ isn't Hermitian!).

Here, we can write $U=e^{iAt}$. There's another subtle point involved here. We do not know the eigenvectors $|u_j\rangle$ of $U$ beforehand (but we do know that for any unitary matrix of size $N\times N$ there exist $N$ orthonormal eigenvectors). Moreover, we need to remind ourselves that if the eigenvalues of $A$ are $\lambda_j$ then the eigenvalues of $e^{iAt}$ will be $e^{i \lambda_j t}$. If we compare this with the form of eigenvalues given in Nielsen and Chuang for $U$ i.e. if $e^{2\pi i \varphi} \equiv e^{ i \lambda_j t}$, we'd find $\varphi = \frac{\lambda_j t}{2\pi}$. In this case, we begin in the state $|b\rangle$ (which can be written as a superposition of the eigenvectors of $U$ i.e. $\sum_{j=1}^{j=N}\beta_j|u_j\rangle$) rather than any particular eigenvector $|u_j\rangle$ of $U$, as far as the second register of qubits is concerned. If we had begun in the state $|u\rangle \otimes (|0\rangle)^{\otimes t}$ we would have ended up with $|u\rangle \otimes |\tilde\varphi\rangle$ i.e. $|u_j\rangle \otimes |\tilde{\frac{\lambda_j t}{2\pi}}\rangle$ (considering that $\lambda_j$ is the eigenvalue associated with the eigenvector $|u_j\rangle$ of $A$). Now, instead if we begin in the superposition of eigenvectors $\sum_{j=1}^{j=N}\beta_j|u_j\rangle$ we should end up with $\sum_{j=1}^{j=N}\beta_j|u_j\rangle\otimes |\tilde{\frac{\lambda_j t}{2\pi}}\rangle$.

## Question:

Part 1: In the HHL09 paper, they wrote about the state of the system after this Phase Estimation step is $\sum_{j=1}^{j=N}\beta_j|u_j\rangle\otimes |\tilde\lambda_j\rangle$. However, from what I wrote above it seems to me that the state of the system should rather be $\sum_{j=1}^{j=N}\beta_j|u_j\rangle\otimes |\tilde{\frac{\lambda_j t}{2\pi}}\rangle$.

What am I missing here? Where did the factor of $\frac{t}{2\pi}$ vanish in their algorithm?

Edit: Part 2 has been asked here to make the individual questions more focused.

I also have several confusions regarding Step 2 and Step 3 of the HHL09 algorithm too, but I decided to post them as separate question threads, as this one is becoming too long. I'll add the links to those question threads, on this post, once they are created.

[$\dagger$]: Homomorphic Encryption Experiments on IBM's Cloud Quantum Computing Platform Huang et al. (2016)

• @Nelimee That $6$ comes from the formula $t = 3 + \lceil { \log_2(2+\frac{1}{2 (0.1)})\rceil} = 3 + 3 = 6$. It denotes the number of qubits in the "first register" needed to represent each $|\lambda_j\rangle$ or $|\frac{\lambda_j t}{2\pi}\rangle$ to $3$-bits of precision and with $90\%$ accuracy. Btw, please note that part of the question has now been shifted here. – Sanchayan Dutta Jun 19 '18 at 8:45

It depends on the papers but I saw 2 approaches:

1. In most of the papers I read about the HHL algorithm and its implementation, the Hamiltonian evolution time $t$ is taken such that this factor disappear, i.e. $t = t_0 = 2\pi$.

2. The approximate eigenvalue is often written $\tilde \lambda$. In some paper this notation really means "the approximation of the true eigenvalue $\lambda$" and in other papers they seems to include $\frac{t}{2\pi}$ in this definition, i.e. "$\tilde \lambda$ is the approximation of the value of $\frac{\lambda t}{2\pi}$".

1. Quantum linear systems algorithms: a primer (Dervovic, Herbster, Mountney, Severini, Usher & Wossnig, 2018): a complete and very good article on HHL algorithm and some improvements that has been discovered. The paper is from the 22th of February, 2018. The value of $t$ you are interested in is first addressed in page 30, in the legend of Figure 5 and is fixed at $2\pi$.

2. Quantum Circuit Design for Solving Linear Systems of Equations (Cao, Daskin, Frankel & Kais, 2013) (take the v2 and not the v3): a detailed implementation of HHL algorithm for a fixed 4x4 matrix. If you plan to use the article let me warn you that there are some mistakes in it. I can provide you the ones I found if you are interested. The value for $t$ (which is denoted as $t_0$ in this paper) is fixed to $2\pi$ in the second page (at the start of the right column).

3. Experimental Quantum Computing to Solve Systems of Linear Equations (Cai, Weedbrook, Su, Chen, Gu, Zhu, Li, Liu, Lu & Pan, 2013): an implementation of HHL algorithm for a 2x2 matrix on an experimental setup. They fix $t = 2\pi$ in the legend of Figure 1.

4. Experimental realization of quantum algorithm for solving linear systems of equations (Pan, Cao, Yao, Li, Ju, Peng, Kais & Du, 2013): implementation of HHL for a 2x2 matrix. The implementation is similar to the one given in the second point above, with the 4x4 matrix. They fix $t_0 = 2\pi$ in page 3, bullet point n°2.

What am I missing here? Where did the factor of $\frac{t}{2\pi}$ vanish in their algorithm?

Remember that in Dirac notation, whatever you write inside the ket is an arbitrary label referring to something more abstract. So, it is true that you are finding the (approximate) eigenvector to $U$, which has eigenvalue $e^{-i\lambda t}$ and therefore what you're extracting is $\lambda t/(2\pi)$, but that is the same as the eigenvector of $A$ with eigenvalue $\lambda$, and it is that which is being referred to in the notation. But if you wanted to be really clear, you could write it as

|approximate eigenvector of $U$ for which eigenvalue is $e^{-i\lambda t}$ and of $A$ for which eigenvalue if $\lambda\rangle$,

but perhaps instead of writing that out every time, we might just write $|\tilde\lambda\rangle$ for brevity!