{'fraction_non_alphanumeric': 0.05296984628568368, 'fraction_numerical': 0.011782978951314874, 'mean_word_length': 4.0451229397460144, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 17, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 14, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
An Augmented Lagrangian Coordination-Decomposition Algorithm for Solving Distributed Non-Convex Programs

A novel augmented Lagrangian method for solving non-convex programs with nonlinear cost and constraint couplings in a distributed framework is presented. The proposed decomposition algorithm is made of two layers: The outer level is a standard multiplier method with penalty on the nonlinear equality constraints, while the inner level consists of a block-coordinate descent (BCD) scheme. Based on standard results on multiplier methods and recent results on proximal regularised BCD techniques, it is proven that the method converges to a KKT point of the non-convex nonlinear program under a semi-algebraicity assumption. Efficacy of the algorithm is demonstrated on a numerical example.


I. INTRODUCTION
When dealing with a large-scale system, such as a power grid for instance, sub-systems are coupled in a complex manner through their dynamics, implying that local control actions may have a major impact throughout the whole network. Implementing NMPC controllers for such systems may result in large-scale non-separable non-convex nonlinear programs (NLP). Solving such programs online in a centralised fashion is a challenging task from a computational point of view. Moreover, autonomy of the agents is likely to be hampered by such a centralised procedure. Therefore, distributed MPC is currently raising much interest, as an advanced control strategy for decentralising computations, thus reducing the computational burden and enabling autonomy of agents, which then only need to solve a local NLP [8]. Potential applications of distributed NMPC abound, from networked systems' control (interconnected power plants) to collaborative control (formation flying). Decomposing a convex NLP is generally performed by applying Lagrangian decomposition methods [4]. A critical requirement for this strategy is strong duality, which is guaranteed by Slater's condition in the convex case. However, such an assumption rarely holds in a non-convex setting. Nevertheless, in an augmented Lagrangian setting, it has been shown that the duality gap can be locally removed by using appropriate penalty functions [10]. Moreover, augmented Lagrangian methods turn out to be much more efficient than standard penalty approaches from a computational point of view, as the risk of ill-conditioning associated with large penalty parameters is reduced by the fast convergence of the dual iteration. However, the resulting quadratic penalty term is not separable even if the coupling constraints are. Several approaches have been explored to remedy this issue [7]. For instance, in [2,12] {jean-hubert.hours, colin.jones}@epfl.ch a local convexification procedure by means of proximal regularisation is analysed. A more recent approach to the decomposition of non-convex programs is the sequential programming scheme presented in [9], which consists in iterative convex approximations and decompositions. Our decomposition strategy is not specifically targeted at separable NLPs, as objective and constraints couplings can be addressed by the proposed BCD scheme. It also differs from [9] in that only local convexification is required to guarantee convergence. Our approach is an augmented Lagrangian method with partial constraint penalisation on the nonlinear equality constraints. It is a two-level optimisation scheme, whose outer layer is a loop on the Lagrange multipliers estimates associated with the nonlinear equality constraints. The inner level consists in approximately solving the primal problem to a critical point. This is performed via an inexact proximal regularised BCD scheme. By means of recent results on non-convex proximal regularised BCD methods [1], we prove convergence of the inner loop to a critical point of the non-convex inner problem under the assumption that the NLP is semi-algebraic. This is the main novelty of our approach, since, until recently, BCD type methods were thought to only have convergence guarantees under the very restrictive assumption that the NLP is convex. Some studies have been conducted in the non-convex case, yet the convergence results are quite weak compared to [1]. The inner loop is the key step of the algorithm, which enables distributed computations. In general, as the augmented Lagrangian is not separable, inner iterations cannot be fully parallelised, but some degree of parallelism can be achieved if the interconnection graph is sparse [4]. Finally, convergence of our two-level optimisation technique to a KKT point of the original NLP is proven under standard assumptions. Eventually, the whole algorithm consists in iteratively solving decoupled Quadratic Programs (QP) and updating multipliers estimates in a coordination step. In Section III, the general framework of our algorithm is presented. In Sections IV and V, the algorithm is described and its convergence analysed. Finally, a numerical example is presented in Section VI.

II. BACKGROUND
Definition 1 (Normal cone to a convex set): Let Ω be a convex set in R n andx ∈ Ω. The normal cone to Ω atx is the set Lemma 1 (Descent Lemma): Let f ∈ C 2 (Ω, R) (twice continuously differentiable in Ω), where Ω is a convex compact set in R n . Let x, y ∈ Ω.
Given M ∈ R n×n , M 2 denotes the induced matrix 2-norm. Definition 2 (Regular point): Definition 3 (Critical point): Let f be a proper lower semicontinuous function. A necessary condition for x * to be a minimiser of f is that where ∂f (x * ) is the sub-differential of f at x * [11]. Points satisfying (4) are called critical points.
The indicator function of a closed subset S of R n is denoted δ S and is defined as Lemma 2 (Sub-differential of indicator function [11]): Given a convex set C, for all x ∈ C, III. PROBLEM FORMULATION We consider NLPs of the following form: where z i ∈ R ni for i ∈ {1, . . . , N } and Z i are polytopes in R ni . The term Q(z 1 , . . . , z N ) is a cost coupling term and G(z 1 , . . . , z N ) ∈ R p a constraint coupling term. They model the different kinds of systems' interactions, which may appear in a distributed NMPC context. For the remainder of the paper, we also define the vector z := z 1 , . . . , z N ∈ R n , with dimension n : where r := m + p. The distributed objective is so that NLP (7) can be rewritten Assumption 1 (Smoothness and semi-algebraicity): and G are twice continuously differentiable and semi-algebraic.
Remark 1: The set of semi-algebraic functions is closed with respect to sum, product and composition. The indicator function of a semi-algebraic set is a semi-algebraic function. Semi-algebraicity is needed in order to apply the results of [1], which are valid for functions satisfying the Kurdyka-Lojasiewicz (KL) property, which is the case for all real semi-algebraic functions [5]. From a control perspective, this means that our results are valid for polynomial systems subject to polynomial constraints and objectives.
Assumption 2: The NLP (10) admits an isolated KKT point satisfies the second order optimality condition, that is where I * ⊆ {1, . . . , q} is the set of indices of active inequality constraints at z * . • For all i ∈ I * , λ * i > 0. IV. OUTER LOOP: PARTIALLY AUGMENTED LAGRANGIAN Given a penalty parameter ρ > 0, the augmented Lagrangian associated with problem (10) is defined as where µ ∈ R r . Only equality constraints H(z) = 0 are penalised. The polytopic inequalities are kept as constraints on the augmented Lagrangian.

A. Algorithm description
The outer loop is similar to a dual ascent on the multiplier estimates µ k , along with iterative updates of the penalty parameter ρ k . At each outer iteration k, the primal problem is solved to a given level of accuracy k > 0. The outer loop is presented in Algorithm 1 below. The notation d(x, S) stands for the distance function between a point x and a set S in R n and is defined by At every outer iteration, a pointz is found, which satisfies inexact optimality conditions. Then, the dual estimate µ is Algorithm 1 Method of multipliers with partial constraint penalisation Input: Objective J, equality constraint H, polytope Z, β > 1 and final tolerance on nonlinear equality constraints η > 0.
• Initial guess for optimiser z 0 ∈ Z, initial guess for multiplier estimates µ 0 ∈ R r , initial penalty parameter ρ 0 > 1, initial tolerance on optimality conditions updated, the penalty parameter increased and the tolerance on the first order optimality conditions is reduced. The main tuning variables are the initial value of the penalty parameter ρ 0 and the growth coefficient β.

B. Convergence analysis
We prove that the augmented Lagrangian iterations are locally convergent to the isolated KKT point satisfying Assumption 2. Our analysis is based on the results of [3]. We start by showing that the inner problem is well-defined, which means that for any optimality tolerance > 0, there exists a pointz, as in Algorithm 1, satisfying d(0, ∇L ρ (z, µ) + N Z (z)) ≤ , under appropriate conditions on ρ and µ.
Lemma 3 (Existence of an inner critical point): There exists δ > 0,ρ > 0 and κ > 0 such that the problem has a unique minimiser. Moreover, for all > 0 and for all µ, ρ ∈ S, there exists z ∈ R n such that Proof: The first part of the statement is a direct consequence of Proposition 2.11 in [3], as the second order optimality condition is satisfied (Assumption 2). From Definition 3, taking µ, ρ ∈ S, this implies that we can findz ∈ Z ∩ B(z * , κ) such that is the open ball of radius κ centered at z * . As ri B(z * , κ) ∩ ri Z = ∅, where ri stands for the relative interior, N Z∩B(z * ,κ) (z) = N Z (z) ∩ N B(z * ,κ) (z), so that The end of the proof follows by taking z :=z.
Lemma 4 (Inexact optimality conditions): Let z ∈ R n , f ∈ C 1 (R n , R) and Z be a convex set in R n . The following equivalence holds: Proof: This directly follows from the definition of the distance as an infimum.
is an isolated KKT point, there exists ν > 0 such that z * is the only critical point of (10) in B(z * , ν). From the statement of Algorithm 1, the penalty parameter ρ is guaranteed to increase to infinity and the optimality tolerance to converge to zero. The following convergence theorem is valid if the iterates z k stay in the ball B(z * , min {ν, κ}) for a large enough k. As noticed in [3], it is generally the case in practice, as warm-starting the inner problem (16) on the previous solution leads the iterates to stay around the same local minimum z * . Thus one can reasonably assume that z k ∈ B(z * , min {ν, κ}).
Theorem 1 (Local convergence to a KKT point): Let {z k } and {µ k } be sequences in R n such that Assume that all limit points of {z k } are regular. We then have that z k → z * andμ k → µ * , whereμ k := µ k + ρ k H(z k ), z * and µ * are defined in Assumption 2.
Proof: As Z ∩ cl B z * , min {κ, ν} is compact, by Weierstrass theorem, there exists an increasing mapping φ : N → N and a point z ∈ Z ∩ cl B z * , min {κ, ν} such that z φ(k) converges to z , which is regular by assumption. We also know by assumption that By definition of the normal cone to Z at z φ(k) , this means that As z is regular, there exists K ∈ N + such that ∇H(z φ(k) ) is full-rank for all k ≥ K. Subsequently, for all k ≥ K, As z φ(k) → z , by continuity of λ φ(k) b − Az φ(k) = 0, there exists λ ≥ 0 such that As d φ(k) → 0, by continuity of ∇H, this implies that µ φ(k) →μ , wherẽ (26) By taking limit in (23), we obtain Feasibility of z immediately follows from the convergence ofμ φ(k) and the fact that ρ φ(k) → +∞. Thus we have H(z ) = 0. This implies that z , μ , λ satisfies the KKT conditions. As z ∈ Z ∩ cl B z * , min {κ, ν} , in which z * is the unique KKT point by assumption, one can claim that z = z * . Taking the KKT conditions on z * and z , it follows that ∇H(z * ) μ * + A I * λ * I * = ∇H(z ) μ + A I λ I , (28) where I * and I are the sets of indices of active constraints at z * and z respectively. Obviously, as z * = z , I * = I . Thus As z * is regular, we finally obtain thatμ =μ * . As all limit points of {z k } converge to z * , one can conclude that z k → z * andμ k → µ * . Remark 2: The convergence result of Theorem 1 is local. Yet convergence can be globalised (meaning that the assumption according to which the iterates lie in a ball centered at z * can be removed) by applying the dual update scheme proposed in [6].

V. INNER LOOP: INEXACT PROXIMAL REGULARISED BCD
The non-convex inner problem (13) is solved in a distributed manner. An inexact BCD scheme is proposed, which is based on the abstract result of [1]. The main idea is to perform local convex approximations of each agent's cost and compute descent updates from it. Contrary to the approach of [9], where splitting is, in some sense, applied after convexification, we show that convexification can be performed after splitting, under the assumption that the KL property is satisfied.

A. Algorithm description
The proposed algorithm is a proximal regularised inexact BCD, which is based on Theorem 6.2 in [1], providing general properties ensuring global convergence to a critical point of the nonsmooth objective, assuming that it satisfied the KL property. The inner minimisation method is presented in Algorithm 2 below. The positive real numbers α l i N i=1 are the regularisation parameters. We assume that Matrices B l i are positive definite and need to be chosen and termination tolerance τ > 0.
carefully in order to guarantee convergence of Algorithm 2, as explained in paragraph V-B. Thus, Algorithm 2 consists in solving convex QPs sequentially. Under appropriate assumptions [4], computations can be partly parallelised. Moreover, computational efficacy can be improved by warm-starting.

B. Convergence analysis
When considering the minimisation of functions f : R n1 × . . . × R n N → R ∪ {+∞} satisfying the KL property [5] and having the following structure where f i are proper lower semicontinuous and P is twice continuously differentiable, two key assumptions need to be satisfied by the iterates z l i of the inexact BCD scheme in order to ensure global convergence [1]: a sufficient decrease property and the relative error condition. The sufficient decrease property asserts that ∀i ∈ {1, . . . , N } , ∀l ∈ N, The relative error condition states that ∀i ∈ 1, . . . , N , ∃b i > 0 such that ∀l ∈ N, (33) This means that at every iteration one can find a vector in the sub-differential of f i computed at the current iterate, which is bounded by the momentum of the sequence z l+1 i − z l i . It is actually a key step in order to apply the KL inequality, as done in [1]. For clarity, we state Theorem 6.2 in [1].
Theorem 2 (Convergence of inexact BCD): Let f be a proper lower semi-continuous function with structure (31) satisfying the KL property and bounded from below. Let α > 0 such that α ≤ α l i for all i ∈ {1, . . . , N } and l ∈ N. Let z l be a sequence satisfying (32) and (33). If z l is bounded, then z l converges to a critical point of f .
Our convergence analysis then mainly consists in verifying that the iterates generated by Algorithm 2 satisfy (32) and (33). The following definitions are necessary for the remainder of the proof: ∀i ∈ {1, . . . , N } , ∀l ∈ N, where µ G ∈ R r is the subvector of multipliers associated with the equality constraint G(z 1 , . . ., z N ) = 0.
Lemma 5 (Bound on hessian norm): where Proof: The proof directly follows from the definitions (34) and the fact that the functions involved in (34) are twice continuously differentiable over compact sets Z i .
Proof: By definition of z l+1 i in Algorithm 2, which by substituting z i with z l i ∈ Z i implies that (39) However, by applying the descent Lemma, one gets Combining this last inequality with inequality (39), one obtains which, by Assumption 3, implies that Lemma 7 (Relative error condition): ∀i ∈ {1, . . . , N }, ∀l ∈ N, where b i := 3C i + α + i and α + i is defined in (30). Proof: Writing the necessary optimality conditions for the i th QP of Algorithm 2, one obtains Thus there exists w l+1 The last equality yields from which it immediately follows However, as Now that the two main ingredients are proven, one can state the theorem guaranteeing convergence of the iterates of Algorithm 2 to a critical point of δ Z + L ρ (., µ). Theorem 3 (Convergence of Algorithm 2): The sequence z l generated by Algorithm 2 converges to a point z * satisfying or equivalently d(0, ∇L ρ (z * , µ) + N Z (z * )) = 0. Proof: The sequence z l is obviously bounded, as the iterates are constrained to stay in the polytope Z. The function δ Z + L ρ (., µ) is bounded from below. Moreover, by Assumption 1, δ Z + L ρ (., µ) satisfies the KL property for all µ ∈ R m and ρ > 0. The sequence z l satisfies the conditions (32) and (33). The parameter α of Theorem 2 can be taken as min α − i i ∈ {1, . . . , N } . Convergence to a critical point z * of δ Z + L ρ (., µ) then follows by applying Theorem 2. We can then guarantee that, given > 0, by taking a sufficiently large number of inner iterations, Algorithm 2 converges to a pointz ∈ Z satisfying d(0, ∇L ρ (z, µ) + N Z (z)) ≤ , as needed by Algorithm 1.

VI. NUMERICAL EXAMPLE
Our algorithm is tested on the following class of nonconvex NLPs: where N is the number of agents and d their dimension. Matrices H i and H i,i+1 are indefinite. When applying Algorithm 2 to minimise the augmented Lagrangian associated with (51), the update of agent i only depends on agents i−1 and i+1, which implies that, after appropriate re-ordering [4], every inner iteration actually consists of two sequential update steps where N/2 QPs are solved in parallel. We fix d = 3, N = 20, a = √ R and b = 0.6R with R = 2. Matrices H i and H i,i+1 are randomly generated. For a fixed number of outer iterations, we gradually increase the number of inner iterations and count the number of random problems, on which the algorithm ends up within a fixed feasibility tolerance with respect to the nonlinear equality constraints x i 2 2 = a 2 . Statistics are reported for 500 random NLPs of the type (51). The initial penalty parameter is set to ρ 0 = 0.1 and the growth coefficient β = 100. The hessian matrices of Algorithm 2 are set to 30ρI d . The algorithm is initialised on random primal and dual initial guesses, feasible with respect to the inequality constraints. Results for a feasibility tolerances of 10 −3 , 10 −4 and 10 −6 are presented in Figure 1. It clearly appears that obtaining an accurate feasibility margin requires a large number of iterations, yet a reasonable feasibility (10 −3 ) can be obtained with approximately hundred iterations.