METHOD AND APPARATUS FOR SOLVING SYSTEM OF NONLINEAR EQUATIONS BASED ON QUANTUM CIRCUIT, AND STORAGE MEDIUM

Information

  • Patent Application
  • 20240354360
  • Publication Number
    20240354360
  • Date Filed
    November 25, 2022
    2 years ago
  • Date Published
    October 24, 2024
    4 months ago
Abstract
A method for solving a system of nonlinear equations on the basis of a quantum circuit includes acquiring a target system of nonlinear equations to be solved, converting the target system to obtain a target system of linear equations, constructing a quantum circuit corresponding to a quantum linear solver used for solving the target system, performing, based on the quantum circuit corresponding to the quantum linear solver, quantum state evolution and measurement on the target system, to solve the target system, and determining, based on an obtained solution of the target system, a solution of the target system to be solved. With the method, the complexity and difficulty in solving a system of nonlinear equations may be reduced, thereby filling in related technical gaps in the field of quantum computation.
Description
BACKGROUND
1. Technical Field

The present disclosure generally relates to the technical field of quantum computation. More specifically, the present disclosure relates to a method, an apparatus, a storage medium and an electronic device for solving a system of nonlinear equations based on a quantum circuit.


2. Background Art

The study of nonlinear equations for the purpose of application, or with the background of physics, mechanics or other subject issues, is not only one of the most important contents in traditional applied mathematics, but also an important part of contemporary mathematics, as well as an important bridge between mathematical theories and practical applications.


On the other hand, nonlinear problems are more common in nature, such as nonlinear finite element analysis, nonlinear dynamics, nonlinear programming, and the like. Therefore, it is important to construct a quantum algorithm for solving nonlinear problems. However, due to the linearity of quantum computation itself, there will be difficulties in constructing the quantum algorithm for solving nonlinear problems, and the study on quantum algorithms for solving a system of nonlinear equations is still relatively scarce.


At present, many significant problems in natural science and engineering technology can be attributed to the study of nonlinear equations. Mathematical models in many fields of real life can be described by nonlinear equations, but in general, it is difficult to obtain efficient solutions of systems of nonlinear equations easily with traditional numerical methods. Therefore, the study on how to accurately and quickly solve a system of nonlinear equations has shown very important theoretical and application values. The quantum computation is a novel computation method, the principle of which is to construct a computing framework based on the theory of quantum mechanics. When solving some problems, quantum computation has an exponential speed-up effect compared with the optimal classical algorithm.


Existing methods for solving a system of nonlinear equations typically consume a lot of computing resources, which may beyond the computing power of traditional computers, and involve a higher computation complexity, a longer time to obtain an accurate solution and more difficult computation. In this context, it is important to develop a more efficient algorithm for solving a system of nonlinear equations.


SUMMARY

In order to solve at least one or more technical problems described in the background section, the present disclosure proposes a method and an apparatus for solving a system of nonlinear ordinary differential equations on the basis of a quantum circuit to address deficiencies in the existing art, which can implement the technology of solving a system of nonlinear ordinary differential equations with a quantum algorithm, and reduce the complexity and difficulty in solving a system of nonlinear ordinary differential equations, thereby filling in related technical gaps in the field of quantum computation. In view of this, the present disclosure provides solutions in the following aspects.


In a first aspect, the present disclosure provides a method for solving a system of nonlinear ordinary differential equations on the basis of a quantum circuit, including: acquiring information about a system of nonlinear ordinary differential equations to be processed; converting the system of nonlinear ordinary differential equations to be processed to obtain a target system of linear ordinary differential equations; constructing a quantum circuit corresponding to a quantum linear solution algorithm, and performing quantum state evolution and measurement on the target system of linear ordinary differential equations, to acquire a solution of the target system of linear ordinary differential equations; calculating, according to the solution of the target system of linear ordinary differential equations, a solution of the system of nonlinear ordinary differential equations to be processed.


In one embodiment, the system of nonlinear ordinary differential equations to be processed is:








du
dt

=



F
1


u

+


F
2



u


2





,


u

(
0
)

=

u
in






where u is a function to be solved of the system of nonlinear ordinary differential equations to be processed, u∈custom-charactern, custom-character is real number space, and F1 and F2 are sparse matrices independent of time with a sparsity s, F1custom-charactern×n, and F2custom-charactern×n2.


In one embodiment, converting the system of nonlinear ordinary differential equations to be processed to obtain the target system of linear ordinary differential equations includes: converting the system of nonlinear ordinary differential equations to be processed into a preset type of system of nonlinear ordinary differential equations in a homotopy perturbation method; and converting the preset type of system of nonlinear ordinary differential equations into the target system of linear ordinary differential equations in a linear embedding method.


In one embodiment, the preset type of system of nonlinear ordinary differential equations is:










dv
0

dt

-


F
1



v
0



=
0

,



v
0

(
0
)

=

u
in












dv
1

dt

-


F
1



v
1


-


F
2




v
0



v
0




=
0

,



v
1

(
0
)

=
0











dv
2

dt

-


F
1



v
2


-


F
2

(



v
0



v
1


+


v
1



v
0



)


=
0

,



v
2

(
0
)

=
0
















dv
c

dt

-


F
1



v
c


-


F
2








j
=
0


c
-
1





v
j



v

c
-
1
-
j





=
0

,



v
c

(
0
)

=
0





where the c is the number of functions to be solved in the preset type of system of nonlinear ordinary differential equations to be processed, vi is a function to be solved of the preset type of system of nonlinear ordinary differential equations be processed, and 0≤i≤c.


In one embodiment, the target system of linear ordinary differential equations is:









d


y



dt

=

A


y




,



y


(
0
)

=

y
in






where custom-character=[{right arrow over (y)}0, {right arrow over (y)}1, {right arrow over (y)}2, . . . , {right arrow over (y)}c], {right arrow over (y)}j satisfies:








y


i

=

{





[


v
0

+

v
1

+

v
2

+


+

v
c


]

,

i
=
0








[



y



i
,
0


,


y



i
,
1


,


y



i
,
2


,


,


y



i
,


β
i

-
1




]

,

1

i

c










yin=[[uin], [uin2, 0, . . . , 0] [uin3, 0 . . . , 0] [uinc]], Bi represents the number of items in {right arrow over (y)}j, and







β
i

=

{





1
,

i
=
0














k
=
i

c



(



k




i



)


,

1

i

c





,


y



i
,
j








represents a jth item in {right arrow over (y)}i, expressed as {right arrow over (y)}i,j=⊗k=0ivai,j,k, and ai,j,k satisfies ai,j,k≥0, i+1≤Σk=0iai,j,k+1≤c+1.


In one embodiment, constructing the quantum circuit corresponding to the quantum linear solution algorithm includes:


constructing a quantum circuit that implements a first functional module V, where the first functional module is defined as










V




"\[LeftBracketingBar]"


0
m





:=


1

α








i




α
i






"\[LeftBracketingBar]"

i





,


α

=


4








j
=
0


j
0





(

-
1

)

j










i
=

j
+
1


b



(




2

b






b
+
j




)



2

2

b





,








j
0

=


blog
(

4

b
/
ϵ

)



,


b
=


κ
2



log

(

κ
/
ϵ

)



;





constructing a quantum circuit that implements a second functional module combination, where the quantum circuit for the second functional module combination includes T, W and T+, the T=Σj∈[N]jcustom-characterj|, W=S·(2T·T+−I4N2), the S is a unitary matrix of an exchange operation module, the T+ is a transposed conjugate of the























"\[RightBracketingBar]"


ψ
j




:=



"\[LeftBracketingBar]"

j






1

d










k




[
N
]

:

A
jk



0





(


A
jk
*







"\[RightBracketingBar]"




k



+


1
-



"\[LeftBracketingBar]"


A
jk



"\[RightBracketingBar]"








"\[RightBracketingBar]"




k

+
N



)

,




and I4N2 is a 4N2-dimensional identity matrix; constructing a quantum circuit that implements a third functional module V+, where the third functional module is a transposed conjugate form of the first functional module; and sequentially inserting the first functional module, the second functional module and the third functional module into the quantum circuit to form the quantum circuit corresponding to the quantum linear solution algorithm.


In a second aspect, an embodiment of the present application further provides an apparatus for solving a system of nonlinear ordinary differential equations on the basis of a quantum circuit, including: an acquisition module configured to acquire information about a system of nonlinear ordinary differential equations to be processed; a conversion module configured to convert the system of nonlinear ordinary differential equations to be processed to obtain a target system of linear ordinary differential equations; a construction module configured to construct a quantum circuit corresponding to a quantum linear solution algorithm, and perform quantum state evolution and measurement on the target system of linear ordinary differential equations, to acquire a solution of the target system of linear ordinary differential equations; a calculation module configured to calculate, according to the solution of the target system of linear ordinary differential equations, a solution of the nonlinear ordinary differential equations to be processed.


In one embodiment, the conversion module includes: a first conversion unit configured to convert the system of nonlinear ordinary differential equations to be processed into a preset type of system of nonlinear ordinary differential equations in a homotopy perturbation method; and a second conversion unit configured to convert the preset type of system of nonlinear ordinary differential equations into the target system of linear ordinary differential equations in a linear embedding method.


In one embodiment, the construction module includes: a first construction unit configured to construct a quantum circuit that implements a first functional module V, where the first functional module is defined as













V


"\[RightBracketingBar]"




0
m




:=


1

α








i



α







"\[RightBracketingBar]"




i



,

α
:=

4







j
=
0


j

0






(

-
1

)

j










i
=

j
+
1


b



(




2

b






b
+
j




)



2

2

b





,


j
0

=


b


log
(

4

b
/
ϵ

)




,


b
=


κ
2



log
(

κ
/
ϵ

)



;





a second construction unit configured to construct a quantum circuit that implements a second functional module combination, where the quantum circuit for the second functional module combination includes T, W and T+, the T=Σj∈[N]jcustom-characterj|, W=S·(2T·T+−I4N2), the S is a unitary matrix of an exchange operation module, the T+ is a transposed conjugate of the T, the





















|

Ψ
j




:=



"\[LeftBracketingBar]"

j






1

d










k




[
N
]

:

A
jk



0





(


A
jk
*







"\[RightBracketingBar]"




k



+


1
-



"\[LeftBracketingBar]"


A
jk



"\[RightBracketingBar]"








"\[RightBracketingBar]"




k

+
N



)

,




and I4N2 is a 4N2-dimensional identity matrix; a third construction unit configured to construct a quantum circuit that implements a third functional module V+, where the third functional module is a transposed conjugate form of the first functional module; and a combination unit configured to sequentially insert the first functional module, the second functional module and the third functional module into the quantum circuit to form the quantum circuit corresponding to the quantum linear solution algorithm.


In a third aspect, the present disclosure further provides a storage medium having a computer program stored thereon, where the computer program is configured to, when executed, cause any one of the methods described above to be implemented.


In a fourth aspect, the present disclosure further provides an electronic device, including a memory and a processor, where the memory has a computer program stored thereon, and the processor is configured to execute the computer program to implement any one of the methods described above.


According to the above solutions of the present disclosure, information about a system of nonlinear ordinary differential equations to be processed is firstly acquired, then the system of nonlinear ordinary differential equations to be processed is converted to obtain a target system of linear ordinary differential equations, a quantum circuit corresponding to a quantum linear solution algorithm is constructed, and quantum state evolution and measurement are performed on the target system of linear ordinary differential equations, to acquire a solution of the target system of linear ordinary differential equations, and according to the solution of the target system of linear ordinary differential equations, a solution of the nonlinear ordinary differential equations to be processed is calculated. Therefore, by means of related properties of quantum, the technology of solving a system of nonlinear ordinary differential equations with a quantum algorithm can be implemented, and the complexity and difficulty in solving the system of nonlinear ordinary differential equations can be reduced, thereby filling in related technical gaps in the field of quantum computation.


In order to solve at least one or more technical problems described in the background section, the present disclosure proposes a method and an apparatus for solving a system of quadratic nonlinear equations on the basis of a quantum circuit to address deficiencies in the existing art, which can implement the technology of solving a system of quadratic nonlinear equations with a quantum algorithm, and reduce the complexity and difficulty in solving a system of quadratic nonlinear equations, thereby filling in related technical gaps in the field of quantum computation. In view of this, the present disclosure provides solutions in the following aspects.


In a first aspect, the present disclosure provides a method for solving a system of quadratic nonlinear equations on the basis of a quantum circuit, including:


acquiring a target system of linear equations, where the target system of linear equations is converted and determined from an initial system of quadratic nonlinear equations; constructing a quantum circuit corresponding to a quantum linear solver, and running the quantum circuit and measuring, to solve the target system of linear equations; determining, based on an obtained solution of the target system of linear equations, a solution of the initial system of quadratic nonlinear equations.


In one embodiment, the initial system of quadratic nonlinear equations is specifically:








F
0

+


F
1


x

+


F
2



x




2





=
0




where x∈Rn, R represents real number space, Fi∈Rn×Rni, and F1 and F2 each have a sparsity s.


In one embodiment, acquiring the target system of linear equations includes: converting the initial system of quadratic nonlinear equations into a preset system of pseudo-linear equations in a homotopy perturbation method; converting the preset system of pseudo-linear equations into a target system of linear equations in a linear embedding method.


In one embodiment, the preset system of pseudo-linear equations is:








v
0

=


-

F
1

-
1





F
0







v
1

=


-

F
1

-
1






F
2

(


v
0



v
0


)







v
2

=


-

F
1

-
1






F
2

(



v
0



v
1


+


v
1



v
0



)










v
c

=


-

F
1

-
1





F
2








j
=
0


j
=

c
-
1






v
j



v

c
-
1
-
j









where F1 is invertible, vi is a variable to be solved in the preset system of pseudo-linear equations, 0≤i≤c, and the c is the number of variables to be solved in the preset system of pseudo-linear equations.


In one embodiment, the target system of linear equations is:








[




A

0
,
0





A

0
,
1




0


0


0


0




0



A

1
,
1





A

1
,
2




0


0


0




0


0



A

2
,
2







0


0




0


0


0






A


c
-
1

,

c
-
2





0




0


0


0


0



A


c
-
1

,

c
-
1






A


c
-
1

,
c






0


0


0


0


0



A

c
,
c





]

[




y
0






y
1






y
2











y

c
-
1







y
c




]

=

[




b
0






b
1






b
2











b

c
-
1







b
c




]





where Aj,i is a







n

i
+
1









j
=
i

c



(



j




i



)

×

n

i
+
1









j
=
i

c



(



j




i



)





dimension matrix, Ai,i+1 is a ni+1βi×ni+2βi+1 dimension matrix, y=y0, y1, . . . , yc, yi satisfies y0=v0+v1+v2+ . . . +vc, yi=[yi,0, yi,1, yi,2, . . . , yi,βi−1], and βi represents the number of items in yi.


In one embodiment, constructing the quantum circuit corresponding to the quantum linear solver includes: constructing a quantum circuit that implements a first functional module V, where the first functional module is defined as













V


"\[RightBracketingBar]"




0
m




:=


1

α








i



α







"\[RightBracketingBar]"




i



,

α
:=

4







j
=
0


j

0






(

-
1

)

j










i
=

j
+
1


b



(




2

b






b
+
i




)



2

2

b





,


j
0

=


b


log
(

4

b
/
ϵ

)




,


b
=


κ
2



log
(

κ
/
ϵ

)



;





constructing a quantum circuit that implements a second functional module combination, where the quantum circuit for the second functional module combination includes T, W and T+, the T=Σj∈[N]jcustom-characterj|, W=S·(2T·T+−I4N2), the S is a unitary matrix of an exchange operation module, the T+ is a transposed conjugate of the T, the























"\[RightBracketingBar]"


ψ
j




:=



"\[LeftBracketingBar]"

j






1

d










k




[
N
]

:

A
jk



0





(


A
jk
*







"\[RightBracketingBar]"




k



+


1
-



"\[LeftBracketingBar]"


A
jk



"\[RightBracketingBar]"








"\[RightBracketingBar]"




k

+
N



)

,




and I4N2 is a 4N2-dimensional identity matrix; constructing a quantum circuit that implements a third functional module V+, where the third functional module is a transposed conjugate form of the first functional module; and sequentially inserting the first functional module, the second functional module and the third functional module into the quantum circuit to form the quantum circuit corresponding to the quantum linear solver.


In a second aspect, the present disclosure provides an apparatus for solving a system of quadratic nonlinear equations on the basis of a quantum circuit, including: an acquisition module configured to acquire a target system of linear equations, where the target system of linear equations is converted and determined from an initial system of quadratic nonlinear equations; a construction module configured to construct a quantum circuit corresponding to a quantum linear solver, and run the quantum circuit and measure, to solve the target system of linear equations; and a determination module configured to determine, based on an obtained solution of the target system of linear equations, a solution of the initial system of quadratic nonlinear equations.


In one embodiment, the acquisition module includes: a first conversion unit configured to convert the initial system of quadratic nonlinear equations into a preset system of pseudo-linear equations in a homotopy perturbation method; and a second conversion unit configured to convert the preset system of pseudo-linear equations into a target system of linear equations in a linear embedding method.


In one embodiment, the construction module includes: a first construction unit configured to construct a quantum circuit that implements a first functional module V, where the first functional module is defined as













V


"\[RightBracketingBar]"




0
m




:=


1

α








i



α







"\[RightBracketingBar]"




i



,

α
:=

4







j
=
0


j

0






(

-
1

)

j










i
=

j
+
1


b



(




2

b






b
+
i




)



2

2

b





,




j0=√{square root over (blog(4b/ϵ))}, b=κ2 log (κ/ϵ); a second construction unit configured to construct a quantum circuit that implements a second functional module combination, where the quantum circuit for the second functional module combination includes T, W and T+, the T=Σj∈[N]jcustom-characterj|, W=S·(2T·T+−I4N2), the S is a unitary matrix of an exchange operation module, the T+ is a transposed conjugate of the T, the























"\[RightBracketingBar]"


ψ
j




:=



"\[LeftBracketingBar]"

j






1

d










k




[
N
]

:

A
jk



0





(


A
jk
*







"\[RightBracketingBar]"




k



+


1
-



"\[LeftBracketingBar]"


A
jk



"\[RightBracketingBar]"








"\[RightBracketingBar]"




k

+
N



)

,




and I4N2 is a 4N2-dimensional identity matrix; a third construction unit configured to construct a quantum circuit that implements a third functional module V+, where the third functional module is a transposed conjugate form of the first functional module; and a combination unit configured to sequentially insert the first functional module, the second functional module and the third functional module into the quantum circuit to form the quantum circuit corresponding to the quantum linear solver.


In a third aspect, the present disclosure provides a storage medium having a computer program stored thereon, where the computer program is configured to, when executed, cause any one of the methods described above to be implemented.


In a fourth aspect, the present disclosure provides an electronic device, including a memory and a processor, where the memory has a computer program stored thereon, and the processor is configured to execute the computer program to implement any one of the methods described above.


According to the above solutions of the present disclosure, a target system of linear equations is firstly converted and determined according to an initial system of quadratic nonlinear equations, then a quantum circuit corresponding to a quantum linear solver is constructed, and the quantum circuit is run and measurement is performed to solve the target system of linear equations, and the solution of the initial system of quadratic nonlinear equations is determined based on the solved solution of the target system of linear equations. Therefore, by means of related properties of quantum, the technology of solving a system of quadratic nonlinear equations with a quantum algorithm can be implemented, and the complexity and difficulty in solving the system of quadratic nonlinear equations can be reduced, thereby filling in related technical gaps in the field of quantum computation.


0 In order to solve at least one or more technical problems described in the background section, the present disclosure proposes a method, an apparatus, a medium and an electronic device for solving a system of linear equations, which aim to reduce the complexity in solving linear system problems, and achieve an acceleration effect in solving linear system problems with quantum. In view of this, the present disclosure provides solutions in the following aspects.


In a first aspect, the present disclosure provides a method for solving a system of linear equations, including: acquiring a matrix A and a vector b in the system of linear equations Ax=b, where the x is an unknown; when a condition number κA of the matrix A is greater than a preset threshold, processing the matrix A and the vector b through a polynomial preprocessor to obtain a matrix A′ and a vector b′, where a condition number κA′ of the matrix A′ is less than the condition number KA of the matrix A; and solving the unknown x through a system of linear equations A′x=b′ constructed from the matrix A′ and the vector b′.


In one embodiment, processing the matrix A and the vector b through the polynomial preprocessor to obtain the matrix A′ and the vector b′ includes: preparing a quantum state |b) of the vector b, and determining a polynomial P(A) based on a polynomial function P(y) and the matrix A, where the y is an independent variable; constructing an operator P based on the polynomial P(A), and multiplying the polynomial P(A) by the matrix A to obtain a matrix A′; and applying the operator P to the quantum state |b) to obtain a quantum state |b′), where the quantum state |b′) is a quantum state corresponding to a vector b′.


In one embodiment, before determining the polynomial P(A) based on the polynomial function P(y) and the matrix A, the method further includes: acquiring an approximate function Km(y) containing parameters, and determining a domain of definition of the approximate function Km(y) containing parameters;


selecting m+2 approximate deviation points from the domain of definition, and substituting the m+2 approximate deviation points into a relational expression composed of the approximate function Km(y) containing parameters and a deviation amplitude E, to obtain a m+2-dimensional system of linear equations Km(yk)−T=(−1)k+1E, where k=0,1,2 . . . m+1, the m is an integer greater than 0, and the Tis any natural number in the domain of definition; solving the system of linear equations Km(yk)−T=(−1)k+1E, to obtain a value of a parameter in the approximate function Km(y) containing parameters; determining a target approximate function based on the value of the parameter, and determining the polynomial function P(y) based on the target approximate function.


In one embodiment, determining the target approximation function based on the value of the parameter includes: substituting the value of the parameter into the approximate function Km(y) containing parameters to obtain an initial approximate function; determining an extreme point of an absolute value of a difference between the initial approximation function and the T; if the extreme point and the m+2 approximation deviation points are equal within an accuracy requirement, determining the initial approximation function as the target approximation function; if the extreme point and the m+2 approximation deviation points are not equal within the accuracy requirement, taking the extreme point as a new point of the m+2 approximate deviation points, and performing the step of substituting the m+2 approximate deviation points into the relational expression composed of the approximate function Km(y) containing parameters and the deviation amplitude E to obtain the m+2-dimensional system of linear equations Km(yk)−T=(−1)k+1E.


In one embodiment, determining the polynomial function P(y) based on the target approximate function includes: determining, based on a relational expression P(y)=K(y)/f(y) of the target approximation function and the polynomial function P(y), the polynomial function P(y), where the f(y)=y.


In one embodiment, if the approximate function Km(y) containing parameters is an even function, the Km(y)=Σi=0mθ2iy2i+1; and if the approximate function Km(y) containing parameters is an odd function, the Km(y)==Σi=0mθ2iy2i+1; where the θ2i, θ2i+1 are parameters, and the i is an integer greater than or equal to 0.


In one embodiment, constructing the operator P based on the polynomial P(A) includes: preparing the operator P corresponding to the polynomial P(A) through quantum signal processing (QSP).


In one embodiment, solving the unknown x through the system of linear equations A′x=b′ constructed from the matrix A′ and the vector b′ includes: determining an operator U for a Harrow-Hassidim-Lloyd algorithm (HHL) algorithm based on the matrix A′; inputting the quantum state |b′) corresponding to the vector b′ and the operator U to the HHL algorithm, to determine a target quantum state including values of the unknown x; and determining a solution result of the unknown x based on the target quantum state.


In one embodiment, determining the operator U for the HHL algorithm based on the matrix A′ includes:


determining the operator U=eiA′t for the HHL algorithm; expanding the operator U=eiA′t according to the Jacobi-Anger expansion to determine an exponential expansion of the operator U, where the exponential expansion of the operator U is:








e

i

γ

t


=



J
0

(
t
)

+

2







k
=
1






(

-
1

)

k




J

2

k


(
t
)




T

2

k


(
γ
)


+

2

i







k
=
0






(

-
1

)

k




J


2

k

+
1


(
t
)




T

2

k


(
γ
)




,




where Jk is a kth order Bessel function of the first kind, and Tk is a kth order Chebyshev polynomial of the first kind; saving the exponential expansion of the operator U in accuracy of a qth order, where the exponential expansion of the operator U in accuracy of the qth order is:









e

i

γ

t


=



f
1

(

γ
,
t

)

+


f
2

(

γ
,
t

)



,

where







f
1

(

γ
,
t

)

=



J
0

(
t
)

+

2







k
=
1

q




(

-
1

)

k




J

2

k


(
t
)




T

2

k


(
γ
)




,




f
2

(

γ
,
t

)

=

2

i







k
=
1

q




(

-
1

)

k




J


2

k

+
1


(
t
)




T

2

k


(
γ
)



,
and





γ
=


P

(
A
)


A


;





preparing an operator f1 corresponding to f1(γ,t) and an operator f2 corresponding to f2(γ, t) through the QSP; and


constructing the operator U through linear combination of unitary operators (LCU), the operator f1 and the operator f2.


In a second aspect, an embodiment of the present application further provides an apparatus for solving a system of linear equations, including: an acquisition unit configured to acquire a matrix A and a vector b in the system of linear equations Ax=b, where the x is an unknown; a processing unit configured to, when a condition number KA of the matrix A is greater than a preset threshold, process the matrix A and the vector b through a polynomial preprocessor to obtain a matrix A′ and a vector b′, where a condition number KA′ of the matrix A′ is less than the condition number KA of the matrix A; and a calculation unit configured to solve the unknown x through a new system of linear equations A′x=b′ constructed from the matrix A′ and the vector b′.


In one embodiment, in the aspect of processing the matrix A and the vector b through the polynomial preprocessor to obtain the matrix A′ and the vector b′, the processing unit is specifically configured to: prepare a quantum state |b) of the vector b, and determine a polynomial P(A) based on a polynomial function P(y) and the matrix A, where the y is an independent variable; construct an operator P based on the polynomial P(A), and multiply the polynomial P(A) by the matrix A to obtain a matrix A′; and apply the operator P to the quantum state |bcustom-character to obtain a quantum state |b′custom-character, where the quantum state |b′) is a quantum state corresponding to a vector b′.


In one embodiment, before determining the polynomial P(A) based on the polynomial function P(y) and the matrix A, the processing unit is further configured to: acquire an approximate function Km(y) containing parameters, and determine a domain of definition of the approximate function Km(y) containing parameters; select m+2 approximate deviation points from the domain of definition, and substitute the m+2 approximate deviation points into a relational expression composed of the approximate function Km(y) containing parameters and a deviation amplitude E, to obtain a m+2-dimensional system of linear equations Km(yk)−T=(−1)k+1E, where k=0,1,2 . . . m+1, the m is an integer greater than 0, and the T is any natural number in the domain of definition; solve the system of linear equations Km(yk)−T=(−1)k+1E, to obtain a value of a parameter in the approximate function Km(y) containing parameters; determine a target approximate function based on the value of the parameter, and determine the polynomial function P(y) based on the target approximate function.


In one embodiment, in the aspect of determining the target approximate function based on the value of the parameter, the processing unit is specifically configured to: substitute the value of the parameter into the approximate function Km(y) containing parameters to obtain an initial approximate function; determine an extreme point of an absolute value of a difference between the initial approximation function and the T; if the extreme point and the m+2 approximation deviation points are equal within an accuracy requirement, determine the initial approximate function as the target approximate function; if the extreme point and the m+2 approximation deviation points are not equal within the accuracy requirement, take the extreme point as a new point of the m+2 approximate deviation points, and perform the step of substituting the m+2 approximate deviation points into the relational expression composed of the approximate function Km(y) containing parameters and the deviation amplitude E to obtain the m+2-dimensional system of linear equations Km(yk)−T=(−1)k+1E.


In one embodiment, in the aspect of determining the polynomial function P(y) based on the target approximate function, the processing unit is specifically configured to: determine, based on a relational expression P(y)=K(y)/f(y) of the target approximation function and the polynomial function P(y), the polynomial function P(y), where the f(y)=y.


In one embodiment, if the approximate function Km(y) containing parameters is an even function, the Km(y)=Σi=0mθ2iy2i+1; and if the approximate function Km(y) containing parameters is an odd function, the Km(y)=Σi=0mθ2iy2i+1; and the θ2i, θ2i+1 are parameters, and the i is an integer greater than or equal to 0.


In one embodiment, in the aspect of constructing the operator P based on the polynomial P(A), the processing unit is specifically configured to: prepare the operator P corresponding to the polynomial P(A) through quantum signal processing (QSP).


In one embodiment, in the aspect of solving the unknown x through the system of linear equations A′x=b′ constructed from the matrix A′ and the vector b′, the calculation unit is specifically configured to: determine an operator U for an HHL algorithm based on the matrix A′;


input the quantum state |b′) corresponding to the vector b′ and the operator U to the HHL algorithm, to determine a target quantum state including a value of the unknown x; and determine a solution result of the unknown x based on the target quantum state.


In one embodiment, in the aspect of determining the operator U for the HHL algorithm based on the matrix A′, the calculation unit is specifically configured to: determine the operator U=eiA′t for the HHL algorithm; expand the operator U=eiA′t according to the Jacobi-Anger expansion to determine an exponential expansion of the operator U, where the exponential expansion of the operator U is:








e

i

γ

t


=



J
0

(
t
)

+

2







k
=
1






(

-
1

)

k




J

2

k


(
t
)




T

2

k


(
γ
)


+

2

i







k
=
0






(

-
1

)

k




J


2

k

+
1


(
t
)




T

2

k


(
γ
)




,




where Jk is a kth order Bessel function of the first kind, and Tk is a kth order Chebyshev polynomial of the first kind; save the exponential expansion of the operator U in precision of a qth order, where the exponential expansion of the operator U in accuracy of the qth order is:









e

i

γ

t


=



f
1

(

γ
,
t

)

+


f
2

(

γ
,
t

)



,

where







f
1

(

γ
,
t

)

=



J
0

(
t
)

+

2







k
=
1

q




(

-
1

)

k




J

2

k


(
t
)




T

2

k


(
γ
)




,




f
2

(

γ
,
t

)

=

2

i







k
=
0

q




(

-
1

)

k




J


2

k

+
1


(
t
)




T

2

k


(
γ
)



,
and





γ
=


P

(
A
)


A


;





prepare an operator f1 corresponding to f1(γ,t) and an operator f2 corresponding to f2(γ, t) through the QSP; and construct the operator U through linear combination of unitary operators (LCU), the operator f1 and the operator f2.


In a third aspect, the present disclosure further provides a storage medium having a computer program stored thereon, where the computer program is configured to, when executed, cause any one of the methods described above to be implemented.


In a fourth aspect, the present disclosure further provides an electronic device, including a memory and a processor, where the memory has a computer program stored thereon, and the processor is configured to execute the computer program to implement any one of the methods described above.


According to the above solutions of the present disclosure, a matrix A and a vector b in an original system of linear equations Ax=b are processed through a polynomial preprocessor to obtain a new system of linear equations A′x=b′, and a condition number κA′ of the matrix A′ in the new system of linear equations is less than the condition number κA of the matrix A in the original system of linear equations, so that the condition number of the input matrix, and thus the complexity of the linear system problem, are reduced. Then, the new system of linear equations is solved to obtain a common unknown x, so that an acceleration effect in solving linear system problems with quantum is achieved.


In order to solve at least one or more technical problems described in the background section, the present disclosure proposes a quantum linear solution method based on a polynomial preprocessor, which can reduce the time, complexity and computation amount in solving linear problems while reducing occupation of hardware resources. In view of this, the present disclosure provides solutions in the following aspects.


In a first aspect, the present disclosure provides a quantum linear solution method based on a polynomial preprocessor, including: acquiring information of a first matrix A and a first vector b in a linear system to be processed; calculating a polynomial p(A) of the first matrix A used for preprocessing the linear system; pre-processing the linear system according to the polynomial p(A) to obtain a second matrix A′ and a second vector b′; and constructing a quantum circuit corresponding to an HHL algorithm, and performing, according to the second matrix A′ and the second vector b′, quantum state evolution and measurement operations to obtain a final quantum state of the evolved quantum circuit.


In one embodiment, obtaining the second matrix A′ and the second vector b′ includes: obtaining a second matrix A′=p(A) A and a second vector b′=p(A)b.


In one embodiment, performing, according to the second matrix A′ and the second vector b′, quantum state evolution and measurement operations to obtain the final quantum state of the evolved quantum circuit includes: inputting, for the quantum circuit, an initial quantum state |b′) and an operator U, and performing quantum state evolution and measurement operations to obtain the final quantum state of the evolved quantum circuit, where the operator U=eiA′t, and the t is a constant.


In one embodiment, calculating the polynomial p(A) of the first matrix A used for preprocessing the linear system includes: constructing a function ƒm(x) and a preset interval S for solving the polynomial, where the function ƒm(x) satisfies fm(x)=pm(x)x, and fm(xi)+(−1)iE=1, the m is a polynomial order, the i=0,1,2, . . . , m+1, the S∈[1/κA, 1], the E is a deviation amplitude, and the κA is a condition number; selecting m+1 points in the preset interval, and acquiring a solution of a system of linear equations fm(xi)+(−1)iE=1, where the m+1 points are respectively x1, x2, . . . , xm+1, and satisfy x1=1/κ, xm+1=1; substituting the acquired solution of the system of linear equations into fm(x) to obtain a set N′ of points corresponding to local maximum values of |1−fm(x)|; and judging whether an absolute value of fm(x)−1 is the same for any x∈N′, and whether a sign of fm(x)−1 changes between plug and minus alternately as x increases, if the absolute value of fm(x)−1 is the same for any x∈N′, and the sign of fm(x)−1 changes between plug and minus alternately as x increases, determining a current fm(x) as an optimal polynomial.


In one embodiment, the method further includes: if the absolute value of fm(x)−1 is not the same for any x∈N′, or the sign of fm(x)−1 does not change between plug and minus alternately as x increases, replacing elements in x1, x2, . . . , xm+1 with elements in the set N′ according to a rule including: replacing x2, . . . , xm with elements in the set N′ while x1, xm+1 remain unchanged, and then returning to perform the step of acquiring a solution of the system of linear equations fm(xi)+(−1)iE=1, until an optimal fm(x) satisfying the judgment condition is obtained.


In one embodiment, constructing the quantum circuit corresponding to the HHL algorithm includes: obtaining several qubits including an auxiliary qubit, a first qubit, and a second qubit, where initial states of the auxiliary qubit and the first qubit are set to |0custom-character, an initial state of the second qubit is set to |b′custom-characteri=1Nb′i|icustom-character, the bi′ is an ith element of the second vector b′, and the N is a dimension number of the second vector; determining a unitary matrix U corresponding to the second matrix A′; constructing a first sub-quantum circuit module for phase estimation, which is used to decompose the |b′custom-character into |b′custom-characterj=1Nβjjcustom-character in eigenspace of the second matrix A′, and transforming the initial state |0custom-character|b′custom-character of the first qubit and the second qubit into Σj=1Nβjjcustom-characterjcustom-character, where the |μjcustom-character is an eigenvector of the second matrix A′, the λj is an eigenvalue of the second matrix A′, and the βj is an eigenvector amplitude of the second matrix A′; constructing a second sub-quantum circuit module for performing a controlled rotation operation, which is used to rotate the auxiliary qubit with |λjcustom-character as a control bit to obtain






















j
=
1

N



β
j


|

λ
j








"\[LeftBracketingBar]"


μ
j







(



1
-


c
2


λ
j
2








"\[LeftBracketingBar]"

0






+


c

λ
j






"\[LeftBracketingBar]"

1





)

,




where the C is a normalization constant; constructing a third sub-quantum circuit module for performing inverse phase estimation, which is used to reset |λjcustom-character to |0custom-character; constructing a measurement operation module for the auxiliary qubit, so that when the measured quantum state of the auxiliary qubit is (1custom-character:












"\[LeftBracketingBar]"


x





=


1








j
=
1

N




(

C


λ
j

-
1




β
j


)

2











j
=
1

N


C


λ
j

-
1




β
j






"\[LeftBracketingBar]"



μ
j






,




and the |x′custom-character∝|xcustom-character; and sequentially forming the first sub-quantum circuit module, the second sub-quantum circuit module, the third sub-quantum circuit module, and the quantum measurement operation module into the quantum circuit corresponding to the HHL algorithm.


In a second aspect, an embodiment of the present application further provides a quantum linear solution apparatus based on a polynomial preprocessor, including: an acquisition module configured to acquire information of a first matrix A and a first vector b in a linear system to be processed; a calculation module configured to calculate a polynomial p(A) of the first matrix A used for preprocessing the linear system; an obtaining module configured to pre-process the linear system according to the polynomial p(A), to obtain a second matrix A′ and a second vector b′; and a construction module configured to construct a quantum circuit corresponding to an HHL algorithm, and perform, according to the second matrix A′ and the second vector b′, quantum state evolution and measurement operations to obtain a final quantum state of the evolved quantum circuit.


In one embodiment, the obtaining module includes: a first obtaining unit configured to obtain a second matrix A′=p(A) A and a second vector b′=p(A)b.


In one embodiment, the construction module includes: an input unit configured to input for the quantum circuit, an initial quantum state |b′) and an operator U, and perform quantum state evolution and measurement operations to obtain the final quantum state of the evolved quantum circuit, where the operator U=eiA′t, and the t is a constant.


In one embodiment, the calculation module includes: a first construction unit configured to construct a function ƒm(x) and a preset interval S for solving the polynomial, where the function ƒm(x) satisfies fm(x)=pm (x)x, and fm(xi)+(−1)iE=1, the m is a polynomial order, the i=0,1,2, . . . , m+1, the S∈[1/κA, 1], the E is a deviation amplitude, and the κA is a condition number; a selection unit configured to select m+1 points in the preset interval, and acquire a solution of a system of linear equations fm(xi)+(−1)iE=1, where the m+1 points are respectively x1, x2, . . . , xm+1, and satisfy x1=1/κ, xm+1=1; a substitution unit configured to substitute the acquired solution of the system of linear equations into fm(x) to obtain a set N′ of points corresponding to local maximum values of |1−fm(x)|; and a judgment unit configured to judge whether an absolute value of fm(x)−1 is the same for any x∈N′, and whether a sign of fm(x)−1 changes between plug and minus alternately as x increases, if the absolute value of fm(x)−1 is the same for any x∈N′, and the sign of fm(x)−1 changes between plug and minus alternately as x increases, determine a current fm(x) as an optimal polynomial.


In one embodiment, the calculation module further includes: a return unit configured to, if the absolute value of fm(x)−1 is not the same for any x∈N′, or the sign of fm(x)−1 does not change between plug and minus alternately as x increases, replace elements in x1, x2, . . . , xm+1 with elements in the set N′ according to a rule including: replace x2, . . . , xm with elements in the set N′ while x1, xm+1 remain unchanged, and then return to perform the step of acquiring a solution of the system of linear equations fm(xi)+(−1)iE=1, until an optimal fm(x) meeting a judgment condition is obtained.


In one embodiment, the construction module includes: a second obtaining unit configured to obtain several qubits including an auxiliary qubit, a first qubit, and a second qubit, where initial states of the auxiliary qubit and the first qubit are set to |0), an initial state of the second qubit is set to |b′custom-characteri=1Nb′i|icustom-character, the bj′ is an ith element of the second vector b′, and the N is a dimension number of the second vector; a determination unit configured to determine a unitary matrix U corresponding to the second matrix A′; a second construction unit configured to construct a first sub-quantum circuit module for phase estimation, which is used to decompose the |b′custom-character into |b′custom-characterj=1Nβjjcustom-character in eigenspace of the second matrix A′, and transform the initial state |0custom-character|b′custom-character of the first qubit and the second qubit into Σj=1Nβjjjcustom-character, where the |μjΣj=1Nβjjjcustom-character is an eigenvector of the second matrix A′, the λj is an eigenvalue of the second matrix A′, and the βj is an eigenvector amplitude of the second matrix A′; a third construction unit configured to construct a second sub-quantum circuit module for performing a controlled rotation operation, which is used to rotate the auxiliary qubit with |λjcustom-character as a control bit to obtain






















j
=
1

N



β
j


|

λ
j








"\[LeftBracketingBar]"


μ
j







(



1
-


c
2


λ
j
2








"\[LeftBracketingBar]"

0






+


c

λ
j






"\[LeftBracketingBar]"

1





)

,




where the C is a normalization constant; a fourth construction unit configured to construct a third sub-quantum circuit module for performing inverse phase estimation, which is used to reset |λjcustom-character to |0custom-character; a fifth construction unit configured to construct a measurement operation module for the auxiliary qubit, so












"\[LeftBracketingBar]"


x





=


1








j
=
1

N




(

C


λ
j

-
1




β
j


)

2











j
=
1

N


C


λ
j

-
1




β
j






"\[LeftBracketingBar]"



μ
j






,




and the |x′custom-character∝|xcustom-character; and a combination unit configured to sequentially form the first sub-quantum circuit module, the second sub-quantum circuit module, the third sub-quantum circuit module, and the quantum measurement operation module into the quantum circuit corresponding to the HHL algorithm.


In a third aspect, the present disclosure further provides a storage medium having a computer program stored thereon, where the computer program is configured to, when executed, cause any one of the methods described above to be implemented.


In a fourth aspect, the present disclosure further provides an electronic device, including a memory and a processor, where the memory has a computer program stored thereon, and the processor is configured to execute the computer program to implement any one of the methods described above.


According to the above solutions of the present disclosure, information of a first matrix A and a first vector b in a linear system to be processed is firstly acquired, a polynomial p(A) of the first matrix A used for preprocessing the linear system is calculated, the linear system is pre-processed according to the polynomial p(A), to obtain a second matrix A′ and a second vector b′, a quantum circuit corresponding to an HHL algorithm is constructed, and according to the second matrix A′ and the second vector b′, quantum state evolution and measurement operations are performed to obtain a final quantum state of the evolved quantum circuit, which can reduce the time, complexity and computation amount in solving linear problems and speed up solution of the quantum linear algorithm, while reducing occupation of hardware resources.





BRIEF DESCRIPTION OF THE DRAWINGS

The description and other objects, features and advantages of the exemplary implementations of the present disclosure will become apparent from the following detailed description with reference to the accompanying drawings. In the accompanying drawings, several implementations of the present disclosure are illustrated by way of example but not limitation, and like or corresponding reference numerals indicate like or corresponding parts, in which:



FIG. 1 is a block diagram showing a hardware structure of a computer terminal for solving a system of nonlinear equations on the basis of a quantum circuit according to an embodiment of the present disclosure;



FIG. 2 is a schematic flowchart of a method for solving a system of nonlinear equations on the basis of a quantum circuit according to an embodiment of the present disclosure;



FIG. 3 is a schematic diagram of a quantum circuit for solving a system of nonlinear equations according to an embodiment of the present disclosure:



FIG. 4 is a schematic diagram of a quantum circuit corresponding to a quantum linear solver according to an embodiment of the present disclosure:



FIG. 5 is a schematic diagram of a quantum circuit for a walk operator W according to an embodiment of the present disclosure;



FIG. 6 is a quantum circuit diagram for constructing an operator U according to an embodiment of the present disclosure;



FIG. 7 is a schematic diagram of a first sub-quantum circuit module corresponding to phase estimation according to an embodiment:



FIG. 8 is a total quantum circuit diagram corresponding to an HHL algorithm according to an embodiment; and



FIG. 9 is a schematic structural diagram of an apparatus for solving a system of nonlinear equations on the basis of a quantum circuit according to an embodiment.





DETAIL DESCRIPTION OF THE INVENTION

The technical solutions provided in the embodiments of the present disclosure will be described clearly and completely below with reference to the accompanying drawings in the embodiments of the present disclosure. Apparently, the described embodiments are some, but not all, of the embodiments of the present disclosure. Based on the embodiments of the present disclosure, all the other embodiments obtained by those of ordinary skill in the art without any creative labor fall into the protection scope of the present disclosure.


It should be understood that the terms “first,” “second,” “third,” and “fourth,” and the like in the claims, description, and drawings of the present disclosure are used to distinguish between different objects, and are not used to describe a particular order. The terms “comprise” and “include,” when used in the description and claims of the present disclosure, specify the presence of stated features, integers, steps, operations, elements, and/or components, but do not preclude the presence or addition of one or more other features, integers, steps, operations, elements, components, and/or groups thereof.


It is also to be understood that the terminology used in the description of the present disclosure here is for the purpose of describing particular embodiments only, and is not intended to be limiting of the present disclosure. As used in the specification and claims of the disclosure, the singular forms “a”, “an” and “the” are intended to include the plural forms as well, unless the context clearly indicates otherwise. It should be further understood that the term “and/or” as used in the description and claims of the disclosure refers to any and all possible combinations of one or more of the associated listed items and includes such combinations.


An embodiment of the present disclosure provides a method for solving a system of nonlinear equations on the basis of a quantum circuit, which may be applied to an electronic device such as computer terminals, specifically general computers, quantum computers, or the like. The following detailed description is made by taking the case of running on a computer terminal as an example.



FIG. 1 is a block diagram showing a hardware structure of a computer terminal for solving a system of nonlinear equations on the basis of a quantum circuit according to an embodiment of the present disclosure. As shown in FIG. 1, the computer terminal may include one or more (only one is shown in FIG. 1) processors 102 (which may include, but are not limited to, processing devices such as microprocessors (MCUs) or programmable logic devices (FPGAs)), and a memory 104 used to store data. Optionally, the computer terminal may further include a transmission device 106 for communication functions and an input/output device 108. It will be understood by those ordinary skilled in the art that the structure shown in FIG. 1 is merely illustrative, and does not form any limitation to the structure of the computer terminal. For example, the computer terminal may include more or fewer components than those shown in FIG. 1, or have a different configuration than that shown in FIG. 1.


The memory 104 may be configured to store software programs and modules of application software, such as the program instructions/modules corresponding to the method for solving a system of nonlinear equations on the basis of a quantum circuit provided in the embodiments the present disclosure. The processor 102, by executing the software programs and modules stored on the memory 104, performs various functional applications and data processing, that is, implement the above method. The memory 104 may include a high-speed random access memory, or include a non-volatile memory such as one or more magnetic storage devices, flash memories, or other non-volatile solid state memories. In some examples, the memory 104 may further include a memory remotely disposed relative to the processor 102, which may be connected to the computer terminal via a network. Examples of such networks include, but are not limited to, the Internet, intranets, local area networks, mobile communication networks, and combinations thereof.


The transmission device 106 is configured to receive or transmit data via a network. Specific examples of such networks may include a wireless network provided by a communication provider of the computer terminal. In one example, the transmission device 106 includes a network interface controller (NIC) that may be connected to another network device through a base station to communicate with the internet. In one example, the transmission device 106 may be a radio frequency (RF) module configured to communicate with the internet wirelessly.


It should be noted that, a true quantum computer is a hybrid structure that consists of two major parts: one is a classical computer responsible for classical computation and control; and the other is a quantum device responsible for running a quantum program and thereby implementing quantum computation. The quantum program is a series of instruction sequences written in a quantum language such as QRunes and able to run on a quantum computer, which supports operations of quantum logic gates and finally implements quantum computation. Specifically, the quantum program is a series of instruction sequences that operate quantum logic gates in a certain time sequence.


In practical applications, limited by the development of quantum device hardware, a quantum computation simulation is typically desired to verify a quantum algorithm, a quantum application, or the like. The quantum computation simulation is a process of simulating running of a quantum program corresponding to a specific problem in a virtual architecture (i.e., a quantum virtual machine) built with resources of a general computer. Generally, the quantum program corresponding to the specific problem has to be constructed. The quantum program, as used in the embodiments of the present disclosure, is a program written in a classical language that characterizes qubits and evolution thereof, in which qubits, quantum logic gates, and the like related to quantum computation are represented by corresponding classical codes.


As an embodied form of the quantum program, the quantum circuit, also called quantum logic circuit, is the most commonly used general-purpose quantum computation model that represents a circuit that operates on qubits under an abstract concept. The quantum circuit includes qubits, circuits (timelines), and various quantum logic gates, and results often need to be read out through quantum measurement operations.


Unlike a traditional circuit connected by metal wires and transmitting voltage or current signals, in the quantum circuit, the circuit may be regarded as being connected by time, which means that the state of qubits evolves naturally with time, and in this process, the circuit keeps being operated according to instruction of a Hamiltonian operator until encountering a logic gate.


One quantum program as a whole corresponds to one total quantum circuit, and the quantum program in the present disclosure refers to the total quantum circuit, where the total number of qubits in the total quantum circuit is the same as the total number of qubits in the quantum program. It may be understood as that: one quantum program may be composed of a quantum circuit, measurement operations for qubits in the quantum circuit, a register for saving measurement results, and control flow nodes (jump instructions). One quantum circuit may contain operations of tens, hundreds, or even thousands or ten thousands of quantum logic gates. An execution process of the quantum program is a process of executing all quantum logic gates according to a certain time sequence. It should be noted that the time sequence refers to a time order of a single quantum logic gate to be executed.


It should be noted that the basic unit in the classical computation is bit, the basic control mode is logic gate, and the circuits may be controlled through combinations of logic gates. Similarly, qubits are handled through quantum logic gates. A quantum logic gate may be used to evolve a quantum state. The quantum logic gate is a basis of the quantum circuit, and may include a single-bit quantum logic gate, such as a Hadamard gate (H gate), a Pauli-X gate (X gate), a Pauli-Y gate (Y gate), a Pauli-Z gate (Z gate), an RX gate, an RY gate, an RZ gate, or the like; or a multi-bit quantum logic gate, such as a CNOT gate, a CR gate, an iSWAP gate, a Toffoli gate, or the like. The quantum logic gate is generally represented in a unitary matrix, which is not only in a matrix form, but also a kind of operation and transformation. Generally, the effect of the quantum logic gate on the quantum state is calculated by multiplying a left side of the unitary matrix by a matrix corresponding to a right vector of the quantum state.


The quantum state, i.e., a logical state of the qubit, is expressed in binary in the quantum algorithm (or quantum program). For example, a group of qubits is q0, q1, and q2, representing the 0th, 1st, and 2nd qubits, which are ranked as q2q1q0 from high to low. The quantum state corresponding to this group of qubits is a superposition of eigenstates corresponding to the qubits in the group. The number of eigenstates corresponding to this group of qubits is 2 to the power of the total number of qubits, that is, 8 eigenstates (determined states): |000>, |001>, |010>, |011>, |100>, |101>, |110>, |111>, each having bits corresponding to the qubits, such as a |000> state, 000, corresponding to q2q1q0 from high to low, where |> is the Dirac symbol.


Taking a single qubit as an example for illustration, the logical state φ of the single qubit may be in a |0custom-character state, a |1custom-character state, and a superposition state (uncertain state) of a |0> state and a |1> state, which may be specifically expressed as φ=c|0custom-character+d|1custom-character, where c and d are complex numbers each representing a quantum state amplitude (probability amplitude), and squares of the amplitude modulus |c|2 and |d|2 respectively represent the probabilities of |0> state and |1> state, where |c|2+|d|2=1. In short, the quantum state is a superposition state of various eigenstates, and when the probability of other eigenstates is 0, it is in the only definite eigenstate.



FIG. 2 is a schematic flowchart of a method for solving a system of nonlinear equations on the basis of a quantum circuit according to an embodiment of the present disclosure. As shown in FIG. 2, the method may include the following steps S201 to S205.


At S201: a target system of nonlinear equations to be solved is acquired. In one embodiment, the target system of nonlinear equations to be solved may include a system of quadratic nonlinear equations or a system of nonlinear ordinary differential equations.


It may be understood that nonlinear problems are very common in nature, such as nonlinear finite element analysis, nonlinear dynamics, nonlinear programming, and the like, and therefore, it is important to construct an algorithm for solving nonlinear problems.


The nonlinear equation is named after a nonlinear relationship between a dependent variable and an independent variable, and is an extension of mathematics after the emergence of various practical problems and the establishment of equations based on problems in real life. The nonlinear equation has gained more and more attentions of people, and become an important research direction of modern mathematics.


The nonlinear equation provides key theoretical support and plays an important role in many fields of science and technology, such as in the fields of mechanics, economics, biotechnology, electronics, and the like. An accurate solution of such equations is hard to obtained, and thus, an approximate solution is usually desired, and corresponding methods for obtaining approximate solutions have gradually gained wide attention. As an important component of nonlinear equations, the system of quadratic nonlinear equations is of great significance both in theory and practice.


In addition, the ordinary differential equation (ODE), a differential equation whose unknown function contains only one independent variable, is gradually developed with calculus, and is an extension of mathematics after the emergence of various practical problems and the establishment of equations based on problems in real life. The ordinary differential equation has gained more and more attentions of people, and become an important research direction of modern mathematics. Similarly, the ordinary differential equation provides key theoretical support and plays an important role in many fields of science and technology, such as in the fields of mechanics, economics, biotechnology, electronics, and the like. These practical problems are finally converted into either obtaining a solution of a differential equation, or studying properties of a solution corresponding to the equation. Most problems in real life can be transformed into obtaining a special solution of a differential equation that satisfies given initial boundary value conditions.


As an important component of ordinary differential equations, the nonlinear dissipative system of ordinary differential equations is of great significance both in theory and practice. The system of nonlinear dissipative ordinary differential equations is far more complex than the system of linear ordinary differential equations, making it almost impossible to solve the system of nonlinear ordinary differential equations in an elementary integral method. Therefore, a method different from that of the linear differential equation theory has to be used to study a solution of a system of nonlinear ordinary differential equations.


In one embodiment, the system of quadratic nonlinear equations may be expressed as:








F
0

+


F
1


x

+


F
2



x




2





=
0




where x∈Rn, R represents real number space, Fi∈Rn×Rni, and F1 and F2 each have a sparsity s.


Exemplarily, it is assumed that the system of quadratic nonlinear equations is:






{






3


x
0


-

x
1

-

0.5

x
0
2


+

0.5

x
0



x
1


+
0.2

=
0








-

x
0


+

3


x
1


-

0.5

x
1
2


+

0.5

x
1



x
0


-
0.2

=
0








then the corresponding F0, F1, F2 are respectively:








F
0

=

(



0.2





-
0.2




)


,


F
1

=

(



3



-
1






-
1



3



)


,


and



F
2


=


(




-
0.5



0.5


0


0




0


0



0
.
5




-

0
.
5





)

.






Give Oracles OF1 and OF2 of non-zero element positions and non-zero element values obtained by querying F1 and F2:


















































O

F
11






"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"

0




=



"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"



f
1

(

j
,
k

)





,

0

j


n
-
1


,

0

k


s
-
1







O

F
12






"\[LeftBracketingBar]"

j









"\[LeftBracketingBar]"

k








"\[LeftBracketingBar]"

z




=



"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"

k








"\[LeftBracketingBar]"


F


1

j

,
k






,

0

j

,

k


n
-
1







O

F
21






"\[LeftBracketingBar]"

j









"\[LeftBracketingBar]"

0




=



"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"



f
2

(

j
,
k

)





,

0

j


n
-
1


,

0

k


s
-
1







O

F
22






"\[LeftBracketingBar]"

j









"\[LeftBracketingBar]"

k








"\[LeftBracketingBar]"

z




=



"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"

k








"\[LeftBracketingBar]"


F


2

j

,
k






,

0

j


n
-
1


,

0

j



n
2

-
1






Give an Oracle OF0 for preparing F0, which implements the function of:










O

F
0






"\[LeftBracketingBar]"

0




=


1



F
0











j
=
1


n
-
1




F

0
,
i






"\[LeftBracketingBar]"

i








In another embodiment, the system of nonlinear ordinary differential equations may be expressed as:









d

u


d

t


=



F
1



u

+


F
2




u


2





,


u

(
0
)

=

u
in






where u is a function to be solved of the system of nonlinear ordinary differential equations, u∈custom-charactern, custom-character is real number space, and F1′ and F2′ are sparse matrices independent of time with a sparsity s, F1′∈custom-charactern×n, and F2′∈custom-charactern×n2, that is, the number of non-zero elements in each row or column of F1′ and F2′ does not exceed s. Assuming that F1′ is a normal matrix, where D=diag (λ1, λ2, . . . , λn) and eigenvalues of F1′ satisfy Re(λn)≤ . . . ≤Re(λ1)<0, given preset parameters OF1′, OF2′ and Ou, OF1′ and OF2′ are respectively used to extract non-zero positions and values of F1′ and F2′, Ou is used to solve |uin), and OF1′, OF2′ and Ou are defined as:


































































O

F
11







"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"

0




=



"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"



f
1

(

j
,
k

)





,

0

j


n
-
1


,


0

k


s
-
1


;






O

F
12







"\[LeftBracketingBar]"

j









"\[LeftBracketingBar]"

k








"\[LeftBracketingBar]"

z




=



"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"

k








"\[LeftBracketingBar]"


F


1

j

,
k






,

0

j

,


k


n
-
1


;






O

F
13







"\[LeftBracketingBar]"

j









"\[LeftBracketingBar]"

0




=



"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"


g

(
j
)





,


0

j


n
-
1


;






O

F
21







"\[LeftBracketingBar]"

j









"\[LeftBracketingBar]"

k




=



"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"



f
2

(

j
,
k

)





,

0

j


n
-
1


,


0

k


s
-
1


;






O

F
22







"\[LeftBracketingBar]"

j









"\[LeftBracketingBar]"

k








"\[LeftBracketingBar]"

z




=



"\[LeftBracketingBar]"

j








"\[LeftBracketingBar]"

k








"\[LeftBracketingBar]"


F


2

j

,
k






,

0

j


n
-
1


,


0

j



n
2

-
1


;
and






O
u





"\[LeftBracketingBar]"

0





=




"\[RightBracketingBar]"





u
in

/



u
in






,




where f1(j,k) and f2(j,k) respectively represent column numbers of a kth non-zero element in row j of F1′ and F2′, g(j) satisfies f1(j,g(j))=j, and a diagonal element of F1′ is regarded as a non-zero element. A preset matrix (Oracle) related to F1′ is constructed with OF13′, and a parameter R is defined to limit









d

u


d

t


=



F
1



u

+


F
2




u


2





,


u

(
0
)

=

u
in


,




which is specifically expressed as:






R
:=


4




u
in







F
2








"\[LeftBracketingBar]"


Re

(

λ
1

)



"\[RightBracketingBar]"







Given R≤∥uin∥, if not satisfied, a suitable constant may be used to scale u to ζu, so that R≤∥uin∥.


At S202: the target system of nonlinear equations to be solved is converted to obtain a target system of linear equations.


In one embodiment, in the case where the target system of nonlinear equations is a system of quadratic nonlinear equations, the system of quadratic nonlinear equations may be converted into a preset system of pseudo-linear equations according to a homotopy perturbation method, and then the preset system of pseudo-linear equations may be converted into a target system of quadratic linear equations in a linear embedding method.


In another embodiment, in the case where the target system of nonlinear equations is a system of nonlinear ordinary differential equations, the system of nonlinear ordinary differential equations may be converted into a preset type of system of nonlinear ordinary differential equations in a homotopy perturbation method, and thus the preset type of system of nonlinear ordinary differential equations may be converted into a target system of linear ordinary differential equations in a linear embedding method.


It may be understood that the homotopy perturbation method is a method that combines the homotopy theory and the perturbation technique, which, differs from the traditional perturbation theory that depends on small parameters, constructs an equation containing embedded parameters by the homotopy technique, and then takes the embedded parameters as small parameters. Therefore, this method not only overcomes shortcomings of the traditional perturbation theory, but also makes full use of various perturbation methods. The essence of the homotopy perturbation method is to transform the nonlinear problem into an infinite number of linear problems to deal with. In this method, an approximate solution of an equation may be written as an addition of a series of infinite series, and this series sum converges to an exact solution thereof. A large number of examples has shown that this method is simple and effective, and a resulting first order approximation solution is often of high accuracy. It should be said that the homotopy perturbation method is a very common method for solving nonlinear problems.


In an implementation scenario, when the target system of nonlinear equations is a system of quadratic nonlinear equations, the system of quadratic nonlinear equations F0+F1x+F2x2=0 is converted into a preset system of pseudo-linear equations related to a series of v0, v1, . . . , vc according to the homotopy perturbation method.


First, a homotopy H is constructed: Rn×[0,1]→Rn, H(v,p)=F0+F1v+F2v2=0, assuming that a solution of H(v,p)=0 is expressed as: v=v0+pv1+p2v2+ . . . +pcvc, it can be obtained by substitution that:










F
1



v
0


+

F
0


=
0







F
1



v
1


+


F
2

(


v
0



v
0


)


=
0








F

1



v
2


+


F
2

(



v
0



v
1


+


v
1



v
0



)


=
0










F
1



v
c


+


F
2








j
=
0


j
=

c
-
1






v
j



v

c
-
1
-
j





=
0





Assuming F1 is invertible, the preset system of pseudo-linear equations is:








v
0

=


-

F
1

-
1







F
0








v
1

=


-

F
1

-
1








F
2


(


v
0



v
0


)







v
2

=


-

F
1

-
1






F
2

(



v
0



v
1


+


v
1



v
0



)










v
c

=


-

F
1

-
1





F
2








j
=
0


j
=

c
-
1






v
j



v

c
-
1
-
j









where vi is a variable to be solved in the preset system of pseudo-linear equations, 0≤i≤c, and the c is the number of variables to be solved in the preset system of pseudo-linear equations.


When p=1, we can get: {tilde over (x)}=v0+v1+ . . . +vc.


Based on the converted preset system of pseudo-linear equations, the preset system of pseudo-linear equations is then converted into a target system of quadratic linear equations in a linear embedding method.


Specifically, the equation F0+F1x+F2x2=0 is converted into a preset system of pseudo-linear equations taking v0, v1, . . . , vc as the variables to be solved according to the homotopy perturbation method. Then, the system of pseudo-linear equations taking v0, v1, . . . , vc as the variables to be solved is embedded into a larger linear system. In other words, the equations are embedded into a higher dimension linear system:






Ay=b


where y=y0, y1, . . . , yc, yi satisfies:








y
0

=


v
0

+

v
1

+

v
2

+

+

v
c







y
i

=

[


y

i
,
0


,

y

i
,
1


,

y

i
,
2


,


,

y

i
,


β
i

-
1




]






where βi represents the number of items in yi. A value of βi is inferred from the expression of yi,j, yi,j is written as:


yi,j=⊗k=0ivai,j,k, and ai,j,k satisfies ai,j,k≥0, and i+1≤Σk=0iai,j,k+1≤c+1. Therefore, it is inferred that yi contains












k
=
i

c



(



k




i



)





items in total, and βi may be expressed as:







β
i

=

{




1
,

i
=
0














k
=
i

c



(



k




i



)


,

1

i

c










Defining {right arrow over (a)}i,j=[ai,j,0, ai,j,1, . . . , ai,j,i], where {right arrow over (a)}i,j corresponds to j one by one, so the following two operations can be constructed:

































0

a
1




"\[RightBracketingBar]"





a



i
,
j








"\[RightBracketingBar]"




0



=



"\[RightBracketingBar]"





a



i
,
j








"\[RightBracketingBar]"




j








0

a
2






"\[RightBracketingBar]"




j






"\[RightBracketingBar]"




0



=




"\[RightBracketingBar]"




j





"\[RightBracketingBar]"





a



i
,
j








Oa1 and Oa2 may be implemented with quantum arithmetic of a O(c) bit quantum circuit with a gate complexity O(ploy(c)).


For yi,j=vai,j,0⊗vai,j,1⊗ . . . ⊗vai,j,i, when {right arrow over (a)}i,j contains a value greater than 0, assuming that a first value greater than 0 is ai,j,k, then











I
n






k






F
1



I
n







i
-
k






)



y

i
,
j



+


v

a

i
,
j
,
0









v

a

i
,
j
,

k
-
1







F
2

(








1
-
0



a

i
,
j
,
k


-
1





v
1



v

a

i
,
j
,
k





-
1
-
1

)







v

a

i
,
j
,
i





=
0




where vai,j,0⊗ . . . ⊗vl⊗vai,j,k−1−l⊗ . . . ⊗vai,j,i∈yi+1, and if we take yi,0=v0⊗i+1, then: F1⊗i+1yi,0=(−F0)⊗i+1.


From the above equations, we obtain bi=[−F0)⊗+1, 0, . . . , 0], and the final equation Ay=b may be expanded as:








[




A

0
,
0





A

0
,
1




0


0


0


0




0



A

1
,
1





A

1
,
2




0


0


0




0


0



A

2
,
2







0


0




0


0


0






A


c
-
1

,

c
-
2





0




0


0


0


0



A


c
-
1

,

c
-
1






A


c
-
1

,
c






0


0


0


0


0



A

c
,
c





]

[




y
0






y
1






y
2











y

c
-
1







y
c




]

=

[




b
0






b
1






b
2











b

c
-
1







b
c




]





where Ai,i is a







n

i
+
1









j
=
i

c



(



j




i



)

×

n

i
+
1









j
=
i

c



(



j




i



)





dimension matrix, Ai,i+1 is a ni+1βi×ni+2βi+1 dimension matrix, y=y0, y1, . . . , yc, yi satisfies y0=v0+v1+v2+ . . . +vc, yi=[yi,0, yi,1, yi,2, . . . , yi,βi−1], and βi represents the number of items in yi.


To optimize the linear system Ay=b, it can be seen from the above equation that A contains a block matrix F1⊗i, i=1, . . . , c, which makes the condition number of A increase with a c exponent, so F1⊗i+1yi,0=(−F0)⊗i+1 is split into:








[





F
1



I
n






i







I





0




0




I
n



F
1



I
n







i
-
1








I


0








0














0













I
n






i






F
1





]

[




v
0







i
+
1











F
0



v
0






i
















F
0







i
-
1







v
0






2











F
0






i






v
0





]

=

[



0




0









0





-

F
0







i
+
1









]





The split linear system has a better condition number and redefines yi,0 as:







y

i
,
0


=

[


v
0







i
+
1





,


F
0



v
0






i





,


,


F
0






i






v
0



]





Continue with the above example, if we take c=2, then components of the optimized y are written as:










y
0

=

[


v
0

+

v
1

+

v
2


]








y
1

=

[



v
0



v
0


,


F
0



v
0


,


v
0



v
1


,


v
1



v
0



]








y
2

=

[


v
0






3




,


F
0



v
0



v
0


,


F
0



F
0



v
0



]








The corresponding matrix A may be written as:






(




F
1




F
2



0



F
2




F
2



0


0


0




0




F
1



I
2





I
4



0


0


0


0


0




0


0




I
2



F
1




0


0


0


0


0




0


0


0




I
2



F
1




0




I
2



F
2




0


0




0


0


0


0




F
1



I
2






F
2



I
2





0





0






0


0


0


0


0




F
1



I
4






I
8






0






0


0


0


0


0


0




I
2



F
1



I
2





I
8





0


0


0


0


0


0


0




I
4



F
1





)







b
=


[


-

F
0


,

0
4

,

-

F
0






2





,

0
4

,

0
4

,

0
8

,

0
8

,

-

F
0






3






]

T


,




where 04=[0,0,0,0], and 08 is similar to 04.


The linear system Ay=b is also adjusted accordingly. In the following, the Ay=b is default as the optimized linear system, and the optimized matrix A has a dimension number of:






N
=








i
=
0

c




n

i
+
1


(


β
i

+
i

)


=




(

n
+
1

)


c
+
1


-
1
-

s

n

+





cn
c

(

n
-
1

)

-

n
c

+
1



(

n
+
1

)

2




n
2







(

n
+
1

)


c
+
1


+


cn

c
+
1


.








For the condition number of the matrix A, some lemmas are given as follows:


Lemma 1: given an n dimension invertible matrix M and i∈custom-character+, and a new matrix is defined:






P
=

[




M


I
n






i







I





0




0




I
n


M


I
n







i
-
1








I


0








0














0













I
n






i





M




]





then P is invertible and satisfies









P

-
1










M

-
1






(

1
-




M

-
1





i
+
1



)



1
-



M

-
1










Lemma 2: ∥A∥ satisfies ∥A∥≤∥F1∥+1+(c+1)∥F2∥.


Lemma 3: given the optimized linear system Ay=b, a truncation order c, when F1, F2 satisfy: |F1-1|<1, and












F
1

-
1





1
-



F
1

-
1








(

c
+
1

)





F
2




<
1

,

A

-
1






satisfies:









A

-
1









F
1

-
1









1
-




F
1

-
1






(

1
+

(

c
+
1

)









F
2




)






Therefore, for a given optimized linear system Ay=b and truncation order c, when F1, F2 satisfy:










F
1

-
1




<
1

,


and






F
1

-
1





1
-



F
1

-
1








(

c
+
1

)





F
2




<
1

,




the condition number κ(A) of the matrix A satisfies:







κ

(
A
)





κ

(

F
1

)

+
1


1
-




F
1

-
1






(

1
+


(

c
+
1

)





F
2





)








In another implementation scenario, when the target system of nonlinear equations is a system of nonlinear ordinary differential equations, the system of nonlinear ordinary differential equations is converted into a preset type of system of nonlinear ordinary differential equations in a homotopy perturbation method.


Exemplarily, a homotopy v(t,p): custom-charactern×[0,1]→custom-charactern is constructed in the homotopy perturbation method, which satisfies:








H

(

v
,
p

)

=



dv
dt

-


F
1



v

-

p


F
2




v


2




=
0


,


v

(

0
,
p

)

=

u

i

n







Assuming that v is expressed as: v=v0+pv1+p2v2+ . . . +pcvc, items with the same power of p are equivalent to obtain the preset type of system of nonlinear ordinary differential equations converted from








du
dt

=



F
1



u

+


F
2




u


2





,


u

(
0
)

=

u

i

n



,




which is specifically:










dv
0

dt

-


F
1




v
0



=
0

,



v
0

(
0
)

=

u

i

n













dv
1

dt

-


F
1




v
1


-


F
2





v
0



v
0




=
0

,



v
1

(
0
)

=
0











dv
2

dt

-


F
1




v
2


-


F
2


(



v
0



v
1


+


v
1



v
0



)


=
0

,



v
2

(
0
)

=
0
















dv
c

dt

-


F
1




v
c


-


F
2









j
=
0


c
-
1





v
j



v

c
-
1
-
j





=
0

,



v
c

(
0
)

=
0





where the c is the number of functions to be solved in the preset type of system of nonlinear ordinary differential equations, and vi is a function to be solved of the preset type of system of nonlinear ordinary differential equations, and 0≤i≤c. The preset type of system of nonlinear ordinary differential equations is a series of ordinary differential equations with nonlinear variables vi.


When p=1, ũ=v0+v1+v2+ . . . +vc.


Based on the converted preset type of system of nonlinear ordinary differential equations, the preset type of system of nonlinear ordinary differential equations may be converted into a target system of linear ordinary differential equations in a linear embedding method.


Specifically, the preset type of system of nonlinear ordinary differential equations is embedded into a system of linear ordinary differential equations with y as the viable, i.e., a target system of linear ordinary differential equations in a linear embedding method:









d


y



dt

=

A


y




,



y


(
0
)

=

y

i

n







where {right arrow over (y)}=[{right arrow over (y)}0, {right arrow over (y)}1, {right arrow over (y)}2, . . . , {right arrow over (y)}c], and {right arrow over (y)}i satisfies:







y


=

{





[


v
0

+

v
1

+

v
2

+

+

v
c


]

,

i
=
0








[



y



i
,

0


,


y



i
,

1


,


y



i
,

2


,


,


y



i
,



β
i

-
1




]

,

1

i

c










βi represents the number of items in {right arrow over (y)}i, {right arrow over (y)}i,j represents a jth item in {right arrow over (y)}i, expressed as {right arrow over (y)}i,j=⊗k=0ivai,j,k, and ai,j,k satisfies: ai,j,k≥0, and i+1≤Σk=0iai,j,k+1≤c+1.


βi satisfies:







β
i

=

{




1
,


i
=
0














k
=
1

c



(



k




i



)


,

1

i

c










yin may be written as:







y

i

n


=

[


[

u

i

n


]

,



[


u

i

n



2


,
0
,


,
0

]

[


u

i

n



3


,
0
,


,
0

]

[

u

i

n



c


]


]





Define {right arrow over (a)}i,j=[ai,j,0, ai,j,1, . . . , ai,j,i] and construct two Oracles:














O

a
1


|


a



i
,

j





|
0



=

|


a



i
,

j







j
















O

a
2


|
j



|
0



=

|
j




|


a



i
,

j








the dimension number of {right arrow over (y)}i is ni+1βi, so the dimension number N of custom-character is:






N
=








i
=
0

c



n

i
+
1




β
i


=




(

n
+
1

)


c
+
1


-
1
-
cn





(

n
+
1

)


c
+
1


.







For initial state preparation of |yincustom-character, we can firstly define:














y

i

n





=



1

M









i
=
0

c



i


,
0



|

u

i

n




i
+
1









where M=Σj=ic∥uin2(i+1), and given the Ou|0custom-character=|uin/∥uincustom-characterdefined above, |yincustom-character may be prepared by querying a frequency of OuO(c).


Specifically,










|
Ψ



=



1

M









i
=
0

c






u

i

n





2


(

i
+
1

)




|
i


,
0






is firstly prepared, then a controlled Ou operation










C
-

O
u


=



1

M









i
=
0

c


|
i


,
0







i
,

0




O
u



i
=
1











is performed to obtain the above initial state of |yincustom-character.


Finally, after constructing the preset matrices (Oracles) of Ai,i and Ai,i+1, the preset matrix (Oracle) of A can be directly constructed by querying the preset matrices (Oracles) of Ai,i and Ai,i+1 once. Therefore, the preset matrix (Oracle) OA may be constructed by querying OF1O(c) and OF2O(1) frequencies.


Exemplarily,









d


y



dt

=

A


y




,



y


(
0
)

=

y

i

n







is solved with a quantum algorithm, and {right arrow over (y)}(t) is written into:









y


(
t
)

=


e

A

t





y


(
0
)



,




define:








T
k

(
z
)

:=







j
=
0

k





z
j


j
!


.






When k is large enough and the evolution time h is short (for example, h≤1/∥A∥), {right arrow over (y)}(h)≈Tk(Ah){right arrow over (y)}(0), and this approximate solution may be used as an initial condition for a next approximation. This process is repeated for m steps, and then the approximation of {right arrow over (y)}(mh) can be obtained.


Let m, k, p∈custom-character+, and define:










C

m
,

k
,

p


(
A
)

:=







j
=
0

d


j








j




I

-







i
=
0


m
-
1









j
=
1

k






i

(

k
+
1

)

+
j










i

(

k
+
1

)

+
j
-
1






A
j


-







i
=
0


m
-
1









j
=
0

k







(

i
+
1

)



(

k
+
1

)


+
j










i

(

k
+
1

)

+
j





I

-






j
=

d
-
p
+
1


d



j









j
-
1




I


,






where d: =m(k+1)+p.


Considering the linear system: Cm,k,p(Ah)|xcustom-character=|0custom-character|yincustom-character where |yincustom-charactercustom-charactern, h∈custom-character+, after evolution, an approximate solution of k order Taylor series is obtained, and the solution at a pth order remain unchanged. In this case, a solution of Cm,k,p(Ah)|xcustom-character=|0custom-character|yincustom-character is expressed as: |xcustom-character=Cm,k,p(Ah)−1|0custom-character|yincustom-character, or may be written as:
















|
x



=








i
=
0


m
-
1









j
=
0

k


|


i

(

k
+
1

)

+
j





|

x

i
,

j





+






j
=
0

p


|


m

(

k
+
1

)

+
j




|

x

m
,

j









and








y

i

n










satisfies
:









|

x

0
,

0





=

|

y

i

n















|

x

i
,

0





=







j
=
0

k

|

x


i
-
1

,

f






,

0

i

m











|

x

i
,

1





=

Ah
|

x

i
,

0






,

0

i

m











|

x

i
,

j





=


(

Ah
j

)

|

x


i
-
1

,

j






,

0

i

m

,

2

j

k











|

x

m
,

j





=

|

x

m
,


j
-
1







,

1

j

p







Then
:









|

x

0
,

0





=

|

y

i

n















|

x

0
,

j





=


(



(
Ah
)

j

/

j
!


)

|

x

0
,

0






,

1

j

k













|

x

1
,

0





=



T
k

(
Ah
)

|

x

0
,

0








exp

(
Ah
)


|

y

i

n














|

x

1
,

j





=


(



(
Ah
)

j

/

j
!


)

|

x

1
,

0






,

1

j

k













|

x

2
,

0





=



T
k

(
Ah
)

|

x

1
,

0








exp

(

2


Ah

)


|

y

i

n





















|

x


m
-
1

,

0





=



T
k

(
Ah
)

|

x


m
-
2

,

0








exp

(

Ah

(

m
-
1

)

)


|

y

i

n














|

x


m
-
1

,

j





=


(



(
Ah
)

j

/

j
!


)

|

x


m
-
1

,

0






,

1

j

k













|

x

m
,

0





=



T
k

(
Ah
)

|

x


m
-
1

,

0








exp

(

2


Ahm

)


|

y

i

n

















|

x

m
,

j





=

|

x

m
,

0








exp

(
Ahm
)


|

y

i

n





,

1

j

p





|xi,0custom-character is an approximate solution of the system at time ih, i=∈{0,1,2, . . . , m}, and |xm,0custom-character=|xm,1custom-character=|xm,2custom-character= . . . =|xm,pcustom-character is an approximate solution of custom-character(t)=eAtyin at t=mh.


At S203: a quantum circuit corresponding to a quantum linear solver used for solving the target system of linear equations is constructed. In one embodiment, a preset matrix Oracle related to the target system of linear equations may be firstly constructed, and then a quantum logic gate function module is constructed based on the Oracle, so as to construct a quantum circuit corresponding to a quantum linear solver based on the quantum logic gate function module, where the quantum logic gate function module includes a first functional module, a second functional module, and a third functional module. In other words, a quantum circuit corresponding to a quantum linear solver including the Oracle and quantum logic gate functional modules is constructed, a quantum state evolution operation is performed on the quantum circuit, and measurement is performed to obtain a quantum state of the evolved quantum circuit. For ease of understanding, construction of the quantum circuit will first be described in detail below in conjunction with FIGS. 3 to 5.



FIG. 3 is a schematic diagram of a quantum circuit for solving a system of nonlinear equations according to an embodiment of the present disclosure. As shown in FIG. 3, a quantum linear quadratic linear solver module and four measurement modules are shown. By constructing a related Oracle of Ay=b, the Oracle may be regarded as an interface for inputting equation information to the quantum circuit, or as input of a quantum linear solver algorithm. Specifically, a quantum state (y) may be output though Oracleinput the Oracle of Ay=b, |ycustom-characteri=0cΣj=0βi|i,jcustom-character|yi,jcustom-character. A first bit register of |ycustom-character is measured, and a normalized approximate solution |y0.0custom-character=|{tilde over (x)}custom-character of the original quadratic nonlinear equation is obtained when |0,0custom-character is measured.


In some embodiments, a system of quantum linear ordinary differential equation solver module and five measurement modules may be included. By constructing a related Oracle of









d


y



dt

=

A


y




,



y


(
0
)

=

y

i

n



,




the Oracle may be regarded as an interface for inputting equation information to the quantum circuit, or as input of a quantum linear solver algorithm for a system of ordinary differential equations. Specifically, the Oracle of








d


y



dt

=

A


y







and evolution time T may be input to output a quantum state |y(T)custom-character at time T. Some bit registers of |y(T)custom-character are measured, and |{right arrow over (y)}0,0(T)custom-character is obtained when |0,0custom-character is measured, where the obtained |{right arrow over (y)}0,0(T)custom-character satisfies: |y0,0(T)custom-character=|ũ(T)custom-character=|u(T)custom-character, so that an output state |ũ(T)custom-character representing an approximate solution of the target system of linear ordinary differential equations can be obtained.


For the construction manner of OracleOA of the matrix A. A may be regarded as a c+1-dimensional block matrix, and an expression and corresponding positions of non-zero block matrix elements of A may be obtained by O(poly(c)) times of classical operations. Then, the Oracles of F1, F2 are used to extract non-zero element positions and non-zero element values of matrix elements in the block matrix, and OA is constructed by implementing arithmetic operations in this process with a quantum circuit, where a query complexity is O(poly(c)) for Oracles of F1, F2 in the construction process. Since only some simple arithmetic operations are involved, a circuit length is also O(poly(c)). Similarly, the operation Ob for preparing |bcustom-character can also be implemented by OF0, where the query complexity is O(poly(c)).


It can be understood that constructing the quantum circuit corresponding to the quantum linear solver mainly includes constructing a sub-quantum circuit corresponding to the quantum linear solver module shown in FIG. 3, including: firstly constructing a quantum circuit that implements a first functional module V, where the first functional module is defined as










V
|

0
m




:=



1

α








i




α
i



|
i




,







α
:=

4








j
=
0


j
0





(

-
1

)

j










i
=

j
+
1


b



(




2

b






b
+
i




)



2

2

b





,








j
0

=


blog
(

4

b
/
ϵ

)



,







b
=


κ
2



log

(

κ
/
ϵ

)



;




next, constructing a quantum circuit that implements a second functional module combination, where the quantum circuit for the second functional module combination includes T, W and T+, the T=Σj∈[N]jcustom-characterj|, W=S·(2T·T+−I4N2), the S is a unitary matrix of an exchange operation module, the T+ is a transposed conjugate of the T, the

















|

Ψ
j




:=

|
j






1

d











k


[
N
]


:


A
jk


0





(



A
jk
*


|
k





+


1
-



"\[LeftBracketingBar]"


A
jk



"\[RightBracketingBar]"





|

k
+
N




)

,




and I4N2 is a 4N2-dimensional identity matrix; further constructing a quantum circuit that implements a third functional module V+, where the third functional module is a transposed conjugate form of the first functional module; and sequentially inserting the first functional module, the second functional module and the third functional module into the quantum circuit to form a quantum circuit corresponding to a quantum linear solver, as shown in FIG. 4.



FIG. 4 is a schematic diagram of a quantum circuit corresponding to a quantum linear solver according to an embodiment of the present disclosure. In the quantum circuit shown in FIG. 4, V, T represent Oracles with different functions, “†” represents transposed conjugate. T represents an overall functional module T of an H gate and a combination of various Oracles, a function of the T module is to transform |jcustom-character into |Ψjcustom-character, the obtained matrix input to the T module is an N order matrix, and the constructed T module may be equivalent to a quantum logic gate in the quantum circuit, the matrix form of which is: Σj∈[N]jcustom-characterj|, where custom-characterj| is a quantum state bra.


It should be noted that to solve the target system of linear equations, a quantum circuit diagram about a walk operator W should be constructed first. It will be understood by those skilled in the art that any simple function can be linearly approximated as a linear combination of other functions, and an inverse function of a matrix can be approximated through a Chebyshev polynomial. The Chebyshev polynomial has to be implemented in the framework of quantum walks.


Since quantum walks are performed on some states |φjcustom-character in space custom-character2Ncustom-character2N define a mapping T:=Σj∈[N]jcustom-characterj|, from custom-characterN to custom-character2Ncustom-character2N
















|

Ψ
j




:=

|
j






1

d











k


[
N
]


:


A
jk


0





(



A

j

k

*


|
k





+


1
-



"\[LeftBracketingBar]"


A
jk



"\[RightBracketingBar]"





|

k
+
N




)




and a walk operator:






W
=

S

(


2

T


T



-
I

)





The operator S performs a flip operation of a product state in |φj), so:

















T




W
n


T


b



=

|

0
m






(
H
)


|
b



+

|


b








custom-character(λ) is a Chebyshev polynomial of the first kind.


It should be noted that as in the above form of |Ψjcustom-character, a combination of vertical lines and angle brackets is used to describe a quantum state, representing that the quantum state is a vector (also called a state vector, a basis vector, or the like), where |Ψjcustom-character represents a ket, and custom-characterΨj| represents a bra.



FIG. 5 is a schematic diagram of a quantum circuit for a walk operator W according to an embodiment of the present disclosure. As shown in FIG. 5, since W=S(2TT−I4N2), S may be constructed through a group of exchange operations (such as a SWAP gate, where the symbol of two connected bold Xes in the qubits of FIG. 5 represent a SWAP gate), and the rest is 2TT−I4N2.


To construct 2TT−I4N2, a unitary operator form of T is desired to be constructed, and the unitary operator Tu is defined to satisfy:












T
u

|
j




0



=

|

Ψ
j








so:

















T
u

(






j

|
j








j
|



(

2
|
0










0
|

-

I

2

N






)

)



T
u



=






j

|

φ
j










φ
j

|



(

2
|
0










0
|

-

I

2

N






)

=


2

T


T



-


I


4


N
2









where















φ
j




=


1

d










k


[
N
]


:


A
jk


0





(



A

j

k

*


|
k






+


1
-



"\[LeftBracketingBar]"


A
jk



"\[RightBracketingBar]"





|

k
+
N




)

,





and:







K
=

2
|
0








0
|

-


I

2

N


.








It should be noted that this schematic diagram merely shows a part of the quantum circuit related to the present application, and the identifiers and connection relationships in the figure are merely examples and do not constitute any limitation to the present disclosure. In addition to the Chebyshev linear solver mentioned above, an HHL algorithm or a variational quantum linear solver can also be used for solution.


Returning to FIG. 2, it further shows that at S204: based on the quantum circuit corresponding to the quantum linear solver, quantum state evolution and measurement are performed on the target system of linear equations to solve the target system of linear equations. In one embodiment, a target matrix in the target system of linear equations is firstly input into the quantum circuit corresponding to the quantum linear solver through the preset matrix Oracle, and based on the input, the quantum linear solver executes the quantum state evolution and measures a quantum register in the corresponding quantum circuit to solve the target system of linear equations. Specifically, according to the constructed matrix A and the corresponding Oracle, a quantum linear solver is used to solve the linear system. The input of the quantum linear solver is the OracleOracle constructed as above, and a solution of the target system of linear equations is obtained by the quantum linear solver.


It should be noted that the quantum OracleOracle is a black box that represents transformation of a certain quantum state. A typical example of the quantum OracleOracle is linear system: O|xcustom-character|0custom-character=|xcustom-character|f(x)custom-character, where in calculation of f(x), the first quantum register is used as input and the second quantum register is used as output. Another example is QRAM, which can be regarded as a kind of OracleOracle. Many quantum algorithms use OracleOracles, but do not consider implementation of the OracleOracle, while it may be decomposed into quantum gates, or implemented as a QRAM. In Qpanda, the “OracleOracle” function can be used to define. OracleOracle is considered to have a name provided by a user.


In quantum applications, an OracleOracle or a combination of Oracles is constructed, and the internal principle of the Oracle or combination is the method flow of the present disclosure. Specifically, the Oracle may be understood as a module (similar to a black box) that completes a specific function in a quantum algorithm, and has specific implementations for specific problems.


At present, the existing quantum circuit can be constructed by means of only the existing single quantum logic gate, double quantum logic gate, and the like, and usually has the following defects:


For a quantum circuit with relatively complex functions, a very large number of qubits are desired. The simulation with a classical computer will consume a huge memory space, involve a very large number of logic gates and a very long time, and some complex algorithms are difficult to implement with the quantum circuit.


On this basis, specific complex functions are realized through Oracle simulation, and controlled and transposed conjugate operations are implemented. Parameters input by a user to the Oracle may include: an Oracle name (used to identify a functional purpose of the Oracle, such as OA1), qubits, matrix elements, and the like.


The advantage of this method is that the Oracle is regarded as a known module as a whole without paying attention to internal implementation details, and the representation of the quantum application scenario such as a quantum circuit, is very simple and clear. Since the Oracle functional modules of classical simulation can be equivalent to quantum logic gates, the constructed quantum circuit is simplified, thereby saving the memory space required for runtime and speeding up simulation verification of the quantum algorithm.


At S205: based on an obtained solution of the target system of linear equations, a solution of the target system of nonlinear equations to be solved is determined.


In an implementation scenario, the quantum state evolution of the quantum state from |bcustom-character to |A−1bcustom-character may be performed by the quantum circuit shown in FIG. 4. Exemplarily, the entire quantum circuit is operated to measure |jcustom-character and |anc), and when both |jcustom-character and |anccustom-character are collapsed to |0custom-character, |A−1bcustom-character can be obtained in the second register.


Similarly, according to the quantum circuit corresponding to the quantum linear solver, a solution of the target system of linear equations can be obtained by measuring some quantum registers on the corresponding quantum circuit, and finally, according to the obtained solution of the target system of linear equations, a solution of a system of quadratic nonlinear equations or of nonlinear ordinary differential equations can be calculated. For a system of quadratic nonlinear equations, since each component in {tilde over (x)} is in the form of a tensor product of vi, the solution of the system of quadratic nonlinear equations is obtained by calculating y0. For a system of nonlinear ordinary differential equations, since each component in {right arrow over (y)} is in the form of a tensor product of vi, the solution of the system of nonlinear ordinary differential equations to be processed is obtained by calculating {right arrow over (y)}(0). For example, some qubit registers of |y(T)) may be measured to obtain a quantum state with a preset accuracy close to a normalized solution of the equation to be solved. The measurement may be divided into two steps: (1) measuring a first qubit register in the defined |xcustom-character, where if the measured value is s, s=|m(k+1)+jcustom-character, j=0, 1, . . . , p, a second quantum register of |xcustom-character has |y(T)custom-character; and (2) measuring a first qubit register of |y(T)custom-character, where if the result is |0,0custom-character, an exact solution with a preset accuracy close to |u(T)/∥u(T)|∥custom-character is obtained in a second qubit register of |y(T)custom-character.


It should be noted that the first qubit register and the second qubit register are a first register and a second register of an output state of the quantum circuit corresponding to the quantum linear solution algorithm, which are subdivisions of the quantum circuit where |A−1bcustom-character is located in FIG. 4. The quantum circuit in FIG. 4 merely shows a part of the quantum circuit related to the present application, and the identifiers and connection relationships in the figure are merely examples and do not constitute any limitation to the present.


Further, for a success rate of the algorithm, |y0custom-character is desired. The success rate is a probability of |0,0custom-character obtained by measuring the first quantum register in |ycustom-characteri=0cΣj=0βi|i,jcustom-character|yi,jcustom-character, and specifically expressed as:







p
=





y
0



2




y


2



,




y=[y0, y1, . . . , yc], assuming ∥y0∥=ηR, where η is a constant. R is defined as: R:=max{4∥F1−12∥F0∥∥F2∥, ∥F0∥}, then |F1−1|<1, and when







R



2

2


:






y
0



2




y


2







η
2

(

1
-

2


R
2



)




η
2

(

1
-

2


R
2



)

+
4


.






Therefore, given the defined system of quadratic nonlinear equations F0+F1x+F2x⊗2=0 and define the parameters: ϵ<0.01, R:=max {4∥F1−12∥F0∥∥F2∥, ∥F0∥}, α:=∥F1−1∥∥F0∥, β:=∥F1−1∥∥F2∥, η=∥x∥/R.







c
=




log

1

4

αβ






4

α


η

R


ϵ

(

1
-

4

αβ


)







,




and G:=∥F1−1∥(1+(c+1)∥F2∥, then it satisfies: when ∥F1−1∥<1, G<1, and







R
<


2

2


,




there is a quantum algorithm with a success rate Ω(1) used to get the normalized quantum state |{tilde over (x)}custom-character, which satisfies ∥{tilde over (x)}−x∥≤ϵ, where x is an exact resolution. A query complexity of the algorithm for Oracles of F0, F1 and F2 is O(s(A)κ(A) poly(log(s(A)κ(A)/ϵ))). Considering the success rate of the algorithm, and through optimization by an amplitude amplification algorithm, the query complexity at the success rate Ω(1) is:







O

(



κ

s


η



1
-

2



R
2

(

1
-
G

)








poly

(

log

(






"\[LeftBracketingBar]"


F
1

-
1









F
0





ϵη

R


)

)


)

.




Following the above example, solve Ay=b, and get the component y0 of y as:







x
~

=


y
0

=


[



-
4.8765625

×

10

-
2



,

5.1265625
×

10

-
2




]

T






A solution obtained in an iterative method is:






x
=

[



-
4.8765625

×

10

-
2



,

5.1265625
×

10

-
2




]






so:








x
-

x
~




=

1.1061184
×


10

-
6


.






In addition, a structure of the matrix A in the nonlinear ordinary differential equation needs to satisfy the form that:








d



y



i
,

j



dt

=



(







j
=
0

i




I
n


j




F
1




I
n



i
-
j





)




y



i
,

j



+







k
=
0

i



(


I
n
k



F
2




I
n



i
-
k




)




ν

a

i
,

j
,

0









v

α

i
,

j
,


k
-
1






(








l
=
0



a

i
,

j
,

k


-
1




v

a

i
,

j
,

k




-
1
-
l

)



ν

a

i
,

j
,


k
+
1










v

a

i
,

j
,

i










where I is an n order identity matrix, vai,j,0⊗ . . . ⊗vai,j,k−1⊗vl⊗vai,j,k−1−l⊗vai,j,k+1⊗ . . . ⊗vai,j,i∈{right arrow over (y)}i+1, so









d


y



dt

=

A


y




,





custom-character(0)=yin, which may be written as:








d
dt

[





y


0







y


1







y


2












y



c
-
1








y


c




]

=


[




A

0
,

0





A

0
,

1




0


0


0


0




0



A

1
,

1





A

1
,

2




0


0


0




0


0



A

2
,

2







0


0




0


0


0






A


c
-
1

,


c
-
2





0




0


0


0


0



A


c
-
1

,


c
-
1






A


c
-
1

,

c






0


0


0


0


0



A

c
,

c





]

[





y


0







y


1







y


2












y



c
-
1








y


c




]





where Ai,i is a ni+1βi-dimensional array, expressed as:







A

i
,

i


=

I



β
i



(







j
=
0

i




I
n
j



F
1




I
n



i
-
j





)







Ai,i+1 is a ni+1βi×ni+2βi+1-dimensional matrix, and |y(t)) is defined to represent custom-character(t):














y

(
t
)




=








i
=
0

c








j
=
0


β
i




i


,
j



|


y

i
,

j


(
t
)







The above quantum linear algorithm is used to solve









d


y



dt

=

A


y




,





custom-character(0)=yin, and custom-character(t) is written as: custom-character(t)=eAtcustom-character(0).


It can be seen that in the present application, a system of quadratic nonlinear equations or a system of nonlinear ordinary differential equations is converted into a matrix and vector information of linear equations in a homotopy perturbation method, and encodes the same to a quantum state, the classical data structure is connected with quantum states in the quantum field, and the evolution operation of encoding the classical data structure into a quantum state is performed to obtain a quantum state of the evolved quantum circuit, which can accelerate solution of the system of quadratic nonlinear equations or system of nonlinear ordinary differential equations of high complexity by means of the superposition property of quantum, thereby expanding the simulation application scenarios of quantum computation.


Compared with the existing art, the present disclosure firstly converts a system of quadratic nonlinear equations or a system of nonlinear ordinary differential equations to determine a target system of linear equations, and then constructs a quantum circuit corresponding to the quantum linear solver, runs the quantum circuit and performs measurement to solve the target system of linear equations, and determines the solution of the system of quadratic nonlinear equations or the system of nonlinear ordinary differential equations based on the solution of the target system of linear equations. Therefore, by means of related properties of quantum, the technology of solving a system of quadratic nonlinear equations or a system of nonlinear ordinary differential equations with a quantum algorithm can be implemented, and the complexity and difficulty in solving a system of quadratic nonlinear equations or a system of nonlinear ordinary differential equations can be reduced, thereby filling in related technical gaps in the field of quantum computation.


According to the foregoing, it is possible to embed the preset system of pseudo-linear equations into the linear system (i.e., form a system of linear equations), and then solve the system of linear equations. In recent years, quantum algorithms for solving a system of linear equations Ax=b have been continuously proposed. For example, Childs et al. have proposed an algorithm for constructing matrix inverses with series expansion, and constructed a new quantum linear solver with a complexity of custom-character(sκ2poly(log(sκ/ϵ))), where performance on ϵ is improved exponentially, and in addition, the complexity is also optimized to O(sκpoly(log(sκ/ϵ))) through VTAA. Recently, the work of using QSVT to construct a quantum linear solver is also proposed, which does not need to be optimized with VTAA, and the complexity is custom-character(sκpoly(log(Sκ/ϵ))).


From the above work, it can be seen that the complexity of the optimal quantum linear solver is linearly related to the system condition number κ, but in practical problems, a linear system with a large condition number may be encountered, such as x is on a polynomial magnitude of n. When solving such problems, the advantage of the quantum linear solver is severely reduced, and there is no quantum advantage anymore. When encountering an ill-conditioned linear system, the classical algorithm for solving a system of linear equations generally pre-processes the linear system to construct a new well-conditioned linear system, and the solution is the same as the original system, thereby effectively solving the ill-conditioned linear system. In view of this, the present application further proposes to add the classical polynomial pre-processing to the HHL algorithm, to construct an HHL algorithm with a polynomial preprocessor, which reduces the condition number of the linear system to be solved, thereby improving the performance of the HHL algorithm and expanding the application scope of the HHL algorithm.


In one embodiment, firstly, the above target system of linear equations may be embedded into a linear system, where the linear system includes a first matrix A and a first vector b, and the linear system includes a condition number κA.


The linear system is a mathematical model that refers to a system composed of linear operators that satisfies both superposition and uniformity (also known as homogeneity). At present, the linear system is the core of many scientific and engineering fields. For a linear system Ax=b to be processed, element information and dimensions of a first matrix A and a first vector b can be obtained respectively, where the first matrix A may be an invertible Hermitian matrix. The Hermitian matrix is a common matrix type that can be simply put as a self-conjugate matrix, in which each element at column j in row i of the matrix is conjugated and equal to an element at column i in row j. Each element on a main diagonal of the Hermitian matrix is a real number with an eigenvalue also being a real number.


Next, a polynomial p(A) for the first matrix A in the linear system is calculated.


Specifically, for the linear system Ax=b to be processed, A is an invertible Hermitian matrix with a condition number being κA=∥A∥∥A−1∥. A preprocessing process of the linear system is to find and calculate the polynomial p(A) for the first matrix A, which may include the following steps:


At step 1: a function ƒm(x) and a preset interval S for solving the polynomial are constructed, where the function ƒm(x) satisfies fm(x)=pm (x)x, and fm(xi)+(−1)iE=1, the m is a polynomial order, the i=0,1,2, . . . , m+1, the S∈[1/κA, 1], and the E is a deviation amplitude.


At step 2: m+1 points in the preset interval are selected, and a solution of a system of linear equations fm(xi)+(−1)iE=1 is acquired, where the m+1 points are respectively x1, x2, . . . , xm+1, and satisfy x1=1/κ, xm+1=1.


At step 3: the acquired solution of the system of linear equations is substituted into fm(x) to obtain a set N′ of points corresponding to local maximum values of |1−fm(x)|.


At step 4: it is judged whether an absolute value of fm(x)−1 is the same for any x∈N′, and whether a sign of fm(x)−1 changes between plug and minus alternately as x increases, and if the absolute value of fm(x)−1 is the same for any x∈N′, and the sign of fm(x)−1 changes between plug and minus alternately as x increases, a current fm(x) is determined as an optimal polynomial.


At step 5: if the absolute value of fm(x)−1 is not the same for any x∈N′, or the sign of fm(x)−1 doesn't change between plug and minus alternately as x increases, elements in x1, x2, . . . , xm+1 are replaced with elements in the set N′ according to a rule including: replacing x2, . . . , xm with elements in the set N′ while x1, xm+1 remain unchanged, and then, the step of acquiring a solution of the system of linear equations fm(xi)+(−1)iE=1 is performed again, until an optimal fm(x) satisfying the judgment condition is obtained.


Specifically, the problem of calculating the polynomial p(A) of the first matrix A used for linear system preprocessing to minimize κA′ may be converted into a minimax approximation problem, which can be then solved by a Remez algorithm to obtain the optimal p(A).


First, it may be defined that: ∥l(x)−p(x)x∥S=maxx∈S|l(x)−p(x)x|, where S=[−1, −1/κ]∪[1/κ, 1], and l(x) is defined as:







l

(
x
)

=

{





1
,

x
>
0








-
1

,

x
>
0





.






Considering that p(x) is an even function, l(x)−p(x)x is an odd function. Therefore, only an interval S=[1/κ, 1] of x>0 is considered, and the optimal p(x) satisfies:








min

p


π
even

2

m








1
-


p

(
x
)


x




S


,




where πeven2m represents a set of 2 m order even functions.


The Remez algorithm is an iterative algorithm for approximating a function through a minimax norm. Assuming that p(x) may be expressed as: pm(x)=Σi=0m+1αix2i+1, it is defined that: fm(x)=pm(x)x=Σi=0m+1αix2i+1. First, m+1 points X={x1, x2, . . . , xm+1} are selected in the interval S∈[1/κ, 1], where x1=1/κ, and xm+1=1. Second, the linear equation system fm(xi)+(−1)iE=1, i=0,1,2, . . . , m+1 is solved to obtain a solution α1, α2, . . . , αm, E. Then, α1, α2, . . . , αm is substituted into fm(x) to obtain a set N′ of points corresponding to local maximum values of g(x)=|1−fm(x)|. Finally, it is judged whether an absolute value of fm(x)−1 is the same for any x∈N′, and whether a sign of fm(x)−1 changes between plug and minus alternately as x increases. If the absolute value of fm(x)−1 is the same for any x∈N′, and the sign of fm(x)−1 changes between plug and minus alternately as x increases, a current fm(x) is determined as an optimal polynomial. If the absolute value of fm(x)−1 is not the same for any x∈N′, or the sign of fm(x)−1 does not change between plug and minus alternately as x increases, elements in X are replaced with elements in the set N′ according to a rule including: x2, . . . , xm are replaced with elements in the set N′ while x1, xm+1 remain unchanged. Then, the step of acquiring a solution of the system of linear equations pm (xi)+(−1)iE=1, i=1, 2, . . . , m+1 is performed again, until an optimal pm(x) meeting the judgment condition is obtained.


Further, according to the above polynomial p(A), the linear system is preprocessed to obtain a second matrix A′ and a second vector b′. Specifically, according to the obtained optimal polynomial p(A), a new linear system A′x=b′ to be processed is constructed, including obtaining a second matrix A′=p(A) A and a second vector b′=p(A)b, so that the new linear system to be processed can be solved more efficiently. In general, a system of linear equations with a smaller condition number is easier to solve. Therefore, p(A) has to make the condition number κA′ of the preprocessed second matrix A′ satisfy κA′<<κA. In one embodiment, according to the above polynomial p(A), preprocessing the linear system to obtain the second matrix A′ and the second vector b′ may further include: constructing an operator P based on the polynomial P(A), and applying the operator P to the quantum state |bcustom-character to obtain a quantum state |b′custom-character, where the quantum state |b′custom-character is a quantum state corresponding to the second vector b′, and multiplying the polynomial P(A) by the first matrix A to obtain the second matrix A′. In an implementation scenario, the operator P corresponding to the polynomial P(A) may be prepared through quantum signal processing (QSP).


It should be noted that p(A) may be selected in various manners, such as a Neumann polynomial, a Chebyshev polynomial, a Least Square polynomial, or the like, which are not limited herein.


Based on the obtained second matrix A′ and second vector b′, a quantum circuit corresponding to an HHL algorithm may be constructed based on the second matrix A′ and the second vector b′, and used to perform quantum state evolution and measurement operations to solve the target system of linear equations.


The problem of solving a system of linear equations is a basic problem in many fields. In recent years, algorithms for solving the system of linear equations Ax=b have been continuously proposed. Harrow et al. proposed a quantum algorithm, i.e., an HHL algorithm, for solving a system of linear equations in 2009, which has a query complexity







𝒪

(



s


κ
2


ϵ



poly

(

log

(


/
ϵ

)

)


)

,




where s represents a sparsity of the first matrix A, and κ=∥A∥∥A−1∥ represents a condition number of the first matrix A.


Specifically, due to an exponential acceleration effect compared with the classical algorithm under certain conditions, the HHL algorithm can be widely applied in data processing, machine learning, numerical calculation, fluid mechanics problem processing and other scenarios in the future. In this embodiment, the HHL algorithm may be used to solve a linear equation, including: inputting, for the quantum circuit, an initial quantum state |b′custom-character and an operator U, and performing quantum state evolution and measurement operations to obtain a final quantum state of the evolved quantum circuit, where the operator U=eiA′t and the t is a constant.


The initial quantum state |b′custom-character and the operator U are input, and a solution x of the linear equation is output, which satisfies A′x=b′, i.e., x=A′−1‘b’. Therefore, the second matrix A′ should be an invertible matrix, and for the second vector b′, since data of the second vector b′ needs to loaded to the quantum circuit as described below, a dimension number of the second vector b′ should be able to be expressed as a positive integer power of 2. If the dimension number cannot be expressed in the form of a positive integer power of 2, zeros are filled into elements of the second vector b′ until the dimension number can be expressed in the form of a positive integer power of 2.


Specifically, constructing the quantum circuit corresponding to the HHL algorithm includes the following steps:


At step 1: several qubits including an auxiliary qubit, a first qubit, and a second qubit are obtained, where initial states of the auxiliary qubit and the first qubit are set to |0custom-character, an initial state of the second qubit is set to |b′custom-characteri=1N b′i|icustom-character, the b′ is an ith element of the second vector b′, and the N is a dimension number of the second vector.


Specifically, the number of the qubits to obtain may be determined by a user according to the needs, or a relatively sufficient number of qubits may be obtained to satisfy the computing needs when the computing resources are sufficient.


The obtained several qubits including the auxiliary qubit, the first qubit, and the second qubit may be specifically represented by qubits. For example, |0custom-character on the initial qubit indicates that the quantum state of the qubit is |0custom-character, and |1custom-character indicates that the initial quantum state is |1custom-character.


It should be noted that, for the convenience of subsequent distinctions, the obtained several qubits may be divided into auxiliary qubits, first qubits, and second qubits. The specific distinction names are not limited herein, and the initial state of each qubit may be prepared by an existing amplitude coding method or quantum state coding method. The initial states of the auxiliary qubit and the first qubit is set to |0), and the initial state of the second qubit is set to |b′custom-characteri=1Nb′i|icustom-character. For example, for a 4-dimensional second vector b′=[b′1, b′2, b′3, b′4], i.e., N=4, data of the second vector b′ is coded to the quantum state amplitude to obtain:



























|

b





=


1




b


1
2

+


b


2
2

+


b


3
2

+


b


4
2











i
=
1

4



b
i

|
i




=


1




b


1
2

+


b


2
2

+


b


3
2

+


b


4
2







b


1

|
0




+



b


2

|

1





+



b


3

|

2





+



b


4

|

3





)

=


1




b


1
2

+


b


2
2

+


b


3
2

+


b


4
2







b


1

|
00




+



b


2

|

01





+



b


3

|

10





+



b


4

|

11





)




thereby loading the data of the second vector b′ to the quantum state amplitudes of two second qubits in the quantum circuit.


At step 2: a unitary matrix U corresponding to the second matrix A′ is determined.


Specifically, since the second matrix A′ is a Hermitian matrix, the unitary matrix U corresponding to the second matrix A′ can be directly determined, or Hamiltonian simulation may be used to complete the transformation from Hermitian matrix to unitary matrix, thereby obtaining the corresponding unitary matrix U=eiA′t; where t is a constant and takes a general value 2π. First, the operator U=eiA′t is expanded according to the Jacobi-Anger expansion to determine an exponential expansion of the operator U, where the exponential expansion of the operator Ü is:








e

i

γ

t


=



J
0

(
t
)

+

2







k
=
1






(

-
1

)

k




J

2

k


(
t
)




T

2

k


(
γ
)


+

2

i







k
=
0






(

-
1

)

k




J


2

k

+
1


(
t
)




T

2

k


(
γ
)




,




where Jk is a kth order Bessel function of the first kind, and Tk is a kth order Chebyshev polynomial of the first kind.


The exponential expansion of the operator U is saved in accuracy of a qth order, where the exponential expansion of the operator U in accuracy of the qth order is:











e

i

γ

t


=



f
1

(

γ
,
t

)

+


f
2

(

γ
,
t

)



,





where








f
1



(

γ
,
t

)


=



J
0

(
t
)

+

2







k
=
1

q




(

-
1

)

k




J

2

k


(
t
)




T

2

k


(
γ
)




,









f
2



(

γ
,
t

)


=

2

i







k
=
0

q




(

-
1

)

k




J


2

k

+
1


(
t
)




T

2

k


(
γ
)



,
and






γ
=


P

(
A
)



A
.









Second, an operator f1 corresponding to f1(γ, t) and an operator f2 corresponding to f2(γ, t) are prepared through quantum signal processing (QSP). Finally, the operator U is constructed through linear combination of unitary operators (LCU), the operator f1 and the operator f2.


When k is a positive integer, Jk(−x)=(−1)kJk(x), Tk(−x)=(−1)kTk(x), and therefore f1(γ, t) is an even function, and f2(γ, t) is an odd function. Therefore, the operators f1, f2 can also be prepared through QSP. After the operators f1, f2 are prepared, an operator U may be constructed through LCU, as shown in FIG. 6. FIG. 6 is a quantum circuit diagram for constructing an operator U according to an embodiment of the present disclosure.


A success rate of constructing U in this way is n, and after amplitude estimation and amplitude amplification, the success rate can be amplified to η≈1. In the subsequent QPE process, O(−log ϵ) controlled U gates will appear, and a success rate of the whole process is







ζ
=

1

ϵ
δ



,




where δ=log η. Therefore, when η approaches 1, the whole algorithm process approaches 1.


At step 3: a first sub-quantum circuit module for phase estimation is constructed and used to decompose the |b′custom-character into |b′custom-characterj=1Nβjj) in eigenspace of the second matrix A′, and the initial state |0custom-character|b′custom-character of the first qubit and the second qubit is converted into Σj=1Nβjjcustom-characterjcustom-character, where the |μjcustom-character is an eigenvector of the second matrix A′, the λj is an eigenvalue of the second matrix A′, and the βj is an eigenvector amplitude of the second matrix A′.


Specifically, a first sub-quantum circuit module for phase estimation is constructed and used to decompose the |b′custom-character into |b′custom-characterj=1Nβjjcustom-character in eigenspace of the second matrix A′, and the initial state |0custom-character|b′custom-character of the first qubit and the second qubit is converted into Σj=1Nβjj|μjcustom-character. It will be understood by those skilled in the art that phase estimation is an important application of quantum Fourier transform (QFT), reflected in that it is the basis of many quantum algorithms.



FIG. 7 is a schematic diagram of a first sub-quantum circuit module corresponding to phase estimation according to an embodiment. The first sub-quantum circuit module shown in FIG. 7 includes: an H gate operation module, a controlled U operator operation module, and a quantum inverse Fourier transform module (module FT in the figure), where the U operator is a unitary matrix U corresponding to the second matrix A′. After passing through the first sub-quantum circuit module, the quantum state of the auxiliary qubit (corresponding to the uppermost time line in FIG. 7) remain unchanged, the initial state |0custom-character of the first qubit (corresponding to the middle time line in FIG. 7) is converted to |λjcustom-character, the initial state |b′custom-character of the second qubit (corresponding to the bottom time line in FIG. 7) is decomposed as |b′custom-characterj=1Nβjjcustom-character.


In fact, the output λj is an estimated value, and the output accuracy of phase estimation can be improved by increasing the number of first qubits. Moreover, in practical applications, an auxiliary quantum register, a first quantum register, and a second quantum register may be provided to respectively store quantum states of the auxiliary qubit, the first qubit, and the second qubit.


At step 4: a second sub-quantum circuit module for performing a controlled rotation operation is constructed to use |λjcustom-character as a control bit to rotate the auxiliary qubit to obtain





















j
=
1

N



β
j

|

λ
j




|

µ
j






(



1
-


C
2



λ
j

2




|
0





+


C

λ
j


|
1




)

,




where the C is a normalization constant.


Specifically, a second sub-quantum circuit module for performing a controlled rotation operation is constructed to use |λjcustom-character as a control bit to rotate the auxiliary qubit to obtain





















j
=
1

N



β
j

|

λ
j




|

µ
j






(



1
-


C
2



λ
j

2




|
0





+


C

λ
j


|
1




)

.




Controlled rotation may also be called “extraction proportion”, because after the phase estimation, the first quantum register will store a series of eigenvalues λj (specifically stored in a ground state |λjcustom-character), while the input state, i.e., initial state |b′custom-character, stored in the second quantum register will be decomposed in the eigenspace of the second matrix A′. Then, through the controlled rotation, a λj value in the ground state will be extracted to the amplitude. The quantum state |0custom-character of the auxiliary qubit is converted to obtain

















j
=
1

N



(



1
-


C
2



λ
j

2




|
0





+


C

λ
j


|
1




)

,




and the quantum state of each qubit is converted from


























j
=
1

N



β
j

|

λ
j




|

µ
j




|
0








to











j
=
1

N



β
j

|

λ
j




|

µ
j






(



1
-


C
2



λ
j

2




|
0





+


C

λ
j


|
1




)




by the second sub-quantum circuit. In order to reduce resource occupation, the auxiliary qubit may be set to 1 bit, and C is a constant generally taking a value 1.


It should be noted that the quantum state |xcustom-character=A′−1|b′custom-character (more specifically, the quantum state close to |xcustom-character) can be obtained from the second quantum register, accompanied by a constant factor C. Later, instead of the process of simple measurement and screening of results, amplitude amplification may be used to increase the probability of success.


At step 5: a third sub-quantum circuit module for performing inverse phase estimation is constructed and used to reset |λjcustom-character to |0custom-character.


Specifically, those skilled in the art will understand that an inverse operation of phase estimation is a restoration process of phase estimation, or a transposed conjugate operation of phase estimation, with a purpose of resetting |λjcustom-character to |0custom-character, specifically converting the quantum state






























j
=
1

N



β
j

|

λ
j




|

µ
j






(



1
-


C
2



λ
j

2




|
0





+


C

λ
j


|
1




)



to
:







j
=
1

N



β
j

|
0



|

µ
j






(



1
-


C
2



λ
j

2




|
0





+


C

λ
j


|
1




)






At step 6: a measurement operation module for the auxiliary qubit is constructed, so that when the measured quantum state of the auxiliary qubit is |1custom-character:














|

x





=


1








j
=
1

N




(

C



λ
j


-
1




β
j


)

2











j
=
1

N


C



λ
j


-
1




β
j

|

µ
j





,


the








|

x








|
x




.




Specifically, the quantum measurement operation is applied to the auxiliary qubit, to measure the auxiliary qubit after the inverse operation of phase estimation. After measurement, the state of the auxiliary qubit will collapse to a definite state, where the probability of collapse to |0> is







1
-


C
2



λ
j

2



,




and the probability of collapse to







|
1

>

is




C

λ
j


.






When the measured quantum state of the auxiliary qubit is |1custom-character and C=1, a definite quantum state can be obtained:










|

x





=


1








j
=
1

N




(

C



λ
j


-
1




β
j


)

2











j
=
1

N


C



λ
j


-
1




β
j

|

µ
j





,




which, as can be seen, is a corresponding result of amplitude normalization of |xcustom-character=A′−1|b′custom-characterj=1Nβjλj−1jcustom-character. In practical applications, |xcustom-character may be obtained correspondingly according to the application scenario required by the user, or |x′custom-character may be directly used for subsequent scenario calculations.


At step 7: the first sub-quantum circuit module, the second sub-quantum circuit module, the third sub-quantum circuit module, and the quantum measurement operation module are sequentially formed into the quantum circuit corresponding to the HHL algorithm.


Specifically, as shown in FIG. 8, where a total quantum circuit diagram corresponding to an HHL algorithm according to an embodiment is shown, the first sub-quantum circuit module, the second sub-quantum circuit module, the third sub-quantum circuit module, and the quantum measurement operation module, based on a time sequence of execution, are sequentially formed into a complete quantum circuit, i.e., a total quantum circuit corresponding to the HHL algorithm.


In an optional implementation, a quantum state |b′custom-character may be obtained by preparing a quantum state |bcustom-character of the first vector b, constructing an operator P based on the polynomial p(A), and then applying the operator P to the quantum state |bcustom-character, where the quantum state |b′custom-character is a quantum state corresponding to the second vector b′. The quantum state |bcustom-characteri=1Nbi|icustom-character and the operator P is applied to the quantum state |bcustom-character to obtain a quantum state |b′custom-character, where the quantum state |b′custom-characteri=1Nb′i|i).


Specifically, the operator P corresponding to the polynomial p(A) is prepared through QSP, where the QSP ensures that a quantum operator mapped from a polynomial defined over a field of real numbers can be prepared if and only if the polynomial is an odd or even function (that is, the polynomial has pure odd or even orders).


It can be seen that an HHL algorithm for polynomial pre-processing is constructed by combining a polynomial preprocessor with a quantum linear solver, which optimizes the complexity of the HHL algorithm with an optimized multiple increasing linearly with the order of the polynomial preprocessing function, thereby improving a capacity of the HHL algorithm in solving ill-conditioned linear systems, and expanding the application range of the HHL algorithm. As the order m increases, the HHL algorithm transitions to a quantum linear solver, and the complexity of the algorithm gradually changes from κ2 to κ.


Compared with the existing art, the present disclosure firstly acquires information of a first matrix A and a first vector b in a linear system to be processed, calculates a polynomial p(A) of the first matrix A used for preprocessing the linear system, pre-processes the linear system according to the polynomial p(A) to obtain a second matrix A′ and a second vector b′, constructs a quantum circuit corresponding to an HHL algorithm, performs quantum state evolution and measurement operations according to the second matrix A′ and the second vector b′, to obtain a final quantum state of the evolved quantum circuit, thereby reducing the time, complexity and computation amount in solving linear problems and speeding up solution of the quantum linear algorithm while reducing occupation of hardware resources.


In some embodiments, there is also a case where the condition number κA of the matrix A is greater than a preset threshold. In this case, a solution of the target system of linear equations may be obtained by constructing a new system of linear equations and then solving the new system of linear equations. Specifically, when the condition number κA of the matrix A is greater than a preset threshold, the first matrix A and the first vector b are processed through a polynomial preprocessor to obtain a third matrix A″ and a third vector b″, where a condition number κA′ of the third matrix A″ is less than the condition number κA of the first matrix A.


In numerical analysis and linear algebra the condition number κA of the first matrix A is defined as:







κ
A

=





"\[LeftBracketingBar]"





"\[LeftBracketingBar]"

A


"\[RightBracketingBar]"





"\[RightBracketingBar]"








"\[LeftBracketingBar]"






"\[LeftBracketingBar]"



A

-
1




"\[RightBracketingBar]"





"\[RightBracketingBar]"








∥A∥, ∥A−1∥ are norms of the first matrix A and an inverse matrix of the first matrix A, respectively.


The condition number of the matrix describes stretching and compressing capabilities of the matrix to vectors. The larger the condition number of the matrix, the higher the complexity in solving the linear system problem, resulting in disappearance of the acceleration effect of using quantum to solve the linear system problem. Therefore, to achieve a solution performance for linear system problems far exceeds that of an ordinary quantum solver, a preprocessor is used to process the first matrix A and the first vector b to reduce the condition number of the matrix. The preprocessor used in the embodiments of the present disclosure is a polynomial preprocessor.


The preset threshold may be 1, 10, 100, 1000, 10000, or any other value, which is not limited herein.


Further, a new system of linear equations A″x=b″ is constructed based on the first matrix A and the first vector b to solve the unknown x.


After processing the matrix A and the vector b through the polynomial preprocessor to obtain the third matrix A″ and the third vector b″, where the condition number κA′ of the third matrix A″ is less than the condition number κA of the matrix A, an ordinary quantum solver may be used to solve the unknown x in the new system of linear equations A″x=b″ constructed from the third matrix A″ and third vector b″.


As mentioned above, the present disclosure provides a polynomial preprocessor to process the matrix A and the vector b in the original system of linear equations Ax=b to obtain a new system of linear equations A″x=b″, where the condition number κA′ of the matrix A′ in the new system of linear equations is less than the condition number κA of the matrix A in the original system of linear equations, so that the condition number of the input matrix, and thus the complexity of the linear system problem, are reduced. Then, the new system of linear equations is solved to obtain a common unknown x, so that an acceleration effect in solving linear system problems with quantum is achieved.


In a specific embodiment of the present disclosure, in the aspect of processing the first matrix A and the first vector b through the polynomial preprocessor to obtain the third matrix A″ and the third vector b″, the method includes the following steps.


At step 1: a quantum state |bcustom-character of the first vector b is prepared, and a polynomial P(A) is determined based on a polynomial function P(y) and the first matrix A, where the y is an independent variable.


At step 2: an operator P is constructed based on the polynomial P(A), and the polynomial P(A) is multiplied by the matrix A to obtain the third matrix A″.


At step 3: the operator P is applied to the quantum state |bcustom-character to obtain a quantum state |b′custom-character, where the quantum state |b′custom-character is a quantum state corresponding to a third vector b″.


The quantum state










|
b



=


1


b










i
=
0


n
-
1




b
i

|
i




,




and we n is a dimension of the first vector b. The operator P is applied to the quantum state (b) to obtain a quantum state |b′custom-character, where the quantum state |b′custom-characteri=0n−1βiicustom-character.


Further, before determining the polynomial P(A) based on the polynomial function P(y) and the matrix A, the method further includes the following steps.


At step 1: an approximate function Km(y) containing parameters is acquired, and a domain of definition of the approximate function Km(y) containing parameters is determined.


At step 2: m+2 approximate deviation points are selected from the domain of definition, and the m+2 approximate deviation points are substituted into a relational expression composed of the approximate function Km(y) containing parameters and a deviation amplitude E, to obtain a m+2-dimensional system of linear equations Km(yk)−T=(−1)k+1E, where k=0,1,2 . . . m+1, the m is an integer greater than 0, and the Tis any natural number in the domain of definition.


At step 3: the system of linear equations Km(yk)−T=(−1)k+1E is solved to obtain a value of a parameter in the approximate function Km(y) containing parameters.


At step 4: a target approximate function is determined based on the value of the parameter.


At step 5: the polynomial function P(y) is determined based on the target approximate function.


If the approximate function Km(y) containing parameters is an even function, the Km(y)=Σi=0mθ2iy2i+1; and if the approximate function Km(y) containing parameters is an odd function, the Km(y)=Σi=0mθ2i+1y2i+2. The θ2i, θ2i+1 are parameters, and the i is an integer greater than or equal to 0. Here, the approximate function Km(y) containing parameters being an even function is taken as an example.


The approximate function Km(y) containing parameters has a domain of definition S=[−|λA|max, −|λA|min]∪[|λA|min, |λA]max], where |λA|max is a maximum absolute value of eigenvalues of the first matrix A, and |λA|min is a minimum absolute value of eigenvalues of the first matrix A. The m+2 approximate deviation points have a size relationship: y0<y1< . . . <ym<ym+1, which may be an arithmetic progression or not and is not limited herein.


The m+2-dimensional system of linear equations Km(yk)−T=(−1)k+1E may be solved to obtain m+2 unknowns E, θ2i, i=0,1,2 . . . m.


The relational expression composed of the approximate function Km(y) containing parameters and the deviation amplitude E may be Km(yk)−T=(−1)k+1E or Km(yk)−T=(−1)kE, which lead to E having opposite signs.


Specifically, in terms of determining the target approximate function based on the value of the parameter, the method includes the following steps.


At step 1: the value of the parameter is substituted into the approximate function Km(y) containing parameters to obtain an initial approximate function.


At step 2: an extreme point of an absolute value of a difference between the initial approximation function and the T is determined.


At step 3: if the extreme point and the m+2 approximation deviation points are equal within an accuracy requirement, the initial approximate function is determined as the target approximate function.


At step 4: if the extreme point and the m+2 approximation deviation points are not equal within the accuracy requirement, the extreme point is taken as a new point of the m+2 approximate deviation points, and the step of substituting the m+2 approximate deviation points into the relational expression composed of the approximate function Km(y) containing parameters and the deviation amplitude E is performed to obtain the m+2-dimensional system of linear equations Km(yk)−T=(−1)k+1E.


The accuracy requirement may include, for example, the same three decimal places. If the extreme point and the corresponding approximate deviation point have the same three decimal places, it is determined that the two are equal within the accuracy requirement. Apparently, the accuracy requirement is not limited to three decimal places, and may include the first digit, two decimal places, five decimal places, seven decimal places, or any other value, which is not limited herein.


Specifically, in the aspect of determining the polynomial function P(y) based on the target approximate function, the method includes: determining, based on a relational expression P(y)=K(y)/f(y) of the target approximation function and the polynomial function P(y), the polynomial function P(y).


Specifically, in the aspect of constructing the operator P based on the polynomial P(A), the method includes: preparing an operator P corresponding to the polynomial P(A) through quantum signal processing (QSP).


The QSP ensures that a quantum operator mapped from a polynomial defined over a field of real numbers can be prepared if and only if the polynomial is an odd or even function (that is, the polynomial has pure odd or even orders).


In a specific embodiment of the present disclosure, solving the unknown x through the system of linear equations A″x=b″ constructed from the third matrix A″ and third vector b″ includes the following steps.


At step 1: an operator U for the HHL algorithm is determined based on the third matrix A″.


At step 2: the quantum state |b′) corresponding to the third vector b″ and the operator U are input to the HHL algorithm, to determine a target quantum state including a value of the unknown x.


At step 3: a solution result of the unknown x is determined based on the target quantum state.


Specifically, in the aspect of determining the operator U for the HHL algorithm based on the third matrix A″, the method includes the following steps.


At step 1: the operator U=eiA′t for HHL is determined.


At step 2: the operator U=eiA′t is expanded according to the Jacobi-Anger expansion to determine an exponential expansion of the operator U, where the exponential expansion of the operator U is:








e

i

γ

t


=



J
0

(
t
)

+

2







k
=
1






(

-
1

)

κ




J

2

k


(
t
)




T

2

k


(
γ
)


+

2

i







k
=
0






(

-
1

)

k




J


2

k

+
1


(
t
)




T

2

k


(
γ
)




,




where Jk is a kth order Bessel function of the first kind, and Tk is a kth order Chebyshev polynomial of the first kind.


At step 3: the exponential expansion of the operator U is saved in accuracy of a qth order, where the exponential expansion of the operator U in accuracy of the qth order is:








e

i

γ

t


=



f
1

(

γ
,

t

)

+


f
2

(

γ
,
t

)



,




where












f
1

(

γ
,

t

)

=



J
0

(
t
)

+

2







k
=
1

q




(

-
1

)

k




J

2

k


(
t
)



T

2

k




(
γ
)




,









f
2



(

γ
,

t

)


=

2

i







k
=
0

q




(

-
1

)

k



J


2

k

+
1




(
t
)




T

2

k


(
y
)



,

and






γ
=

P


(
A
)



A
.









At step 4: an operator f1 corresponding to f1(γ,t) and an operator f2 corresponding to f2(γ, t) are prepared through the QSP.


At step 5: the operator U is constructed through linear combination of unitary operators (LCU), the operator f1 and the operator f2.


When k is a positive integer, Jk(−x)=(−1)kJk(x), Tk(−x)=(−1)kTk(x), and therefore f1(γ,t) is an even function, and f2(γ,t) is an odd function. Therefore, operators f1, f2 can also be prepared through QSP. After the operators f1, f2 are prepared, the operator U may be constructed through LCU, as shown in FIG. 6.


A success rate of constructing U in this way is η, and after amplitude estimation and amplitude amplification, the success rate can be amplified to n=1. In the subsequent QPE process, 0(−log ϵ) controlled U gates will appear, and a success rate of the whole process is







ζ
=

1

ϵ
δ



,




where δ=log η. Therefore, when η approaches 1, the whole algorithm process approaches 1.


Specifically, in the aspect of inputting the quantum state |b′) corresponding to the vector b′ and the operator U to the HHL algorithm to determine the target quantum state including a value of the unknown x, the method includes the following steps.


At step 1: a quantum state |b′) is prepared.


At step 2: the initial quantum state is evolved to a target quantum state containing a value of the unknown x through quantum phase estimation (QPE) and controlled rotation operations, where a quantum logic gate required in the QPE is the U. Specifically, the HHL algorithm needs three registers, each containing at least one qubit. Before executing the HHL algorithm, the three registers need to be initialized. In other words, the qubits in the three registers are all initialized to |0custom-character. Then, the |0custom-character in the third register are converted into a quantum state |bcustom-character of the vector b in a general method. Then, a step 2023 of applying the operator P to the quantum state |bcustom-character to obtain a quantum state |b′custom-character, where the quantum state |b′custom-characteri=0n−1βi|uicustom-character, is performed.


In the HHL algorithm, the QPE specifically includes: sequentially acting on the second register through an H gate, acting on the second register and the third register through a controlled U gate, and acting on the second register through a FT gate, where the FT gate is a quantum logic gate corresponding to the Fourier transform. Through the QPE, the quantum state Σi=0n−1βi|uicustom-character is converted into the quantum state |ψ1custom-characteri=0n−1βi|custom-character|uicustom-character, where custom-character is an estimated value of the eigenvalue λi obtained by QPE.


The controlled operation specifically includes applying a controlled R gate to the first register and the second register, where the second register is a control bit and the first register is a controlled bit, and performing controlled rotation on the quantum state of the auxiliary bit through a sub-quantum state |custom-charactercustom-character of the quantum state |ψ1custom-characteri=0n−1βi|custom-character|uicustom-character, to obtain a quantum state



















|


ψ
2




=







i
=
0


n
-
1




β
i

|




|

u
i






(



1
-


C
2


2




|
0





+


C

|
1




)

.




The first register is measured, and when the quantum state of the first register is |1custom-character, a quantum state


















"\[LeftBracketingBar]"


ψ
3




=







i
=
0


n
-
1





C


β
i




λ
~

ι







"\[RightBracketingBar]"






λ
~

ι








"\[LeftBracketingBar]"


u
i








"\[RightBracketingBar]"




1






is obtained.


The quantum state


















"\[LeftBracketingBar]"


ψ
3




=







i
=
0


n
-
1





C


β
i




λ
~

ι







"\[RightBracketingBar]"






λ
~

ι








"\[LeftBracketingBar]"


u
i








"\[RightBracketingBar]"




1






is evolved to a target quantum state











"\[LeftBracketingBar]"


ψ
4




=







i
=
0


n
-
1







β
i




λ
~

ι






"\[LeftBracketingBar]"


u
i









containing a value of the unknown x through an inverse operation of QPE.


Specifically, a solution result of the unknown x is determined based on the target quantum state, that is, the unknown x takes a value












β
i




λ
~

ι






"\[LeftBracketingBar]"


u
i





.




After practical verification, an optimization multiple of the condition number κA of the matrix A achieved by the method for solving a system of linear equations including a polynomial preprocessor provided in the embodiments of the present disclosure is χ=









κ
A


κ

A





l

,




which means that the optimization multiple of the condition number has a linear relationship with the order of the input polynomial Pm. The overall solution complexity is O(ls2κA′2 log N), which is optimized by







ω
_

=




κ
A

2


l



κ

A



2




l





compared with the solution complexity O(s2κA2 log N) of HHL. It can be seen that the degree of optimization increases linearly with l. When l→O(κA), κA′→1, the solution complexity of the method for solving a system of linear equations including a polynomial preprocessor provided in the embodiments of the present disclosure is reduced to the minimum O(s2κA log N).



FIG. 9 is a schematic structural diagram of an apparatus for solving a system of nonlinear equations on the basis of a quantum circuit according to an embodiment. As shown in FIG. 9, the apparatus includes:

    • an acquisition module 901 configured to acquire a target system of nonlinear equations to be solved;
    • a conversion module 902 configured to convert the target system of nonlinear equations to obtain a target system of linear equations;
    • a construction module 903 configured to construct a quantum circuit corresponding to a quantum linear solver used for solving the target system of linear equations;
    • a determination module 904 configured to perform, based on the quantum circuit corresponding to the quantum linear solver, quantum state evolution and measurement on the target system of linear equations, to solve the target system of linear equations, and determine, based on an obtained solution of the target system of linear equations, a solution of the target system of nonlinear equations to be solved. The specific details may be referred to the above description in conjunction with FIG. 2, and will not be repeated here.


An embodiment of the present disclosure further provides a storage medium having a computer program stored thereon, where the computer program is configured to, when executed, cause the method for solving a system of nonlinear equations on the basis of a quantum circuit described above with reference to FIG. 2 to be implemented.


In addition, an embodiment of the present disclosure further provides an electronic device, including a memory and a processor, where the memory has a computer program stored thereon, and the processor is configured to execute the computer program to implement the method for solving a system of nonlinear equations on the basis of a quantum circuit described above with reference to FIG. 2.


As used in the description and the claims, the term “if” may be interpreted contextually as “when . . . ” or “once” or “in response to determining” or “in response to detecting.” Similarly, the phrase “if it is determined” or “if [the described condition or event] is detected” may be interpreted contextually as meaning “upon determining” or “in response to determining” or “upon detecting [the described condition or event]” or “in response to detecting [the described condition or event].”


Although the implementations of the present disclosure have been described above, they are merely the embodiments used for facilitating understanding the present disclosure, and are not intended to limit the scope and application scenarios of the present disclosure. Any skilled person in the art of the present disclosure can make any modification and variation in implementation forms and details without departing from the spirit and scope revealed in the present disclosure, but the patent protection scope of the present disclosure shall still be subject to the scope defined in the attached claims.

Claims
  • 1. A method for solving a system of nonlinear equations on the basis of a quantum circuit, the method comprising: acquiring a target system of nonlinear equations to be solved;converting the target system of nonlinear equations to be solved to obtain a target system of linear equations;constructing a quantum circuit corresponding to a quantum linear solver used for solving the target system of linear equations;performing, based on the quantum circuit corresponding to the quantum linear solver, quantum state evolution and measurement on the target system of linear equations, to solve the target system of linear equations; anddetermining, based on an obtained solution of the target system of linear equations, a solution of the target system of nonlinear equations to be solved.
  • 2-13. (canceled)
  • 14. The method of claim 1, further comprising: embedding the target system of linear equations into a linear system, wherein the linear system includes a first matrix A and a first vector b, and the linear system includes a condition number κA;calculating a polynomial p(A) of the first matrix A in the linear system; andpre-processing the linear system according to the polynomial p(A) to obtain a second matrix A′ and a second vector b′.
  • 15. The method of claim 14, wherein the pre-processing of the linear system includes: the second matrix A′=p(A)A, and the second vector b′=p(A)b.
  • 16. The method of claim 14, wherein the pre-processing of the linear system includes: constructing an operator P based on the polynomial P(A), and applying the operator P to the quantum state |b) to obtain a quantum state |b′, wherein the quantum state |b′ is a quantum state corresponding to the second vector b′; andmultiplying the polynomial P(A) by the first matrix A to obtain the second matrix A′.
  • 17. The method of claim 16, wherein the constructing of the operator P based on the polynomial P(A) includes: preparing the operator P corresponding to the polynomial P(A) through quantum signal processing (QSP).
  • 18. The method of claim 14, wherein the constructing of the quantum circuit further includes constructing a quantum circuit corresponding to a Harrow-Hassidim-Lloyd algorithm (HHL) algorithm based on the second matrix A′ and the second vector b′; and the performing of the quantum state evolution and measurement further includes performing quantum state evolution and measurement operations by the quantum circuit corresponding to the HHL algorithm to solve the target system of linear equations.
  • 19. The method of claim 18, wherein the constructing of the quantum circuit includes: obtaining several qubits including an auxiliary qubit, a first qubit, and a second qubit, wherein initial states of the auxiliary qubit and the first qubit are set to |0, an initial state of the second qubit is set to |b′=Σi=1Nbi′|i, the bi′ is an ith element of the second vector b′, and the N is a dimension number of the second vector;determining a unitary matrix U corresponding to the second matrix A′;constructing a first sub-quantum circuit module for phase estimation, which is used to decompose the |b′ into |b′=Σj=1Nβj|μj in eigenspace of the second matrix A′, and transforming the initial state |0|b′ of the first qubit and the second qubit into Σj=1Nβj|λj|μj, wherein the |μj is an eigenvector of the second matrix A′, the λj is an eigenvalue of the second matrix A′, and the βj is an eigenvector amplitude of the second matrix A′;constructing a second sub-quantum circuit module for performing a controlled rotation operation, which uses |λj as a control bit to rotate the auxiliary qubit to obtain
  • 20. The method of claim 18, wherein the performing of the quantum state evolution and measurement operations includes: inputting the initial quantum state |b′ and the unitary matrix U to the quantum circuit corresponding to the HHL algorithm, and performing quantum state evolution and measurement operations by the quantum circuit corresponding to the HHL algorithm to obtain a final quantum state of the evolved quantum circuit; andsolving the target system of linear equations based on the final quantum state.
  • 21. The method of claim 18, wherein the performing of quantum state evolution and measurement operations includes: in response to the condition number κA greater than a preset threshold, constructing a third matrix A″ and a third vector b″ based on the linear system;determining an operator U for the HHL algorithm based on the third matrix A″;inputting the quantum state |b′) corresponding to the third vector b″ and the operator U to the HHL algorithm, to determine a target quantum state including a value of the unknown x; anddetermining a solution result of the unknown x based on the target quantum state, to solve the target system of linear equations.
  • 22. The method of claim 21, wherein the determining of the operator U for the HHL algorithm based on the third matrix A″ includes: determining the operator U=eiA′t for the HHL algorithm;expanding the operator U=eiA′t according to the Jacobi-Anger expansion to determine an exponential expansion of the operator U, wherein the exponential expansion of the operator U is:
  • 23. The method of claim 14, wherein the calculating of the polynomial p(A) of the first matrix A in the linear system includes: constructing a function ƒm(x) and a preset interval S for solving the polynomial, wherein the function ƒm(x) satisfies fm(x)=pm(x)x, and fm(x2)+(−1)iE=1, the m is a polynomial order, the i=0,1,2, . . . , m+1, the S∈[1/κA, 1], the E is a deviation amplitude, and the κA is a condition number;selecting m+1 points in the preset interval, and acquiring a solution of a system of intermediate linear equations fm(x2)+(−1)iE=1, wherein the m+1 points are respectively x1, x2, . . . , xm+1, and satisfy x1=1/κ, xm+1=1;substituting the acquired solution of the system of intermediate linear equations into fm(x) to obtain a set N′ of points corresponding to local maximum values of |1−fm(x)|;judging whether an absolute value of fm(x)−1 is the same for any x∈N′, and whether a sign of fm(x)−1 changes between plug and minus alternately as x increases; anddetermining, according to a judgment result, the polynomial p(A) of the first matrix A.
  • 24. The method of claim 23, wherein the determining of the polynomial p(A) of the first matrix A includes: in response to that the sign of fm(x)−1 changes between plug and minus alternately, determining a current fm(x) as an optimal polynomial and as the polynomial p(A) of the first matrix A; orin response to that the sign of fm(x)−1 does not change between plug and minus alternately, replacing elements in x1, x2, . . . , xm+1 with elements in the set N′ according to a rule including: replacing x2, . . . , xm with elements in the set N′ while x1, xm+1 remain unchanged, and then returning to perform the step of acquiring a solution of the system of intermediate linear equations fm(xi)+(−1)iE=1, until an optimal fm(x) satisfying the judgment condition is obtained, so as to determine the polynomial p(A) of the first matrix A.
  • 25. The method of claim 23, further comprising, before the calculating of the polynomial p(A) of the first matrix A in the linear system: acquiring an approximate function Km(y) containing parameters, and determining a domain of definition of the approximate function Km(y) containing parameters;selecting m+2 approximate deviation points from the domain of definition, and substituting the m+2 approximate deviation points into a relational expression composed of the approximate function Km(y) containing parameters and a deviation amplitude E, to obtain a m+2-dimensional system of intermediate linear equations Km(yk)−T=(−1)k+1E, where k=0,1,2 . . . m+1, the m is an integer greater than 0, and the Tis any natural number in the domain of definition;solving the system of intermediate linear equations Km(yk)−T=(−1)k+1E, to obtain a value of a parameter in the approximate function Km(y) containing parameters;determining a target approximate function based on the value of the parameter; anddetermining, based on the target approximate function, a polynomial function P(y) for determining the polynomial P(A).
  • 26. The method of claim 25, wherein the determining of the target approximate function based on the value of the parameter includes: substituting the value of the parameter into the approximate function Km(y) containing parameters to obtain an initial approximate function;determining an extreme point of an absolute value of a difference between the initial approximation function and the T;if the extreme point and the m+2 approximation deviation points are equal within an accuracy requirement, determining the initial approximate function as the target approximate function; orif the extreme point and the m+2 approximation deviation points are not equal within the accuracy requirement, taking the extreme point as a new point of the m+2 approximate deviation points, and performing the step of substituting the m+2 approximate deviation points into the relational expression composed of the approximate function Km(y) containing parameters and the deviation amplitude E to obtain the m+2-dimensional system of intermediate linear equations Km(γk)−T=(−1)k+1E.
  • 27. The method of claim 26, wherein the determining of the polynomial function P(y) includes: determining, based on a relational expression P(y)=K(y)/f(y) of the target approximation function and the polynomial function P(y), the polynomial function P(y), wherein the f(y)=y.
  • 28. The method of claim 22, wherein, if the approximate function Km(y) containing parameters is an even function, the Km(y)=Σi=0mθ2iy2i+1; and if the approximate function Km(y) containing parameters is an odd function, the Km(y)=Σi=0mθ2i+1y2i+2; wherein the θ2i, θ2i+1 are parameters, and the i is an integer greater than or equal to 0.
  • 29. (canceled)
  • 30. A non-transitory storage medium storing a computer program configured to, when executed, cause a method for solving a system of nonlinear equations on the basis of a quantum circuit to be implemented, the method comprising:acquiring a target system of nonlinear equations to be solved;converting the target system of nonlinear equations to be solved to obtain a target system of linear equations;constructing a quantum circuit corresponding to a quantum linear solver used for solving the target system of linear equations;performing based on the quantum circuit corresponding to the quantum linear solver, quantum state evolution and measurement on the target system of linear equations, to solve the target system of linear equations; anddetermining, based on an obtained solution of the target system of linear equations, a solution of the target system of nonlinear equations to be solved.
  • 31. An electronic device; comprising: a memory having a computer program stored thereon; anda processor configured to execute the computer program to implement a method for solving a system of nonlinear equations on the basis of a quantum circuit, the method comprising:acquiring a target system of nonlinear equations to be solved;converting the target system of nonlinear equations to be solved to obtain a target system of linear equations;constructing a quantum circuit corresponding to a quantum linear solver used for solving the target system of linear equations;performing, based on the quantum circuit corresponding to the quantum linear solver, quantum state evolution and measurement on the target system of linear equations, to solve the target system of linear equations; anddetermining, based on an obtained solution of the target system of linear equations, a solution of the target system of nonlinear equations to be solved.
  • 32. The non-transitory storage medium of claim 30, wherein the method further comprises: embedding the target system of linear equations into a linear system, wherein the linear system includes a first matrix A and a first vector b, and the linear system includes a condition number κA;calculating a polynomial p(A) of the first matrix A in the linear system; andpre-processing the linear system according to the polynomial p(A) to obtain a second matrix A′ and a second vector b′.
  • 33. The electronic device of claim 31, wherein the method further comprises: embedding the target system of linear equations into a linear system, wherein the linear system includes a first matrix A and a first vector b, and the linear system includes a condition number κA;calculating a polynomial p(A) of the first matrix A in the linear system; andpre-processing the linear system according to the polynomial p(A) to obtain a second matrix A′ and a second vector b′.
Priority Claims (4)
Number Date Country Kind
202111421916.0 Nov 2021 CN national
202111425744.4 Nov 2021 CN national
202111470416.6 Dec 2021 CN national
202111527379.8 Dec 2021 CN national
CROSS REFERENCE TO RELATED APPLICATIONS AND CLAIM OF PRIORITY

This application claims benefit under 35 U.S.C. 119, 120, 121, or 365 (c), and is a National Stage entry from International Application No. PCT/CN2022/134387, filed Nov. 25, 2022, which claims priority to the benefit of Chinese Patent Application Nos. 202111421916.0 filed on Nov. 26, 2021, 202111425744.4 filed on Nov. 26, 2021, 202111470416.6 filed on Dec. 3, 2021 and 202111527379.8 filed on Dec. 14, 2021 in the China Intellectual Property Office, the entire contents of which are incorporated herein by reference.

PCT Information
Filing Document Filing Date Country Kind
PCT/CN2022/134387 11/25/2022 WO