Solving linear equation systems on noisy intermediate‐scale quantum computers

Quantum computing promises classically unparalleled benefits for various applications. Its properties are exploited in the Harrow‐Hassidim‐Lloyd (HHL) algorithm that, in conjunction with quantum phase estimation, is capable of constructing quantum states that are proportional to the solution of linear equation systems and does so exponentially faster than the fastest known classical algorithms. We explore this capability by computing the nodal displacements of a 1‐dimensional loaded cantilever, discretized by using the finite element method (FEM).


.1 Introduction
A general qubit state can be written as |Ψ = α 1 |0 +α 2 |1 [1], where the coefficients α i are the amplitudes of the qubit state. Akin to classical linear algebra, |0 and |1 are basis vectors, representing classical computing states. For multi-qubit register states, the tensor product of each qubit's basis vectors is taken to get the basis vectors of the full system, resulting in 2 n basis vectors for an n-qubit system. Hence, the available space for performing computations grows exponentially with each added qubit. A computation itself is a rotation U in the Hilbert space spanned by the basis vectors. Applying this rotation to a general n-qubit state can be expressed by applying this rotation to each basis vector [1]. This gives rise to quantum parallelism. If each basis vector encodes one datum of a set of inputs to an algorithm, applying U to |Ψ will result in a linear combination of each individual solution.

The Harrow-Hassidim-Lloyd algorithm
The exponentially increasing computing space of quantum registers can be used to solve sparse linear equation systems faster and taking up less space than classical approaches by using analog-encoding [2], where the information about the original matrices and vectors is encoded in the amplitudes of multi-qubit states. A requirement for HHL is that the matrix must be hermitian, which is automatically given for orthogonal matrices, and of dimension that is some power of 2. Both can be ensured by appropriate transformations [3]. The problem to solve is Kū =f , and the goal is to construct |ū = K −1 f . The trick is that if K is Hermitian, it can be implemented on a quantum computer as a Hamiltonian, creating the time evolution |ū(t) = exp(−itK −1 ) f , which by itself would suffice to gather some global information about the result, such as the largest displacement [3].
To get the full solution |ū , quantum phase estimation [4] is employed. It introduces another qubit register of the same size as f , initialized as |0 . The Hamiltonian here is constructed by tensoring K with the momentum operatorp, which induces a translation in the second qubit register proportional to the eigenvalues λ i of K. The time evolution is where β i are some unknown, irrelevant amplitudes, and |ψ i are the eigenvectors corresponding to the eigenvalues λ i of K. The quantum computer can now use the eigenvalues λ i saved in the second quantum register to apply the phase U φ = exp(idλ −1 i ) to the full state, after which the first computation from equation 1 is undone by applying the inverse transformation to get where the linearization on the right side is sound if d is sufficiently small. The relevant part idK −1 f can be extracted in logarithmic time [4]. This procedure can be modified to get the extraction in poly( κ 2 s 2 log ε log N ) [6]. A classical conjugated gradient solver for sparse matrices has complexity O( κs log ε N log N ) [5], where s N is the sparsity of the stiffness matrix, κ = |λmax| |λmin| its condition number, and ε the precision. The modified HHL is much more sensitive to the sparsity and condition number of the system matrix, but exponentially faster in the input size N .  Fig. 1: a: The discretized cantilever, with applied boundary condition and load, and b: comparison of reference configuration, analytical and quantum-computed deformation.

Cantilever example problem
Since current hardware is limited to a small number of qubits, the problem size must be small enough to fit on available devices. Here, we used IBM's "qasm simulator" and qiskit [7] to test the problem setup. This simulator uses a noise model reflecting the noisy environment on real hardware. The results are sampled from noisy versions of the quantum circuit. We analyze a cantilever discretized into 4 truss elements, bounded by 5 nodes, fixed at node n 1 and loaded at node n 5 (cf. figure 1a). After applying the boundary condition, the reduced system is described by K ∈ R 4×4 andf ∈ R 4 , so log 2 (4·4)+log 2 (4)+log 2 (4) = 8 qubits are needed. The whole circuit has a depth (total number of quantum gate operations) of 110, of which 58 gates are CNOTs, which are especially noisy. We use a simple Jacobi preconditioner P −1 ij = δij Kij to further improve the condition number of K. The deformed cantilever is depicted in figure 1b, comparing the analytical and the mean quantum-computed displacements along with the standard deviation for each nodal displacement. 100 noisy runs of the circuit were performed, summarized in table 1.

Summary
Figure 1b seems to show convincing agreement of the analytical and quantum-computed solutions, table 1 shows a deviation of the quantum-computed mean from the analytical solution of ∆u 2 = 13.3%. This is a significant error, although the analytical result lies inside the standard deviation. The standard deviation is high in all cases, ranging from 6.8% up to 29.8%. In a simple system like a cantilever, these errors do not seem catastrophic for the final result as long as the displacements are small, but they will accumulate for larger systems and cause noticeable deviations from the analytical solution. The main source of these errors is likely the high number of CNOT gates. A better approach might be to use a variational quantum eigensolver [8], which is only quadratically faster than classical algorithms, but less prone to noise by depending on fewer CNOT gates.