Solving Lyapunov equation by quantum algorithm

He SUN,Jing ZHANG

Department of Automation,Tsinghua University,Beijing 100084,China

Solving Lyapunov equation by quantum algorithm

He SUN,Jing ZHANG?

Department of Automation,Tsinghua University,Beijing 100084,China

Lyapunov equation is one of the most basic and important equations in control theory,which has various applications in,e.g.,stability analysis and robust analysis of linear control systems.Inspired by the recent progresses of quantum algorithms,we find that solving Lyapunov equation can be exponentially accelerated by quantum algorithms rather than traditional classical algorithms.Our algorithm is more efficient especially when the system matrix is sparse and has a low condition number.The results presented in this paper open up new dimensions of research in controlling classical system by quantum information processors,which has rarely been considered in the existing literature.

Quantum control,linear Lyapunov equation,large-scale complex system

1 Introduction

Linear systems are a class of physical systems which satisfy the superposition principle and homogeneity between the input and output.Using a set of first-order differential equations,we can describe the dynamics of the linear systems by the state-space representation,analyze the qualitative properties of the system such as stability,controllability and observability,and solve the control laws to steer the system behaviors by digital signal processors such as computers.A considerable number of studies have been focused on the improvement of the numerical algorithms for solving the control laws to control linear systems in the past few decades[1–3].

For designing the control to steer practical linear systems,it is crucial to find an efficient way considering the high dimension of their dynamical equations.Such large-scale linear systems can be easily seen in various systems,such as power networks,urban traffic networks,digital communications networks,flexible manufacturing networks,ecological systems and economic systems,etc.[4].Another significant source of largescale linear control systems comes from discretizing the partial derivatives of spatial coordinates for a system dominated by a partial differential equation[5].These have posed a challenge to the classical computers,with increasing scale of the practical system in contrast to the limited computational capacity of the classical computers.When the miniaturization of integrated circuits reaches a certain degree,the heat dissipation problem and the quantum effects will make the increasing of the computational capacity of the classical computers difficult to continue as declared by the Moore’s law.

The main strategy dealing with large-scale linear systems currently is to exploit the special structure,such as sparseness,of the system matrices.Some approaches based on the sparse-matrix algorithms can reduce the computational complexity to some extent[5,6].In recent few decades,a new area called“quantum computation”arose and soon came into the researchers’attention[7].Pointed out by Feynman in 1982,quantum computation is a kind of computational theory based on the superposition principle of quantum mechanics,which totally differs from classical computation[8].Due to the quantum parallelism induced by quantum superposition principle and quantum entanglement,various quantum algorithms are proposed and proven to be much more efficient than the classical algorithms in dealing with some specific problems.Even some problems with unaccessible computational complexity on classical computers can be effectively solved using quantum algorithms.One possible example is the Shor’s algorithm[9,10],which shows that factoring an integer into the product of the prime factors can be done in polynomial time by quantum computers,rather than the subexponential time by classical computers.Another possible example is the Grover’s algorithm,which achieves a quadratic acceleration compared with the classical algorithm in searching data in an unordered database[11].Hitherto,there have been various quantum algorithms applied in algebraic and number theory[12,13],graph theory[14–16],pattern matching[17–19],quantum simulation[20–22],and many other related fields.

A more recent progress of quantum algorithm in 2009 shows that solving linear algebraic equation Ax=b can be exponentially speeded up by quantum algorithm[23],especially when the matrix A is sparse and the condition number(the ratio between the largest and smallest eigenvalue)of A is low.So far,this algorithm has been applied to various related problems such as solving ordinary differential equations[24–26],data fitting[27],and support vector machine[28].Several experiments have also been carried out to verify these quantum algorithms[29,30]in optical and nuclear magnetic resonance(NMR)systems.Borrowing idea from these quantum algorithms,we here propose a quantum algorithm to solve the Lyapunov equation,for which the classical numerical methods need a computing time no less than O(n3)and a memory space of O(n2)[31],where n is the order of the Lyapunov equation.It is shown that an norder Lyapunov equation can be solved by our quantum algorithm with computational complexity of O(log(n))and thus is exponentially speeded up.The paper is organized as follows:in Section2,we convert the Lyapunov equation from the matrix form into vector form by tensor product.Our quantum algorithm and the corresponding quantum circuits are presented in Section3.In Section4,we give a numerical example of our algorithm.The conclusion and forecast of the future work are drawn in Section5.

NotationIn this paper,we use standard notation of quantum mechanics.|φ〉denotes the state of a quantum system,which is a normalized column vector defined on the complex Hilbert space,?indicates a tensor product.We use log(x)to denote a logarithm of base 2 throughout the paper.

2 Problem formulation

The n-order continuous-time Lyapunov equation we are interested in can be expressed as

where A is an n×n square matrix and ATis the transpose of A.Q=QTis a non-negative definite matrix.We use k1,k2,...,knto denote the eigenvalues of A,then the Lyapunov equation has a unique solution of Hermitian n×n matrix X if and only if ki+kj≠0 for all i and j.If the solution X is non-negative definite,then all the eigenvalues of A must have negative real part.Lyapunov equation plays a very important role in control theory,communications and dynamic systems.

The Lyapunov equation can be transformed into a linear algebraic equation in the Kronecker product form.Let A?B denote the Kronecker product of two n-order matrices A and B,i.e.,

and define vec(C)as

where ciis the ith column vector of the n-order matrix C.Thus,the Lyapunov equation(1)can be converted into the following linear algebraic equation:

where I is the n×n identity matrix.Here,the solution of equation vec(X)contains all the elements of the matrix X.The coefficient matrix of equation(4),i.e.,H=I?AT+A?I,can then be expressed as

3 Quantum algorithm and circuit for Lyapunov equation

Based on the preparations given in Section2 and borrowing ideas from[23],we give our quantum algorithm for solving Lyapunov equation as follows:

1)To solve equation(4),we first need to encode the vector?vec(Q)into the following quantum state

where qiis the ith element of the vector?vec(Q)and|i〉is the ith basis vector of the quantum Hilbert space.Usually,we set Q to be identity matrix I,then we have

The above quantum state can be efficiently prepared using the method in[32].|?vec(Q)〉can also be expressed as

where|uj〉is the jth eigenvector of the matrix H and βjis the corresponding coefficient.

2)We use the technique of phase estimation(see,e.g.,[33])to find the eigenvalues λ1,λ2,...,λn2of the matrix H and store them in another quantum register(labelled as M here).The quantum circuit[7]for phase estimation is given in Fig.1.

This step is first done by conducting the Walsh-Hadamard transformation(given by the “H”gates in Fig.1)onto the register M which composed of m qubits all initialized in|0〉to obtain a superposition state with equal weights of each computational basis states{|τ〉}

where m must be large enough to guarantee the computational accuracy of the phase estimation algorithm.Then we apply the operationto the quantum state|Ψ〉? |? vec(Q)〉,and in addition an inverse Fourier transformation to the state on the register M.After that,we can obtain the quantum state

The most vital part of phase estimation is the quantum simulation of the controlled Hamiltonian which produces the transformation eiHtfrom an Hermitian matrix H with t conditioned by the eigenvalue λjs.It is pointed out that eiHtcan be simulated in time O(log(n)s4t)when the Hermitian matrix H is sparse[34].Here,the word“sparse”means that H has at most s(? n)nonzero entries in each row which can be read with time in scale of O(s).Due to this reason,we here assume A to be symmetric and sparse.As a result,H=I?AT+A?I is also symmetric and sparse,thus eiHtcan be efficiently simulated.

Fig.1 Quantum circuit for phase estimation.Two quantum registers are needed in this circuit:one is used to store the vector?vec(Q),the other,which is labelled as M,is initialized in|0〉?m,and is used to store the eigenvalues.The phase estimation is conducted by Walsh-Hadamard transformation,controlled Hamiltonian simulations and inverse Fourier transformation.

3)Add an ancilla qubit|0〉and apply rotation gates to it.The rotation gates are controlled by register M.Conditioned on|λj〉,the rotation gates transform the ancilla qubit from|0〉towhere C is a normalization constant.After this step,the state of the system can be expressed as

This step can be realized via the quantum circuit given in Fig.2.

Fig.2 Quantum circuit for the controlled rotation.2msingle-qubit gates which are controlled by register M are needed.When the state of register M is|i〉,the gate Riwill act on the ancilla qubit and change the state of ancilla qubit from|0〉to

4)Finally,we apply the inverse phase estimation to the system,which change the state of the register M into the initial states|0〉?m,and we get the quantum state

Then we measure the ancillary qubit.If the state of ancilla qubit after measurement is|1〉,then the quantum state of system will collapse to

That is the final result we need.Otherwise,the computation process fails,we need to repeat the above steps until we get a successful result.A brief quantum circuit for solving Lyapunov equation is shown in Fig.3.

Fig.3 Quantum circuit for solving Lyapunov equation.Two multi-qubit quantum registers initialized in|? vec(Q)〉and|0〉?mrespectively,and an ancilla qubit initialized in|0〉is needed.The whole process can be divided into four steps:phase estimation,controlled rotation,inverse phase estimation,and measurement.

4 Numerical example

We give a numerical example here to illustrate our quantum algorithm.The system matrix A is chosen to be

a 2×2 symmetric matrix with eigenvalues 1 and 2.Then the Lyapunov equation(1)can be equivalently transformed into a linear algebraic equation with 4 unknowns

The eigenvalues of H are λ1=2, λ2=3, λ3=3,λ4=4,which can be encoded in binaries of three bits,thus three qubits are enough for the storage of these eigenvalues.What’s more,two qubits for the storage of vector b and an ancilla qubit are needed.To sum up,we altogether need six qubits to solve equation(15).The quantum circuit designed to complete the computation is given in Fig.4.

Fig.4 Quantum circuit for the numerical example.Six qubits initialized in|0〉|0〉?3|b〉are needed.The Walsh-Hadamard transformation,controlled Hamiltonians,inverse Fourier transformation,controlled rotation and measurement are given in this circuit.The inverse phase estimation is represented by the“INV”part.

We first initialize the system into the quantum state

and conduct the Walsh-Hadamard transformation.ujand βjdenote the eigenvector and coefficient correspond to λjrespectively.According to the clearly established equation

we can convert the controlled-e2πiHtinto the tensor product of two controlled-e2πiAtgates,as U1,U2,U3shown in Fig.4.Here,

After this,we obtain the quantum state:

The next part,which consists of one swap gate,three Hadamard gates and three controlled single-qubit gates,is for the inverse Fourier transformation.Here,

The quantum state of the system then becomes

Subsequently,the controlled rotation part computes the reciprocal of eigenvalues and transfers that information into the state of ancilla qubit.To achieve this,we need three controlled single-qubit gates

each conditioned by three qubits and act on the eigenvalue 2,3,4,respectively.Then we get the quantum state

After that,we apply the inverse phase estimation to the system(represented as“INV”in Fig.4),namely,we apply the inverse of all the previous procedure except the controlled rotation.As a result,we get

Finally,we take a measurement to the ancilla qubit.If the measure result is|1〉,then the two qubits initially store|b〉would collapse to the state

which is the final result after normalization.

5 Conclusions

In this paper,we present a quantum algorithm for solving linear continuous-time Lyapunov equation.This algorithm can be finished with time scale of O(logn)for an n-order Lyapunov equation,which means that an exponential acceleration can be realized by our quantum algorithms in comparison to the existing numerical algorithms on classical computers.The method in this paper can also be naturally extended to solve more general Sylvester equations and discrete-time Lyapunov equations.Our study opens up a new dimension of research by using quantum algorithms to speed up the process for solving the control laws in classical control systems.

[1]P.H.Petkov,N.D.Christov,M.M.Konstantinov.Computational Methods for Linear Control Systems.Hertfordshire:Prentice-Hall,1991.

[2]B.N.Datta.Numerical Methods for Linear Control Systems.Boston:Elsevier Acadamic Press,2004.

[3]V.Sima.Algorithms for Linear-quadratic Optimization.New York:Marcel Dekker,1996.

[4]N.R.Sandell,P.Varaiya,M.Athans,et al.Survey of decentralized control methods for large scale systems.IEEE Transactions on Automatic Control,1978,23(2):108–128.

[5]P.Benner.Solving large-scale control problems.IEEE Control Systems,2004,24(1):44–59.

[6]P.Benner,J.R.Li,T.Penzl.Numerical solution of largescale Lyapunov equations,Riccati equations,and linear-quadratic optimal control problems.Numerical Linear Algebra with Applications,2008,15(9):755–777.

[7]M.A.Nielsen,I.L.Chuang.Quantum Computation and Quantum Information.Cambrige:Cambridge University Press,2000.

[8]R.P.Feynman.Simulating physics with computers.International Journal of Theoretical Physics,1982,21(6):467–488.

[9]P.W.Shor.Algorithms for quantum computation:discrete logarithms and factoring.Proceedings of the 35th Annual Symposium on Foundations of Computer Science.Piscataway:IEEE,1994:124–134.

[10]P.W.Shor.Polynomial-time algorithms for prime factorization and discrete logarithms on a quantum computer.SIAM Review,1999,41(2):303–332.

[11]L.K.Grover.A fast quantum mechanical algorithm for database search.Proceedings of the 28th ACM Symposium on Theory of Computing.Philadelphia:ACM,1996:212–219.

[12]S.Hallgren.Polynomial-time quantum algorithms for Pell’s equation and the principal ideal problem.Journal of the ACM,2007,54(1):1–19.

[13]S.Hallgren.Fast quantum algorithms for computing the unit group and class group of a number field.Proceedings of the 37th ACM Symposium on Theory of Computing,New York:ACM,2005:468–474.

[14]C.D¨urr,M.Heiligman,P.Hoyer,et al.Quantum query complexity of some graph problems.SIAM Journal on Computing,2004,35(6):1310–1328.

[15]A.Ambainis,A.M.Childs,Y.K.Liu.Quantum property testing for bounded-degree graphs.Proceedings of RANDOM 2011:Lecture Notes in Computer Science 6845.Princeton:Springer,2011:365–376.

[16]F.Magniez,M.Santha,M.Szegedy.Quantum algorithms for the triangle problem.SIAM Journal on Computing,2007,37(2):413–424.

[17]N.Wiebe,A.Kapoor,K.M.Svore.Quantum algorithms for nearest-neighbor methods for supervised and unsupervised learning.Quantum Information and Computation,2014,15(3):318–358.

[18]S.Lloyd,S.Garnerone,P.Zanardi.Quantum algorithms for topological and geometric analysis of data.Nature Communications,2016,7:DOI 10.1038/ncomms10138.

[19]S.Lloyd,M.Mohseni,P.Rebentrost.Quantum principal component analysis.Nature Physics,2014,10(9):631–633.

[20]D.S.Abrams,S.Lloyd.Simulation of many-body Fermi systems on a universal quantum computer.Physical Review Letters,1997,79(13):2586–2589.

[21]I.Kassal,S.P.Jordan,P.J.Love,et al.Quantum algorithms for the simulation of chemical dynamics.Proceedings of the National Academy of Sciences of the United States of America,2008,105(48):18681–18686.

[22]L.A.Wu,M.S.Byrd,D.A.Lidar.Polynomial-time simulation of pairing models on a quantum computer.Physical Review Letters,2002,89(5):DOI 10.1103/PhysRevLett.89.057904.

[23]A.W.Harrow,A.Hassidim,S.Lloyd.Quantum algorithm for solving linear systems of equations.Physical Review Letters,2009,15(103):DOI 10.1103/PhysRevLett.103.150502.

[24]D.W.Berry.High-order quantum algorithm for solving linear differential equations.Journal of Physics A:Mathematical and Theoretical,2014,47(10):DOI 10.1088/1751-8113/47/10/105301.

[25]S.K.Leyton,T.J.Osborne.A quantum algorithm to solve nonlinear differential equations.arXiv,2008:arXiv:0812.4423.

[26]Y.D.Cao,A.Papageorgiou,I.Petras,et al.Quantum algorithm and circuit design solving the Poisson equation.New Journal of Physics,2013,15(1):DOI 10.1088/1367-2630/15/1/013021.

[27]N.Wiebe,D.Braun,S.Lloyd.Quantum algorithm for data fitting.Physical Review Letters,2012,109(5):DOI 10.1103/PhysRevLett.109.05050.

[28]P.Rebentrost,M.Mohseni,S.Lloyd.Quantum support vector machine for big data classification.Physical Review Letters,2014,113(13):DOI 10.1103/PhysRevLett.113.130503.

[29]X.D.Cai,C.Weedbrook,Z.E.Su,et al.Experimental quantum computing to solve systems of linear equations.Physical Review Letters,2013,110(23):DOI 10.1103/PhysRevLett.110.230501.

[30]J.Pan,Y.D.Cao,X.W.Yao,et al.Experimental realization of quantum algorithm for solving linear systems of equations.Physical Review A,2014,89(2):DOI 10.1103/PhysRevA.89.022313.

[31]S.J.Hammarling.Numerical solution of the stable,non-negative definite Lyapunov equation.IMA Journal of Numerical Analysis,1982,2(3):303–323.

[32]L.Grover,T.Rudolph.Creating superpositions that correspond to efficiently integrable probability distributions.arXiv,2002:arXiv:quant-ph/0208112.

[33]A.Luis,J.Pe?rina.Optimum phase-shift estimation and the quantum description of the phase difference.Physical Review A,1996,54(5):4564–4570.

[34]D.W.Berry,G.Ahokas,R.Cleve,et al.Efficient quantum algorithms for simulating sparse Hamiltonians.Communications in Mathematical Physics,2007,270(2):359–371.

28 July 2017;revised 11 September 2017;accepted 11 September 2017

DOIhttps://doi.org/10.1007/s11768-017-7091-0

?Corresponding author.

E-mail:jing-zhang@mail.tsinghua.edu.cn.Tel:+86-13661167656.

This paper is dedicated to Professor T.J.Tarn on the occasion of his 80th birthday.

?2017 South China University of Technology,Academy of Mathematics and Systems Science,CAS,and Springer-Verlag GmbH Germany

He SUNreceived his B.Sc.and M.Sc.degree from School of Information Science and Technology,Nankai University,Tianjin,China,in 2013,and Department of Automation,Tsinghua University,Beijing,China,in 2017,respectively.E-mail:SunHe1992@yeah.net.

Jing ZHANGreceived his B.Sc.degree from Department of Mathematical Science and Ph.D.degree from Department of Automation,Tsinghua University,Beijing,China,in 2001 and 2006,respectively.From 2006 to 2008,he was a Postdoctoral Fellow at the Department of Computer Science and Technology,Tsinghua University,Beijing,China,and a Visiting Researcher from2008to2009 at the Advanced Science Institute,the Institute of Physical and Chemical Research(RIKEN),Japan.In 2010,he worked as a Visiting Assistant Professor at Department of Physics and National Center for Theoretical Sciences,National Cheng Kung University,Taiwan.He is now an Associate Professor at the Department of Automation,Tsinghua University,Beijing,China.His research interests include quantum control and micro/nano photonics.E-mail:jing-zhang@mail.tsinghua.edu.cn.