Quantum-Solving Algorithm for d’Alembert Solutions of the Wave Equation

When faced with a quantum-solving problem for partial differential equations, people usually transform such problems into Hamiltonian simulation problems or quantum-solving problems for linear equation systems. In this paper, we propose a third approach to solving partial differential equations that differs from the two approaches. By using the duality quantum algorithm, we construct a quantum-solving algorithm for solving the first-order wave equation, which represents a typical class of partial differential equations. Numerical results of the quantum circuit have high precision consistency with the theoretical d’Alembert solution. Then the routine is applied to the wave equation with either a dissipation or dispersion term. As shown by complexity analysis for all these cases of the wave equation, our algorithm has a quadratic acceleration for each iteration compared to the classical algorithm.


Introduction
Most scientific problems can be solved by studying the laws governing the evolution of physical quantities in space and time. Therefore, partial differential equations undoubtedly play an extremely important role in the field of natural sciences. However, the problem of solving partial differential equations is extremely difficult. While if quantum algorithms are introduced and the problems of partial differential equations are solved on a quantum computer, it can achieve accelerated characteristics compared to classical algorithms.
The usual quantum algorithm for solving partial differential equations proceeds as follows. First, discretize the space so that the function f (x, t) becomes a vector f(t) and map its normalized components to the quantum state components, i.e., |x(t) = ∑ i f i (t)|x i , where f i (t) is the i-th component of the vector f(t) after normalization. Next, the vectors encoded onto the quantum states are mapped onto a fixed model. Most quantum algorithms for solving partial differential equations rely on Hamiltonian simulations [1][2][3][4] or a linear equation system-solving algorithm (HHL algorithm) [5].
In the following, the main ideas of the above two solution methods will be briefly reviewed with examples. The solution method based on Hamiltonian simulation [6][7][8] that maps partial differential equations to the Schrödinger equation will be introduced first. This method maps equations with a similar structure to the Schrödinger equation to the Schrödinger equation and transforms the equation solving problem into a Hamiltonian simulation problem. For example, solving the Black-Scholes equation [9] ∂ f ∂t (1) The equation can be written in the following form It is obvious that the equation is formally similar to the Schrödinger equation. Thus the A operator can be mapped to the Hamiltonian in Schrödinger's equation in such a way that A = ibp + aI + cp 2 ,p = −i∂ x . One can split A into Hermitian and anti-Hermitian parts, i.e., A = A H + A aH , where The vector f(t), obtained by discretizing the function, is encoded onto the state vector |x( ) , using the Trotter product formula |x( ) = e A |x 0 ≈ e A H e A aH |x 0 .
The problem of solving the partial differential equation is transformed into the problem of a Hamiltonian simulation. The process of simulating the action of the above Hamiltonian, i.e., the quantum state |x 0 , evolves under the designed Hamiltonian to obtain the final state. The solution of the original equation at different moments can be obtained by measuring the final state for different iterations.
In fact, it is efficient to use a Hamiltonian simulation to construct quantum algorithms for solving partial differential equations, which can solve first-order partial differential equations (requiring that the Hermitian and anti-Hermitian parts of the matrix A decomposition commute with each other) and second-order partial differential equations such as the wave equation. However, not all partial differential equations have the algebraic structure of Schrödinger's equation. The quantum algorithm [7,[10][11][12] for solving partial differential equations is presented below using the HHL algorithm. For the partial differential equation with the following structure after spatial discretizatioṅ Using Euler's method to discretize time gives Let x j = x t j , the partial differential equation can be transformed into the following linear equation system; as an example, only the result of j ≤ 2 is given here, This system of equations is then solved using the HHL algorithm to obtain the following quantum states The quantum state contains the solution of the partial differential equation t 0 to t j moments. Therefore, are there any other efficient algorithms for solving partial differential equations besides the above two methods? In this article, a third method different from the above two, the duality quantum algorithm [13][14][15][16] with amplitude amplification [17][18][19][20], is used to construct a quantum algorithm for solving the partial differential equation.
The duality quantum algorithm also brings a speed-up effect compared to the classical algorithm [21][22][23].
This paper is organized as follows. First of all, the duality quantum algorithm will be used to solve the first-order wave equation with the d'Alembert solution. In the second and third parts of this paper, we will use the duality quantum algorithm to construct a solution algorithm to the wave equations with dissipation and dispersion terms. In these three parts of the paper, for these three problems, we use numerical simulations and present the results of the quantum algorithm solutions in the form of pictures for comparison with the theoretical values. At the end of this paper, we will analyze the complexity of our algorithm for solving the wave equations.

Duality Quantum Algorithm for Solving the First-Order Wave Equation
When talking about the wave equation, people must first think of the second-order linear hyperbolic type equation The general solution can be written as u(x, t) = f (x − kt) + g(x + kt), where f , g are two arbitrary functions. f (x − kt) and g(x + kt) represent waves passing along the x-axis with constant velocity to the right and to the left. Since Equation (9) is a linear homogeneous equations, its solutions are superposed. Therefore, f and g are two traveling waves that propagate independently without interfering with each other. If one focuses on only one of these two waves, Equation (9) degenerates to a linear hyperbolic equation of the first-order In this paper, we discretize the continuous independent variable x into N points, i.e., x = (x 0 , x 1 , · · · , x N−1 ). Then the spatial part of the function u x, t j at the moment t j is discretized into the vector u x, t j = u x 0 , t j , u x 1 , t j , · · · , u x N−1 , t j , encode it onto the computational basis and define the quantum state |ψ j as In the following, we will give the quantum algorithm for solving Equation (10) based on the non-unitary evolution of the quantum system. First, the Taylor expansion for each order partial differential term of Equation (10) is Pluging Equation (13) into Equation (10), the difference equation form of Equation (10) is obtained as Its local truncation error is o(τ + h). When τ, h → 0, Equation (14) approximates the original Equation (10). Organizing Equation (14) leads to Let ∆ = τ h , then the following iterative relation can be obtained from Equation (15).
Taking the periodic boundary condition that u(x N , t) = (x 0 , t), the equation describing the whole system can be written in the following form Then the state |ψ j+1 of the system at the next moment, i.e., the moment t j+1 , can be expressed as A|ψ j . It is obvious that the A-matrix is not an unitary matrix, so there is no way to achieve it directly by the product of quantum logic gates. Instead, the A-matrix has to be split into linear combinations of the unitary operators by the duality model of quantum computation, i.e., By introducing an auxiliary qubit, the operation of the linear combination of unitary operators can be realized and thus equivalently the non-unitary evolution, i.e., A|ψ j . Its quantum circuit is shown in Figure 1. Where the matrix A 1 can be decomposed into C n (X) gates as well as X gates with O(log 2 N). The specific quantum circuit that implements the A 1 operation is shown in According to Lemma 5.5 and Lemma 7.1 in article [24], a total of 4N − 5 log 2 N − 4 CNOT gates and single-qubit rotating gates are needed if the A 1 operation continues to be disassembled. The following will explain the duality quantum algorithm for the solution of the first-order wave equation according to Figure 1, where first the auxiliary qubit passes through the W 0 gate, which has the following effect

… …
Next, after two controlled quantum gates |0 0| ⊗ A 0 and |1 1| ⊗ A 1 , the quantum state evolves as Then the quantum state after the Hadamard transformation is Finally, after the measurement to select the state of the auxiliary qubit as 0, the state of the working qubits at this time is |ψ j+1 , that is Define the coefficients C j as The amplitude under the computational basis of the quantum state |ψ j+1 is enlarged by a factor of C j to obtain the column vector u(x i , t j+1 ), which is the state of the system at the moment t j+1 . The analysis yields that the computational complexity of this algorithm is O(N) per iteration, while the complexity of the classical algorithm is O N 2 . The specific calculation of the complexity is presented at the end of this paper.
The following equation will be used as an example (25) to show the duality quantum algorithm for the solution of the first-order wave equation. First of all, one period of the function, i.e., [0, 1], is chosen, and this interval is discretized into 32 points, i.e., h = 0.03125 in Equation (15). Thus, the function value of 32 discrete points can be encoded using 5 qubits. Choose ∆ = 0.1 in Equation (16). Then τ = 0.003125, which represents the time interval for each iteration of the system evolution. According to Equation (20), the effect of the action of W 0 can be determined as This gives W 0 as R y (−0.33). Thus far, the quantum circuit for each iteration of the solution to Equation (25) can be given in Figure 3. Numerical simulation of the first 10 iterations of this quantum circuit, whose results are shown in Figure 4. Where the orange curve represents the theoretical value. The blue points represent the results given by the numerical simulation of the quantum solution algorithm.

Duality Quantum Algorithm for the Solution of the Traveling Wave Dissipation Problem
By adding the dissipation term to Equation (10), the traveling wave equation with dissipation is obtained as The Taylor expansion for each term of Equation (27) is Differentiating Equation (27) yields The collation leads to Let ∆ = τ h and ∆ 1 = α ∆ h , then the following iterative relation can be obtained Take the periodic boundary condition that u(x N , t) = u(x 0 , t). Then, the equation describing the whole system can be written in the following form In Equation (32) in which Following the encoding method of Equation (12), then the state |ψ j+1 of the system at the next moment, i.e., the moment t j+1 , can be expressed as A|ψ j . The matrix A = (k∆ − 2∆ 1 Thus the operation of a linear combination of unitary operators can be equivalently implemented by introducing two auxiliary qubits, whose quantum circuit is shown in Figure 5. … … In the following, the duality quantum algorithm for solving the dissipation problem of the first-order wave equation is explained in conjunction with the quantum circuit ( Figure 5), where the first the auxiliary qubits pass through the W 1 gate, which has the following effect Next, after three controlled quantum gates |00 00| ⊗ A 0 , |01 01| ⊗ A 1 and |10 10| ⊗ A 2 the quantum state evolves as Then, the quantum state evolves after the Hadamard transformation of two auxiliary qubits as Finally, the auxiliary qubits are measured to select the state of the auxiliary qubit as 00; then the state of the working qubit at this time is |ψ j+1 , that is Define the coefficient C j as The column vector u(x i , t j+1 ), which is the state of the system at t j+1 moments, is obtained by amplifying the amplitude under the computational basis of the quantum state |ψ j+1 by a factor of C j . According to Lemma 5.5 and Lemma 7.1 in article [24] and combined with Figure 5, the computational complexity of this algorithm per iteration is O(N), while the complexity of the classical algorithm is O N 2 . Thus, the present algorithm has the property of speeding up in each iteration compared to the classical algorithm. The specific calculation of the complexity is presented at the end of this paper.
The following equation will be used as an example to show the duality quantum algorithm for the solution of the dissipation problem of the first-order wave equation. First, one period of the function is chosen, i.e., [0, 1], and this interval is discretized into 32 points, i.e., h = 0.03125 in Equation (28). Thus, the function value of 32 discrete points can be encoded using 5 qubits. Choose ∆ = 0.2 and ∆ 1 = 1.28 in Equation (31). Then τ = 0.00625, which represents the time interval for each iteration of system evolution. According to Equation (35), the effect of the action of W 1 can be determined as It is constructed as shown in Figure A1 with Equations (A1) and (A2). The revolving gate R n (θ) of the first auxiliary qubit is constructed according to Equation (A1), such that It is obtained that R n (θ) is R y (1.442). According to Equation (A3), the controlled operator U 1 is to achieve the following action The controlled operator U 1 can be obtained as R y (1.298). At this point, we can give the quantum circuit for each iteration of the solution Equation (40), as shown in Figure 6

Duality Quantum Algorithm for Solving Traveling Wave Dispersion Problems
By adding the dispersion term to Equation (10), the traveling wave equation with dispersion is obtained as ∂u ∂t Taylor expansion of the terms of Equation (44) ∂u ∂t Differentiating Equation (44) yields Let ∆ = τ h , ∆ 2 = β ∆ h 2 , and the following iterative relation can be obtained Take the periodic boundary condition that u(x N , t) = u(x 0 , t). Then the equation describing the whole system can be written in the following form where in which a = k∆ + 3∆ 2 It can be seen that A = (k∆ Therefore, the operation of the linear combination of unitary operators can be equivalently implemented by introducing two auxiliary qubits, whose quantum circuit is shown in Figure 8. It is not difficult to find A 3 = A 2 2 , and the following will be combined with the quantum circuit ( Figure 8) to explain the duality quantum algorithm for the solution of the dispersion problem of the first-order wave equation. First, the auxiliary qubits pass through the W 2 gate, the effect of which is as follows

… …
Next, after four controlled quantum gates, |00 00| ⊗ A 0 , |01 01| ⊗ A 1 , |10 10| ⊗ A 2 and |11 11| ⊗ A 3 , the quantum state evolves as Then, the Hadamard transform is performed for the two auxiliary qubits, and the quantum state evolves as Finally, the state of the auxiliary qubit is selected as 00 after measurement, then the state of the working qubits at this time is |ψ j+1 , that is Define the coefficient C j as The amplitude under the computational basis of the quantum state |ψ j+1 is enlarged C j times to obtain the column vector u x i , t j+1 , which is the state of the system at the moment t j+1 . According to Lemma 5.5 and Lemma 7.1 in article [24], and combined with the analysis of Figure 8, we can obtain that the computational complexity of this algorithm for each iteration is O(N), while the complexity of the classical algorithm is O N 2 . The specific calculation of the complexity is presented at the end of this paper.
The following equation will be used as an example to show the duality quantum algorithm for the solution of the dispersion problem of the first-order wave equation. Firstly, one period of the function is chosen, i.e., [0, 1], and this interval is discretized into 32 points, i.e., h = 0.03125 in Equation (47). Thus, the function value of 32 discrete points can be encoded using 5 qubits. Choose ∆ = 0.05 and ∆ 2 = 0.512 in Equation (47). Then τ = 0.0015625, which represents the time interval for each iteration of the system evolution. According to Equation (52), the effect of the action of W 2 can be determined as It is constructed as shown in Figure A1 with Equations (A1) and (A2). Construct the revolving gate R n (θ) of the first auxiliary qubit according to Equation (A1), such that It is obtained that R n (θ) is R y (1.078). According to Equation (A3), two controlled operators U 1 , U 2 are to be realized as follows The controlled operator U 1 can be obtained as R y (0.4563) and U 2 as R y (5.639). The quantum circuit for each iteration of the solution Equation (57) so far is given in Figure 9. And the specific quantum circuit of the R-operation in Figure 9 is shown in Figure 10. The result of numerically simulating the first 10 iterations of this quantum circuit is shown in Figure 11. The orange curve represents the theoretical value. The blue points represent the values solved by the numerical simulation quantum algorithm.

Discussion
For a quantum algorithm that solves a d-dimensional partial differential equation (meaning that there are d spatial variables) with a spatial discretization number of N, the output is an approximation C( f ) of the function f with an error . In fact, for the problem of quantum algorithms solving partial differential equations, the number of discrete points N and the error are interrelated [12,25]. The correlations are as follows For the preparation of the initial state, its complexity is O(poly log(N)). The complexity of each iteration of the algorithm in this paper will be given below. First, according to Lemma 5.5 and Lemma 7.1 in article [24], it can be obtained that the controlled gate C n−1 (U) for n qubits having n − 1 control qubits can be split into CNOT gates and single-qubit gates for a total of 2 n+1 − 5, where n ≥ 3. For the quantum circuit in Figure 1, the total number of elementary quantum gates required is For the quantum circuit in Figure 5, the total number of elementary quantum gates required is For the quantum circuit in Figure 8, the total number of elementary quantum gates required is It can be found that the quantum-solving algorithm given in this paper has a quadratic acceleration for each iteration compared to the classical algorithm. However, the state of the auxiliary qubits needs to be selected after measurement at the end of each iteration. This result is probabilistic, and the overall success rate of the algorithm decreases exponentially as the number of iterations increases if the selection is made after each iteration. Therefore, ensuring an overall higher success rate requires the use of the quantum search algorithm [17][18][19][20] to amplify the amplitude of the target state before measurement. Under ideal conditions of the device, it is proved that the Grover-Long algorithm can achieve a 100% success rate in all cases [17,26,27]. Thus, using the Grover-Long algorithm under ideal conditions to amplify the amplitude, it is possible to obtain a 100% success rate every time. If the complexity of each iteration step is considered synthetically, then the complexity of each iteration step is O(N) + O( √ M), where M is the dimension of the auxiliary qubits space.
In the future, our quantum algorithms are expected to be combined with finite element methods to solve complex practical problems, such as those related to fluid dynamics [28].

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A
The quantum circuit for the preparation of two-qubit arbitrary quantum states is as follows in Figure A1. First, the rotation operator R n (θ) is constructed, such that