# Solving Linear Equations with Separable Problem Data over Directed
Networks††thanks: This work was partially supported by NSF Awards 1917177 and
1947050.

Priyank Srivastava Jorge Cortés The authors are with the Department of
Mechanical and Aerospace Engineering, UC San Diego<EMAIL_ADDRESS>
###### Abstract

This paper deals with linear algebraic equations where the global coefficient
matrix and constant vector are given respectively, by the summation of the
coefficient matrices and constant vectors of the individual agents. Our
approach is based on reformulating the original problem as an unconstrained
optimization. Based on this exact reformulation, we first provide a gradient-
based, centralized algorithm which serves as a reference for the ensuing
design of distributed algorithms. We propose two sets of exponentially stable
continuous-time distributed algorithms that do not require the individual
agent matrices to be invertible, and are based on estimating non-distributed
terms in the centralized algorithm using dynamic average consensus. The first
algorithm works for time-varying weight-balanced directed networks, and the
second algorithm works for general directed networks for which the
communication graphs might not be balanced. Numerical simulations illustrate
our results.

###### Index Terms:

Linear algebraic equations; Distributed algorithms; Directed graphs.

##  I Introduction

The importance of solving linear algebraic equations is paramount. They appear
frequently in core mathematics as well as in applications, in physics and
engineering. Nonlinear systems can often be well understood by their linear
approximation. Due to the recent development of large-scale networks coupled
with parallel processing power and fast communication capabilities, there is a
growing effort aimed at developing distributed algorithms to solve systems of
linear equations. Distributed algorithms preserve the privacy of the agents,
are robust against single point of failures, and scale well with the network
size. Keeping these considerations in mind, this paper is a contribution to
the growing body of distributed algorithms to solve linear algebraic
equations.

_Literature Review:_ Justifying the ubiquity of linear equations, there is a
vast and expanding literature to solve them efficiently, cf. [1, 2, 3] and
references therein. However, most of the works consider the information
structure where each agent knows some rows of the coefficient matrix and the
constant vector. In those cases, the collective problem has a solution if and
only if the individual equations are solvable. Instead, the problem structure
considered here is different, and assumes that each agent has a full
coefficient matrix and constant vector of its own. This setting appears
frequently in distributed sensor fusion, where sensors are spatially
distributed and they seek to build a global state estimate (e.g., about the
location of a source or the position of a target) from local measurements, cf.
[4, 5]. All the works in this category rely on the communication graph being
_undirected_. The work [4] relies on the positive definiteness of the
individual matrices to compute the updates and prove stability. [5] uses
element-wise average consensus for the coefficient matrix as well as the
constant vector, which does not scale with either the problem dimension or the
network size, and is not desirable from a privacy standpoint. [6] also
exploits the positive definite property of the individual matrices and
requires the agents to know the state as well as the matrices of the
neighbors. [7] proposes a distributed algorithm without any positive
definiteness condition, but agents are allowed to converge to different
solutions. Our approach here uses dynamic average consensus [8, 9] to estimate
certain non-distributed terms in a gradient-based algorithm for the
reformulated optimization problem. We also draw inspiration from [10, 11] on
distributed optimization to extend our treatment to deal with unbalanced
networks. However, unlike the aforementioned works where the desired solution
is not an equilibrium of the dynamics, requiring a diminishing time-varying
stepsize-like parameter to ensure convergence, here we make sure that any
solution of the linear equation is indeed an equilibrium of the proposed
dynamics. This enables us to employ Lyapunov stability analysis to establish
convergence and offers a framework to study robustness against disturbances
and errors. Our work [12] requires bidirectional 2-hop communication. In
contrast, the distributed algorithms here require information exchange only
with immediate neighbors and work for arbitrary directed graphs.

00footnotetext: We employ the following notation. , $\mathbb{R}_{>0}$ and
$\mathbb{Z}$ denote the set of real numbers, positive real numbers, and
integers, resp. $|\mathcal{X}|$ denotes the cardinality of a set
$\mathcal{X}$. $\operatorname{\mathbf{1}}$, $\operatorname{\mathbf{0}}$ and
$I$ denote a vector or matrix of all ones and zeros, and an identity matrix of
appropriate dimension, resp. We let lowercase letters to denote vectors and
uppercase letters to denote matrices. $\|x\|$ and $\|A\|$ denote the 2-norm of
a vector $x$ and the induced 2-norm of a matrix $A$, resp.
$\operatorname{diag}(x)$ denotes the diagonal matrix obtained after arranging
the entries of the vector $x$ along the principal diagonal. $A_{ij}$ denotes
the $ij$th element of a matrix $A$, $A^{\top}$ its transpose, $A^{-1}$ its
inverse (if it exists) and $\operatorname{\operatorname{null}}(A)$ its null
space. $A\otimes B$ denotes the Kronecker product between two matrices $A$ and
$B$. Unless otherwise stated, $\operatorname{\mathbf{x}}\in^{mn}$ denotes the
concatenated vector obtained after stacking the vectors
$\\{x_{i}\\}_{i=1}^{n}\in^{m}$. $A\succ\operatorname{\mathbf{0}}$ and
$A\succeq\operatorname{\mathbf{0}}$ imply that a matrix $A$ is positive
definite and semidefinite, resp. For a symmetric matrix $A$,
$\lambda_{\max}(A)$ and $\lambda_{\min}(A)$ denote its maximum and minimum
eigenvalue, resp. Regardless of the multiplicity of eigenvalue 0,
$\lambda_{2}(A)$ denotes the minimum non-zero eigenvalue of a positive
semidefinite matrix $A$. For two vectors $x,y\in^{n}$, $[x;y]$ denotes the
concatenated vector containing the entries of $x$ and $y$, in that order, and
$x>y$ means that the inequality holds elementwise.

_Statement of Contributions:_ We consider linear algebraic equations where the
coefficient matrices and constant vector for the overall problem are given,
respectively, by the summation of the individual agents’ coefficient matrices
and constant vectors. Our starting point is the exact reformulation of this
problem as a constrained optimization problem. Using the observation that the
optimal value of this optimization is zero, we reformulate it as optimization
of an unconstrained function, and propose a centralized algorithm which works
for weight-balanced networks and serves as a reference for the design of
distributed algorithms. Using dynamic average consensus, we then propose a
distributed algorithm that does not require the agent matrices to be positive
definite, works for time-varying weight-balanced networks and is guaranteed to
converge to a solution of the original problem exponentially fast. Building on
the insights gained in establishing these results, we propose a distributed
algorithm that is not limited to weight-balanced networks and is also
guaranteed to converge to a solution of the linear equation exponentially
fast.

##  II Preliminaries

Here we review basic notions from graph theory [13, 10, 14] and dynamic
average consensus [8, 9].

_Graph Theory:_ Let
$\operatorname{\mathcal{G}}=(\operatorname{\mathcal{V}},\operatorname{\mathcal{E}},\operatorname{\mathsf{A}})$
be a weighted directed graph (or digraph), with $\operatorname{\mathcal{V}}$
as the set of vertices (or nodes) and
$\operatorname{\mathcal{E}}\subseteq\operatorname{\mathcal{V}}\times\operatorname{\mathcal{V}}$
as the set of edges: $(v_{i},v_{j})\in\operatorname{\mathcal{E}}$ iff there is
an edge from node $v_{i}$ to node $v_{j}$. With
$|\operatorname{\mathcal{V}}|=n$, the _adjacency matrix_
$\operatorname{\mathsf{A}}\in^{n\times n}$ of $\operatorname{\mathcal{G}}$ is
such that $\operatorname{\mathsf{A}}_{ij}>0$ if
$(v_{i},v_{j})\in\operatorname{\mathcal{E}}$ and
$\operatorname{\mathsf{A}}_{ij}=0$, otherwise. A directed path is an ordered
sequence of vertices such that any pair of consecutive vertices is an edge. A
digraph is strongly connected if there is a directed path between any two
distinct vertices. The out- and in-degree of a node are, resp., the number of
outgoing edges from and incoming edges to it. The weighted out-degree and
weighted in-degree of a node $v_{i}$ are
$d^{\text{out}}(v_{i})=\sum_{j=1}^{n}\operatorname{\mathsf{A}}_{ij}$ and
$d^{\text{in}}(v_{i})=\sum_{j=1}^{n}\operatorname{\mathsf{A}}_{ji}$, resp. The
_out-degree matrix_ $\operatorname{\mathsf{D}}^{\text{out}}\in^{n\times n}$
and _in-degree matrix_ $\operatorname{\mathsf{D}}^{\text{in}}\in^{n\times n}$
are diagonal matrices defined as
$\operatorname{\mathsf{D}}^{\text{out}}_{ii}=d^{\text{out}}(v_{i})$ and
$\operatorname{\mathsf{D}}^{\text{in}}_{ii}=d^{\text{in}}(v_{i})$, resp. A
graph is weight-balanced if
$\operatorname{\mathsf{D}}^{\text{out}}=\operatorname{\mathsf{D}}^{\text{in}}$.
The _Laplacian_ $\operatorname{\mathsf{L}}\in^{n\times n}$ is
$\operatorname{\mathsf{L}}=\operatorname{\mathsf{D}}^{\text{in}}-\operatorname{\mathsf{A}}$.
All eigenvalues of $\operatorname{\mathsf{L}}$ have nonnegative real parts,
$0$ is simple with left eigenvector $\operatorname{\mathbf{1}}$ iff
$\operatorname{\mathcal{G}}$ is strongly connected, and
$\operatorname{\mathsf{L}}\operatorname{\mathbf{1}}=\operatorname{\mathbf{0}}$
iff $\operatorname{\mathcal{G}}$ is weight-balanced iff
$\operatorname{\mathsf{L}}+\operatorname{\mathsf{L}}^{\top}$ is positive
semidefinite, cf. [13, Theorem 1.37]. If $\operatorname{\mathcal{G}}$ is
strongly connected, it follows from [14, Lemma 3] that there exists a positive
right eigenvector $\bar{v}\in^{n}$ associated to $0$.

_Dynamic Average Consensus:_ Consider a group of $n\in~{}\mathbb{Z}_{>1}$
agents communicating over a weight-balanced digraph
$\operatorname{\mathcal{G}}$ with Laplacian $\operatorname{\mathsf{L}}$. Each
agent $i\in\\{1,\dots,n\\}$ has a state $x_{i}\in$ and an input $z_{i}\in$.
The dynamic average consensus algorithm aims at making all the agents track
the average $\frac{1}{n}\sum_{i=1}^{n}z_{i}$ asymptotically. Here we present
the algorithm following [8], where it was introduced for undirected graphs.
Consider

$\displaystyle\dot{\operatorname{\mathbf{x}}}=-\operatorname{\mathsf{L}}\operatorname{\mathbf{x}}+\dot{\operatorname{\mathbf{z}}}.$

If $\sum_{i=1}^{n}x_{i}(0)=\sum_{i=1}^{n}z_{i}(0)$ and the input
$\operatorname{\mathbf{z}}$ is bounded, then
$x_{i}(t)\to\frac{1}{n}\sum_{i=1}^{n}z_{i}(t)$ as $t\rightarrow\infty$ for
$i\in\\{1,\dots,n\\}$, cf. [8].

##  III Problem Formulation

Consider a group of $n$ agents interacting over a digraph that seek to solve
in a distributed way the linear algebraic equation

$\displaystyle\underbrace{\left(\sum\limits_{i=1}^{n}A_{i}\right)}_{A}x=\underbrace{\left(\sum\limits_{i=1}^{n}b_{i}\right)}_{b},$
(1)

where $x\in^{m}$ is the unknown solution vector, and $A_{i}\in^{m\times m}$
and $b_{i}\in^{m}$ are the coefficient matrix and constant vector
corresponding to agent $i\in\\{1,\dots,n\\}$. We assume that (1) has at least
one solution. The formulation (1) includes, as a particular case, scenarios
where each agent $i$ knows only some rows of the coefficient matrix $A$ and
constant vector $b$. Our approach consists of first formulating (1) as a
system involving $n$ unknown solution vectors, one per agent, and then
reformulating it as a convex optimization problem. Based on this
reformulation, we propose two sets of (out-)distributed algorithms (where each
agent only needs information from its out-neighbors) to find the solutions of
(1). We start by endowing each agent with its own version $x_{i}\in^{m}$ of
$x$. Then (1) can be equivalently written as

$\displaystyle\sum\limits_{i=1}^{n}A_{i}x_{i}$
$\displaystyle=\sum\limits_{i=1}^{n}b_{i},$ (2a) $\displaystyle x_{i}$
$\displaystyle=x_{j}\quad\forall i,j.$ (2b)

Equation (2b) ensures that $x_{i}=x$ for all the agents. Clearly the set of
equations (2) and the original problem (1) are equivalent. Next we formulate
(2) as a convex optimization problem. Consider the quadratic function
$f:^{mn}\to$

$\displaystyle
f(\operatorname{\mathbf{x}})=\Big{(}\sum\limits_{i=1}^{n}(A_{i}x_{i}-b_{i})\Big{)}^{\top}\Big{(}\sum\limits_{i=1}^{n}(A_{i}x_{i}-b_{i})\Big{)},$

which is convex and attains its minimum over the solution set of
$\eqref{eq:xcopy}$. For convenience, we use
$\operatorname{\mathbf{L}}=\operatorname{\mathsf{L}}\otimes I$ and
$f(\operatorname{\mathbf{x}})=(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}})^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}})$,
where $\operatorname{\mathds{1}}=\operatorname{\mathbf{1}}\otimes I$,
$\operatorname{\mathbf{A}}\in^{mn\times mn}$ denotes the block-diagonal matrix
obtained after putting the matrices $\\{A_{i}\\}_{i=1}^{n}$ along the
principal diagonal, and
$\operatorname{\mathbf{b}}=[b_{1};\ldots;b_{n}]\in^{mn}$. If
$\operatorname{\mathcal{G}}$ is strongly connected, the solutions of (2) are
the same as the optimizers of

$\displaystyle\min\limits_{\operatorname{\mathbf{x}}}$ $\displaystyle
f(\operatorname{\mathbf{x}})$ (3) s.t.
$\displaystyle\operatorname{\mathbf{L}}^{\top}\operatorname{\mathbf{x}}=\operatorname{\mathbf{0}}.$

######  Remark 1.

_(Distributed algorithmic solutions to optimization problem):_ The problem (3)
can be solved over an undirected graph by reformulating it using the
techniques in [15] and employing the saddle-point dynamics, cf. [16, 17].
These dynamics involve terms of the form $\operatorname{\mathbf{L}}^{\top}$
and, to be implemented over a digraph, would need information from in- as well
as out-neighbors and hence are not suitable for our setup. Works that deal
with distributed optimization under consensus constraints over digraphs, see
e.g. [18, 10] and references therein, require the objective function to be
separable, and therefore are not applicable either here. $\bullet$

##  IV Distributed Algorithms Over Weight-Balanced Networks

We present distributed algorithms to solve (1) over weight-balanced networks.

###  IV-A Centralized Algorithm

We first introduce a centralized algorithm using the fact that the objective
function vanishes at the optimizers of (3). Let

$\displaystyle\min\limits_{\operatorname{\mathbf{x}}}\quad\frac{1}{2}\alpha\operatorname{\mathbf{x}}^{\top}(\operatorname{\mathbf{L}}+\operatorname{\mathbf{L}}^{\top})\operatorname{\mathbf{x}}+\beta
f(\operatorname{\mathbf{x}}),$ (4)

where $\alpha,\beta>0$. Clearly, (3) and (4) have the same set of solutions if
$\operatorname{\mathcal{G}}$ is strongly connected and weight-balanced. Since
problem (4) is unconstrained, one can use gradient descent to find its
optimizers. However, the gradient
$-\alpha(\operatorname{\mathbf{L}}+\operatorname{\mathbf{L}}^{\top})\operatorname{\mathbf{x}}-\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}})$
of the objective function in (4) involves terms with
$\operatorname{\mathbf{L}}^{\top}$, whose computation would require
information from in-neighbors. Instead, we consider the following gradient-
based dynamics

$\displaystyle\dot{\operatorname{\mathbf{x}}}=-\alpha\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}-\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}}).$
(5)

Whenever convenient, we refer to (5) as
$\operatorname{\psi_{\textnormal{grad}}}$. Note that the first term in the
dynamics (5) is distributed, meaning that each agent can implement it with
information from its out-neighbors. The second term, however, requires
collective information from all the agents because of the summation across the
network. Nevertheless, this algorithm serves as the basis for our distributed
algorithm design in the next section.

The next result formally characterizes the equivalence between the equilibria
of (5) and the solutions of (1).

######  Lemma IV.1.

_(Equivalence between ( 5) and (1)):_ Let $\operatorname{\mathcal{G}}$ be a
strongly connected and weight-balanced digraph. Then for all
$\alpha,\beta\in\mathbb{R}_{>0}$, $\operatorname{\mathbf{x}}^{*}$ is an
equilibrium of (5) if and only if
$\operatorname{\mathbf{x}}^{*}=\operatorname{\mathbf{1}}\otimes x^{*}$, where
$x^{*}\in^{m}$ solves (1).

###### Proof:

The implication from right to left is immediate. To prove the implication in
the other direction, let $\bar{x}\in^{m}$ be a solution of (1) and consider
$\bar{\operatorname{\mathbf{x}}}=\operatorname{\mathbf{1}}\otimes\bar{x}$.
Since $\operatorname{\mathbf{x}}^{*}$ and $\bar{\operatorname{\mathbf{x}}}$
are equilibria of (5),

$\displaystyle\alpha\operatorname{\mathbf{L}}(\operatorname{\mathbf{x}}^{*}-\bar{\operatorname{\mathbf{x}}})+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}(\operatorname{\mathbf{x}}^{*}-\bar{\operatorname{\mathbf{x}}})=\operatorname{\mathbf{0}}.$
(6)

Let
$\mathbf{Q}_{11}=\frac{1}{2}\alpha(\operatorname{\mathbf{L}}+\operatorname{\mathbf{L}}^{\top})+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}$.
Then (6) implies

$\displaystyle(\operatorname{\mathbf{x}}^{*}-\bar{\operatorname{\mathbf{x}}})^{\top}\mathbf{Q}_{11}(\operatorname{\mathbf{x}}^{*}-\bar{\operatorname{\mathbf{x}}})=0.$

Since $\operatorname{\mathcal{G}}$ is weight-balanced,
$(\operatorname{\mathbf{L}}+\operatorname{\mathbf{L}}^{\top})\succeq\operatorname{\mathbf{0}}$.
This along with the fact that
$\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}\succeq\operatorname{\mathbf{0}}$
implies
$\operatorname{\mathbf{L}}^{\top}(\operatorname{\mathbf{x}}^{*}-\bar{\operatorname{\mathbf{x}}})=\operatorname{\mathbf{0}}$
and
$\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}(\operatorname{\mathbf{x}}^{*}-\bar{\operatorname{\mathbf{x}}})=\operatorname{\mathbf{0}}$.
Therefore, $\operatorname{\mathbf{x}}^{*}=\operatorname{\mathbf{1}}\otimes
x^{*}$, for some $x^{*}\in~{}^{m}$ which satisfies $Ax^{*}=A\bar{x}=b$, as
claimed. ∎

The next result characterizes the convergence of (5).

######  Proposition IV.2.

_(Exponential stability of ( 5)):_ Let $\operatorname{\mathcal{G}}$ be a
strongly connected and weight-balanced digraph. Then for all
$\alpha,\beta\in\mathbb{R}_{>0}$, any trajectory of (5) converges
exponentially to a point of the form
$\operatorname{\mathbf{x}}^{*}=\operatorname{\mathbf{1}}\otimes x^{*}$, where
$x^{*}\in^{m}$ solves (1).

###### Proof:

Consider a vector $\mathbf{w}\in^{mn}$ in the null space of $\mathbf{Q}_{11}$.
Using the same line of arguments as in the proof of Lemma IV.1, this implies
that $\operatorname{\mathbf{L}}^{\top}\mathbf{w}=\operatorname{\mathbf{0}}$
and
$\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}\mathbf{w}=\operatorname{\mathbf{0}}$.
Therefore, along (5),

$\displaystyle\dot{\operatorname{\mathbf{x}}}^{\top}\mathbf{w}=-(\alpha\operatorname{\mathbf{x}}^{\top}\operatorname{\mathbf{L}}^{\top}+\beta(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}})^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}})\mathbf{w}=0.$

This means that the dynamics (5) are orthogonal to the null space of
$\mathbf{Q}_{11}$ and hence the component of $\operatorname{\mathbf{x}}$ in
the null space of $\mathbf{Q}_{11}$, say
$\operatorname{\mathbf{x}}_{\operatorname{\operatorname{null}}}$, remains
constant. Given the initial condition $\operatorname{\mathbf{x}}(0)$, consider
the particular equilibrium $\operatorname{\mathbf{x}}^{*}$ of (5) satisfying
$\operatorname{\mathbf{x}}^{*}_{\operatorname{\operatorname{null}}}=\operatorname{\mathbf{x}}(0)_{\operatorname{\operatorname{null}}}$.
Since different equilibria differ only in their null space component,
$\operatorname{\mathbf{x}}^{*}$ defined this way is unique. Consider the
Lyapunov function candidate $V:^{mn}\to$

$\displaystyle
V(\operatorname{\mathbf{x}})=\frac{1}{2}(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}).$

The Lie derivative of $V$ along the dynamics (5) is given by

$\displaystyle\operatorname{\mathcal{L}}_{\operatorname{\psi_{\textnormal{grad}}}}V=$
$\displaystyle-(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}(\alpha\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}}))$
$\displaystyle=$
$\displaystyle-(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}\mathbf{Q}_{11}(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})\leq-2\lambda_{2}(\mathbf{Q}_{11})V.$

The last inequality follows from applying the Courant-Fischer theorem [19,
Theorem 4.2.11] together with the fact that
$(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}\mathbf{w}=0$
as $\operatorname{\mathbf{x}}_{\operatorname{\operatorname{null}}}$ is
constant. Using the monotonicity theorem [19, Corollary 4.3.3], we further
have

$\displaystyle\operatorname{\mathcal{L}}_{\operatorname{\psi_{\textnormal{grad}}}}V\leq-2\min\left\\{\frac{1}{2}\alpha\lambda_{2}(\operatorname{\mathsf{L}}+\operatorname{\mathsf{L}}^{\top}),\beta\lambda_{2}(\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}})\right\\}V.$

Hence, the dynamics (5) is exponentially stable with a rate depending on
$\alpha,\beta,\operatorname{\mathsf{L}}$ and $\\{A_{i}\\}_{i=1}^{n}$. ∎

###  IV-B Distributed Algorithm

We present a distributed algorithm to find a solution of (1), which is based
on the centralized algorithm (5) and involves employing dynamic average
consensus (cf. Section II) to estimate the aggregate
$\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}})$.
Formally,

$\displaystyle\dot{\operatorname{\mathbf{x}}}=$
$\displaystyle-\alpha\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}-n\beta\operatorname{\mathbf{A}}^{\top}\mathbf{y},$
(7a) $\displaystyle\dot{\mathbf{y}}=$
$\displaystyle-\alpha\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}-n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\mathbf{y}-\gamma\operatorname{\mathbf{L}}\mathbf{y},$
(7b)

with design parameter $\gamma>0$. Here, each agent $i\in\\{1,\dots,n\\}$
updates $y_{i}\in^{m}$ which estimates the average mismatch
$\frac{1}{n}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}})$.
The dynamics (7) is distributed as each agent just needs to know its state and
that of its out-neighbors. Whenever convenient, we refer to it as
$\operatorname{\psi_{\textnormal{gdac}}}$. The following result characterizes
the equilibria of (7) and shows that the total deviation from the average
mismatch is conserved.

######  Lemma IV.3.

_(Equilibria of ( 7) and invariance of total deviation):_ Let
$\operatorname{\mathcal{G}}$ be a strongly connected and weight-balanced
digraph. Then, if $(\operatorname{\mathbf{x}}^{*},\operatorname{\mathbf{0}})$
is an equilibrium of (7) then
$\operatorname{\mathbf{x}}^{*}=\operatorname{\mathbf{1}}\otimes x^{*}$, where
$x^{*}\in^{m}$. Moreover, for all $\alpha,\beta,\gamma\in\mathbb{R}_{>0}$,
$\operatorname{\mathds{1}}^{\top}(\mathbf{y}-\operatorname{\mathbf{A}}\operatorname{\mathbf{x}})$
remains constant along the evolution of (7).

###### Proof:

Let $(\operatorname{\mathbf{x}}^{*},\operatorname{\mathbf{0}})$ be an
equilibrium of (7). From (7a), it follows that
$\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}=\operatorname{\mathbf{0}}$,
and hence $\operatorname{\mathbf{x}}^{*}=\operatorname{\mathbf{1}}\otimes
x^{*}$ for some $x^{*}\in^{m}$, establishing the first statement. Now,
consider the derivative
$\operatorname{\mathds{1}}^{\top}(\dot{\mathbf{y}}-\operatorname{\mathbf{A}}\dot{\operatorname{\mathbf{x}}})=-\gamma\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{L}}\mathbf{y}=\operatorname{\mathbf{0}}$.
Hence,
$\operatorname{\mathds{1}}^{\top}(\mathbf{y}-\operatorname{\mathbf{A}}\operatorname{\mathbf{x}})$
is conserved along the evolution of (7). ∎

######  Remark 2.

_(Distributed initialization of the $\operatorname{\psi_{\textnormal{gdac}}}$
algorithm):_ From Lemma IV.3, we observe that in order for a trajectory of (7)
to converge to an equilibrium of the form
$(\operatorname{\mathbf{x}}^{*},\mathbf{y}^{*})=(\operatorname{\mathbf{1}}\otimes
x^{*},\operatorname{\mathbf{0}})$, where $x^{*}\in^{m}$ solves (1), its
initial condition must satisfy
$\operatorname{\mathds{1}}^{\top}\mathbf{y}(0)=\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}(0)-\operatorname{\mathbf{b}}).$
This could be implemented in a distributed way if each agent
$i\in\\{1,\dots,n\\}$ chooses its initial states satisfying
$y_{i}(0)=A_{i}x_{i}(0)-b_{i}$. One trivial selection, for example, is
$\operatorname{\mathbf{x}}(0)=\operatorname{\mathbf{0}}$ and
$\mathbf{y}(0)=-\operatorname{\mathbf{b}}$. $\bullet$

The next result characterizes the convergence of (7).

######  Theorem IV.4.

_(Exponential stability of ( 7) over balanced networks):_ Let
$\operatorname{\mathcal{G}}$ be a strongly connected and weight-balanced
digraph and assume
$\operatorname{\operatorname{null}}(A)\subseteq\operatorname{\operatorname{null}}(A_{i})$,
for all $i\in\\{1,\dots,n\\}$. Let $\alpha,\beta\in\mathbb{R}_{>0}$ and define

$\displaystyle\bar{\gamma}=\max\left\\{\frac{2}{\lambda_{2}{(\operatorname{\mathsf{L}}+\operatorname{\mathsf{L}}^{\top})}}\lambda_{\max}\left(\frac{\mathbf{Q}_{12}^{\top}\mathbf{Q}_{12}}{\lambda_{2}(\mathbf{Q}_{11})}-n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\right),0\right\\},$

where
$\mathbf{Q}_{11}=\frac{1}{2}\alpha(\operatorname{\mathbf{L}}+\operatorname{\mathbf{L}}^{\top})+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}$
and
$\mathbf{Q}_{12}=\frac{1}{2}(n\beta\operatorname{\mathbf{A}}^{\top}+\alpha\operatorname{\mathbf{L}}^{\top}\operatorname{\mathbf{A}}^{\top}+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top})$.
Then, for all $\gamma\in(\bar{\gamma},\infty)$, any trajectory of (7) with
initial condition satisfying
$\operatorname{\mathds{1}}^{\top}\mathbf{y}(0)=\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}(0)-\operatorname{\mathbf{b}})$
converges exponentially to
$(\operatorname{\mathbf{x}}^{*},\operatorname{\mathbf{0}})$, where
$\operatorname{\mathbf{x}}^{*}=\operatorname{\mathbf{1}}\otimes x^{*}$ and
$x^{*}\in^{m}$ solves (1).

###### Proof:

Define the error variable

$\displaystyle\mathbf{e}=\mathbf{y}-\frac{1}{n}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}}),$
(8)

measuring the difference between the agents’ estimates and the actual value of
average mismatch. Note that

$\displaystyle\dot{\mathbf{e}}=-\alpha\mathbf{\Pi}\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}-n\beta\mathbf{\Pi}\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\mathbf{y}-\gamma\operatorname{\mathbf{L}}\mathbf{y},$

where
$\mathbf{\Pi}=I-\frac{1}{n}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}$.
Rewriting (7) in terms of $\operatorname{\mathbf{x}}$ and $\mathbf{e}$,

$\displaystyle\dot{\operatorname{\mathbf{x}}}$
$\displaystyle=-\alpha\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}-\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}})-n\beta\operatorname{\mathbf{A}}^{\top}\mathbf{e},$
(9a) $\displaystyle\dot{\mathbf{e}}$
$\displaystyle=-\alpha\mathbf{\Pi}\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}-\beta\mathbf{\Pi}\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}})$
(9b)
$\displaystyle\quad-n\beta\mathbf{\Pi}\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\mathbf{e}-\gamma\operatorname{\mathbf{L}}\mathbf{e}.$

From the proof of Proposition IV.2, we know that if $\mathbf{w}\in^{mn}$ is in
the null space of $\mathbf{Q}_{11}$, then
$\operatorname{\mathbf{L}}^{\top}\mathbf{w}=\operatorname{\mathbf{0}}$ and
$\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}\mathbf{w}=~{}\operatorname{\mathbf{0}}$.
Therefore, $\mathbf{w}=\operatorname{\mathbf{1}}\otimes w$, where $w\in^{m}$
belongs to $w\in\operatorname{\operatorname{null}}(A)$. By hypothesis,
$A_{i}w=\operatorname{\mathbf{0}}$ for all $i\in\\{1,\dots,n\\}$. Therefore,
from (9a), $\dot{\operatorname{\mathbf{x}}}^{\top}\mathbf{w}=0$, and the
$\operatorname{\mathbf{x}}$ component of the equilibrium
$(\operatorname{\mathbf{x}}^{*},\mathbf{y}^{*})$ of (7) satisfies
$\operatorname{\mathbf{x}}^{*}_{\operatorname{\operatorname{null}}}=\operatorname{\mathbf{x}}(0)_{\operatorname{\operatorname{null}}}$
and is unique. With the initialization of the statement, it follows from Lemma
IV.3 that
$\mathbf{y}^{*}=\operatorname{\mathbf{1}}\otimes\frac{1}{n}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}^{*}-\operatorname{\mathbf{b}})$.
Substituting this value of $\mathbf{y}^{*}$ in (7a) and following the proof of
Lemma IV.1, one can establish that the corresponding equilibrium is of the
form $(\operatorname{\mathbf{1}}\otimes x^{*},\operatorname{\mathbf{0}})$,
where $x^{*}\in^{m}$ is a solution of (1). Consider the Lyapunov function
candidate $V_{2}:^{2mn}\to$

$\displaystyle
V_{2}(\operatorname{\mathbf{x}},\mathbf{e})=\frac{1}{2}(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})+\frac{1}{2}\mathbf{e}^{\top}\mathbf{e}.$

The Lie derivative of $V_{2}$ along $\eqref{eq:transformed}$ is given by

$\displaystyle\operatorname{\mathcal{L}}_{\operatorname{\psi_{\textnormal{gdac}}}}V_{2}=$
$\displaystyle-(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}(\alpha\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}}))$
$\displaystyle-n\beta(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}\operatorname{\mathbf{A}}^{\top}\mathbf{e}-\mathbf{e}^{\top}\mathbf{\Pi}\operatorname{\mathbf{A}}(\alpha\operatorname{\mathbf{L}}\operatorname{\mathbf{x}}+n\beta\operatorname{\mathbf{A}}^{\top}\mathbf{e})$
$\displaystyle-\mathbf{e}^{\top}(\beta\mathbf{\Pi}\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}-\operatorname{\mathbf{b}})+\gamma\operatorname{\mathbf{L}}\mathbf{e})$
$\displaystyle=$
$\displaystyle-\begin{bmatrix}\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\\\
\mathbf{e}\end{bmatrix}^{\top}\begin{bmatrix}\mathbf{Q}_{11}&\mathbf{Q}_{12}\\\
\mathbf{Q}_{12}^{\top}&\mathbf{Q}_{22}\end{bmatrix}\begin{bmatrix}\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\\\
\mathbf{e}\end{bmatrix},$

where
$\mathbf{Q}_{22}=\frac{1}{2}\gamma(\operatorname{\mathbf{L}}+\operatorname{\mathbf{L}}^{\top})+n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}$
and we have used the fact that due to the mentioned initialization,
$\operatorname{\mathds{1}}^{\top}\mathbf{e}=\operatorname{\mathbf{0}}$ from
Lemma IV.3. Since
$\operatorname{\mathbf{x}}_{\operatorname{\operatorname{null}}}$ is constant,
$(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}\mathbf{w}=0$
and from the Courant-Fischer theorem [19, Theorem 4.2.11],

$\displaystyle-(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}\mathbf{Q}_{11}(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})\leq-\lambda_{2}(\mathbf{Q}_{11})(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}).$

Also, since
$\operatorname{\mathds{1}}^{\top}\mathbf{e}=\operatorname{\mathbf{0}}$ and
$\operatorname{\mathcal{G}}$ is weight-balanced, it again follows from the
Courant-Fischer theorem that

$\displaystyle-\mathbf{e}^{\top}\mathbf{Q}_{22}\mathbf{e}\leq-\frac{1}{2}\gamma\lambda_{2}(\operatorname{\mathsf{L}}+\operatorname{\mathsf{L}}^{\top})\mathbf{e}^{\top}\mathbf{e}-n\beta\mathbf{e}^{\top}\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\mathbf{e}.$

Therefore, we can upper bound the Lie derivative as

$\displaystyle\operatorname{\mathcal{L}}_{\operatorname{\psi_{\textnormal{gdac}}}}V_{2}\leq-\begin{bmatrix}\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\\\
\mathbf{e}\end{bmatrix}^{\top}\underbrace{\begin{bmatrix}\lambda_{2}(\mathbf{Q}_{11})I&\mathbf{Q}_{12}\\\
\mathbf{Q}_{12}^{\top}&\bar{\mathbf{Q}}_{22}\end{bmatrix}}_{\bar{\mathbf{Q}}}\begin{bmatrix}\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\\\
\mathbf{e}\end{bmatrix},$

where
$\bar{\mathbf{Q}}_{22}=\frac{1}{2}\gamma\lambda_{2}(\operatorname{\mathsf{L}}+\operatorname{\mathsf{L}}^{\top})I+n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}$.
Next, we examine the positive definiteness of $\bar{\mathbf{Q}}$. Using the
Schur complement [20], $\bar{\mathbf{Q}}\succ\operatorname{\mathbf{0}}$ iff

$\displaystyle\frac{1}{2}\gamma\lambda_{2}(\operatorname{\mathsf{L}}+\operatorname{\mathsf{L}}^{\top})I+n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}-\frac{1}{\lambda_{2}(\mathbf{Q}_{11})}\mathbf{Q}_{12}^{\top}\mathbf{Q}_{12}\succ\operatorname{\mathbf{0}}.$

Hence, $\bar{\mathbf{Q}}\succ\operatorname{\mathbf{0}}$ if
$\gamma>\bar{\gamma}$, and
$\operatorname{\mathcal{L}}_{\operatorname{\psi_{\textnormal{gdac}}}}V_{2}\leq-2\lambda_{\min}(\bar{\mathbf{Q}})V_{2}$.
∎

The null space condition in Theorem IV.4 makes sure that
$\operatorname{\mathbf{x}}^{*}_{\operatorname{\operatorname{null}}}$ remains
invariant along the evolution of (7) and all the agents approach the solution
of (1) closest to $\operatorname{\mathbf{x}}(0)$. This condition is
automatically satisfied if the matrix $A$ is full rank, or in other words,
equation (1) has a unique solution. We believe (and simulations also suggest)
that if this condition is not satisfied, the $\operatorname{\mathbf{x}}$
component of the dynamics still converges to a solution of (1).

######  Remark 3.

_(Lower bound on $\gamma$):_ The lower bound $\bar{\gamma}$ in Theorem IV.4 is
conservative in general. In fact, the algorithm may converge even if this
condition is not satisfied, something that we have observed in simulation.
Note also that although $\alpha$ and $\beta$ are free parameters, they should
still be carefully chosen as $\bar{\gamma}$ depends on them. $\bullet$

The result above can be extended to time-varying networks. In case
$\operatorname{\mathcal{G}}(t)$ is time-varying, the algorithm in (7) reads as

$\displaystyle\dot{\operatorname{\mathbf{x}}}$
$\displaystyle=-\alpha\operatorname{\mathbf{L}}(t)\operatorname{\mathbf{x}}-n\beta\operatorname{\mathbf{A}}^{\top}\mathbf{y},$
(10a) $\displaystyle\dot{\mathbf{y}}$
$\displaystyle=-\alpha\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}(t)\operatorname{\mathbf{x}}-n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\mathbf{y}-\gamma\operatorname{\mathbf{L}}(t)\mathbf{y}.$
(10b)

The next result formally characterizes the convergence of (10). Its proof is
similar to that of Theorem IV.4 and hence omitted.

######  Theorem IV.5.

_(Exponential stability of ( 10) over time-varying balanced networks):_ Let
$\\{\operatorname{\mathcal{G}}(t)\\}_{t=0}^{\infty}$ be a sequence of strongly
connected and weight-balanced digraphs with uniformly bounded edge weights
(i.e., there exists $a\in(0,\infty)$ such that
$\operatorname{\mathsf{A}}_{ij}(t)<a$ for all $(i,j)$ and $t\geq 0$), and
assume
$\operatorname{\operatorname{null}}(A)\subseteq\operatorname{\operatorname{null}}(A_{i})$,
for all $i\in\\{1,\dots,n\\}$. Let $\alpha,\beta\in\mathbb{R}_{>0}$ and define
$\bar{\gamma}(t)$ as

$\displaystyle\max\left\\{\frac{2}{\lambda_{2}(\operatorname{\mathsf{L}}(t)+\operatorname{\mathsf{L}}(t)^{\top})}\lambda_{\max}\left(\frac{\mathbf{Q}_{12}(t)^{\top}\mathbf{Q}_{12}(t)}{\lambda_{2}(\mathbf{Q}_{11}(t))}-n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\right),0\right\\},$

where
$\mathbf{Q}_{11}(t)=\frac{1}{2}\alpha(\operatorname{\mathbf{L}}(t)+\operatorname{\mathbf{L}}(t)^{\top})+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}$
and
$\mathbf{Q}_{12}(t)=\frac{1}{2}(n\beta\operatorname{\mathbf{A}}^{\top}+\alpha\operatorname{\mathbf{L}}(t)^{\top}\operatorname{\mathbf{A}}^{\top}+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top})$.
Then for all $\gamma\in(\hat{\gamma},\infty)$, where
$\hat{\gamma}=\sup\limits_{t\geq 0}\bar{\gamma}(t)$, any trajectory of (10)
with initial conditions
$\operatorname{\mathds{1}}^{\top}\mathbf{y}(0)=\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}(0)-\operatorname{\mathbf{b}})$
converges exponentially to
$(\operatorname{\mathbf{x}}^{*},\operatorname{\mathbf{0}})$, where
$\operatorname{\mathbf{x}}^{*}=\operatorname{\mathbf{1}}\otimes x^{*}$ and
$x^{*}\in^{m}$ solves (1).

##  V Distributed Algorithm Over Unbalanced Networks

In this section, we extend our approach to solve problem (1) over graphs that
are not necessarily balanced. In those scenarios, since
$\operatorname{\mathsf{L}}\operatorname{\mathbf{1}}\neq\operatorname{\mathbf{0}}$,
the one-to-one correspondence between the desired equilibria of (5) or (7) and
the solutions of (1) does not hold anymore. To overcome this, we propose

$\displaystyle\dot{\operatorname{\mathbf{x}}}$
$\displaystyle=-\alpha\operatorname{\mathbf{L}}\bar{\mathbf{V}}\operatorname{\mathbf{x}}-n\beta\operatorname{\mathbf{A}}^{\top}\mathbf{y},$
(11a) $\displaystyle\dot{\mathbf{y}}$
$\displaystyle=-\alpha\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}\bar{\mathbf{V}}\operatorname{\mathbf{x}}-n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\mathbf{y}-\gamma\operatorname{\mathbf{L}}\bar{\mathbf{V}}\mathbf{y},$
(11b)

where $\bar{\mathbf{V}}=\operatorname{diag}(\bar{\mathbf{v}})$,
$\bar{\mathbf{v}}=\operatorname{\mathbf{1}}\otimes\bar{v}$, and $\bar{v}$ is a
positive right eigenvector with eigenvalue $0$ of $\operatorname{\mathsf{L}}$.
Exponential stability of (11) can be established by interpreting
$\operatorname{\mathsf{L}}\cdot\operatorname{diag}(\bar{v})$ as the Laplacian
of a weight-balanced graph and then following the same steps as in the proof
of Theorem IV.4, but we omit it here for reasons of space. Although (11) is
distributed, it assumes that agents have a priori knowledge of the
corresponding entries of $\bar{v}$ which might be limiting in practice. To
deal with this limitation, we propose an algorithm that does not require such
knowledge by augmenting (11) with an additional dynamics converging to
$\bar{\mathbf{v}}$,

$\displaystyle\dot{\operatorname{\mathbf{x}}}$
$\displaystyle=-\alpha\operatorname{\mathbf{L}}\mathbf{V}\operatorname{\mathbf{x}}-n\beta\operatorname{\mathbf{A}}^{\top}\mathbf{y},$
(12a) $\displaystyle\dot{\mathbf{y}}$
$\displaystyle=-\alpha\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}\mathbf{V}\operatorname{\mathbf{x}}-n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\mathbf{y}-\gamma\operatorname{\mathbf{L}}\mathbf{V}\mathbf{y},$
(12b) $\displaystyle\dot{\mathbf{v}}$
$\displaystyle=-\operatorname{\mathbf{L}}\mathbf{v},$ (12c)

where $\mathbf{V}=\operatorname{diag}(\mathbf{v})$. Whenever convenient, we
refer to dynamics (12) as $\operatorname{\psi_{\textnormal{dist}}}$. Note
that, unlike all the dynamics discussed so far,
$\operatorname{\psi_{\textnormal{dist}}}$ is nonlinear.

######  Remark 4.

_(Distributed nature of ( 12)):_ The dynamics (12) is out-distributed, but
requires each agent $i\in\\{1,\dots,n\\}$ to have knowledge of its in-degree
because
$\operatorname{\mathsf{L}}=\operatorname{\mathsf{D}}^{\textnormal{in}}-\operatorname{\mathsf{A}}$
and the graph is not weight-balanced. If we use instead the out-Laplacian
$\operatorname{\mathsf{L}}=\operatorname{\mathsf{D}}^{\textnormal{out}}-\operatorname{\mathsf{A}}$,
then one could still define an equivalent algorithm for (11) with
$\operatorname{\mathbf{L}}\bar{\mathbf{V}}$ replaced by
$\bar{\mathbf{V}}\operatorname{\mathbf{L}}$, but (12c) would look like
$\dot{\mathbf{v}}=-\operatorname{\mathbf{L}}^{\top}\mathbf{v}$, which would
require state information from in-neighbors too. $\bullet$

The next result characterizes the convergence of (12).

######  Theorem V.1.

_(Exponential stability of ( 12) over unbalanced networks):_ Let
$\operatorname{\mathcal{G}}$ be a strongly connected digraph and assume
$\operatorname{\operatorname{null}}(A)\subseteq\operatorname{\operatorname{null}}(A_{i})$,
for all $i\in\\{1,\dots,n\\}$. Let $\alpha,\beta\in\mathbb{R}_{>0}$ and define

$\displaystyle\bar{\gamma}=\max\left\\{\frac{2}{\lambda_{2}{(\operatorname{\mathsf{L}}\bar{V}+\bar{V}\operatorname{\mathsf{L}}^{\top})}}\lambda_{\max}\left(\frac{\mathbf{Q}_{12}^{\top}\mathbf{Q}_{12}}{\lambda_{2}(\mathbf{Q}_{11})}-n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}\right),0\right\\},$

where
$\mathbf{Q}_{11}=\frac{1}{2}(\alpha\operatorname{\mathbf{L}}\bar{\mathbf{V}}+\bar{\mathbf{V}}\operatorname{\mathbf{L}}^{\top})+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}$,
$\mathbf{Q}_{12}=\frac{1}{2}(n\beta\operatorname{\mathbf{A}}^{\top}+\alpha\bar{\mathbf{V}}\operatorname{\mathbf{L}}^{\top}\operatorname{\mathbf{A}}^{\top}+\beta\operatorname{\mathbf{A}}^{\top}\operatorname{\mathds{1}}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top})$,
$\bar{v}$ is the positive eigenvector with eigenvalue $0$ of
$\operatorname{\mathsf{L}}$ satisfying
$\operatorname{\mathbf{1}}^{\top}\bar{v}=1$, and
$\bar{V}=\operatorname{diag}(\bar{v})$. Then, for all
$\gamma\in(\bar{\gamma},\infty)$, any trajectory of (12) with initial
condition satisfying
$\operatorname{\mathds{1}}^{\top}\mathbf{y}(0)=\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}(0)-\operatorname{\mathbf{b}})$
and $\mathbf{v}(0)=\frac{1}{n}\operatorname{\mathbf{1}}$, converges
exponentially to
$(\operatorname{\mathbf{x}}^{*},\operatorname{\mathbf{0}},\bar{\mathbf{v}})$,
where $\operatorname{\mathbf{x}}^{*}=\operatorname{\mathbf{1}}\otimes x^{*}$
and $x^{*}\in^{m}$ solves (1), and
$\bar{\mathbf{v}}=\operatorname{\mathbf{1}}\otimes\bar{v}$.

###### Proof:

From [11, Proposition 2.2], we have that
$\mathbf{v}(t)>\operatorname{\mathbf{0}}$ for all $t\geq 0$. Also, since
$\operatorname{\mathbf{1}}^{\top}\operatorname{\mathsf{L}}=\operatorname{\mathbf{0}}$,
$\operatorname{\mathbf{1}}^{\top}\mathbf{v}$ is conserved along the evolution
of (12c). Hence $\mathbf{v}(t)\to\bar{\mathbf{v}}$ exponentially fast with a
rate determined by the non-zero eigenvalue of $\operatorname{\mathsf{L}}$ with
the smallest real part. Let us interpret the dynamics (12a)-(12b) as the
dynamics (11) with some disturbance $\mathbf{d}(t)$ defined by

$\displaystyle\mathbf{d}=\begin{bmatrix}\mathbf{d}^{\operatorname{\mathbf{x}}}\\\
\mathbf{d}^{\mathbf{y}}\end{bmatrix}=\begin{bmatrix}-\alpha\operatorname{\mathbf{L}}(\mathbf{V}-\bar{\mathbf{V}})\operatorname{\mathbf{x}}\\\
-\alpha\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}(\mathbf{V}-\bar{\mathbf{V}})\operatorname{\mathbf{x}}-\gamma\operatorname{\mathbf{L}}(\mathbf{V}-\bar{\mathbf{V}})\mathbf{y}\end{bmatrix},$

which goes to $\operatorname{\mathbf{0}}$ as $t\to\infty$. Consider a vector
$\mathbf{w}\in\operatorname{\operatorname{null}}(\mathbf{Q}_{11})$. Then as in
the proof of Theorem IV.4, $\mathbf{w}=\operatorname{\mathbf{1}}\otimes w$,
where $w\in\operatorname{\operatorname{null}}(A)$ and by hypothesis,
$A_{i}w=\operatorname{\mathbf{0}}$ for all $i\in\\{1,\dots,n\\}$. Since
$\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{L}}=\operatorname{\mathbf{0}}$,
therefore, $\mathbf{w}^{\top}\mathbf{d}^{\operatorname{\mathbf{x}}}=0$ and we
still have $\mathbf{w}^{\top}\dot{\operatorname{\mathbf{x}}}=0$, and the
$\operatorname{\mathbf{x}}$ component of the equilibrium
$(\operatorname{\mathbf{x}}^{*},\mathbf{y}^{*},\bar{\mathbf{v}})$ of (12)
satisfies
$\operatorname{\mathbf{x}}^{*}_{\operatorname{\operatorname{null}}}=\operatorname{\mathbf{x}}(0)_{\operatorname{\operatorname{null}}}$
and is unique. With the initialization of the statement and following the same
steps as in the proof of Lemma IV.3, one can establish that
$\mathbf{y}^{*}=\operatorname{\mathbf{1}}\otimes\frac{1}{n}\operatorname{\mathds{1}}^{\top}(\operatorname{\mathbf{A}}\operatorname{\mathbf{x}}^{*}-\operatorname{\mathbf{b}})$.
Substituting this value of $\mathbf{y}^{*}$ in (12a) and following the proof
of Lemma IV.1, one can establish that the corresponding equilibrium is of the
form $(\operatorname{\mathbf{1}}\otimes
x^{*},\operatorname{\mathbf{0}},\bar{\mathbf{v}})$, where $x^{*}\in^{m}$ is a
solution of (1). Consider now the Lyapunov function candidate
$V_{3}:^{3mn}\to$

$\displaystyle
V_{3}(\operatorname{\mathbf{x}},\mathbf{e},\mathbf{v})=V_{2}(\operatorname{\mathbf{x}},\mathbf{e})+\frac{\delta}{2}(\mathbf{v}-\bar{\mathbf{v}})^{\top}\operatorname{\mathbf{P}}(\mathbf{v}-\bar{\mathbf{v}}),$

where $\delta>0$, $\operatorname{\mathbf{P}}=\bar{\mathbf{V}}^{-1}$,
$\mathbf{e}$ is defined as in (8), and $V_{2}$ is the same function as in the
proof of Theorem IV.4. The Lie derivative of $V_{3}$ along
$\eqref{eq:distributed_v}$ is given by

$\displaystyle\operatorname{\mathcal{L}}_{\operatorname{\psi_{\textnormal{dist}}}}V_{3}=$
$\displaystyle-\begin{bmatrix}\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\\\
\mathbf{e}\end{bmatrix}^{\top}\begin{bmatrix}\mathbf{Q}_{11}&\mathbf{Q}_{12}\\\
\mathbf{Q}_{12}^{\top}&\mathbf{Q}_{22}\end{bmatrix}\begin{bmatrix}\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\\\
\mathbf{e}\end{bmatrix}+(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}\mathbf{d}^{\operatorname{\mathbf{x}}}$
$\displaystyle+\mathbf{e}^{\top}\mathbf{d}^{\mathbf{e}}-\delta(\mathbf{v}-\bar{\mathbf{v}})^{\top}(\operatorname{\mathbf{L}}^{\top}\operatorname{\mathbf{P}}+\operatorname{\mathbf{P}}\operatorname{\mathbf{L}})(\mathbf{v}-\bar{\mathbf{v}}),$

where
$\mathbf{d}^{\mathbf{e}}=-\alpha\mathbf{\Pi}\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}(\mathbf{V}-\bar{\mathbf{V}})\operatorname{\mathbf{x}}-\gamma\operatorname{\mathbf{L}}(\mathbf{V}-\bar{\mathbf{V}})\mathbf{e}$,
and
$\mathbf{Q}_{22}=\frac{1}{2}\gamma(\operatorname{\mathbf{L}}\bar{\mathbf{V}}+\bar{\mathbf{V}}\operatorname{\mathbf{L}}^{\top})+n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}$.
Interestingly, $\operatorname{\mathsf{L}}\bar{V}$ can be interpreted as the
Laplacian of a weight-balanced graph and as a result,
$\operatorname{\mathbf{L}}\bar{\mathbf{V}}+\bar{\mathbf{V}}\operatorname{\mathbf{L}}^{\top}\succeq\operatorname{\mathbf{0}}$
implying that
$\operatorname{\mathbf{L}}^{\top}\operatorname{\mathbf{P}}+\operatorname{\mathbf{P}}\operatorname{\mathbf{L}}\succeq\operatorname{\mathbf{0}}$.
Once again, following Lemma IV.3, one can establish that with the
initialization of the statement,
$\operatorname{\mathds{1}}^{\top}\mathbf{e}=\operatorname{\mathbf{0}}$ and
therefore using the Courant-Fischer theorem [19, Theorem 4.2.11] together with
the fact that
$(\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*})^{\top}\mathbf{w}=0$
due to invariance of
$\operatorname{\mathbf{x}}_{\operatorname{\operatorname{null}}}$, we can upper
bound the Lie derivative as

$\displaystyle\operatorname{\mathcal{L}}_{\operatorname{\psi_{\textnormal{dist}}}}V_{3}\leq$
$\displaystyle-\begin{bmatrix}\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\\\
\mathbf{e}\end{bmatrix}^{\top}\underbrace{\begin{bmatrix}\lambda_{2}(\mathbf{Q}_{11})I&\mathbf{Q}_{12}\\\
\mathbf{Q}_{12}^{\top}&\bar{\mathbf{Q}}_{22}\end{bmatrix}}_{\bar{\mathbf{Q}}}\begin{bmatrix}\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\\\
\mathbf{e}\end{bmatrix}$
$\displaystyle+\alpha\|\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\|\|\operatorname{\mathbf{L}}\|\|\mathbf{v}-\bar{\mathbf{v}}\|(\|\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\|+\|\operatorname{\mathbf{x}}^{*}\|)$
$\displaystyle+\alpha\|\mathbf{e}\|\|\mathbf{\Pi}\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}\|\|\mathbf{v}-\bar{\mathbf{v}}\|(\|\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\|+\|\operatorname{\mathbf{x}}^{*}\|)$
$\displaystyle+\gamma\|\mathbf{e}\|\|\operatorname{\mathbf{L}}\|\|\mathbf{v}-\bar{\mathbf{v}}\|\|\mathbf{e}\|-\delta\lambda_{2}(\operatorname{\mathbf{L}}^{\top}\operatorname{\mathbf{P}}+\operatorname{\mathbf{P}}\operatorname{\mathbf{L}})\|\mathbf{v}-\bar{\mathbf{v}}\|^{2},$

where
$\bar{\mathbf{Q}}_{22}=\frac{1}{2}\gamma\lambda_{2}(\operatorname{\mathsf{L}}\bar{V}+\bar{V}\operatorname{\mathsf{L}}^{\top})I+n\beta\operatorname{\mathbf{A}}\operatorname{\mathbf{A}}^{\top}$.
Define
$\operatorname{\mathbf{z}}=[\|\operatorname{\mathbf{x}}-\operatorname{\mathbf{x}}^{*}\|;\|\mathbf{e}\|;\|\mathbf{v}-\bar{\mathbf{v}}\|]$.
If $\gamma>\bar{\gamma}$, then
$\bar{\mathbf{Q}}\succ\operatorname{\mathbf{0}}$ and from the Courant-Fischer
theorem, we have

$\displaystyle\operatorname{\mathcal{L}}_{\operatorname{\psi_{\textnormal{dist}}}}V_{3}\leq-\operatorname{\mathbf{z}}^{\top}\underbrace{\begin{bmatrix}\lambda_{\min}(\bar{\mathbf{Q}})&0&\hat{\mathbf{Q}}_{13}(\operatorname{\mathbf{z}})\\\
0&\lambda_{\min}(\bar{\mathbf{Q}})&\hat{\mathbf{Q}}_{23}(\operatorname{\mathbf{z}})\\\
\hat{\mathbf{Q}}_{13}(\operatorname{\mathbf{z}})&\hat{\mathbf{Q}}_{23}(\operatorname{\mathbf{z}})&\delta\lambda_{2}(\operatorname{\mathbf{L}}^{\top}\operatorname{\mathbf{P}}+\operatorname{\mathbf{P}}\operatorname{\mathbf{L}})\end{bmatrix}}_{\hat{\mathbf{Q}}(\operatorname{\mathbf{z}})}\operatorname{\mathbf{z}},$

where
$\hat{\mathbf{Q}}_{23}(\operatorname{\mathbf{z}})=-\frac{1}{2}\alpha\|\mathbf{\Pi}\operatorname{\mathbf{A}}\operatorname{\mathbf{L}}\|(\operatorname{\mathbf{z}}+\|\operatorname{\mathbf{x}}^{*}\|)-\frac{1}{2}\gamma\|\operatorname{\mathbf{L}}\|\operatorname{\mathbf{z}}$
and
$\hat{\mathbf{Q}}_{13}(\operatorname{\mathbf{z}})=-\frac{1}{2}\alpha\|\operatorname{\mathbf{L}}\|(\operatorname{\mathbf{z}}+\|\operatorname{\mathbf{x}}^{*}\|)$.
Using the Schur complement, one can verify that for a given value of
$\operatorname{\mathbf{z}}$,
$\hat{\mathbf{Q}}(\operatorname{\mathbf{z}})\succ\operatorname{\mathbf{0}}$
iff
$\delta>\bar{\delta}(\operatorname{\mathbf{z}})=\dfrac{1}{\lambda_{\min}(\bar{\mathbf{Q}})\lambda_{2}(\operatorname{\mathbf{L}}^{\top}\operatorname{\mathbf{P}}+\operatorname{\mathbf{P}}\operatorname{\mathbf{L}})}(\hat{\mathbf{Q}}_{13}(\operatorname{\mathbf{z}})^{2}+\hat{\mathbf{Q}}_{23}(\operatorname{\mathbf{z}})^{2}).$
Hence, if $\delta>\bar{\delta}(\operatorname{\mathbf{z}}(0))$, then
$\operatorname{\mathcal{L}}_{\operatorname{\psi_{\textnormal{dist}}}}V_{3}\leq-\lambda_{\min}(\hat{\mathbf{Q}}(\operatorname{\mathbf{z}}(0)))\operatorname{\mathbf{z}}^{\top}\operatorname{\mathbf{z}}$.
This along with the fact that
$\frac{1}{2}\min\\{1,\delta\lambda_{\min}(\operatorname{\mathbf{P}})\\}\|\operatorname{\mathbf{z}}\|^{2}\leq
V_{3}\leq\frac{1}{2}\max\\{1,\delta\lambda_{\max}(\operatorname{\mathbf{P}})\\}\|\operatorname{\mathbf{z}}\|^{2}$,
implies that $V_{3}$ satisfies the hypotheses of [21, Theorem 4.10] for
exponential stability. ∎

The exponential convergence of algorithms (5) and (7) for weight-balanced
graphs, and (11) for unbalanced graphs follows from their linear nature. For
algorithm (12), exponential convergence could be attributed to the fact that
the dynamics (12c) converge exponentially and hence, after some time,
(12a)-(12b) and (11) are essentially the same.

12345678910 (a) $\operatorname{\mathcal{G}}_{1}$

12345678910 (b) $\operatorname{\mathcal{G}}_{2}$

Figure 1: Communication topologies among the agents. The edge weights are
adjusted to make the graphs either weight-balanced or unbalanced, as needed.
Figure 2: Evolution of the error between the actual solution and the average
state using the proposed algorithms from initial condition
$\operatorname{\mathbf{x}}(0)=\operatorname{\mathbf{0}}$,
$\mathbf{y}(0)=-\operatorname{\mathbf{b}}$, over the graphs shown in Fig. 1.
The algorithms are implemented in discrete time with a stepsize of $2.5\times
10^{-3}$, and the values of $\alpha=2$, $\beta=0.1$ and $\gamma=20$. Straight
lines correspond to exponential convergence.

##  VI Simulations

We consider 10 agents communicating over the digraphs shown in Fig. 1, seeking
to solve problem (1) with $\\{A_{i}\\}_{i=1}^{10}\in^{5\times 5}$ and
$\\{b_{i}\\}_{i=1}^{10}\in^{5}$. Since the proposed dynamics are in continuous
time, we use a first-order Euler discretization with stepsize $2.5\times
10^{-3}$ for the MATLAB implementation. The edge weights for various cases are
adjusted to make the graphs weight-balanced and unbalanced, resp. For the
time-varying case, at every iteration, the communication graph is switched
randomly between $\operatorname{\mathcal{G}}_{1}$ and
$\operatorname{\mathcal{G}}_{2}$. In Fig. 2, we plot the evolution of the
error between the actual solution of (1) and the average state
$\bar{x}=\frac{1}{n}\operatorname{\mathds{1}}^{\top}\operatorname{\mathbf{x}}$
using (7), (10) and (12). The initial conditions for all the algorithms are
chosen according to Remark 2. Even though $\operatorname{\mathcal{G}}_{2}$
(with $4.6$ as the minimum of the real parts of non-zero eigenvalues of
$\operatorname{\mathsf{L}}$ and
$\lambda_{2}(\operatorname{\mathsf{L}}+\operatorname{\mathsf{L}}^{\top})=7.6$,
for the weight-balanced case) is more connected than
$\operatorname{\mathcal{G}}_{1}$ (with $1.9$ as the minimum of the real parts
of non-zero eigenvalues of $\operatorname{\mathsf{L}}$ and
$\lambda_{2}(\operatorname{\mathsf{L}}+\operatorname{\mathsf{L}}^{\top})=3.8$,
for the weight-balanced case), convergence is slower. The error in the time-
varying case is lower and upper bounded by the error for
$\operatorname{\mathcal{G}}_{1}$ and $\operatorname{\mathcal{G}}_{2}$, resp.

##  VII Conclusions and Future Work

We have presented continuous-time algorithms to solve linear algebraic
equations whose problem data is represented as the summation of the data of
individual agents. The proposed algorithms are distributed over general
directed networks, do not require the individual agent matrices to be positive
definite, and are guaranteed to converge to a solution of the linear equation
exponentially fast. Future work will involve formally characterizing the
convergence when the null space condition is not satisfied, and explore the
design of distributed algorithms for finding least-square solutions when exact
ones do not exist, extension to cases where the problem data is time-varying,
and the communication graph is unbalanced and time-varying.

## References

  * [1] D. P. Bertsekas and J. N. Tsitsiklis, _Parallel and Distributed Computation: Numerical Methods_. Athena Scientific, 1997.
  * [2] S. Mou, J. Liu, and A. S. Morse, “A distributed algorithm for solving a linear algebraic equation,” _IEEE Transactions on Automatic Control_ , vol. 60, no. 11, pp. 2863–2878, 2015.
  * [3] B. D. O. Anderson, S. Mou, A. S. Morse, and U. Helmke, “Decentralized gradient algorithm for solution of a linear equation,” _Numerical Algebra, Control and Optimization_ , vol. 6, no. 3, pp. 319–328, 2016.
  * [4] D. P. Spanos, R. Olfati-Saber, and R. M. Murray, “Distributed sensor fusion using dynamic consensus,” in _IFAC World Congress_ , Prague, CZ, July 2005, electronic proceedings.
  * [5] L. Xiao, S. Boyd, and S. Lall, “A scheme for robust distributed sensor fusion based on average consensus,” in _Symposium on Information Processing of Sensor Networks_ , Los Angeles, CA, Apr. 2005, pp. 63–70.
  * [6] J. Lu and C. Y. Tang, “A distributed algorithm for solving positive definite linear equations over networks with membership dynamics,” _IEEE Transactions on Control of Network Systems_ , vol. 5, no. 1, pp. 215–227, 2018\.
  * [7] X. Wang and S. Mou, “A distributed algorithm for achieving the conservation principle,” in _American Control Conference_ , Milwaukee, WI, June 2018, pp. 5863–5867.
  * [8] D. P. Spanos, R. Olfati-Saber, and R. M. Murray, “Dynamic consensus on mobile networks,” in _IFAC World Congress_ , Prague, Czech Republic, July 2005\.
  * [9] S. S. Kia, B. V. Scoy, J. Cortés, R. A. Freeman, K. M. Lynch, and S. Martinez, “Tutorial on dynamic average consensus: The problem, its applications, and the algorithms,” _IEEE Control Systems_ , vol. 39, no. 3, pp. 40–72, 2019.
  * [10] B. Touri and B. Gharesifard, “Continuous-time distributed convex optimization on time-varying directed networks,” in _IEEE Conf. on Decision and Control_ , Osaka, Japan, 2015, pp. 724–729.
  * [11] ——, “A modified saddle-point dynamics for distributed convex optimization on general directed graphs,” _IEEE Transactions on Automatic Control_ , vol. 65, no. 7, pp. 3098–3103, 2020.
  * [12] P. Srivastava and J. Cortés, “Network optimization via smooth exact penalty functions enabled by distributed gradient computation,” _IEEE Transactions on Control of Network Systems_ , 2021, to appear.
  * [13] F. Bullo, J. Cortés, and S. Martinez, _Distributed Control of Robotic Networks_ , ser. Applied Mathematics Series. Princeton University Press, 2009.
  * [14] Z. Li and D. Zhisheng, _Cooperative Control of Multi-Agent Systems: A Consensus Region Approach_. CRC Press, 2015\.
  * [15] A. Cherukuri and J. Cortés, “Distributed algorithms for convex network optimization under non-sparse equality constraints,” in _Allerton Conf. on Communications, Control and Computing_ , Monticello, IL, Sept. 2016, pp. 452–459.
  * [16] K. Arrow, L. Hurwitz, and H. Uzawa, _Studies in Linear and Non-Linear Programming_. Stanford, CA: Stanford University Press, 1958.
  * [17] A. Cherukuri, B. Gharesifard, and J. Cortés, “Saddle-point dynamics: conditions for asymptotic stability of saddle points,” _SIAM Journal on Control and Optimization_ , vol. 55, no. 1, pp. 486–511, 2017.
  * [18] B. Gharesifard and J. Cortés, “Distributed continuous-time convex optimization on weight-balanced digraphs,” _IEEE Transactions on Automatic Control_ , vol. 59, no. 3, pp. 781–786, 2014.
  * [19] R. A. Horn and C. R. Johnson, _Matrix Analysis_. Cambridge University Press, 1985.
  * [20] S. Boyd and L. Vandenberghe, _Convex Optimization_. Cambridge University Press, 2009.
  * [21] H. K. Khalil, _Nonlinear Systems_ , 3rd ed. Prentice Hall, 2002.