Amplitude Amplification for Optimization via Subdivided Phase Oracle

We propose an algorithm using a modified variant of amplitude amplification to solve combinatorial optimization problems via the use of a subdivided phase oracle. Instead of dividing input states into two groups and shifting the phase equally for all states within the same group, the subdivided phase oracle changes the phase of each input state uniquely in proportion to their objective value. We provide visualization of how amplitudes change after each iteration of applying the subdivided phase oracle followed by conventional Grover diffusion in the complex plane. We then show via numerical simulation that for normal, skew normal, and exponential distribution of objective values, the algorithm can be used to amplify the probability of measuring the optimal solution to a significant degree independent of the search space size. In the case of skew normal and exponential distributions, this probability can be amplified to be close to unity, making our algorithm near deterministic. We then modify our algorithm in order to demonstrate how it can be extended to a broader set of objective value distributions. Finally, we discuss the speedup compared to classical schemes using the query complexity model, and show that our algorithm offers a significant advantage over these classical approaches.


I. INTRODUCTION
Grover's search algorithm [1] is one of the most celebrated quantum algorithms which demonstrated the computational power of quantum computers that classical computers cannot reach. The problem that Grover's algorithm solves can be explained as finding one of the t inputs among N total inputs that matches some certain constraints, f (x) = 1, given by an oracle function f (x). The assumption about the oracle is that nothing is known about its construction or whether the problem has any structure we can exploit to infer the solution from previously seen inputs and outputs. Grover's algorithm can accomplish this task in just O( N/t) queries, compared to the classical approach which requires O(N/t) queries as we need to exhaustively check all inputs.
Since the original algorithm proposed by Grover in 1996, there have been multiple advancements to the algorithm. It is shown that Grover's algorithm is a special case, and is optimal [2]- [4], of a more general approach, later known as "amplitude amplification" technique [5]- [9]. The use of a non-π phase shift in place of the usual π-flip of Grover oracle was investigated in [10]- [13] and later generalized into a technique to increase the success probability of Grover to unity [14], [15]. It also finds applications in solving the souffle problem, too few or too many iterations deteriorates the success probability, which is known as the fixed-point searching [16]- [19].
In recent years, there have been investigations on the viability of using amplitude amplification technique with a new variant of non-trivial phase shift oracle, one that does not simply recognize two classes of inputs but rather acts on each input uniquely depending on some cost functions [20]- [22]. In this paper, we aim to extend the use of this non-trivial phase shift oracle to solve combinatorial optimization problems, instead of the usual adaptive Grover's algorithm [23]- [25] where we try to find better solutions by running the algorithm repeatedly until we cannot find a better solution.
The contribution of this paper is to introduce the use of a non-trivial phase shift oracle which we refer to as a subdivided phase oracle (SPO) [20], replacing the canonical Grover oracle in the amplitude amplification process, to solve combinatorial optimization problems. We show that this approach can be used to efficiently find the solutions of a large class of optimization problems, characterized by a number of objective functions.
We describe the setting of maximization problems (section II) which we then discuss one approach to construct a subdivided phase oracle operator from the given objective function and how the objective values can be embedded into the phase of each basis state (section III). We define our algorithm and provide visualization of how amplitudes change after each iteration of applying the subdivided phase oracle followed by conventional Grover diffusion (inverse around the average). This modified amplitude amplification process is complex and non-trivial to mathematically derive closed form expressions. Instead, we give intuitions on how it can amplify the probability of measuring the solution states and discuss the relationship between the distribution of objective values, oracle construction, and performance of the algorithm (section IV). We then analyze the performance of the algorithm and show results via simulations of applying it to several classes of objective functions, plotting the probability of measuring the solutions with respect to the number of iterations (section V). We also discuss a dynamic approach where the oracle operator changes at each iteration depending on current amplitude values, although this approach is purely theoretical since it is impossible to keep track of all amplitudes in large space, it suggests that a heuristic or variational [26] approach may exist in selecting parameters to construct the oracle. We then propose a sampling scheme based on the modified amplitude amplification process, plot the trade-off between the number of sampling trials and iterations, and derive conditions to which the sampling scheme gives a better expected run time than classical sampling or exhaustive search (section VI). We end the paper by discussing possible ways of extending the ideas we present here, pointers to efficiently realize the subdivided phase oracle in circuit models, connections to QAOA [27]- [29] and pointers to make it variational, and utilizing the visualization we developed to create a more suitable diffusion operator to increase the performance of our algorithm (section VII).

II. OPTIMIZATION PROBLEM SETTING
The problem we consider is the combinatorial optimization problem. For example, "what is the shortest route starting from city A to city B given a map of roads currently available?". An optimization problem can be described by three elements; a set of all possible inputs X, an objective function f (x) defined over part of the input space, and a Boolean function g(x) which checks for the feasibility of the input. Using the above example, X refers to all possible routes connecting A and B, f (x) refers to the distance of each route, and since some roads might be closed due to constructions or accidents, g(x) refers to the state of whether a route can be used at the time of travel. The task is to find an input x ∈ X such that f (x) is the maximum or the minimum while g(x) = 1, arg max x∈X f (x) such that g(x) = 1. (1) In this study, we will focus on the combinatorial optimization problem, where the input is defined over the set of binary string x = {0, 1} n of length n, which gives N = 2 n possible inputs. The properties of problems we consider in this work and their objective functions are as follow: • We will only study the task of maximizing f (x). Since the minimization can be transformed into maximization by multiplying every objective value by −1. • We assume that all N inputs are feasible. Since the constraint function g(x) can be absorbed into where C is a positive constant chosen such that any feasible input has higher objective value than all unfeasible inputs (f (x) > f (y) for all x, y when g(x) = 1 and g(y) = 0). • We assume the range of f (x) are non-negative real values. Since we are considering finite search space, the objective values are bounded and minimum value exists.

III. SUBDIVIDED PHASE ORACLE
Here in place of the two-class recognition oracle O G which is used in amplitude amplification for canonical Grover, we will define O φ which changes the phase of each input in proportion to its objective value given by the objective function f (x). Even though the word "oracle" normally refers to a black-box Boolean function which recognizes solutions to decision or search problems, we will keep using it to refer to this phase shift operator to reflect the modification from amplitude amplification technique, and to be consistent with [20]- [22]. We will refer to the oracle used in this paper as the subdivided phase oracle or just the phase oracle, as taken from [20], since we subdivided the objective value into the phase of each state with equal proportion. In this section we will explain how one could construct a subdivided phase oracle from a given objective function.
The action of the phase oracle O φ is to embed the objective value into the phase of each computational basis state, where φ(x) is a function that maps the set of possible inputs {0, 1, 2, · · · , N − 1} to the phase angle with which each state is shifted in the range [θ, 2π + θ]. The mapping from objective where k is a parameter we can adjust and θ is the minimum objective value. The action of θ is to introduce a global phase which can be ignored, meaning the action of the phase oracle can be written as follows, We make a distinction between f (x) and φ(x) so that f (x) is always referred to the objective function while φ(x) is the mapped function that will be used to perform the phase shift to the computational bases. When the context is clear which k is used or when it is not important to denote the selected k value, we will omit the k and write the oracle as just O φ . From Eq. (4), it can be seen that the action of the oracle is invariant under shift of the objective values. Furthermore, the nature of the oracle does not depend nor exploit the structure of the problem, which leads to the following two simplifications: can be taken to be a non-decreasing function, f (x) ≤ f (y) for all x < y, since the action of oracle is invariant under input permutation. • We take the minimum objective value to be 0, f (x min ) = 0. As we have shown in Eq. (3) and Eq. (4) that the subdivided phase oracle is invariant under shift to the objective values.

IV. AMPLITUDE AMPLIFICATION USING PHASE ORACLE
Having established the problem definition, we now present our algorithm as below. Apply O φ (k) (Subdivided Phase Oracle) 5: Apply D (Diffusion Operator) 6: It has the same structure as Grover's algorithm but uses the subdivided phase oracle O φ (k) in place of the binary oracle operator. The diffusion operator is given by where I is the identity operator.
A. Visualizing the effect of subdivided phase oracle Before applying our algorithm to various optimization problems, it is in order to gain visual intuition about how it leads to amplitude amplification. Unlike the Grover's algorithm, where all vector amplitudes are real numbers, the subdivided phase oracle produces complex amplitudes. Therefore each amplitude corresponding to a computational basis state |x can be represented as a vector in the complex plane. Figure 1(a) represents the amplitudes of the equal superposition state |Ψ 0 , where all vectors are real and have the same length given by 1/ √ N , where N = 2 n . The action of the phase oracle O φ (k) is to rotate each vector corresponding to all the computational basis states around the origin as seen in Figure 1(b). This introduction of relative phases does not change the probabilities of measuring the individual computational basis state. The amplitude amplification occurs after the application of the diffusion operator D, which can be visualized in the following three steps: 1) Find the mean amplitude vector and scale its length by two, represented by a yellow arrow in Figure 1(c). 2) Compute the vector differences between the amplitude vectors and the scaled mean amplitude vector, depicted by purple arrows in Figures 1(d) and (e). 3) Finally, replace the old amplitude vectors with the newly computed vectors from Step 2. This effectively shifts these updated amplitude vectors to the origin, as shown in Figure 1(f). The vector corresponding to the probability amplitude of the optimal computational basis state representing the solution to the optimization problem is given by the most saturated red color in Fig. 1. We can see that application of the subdivided phase oracle followed by the diffusion operator increases its magnitude, amplifying the probability of measuring the optimal solution.

B. Amplifiability via subdivided phase oracle
Grover oracle, when viewed in the same picture can be tracked easily, as the π phase flip would make the group of marked state's amplitude vectors point to the reverse direction. When the diffusion operator is applied, the mean would lie on the real number line and after the difference vectors are created, all amplitude vectors would still lie on the real number line. It is noteworthy that when everything stays on the real number line, the amplitude of the marked states only gets larger when the mean vector is pointing in the opposite direction. The mean vector would get smaller and smaller after each iteration until it stops pointing in the opposite direction and starts to point in the same direction as vectors of marked states.
Sequences of applying O φ (k) followed by the diffusion operator lead to a more complex scenario than the Grover oracle, as both the directions and norm of the amplitude vectors and the mean vectors largely depend on the choice of O φ (k). The visualization that we discussed in the previous subsection however offers an intuitive picture when and how the subdivided phase oracles leads to amplitude amplification. Since we are interested in using this algorithm to solve optimization problems, we are only interested in the optimal state vector. The optimal state solution will grow larger up to a certain point and then shrink back and repeat this process in cycles. The condition when the vector gets larger after diffusion is only when ||2 where − → α denotes the amplitude vector. This means the angle between the mean vector and the amplitude vector of the optimal state θ( − → α best , − → α mean ) should be between π/2 and 3π/2 in order to observe significant amplitude amplification after every iteration of our algorithm. For angles outside this range, amplitude amplification occurs but its magnitude is only marginal. As we will see later, the choice of k affects this angle and hence the level of amplitude amplification.

V. RESULTS
Before we dive into the simulations and their results, we first need to discuss the connections between the objective function f (x) and the distribution of its mapped objective values. Recall that we want to have the angle between the best state amplitude vector and the mean vector be as close as opposite direction as possible, this success condition of the amplification relies heavily on the distribution of the objective values. Hence we will categorize the sets of simulations based on the density of the objective function.
In addition to the normal, skew normal, and exponential distributions, we consider the case where the objective function follows basic classes of polynomials and exponential. Although these objective functions are artificially generated and probably will not come up from natural real-world problems, we include them in our study to look at a more general behavior of the algorithm.

A. Normally distributed objective values
We begin by considering the objective values to be normally distributed according to the probability density function, where µ denotes the mean and σ the standard deviation of the distribution, as shown in Fig. 2(a). This is a natural assumption as many instances of NP-hard optimization problems can be transformed into the Ising problem [30]. Such problems can be described by using local interactions between a small number of bits of the entire bit string. This property in turn allows for the distribution of objective values of a randomly generated large instance of such NP-hard problems with local cost term sampled from a single distribution to be well approximated by a normal distribution. This approximation becomes more accurate with increasing size of the problem due to the central limit theorem and the law of large numbers.
Normal distribution is symmetric, which makes it possible for the phase difference between the mean amplitude vector and the optimal solution to be π, leading to excellent amplification results. In fact, using the visualization technique discussed in subsection IV-A, it is straightforward to observe that the mean amplitude vector and the vector corresponding to the optimal solution |N − 1 can be made to always lie on the real axis for appropriately chosen k.
Choosing k such that the mean amplitude vector and the solution vector to lie on the real axis is actually the optimal k value. Since the normal distribution is symmetric, it is the case that applying the algorithm does not only amplify the best state |N − 1 but it also amplify the worst state |0 to the same probability.

B. Skew normal distribution
We have demonstrated that optimization problems with normally distributed objective values can be tackled easily by our algorithm. The key property that made this possible was the symmetry of the normal distribution. Real-world optimization problems on the other hand are unlikely to possess such nice property. In this and the following subsections, we evaluate the performance of our algorithm for optimization problems with objective value distributions which are not symmetric. We show that not only do they not pose an obstacle to our algorithm, asymmetry in the distributions works in our favor.
To capture the deviation from normally distributed objective values we now turn to skew normal distributions, where φ norm (x) is the probability density function of the normal distribution given in Eq. (6), and The parameter α controls the shape of the distribution, as seen in Fig. 2(b)-(c). For α = 0, φ skew (x; α) becomes the normal distribution of Eq. (6), while for |α| > 0, the distribution in Eq. (7) becomes asymmetric. In the limit of α → ∞, φ skew (x; α) becomes a half-normal distribution. We observe that as for the case of a normal distribution, there is always an optimal parameter k which amplifies the amplitude of the optimal solution to a significant degree as depicted in Fig 2(e).

C. Exponential distribution
We have seen that our algorithm can be applied to optimization problems involving normal and skew normal distributions of objective values. In order to assess the applicability of the algorithm further, we expand our discussion to exponential distribution, where λ is the rate parameter, as depicted in Fig. 2(d).
Despite their rather different definitions, all three objective distribution functions described above share a number of important features when it comes to amplitude amplification, which we discuss in the following subsection.

D. Common observations
We have seen how Algorithm 1 can be applied to optimization problems with normal, skew normal, and exponential distribution of objective values. In Fig. 2(e), we display how the probability of measuring the optimal value, P solution , gets amplified with every successive application of the subdivided phase oracle O φ (k) and the diffusion operator D. Initially, P solution = 1/N ≈ 0, as the input size of the optimization problem is set to N = 2 25 . We observe that the type of distribution affects both how many iterations are needed to achieve the maximum amplitude amplification, given by P * solution , as well as the value of P * solution . Interestingly, when the objective values are distributed according to the normal distribution, Algorithm 1 yields P * solution ≈ 0.5. For the case of the two asymmetric objective value distributions on the other hand, the probability of measuring the optimal value gets amplified more strongly. P * solution = 0.94 for the skew normal distribution, and P * solution = 0.97 for the exponential distribution. Figure 2(e) displays the number of iterations needed to achieve maximum amplification of the optimal solution. Extending the number of iterations further reveals that P solution begins to oscillate as depicted in Fig. 3. This is reminiscent of the behavior observed in amplitude amplification using the conventional binary oracle. Unlike the binary oracle, the oscillations of P solution display amplitude modulation not observed previously, leading to a number of oscillations of P solution before the value of P * solution is reached again. This beating behavior is a consequence of the fact that the subdivided phase oracle introduces complex probability amplitudes into the state vector |Ψ j . This basic pattern is observed for all distributions and tested objective functions. Increasing k in the subdivided phase oracle O φ (k) leads to decreased frequency of these oscillations.
As discussed in Section IV, in order to obtain the maximum probability P * solution , we have to choose the optimal k in the subdivided phase oracle O φ (k). Algorithm 1 is sensitive to this choice and deviating from the optimal k rapidly decreases the achievable amplitude amplification. In Fig. 4, we show how the maximum amplified probability P * solution depends on the deviation of k from its optimal value for normal and skew normal distributions.
In Grover's algorithm or the usual amplitude amplification algorithm, if the diffusion operator is kept the same, being the inversion around the mean operator 2 |Ψ 0 Ψ 0 | − I, it is known that in large search space, only π phase shift on the solution states are acceptable. And if one were to change the oracle to use phase angle other than π, one would need to change the diffusion operator to match the angle of the oracle under some constraints as explored in [10], [12]. The same effect can also be seen in our results in Fig 3. The sensitivity to the angle or the k value is not just for the efficiency in terms of number of iterations but is crucial in order to make algorithm 1 work.
An important observation to make at this point is that the effectiveness of our algorithm is not affected by the size of the input space. This is a desired feature and signifies the applicability of our subdivided phase oracle approach to optimization even for extremely large data sets. Below, we discuss optimization problems which do not share this property and the degree of amplitude amplification diminishes with increasing input space. This will lead to a modified algorithm which overcomes this issue.

E. Overcoming size-induced limits to amplifiability
We have seen above that Algorithm 1 possesses a number of desired properties when applied to common distributions of objective values. One such property is that its effectivity does not diminish with increasing size of the input space of the optimization problem. In this subsection we demonstrate that this is not always the case and introduce a modification to Algorithm 1 that counteracts this issue.
We focus on the following group of objective functions, They map each state to a unique object value, therefore we refer to them as injective objective functions. These injective objective functions are not likely to appear naturally in optimization problems in this simple form. However, they serve the purpose of demonstrating an important property of Algorithm 1 that may be encountered when using objective functions not explicitly covered in our manuscript. The primary quantity that we are concerned about in this subsection is the relative amplified probability of the optimal solution, P solution , with its initial value, P initial = 1/N , before the first iteration of Algorithm 1. In particular, we are interested in determining how the ratio P solution /P initial behaves for increasing input size space N . Figure 5 displays the scaling of this probability ratio as a function of the number of iterations of Algorithm 1 for increasing N . Interestingly, we can observe that the ratio P solution /P initial converges to a single curve for increasing N . This is in stark contrast to behavior observed for normal, skew normal, and exponential distributions where the optimal solution amplitude could be amplified easily independently of the input space size N .
One might suspect that we would see a different picture if we fix the objective value ranges but increase the possible data points in between while keeping the density fixed (sampling more from the same density). Instead of mapping |x to x, we map |x to x/m where m is some constant and our total states then become mN instead of N . We would still see the same picture, the highest achievable probability relative to initial 1/N approach a certain limit. Figure 5 shows the behavior of f q (x) = x 2 but identical observations were made for the remaining injective objective functions as well. In the remainder of this section, we present a modification that overcomes this size-induced limit to amplification.
So far we have only considered applying O φ (k) with constant k for all iterations of the algorithm. There are no real restrictions that this must be the case, and in fact, we can show that by changing the parameter k to the operator O φ (k) for each iteration, we can achieve better results and overcome the above mentioned issue. Recall that in Fig 3, the general trend is that the smaller the phase angle we choose, the shorter the period of the oscillating probability. This gives us some intuition that a smaller angle can be used in early iterations to quickly boost the amplitude of the solution state then later use a larger angle which allows it to get past the amplification limit of the smaller angle. An interesting effect of varying the k at each iteration is that it allows the injective objective functions case to get past the negligible boosted probability limit by a single fixed k. We show in Fig 6 for injective families the probability of measuring the solution state at each iteration when using a greedy approach. O φ at each iteration is chosen such that the P (|N − 1 ) is highest by scanning the range of k.  26 . O φ (k) for each line is constructed using the same k which is optimal or very close to optimal. Different colors denote different size in state space. The highest probability relative to initial probability approach a certain value when N grows. Noted that the optimal state |N − 1 for each line is different but the objective functions are the same.
Although the greedy approach to dynamically choose k is not the optimal strategy, it serves as one point in favor of evidence that in general even the injective families of objective functions which cannot be boosted via a fixed choice of O φ (k) are boostable with this amplitude amplification technique.
One way to avoid dealing with the complexity when using the subdivided phase oracle O φ is to alternate the direction of the rotation for each iteration which can be done by applying O φ (k) on the even iteration and O φ (−k) on the odd iteration. This way the mean vector can be tracked easily as it always stays on the same axis throughout all iterations. Although this approach solves the problem of choosing the optimal k since the algorithm performs best when k approaches 0, the number of iterations requires to amplify the probability grows much quicker than the fixed k method, and the magnitude of the solution state is smaller. One interesting side effect of this approach is that the states with objective value f (x) higher than the mean get larger while those lower than the mean get smaller. This approach of alternating between O φ (k) and O φ (−k) was explored in detail in [21] with the help of an interference pattern created by adding one ancilla qubit to select the rotation direction in superposition which can be thought of as doubling the number of amplitude vectors and splitting them into two groups. Having them rotate in opposite directions forces the mean vector to always stay on the real number line.

VI. QUANTUM ADVANTAGE
We have shown that using our modified amplitude amplification algorithm, it is possible to raise the probability of measuring the solution to a significant degree depending on the distribution of objective values. Since the amplified probability does not reach unity, we still need to resort to classical random sampling after applying the amplitude amplification algorithm. If we look at Fig 2(e), the increase in probability P solution is faster for the first half than the later half of the iterations required to reach the peak which means that there is an optimal stopping point in which would give us smaller total overall complexity.
In order to properly quantify the speedups against classical scheme, we will use 'expected number of trials' to compare quantum and classical schemes. Let us define: where P C (P Q ) is the probability of measuring the solution classically (quantumly using Algorithm 1), and E C (E Q ) is the expected number of queries to get (measure) the solution classically (quantumly), t is the number of iterations or time step count for each trial in the quantum case, and P(|N − 1 ) is the probability of measuring the optimal solution state |N − 1 .
We make a reasonable simplification that we treat one iteration of the amplitude amplification algorithm as a single query. Query complexity is the most used metric to quantify the speedups of Grover's algorithm and amplitude amplification thus we will use this in the same manner.  Fig. 2. Location of the minimum of E Q indicates the optimal stopping point for the Algorithm 1 after which sampling the state vector |Ψ j leads to the best result.
As we have seen in the previous sections, the increase in probability does not scale linearly with number of iterations, hence it is not always the optimal strategy to sample after boosting the probability of the solution state to the highest but rather boost it enough such that the expected number of trials is minimum. Using the simulations and objective functions mentioned in the previous section, Fig 7 shows the expected number of queries and the trade-off between further amplifying the probability and the number of sampling. Now that we have defined the expected number of trials to get the optimal solution, we will look at another metric, the success probability of finding the solution. Since one trial has the success probability P (|N − 1 ), the success probability when we sample for p trials is Using the values obtained in Section V, we see that even a very small number of trials p yields probability of success to be very close to unity.

VII. CONCLUSION
Even though Amplitude amplification has been known since its discovery decades ago as one of the key quantum tools for achieving an advantage over classical computers, there are still many things to be learned from it. In this paper, we showed evidence that the amplitude amplification algorithm in combination with subdivided phase oracle can be used to solve optimization problems, and depending on the structure and distribution of their objective values, an advantage can be gained over classical schemes. Using a fixed angle oracle O φ (k), normal distribution of objective values can be amplified such that the probability of observing the solution state is of a significant degree. In the case of skew normal and exponential distributions, this probability becomes near unity. If changing k at each time-step is allowed, the algorithm is general enough to amplify other distributions of objective values, such as those that are uniform or the injective objective functions. The process of choosing an optimal k is non-trivial and needs further study but it has shown promising path in the next step to the more general amplitude amplification of more than two classes of inputs.
There are still so many things to be learned from using the subdivided phase oracle with the diffusion operator, such as the mathematical derivation of the process explained in section III or the methods of choosing k when little is known about the problem. The visualization technique that we developed and propose in this paper, albeit a simple approach, should be proved useful in designing a better suited diffusion operator to be used with respect to the choice of O φ (k) and the distribution. As with past research on Grover's algorithm, it is known that there is a relationship between the choice of phase shift operator and the angle of diffusion operator for the unstructured search. It will be interesting to see if the same thing can be derived when the oracle does not separate two groups of states but each of them uniquely.
From another perspective, the whole process can be viewed as a form of QAOA [27], [29], replacing the mixer layer with the diffusion operator [28]. Applying the low-depth QAOA might be a good starting point for selecting k and from there we can vary k using some heuristic algorithm.
Another approach to ease the sensitivity of the choice of k and to improve the potency of the algorithm is to transform the objective function f (x) into another function that exhibits better characteristics. So far we have only discussed the theoretical aspect of using this algorithm without mentioning the implementation. The embedding of objective function could be referred to the ideas presented in [24], [25] where techniques of directly embedding polynomials can be done without the need for ancilla qubits which might enable the idea of transforming the objective function without an excessively complex circuit design.