Constrained Optimal Consensus in Dynamical Networks

In this paper, an optimal consensus problem with local inequality constraints is studied for a network of single-integrator agents. The goal is that a group of single-integrator a gents rendezvous at the optimal point of the sum of local convex objective functions. The local objective functions are only available to the corresponding agents that only need to know their relative distances from their neighbors in order to seek the final optimal point. This point is supposed to be confined by some local inequality constraints. To tackle this problem, we integrate the primal dual gradient-based optimization algorithm with a consensus protocol to drive the agents toward the agreed point that satisfies KKT conditions. The asymptotic convergence of the solution of the optimization problem is proven with the help of LaSalle's invariance principle for hybrid systems. A numerical example is presented to show the effectiveness of our protocol.


1.
Introduction. Over the last decade, cooperative control in a network of autonomous agents have been considered in scientific communities by virtue of big breakthroughs in wireless communication technology. Among these problems, consensus in dynamical networks is a central problem that has been studied from many aspects [15,2,3,16]. In particular, the problem of optimal consensus among networked agents has recently gained considerable attention. In this setup, the final consensus value is required to minimize the sum of individual uncoupled convex functions. For instance, the paper [14] resolved the optimal consensus problem over a network of single-integrator agents with time-varying objective function under the confining condition that Hessians associated with all local convex functions being The results of the current paper provides further developments compared to the existing literature in this area. i) Compared to [18,5,20], in the present approach the agents do not need to exchange the information of their dual variables and can reach optimal consensus by only knowing their relative positions with respect to their neighbors. ii) From design perspective, the penalty-based protocol studied in [18,5,8,21,20,13] only admits linear consensus paradigm. This restricts the protocol illustrated in these references from adopting nonlinear consensus strategies that can in turn deliver fast convergence outcomes, see e.g. [14]. Besides, in the case of high order dynamics, this approach does not work, see e.g. [14,19]. The algorithm introduced here does not have such limitation. iii) Even though the problem studied here is closely related to that of in [14,19,10,8], unlike the current paper, these references only addressed unconstrained optimization. iv) While the references [12,9] explored constrained optimization problems with convex set constraints, the projection operator utilized therein is difficult to handle in real-time specially when a large number of constraints are involved. Since a closed convex set can be approximated by a polyhedron set that is constituted by a set of linear equalities and inequalities, one can cast the optimization problem of [12,9] into the present formulation and adopt an easy-to-handle gradient-based primal-dual method discussed here to resolve it. v) The proposed algorithm does not achieve spiral trajectories toward the final point compared to the existing penalty-based algorithms (see Section 4). This paper is organized as follows. The problem formulation is given Section 2. Then, our proposed solution is presented in Section 3. A numerical example is presented in Section 4. Finally, the concluding remarks and suggestions for future studies are given in Section 5.
2. Problem Formulation. Consider N physical agents over a network with timeinvariant undirected graph G = (N , E, A), where N = {1, . . . , N } is the node set, E ⊆ N × N is the edge set , and A = [a ij ] is a weighted adjacency matrix. Each pair e = (i, j) ∈ E indicates link between the node i and the node j in an undirected graph. Suppose that each agent is described by the continuous-time single-integrator dynamicsẋ where x i (t) ∈ R represents the position of agent i, and u i (t) is the control input to agent i. We shall drop the argument t throughout this paper unless it is necessary. It is worthwhile noting that here we consider only one dimensional agents for the sake of simplicity in notations. However, it is straightforward to show that our algorithm can be extended to higher dimensional dynamics as each dimension is decoupled from others and can be treated independently. The agents are supposed to reach at an agreed point that shall minimize a convex optimization problem as in which f i (·) : R → R is the local cost function associated with node i in the network. Furthermore, g i (·) : R → R represents a constraint on the optimal position and is associated with node i. It is supposed that each agent knows only its associated cost function and inequality constraint function. We assumed only one inequality constraint per node for simplicity; our algorithm can solve the same problem with a desired number of inequality constraints. We consider the following assumptions in relation to the problem (2). Assumption 1.
(i) The objective functions f i (·), i ∈ N , and g i (·), i ∈ N , are convex and continuously differentiable on R.
The Assumption 1 and 2 fulfill the solution existence conditions for the optimization problem (2). Note that the constrained optimal consensus problem that was defined above is equivalent to the following distributed convex optimization problem In the problem (3), the consensus constraint, i.e. x i = x j , i, j = 1, . . . , N , is imposed to guarantee the same decision variable is achieved eventually. Here, the agents with dynamics as in (1) shall seek the optimal point x * i.e. x i = x * , i ∈ N , which minimize the collective cost function N i=1 f i (x i ) in a distributed fashion, given inequality constraints g i (x i ) ≤ 0, i ∈ N . To this end, each agent searches for the minimum of its associated cost function, f i (x i ), with regards to its local inequality constraint, g i (x i ), not knowing other local cost functions and constraint inequalities. Furthermore, all agents shall reach an agreement on their positions through only knowing the relative distances from their neighboring agents. Now, we shall design the control input u i to fulfill these requirements.
One can say that the problem (3) consists of a minimization sub-problem, with inequality constraints, and a consensus sub-problem. This splitting is the cornerstone of our approach to resolve the problem (2). The minimization sub-problem can be defined as and the consensus sub-problem is Before proceeding to solving the above mentioned sub-problems, i.e. the minimization sub-problem (4) and the consensus sub-problem (5), we present some optimality conditions for the optimization sub-problem through the following lemma. Later in Section 3, we will use these conditions to show the convergence of our algorithm.
is the optimal solution of the problem (4) if and only if there exist Lagrangian multipliers, λ * i > 0 , i = 1, . . . , N , such that the following conditions are satisfied To solve the minimization sub-problem (4), we focus on the primal-dual method that seeks the saddle point of the Lagrangian associated with convex optimization sub-problem (4). The Lagrangian is defined by is convex inx and concave in λ. We have the following properties for L(x,λ), where x * ,λ * is said to be the saddle point of L(x,λ) [1, p. 238]. The following inequalities hold for all (x,λ) ∈ dom(L) (8), one can define the Lagrangian function for node i as In the sequel, we will use L and L i to denote the aggregate Lagrangian (8) and the Lagrangian corresponding to node i, i.e. (12), respectively. Hence, the main task of this paper is to find the saddle point of (8) while consensus on the agents' states is also achieved.
3. Main Results. We propose the following dynamics to find the saddle point of (8) and satisfy the consensus constraint (5) where α > 0 and is the set of neighbors corresponding to node i. The operators ∇ xi and ∇ λi are the partial derivative with respect to x i and λ i , respectively. Note that −α∇ xi L i + h i acts as the control input for agent i, i.e.
In (14), a positive projection is used to ensure that Lagrangian multipliers remain non-negative. For scalars, [p] + q = p if p > 0 or q > 0, and [p] + q = 0 otherwise. When [p] + q = 0, the projection is said to be active. Therefore, in (14) when λ i > 0 and g i (x i ) < 0,λ i < 0 and λ i decreases until it reaches 0 where the projection becomes active and it remains 0 until the sign of g i (x i ) turns. Note that we start with λ i (0) > 0; therefore, λ i ≥ 0 for all t > 0. One can define the set of active projection by σ = {i : λ i = 0, g i (x i ) < 0}. Note that the control command (15), consists of two parts. The first part is to minimize the local cost function and the second part is associated with the consensus error. The following lemma is instrumental to some of the results presented in this paper.
Before proving that the algorithm in (13) and (14) yields the saddle point of (8), we show that the positions of agents, x i , i ∈ N , reach consensus, when taking control input as u i = −α∇ xi L i + h i , i ∈ N . This is established in the next proposition.
Proposition 1. Suppose that the graph G is connected and undirected. Then, there exists some finite t 1 such that the agents (1) under the protocol (15) Proof. The compact form of dynamics of agents (1) with (15) Let the network's consensus error be defined byē x = Πx, where Π = I N − 1 N 1 N 1 N andx denotes the all states of the whole network, that is defined byx = [x 1 . . . x N ] . Note that 1 Π = 0 and Π1 = 0. Thus, one can writeė where D = [d ik ] ∈ R N ×|E| is the incidence matrix associated with the topology G. And, its entries i.e. d ik , are obtained by assigning an arbitrary orientation for the edges in G. For instance, if one considers the k th edge i.e. e k = (i, j), then d ik = −1 if the edge e k leaves node i, d ik = 1 if it enters node i, and d ik = 0 otherwise. We choose the Lyapunov candidate function V (ē x ) = 1 2ē xēx . By taking time derivative from V (ē x ) along the trajectories ofē x , one can writė where v 2 (DD ) denotes the smallest non-zero eigenvalue of DD . In the above, the first inequality is resulted from Lemma 3.1, and the second inequality is resulted from the assumption ∇ xi L i − ∇ xj L j < ω 0 given in the statement of the proposition. From (17), one can say thaṫ where 0 < θ < 1. For ē x ≥ ω 0 α βv 2 (DD ) − θ , we obtainV (ē x ) < 0. Now, we are ready to invoke Theorem 5.1 from [7] that guarantees that by choosing β large enough, one can make the consensus error, δ 0 , as small as desired.
Remark 1. Assumption ∇ xi L i − ∇ xj L j < ω 0 in Proposition 1 seems to be unreasonable at the first glance as it assumes that the primal and dual variables x i and λ i , i = 1, . . . , N , must remain bounded. However, we will show by the following lemma that this requirement always holds. It is worthwhile mentioning that by choosing a conservative bound on ω 0 one can adjust the protocol's parameters to reach consensus with any desired accuracy.
We now assert that the trajectories generated by the dynamics (13) and (14) are globally bounded.
Lemma 3.2. Given that the graph G is connected and undirected, the solutions of (13) and (14) are globally bounded.
Proof. We study boundedness of the solutions of (13) and (14) by Lyapunov stability analysis. Let us define a quadratic Lyapunov function as In the above equation, x * ,λ * represents a saddle point equilibrium associated with L(x,λ). By taking derivative from both sides of (19) along the trajectories (13) and (14), with respect to time, we will havė where Suppose that for some index i, the projection becomes active i.e. i ∈ σ. In this case λ i = 0 and ∇ λi L i = g i (x i ) < 0. It is worthwhile noting that λ i < 0 never holds when parameters are initialized by positive values. Thus, in this case one can conclude that (λ i − λ * i )∇ λi L i ≥ 0 due to the fact that ∇ λi L i < 0 and λ * i ≥ 0. On the other hand, for the agents the projection is not active, Thus, we can assert that the following inequality holds.
Then, from (9) and (10), we havė The inequality (22) is due to (11). Furthermore, the last equality results from the fact that j∈Ni (x i − x j ) = 0 in a network with the undirected graph G. It is easy to show that − N i=1 x i j∈Ni (x i − x j ) ≤ 0 in an undirected graph. Hence, W (x,λ) ≤ 0, and, thus, the proof is concluded.
The dynamics (13) and (14) can be regarded as a hybrid system due to switching projection operator on the right side of the relation (14). Thus, before proceeding to the main result of this section, we introduce the LaSalle's invariance principle for hybrid systems through a lemma first given in [11] and later summarized in [4]. Next, in the light of the above lemma, we express the main result of this section. Theorem 3.4. Assume that f i (x i ) and g i (x i ), i ∈ N , are twice continuously differentiable on R. Given Assumption 1 and 2, the dynamics (13) and (14) will converge to x * ,λ * that is the solution to the optimization sub-problem (4).
Proof. To prove the theorem, it suffices to show that dynamics (13) and (14) will converge to a saddle point associated with the Lagrangian function (8).
To this end, we split the proof into two parts. We first illustrate that the Lyapunov function is always decreasing. Then, in the second part, we appeal to Lemma 3.3 to establish that the optimality conditions in Lemma 2.1 hold. To examine the above Lyapunov function, we only need to consider two scenarios, namely, the one in which the index set σ changes and the other one where this set is fixed. One should note that in the former case the Lyapunov function (24) might be discontinuous asλ i switches when σ changes according to (14). However, in the latter, the Lyapunov function (24) is continuous. In the following, we establish that in both cases the positive function (24) is always non-increasing. We first assume that σ is fixed. Taking derivative of V (ẋ,λ; σ) along the trajectories (13) and (14) with respect to time, we obtaiṅ The above equations can be simplified by expanding some of its terms into two cases, namely, i ∈ σ and i / ∈ σ. Note that when i ∈ σ, λ i = 0,λ i = 0. Thus, we can Then after a simple algebraic simplification, it is easy to verify thaṫ From the definition of h i , we attain the following equality: Now, with substituting (28) in (27), we obtaiṅ From Assumption 1 and that β α > 0, it is attained thaṫ V (ẋ,λ; σ) ≤ 0.
In the following, we will show that the same property holds even when the set σ changes. Consider conditions under which the index set σ varies: (1) Consider the case at given time index, say t 0 , the index set σ is enlarged. This happens when there is a larger number of constraints with g i (x i (t + 0 )) < 0 compared to those with Here t − 0 and t + 0 stand for the moment just before and after t 0 , respectively. (2) Now suppose that the index set σ shrinks. This case occurs when the set loses a constraint i at time t 0 and g i (x i (t + 0 )) becomes positive. Since g i (·) is a continuous function and x i is continuous as well, it can be said that this function has passed through zero to become positive. The latter supports that the new termλ 2 i is added to V (ẋ,λ; σ) but since g i (x i (t + 0 )) = g i (x i (t − 0 )), no discontinuity happens. Therefore, one can say V (ẋ,λ; σ) does not change in this case and, therefore, remains non-increasing according to (30). Now, we invoke Lemma 3.3 that presents LaSalles invariance principle for hybrid systems. From Lemma 3.2, we conclude that whole space R 2N represents an invariant set for the hybrid dynamics (13) and (14). On the other hand, in the first part of the proof, we showed that the Lyapunov function (24) decreases along the trajectories produced by (13) and (14). According to the statement of Lemma 3.3 there should exist maximal invariant set, say M , that satisfies conditions (a) and (b) stated in Lemma 3.3. In the sequel, we will show that (13) and (14) will stabilize at the point in which conditions (a) and (b) are met; moreover, the KKT conditions (6) and (7) are also fulfilled.
We first attend to part (a). From the equation (29), we attainẋ i = 0, i ∈ N , i.e.x ≡x * since one can derive from (13) thatx is continuous. Also, (7) is satisfied. As forλ, assume that g i (x * i ) > 0, then, λ i will grow unboundedly that it contradicts its boundedness shown earlier in Lemma 3.2. Therefore, g i (x * i ) ≤ 0, then two possible cases happen: i) λ i would decrease until it reaches at zero, producing a discontinuity once the projection becomes active. This would contradict with part (b) of Lemma 3.3. ii) λ i = 0; the corresponding projection is active for some i. Thus, g i (x * i ) ≤ 0 and λ * i = 0 always hold, and, (6) is met. In the above, we showed that the equilibrium point of the dynamics (13) and (14) is a saddle point of the Lagrangian function (8), and in the light of Saddle Point Theorem [17,Theorem 4.7], it is the optimal solution to (4).
One should note that through Proposition 1, we showed consensus on states, i.e. x i = x j , i, j = 1, . . . , N . Furthermore, by Theorem 3.4, we proved that the control inputs (15) drive the agents towards the saddle point of the Lagrangian associated with (4). Hence, the optimal consensus problem (3) associated with the network of single-integrator agents (1) is resolved.

Remark 2.
There exists a trade-off between size of the control command and permitted consensus error when selecting parameters α and β. As β increases, according to Proposition 1, the consensus error becomes smaller while the control input size attains a larger value. On the other hand, with small α, the consensus error decreases; however, this decelerates the optimization process. 4. Simulation Example. As mentioned earlier, the results of current paper also hold when agents are modelled by several integrators i.e. x i ∈ R m . We exploit this fact and consider the following scenario that clearly exhibits the results of this paper through a numerical simulation. Consider four agents that move in a 2-D space and are connected under a ring topology. Assume that each agent is modelled by one single-integrator dynamics per coordinate. Their local objective functions are as f 1 (x 11 , has a local constraint as g 1 (x 11 , x 12 ) = −x 11 − x 12 + 1 ≤ 0. Agent 2 suffers the constraint g 2 (x 21 , x 22 ) = x 2 21 + x 2 22 − 2 ≤ 0. Agent 3 has the local constraint g 3 (x 31 , x 32 ) = x 2 31 + x 2 32 − 1 ≤ 0, while agent 4 has no constraint. Let α = 0.1 and β = 10 be the control law's coefficients as in (15). Under the control law (15), the trajectories of agents' positions are shown in Fig. 1 when the initial positions of the agents 1, 2, 3, and 4 are set as x 1 = 2 3 , x 2 = 1 4 , x 3 = 3 4 , and x 4 = 5 0 , respectively. We set the initial values for the Lagrangian multipliers as zero. The optimal solution to the problem is 0.85 0.53 . Among many existing penalty-based algorithm, due to the page limitation, we only compare our result with that of the algorithm proposed by [21] on the above  Figure 2. Simulation results of [21]: States' trajectories for a ring network of single-integrator agents example (see Fig. 2). As it is observed, with the primal-dual dynamics proposed in [21], the agents spiral around the optimal point in a preplexed way to reach the optimal point. Such trajectories towards the final point will impose too much energy consumption and practically are not feasible to achieve.

5.
Conclusion. We studied constrained optimal consensus problem for an undirected network of single-integrator agents. We proposed a fusion algorithm in which: i) a primal-dual gradient method was used to satisfy KKT conditions for constrained convex optimization problems, and ii) a consensus protocol was adopted to make all agents reach the agreed optimal value. Then, through the theory of stability of perturbed systems, we showed that this algorithm delivers consensus. Moreover, we proved that the equilibrium point of the network's dynamics coincides with the optimal solution to the optimization problem, adopting the LaSalle's invariance principle for hybrid systems. Finally, we illustrated the performance of our proposed algorithm through a numerical example.