# Vibration Analysis of Cyclic Symmetrical Systems by Quantum Algorithms

Authors: Anuradha Mahasinghe
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

Quantum computers have provided exponentially faster solutions to several physical and engineering problems over existing classical solutions. In this paper we present two quantum algorithms to analyze the forced vibration of a mechanical system with cyclic symmetry. Our main algorithm solves the equation of motion of an undamped and nonelastic rotating system with cyclic symmetry consisting ofn sectors, by encoding the displacements of each sector at time t in a quantum state. The runtime of this algorithm is polylog in both n and t, thus exponentially faster than the analogous classical algorithm. Also we consider damped and elastic systems with cyclic symmetry and present another quantum algorithm to solve it in runtime polylog in n and polynomial in t.


---

## Body

## 1. Introduction

Several quantum algorithms were proven to run exponentially faster than the best-known classical algorithms. These algorithms were built to solve specific mathematical problems, for which the classical techniques seemed inadequate, with respect to resource cost and efficiency. Quantum algorithms designed for tasks such as searching databases [1], prime factoring [2], and phase estimation [3] are of fundamental importance in quantum computing. Many later developments in quantum algorithms over the circuit-gate framework are based on these primary algorithms.A remarkable quantum algorithm designed by Harrow, Hassidim, and Loyd in 2009 [4], which is widely known as theHHL algorithm, was of particular importance in the field in the last decade. Based upon the more primary quantum eigenvalue estimation algorithm, its objective is to solve a sparse linear system of equations efficiently, which is unimaginable in the realm of classical algorithms. More precisely, the HHL algorithm is capable of solving a sparse linear system of n equations and n variables in a runtime and resource cost O(polylog n). It is thus exponentially faster than the most efficient classical algorithm for solving a linear system. Experimental verification of the HHL algorithm can be found in [5].The HHL is one of the powerful quantum algorithms, with promising applications to different fields. Also it has been the key ingredient in several new quantum algorithms, yielding an exponential efficiency gain over best-known classical techniques. Examples include the quantum algorithm to solve thed–dimensional Poisson equation with complexity linear in d [6], one for solving a system of first-order differential equations [7] and another algorithm for curve fitting and calculating the effective resistance of electric networks [8]. It is interesting to see that most of these applications solve specific problems in the physical world, with straightforward applications to engineering. Therefore, it is an interesting task to find more physical and engineering problems, for which moderated versions of the HHL could be used to provide exponentially faster solutions over existing classical algorithms.In this paper, we consider one of the important problems in mechanics and engineering: the forced vibration of a mechanical system with cyclic symmetry, for which we develop efficient quantum algorithms, based on the HHL approach. Vibration analysis of cyclic symmetrical systems had been a subject of interest in different disciplines such as acoustics and mechanical engineering. It is customary to describe the dynamics of such systems using ordinary differential equations, in which the coefficients are distributed in a cyclic symmetrical way. In many situations, large-scale computational models are required to solve them. Solving them classically requires a large resource cost as well as a considerable time, since both scale rapidly with the number of components in the system. This is the context in which we develop our quantum algorithms.The paper is organized as follows: In Section2 we present the equation of motion for the forced vibration of a mechanical system with cyclic symmetry. In Section 3 we describe the relevant mathematical principles and transformations and express the solution in a particular form, which enables efficient quantum implementation. In Section 4 we present the corresponding quantum algorithms, describing why our efficiency gain is impossible with classical or other potential quantum methods. Section 5 contains our discussion of the methods.


## 2. Vibrating Systems with Cyclic Symmetry

A vibrating system with cyclic symmetry, which is also called arotationally periodic system, is a system with n sectors, each of which contains an identical mechanical structure with identical coupling to forward nearest neighbors and to ground [9]. For simplicity, we consider a system in which each sector has a single degree of freedom. It can be considered as a cyclic chain of n identical and identically coupled oscillators with the dynamics of which are captured by the dimensionless transverse displacements.The nature of the system imposes a cyclic structure on the factors effecting these displacements. For example, the acceleration of the system couples the masses to be arranged in a cyclic symmetrical way, and the product of the coupled masses and accelerations will be equal to the force applied to the system, if the other factors effecting the acceleration of the system were negligible. If the factors such as damping and elasticity were to be considered, they also act in a cyclic symmetrical way on each sector.Thus, the main factor which distinguishes such a system from a general mechanical system is what is calledcyclic symmetry. Mathematically, this is characterized by the notion ofcirculant matrices. A circulant matrix C is a square matrix in which the elements of each row of C are identical to those of the previous row but are moved one position to the right and wrapped around [10]. More specifically, the circulant matrix C, which is symbolized by circ(c1,c2,…,cn) is the following square matrix.(1)C=c1c2⋯cncnc1⋯cn-1⋮⋮⋱⋮c2c3⋯c1.In the forthcoming subsection we discuss how the dynamics of vibrating systems are governed by circulant matrices.


### 2.1. Simple Vibrating Systems with Cyclic Symmetry

The motion of the simplest vibrating system with cyclic symmetry can be described by Newton’s second law. We will refer to this as the simple vibrating system with cyclic symmetry. In such a system, the accelerations of then sectors are coupled by a circulant matrix which contains the individual masses of each sector, in which the product with the acceleration vector is equal to the coupled forcing vector, obeying Newton’s second law of motion. Therefore, the motion of a simple vibrating system with cyclic symmetry is described by the equation(2)Mu¨=f,where(3)M=circm1,m2,…,mn,m1,m2,…,mn being the individual masses of the n sectors.Let the forces acting on the coupled system be denoted byf1,f2,…,fn. In a practical perspective, the masses will be distributed in a diagonal matrix. However, as an inspiration to the damped elastic vibration which we will discuss later, we may treat the matrix M as a circulant. The system represents much complex situation in which the masses are supposed to be arranged with circulant symmetry. It should be noted that each fj is time-varied and hence a function of time. Thus, we consider the time-variation of the forcing vector described by(4)fjt=fj0et,where e(t) is a function of time and fj(0) is the jth initial force applied on the coupled system. In most of the practical situations, e(t) is a trigonometric function, in particular either sine or cosine, or a combination of them. The function e(t) in our case is supposed to be a function that can be encoded and efficiently implemented in a quantum circuit. This assumption is due to the fact that forcing functions normally take sine or cosine forms, which can be implemented efficiently in a quantum circuit.


### 2.2. Damped Elastic Vibrating Systems with Cyclic Symmetry

Consider a simple vibrating system with cyclic symmetry, with additional effects due to damping and elasticity. Such a system hasn masses, distributed in a circulant matrix, coupled with accelerations, as in a simple vibrating system with cyclic symmetry, with additional effects (decelerations) due to damping and stiffness. Adjacent masses are supposed to be elastically coupled through linear springs, each with nondimensional stiffness v. An individual oscillator, together with the forward nearest neighbor elastic coupling, forms one fundamental sector. We suppose there are n such sectors in the overall system.This vibrating system is modelled by a system of second-order ordinary differential equations, in which the mass, damping, and stiffness matrices are circulants denoted byM=circ(m1,m2,…,mn), D=circ(d1,d2,…,dn), and K=circ(k1,k2,…,kn), respectively, and u is the displacement vector of the coupled system [9]. Thus, the dynamics of the this system is described by(5)Mu¨+Du˙+Ku=f.In the forthcoming section, we describe the approximated solutions for (2) and (5) in detail. Our quantum algorithms will be motivated by these solutions.


## 2.1. Simple Vibrating Systems with Cyclic Symmetry

The motion of the simplest vibrating system with cyclic symmetry can be described by Newton’s second law. We will refer to this as the simple vibrating system with cyclic symmetry. In such a system, the accelerations of then sectors are coupled by a circulant matrix which contains the individual masses of each sector, in which the product with the acceleration vector is equal to the coupled forcing vector, obeying Newton’s second law of motion. Therefore, the motion of a simple vibrating system with cyclic symmetry is described by the equation(2)Mu¨=f,where(3)M=circm1,m2,…,mn,m1,m2,…,mn being the individual masses of the n sectors.Let the forces acting on the coupled system be denoted byf1,f2,…,fn. In a practical perspective, the masses will be distributed in a diagonal matrix. However, as an inspiration to the damped elastic vibration which we will discuss later, we may treat the matrix M as a circulant. The system represents much complex situation in which the masses are supposed to be arranged with circulant symmetry. It should be noted that each fj is time-varied and hence a function of time. Thus, we consider the time-variation of the forcing vector described by(4)fjt=fj0et,where e(t) is a function of time and fj(0) is the jth initial force applied on the coupled system. In most of the practical situations, e(t) is a trigonometric function, in particular either sine or cosine, or a combination of them. The function e(t) in our case is supposed to be a function that can be encoded and efficiently implemented in a quantum circuit. This assumption is due to the fact that forcing functions normally take sine or cosine forms, which can be implemented efficiently in a quantum circuit.


## 2.2. Damped Elastic Vibrating Systems with Cyclic Symmetry

Consider a simple vibrating system with cyclic symmetry, with additional effects due to damping and elasticity. Such a system hasn masses, distributed in a circulant matrix, coupled with accelerations, as in a simple vibrating system with cyclic symmetry, with additional effects (decelerations) due to damping and stiffness. Adjacent masses are supposed to be elastically coupled through linear springs, each with nondimensional stiffness v. An individual oscillator, together with the forward nearest neighbor elastic coupling, forms one fundamental sector. We suppose there are n such sectors in the overall system.This vibrating system is modelled by a system of second-order ordinary differential equations, in which the mass, damping, and stiffness matrices are circulants denoted byM=circ(m1,m2,…,mn), D=circ(d1,d2,…,dn), and K=circ(k1,k2,…,kn), respectively, and u is the displacement vector of the coupled system [9]. Thus, the dynamics of the this system is described by(5)Mu¨+Du˙+Ku=f.In the forthcoming section, we describe the approximated solutions for (2) and (5) in detail. Our quantum algorithms will be motivated by these solutions.


## 3. Solving Vibrating Systems

### 3.1. Diagonalization of Circulant Matrices

Due to the symmetry in the distribution of elements in a circulant matrix, it has interesting spectral properties. One may find a descriptive analysis of the spectral properties of circulants in [10]. Considering the circulant C=circ(c1,c2,…,cn), we define therepresenter of C as the complex polynomial,(6)pz=c1+c2z+c3z2+⋯+cnzn-1.Then, it turns out that the eigenspectrum of C is the list of complex numbers p(1), p(ω), p(ω2),…,p(ωn-1), where ω is the nth root of unity given by ω=exp(2πi/n). Considering the special circulant matrix of dimensions n×n, named thecyclic forward shift, denoted by Π, defined as,(7)Π=circ0,1,0,0,…,0,it is possible to find an alternative representation of a circulant matrix C, which restates the circulant as a linear combination of powers of the cyclic forward shift, with the same coefficients as in the representer of C, as follows:(8)C=c1I+c2Π+c3Π2+⋯+cnΠn-1.The cyclic forward shift matrix is readily diagonalized by the Fourier matrix of order n. More precisely, it is possible to express the cyclic forward shift matrix as(9)Π=F†ΩF,where F is the n×n Fourier matrix defined by(10)F=111111ωω2⋯ωn-11ω2ω4⋯ω2n-1⋮⋮⋮⋱⋮1ωn-1ω2n-1⋯ωn-1n-1,and(11)Ω=diag1,ω,ω2,…,ωn-1.Then it turns out that also the general circulant matrix C is diagonalized by the Fourier matrix as follows:(12)C=F†diagp1,pω,…,pωn-1F.Here p1,pω,…, and pωn-1 are the eigenvalues of the circulant C. Therefore, the eigenvalues of a circulant are derivable simply by Fourier transforming a single row or a column.


### 3.2. Central Differences

The first- and second-order central differences of a variableu are defined as follows:(13)u˙≈ui+1-ui-12h,(14)u¨≈ui+1-2ui+ui-1h2.In a domain with p+2 points, with the assumption that the value of u would be zero at the two extremes of the domain, that is,(15)u0=up+1=0,and the central differences are expressible in terms of two tridiagonal matrices L1 and L2. We write down the discretized approximations of the derivatives using these matrices,(16)u˙=L1u,(17)u¨=L2u,for u being a vector u1,u2,⋯,upT, while u˙ and u¨ are defined similarly. The matrices L1 and L2 correspond to the first- and second-order central differences and then take the following form:(18)L1=010⋯000-101⋯0000-10⋯000⋮⋮⋮⋱⋮⋮⋮000⋯010000⋯-101000⋯0-10,(19)L2=-210⋯0001-21⋯00001-2⋯000⋮⋮⋮⋱⋮⋮⋮000⋯-210000⋯1-21000⋯01-2.Although the two tridiagonal matrices in (18) and (19) are not circulants, they exhibit a certain degree of symmetry, close to circulant. Such symmetry is studied in mathematics underToeplitz operators. More precisely, a matrix is said to be Toeplitz, if its elements are constant along all diagonals parallel to the main diagonal. An immediate consequence of this definition is that, Toeplitz matrices are closed under matrix addition and scalar multiplication (we will use this fact in approximating the displacements of damped elastic vibrating systems with cyclic symmetry). Further, it follows from this definition that a circulant is indeed a Toeplitz. Also the two noncirculants L1 and L2 are Toeplitz.


### 3.3. Solution of a Simple Vibrating System with Cyclic Symmetry

First we consider the simplest forced vibration in which the damping and elasticity effects can simply be ignored. Recall the motion of such a system is governed by (2), where M is the circulant matrix containing the individual masses; the system is in a coupled state in its present form. Since this condition would not support further analysis, our first step would be the decoupling the system, using diagonalization results of circulants described in Section 3.1.Accordingly, it is straightforward to see that the diagonalized version of (2) takes the form(20)F†ΛMFu¨=f,where ΛM is the diagonal matrix consisting of the eigenvalues λ1,λ2,⋯,λn of the mass circulant M,(21)ΛM=diagλ1,λ2,⋯,λn.Left-multiplying both sides of this equation by F,(22)ΛMFu¨=Ff.For simplicity, we assume the circulant M to be invertible. That is, all eigenvalues of M are non-zero, and the diagonal matrix ΛM is hence invertible, in which the inverse consists of eigenvalue reciprocals of M;(23)ΛM-1=diag1λ1,1λ2,⋯,1λn.However, this assumption can be lifted, using M÷, theMoore–Penrose pseudoinverse of ΛM, in case M is singular. That is equivalent to multiplying both sides of (20) by the pseudoinverse of M, as it can be seen that the pseudoinverse of a circulant C with eigenvalues ξ1,ξ2,⋯,ξn is given by(24)C÷=F†Λ÷F,where(25)Λ÷=diagξ1÷,ξ1÷,⋯,ξ1÷,and(26)ξj÷=1ξj,ifξj≠00,otherwise.Thus, we arrive at the following (in which the eigenvalue reciprocal matrixΛM-1 would be replaced by ΛM÷ if M is singular):(27)u¨=F†ΛM-1Ff.By the spectral properties of circulant matrices, it can be easily seen that the decoupled force of the jth sector at time t (the jth element of Ff) is a function of j, ω, t and f1,f2,⋯,fn. Let this be denoted by G(j,ω,f(t)). Since each entry of the Fourier transformed forcing vector is a linear combination of fj’s, G(j,ω,f(t)) can be alternatively expressed as follows:(28)Gj,ω,ft=Gj,ω,f10,f20,⋯,fn0ht.Since it is possible to derive an expression for the acceleration of each sector at a given time, now we discretize the domain with respect to time and approximate the displacements by central differences. Assuming the initial displacement of each sector is zero, we apply the central difference formula in (17).Then the discretized version of the problem is the linear system:(29)Ax=b,where (i)
A is the real symmetric 3–sparse square matrix of order np, given by A=In⊗L2.(ii)
x is the np×1 vector consisting of n blocks, where each block has p entries. The kth element of the jth block gives u(j,k); the displacement of the jth sector at time t=k.(iii)
b is the np×1 vector consisting of n blocks, where each block has p entries. The kth element of the jth block gives G(j,ω,f1(0),f2(0),⋯,fn(0))h(k).


### 3.4. Solution of a Damped Elastic Vibrating System with Cyclic Symmetry

Recall that the forced vibration of ann–sector damped elastic vibrating system with cyclic symmetry is governed by the system of ordinary differential equations Mu¨+Du˙+Ku=f, where M, D, and K are n×n circulant matrices corresponding to the mass, damping, and stiffness, and f is the forcing vector. Using the diagonalization of circulants we described previously, this can be restated as(30)F†ΛMFu¨+F†ΛDFu˙+F†ΛKFu=f,where ΛM, ΛD, and ΛK are diagonal matrices consisting of the eigenvalues of M, D, and K, respectively. Using the transformation q=Fu and left-multiplying by F, we arrive at the decoupled system:(31)ΛMq¨+ΛDq˙+ΛKq=Ff.In other words, we have a system of n equations, in which the jth equation is expressed as(32)λjq¨+μjq˙+νjq=Ffj,where λj, μj, and νj denote the jth eigenvalue of M, D, and K respectively, while the jth term of the Fourier transformed forcing vector is denoted by Ffj. We keep the same notation as in the case of simple vibrating system:(33)λjq¨+μjq˙+νjq=Gj,ω,ftDiscretizing the domain with p points, again we make use of the central differences given by (16) and (17). Then the numerical solution for the jth system is given by(34)λjL2+μjL1+νjIqj=Ffj,where qj is the p×1 vector given by qj1,qj2,⋯,qjmT, and I is the p×p identity matrix.RecallL1, L2, and I are Toeplitz matrices; also a linear combination of them is Toeplitz. Now, the jth Toeplitz system can be expressed as(35)Tjqj=Gj,ω,ft,where(36)Tj=2λj+μj+νj-λj0⋯000-λj-μj2λj+μj+νj-λj⋯0000-λj-μj2λj+μj+νj⋯000⋮⋮⋮⋱⋮⋮⋮000⋯λj+μj+νj-λj0000⋯-λj-μj2λj+μj+νj-λj000⋯0-λj-μj2λj+μj+νj,Thus, the original problem boils down to solve the linear system:(37)A^x^=b^,where (i)
A^ is the 3–sparse block diagonal matrix of dimensions np×np, in which the jth block is the Toeplitz matrix Tj.(ii)
x^ is the np×1 vector consisting of n blocks, where each block has p entries. The kth element in the jth block gives q(j,k), the Fourier transformed displacement of the jth sector at time t=k.(iii)
b^ is the np×1 vector consisting of n blocks, where each block has p entries. The kth element in the jth block gives G(j,ω,f1(0),f2(0),⋯,fn(0))h(k).Thus, we observe that in both cases, the original problem can be transformed into a sparse linear system of equations. Recall [4] describes an efficient quantum algorithm for this task; we make use of a similar technique as a part of our algorithms in solving the vibration problem.


## 3.1. Diagonalization of Circulant Matrices

Due to the symmetry in the distribution of elements in a circulant matrix, it has interesting spectral properties. One may find a descriptive analysis of the spectral properties of circulants in [10]. Considering the circulant C=circ(c1,c2,…,cn), we define therepresenter of C as the complex polynomial,(6)pz=c1+c2z+c3z2+⋯+cnzn-1.Then, it turns out that the eigenspectrum of C is the list of complex numbers p(1), p(ω), p(ω2),…,p(ωn-1), where ω is the nth root of unity given by ω=exp(2πi/n). Considering the special circulant matrix of dimensions n×n, named thecyclic forward shift, denoted by Π, defined as,(7)Π=circ0,1,0,0,…,0,it is possible to find an alternative representation of a circulant matrix C, which restates the circulant as a linear combination of powers of the cyclic forward shift, with the same coefficients as in the representer of C, as follows:(8)C=c1I+c2Π+c3Π2+⋯+cnΠn-1.The cyclic forward shift matrix is readily diagonalized by the Fourier matrix of order n. More precisely, it is possible to express the cyclic forward shift matrix as(9)Π=F†ΩF,where F is the n×n Fourier matrix defined by(10)F=111111ωω2⋯ωn-11ω2ω4⋯ω2n-1⋮⋮⋮⋱⋮1ωn-1ω2n-1⋯ωn-1n-1,and(11)Ω=diag1,ω,ω2,…,ωn-1.Then it turns out that also the general circulant matrix C is diagonalized by the Fourier matrix as follows:(12)C=F†diagp1,pω,…,pωn-1F.Here p1,pω,…, and pωn-1 are the eigenvalues of the circulant C. Therefore, the eigenvalues of a circulant are derivable simply by Fourier transforming a single row or a column.


## 3.2. Central Differences

The first- and second-order central differences of a variableu are defined as follows:(13)u˙≈ui+1-ui-12h,(14)u¨≈ui+1-2ui+ui-1h2.In a domain with p+2 points, with the assumption that the value of u would be zero at the two extremes of the domain, that is,(15)u0=up+1=0,and the central differences are expressible in terms of two tridiagonal matrices L1 and L2. We write down the discretized approximations of the derivatives using these matrices,(16)u˙=L1u,(17)u¨=L2u,for u being a vector u1,u2,⋯,upT, while u˙ and u¨ are defined similarly. The matrices L1 and L2 correspond to the first- and second-order central differences and then take the following form:(18)L1=010⋯000-101⋯0000-10⋯000⋮⋮⋮⋱⋮⋮⋮000⋯010000⋯-101000⋯0-10,(19)L2=-210⋯0001-21⋯00001-2⋯000⋮⋮⋮⋱⋮⋮⋮000⋯-210000⋯1-21000⋯01-2.Although the two tridiagonal matrices in (18) and (19) are not circulants, they exhibit a certain degree of symmetry, close to circulant. Such symmetry is studied in mathematics underToeplitz operators. More precisely, a matrix is said to be Toeplitz, if its elements are constant along all diagonals parallel to the main diagonal. An immediate consequence of this definition is that, Toeplitz matrices are closed under matrix addition and scalar multiplication (we will use this fact in approximating the displacements of damped elastic vibrating systems with cyclic symmetry). Further, it follows from this definition that a circulant is indeed a Toeplitz. Also the two noncirculants L1 and L2 are Toeplitz.


## 3.3. Solution of a Simple Vibrating System with Cyclic Symmetry

First we consider the simplest forced vibration in which the damping and elasticity effects can simply be ignored. Recall the motion of such a system is governed by (2), where M is the circulant matrix containing the individual masses; the system is in a coupled state in its present form. Since this condition would not support further analysis, our first step would be the decoupling the system, using diagonalization results of circulants described in Section 3.1.Accordingly, it is straightforward to see that the diagonalized version of (2) takes the form(20)F†ΛMFu¨=f,where ΛM is the diagonal matrix consisting of the eigenvalues λ1,λ2,⋯,λn of the mass circulant M,(21)ΛM=diagλ1,λ2,⋯,λn.Left-multiplying both sides of this equation by F,(22)ΛMFu¨=Ff.For simplicity, we assume the circulant M to be invertible. That is, all eigenvalues of M are non-zero, and the diagonal matrix ΛM is hence invertible, in which the inverse consists of eigenvalue reciprocals of M;(23)ΛM-1=diag1λ1,1λ2,⋯,1λn.However, this assumption can be lifted, using M÷, theMoore–Penrose pseudoinverse of ΛM, in case M is singular. That is equivalent to multiplying both sides of (20) by the pseudoinverse of M, as it can be seen that the pseudoinverse of a circulant C with eigenvalues ξ1,ξ2,⋯,ξn is given by(24)C÷=F†Λ÷F,where(25)Λ÷=diagξ1÷,ξ1÷,⋯,ξ1÷,and(26)ξj÷=1ξj,ifξj≠00,otherwise.Thus, we arrive at the following (in which the eigenvalue reciprocal matrixΛM-1 would be replaced by ΛM÷ if M is singular):(27)u¨=F†ΛM-1Ff.By the spectral properties of circulant matrices, it can be easily seen that the decoupled force of the jth sector at time t (the jth element of Ff) is a function of j, ω, t and f1,f2,⋯,fn. Let this be denoted by G(j,ω,f(t)). Since each entry of the Fourier transformed forcing vector is a linear combination of fj’s, G(j,ω,f(t)) can be alternatively expressed as follows:(28)Gj,ω,ft=Gj,ω,f10,f20,⋯,fn0ht.Since it is possible to derive an expression for the acceleration of each sector at a given time, now we discretize the domain with respect to time and approximate the displacements by central differences. Assuming the initial displacement of each sector is zero, we apply the central difference formula in (17).Then the discretized version of the problem is the linear system:(29)Ax=b,where (i)
A is the real symmetric 3–sparse square matrix of order np, given by A=In⊗L2.(ii)
x is the np×1 vector consisting of n blocks, where each block has p entries. The kth element of the jth block gives u(j,k); the displacement of the jth sector at time t=k.(iii)
b is the np×1 vector consisting of n blocks, where each block has p entries. The kth element of the jth block gives G(j,ω,f1(0),f2(0),⋯,fn(0))h(k).


## 3.4. Solution of a Damped Elastic Vibrating System with Cyclic Symmetry

Recall that the forced vibration of ann–sector damped elastic vibrating system with cyclic symmetry is governed by the system of ordinary differential equations Mu¨+Du˙+Ku=f, where M, D, and K are n×n circulant matrices corresponding to the mass, damping, and stiffness, and f is the forcing vector. Using the diagonalization of circulants we described previously, this can be restated as(30)F†ΛMFu¨+F†ΛDFu˙+F†ΛKFu=f,where ΛM, ΛD, and ΛK are diagonal matrices consisting of the eigenvalues of M, D, and K, respectively. Using the transformation q=Fu and left-multiplying by F, we arrive at the decoupled system:(31)ΛMq¨+ΛDq˙+ΛKq=Ff.In other words, we have a system of n equations, in which the jth equation is expressed as(32)λjq¨+μjq˙+νjq=Ffj,where λj, μj, and νj denote the jth eigenvalue of M, D, and K respectively, while the jth term of the Fourier transformed forcing vector is denoted by Ffj. We keep the same notation as in the case of simple vibrating system:(33)λjq¨+μjq˙+νjq=Gj,ω,ftDiscretizing the domain with p points, again we make use of the central differences given by (16) and (17). Then the numerical solution for the jth system is given by(34)λjL2+μjL1+νjIqj=Ffj,where qj is the p×1 vector given by qj1,qj2,⋯,qjmT, and I is the p×p identity matrix.RecallL1, L2, and I are Toeplitz matrices; also a linear combination of them is Toeplitz. Now, the jth Toeplitz system can be expressed as(35)Tjqj=Gj,ω,ft,where(36)Tj=2λj+μj+νj-λj0⋯000-λj-μj2λj+μj+νj-λj⋯0000-λj-μj2λj+μj+νj⋯000⋮⋮⋮⋱⋮⋮⋮000⋯λj+μj+νj-λj0000⋯-λj-μj2λj+μj+νj-λj000⋯0-λj-μj2λj+μj+νj,Thus, the original problem boils down to solve the linear system:(37)A^x^=b^,where (i)
A^ is the 3–sparse block diagonal matrix of dimensions np×np, in which the jth block is the Toeplitz matrix Tj.(ii)
x^ is the np×1 vector consisting of n blocks, where each block has p entries. The kth element in the jth block gives q(j,k), the Fourier transformed displacement of the jth sector at time t=k.(iii)
b^ is the np×1 vector consisting of n blocks, where each block has p entries. The kth element in the jth block gives G(j,ω,f1(0),f2(0),⋯,fn(0))h(k).Thus, we observe that in both cases, the original problem can be transformed into a sparse linear system of equations. Recall [4] describes an efficient quantum algorithm for this task; we make use of a similar technique as a part of our algorithms in solving the vibration problem.


## 4. Quantum Algorithms

### 4.1. QFT and HHL

Quantum Fourier transform (QFT) is the key ingredient of many interesting quantum algorithms, including quantum prime factoring [2], quantum period finding [3], and the Simon’s problem [11]. It is also regarded the quantum version of discrete Fourier transform, which is classically achieved in engineering applications via fast Fourier transform (FFT) [12–14]. Although FFT has runtime O(nlogn), QFT has the significant advantage of having a runtime O(polylog n). Also it can be implemented with O(log⁡n2) quantum gates, consisting of the Hadamard, controlled phase, and SWAP operators.Quantum eigenvalue estimation, which is based on the QFT, is the main tool in the HHL algorithm. Given a unitary matrixU embedded in a quantum circuit, the quantum eigenvalue estimation generates the eigendecomposition of an input vector as the output. Since the HHL algorithm has been used as a part of our algorithms, we present a brief outline of the HHL steps.The objective of the HHL algorithm is to solve a sparse linear systemAx=b, where A is a Hermitian matrix (the condition of Hermitianity can be relaxed by embedding the matrix and its conjugate transpose in a block matrix). It encodes b in a quantum states b and generates a quantum state x, which is proportional to the solution vector x.The key idea is the eigendecomposition of a Hermitian matrix. Considering this fact,b can be expressed as a linear combination of the eigenvectors of A:(38)b=∑lβlvl.When this vector b is made the input to the quantum eigenvalue estimation, the output is the state; ∑lβlξlvl. The next step of the HHL is to append a qubit 0 to this state and make ∑l‍βlξlvl0. Then a controlled rotation is applied on |ξl〉, and the result is(39)∑lβlξlvl1-c2ξl20+cξl1,where c is the normalizing constant, which is the ratio of the minimum and maximum eigenvalues of A. The next step is to measure the last qubit in the standard basis, conditioned on seeing 1. Then the system collapses to the state ∑lβl/ξl‍ξlvl.Since the gates used in the quantum algorithm are unitary operations, it is possible to reverse them, that is, applying the conjugate transposes of the sequence of gates. This procedure is nameduncomputing. The quantum eigenvalue estimation is a procedure in which all the steps are defined by unitary gates, making the uncomputing of eigenvalue estimation possible. After applying the relevant conjugate transposes, the system is in the state ∑l‍βl/ξlvl, which gives the desired outcome x of the algorithm.Consider the Halmos dilationUS of the nonunitary matrix S,(40)US=SI-SS†I-S†S-S.We apply this unitary dilation to ΛM-1, which can be constructed through a mapping, corresponding to the eigenvalues of M. Note that UΛM-1 is a 2–sparse matrix which can be considered a block matrix with four blocks, in which each block is an n×n diagonal matrix. Thus, it can be prepared and applied efficiently on a state. It can be easily seen that(41)F†OOF†UΛM-1FOOFf0=F†ΛM-1FfF†I-ΛM-1†ΛM-1Ff,where f0 is the 2n×1 vector consisting of n× 1 blocks f and zero.An efficient quantum circuit implementation for the unitary dilation can be done by making a Hermitian matrix from that and through sparse Hamiltonian simulation. Similar methods have been adopted in [15], for implementing arbitrary sparse unitaries efficiently. We consider the matrix HΛM-1 in which unitary dilation of ΛM-1 and its conjugate transpose are embedded in.(42)HΛM-1=OUΛM-1UΛM-1†O.It can be easily seen that HΛM-1 is a Hermitian 2–sparse matrix. Also it follows that HΛM-1=1. Therefore,(43)e-iHΛM-1θ=cos⁡θI-isinθHΛM-1.Letting θ=π/2,(44)e-iπ/2HΛM-1=-isinθHΛM-1.Therefore,(45)HΛM-1=ie-iπ/2HΛM-1.Then it is easily verified that, given the input with three registers 10f to e-i(π/2)HM-1, the output would be (46)10⊗F†ΛM-1FfF†I-ΛM-1†ΛM-1Ffin which the desired state is derived upon measurement of the second register, conditioned on seeing 0.The time-varied force applied on the system is an efficiently implementable function in a quantum circuit. This has to be plugged into the system after discretizing the time domain. This can be accomplished through the discretization of the functione(t), encoded in a matrix E. Since e(t) is usually sine or cosine, the corresponding operator E is assumed to be unitary and efficient.


### 4.2. Simple Vibrating Systems with Cyclic Symmetry

Now we present the algorithmic steps in approximating the solutionu(j,t) of the circulant system of ordinary differential equations Mu¨=f. Since a circulant matrix is determined by its eigenvalues, we assume the input which corresponds to the circulant M to be in terms of its eigenvalues. By sparse Hamiltonian simulation, we construct HΛM-1 (in (42), (45)) which will be implemented efficiently. The second input is the original forcing vector, acting on the coupled components of the vibrating system, and it is assumed to be encoded in a quantum state. Whenever a qubit is measured in the algorithm, that is performed with respect to the standard basis.Algorithm1 requires log⁡n qubits to encode f, O(log⁡p) ancillary qubits to perform the Hadamard transformations and to append at appropriate stages. To solve the linear system through sparse Hamiltonian simulation requires O(polylog np) qubits. Recall that all the operations used in the algorithm are efficient, scaling down to O( polylog np). A classical algorithm with analogous steps would require a space of polynomial np to store the information, and also it would run with time complexity O(poly np). Therefore, the above algorithm has an exponential saving, compared to a classical algorithm designed for the same task (more on complexity and the selection of HHL for this task will be described in Section 4.4). In particular, if the classical algorithm is based upon approximating the displacement of each sector by finite differences, it would require a linear system solver based method, in which the HHL approach is exponentially advantageous over it. Also the classical fast Fourier transform (FFT) has complexity nlogn, while the QFT is again exponentially advantageous over the classical FFT, guaranteeing the exponential saving in our algorithm.Algorithm 1:Simple vibrating system.
IN:f, HΛM-1, A=In⊗L2 OUT:u(j,t) (i) Append f to two qubits in state 10 and prepare the state 10f. (ii) Apply the sequence of operations I4⊗Fn, HΛM-1 and I4⊗Fn† respectively to the prepared state, and produce 0 ⊗ F † Λ M - 1 F f 0 + F † I - Λ M - 1 † Λ M - 1 F f 1 . (iii) Measure the first two qubits in the standard basis, conditioned on seeing 00. The outcome is proportional to F † Λ M - 1 F f . (iv) Applying H to log log⁡p qubits in state 0, prepare H⊗logp0⊗log⁡p. (v) Append H⊗log⁡p0⊗log⁡p to the outcome F†ΛM-1Ff of the measurement. (vi) Apply In⊗E to the new state and prepare b. (vii) Run the eigenvalue estimation of eiA, letting b be the input, and get the output ∑‍βlξlvl. (viii) Append a qubit 0 to the output to get ∑‍βlξlvl0. (ix) Perform controlled rotation on ξl. (x) Measure the last qubit in the standard basis, conditioned on seeing 1.The outcome is proportional to ∑‍βl/ξlξlvl. (xi) Uncompute the eigenvalue estimation to get ∑‍βl/ξlvl, which is proportional to the state x encoding the solution.


### 4.3. Damped Elastic Vibrating Systems with Cyclic Symmetry

Given the matrixA^, Algorithm 2 requires almost the same resource cost as Algorithm 1 and runtime. If the entries in A^ are assumed to be in the input then it would be the cost for solving the linear system. In particular, it will be polylog in n and polynomial in t.Algorithm 2:Damped elastic vibrating system.
IN:f, A^ OUT:u(j,t) (i) Applying H to log⁡p qubits in state 0, prepare H⊗log⁡p0⊗log⁡p. (ii) Apply F to f. (iii) Append Ff to H⊗log⁡p0⊗log⁡p. (iv) Apply In⊗E to the new state and prepare b^. (v) Run the eigenvalue estimation of eiA^, letting b^ be the input, and get the output Σβ^lξ^lv^l. (vi) Append a qubit 0 to the output to get ∑β^lξ^lv^l0. (vii) Perform controlled rotation on ξ^l. (viii) Measure the last qubit in the standard basis, conditioned on seeing 1.The outcome is proportional to ∑‍β^l/ξ^lξ^lv^l. (ix) Uncompute the eigenvalue estimation to get ∑‍β^l/ξ^lv^l, which is proportional to the state x^. (x) Apply the permutation P on the outcome of the previous step. (xi) Apply In⊗Fp†. This produces a state which is proportional to an np×1 vector consisting of p blocks, each block has n entries, where the k th element in the j th block gives u(k,j).


### 4.4. Achieving Quantum Speed-Up with HHL

It is worthwhile to discuss further the complexity of our algorithms, describing why these solutions are derivable with the same efficiency by neither classical nor other quantum algorithms. First we compare our methods with existing or possible counterpart classical algorithms. Recall the efficiency of HHL is polylog ofn for a sparse system with n entries, which is scarcely imaginable in the classical realm, also the efficiency of non-HHL techniques we have used here are not achievable classically. We will explain this further by considering some concrete, easily realizable examples and standard facts from the circuit-gate framework of quantum computation.The major techniques we have used here are related to the properties of circulant matrices and unitary-Hermitian embedding of general complex matrices. For instance let us take a cyclic symmetrical system with three components, having masses 1, 2, and 3, respectively. Then the corresponding circulant is(47)M=123312231.It is immediately verified that the circulant in (47) is invertible. Further, it can be seen that the eigenvalues of M-1 are 1/6,-0.5+0.29i, and -0.5-0.29i. Therefore,(48)ΛM-1=16000-0.5+0.29i000-0.5-0.29i.It is clear that any classical implementation of the matrix in (48) on some classical vector requires a minimum of three registers. Accordingly, any classical implementation of this matrix on a classical state when there are n components in the cyclic symmetrical system cannot be done unless we use at least n registers. However, this can be implemented on a quantum state in complexity O(polylogn) [16]. Thus, we implement the matrix-vector product in O(polylogn), which cannot be achieved classically.The exact matrix-vector product in our method is application of the matrixΛM-1 on the already prepared quantum state F|f〉. Recall the time-varied force which governs the system in most real-world applications is a trigonometric function, in particular either sine or cosine, or a combination of them; the function e(t) can be encoded and efficiently implemented in a quantum circuit. For instance, suppose the forces Ff0 that act on the system represented by (47) are given by(49)E=sin⁡t000cos⁡t000sin⁡t,f0=12-1.It is possible to implement the operation Ef0 in O(polylogn) using standard quantum circuits in literature [15, 16].Recall we apply the matrixΛM-1 to the already prepared quantum state Ff and it will undergo multiplication by F† consequently, it is important we discuss on this operation both classically and quantum mechanically. For instance, when we analyze the system given by (47) and (49), the corresponding two Fourier matrices in conjugate pairs are as follows:(50)F=131111-12+32i-12-32i1-12-32i-12+32i,(51)F†=131111-12-32i-12+32i1-12+32i-12-32i. These operations are not implementable on a classical computer unless at least n registers are used. Considering the circulants of order n, we use n×n Fourier transform, for which the fastest possible classical methods are based on fast Fourier transforms (FFT) [17, 18]. Thus, the classical complexity boils down eventually at least to linear of n as the classical bound O(nlogn) cannot be overcome. Contrasting to this, the n×n quantum Fourier transform or theQFT, the quantum counterpart ofFFT, can be implemented in a quantum circuit in O(polylogn) [19].The next instrumental non-HHL operation in our algorithms is the implementation of the Hermitian matrix, in which the Halmos dilation of a diagonal matrix is embedded. It is this Hermitian matrix that we claim to implement with an exponential efficiency gain over classical algorithms. Therefore, it is worthwhile to consider a concrete example in order to gain some intuition on its implementation in classical and quantum realms. Let us consider the system given by (47) and (49) again. The Halmos dilation of ΛM-1 in (48) is given by(52)UΛM-1=0.167000.986000-0.5+0.29i000.816000-0.5-0.29i000.8160.98600-0.1670000.816000.5-0.29i0000.816000.5+0.29i.The Hermitian embedding of the matrix in (52) is given by (53)HΛM-1=0000000.167000.986000000000-0.5+0.29i000.816000000000-0.5-0.29i000.8160000000.98600-0.1670000000000.816000.5-0.29i0000000000.816000.5+0.29i0.167000.98600000000-0.5-0.29i000.816000000000-0.5+0.29i000.8160000000.98600-0.1670000000000.816000.5+0.29i0000000000.816000.5-0.29i000000,Clearly, any classical implementation of the matrix in (53) upon a classical vector can be done in complexity not less than 12. Considering an n×n system, one may require O(polyn) resources and hence the same time complexity. It can be seen from (45) that our Hermitian matrix is expressible as the complex exponent of itself, with some global phase factor which can be ignored. This observation enables us using thesparse Hamiltonian lemma, a well-known remarkable result from quantum computation [15]. It guarantees the Hermitian matrix can be implemented on a quantum state in complexity O(polylogn). Obviously, this can never be done using a classical computer. Other operations such as tensor products with sparse matrices are standard efficient constructions in quantum computation in the circuit+gate framework.Once the efficiency of non–HHL methods are discussed, it is worthwhile to look into the direct application of HHL in our methods. Our main application is solving system in (37), in which both matrix and known vector were already proven by us to be constructible in O(polylogn) time. What is left is solving a sparse linear system, for which we can use the HHL in O(polylogn) time. It is noteworthy that no classical computer can achieve this, as even theconjugate gradient algorithm [20]; the fastest classical method known so far for solving linear systems would not yield the efficiency of the HHL, as the latter is exponentially faster when coming to sparse systems.It is worthwhile considering other available quantum techniques which might provide some insight on analyzing cyclic symmetrical systems. One recent advance is thequantum principal component analysis (QPCA) algorithm [21], which shares several common features with HHL, inclusive of being a quantum protocol that takes time polylog of the dimension such that its classical counterpart takes time polynomial in the dimension of the system. Despite the fact that several steps in our algorithms can still be applicable when the HHL is replaced by QPCA, the key ingredient in the latter is efficient exponentiation of a low-rank dense matrix. The interesting characteristic of cyclic symmetrical systems, in particular when considering the case of a simple vibration is the sparse linear system it eventually reduces down to, with the aid of finite differences. Recalling the HHL particularly makes use of sparsity, it suits better for our purpose over QPCA which considers rank of the matrix over sparsity.Another possible quantum algorithm isquantum neural networks (QNN) [22] which has proven to be applicable in several fields [23]. This too has several similarities including the application of eigenvalue estimation as the key ingredient. However, the reduction of a cyclic symmetrical system by central differences to a sparse linear system is still of high significance, which can only be fully appreciated if either HHL or any other exponentially faster linear system solver was applied. Quantum-based learning algorithms such asquantum reinforcement learning (QRA) too have provided significantly faster results over classical algorithms [24], which could potentially be helpful in vibration analysis as well. Despite the probabilistic component, the HHL is a direct algorithm built over the circuit-gate framework of quantum computation, for which the complexity is explicitly expressible in terms of the size of the system. On the other hand, our HHL-based approaches have analytically guaranteed speed-up over their classical counterparts. This explains our preference for HHL over quantum learning algorithms.The latest quantum computing paradigm which has shown much progress in recent years isquantum annealing, which is also used in the D-Wave machine, the first-ever commercially available quantum computer [25]. This is a hardware meta-heuristic which has been modelled upon the adiabatic framework of quantum computation, aimed at solving NP-hard problem instances efficiently. As for many other heuristics, the complexity of quantum annealing is still not explored completely, although several bounds have been suggested [26] and several numerical results are available in literature [27]. Indeed, many heuristics including certain learning algorithms were designed with the objective of solving NP problem instances in a reasonable time. It is important to mention that the equation of motion of a cyclic symmetrical system can be classically solved in polynomial time; hence our aim is to make it polylog. The same complexity reduction was achieved by HHL as well as its applications [6–8] which appeared consequently. Therefore, our selection of HHL as the base algorithm in this work is justifiable even in the context of quantum algorithms.


## 4.1. QFT and HHL

Quantum Fourier transform (QFT) is the key ingredient of many interesting quantum algorithms, including quantum prime factoring [2], quantum period finding [3], and the Simon’s problem [11]. It is also regarded the quantum version of discrete Fourier transform, which is classically achieved in engineering applications via fast Fourier transform (FFT) [12–14]. Although FFT has runtime O(nlogn), QFT has the significant advantage of having a runtime O(polylog n). Also it can be implemented with O(log⁡n2) quantum gates, consisting of the Hadamard, controlled phase, and SWAP operators.Quantum eigenvalue estimation, which is based on the QFT, is the main tool in the HHL algorithm. Given a unitary matrixU embedded in a quantum circuit, the quantum eigenvalue estimation generates the eigendecomposition of an input vector as the output. Since the HHL algorithm has been used as a part of our algorithms, we present a brief outline of the HHL steps.The objective of the HHL algorithm is to solve a sparse linear systemAx=b, where A is a Hermitian matrix (the condition of Hermitianity can be relaxed by embedding the matrix and its conjugate transpose in a block matrix). It encodes b in a quantum states b and generates a quantum state x, which is proportional to the solution vector x.The key idea is the eigendecomposition of a Hermitian matrix. Considering this fact,b can be expressed as a linear combination of the eigenvectors of A:(38)b=∑lβlvl.When this vector b is made the input to the quantum eigenvalue estimation, the output is the state; ∑lβlξlvl. The next step of the HHL is to append a qubit 0 to this state and make ∑l‍βlξlvl0. Then a controlled rotation is applied on |ξl〉, and the result is(39)∑lβlξlvl1-c2ξl20+cξl1,where c is the normalizing constant, which is the ratio of the minimum and maximum eigenvalues of A. The next step is to measure the last qubit in the standard basis, conditioned on seeing 1. Then the system collapses to the state ∑lβl/ξl‍ξlvl.Since the gates used in the quantum algorithm are unitary operations, it is possible to reverse them, that is, applying the conjugate transposes of the sequence of gates. This procedure is nameduncomputing. The quantum eigenvalue estimation is a procedure in which all the steps are defined by unitary gates, making the uncomputing of eigenvalue estimation possible. After applying the relevant conjugate transposes, the system is in the state ∑l‍βl/ξlvl, which gives the desired outcome x of the algorithm.Consider the Halmos dilationUS of the nonunitary matrix S,(40)US=SI-SS†I-S†S-S.We apply this unitary dilation to ΛM-1, which can be constructed through a mapping, corresponding to the eigenvalues of M. Note that UΛM-1 is a 2–sparse matrix which can be considered a block matrix with four blocks, in which each block is an n×n diagonal matrix. Thus, it can be prepared and applied efficiently on a state. It can be easily seen that(41)F†OOF†UΛM-1FOOFf0=F†ΛM-1FfF†I-ΛM-1†ΛM-1Ff,where f0 is the 2n×1 vector consisting of n× 1 blocks f and zero.An efficient quantum circuit implementation for the unitary dilation can be done by making a Hermitian matrix from that and through sparse Hamiltonian simulation. Similar methods have been adopted in [15], for implementing arbitrary sparse unitaries efficiently. We consider the matrix HΛM-1 in which unitary dilation of ΛM-1 and its conjugate transpose are embedded in.(42)HΛM-1=OUΛM-1UΛM-1†O.It can be easily seen that HΛM-1 is a Hermitian 2–sparse matrix. Also it follows that HΛM-1=1. Therefore,(43)e-iHΛM-1θ=cos⁡θI-isinθHΛM-1.Letting θ=π/2,(44)e-iπ/2HΛM-1=-isinθHΛM-1.Therefore,(45)HΛM-1=ie-iπ/2HΛM-1.Then it is easily verified that, given the input with three registers 10f to e-i(π/2)HM-1, the output would be (46)10⊗F†ΛM-1FfF†I-ΛM-1†ΛM-1Ffin which the desired state is derived upon measurement of the second register, conditioned on seeing 0.The time-varied force applied on the system is an efficiently implementable function in a quantum circuit. This has to be plugged into the system after discretizing the time domain. This can be accomplished through the discretization of the functione(t), encoded in a matrix E. Since e(t) is usually sine or cosine, the corresponding operator E is assumed to be unitary and efficient.


## 4.2. Simple Vibrating Systems with Cyclic Symmetry

Now we present the algorithmic steps in approximating the solutionu(j,t) of the circulant system of ordinary differential equations Mu¨=f. Since a circulant matrix is determined by its eigenvalues, we assume the input which corresponds to the circulant M to be in terms of its eigenvalues. By sparse Hamiltonian simulation, we construct HΛM-1 (in (42), (45)) which will be implemented efficiently. The second input is the original forcing vector, acting on the coupled components of the vibrating system, and it is assumed to be encoded in a quantum state. Whenever a qubit is measured in the algorithm, that is performed with respect to the standard basis.Algorithm1 requires log⁡n qubits to encode f, O(log⁡p) ancillary qubits to perform the Hadamard transformations and to append at appropriate stages. To solve the linear system through sparse Hamiltonian simulation requires O(polylog np) qubits. Recall that all the operations used in the algorithm are efficient, scaling down to O( polylog np). A classical algorithm with analogous steps would require a space of polynomial np to store the information, and also it would run with time complexity O(poly np). Therefore, the above algorithm has an exponential saving, compared to a classical algorithm designed for the same task (more on complexity and the selection of HHL for this task will be described in Section 4.4). In particular, if the classical algorithm is based upon approximating the displacement of each sector by finite differences, it would require a linear system solver based method, in which the HHL approach is exponentially advantageous over it. Also the classical fast Fourier transform (FFT) has complexity nlogn, while the QFT is again exponentially advantageous over the classical FFT, guaranteeing the exponential saving in our algorithm.Algorithm 1:Simple vibrating system.
IN:f, HΛM-1, A=In⊗L2 OUT:u(j,t) (i) Append f to two qubits in state 10 and prepare the state 10f. (ii) Apply the sequence of operations I4⊗Fn, HΛM-1 and I4⊗Fn† respectively to the prepared state, and produce 0 ⊗ F † Λ M - 1 F f 0 + F † I - Λ M - 1 † Λ M - 1 F f 1 . (iii) Measure the first two qubits in the standard basis, conditioned on seeing 00. The outcome is proportional to F † Λ M - 1 F f . (iv) Applying H to log log⁡p qubits in state 0, prepare H⊗logp0⊗log⁡p. (v) Append H⊗log⁡p0⊗log⁡p to the outcome F†ΛM-1Ff of the measurement. (vi) Apply In⊗E to the new state and prepare b. (vii) Run the eigenvalue estimation of eiA, letting b be the input, and get the output ∑‍βlξlvl. (viii) Append a qubit 0 to the output to get ∑‍βlξlvl0. (ix) Perform controlled rotation on ξl. (x) Measure the last qubit in the standard basis, conditioned on seeing 1.The outcome is proportional to ∑‍βl/ξlξlvl. (xi) Uncompute the eigenvalue estimation to get ∑‍βl/ξlvl, which is proportional to the state x encoding the solution.


## 4.3. Damped Elastic Vibrating Systems with Cyclic Symmetry

Given the matrixA^, Algorithm 2 requires almost the same resource cost as Algorithm 1 and runtime. If the entries in A^ are assumed to be in the input then it would be the cost for solving the linear system. In particular, it will be polylog in n and polynomial in t.Algorithm 2:Damped elastic vibrating system.
IN:f, A^ OUT:u(j,t) (i) Applying H to log⁡p qubits in state 0, prepare H⊗log⁡p0⊗log⁡p. (ii) Apply F to f. (iii) Append Ff to H⊗log⁡p0⊗log⁡p. (iv) Apply In⊗E to the new state and prepare b^. (v) Run the eigenvalue estimation of eiA^, letting b^ be the input, and get the output Σβ^lξ^lv^l. (vi) Append a qubit 0 to the output to get ∑β^lξ^lv^l0. (vii) Perform controlled rotation on ξ^l. (viii) Measure the last qubit in the standard basis, conditioned on seeing 1.The outcome is proportional to ∑‍β^l/ξ^lξ^lv^l. (ix) Uncompute the eigenvalue estimation to get ∑‍β^l/ξ^lv^l, which is proportional to the state x^. (x) Apply the permutation P on the outcome of the previous step. (xi) Apply In⊗Fp†. This produces a state which is proportional to an np×1 vector consisting of p blocks, each block has n entries, where the k th element in the j th block gives u(k,j).


## 4.4. Achieving Quantum Speed-Up with HHL

It is worthwhile to discuss further the complexity of our algorithms, describing why these solutions are derivable with the same efficiency by neither classical nor other quantum algorithms. First we compare our methods with existing or possible counterpart classical algorithms. Recall the efficiency of HHL is polylog ofn for a sparse system with n entries, which is scarcely imaginable in the classical realm, also the efficiency of non-HHL techniques we have used here are not achievable classically. We will explain this further by considering some concrete, easily realizable examples and standard facts from the circuit-gate framework of quantum computation.The major techniques we have used here are related to the properties of circulant matrices and unitary-Hermitian embedding of general complex matrices. For instance let us take a cyclic symmetrical system with three components, having masses 1, 2, and 3, respectively. Then the corresponding circulant is(47)M=123312231.It is immediately verified that the circulant in (47) is invertible. Further, it can be seen that the eigenvalues of M-1 are 1/6,-0.5+0.29i, and -0.5-0.29i. Therefore,(48)ΛM-1=16000-0.5+0.29i000-0.5-0.29i.It is clear that any classical implementation of the matrix in (48) on some classical vector requires a minimum of three registers. Accordingly, any classical implementation of this matrix on a classical state when there are n components in the cyclic symmetrical system cannot be done unless we use at least n registers. However, this can be implemented on a quantum state in complexity O(polylogn) [16]. Thus, we implement the matrix-vector product in O(polylogn), which cannot be achieved classically.The exact matrix-vector product in our method is application of the matrixΛM-1 on the already prepared quantum state F|f〉. Recall the time-varied force which governs the system in most real-world applications is a trigonometric function, in particular either sine or cosine, or a combination of them; the function e(t) can be encoded and efficiently implemented in a quantum circuit. For instance, suppose the forces Ff0 that act on the system represented by (47) are given by(49)E=sin⁡t000cos⁡t000sin⁡t,f0=12-1.It is possible to implement the operation Ef0 in O(polylogn) using standard quantum circuits in literature [15, 16].Recall we apply the matrixΛM-1 to the already prepared quantum state Ff and it will undergo multiplication by F† consequently, it is important we discuss on this operation both classically and quantum mechanically. For instance, when we analyze the system given by (47) and (49), the corresponding two Fourier matrices in conjugate pairs are as follows:(50)F=131111-12+32i-12-32i1-12-32i-12+32i,(51)F†=131111-12-32i-12+32i1-12+32i-12-32i. These operations are not implementable on a classical computer unless at least n registers are used. Considering the circulants of order n, we use n×n Fourier transform, for which the fastest possible classical methods are based on fast Fourier transforms (FFT) [17, 18]. Thus, the classical complexity boils down eventually at least to linear of n as the classical bound O(nlogn) cannot be overcome. Contrasting to this, the n×n quantum Fourier transform or theQFT, the quantum counterpart ofFFT, can be implemented in a quantum circuit in O(polylogn) [19].The next instrumental non-HHL operation in our algorithms is the implementation of the Hermitian matrix, in which the Halmos dilation of a diagonal matrix is embedded. It is this Hermitian matrix that we claim to implement with an exponential efficiency gain over classical algorithms. Therefore, it is worthwhile to consider a concrete example in order to gain some intuition on its implementation in classical and quantum realms. Let us consider the system given by (47) and (49) again. The Halmos dilation of ΛM-1 in (48) is given by(52)UΛM-1=0.167000.986000-0.5+0.29i000.816000-0.5-0.29i000.8160.98600-0.1670000.816000.5-0.29i0000.816000.5+0.29i.The Hermitian embedding of the matrix in (52) is given by (53)HΛM-1=0000000.167000.986000000000-0.5+0.29i000.816000000000-0.5-0.29i000.8160000000.98600-0.1670000000000.816000.5-0.29i0000000000.816000.5+0.29i0.167000.98600000000-0.5-0.29i000.816000000000-0.5+0.29i000.8160000000.98600-0.1670000000000.816000.5+0.29i0000000000.816000.5-0.29i000000,Clearly, any classical implementation of the matrix in (53) upon a classical vector can be done in complexity not less than 12. Considering an n×n system, one may require O(polyn) resources and hence the same time complexity. It can be seen from (45) that our Hermitian matrix is expressible as the complex exponent of itself, with some global phase factor which can be ignored. This observation enables us using thesparse Hamiltonian lemma, a well-known remarkable result from quantum computation [15]. It guarantees the Hermitian matrix can be implemented on a quantum state in complexity O(polylogn). Obviously, this can never be done using a classical computer. Other operations such as tensor products with sparse matrices are standard efficient constructions in quantum computation in the circuit+gate framework.Once the efficiency of non–HHL methods are discussed, it is worthwhile to look into the direct application of HHL in our methods. Our main application is solving system in (37), in which both matrix and known vector were already proven by us to be constructible in O(polylogn) time. What is left is solving a sparse linear system, for which we can use the HHL in O(polylogn) time. It is noteworthy that no classical computer can achieve this, as even theconjugate gradient algorithm [20]; the fastest classical method known so far for solving linear systems would not yield the efficiency of the HHL, as the latter is exponentially faster when coming to sparse systems.It is worthwhile considering other available quantum techniques which might provide some insight on analyzing cyclic symmetrical systems. One recent advance is thequantum principal component analysis (QPCA) algorithm [21], which shares several common features with HHL, inclusive of being a quantum protocol that takes time polylog of the dimension such that its classical counterpart takes time polynomial in the dimension of the system. Despite the fact that several steps in our algorithms can still be applicable when the HHL is replaced by QPCA, the key ingredient in the latter is efficient exponentiation of a low-rank dense matrix. The interesting characteristic of cyclic symmetrical systems, in particular when considering the case of a simple vibration is the sparse linear system it eventually reduces down to, with the aid of finite differences. Recalling the HHL particularly makes use of sparsity, it suits better for our purpose over QPCA which considers rank of the matrix over sparsity.Another possible quantum algorithm isquantum neural networks (QNN) [22] which has proven to be applicable in several fields [23]. This too has several similarities including the application of eigenvalue estimation as the key ingredient. However, the reduction of a cyclic symmetrical system by central differences to a sparse linear system is still of high significance, which can only be fully appreciated if either HHL or any other exponentially faster linear system solver was applied. Quantum-based learning algorithms such asquantum reinforcement learning (QRA) too have provided significantly faster results over classical algorithms [24], which could potentially be helpful in vibration analysis as well. Despite the probabilistic component, the HHL is a direct algorithm built over the circuit-gate framework of quantum computation, for which the complexity is explicitly expressible in terms of the size of the system. On the other hand, our HHL-based approaches have analytically guaranteed speed-up over their classical counterparts. This explains our preference for HHL over quantum learning algorithms.The latest quantum computing paradigm which has shown much progress in recent years isquantum annealing, which is also used in the D-Wave machine, the first-ever commercially available quantum computer [25]. This is a hardware meta-heuristic which has been modelled upon the adiabatic framework of quantum computation, aimed at solving NP-hard problem instances efficiently. As for many other heuristics, the complexity of quantum annealing is still not explored completely, although several bounds have been suggested [26] and several numerical results are available in literature [27]. Indeed, many heuristics including certain learning algorithms were designed with the objective of solving NP problem instances in a reasonable time. It is important to mention that the equation of motion of a cyclic symmetrical system can be classically solved in polynomial time; hence our aim is to make it polylog. The same complexity reduction was achieved by HHL as well as its applications [6–8] which appeared consequently. Therefore, our selection of HHL as the base algorithm in this work is justifiable even in the context of quantum algorithms.


## 5. Discussion

Considering the real-world engineering applications of vibrating systems with cyclic symmetry, large-scale computational models which require significant resource cost and time must be built in order to analyze such systems in a classical setting. We proposed quantum algorithms for this task, which could earn an exponential saving compared to their classical counterparts. Thus developed quantum methods were based on the well-known HHL algorithm for solving linear systems. The speed-up suggested by the HHL-based methods for solving the equations of motion of vibrating systems with cyclic symmetry is hardly imaginable in a classical setting.In order to make our algorithms more accurate, one might increase the dimension of discretized domain. Further, lettingp very large does not significantly affect the efficiency of the algorithm. However, it will be more advantageous if the optimum p could be determined, considering the practical needs of accuracy and efficiency.



---

*Source: 4969351-2019-06-11.xml*