Controller-based Energy-Aware Wireless Sensor Network Routing using Quantum Algorithms

Energy efficient routing in wireless sensor networks has attracted attention from researchers in both academia and industry, most recently motivated by the opportunity to use SDN (software defined network)-inspired approaches. These problems are NP-hard, with algorithms needing computation time which scales faster than polynomial in the problem size. Consequently, heuristic algorithms are used in practice, which are unable to guarantee optimally. In this short paper, we show proof-of-principle for the use of a quantum annealing processor instead of a classical processor, to find optimal or near-optimal solutions very quickly. Our preliminary results for small networks show that this approach using quantum computing has great promise and may open the door for other significant improvements in the efficacy of network algorithms.

spaces and can perform computations in exponentially large dimensions to solve complex problems. Quantum algorithms are known to outperform classical algorithms on many challenging problems such as integer factoring [3], search [4], Fourier transform [5] and training machine learning models [6,7]. Several quantum and quantum-classical hybrid algorithms have also been proposed to address different types of optimization problems. For instance, the HHL algorithm for least square fitting [8], quantum semidefinite programming algorithm for semidefinite programming [9], quantum approximate optimization algorithm for combinatorial optimization [10], adiabatic quantum computing for quadratic programming as well as NP-complete problems [11] and variational quantum eigensolver for nonlinear optimization [12].
In this work, we focus on the D-Wave adiabatic quantum computers, which were designed to approximately solve the quadratic unconstrained binary optimization (QUBO) problem, which is NP-hard. As such, these machines have been used to address many NP-hard problems in the literature in purely quantum as well as quantum-classical hybrid approaches. D-Wave machines have been used to address the graph partitioning [13,14] and graph coloring problems [15]. Warren addresses the traveling salesman problem using the D-Wave quantum annealing computers [16]. Date et al. propose a quantum-classical hybrid algorithm to train restricted Boltzmann machines and deep belief networks [17]. One of the biggest challenges pertaining to adiabatic quantum computing is the embedding a given QUBO problem onto the hardware of the quantum computer [18]. In this work, we first convert the energy efficient network routing problem as a QUBO problem and then leverage two D-Wave processors, both 2000Q and Advantage to solve it, we find that both can provide faster (in terms of raw processing time) solutions than state-of-the-art-classical solvers for some small network problems. While these problems are too small to directly show a quantum advantage, these are still hopeful signs in terms of the ability of these devices to find a "good solution quickly".

B. Novelty and Contribution
This work has made three breakthroughs in the field of quantum computing, being the first to: • Engage the computation power of a QPU to network design, in this instance focused on energy management. • Compare the 2000Q and Advantage System1. 1 processors in solving a network design problem.
• Apply the Domain Wall Encoding scheme [19] for QPU in a practical engineering problem. The work demonstrates the merits of applying a quantum annealing QPU in network design, but also provides a simulation platform that can be used by other researchers for future QPU-facilitated design and test problems.
We start by discussing the current state of the art and related work which we build upon. We then discuss the general formulation of the problem by conventional methods and how to translate it to the recently proposed domain-wall encoding, and propose an algorithm for mapping the problem. We then describe the details of our experimental methods and report our results. Finally we discuss the consequences of the results and the longer term outlook.

II. RELATED WORK
Energy consumption in network routing is caused by neighbourhood discovery, communication and computation. Energy efficient routing with quality of service (QoS) guarantee in different applications or diverse wireless sensor networks (WSNs) can be viewed as an interesting area for future investigation [1]. Unbalanced energy consumption among the nodes causes network partition and node failures, where transmission from some nodes to the sink becomes blocked [20]. Energy efficiency can be improved at various layers of the communication protocol stack of WSN. For hardware-related energy efficiency, topics have been focusing on lower power electronics, power-off mode and energy efficient modulation. For network-layer related energy Lee et. al. [21] proposed the routing schemes in consideration of both the link quality and the residual energy level. It discusses the mechanism for forwarding route requesting packets by calculation of the probability to forward or not by taking into account the link quality and residual energy level two metrics. The scheme is evaluated against the plain AODV (Ad hoc On-Demand Distance Vector) algorithm. In [22], Liu, Su, and Chou presented a simple and highly efficient strategy to form the energy aware path from source to sink node in wireless sensor network. When the event path and query path intersects, an anchor node is discovered or it can be found within the candidate region around the intersection point. It solves the spiral problem in rumour routing by keeping the event and query paths as straight as possible and proves to outperform rumor routing and achieves higher successful path discovery ratios and lower hop counts and saves more energy. In [23] Zhao et. al. solve the optimal routing path issue of wireless sensor network by formulating the problem using the optimal path set which consists of possible combination of nodes that falls within the optimal transmission range of the source node and have minimum cost value corresponding to the cost function with energy and loss rate as the input parameters. The author firstly derives the metrics to evaluate the cost function and then proposes a brute force scheme to scan all the possible combination of nodes to determine the optimal path set. Experiments have shown that the scheme can provide robust connectivity and prolong the lifetime of the network compared with benchmarks across different scenarios. Yao, Cao, and Vasilakos [24] solve the energy efficiency problem within wireless sensor network while not violating QoS metrics by formulating the problem in the framework of the open vehicle routing problem in operation research. As the OVR (open vehicle routing) problem proves to be NP-hard, the authors subsequently proposed two different heuristic algorithms to approximate the design outcome. The outputs generated demonstrated to outperform baseline protocols in achieving longation of the network lifetime within the expected delivery latency bounds.Due to complexity all today's network routing solutions for WSNs rely on heuristic algorithms In this work, the computation for the optimal set of paths is achieved by formulating the problem as a quadratic unconstrained binary optimization (QUBO) problem that can be easily mapped to the Ising model.QUBO has been used to solve network routing and similar optimisation problems in literature [25]- [27] and it is just one of possible formulation selected to demonstrate the power of quantum. The optimization can be solved by applying the quantum annealing technique used in physics to attain the ground state of the final state of the Hamiltonian system specified by the Ising model co-efficients. We are using the Quantum Processor Unit (QPU) by D-Wave Systems Inc in this work.

III. PROBLEM FORMULATION
For this problem, we assume the network controller monitors the data traffic at a regular basis ∆t, starting from t org = 0. Multiple source nodes exist in the network s i while only one destination d,i.e. the WSN sink node is present. For each packet stream p j with average data rate r j , there is associated a pair of destination node and source node (s i,j , d).It describes the packet stream p j routed to the sink node d from source node s i . Let the number of available routes to choose be a constant number K for all the packet streams, an approach similar to [28]. For each route route j,k , there exists a set of edges edge m that connect with each other to compose the route.route j,k describes the k th route for packet stream j.K is the maximum number of routes available for each packet stream to be routed from its source node to the destination node.
At a given moment t ∈ {n∆t, (n + 1)∆t}, let there be L packet streams. Consequently there would be a binary set x i of size K * L for x i represents the i (mod K)th optional route for path n= i k . x i ∈ {0, 1}. The edges associated with x i is edge i,j . For each edge j, the length is l j . the j th packet stream r j the average data rate for the j th packet stream (s i,j , d) the pair of the j th packet stream starting from the i th source node to the sink route j,k the k th path for the j th packet stream edgem the m th edge K the maximum number of paths available for every packet stream l j the length of the j th edge p i,j the j th path for packet stream starting from the i th source node edg i the i th edgē X i the indicator vector to describe which path is selected for packet stream starting from source node i xp i,j the j th path for the packet stream starting from source node i. 1: the path is selected;0: the path is not selected e i,j the energy consumption for packet stream starting from source node i and travel through path j e i,j,k the energy consumption at k th edge with source number i whose path index is j Cmax maximum link rate per edge

A. Example Illustration
the transmitter energy consumption of l bits of data across distance d E Rx (l, d) the receiver energy consumption of l bits of data across distance d f (d j ) the transmission energy consumption for the j th edge with distance d L the number of packet streams I i,j the indicator function to indicate whether the j th edge is selected for packet stream starting from node i or not There are five source nodes and one sink in Figure 1. Suppose only node 1 and node 3 are transmitting at the moment of t, each with an average data rate of r 1 and r 3 . For node 1 to sink 6, there are two paths p 1,1 and p 1,2 . For node 3 to sink 6, there are also two paths p 3,1 and p 3,2 . From the graph we realise that indicates whether the jth path for traffic stream starting from node i is selected or not. For each path, we assume the energy consumption per interval t is e i,j = edg k e i,j,k , where e i,j,k is the energy consumption at k th edge with source number i whose path index is j.
The conditions indicate for each traffic stream, only one path can be selected.
There exists another condition, which is the edge load should not exceed the maximum capacity. Assume for the edges within the network, the maximum capacity is uniform as C max . As a simple combinatorial problem, it is easy to deduce that there will be 4 combinations. Take one combination [X 1,1 ,X 3,1 ], for example, there are three edges in use.For edg 1 and edg 3 , the edge load is r 1 ≤ C max and r 3 ≤ C max respectively. For edg 2 , the edge load is

B. Energy Model
In this work, we use the energy model from [29].
The second item in the formula is the transmission energy for l bits and the first item is the device holding energy for l bits.d is the distance of the edge connecting the transmitter and the receiver.

C. Objective Function
The objective function is the summation of energy consumption per path per edge per computation interval subject to bandwidth capacity and the encoding format according to domain wall encoding [19].
Suppose f (d j ) is the transmission power consumption on edge j with length d j .SupposeX = x i |i ≤ K indicates whether a path has been selected or not. K is the total number of available paths. We go over all the edges and within each edge j, we go through all the paths.And within each valid path i that is indicated by x i , we add the edge power consumption which includes both the transmitter and receiver power consumption.
We apply slack variable technique to mitigate the inequality and equality constraint. As x i is either 0 or 1, it can be equivalently transferred to x 2 i such that the objective function becomes a quadratic function with a constant term, which will be omitted in computation.
As the second constraint sum xi∈n x i = 1∀n of the optimisation problem falls into the one hot encoding format. In order to save the computational resources (number of logical/physical qubits) so as to ease the actual computational process, we converted the encoding of the whole problem into domain wall encoding, the explanation of which is as follows.

IV. TRANSLATION TO DOMAIN-WALL ENCODING
Motivated by the enhanced performance seen in [30,31] (although not yet shown on a real engineering problem), we employ the domain-wall encoding scheme first proposed in [19], translate to this encoding we first replace one-hot constraints with: Where Z i = −2x i +1 are Ising variables used to construct that encoding (note that thex is used to distinguish these variables from the original variables x), since we wish to work in a QUBO formulation , we substitute to obtain, which can directly be used to replace variables under a one hot constraint such that since this translation takes linear terms to linear terms, a quadratic formula in x will also be be quadratic inx Algorithm 1: Mapping Algorithm Sub procedure getEdgeM in line 2 of algorithm 2 is to form a virtual three dimensional matrix that pID = M i,j,k . pID is the path ID assigned consecutively when running the P athCollector algorithm. In M i,j,k ,i, j is the node 1: Call the sub procedure to collect all feasible paths -PathCollector 2: Call the sub procedure to assign paths to respective edges -getEdgeM 3: Call the sub procedure to formulate QUBO problem -makeEffArray 4: Call the sub procedure to encode the QUBO problem -makeEncoding 5: Call the QPU API Solver Algorithm 2: Hybrid Algorithm Procedure Require: adjacency matrix, destination ID, source ID array,maxflownum for all destination and source pair do 2: no node has been assigned the relay role yet while path amount is less than maxflownum do 4: while the last relay node is not the destination node do Find the next relay node 6: Tick this node as assigned end while 8: end while end for Algorithm 3: Path Collector ID and (i, j) indicates the edge that connects node i and node j. k indicates the k th path that goes through the edge (i, j). In implementation, a two dimension matrix of size N * N by P ath Amount is created instead for manipulation convenience. In this case, (i, j, k) in the virtual three dimensional matrix corresponds to the ((i − 1) * N + j) th position in row and k th position in column in the real two dimensional matrix.

VI. EXPERIMENT
The goal of the experiment is to evaluate the performance of the QPU against classical solvers (Cplex and Gurobi) in multi-objective routing problem that has been formulated into a QUBO. We expect that QPU given current hardware maturity can guarantee a solution quaility as good as the classical solvers whilst at a faster speed.

A. Configuration and Set-up
Two random number generators are used. One is to generate the probability following uniform distribution. That is p. If p > 0.5, a value is generated by the second random number generator uniformly distributed over 1 to 5 and it is assigned to be the corresponding flow rate.
Since we are interested in the ability of the annealers to provide samples very quickly, we have used 10 samples for every anneal, except where stated otherwise. For the classical solvers, we deployed timers to record the time before the solver call and after the solver call and calculate the lap between them. For the QPU, we use qpu sampling time within the 'timing' info as the processing time per run. We are using the default anneal time of 20 µs for each run.
We further use the fixed variable technique [32] to slim the effective QUBO size submitted to the QPU.
There are two target experiment types that we have done. The first is that we go over all the possible combinations of graph size up to size 12 and source number excluding the randomness of the flow rate per source and the second is that we apply Erdos-Renyí graph generation algorithm (assigning different edge existence probabilities) and generate 20 problem samples each graph size from (4 to 12) and run the statistical analysis.

B. Experimental methods and data reporting
There are several quantities which we used as axes on our plots, for the readers convenience we define them here (see section III for details of problem formulation): • Source number: the number of nodes on the network graph which act as sources within the network. • Graph size: the size of the graph used in the network problem, larger graphs will typically lead to larger QUBOs. • Edge probability: the probability of an edge appearing within an Erdös Renyi random graph used to construct the routing problems. • QUBO size: the number of binary variables used in the problem which is passed to the annealer or classical solver, this number is always quoted before minor embedding is performed to allow a fair comparison across solvers. The QUBO size we report is after variable fixing has been performed. • Performance degradation graph size: the point along axis one step ahead of the cross point where Correct Rate intersects with Incorrect Rate or Embedding Error Rate. • Correctness rate: the number of solutions from QPU that reach the minimum to the overall number of problem instances (the fraction of the samples which returned an optimal solution). • Processing time: the time taken to attain a feasible solution • Embedding error rate: the number of problem instances that fail to be embedded onto the QPU to the overall number of problem instances. • Incorrect rate: the number of solutions from the QPU where none of the returned solutions reach the minimum energy found for all solvers. Note that Correctness rate +Incorrect rate+ Embedding error rate sum to one. We often plot the correctness rate for different solvers and we use this as a key metric to compare performance between the QPUs and classical solver over different quantities such as QUBO size, source number or the graph size.
We used macOS Sierra version 10.12.6 to run the classical algorithms. The processor is 3.4GHz Intel Core i7 and the memory is 16GB 1333MHz DDR3.
For all experiments reported here we performed 10 reads on the quantum annealer, we have chosen this relatively small number of reads to assess the ability to attain a "good solution quickly" as it is likely that solving network problems like those described here will be very time-sensitive in most real applications. It is likely that some quantities, such as the probably that a valid solution is ever found for a given problem, would be substantially improved by taking more reads.

VII. RESULTS
Before getting into discussion of how the plots depend on the properties of the graph problem we are solving, it is worth briefly stepping back and seeing how they depend on lower level properties, from figure 2, we see that for a broad sampling of the data, the correctness rate correlates strongly with QUBO size with a relatively narrow window where success probability drops off from approximately 1 to approximately 0. We also observe that the the success probability is increased by taking more reads, however in this paper we are interested in being able to get a good solution very quickly, so we base the data in the remainder of the paper on 10 reads.
A. Correctness 1) Correctness Rate based on brute force graph generation algorithm: We ran over all the possible combinations of graph topology and traffic source for graph size 4 and 5. The plots in figs. 3,4 show the correctness rate (defined as the fraction of feasible solutions obtained in a given experiment) against the source number. There are three categories that the result data from QPU can fall into: 1. QPU reports embedding error and hence can not solve the problem; 2. QPU solves the problem but it is not the most optimal; 3. QPU provides the most optimal solution. We note that for all cases at graph size 5, the the problem can be embedded and solved to optimality by the QPU.
We also note from fig. 19 that if more samples were taken the range in which this transition occurs increases, but the qualitative shape remains the same. In this paper we are restricting ourselves to studying how the annealer performs when restricted to running for a very short time and therefore a small number of samples, but it is worth remarking that many of the problems tested in this paper could be solved eventually with more anneals.
From these plots, we can tell from the data that advantage sys1.1 has an absolute advantage in terms of speed as QUBO size increases to around 10. Up to around size 20, the QPUs are typically able to solve the problems within the 10 reads we take. By a one-by-one eye-check of the data, we can tell that in these experiments 2000Q also outperforms classical solvers in terms of speed. As the solution is not found as quickly as the advantage sys1.1, in this plot, the speediness correctness rate (the ratio of the fraction of cases where a solution is found faster) for 2000Q is always 0. Figures 3,4 depict results for problems generated exhaustively for graph size 5 with source numbers 1 − 4 and all possible connected graphs. The measures are explained in section VI-B, the correct and incorrect rate refer to cases where the most optimal solution was or was not found respectively. The embedding error rate refers to cases where the problem could not be embedded successfully. We can tell from these figures that the correctness rate decreases as the source number increases for both advantage system and 2000Q. This is probably because the QUBO size increases when the number of sources increases. Furthermore, advantage system keeps a more than 60% faster rate than both the classical solvers and 2000Q across all possible number of source nodes. Even further, at the highest number of source nodes, the faster rate for advantage system reaches the highest value among all. We suspect it is because given highest number of source nodes, the practical QUBOs submitted for all the four solvers become more complicated. While classical solvers can tackle smaller size of QUBO problem at a faster rate (we suspect if the QPU processing time might decrease at the same solution quality if we reduce the sample number per run to 5 or 3), they are losing ground of larger QUBO size to QPU solvers possibly because of the parallel solution searching mechanism facilitated by quantum mechanics.   From the analysis for data collected from all the problem instances for graph size 4 that are feasible for submission to the QPU, we can tell that the QPU demonstrates an absolute advantage over the classical solvers we tested in terms of solution quality across all the problem space. It doesn't show an advantage in the speed with number of reads equal to 3000. It is because the overall QUBO size is small (= 5) so it is fast enough for classical solvers to attain a correct solution and when the number of runs is to be decreased to 5, the QPU speed will be overriding those by the classical solvers without degrading the solution quality.
2) Correctness Rate based on probabilistic graph generation algorithm: Figures 5,6,7,8,10,9,11 depict relevant quantities for data collected by using the Erdos-Renyi graph generator with the edge probability set to 0.6,0.7 and 0.9 separately. We generated 20 graph samples per edge probability and do the average in the analysis. For figures 5,6,7,8,10,9, the correct and incorrect rate refer to cases where the most optimal solution was or was not found respectively. The embedding error rate refers to cases where the problem could not be embedded successfully, and the measures are explained in section VI-B. These plots differ only in the edge probability and the QPU used.
We can tell that for edge probability 0.6, advantage system's performance starts to degrade at graph size 10 while 2000Q system's performance starts to degrade at graph size 9. Embedding errors are reported occasionally for graph size equal to 12.
As we increase the edge probability to 0.7 ,intuitively we suspect the graph becomes more complex and 11 shows that as the edge probability increases, the trimmed QUBO size increases. A larger share of embedding error cases emerges.For advantage system, the performance starts to degrade at graph size 9 while for 2000Q system, it is at graph size 8. Furthermore, more we notice that for advantage system, embedding errors occur by around 20% for graph size equal to 11 and this figure increases to 80% when graph size reaches 12. For 2000Q, embedding errors appear at graph size 8 at around 20%and reaches 100% at graph size 11 and onward.
At the edge probability to 0.9, the degradation graph size is 8 for advantage system while 7 for 2000Q. Furthermore, advantage system starts to report embedding error at higher rate when graph size reaches 10 and on wards while for 2000Q, the graph size value is 8 with around 80% embedding error rate. B. Speediness 1) Speediness based on brute force graph generation algorithm: For network size of 5, we compare the processing time between the classical solvers and the QPUs, we ignore cases where embedding errors were encountered for this analysis. The percentage of problems that reports embedding error is 57%. We observe from figures 12,13,14,15 that QPUs have a constant value of processing time across almost all the valid problem instances. For 2000Q the processing time is 0.0025 seconds while for Advantage system the    processing time is nearly half that value. There exists a larger portion of valid problem instances for 2000Q that consume longer processing time that for Advantage system against both Gurobi and Cplex. We further notice that Cplex works noticeably longer than Gurobi for all valid problems.    16,17,18,19 show the average processing time (in seconds) per graph size for QPUs and classical solvers respectively. We excluded problems where embedding errors occurred on either QPU. Hence readers will notice that for figures 16,19 and 17,18, the plots of the processing time for classical solvers have similar patterns respectively. To be more specific, from previous figures, we can tell that 2000Q is running faster to encounter embedding errors as the graph size increases, compared with the Advantage system.
2) Speediness based on probabilistic graph generation algorithm: As the same with what has been observed from 12, 2000Q shows a constant processing time on average. Whilst for Advantage system, the average processing time goes up the graph size increases but stops increasing certain graph size and beyond. For Edge probability equal to 0.9, this occurs at size 9, for edge probability equal to 0.7, this happens at size 11 while for edge probability equal to 0.6, the we cannot see this effect and hypothosize this may happen beyond size 12 -our experiment graph size upper limit. From previous figures, we have come to agree that for edge probability equal to 0.9 with graph size larger than 9, only problem instances below graph size 9 can consistently be embedded well onto the QPUs and return valid solutions. The similar conjecture can be applied to cases where the edge probability is less with larger turning points. For the phenomenon that average processing time is lessened for even larger graph size, we suspect the graph size alone can not determine the complexity of the network. Other factors take effect such as the structure of the network.
The Cplex average processing time seems to correlate less with the graph size and the edge probability whilst for Gurobi, 16,19 show that the average processing time goes up as the graph size increases in general and with larger edge probability. If this trend continues, it suggests that for larger problems there may be a crossover where Cplex becomes more effective.

VIII. CONCLUSION
We have studied the performance of two D-Wave QPUs on small network routing problems with a limited number of reads from the processor. By comparing this to the well used Gurobi and Cplex classical solvers on a standard workstation, we find that both the 2000Q and Advantage yield superior performance in terms of absolute runtime. This is an encouraging proof-of-concept result for the application of similar QPUs to this kind of problem. We find that the most relevant quantity to determine QPU performance is the overall size of the QUBO which we apply it to, although we did find that both the size of the underlying network graph, and the number of sources had a significant effect as well. While most of the problems here involved small QUBOs, larger ones were occasionally also generated. Even within a few reads the QPUs were able to solve most QUBOs below size about 20 and were not able to beyond this size. While this range is still accessible by exhaustive search methods, this study still provides useful proof-of-concept, as the problems could often be solved very quickly. This work suggests a route to practical quantum advantage where problems where problems which can be solved classically still yield an advantage by being solved much more quickly.

A. Discussion
Some preliminary experiments (not presented in this paper) are conducted in ns3, the codes of which has been made public via https://github.com/cjie3331/quantumrouting. It demonstrated the feasibility of applying a QPU in the network routing design softwaredly. By experiments (presented in this paper) conducted in python, we have shown the advantage of QPU over classical solvers towards the same QUBO problem tailored to optimal route selection. It is our interest to run the performance comparison in ns3 between QPU and classical solvers in the near future and then compare the overall design to current state-of-art heuristic algorithm such as various localised learning algorithms. Due to the limited number of qubits supported by current QPU hardware, it is our plan to employ the clustering concept in sensor network to assign controller to each micro-network in a hierarchical manner. Last but not least, we intend to implement the overall hierarchical design into automated vehicular communication based on the justification that faster computation turn-around time can more seamlessly monitor/manage the communication process to its best effort.