Local Synchronization of Sampled-Data Systems on Lie Groups
Philip James McCarthy Christopher Nielsen


Analysis of Low Excitation HDO Transitions Toward the High-Mass Star-forming Regions
                            G34.26+0.15, W51e_1/e_2, and W49N
    Magda  Kulczak-Jastrzȩbska
    Received: date / Accepted: date
==================================================================================================================================================

empty
empty










We present a smooth distributed nonlinear control law for local synchronization of identical driftless kinematic agents on a Cartesian product of matrix Lie groups with a connected communication graph. If the agents are initialized sufficiently close to one another, then synchronization is achieved exponentially fast. We first analyze the special case of commutative Lie groups and show that in exponential coordinates, the closed-loop dynamics are linear. We characterize all equilibria of the network and, in the case of an unweighted, complete graph, characterize the settling time and conditions for deadbeat performance. Using the Baker-Campbell-Hausdorff theorem, we show that, in a neighbourhood of the identity element, all results generalize to arbitrary matrix Lie groups.






§ INTRODUCTION



The sampled-data setup is ubiquitous in applied control. In the LTI case, the plant may be exactly discretized, and a discrete-time controller can be designed such that closed-loop stability is achieved for non-pathological sampling-periods. Such stability guarantees cannot generally be enforced for nonlinear plants, as nonlinear ODEs generally do not have closed-form solutions. The standard approach to nonlinear sampled-data control design, emulation, is to approximate the discretized plant dynamics. If the sampling-period is sufficiently small, then the actual closed-loop system is stable. This technique has two key shortcomings <cit.>: 1) it may not be possible for a given approximate discretization method, e.g., Euler's method; 2) it relies on fast sampling, which may not be possible due to hardware limitations.

The limitations of emulation do not necessarily pose a problem for the class of systems on matrix Lie groups, which are nonlinear, yet have dynamics that yield exact closed-form solutions <cit.>, thereby enabling direct design. To our knowledge, sampled-data control of systems on Lie groups has not yet been explored in the literature. However, the closely related class of bilinear systems has been studied in the discrete-time <cit.> and sampled-data settings <cit.>.

Many engineering systems are modelled on Lie groups. The motion of robots in a plane is modelled on 2 <cit.>, and their motion in space, such as that of UAVs', is modelled on 3 <cit.>. Quantum systems evolve on the unitary groups 𝖴(n) <cit.> and 𝖲𝖴(n) <cit.>. Some circuits can be modelled using Lie groups <cit.> while oscillator networks <cit.> evolve on 2.

Synchronization of networks on 3 was achieved using passivity in <cit.>. Synchronization under sampling was studied for a network of Kuramoto-like oscillators in <cit.> and harmonic oscillators with a time-varying period in <cit.>, and path following nonlinear agents in <cit.>, but the analyses in these works were not conducted from a Lie-theoretic perspective. The Kuramoto network model was extended from 2 to n in <cit.>. A framework for coordinated motion on Lie groups was developed in <cit.>, where the synchronization problem that we consider is a special case of what the authors call bi-invariant coordination. In <cit.>, linear consensus algorithms were applied to systems on Lie groups in continuous-time. The most salient difference between the current paper and <cit.> is the consideration of the sampled-data setup. Further, in contrast to <cit.>, we take a global perspective, and explore the geometry of the problem to much greater depth.

Lie groups are not vector spaces, but their structure facilitates analysis and control design in global coordinates. The Lie structure has been leveraged, for example, for motion tracking on 3 <cit.>, and the control of UAV <cit.> and spacecraft <cit.> orientation on 3. We too take a global perspective in our control design and, when possible, analysis. 

We present a control law that achieves synchronization for a network of identical agents on any matrix Lie group with driftless dynamics with a connected communication graph. The current paper generalizes and extends the results of <cit.>, which considered only unweighted graphs with agents on one-parameter Lie subgroups. The controller requires that each agent have access to its relative state with respect to each of its neighbours. For example, on 3, relative position and orientation can be attained using machine vision <cit.>. We examine the special case where the error dynamics evolve on a Cartesian product of one-parameter subgroups – a “generalized cylinder” – and the general case. We prove that in both cases, that if the agents are initialized sufficiently close to one another, that synchronization is achieved exponentially fast. For a generalized cylinder, we characterize the performance in the case of an unweighted, complete graph.




 §.§ Notation and Terminology


If N ∈, then _N 1, …, N. Given a matrix M ∈^n × n, M^⊤ is its (non-Hermitian) transpose, and λ_max(M) and λ_min(M) are its its eigenvalues of greatest and least magnitude, respectively. If x ∈^n, then x is its Euclidean norm; if M ∈^n × n, then M is its induced Euclidean norm. Let 1_n ∈^n and 0_n ∈^n denote the column vector of ones and zeros, respectively. Let 0_m × n∈^m × n denote the matrix of zeros. Let ^- denote the set of nonpositive real numbers. Given an equivalence relation ∼ on a set R, and an element x ∈ R, let [x] ∈ R/∼ be the coset containing x.

Weighted, directed graphs are used to model communication constraints between agents. A graph 𝒢 is a triple (𝒱, ℰ,w) consisting of a finite set of vertices 𝒱 = _N, a set of edges ℰ⊆𝒱×𝒱, and a weight function w : ℰ→ [0, 1] ⊂. The weight w_ij w((i,j)) is nonzero only if (i, j) ∈ℰ. If agent i has access to its relative state with respect to agent j, then (i,j) ∈ℰ. Define vertex i's neighbour set as 𝒩_i {j ∈_N: (i,j) ∈ℰ}. We assume that 𝒢 has no self-loops. If 𝒢 is unweighted, then for all i ≠ j ∈_N, w_ij∈{0,1}. Associated with 𝒢 is the Laplacian L ∈^N × N, defined elementwise as


    L_ij = 
    		{[        -w_ij,        i ≠ j,; ∑_j ∈𝒩_iw_ij,        i = j. ].


The ith row of the Laplacian L is denoted ℓ_i.




§ SAMPLED-DATA SYNCHRONIZATION PROBLEM



We consider a network of N controlled agents, each modelled by the differential equation


    Ẋ_i = X_i(∑_j=1^mB_i,ju_i,j),      i ∈_N.


Here X_i ∈ where ⊂n, is an m-dimensional connected matrix Lie group over the complex field  which includes, as a special case, real matrix Lie groups. The matrices B_i,j are elements of the Lie algebra , which is a vector space over a field  equal to either  or , associated with , and u_i  (u_i,1, …, u_i,m) ∈𝔽^m is the control input. Note that the Lie algebra of a complex Lie group  may in fact be a real vector space. For example, the Lie algebra 𝔰𝔲(2) of the complex Lie group 2 is a vector space over the field of reals despite its vectors being matrices with possibly complex entries. Equation (<ref>) is a kinematic model of a system evolving on a matrix Lie group . Each agent is assumed to be fully actuated in the sense that


    (∀ i ∈_N)  _𝔽B_i,1, …, B_i,m = .


Under this assumption, without loss of generality, we take the
system (<ref>) to be driftless since the inputs u_i,
i ∈_N, can be chosen to cancel any drift vector field. We are
interested in the sampled-data control of this multi-agent system in
which each agent's control law is implemented on an embedded computer,
which we explicitly model using the setup in Figure <ref>. The
blocks H and S in Figure <ref> are, respectively, the ideal
hold and sample operators. Sample and hold are, respectively,
idealized models of A/D and D/A conversion.
 


The following assumption is made throughout this paper.


All sample and hold blocks operate at the same period T > 0 and the blocks are synchronized for the multi-agent system (<ref>). ◂


Under Assumption <ref>, letting X_i[k]  X_i(kT) and u_i[k]  u_i(kT), the discretized dynamics of each agent are given by


    X_i^+ = X_iexp(T∑_j=1^mB_i,ju_i,j),      i ∈_N


which is an exact discretization of (<ref>). For each i ∈_N, define Ω_ i ∑_j=1^mB_i,ju_i,j∈. Then the discrete-time dynamics can be compactly expressed as


    X_i^+ = X_iexp(TΩ_i),      i ∈_N.







 §.§ The Synchronization Problem



Given a network of N agents with kinematic dynamics (<ref>), we define the error quantities E_ij X_i^-1X_j, i,j ∈_N. Observe that E_ij = I if, and only if, X_i = X_j. The error matrix E_ij is called left-invariant <cit.>, since for all X ∈, (XX_i)^-1(XX_j) = X_i^-1X_j. The class of systems considered has  for its state space, which is generally not a vector space, so we do not use X_i - X_j as a measure of error.

Local Synchronization on Matrix Lie Groups
Given a network of N agents with continuous-time dynamics (<ref>), sampling period T > 0 and an unweighted, connected communication graph 𝒢 = (𝒱, ℰ), find, if possible, distributed control laws Ω_i, i ∈_N, such that for all initial errors in a neighbourhood of the identity in ^N, for all i,j ∈_N, E_ij→ I as t →∞.


By a distributed control law we mean that for each agent i, the control signal Ω_i can depend on E_ij only if (i, j) ∈ℰ. In this paper we propose the distributed feedback control law


    Ω_i 1/Tlog((∏_j ∈𝒩_iE_ij^_ij)^1/K)


where K ∈ is a gain and the matrix logarithm need not be the
principal logarithm. The control law (<ref>) does not
require agent i to know agent j's state X_j, nor its own state
X_i, but instead requires knowledge of the relative state
E_ij. The expression (<ref>) is well-defined so long
as the product ∏_j ∈𝒩_iE_ij has no eigenvalues
in ^-, as discussed in Section <ref>. This
control law is expressed in global coordinates but we only prove local
exponential stability of the synchronized state. When the control
law (<ref>) is well-defined, the closed-loop
discrete-time dynamics are


    X_i^+ = X_i(∏_j ∈𝒩_iE_ij^_ij)^1/K,      i ∈_N


and the synchronization error dynamics are



    E_ij^+    = (X_i^+)^-1X_j^+ 
       = (∏_p ∈𝒩_iE_ip^_ip)^-1/KX_i^-1X_j(∏_q ∈𝒩_jE_jq^_jq)^1/K
       = (∏_p ∈𝒩_iE_ip^_ip)^-1/KE_ij(∏_q ∈𝒩_jE_jq^_jq)^1/K.




The order of multiplication in (<ref>) need not be common to all agents or even constant.
▴








The control law (<ref>) is motivated by exponential coordinates for Lie groups, classical consensus algorithms in ^n, and the notion of Riemannian mean of rotations on 3, which on a one-parameter subgroup thereof can be explicitly computed as ∏_i=1^NR_i^1/N <cit.>.
 
A key advantage of direct design over emulation, is that stability can be guaranteed at the sampling instants. As mentioned in the Introduction, on 3, the relative error E_ij can be computed using machine vision, where the speed of sampling is limited by the frame rate of the camera, for example, 25 Hz <cit.>. This limits the feasibility of emulation-based design. However, direct design does not guarantee good performance between sampling instants. But in the specific case of the plant and problem discussed in this paper, achieving synchronization at the sampling instants implies synchronization between the sampling instants.


If E_ij[k] = X_i[k]^-1X_j[k] asymptotically approaches I as k →∞, then E_ij(t) = X_i(t)^-1 X_j(t) asymptotically approaches I as t →∞, where X_i(t) and X_j(t) evolve according to (<ref>).



If E_ij[k] → I, then the proposed control law (<ref>) satisfies Ω_i[k] →0_n × n. Let 0 < δ < T. Then


    lim_k →∞E_ij(kT + δ) = lim_k →∞exp(δΩ_i[k])^-1E_ij[k]exp(δΩ_j[k])
           = lim_k →∞exp(δΩ_i[k])^-1lim_k →∞E_ij[k]lim_k →∞exp(δΩ_j[k])
           = I^3 = I.


Since δ is arbitrary, this implies that E_ij(t) → I.


Proposition <ref> means that asymptotically stabilizing the set where E_ij = I, for all i,j ∈_N, at the sample instances is sufficient for solving the synchronization problem. Thus, we can conduct all analysis in the discrete-time setting and do not rely on T being sufficiently small.

Our main result is the following theorem, which we prove in Section <ref>.


For any Lie group  with connected communication graph 𝒢, if the gain K of each agent's controller (<ref>) satisfies (<ref>), then the equilibrium {E_ij = I : i,j ∈_N} is locally uniformly exponentially stable.





§ PRELIMINARIES







 §.§ Functions of matrices


For every nonsingular matrix X ∈^n × n there are (infinitely many) A ∈^n × n such that exp(A) = X, see <cit.>. Every such matrix A is a non-primary logarithm of X, which we denote by log(X). If, in addition to being nonsingular, the matrix X has no eigenvalues in ^-, then it has a (unique) principal logarithm.


Let X ∈^n × n have no eigenvalues in ^-. There is a unique logarithm A ∈^n × n of X, all of whose eigenvalues lie in the strip {z : -π < Im(z) < π}. If X ∈^n × n, then A ∈^n × n.


The unique matrix A from Theorem <ref> is called the principal logarithm of X and is denoted (X). Unlike complex numbers, it is not possible to express log(X) as a function of (X) for arbitrary non-singular matrices. If X - I < 1, then

    (X) = ∑_k=1^∞(-1)^k-1/k (X - I)^k.


Any matrix logarithm is a right inverse of the matrix exponential, but not necessarily a left inverse. On a matrix Lie group, the principal logarithm  is a left inverse, but only in a neighbourhood of the identity. Choose r > 0 such that (<ref>) converges on X ∈ : X = exp(A), A < r, e.g., r = (2) is a valid choice with any Lie group. Larger values of r may be possible for specific Lie groups. The set


    U X ∈: X = exp(A), A ∈, A < r


is an open neighbourhood of I in  in the group topology in which : U → provides an inverse.

Borrowing from the definition of complex powers of scalars <cit.> and the form of the square root of a matrix on a Lie group <cit.>, we define the Kth root of a matrix in the following way.


Let X ∈^n × n have no eigenvalues in ^-. Given K ∈, the principal Kth root of X is


    X^1/Kexp(1/K(X)).



If X ∈, then X^1/K∈, due to the Lie correspondence : →.


If X^1/K is well-defined, then for K ∈,


    (X^1/K)^K = exp(∑_i = 1^K1/K(X)) = exp((X)) = A.


Thus, in this case, X = X^1/KX^1/K⋯ X^1/K (K times), which is the intuitive notion of a Kth root. The somewhat indirect definition (<ref>) allows for Kth roots for K ∈.
▴


Throughout this paper, we use an important algebraic property of the logarithm of a matrix power.


  If X ∈^n × n has no eigenvalues in ^-,
  then for α∈ [-1,1], we have
  (X^α) = α(X).





 §.§ Exponential Coordinates and One-Parameter Subgroups





Given a Lie group , a one-parameter subgroup is a continuous morphism of groups ϕ : →.


Although this terminology is standard, it is technically the image of the map ϕ that is a subgroup of . The subgroup ϕ() ⊂ is a one-dimensional manifold and there exists a unique H ∈ such that ϕ(θ) = exp(θ H) for all θ∈ <cit.>.

To generalize the concept of one-parameter subgroups to higher dimensional manifolds, we consider generalized cylinders. A generalized cylinder is an m-dimensional manifold that is diffeomorphic to 𝕋^k ×^m - k. Such a diffeomorphism exists if and only if there exist m commutative and everywhere-linearly-independent vector fields on the manifold <cit.>. If the manifold is a Lie group , then this simplifies to its Lie algebra  having a commutative basis.

Let  be such a manifold and fix a commutative basis ℋ{H_1,…,H_m} for its Lie algebra . Consider the one-parameter groups ϕ_i : → associated with each H_i. The image of ϕ(t_1,…,t_m) ϕ_1(t_1)ϕ_2(t_2)⋯ϕ_m(t_m) is . Without loss of generality, let ϕ_i, i ∈_k, 0 ≤ k ≤ m have nonzero kernel, and let ϕ_i, i ∈{k + 1,…,m} have zero kernel.

Fixing such a basis ℋ, the  map induces local coordinates on ∩ U. Given X ∈, by commutativity of H_1,…,H_m,


    X    = exp(t_1 H_1)⋯exp(t_m H_m) 
       = exp(t_1 H_1 + ⋯ + t_m H_m).


If X ∈∩ U, then (X) = t_1 H_1 + ⋯ + t_m
H_m.
Then, by linear independence of H_1,…,H_m, t_1,…,t_m can
be uniquely determined, yielding local coordinates
(t_1,…,t_m) ∈^m. Thus, a Lie group  can be locally
identified with an open subset of the vector space ^m
containing the origin. Note that, by commutativity of , these
local coordinates coincide with the familiar exponential coordinates
of both the first and second kind.



 §.§ Properties of the composed flow

The map ϕ : ^m →, defined in the previous section, is critical to our analysis throughout this paper. In this section, we establish important properties of ϕ when  is a generalized cylinder, we then show that these properties hold approximately for any Lie group in a neighbourhood of the identity.

By definition, ϕ is surjective onto its image, but it is not necessarily injective. Let p : ^m →^m/(ϕ) be the projection of ^m onto the quotient space ^m/(ϕ). There exists a unique isomorphism of groups ϕ' such that the following diagram commutes <cit.>.


    ^m [drr]_ϕ[rr]^p 	      ^m/(ϕ) [d]^ϕ'


The bijection ϕ' yields alternative global coordinates on the quotient group ^m/(ϕ); it will be used in our characterization of equilibria. When  is a generalized cylinder, the map ϕ has several important properties.


If  is a generalized cylinder, then ϕ : ^m → is a morphism of groups.



It is clear that ϕ(0_m) = exp(0)⋯exp(0) = I.

Let _i = (_i^(1),…,_i^(m)) ∈^m and _j = (_j^(1),…,_j^(m)) ∈^m, where ϕ(_i) = X_i and ϕ(_j) = X_j. By commutativity of H_1,…,H_m,


    ϕ(_i + _j)    = exp((t_i^(1) + t_j^(1))H_1)⋯exp((t_i^(m) + t_j^(m))H_m) 
       = exp(t_i^(1)H_1)⋯exp(t_i^(m)H_m)exp(t_j^(1)H_1)⋯exp(t_j^(m)H_m) 
       = ϕ(_i)ϕ(_j).




Let  be a generalized cylinder. If K > 0 and ϕ() ∈ U, then ϕ(/K) = ϕ()^1/K.



By the commutativity of H_1,…,H_m,


    ϕ(/K)    = exp(t_1/KH_1)⋯exp(t_m/KH_m) 
       = exp(1/K(t_1H_1 + ⋯ + t_mH_m)) 
       = exp(1/K(ϕ())) 
       = ϕ()^1/K.






§ EQUILIBRIA ON GENERALIZED CYLINDERS


Since we consider driftless kinematic models, the system is at equilibrium if, and only if, every agent's input is zero, i.e., for all i ∈_N, Ω_i = 0_n × n. We show that all equilibria are isolated and exhibit the same stability properties.

Hereinafter, we use the notation _ijϕ^-1(E_ij) and [ _11^⊤ _12^⊤     ⋯ _1N^⊤ ]^⊤∈^Nm. Define _𝕋 and _ as the projections under p of  onto 𝕋^Nk and ^N(m - k), respectively.


If the controller (<ref>) is well-defined, then the equilibria of (<ref>) on a generalized cylinder are characterized by [1/K(L ⊗ I_k)_𝕋] = [0_Nk], (L ⊗ I_m-k)_ = 0_N(m-k), [_11] = [0_m].



The sampled dynamics of each agent are given by (<ref>). Therefore, the system is at equilibrium if, and only if, for all i ∈_N, exp(TΩ_i) = I. If  is a generalized cylinder, then, by commutativity and Definition <ref>, this condition becomes


    I = exp(TΩ_i) = (∏_j ∈𝒩_iE_ij^_ij)^1/K = ∏_j ∈𝒩_iE_1i^-_ij/KE_1j^_ij/K.


In the global coordinates admitted by ϕ' we have


    ϕ'^-1(I)    = ∑_j ∈𝒩_i(ϕ'^-1(E_1i^-_ij/K) + ϕ'^-1(E_1j^_ij/K)) 
    
    	[0_m]    = ∑_j ∈𝒩_i([_ij/K_1j] - [_ij/K_1i])
    	= -[1/K(ℓ_i ⊗ I_m)].




We “stack" the inputs for all agents i, yielding the equation


    [1/K(L ⊗ I_m)] = [0_Nm],


which can be rewritten as


    [1/K[         L ⊗ I_k 0_Nk × N(m - k); 0_N(m - k) × Nk     L ⊗ I_m - k ][ _𝕋;  _ ]] = 
    						[[       0_Nk; 0_N(m - k) ]].


By assumption, (ϕ_i) = {0} for i ∈{k+1,…,m}, so the condition [1/K(L ⊗ I_m - k)] = [0_N(m - k)] simplifies to equality on ^N(m - k), rather than congruence on a quotient space. Lastly, since _11 is the error of agent 1 with itself, [_11] = [0_m].



On a generalized cylinder, all equilibria are isolated.



By assumption, (ϕ_i) = {0}, for all i ∈{k+1,…,m}. Thus /(ϕ_i) ≅. Thus, for j ∈_N, [_1j^(i)] = [0] simplifies to _1j^(i) = 0.

By assumption, (ϕ_i) ≠{0}, for all i ∈_k. The
map ϕ_i can be viewed as a flow, thus, by <cit.>, there exists a d_i > 0 such that for every
r_i ∈ [0], we have r_i = q_id_i for some q_i ∈. Thus,
for all j ∈_N, if _1j^(i), _1j^(i)∈ [0],
_1j^(i)≠_1j^(i), then
|_1j^(i) - _1j^(i)| ≥ d_i.



On a generalized cylinder, every equilibrium has the same stability properties as the identity.



Let {Ξ_i1,…,Ξ_iN}∈^N be an equilibrium. Define E̅_ijΞ_ij^-1E_ij. Then E̅_ij = I if and only if E̅_ij = Ξ_ij. The error dynamics (<ref>) can be expressed in terms of E̅_ij:


    Ξ_ijE̅_ij^+    = (∏_p ∈𝒩_i(Ξ_1iE̅_1i)^-_ip(Ξ_1pE̅_1p)^_ip)^-1/KΞ_ijE̅_ij×(∏_q ∈𝒩_j(Ξ_1jE̅_1j)^-_jq(Ξ_1qE̅_1q)^_jq)^1/K
    E̅_ij^+    = (∏_p ∈𝒩_iE̅_1i^-_ipE̅_1p^_ip)^-1/KE̅_ij(∏_q ∈𝒩_jE̅_1j^-_jqE̅_1q^_jq)^1/K(∏_p ∈𝒩_iΞ_1i^-_ipΞ_1p^_ip)^-1/K(∏_q ∈𝒩_jΞ_1j^-_jqΞ_1q^_jq)^1/K
       = (∏_p ∈𝒩_iE̅_ip^_ip)^-1/KE̅_ij(∏_q ∈𝒩_jE̅_jq^_jq)^1/K(∏_p ∈𝒩_iΞ_ip^_ip)^-1/K(∏_q ∈𝒩_jΞ_jq^_jq)^1/K.


By (<ref>), the Ξ_ij product terms equal
identity, thus the E̅_ij dynamics have the same form as the
dynamics (<ref>) of E_ij, and therefore have the same
qualitative behaviour.


Proposition <ref> says that the dynamics near every equilibrium “look the same”. Thus, by analyzing only the equilibrium at identity, we characterize the behaviour near all equilibria.




§ SYNCHRONIZATION ON GENERALIZED CYLINDERS



In this section, we consider the case where  is a generalized
cylinder. This means that X_i ∈, i ∈_N, which implies
E_ij∈, for all i,j ∈_N. This is done only to
simplify discussion. The results of this section hold under the weaker
assumption that the errors E_ij lie on a generalized cylinder.


Any generalized cylinder , on which (<ref>) is well-defined for all forward time, is positively invariant for (<ref>).



Let k ∈ be arbitrary and suppose that for all i, j ∈_N, E_ij[k] ∈. Let q ∈_N be arbitrary. By elementary group theory ∏_j ∈𝒩_q E_qj[k] ∈ and therefore, by hypothesis, its Kth root is well-defined. Thus, by definition of log, exp(TΩ_q[k])∈. Since E_ij^+ = exp(-TΩ_i)E_ijexp(TΩ_j),  we have E_ij[k+1] ∈. Induction on the time index proves positive invariance of .
















Using the exponential coordinates from Section <ref> we identify each relative error E_ij∈∩ U with its exponential coordinates _ij∈^m. We henceforth impose that the synchronization errors and their products over neighbour sets be close to the identity.


For all i,j ∈_N, we have E_ij∈ U and ∏_j ∈𝒩_iE_ij^w_ij∈ U. ◂


Let  be a generalized cylinder and suppose E_ij[0] ∈, i,j ∈_N.  By Proposition <ref>, for all k ≥ 0, E_ij[k] ∈.

It follows from its definition that ϕ is a local diffeomorphism in a neighbourhood of the identity element. We first apply the identity that for all i,j ∈_N, E_ij = E_1i^-1E_1j to (<ref>):


    E_ij^+ = (∏_p ∈𝒩_iE_1i^-w_ipE_1p^w_ip)^-1/KE_ij(∏_q ∈𝒩_jE_1j^-w_jqE_1q^w_jq)^1/K.


Applying Proposition <ref> and Lemma <ref> to (<ref>), we have








    _ij^+    = _ij - 1/K∑_p ∈𝒩_i_ip(_1p - _1i) + 1/K∑_q ∈𝒩_j_jq(_1q - _1j) 
       = _ij - 1/K(∑_p ∈𝒩_i_ip_1p - (∑_p ∈𝒩_i_ip)_1i) 
    		+ 1/K(∑_q ∈𝒩_j_jq_1q - (∑_q ∈𝒩_j_jq)_1j) 
       = _ij + 1/K(ℓ_i ⊗ I_m) - 1/K(ℓ_j ⊗ I_m)
       = _ij + 1/K((ℓ_i - ℓ_j) ⊗ I_m).


Setting i = 1 and “stacking" the last line for all j, we obtain


    ^+ = ((I_N + 1/K(1_Nℓ_1 - L)) ⊗ I_m) .


Thus, the local error dynamics are linear. It is interesting to note that the form of (<ref>) implies that the dynamics on each one-parameter subgroup are decoupled. The eigenvalues of the state matrix in (<ref>) are the m-times-repeated eigenvalues of I + (1_Nℓ_1 - L)/K <cit.>.The linear dynamics (<ref>) are (exponentially) stable if and only if the matrix I + (1_Nℓ_1 - L)/K is Schur. We must therefore establish conditions on the gain K such that all eigenvalues of I + (1_Nℓ_1 - L)/K are in the open unit disc.

The Laplacian L of the graph 𝒢 is positive semidefinite, with a zero eigenvalue of algebraic multiplicity equal to the number of connected components in 𝒢 <cit.>; the eigenvector associated with the 0 eigenvalue is 1_N.



The spectrum of 1_Nℓ_1 - L equals σ(-L).



Let J be the Jordan form of L and let V ∈^N × N be the nonsingular matrix such that J = V^-1LV, where the first column V_1 is in the span of 1_N. We have


    V^-1(1_Nℓ_1 - L)V = V^-11_Nℓ_1 V - J.


Since V_1 is in the span of 1_N and V^-1V = I, we have (V^-11_N)_i = 0 for all i ≠ 1. Also because V_1 is in the span of 1_N, we have (ℓ_1V)_1 = 0. Therefore, V^-11_Nℓ_1 V is strictly upper triangular. Therefore, the eigenvalues of (<ref>) are its diagonal elements, which are the diagonal elements of -J, which are the negatives of the eigenvalues of L.



The spectrum of I + (1_Nℓ_1 - L)/K is the image of 1 - σ(L)/K.



The result follows from Lemma <ref> and applying the Spectral Mapping Theorem <cit.> using the function f : →, f(x) = 1 - x/K.


Since the graph is assumed to be connected, L has a simple eigenvalue at 0. By Lemma <ref>, this eigenvalue gets mapped to 1 in the spectrum of I + (1_Nℓ_1 - L)/K.

Let λ be an eigenvalue of L and define the function f(x) = 1 - x/K as in the proof of Lemma <ref>. Applying this function to λ we have


    f(λ) = 1 - |λ|/K^j∠λ 
    	= (1 - |λ|/Kcos(∠λ)) -	j|λ|/Ksin(∠λ)


For stability, we require f(λ) to be in the open unit disc. The squared magnitude of f(λ) is


    |f(λ)|^2    = (1 - |λ|/Kcos(∠λ))^2 + (|λ|/Ksin(∠λ))^2 
    
    
    	   = (|λ|/K)^2 - 2|λ|/Kcos(∠λ) + 1


Then |f(λ)|^2 < 1 if, and only if


    (|λ|/K)^2 - 2|λ|/Kcos(∠λ) < 0.


Since we have already addressed the simple eigenvalue at 0, we assume that λ≠ 0. Therefore, dividing by |λ| we seek conditions on K such that, for all λ∈σ(L)\0, 2cos(∠λ) > |λ|/K. This is equivalent to the condition


    (∀λ∈σ(L)\0)    K > |λ|/2cos(∠λ) =|λ|^2/2Re(λ).


If (<ref>) holds, then all eigenvalues of I + (1_Nℓ_1 - L)/K, except the single eigenvalue at 1, are in the open unit disc.

The next result provides a lower bound on the controller gain K, as a function of the number of agents N, using the properties of the eigenvalues of the Laplacian of a directed graph <cit.>.


If K > K_min(N), where


    K_min(N) {[              N/2           N ≤ 9,; 1/8^2(π/2N)(π/N)     10 ≤ N ≤ 18,;            N - 1          N ≥ 19, ].


then I + (1_Nℓ_1 - L)/K has a single eigenvalue at 1 and all others in the open unit disc.



See Appendix.


The results of <cit.> allow us to find a tighter bound on K than the Gershgorin Disc Theorem, which is used, for example, in <cit.>.


If 𝒢 is symmetric, then σ(L) ⊂ [0,N]. Thus (<ref>) in Lemma <ref> simplifies to K_min(N) = N/2.
▴


By Lemma <ref>, there is no K for which I + (1_Nℓ_1 - L)/K is Schur. However, this does not preclude stability of (<ref>), because the
eigenvalue of 1 corresponds to the dynamics of _11, the error of agent 1 with itself, which is identically zero.


Let  be a generalized cylinder. If the gain K of each agent's controller (<ref>) satisfies (<ref>), then the equilibrium = 0_Nm of (<ref>) is locally exponentially stable
  
  



Since E_11(t) = X_1^-1(t)X_1(t) ≡ I, it follows immediately that _11(t) ≡0_m. Therefore the (N - 1)m dimensional subspace 𝒱∈^Nm: _11 = 0_m is invariant under the dynamics (<ref>). As a result, we have


    σ((I_N + (1_Nℓ_1 - L)/K) ⊗ I_m) 
    	 = σ((I_N + (1_Nℓ_1 - L)/K) ⊗ I_m | 𝒱) ⊔{1,…,1_m times},


where (I_N + (1_Nℓ_1 - L)/K) ⊗ I_m | 𝒱 is the restriction to the subspace 𝒱. If the gain K of each agent's controller (<ref>) satisfies (<ref>), then by Lemma <ref> (I_N + (1_Nℓ_1 - L)/K) ⊗ I_m | 𝒱 is Schur.
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  


By Proposition <ref>, analogous results hold for all equilibria. We emphasize that Theorem <ref> does
not rely on Jacobian linearization of the nonlinear dynamics E_ij^+. The system in exponential coordinates evolves according to linear dynamics.




§ PERFORMANCE WITH AN UNWEIGHTED COMPLETE GRAPH ON A GENERALIZED CYLINDER



We define the ε settling time of error E_ij to be the smallest k∈ such that, for all k ≥k, E_ij[k] = E_ij[0]^α, where |α| ≤ε, 0 < ε < 1.


If 𝒢 is complete, then the ε settling time, where ε∈ (0,1), is


    T_s = ⌈ε/(|K - N|/K)⌉.





    E_ij^+    = (E_ij∏_p ∈ℕ_N ∖{i,j}E_ip)^-1/KE_ij(E_ij^-1∏_p ∈ℕ_N ∖{i,j}E_jp)^1/K
    
    	   = E_ij^K - 2/K(∏_p ∈ℕ_N ∖{i,j}E_ij^-1)^1/K
       = E_ij^K - 2/KE_ij^-N - 2/K
       = E_ij^K - N/K


Thus E_ij[k] = E_ij[0]^(K - N/K)^k.






Therefore, the ε settling time is computed thus


    |K - N/K|^k = εk = ε/(|K - N|/K)


Since k is a time-step, we round up to the nearest integer.


The derivative of the settling time with respect to K is


    ∂ T_s/∂ K   = (ε)(|K - N|^2 + K(N - K))/K|K - N|^2((|K - N|/K))^2
    
    
    	   = (ε)N(N - K)/K|K - N|^2((|K - N|/K))^2.


If K > N, then (<ref>) is positive, so increasing K, i.e., reducing the gain 1/K, delays synchronization, which agrees with intuition. But, interestingly, if K < N, then (<ref>) is negative, so increasing K hastens synchronization. Although (<ref>) is undefined at K = N, these observations suggest that K = N is the minimizer of T_s.


If 𝒢 is complete and K = N, then synchronization is achieved at time-step k = 1.



Setting K = N in (<ref>), we have E_ij^+ = I.





§ GENERAL LIE GROUPS




For our purposes, the only difference between a generalized cylinder and any other Lie group is commutativity. Commutativity is the key property yielding Proposition <ref> and Lemma <ref>, from which all subsequent results follow. We now show that in a neighbourhood of the identity of any Lie group , commutativity holds approximately. Which has the very important implication that all our results for generalized cylinders hold mutatis mutandis on any Lie group in a neighbourhood of the identity. In particular, we obtain Theorem <ref>.

The Baker-Campbell-Hausdorff (BCH) formula relates the product of two elements on the Lie group  to an analytic function of their principal logarithms. If A,B ∈, then the BCH formula has the series representation:


    (exp(A)exp(B))    = A + B + 1/2[A,B] + 1/12[A,[A,B]] 
    	 - 1/12[B,[A,B]] + ⋯,


where the remaining terms are nested brackets of increasing order <cit.>. We will use (<ref>) to derive a linear approximation of the error dynamics on an arbitrary Lie group  near the identity, or equivalently, near the origin on the associated Lie algebra .


The linearization of the BCH formula at the origin of  is (exp(A)exp(B)) ≈ A + B.



All nonlinear terms in (<ref>) are of the form [A,ad_A^k(B)] and [B,ad_B^k(A)], k ∈_≥ 0. Direct computation verifies


    ∂/∂ A[A,ad_A^k(B)] = [A,∂ad_A^k(B)/∂ A].


Thus, if (<ref>) is linearized at the origin, then all nonlinear terms vanish, and (exp(A)exp(B)) ≈ A + B.



Near the identity, we have exp(A + B) ≈exp(A)exp(B) ≈exp(B)exp(A).


Thus, commutativity is satisfied approximately in a neighbourhood of the identity of any Lie group . Therefore, all our results for generalized cylinders apply mutatis mutandis to arbitrary matrix Lie groups.


Given _i and _j in a sufficiently small neighbourhood of zero, we have ϕ(_i + _j) ≈ϕ(_i)ϕ(_j).



For  sufficiently small, ϕ(/K) ≈ϕ(t)^1/K.


We state the analogues of key results for generalized cylinders. Their proofs, as well as the proof of Theorem <ref>, are straightforward applications of Corollaries <ref> and <ref>.


The equilibrium {E_ij = I : i,j ∈_N} is isolated.



Every equilibrium has the same stability properties as the identity.


To illustrate how the proofs of these analogues differ, we present the difference between the proofs of Propositions <ref> and <ref>.


The proof differs from that of Proposition <ref> in only one line of the arithmetic. By Corollaries <ref> and <ref>,


    E̅_ij^+    ≈(∏_p ∈𝒩_iE̅_1i^-_ipE̅_1p^_ip)^-1/KE̅_ij(∏_q ∈𝒩_jE̅_1j^-_jqE̅_1q^_jq)^1/K(∏_p ∈𝒩_iΞ_1i^-_ipΞ_1p^_ip)^-1/K(∏_q ∈𝒩_jΞ_1j^-_jqΞ_1q^_jq)^1/K.


The rest of the proof is identical.





§ SIMULATIONS







 §.§ Comparison with Kuramoto Network on
  SO(2)



The Lie group 2 is one dimensional, thus, it is a one-parameter subgroup of n for any n ≥ 2. 2 is the group of rotations in the plane, which can be interpreted locally as a position on the circumference of a circle. Given an element R ∈2, its local coordinate ∈ is often called the “phase” or “angle”. The Kuramoto oscillator is a popular model of synchronization of networks of oscillators. We can view a Kuramoto network of N agents as a control system, where agent i has phase θ_i ∈ with dynamics


    θ̇_i = u_i,      u_i = -∑_j ∈𝒩_ia_ijsin(θ_i - θ_j),


where a_ij∈ is the coupling strength between agents i and j. System (<ref>) can be modelled as a system on a Lie group in the form of (<ref>), where


    R_i = ϕ(θ) = [  cos(θ) -sin(θ);  sin(θ)  cos(θ) ],   Ṙ_i = R_i[  0 -1;  1  0 ]u_i.


We simulate using N = 3 and a_ij = 1 for all i,j ∈_N. It can be shown that with this choice of parameters, that (<ref>) achieves phase synchronization <cit.>. Sampling with period T = 0.1, we see in Figure <ref> that synchronization is preserved under sampling. But in Figure <ref>, we see that sampling with period T = 0.8, that sampling destroys synchronization.

 
 


 
 


We simulate this network again using the proposed controller with K = 2 and T = 0.8.





As seen in Figure <ref>, synchronization is achieved at T = 0.8, whereas it was lost using the naïvely discretized Kuramoto coupling.




 §.§ Deadbeat Performance on SO(2)


To illustrate Proposition <ref>, we simulate a network with a complete connectivity graph on 2 with K = N = 40, T = 1 and initial phases evenly spaced from -π/(N + 1) to π/(N + 1):


    θ_i  -π/41 + i2π/1599.


As seen in Figure <ref>, all error phases are driven to zero in a single time-step.








 §.§ Network on SU(2)


We simulate a network on 2 to demonstrate Theorem <ref> on a complex, non-commutative Lie group. We simulate a network with N = 6, K = 3.5, and graph Laplacian












    L = [  0.5 -0.1 -0.1 -0.1 -0.1 -0.1;    0  0.8 -0.2 -0.2 -0.2 -0.2;    0    0  0.9 -0.3 -0.3 -0.3;    0    0    0  0.8 -0.4 -0.4;    0    0    0    0  0.5 -0.5;    0    0    0    0    0    0;      ]


The Pauli matrices constitute the canonical basis of 𝔰𝔲(2):


    σ_1 = [ 0 j; j 0 ],   σ_2 = [  0 -1;  1  0 ],   σ_3 = [  j  0;  0 -j ].


We use the Pauli matrices to generate the initial conditions:


    U_i(0) exp(a_iσ_1 + b_iσ_2 + c_iσ_3),


where a_i  -0.32 + i0.6/(N-1), b_i  -0.06 +  i0.3/(N - 1), c_i  -0.42 + i0.6/(N - 1).

For visualization, we plot the Euclidean norms of E_1j - I, j ∈{2,…,N}. As seen in Figure <ref>, the errors tend to identity, thus synchronization is achieved.








§ FUTURE RESEARCH

Future work includes extending our results to agents with dynamic
models and relaxing the assumption that the agents are fully
actuated. The latter could first be addressed by assuming that the Lie
algebra generated by the input vector fields equals . It would
also be of interest to extend our results to time-varying connectivity
graphs.





IEEEtran







§ APPENDIX


Write λ∈ in Cartesian form λσ + jω, σ, ω∈. The eigenvalues of the Laplacian of a directed graph lie in the closed interior of the region in , whose boundary is defined the parametrized curves: c_i(σ) σ + jω_i(σ), i ∈_5, and their complex conjugates c̅_i, i ∈{2,3,4} <cit.>, where the ω_i are defined by the loci:



  * (σ - 1)^2 + ω_1(σ)^2 = (N - 1)^2,

  * ω_2(σ) = (π/N)σ,

  * ω_3(σ) = 1/2(π/2N),

  * ω_4(σ) = (π/N)(N - σ),

  * (σ + 1 - N)^2 + ω_5(σ)^2 = (N - 1)^2.


If N = 2, then this region reduces to the interval [0, N], so K > N/2 implies |f(λ)| < 1 for all λ∈σ(L)∖{0}. If N = 3, then this region reduces to the rhombus with vertices 0, N, and ± jN/2√(3). For 4 ≤ N ≤ 18, the region is a hexagon, as illustrated in Figure <ref>. If N ≥ 19, then the region appears as in Figure <ref>.



 
 


 
 


To lower bound K using only the number of agents N, we maximize the lower bound on K in (<ref>), denoted by g = 0.5|λ|^2/Re(λ) = 0.5(σ^2 + ω^2)/σ, over this region. Since the region is a compact set, the maximum of g is attained at either a critical point or at a point on the boundary of this set. The differential of g is


    dg    = [ ∂ g/∂σ ∂ g/∂ω ]
       = [ 4σ^2 - 2(σ^2 + ω^2)/4σ^2                      ω/σ ]
    
    
       = [ 1/2(1 - ω^2/σ^2)              ω/σ ],


which vanishes nowhere, thus g has no critical points. Therefore, g attains its maximum at a point on the boundary. We parametrize the boundary of the region by σ, and maximize g on this compact, one-dimensional set. Since the Laplacian is a real matrix, its eigenvalues appear in complex conjugate pairs, so it suffices to consider the upper half complex plane.

Let σ_ij denote the value of σ at which locus i intersects locus j. Solving ω_i(σ) = ω_j(σ) for σ, we find:



  * σ_35 = N - 1 - √((N - 1)^2 - (1/2(π/2N))^2),

  * σ_13 = 1 + √((N - 1)^2 - (1/2(π/2N))^2),

  * σ_14 = (N - 1)cos(2π/N) + 1 or N,

  * σ_25 = (N - 1)(1 - cos(2π/N)) or 0,



  * σ_23 = 1/2(1 + (π/N)),

  * σ_34 = N - 1/2(1 + (π/N)),

  * σ_24 = N/2.


This boundary is illustrated for 4 ≤ N ≤ 18 in Figure <ref>, and N ≥ 19 in Figure <ref>.

 



 



The Lie derivatives of g in the direction of c_i are:



  * L_c_1g(σ) = -N(N - 2)/2σ^2,

  * L_c_2g(σ) = 1/2(1 + ^2(π/N)),

  * L_c_3g(σ) = 1/2 - ^2(π/2N)/8σ^2),

  * L_c_4g(σ) = ^2(π/N)/2 - N^2^2(π/N)/2 σ^2,

  * L_c_5g(σ) = 0.


Let σ_i^⋆ denote the value of σ at which L_c_ig vanishes, which are the critical points of the restriction of g to the boundary. We have



  * L_c_1g(σ) = 0 if and only if N = 0 or N = 2,

  * L_c_2g(σ) ≠ 0 for all σ∈, N ∈,

  * σ_3^⋆ = 1/2(π/2N),

  * σ_4^⋆ = Ncos(π/N),

  * L_c_5g(σ) =0 identically.


We determine whether σ_i^⋆ is a local maximum or minimum by examining the second Lie derivative of g evaluated at σ_i^⋆:


    L^2_c_ig(σ) = dL_c_ig(σ)dc_i/dσ
         = [ ∂^2g/∂σ^2 + ∂^2g/∂ω∂σdω_i/dσ + ∂ g/∂ωd^2ω_i/dσ^2                     ∂^2g/∂σ∂ω + ∂^2g/∂ω^2dω_i/dσ ][       1; dω_i/dσ ]
         = ω^2/σ^3 - 2ω/σ^2dω_i/dσ + ω/σd^2ω_i/dσ^2 + 1/σ(dω_i/dσ)^2.


For i = 3 and i = 4 we have


	
  * L^2_c_3g(σ) = ^2(π/2N)/(4σ^3) > 0 for all σ≥ 0, N ≥ 2,
	
  * L^2_c_4g(σ) = N^2^2(π/2)/σ^3 >0 for all σ≥ 0, N ≥ 2.


Thus, the critical points on loci 3 and 4 are minima, thus the maxima of g on these loci restricted to the boundary are attained at the intersection points, which we now characterize.

Let g_i ∈_≥ 0 be the value of g at c_i(σ_i^⋆), let g_ij∈_≥ 0 be the value of g at c_i(σ_ij) = c_j(σ_ij), and let g_N be the value of g at σ = N, ω = 0. Since the value of g is constant on locus 5, we do not consider the intersection points of locus 5 with loci 2 or 3. We have:



  * g_5 = N - 1,

  * g_13 = 1 + N(N - 2)/2(1 + √((N - 1)^2 - (1/2(π/N))^2)),

  * g_14 = 1 + N(N - 2)/2(1 + (N - 1)cos(2π/N)),

  * g_23 = 1/8^2(π/2N)(π/N),

  * g_34 = ^2(π/2 N) + (2N - 1 - (π/N))^2/4(2N - 1 - (π/N)),

  * g_24 = N(^2(π/N)+ 1)/4,

  * g_N = N/2.


Finally, we identify the maximum value of g on the boundary. For N = 3, the boundary is defined by only loci 2 and 4. N = 3 is also the only case in which loci 2 and 4 intersect. In can be shown numerically that for N = 3, that g_N maximizes g. It can be verified numerically that if 4 ≤ N ≤ 9, then g_N is the maximum of g, and if 10 ≤ N ≤ 18, then g_23 is the maximum of g, which proves the first two cases in (<ref>).


 



 


For N ≥ 19, we now establish that max{g_5, g_13, g_14} = g_5. Notice that g_13 and g_14 can be expressed:


    g_13 = 1 + N(N - 2)/2σ_13,      g_14 = 1 + N(N - 2)/2σ_14.


From their definitions, g_5 ≥ g_13 if and only if


    N - 1 ≥ 1 + N(N - 2)/2σ_13σ_13≥N/2.


Similarly, we find that g_5 ≥ g_14 if and only if σ_14 > N/2. By the geometry of the region as discussed in <cit.> and illustrated in Figure <ref>, these inequalities hold for all N ≥ 19.













































