1. Introduction

Our work is focused on a distributed algorithm for data flow routing through multihop ad hoc and sensor networks, where all data has to be delivered to the destinations in time. An example of a target application is a network periodically sensing some consumption variables (like gas consumption, water consumption, etc.) in large objects. Each sensing device produces a data flow of a particular volume, which is supposed to be routed through the network. The objective is to optimize the energy consumption for the data transfer (minimal possible consumption in the whole network), while constrained by the communication capacities (maximum data volume which can be transferred by a device per time unit) and by maximum communication delay.

We assume a time division multiple access (TDMA) protocol (e.g., GTS allocation in IEEE 802.15.4 [1]) which ensures collision-free communication and causes constant communication delay. This approach is used quite often in industry and time-critical systems. Due to the TDMA mechanism assumed, the worst-case delay from the source node to the sink node is the sum of the particular delays for each of the hops, assumed to be an integer (derived from the parameters like TDMA period, worst-case execution time of the communication stack…). In a particular setting, we may assume a unit hop delay (e.g., the same TDMA period). In this paper we assume the unit hop delay (the deadlines are expressed as the number of communication hops between devices) which is very transparent for the reader.

There are many communication protocols designed to find the exact energy optimal data routing in wireless sensor networks. However, to comply with the communication capacities, they usually need a central computational point with the knowledge of the actual network structure and parameters (e.g., [2]). The existence of such a computational point decreases the robustness of the system against the network damage. Furthermore, the routing of information about the actual network structure and parameters has to be solved in the case of the centralized algorithm.

In this paper, we propose a distributed algorithm, which computes the energy optimal real-time data routing without the need of any central computational or data point. The algorithm supposes that each node knows only the cost (energy consumption per transmitted data unit) of its outgoing communication links and the data which is supposed to send and receive. The whole algorithm is mathematically derived using the convex programming theory. The main purpose of this paper is to present the principle of new distributed routing algorithms rather than to present an application-ready algorithm. We believe that the presented approach can lead to new, efficient, and highly adaptive routing algorithms for ad hoc and sensor networks. Moreover, the approach used in this work can be adapted for general distribution of convex optimization problems into the network.

2. Related Works

Traditionally, the routing problems for data networks are formulated as linear or convex multicommodity network flow routing problems, for example [2–4]. One of the advantages of this method is that several cost functions and constraints can be put together (e.g., different types of capacity and energy consumption constraints and real-time constraints). Using the same underlying model, we can easily combine the solution of different works focused on partial problems.

Several papers have been performed in the area of real-time routing in multihop wireless sensor networks. In [5], a well-known soft real-time communication protocol SPEED is presented. Several works use the relation between the message propagation speed and transmitting energy to balance the trade-off between the energy consumption and communication delay [6–8]. In [9] the authors deal with real-time communications over cluster-tree sensor networks, where they evaluate the end-to-end communication delay. However, none of these algorithms can ensure real-time and energy optimal data routing especially in high-loaded networks.

There are several works, which focus on the decomposition of network problems described by strictly convex optimization. A systematic presentation and classification of the decomposition techniques for network utility maximization (NUM) is presented in [10–12].

In [13–15], the authors use dual decomposition to decompose cross-layer optimization problems into the optimization of separated layers. The presented approaches lead to structural decomposition (e.g., to the routing layer, the capacity layer…) which is not suitable for the derivation of the in-network distributed algorithm. In [16], a general distributed algorithm for a strictly convex optimization problem with a common parameter for all nodes is presented. The decomposition of an optimal routing problem is presented, for example, in [17, 18], where the authors focus on the node-path formulation and based on dual decomposition derive an algorithm suitable for networks with a small number of communication paths.

In [19, 20], the authors use the duality theorem and node-link problem formulation to analytically derive the maximum network lifetime for linear networks. In [21], the authors use the node-link problem formulation and derive the distributed routing algorithm with an extension for the queuing delay. However, all these algorithms are limited to strictly convex cost functions and fail in the case of linear cost functions.

This paper also continues our previous work. It joins and extends the results from [22, 23]. In [22], we have derived a distributed in-network routing algorithm based on dual decomposition. The algorithm consists of one iteration loop and is even able to solve problems with a linear objective function. In [23], we have presented a centralized algorithm for energy optimal real-time routing, which is represented as a min-cost multicommodity network flow routing problem with side constraints.

2.1. Contribution and Outline

The main contributions of this paper are(1)introduction of a new mathematically derived, distributed algorithm for energy optimal real-time routing based on network replication and dual decomposition;(2)introduction of a new approach for an in-network distribution of real-time routing problems with constant communication delay and linear objective function;(3)application of a distribution procedure for routing problems presented in [22] on a real-time routing problem [23].

The paper is organized as follows. Section 3 briefly describes the multicommodity network flow model. In Section 4 we define the real-time routing model. In Section 5, the distributed algorithm and its derivation are presented. An example and computational complexity experiments are given in Section 6. Section 7 concludes the paper and mentions the future work. The proof of the algorithm convergence is presented in the Appendix.

3. Multicommodity Network Flow Model

In this section, we briefly summarize the basic terminology and specify the multicommodity network flow model. For more details see, for example, [3, 4].

The network is represented by an oriented graph, where for each device able to send or receive data, a node of the graph exists. The nodes are labeled as . Directed communication links are represented as ordered pairs of distinct nodes. The links are labeled as . We define the set of the links leaving the node as and the set of the links incoming to node as . The network structure is described with two incidence matrices in node-link form. The matrix of the incoming links is denoted as , and the matrix of the outgoing links is denoted as :

By we denote an index of the communication demand and by we denote a set of all communication demands. The communication demands can be seen as the flow of various commodities incoming/leaving the network in some nodes. The flow of each communication demand has to satisfy the flow conservation law at each node (for a given commodity, the sum of the flow incoming to the node is equal to the sum of the flow leaving the node) where the column vector denotes the flow coming into the network, the denotes the flow leaving the network, and the denotes the flow routed through the network for demand . Notice that a multisource multisink problem can be described in this way (e.g., the data gathering problem).

In this work, the node capacity constraints are used. They are the most common capacity constraint used for wireless networks. The total volume of the flow leaving one node has to satisfy the capacity constraint , where is the column vector of the node capacities and matrix describes the constraint structure (i.e., ).

In summary, the network flow model imposes the following constraints on the network flow variables :

The task of the total energy minimization is to minimize the cost function by setting the flow vector for all subject to the system of inequalities (3). Vector is a column vector of the energy consumption per data unit transmitted. The components of vector correspond with the energy consumption for the individual links in the network. Their values are usually determined directly from the transmission energy level in the nodes, which is needed for a reliable connection. This is done by the MAC layer of the communication protocol.

Another target application of this approach is a network, where the energy supplies of the nodes can be recharged with different maintenance costs (e.g., depending on the node reachability). In such an application, the recharge price is included in vector , and we minimize the network long-time maintenance cost.

4. Real-Time Multicommodity Network Flow Model

In this section, we extend the multicommodity flow model by real-time constraints, which guarantee a controllable routing delay through the network. Each communication demand has its own deadline, and the communication delay of each demand has to be shorter than its deadline. We model the hop delay as an integer value associated with each communication link. For transparent model derivation, we assume the same communication delay over the entire network (i.e., each communication hop causes a delay equal to one). However, the model can be extended to a more general form, where the communication delays are integers.

4.1. Mathematical Model of Real-Time Routing

Let vector denote the flow of the communication demand with an integer communication delay in the network. Then the flow vector , independent of the flow delay of demand , is equal to the sum of the flow vectors over all acceptable delays: , where denotes the deadline of the communication demand . Using this equation, we can rewrite the capacity constraint from (3) into a new form:

Vector stands for the flow of the demand leaving the network with the communication delay , and vector denotes the flow of demand coming into the network with the initial delay . The flow of each demand may come into the network and leave it in more nodes. If all flow fragments of one demand coming into the network in different nodes have the same initial delay, then , and for . The flow of demand leaving the network prior to the deadline is

Through (4) and (5), we have converted the real-time constraint (i.e., the delay has to be shorter than the deadline) to the structural constraint. Only the flow, whose delay is shorter than the deadline, is represented. The flow, which does not meet the deadline, violates the flow conservation law, and then the network flow constraints are not satisfied; that is, this solution is not feasible.

If the flow is sent through the network, the flow delay is increased by each communication hop. The flow of demand coming into node with communication delay has to either leave the network in node with the same delay or reach the neighbor node with delay . The flow conservation law from (2) can be rewritten in the delay awareness form as

In summary, the real-time energy optimal multicommodity flow routing problem can be written as

Problem (7) describes the problem in a form, which can be solved in centralized way by any linear programming solver.

4.2. Intuitive Presentation of Extended Graph

In this section, we illustrate, in an intuitive way, the graph transformation, which has been discussed in Section 4.1 by mathematical equations. New variables have appeared for each communication link as well as new constraining equations for each node. These variables and constraints can be seen as virtual layers of the network where each layer represents a different communication delay . The number of the network layers is equal to the integer deadline of the demand plus one (the number of allowed communication hops plus a zero layer). As consistent with the structure of (6), all communication links are redirected to the nodes in the higher layer, which means that the flow is routed not only in node space but also in delay space. Because the number of layers is limited by the deadline and the flow can leave the network only in virtual nodes of the sink nodes, all possible routings through this transformed network hold the deadlines. An easy example of the graph replication is shown in Figure 1.

Figure 1: Example for the intuitive presentation.

5. Decomposition of the Routing Problem

To decompose the routing algorithm, we use the gradient optimization method to solve its dual problem. However, the linearity of the cost function of the problem (7) would cause oscillations in the gradient algorithm and prevents us finding the optimal solution. Therefore, we use the proximal-point method (for details see [4]) to modify the problem into a strictly convex form, which allows the usage of the gradient method. Moreover, we rewrite the problem into an equality form for a more transparent presentation: where the objective function is

We have added slack variables into problem (8) to convert the problem into the equality form. The variables , , and in (9) are the proximal-point variables, which have been added to convert the objective function from a linear form to a strictly convex one. Please notice that the set of optimal solutions for problem (8) is the same as for the problem (7). For the optimal solution of problem (8) it holds , , . In this way the routing problem has been separated into two nested subproblems. The internal subproblem is the minimization over the variables , , , and it is strictly convex. The outer subproblem minimizes the internal one by the proximal-point variables , , .

5.1. Dual Problem

To solve the internal subproblem of (8) (minimization over the variables , , ) we present its dual problem, which allows us to derive the distributable gradient algorithm. According to Slater’s conditions (see e.g., [24]) the optimal solutions of the dual and primal problems have the same optimal values in this case.

The Lagrangian function of the problem (8) is where , , and are primal variables and , and are dual variables. The dual function is The minimizers of the dual function (11) are where symbol denotes a positive or zero value in each component of vector and .

The dual problem of the internal subproblem of (8) is The gradients of the dual function (11) are The gradients of the dual problem (13) are

5.2. Dual Gradient Algorithm

Using the dual problem (13) and the dual function (11) we rewrite the routing problem (8) into the form: A gradient algorithm created from problem (16) consists of 2 nested loops, and it is described as (17). The internal loop solves the dual problem (13) using the gradients (14). The outer loop minimizes problem (16) over the proximal-point variables using the gradients (15): The is a constant step size of the algorithm.

The important property of algorithm (17) is that it can be distributed as an in-network algorithm. Each node is responsible for the computation of the flow volume of the links leaving the node and for computation of the other corresponding variables. All the components of the variable vectors are a function of the node local variables and of the variables of the neighboring nodes that are within one hop communication distance. Only peer-to-peer communication during the algorithm is needed.

However, the distributed version of such an algorithm would have problems with the termination of the internal loop and with the synchronization of the loops between the nodes. Moreover the nested loops would cause an increase in the iterations.

To overcome these problems, we join both loops of the gradient algorithm into one loop, where we update all the variables simultaneously. Using (12), (14), and (15) we derive the iterative algorithm which is presented in Table 1, where denotes the iteration number.

Table 1: Distributed routing algorithm.

The correctness of such an algorithm is not seen directly from its derivation and has to be proven. The proof of the algorithm convergence is not a trivial problem, and its simplified version is presented in the Appendix. A necessary condition for the algorithm convergence assumed in the proof is . Moreover, we have performed several simulation experiments to test the algorithm convergence in Section 6.

The initial variables , , , , , are set to arbitrary initial values. The closer the values are to the final solution, the faster the algorithm converges. This property can be used in the case of minor changes of the network structure during its operation or in case of a precomputed routing, for example, based on Dijkstra’s algorithm.

5.3. Variables Initialization

As we have mentioned in Section 5.2, the number of iterations depends on the initial variables , , , , , . The closer the initial variables are to the final solution the less number of iterations are needed. In this work, we focus only on the dual variables. The proximal-point variables , , are set to zero. The variables , cannot be estimated without the knowledge of the final routing. We set them to zero.

Let us suppose some initial variable which satisfies a condition: According to (12), such initial variables do not cause any changes of the algorithm variables on its own (i.e., no variable changes if ).

According to the complementary slackness condition (for details see, e.g., [4]) for the optimal solution holds for two neighboring nodes: if then , if then , and the link is a part of the shortest path for the demand . (We use the index to denote the vector component corresponding to the link , to denote end node of the link, and to denote the start node of the link.) This fact leads us directly to Dijkstra’s algorithm which can be used in a distributed way. If for all sink nodes we set for all and for the other nodes we set , we get an initial setting, which satisfies condition (18). Moreover, this initial setting is much closer to the optimal solution than the setting for all . In Section 6.3, we present several experiments to evaluate the heuristic behavior in comparison with the zero initial setting. Please notice that due to the capacity constraints, the heuristic solution does not need to be equal to the final optimal solution. According to our experiments it rapidly decreases the number of iterations.

6. Experiments

To demonstrate the behavior and the correctness of the distributed routing algorithm, we have performed several experiments in Matlab. We have focused on a basic problem, where for each communication demand one node sends data flow to one sink node (i.e., a multicommodity mono-source, mono-sink problem). Therefore, for the source node and for the sink node of the communication demand .

The random networks for the experiments have been constructed as follows. We consider a square field of size [size × size], where the size is changing during the experiments. The field is divided into subsquares of size . One node is randomly placed into each subsquare, and the communication distance is set to 1.7 (i.e., node can communicate with node , if and only if their Euclidean distance is less than 1.7). Please notice that our network is close to the “unit-disk network" [2]. The communication costs per transmitted data flow unit have been set as the square power of the distance between the nodes. The link capacities have been set to two , and the maximum number of communication hops to for all . The constants of the algorithm have been set as and . The initial values , , , , have been set to 0 and for all experiments except for those in Sections 6.2 and 6.4. Only feasible problems are used.

During the experiments we evaluate as a number of iterations needed to achieve less than 1% deviation of the cost function from the optimal value, less than 1% capacity violation, and less than 1% flow conservation law violation during the last 500 iterations. The 500 iterations are included in the statistics. (the optimal value was computed separately by a centralized algorithm for evaluation purposes only).

Unfortunately, as the presented approach of dual decomposition of the routing problems in node-link form is a new technique, there are only few works in this area at this moment. The most similar work can be found in [21]. However, it is focused on slightly different problems. It uses different energy consumption and communication delay models, and it is limited only on problems with strictly convex objective functions. Our work is focused on the problems with linear objective functions and presents different qualities than [21]. Due to the differences in the used models and algorithms derivations the experimental comparison would have a disputable contribution and would be strongly dependent on the chosen problems.

6.1. Example

To present the resulting optimal data flow routing in the network we have performed an experiment based on the network described above. The field size has been set to 9 (i.e., 81 nodes in the network), and the number of communication demands has been set to 8. The link capacities have been set to for this problem. The initial variables have been set according to Section 5.3.

The optimal data flow routing in the network is shown in Figure 2. The 3D routing model in the node-delay space is presented in Figure 3 where the vertical axis represents the communication delay. No communication demand routing has more than 6 communication hops.

Figure 2: Optimal data flow routing (multicommodity, mono-source, mono-sink problem).

Figure 3: Optimal data flow routing in node-delay space.

The algorithm behavior can be seen in [25] on video. The flow routing in the network and in the node-delay space and the progress of the flow conservation law are presented.

6.2. Algorithm Convergence

We have performed a set of experiments with random initial variables on random networks, and we have evaluated the algorithm convergence. The field size has been set to 7. There are 8 communication demands in the network. The initial values have been set randomly from intervals: for proximal-point variables, and for dual variables. The intervals have been chosen as the double value of typical optimal values.

The algorithm has been run 2000 times on random networks. The results are presented in Figure 4. There, the number of iterations is placed on the horizontal axis, and the number of experiments, which has been finished before the number of iterations, is on the vertical axis.

Figure 4: Algorithm convergence with random initial variables.

This experiment provides an important practical verification of the theoretical proof of the algorithm convergence. 95% of the experiments have been finished before 14300 iterations.

6.3. Number of Iterations

To demonstrate the statistical behavior of the algorithm, we have gradually increased the number of nodes from 16 to 100 (i.e., the field size from 4 to 10) for 8 communication demands. The computation has been repeated, on random networks 1000 times for each number of nodes.

The results have been evaluated as a maximum, average, and minimum number of iterations needed to achieve the optimal value, and it is presented in Figure 5. There, the experiment progress is presented for the basic algorithm without the initial heuristic in black (the initial variables have been set to zero) and for the algorithm with the initial heuristic according to Section 5.3 in gray.

Figure 5: The number of iterations in relation to the number of nodes. In black, for zero initial variables. In gray, for initialization according to Section 5.3.

Similar experiments have been performed for the iterations dependence on the number of communication demands. We have gradually increased the number of communication demands from 1 to 16 in networks with 49 nodes. The computation has been repeated, on random networks 1000 times for each number of communication demands. The results are presented in Figure 6 for both algorithms with and without the initial heuristic.

Figure 6: The number of iterations in relation to the number of communication demands. In black, for zero initial variables. In gray, for initialization according to Section 5.3.

In Table 2, we present the percentages of infeasible problems for each number of communication demands. The values say how much percent of the generated problems have not been feasible due to the capacity and real-time constraints. From Table 2, it follows that the network has been close to saturation and that the constraints affect the final routings.

Table 2: Percentage of infeasible problems in the experiment for the number of communication demands.

The important outcome of these experiments is the observation that the number of iterations is approximately linear. It follows that the algorithm is easily applied to networks with many nodes. We can see a significant improvement of the initial heuristic in the Figures.

6.4. Network Change

The advantage of the one-loop algorithm presented in this work is that it can automatically adjust the data routing in case of network changes. In order to evaluate the algorithm behavior in this case, we have simulated a dying node in the network as follow.

(1) We generated a random network with communication demands and found the optimal solution. (2) In the original problem from step 1, we removed one node and measured the number of iterations to find the new optimal solution. (3) We repeated step 2 for 15 different nodes with the largest data flow.

The simulation has been performed for 942 random original networks with field size set to 7 and with 8 communication demands (i.e., ).

The results are presented in Figure 7. There, the number of iterations is placed on the horizontal axis, and the number of experiments, which has been finished before the number of iterations, is on the vertical axis.

Figure 7: Algorithm convergence for network change simulation.

of the experiments have been finished before 4466 iterations. Moreover, less than of the dual variables and of the primal variables have been changed during the experiments, which significantly decrease the amount of transmitted data.

Video presentation can be seen in [25].

7. Conclusion

We have presented a distributed algorithm for the real-time energy optimal routing in ad hoc and sensor networks for systems with linear cost functions and constant communication delays. Such networks are used very often in the industrial environment, where TDMA-derived mechanisms or mechanisms with a significant stack computation delay are used. We have described the routing problem as a multicommodity network flow optimization problem, modified it into a time aware form, and used the dual decomposition method to derive the distributed algorithm. The algorithm does not need any central computational node with knowledge about the whole network structure, and it only uses local communication between the neighboring nodes. This rapidly increases the robustness of the algorithm in the case of partial network damage. We have performed several simulations to evaluate the algorithm behavior and to test its convergence. The mathematical proof of the algorithm convergence is available in the Appendix.

The main purpose of this paper was to present the basic concept of a new in-network distributed routing algorithms for real-time data and its derivation. From the experimental section it is seen that the algorithm is not application-ready because of the high number of iterations, which leads to high number of communications. However, the main strength of the algorithm is not to find the whole optimal routing in an unknown network, but to adapt an existing routing in case of local network changes (dead/new node, loss of connection, communication costs, or capacities change, etc.) where the number of data communications could be significantly decreased. Moreover, the algorithm can easily adapt to slow network changes. The algorithm is suitable for static networks with sporadic changes or for networks with slow continuous changes, so the routing can be adjusted.

According to our preliminary experiments the number of iterations can be significantly decreased in future work. The algorithm can be extended by heuristics based on the partial knowledge about the network structure (e.g., node geographical position) and heuristics based on Newton’s method. The results of our preliminary testing indicate that Newton’s method-based heuristics can decrease the number of iterations more than 3 times.

Considering the fact that the algorithm is based on Linear programming formulation, we believe that the principle of the algorithm and the approach used to its derivation can be used to solve many different problems in the sensor networks area, like resource sharing, network localization, object tracking, and so forth.

Appendix A

Proof of the Algorithm Convergence

We prove the convergence of the algorithm presented in Table 1 of the paper for and for as follows. First, we rewrite the equations of the algorithm into a more transparent form, which significantly simplifies the proof presentation. Then we define a merit function such that and for the optimal solution, and we show that is nonincreasing during the algorithm computation. Next, we assume the merit function to be nondecreasing for all and show, for feasible problems, that for all we have the optimal solution. We use the marking , , to denote the th component of the vectors. We simplify the notation of , , , and so forth, instead of , , and so forth, for a more compact and transparent description. Let us remind the reader that , , and that according to Slater’s conditions [24] the optimal solutions of the dual and primal problems have the same optimal values in this case.

We rewrite the problem (8) into a more transparent form: An example of the matrix and vectors transformation is where is an identity matrix. If the sizes of new vectors and matrices are The Lagrangian function of the problem (A.1) is and the algorithm equations from Table 1 can be expressed as It can be easily verified, that (A.5) are identical with the equations from Table 1. We define two diagonal matrices and as Now we can write . Then by differentiation of the Lagrangian function (A.4) we get We define the merit function as According to the Karush-Kuhn-Tucker conditions (see e.g., [24, 26]) if the solution in the th iteration of the algorithm is the optimal solution of the problem (A.1). The gradient of the merit function is: We can define a column change vector as . For we can express one iteration step of the merit function (A.10): From (A.12), it results that , and then the merit function is nonincreasing during the algorithm. Let us assume that there is a , such that the merit function is nondecreasing for all (i.e., for all ).

From (A.12), we can write for all : And it follows that Let us suppose for some that and for some index , it holds: Then we can write from  (A.15) Under the condition from (A.16) follows Under the condition and (A.17) we can write from (A.15), (A.7), and(A.5) It means which is in contradiction with the conditions (A.14) and (A.13). It follows for all that And according to (A.7) for all , it holds: From (A.18), (A.14), and (A.15) follows The result of (A.21) is that there are a finite number of changes of the matrices and such as and for . So for some , it holds for all . Based on this fact we can express from (A.5) using (A.13) for all: and from (A.13) and (A.22) follows for all From (A.23) we get a condition, which says that for all the change of variable is a circulation in the network: Based on (A.24) and (A.13) we write: And the result of (A.25) for all is To express we use the problem (A.1). According to the Karush-Kuhn-Tucker conditions and (A.26) we can write a dual function of this problem. Please notice that for all according to (A.26) and(A.7) Then for iteration it holds: From (A.28) it follows that for it holds: if then . According to the duality gap theorem (see, e.g., [26]) it holds: , where represents the primal function and a set of feasible solutions. It follows that if and the merit function (A.10) is nondecreasing then the original problem is not feasible. It follows for feasible problems We have presented a merit function which is nonincreasing during the algorithm for some and which is equal to zero for an optimal feasible solution. Next, we have shown for feasible problems that if the merit function is nondecreasing for all then according to (A.26) and (A.29) the merit function , and the solution is an optimal solution of the original problem (A.1).