# Composite Distributed Learning and Synchronization of Nonlinear Multi-Agent
Systems with Complete Uncertain Dynamics

Emadodin Jandaghi, Dalton L. Stein, Adam Hoburg, Mingxi Zhou and Chengzhi Yuan
Emadodin Jandaghi, Dalton L. Stein and Chengzhi Yuan are with the Department
of Mechanical, Industrial and Systems Engineering, Kingston, RI 02881, USA
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>Huborg is with
the Department of Electrical, Computer and Biomedical Engineering, Kingston,
RI 02881, USA<EMAIL_ADDRESS>Zhou is with the Graduate School of
Oceanography, , Kingston, RI 02881, USA<EMAIL_ADDRESS>
###### Abstract

This paper addresses the challenging problem of composite synchronization and
learning control in a network of multi-agent robotic manipulator systems
operating under heterogeneous nonlinear uncertainties within a leader-follower
framework. A novel two-layer distributed adaptive learning control strategy is
introduced, comprising a first-layer distributed cooperative estimator and a
second-layer decentralized deterministic learning controller. The primary
objective of the first layer is to facilitate each robotic agent’s estimation
of the leader’s information. The second layer is responsible for both enabling
individual robot agents to track desired reference trajectories and accurately
identifying and learning their nonlinear uncertain dynamics. The proposed
distributed learning control scheme represents an advancement in the existing
literature due to its ability to manage robotic agents with completely
uncertain dynamics including uncertain mass matrices. This framework allows
the robotic control to be environment-independent which can be used in various
settings, from underwater to space where identifying system dynamics
parameters is challenging. The stability and parameter convergence of the
closed-loop system are rigorously analyzed using the Lyapunov method.
Numerical simulations conducted on multi-agent robot manipulators validate the
effectiveness of the proposed scheme. The identified nonlinear dynamics can be
saved and reused whenever the system restarts.

##  I INTRODUCTION

Robotics has many applications from manufacturing to surgical procedures [1,
2, 3]. However, controlling robots in space and underwater is challenging due
to unpredictable robot dynamics in such environments. Also, the increasing
demand for high precision and operational complexity has shifted the focus
towards cooperatively utilizing multiple standard robots. This approach
benefited from improved efficiency, cost reduction, and redundancy [4, 5].
Given these advantages, numerous studies have been done to formulate diverse
decentralized control strategies for coordinating multiple robotic arms
effectively [6].

Despite extensive research, there are several outstanding challenges. One
persistent challenge is the management of model uncertainties which can impair
the performance of distributed control systems [7]. Some studies have examined
the role of uncertainties in robotic arm control [8], however, these did not
extend their findings to the synchronization of multiple robots. [9] explored
a virtual leader-follower strategy for robot manipulators under uncertainties
and disturbances. While a high-gain observer was used for velocity estimation,
this approach risks exciting unmodeled high-frequency dynamics and magnifying
measurement noise. Such issues make controller implementation challenging and
necessitate meticulous parameter tuning. Another drawback in the existing
literature is the assumption of homogeneous system dynamics across all robots,
which is often not the case in real-world applications. Meanwhile, [10]
ignored system uncertainties and assumed full model knowledge. Furthermore,
while some research has considered non-identical robotic systems [9], others
[11] mostly focused solely on achieving adaptive tracking control, without
consideration for the convergence of controller parameters to their optimal
states. Neural networks have been employed to tackle uncertainties, but they
also suffer from the unrealistic assumption of uniform robot dynamics.

Previous literature has not fully tackled the identification of all system
nonlinearities without certain assumptions. These include having a known Mass
or Inertia matrix [12], using large gains to suppress errors caused by the
Coriolis and Centripetal force matrix [8], and assuming a structurally
balanced, signed switching network for multiple systems in their adaptive
control law [13]. In contrast, our work uniquely identifies the nonlinear,
uncertain dynamics of each robot without making any assumptions about the
certainty or structure of each system.

In this study, we make a significant step forward compared to previous
research. This framework is completely independent of any multi-agent system
attributes, rendering this method universally applicable to any nonlinear
Euler-Lagrange (EL) system dynamics. We tackle the challenges of achieving
synchronization control and integrating learning capabilities into multi-robot
systems with heterogeneous nonlinear uncertain dynamics without any knowledge
of each robotic system’s nonlinearities, which operate under a virtual leader-
following framework. These robots have the same kinematics and the same degree
of freedom but have distinct dynamic system parameters like masses, inertia,
and stiffness.

Our control architecture employs a fixed directed graph communication, with a
virtual leader modeled as a Linear Time-Invariant (LTI) system. Only select
agents access the leader’s data directly. The control strategy is dual-
layered: The first layer focuses on cooperative estimation, allowing agents to
have inter-agent communication and share leader-estimated states and system
matrices while keeping their own plant data private (i.e., robot’s
position/velocity). The second layer involves a decentralized adaptive
learning controller, designed to regulate each robot’s state and pinpoint its
unique dynamics using first-layer estimates. The second layer operates without
any data sharing, allowing each local robot to implement its adaptive learning
controller in a completely decentralized way.

Our adaptive learning control law uses precise function approximation with
Radial Basis Function (RBF) Neural Networks (NN) for the identification of
systems with completely uncertain dynamics. This enhances control over
robotics in space and underwater environments, where the system dynamics are
often unpredictable. We confirm the efficacy of the approach through
mathematically rigorous validation and comprehensive simulation studies.
Presenting a thorough convergence analysis for the neural network weights
guarantees that the learned nonlinearities are only stored once and can be
used in the controller after the system is turned off, eliminating the need
for the controller to adjust every time the system is restarted.

The rest of the paper contains the following sections: Section II provides an
initial overview of graph theory and Radial basis function NNs. Distributed
cooperative estimator and decentralized deterministic learning control law are
discussed in Section III. The simulation studies is provided in Section IV.
Finally, Section V concludes the paper.

##  II PRELIMINARIES AND PROBLEM STATEMENT

###  II-A Notation and Graph Theory

We denote the sets of real numbers as  $\mathbb{R}$ and  $\mathbb{R}^{+}$ .
$\mathbb{R}^{m\times n}$ represents the set of real  $m\times n$ matrices, and
$\mathbb{R}^{n}$ is the set of real  $n\times 1$ vectors. The identity matrix
is symbolized as  $\mathbf{I}$ . The vector with all elements being 1 in an
$n$-dimensional space is represented as  $\mathbf{1}_{n}$ . The sets
$\mathbb{S}_{+}^{n}$ and  $\mathbb{S}_{-}^{n+}$ stand for real symmetric
$n\times n$ and positive definite matrices, respectively. A block diagonal
matrix with matrices  $X_{1},X_{2},\dots,X_{p}$ on its main diagonal is
denoted by  $\text{diag}\\{X_{1},X_{2},\ldots,X_{p}\\}$ . The notation
$A\otimes B$ signifies the Kronecker product of matrices $A$ and $B$. For a
matrix  $A$, $\vec{A}$ is the vectorization of $A$ by stacking its columns on
top of each other. Given two integers $k_{1}$ and  $k_{2}$ with $k_{1}<k_{2}$,
$I[k_{1},k_{2}]=\\{k_{1},k_{1}+1,\dots,k_{2}\\}$ . For a vector
$x\in\mathbb{R}^{n}$ , its norm is defined as
$|x|:=(\mathbf{x}^{T}\mathbf{x})^{1/2}$ . For a square matrix  $A$,
$\lambda_{i}(A)$ denotes its $i$ -th eigenvalue, while
$\lambda_{\text{min}}(A)$ and  $\lambda_{\text{max}}(A)$ represent its maximum
and minimum eigenvalues, respectively. A directed graph denoted as  $G=(V,E)$
comprises nodes in the set  $V=\\{1,2,\dots,N\\}$ and edges in  $E\subseteq
V\times V$ . An edge from node $i$ to node $j$ is represented as  $(i,j)$ ,
with $i$ as the parent node and $j$ as the child node. Node $i$ is also termed
a neighbor of node $j$. $N_{i}$ considered as the subset of $V$ consisting of
the neighbors of node $i$. A sequence of edges in $G$,
$(i_{1},i_{2}),(i_{2},i_{3}),\dots,(i_{k},i_{k+1})$ , is called a path from
node $i_{1}$ to node $i_{k+1}$. Node $i_{k+1}$ is reachable from node $i_{1}$.
A directed tree is a graph where each node, except for a root node, has
exactly one parent. The root node is reachable from all other nodes. A
directed graph $G$ contains a directed spanning tree if at least one node can
reach all other nodes. The weighted adjacency matrix of $G$ is a non-negative
matrix  $A=[a_{ij}]\in\mathbb{R}^{N\times N}$ , where $a_{ii}=0$ and
$a_{ij}>0\implies(j,i)\in E$ . The Laplacian of $G$ is denoted as
$L=[l_{ij}]\in\mathbb{R}^{N\times N}$ , where  $l_{ii}=\sum_{j=1}^{N}a_{ij}$
and  $l_{ij}=-a_{ij}$ if $i\neq j$ . It is established that $L$ has at least
one eigenvalue at the origin, and all nonzero eigenvalues of $L$ have positive
real parts. Furthermore, as stated in [14], Lemma 3.3, $L$ has one eigenvalue
at the origin and  $(N-1)$ eigenvalues with positive real parts if and only if
$G$ includes a directed spanning tree.

###  II-B Radial Basis Function NNs

The RBF Neural Networks (NN) can be described as
$f_{nn}(Z)=\sum_{i=1}^{N}w_{i}s_{i}(Z)=W^{T}S(Z)$ , where
$Z\in\Omega_{Z}\subseteq\mathbb{R}^{q}$ and
$W=w_{1},...,w_{N}^{T}\in\mathbb{R}^{N}$ as input and weight vectors
respectively.  $N$ indicates the number of NN nodes,
$S(Z)=[s_{1}(||Z-\mu_{i}||),...,s_{N}(||Z-\mu_{i}||)]^{T}$ with
$s_{i}(\cdot)$ is a RBF, and  $\mu_{i}(i=1,...,N)$ is distinct points in the
state space. The Gaussian function
$s_{i}(||Z-\mu_{i}||)=\operatorname{exp}\left[-\frac{(Z-\mu_{i})^{T}(Z-\mu_{i})}{\eta_{i}^{2}}\right]$
is generally used for RBF, where
$\mu_{i}=[\mu_{i1},\mu_{i2},...,\mu_{iN}]^{T}$ is the center and  $\eta_{i}$
is the width of the receptive field. The Gaussian function categorized by
localized radial basis function  $s$ in the sense that
$s_{i}(||Z-\mu_{i}||)\rightarrow 0$ as  $||Z||\rightarrow\infty$ . It has been
shown in [15, 16] for any continuous function  $f(Z):\Omega_{Z}\rightarrow R$
where $\Omega_{Z}\subset R^{p}$ is a compact set. For the NN approximator,
where the node number N is sufficiently large, there exists an ideal constant
weight vector  $W^{*}$ , such that for each  $\epsilon^{*}>0$ ,
$f(Z)=W^{*T}S(Z)+\epsilon(Z)$ , $\forall Z\in\Omega_{Z}$, where $\epsilon(Z)$
is the approximation error. The following lemma describes the PE (Persistent
Excitation) condition for RBF NN from [17].

######  Lemma 1.

Consider any continuous recurrent trajectory  $Z(t):[0,\infty)\rightarrow
R^{q}$. $Z(t)$ remains in a bounded compact set
$\Omega_{Z}\subset\mathbb{R}^{q}$ . For RBF NN  $W^{T}S(Z)$ with centers
placed on a regular lattice (large enough to cover compact set  $\Omega_{Z}$
), the regressor subvector  $S_{\zeta}(Z)$ consisting of RBFs with centers
located in a small neighborhood of $Z(t)$ is persistently exciting.

###  II-C Problem Statement

In this paper, we consider a multi-robot manipulator system of $N$ robots with
heterogeneous uncertain nonlinear dynamics described with Euler-Lagrange (EL)
mechanical systems as:

$M_{i}(\theta_{i})\ddot{\theta}_{i}+C_{i}(\theta_{i},\dot{\theta}_{i})\dot{\theta}_{i}+g_{i}(\theta_{i})=\tau_{i},\quad
i\in\mathbf{I}[1,N],$ (1)

where  $\tau_{i}\in\mathbb{R}^{n}$ is the vector of input signals. The
subscript $i$ denotes the $i^{\textrm{th}}$ robotic agent. For each
$i\in\mathbf{I}[1,N]$ ,
$\theta_{i}=[\theta_{i1},\theta_{i2},...,\theta_{in}]^{T}\in\mathbb{R}^{n}$
are the joint positions,  $\dot{\theta}_{i},\ddot{\theta}_{i}$ are the joint
velocities and accelerations, respectively.  $M\in\mathbb{S}_{+}^{n}$ is a
positive definite mass matrix,  $C\in\mathbb{R}^{n\times n}$ is the matrix
containing Coriolis and Centripetal forces,  $g\in\mathbb{R}^{n}$ is the
gravity term. It is crucial to note that as opposed to our previous work [12],
we make no assumptions about the certainty of any system matrices ($M,C,$ and
$g$) throughout designing our distributed controller and consider all system
matrices to be uncertain. We can rewrite the system dynamics as:

$\displaystyle\dot{x}_{i1}=x_{i2}$ (2)
$\displaystyle\dot{x}_{i2}=M_{i}^{-1}(x_{i1})[\tau_{i}-C_{i}(x_{i})x_{i2}-g_{i}(x_{i1})]$

where $x_{i1}=\theta_{i}$, $x_{i2}=\dot{\theta}_{i}$, and let $x_{i}$ be a
column vector of $x_{i1}$ and $x_{i2}$. In a leader-follower setting, the
leader’s dynamics is:

$\displaystyle\dot{\chi}_{0}=A_{0}\chi_{0}$ (3)

Here, ”0” marks the leader node. $\chi_{0}$ is a column vector of $x_{01}$ and
$x_{02}$, which are the leader’s states. $A_{0}$ is a constant system matrix.
Only robots directly linked to the leader have access to $\chi_{0}$ and
$A_{0}$. For a multi-robot system with $N$ follower robots and one leader, we
describe an adjacency matrix $A$ and set $V=\\{0,1,\ldots,N\\}$, with ”0”
being the leader. We make two assumptions, which do not sacrifice generality:

######  Assumption 1.

All eigenvalues of $A_{0}$ are imaginary.

######  Assumption 2.

The digraph $G$ has a directed spanning tree rooted at node 0.

Assumption 1 ensures all the states of the leader dynamics remain periodic and
uniformly bounded. The Laplacian $L$ of the graph can be divided as follows
using the Assumption 2:

$L=\begin{bmatrix}\sum_{j=1}^{N}a_{0j}&-[a_{01},\ldots,a_{0N}]\\\
-\Phi_{1N}&H\end{bmatrix}$

Here, $\Phi$ is a diagonal matrix and $H$ has only positive real parts for its
non-zero eigenvalues. The main problem to be addressed in this paper:

######  Problem 1.

Given the robotic system with $N$ agents as described in (2) and a leader as
(3), along with Assumption 1 and Assumption 2, we aim to create a control plan
that uses only local data for each robot. The goals are:

  1. 1.

Cooperative Synchronization: All robots should synchronize with the common
path set by the leader, meaning
$\lim_{{t\to\infty}}(x_{i1}(t)-x_{01}(t))=0\quad\forall i\in\mathbf{I}[1,N]$ .

  2. 2.

Decentralized Learning: Each robot should learn its own nonlinear uncertain
dynamics via its local adaptive controller using RBF NNs.

######  Remark 1.

Unlike existing works that only focus on synchronization, our problem also
aims for both control and learning of robotic dynamics. This involves
challenges like diverse robot behaviors, limited information from the leader,
and complete uncertain nonlinear dynamics.

A two-layer framework is proposed as follows:

  * •

The first layer uses a distributed estimator to figure out the leader’s state
and system information.

  * •

The second layer uses a decentralized controller for synchronization and
learning the local robot dynamics.

Only the first layer requires data sharing between nearby robots and the
second layer just works on local data. The next section will explore each
layer. In addition to a Theorem 1 on the stability and convergence of the
first layer, we will present two crucial theorems. Theorem 2 will address the
second layer’s learning and control performance, and Theorem 3 will combine
insights from both layers to provide conclusions.

##  III DISTRIBUTED ADAPTIVE CONTROL FRAMEWORK

###  III-A First Layer: Distributed Cooperative Estimator

In a distributed control setting, each robotic agent can access only its own
data and the data of its immediate neighbors. Consequently, the leader’s
parameters, such as $\chi_{0}$ and $A_{0}$, may not be available to all
agents. This constraint leads us to create a distributed cooperative
estimator, enabling each agent to estimate the leader’s state and dynamics
through inter-agent collaboration:

$\small\dot{\hat{\chi}}_{i}(t)=A_{i}(t)\hat{\chi}_{i}(t)+\beta_{i1}\sum_{j=0}^{N}a_{ij}(\hat{\chi}_{j}(t)-\hat{\chi}_{i}(t)),\
\forall i\in I[1,N]\ $ (4)

The observer states for each robot, represented by
$\hat{\chi}_{i}=[\hat{x}_{i1},\hat{x}_{i2}]^{T}$, are used to estimate the
leader’s state $\chi_{0}=[x_{01},x_{02}]^{T}$. $\beta_{i1}$ are positive
design constant numbers. The time-varying system parameters $\hat{A}_{i}(t)$
are updated by the following equation for all $i$:

$\small\dot{\hat{A}}_{i}(t)=\beta_{i2}\sum_{j=0}^{N}a_{ij}(\hat{A}_{j}(t)-\hat{A}_{i}(t)),\quad\forall
i\in I[1,N]$ (5)

The matrices $\hat{A}_{i}$ are used to estimate the leader’s system matrix
$A_{0}$ and have dimensions $n\times n$. The constants $\beta_{i2}$ are all
positive numbers. We observe that only local information, such as the observer
states $\hat{\chi}_{i}$ and $\hat{A}_{i}$, need to be shared among neighboring
robotic agents. For each agent $i$ in the range $[1,N]$, we define the local
errors as $\tilde{\chi}_{i}=\hat{\chi}_{i}-\chi_{0}$ and
$\tilde{A}_{i}=\hat{A}_{i}-A_{0}$. This allows us to describe the local error
dynamics for agent $i$:

$\displaystyle\dot{\tilde{\chi}}_{i}(t)$
$\displaystyle=A_{0}\tilde{\chi}_{i}(t)+\tilde{A}_{i}(t)\tilde{\chi}_{i}(t)+\tilde{A}_{i}(t)\chi_{0}(t)$
(6)
$\displaystyle+\beta_{i1}\sum_{j=0}^{N}a_{ij}(\tilde{\chi}_{j}(t)-\tilde{\chi}_{i}(t)),\quad\forall
i\in I[1,N]$ $\displaystyle\dot{\tilde{A}}_{i}(t)$
$\displaystyle=\beta_{i2}\sum_{j=0}^{N}a_{ij}-(\tilde{A}_{j}(t)-\tilde{A}_{i}(t)),\quad\forall
i\in I[1,N].$

Now consider
$\tilde{\chi}=\text{col}\\{\tilde{\chi}_{1},\ldots,\tilde{\chi}_{N}\\}$,
$\tilde{A}=\text{col}\\{\tilde{A}_{1},\ldots,\tilde{A}_{N}\\}$,
$\tilde{A}_{b}=\text{diag}\\{\tilde{A}_{1},\ldots,\tilde{A}_{N}\\}$,
$B_{\beta_{1}}=\text{diag}\\{\beta_{11},\ldots,\beta_{N1}\\}$, and
$B_{\beta_{2}}=\text{diag}\\{\beta_{12},\ldots,\beta_{N2}\\}$. Using these
definitions, the error dynamics for the entire network system are:

$\displaystyle\dot{\tilde{\chi}}(t)$ $\displaystyle=((I_{N}\otimes
A_{0})-B_{\beta_{1}}(H\otimes I_{2n}))\tilde{\chi}(t)$ (7)
$\displaystyle+\tilde{A}_{b}(t)\otimes\tilde{\chi}(t)+\tilde{A}_{b}(t)(1_{N}\otimes\chi_{0}(t)),$
$\displaystyle\dot{\tilde{A}}(t)$ $\displaystyle=-B_{\beta_{2}}(H\otimes
I_{n})\tilde{A}(t).$

Here, $H$ is specified in II-C. We rely on Assumptions 1 and 2 to support the
following Theorem:

######  Theorem 1.

Given the error system (7) and under Assumptions 1 and 2, for all $i\in
I[1,N]$ and any initial conditions
$\chi_{0}(0),\hat{\chi}_{i}(0),\hat{A}_{i}(0)$, we have
$\lim_{t\to\infty}\tilde{A}_{i}(t)=0$ and
$\lim_{t\to\infty}\tilde{\chi}_{i}(t)=0$ exponentially.

The proof of Theorem 1 is similar to that of [12], which will be omitted here
to save some space.

###  III-B Second Layer: Decentralized Deterministic Learning Control

This section introduces our new decentralized deterministic learning control
law designed to improve synchronization performance. Given that the leader’s
state information $\chi_{0}$ is unavailable to all robotic agents,
$\hat{\chi}_{i}$ will serve as the tracking reference signal for tracking
control design and implementation for each agent. Additionally, all system
matrices for each agent are considered to be unstructurally uncertain.
Meaning, unlike in previous studies, we make no assumptions about agents’
system structures, and the tracking controller is implemented independently on
each agent. Therefore, the learning control objective is met when each agent’s
system states $x_{i}$ converge to the tracking reference signal
$\hat{\chi}_{i}$ in a decentralized manner. For the stability analysis, we
implement the EL system property below:

######  Property 1.

The Matrix  $\dot{M}_{i}(\theta_{i})-2C_{i}(\theta_{i},\dot{\theta}_{i})$ is
Skew-symmetric which means
$\frac{1}{2}\dot{M}_{i}(\theta_{i})=2C_{i}(\theta_{i},\dot{\theta}_{i})$ .

To introduce the decentralized deterministic learning control law, consider
the $i^{th}$ robotic agent with a reference signal in its controller as
$\dot{x}_{ri}$ and $\ddot{x}_{ri}$, while assuming a filtered output signal
$r_{i}$ be:

$\displaystyle r_{i}$ $\displaystyle=\dot{x}_{i1}-\dot{x}_{ri}$ (8)
$\displaystyle=\dot{e}_{i}+\lambda_{i}e_{i}\quad\forall i\in I[1,N]$

where $\lambda>0$ and $e_{i}\in\mathbb{R}^{n}$ is the position tracking error:

$\displaystyle e_{i}=x_{i1}-\hat{x}_{i1},\quad\forall i\in I[1,N]$ (9)

Therefore:

$\displaystyle\dot{x}_{ri}$
$\displaystyle=\dot{\hat{x}}_{i1}-\lambda({e_{i}})\quad\forall i\in I[1,N]$
(10)

where $\dot{x}_{ri2}$ is the second derivative of reference signal with:

$\small
M_{i}(x_{i1})\ddot{x}_{ri}+C_{i}(x_{i},\dot{x}_{i})\dot{x}_{ri}+g_{i}(x_{i})=H_{i}(\chi_{i})$
(11)

where  $\chi_{i}=col\\{x,\dot{x},\dot{x}_{r},\ddot{x}_{r}\\}$ . The unknown
nonlinear function $H_{i}(\chi_{i})$ is aimed to be approximated using the
universal approximation capability of RBF NN. Specifically, according to II-B,
there exist RBF NNs ${W^{*}_{i}}^{T}S_{i}(\chi_{i})$ such that

$\displaystyle
H_{i}(\chi_{i})=W_{i}^{*T}S_{i}(\chi_{i})+\epsilon_{i}(\chi_{i})\quad\forall
i\in I[1,N]$ (12)

with $W^{*}_{i}$ as the ideal constant weights, and
$|\epsilon_{i}(\chi_{i})|\leq\epsilon^{*}_{i}$ is the ideal approximation
errors which can be made arbitrarily small given the sufficiently large number
of neurons. Assuming $\hat{W}_{i}$ as the estimate of $W^{*}_{i}$, we
construct the implementable decentralized deterministic learning control law
as:

$\tau_{i}=\hat{W}_{i}^{T}S_{i}(\chi_{i})-K_{i}r_{i},\quad\forall i\in I[1,N].$
(13)

where
$\hat{W}_{i}^{T}S_{i}(\chi_{i})=[\hat{W}_{i1}^{T}S_{i1}(\chi_{i}),\ldots,\hat{W}_{in}^{T}S_{in}(\chi_{i})]^{T}$
and $K_{i}\in\mathbb{R}^{n\times n}$ are used to approximate the unknown
nonlinear function vector $H_{i}(\chi_{i})$. A robust self-adaptation law for
online updating $\hat{W}_{i}$ is constructed as:

$\displaystyle\dot{\hat{W}}_{ij}=-\Gamma_{i}[S_{ij}(\chi_{i})r_{ij}+\sigma_{i}\hat{W}_{ij}]$
(14) $\displaystyle\quad\forall i\in I[1,N]\text{ and }\forall j\in I[1,n].$

where $\Gamma_{i}$ and $\sigma_{i}$ are positive constant numbers with
$\sigma_{i}$ being very small, $r_{ij}$ represents the $j^{th}$ element of
$r_{i}$. Substituting (11), (12) and (13) into (8) for all $i\in I[1,N]$
yields:

$\displaystyle\dot{r}_{i}=M_{i}^{-1}(x_{i})(\tilde{W}_{i}^{T}S_{i}(\chi_{i})-\epsilon_{i}(\chi)-K_{i}r_{i}-C_{i}(x_{i},\dot{x}_{i})r_{i})$
(15)

where $\tilde{W}_{i}:=\hat{W}_{i}-W^{*}_{i}$. Based on (14) and (15), we have
the following theorem.

######  Theorem 2.

Given systems (14) and (15).If there exists a sufficiently large compact set
$\Omega_{\chi_{i}}$ such that  $\chi_{i}\in\Omega_{\chi_{i}}\ \forall
i\in\mathbf{I}[1,N]$ , then for any bounded initial conditions with
$\hat{W}_{i}(0)=0(\forall i\in\mathbf{I}[1,N])$ we have: (i) all the signals
in the system remain uniformly bounded; (ii) the position tracking error
$x_{i1}-\hat{x}_{i1}$ converges exponentially to a small neighborhood around
the origin, by choosing the design parameters with
$K_{i}\in\mathbb{S}_{+}^{n}\ \forall i\in\mathbf{I}[1,N]$ . (iii) along the
system trajectory denoted by  $\varphi(\chi_{i}(t))|_{t\geq T_{i}}$ starting
from $T_{i}$ which represents the settling time of tracking control, the local
estimated neural weights  $\hat{W}_{i\psi}$ converge to small neighborhoods
close to the corresponding ideal values  $W_{i\phi}^{*}$ , and locally-
accurate identification of nonlinear uncertain dynamics defined in (11) can be
obtained by  $\hat{W}^{T}_{i}S_{i}(\chi_{i})$ as well as
$\bar{W}_{i}^{T}S_{i}(\chi_{i})$ along the system trajectory
$\psi(\chi_{i}(t))|_{t\geq T_{i}}$ , where

$\displaystyle\bar{W}_{i}^{T}=mean_{t\in[t_{ia},t_{ib}]}\hat{W}_{i}(t),\quad\forall
i\in\mathbf{I}[1,N],$ (16)

with  $[t_{ia},t_{ib}](t_{ib}>t_{ia}>T_{i})$ being a time segment after the
transient period of tracking control.

Proof. (i) Given the systems (14) and (15) consider the following Lyaponuv
function candidate,

$\small
V_{i}=\frac{1}{2}r_{i}^{T}M_{i}(x_{i})r_{i}+\frac{1}{2}\sum_{j=1}^{n}\tilde{W}_{ij}^{T}\Gamma^{-1}\tilde{W}_{ij}.$
(17)

Temporal differentiation of (17) yields,

$\small\dot{V}_{i}=r_{i}^{T}M_{i}(x)\dot{r}_{i}+\frac{1}{2}r_{i}^{T}\dot{M}_{i}(x)r_{i}+\sum_{j=1}^{n}\dot{\tilde{W}}_{ij}^{T}\Gamma_{i}^{-1}\tilde{W}_{ij}.$
(18)

Exploiting (15) and the assumption in Property 1 in (18) yields

$\small\begin{split}\begin{aligned}
\dot{V}_{i}&=-r_{i}^{T}K_{i}r_{i}-r_{i}^{T}\epsilon(\chi_{i})+r_{i}^{T}\tilde{W}_{i}^{T}S_{i}(\chi)\\\
&+\sum_{j=1}^{n}\dot{\tilde{W}}_{ij}^{T}\Gamma_{i}^{-1}\tilde{W}_{ij}.\end{aligned}\end{split}$
(19)

Utilizing the fact that  $\dot{\tilde{W}}=\dot{\hat{W}}$ in (14), and
substituting (14) into (19),

$\small\begin{split}\begin{aligned}
\dot{V}_{i}&=-r_{i}^{T}K_{i}r_{i}-r_{i}^{T}\epsilon(\chi_{i})+r_{i}^{T}\tilde{W}_{i}^{T}S_{i}(\chi_{i})\\\
&+\sum_{j=1}^{n}\left(-\Gamma_{i}[S_{ij}(\chi_{i})r_{ij}+\sigma_{i}\hat{W}_{ij}]\right)^{T}\Gamma_{i}^{-1}\tilde{W}_{ij}\\\
&=-r_{i}^{T}K_{i}r_{i}-r_{i}^{T}\epsilon(\chi_{i})-\sum_{j=1}^{n}\sigma_{i}\hat{W}_{ij}^{T}\tilde{W}_{ij}\end{aligned}\end{split}$
(20)

now, we select  $K_{i}=K_{i1}+K_{i2}$ with $K_{i1}$ and $K_{i2}$ being
positive definite to yield:

$\displaystyle\dot{V}_{i}$
$\displaystyle=-r_{i}^{T}K_{i1}r_{i}-r_{i}^{T}K_{i2}r_{i}-r_{i}^{T}\epsilon(\chi_{i})-$
(21) $\displaystyle\sum_{j=1}^{n}\sigma_{i}\hat{W}_{i}^{T}j\tilde{W}_{ij}.$

As a result, we have for all  $i\in\mathbf{I}[1,N]$ and  $j\in\mathbf{I}[1,n]$
:

$\displaystyle-\sigma_{i}\hat{W}_{i}^{T}j\tilde{W}_{ij}$
$\displaystyle=-\sigma_{i}(\tilde{W}_{ij}+W^{*}_{ij})^{T}\tilde{W}_{ij}$ (22)
$\displaystyle\leq-\sigma_{i}\tilde{W}_{ij}^{2}-\sigma_{i}\tilde{W}_{ij}W^{*}_{ij}$
$\displaystyle\leq-\frac{\sigma_{i}}{2}\|\tilde{W}_{ij}\|^{2}+\frac{\sigma_{i}}{2}\|W^{*}_{ij}\|^{2}.$

Using the same approach we have:

$\small-
r_{i}^{T}K_{i2}r_{i}-r_{i}^{T}\epsilon_{i}\leq\frac{\epsilon_{i}^{T}\epsilon_{i}}{4\lambda_{\text{min}}(K_{i2})}\leq\frac{||\epsilon_{i}^{*}||^{2}}{4\lambda_{\text{min}}(K_{i2})}.$
(23)

Based on (22) and (23), we conclude:

$\displaystyle\dot{V}_{i}\leq$ $\displaystyle-
r_{i}^{T}K_{i1}r_{i}+\frac{||\epsilon_{i}^{*}||^{2}}{4\lambda_{\text{min}}(K_{i2})}$
(24)
$\displaystyle-\frac{1}{2}\sum_{j=1}^{n}\sigma_{i}\|\tilde{W}_{ij}\|^{2}+\frac{1}{2}\sum_{j=1}^{n}\sigma_{i}\|W^{*}_{ij}\|^{2}.\quad\forall
i\in I[1,N]$

Based on (24), and given that $K_{i1}$ is positive definite, we can conclude
that $\dot{V}_{i}$ is negative definite under the condition that:

$\displaystyle\|r_{i}\|$
$\displaystyle>\frac{\|\epsilon^{*}_{i}\|}{\sqrt{2\lambda_{\text{min}}(K_{i1})\lambda_{\text{min}}(K_{i2})}}+$
(25)
$\displaystyle\sqrt{\frac{2\sigma_{i}}{\lambda_{\text{min}}(K_{i1})}}\sum_{j=1}^{n}\|W^{*}_{ij}\|$
$\displaystyle\text{or}\quad\sum_{j=1}^{n}\|\tilde{W}_{ij}\|$
$\displaystyle>\frac{\|\epsilon_{i}\|}{\sqrt{2\sigma_{i}\lambda_{\text{min}}(K_{i2})}}+\sum_{j=1}^{n}\|W^{*}_{ij}\|.\quad\forall
i\in I[1,N]$

The signals $r_{i}$ and  $\tilde{W}_{ij}$ are bounded, leading to boundedness
of  $\tilde{W}_{ij}$ . As a result,  $\hat{W}_{ij}$ is also bounded. Given
that the regressor vector  $S_{ij}(\chi_{i})$ is bounded according to [18],
the feedback control law $\tau_{i}$ from (13) is bounded as well. This ensures
boundedness for all signals in the closed-loop system, thereby proving the
first part.

(ii) For the second part of the proof, we examine the Lyapunov function for
the dynamics of $r_{i}$ given by (15) as,

$\small V_{ri}=\frac{1}{2}r_{i}^{T}M_{i}(x_{i1})r_{i}.$ (26)

Temporal differentiation of (26) and substituting in (15) for
$M_{i}(\theta)\dot{r}_{i}$ yields,

$\displaystyle\dot{V}_{ri}=-r_{i}^{T}K_{i}r_{i}+r_{i}^{T}\tilde{W}_{i}^{T}S_{i}(\chi_{i})-\epsilon_{i}^{T}r_{i}.\quad\forall
i\in I[1,N]$ (27)

Using an approach similar to the one used for inequality (23), and by setting
$K_{i}=K_{i1}+2K_{i2}$ , we can demonstrate:

$\displaystyle-r_{i}^{T}K_{i2}r_{i}-r_{i}^{T}\epsilon_{i}$
$\displaystyle\leq\frac{\epsilon_{i}^{2}}{4\lambda_{\min}(K_{i2})}$ (28)
$\displaystyle\leq\frac{\epsilon_{i}^{*2}}{4\lambda_{\min}(K_{i2})},\quad\forall
i\in I[1,N]$ $\displaystyle-r_{i}^{T}K_{i2}r_{i}$
$\displaystyle+r_{i}\tilde{W}_{i}^{T}S_{i}(\chi_{i})$ (29)
$\displaystyle\leq\frac{\tilde{W}_{i}^{*2}s_{i}^{*2}}{4\lambda_{\min}(K_{i2})},\quad\forall
i\in I[1,N].$

In this case, $s_{i}^{*}$ is a fixed positive number so that
$||S_{i}(\chi_{i})||\leq s_{i}^{*}$ for all  $i\in I[1,N]$ . The existence of
such a $s_{i}^{*}$ is confirmed by [18]. Therefore, we arrive at:

$\displaystyle\dot{V}_{ri}$ $\displaystyle\leq-
r_{i}^{T}\lambda_{\min}(K_{i1})r_{i}+\delta_{i}$ (30)
$\displaystyle\leq-\rho_{i}V_{ri}+\delta_{i},\quad\forall i\in I[1,N]$

Where

$\displaystyle\rho_{i}=\min\\{2\lambda_{\min}(K_{i1}),\frac{2\lambda_{\min}(K_{i1})}{\lambda_{\max}(M_{i})}\\}$
$\displaystyle\delta_{i}=\frac{\delta_{i}^{*2}}{4\lambda_{\min}(K_{i2})}+\frac{\tilde{W}_{i}^{*2}s_{i}^{*2}}{4\lambda_{\min}(K_{i2})}.$

Solving 30 gives us:

$\displaystyle 0\leq V_{ri}(t)$ $\displaystyle\leq V_{ri}(0)\exp(-\rho_{i}t)$
(31) $\displaystyle+\frac{\delta_{i}}{\rho_{i}},\quad\forall t\geq 0,\forall
i\in I[1,N]$

from the above solutions, it becomes clear that there exists a finite time
$T_{i}>0$ such that $r_{i}$ will exponentially converge to a small area near
zero. This further confirms that the tracking errors $e_{i}$ will also reach a
small zone around zero. The size of this zone can be minimized by carefully
choosing the design parameters $K_{i1}$ with  $\lambda_{\min}(K_{i1})>0$ .
This completes the proof for the second part of the Theorem.

(iii): From the above proof of the second part, it has been established that
for all  $i\in I[1,N]$ , there exists a finite time $T_{i}>0$ such that for
all  $t\geq T_{i}$ , the tracking errors  $e_{i}$ as well as  $r_{i}$ will
converge to a vicinity near zero. Because  $e_{i}=x_{i1}-\hat{x}_{i1}$ ,
$\hat{x}_{i1}$ will eventually converge to  $x_{01}$ as per Theorem 1, and
$x_{01}$ is a periodic signal under Assumption 1,  $x_{i1}$ will also become a
periodic signal after a finite time  $T_{i}$ . Furthermore, while  $e_{i}$
converges to zero,  $\dot{e_{i}}$ converges to zero, leading  $x_{i2}$
converges to  ${x}_{02}$ . Since the leader dynamics is a smooth continuous
LTI system, periodicity of  $x_{01}$ implies that  $\dot{x}_{01}$ is also
periodic and thus  $x_{i2}$ is periodic after a finite time  $T_{i}$ .
Consequently, the inputs of RBF NNs  $(\chi_{i})$ are made as periodic signals
for all  $t\geq T_{i}$ . According to Lemma 1, the partial PE condition of the
localized RBF NN regression subvector  $S_{i\phi}(\chi_{i})$ along the system
trajectory  $\phi_{i}(\chi_{i}(t)|_{t\geq T_{i}}$ is guaranteed.

$\displaystyle\hat{W}_{i}^{T}S_{i}(\chi_{i})$
$\displaystyle=\hat{W}_{i\phi}^{T}S_{i\phi}(\chi_{i})$
$\displaystyle+\hat{W}_{i\phi^{\bar{}}}^{T}S_{i\phi^{\bar{}}}(\chi_{i}),\quad\forall
i\in I[1,N]$

Accordingly, system dynamics of (19) can be expressed as:

$\displaystyle\dot{r}_{i}$
$\displaystyle=M_{i}^{-1}(x_{i1})(\widetilde{W}_{i\phi}^{T}S_{i\phi}(\chi_{i})-K_{i}r_{i}-C_{i}(\chi_{i})r_{i}-\epsilon^{\prime}_{i\phi}),$
$\displaystyle\forall i\in I[1,N]$

where
$\epsilon^{\prime}_{i\phi}=\epsilon_{i}-\hat{W}_{i\phi^{\bar{}}}^{T}S_{i\phi^{\bar{}}}(\chi_{i})$
is the localized ideal NN approximation error along the tracking trajectory.
Thus, the overall closed-loop adaptive learning system can be described by:

$\displaystyle\scalebox{0.63}{$\begin{pmatrix}\dot{r}_{i}\\\
\dot{\widetilde{W}}_{i\phi,1}\\\ \dot{\widetilde{W}}_{i\phi,2}\\\ \vdots\\\
\dot{\widetilde{W}}_{i\phi,n}\end{pmatrix}=\begin{pmatrix}-M_{i}^{-1}(x_{i1})K_{i}&M_{i}^{-1}(x_{i1})\begin{bmatrix}S_{i\phi,1}^{T}(\chi_{i})&0&0&0\\\
0&S_{i\phi,2}^{T}(\chi_{i})&0&0\\\ \vdots&\vdots&\ddots&\vdots\\\
0&0&0&S_{i\phi,n}^{T}(\chi_{i})\end{bmatrix}\\\
-\Gamma_{i}S_{i\phi,1}(\chi_{i})\\\ -\Gamma_{i}S_{i\phi,2}(\chi_{i})\\\
\vdots&0\\\ -\Gamma_{i}S_{i\phi,n}(\chi_{i})\end{pmatrix}$}$
$\displaystyle\scalebox{0.63}{$\times\begin{pmatrix}r_{i}\\\
\widetilde{W}_{i\phi,1}\\\ \widetilde{W}_{i\phi,2}\\\ \vdots\\\
\widetilde{W}_{i\phi,n}\end{pmatrix}+\begin{pmatrix}-M_{i}^{-1}(x_{i1})\epsilon_{i\phi}\\\
-\Gamma_{i}\sigma_{i}\widehat{W}_{i\phi,1}\\\
-\Gamma_{i}\sigma_{i}\widehat{W}_{i\phi,2}\\\ \vdots\\\
-\Gamma_{i}\sigma_{i}\widehat{W}_{i\phi,n}\end{pmatrix}$},\forall i\in I[1,N]$
(32)

$\displaystyle\scalebox{0.63}{$\begin{pmatrix}\dot{\widetilde{W}}_{i\phi\bar{},1}\\\
\dot{\widetilde{W}}_{i\phi\bar{},2}\\\ \vdots\\\
\dot{\widetilde{W}}_{i\phi\bar{},n}\end{pmatrix}=\begin{pmatrix}-\Gamma_{i}(S_{i\phi\bar{},1}(\chi_{i})r_{i}+\sigma_{i}\hat{W}_{i\phi\bar{},1})\\\
-\Gamma_{i}(S_{i\phi\bar{},2}(\chi_{i})r_{i}+\sigma_{i}\hat{W}_{i\phi\bar{},2})\\\
\vdots\\\
-\Gamma_{i}(S_{i\phi\bar{},n}(\chi_{i})r_{i}+\sigma_{i}\hat{W}_{i\phi\bar{},n})\end{pmatrix}$},\forall
i\in I[1,N]$ (33)

Based on [18], the local approximation error $\epsilon^{\prime}_{i\phi}$ and
$\hat{W}_{i\phi^{\bar{}}}^{T}S_{i\phi^{\bar{}}}(\chi_{i})$ are both small.
Moreover, $\epsilon^{\prime}_{i\phi}$ is proportional to $\epsilon_{i}$.
Studies [18] [19] have thoroughly examined the stability and convergence of
the closed-loop system given by (27). Specifically, it’s established that the
PE condition of $S_{i\phi}(\chi_{i})$ leads to exponential convergence of
$(r_{i},\tilde{W}_{i\phi})$ to zero for all  $i\in I[1,N]$ .

Further, since  $\epsilon_{i\phi}$ is proportional to  $\epsilon_{i}$ , and
$\sigma_{i}$ can be made as small as desired,  $\tilde{W}_{i\phi}$ will also
converge to an arbitrarily small vicinity of zero. This size is dictated by
$\epsilon_{i}$ and  $\sigma_{i}$ .

$\displaystyle H_{i}(\chi_{i})$
$\displaystyle=\hat{W}_{i\phi}^{T}S_{i\phi}(\chi_{i})+\epsilon_{i\phi,1}$ (34)
$\displaystyle=\bar{W}_{i\phi}^{T}S_{i\phi}(\chi_{i})+\epsilon_{i\phi,2},\quad\forall
i\in I[1,N]$

Here, $\epsilon_{i\phi,1}$ and $\epsilon_{i\phi,2}$ are approximation errors.
They are proportional to $\epsilon_{i\phi}$ due to the proven convergence of
$\tilde{W}_{i\phi}$ to zero. As for neurons far from the trajectory,
$S_{i\phi^{\bar{}}}(\chi_{i})$ is minimal, impacting the neural weight
adaptation only slightly. Therefore, the full RBF NN can still accurately
approximate the unknown function $H_{i}(\chi_{i})$ along the trajectory for
$t\geq T_{i}$.

$\displaystyle Hi(\chi_{i})$
$\displaystyle=\hat{W}_{i}^{T}S_{i}(\chi_{i})+\epsilon_{i1}$ (35)
$\displaystyle=\bar{W}_{i}^{T}S_{i}(\chi_{i})+\epsilon_{i2},\quad\forall i\in
I[1,N]$

The approximation errors $\epsilon_{i1}$ and $\epsilon_{i2}$ are proportional
to $\epsilon_{i\phi,1}$ and $\epsilon_{i\phi,2}$ respectively. This concludes
the proof.

######  Remark 2.

From (13) and (14), we observe that the second-layer deterministic learning
control is fully decentralized. Each local control law functions
independently, without the need for sharing information with neighboring
agents. This is in contrast to the first-layer’s cooperative estimation.

######  Theorem 3.

Given a multi-robot manipulator system as described by systems (2) and virtual
leader dynamics (3), along with the network communication topology $G$. Under
Assumption 1 and 2, Problem 1 can be addressed using a distributed control
protocol comprised of (13) (14) and (4) (5). All design parameters should
satisfy the requirements in Theorem 1 and 2.

######  Remark 3.

The proposed two-layer distributed learning control scheme consists of a
cooperative estimation law, given by (4) and (5), and a decentralized
deterministic learning control law, detailed in (13) and (14) for the second
layer. More specifically, the first layer requires inter-agent communication
for sharing estimated information  $(\hat{\chi}_{i},\hat{A}_{i})$ . The second
layer only demands local robot plant states  $\chi_{i}$ and local estimated
information  $\hat{\chi}_{i}$ for feedback.

######  Remark 4.

Despite extensive research in distributed control for multi-robot systems
(e.g., [19]), most existing techniques focus solely on control, largely
overlooking the learning aspects. The distributed learning control scheme
proposed here fills this gap. It provides a capability currently missing in
multi-agent/multi-robot literature: achieving both synchronization and
accurate identification/learning of completely nonlinear uncertain dynamics
for a team of heterogeneous robot manipulators.

##  IV SIMULATION

The multiple 2-DOF robot manipulator system, as described in (2), is
considered with parameters:

$\displaystyle\begin{aligned} M_{i}(q_{i})&=\begin{bmatrix}M_{i11}&M_{i12}\\\
M_{i21}&M_{i22}\end{bmatrix},&C_{i}(q_{i},\dot{q}_{i})&=\begin{bmatrix}C_{i11}&C_{i12}\\\
C_{i21}&C_{i22}\end{bmatrix},\\\ F_{i}(\dot{q}_{i})&=\begin{bmatrix}F_{i11}\\\
F_{i21}\end{bmatrix},&g_{i}(q_{i})&=\begin{bmatrix}g_{i11}\\\
g_{i21}\end{bmatrix}.\end{aligned}$

with

$\displaystyle\begin{aligned}
M_{i11}&=m_{i1}l_{ic1}^{2}+m_{i2}(l_{i1}^{2}+l_{ic2}^{2}+2l_{i1}l_{ic2}\cos(qi2))\\\
&\quad+I_{i1}+I_{i2},\\\
M_{i12}&=m_{i2}(l_{ic2}^{2}+l_{i1}l_{ic2}\cos(q_{i2}))+I_{i2},\\\
M_{i21}&=m_{i2}(l_{ic2}+l_{i1}l_{ic2}\cos(qi2))+I_{i2},\\\
M_{i22}&=m_{i2}l_{ic2}^{2}+I_{i2},\\\
C_{i11}&=-m_{i2}l_{i1}l_{ic2}\dot{q}_{i2}\sin(q_{i2}),\\\
C_{i12}&=-m_{i2}l_{i1}l_{ic2}(\dot{q}_{i1}+\dot{q}_{i2})\sin(q_{i2}),\\\
C_{i21}&=m_{i2}l_{i1}l_{ic2}\dot{q}_{i1}\sin(q_{i2}),\\\ C_{i22}&=0,\\\
g_{i11}&=(m_{i1}l_{ic2}+m_{i2}l_{i1})g\cos(q_{i1})\\\
&\quad+m_{i2}l_{ic2}g\cos(q_{i1}+q_{i2}),\\\
G_{i22}&=m_{i2}l_{ic2}g\cos(q_{i1}+q_{i2}).\end{aligned}$ (36)

Figure 1: Schematic for the inverted double pendulum system..

Figure 2: Network Topology with Agent 0 Serving as the Virtual Leader

In Fig. 1, we identify the core components and their specifications for the
$i$-th robotic agent.  $l_{ic1},l_{ic2}$ represent half of these lengths, and
$F_{i11},F_{i21}$ are defined as constants. The inertia of the links is given
by  $I_{i1}$ and  $I_{i2}$ , and the detailed values are provided in Table I.
By employing  $N=5$ manipulators, the robotic agents are set to follow the
reference trajectory provided by the virtual leader. A leader dynamic is
formulated to generate a periodic signal for the synchronization control:

$\begin{bmatrix}\dot{x}_{01}\\\
\dot{x}_{02}\end{bmatrix},=\begin{bmatrix}0&0&1&0\\\ 0&0&0&1\\\ -1&0&0&0\\\
0&-1&0&0\end{bmatrix}\begin{bmatrix}x_{01}\\\
x_{02}\end{bmatrix},\begin{bmatrix}x_{01}(0)\\\
x_{02}(0)\end{bmatrix}=\begin{bmatrix}0\\\ 0.8\\\ 0.8\\\ 0\end{bmatrix}$

The network topology $G$ amongst the five robots is depicted in Fig 2. This
topology includes a spanning tree with the virtual leader, referred to as
agent 0, positioned at the center. The system states encompass both measured
signals and reference signals, totaling eight signals. Originally, if we use a
neural network with 4 nodes for each signal, we end up with a massive number
of nodes $4^{8}=65,536$ which is computationally expensive to train. However,
the reference signals are expressed in terms of the measured signals, as shown
in (8), which accentuates the significance of the measured signals.
Consequently, we construct the Gaussian RBF NN using the dominant four
dimensions of the system state. This strategy reduces the size of the NN to
$4^{4}=256$ which not only saves computational resources but also allows us to
build a more precise model by dedicating more nodes to each of these four key
dimensions. The range of each dimension lies within $[-1.2,1.2]$, with a width
parameter  $\gamma_{i}=0.8$ . The observer and controller parameters are set
to  $\beta_{1}=\beta_{2}=1$, $\gamma_{i}=10$ ,  $K=10$ , and
$\sigma_{i}=0.001$ for all  $i\in[1,5]$ . The initial conditions are defined
as  $x_{11}(0)=[0.2\ 0.1]^{T}$ ,  $x_{21}(0)=[0.3\ 0.5]^{T}$ ,
$x_{31}(0)=[0.4\ 0.1]^{T}$, $x_{41}(0)=[0.2\ 0.5]^{T}$ ,  $x_{51}(0)=[0.4\
0.1]^{T}$ , and  $x_{i2}=[0\ 0]^{T}$, $\forall i\in[1,5]$ . Initial conditions
for all distributed observer states  $\hat{\chi}_{i}$ and NN weights
$\hat{W}_{i}$ are uniformly initialized to zero. We conducted an assessment of
two key components within our system setup: the distributed cooperative
estimator, described by (4) and (5), and the decentralized deterministic
learning control law, encapsulated by (13) and (14). Despite varying
complexities and nonlinear uncertainties among multiple robotic agents, the
results indicate satisfactory tracking performance. This is further confirmed
by Fig. 3, which shows a rapid convergence of tracking errors to zero, as
discussed in Theorem 2.

TABLE I: Parameters of the robot. Parameter | Robot number
---|---
|  | 1 | 2 | 3 | 4 | 5
$m1$ (kg) |  | 2 | 2.2 | 2.3 | 1.9 | 2.4
$m2$ (kg) |  | 0.85 | 0.9 | 1 | 0.9 | 1.5
$l1$ (m) |  | 0.35 | 0.5 | 0.6 | 0.52 | 0.57
$l2$ (m) |  | 0.31 | 0.4 | 0.5 | 0.48 | 0.53
$I1\times 10^{-3}$ (kg m2) |  | 61.25 | 70 | 72.14 | 67.21 | 73.42
$I2\times 10^{-3}$ (kg m2) |  | 20.42 | 25.21 | 27.1 | 25.4 | 22.63

Figure 3: Tracking error for each agent.

Figure 4: Joint angle tracking response for each agent

The tracking response of joint angles from the second layer in comparison to
the desired signals originating from the first layer and the leader’s signal
is represented in Fig. 4. It is observed that all signals rapidly converge to
align with the leader’s signal within the initial seconds. Fig. 5 shows how
quickly the Neural Network weights settle for all robotic agents. This fast
convergence is in accordance with Theorem 2.

Figure 5: NN weights Convergence for all agents.

Moreover, we graphed the NN approximation associated with the unknown dynamic
variables  $H_{i}(\chi_{i})$ for each agent in Fig. 6, employing RBF NN
approximation represented as  $\hat{W}_{i}^{T}S_{i}(\chi_{i})$ . In our
controller, all system parameters are treated as having unstructured
uncertainty, indicating that the control algorithm is environment-independent.
Whether the system operates underwater, where buoyancy forces vary with the
depth of a robot’s arm, or in space, where the inertia matrix is unknown, the
algorithm adapts accordingly. Regardless of these conditions, after completing
the training process, the stored weights for dynamic parameter identification
can be utilized in the controller, even after the system has been turned off.

Figure 6: NN approximation

##  V CONCLUSION

In this study, The challenge of implementing composite synchronization and
adaptive learning control for a network of multi-robot manipulators, with
complete nonlinear uncertain dynamics was successfully met. The novel approach
contains a distributed cooperative estimator in the first layer for estimating
the virtual leader’s states, and a decentralized deterministic learning
controller in the second layer to track the leader’s states while identifying
each robot’s distinct nonlinear uncertain dynamics. The key strengths of This
method include simultaneously achieving synchronization control and learning
robots’ complete nonlinear uncertain dynamics in a decentralized manner that
eliminates the need for global connectivity. This improves robotics control in
space and underwater settings, where system dynamics are typically uncertain.
The method’s effectiveness is verified through detailed mathematical
validation and simulations. Stored identified nonlinearities can be used in
the controller once, eliminating frequent adjustments when the system
restarts. We aim to extend this framework to other robots with Euler-
Lagrangian dynamic systems.

## References

  * [1] R. Cui and W. Yan, “Mutual synchronization of multiple robot manipulators with unknown dynamics,” _Journal of Intelligent & Robotic Systems_, vol. 68, pp. 105–119, 2012.
  * [2] N. Seenu, K. C. RM, M. Ramya, and M. N. Janardhanan, “Review on state-of-the-art dynamic task allocation strategies for multiple-robot systems,” _Industrial Robot: the international journal of robotics research and application_ , vol. 47, no. 6, pp. 929–942, 2020.
  * [3] E. Jandaghi, X. Chen, and C. Yuan, “Motion dynamics modeling and fault detection of a soft trunk robot,” in _2023 IEEE/ASME International Conference on Advanced Intelligent Mechatronics (AIM)_. IEEE, 2023, pp. 1324–1329.
  * [4] C. Yuan, S. Licht, and H. He, “Formation learning control of multiple autonomous underwater vehicles with heterogeneous nonlinear uncertain dynamics,” _IEEE transactions on cybernetics_ , vol. 48, no. 10, pp. 2920–2934, 2017.
  * [5] N. Hazon and G. A. Kaminka, “On redundancy, efficiency, and robustness in coverage for multiple robots,” _Robotics and Autonomous Systems_ , vol. 56, no. 12, pp. 1102–1114, 2008.
  * [6] J. Huang, C. Wen, W. Wang, and Y.-D. Song, “Adaptive finite-time consensus control of a group of uncertain nonlinear mechanical systems,” _Automatica_ , vol. 51, pp. 292–301, 2015.
  * [7] M. Abdelatti, C. Yuan, W. Zeng, and C. Wang, “Cooperative deterministic learning control for a group of homogeneous nonlinear uncertain robot manipulators,” _Science China Information Sciences_ , vol. 61, pp. 1–19, 2018.
  * [8] Q. Liu, D. Li, S. S. Ge, R. Ji, Z. Ouyang, and K. P. Tee, “Adaptive bias rbf neural network control for a robotic manipulator,” _Neurocomputing_ , vol. 447, pp. 213–223, 2021.
  * [9] A. Rodriguez-Angeles and H. Nijmeijer, “Mutual synchronization of robots via estimated state feedback: a cooperative approach,” _IEEE Transactions on control systems technology_ , vol. 12, no. 4, pp. 542–554, 2004.
  * [10] S.-J. Chung and J.-J. E. Slotine, “Cooperative robot control and concurrent synchronization of lagrangian systems,” _IEEE transactions on Robotics_ , vol. 25, no. 3, pp. 686–700, 2009.
  * [11] H. Wang, “Flocking of networked uncertain euler–lagrange systems on directed graphs,” _Automatica_ , vol. 49, no. 9, pp. 2774–2779, 2013.
  * [12] X. Dong, C. Yuan, P. Stegagno, W. Zeng, and C. Wang, “Composite cooperative synchronization and decentralized learning of multi-robot manipulators with heterogeneous nonlinear uncertain dynamics,” _Journal of the Franklin Institute_ , vol. 356, no. 10, pp. 5049–5072, 2019.
  * [13] D. Liang and J. Huang, “Leader-following bipartite consensus of multiple uncertain euler-lagrange systems over signed switching digraphs,” _Neurocomputing_ , vol. 405, pp. 96–102, 2020.
  * [14] W. Ren and R. W. Beard, “Consensus seeking in multiagent systems under dynamically changing interaction topologies,” _IEEE Transactions on automatic control_ , vol. 50, no. 5, pp. 655–661, 2005.
  * [15] J. Park and I. W. Sandberg, “Universal approximation using radial-basis-function networks,” _Neural computation_ , vol. 3, no. 2, pp. 246–257, 1991.
  * [16] ——, “Approximation and radial-basis-function networks,” _Neural computation_ , vol. 5, no. 2, pp. 305–316, 1993.
  * [17] C. Wang and D. J. Hill, _Deterministic learning theory for identification, recognition, and control_. CRC Press, 2018.
  * [18] ——, _Deterministic learning theory for identification, recognition, and control_. CRC Press, 2009.
  * [19] C. Yuan and C. Wang, “Persistency of excitation and performance of deterministic learning,” _Systems & control letters_, vol. 60, no. 12, pp. 952–959, 2011.


# Probing non-unitarity of neutrino mixing in the scenario of Lorentz
violation and dark nonstandard interaction

Trisha Sarkar<EMAIL_ADDRESS>Indian Institute of Technology Jodhpur,
Jodhpur 342037, India

###### Abstract

Neutrino flavour oscillation is one of the primary indication of the existence
of physics beyond standard model. The presence of small neutrino mass is
indispensable to explicate the oscillation among different flavours of
neutrino. By the addition of a right handed neutral lepton with the standard
model fermions, it is possible to generate tiny neutrino mass. Such additional
fermion may induce non-unitarity to the $3\times 3$ PMNS mixing matrix which
influences the propagation of neutrino in space-time. In this work the effect
of both unitary and non-unitary mixing matrix is analyzed in neutrino
oscillation in presence of two new physics scenarios, Lorentz violation and
dark non-standard interaction. Lorentz symmetry violation mainly appears at
the Planck scale, which may also be manifested at a lower energy level. On the
other hand, dark non standard interaction arises due to the interaction of
neutrino with the environmental dark matter which contributes as a
perturbative correction to the neutrino mass. The analysis is carried out in
the context of long baseline DUNE and short baseline Daya Bay experimental set
up. The signature of dark nonstandard interaction is observable in both DUNE
and Daya Bay set up in terms of large value of neutrino survival and
oscillation probability respectively and is a possible explanation for the
excess flux observed at $\sim 5$ MeV in Daya Bay experiment. The signature of
Lorentz violation is also possible to be observed in the Daya Bay experiment.

##  I Introduction

Neutrino oscillation is a well established phenomenon at present in the domain
of particle physics. The presence of small neutrino mass is essentially
required for the neutrino flavour oscillation to occur which is one of the
most primary necessity to look for the physics beyond standard model (SM).
While in SM neutrinos are considered to be left handed (LH) massless Weyl
fermions, in several extensions of SM neutrino mass may appear due to some
radiative correction Zee:1980ai ; Frigerio:2002in ; Ahriche:2015loa  or by the
inclusion of right handed (RH) chiral neutrinos Kersten:2007vk ; Tang:2014hna
; King:2003jb . In generally, the non zero neutrino mass is expected to appear
due to breaking of some hidden symmetry MohseniSadjadi:2017jne ;
Lindner:2014oea . The smallness of neutrino mass can be explained from type-I
Seesaw mechanism which involves the addition of heavy neutral leptons
Branco:2020yvs .

Due to the presence of these additional leptons the unitary neutrino mixing
matrix, which is a subset of a general lepton mixing matrix, acquires non-
unitary properties. Such non-unitarity condition affects the neutrino
oscillation and its propagation through the medium. Non unitraity (NU) of
mixing matrix also gives rise to lepton flavour violation (LFV) which is
another signature of beyond standard model (BSM) phenomena Forero:2011pc . The
Seesaw scale may exist at any point within electroweak (EW) and Planck scale.
Although high scale Seesaw mechanism are unable to provide momentous deviation
from unitarity, several low scale mechanisms are capable of doing so
Malinsky:2009df ; Forero:2021azc . The effect of non-unitary mixing matrix has
been discussed in several previous works in context of different aspects of
neutrino oscillation Berryman:2014yoa ; Cai:2017mow ; Miranda:2019ynh . In
addition to the feature of non-unitarity of the mixing matrix, two of the
other most profound BSM scenarios in the neutrino sector are the violation of
Lorentz symmetry and the presence of nonstandard interaction (NSI).

Lorentz symmetry is conserved in two of the foremost theories of nature, SM of
particle physics and general theory of relativity (GR). Various models have
been developed allowing the spontaneous violation of Lorentz symmetry at the
Planck scale ($\sim 10^{19}$ GeV) which unifies SM with gravity. Although it
is impossible to test any theory at this energy scale, it is expected that
such new physics (NP) signals are manifested at lower energy scale which can
be probed with current experimental facilities. Mathematically, SM is
portrayed as the renormalizable effective field theory (EFT) at lower energy
scale consisted of the operators with dimension $D\leq 4$, while the NP
signals are represented in terms of non-renormalizable higher dimensional
operators ($D>4$) suppressed by the scale of NP. The technique of
interferometry has been widely utilized to investigate the Lorentz symmetry
Michelson:1887zz ; PhysRev.42.400 ; Kobakhidze:2007iz ; 2017arXiv170709016P .
Neutrino oscillation provides a natural source of interferometer
CaboMontesdeOca:2020xeq  which is an excellent tool to search for Lorentz
violation (LV). Currently ongoing oscillation experiments are able to search
for the violation of such a fundamental symmetry. Several experimental
facilities have constrained the strength of LV Barenboim:2018ctx ; Super-
Kamiokande:2014exs ; Glashow:2004im .

Non standard interaction is another elegant way to incorporate NP in the
neutrino sector. Introduced by Wolfenstein in 1978 in terms of flavour
changing neutral current (FCNC) when the neutrino undergoes the coherent
forward elastic scattering in the matter PhysRevD.17.2369 , NSI has been
widely studied in many literatures Farzan:2017xzy ; Proceedings:2019qno ;
Ohlsson:2012kf ; Antusch:2008tz . It should be mentioned that the discussions
in these works involve neutrino interaction mediated by a vector boson, which
directly contributes to the matter potential, while the neutrinos may also
interact via a scalar mediator, giving rise to scalar NSI Ge:2018uhz ;
Medhi:2021wxj . Since the universe contains a large abundance of dark matter
(DM), it is also possible that the scalar particle may belong to the dark
sector, generating dark NSI article . Unlike the case of vector NSI, which is
possible to be represented in terms of dimension-$6$ four fermion operator and
introduces modification to the matter potential, the neutrino-scalar
interaction induces correction to the neutrino mass term. The bounds on the
NSI parameters corresponding to the scalar mediator is much weaker, as
compared to the case of vector NSI. Borexino collaboration has been able to
put constraints on the scalar NSI parameters Borexino:2017rsf ; Khan:2019jvr .

In this work, the comparative analysis of the unitary and non-unitary mixing
matrix are presented, in the sector of Lorentz violation (LV) and dark NSI
within the framework of neutrino oscillation. The results are determined in
terms of neutrino survival and oscillation probability in the context of long
baseline (LBL) accelerator experiment, DUNE and short baseline (SBL) reactor
experiment, Daya Bay. The signature of dark NSI with the consideration of NU
scenario is observable in DUNE experiment within the energy range $\sim
2.5-2.8$ GeV at which the maximum flux is obtained. It is also identified that
the same combined scenario of dark NSI with NU mixing could be a possible
explanation for the excess flux observed at $\sim 5$ MeV in Daya Bay, known as
the ”$5$ MeV bump” Jaffe:2021wgf ; Dentler:2017tkw ; Berryman:2018jxt .

The plan of this work is as follows. In section II the formalism of this work
is described including the different aspects of NP scenarios. In section III,
the results of this analysis are discussed and finally the conclusion is
mention in section IV.

##  II Formalism

In this section the formalism of this work is presented. In section II.1, the
facet of non-unitarity of neutrino mixing matrix is discussed, while the
scenarios of LV and dark NSI are described in section II.2 and II.3
respectively. In section II.4, the analytical framework to determine the
neutrino oscillation probability is presented.

###  II.1 Non unitarity of mixing matrix

The origin of non-unitary neutrino mixing lies in the existence of non zero
neutrino mass. Recent measurement by KATRIN collaboration estimates the
neutrino mass to be $m_{\nu}<0.8$ eV KATRIN:2022ayy . Neutrino mass can be
either of Dirac or Majorana in nature. However, it is not favourable for the
neutrino to acquire pure Dirac mass term, as the corresponding Yukawa
couplings are extremely small $\sim 10^{-12}$ for the generation of neutrino
mass $m_{\nu}\sim 0.1$ eV deGouvea:2016qpx  which is much smaller as compared
to the other SM fermions. This also demands the lepton number to be a
fundamental symmetry of SM, which originally appears accidentally in SM. Such
problem is possible to overcome if the neutrinos are considered to be of
Majorana nature which beings out the concept of the well known Seesaw
mechanism Lindner:2001hr . This can be achieved by the addition of non-
renormalizable and lepton number violating dimension-$5$ Weinberg operator
with the SM effective Lagrangian, which generates tiny Majorana neutrino mass
after electroweak symmetry breaking (EWSB). If the mass of the RH neutrino is
much larger than the vacuum expectation value of the SM Higgs boson ($\sim
246$ GeV), it is possible to generate tiny neutrino mass, without undergoing
the problem of small Yukawa coupling Kanemura:2011jj .

However, the addition of the RH neutrinos affects the flavour oscillation and
the interaction nature of the neutrinos. While incorporating the RH neutrino
with its LH counterpart, the contribution of RH neutrino to the weak CC
reaction gives rise to NU property of the $3\times 3$ PMNS mixing matrix. In
the presence of $n$ number of neutrinos, the generic mixing matrix takes the
following form

$U^{n\times n}=\begin{pmatrix}T^{3\times 3}&Q^{3\times(n-3)}\\\ V^{(n-3)\times
3}&U^{(n-3)\times(n-3)},\end{pmatrix}$ (1)

where $T^{3\times 3}$ and $U^{(n-3)\times(n-3)}$ represent the mixing in the
sector of LH and RH neutrinos respectively, provided there are three active LH
neutrinos and $(n-3)$ number of RH neutrinos. $Q^{3\times(n-3)}$ and
$V^{(n-3)\times 3}$ contain the coupling parameters in the RH sector. In case
of unitary mixing, $T^{3\times 3}$ is represented as the PMNS mixing matrix,
given by

$\displaystyle U_{PMNS}=$
$\displaystyle\begin{pmatrix}c_{12}c_{13}&s_{12}c_{13}&s_{13}e^{-i\delta}\\\
-s_{12}c_{23}-c_{12}s_{13}s_{23}e^{i\delta}&c_{12}c_{23}-s_{12}s_{13}s_{23}e^{i\delta}&c_{13}s_{23}\\\
s_{12}s_{23}-c_{12}s_{13}c_{23}e^{i\delta}&-c_{12}s_{23}-s_{12}s_{13}c_{23}e^{i\delta}&c_{13}c_{23}\end{pmatrix}.$
(2)

Here $c_{ij}=cos\theta_{ij}$, $s_{ij}=sin\theta_{ij}$ and $\delta$ is the CP
violating phase. In presence of NU, $T^{3\times 3}$ can be parametrized in a
model independent way, assuming that the additional neutral heavy leptons do
not contribute to neutrino oscillation Forero:2021azc . One of the most
convenient parametrization of $T^{3\times 3}$ is expressed as follows
Fernandez-Martinez:2007iaa ; Escrihuela:2015wra

$T^{3\times 3}=\begin{pmatrix}\alpha_{00}&0&0\\\ \alpha_{10}&\alpha_{11}&0\\\
\alpha_{20}&\alpha_{21}&\alpha_{22}\end{pmatrix}U_{PMNS},$ (3)

where $\alpha_{ij}$’s are the parameters that invoke NU and $U_{PMNS}$ is
given by eqn. (II.1). In generally, the diagonal parameters are real , while
the off-diagonal ones are complex
($\alpha_{ij}=|\alpha_{ij}|e^{\phi_{ij}},\leavevmode\nobreak\ i\neq j$) that
can provide an additional source of CP violation. The constraints on NU
parameters in eqn. (3) are obtained from several decay processes
Antusch:2006vwa ; Fernandez-Martinez:2007iaa ; Escrihuela:2019mot  and
oscillation experiments Escrihuela:2015wra ; Miranda:2019ynh .

###  II.2 Violation of Lorentz symmetry

It is widely known that SM is invariant under Lorentz transformation, which
leads to the conservation of CPT symmetry in the universe. It is assumed in
several extensions of SM that the Lorentz symmetry is possible to be violated
at some higher energy scale, possibly at Planck scale. The origin of Lorentz
violation (LV) is discussed extensively in various theoretical models of
string theory Mavromatos:2007xe  and quantum gravity Cuadros-Melgar:2011whf ;
Moffat:2009ks ; Eichhorn:2019ybe . In string theory, the violation of Lorentz
symmetry can occur spontaneously in presence of a non-perturbative vacuum
generating tensor fields acquiring non-zero vacuum expectation value
Bluhm:2013mu ; 2014IJMPS..3060274H . In non-commutating field theory, the
phenomenon of Lorentz violation occurs naturally. In quantum gravity theory,
the Lorentz violation may appear by introducing the concept of anisotropic
scaling of the space and time Horava:2009if .

The observation of Lorentz violation is possible at Planck scale ($E_{Pl}\sim
10^{19}$ GeV) which is not accessible by any current terrestrial experimental
facility. At electroweak (EW) scale ($E_{EW}\sim 100$ GeV) the effect of
Lorentz violation is suppressed by a factor of $E_{EW}/E_{Pl}\sim 10^{-17}$.
However, several experimental techniques have been developed over the years to
test the existence of LV. The system of three flavour neutrino oscillation is
one of the most commendatory sector to study such effects which can also be
understood as interference between different distinguishable states
Cabo:2021exm  . The analysis of LV in the framework of neutrino oscillation
has been done previously in several works Bahcall:2002ia ; Antonelli:2018fbv ;
Torri:2020dec ; Super-Kamiokande:2014exs ; Kostelecky:2004hg ;
KumarAgarwalla:2019gdj ; Sahoo:2021dit ; Arias:2006vgq . The Lagrangian
density corresponding to LV considering the effect of CPT violation is
expressed as Diaz:2015dxa

$\mathcal{L}_{LV}\subset
a_{\alpha\beta}^{\lambda}(\bar{\nu}_{\alpha}\gamma_{\lambda}P_{L}\nu_{\beta}),$
(4)

where $P_{L}=(1-\gamma_{5})/2$ and $a_{\alpha\beta}^{\lambda}$ is the
parameter exhibiting the strength of LV. Here $\lambda$ is the spacetime index
($\lambda=0,1,2,3$) and $\alpha,\beta$ are the flavour indices,
$\alpha,\beta\in\\{e,\mu,\tau\\}$. In this analysis, only the isotropic
component of LV parameters ($\lambda=0$) are taken into account. The
interaction Hamiltonian corresponding to the phenomenon of LV is given by

$\mathcal{H}_{LV}=\begin{pmatrix}a_{ee}&a_{e\mu}&a_{e\tau}\\\ a_{\mu
e}&a_{\mu\mu}&a_{\mu\tau}\\\ a_{\tau e}&a_{\tau\mu}&a_{\tau\tau}\end{pmatrix}$
(5)

Here $a_{\alpha\beta}^{\lambda=0}$ is simply denoted by $a_{\alpha\beta}$. As
$\mathcal{H}_{LV}$ is Hermitian, $a_{\alpha\beta}^{*}=a_{\beta\alpha}$. Due to
the consideration of CPT-violating scenario,
$a^{\lambda}_{\alpha\beta}\xrightarrow[\text{}]{\text{CPT}}-a^{\lambda}_{\alpha\beta}$.
The off-diagonal elements are generally complex,
$a_{\alpha\beta}=|a_{\alpha\beta}|e^{i\phi_{\alpha\beta}}$. The constraints on
these parameters are obtained from various experiments such as IceCube
PhysRevD.82.112003 , MINOS MINOS:2010kat , Double Chooz PhysRevD.86.112009 .

###  II.3 Dark non standard interaction

In neutrino sector another way to introduce the BSM physics is via the
inclusion of the non standard interaction (NSI) which appears in several mass
models. Considering SM as a lower energy approximation of a full theory
existing at a much higher energy scale ($\Lambda$), the effective Lagrangian
is possible to be represented in terms of a number of higher dimensional
Lorentz structures. Each of the operator is non-renormalizable, in general,
and is restricted by the NP energy scale, $\Lambda$ at which the phenomenon of
LV might also occur.

NSI is exchanged either by a vector or a scalar mediator. The vector mediator
is either charged or neutral in nature. The charged current (CC) vector NSI is
tightly constrained, while the constraints on the neutral current (NC) vector
NSI is much weaker. The importance of NC NSI has been discussed in several
literatures Farzan:2017xzy ; Antusch:2008tz ; Proceedings:2019qno . The NC NSI
contributes a subleading effect to the neutrino oscillation and the
Hamiltonian corresponding to it is given by

$\mathcal{H}_{NC-
NSI}^{vec}=A\begin{pmatrix}\epsilon_{ee}&\epsilon_{e\mu}&\epsilon_{e\tau}\\\
\epsilon_{\mu e}&\epsilon_{\mu\mu}&\epsilon_{\mu\tau}\\\ \epsilon_{\tau
e}&\epsilon_{\tau\mu}&\epsilon_{\tau\tau}\end{pmatrix}$ (6)

The parameters, $\epsilon_{\alpha\beta}$ determine the strength of NSI which
corresponds to the dimension-$6$ four fermion lepton number conserving
operator. The off-diagonal elements ($\alpha\neq\beta$) are complex, in
general. $A=\pm G_{F}n_{e}$, with $G_{F}$ and $n_{e}$ representing Fermi
constant and the number density of electron in the medium respectively. The
${}^{\prime}+^{\prime}$ and ${}^{\prime}-^{\prime}$ signs correspond to
neutrino and anti-neutrino respectively. In a recent analysis, the degeneracy
existing between the scenario of LV and vector NSI is discussed along with the
role of atmospheric neutrino to lift up this degeneracy sahoo2022core .

While the vector NSI affects the matter potential during neutrino propagation
through a medium, the NSI mediated by a scalar particle has completely
different consequences. Instead of modifying the SM matter potential, the
scalar NSI directly influences the neutrino mass term. The effect of scalar
mediator is added up to the mass term as a perturbation. According to the
modern cosmological model, the universe is filled with dark matter (DM)
particles ($\sim 27\%$). From several astrophysical constraints the local dark
matter density is estimated to be $\rho_{\chi}\sim 0.47$ GeV/cm3. Since the
number density of the DM particles is inversely proportional to its mass,
there can be a large abundance of DM present in the universe if their mass is
sufficiently small. These light DM particles ($m_{\chi}\leq 100$ eV) can be
either scalar or vector, and they may also interact with neutrinos. Such
interaction generates another kind of NSI, named as dark NSI. If the DM
mediator is a scalar particle ($\phi$) of mass $m_{\phi}$, the neutrino-DM
interaction Lagrangian is given by Ge:2019tdi ; article

$-\mathcal{L}_{NSI}^{dark}\subset\frac{1}{2}M_{\alpha\beta}\bar{\nu}_{\alpha}\nu_{\beta}+y_{\alpha\beta}\bar{\nu}_{\alpha}\phi\nu_{\beta}+\frac{1}{2}m_{\phi}^{2}\phi^{2}$
(7)

Here $y_{\alpha\beta}$ is the Yukawa coupling corresponding to neutrino-DM
interaction and $M_{\alpha\beta}$ is the neutrino mass matrix. The
corresponding interaction Hamiltonian is given by Ge:2019tdi ; article

$\mathcal{H}_{NSI}^{dark}=\frac{1}{2E_{\nu}}(MM^{\dagger}+\delta M^{2})$ (8)

Here the correction to the mass term is represented as, $\delta
M^{2}_{\alpha\beta}=\frac{2\rho_{\phi}}{m_{\phi}^{2}}y_{\alpha
j}y^{*}_{j\beta}$. Such coupling of neutrino with DM is shown to violate the
CPT symmetry Ge:2019tdi . The correction to the mass term due to the presence
of dark NSI is parametrized as follows Ge:2018uhz ; Medhi:2021wxj

$\delta M\equiv\sqrt{\Delta
m_{31}^{2}}\begin{pmatrix}\eta_{ee}&\eta_{e\mu}&\eta_{e\tau}\\\ \eta_{\mu
e}&\eta_{\mu\mu}&\eta_{\mu\tau}\\\ \eta_{\tau
e}&\eta_{\tau\mu}&\eta_{\tau\tau}.\end{pmatrix}$ (9)

Here $\eta_{\alpha\beta}$ determines the strength of dark NSI. From the recent
data obtained from Borexino $2017$, the constraint on only one diagonal NSI
parameter is estimated Ge:2018uhz .

In the next section, the analytical description to obtain neutrino oscillation
probability is discussed in the scenario of LV and dark NSI with the
consideration of unitary and non-unitary neutrino mixing matrix.

###  II.4 Neutrino oscillation in matter

The phenomenon of neutrino oscillation is based on the fact that neutrino
flavour eigen states are represented as a superposition of the mass eigen
states. Considering an initial neutrino beam $\nu_{\alpha}$
($\alpha=e,\mu,\tau$), the relation between the flavour and the mass eigen
states is represented as

$\ket{\nu_{\alpha}(t=0)}=\sum_{i=1}^{3}U^{3\times 3*}_{\alpha
i}\ket{\nu_{i}(t=0)}$ (10)

The matrix $U^{3\times 3}$ might be unitary ($U_{PMNS}$) or non-unitary
($T^{3\times 3}$), represented by eqn. (2) and (3) respectively. Now, the time
evolution of the mass eigen state is expressed as,
$\ket{\nu_{i}(t)}=e^{-iH_{m}t}\ket{\nu_{i}(t=0)}=e^{-iE_{i}t}$, where $H_{m}$
is the interaction Hamiltonian in mass basis and $E_{i}$ is the corresponding
energy eigen value. In mass basis the effective Hamiltonian in presence of
Lorentz violation (LV) is represented as

$\displaystyle H_{m}=\frac{1}{2E_{\nu}}\begin{pmatrix}0&0&0\\\ 0&\Delta
m_{21}^{2}&0\\\ 0&0&\Delta m_{31}^{2}\end{pmatrix}+U^{3\times
3}A\begin{pmatrix}1&0&0\\\ 0&0&0\\\ 0&0&0\end{pmatrix}U^{3\times 3\dagger}$
$\displaystyle+U^{3\times 3}\begin{pmatrix}a_{ee}&a_{e\mu}&a_{e\tau}\\\ a_{\mu
e}&a_{\mu\mu}&a_{\mu\tau}\\\ a_{\tau
e}&a_{\tau\mu}&a_{\tau\tau}\end{pmatrix}U^{3\times 3\dagger}$ (11)

while considering the effect of dark NSI, the effective Hamiltonian is
modified as

$\displaystyle H_{m}=\frac{1}{2E_{\nu}}\begin{pmatrix}0&0&0\\\ 0&\Delta
m_{21}^{2}&0\\\ 0&0&\Delta m_{31}^{2}\end{pmatrix}+U^{3\times
3}A\begin{pmatrix}1&0&0\\\ 0&0&0\\\ 0&0&0\end{pmatrix}U^{3\times 3\dagger}$
$\displaystyle+\frac{1}{2E_{\nu}}(MM^{\dagger}+\delta M^{2})$ (12)

Figure 1: Variation of oscillation probability $P_{\mu e}$ in
($E_{\nu},\delta$) plane in the scenario of Lorentz violation (LV), in the
context of DUNE experiment ($E_{\nu}=1-10$ GeV, $L=1300$ km) for NO,
considering unitary (left panel) and non-unitary (right panel) neutrino
mixing.

Figure 2: Same as fig. 1 in the scenario of dark NSI.

Figure 3: Variation of oscillation probability $P_{\mu e}$ in
($E_{\nu},\delta$) plane in the scenario of Lorentz violation (LV), in the
context of DUNE experiment ($E_{\nu}=1-10$ GeV, $L=1300$ km) for IO,
considering unitary (left panel) and non-unitary (right panel) neutrino
mixing.

Figure 4: Same as fig. 3 in case of dark NSI.

The unperturbed part of the mass matrix $M$ is given by, $M\equiv U^{3\times
3}.Diag(m_{1},m_{2},m_{3}).U^{3\times 3\dagger}$. The time evolution of the
flavour state is expressed as,
$\ket{\nu_{\alpha}(t)}=e^{-iH_{f}t}\ket{\nu_{\alpha}(t=0)}$, $H_{f}$ is the
Hamiltonian in the flavour basis, $H_{f}=U^{3\times 3}H_{m}U^{3\times
3\dagger}$. To calculate the probability of neutrino oscillation, the
framework of ref. Ohlsson:1999xb  is taken into consideration, according to
which the flavour evolution operator at ultra relativistic limit is expressed
as

$\displaystyle
U_{f}(L)=e^{-iH_{f}L}=\phi\sum_{i=1}^{3}e^{-iL\lambda_{a}}\frac{1}{3\lambda_{a}^{2}+C_{1}}\times$
$\displaystyle\left[(\lambda_{a}^{2}+C_{1})I_{3}+\lambda_{a}\tilde{T}+\tilde{T}^{2}\right]$
(13)

Here $T=H_{m}-tr(H_{m})I_{3}/3$, $I_{3}$ being $3\times 3$ identity matrix.
The matrix, $T$ is $3\times 3$ traceless Hermitian matrix with eigenvalues
$\lambda_{a}$. $\phi=exp[-iLtr(H_{m})I_{3}/3]$, $\tilde{T}=U^{3\times
3}TU^{3\times 3\dagger}$ and $C_{1}=Det(T)tr(T^{-1})$.

Parameters | Best fit$\pm\leavevmode\nobreak\ 1\sigma$ (NO) | Best fit$\pm\leavevmode\nobreak\ 1\sigma$ (IO)
---|---|---
$\theta_{12}^{o}$ | $33.82^{+0.78}_{-0.76}$ | $33.82^{+0.78}_{-0.75}$
$\theta_{23}^{o}$ | $49.7^{+0.9}_{-1.1}$ | $49.7^{+0.9}_{-1.0}$
$\theta_{13}^{o}$ | $8.61^{+0.12}_{-0.13}$ | $8.65^{+0.12}_{-0.13}$
$\delta$ | $217^{+40}_{-28}$ | $280^{+25}_{-28}$
$\Delta m_{21}^{2}$ (eV2) | $7.39^{+0.21}_{-0.20}\times 10^{-5}$ | $7.39^{+0.21}_{-0.20}\times 10^{-5}$
$\Delta m_{3l}^{2}$ (eV2) | $+2.525^{+0.033}_{-0.031}\times 10^{-3}$ | $-2.512^{+0.034}_{-0.031}\times 10^{-3}$
Table 1: Values of neutrino oscillation parameters Esteban:2018azc  with
$1\sigma$ interval

Parameters | Best fit with $3\sigma$ interval
---|---
$\alpha_{00}$ | $0.93$
$\alpha_{11}$ | $0.95$
$\alpha_{22}$ | $0.95$
$|\alpha_{10}|$ | $3.6\times 10^{-2}$
$|\alpha_{20}|$ | $0.13$
$|\alpha_{21}|$ | $0.021$
Table 2: Non-unitarity parameters Miranda:2019ynh  with $3\sigma$ interval.
The phase factors are considered to be zero.

##  III Results and Discussion

Figure 5: Variation of $\bar{\nu}_{ee}$ survival probability $\bar{P}_{ee}$ in
($E_{\nu},\delta$) plane in the scenario of Lorentz violation (LV), in the
context of Daya Bay experiment ($E_{\nu}=0.8-6$ MeV, $L=2$ km) for NO,
considering unitary (left panel) and non-unitary (right panel) neutrino
mixing.

Figure 6: Same as fig. 5 in case of dark NSI

In this section the effects of unitary and non unitary mixing matrix in the
neutrino oscillation and survival probability are analysed and compared for
both normal (NO) and inverted mass ordering (IO), considering two of the new
physics scenarios, LV and dark NSI, in the context of DUNE and Daya Bay
experimental set ups respectively.

($1$) DUNE is an LBL accelerator experiment having the baseline length
$L=1300$ km, which operates within the energy range $E_{\nu}=1-10$ GeV. The
oscillation channel mainly observed in DUNE is $\nu_{\mu}\rightarrow\nu_{e}$
and the maximum flux is observed around the energy, $E_{\nu}\sim 2.5-3$ GeV
Brailsford:2018dzn .

($2$) Daya Bay is an SBL reactor experiment with baseline length of $L=2$ km
operating within the energy range $E_{\nu}=0.8-6$ MeV which observes the
$\bar{\nu}_{e}$ beam and measures its survival probability DayaBay:2012fng .
It should be mentioned that at $E_{\nu}\approx 4-5$ MeV an excess of flux is
observed which is larger than the expected value predicted from the theory
an2016measurement . Such discrepancy is also observed at $E_{\nu}\approx 1$
MeV in a much lesser extent an2016measurement . Further, it should be be
mentioned that at $E_{\nu}\sim 1$ MeV the energy resolution of the detector is
$\sim 8.5\%$ an2021antineutrino .

The Earth matter density considered in this work is $\rho=2.8$ gm/c.c. which
generates the matter potential $A=1.01\times 10^{-13}$ eV. The mass eigen
values of the three neutrinos are considered as, $m_{1}=10^{-5}$ eV,
$m_{2}=\sqrt{\Delta m_{21}^{2}+m_{1}^{2}},\leavevmode\nobreak\
m_{3}=\sqrt{\Delta m_{31}^{2}+m_{2}^{2}}$. In this analysis only one diagonal
NSI parameter is considered to be non zero, given by $\eta_{ee}\in[-0.3,0.3]$
Medhi:2021wxj . The LV parameter is chosen as, $a_{e\mu}=10^{-24}$ GeV with
the corresponding complex phase, $\phi_{e\mu}=0$, while the rest of the
parameters are taken to be zero. The standard neutrino oscillation parameters
in case of normal (NO) and inverted mass ordering (IO) are illustrated in
table 1 Esteban:2018azc . The bounds on NU parameters ($\alpha_{ij}$) are
presented in table 2.

DUNE
---
Mass ordering | ($E_{\nu}$,$\delta$) | Unitary (U) | Non-unitary (NU)
| LV | NSI | LV | NSI
NO |  $E_{\nu}$ (GeV) | $1.5-2.5$ | $2-2.5$ | $2-2.5$ | $2-2.75$
$\delta$ (o) | $189-215$ | $189-200$ | $230-257$ | $189-230$
IO |  $E_{\nu}$ (GeV) |  $1$, $3.5-4$ | $1$ | $1$ | $1$
$\delta$ (o) | $252$ | $252$ | $252$ | $252$
Table 3: The specific ranges of the neutrino energy $E_{\nu}$ and CP violating
phase $\delta$ corresponding to $P_{\mu e,max}\gtrsim 0.225$ in case of NO and
$P_{\mu e,max}\gtrsim 0.09$ in case of IO, for the two NP scenarios of LV and
dark NSI, under the assumption of unitary and non-unitary neutrino mixing, in
context of LBL DUNE experimental set up ($E_{\nu}=1-10$ GeV, $L=1300$ km).

In fig. 1 and 2, the variation of the transition probability in the
$\nu_{\mu}\rightarrow\nu_{e}$ channel is depicted in the plane of
($E_{\nu},\delta$) in the scenarios of LV and dark NSI respectively, for NO in
the context of DUNE experimental set up. It can be seen from the figure that
$P_{\mu e}$ attains large value ($\gtrsim 0.225$), with the consideration of
NU mixing only. In the LV scenario, $P_{\mu e}\gtrsim 0.225$ spans the region
confined within $189^{o}\lesssim\delta\lesssim 200^{o}$, corresponding to the
energy range $E_{\nu}\approx 2-2.5$ GeV, while in case of the dark NSI, this
is further extended to a wider region of $(E_{\nu},\delta)$ plane $i.e.$
$189^{o}\lesssim\delta\lesssim 230^{o}$ for the energy interval
$E_{\nu}\approx 2-2.8$ GeV. However, for the unitary neutrino mixing $P_{\mu
e}$ does not attain such a large value and exists below the limit $\sim 0.125$
throughout the entire plane.

Figure 7: Variation of $\bar{\nu}_{ee}$ survival probability $\bar{P}_{ee}$ in
($E_{\nu},\delta$) plane in the scenario of Lorentz violation (LV), in the
context of Daya Bay experiment ($E_{\nu}=0.8-6$ MeV, $L=2$ km) for IO,
considering unitary (left panel) and non-unitary (right panel) neutrino
mixing.

Figure 8: Same as fig. 7 for dark NSI.

In fig. 3 and 4 similar results are presented for DUNE set up for IO. Here the
maximum expanse of $P_{\mu e}$ is much lowered $i.e.\leavevmode\nobreak\
P_{\mu e,max}\gtrsim 0.09$, as compared to the case of NO. From the two
figures, it is observed that $P_{\mu e}$ is very low in the major portion of
the ($E_{\nu},\delta$) plane, for both of the NP scenarios of LV and dark NSI.
Also, $P_{\mu e}$ is comparatively much lower in case of unitary mixing. For
NU mixing in case of LV scenario, $P_{\mu e,max}$ is confined in a very tiny
region at $(E_{\nu},\delta)\sim(1\leavevmode\nobreak\ \text{GeV},252^{o})$,
while for dark NSI, $P_{\mu e,max}$ corresponds to the range
$252^{o}\lesssim\delta\lesssim 260^{o}$ around $E_{\nu}\sim 1$ GeV. The
results for fig. 1, 2, 3 and 4 are illustrated in table 3.

One distinguishable feature to be noticed from the above figures is that DUNE
experiment is capable of exhibiting the signature of dark NSI in terms of
large value of $P_{\mu e}(\gtrsim 0.225)$, when NO is taken into
consideration, around the energy $E_{\nu}\approx 2.5-2.8$ GeV , which
corresponds to the energy of the maximum flux.

Similar to DUNE, the estimation of the survival probability in
$\bar{\nu}_{e}\rightarrow\bar{\nu}_{e}$ channel in context of Daya Bay
experiment are presented in fig. 5 and 6 , for LV and dark NSI scenario
considering unitarity and non-unitarity of neutrino mixing. Here the survival
probability of electron anti-neutrino ($\bar{P}_{ee}$) is plotted in
($E_{\nu},\delta$) plane.

Daya Bay
---
Mass ordering | ($E_{\nu}$,$\delta$) | Unitary (U) | Non-unitary (NU)
| LV | NSI | LV | NSI
NO |  $E_{\nu}$ (MeV) | $1.8-2.3$ | - |  $1.8-2.8$, $1$ | $1.8-2.7$
$\delta$ (o) | $189-257$ | - |  $189-257$, $223-257$ | $189-257$
IO |  $E_{\nu}$ (MeV) | $1.7-2.3$ | - |  1.8-2.2, 1 | $3.2-6$
$\delta$ (o) | $252-307$ | - |  252-307, 265-307 | $252-307$
Table 4: The specific regions of ($E_{\nu},\delta$) plane corresponding to
$\bar{P}_{ee,max}\gtrsim 0.96$ (for both NO and IO), in the two NP scenarios
of LV and dark NSI, under the assumption of unitary and non-unitary neutrino
mixing, in the context of SBL Daya Bay experimental set up ($E_{\nu}=0.8-6$
MeV, $L=2$ km).

In the scenario of LV, $\bar{P}_{ee}$ attains its maximum value
($\bar{P}_{ee,max}\gtrsim 0.96$) in a broad strip around $E_{\nu}\sim 2$ MeV
within the entire allowed range of $\delta$ $i.e.$
$189^{o}\lesssim\delta\lesssim 257^{o}$ for both unitary and non-unitary
mixing. Additionally, in the NU scenario, there is a narrow strip confined
within $223^{o}\lesssim\delta\lesssim 257^{o}$ at $E_{\nu}\sim 1$ MeV. For
dark NSI $\bar{P}_{ee,max}$ corresponds to a broad region around
$E_{\nu}\approx 2$ MeV along the total range of $\delta$ for non-unitary
mixing, while for unitary mixing $\bar{P}_{ee}$ does not even surpass $\sim
0.96$, as observed from the two panels of fig. 6.

Similar plots are presented for IO in the two panels of fig. 7 and 8. It is
observed that the maximum value attained by $\bar{P}_{ee}$ remains the same as
NO $i.e.$ $\bar{P}_{ee,max}\gtrsim 0.96$. In fig. 7, $\bar{P}_{ee,max}$
corresponds to a wide band spread along the full allowed range of $\delta$
$i.e.\leavevmode\nobreak\ 252^{o}\lesssim\delta\lesssim 307^{o}$ for IO in
case of unitary mixing, within $1.8\leavevmode\nobreak\ \text{MeV}\lesssim
E_{\nu}\lesssim 2.3$ MeV. The region is slightly modified and becomes confined
to a region within the energy range, $1.7\leavevmode\nobreak\
\text{MeV}\lesssim E_{\nu}\lesssim 2.1$ MeV, for NU scenario. Moreover,
$\bar{P}_{ee}$ also attains the large value ($\gtrsim 0.96$) around
$E_{\nu}\approx 1$ MeV in a very thin strip , within
$265^{o}\lesssim\delta\lesssim 307^{o}$.

In case of dark NSI, as shown in fig. 8, $\bar{P}_{ee,max}$ corresponds to a
very wide region in the energy interval $4.8\leavevmode\nobreak\
\text{MeV}\lesssim E_{\nu}\lesssim 6$ MeV at $\delta\sim 252^{o}$ for NU
mixing which continues to enlarge further with increase in $\delta$. For
unitary mixing $\bar{P}_{ee}$ remains well below $0.96$ in the entire
$(E_{\nu},\delta)$ plane.

In case of Daya Bay, it is noticed that the large value of $\bar{P}_{ee}$
($\gtrsim 0.96$) appears only in the scenario of dark NSI with NU neutrino
mixing, in the energy range $E_{\nu}\approx 4-5$ MeV at which the excess
$\bar{\nu}_{e}$ flux is observed than that predicted from the theory. The
signature of LV is also apparent in Daya Bay set up for both NO and IO, at
$E_{\nu}\approx 1$ MeV where $\bar{P}_{ee}\gtrsim 0.96$ and the energy
resolution is quite good ($\sim 8.5\%$), although it is more perceptible in
case of NO.

Further, it can be pointed out from all the figures that the nature of $P_{\mu
e}$ and $\bar{P}_{ee}$ in case of LV and dark NSI are quite different in case
of both unitary and non-unitary neutrino mixing scenario in both the
experimental set ups having different baseline lengths, which implies that
these two BSM effects are completely non-degenerate, unlike the existent
degeneracy between vector NSI and LV scenario sahoo2022core .

##  IV Conclusion

The effect of unitary and non-unitary mixing matrix on neutrino oscillation
are analyzed and compared in the new physics scenarios of LV and dark NSI. The
results are determined in the context of LBL and SBL experimental set up, DUNE
and Daya Bay respectively, considering both kinds of neutrino mass orderings.
In DUNE experiment the oscillation probability $P_{\mu e}$ is observed to be
much higher for NO as compared to IO, in both the scenarios of LV and dark
NSI. Also in case of NO, the incorporation of non-unitary mixing matrix is
able to provide relatively higher value of $P_{\mu e}$ in certain regions of
($E_{\nu},\delta$) plane. An interesting feature observed here is that at
$E_{\nu}\approx 2.5-2.8$ GeV large value of $P_{\mu e}(\gtrsim 0.225)$ is
noticed only in case of dark NSI for non-unitary mixing scenario. This points
towards the existence of the new physics in the sector of dark NSI and shows
the potential of DUNE to probe this sector at the energy corresponding to the
maximum flux. However, no such distinguishable feature for LV is observed at
DUNE.

The distinct signature of dark NSI is also observed in case of SBL Daya Bay
set up. At $E_{\nu}\approx 4-5$ MeV, $\bar{P}_{ee}$ possesses large value
($\gtrsim 0.96$) only in the scenario of dark NSI for IO, which could also be
a possible explanation for the observed discrepancy in $\bar{\nu}_{e}$ flux at
this energy interval. Further, $\bar{P}_{ee}$ is also large $(\gtrsim 0.96)$
at $E_{\nu}\approx 1$ MeV at which the energy resolution of the detector is
$\sim 8.5\%$ . This is only observed in the LV scenario for non-unitary
neutrino mixing. Although such feature is noticed for both NO and IO, it is
more prominently discernible in case of NO only.

## References

  * [1] A. Zee.  A Theory of Lepton Number Violation, Neutrino Majorana Mass, and Oscillation.  Phys. Lett. B, 93:389, 1980.  [Erratum: Phys.Lett.B 95, 461 (1980)].
  * [2] Michele Frigerio and Alexei Yu. Smirnov.  Radiative corrections to neutrino mass matrix in the standard model and beyond.  JHEP, 02:004, 2003.
  * [3] Amine Ahriche, Kristian L. McDonald, and Salah Nasri.  A Radiative Model for the Weak Scale and Neutrino Mass via Dark Matter.  JHEP, 02:038, 2016.
  * [4] Jörn Kersten and Alexei Yu. Smirnov.  Right-Handed Neutrinos at CERN LHC and the Mechanism of Neutrino Mass Generation.  Phys. Rev. D, 76:073005, 2007.
  * [5] Yi-Lei Tang.  One Right-handed Neutrino to Generate Complete Neutrino Mass Spectrum in the Framework of NMSSM.  Nucl. Phys. B, 890:263–278, 2014.
  * [6] S. F. King.  Neutrino mass models.  Rept. Prog. Phys., 67:107–158, 2004.
  * [7] H. Mohseni Sadjadi and A. P. Khosravi.  Symmetry breaking, and the effect of matter density on neutrino oscillation.  JCAP, 04:008, 2018.
  * [8] Manfred Lindner, Steffen Schmidt, and Juri Smirnov.  Neutrino Masses and Conformal Electro-Weak Symmetry Breaking.  JHEP, 10:177, 2014.
  * [9] G. C. Branco, J. T. Penedo, Pedro M. F. Pereira, M. N. Rebelo, and J. I. Silva-Marcos.  Type-I Seesaw with eV-Scale Neutrinos.  JHEP, 07:164, 2020.
  * [10] D. V. Forero, S. Morisi, M. Tortola, and J. W. F. Valle.  Lepton flavor violation and non-unitary lepton mixing in low-scale type-I seesaw.  JHEP, 09:142, 2011.
  * [11] Michal Malinsky, Tommy Ohlsson, Zhi-zhong Xing, and He Zhang.  Non-unitary neutrino mixing and CP violation in the minimal inverse seesaw model.  Phys. Lett. B, 679:242–248, 2009.
  * [12] D. V. Forero, C. Giunti, C. A. Ternes, and M. Tortola.  Nonunitary neutrino mixing in short and long-baseline experiments.  Phys. Rev. D, 104(7):075030, 2021.
  * [13] Jeffrey M. Berryman, André de Gouvêa, Daniel Hernández, and Roberto L. N. Oliveira.  Non-Unitary Neutrino Propagation From Neutrino Decay.  Phys. Lett. B, 742:74–79, 2015.
  * [14] Yi Cai, Tao Han, Tong Li, and Richard Ruiz.  Lepton Number Violation: Seesaw Models and Their Collider Tests.  Front. in Phys., 6:40, 2018.
  * [15] Luis Salvador Miranda, Pedro Pasquini, Ushak Rahaman, and Soebur Razzaque.  Searching for non-unitary neutrino oscillations in the present T2K and NO$\nu$A data.  Eur. Phys. J. C, 81(5):444, 2021.
  * [16] Albert Abraham Michelson and Edward Williams Morley.  On the Relative Motion of the Earth and the Luminiferous Ether.  Am. J. Sci., 34:333–345, 1887.
  * [17] Roy J. Kennedy and Edward M. Thorndike.  Experimental establishment of the relativity of time.  Phys. Rev., 42:400–418, Nov 1932.
  * [18] Archil Kobakhidze and Bruce H. J. McKellar.  Particle interference as a test of Lorentz-violating electrodynamics.  Phys. Rev. D, 76:093004, 2007.
  * [19] Andrea Petrucci.  Clues to detect the Pilot Wave in a photon double-slit interference experiment.  arXiv e-prints, page arXiv:1707.09016, July 2017.
  * [20] Alejandro Cabo Montes de Oca and Nana Geraldine Cabo Bizet.  Neutrino oscillations as many-particle induced interference between distinguishable particles.  5 2020.
  * [21] Gabriela Barenboim, Mehedi Masud, Christoph A. Ternes, and Mariam Tórtola.  Exploring the intrinsic Lorentz-violating parameters at DUNE.  Phys. Lett. B, 788:308–315, 2019.
  * [22] K. Abe et al.  Test of Lorentz invariance with atmospheric neutrinos.  Phys. Rev. D, 91(5):052003, 2015.
  * [23] Sheldon Lee Glashow.  Atmospheric neutrino constraints on Lorentz violation.  7 2004.
  * [24] L. Wolfenstein.  Neutrino oscillations in matter.  Phys. Rev. D, 17:2369–2374, May 1978.
  * [25] Y. Farzan and M. Tortola.  Neutrino oscillations and Non-Standard Interactions.  Front. in Phys., 6:10, 2018.
  * [26] Neutrino Non-Standard Interactions: A Status Report, volume 2, 2019.
  * [27] Tommy Ohlsson.  Status of non-standard neutrino interactions.  Rept. Prog. Phys., 76:044201, 2013.
  * [28] Stefan Antusch, Jochen P. Baumann, and Enrique Fernandez-Martinez.  Non-Standard Neutrino Interactions with Matter from Physics Beyond the Standard Model.  Nucl. Phys. B, 810:369–388, 2009.
  * [29] Shao-Feng Ge and Stephen J. Parke.  Scalar Nonstandard Interactions in Neutrino Oscillation.  Phys. Rev. Lett., 122(21):211801, 2019.
  * [30] Abinash Medhi, Debajyoti Dutta, and Moon Moon Devi.  Exploring the effects of scalar non standard interactions on the CP violation sensitivity at DUNE.  JHEP, 06:129, 2022.
  * [31] Shao-Feng Ge.  New physics with scalar and dark non-standard interactions in neutrino oscillation.  Journal of Physics: Conference Series, 1468:012125, 02 2020.
  * [32] M. Agostini et al.  First Simultaneous Precision Spectroscopy of $pp$, 7Be, and $pep$ Solar Neutrinos with Borexino Phase-II.  Phys. Rev. D, 100(8):082004, 2019.
  * [33] Amir N. Khan, Werner Rodejohann, and Xun-Jie Xu.  Borexino and general neutrino interactions.  Phys. Rev. D, 101(5):055047, 2020.
  * [34] David E. Jaffe.  A decade of discoveries by the Daya Bay reactor neutrino experiment.  Mod. Phys. Lett. A, 36(21):2130021, 2021.
  * [35] Mona Dentler, Álvaro Hernández-Cabezudo, Joachim Kopp, Michele Maltoni, and Thomas Schwetz.  Sterile neutrinos or flux uncertainties? — Status of the reactor anti-neutrino anomaly.  JHEP, 11:099, 2017.
  * [36] Jeffrey M. Berryman, Vedran Brdar, and Patrick Huber.  Particle physics origin of the 5 MeV bump in the reactor antineutrino spectrum?  Phys. Rev. D, 99(5):055045, 2019.
  * [37] M. Aker et al.  KATRIN: status and prospects for the neutrino mass and beyond.  J. Phys. G, 49(10):100501, 2022.
  * [38] Andre de Gouvea.  Neutrino Mass Models.  Ann. Rev. Nucl. Part. Sci., 66:197–217, 2016.
  * [39] Manfred Lindner, Tommy Ohlsson, and Gerhart Seidl.  Seesaw mechanisms for Dirac and Majorana neutrino masses.  Phys. Rev. D, 65:053014, 2002.
  * [40] Shinya Kanemura, Takehiro Nabeshima, and Hiroaki Sugiyama.  Neutrino Masses from Loop-Induced Dirac Yukawa Couplings.  Phys. Lett. B, 703:66–70, 2011.
  * [41] E. Fernandez-Martinez, M. B. Gavela, J. Lopez-Pavon, and O. Yasuda.  CP-violation from non-unitary leptonic mixing.  Phys. Lett. B, 649:427–435, 2007.
  * [42] F. J. Escrihuela, D. V. Forero, O. G. Miranda, M. Tortola, and J. W. F. Valle.  On the description of nonunitary neutrino mixing.  Phys. Rev. D, 92(5):053009, 2015.  [Erratum: Phys.Rev.D 93, 119905 (2016)].
  * [43] S. Antusch, C. Biggio, E. Fernandez-Martinez, M. B. Gavela, and J. Lopez-Pavon.  Unitarity of the Leptonic Mixing Matrix.  JHEP, 10:084, 2006.
  * [44] F. J. Escrihuela, L. J. Flores, and O. G. Miranda.  Neutrino counting experiments and non-unitarity from LEP and future experiments.  Phys. Lett. B, 802:135241, 2020.
  * [45] Nikolaos E. Mavromatos.  Lorentz Invariance Violation from String Theory.  PoS, QG-PH:027, 2007.
  * [46] M. Tsoukalas Cuadros-Melgar:2011whf B. Cuadros-Melgar, E. Papantonopoulos and V. Zamarias.  Massive Gravity with Anisotropic Scaling.  Phys. Rev. D, 85:124035, 2012.
  * [47] J. W. Moffat.  Lorentz Violation of Quantum Gravity.  Class. Quant. Grav., 27:135016, 2010.
  * [48] Astrid Eichhorn, Alessia Platania, and Marc Schiffer.  Lorentz invariance violations in the interplay of quantum gravity with matter.  Phys. Rev. D, 102(2):026007, 2020.
  * [49] Robert Bluhm.  Observational Constraints on Local Lorentz Invariance, pages 485–507.  2014\.
  * [50] Koji Hashimoto.  Lorentz violation and vacuum structure in string theories.  In International Journal of Modern Physics Conference Series, volume 30 of International Journal of Modern Physics Conference Series, page 1460274, May 2014.
  * [51] Petr Horava.  Spectral Dimension of the Universe in Quantum Gravity at a Lifshitz Point.  Phys. Rev. Lett., 102:161301, 2009.
  * [52] Alejandro Cabo and Nana Geraldine Cabo Bizet.  About the neutrino oscillation-like effects in general physical systems: On interference between distinguishable particles…  Eur. Phys. J. Plus, 136(10):1042, 2021.
  * [53] John N. Bahcall, V. Barger, and Danny Marfatia.  How accurately can one test CPT conservation with reactor and solar neutrino experiments?  Phys. Lett. B, 534:120–123, 2002.
  * [54] Vito Antonelli, L. Miramonti, and M. D. C. Torri.  Neutrino oscillations and Lorentz Invariance Violation in a Finslerian Geometrical model.  Eur. Phys. J. C, 78(8):667, 2018.
  * [55] Marco Danilo Claudio Torri.  Neutrino Oscillations and Lorentz Invariance Violation.  Universe, 6(3):37, 2020.
  * [56] V. Alan Kostelecky and Matthew Mewes.  Lorentz violation and short-baseline neutrino experiments.  Phys. Rev. D, 70:076002, 2004.
  * [57] Sanjib Kumar Agarwalla and Mehedi Masud.  Can Lorentz invariance violation affect the sensitivity of deep underground neutrino experiment?  Eur. Phys. J. C, 80(8):716, 2020.
  * [58] Sadashiv Sahoo, Anil Kumar, and Sanjib Kumar Agarwalla.  Probing Lorentz Invariance Violation with atmospheric neutrinos at INO-ICAL.  JHEP, 03:050, 2022.
  * [59] P. Arias, J. Gamboa, J. Lopez-Sarrion, F. Mendez, and Ashok K. Das.  CPT / Lorentz invariance violation and neutrino oscillation.  Phys. Lett. B, 650:401–406, 2007.
  * [60] Jorge S. Diaz.  Correspondence between nonstandard interactions and CPT violation in neutrino oscillations.  6 2015.
  * [61] R. et al Abbasi.  Search for a lorentz-violating sidereal signal with atmospheric neutrinos in icecube.  Phys. Rev. D, 82:112003, Dec 2010.
  * [62] P. Adamson et al.  A Search for Lorentz Invariance and CPT Violation with the MINOS Far Detector.  Phys. Rev. Lett., 105:151601, 2010.
  * [63] Y. Abe et al.  First test of lorentz violation with a reactor-based antineutrino experiment.  Phys. Rev. D, 86:112009, Dec 2012.
  * [64] Sadashiv Sahoo, Anil Kumar, Sanjib Kumar Agarwalla, and Amol Dighe.  Core-passing atmospheric neutrinos: a unique probe to discriminate between lorentz violation and non-standard interactions.  arXiv preprint arXiv:2205.05134, 2022.
  * [65] Shao-Feng Ge and Hitoshi Murayama.  Apparent CPT Violation in Neutrino Oscillation from Dark Non-Standard Interactions.  4 2019.
  * [66] Tommy Ohlsson and Hakan Snellman.  Three flavor neutrino oscillations in matter.  J. Math. Phys., 41:2768–2788, 2000.  [Erratum: J.Math.Phys. 42, 2345 (2001)].
  * [67] Ivan Esteban, M. C. Gonzalez-Garcia, Alvaro Hernandez-Cabezudo, Michele Maltoni, and Thomas Schwetz.  Global analysis of three-flavour neutrino oscillations: synergies and tensions in the determination of $\theta_{23}$, $\delta_{CP}$, and the mass ordering.  JHEP, 01:106, 2019.
  * [68] Dominic Brailsford.  DUNE: Status and Perspectives.  In Prospects in Neutrino Physics, pages 4–12, 4 2018.
  * [69] F. P. An et al.  Observation of electron-antineutrino disappearance at Daya Bay.  Phys. Rev. Lett., 108:171803, 2012.
  * [70] Feng Peng An, AB Balantekin, HR Band, M Bishai, S Blyth, I Butorov, D Cao, GF Cao, J Cao, WR Cen, et al.  Measurement of the reactor antineutrino flux and spectrum at daya bay.  Physical review letters, 116(6):061801, 2016.
  * [71] FP An, AB Balantekin, M Bishai, S Blyth, GF Cao, J Cao, JF Chang, Y Chang, HS Chen, SM Chen, et al.  Antineutrino energy spectrum unfolding based on the daya bay measurement and its applications.  Chinese Physics C, 45(7):073001, 2021.