# Byzantine-Resilient Distributed Optimization of
Multi-Dimensional Functions

Kananart Kuwaranancharoen, Lei Xin, Shreyas Sundaram This research was
supported by NSF CAREER award 1653648. The authors are with the School of
Electrical and Computer Engineering at Purdue University. Email:
<EMAIL_ADDRESS>Both of the first two authors contributed
equally to this paper.

###### Abstract

The problem of distributed optimization requires a group of agents to reach
agreement on a parameter that minimizes the average of their local cost
functions using information received from their neighbors. While there are a
variety of distributed optimization algorithms that can solve this problem,
they are typically vulnerable to malicious (or “Byzantine”) agents that do not
follow the algorithm. Recent attempts to address this issue focus on single
dimensional functions, or provide analysis under certain assumptions on the
statistical properties of the functions at the agents. In this paper, we
propose a resilient distributed optimization algorithm for multi-dimensional
convex functions. Our scheme involves two filtering steps at each iteration of
the algorithm: (1) distance-based and (2) component-wise removal of extreme
states. We show that this algorithm can mitigate the impact of up to $F$
Byzantine agents in the neighborhood of each regular node, without knowing the
identities of the Byzantine agents in advance. In particular, we show that if
the network topology satisfies certain conditions, all of the regular states
are guaranteed to asymptotically converge to a bounded region that contains
the global minimizer.

##  I Introduction

The problem of distributed optimization requires a network of agents to reach
agreement on a parameter that minimizes the average of their objective
functions using local information received from their neighbors. This
framework is motivated by various applications including machine learning,
power systems, and robotic networks, and there are a variety of approaches to
solve this problem [1, 2, 3, 4]. However, these existing works typically make
the assumption that all agents follow the prescribed protocol; indeed, such
protocols fail if even a single agent behaves in a malicious or incorrect
manner [5].

A handful of recent papers have considered this problem for the case where
agent misbehavior follows prescribed patterns [6, 7]. A more general (and
serious) form of misbehavior is captured by the Byzantine adversary model from
computer science, where misbehaving agents can send arbitrary (and
conflicting) values to their neighbors at each iteration of the algorithm.
Under such Byzantine behavior, it has been shown that it is impossible to
guarantee computation of the true optimal point [5, 8]. Thus, some papers have
formulated distributed optimization algorithms that allow the non-adversarial
nodes to converge to a certain region surrounding the the true minimizer,
regardless of the adversaries’ actions [5, 9, 10]. However, the above works
focus on single dimensional functions. The extension to general multi-
dimensional functions remains largely open, however, since even the region
containing the true minimizer of the functions is challenging to characterize
in such cases [11]. The recent papers [12, 13] consider a vector version of
the resilient decentralized machine learning problem by utilizing block
coordinate descent. Those papers show that the states of regular nodes will
converge to the statistical minimizer with high probability, but the analysis
is restricted to i.i.d training data across the network.111We note that there
is also a branch of literature pertaining to optimization in a client-server
architecture [14] [15]; this differs from our setting where we consider a
fully distributed network of agents.

In this paper, we propose an algorithm that extends the “local-filtering”
dynamics proposed in [5, 9] (for single-dimensional convex functions) to the
multi-dimensional case. However, this extension is non-trivial, as simply
applying the filtering operations proposed in those papers to each coordinate
of the parameter vector does not appear to yield clear guarantees. Instead, we
show that by having each node apply an additional filtering step on the
parameter vectors that it receives from its neighbors at each step (based on
the distance of those vectors from a commonly chosen reference point), one can
recover certain performance guarantees in the face of Byzantine adversaries.

##  II Notation and Terminology

###  II-A General Notation

Let $\mathbb{R}$ and $\mathbb{N}$ denote the set of real and natural numbers,
respectively. For $N\in\mathbb{N}$, let $[N]$ denote the set
$\\{1,2,\ldots,N\\}$. Vectors are taken to be column vectors, unless otherwise
noted. We use $[x]_{p}$ to represent the $p$-th component of a vector $x$. The
cardinality of a set is denoted by $|\cdot|$, and the Euclidean norm on
$\mathbb{R}^{d}$ is denoted by $\|\cdot\|_{2}$. We denote by $\langle
u,v\rangle$ the Euclidean inner product of $u$ and $v$ i.e., $\langle
u,v\rangle=u^{T}v$ and by $\angle(u,v)$ the angle between vectors $u$ and $v$
i.e., $\angle(u,v)=\arccos\big{(}\frac{\langle
u,v\rangle}{\|u\|_{2}\|v\|_{2}}\big{)}$. The Euclidean ball in $d$-dimensional
space with center at $x_{0}$ and radius $r$ is denoted by
$\mathcal{B}(x_{0},r)\triangleq\\{x\in\mathbb{R}^{d}:\|x-x_{0}\|_{2}\leq
r\\}$. Given a convex function $f:\mathbb{R}^{d}\rightarrow\mathbb{R}$, the
set of subgradients of $f$ at any point $x\in\mathbb{R}^{d}$ is denoted by
$\partial f(x)$.

###  II-B Graph Theory

We denote a network by a directed graph
$\mathcal{G}=(\mathcal{V},\mathcal{E})$, which consists of the set of nodes
$\mathcal{V}=\\{v_{1},v_{2},\ldots,v_{N}\\}$ and the set of edges
$\mathcal{E}\subseteq\mathcal{V}\times\mathcal{V}$. If
$(v_{i},v_{j})\in\mathcal{E}$, then node $v_{j}$ can receive information from
node $v_{i}$. The in-neighbor and out-neighbor sets are denoted by
$\mathcal{N}^{-}_{i}=\\{v_{j}\in\mathcal{V}:\;(v_{j},v_{i})\in\mathcal{E}\\}$
and
$\mathcal{N}^{+}_{i}=\\{v_{j}\in\mathcal{V}:\;(v_{i},v_{j})\in\mathcal{E}\\}$,
respectively. A path from node $v_{i}\in\mathcal{V}$ to node
$v_{j}\in\mathcal{V}$ is a sequence of nodes
$v_{k_{1}},v_{k_{2}},\dots,v_{k_{l}}$ such that $v_{k_{1}}=v_{i}$,
$v_{k_{l}}=v_{j}$ and $(v_{k_{r}},v_{k_{r+1}})\in\mathcal{E}$ for $1\leq r\leq
l-1$. Throughout the paper, the terms nodes and agents will be used
interchangeably.

###### Definition 1

A graph $\mathcal{G}=(\mathcal{V},\mathcal{E})$ is said to be rooted at node
$v_{i}\in\mathcal{V}$ if for for all nodes
$v_{j}\in\mathcal{V}\setminus\\{v_{i}\\}$, there is a path from $v_{i}$ to
$v_{j}$. A graph is said to be rooted if it is rooted at some node
$v_{i}\in\mathcal{V}$.

We will rely on the following definitions from [16].

######  Definition 2 ($r$-reachable set)

For any given $r\in\mathbb{N}$, a subset of nodes
$\mathcal{S}\subseteq\mathcal{V}$ is said to be $r$-reachable if there exists
a node $v_{i}\in\mathcal{S}$ such that
$|\mathcal{N}^{-}_{i}\setminus\mathcal{S}|\geq r$.

######  Definition 3 ($r$-robust graphs)

For $r\in\mathbb{N}$, a graph $\mathcal{G}$ is said to be $r$-robust if for
all pairs of disjoint nonempty subsets
$\mathcal{S}_{1},\mathcal{S}_{2}\subset\mathcal{V}$, at least one of
$\mathcal{S}_{1}$ or $\mathcal{S}_{2}$ is $r$-reachable.

The above definitions capture the idea that sets of nodes should contain a
node that has a sufficient number of neighbors outside that set. This will be
important for the local decisions made by each node in the network under our
algorithm, and will allow information from the rest of the network to
penetrate into different sets of nodes.

###  II-C Adversarial Behavior

###### Definition 4

A node $v_{j}\in\mathcal{V}$ is said to be Byzantine if during each iteration
of the prescribed algorithm, it is capable of sending arbitrary (and perhaps
conflicting) values to different neighbors.

The set of Byzantine nodes is denoted by $\mathcal{A}\subset\mathcal{V}$. The
set of regular nodes is denoted by
$\mathcal{R}=\mathcal{V}\setminus\mathcal{A}$.

The identities of the Byzantine agents are unknown to regular agents in
advance. Furthermore, we allow the Byzantine agents know the entire topology
of the network and functions equipped by the regular nodes (such worst case
behavior is typical in the study of such adversarial models [5, 8, 12]).

######  Definition 5 ($F$-local model)

For $F\in\mathbb{N}$, we say that the set of adversaries $\mathcal{A}$ is an
$F$-local set if $|\mathcal{N}^{-}_{i}\cap\mathcal{A}|\leq F$, for all
$v_{i}\in\mathcal{R}$.

Thus, the $F$-local model captures the idea that each regular node has at most
$F$ Byzantine in-neighbors.

##  III Problem formulation

Consider a group of $N$ agents $\mathcal{V}$ interconnected over a graph
$\mathcal{G}=(\mathcal{V},\mathcal{E})$. Each agent $v_{i}\in\mathcal{V}$ has
a local convex cost function $f_{i}:\mathbb{R}^{d}\rightarrow\mathbb{R}$. The
objective is to collaboratively solve the following minimization problem:

$\displaystyle\mathop{\min}_{x\in\mathbb{R}^{d}}\frac{1}{N}\sum_{v_{i}\in\mathcal{V}}f_{i}(x)$
(1)

where $x$ is the common decision variable. We assume that nodes can only
communicate with their immediate neighbors to solve the above problem.
However, since Byzantine nodes are allowed to send arbitrary values to their
neighbors at each iteration of any algorithm, it is not possible to solve
Problem (1) under such misbehavior (since one is not guaranteed to infer any
information about the true functions of the Byzantine agents) [5, 8]. Thus,
the optimization problem is recast into the following form:

$\displaystyle\mathop{\min}_{x\in\mathbb{R}^{d}}\frac{1}{|\mathcal{R}|}\sum_{v_{i}\in\mathcal{R}}f_{i}(x).$
(2)

We will now propose an algorithm that allows the regular nodes to
approximately solve the above problem (as characterized later in the paper).
We will make the following assumption throughout.

###### Assumption 1

For all $v_{i}\in\mathcal{V}$, the functions $f_{i}(x)$ are convex, and the
sets $\operatorname*{argmin}f_{i}(x)$ are non-empty and bounded.

##  IV A Resilient Distributed Optimization Algorithm

The algorithm that we propose is stated as Algorithm 1. At each time-step $k$,
each node $v_{i}\in\mathcal{V}$ maintains and updates a vector $x_{i}[k]$,
which is its estimate of the solution to Problem (2). After presenting the
algorithm, we describe each of the steps and the update rule.

Algorithm 1 Distance-MinMax Filtering Dynamics

Input Network $\mathcal{G}$, functions $\\{f_{i}\\}_{i=1}^{N}$, the parameter
$F$

1:Each $v_{i}\in\mathcal{R}$ sets $x_{i}^{*}\leftarrow$ optimize($f_{i}$)

2:$\hat{x}\leftarrow$ resilient_consensus ($F,\;\\{x_{i}^{*}\\}$)

3:Each $v_{i}\in\mathcal{R}$ sets $x_{i}[0]=x_{i}^{*}$

4:for $k\in\mathbb{N}$ do

5: for $v_{i}\in\mathcal{R}$ do $\triangleright$ Implement in parallel

6: broadcast($\mathcal{N}_{i}^{+},x_{i}[k]$)

7: $\mathcal{S}_{i}[k]\leftarrow$ receive($\mathcal{N}_{i}^{-}$)

8: $\mathcal{S}_{i}^{\text{dist}}[k]\leftarrow$
dist_filter($F,\;\hat{x},\;\mathcal{S}_{i}[k]$)

9: $\mathcal{S}_{i}^{\text{mm}}[k]\leftarrow$
minmax_filter($F,\;\mathcal{S}_{i}^{\text{dist}}[k]$)

10: $z_{i}[k]\leftarrow$ average($\mathcal{S}_{i}^{\text{mm}}[k]$)

11: $x_{i}[k+1]\leftarrow$ gradient($f_{i},z_{i}[k]$)

12: end for

13:end for

Note that Byzantine nodes do not necessarily need to follow the above
algorithm, and can update their states however they wish. We now explain each
step used in Algorithm 1.

  1. 1.

$x_{i}^{*}\leftarrow$ optimize ($f_{i}$)
Each node $v_{i}\in\mathcal{R}$ finds the minimizer $x_{i}^{*}$ of its local
function $f_{i}$ (using any appropriate algorithm).

  2. 2.

$\hat{x}\leftarrow$ resilient_consensus ($F,\;\\{x_{i}^{*}\\}$)
The nodes run a resilient consensus algorithm to calculate a consensus point
(which we term an auxiliary point) $\hat{x}\in\mathbb{R}^{d}$, with each node
$v_{i}\in\mathcal{R}$ setting its initial vector to be its individual
minimizer $x_{i}^{*}$. For example, $d$ parallel versions of the the resilient
scalar consensus algorithm from [16] can be applied (one for each component of
the parameter vector). This is guaranteed to return a consensus value
$\hat{x}$ that is in the smallest hyperrectangle containing all of the
minimizers of the regular nodes’ functions, regardless of the actions of any
$F$-local set of adversaries (under the network conditions provided in the
next section) [16].

  3. 3.

broadcast ($\mathcal{N}_{i}^{+},x_{i}[k]$)
Node $v_{i}\in\mathcal{R}$ broadcasts its current state $x_{i}[k]$ to its out-
neighbors $\mathcal{N}_{i}^{+}$.

  4. 4.

$\mathcal{S}_{i}[k]\leftarrow$ receive($\mathcal{N}_{i}^{-}$)
Node $v_{i}\in\mathcal{R}$ receives the current states from its in-neighbors
$\mathcal{N}_{i}^{-}$. So, at time step $k$, node $v_{i}$ possesses the set of
states
$\mathcal{S}_{i}[k]\triangleq\\{x_{j}[k]:j\in\mathcal{N}_{i}^{-}\\}\cup\\{x_{i}[k]\\}$.

  5. 5.

$\mathcal{S}_{i}^{\text{dist}}[k]\leftarrow$ dist_filter
($F,\;\hat{x},\;\mathcal{S}_{i}[k]$)
Node $v_{i}\in\mathcal{R}$ computes

$\displaystyle
D_{ij}[k]\triangleq\|x_{j}[k]-\hat{x}\|_{2}\;\;\text{for}\;\;x_{j}[k]\in\mathcal{S}_{i}[k]\setminus\\{x_{i}[k]\\}.$
(3)

Then, node $v_{i}$ removes the states
$x_{j}[k]\in\mathcal{S}_{i}[k]\setminus\\{x_{i}[k]\\}$ that produce the
$F$-highest values of $D_{ij}[k]$. The remaining states $x_{j}[k]$ are stored
in $\mathcal{S}_{i}^{\text{dist}}[k]$.

  6. 6.

$\mathcal{S}_{i}^{\text{mm}}[k]\leftarrow$ minmax_filter
($F,\;\mathcal{S}_{i}^{\text{dist}}[k]$)
Node $v_{i}$ further removes states that have extreme values in any of their
components. More specifically, node $v_{i}\in\mathcal{R}$ removes the state
$x_{j}[k]\in\mathcal{S}_{i}^{\text{dist}}[k]\setminus\\{x_{i}[k]\\}$ if there
exists $p\in[d]$ such that $[x_{j}[k]]_{p}$ is in the $F$-highest or
$F$-lowest values from the set of scalars
$\big{\\{}[x_{l}[k]]_{p}:x_{l}[k]\in\mathcal{S}_{i}^{\text{dist}}[k]\setminus\\{x_{i}[k]\\}\big{\\}}$.
The remaining states are stored in $\mathcal{S}_{i}^{\text{mm}}[k]$.

  7. 7.

$z_{i}[k]\leftarrow$ average ($\mathcal{S}_{i}^{\text{mm}}[k]$)
Each node $v_{i}\in\mathcal{R}$ computes

$z_{i}[k]=\sum_{v_{j}\in\mathcal{N}_{i}^{\text{mm}}[k]}\frac{1}{|\mathcal{S}_{i}^{\text{mm}}[k]|}x_{j}[k]$
(4)

where
$\mathcal{N}_{i}^{\text{mm}}[k]\triangleq\\{v_{s}:x_{s}[k]\in\mathcal{S}_{i}^{\text{mm}}[k]\\}$.

  8. 8.

$x_{i}[k+1]\leftarrow$ gradient ($f_{i},z_{i}[k]$)
Node $v_{i}\in\mathcal{R}$ computes the gradient update as follows:

$\displaystyle x_{i}[k+1]=z_{i}[k]-\eta[k]g_{i}[k]$ (5)

where $g_{i}[k]\in\partial f_{i}(z_{i}[k])$ and $\eta[k]$ is the step-size at
time-step $k$.

###### Remark 1

The role of the auxiliary point $\hat{x}$ is to give the states $x_{i}[k]$ a
“sense of good direction”. In other words, since adversarial nodes can try to
pull the regular nodes away from the true minimizer, the auxiliary point
provides a common reference point for each regular node with which to evaluate
the states in its neighborhood. This motivates the distance-based filtering
step of Algorithm 1 (in line 8), which removes the $F$ states that are
furthest away from the auxiliary point at each time-step. Note that the
auxiliary point will, in general, be different from the true minimizer of
interest. We also note that the resilient consensus algorithm from [16] for
computing the auxiliary point will only provide asymptotic consensus in
general. The regular agents will converge to consensus exponentially fast
under that algorithm, however, and thus we expect that running the resilient
consensus algorithm to compute the auxiliary point simultaneously with the
other update steps (given by lines 6-11) will also lead to the same
guarantees. In the interest of space, we do not provide that analysis here.

##  V Assumptions and Main Results

We will make the following assumptions in our analysis.

###### Assumption 2

There exists $L>0$ such that $\|g_{i}(x)\|_{2}\leq L$ for all
$x\in\mathbb{R}^{d}$ and $v_{i}\in\mathcal{V}$, where $g_{i}(x)\in\partial
f_{i}(x)$.

###### Assumption 3

The step-sizes used in line 11 of Algorithm 1 satisfy
$\lim_{k\to\infty}\eta[k]=0$, $\eta[k+1]<\eta[k]$ for all $k$, and
$\sum_{k=0}^{\infty}\eta[k]=\infty$.

###### Assumption 4

The underlying communication graph $\mathcal{G}$ is
$\big{(}(2d+1)F+1\big{)}$-robust, and the Byzantine agents form a $F$-local
set.

We will also use the following lemma from [16].

###### Lemma 1

Suppose a graph $\mathcal{G}$ satisfies Assumption 4. Let
$\mathcal{G}^{\prime}$ be a graph obtained by removing $(2d+1)F$ or fewer
incoming edges from each node in $\mathcal{G}$. Then $\mathcal{G}^{\prime}$ is
rooted.

###  V-A Convergence to Consensus

We first show that the states $x_{i}[k]$ of all regular nodes
$v_{i}\in\mathcal{R}$ reach consensus under Algorithm 1.

######  Theorem 1 (Consensus)

Under Assumptions 2, 3, and 4,
$\lim_{k\rightarrow\infty}\|x_{i}[k]-x_{j}[k]\|=0$ for all
$v_{i},v_{j}\in\mathcal{R}$.

###### Proof:

We will argue that all regular nodes $v_{i}\in\mathcal{R}$ reach consensus on
each component of their vectors $x_{i}[k]$, which will then prove the result.
For all $p\in[d]$ and for all $v_{i}\in\mathcal{R}$, from (5), the $p$-th
component of the vector $x_{i}[k]$ evolves as

$[x_{i}[k+1]]_{p}=[z_{i}[k]]_{p}-\eta[k][g_{i}[k]]_{p}.$ (6)

From (4), the quantity $z_{i}[k]$ is an average of a subset of the parameter
vectors from node $v_{i}$’s neighborhood. In particular, the set
$\mathcal{S}_{i}^{mm}[k]$ is obtained by removing at most $(2d+1)F$ of the
vectors received from $v_{i}$’s neighbors ($F$ vectors removed by the distance
based filtering in line 8, and up to $2F$ additional vectors removed by the
minmax filtering step on each of the $d$ components in line 9 of the
algorithm). Thus, at each time-step $k$, component $[z_{i}[k]]_{p}$ is an
average of at least $|\mathcal{N}_{i}^{-}|-(2d+1)F$ of $v_{i}$’s neighbors
values on that component. Since the graph is $((2d+1)F+1)$-robust and the
Byzantine agents form an $F$-local set (by Assumption 4), and leveraging Lemma
1 and the fact that the term $\eta[k]g_{i}[k]$ asymptotically goes to zero (by
Assumptions 2 and 3), we can use an identical argument as in Theorem 6.1 from
[5] to show that the scalar dynamics (6) converge to consensus. ∎

###  V-B What Region Do the States Converge To?

We now analyze the trajectories of the states of the agents under Algorithm 1.
We start with the following result about the quantity $z_{i}[k]$ calculated in
line 10 of the algorithm (and given by equation (4)).

###### Proposition 1

For all $k\in\mathbb{N}$ and $v_{i}\in\mathcal{R}$, if there exists
$R_{i}[k]\in\mathbb{R}_{\geq 0}$ such that $\|x_{j}[k]-\hat{x}\|_{2}\leq
R_{i}[k]$ for all
$v_{j}\in(\mathcal{N}_{i}^{-}\cap\mathcal{R})\cup\\{v_{i}\\}$ then
$\|z_{i}[k]-\hat{x}\|_{2}\leq R_{i}[k]$.

###### Proof:

Consider the step $\mathcal{S}_{i}^{\text{dist}}[k]\leftarrow$ dist_filter
($F,\;\hat{x},\;\mathcal{S}_{i}[k]$) in Algorithm 1. We will first prove the
following claim. For each $v_{i}\in\mathcal{R}$, there exists
$v_{r}\in(\mathcal{N}_{i}^{-}\cap\mathcal{R})\cup\\{v_{i}\\}$ such that
$\|x_{j}[k]-\hat{x}\|_{2}\leq\|x_{r}[k]-\hat{x}\|_{2}$ for all
$v_{j}\in\\{v_{s}:x_{s}[k]\in\mathcal{S}_{i}^{\text{dist}}[k]\\}$.

There are two possible cases. First, if the set
$\mathcal{S}_{i}^{\text{dist}}[k]$ contains only regular nodes, we can simply
choose $v_{r}\in(\mathcal{N}_{i}^{-}\cap\mathcal{R})\cup\\{v_{i}\\}$ to be the
node whose state $x_{r}[k]$ is furthest away from $\hat{x}$. Next, consider
the case where $\mathcal{S}_{i}^{\text{dist}}[k]$ contains the states of one
or more Byzantine nodes. Since node $v_{i}\in\mathcal{R}$ removes the $F$
states from $\mathcal{N}_{i}^{-}$ that are furthest away from $\hat{x}$ (in
line 8 of the algorithm), and there are at most $F$ Byzantine nodes in
$\mathcal{N}_{i}^{-}$, there is at least one regular state removed by the node
$v_{i}$. Let $v_{r}$ be one of the regular nodes whose state is removed. We
then have $D_{ir}[k]\geq D_{ij}[k]$, for all
$v_{j}\in\\{v_{s}:x_{s}[k]\in\mathcal{S}_{i}^{\text{dist}}[k]\\}$ which proves
the claim.

Consider the step $\mathcal{S}_{i}^{\text{mm}}[k]\leftarrow$ minmax_filter
($F,\;\mathcal{S}_{i}^{\text{dist}}[k]$) in Algorithm 1. We have that
$\mathcal{S}_{i}^{\text{mm}}[k]\subset\mathcal{S}_{i}^{\text{dist}}[k]$. Then,
consider the step $z_{i}[k]\leftarrow$ average
($\mathcal{S}_{i}^{\text{mm}}[k]$). We have

$z_{i}[k]-\hat{x}=\sum_{v_{j}\in\mathcal{N}_{i}^{\text{mm}}[k]}\frac{1}{|\mathcal{S}_{i}^{\text{mm}}[k]|}\big{(}x_{j}[k]-\hat{x}\big{)}.$

Since $\|x_{j}[k]-\hat{x}\|_{2}\leq\|x_{r}[k]-\hat{x}\|_{2}$ for all
$v_{j}\in\mathcal{N}_{i}^{\text{mm}}[k]$ (where $v_{r}$ is the node identified
in the claim at the start of the proof), we obtain

$\|z_{i}[k]-\hat{x}\|_{2}\leq\sum_{v_{j}\in\mathcal{N}_{i}^{\text{mm}}[k]}\frac{1}{|\mathcal{S}_{i}^{\text{mm}}[k]|}\|x_{j}[k]-\hat{x}\|_{2}\leq
R_{i}[k].$

∎

Next, we will establish certain quantities that will be useful for our
analysis of the convergence region. Since the set
$\operatorname*{argmin}f_{i}(x)$ is non-empty for all $v_{i}\in\mathcal{V}$
(by Assumption 1), we define $x_{i}^{*}\in\mathbb{R}^{d}$ to be a minimizer of
$f_{i}(x)$ for all $v_{i}\in\mathcal{V}$. For $\epsilon>0$, define

$\displaystyle\mathcal{C}_{i}(\epsilon)\triangleq\\{x\in\mathbb{R}^{d}:f_{i}(x)\leq
f_{i}(x_{i}^{*})+\epsilon\\}.$ (7)

Since the set $\operatorname*{argmin}f_{i}(x)$ is bounded for all
$v_{i}\in\mathcal{V}$, there exists $\delta_{i}(\epsilon)\in(0,\infty)$ such
that
$\mathcal{C}_{i}(\epsilon)\subseteq\mathcal{B}(x_{i}^{*},\delta_{i}(\epsilon))$
for all $v_{i}\in\mathcal{V}$.

###### Proposition 2

Consider a convex function $f$ that has bounded subgradients, and suppose the
set $\operatorname*{argmin}f(x)$ is non-empty and bounded. Then for all
$\epsilon>0$ there exists $\theta(\epsilon)\geq 0$ such that
$\angle(-g(x),\;x^{*}-x)\leq\theta(\epsilon)<\frac{\pi}{2}$ for all
$x\notin\mathcal{C}(\epsilon)$, $g(x)\in\partial f(x)$, and
$x^{*}\in\operatorname*{argmin}f(x)$, where $\mathcal{C}(\epsilon)$ is defined
in the same way as (7).

###### Proof:

From the definition of convex functions, for any $x$, $y\in\mathbb{R}^{d}$, we
have $f(y)\geq f(x)+\langle g(x),\;y-x\rangle$, where $g(x)\in\partial f(x)$.
Substitute a minimizer $x^{*}$ of the function $f$ into the variable $y$ to
get

$\displaystyle-\langle g(x),\;x^{*}-x\rangle\geq f(x)-f(x^{*}).$ (8)

Let $\hat{\theta}(x)\triangleq\angle(g(x),x-x^{*})$. The inequality (8)
becomes

$\displaystyle\|g(x)\|_{2}\;\|x^{*}-x\|_{2}\cos\hat{\theta}(x)$
$\displaystyle\geq f(x)-f(x^{*}).$

Fix $\epsilon>0$, and suppose that $x\notin\mathcal{C}(\epsilon)$. Applying
$\|g(x)\|_{2}\leq L$, we have

$\displaystyle\cos\hat{\theta}(x)\geq\frac{f(x)-f(x^{*})}{\|g(x)\|_{2}\;\|x^{*}-x\|_{2}}\geq\frac{f(x)-f(x^{*})}{L\;\|x^{*}-x\|_{2}}.$
(9)

Let $\tilde{x}$ be the point on the line connecting $x^{*}$ and $x$ such that
$f(\tilde{x})=f(x^{*})+\epsilon$. We can rewrite the point $x$ as

$\displaystyle x=x^{*}+t(\tilde{x}-x^{*})\quad\text{where}\quad
t=\frac{\|x-x^{*}\|_{2}}{\|\tilde{x}-x^{*}\|_{2}}\geq 1.$

Suppose $\mathcal{C}(\epsilon)\subseteq\mathcal{B}(x^{*},\delta(\epsilon))$.
Consider the term on the RHS of (9). We have

$\displaystyle\frac{f(x)-f(x^{*})}{\|x-x^{*}\|_{2}}$
$\displaystyle=\frac{f(x^{*}+t(\tilde{x}-x^{*}))-f(x^{*})}{t\|\tilde{x}-x^{*}\|_{2}}$
$\displaystyle\geq\frac{f(x^{*}+t(\tilde{x}-x^{*}))-f(x^{*})}{t\max_{x\in\mathcal{C}(\epsilon)}\|x-x^{*}\|_{2}}$
$\displaystyle\geq\frac{f(x^{*}+t(\tilde{x}-x^{*}))-f(x^{*})}{t\delta(\epsilon)}.$
(10)

Since the quantity $\frac{f(x^{*}+t(\tilde{x}-x^{*}))-f(x^{*})}{t}$ is non-
decreasing in $t\in(0,\infty)$ [17, Lemma 2.80], the inequality (10) becomes

$\displaystyle\frac{f(x)-f(x^{*})}{\|x-x^{*}\|_{2}}\geq\frac{f(\tilde{x})-f(x^{*})}{\delta(\epsilon)}=\frac{\epsilon}{\delta(\epsilon)}.$
(11)

Therefore, combining (9) and (11), we obtain

$\displaystyle\cos\hat{\theta}(x)\geq\frac{\epsilon}{L\delta(\epsilon)}.$ (12)

However, from the definition of convex functions, we have

$\displaystyle f(x^{*})\geq f(\tilde{x})+\langle
g(\tilde{x}),\;x^{*}-\tilde{x}\rangle.$

Since $\|g(\tilde{x})\|_{2}\leq L$ and
$\|x^{*}-\tilde{x}\|_{2}\leq\delta(\epsilon)$, we get

$\displaystyle\epsilon=f(\tilde{x})-f(x^{*})\leq-\langle
g(\tilde{x}),\;x^{*}-\tilde{x}\rangle\leq L\delta(\epsilon).$

Since $\epsilon>0$ and $\frac{\epsilon}{L\delta(\epsilon)}\leq 1$ from the
above inequality, the inequality (12) becomes

$\displaystyle\hat{\theta}(x)\leq\arccos\Big{(}\frac{\epsilon}{L\delta(\epsilon)}\Big{)}\triangleq\theta(\epsilon)<\frac{\pi}{2}.$

∎

From Proposition 2, if $f_{i}$ satisfies Assumptions 1 and 2 for all
$v_{i}\in\mathcal{R}$, then for all $\epsilon>0$, we have
$\angle(-g_{i}(x),\;x_{i}^{*}-x)\leq\theta_{i}(\epsilon)<\frac{\pi}{2}$ for
all $x\notin\mathcal{C}_{i}(\epsilon)$ and $v_{i}\in\mathcal{R}$.

Define $\tilde{R}_{i}\triangleq\|x_{i}^{*}-\hat{x}\|_{2}$ and

$\displaystyle
R^{*}\triangleq\inf_{\epsilon>0}\Big{\\{}\max_{v_{i}\in\mathcal{R}}\big{\\{}\max\\{\tilde{R}_{i}\sec\theta_{i}(\epsilon),\;\tilde{R}_{i}+\delta_{i}(\epsilon)\\}\big{\\}}\Big{\\}}.$
(13)

We now come to the main result of this paper, showing that the states of all
the regular nodes will asymptotically converge to a ball of radius $R^{*}$
around the auxiliary point $\hat{x}$ under Algorithm 1.

######  Theorem 2 (Convergence)

If the states are updated using Algorithm 1, and Assumptions 1, 2, 3 and 4
hold, then for all $v_{i}\in\mathcal{R}$,
$\limsup_{k}\|x_{i}[k]-\hat{x}\|_{2}\leq R^{*}$, regardless of the actions of
any $F$-local set of Byzantine adversaries.

The full proof of the theorem is somewhat long, and is thus provided in the
Appendix. We present here a sketch of the proof. We work towards the proof of
Theorem 2 in several steps. For a fixed $\epsilon>0$ and for all
$\xi\in\mathbb{R}_{>0}$, define the convergence radius

$\displaystyle
s^{*}(\xi,\epsilon)\triangleq\max_{v_{i}\in\mathcal{R}}\big{\\{}\max\\{\tilde{R}_{i}\sec\theta_{i}(\epsilon),\tilde{R}_{i}+\delta_{i}(\epsilon)\\}\big{\\}}+\xi.$

Throughout the paper, we will omit the dependence of $\epsilon$ in
$s^{*}(\xi,\epsilon)$ for notational simplicity.

The key idea is that for a fixed $\xi>0$, we partition the space
$\mathbb{R}^{d}$ into $3$ regions:

  1. 1.

$\mathcal{B}(\hat{x},\;\max_{v_{i}\in\mathcal{R}}\\{\tilde{R}_{i}+\delta_{i}\\})$,

  2. 2.

$\mathcal{B}(\hat{x},\;s^{*}(\xi))\setminus\mathcal{B}(\hat{x},\;\max_{v_{i}\in\mathcal{R}}\\{\tilde{R}_{i}+\delta_{i}\\})$,
and

  3. 3.

$\mathbb{R}^{d}\setminus\mathcal{B}(\hat{x},\;s^{*}(\xi))$.

We refer to $\mathcal{B}(\hat{x},\;s^{*}(\xi))$ as the convergence region and
note that $R^{*}=\inf_{\epsilon>0}s^{*}(0)$. In Steps 1, 2 and 3, we analyze
the gradient update (5) for each regular agent $v_{i}\in\mathcal{R}$, and in
Step 4, we consider the sequence of distance updates of the regular agent
furthest from the auxiliary point. Finally, we will take $\xi$ to go to zero
and take the infimum of the distance bound over $\epsilon>0$.

Step 1: First, we show that if $k$ is sufficiently large and the state
$z_{i}[k]$ is in Region 1 i.e.,

$\displaystyle
z_{i}[k]\in\mathcal{B}(\hat{x},\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\})\subset\mathcal{B}(\hat{x},s^{*}(\xi)),$

then by applying the gradient update (5), the state $x_{i}[k+1]$ is still in
the convergence region i.e., $x_{i}[k+1]\in\mathcal{B}(\hat{x},s^{*}(\xi))$.
To do this, we leverage the fact that the magnitude of the gradient step
satisfies $\eta[k]\|g_{i}[k]\|_{2}\leq\eta[k]L$ by Assumption 2 and $\eta[k]L$
decreases as $k$ increases by Assumption 3.

Step 2: We find the relationship between the terms
$\|x_{i}[k+1]-\hat{x}\|_{2}$ and $\|z_{i}[k]-\hat{x}\|_{2}$ which will be used
in the subsequent step. Specifically, for $v_{i}\in\mathcal{R}$, if the state
$z_{i}[k]$ is in region 2 or 3, i.e.,
$\|z_{i}[k]-\hat{x}\|_{2}>\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\}$,
then we show

$\|x_{i}[k+1]-\hat{x}\|_{2}^{2}\leq\|z_{i}[k]-\hat{x}\|_{2}^{2}\\\
-\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\;\eta[k]\;\|g_{i}[k]\|_{2})$ (14)

where we define
$\Delta_{i}:[\tilde{R}_{i},\infty)\times\mathbb{R}_{+}\to\mathbb{R}$ to be the
function

$\displaystyle\Delta_{i}(p,l)\triangleq
2l\Big{(}\sqrt{p^{2}-\tilde{R}_{i}^{2}}\cos\theta_{i}-\tilde{R}_{i}\sin\theta_{i}\Big{)}-l^{2}.$

Step 3: We show that if $k$ is sufficiently large and the state $z_{i}[k]$ is
in Region 2 i.e.,
$\|z_{i}[k]-\hat{x}\|_{2}\in\big{(}\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\},s^{*}(\xi)\big{]}$
then by applying the gradient update (5), we have that the state $x_{i}[k+1]$
is still in the convergence region i.e.,
$x_{i}[k+1]\in\mathcal{B}(\hat{x},s^{*}(\xi))$. This is done by showing that
small $\eta[k]L$ makes the RHS of (14) bounded above by
$\big{(}s^{*}(\xi)\big{)}^{2}$.

From Step 1 and 3, and Proposition 1, we can conclude that for sufficiently
large $k$, for each regular node $v_{i}\in\mathcal{R}$, if the state
$x_{i}[k]$ is in the convergence region, then the state $x_{i}[k+1]$ is still
in the convergence region, i.e., the regular states cannot leave the
convergence region.

Step 4: We show that the states of all regular nodes eventually enter the
convergence region. Specifically, from Proposition 1, we have that

$\displaystyle\max_{v_{i}\in\mathcal{R}}\|z_{i}[k]-\hat{x}\|_{2}\leq\max_{v_{i}\in\mathcal{R}}\|x_{i}[k]-\hat{x}\|_{2}.$
(15)

On the other hand, for the agents whose state $z_{i}[k]$ is in Region 3, we
show that the term
$\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\;\eta[k]\;\|g_{i}[k]\|_{2})$ in (14) is
bounded below as

$\displaystyle\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\eta[k]\;\|g_{i}[k]\|_{2})>c_{i}\eta[k]$
(16)

for sufficiently large $k$, where $c_{i}$ is a positive constant.
Incorporating (15) and (16) into (14), we can conclude that the quantity
$\max_{v_{i}\in\mathcal{R}}\|x_{i}[k]-\hat{x}\|_{2}^{2}$ strictly decreases if
it is greater than $s^{*}(\xi)$. Eventually, every state will enter the
convergence region, which completes the proof.

To gain insight into the convergence region, we provide the following result.

###### Proposition 3

Let $x^{*}$ be a solution of Problem (2). If Assumptions 1 and 2 hold, then
$x^{*}\in\mathcal{B}(\hat{x},R^{*})$ where $R^{*}$ is defined in (13).

###### Proof:

We will show that the summation of any subgradients of the regular nodes’
functions at any point outside the region $\mathcal{B}(\hat{x},R^{*})$ cannot
be zero.

Let $x_{0}$ be a point outside $\mathcal{B}(\hat{x},R^{*})$. Since
$\|x_{0}-\hat{x}\|_{2}>\max_{v_{i}\in\mathcal{R}}\\{\tilde{R}_{i}+\delta_{i}(\epsilon)\\}$
for some $\epsilon>0$, we have that $x_{0}\notin\mathcal{C}_{i}(\epsilon)$ for
all $v_{i}\in\mathcal{R}$. By the definition of $\mathcal{C}_{i}(\epsilon)$ in
(7), we have $f_{i}(x_{0})>f_{i}(x_{i}^{*})+\epsilon$ for all
$v_{i}\in\mathcal{R}$. Since the functions $f_{i}$ are convex, we obtain
$g_{i}(x_{0})\neq\mathbf{0}$ for all $v_{i}\in\mathcal{R}$ where
$g_{i}(x_{0})\in\partial f_{i}(x_{0})$.

Consider the angle between the vectors $x_{0}-x_{i}^{*}$ and $x_{0}-\hat{x}$.
Suppose $\tilde{R}_{i}>0$; otherwise, we have
$\angle(x_{0}-x_{i}^{*},\;x_{0}-\hat{x})=0$. Using Lemma 2 in the Appendix, we
can bound the angle as follows:

$\displaystyle\angle(x_{0}-x_{i}^{*},\;x_{0}-\hat{x})$
$\displaystyle\leq\max_{y\in\mathcal{B}(\hat{x},\tilde{R}_{i})}\angle(x_{0}-y,\;x_{0}-\hat{x})$
$\displaystyle=\arcsin\Big{(}\frac{\tilde{R}_{i}}{\|x_{0}-\hat{x}\|_{2}}\Big{)}.$

Since
$\|x_{0}-\hat{x}\|_{2}>\max_{v_{i}\in\mathcal{R}}\\{\tilde{R}_{i}\sec\theta_{i}\\}$
and $\arcsin(x)$ is an increasing function in $x\in[-1,1]$, we have

$\displaystyle\angle(x_{0}-x_{i}^{*},\;x_{0}-\hat{x})$
$\displaystyle<\arcsin\Big{(}\frac{\tilde{R}_{i}}{\tilde{R}_{i}\sec\theta_{i}}\Big{)}$
$\displaystyle=\arcsin(\cos\theta_{i})$
$\displaystyle=\frac{\pi}{2}-\theta_{i}.$

Using Proposition 2 and the inequality above, we can bound the angle between
the vectors $g_{i}(x_{0})$ and $x_{0}-\hat{x}$ as follows:

$\displaystyle\quad\angle(g_{i}(x_{0}),\;x_{0}-\hat{x})$
$\displaystyle\leq\angle(g_{i}(x_{0}),\;x_{0}-x_{i}^{*})+\angle(x_{0}-x_{i}^{*},\;x_{0}-\hat{x})$
$\displaystyle<\theta_{i}+\Big{(}\frac{\pi}{2}-\theta_{i}\Big{)}=\frac{\pi}{2}.$

Note that the first inequality is obtained from [18, Corollary 12]. Let
$u=\frac{x_{0}-\hat{x}}{\|x_{0}-\hat{x}\|_{2}}$. Compute the inner product

$\displaystyle\Big{\langle}\sum_{v_{i}\in\mathcal{R}}g_{i}(x_{0}),\;u\Big{\rangle}$
$\displaystyle=\sum_{v_{i}\in\mathcal{R}}\langle g_{i}(x_{0}),\;u\rangle$
$\displaystyle=\sum_{v_{i}\in\mathcal{R}}\|g_{i}(x_{0})\|_{2}\cos\angle(g_{i}(x_{0}),\;x_{0}-\hat{x})$
$\displaystyle>0$

since $\|g_{i}(x_{0})\|_{2}>0$ and
$\cos\angle(g_{i}(x_{0}),\;x_{0}-\hat{x})>0$ for any $v_{i}\in\mathcal{R}$.
This implies that $\sum_{v_{i}\in\mathcal{R}}g_{i}(x_{0})\neq\mathbf{0}$.
Since we can arbitrarily choose $g_{i}(x_{0})$ from the set $\partial
f_{i}(x_{0})$, we have $\mathbf{0}\notin\partial f(x_{0})$ where
$f(x)=\frac{1}{|\mathcal{R}|}\sum_{v_{i}\in\mathcal{R}}f_{i}(x)$. ∎

Thus, Theorem 2 and Proposition 3 show that Algorithm 1 causes all regular
nodes to converge to a region that also contains the true solution, regardless
of the actions of any $F$-local set of Byzantine adversaries. The size of this
region scales with the quantity $R^{*}$. Loosely speaking, this quantity
becomes smaller as the minimizers of the local functions of the regular agents
get closer together. More specifically, consider a fixed $\epsilon>0$. If the
functions $f_{i}(x)$ are translated so that the minimizers $x_{i}^{*}$ get
closer together (i.e., $\tilde{R}_{i}$ is smaller and $\theta(\epsilon)$ is
fixed), and the auxiliary point $\hat{x}$ is in the hyperrectangle containing
the minimizers (which is the case when we run a resilient consensus algorithm
such as the one in [16]) then $R^{*}$ also decreases, and the state $x_{i}[k]$
is guaranteed to become closer to the true minimizer as $k$ goes to infinity.

##  VI Numerical Experiment

Figure 1: Function value $f(x)$ evaluated at the average of regular agents’
states (blue solid line), at the regular agent’s state that gives the lowest
and highest value (orange and yellow solid line respectively), at the
auxiliary point (blue dashed line), and at the minimizer (red dashed line)

In our numerical experiment, we set the total number of agents to be $n=100$
and the number of dimensions to be $d=3$. We construct a $15$-robust network
in order to tolerate up to $2$ Byzantine agents in the neighborhood of any
regular nodes i.e., $F=2$. Each regular node possesses a quadratic function
$f_{i}(x)=\frac{1}{2}x^{T}Q_{i}x+b_{i}^{T}x$ where $Q_{i}$ and $b_{i}$ are
randomly chosen with $Q_{i}$ guaranteed to be positive definite.222To satisfy
Assumption 2 (bounded gradients), we saturate the gradients during the updates
when the norms are sufficiently large. Each Byzantine agent selects the
transmitted vector based on the target regular node. Specifically, the
transmitted vector is picked uniformly at random inside the hyper-rectangle
that guarantees that its state is not discarded by the target node.

Let $f(x)\triangleq\frac{1}{|\mathcal{R}|}\sum_{v_{i}\in\mathcal{R}}f_{i}(x)$
be the objective function (Problem 2) evaluated at $x$,
$f^{*}\triangleq\min_{x\in\mathbb{R}^{d}}f(x)$ be the optimal value of the
objective function and
$\bar{x}[k]\triangleq\frac{1}{|\mathcal{R}|}\sum_{v_{i}\in\mathcal{R}}x_{i}[k]$
be the average of the regular nodes’ states at time-step $k$. In Figure 1, the
objective value evaluated at the average of all regular states $f(\bar{x}[k])$
is near the optimal value $f^{*}$ after $40$ iterations. Furthermore, the
maximum of the objective values among all regular agents
$\max_{v_{i}\in\mathcal{R}}f(x_{i}[k])$ and the minimum of the objective
values among all regular agents $\min_{v_{i}\in\mathcal{R}}f(x_{i}[k])$
converge to the same value, which is an implication of the fact that the
regular agents reach consensus (from Theorem 1).

###### Remark 2

The bound in Theorem 2, for large $k$, does not preclude
$\|\bar{x}[k]-x^{*}\|_{2}>\|\hat{x}-x^{*}\|_{2}$, where $x^{*}$ is a minimizer
of $f$, and $f(\bar{x}[k])>f(\hat{x})$. However, our experiments (with various
settings), including the one above, show that the objective function value
evaluated at the average of the states is much closer to the optimal value
than that of the auxiliary point in practice, i.e., $f(\bar{x}[k])-f^{*}\ll
f(\hat{x})-f^{*}$.

##  VII Conclusion and Future work

In this paper, we developed a resilient distributed optimization algorithm for
multi-dimensional functions. Our results guarantee that the regular states
asymptotically reach consensus and enter a bounded region that contains the
global minimizer, irrespective of the actions of Byzantine agents, if the
network topology satisfies certain conditions. We characterized the size of
this region. A promising avenue for future research would be to further refine
the size of the convergence region, and to relax the conditions on the network
topology.

###  -A Additional Lemma

We provide a lemma which will be utilized in Step 2 of the proof of Theorem 2.

###### Lemma 2

If $x\notin\mathcal{B}(\hat{x},R)$ then

$\displaystyle\max_{y\in\mathcal{B}(\hat{x},R)}\angle(x-y,\;x-\hat{x})=\arcsin\Big{(}\frac{R}{\|x-\hat{x}\|_{2}}\Big{)}.$

###### Proof:

Since the angle is measured with respect to the vector $x-\hat{x}$, consider
any 2-D planes passed through the center $\hat{x}$ and the point $x$. Since
the planes pass through $\hat{x}$, the intersections between of the ball
$\mathcal{B}(\hat{x},R)$ and the planes are great circles of radius $R$. Thus,
all of the intersections generated from each plane are identical and we can
consider the angle using a great circle instead of the ball. From geometry,
the maximum angle $\phi=\angle(x-y^{*},\;x-\hat{x})$ only occurs when the ray
starting from the point $x$ touches the circle at point $y^{*}$. Therefore,
$\angle(\hat{x}-y^{*},x-y^{*})=\frac{\pi}{2}$ and $\|\hat{x}-y^{*}\|_{2}=R$.
We have

$\displaystyle\sin\phi=\frac{\|\hat{x}-y^{*}\|_{2}}{\|\hat{x}-x\|_{2}}=\frac{R}{\|\hat{x}-x\|_{2}}$

and the result follows. ∎

###  -B Proof of Theorem 2

We work towards the proof of Theorem 2 in several steps. For a fixed
$\epsilon>0$ and for all $\xi\in\mathbb{R}_{>0}$, define the convergence
radius

$\displaystyle
s^{*}(\xi)\triangleq\max_{v_{i}\in\mathcal{R}}\big{\\{}\max\\{\tilde{R}_{i}\sec\theta_{i},\tilde{R}_{i}+\delta_{i}\\}\big{\\}}+\xi.$

The key idea is that we partition the space $\mathbb{R}^{d}$ into $3$ regions:

  1. 1.

$\mathcal{B}(\hat{x},\;\max_{v_{i}\in\mathcal{R}}\\{\tilde{R}_{i}+\delta_{i}\\})$,

  2. 2.

$\mathcal{B}(\hat{x},\;s^{*}(\xi))\setminus\mathcal{B}(\hat{x},\;\max_{v_{i}\in\mathcal{R}}\\{\tilde{R}_{i}+\delta_{i}\\})$,
and

  3. 3.

$\mathbb{R}^{d}\setminus\mathcal{B}(\hat{x},\;s^{*}(\xi))$.

We refer to $\mathcal{B}(\hat{x},\;s^{*}(\xi))$ as the convergence region.
First, for sufficiently large $k$, when every state $x_{i}[k]$ is in the
convergence region (region 1 or 2), by using Proposition 1 and the gradient
update (5), we show that the regular states cannot leave the convergence
region (Steps 1 to 3 below). Next, we show that if some state $x_{i}[k]$ is
outside of the convergence region (i.e., in region 3), then by applying
Proposition 1 and the gradient update (5) repeatedly, every state will
eventually enter the convergence region (Step 4). These steps will lead to the
proof of Theorem 2. Throughout the following exposition, we fix $\xi$ to be a
strictly positive constant.

Step 1: We want to show that if
$z_{i}[k]\in\mathcal{B}(\hat{x},\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\})\subset\mathcal{B}(\hat{x},s^{*}(\xi))$
and $k$ is large enough, then by applying the gradient update (5), the state
$x_{i}[k+1]\in\mathcal{B}(\hat{x},s^{*}(\xi))$ i.e., $x_{i}[k+1]$ is still in
the convergence region. Let $k^{*}_{1}\in\mathbb{N}$ be such that
$\eta[k^{*}_{1}]\leq\frac{\xi}{L}$.

###### Lemma 3

For all $v_{i}\in\mathcal{R}$ and $k\geq k^{*}_{1}$, if
$\|z_{i}[k]-\hat{x}\|_{2}\leq\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\}$
then $\|x_{i}[k+1]-\hat{x}\|_{2}\leq s^{*}(\xi)$.

###### Proof:

From the update (5), we have

$\displaystyle\|x_{i}[k+1]-\hat{x}\|_{2}$
$\displaystyle=\|z_{i}[k]-\hat{x}-\eta[k]g_{i}[k]\|_{2}$
$\displaystyle\leq\|z_{i}[k]-\hat{x}\|_{2}+\eta[k]\|g_{i}[k]\|_{2}.$

Using Assumption 2 and 3, and
$\|z_{i}[k]-\hat{x}\|_{2}\leq\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\}$,
we obtain

$\displaystyle\|x_{i}[k+1]-\hat{x}\|_{2}$
$\displaystyle\leq\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\}+\eta[k]L$
$\displaystyle\leq\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\}+\eta[k_{1}^{*}]L.$

By the definition of $k_{1}^{*}$ and $s^{*}(\xi)$, the above inequality
becomes

$\displaystyle\|x_{i}[k+1]-\hat{x}\|_{2}\leq\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\}+\xi\leq
s^{*}(\xi).$

∎

Step 2: We want to analyze the relationship between the terms
$\|x_{i}[k+1]-\hat{x}\|_{2}$ and $\|z_{i}[k]-\hat{x}\|_{2}$ which will be used
in the subsequent step.

For $v_{i}\in\mathcal{R}$, define
$\Delta_{i}:[\tilde{R}_{i},\infty)\times\mathbb{R}_{+}\to\mathbb{R}$ to be the
function

$\displaystyle\Delta_{i}(p,l)\triangleq
2l(\sqrt{p^{2}-\tilde{R}_{i}^{2}}\cos\theta_{i}-\tilde{R}_{i}\sin\theta_{i})-l^{2},$
(17)

and $\Gamma_{i}:[\tilde{R}_{i},\infty)\times\mathbb{R}_{+}\to\mathbb{R}$ to be
the function

$\displaystyle\Gamma_{i}(p,l)\triangleq p^{2}-\Delta_{i}(p,l).$ (18)

###### Lemma 4

For all $v_{i}\in\mathcal{R}$ and $k\in\mathbb{N}$, if
$\|z_{i}[k]-\hat{x}\|_{2}>\tilde{R}_{i}+\delta_{i}$ then

$\displaystyle\|x_{i}[k+1]-\hat{x}\|_{2}^{2}$
$\displaystyle\leq\|z_{i}[k]-\hat{x}\|_{2}^{2}-\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\;\eta[k]\;\|g_{i}[k]\|_{2})$
$\displaystyle=\Gamma_{i}(\|z_{i}[k]-\hat{x}\|_{2},\;\eta[k]\;\|g_{i}[k]\|_{2}).$

###### Proof:

Consider $v_{i}\in\mathcal{R}$. Define $\|z_{i}[k]-\hat{x}\|_{2}\triangleq
s_{i}[k]$ and $z_{i}[k]\notin\mathcal{B}(\hat{x},\tilde{R}_{i}+\delta_{i})$.
By Lemma 2, we have

$\displaystyle\qquad\angle(x_{i}^{*}-z_{i}[k],\hat{x}-z_{i}[k])$
$\displaystyle\leq\max_{y\in\mathcal{B}(\hat{x},\tilde{R}_{i})}\angle(y-z_{i}[k],\hat{x}-z_{i}[k])$
$\displaystyle=\arcsin\frac{\tilde{R}_{i}}{s_{i}[k]}\triangleq\phi_{i}(s_{i}[k]).$
(19)

Let
$\psi_{i}[k]\triangleq\angle(x_{i}[k+1]-z_{i}[k],\;\hat{x}-z_{i}[k])=\angle(-\eta[k]g_{i}[k],\;\hat{x}-z_{i}[k])$.
Since for all $v_{i}\in\mathcal{R}$ and $k$,

$\angle(-\eta[k]g_{i}[k],\;\hat{x}-z_{i}[k])\leq\angle(-\eta[k]g_{i}[k],\;x_{i}^{*}-z_{i}[k])\\\
+\angle(x_{i}^{*}-z_{i}[k],\;\hat{x}-z_{i}[k])$

by [18, Corollary 12], applying Proposition 2 and inequality (19), we have

$\displaystyle\psi_{i}[k]\leq\theta_{i}+\phi_{i}(s_{i}[k])\triangleq\tilde{\psi}_{i}(s_{i}[k]).$
(20)

Note that $\tilde{\psi}_{i}(s_{i}[k])<\pi$ since $\theta_{i}<\frac{\pi}{2}$
and $\phi_{i}(s_{i}[k])\leq\frac{\pi}{2}$. Consider the triangle which has the
vertices at $x_{i}[k+1]$, $z_{i}[k]$, and $\hat{x}$. We can calculate the
square of the distance by using the law of cosines:

$\|x_{i}[k+1]-\hat{x}\|_{2}^{2}=\;\|x_{i}[k+1]-z_{i}[k]\|^{2}_{2}+\|\hat{x}-z_{i}[k]\|^{2}_{2}\\\
-2\|x_{i}[k+1]-z_{i}[k]\|_{2}\cdot\|z_{i}[k]-\hat{x}\|_{2}\\\
\times\cos\angle(x_{i}[k+1]-z_{i}[k],\;\hat{x}-z_{i}[k]).$

Using the update (5) and $\|z_{i}[k]-\hat{x}\|_{2}=s_{i}[k]$, we get

$\displaystyle\|x_{i}[k+1]-\hat{x}\|_{2}^{2}$
$\displaystyle=\eta[k]^{2}\;\|g_{i}[k]\|_{2}^{2}+s_{i}[k]^{2}$
$\displaystyle\quad-2s_{i}[k]\;\eta[k]\;\|g_{i}[k]\|_{2}\cos\psi_{i}[k]$
$\displaystyle\leq\eta[k]^{2}\;\|g_{i}[k]\|_{2}^{2}+s_{i}[k]^{2}$
$\displaystyle\quad-2s_{i}[k]\;\eta[k]\;\|g_{i}[k]\|_{2}\cos\tilde{\psi}_{i}(s_{i}[k]).$
(21)

The inequality above is obtained by using the inequality (20). In addition, we
can simplify the term $s_{i}[k]\cos\tilde{\psi}_{i}(s_{i}[k])$ as follows:

$\displaystyle s_{i}[k]\cos\tilde{\psi}_{i}(s_{i}[k])$
$\displaystyle=s_{i}[k]\cos\Big{(}\arcsin\frac{\tilde{R}_{i}}{s_{i}[k]}+\theta_{i}\Big{)}$
$\displaystyle=s_{i}[k]\cos\Big{(}\arcsin\frac{\tilde{R}_{i}}{s_{i}[k]}\Big{)}\cos\theta_{i}$
$\displaystyle\quad-
s_{i}[k]\sin\Big{(}\arcsin\frac{\tilde{R}_{i}}{s_{i}[k]}\Big{)}\sin\theta_{i}$
$\displaystyle=\sqrt{s_{i}[k]^{2}-\tilde{R}_{i}^{2}}\cos\theta_{i}-\tilde{R}_{i}\sin\theta_{i}.$

Therefore, we can rewrite (21) as

$\|x_{i}[k+1]-\hat{x}\|_{2}^{2}\leq\eta[k]^{2}\;\|g_{i}[k]\|_{2}^{2}+s_{i}[k]^{2}\\\
-2\;\eta[k]\;\|g_{i}[k]\|_{2}\;(\sqrt{s_{i}[k]^{2}-\tilde{R}_{i}^{2}}\cos\theta_{i}-\tilde{R}_{i}\sin\theta_{i}).$

∎

Step 3: We want to show that if $z_{i}[k]\in\mathcal{B}(\hat{x},s^{*}(\xi))$,
$z_{i}[k]\notin\mathcal{B}(\hat{x},\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\})$
and $k$ is large enough, then by applying the gradient update (5), the state
$x_{i}[k+1]\in\mathcal{B}(\hat{x},s^{*}(\xi))$ i.e., $x_{i}[k+1]$ is still in
the convergence region.

Let

$\displaystyle a_{i}^{+}$
$\displaystyle\triangleq-\tilde{R}_{i}\sin\theta_{i}+\sqrt{\big{(}s^{*}(\xi)\big{)}^{2}-\tilde{R}_{i}^{2}\cos^{2}\theta_{i}}\;,$
$\displaystyle a_{i}^{-}$
$\displaystyle\triangleq-\tilde{R}_{i}\sin\theta_{i}-\sqrt{\big{(}s^{*}(\xi)\big{)}^{2}-\tilde{R}_{i}^{2}\cos^{2}\theta_{i}}\;,\quad\text{and}$
$\displaystyle b_{i}$ $\displaystyle\triangleq
2\big{(}\sqrt{\big{(}s^{*}(\xi)\big{)}^{2}-\tilde{R}_{i}^{2}}\cos\theta_{i}-\tilde{R}_{i}\sin\theta_{i}\big{)}.$

Let $k_{2}^{*}$ be such that

$\displaystyle\eta[k_{2}^{*}]\leq\frac{1}{L}\min_{v_{i}\in\mathcal{R}}\big{\\{}\min\\{a_{i}^{+},\;b_{i}\\}\big{\\}}.$

###### Lemma 5

For all $i\in\mathcal{R}$ and $k\geq k_{2}^{*}$, if
$\|z_{i}[k]-\hat{x}\|_{2}\in\big{(}\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\},s^{*}(\xi)\big{]}$
then $\|x_{i}[k+1]-\hat{x}\|_{2}\leq s^{*}(\xi)$.

###### Proof:

Consider the function $\Gamma_{i}(p,l)$ defined in (18). Compute the second
derivative with respect to $p$:

$\displaystyle\frac{\partial\Gamma_{i}}{\partial
p}=2p-2lp(p^{2}-\tilde{R}_{i}^{2})^{-\frac{1}{2}}\cos\theta_{i},$
$\displaystyle\frac{\partial^{2}\Gamma_{i}}{\partial
p^{2}}=2+2l\tilde{R}_{i}^{2}(p^{2}-\tilde{R}_{i}^{2})^{-\frac{3}{2}}\cos\theta_{i}.$

Note that $\frac{\partial^{2}\Gamma_{i}}{\partial p^{2}}>0$ for all
$p\in(\tilde{R}_{i},\infty)$. This implies that

$\sup_{p\in(\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\},\;s^{*}(\xi)]}\Gamma_{i}(p,l)\leq\max_{p\in[\tilde{R}_{i},\;s^{*}(\xi)]}\Gamma_{i}(p,l)\\\
=\max\big{\\{}\Gamma_{i}(\tilde{R}_{i},l),\;\Gamma_{i}(s^{*}(\xi),l)\big{\\}}.$
(22)

First, let consider $\Gamma_{i}(\tilde{R}_{i},\;l)$. From the definition of
$\Gamma_{i}$ in (18), we have that

$\displaystyle\Gamma_{i}(\tilde{R}_{i},l)\leq\big{(}s^{*}(\xi)\big{)}^{2}\quad\Longleftrightarrow\quad
l\in[a_{i}^{-},\;a_{i}^{+}].$ (23)

Note that $a_{i}^{+}>0$ and $a_{i}^{-}<0$ since
$s^{*}(\xi)\geq\tilde{R}_{i}+\xi$. Using Assumption 2, Assumption 3, and the
definition of $k_{2}^{*}$, we have that for all $v_{i}\in\mathcal{R}$ and
$k\geq k_{2}^{*}$,

$\displaystyle\eta[k]\;\|g_{i}[k]\|_{2}\leq\eta[k]L\leq\min_{v_{j}\in\mathcal{R}}\big{\\{}\min\\{a_{j}^{+},\;b_{j}\\}\big{\\}}\leq
a_{i}^{+}.$

By (23), we obtain that for all $v_{i}\in\mathcal{R}$ and $k\geq k_{2}^{*}$,

$\displaystyle\Gamma_{i}(\tilde{R}_{i},\;\eta[k]\;\|g_{i}[k]\|_{2})\leq\big{(}s^{*}(\xi)\big{)}^{2}.$
(24)

Now, let consider $\Gamma_{i}(s^{*}(\xi),\;l)$. From the definition of
$\Gamma_{i}$ in (18), we have that

$\displaystyle\Gamma_{i}(s^{*}(\xi),l)\leq\big{(}s^{*}(\xi)\big{)}^{2}\quad\Longleftrightarrow\quad
l\in[0,\;b_{i}].$ (25)

Note that $b_{i}>0$ since $s^{*}(\xi)\geq\tilde{R}_{i}\sec\theta_{i}+\xi$.
Using Assumption 2, Assumption 3, and the definition of $k_{2}^{*}$, we have
that for all $v_{i}\in\mathcal{R}$ and $k\geq k_{2}^{*}$,

$\displaystyle\eta[k]\;\|g_{i}[k]\|_{2}\leq\eta[k]L\leq\min_{v_{j}\in\mathcal{R}}\big{\\{}\min\\{a_{j}^{+},\;b_{j}\\}\big{\\}}\leq
b_{i}.$

By (25), we obtain that for all $v_{i}\in\mathcal{R}$ and $k\geq k_{2}^{*}$,

$\displaystyle\Gamma_{i}(s^{*}(\xi),\;\eta[k]\;\|g_{i}[k]\|_{2})\leq\big{(}s^{*}(\xi)\big{)}^{2}.$
(26)

Combine (24) and (26) to get that for all $v_{i}\in\mathcal{R}$ and $k\geq
k_{2}^{*}$,

$\max\big{\\{}\Gamma_{i}(\tilde{R}_{i},\;\eta[k]\;\|g_{i}[k]\|_{2}),\;\Gamma_{i}(s^{*}(\xi),\;\eta[k]\;\|g_{i}[k]\|_{2})\big{\\}}\\\
\leq\big{(}s^{*}(\xi)\big{)}^{2}.$ (27)

Using Lemma 4, (22) and (27), respectively, we obtain that for all
$v_{i}\in\mathcal{R}$ and $k\geq k_{2}^{*}$,

$\displaystyle\quad\|x_{i}[k+1]-\hat{x}\|_{2}^{2}$
$\displaystyle\leq\Gamma_{i}(\|z_{i}[k]-\hat{x}\|_{2},\;\eta[k]\;\|g_{i}[k]\|_{2})$
$\displaystyle\leq\sup_{p\in(\max_{v_{j}\in\mathcal{R}}\\{\tilde{R}_{j}+\delta_{j}\\},\;s^{*}(\xi)]}\Gamma_{i}(p,\;\eta[k]\;\|g_{i}[k]\|_{2})$
$\displaystyle\leq\max\big{\\{}\Gamma_{i}(\tilde{R}_{i},\;\eta[k]\;\|g_{i}[k]\|_{2}),\;\Gamma_{i}(s^{*}(\xi),\;\eta[k]\;\|g_{i}[k]\|_{2})\big{\\}}$
$\displaystyle\leq\big{(}s^{*}(\xi)\big{)}^{2}.$

∎

Step 4: We now show that the states of all regular nodes eventually enter the
convergence region.

###### Lemma 6

There exists $k\geq\max\\{k_{1}^{*},k_{2}^{*}\\}$ such that for all
$v_{i}\in\mathcal{R}$, $\|x_{i}[k]-\hat{x}\|_{2}\leq s^{*}(\xi)$.

###### Proof:

Assume for the purpose of contradiction that for all
$k\geq\max\\{k_{1}^{*},k_{2}^{*}\\}$ there exists $v_{i}\in\mathcal{R}$ such
that $\|x_{i}[k]-\hat{x}\|_{2}>s^{*}(\xi)$. Define the sets of nodes

$\displaystyle\mathcal{I}_{x}[k]\triangleq\\{v_{i}\in\mathcal{R}:\|x_{i}[k]-\hat{x}\|_{2}>s^{*}(\xi)\\}$

and

$\displaystyle\mathcal{I}_{z}[k]\triangleq\\{v_{i}\in\mathcal{R}:\|z_{i}[k]-\hat{x}\|_{2}>s^{*}(\xi)\\}.$

Note that $\mathcal{I}_{x}[k]\neq\emptyset$ for all
$k\geq\max\\{k_{1}^{*},k_{2}^{*}\\}$ by the assumption. In addition, from
Lemma 3 and Lemma 5, we have that for all $v_{i}\in\mathcal{R}$ and
$k\geq\max\\{k_{1}^{*},k_{2}^{*}\\}$, if $\|z_{i}[k]-\hat{x}\|_{2}\leq
s^{*}(\xi)$ then $\|x_{i}[k+1]-\hat{x}\|_{2}\leq s^{*}(\xi)$. This implies
that for all $k\geq\max\\{k_{1}^{*},k_{2}^{*}\\}$,
$\emptyset\neq\mathcal{I}_{x}[k+1]\subseteq\mathcal{I}_{z}[k]$.

From Lemma 4, for $k\geq\max\\{k_{1}^{*},k_{2}^{*}\\}$ and
$v_{i}\in\mathcal{I}_{z}[k]$,

$\|x_{i}[k+1]-\hat{x}\|_{2}^{2}\leq\|z_{i}[k]-\hat{x}\|_{2}^{2}\\\
-\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\;\eta[k]\;\|g_{i}[k]\|_{2}).$

For any $k\geq\max\\{k_{1}^{*},k_{2}^{*}\\}$, we have

$\displaystyle\quad\max_{v_{i}\in\mathcal{I}_{z}[k]}\|x_{i}[k+1]-\hat{x}\|_{2}^{2}$
$\displaystyle\leq\max_{v_{i}\in\mathcal{I}_{z}[k]}\Big{(}\|z_{i}[k]-\hat{x}\|_{2}^{2}$
$\displaystyle\qquad\qquad-\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\eta[k]\;\|g_{i}[k]\|_{2})\Big{)}$
$\displaystyle\leq\max_{v_{i}\in\mathcal{I}_{z}[k]}\|z_{i}[k]-\hat{x}\|_{2}^{2}$
$\displaystyle\qquad\qquad-\min_{v_{i}\in\mathcal{I}_{z}[k]}\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\;\eta[k]\;\|g_{i}[k]\|_{2}).$
(28)

Consider the first term on the RHS of (28). By the definition of
$\mathcal{I}_{z}[k]$, we have

$\displaystyle\max_{v_{i}\in\mathcal{I}_{z}[k]}\|z_{i}[k]-\hat{x}\|_{2}^{2}=\max_{v_{i}\in\mathcal{R}}\|z_{i}[k]-\hat{x}\|_{2}^{2}.$

Consider the term on the LHS of (28). Since
$\emptyset\neq\mathcal{I}_{x}[k+1]\subseteq\mathcal{I}_{z}[k]$ and the
definition of $\mathcal{I}_{x}[k+1]$, we have

$\displaystyle\max_{v_{i}\in\mathcal{I}_{z}[k]}\|x_{i}[k+1]-\hat{x}\|_{2}^{2}$
$\displaystyle\geq\max_{v_{i}\in\mathcal{I}_{x}[k+1]}\|x_{i}[k+1]-\hat{x}\|_{2}^{2}$
$\displaystyle=\max_{v_{i}\in\mathcal{R}}\|x_{i}[k+1]-\hat{x}\|_{2}^{2}.$

Therefore, the inequality (28) becomes

$\max_{v_{i}\in\mathcal{R}}\|x_{i}[k+1]-\hat{x}\|_{2}^{2}\leq\max_{v_{i}\in\mathcal{R}}\|z_{i}[k]-\hat{x}\|_{2}^{2}\\\
-\min_{v_{i}\in\mathcal{I}_{z}[k]}\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\;\eta[k]\;\|g_{i}[k]\|_{2}).$
(29)

Consider the definition of the function $\Delta_{i}$ in (17). It is clear that
if $p_{1}>p_{2}\geq\tilde{R}_{i}$ then
$\Delta_{i}(p_{1},l)>\Delta_{i}(p_{2},l)$. Then, for all
$k\geq\max\\{k_{1}^{*},k_{2}^{*}\\}$ and $v_{i}\in\mathcal{I}_{z}[k]$, we get

$\displaystyle\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\eta[k]\;\|g_{i}[k]\|_{2})>\Delta_{i}(s^{*}(\xi),\eta[k]\;\|g_{i}[k]\|_{2}).$
(30)

Furthermore, the function $\Delta_{i}$ satisfies

$\Delta_{i}(p,l)\geq\Big{(}\sqrt{p^{2}-\tilde{R}_{i}^{2}}\cos\theta_{i}-\tilde{R}_{i}\sin\theta_{i}\Big{)}l\\\
\Longleftrightarrow\quad
l\in\Big{[}0,\;\sqrt{p^{2}-\tilde{R}_{i}^{2}}\cos\theta_{i}-\tilde{R}_{i}\sin\theta_{i}\Big{]}.$
(31)

For $x\notin\mathcal{C}_{i}(\epsilon)$, from the definition of convex
functions, we have $-\langle g_{i}(x),\;x_{i}^{*}-x\rangle\geq
f_{i}(x)-f_{i}(x_{i}^{*})$. Using the inequality (11), we obtain

$\displaystyle\|g_{i}(x)\|_{2}\geq\frac{f_{i}(x)-f_{i}(x_{i}^{*})}{\|x-x_{i}^{*}\|_{2}}\geq\frac{\epsilon}{\delta_{i}(\epsilon)}\triangleq\kappa_{i}(\epsilon)>0,$
(32)

where $\mathcal{C}_{i}(\epsilon)$ is defined in (7). Let $k_{3}^{*}$ be such
that

$\displaystyle\eta[k_{3}^{*}]\leq\frac{1}{2L}\min_{v_{i}\in\mathcal{R}}b_{i}.$

Using (32), Assumption 2, and the definition of $k_{3}^{*}$, we have that for
all $k\geq k_{3}^{*}$ and $v_{i}\in\mathcal{I}_{z}[k]$,

$\displaystyle\eta[k]\kappa_{i}\leq\eta[k]\;\|g_{i}[k]\|_{2}\leq\eta[k]L\leq\frac{b_{i}}{2}.$
(33)

Since $\eta[k]\;\|g_{i}[k]\|_{2}\in[0,\;\frac{b_{i}}{2}]$, we can apply (31)
to get

$\displaystyle\Delta_{i}\big{(}s^{*}(\xi),\;\eta[k]\;\|g_{i}[k]\|_{2}\big{)}\geq\frac{b_{i}}{2}\eta[k]\;\|g_{i}[k]\|_{2}.$
(34)

Let $k^{*}=\max\\{k_{1}^{*},k_{2}^{*},k_{3}^{*}\\}$. Combine (30), (34), and
the first inequality of (33) to get that for all $k\geq k^{*}$ and
$v_{i}\in\mathcal{I}_{z}[k]$,

$\displaystyle\Delta_{i}(\|z_{i}[k]-\hat{x}\|_{2},\eta[k]\;\|g_{i}[k]\|_{2})>\frac{1}{2}b_{i}\kappa_{i}\eta[k].$
(35)

On the other hand, from Proposition 1, we have that for all
$v_{i}\in\mathcal{R}$,

$\displaystyle\|z_{i}[k]-\hat{x}\|_{2}$
$\displaystyle\leq\max_{v_{j}\in(\mathcal{N}_{i}^{-}\cap\mathcal{R})\cup\\{v_{i}\\}}\|x_{j}[k]-\hat{x}\|_{2}$
$\displaystyle\leq\max_{v_{j}\in\mathcal{R}}\|x_{j}[k]-\hat{x}\|_{2},$

which implies that

$\displaystyle\max_{v_{i}\in\mathcal{R}}\|z_{i}[k]-\hat{x}\|_{2}\leq\max_{v_{i}\in\mathcal{R}}\|x_{i}[k]-\hat{x}\|_{2}.$

Apply the above inequality and the inequality (35) to (29) to get that for all
$k\geq k^{*}$,

$\displaystyle\quad\max_{v_{i}\in\mathcal{R}}\|x_{i}[k+1]-\hat{x}\|_{2}^{2}$
$\displaystyle<\max_{v_{i}\in\mathcal{R}}\|x_{i}[k]-\hat{x}\|_{2}^{2}-\frac{1}{2}\min_{v_{i}\in\mathcal{I}_{z}[k]}b_{i}\kappa_{i}\eta[k]$
$\displaystyle\leq\max_{v_{i}\in\mathcal{R}}\|x_{i}[k]-\hat{x}\|_{2}^{2}-\frac{1}{2}\min_{v_{i}\in\mathcal{R}}b_{i}\kappa_{i}\eta[k].$
(36)

Apply the inequality (36) recursively to get that for $K\geq 1$,

$\max_{v_{i}\in\mathcal{R}}\|x_{i}[k^{*}+K]-\hat{x}\|_{2}^{2}<\max_{v_{i}\in\mathcal{R}}\|x_{i}[k^{*}]-\hat{x}\|_{2}^{2}\\\
-\frac{1}{2}\min_{v_{i}\in\mathcal{R}}b_{i}\kappa_{i}\sum_{k=k^{*}}^{k^{*}+K-1}\eta[k].$

Since $\sum_{k=1}^{\infty}\eta[k]=\infty$ from Assumption 3, we have that
there exists $K^{*}$ such that

$\displaystyle\max_{v_{i}\in\mathcal{R}}\|x_{i}[k^{*}+K^{*}]-\hat{x}\|_{2}^{2}\leq\big{(}s^{*}(\xi)\big{)}^{2}.$

This implies that for all $v_{i}\in\mathcal{R}$,

$\displaystyle\|x_{i}[k^{*}+K^{*}]-\hat{x}\|_{2}^{2}\leq\big{(}s^{*}(\xi)\big{)}^{2}$

which contradicts our assumption. ∎

We are now in place to prove Theorem 2 by combining the results from Lemma 3,
5, and 6. Intuitively, Lemma 6 states that there is a time-step $k$ that
$x_{i}[k]\in\mathcal{B}(\hat{x},s^{*}(\xi))$ for all $v_{i}\in\mathcal{R}$.
Then, by applying Lemma 3 and 5, and Proposition 1, we can conclude that the
state $x_{i}[k]$ will remain inside $\mathcal{B}(\hat{x},s^{*}(\xi))$ for all
subsequent time-steps.

###### Proof:

Let $k^{*}\geq\max\big{\\{}k_{1}^{*},k_{2}^{*}\big{\\}}$ be the time index
that satisfies Lemma 6, i.e.,

$\displaystyle\max_{v_{i}\in\mathcal{R}}\|x_{i}[k^{*}]-\hat{x}\|_{2}\leq
s^{*}(\xi).$

Applying Proposition 1, we get

$\displaystyle\max_{v_{i}\in\mathcal{R}}\|z_{i}[k^{*}]-\hat{x}\|_{2}\leq\max_{v_{i}\in\mathcal{R}}\|x_{i}[k^{*}]-\hat{x}\|_{2}\leq
s^{*}(\xi).$ (37)

On the other hand, from Lemma 3 and 5, we have that for
$k\geq\max\\{k_{1}^{*},k_{2}^{*}\\}$, if $\|z_{i}[k]-\hat{x}\|_{2}\leq
s^{*}(\xi)$ then $\|x_{i}[k+1]-\hat{x}\|_{2}\leq s^{*}(\xi)$. Therefore, the
inequality (37) implies that

$\displaystyle\max_{v_{i}\in\mathcal{R}}\|x_{i}[k^{*}+1]-\hat{x}\|_{2}\leq
s^{*}(\xi).$

Applying the above procedure recursively, we get

$\displaystyle\|x_{i}[k]-\hat{x}\|_{2}\leq s^{*}(\xi)\quad\text{for
all}\;\;k\geq k^{*}\quad\text{and}\;\;v_{i}\in\mathcal{R}.$

However, since we can choose $\xi$ to be arbitrary small, we have that for all
$v_{i}\in\mathcal{R}$,

$\displaystyle\limsup_{k}\|x_{i}[k]-\hat{x}\|_{2}$ $\displaystyle\leq
s^{*}(0)$
$\displaystyle=\max_{v_{j}\in\mathcal{R}}\big{\\{}\max\\{\tilde{R}_{j}\sec\theta_{j},\tilde{R}_{j}+\delta_{j}\\}\big{\\}}.$

Furthermore, the variables $\theta_{i}(\epsilon)$ and $\delta_{i}(\epsilon)$
depend on $\epsilon$ which is defined in (7) and the above inequality holds
for any $\epsilon>0$. Hence, for all $v_{i}\in\mathcal{R}$,

$\limsup_{k}\|x_{i}[k]-\hat{x}\|_{2}\\\
\leq\inf_{\epsilon>0}\Big{\\{}\max_{v_{j}\in\mathcal{R}}\big{\\{}\max\\{\tilde{R}_{j}\sec\theta_{j}(\epsilon),\;\tilde{R}_{j}+\delta_{j}(\epsilon)\\}\big{\\}}\Big{\\}}.$

∎

## References

  * [1] A. Nedić and A. Olshevsky, “Distributed optimization over time-varying directed graphs,” _IEEE Transactions on Automatic Control_ , vol. 60, no. 3, pp. 601–615, 2014.
  * [2] R. Xin, C. Xi, and U. A. Khan, “Frost—fast row-stochastic optimization with uncoordinated step-sizes,” _EURASIP Journal on Advances in Signal Processing_ , vol. 2019, no. 1, p. 1, 2019.
  * [3] 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, 2013.
  * [4] M. Zhu and S. Martínez, “On distributed convex optimization under inequality and equality constraints,” _IEEE Transactions on Automatic Control_ , vol. 57, no. 1, pp. 151–164, 2011.
  * [5] S. Sundaram and B. Gharesifard, “Distributed optimization under adversarial nodes,” _IEEE Transactions on Automatic Control_ , vol. 64, no. 3, pp. 1063–1076, 2018.
  * [6] N. Ravi, A. Scaglione, and A. Nedić, “A case of distributed optimization in adversarial environment,” in _IEEE International Conference on Acoustics, Speech and Signal Processing_ , 2019, pp. 5252–5256.
  * [7] S. X. Wu, H.-T. Wai, A. Scaglione, A. Nedić, and A. Leshem, “Data injection attack on decentralized optimization,” in _IEEE Int. Conf. on Acoustics, Speech and Signal Processing_ , 2018, pp. 3644–3648.
  * [8] L. Su and N. Vaidya, “Byzantine multi-agent optimization: Part i,” _arXiv preprint arXiv:1506.04681_ , 2015.
  * [9] L. Su and N. H. Vaidya, “Fault-tolerant multi-agent optimization: optimal iterative distributed algorithms,” in _ACM Symposium on Principles of Distributed Computing_ , 2016, pp. 425–434.
  * [10] C. Zhao, J. He, and Q.-G. Wang, “Resilient distributed optimization algorithm against adversary attacks,” in _IEEE International Conference on Control & Automation (ICCA)_, 2017, pp. 473–478.
  * [11] K. Kuwaranancharoen and S. Sundaram, “On the location of the minimizer of the sum of two strongly convex functions,” in _IEEE Conference on Decision and Control (CDC)_ , 2018, pp. 1769–1774.
  * [12] Z. Yang and W. U. Bajwa, “Byrdie: Byzantine-resilient distributed coordinate descent for decentralized learning,” _IEEE Transactions on Signal and Information Processing over Networks_ , 2019.
  * [13] ——, “Bridge: Byzantine-resilient decentralized gradient descent,” _arXiv preprint arXiv:1908.08098_ , 2019.
  * [14] N. Gupta and N. H. Vaidya, “Byzantine fault tolerant distributed linear regression,” _arXiv preprint arXiv:1903.08752_ , 2019.
  * [15] P. Blanchard, R. Guerraoui, J. Stainer _et al._ , “Machine learning with adversaries: Byzantine tolerant gradient descent,” in _Advances in Neural Information Processing Systems_ , 2017, pp. 119–129.
  * [16] H. J. LeBlanc, H. Zhang, X. Koutsoukos, and S. Sundaram, “Resilient asymptotic consensus in robust networks,” _IEEE Journal on Selected Areas in Communications_ , vol. 31, no. 4, pp. 766–781, 2013.
  * [17] B. S. Mordukhovich and N. M. Nam, “An easy path to convex analysis and applications,” _Synthesis Lectures on Mathematics and Statistics_ , vol. 6, no. 2, pp. 1–218, 2013.
  * [18] D. Castano, V. E. Paksoy, and F. Zhang, “Angles, triangle inequalities, correlation matrices and metric-preserving and subadditive functions,” _Linear Algebra and its Applications_ , vol. 491, pp. 15–29, 2016.


# Analyzing stability of a delay differential equation involving two delays

Sachin Bhalekar
Department of Mathematics,
Shivaji University, Kolhapur - 416004, India
Email<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Analysis of the systems involving delay is a popular topic among applied
scientists. In the present work, we analyze the generalized equation
$D^{\alpha}x(t)=g\left(x(t-\tau_{1}),x(t-\tau_{2})\right)$ involving two
delays viz. $\tau_{1}\geq 0$ and $\tau_{2}\geq 0$. We use the the stability
conditions to propose the critical values of delays. Using examples, we show
that the chaotic oscillations are observed in the unstable region only. We
also propose a numerical scheme to solve such equations.

##  1 Introduction

Equations involving nonlocal operators are crucial in modeling natural
systems. The memory involved in such systems cannot be modeled properly by
using local operators such as ordinary integer order derivative. The
fractional order derivative and the time delay (lag) are proved suitable for
this task.

If the order of the derivative involved in the model is non-integer then it
can be called as a fractional derivative. Analysis and applications of
fractional derivatives can be found in [2, 1, 3, 4, 5, 6].

If the modeling differential equation includes the past values of state
variable then it is called as a delay differential equation (DDE). Basic
analysis and various applications of DDE are discussed in [7, 8, 9, 10].

Stability of fractional order delay differential equations (FDDEs) is
discussed by Matignon in [11]. BIBO stability of linear time invariant FDDEs
is discussed in [12, 13]. A numerical algorithm is proposed by Hwang and Cheng
[14] to test the stability of FDDEs. Stability of a class of FDDE is studied
in [15]. In [16], the present author proposed stability and bifurcation
analysis of generalized FDDE. Chaos in various FDDEs is analyzed by Bhalekar
and Daftardar-Gejji [17, 18, 19]. Application of FDDE in NMR is discussed by
Bhalekar et al. [20, 21].

Equations involving single delay are extensively studied in the literature.
However, the analysis of the systems with multiple delay is complicated and
hence the work is relatively rare. Geometry of stable regions in differential
equations with two delays is presented in [22] and references cited therein.
Hopf bifurcation in these systems is studied in [23, 24, 25] and references
cited therein. Applications of systems involving multiple delay can be found
in various natural systems. Bi and Ruan [26] discussed these systems in tumor
and immune system interaction models. Prey-predator system with multiple delay
is analyzed by Gakkhar and Singh [27].

In this work, we extend our analysis in [16] to consider two delay
differential equation with fractional order. We present stability and
bifurcation analysis of these equations, propose numerical method and solve
examples.

The paper is organized as follows: Basic definitions are listed in Section 2.
The expressions for critical curves are derived in Section 3. Section 4 deals
with the numerical method for solving general two delay models. Illustrated
examples are presented in Section 5. Conclusions are summarized in Section 6.

##  2 Preliminaries

In this section, we discuss some basic definitions [1, 2, YLG].

###### Definition 2.1

A real function $f(t),\,t>0$ is said to be in space
$C_{\alpha},\,\alpha\in\Re$ if there exists a real number $p$ ($>\alpha$),
such that $f(t)=t^{p}f_{1}(t)$ where $f_{1}(t)\in C[0,\infty)$.

###### Definition 2.2

A real function $f(t),$ $t>0$ is said to be in space $C_{\alpha}^{m},\,m\in
I\\!\\!N\bigcup\left\\{0\right\\}$ if $f^{(m)}\in C_{\alpha}.$

###### Definition 2.3

Let $f\in C_{\alpha}$ and $\alpha\geq-1,$ then the (left-sided)
Riemann–Liouville integral of order $\mu,\,\mu>0$ is given by

$I^{\mu}f(t)=\frac{1}{\Gamma(\mu)}\int_{0}^{t}(t-\tau)^{\mu-1}f(\tau)\,d\tau,\quad
t>0.$ (1)

###### Definition 2.4

The Caputo fractional derivative of $f$, $f\in C_{-1}^{m},m\in I\\!\\!N,$ is
defined as:

$\displaystyle D^{\mu}f(t)$
$\displaystyle=\frac{d^{m}}{dt^{m}}f(t),\quad\quad\quad\mu=m$ (2)
$\displaystyle=I^{m-\mu}\frac{d^{m}f(t)}{dt^{m}},\quad m-1<\mu<m.$

Note that for $m-1<\mu\leq m,\,\,m\in I\\!\\!N,$

$\displaystyle I^{\mu}D^{\mu}f(t)$
$\displaystyle=f(t)-\sum^{m-1}_{k=0}\frac{d^{k}f}{dt^{k}}(0)\frac{t^{k}}{k!},$
(3) $\displaystyle I^{\mu}t^{\nu}$
$\displaystyle=\frac{\Gamma(\nu+1)}{\Gamma(\mu+\nu+1)}t^{\mu+\nu}.$ (4)

##  3 Main Results

We consider the following generalized delay differential equation (GDDE)

$D^{\alpha}x(t)=g\left(x(t-\tau_{1}),x(t-\tau_{2})\right),$ (5)

where $D^{\alpha}$ is a Caputo fractional derivative of order
$\alpha\in(0,1]$, $g$ is continuously differentiable function and
$\tau_{1}\geq 0$, $\tau_{2}\geq 0$ are delays.

A number $x^{*}\in\mathbb{R}$ is called an equilibrium point of (5) if

$g\left(x^{*},x^{*}\right)=0.$ (6)

Linearization of the system (5) in the neighborhood of equilibrium point
$x^{*}$ is given by

$\displaystyle D^{\alpha}\xi=a\xi_{\tau_{1}}+b\xi_{\tau_{2}},$ (7)

where $x_{\tau}(t)=x(t-\tau)$ and $a$, $b$ are partial derivatives of $g$ with
respect to first and second variables evaluated at $\left(x^{*},x^{*}\right)$
respectively. The characteristic equation of the system is

$\lambda^{\alpha}=a\exp(-\lambda\tau_{1})+b\exp(-\lambda\tau_{2}).$ (8)

###  3.1 Stability of equilibrium

We present the stability analysis by using the same strategy as in [16]. If
all the eigenvalues $\lambda_{i}$ of characteristic equation (8) satisfy

$Re(\lambda_{i})<0,\,\forall i.$ (9)

then the equilibrium $x^{*}$ is asymptotically stable.

Hence, the equation without delay (i.e. with $\tau_{1}=\tau_{2}=0$) is
asymptotically stable if

$a+b<0.$ (10)

Now, we assmue that $\tau_{1}>0$, $\tau_{2}>0$ and $\lambda=u+\imath v$,
$u,v\in\Re$. The stability properties of equilibrium will change at
$\lambda=\imath v$. In this case, the characteristic equation takes the form

$\left(\imath v\right)^{\alpha}=a\exp(-\imath v\tau_{1})+b\exp(-\imath
v\tau_{2}).$ (11)

Simplifying, we get

$\displaystyle
v^{\alpha}\cos\left(\frac{\alpha\pi}{2}\right)-a\cos(v\tau_{1})$
$\displaystyle=$ $\displaystyle b\cos(v\tau_{2})$ (12) $\displaystyle
v^{\alpha}\sin\left(\frac{\alpha\pi}{2}\right)+a\sin(v\tau_{1})$
$\displaystyle=$ $\displaystyle-b\sin(v\tau_{2}).$ (13)

Squaring and adding

$v^{2\alpha}+a^{2}-2av^{\alpha}\cos\left(\frac{\alpha\pi}{2}+v\tau_{1}\right)=b^{2}.$
(14)

This gives

$\tau_{1}=\frac{1}{v}\left(-\frac{\alpha\pi}{2}+\arccos\left[\frac{v^{2\alpha}+a^{2}-b^{2}}{2av^{\alpha}}\right]\right).$
(15)

Using similar arguments and rewriting system (12)–(13), we get

$\tau_{2}=\frac{1}{v}\left(-\frac{\alpha\pi}{2}+\arccos\left[\frac{v^{2\alpha}-a^{2}+b^{2}}{2bv^{\alpha}}\right]\right).$
(16)

Equations (15) and (16) are parametric expressions of $\tau_{1}$ and
$\tau_{2}$ in parameter $v$. These can be used to plot critical curves.

##  4 Numerical Method

In this section, we present a numerical scheme based on [28] to solve GDDE (5)
with initial function $x(t)=\phi(t),$ $t\leq 0$. This initial value problem
(IVP) is equivalent to the integral equation

$\displaystyle x(t)$ $\displaystyle=$
$\displaystyle\int_{0}^{t}\frac{\left(t-\zeta\right)^{\alpha-1}}{\Gamma(\alpha)}f\left(x(\zeta-\tau_{1}),x(\zeta-\tau_{2})\right)d\zeta$
(17) $\displaystyle+\phi(0),\quad t\in[0,T].$

We consider the uniform grid
$\left\\{t_{n}=nh:n=-k,-k+1,\cdots,-1,0,1,\cdots,N\right\\}$ where $k$ and $N$
are the integers chosen so that $h=T/N=\tau_{1}/k_{1}=\tau_{2}/k_{2}$ and
$k=max\left\\{k_{1},k_{2}\right\\}$. Let $x_{j}$ be the approximation to
$x(t_{j})$. Note that $x_{j}=\phi(t_{j}),$ if $j\leq 0$ and
$x(t_{j}-\tau_{1})=x_{j-k_{1}}$, $x(t_{j}-\tau_{2})=x_{j-k_{2}}$.

Evaluating (17) at $t=t_{n+1}$ and using product trapezoidal formula, we get

$x_{n+1}=\phi(0)+\frac{h^{\alpha}}{\Gamma(\alpha+2)}\sum_{j=0}^{n+1}a_{j,n+1}f\left(x_{j-k_{1}},x_{j-k_{2}}\right),$
(18)

where

$a_{j,n+1}=\left\\{\begin{array}[]{ccc}n^{\alpha+1}-(n-\alpha)(n+1)^{\alpha}&\mbox{
if}&j=0;\\\ (n-j+2)^{\alpha+1}+(n-j)^{\alpha+1}&\\\
-2(n-j+1)^{\alpha+1}&\mbox{if}&1\leq j\leq n;\\\ 1&\mbox{if}&j=n+1.\\\
\end{array}\right.$

##  5 Illustrative Examples

In this section, we analyze stability and solve some examples using the
numerical scheme presented in preceding section.

###### Example 5.1

Consider the fractional order generalization of Uçar system [29]

$D^{\alpha}x(t)=\delta x(t-\tau_{1})-\epsilon\left[x(t-\tau_{2})\right]^{3},$
(19)

where $\delta$ and $\epsilon$ are positive parameters.

There are three equilibrium points of the system (19) viz.
$0,\pm\sqrt{\delta/\epsilon}$. If $\tau_{1}=\tau_{2}=0$ then $0$ is unstable
whereas $\pm\sqrt{\delta/\epsilon}$ are stable. We take $\delta=1$.

The characteristic equation of the system (19) at
$x^{*}=\pm\sqrt{\delta/\epsilon}$ is

$\lambda^{\alpha}-\exp(-\lambda\tau_{1})+3\exp(-\lambda\tau_{2})=0.$ (20)

Note that the characteristic equation and hence stability does not depend on
$\epsilon$. We take $\epsilon=1$.

Using (15) and (16) we obtain the critical values of delay as

$\tau_{1}=\frac{1}{v}\left(-\frac{\alpha\pi}{2}+\arccos\left[\frac{v^{2\alpha}-8}{2v^{\alpha}}\right]\right)$
(21)

and

$\tau_{2}=\frac{1}{v}\left(-\frac{\alpha\pi}{2}+\arccos\left[\frac{v^{2\alpha}+8}{-6v^{\alpha}}\right]\right).$
(22)

The following stability result based on [24] for integer order case $\alpha=1$
is proposed in [29].

###### Theorem 5.1

[29] The system $\dot{x}(t)=\delta
x(t-\tau_{1})-\epsilon\left[x(t-\tau_{2})\right]^{3}$ is locally stable if one
of the following conditions hold:
(i) $\tau_{2}\in[0,0.366667/\delta)$
(ii) $\tau_{2}\in[0.366667/\delta,0.43521/\delta)$ and
$\tau_{1}\in[0,0.0697127/\delta)$.

It can be verified that the stability region described in this Theorem can be
readily obtained by using expressions (21) and (22). Further, as shown in Fig.
1 we can extend this region for higher values of delay also.

Fig. 2 shows the critical curves for different values of $\alpha$.

Figure 1: Extended stability region of Uçar system for $\alpha=1$ Figure 2:
Critical curves of Uçar system for different values of $\alpha$

Now, we verify these stability results using numerical computations. In Fig.
3, we have taken $\tau_{1}=0.4$, $\tau_{2}=0.5$ and $\alpha=0.9$. These values
are in stable region and the figure shows the orbit converging to an
equilibrium state. In Fig. 4, we take the parameter values $\tau_{1}=1.6$,
$\tau_{2}=1.4$ and $\alpha=0.9$ in unstable region. In this case, the system
exhibit chaotic oscillations.

Figure 3: Stable orbit of Uçar system for $\alpha=0.9$, $\tau_{1}=0.4$ and
$\tau_{2}=0.5$ Figure 4: Chaotic attractor of Uçar system for $\alpha=0.9$,
$\tau_{1}=1.6$ and $\tau_{2}=1.4$

###### Example 5.2

Consider the generalization of fractional order Ikeda equation [30] with two
delays

$D^{\alpha}x(t)=-3x\left(t-\tau_{1}\right)+24\sin\left(x\left(t-\tau_{2}\right)\right),\quad
0<\alpha\leq 1.$ (23)

As discussed in [16], there are seven equilibrium points. The points $x_{0}=0$
and $x_{3,4}=\pm 7.49775$ are unstable whereas $x_{1,2}=\pm 2.7859$ and
$x_{5,6}=\pm 7.95732$ are stable at $\tau_{1}=\tau_{2}=0$. The characteristic
equation of the system (23) at the equilibrium point $x^{*}$ is

$\lambda^{\alpha}=-3\exp(-\lambda\tau_{1})+24\cos(x^{*})\exp(-\lambda\tau_{2}).$
(24)

If we take $x^{*}=x_{1,2}$ then the characteristic equation becomes

$\lambda^{\alpha}=-3\exp(-\lambda\tau_{1})-22.4977\exp(-\lambda\tau_{2}).$
(25)

This gives $a=-3$ and $b=-22.4977$. The critical curve in this case is given
in Fig. 5 for $\alpha=0.7$.

Figure 5: Critical curve of Ikeda system for $\alpha=0.7$

In Fig. 6, we have presented stable solution of system (23) with
$\tau_{1}=0.02$ and $\tau_{2}=0.01$ in stable region. If we consider the
values $\tau_{1}=0.01$ and $\tau_{2}=0.1$ in unstable region then we get
chaotic attractor (cf. Fig. 7).

Figure 6: Stable orbit of Ikeda system for $\alpha=0.7$, $\tau_{1}=0.02$ and
$\tau_{2}=0.01$ Figure 7: Chaotic attractor shown by Ikeda system for
$\alpha=0.7$, $\tau_{1}=0.01$ and $\tau_{2}=0.1$

##  6 Conclusion

In [16], we have presented the complete analysis of the generalized equation
involving delay. However, the analysis of the systems with multiple delay is
not as simple as that of single delay. In this article, we have presented the
expressions for critical values of generalized two delay system. We obtain the
parametric relation between the delays $\tau_{1}$ and $\tau_{2}$. The region
bounded by horizontal axis and the critical curve is stable region if the
eigenvalues are stable at $\tau_{1}=\tau_{2}=0$. The nonlinear system may
exhibit chaotic oscillations in unstable region. We also have presented
numerical scheme to solve these systems. Further, the results are illustrated
with two examples.

## Acknowledgements

Author acknowledges the Science and Engineering Research Board (SERB), New
Delhi, India for the Research Grant (Ref. MTR/2017/000068) under Mathematical
Research Impact Centric Support (MATRICS) Scheme.

## References

  * [1] I. Podlubny, Fractional Differential Equations (Academic Press, New York, 1999).
  * [2] S. G. Samko, A. A. Kilbas and O. I. Marichev, Fractional Integrals and Derivatives: Theory and Applications (Gordon and Breach, Yverdon, 1993).
  * [3] A. A. Kilbas, H. M. Srivastava, J. J. Trujillo, Theory and applications of fractional differential equations (Elsevier, Amsterdam, 2006).
  * [4] F. Mainardi, Fractional calculus and waves in linear viscoelasticity (Imperial College Press, London, 2010).
  * [5] R. L. Magin, Fractional calculus in bioengineering (Begll House Publishers, USA, 2006).
  * [6] X. Yang, D. Baleanu and H. M. Srivastava, Local Fractional Integral Transforms and Their Applications (Academic Press, London, 2015).
  * [7] H. Smith, An introduction to delay differential equations with applications to the life sciences (Springer, New York, 2010).
  * [8] M. Lakshmanan, D. V. Senthilkumar, Dynamics of nonlinear time-delay systems (Springer, Heidelberg, 2010).
  * [9] C. Lainscsek, P. Rowat, L. Schettino, D. Lee, D. Song, C. Letellier, and H. Poizner, “Finger tapping movements of Parkinson’s disease patients automatically rated using nonlinear delay differential equations,” Chaos 22, 013119 (2012).
  * [10] C. Lainscsek, and T. J. Sejnowski, “Electrocardiogram classification using delay differential equations,” Chaos 23, 023132 (2013)
  * [11] D. Matignon, “Stability results for fractional differential equations with applications to control processing.” In: Computational Engineering in Systems and Application multiconference, vol. 2, pp. 963–968, IMACS, IEEE-SMC Proceedings, Lille, France, July (1996).
  * [12] M. A. Pakzad and S. Pakzad, “Stability map of fractional order time-delay systems,” WSEAS Trans. Syst. 11(10), 541–550 (2012).
  * [13] C. Bonnet and J. R. Partington, “Analysis of fractional delay systems of retarded and neutral type,” Automatica 38(7), 1133–1138 (2002).
  * [14] C. Hwang, Y. C. Cheng, “A numerical algorithm for stability testing of fractional delay systems,” Automatica 42, 825–831 (2006).
  * [15] S. Bhalekar, “Stability analysis of a class of fractional delay differential equations,” Pramana- Journal of Physics 81(2), 215–224 (2013).
  * [16] S. Bhalekar, “Stability and bifurcation analysis of a generalized scalar delay differential equation,” Chaos 26(8), 084306 (2017).
  * [17] S. Bhalekar and V. Daftardar-Gejji, “Fractional ordered Liu system with time-delay,” Commun. Nonlinear Sci. Numer. Simulat. 15(8), 2178–2191 (2010).
  * [18] V. Daftardar-Gejji, S. Bhalekar and P. Gade, “Dynamics of fractional ordered Chen system with delay,” Pramana- Journal of Physics 79(1), 61–69 (2012).
  * [19] S. Bhalekar, “Dynamical analysis of fractional order Ucar prototype delayed system,” Signals Image and Video Processing 6(3), 513–519 (2012).
  * [20] S. Bhalekar, V. Daftardar-Gejji, D. Baleanu and R. Magin, “Fractional Bloch equation with delay,” Comput. Math. Appl. 61(5), 1355–1365 (2011).
  * [21] S. Bhalekar, V. Daftardar-Gejji, D. Baleanu and R. Magin, “Generalized fractional order Bloch equation with extended delay,” Int. J. Bifurc. Chaos 22(4), 1250071 (2012).
  * [22] J. K. Hale and W. Huang, “Global geometry of the stable regions for two delay differential equations,” J. Math. Anal. Appl. 178, 344–362 (1993).
  * [23] J. Belair and S. A. Campbell, “Stability and bifurcations of equilibria in a multiple-delayed differential equation,” SIAM J. Appl. Math. 54, 1402–1424 (1994).
  * [24] X. Li and S. Ruan, “Stability and bifurcation in delay differential equations with two delays,” J. Math. Anal. Appl. 236, 254–280 (1999).
  * [25] X. P. Wu and L. Wang, “Zero-Hopf bifurcation analysis in delayed differential equations with two delays,” J. Franklin Inst. 354, 1484–1513 (2017).
  * [26] P. Bi and S. Ruan, “Bifurcations in delay differential equations and applications to tumor and immune system interaction models,” SIAM J. Appl. Dyn. Sys. 12(4), 1847–1888 (2013).
  * [27] S. Gakkhar and A. Singh, “Complex dynamics in a prey predator system with multiple delays,” Commun. Nonlinear Sci. Numer. Simulat. 17, 914–929 (2012).
  * [28] V. Daftardar-Gejji, Y. Sukale and S. Bhalekar, “Solving fractional delay differential equations: A new approach,” Fract. Calc. Appl. Anal. 18, 400–418 (2015).
  * [29] S. Bhalekar, “On the Uçar prototype model with incommensurate delays,” Signal, Image and Video Processing 8(4), 635–639 (2014).
  * [30] J. G. Lu, “Chaotic dynamics of the fractional-order Ikeda delay system and its synchronization,” Chinese Physics 15(2) 301–305 (2006).