375386

# A criterion for continuity in a parameter of solutions to generic boundary-
value problems for higher-order differential systems

Vladimir Mikhailets  Institute of Mathematics, National Academy of Sciences of
Ukraine, 3 Tereshchenkivs’ka, Kyiv, 01601, Ukraine; National Technical
University of Ukraine ”Kyiv Polytechnic Institute”, 37 Prospect Peremogy,
Kyiv, 03056, Ukraine<EMAIL_ADDRESS>Aleksandr Murach  Institute of
Mathematics, National Academy of Sciences of Ukraine, 3 Tereshchenkivs’ka,
Kyiv, 01601, Ukraine; Chernihiv National Pedagogical University, 53 Het’mana
Polubotka, Chernihiv, 14013, Ukraine<EMAIL_ADDRESS>Vitalii Soldatov
Institute of Mathematics, National Academy of Sciences of Ukraine, 3
Tereshchenkivs’ka, Kyiv, 01601, Ukraine<EMAIL_ADDRESS>
(Date: 15/08/2016)

###### Abstract.

We consider the most general class of linear boundary-value problems for
ordinary differential systems, of order $r\geq 1$, whose solutions belong to
the complex space $C^{(n+r)}$, with $0\leq n\in\mathbb{Z}$. The boundary
conditions can contain derivatives of order $l$, with $r\leq l\leq n+r$, of
the solutions. We obtain a constructive criterion under which the solutions to
these problems are continuous with respect to the parameter in the normed
space $C^{(n+r)}$. We also obtain a two-sided estimate for the degree of
convergence of these solutions.

###### Key words and phrases:

Differential system, boundary-value problem, continuity in parameter.

###### 2010 Mathematics Subject Classification:

34B08

##  1\. Introduction

Questions concerning the validity of passage to the limit in parameter-
dependent Cauchy problems and boundary-value problems are very important in
the theory of ordinary differential equations. Fundamental results on the
continuous in the parameter of solutions to the Cauchy problem for nonlinear
systems were obtained by Gikhman [3], Krasnosel’skii and S. Krein [16],
Kurzweil and Vorel [17]. These results were refined and supplemented for
linear systems by Levin [18], Opial [24], Reid [25], and Nguyen [23].

Parameter-dependent boundary-value problems are far less studied than the
Cauchy problem, which is connected with great diversity of boundary
conditions. Kiguradze [10, 11, 12] and Ashordia [1] introduced and studied a
class of general linear boundary-value problems for systems of first-order
ordinary differential equations. The solutions $y$ to these problems are
supposed to be absolutely continuous on a compact interval $[a,b]$, and the
boundary condition is given in the form $By=q$, where
$B:C([a,b],\mathbb{R}^{m})\to\mathbb{R}^{m}$ is an arbitrary continuous linear
operator ($m$ is the number of differential equations in the system).
Kiguradze and Ashordia obtained conditions under which the solutions to these
problems are continuous in the normed space $C([a,b],\mathbb{R}^{m})$ with
respect to the parameter. Recently [15, 20] these results were refined and
extended to complex-valued functions and systems of higher-order differential
equations.

Other broad classes of linear boundary-value problems were introduced and
investigated in [22, 14, 4] and [19, 27]. These classes relate closely to the
classical scales of the function Sobolev spaces and the spaces $C^{(l)}$ of
continuously differential functions, respectively. The boundary conditions for
the problems from these classes are given in the form $By=q$, where $B$ is an
arbitrary continuous linear operator that acts from the corresponding function
space to the finite-dimensional complex space. It is naturally to call these
boundary-value problems generic with respect to the indicated function space.
Constructive sufficient conditions for continuity of their solutions in the
parameter in the corresponding spaces are obtained. These results were applied
to the investigation of multipoint boundary-value problems [13], Green’s
matrices of boundary-value problems [15, 20], in the spectral theory of
differential operators with distributional coefficients [6, 5, 7].

In the present paper, we investigate the class of generic boundary-value
problems for differential systems of order $r\geq 1$ with respect to the
complex space $C^{(n+r)}$, where the integer $n\geq 0$. This class was
introduced in [19] for $r=1$ and in [27] for $r\geq 2$. We assume that all
coefficients and right-hand sides of these systems belong to the space
$C^{(n)}$. Since the solution $z$ to such a system runs through the whole
space $(C^{(n+r)})^{m}$ provided that its right-hand side runs through the
whole space $(C^{(n)})^{m}$, the boundary condition is given in the most
general form $Bz=q$ where $B:(C^{(n+r)})^{m}\to\mathbb{C}^{rm}$ is an
arbitrary continuous linear operator and $q\in\mathbb{C}^{rm}$. This condition
covers all the classical types of boundary conditions such as initial
conditions in the Cauchy problem, multipoint conditions, integral conditions,
and also nonclassical conditions containing the derivatives $z^{(l)}$ with
$r\leq l\leq n+r$ of the unknown function.

The sited papers [19, 27] provide a constructive sufficient conditions under
which the solution to an arbitrary parameter-dependent boundary-value problem
from this class depends continuously on the parameter in the normed space
$(C^{(n+r)})^{m}$. The goal of the present paper is to show that these
conditions are also necessary. Thus, we will prove a criterion for continuity
in the parameter of the solution to this problem in the indicated space.
Besides, we will obtain a two-sided estimate for the degree of convergence of
the solutions.

The approach developed in this paper can be applied to generic boundary-value
problems with respect to other function spaces [8, 21].

##  2\. Main results

We arbitrarily choose a compact interval $[a,b]\subset\mathbb{R}$ and integers
$n\geq 0$, $r\geq 1$, and $m\geq 1$. We use the Banach spaces

$(C^{(l)})^{m}:=C^{(l)}([a,b],\mathbb{C}^{m})\quad\mbox{and}\quad(C^{(l)})^{m\times
m}:=C^{(l)}([a,b],\mathbb{C}^{m\times m}),$ (1)

with $0\leq l\in\mathbb{Z}$. They consist respectively of all vector-valued
functions and $(m\times m)$-matrix-valued functions whose components belong to
the space $C^{(l)}:=C^{(l)}([a,b],\mathbb{C})$ of all $l$ times continuously
differentiable functions $x:[a,b]\rightarrow\mathbb{C}$. This space is endowed
with the norm

$\|x\|_{(l)}:=\sum_{j=0}^{l}\max\bigl{\\{}|x^{(j)}(t)|:\,t\in[a,b]\bigr{\\}}.$

The norms in the spaces (1) are equal to the sum of the norms in $C^{(l)}$ of
all components of the vector-valued or matrix-valued functions and are also
denoted by $\|\cdot\|_{(l)}$. It will be always clear from the context to
which space (scalar or vector-valued or matrix-valued functions) the norm
$\|\cdot\|_{(l)}$ relates.

Let $\varepsilon_{0}>0$. We consider the following linear boundary-value
problem depending on the parameter $\varepsilon\in[0,\varepsilon_{0})$:

$\displaystyle L(\varepsilon)z(t,\varepsilon)\equiv
z^{(r)}(t,\varepsilon)+\sum_{j=1}^{r}K_{r-j}(t,\varepsilon)z^{(r-j)}(t,\varepsilon)=f(t,\varepsilon),\quad
a\leq t\leq b,$ (2) $\displaystyle
B(\varepsilon)z(\cdot,\varepsilon)=q(\varepsilon).$ (3)

We suppose for every $\varepsilon\in[0,\varepsilon_{0})$ that

$z(\cdot,\varepsilon)\in(C^{(n+r)})^{m}$

is an unknown vector-valued function, whereas the matrix-valued functions
$K_{r-j}(\cdot,\varepsilon)\in(C^{(n)})^{m\times m}$, with
$j\in\\{1,\ldots,r\\}$, vector-valued function
$f(\cdot,\varepsilon)\in(C^{(n)})^{m}$, continuous linear operator

$B(\varepsilon):(C^{(n+r)})^{m}\to\mathbb{C}^{rm}$ (4)

and vector $q(\varepsilon)\in\mathbb{C}^{rm}$ are arbitrarily given. In the
paper, we interpret vectors as columns.

With this problem, we associate the continuous linear operator

$(L(\varepsilon),B(\varepsilon)):(C^{(n+r)})^{m}\rightarrow(C^{(n)})^{m}\times\mathbb{C}^{rm}.$
(5)

The operator (5) is Fredholm with index zero for every
$\varepsilon\in[0,\varepsilon_{0})$ because this operator is a finite-
dimensional perturbation of the isomorphism that acts between the spaces
$(C^{(n+r)})^{m}$ and $(C^{(n)})^{m}\times\mathbb{C}^{rm}$ and corresponds to
the Cauchy problem for the differential system (2) (see, e.g., [9]).

For the boundary-value problem (2), (3), we consider the following three

Limit Conditions as $\varepsilon\to 0+$:

  * (I)

$K_{r-j}(\cdot,\varepsilon)\to K_{r-j}(\cdot,0)$ in $(C^{(n)})^{m\times m}$
for each $j\in\\{1,\ldots,r\\}$;

  * (II)

$B(\varepsilon)z\to B(0)z$ in $\mathbb{C}^{rm}$ for every
$z\in(C^{(n+r)})^{m}$;

  * (III)

$f(\cdot,\varepsilon)\to f(\cdot,0)$ in $(C^{(n)})^{m}$, and
$q(\varepsilon)\to q(0)$ in $\mathbb{C}^{rm}$.

We also consider

Condition (0). The limiting homogeneous boundary-value problem

$L(0)\,z(t,0)=0,\quad a\leq t\leq b,\quad\mbox{and}\quad B(0)\,z(\cdot,0)=0$
(6)

has only the trivial solution.

Note that a criterion for the problem (6) satisfy Condition (0) is given in
[27, Theorem 2]. This criterion is formulated in terms of the matrices
$K_{r-j}(\cdot,0)$ and the operator $B(0)$.

Now we introduce our

Basic Definition. We say that the solution to the boundary-value problem (2),
(3) is continuous in the parameter $\varepsilon$ at $\varepsilon=0$ if the
following two conditions are fulfilled:

  * $(\ast)$

There exists a positive number $\varepsilon_{1}<\varepsilon_{0}$ such that for
an arbitrary number $\varepsilon\in\nobreak[0,\varepsilon_{1})$, function
$f(\cdot,\varepsilon)\in(C^{(n)})^{m}$, and vector
$q(\varepsilon)\in\mathbb{C}^{rm}$ this problem has a unique solution
$z(\cdot,\varepsilon)\in(C^{(n+r)})^{m}$.

  * $(\ast\ast)$

Limit Condition (III) implies the convergence $z(\cdot,\varepsilon)\to
z(\cdot,0)$ in $(C^{(n+r)})^{m}$ as $\varepsilon\to 0+$.

Main Theorem. The solution to the boundary-value problem (2), (3) is
continuous in the parameter $\varepsilon$ at $\varepsilon=0$ if and only if
this problem satisfies Condition (0) and Limit Conditions (I) and (II).

We supplement this result by the following theorem. Let

$d_{n}(\varepsilon):=\|L(\varepsilon)z(\cdot,0)-f(\cdot,\varepsilon)\|_{(n)}+\|B(\varepsilon)z(\cdot,0)-q(\varepsilon)\|_{\mathbb{C}^{rm}}.$

Estimation Theorem. Suppose that the boundary-value problem (2), (3) satisfies
Condition (0) and Limit Conditions (I) and (II). Then there exist positive
numbers $\varepsilon_{2}<\varepsilon_{1}$ and $\varkappa_{1}$, $\varkappa_{2}$
such that for every $\varepsilon\in(0,\varepsilon_{2})$ the following two-
sided estimate holds:

$\varkappa_{1}\,d_{n}(\varepsilon)\leq\|z(\cdot,0)-z(\cdot,\varepsilon)\|_{(n+r)}\leq\varkappa_{2}\,d_{n}(\varepsilon).$
(7)

Here, the numbers $\varepsilon_{2}$, $\varkappa_{1}$, and $\varkappa_{2}$ do
not depend on $z(\cdot,0)$ and $z(\cdot,\varepsilon)$.

According to (7), the error and discrepancy of the solution
$z(\cdot,\varepsilon)$ to the boundary-value problem (2), (3) are of the same
degree provided that we consider $z(\cdot,0)$ as an approximate solution to
this problem.

We will prove these theorems in the next section.

Note that every linear continuous operator (4) can be uniquely represented in
the form

$B(\varepsilon)z=\sum_{k=1}^{n+r}\beta_{k}(\varepsilon)\,z^{(k-1)}(a)+\int_{a}^{b}(d\Phi(t,\varepsilon))z^{(n+r)}(t),\quad\mbox{with}\quad
z\in(C^{(n+r)})^{m}.$ (8)

Here, each $\beta_{k}(\varepsilon)$ is a number $rm\times m$-matrix, and
$\Phi(\cdot,\varepsilon)$ is an $rm\times m$-matrix-valued function formed by
scalar functions that are of bounded variation on $[a,b]$, right-continuous on
$(a,b)$, and equal to zero at $t=a$. (Certainly, the integral is understood in
the Riemann-Stieltjes sense.) Representation (8) follows from the known
description of the dual of $C^{(n+r)}$; see, e.g., [2, p. 344]. Using this
representation, we can reformulate Limit Condition (II) in an explicit form.
Namely, Limit Condition (II) is equivalent to that the following four
conditions are fulfilled as $\varepsilon\to 0+$:

  * (2a)

$\beta_{k}(\varepsilon)\to\beta_{k}(0)$ for every $k\in\\{1,\ldots,n+r\\}$;

  * (2b)

$\|V_{a}^{b}\Phi(\cdot,\varepsilon)\|_{\mathbb{C}^{rm\times m}}=O(1)$;

  * (2c)

$\Phi(b,\varepsilon)\to\Phi(b,0)$;

  * (2d)

$\int_{a}^{t}\Phi(s,\varepsilon)ds\to\int_{a}^{t}\Phi(s,0)ds$ for every
$t\in(a,b]$.

(Of course, the convergence is considered in $\mathbb{C}^{rm\times m}$.) This
equivalence follows from the F. Riesz criterion of the weak convergence of
linear continuous functionals on $C^{(0)}$ (see, e.g., [26, Ch. III, Sect.
55]). It is useful to compare these conditions with the criterion of the norm
convergence of the operators $B(\varepsilon)\to B(0)$ as $\varepsilon\to 0+$.
It asserts that this norm convergence is equivalent to the fulfilment of
condition (2a) and

$V_{a}^{b}\bigl{(}\Phi(\cdot,\varepsilon)-\Phi(\cdot,0)\bigr{)}\to 0$ (9)

as $\varepsilon\to 0+$. Condition (9) is much stronger than the system of
conditions (2b)–(2d). Indeed, condition (9) implies the uniform convergence of
the functions $\Phi(t,\varepsilon)$ to $\Phi(t,0)$ on $[a,b]$ as
$\varepsilon\to 0+$ whereas conditions (2b)–(2d) do not entail the pointwise
convergence of these functions at least in one point of $(a,b)$.

##  3\. Proofs

###### Proof of Main Theorem..

The sufficiency of Conditions (0), (I) and (II) for problem (2), (3) to
satisfy Basic definition is proved in [27, Theorem 3]. Let us prove the
necessity. Assume that this problem satisfies Basic definition. Then Condition
(0) is fulfilled. It remains to prove that this problem satisfies Conditions
(I) and (II). We divide our reasoning into three steps.

_Step $1$._ We will prove here that the boundary-value problem (2), (3)
satisfies Limit Condition (I). In the $r\geq 2$ case, we previously reduce
this problem to a boundary-value problem for first-order differential system.
To this end we put, as usual,

$\displaystyle
y(\cdot,\varepsilon):=\mathrm{col}\bigl{(}z(\cdot,\varepsilon),z^{\prime}(\cdot,\varepsilon),\ldots,z^{(r-1)}(\cdot,\varepsilon)\bigr{)}\in(C^{(n+1)})^{rm},$
$\displaystyle
g(\cdot,\varepsilon):=\mathrm{col}\bigl{(}0,f(\cdot,\varepsilon)\bigr{)}\in(C^{(n)})^{rm},$

and

$A(\cdot,\varepsilon):=\left(\begin{array}[]{ccccc}O_{m}&I_{m}&O_{m}&\ldots&O_{m}\\\
O_{m}&O_{m}&I_{m}&\ldots&O_{m}\\\ \vdots&\vdots&\vdots&\ddots&\vdots\\\
O_{m}&O_{m}&O_{m}&\ldots&I_{m}\\\
K_{0}(\cdot,\varepsilon)&K_{1}(\cdot,\varepsilon)&K_{2}(\cdot,\varepsilon)&\ldots&K_{r-1}(\cdot,\varepsilon)\\\
\end{array}\right)\in(C^{(n)})^{rm\times rm},$

with $O_{m}$ and $I_{m}$ denoting zero and identity $(m\times m)$-matrices,
respectively. In view of representation (8), we also put

$N(\varepsilon)y:=\sum_{k=1}^{r-1}\beta_{k}(\varepsilon)y_{k}(a)+\sum_{k=r}^{n+r}\beta_{k}(\varepsilon)y_{r}^{(k-r)}(a)+\int_{a}^{b}(d\Phi(t,\varepsilon))y_{r}^{(n+1)}(t)$
(10)

for an arbitrary vector-valued function $y=\mathrm{col}(y_{1},\ldots,y_{r})$
with $y_{1},\ldots,y_{r}\in(C^{(n+1)})^{m}$. The linear mapping $y\mapsto Ny$
acts continuously from $(C^{(n+1)})^{rm}$ to $\mathbb{C}^{rm}$.

Evidently, a function $z(\cdot,\varepsilon)\in(C^{(n+r)})^{m}$ is a solution
to the boundary-value problem (2), (3) if and only if the function
$y(\cdot,\varepsilon)$ is a solution to the boundary-value problem

$\displaystyle
y^{\prime}(t,\varepsilon)+A(t,\varepsilon)y(t,\varepsilon)=g(t,\varepsilon),\quad
a\leq t\leq b,$ (11) $\displaystyle
N(\varepsilon)y(\cdot,\varepsilon)=q(\varepsilon).$ (12)

In the $r=1$ case, we put $y(\cdot,\varepsilon):=z(\cdot,\varepsilon)$,
$g(\cdot,\varepsilon):=f(\cdot,\varepsilon)$,
$A(\cdot,\varepsilon):=A_{0}(\cdot,\varepsilon)$, and
$N(\varepsilon):=B(\varepsilon)$ for the sake of uniformity in notation on
Step 1, so the problem (2), (3) coincides with the problem (11), (12).

According to condition ($\ast$) of Basic Definition, the latter problem has a
unique solution $y(t,\varepsilon)$ for every
$\varepsilon\in[0,\varepsilon_{1})$. Limit Condition (I) is equivalent to that
$A(\cdot,\varepsilon)\to A(\cdot,0)$ in $(C^{(n)})^{rm\times rm}$ as
$\varepsilon\to 0+$. Let us prove this convergence.

We previously note the following. If $g(\cdot,\varepsilon)$ and
$q(\varepsilon)$ are independent of $\varepsilon\in[0,\varepsilon_{1})$, then
$z(\cdot,\varepsilon)\to z(\cdot,0)$ in $(C^{(n+r)})^{m}$ as $\varepsilon\to
0+$ by condition ($\ast\ast$) of Basic Definition. The latter convergence is
equivalent to that $y(\cdot,\varepsilon)\to y(\cdot,0)$ in $(C^{(n+1)})^{rm}$
as $\varepsilon\to 0+$.

Given $\varepsilon\in[0,\varepsilon_{1})$, we now consider the matrix
boundary-value problem

$\displaystyle
Y^{\prime}(t,\varepsilon)+A(t,\varepsilon)Y(t,\varepsilon)=0_{rm},\quad a\leq
t\leq b,$ (13) $\displaystyle[N(\varepsilon)Y(\cdot,\varepsilon)]=I_{rm}.$
(14)

Here, the unknown matrix-valued function
$Y(\cdot,\varepsilon):=(y_{j,k}(\cdot,\varepsilon))_{j,k=1}^{rm}$ belongs to
the space $(C^{(n+1)})^{rm\times rm}$, and we put

$[N(\varepsilon)Y(\cdot,\varepsilon)]:=\left(N(\varepsilon)\begin{pmatrix}y_{1,1}(\cdot,\varepsilon)\\\
\vdots\\\ y_{rm,1}(\cdot,\varepsilon)\\\
\end{pmatrix}\;\ldots\;N(\varepsilon)\begin{pmatrix}y_{1,rm}(\cdot,\varepsilon)\\\
\vdots\\\ y_{rm,rm}(\cdot,\varepsilon)\\\ \end{pmatrix}\right).$

This problem is a union of $rm$ boundary-value problems (11), (12) whose
right-hand sides are independent of $\varepsilon$. So, it has a unique
solution $Y(\cdot,\varepsilon)\in(C^{(n+1)})^{rm\times rm}$, and
$Y(\cdot,\varepsilon)\to Y(\cdot,0)$ in $(C^{(n+1)})^{rm\times rm}$ as
$\varepsilon\to 0+$. Besides, $\det Y(t,\varepsilon)\neq 0$ for every
$t\in[a,b]$. Indeed, if this were wrong, the function columns of the matrix
$Y(\cdot,\varepsilon)$ would be linearly dependent, which would contradict the
boundary condition (14). Thus, we obtain the required convergence

$A(\cdot,\varepsilon)=-Y^{\prime}(\cdot,\varepsilon)(Y(\cdot,\varepsilon))^{-1}\to-Y^{\prime}(\cdot,0)(Y(\cdot,0))^{-1}=A(\cdot,0)$

in $(C^{(n)})^{rm\times rm}$ as $\varepsilon\to 0+$. It implies that

$\|K_{r-j}(\varepsilon)\|_{(n)}=O(1)\quad\mbox{as}\quad\varepsilon\to
0+\quad\mbox{for each}\quad j\in\\{1,\ldots,r\\}.$ (15)

_Step $2$._ Let us prove that $\|B(\varepsilon)\|=O(1)$ as $\varepsilon\to
0+$; here, $\|\cdot\|$ denotes the norm of a bounded operator from
$(C^{(n+r)})^{m}$ to $\mathbb{C}^{rm}$. Suppose the contrary; then there
exists a sequence
$(\varepsilon^{(k)})_{k=1}^{\infty}\subset(0,\varepsilon_{1})$ such that
$\varepsilon^{(k)}\to 0$ and $0<\|B(\varepsilon^{(k)})\|\to\infty$ as
$\varepsilon\to 0+$. For every integer $k\geq 1$, we can choose a function
$w_{k}\in(C^{(n+r)})^{m}$ such that $\|w_{k}\|_{(n+r)}=1$ and
$\|B(\varepsilon^{(k)})w_{k}\|_{\mathbb{C}^{rm}}\geq\|B(\varepsilon^{(k)})\|/2$.
We now let $z(\cdot,\varepsilon^{(k)}):=\|B(\varepsilon^{(k)})\|^{-1}\,w_{k}$,
$f(\cdot,\varepsilon^{(k)}):=L(\varepsilon^{(k)})\,z(\cdot,\varepsilon^{(k)})$,
and $q(\varepsilon^{(k)}):=B(\varepsilon^{(k)})\,z(\cdot,\varepsilon^{(k)})$.
Since $z(\cdot,\varepsilon^{(k)})\to 0$ in $(C^{(n+r)})^{m}$ as $k\to\infty$,
we get $f(\cdot,\varepsilon^{(k)})\to 0$ in $(C^{(n)})^{m}$ because the
boundary-value problem (2), (3) satisfies Limit Condition (I) as we have
proved on step 1. Since $1/2\leq\|q(\varepsilon^{(k)})\|_{\mathbb{C}^{rm}}\leq
1$ for every $k$, there exists a subsequence
$(q(\varepsilon^{(k_{p})}))_{p=1}^{\infty}\subset(q(\varepsilon^{(k)}))_{k=1}^{\infty}$
and a nonzero vector $q(0)\in\mathbb{C}^{rm}$ such that
$q(\varepsilon^{(k_{p})})\to q(0)$ in $\mathbb{C}^{rm}$ as $p\to\infty$.

Thus, for every integer $p\geq 1$, the function
$z(\cdot,\varepsilon^{(k_{p})})\in(C^{(n+r)})^{m}$ is a unique solution to the
boundary-value problem

$\displaystyle
L(\varepsilon^{(k_{p})})\,z(t,\varepsilon^{(k_{p})})=f(t,\varepsilon^{(k_{p})}),\quad
a\leq t\leq b,$ $\displaystyle
B(\varepsilon^{(k_{p})})\,z(\cdot,\varepsilon^{(k_{p})})=q(\varepsilon^{(k_{p})}).$

Since $f(\cdot,\varepsilon^{(k_{p})})\to 0$ in $(C^{(n)})^{m}$ and
$q(\varepsilon^{(k_{p})})\to q(0)$ in $\mathbb{C}^{rm}$ as $p\to\infty$, it
follows from condition $(\ast\ast)$ of Basic Definition that the function
$z(\cdot,\varepsilon^{(k_{p})})$ converges to the unique solution $z(\cdot,0)$
of the limiting boundary-value problem

$L(0)z(t,0)=0,\quad a\leq t\leq b,\qquad\mbox{and}\qquad B(0)z(\cdot,0)=q(0).$

This convergence holds in the space $(C^{(n+r)})^{m}$ as $\varepsilon\to 0+$.
But recall that $z(\cdot,\varepsilon^{(k_{p})})\to 0$ in the same space.
Therefore $z(\cdot,0)\equiv 0$, which contradicts the boundary condition
$B(0)z(\cdot,0)=q(0)$ with $q(0)\neq 0$. Thus, our assumption is wrong, and we
have thereby proved that $\|B(\varepsilon)\|=O(1)$ as $\varepsilon\to 0+$.

_Step $3$._ We can now prove that the boundary-value problem (2), (3)
satisfies Limit Condition (II). Owing to (15) and Step 2, there exist numbers
$\varkappa^{\prime}>0$ and $\varepsilon^{\prime}\in(0,\varepsilon_{1})$ such
that $\|(L(\varepsilon),B(\varepsilon))\|\leq\varkappa^{\prime}$ for each
$\varepsilon\in[0,\varepsilon^{\prime})$. Here,
$\|(L(\varepsilon),B(\varepsilon))\|$ is the norm of the bounded operator (5).
We choose a function $z\in(C^{(n+r)})^{m}$ arbitrarily and then put
$f(\cdot,\varepsilon):=L(\varepsilon)z$ and $q(\varepsilon):=B(\varepsilon)z$
for each $\varepsilon\in[0,\varepsilon^{\prime})$. Thus,

$z=(L(\varepsilon),B(\varepsilon))^{-1}(f(\cdot,\varepsilon),q(\varepsilon))\quad\mbox{for
every}\quad\varepsilon\in[0,\varepsilon^{\prime}).$

Here, of course, $(L(\varepsilon),B(\varepsilon))^{-1}$ denotes the inverse of
(5). (The operator (5) is invertible due to condition ($\ast$) of Basic
Definition.) Then, whenever $0<\varepsilon<\varepsilon^{\prime}$, we have

$\displaystyle\bigl{\|}B(\varepsilon)z-B(0)z\bigr{\|}_{\mathbb{C}^{rm}}\leq\bigl{\|}(f(\cdot,\varepsilon),q(\varepsilon))-(f(\cdot,0),q(0))\bigr{\|}_{(C^{n})^{m}\times\mathbb{C}^{rm}}$
$\displaystyle\qquad=\bigl{\|}(L(\varepsilon),B(\varepsilon))(L(\varepsilon),B(\varepsilon))^{-1}\bigl{(}(f(\cdot,\varepsilon),q(\varepsilon))-(f(\cdot,0),q(0))\bigr{)}\bigr{\|}_{(C^{(n)})^{m}\times\mathbb{C}^{rm}}$
$\displaystyle\qquad\leq\varkappa^{\prime}\,\bigl{\|}(L(\varepsilon),B(\varepsilon))^{-1}\bigl{(}(f(\cdot,\varepsilon),q(\varepsilon))-(f(\cdot,0),q(0))\bigr{)}\bigr{\|}_{(n+r)}$
$\displaystyle\qquad=\varkappa^{\prime}\,\bigl{\|}(L(0),B(0))^{-1}(f(\cdot,0),q(0))-(L(\varepsilon),B(\varepsilon))^{-1}(f(\cdot,0),q(0))\bigr{\|}_{(n+r)}\to
0$

as $\varepsilon\to 0+$ due to condition ($\ast\ast$) of Basic Definition.
Since $y\in(C^{(n+r)})^{m}$ is arbitrary, we have proved that the boundary-
value problem (2), (3) satisfies Limit Condition (II). ∎

###### Proof of Estimation Theorem.

We will first prove the left-hand part of (7). Limit Conditions (I) and (II)
yield the strong convergence of $(L(\varepsilon),B(\varepsilon))$ to
$(L(0),B(0))$ as $\varepsilon\to 0+$. Recall that we consider
$(L(\varepsilon),B(\varepsilon))$, with $0\leq\varepsilon<\varepsilon_{0}$, as
a bounded operator from $(C^{(n+r)})^{m}$ to
$(C^{(n)})^{m}\times\mathbb{C}^{rm}$. Hence, there exist numbers
$\varkappa^{\prime}>0$ and $\varepsilon^{\prime}\in(0,\varepsilon_{1})$ such
that the operator norm
$\|(L(\varepsilon),B(\varepsilon))\|\leq\varkappa^{\prime}$ for each
$\varepsilon\in[0,\varepsilon^{\prime})$. Indeed, if the contrary were true,
we would obtain a sequence of positive numbers
$(\varepsilon^{(k)})_{k=1}^{\infty}$ such that $\varepsilon^{(k)}\to 0$ and
$\|(L(\varepsilon^{(k)}),B(\varepsilon^{(k)}))\|\to\infty$ as $k\to\infty$,
which would contradict the above strong convergence in view of the Banach-
Steinhaus theorem. We see now that the left-hand part of the two-sided
estimate (7) holds for each $\varepsilon\in(0,\varepsilon^{\prime})$, with
$\varkappa_{1}:=1/\varkappa^{\prime}$.

Let us prove the right-hand part of this estimate. According to Main Theorem,
the boundary-value problem (2), (3) satisfies Basic Definition. Therefore
operator (5) is invertible for every $\varepsilon\in[0,\varepsilon_{1})$, and,
moreover, its inverse $(L(\varepsilon),B(\varepsilon))^{-1}$ converges
strongly to $(L(0),B(0))^{-1}$ as $\varepsilon\to 0+$. Indeed, for arbitrary
$f\in(C^{n})^{m}$ and $q\in\mathbb{C}^{rm}$, we conclude by condition
$(\ast\ast)$ of Basic Definition that

$(L(\varepsilon),B(\varepsilon))^{-1}(f,q)=:z(\cdot,\varepsilon)\to
z(\cdot,0):=(L(0),B(0))^{-1}(f,q)$

in $(C^{(n+r)})^{m}$ as $\varepsilon\to 0+$. Hence, there exists positive
numbers $\varepsilon_{2}<\varepsilon^{\prime}$ and $\varkappa_{2}$ such that
the norm of the inverse operator
$\|(L(\varepsilon),B(\varepsilon))^{-1}\|\leq\varkappa_{2}$ for each
$\varepsilon\in[0,\varepsilon_{2})$. This follows from the Banach-Steinhaus
theorem in a way analogous to that used above. We now see that the right-hand
part of (7) is true for $\varepsilon\in(0,\varepsilon_{2})$. ∎

##  4\. Application

Let us give an application of Main Theorem to multipoint boundary-value
problems. For the parameter-dependent differential system (2), we pose the
multipoint boundary condition

$B(\varepsilon)z(\cdot,\varepsilon)\equiv\sum_{j=0}^{p}\sum_{k=1}^{\omega_{j}}\sum_{l=0}^{n+r}{\alpha_{j,k}^{(l)}(\varepsilon)z^{(l)}(t_{j,k}(\varepsilon),\varepsilon)}=q(\varepsilon).$
(16)

Here, the integers $p\geq 1$ and $\omega_{j}\geq 1$, matrixes
$\alpha_{j,k}^{(l)}(\varepsilon)\in\mathbb{C}^{rm\times m}$, and points
$t_{j,k}(\varepsilon)\in[a,b]$ are arbitrarily chosen for all admissible
values of $j$, $k$, $l$, and $\varepsilon$. As above,
$q(\varepsilon)\in\mathbb{C}^{rm}$. In (16), we use the summation over the two
indexes $j$ and $k$ instead of one index because the points
$t_{j,k}(\varepsilon)$, with $j\geq 1$, will be supposed to converge to a
certain point $t_{j}$ as $\varepsilon\to 0+$, whereas the points
$t_{0,k}(\varepsilon)$ will not be supposed to have a limit as $\varepsilon\to
0+$. So, we divide the set of all points $t_{j,k}(\varepsilon)$ into $p+1$
sets $\\{t_{j,k}(\varepsilon):k=1,\ldots,\omega_{j}\\}$, with
$j=0,1,\ldots,p$, in accordance with their behavior as $\varepsilon\to 0+$.

In this connection, we consider the following limiting boundary-value problem:

$\displaystyle L(0)z(t,0)=f(\cdot,0),\quad a\leq t\leq b,$ (17) $\displaystyle
B(0)z(\cdot,0)\equiv\sum_{j=1}^{p}\sum_{l=0}^{n+r}{\alpha_{j}^{(l)}z^{(l)}(t_{j},0)}=q(0).$
(18)

Here, $\alpha_{j}^{(l)}\in\mathbb{C}^{rm\times m}$, $t_{j}\in[a,b]$ for all
admissible values of $j$ and $l$, and as above, $q(0)\in\mathbb{C}^{rm}$.

It is evident that, for every $\varepsilon\in[0,\varepsilon_{0})$, the linear
mapping $z(\cdot,\varepsilon)\mapsto B(\varepsilon)z(\cdot,\varepsilon)$ is a
continuous operator from $(C^{(n+r)})^{m}$ to $\mathbb{C}^{rm}$. Thus, the
multipoint boundary-value problem (2), (16) is generic with respect to the
space $(C^{(n+r)})^{m}$ for every $\varepsilon\in(0,\varepsilon_{0})$ as well
as the the limiting problem (17), (18).

For these boundary-value problems, we can give constructive sufficient
conditions for the convergence $B(\varepsilon)z\to B(0)z$ in $\mathbb{C}^{rm}$
for every $z\in(C^{(n+r)})^{m}$. Note that it is scarcely possible to use the
system of conditions (2a)–(2d) for the verification of this convergence
because it is difficult in practice to find the function
$\Phi(\cdot,\varepsilon)$ in the canonical representation (8) of the
multipoint boundary operator $B(\varepsilon)$.

Theorem A1. Suppose that the boundary-value problem (2), (16) satisfies the
following conditions as $\varepsilon\to 0+$:

  * (d1)

$t_{j,k}(\varepsilon)\to t_{j}$ for all $j\in\\{1,\ldots,p\\}$ and
$k\in\\{1,\ldots,\omega_{j}\\}$;

  * (d2)

$\sum_{k=1}^{\omega_{j}}\alpha_{j,k}^{(l)}(\varepsilon)\to\alpha_{j}^{(l)}$
for all $j\in\\{1,\ldots,p\\}$ and $l\in\\{0,\ldots,n+r\\}$;

  * (d3)

$\alpha_{j,k}^{(n+r)}(\varepsilon)=O(1)$ for all $j\in\\{1,\ldots,p\\}$ and
$k\in\\{1,\ldots,\omega_{j}\\}$;

  * (d4)

$\|\alpha_{j,k}^{(l)}(\varepsilon)\|_{\mathbb{C}^{rm\times
m}}\\!\cdot\\!|t_{j,k}(\varepsilon)-t_{j}|\to 0$ for all
$j\in\\{1,\ldots,p\\}$, $k\in\\{1,\ldots,\omega_{j}\\}$, and
$l\in\\{0,\ldots,n+r-1\\}$;

  * (d5)

$\alpha_{0,k}^{(l)}(\varepsilon)\to 0$ for all $k\in\\{1,\ldots,\omega_{0}\\}$
and $l\in\\{0,\cdots,n+r\\}$.

Then this problem satisfies Limit Condition (II).

Owing to Main Theorem and Theorem A1, we immediately obtain the following
result:

Theorem A2. Assume that the multipoint boundary-value problem (2), (16)
satisfies Limit Condition (I) and conditions (d1)–(d5) and that the limiting
problem (17), (18) satisfies Condition $(0)$. Then the solution to the problem
(2), (16) is continuous in the parameter $\varepsilon$ at $\varepsilon=0$.

###### Proof of Theorem A.1.

Let $\|\cdot\|$ denote the norm of a number vector or matrix, this norm being
the sum of the absolute values of all their elements. For an arbitrary
function $z\in(C^{(n+r)})^{m}$ and sufficiently small $\varepsilon>0$, we have
the inequality

$\displaystyle\|B(\varepsilon)z-B(0)z\|\leq$
$\displaystyle\sum_{k=1}^{\omega_{0}}\sum_{l=0}^{n+r}\|\alpha_{0,k}^{(l)}(\varepsilon)\|\cdot\|z^{(l)}(t_{0,k}(\varepsilon))\|$
(19)
$\displaystyle+\sum_{j=1}^{p}\sum_{l=0}^{n+r}\,\left\|\sum_{k=1}^{\omega_{j}}\alpha_{j,k}^{(l)}(\varepsilon)z^{(l)}(t_{j,k}(\varepsilon))-\alpha_{j}^{(l)}z^{(l)}(t_{j})\right\|.$

Here, by condition (d5), we can write

$\begin{gathered}\|\alpha_{0,k}^{(l)}(\varepsilon)\|\cdot\|z^{(l)}(t_{0,k}(\varepsilon))\|\leq\|\alpha_{0,k}^{(l)}(\varepsilon)\|\cdot\|z\|_{(n+r)}\to
0\end{gathered}$ (20)

for all admissible values of $k$ and $l$. This and the other limits in the
proof are considered as $\varepsilon\to 0+$. Besides,

$\displaystyle\left\|\sum_{k=1}^{\omega_{j}}\alpha_{j,k}^{(l)}(\varepsilon)z^{(l)}(t_{j,k}(\varepsilon))-\alpha_{j}^{(l)}z^{(l)}(t_{j})\right\|$
(21)
$\displaystyle\leq\left\|\sum_{k=1}^{\omega_{j}}\alpha_{j,k}^{(l)}(\varepsilon)\cdot\left(z^{(l)}(t_{j,k}(\varepsilon))-z^{(l)}(t_{j})\right)\right\|+\left\|\left(\sum_{k=1}^{\omega_{j}}\alpha_{j,k}^{(l)}(\varepsilon)-\alpha_{j}^{(l)}\right)\cdot
z^{(l)}(t_{j})\right\|$
$\displaystyle\leq\sum_{k=1}^{\omega_{j}}\|\alpha_{j,k}^{(l)}(\varepsilon)\|\cdot\|z^{(l)}(t_{j,k}(\varepsilon))-z^{(l)}(t_{j})\|+\left\|\sum_{k=1}^{\omega_{j}}\alpha_{j,k}^{(l)}(\varepsilon)-\alpha_{j}^{(l)}\right\|\cdot\|z\|_{(n+r)}.$

Here, according to condition (d2), we have

$\left\|\sum_{k=1}^{\omega_{j}}\alpha_{j,k}^{(l)}(\varepsilon)-\alpha_{j}^{(l)}\right\|\cdot\|z\|_{(n+r)}\to
0.$ (22)

Moreover,

$\sum_{k=1}^{\omega_{j}}\|\alpha_{j,k}^{(l)}(\varepsilon)\|\cdot\|z^{(l)}(t_{j,k}(\varepsilon))-z^{(l)}(t_{j})\,\|\to
0.$ (23)

Indeed, if $l=n+r$, then formula (23) is a direct consequence of conditions
(d1), (d3) and the continuity of $z^{(l)}$. If $l\leq n+r-1$, then this
formula is a consequence of the Lagrange theorem of the mean and condition
(d4), namely

$\sum_{k=1}^{\omega_{j}}\|\alpha_{j,k}^{(l)}(\varepsilon)\|\cdot\|z^{(l)}(t_{j,k}(\varepsilon))-z^{(l)}(t_{j})\,\|\leq\|z\|_{(n+r)}\sum_{k=1}^{\omega_{j}}\|\alpha_{j,k}^{(l)}(\varepsilon)\|\cdot|t_{j,k}(\varepsilon)-t_{j}|\to
0.$

According to (21)–(23), we have the convergence

$\left\|\sum_{k=1}^{\omega_{j}}\alpha_{j,k}^{(l)}(\varepsilon)z^{(l)}(t_{j,k}(\varepsilon))-\alpha_{j}^{(l)}z^{(l)}(t_{j})\right\|\to
0.$ (24)

Now formulas (19), (20), and (24) imply that $\|B(\varepsilon)z-B(0)z\|\to 0$,
i.e. Limit Condition (II) is satisfied. ∎

##  5\. Concluding remarks

It is easy to check that conditions $(\ast)$ and $(\ast\ast)$ in Basic
Definition are separately equivalent to the corresponding ones:

  * $(\star)$

there exists a positive number $\varepsilon_{1}<\varepsilon_{0}$ such that the
operator (5) is invertible for every $\varepsilon\in[0,\varepsilon_{1})$;

  * $(\star\star)$

the operator $(L(\varepsilon),B(\varepsilon))^{-1}$, which is inverse of (5),
converges to $(L(0),B(0))^{-1}$ as $\varepsilon\to\infty$ in the strong
operator topology.

It is also not difficult to prove that the system of Limit Conditions (I) and
(II) is equivalent to

Limit Condition (IV). The operator $(L(\varepsilon),B(\varepsilon))$ converges
to $(L(0),B(0))$ as $\varepsilon\to\infty$ in the strong operator topology.

So, Main Theorem implies that, the equivalence

$\bigl{(}\mbox{Limit Condition
(IV)}\bigr{)}\;\Leftrightarrow\;\bigl{(}(\star)\wedge(\star\star)\bigr{)}$

holds true under Condition (0).

This is a surprising result at first sight. Indeed, if $X$ and $Y$ are
arbitrary infinite-dimensional Banach spaces with Schauder basis, then the set
of all irreversible operators is sequentially dense in $\mathcal{L}(X,Y)$ in
the strong operator topology; here, $\mathcal{L}(X,Y)$ stands for the linear
space of all continuous linear operators from $X$ to $Y$. Moreover, the
mapping $\mathrm{inv}:T\mapsto T^{-1}$ given on the set of all invertible
operators $T\in\mathcal{L}(X,Y)$ is everywhere discontinuous in this topology.
We will prove these propositions in Appendix.

## Appendix

Let $X$ and $Y$ are complex or real Banach spaces. Assume that at least one of
the spaces $X$ and $Y$ has Schauder basis. Specifically, we can put
$X:=(C^{(n+r)})^{m}$ and $Y:=(C^{(n)})^{m}\times\mathbb{C}^{rm}$.

###### Proposition 1

The set of all finite-dimensional operators from $\mathcal{L}(X,Y)$ is
sequentially dense in $\mathcal{L}(X,Y)$ in the strong operator topology.

Since every finite-dimensional operator from $\mathcal{L}(X,Y)$ is
irreversible, this proposition implies that the set of all irreversible
operators from $\mathcal{L}(X,Y)$ is sequentially dense in $\mathcal{L}(X,Y)$
in the strong operator topology.

###### Proof of Proposition 1.

We apart examine the cases where $X$ or $Y$ has Schauder basis.

Consider first the case where $X$ has a certain Schauder basis
$(e_{k})_{k=1}^{\infty}$. Then every vector $x\in X$ admits the unique
representation

$x=\sum_{k=1}^{\infty}e_{k}^{\ast}(x)\,e_{k},$ (25)

where each $e_{k}^{\ast}$ is a certain continuous linear functional on $X$.
Using this basis, we define the operator of the $n$-th partial sum

$S_{n}:x\mapsto\sum_{k=1}^{n}e_{k}^{\ast}(x)\,e_{k},\quad\mbox{with}\quad x\in
X.$

It is evident that $S_{n}\rightarrow I_{X}$ as $n\to\infty$ in the strong
operator topology in $\mathcal{L}(X,X)$; here, as usual, $I_{X}$ denotes the
identical operator on $X$. Now, given an operator $T\in\mathcal{L}(X,Y)$, we
put $T_{n}:=TS_{n}$ for each integer $n\geq 1$. Then every operator $T_{n}$ is
finite-dimensional, and $T_{n}\to T$ as $n\to\infty$ in the strong operator
topology in $\mathcal{L}(X,Y)$.

The case where $Y$ has Schauder basis is examined similarly. Namely, let
$(\tilde{e}_{k})_{k=1}^{\infty}$ be a Schauder basis of $Y$, and let
$\tilde{e}_{k}^{\ast}(y)$ be the $k$-th coordinate of $y\in Y$ in this basis.
For every integer $n\geq 1$, we define the operator

$\widetilde{S}_{n}:y\mapsto\sum_{k=1}^{n}\tilde{e}_{k}^{\ast}(y)\,\tilde{e}_{k},\quad\mbox{with}\quad
y\in Y.$

Given an operator $T\in\mathcal{L}(X,Y)$, we put
$\widetilde{T}_{n}:=\widetilde{S}_{n}T$ and conclude that $\widetilde{T}_{n}$
is finite-dimensional and that $\widetilde{T}_{n}\to T$ as $n\to\infty$ in the
strong operator topology in $\mathcal{L}(X,Y)$. ∎

###### Proposition 2

Assume that the spaces $X$ and $Y$ are isomorphic. Then the mapping
$\mathrm{inv}:T\mapsto T^{-1}$ given on the set of all invertible operators
$T\in\mathcal{L}(X,Y)$ is everywhere discontinuous in the strong operator
topology.

The assumption made in this proposition is natural; indeed, if $X$ and $Y$ are
not isomorphic, then there are not any invertible operators in
$\mathcal{L}(X,Y)$.

###### Proof of Proposition 2.

Owing to the Banach–Steinhaus theorem, it is sufficient to prove that, for an
arbitrary invertible operator $T\in\mathcal{L}(X,Y)$, there exists a sequence
of invertible operators $(T_{n})_{n=1}^{\infty}\subset\mathcal{L}(X,Y)$ such
that $T_{n}\to T$ as $n\to\infty$ in the strong operator topology in
$\mathcal{L}(X,Y)$ and that the sequence $(\|T_{n}^{-1}\|)_{n=1}^{\infty}$ is
not bounded. Here, of course, $\|T_{n}^{-1}\|$ stands for the norm of the
inverse operator $T_{n}^{-1}:Y\to X$. Without loss of generality, we can
assume that $Y=X$ and $T=I_{X}$. Indeed, since the Banach spaces $X$ and $Y$
are isomorphic, the general case is trivially reduced to the case assumed.

Let us use the reasoning given in the proof of Proposition 1. For every
integer $n\geq 1$, we consider the operator

$I_{n}:x\mapsto\sum_{k=1}^{n-1}e_{k}^{\ast}(x)\,e_{k}+\frac{1}{n}\,e_{n}^{\ast}(x)\,e_{n}+\sum_{k=n+1}^{\infty}e_{k}^{\ast}(x)\,e_{k},\quad\mbox{with}\quad
x\in X.$ (26)

Every operator $I_{n}$ belongs to $\mathcal{L}(X,Y)$ and is invertible;
moreover,

$I_{n}^{-1}:x\mapsto\sum_{k=1}^{n-1}e_{k}^{\ast}(x)\,e_{k}+n\,e_{n}^{\ast}(x)\,e_{n}+\sum_{k=n+1}^{\infty}e_{k}^{\ast}(x)\,e_{k},\quad\mbox{with}\quad
x\in X.$ (27)

It follows from (25) and (26) that $I_{n}\to I_{X}$ as $n\to\infty$ in the
strong operator topology in $\mathcal{L}(X,Y)$. However, the sequence
$(\|I_{n}^{-1}\|)_{n=1}^{\infty}$ is not bounded because

$\|I_{n}^{-1}\|\geq\frac{\|I_{n}^{-1}e_{n}\|}{\|e_{n}\|_{X}}=n\quad\mbox{for
each integer}\quad n\geq 1$

due to (27). ∎

Note that we cannot replace the strong operator topology with the uniform one
in Propositions 1 and 2. As is known, the mapping $\mathrm{inv}$ is everywhere
continuous in the uniform operator topology, and the closure of the set of all
finite-dimensional operators in this topology contains only compact (hence,
irreversible) operators.

## References

  * [1] Malkhaz Ashordia, _Criteria of correctness of linear boundary value problems for systems of generalized ordinary differential equations_ , Czechoslovak Math. J. 46 (1996), no. 3, 385–404.
  * [2] Nelson Dunford and Jacob T. Schwartz, _Linear operators. part I: General theory_ , Interscience, New York, 1958.
  * [3] I. I. Gihman, _Concerning a theorem of N. N. Bogolyubov_ , Ukr. Mat. Zh. 4 (1952), 215–219 (Russian).
  * [4] E. V. Gnyp, T. I. Kodlyuk, and V. A. Mikhailets, _Fredholm boundary-value problems with parameter in Sobolev spaces_ , Ukrainian Math. J. 67 (2015), no. 5, 658–667.
  * [5] Andrii Goriunov, Vladimir Mikhailets, and Konstantin Pankrashkin, _Formally self-adjoint quasi-differential operators and boundary-value problems_ , Electron. J. Differential Equations (2013), no. 101, 1–16.
  * [6] A. S. Goryunov and V. A. Mikhaĭlets, _Resolvent convergence of sturm-liouville operators with singular potentials_ , Math. Notes 87 (2010), no. 1-2, 287–292.
  * [7] A. S. Goryunov and V. A. Mikhailets, _Regularization of two-term differential equations with singular coefficients by quasiderivatives_ , Ukrainian Math. J. 63 (2012), no. 9, 1361–1378.
  * [8] E. Hnyp, V. Mikhailets, and A. Murach, _A criterion for continuity in a parameter of solutions to generic boundary-value problems for differential systems in Sobolev spaces (to appear)_.
  * [9] Lars Hörmander, _The analysis of linear partial differential operators. III: Pseudo-differential operators_ , Springer-Verlag, Berlin, 1985\.
  * [10] I. T. Kiguradze, _Some singular boundary value problems for ordinary differential equations_ , Tbilisi University, Tbilisi, 1975 (Russian).
  * [11] I. T. Kiguradze, _Boundary value problems for systems of ordinary differential equations_ , J. Soviet Math. 43 (1988), no. 2, 2259–2339.
  * [12] I. T. Kiguradze, _On boundary value problems for linear differential systems with singularities_ , Differential Equations 39 (2003), no. 2, 212–225.
  * [13] T. I. Kodliuk and V. A. Mikhailets, _Multipoint boundary-value problems with parameter in Sobolev spaces_ , Dopov. Nats. Akad. Nauk Ukr. Mat. Prirodozn. Tekh. Nauki (2012), no. 11, 15–19 (Russian).
  * [14] Tatiana I. Kodliuk and Vladimir A. Mikhailets, _Solutions of one-dimensional boundary-value problems with a parameter in Sobolev spaces_ , J. Math. Sci. (N. Y.) 190 (2013), no. 4, 589–599.
  * [15] T. I. Kodlyuk, V. A. Mikhailets, and N. V. Reva, _Limit theorems for one-dimensional boundary-value problems_ , Ukrainian Math. J. 65 (2013), no. 1, 77–90.
  * [16] M. A. Krasnoselskiĭ and S. G. Kreĭn, _On the principle of averaging in nonlinear mechanics_ , Uspehi Mat. Nauk (N.S.) 10 (1955), no. 3, 147–152 (Russian).
  * [17] Yaroslav Kurcveĭl and Zdeněk Vorel, _Continuous dependence of solutions of differential equations on a parameter_ , Czechoslovak Math. J. 7 (82) (1957), 568–583 (Russian).
  * [18] A. Ju. Levin, _Passage to the limit for nonsingular systems $\dot{X}=A_{n}(t)X.$_, Dokl. Akad. Nauk SSSR 176 (1967), 774–777 (Russian).
  * [19] V. A. Mikhailets and G. A. Chekhanova, _Fredholm boundary-value problems with parameter on the spaces $C^{(n)}[a;b]$_, Dopov. Nats. Akad. Nauk Ukr. Mat. Prirodozn. Tekh. Nauki. (2014), no. 7, 24–28 (Russian).
  * [20] V. A. Mikhailets and G. A. Chekhanova, _Limit theorems for general one-dimensional boundary-value problems_ , J. Math. Sci. (N. Y.) 204 (2015), no. 3, 333–342.
  * [21] V. A. Mikhailets, A. A. Murach, and V. Soldatov, _Continuity in a parameter of solutions to generic boundary-value problems_ , Electron. J. Qual. Theory Differ. Equ. (2016), no. 87, 1–16.
  * [22] V. A. Mikhaĭlets and N. V. Reva, _Passage to the limit in systems of linear differential equations_ , Dopov. Nats. Akad. Nauk Ukr. Mat. Prirodozn. Tekh. Nauki (2008), no. 8, 28–30 (Russian).
  * [23] Tkhe Khoan Nguen, _Dependence of the solutions of a linear system of differential equations on a parameter_ , Differential Equations 29 (1993), no. 6, 830–835.
  * [24] Zdzislaw Opial, _Continuous parameter dependence in linear systems of differential equations_ , J. Differential Equations 3 (1967), 571–579.
  * [25] William T. Reid, _Some limit theorems for ordinary differential systems_ , J. Differential Equations 3 (1967), 423–439.
  * [26] Frigyes Riesz and Béla Sz.-Nagy, _Functional analysis_ , Frederick Ungar Publishing Co., New York, 1955.
  * [27] V. O. Soldatov, _On the continuity in a parameter for the solutions of boundary-value problems total with respect to the spaces $C^{(n+r)}[a,b]$_, Ukrainian Math. J. 67 (2015), no. 5, 785–794.


# Entanglement as a sufficient condition for quantum synchronization between
two mechanical oscillators

Manju<EMAIL_ADDRESS>Shubhrangshu Dasgupta, and Asoka Biswas
Department of Physics, Indian Institute of Technology Ropar, Rupnagar, Punjab
140001, India

###### Abstract

We present an optomechanical model to show that entanglement can be a
sufficient condition for quantum synchronization of two mechanical
oscillators. As both these entities can be characterized in terms of variances
of a set of EPR-like conjugate quadratures, we investigate whether this leads
to a specific condition for simultaneous existence of the both. In our model,
one of the oscillators makes the cavity, while the other is kept suspended
inside the cavity, and the always-on coupling between the two is mediated via
the same cavity mode. We show that in presence of amplitude modulation with
the same frequency as that of the oscillators, these oscillators get nearly
complete quantum synchronized and entangled simultaneously in the steady
state. We also show that entanglement always becomes accompanied by quantum
synchronization, though the reverse is not necessarily true. Thus,
entanglement becomes a sufficient condition for the quantum synchronization.
This behaviour can be observed for a large range of system parameters.

##  I INTRODUCTION

Spontaneous synchronization is a natural phenomenon that can be often
experienced in avian flight and flashing of fireflies strogatz . It was also
observed by Huygens in a classical clock pendulum with a common support
huygens1897oeuvres . It has been explored in many different areas, namely,
neuron networks leone2015synchronization ; tang2011synchronization ;
cao2017fixed , chemical reaction vaidyanathan2015dynamics , nonlinear dynamics
lee2013quantum1 ; walter2014quantum ; shirasaka2017optimizing , and electrical
circuits and communications van1927frequency ; pecora1990synchronization .

Classical synchronization between two nonlinearly coupled oscillators occurs
due to a rephasing and energy redistribution of their motions. This refers to
a limit cycle in phase space, where their respective generalized positions
$q_{j}(t)$ and generalized linear momenta $p_{j}(t)$ ($j\in 1,2$) become
equal. In quantum regime, however, the one cannot measure these quadrataures
with certainty at the same time $t$, according to Heisenberg’s uncertainty
principle. The best estimate of the quantum synchronization (QS) therefore
corresponds to a state of two oscillators, for which the total uncertainty of
their joint quadratures becomes minimum. In this regard, Mari et al. proposed
the following figure of merit mari2013measures :

$\displaystyle S_{qm}=\left\langle
q_{-}^{2}(t)+p_{-}^{2}(t)\right\rangle^{-1}$ (1)

in terms of the synchronization errors
$q_{-}(t)=\frac{1}{\sqrt{2}}\left[q_{1}(t)-q_{2}(t)\right]$ and
$p_{-}(t)=\frac{1}{\sqrt{2}}\left[p_{1}(t)-p_{2}(t)\right]$, where
$\langle\cdots\rangle$ denotes the expectation (or mean) values of the
corresponding operators. Expanding these operators around their mean values,
$q_{-}(t)=Q_{-}(t)+\delta q_{-}(t)$, $p_{-}(t)=P_{-}(t)+\delta p_{-}(t)$, and
using the limit of zero mean of these quadrature differences, a modified form
of the measure of QS can be obtained as

$S_{q}=\frac{1}{\left\langle\delta q_{-}(t)^{2}+\delta
p_{-}(t)^{2}\right\rangle}\leq 1\;.\\\ $

On the other hand, quantum correlation between the interacting systems can
also manifest itself as entanglement. In the context of two oscillators (two
bosonic modes), entanglement has been characterized in terms of uncertainties
of the joint quadratures. The criterion for entanglement can be derived
starting from the uncertainty relation of $q_{-}$ and $p_{-}$ and then
applying the Peres’s separability criterion peres1996separability  based on
the partial transposition. Any bimodal state can be said to be entangled if it
violates the partially transposed uncertainty relation, as given by
duan2000inseparability

$\langle(\delta q_{-})^{2}\rangle+\langle(\delta p_{+})^{2}\rangle\geq 1\;,$
(2)

where the partial transposition on the second subsystem changes $p_{-}$ into
$p_{+}=\frac{1}{\sqrt{2}}\left[p_{1}+p_{2}\right]$. A stronger version of the
above inequality, as follows, can be obtained by using a simple algebraic
identity mancini2002entangling ; biswas2005NJP :

$E_{D}=\langle(\delta q_{-})^{2}\rangle\langle(\delta
p_{+})^{2}\rangle\geq\frac{1}{4}\;.$ (3)

It can be observed from the above, both the QS measure and the entanglement
criterion are derived from the Heisenberg uncertainty principle for a pair of
EPR-like variables. This suggests that both these seemingly different features
originate purely from the same class of quantum correlation (namely, the EPR-
type) and therefore can appear in the same system at a similar time-scale. In
fact, in our recent work garg2022quantum , we have shown that it is indeed
possible to both quantum-synchronize and entangle two oscillators in a certain
parameter regime, by coupling them with two different cavities. In this paper,
we will investigate if such possibility is generic in nature. With an always-
on nonlinear coupling between two oscillators, we will show that they remain
quantum-synchronized as long as entangled. More importantly, without the
entanglement, they loose synchronization. Thus, the entanglement can be
considered as a sufficient condition for achieving QS.

We note that whether two coupled systems can be simultaneously synchronized
and entangled has been investigated in several systems. Manzano et al. showed
that any two coupled oscillators in a network, initially prepared in a
separable state, can be both synchronized and entangled
manzano2013synchronization , for a suitable choice of the interaction strength
to the other oscillators in the network. They further showed that in presence
of synchronization, the entanglement is retained despite the decoherence. It
was further shown in witthaut2017classical ; lee2013quantum  that many-body
system can be entangled in presence of synchronization. In the classical
(quantum) limit, two single-mode cavities get synchronized (entangled).
However, all the above works dealt only with classical synchronization. On the
contrary, it is more meaningful to focus on QS, when one looks for its
relationship with entanglement, which has no classical analogy. In this paper,
we specifically investigate if there is any inherent quantum correlation that
can lead to both QS and entanglement. Previous attempts in this regard include
study of quantum discord giorgi2012 ; manzano2013synchronization , linear
entropy chimera , mutual information based on von Neumann entropy
ameri2015mutual  and Renyi entropy chimera2015 . In these works, similar
steady state behavior of these measures of quantum correlation and the QS was
reported. But their numerical studies could not essentially relate them as
arising from the same physical origin. The present work addresses this exact
loophole of the earlier reports.

To do the relevant analysis, optomechanical system poses as a suitable
platform. In our model, two mechanical oscillators indirectly interact with
each other, via their common coupling to the same cavity mode. One of the
oscillators is suspended inside the cavity, while the other is usual movable
mirror of the same cavity. Clearly, their interaction is always on and weakly
nonlinear, but still leads to near-complete QS and entanglement between them.
We will specifically show that the quantum correlation manifests itself in the
simultaneous existence of the two. This happens at the modulation frequency of
the driving field. It is recently reported that the continuous measurement of
a single bosonic system can enhance its phase synchronization with the
external drive kwek , when it is subjected to squeezed driving field, and both
negative and nonlinear damping. A feedback control of the system based on the
continuous measurement of an additional bath can also enhance synchronization
nakao . In this paper, we however consider a coupled system (not a single
system entrained with a drive) and our result suggests an opposite effect upon
measurement. If one of the oscillators is measured, it destroys the
entanglement and hence the complete QS (note that we do not study the phase
synchronization in this paper). In addition, the dynamics of the oscillators
does not involve any two-photon processes, either driving or damping, unlike
the van der Pol oscillators considered in kwek ; nakao .

We organize the paper in the following way. In the following Section, we
introduce the model and the Hamiltonian, and derive all the equations of
motion for the fluctuation operators. We obtain the analytic expressions of
the fluctuations and its spectrum in Sec. III. In Sec. IV, we numerically
solve the dynamics of the system and will study how both entanglement and QS
between the oscillators can arise at similar time-scale. Finally, we conclude
the paper in Sec. V.

Figure 1: Schematic diagram of a membrane-in-middle optomechanical system
(MIMOS), in which a membrane is suspended in the middle of a driven cavity.

##  II THEORETICAL MODEL

The optomechanical system under consideration is shown in Fig. 1. This system
consists of an optical cavity with one mirror fixed, while the other is
movable, acting as a mechanical oscillator (with resonance frequency
$\omega_{m1}$). Another oscillating membrane is also suspended inside the
cavity. The cavity is driven by a laser with frequency $\omega_{l}$ and time-
modulated amplitude $E(t)$. This model is quite different from the other
optomechanical setups in which two membranes (or mirrors) are optically
coupled to two different cavities, while an additional mechanical coupling
between them has to be externally introduced to facilitate the
synchronization. In our model, the coupling between the mirror and the
membrane is always on and there is no need to introduce any additional
coupling. Their interaction with the same cavity mode generates an effective
photon-number-dependent mirror-membrane coupling, which leads to their
synchronization and entanglement as well.

The coupling between the optical cavity and membrane depends upon the position
of membrane, placed relative to the nodes and antinodes of cavity mode. An
one-dimensional calculation gives the frequency of the cavity as a function of
the mirror displacement $q_{1}$ and the membrane displacement $q_{2}$, as
follows thompson2008strong :

$\displaystyle\centering\omega_{cav}(q_{1},{q_{2}})=\left(\frac{c}{L+q_{1}}\right)\cos^{-1}\left(r_{c}\cos{\frac{4\pi
q_{2}}{\lambda}}\right)\;.\@add@centering$ (4)

Here $L$ is the length of cavity, $r_{c}$ is the reflectivity of the membrane
and $\lambda$ is the wavelength of laser field. Expanding
$\omega_{cav}(q_{1},q_{2})$ about the equilibrium positions $q_{10}=0$ and
$q_{20}$ of the respective cavity, up to second order in Taylor series, we
have:

$\centering\omega_{cav}(q_{1},q_{2})=\omega_{c}-g_{1}^{(1)}q_{1}-g_{1}^{(2)}q_{2}+g_{2}^{(1)}q_{1}^{2}+g_{2}^{(2)}q_{2}^{2}-g_{3}q_{1}q_{2}\;,\@add@centering$
(5)

where we can identify

$\omega_{c}=\omega_{cav}(0,q_{20})-\left.\frac{\partial\omega_{cav}}{\partial
q_{2}}\right|_{(0,q_{20})}q_{20}+\frac{1}{2}\left.\frac{\partial^{2}\omega_{cav}}{\partial
q_{2}^{2}}\right|_{(0,q_{20})}q_{20}^{2}\;,$

and

$\displaystyle g_{1}^{(1)}$ $\displaystyle=$
$\displaystyle-\left.\frac{\partial\omega_{cav}}{\partial
q_{1}}\right|_{(0,q_{20})}+\frac{1}{2}\left.\frac{\partial^{2}\omega_{cav}}{\partial
q_{1}\partial q_{2}}\right|_{(0,q_{20})}q_{20}\;,$ $\displaystyle g_{1}^{(2)}$
$\displaystyle=$ $\displaystyle-\left.\frac{\partial\omega_{cav}}{\partial
q_{2}}\right|_{(0,q_{20})}+\left.\frac{\partial^{2}\omega_{cav}}{\partial
q_{2}^{2}}\right|_{(0,q_{20})}q_{20}\;,$ $\displaystyle g_{2}^{(j)}$
$\displaystyle=$
$\displaystyle\frac{1}{2}\left.\frac{\partial^{2}\omega_{cav}}{\partial
q_{j}^{2}}\right|_{(0,q_{20})}\;\;,\;\;j=1,2$ $\displaystyle g_{3}$
$\displaystyle=$
$\displaystyle-\frac{1}{2}\left.\frac{\partial^{2}\omega_{cav}}{\partial
q_{1}\partial q_{2}}\right|_{(0,q_{20})}\;.$ (6)

The Hamiltonian of the system then takes the following form ($\hbar=1$):

$H=\omega_{cav}(q_{1},q_{2})a^{\dagger}a+\sum_{j=1,2}\frac{\omega_{mj}}{2}\left(q_{j}^{2}+p_{j}^{2}\right)+\iota
E\left[1+\eta_{D}\cos(\Omega_{D}t)\right]\left(a^{\dagger}-a\right)\;,$ (7)

where $\omega_{mj}$ $q_{j}$, and $p_{j}$ are the frequency, dimensionless
position and momentum operators, respectively, of the $j$th oscillator,
$a^{\dagger}(a)$ is the creation (annihilation) operator of cavity mode,
satisfying the commutation relation $\left[a,a^{\dagger}\right]=1$, and
$\eta_{D}$ and $\Omega_{D}$ are the amplitude and the frequency of the
modulating field. The quadratures satisfy the following commutation relation:
$\left[q_{j},p_{j^{\prime}}\right]=\iota\delta_{jj^{\prime}}.$

On putting the expression (5) of $\omega_{cav}(q_{1},q_{2})$, the Hamiltonian
in the frame rotating with the laser frequency $\omega_{l}$ can be written as

$\displaystyle H=\Delta
a^{\dagger}a+\sum_{j=1,2}\left[\frac{\omega_{mj}}{2}\left(q_{j}^{2}+p_{j}^{2}\right)+\left(-g_{1}^{(j)}q_{j}+g_{2}^{(j)}q_{j}^{2}\right)a^{\dagger}a\right]-g_{3}q_{1}q_{2}a^{\dagger}a+\iota
E\left[1+\eta_{D}\cos(\Omega_{D}t)\right]\left(a^{\dagger}-a\right)\;.$ (8)

Here $\Delta=\omega_{c}-\omega_{l}$ denotes the cavity mode detuning.
Considering the dissipation of the system, we can now obtain the following
quantum Langevin equations for the relevant operators

$\displaystyle\frac{dq_{j}}{dt}$ $\displaystyle=$
$\displaystyle\omega_{mj}p_{j}\;,$ $\displaystyle\frac{dp_{j}}{dt}$
$\displaystyle=$
$\displaystyle-\omega_{mj}q_{j}+g_{1}^{(j)}a^{\dagger}a-2g_{2}^{(j)}a^{\dagger}aq_{j}+g_{3}a^{\dagger}aq_{3-j}-\gamma_{mj}p_{j}+\xi_{j}(t)\;,$
$\displaystyle\frac{da}{dt}$ $\displaystyle=$
$\displaystyle-(\kappa+\iota\Delta)a+\iota\sum_{j=1,2}\left(g_{1}^{(j)}q_{j}-g_{2}^{(j)}q_{j}^{2}\right)a+\iota
g_{3}q_{1}q_{2}a+E\left[1+\eta_{D}\cos\left(\Omega_{D}t\right)\right]+\sqrt{2\kappa}a_{in}\;,$
(9)

where we have used the Hamiltonian (8). Here $\kappa$ is the decay rate of
cavity mode, $\gamma_{mj}$ is the damping rate of $j$th mechanical oscillator,
and $a_{in}$ is the input noise operator with the following two-time
correlation functions: walls2008input :

$\left\langle
a_{in}(t)a_{in}^{{\dagger}}(t^{\prime})\right\rangle=\delta(t-t^{\prime})\;,\;\left\langle\
a_{in}^{{\dagger}}(t)a_{in}(t^{\prime})\right\rangle=0\;.$ (10)

On the other hand, the thermal bath at an equilibrium temperature $T$ is
described by a Brownian noise $\xi_{j}(t)$ with zero mean and the following
two-time correlation function:

$\left\langle\xi_{j}(t)\xi_{j}\left(t^{\prime}\right)\right\rangle=\frac{\gamma_{mj}}{2\pi\omega_{mj}}\int\omega
e^{-\iota\omega\left(t-t^{\prime}\right)}\left[1+\operatorname{coth}\left(\frac{\hbar\omega}{2k_{B}T}\right)\right]d\omega\;,$
(11)

where $k_{B}$ is the Boltzmann constant. For the case of large quality factor
of the mechanical oscillator, the expression above reduces to the following
Markovian-approximated form: giovannetti2001phase :

$\displaystyle\left\langle\xi_{j}(t)\xi_{j}\left(t^{\prime}\right)\right\rangle=\gamma_{mj}\left(2\bar{n}_{mj}+1\right)\delta\left(t-t^{\prime}\right)\;,$
(12)

where $\bar{n}_{mj}=1/\left[\exp\left(\hbar\omega_{mj}/k_{B}T\right)-1\right]$
is the mean phonon number of the $j$th mechanical oscillator.

###  II.1 Solution in mean-field approximation

Since it is difficult to solve equation (9) analytically, we use the mean-
field approximation to simplify the calculation by rewriting the operators as
the sum of their mean values and quantum fluctuation near mean value, i.e.,
$a\rightarrow A+\delta a$, $q_{j}\rightarrow{Q_{j}}+\delta q_{j}$, and
$p_{j}\rightarrow{P_{j}}+\delta p_{j}$, where $A=\left\langle a\right\rangle$,
$Q_{j}=\left\langle q_{j}\right\rangle$. and $P_{j}=\left\langle
p_{j}\right\rangle$. Thus, the quantum Langevin equations can be split into
two sets of equations, namely, (i) the nonlinear equations for the mean
values:

$\displaystyle\frac{dQ_{j}}{dt}$ $\displaystyle=$
$\displaystyle\omega_{mj}P_{j}\;,$ $\displaystyle\frac{dP_{j}}{dt}$
$\displaystyle=$
$\displaystyle-\omega_{mj}Q_{j}+\left(g_{1}^{(j)}+g_{3}Q_{3-j}\right)|\alpha|^{2}-2g_{2}^{(j)}Q_{j}|\alpha|^{2}-\gamma_{mj}P_{j}\;,$
$\displaystyle\frac{dA}{dt}$ $\displaystyle=$
$\displaystyle-(\kappa+\iota\Delta)A+\iota\sum_{j=1,2}\left(g_{1}^{(j)}Q_{j}-g_{2}^{(j)}Q_{j}^{2}\right)A+\iota
g_{3}Q_{1}Q_{2}A+E\left[1+\eta_{D}\cos\left(\Omega_{D}t\right)\right]\;,$ (13)

and (ii) the linearized equations for the quantum fluctuations:

$\displaystyle\frac{d}{dt}\delta q_{j}$ $\displaystyle=$
$\displaystyle\omega_{j}\delta p_{j}\;,$ $\displaystyle\frac{d}{dt}\delta
p_{j}$ $\displaystyle=$
$\displaystyle-\left(\omega_{mj}+2g_{2}^{(j)}|A|^{2}\right)\delta
q_{j}+g_{3}|A|^{2}\delta q_{3-j}+G_{j}\left(A\delta a^{\dagger}+A^{*}\delta
a\right)-\gamma_{mj}\delta p_{j}+\xi_{j}(t)\;,$
$\displaystyle\frac{d}{dt}\delta a$ $\displaystyle=$ $\displaystyle-\iota
F\delta a+\iota G_{1}A\delta q_{1}+\iota G_{2}A\delta q_{2}-\kappa\delta
a+\sqrt{2\kappa}a_{in}\;,$ (14)

where

$\displaystyle G_{j}$ $\displaystyle=$ $\displaystyle
g_{1}^{(j)}-2g_{2}^{(j)}Q_{j}+g_{3}Q_{3-j}\;,$ $\displaystyle F$
$\displaystyle=$ $\displaystyle\Delta-
g_{1}^{(1)}Q_{1}-g_{1}^{(2)}Q_{2}+g_{2}^{(1)}Q_{1}^{2}+g_{2}^{(2)}Q_{2}^{2}-g_{3}Q_{1}Q_{2}\;.$
(15)

Here, we have ignored the second and higher-order terms in fluctuations.

In order to calculate the desired markers, $S_{q}$ and $E_{D}$, respectively
for QS and entanglement, we need to obtain the solutions for the quadrature
fluctuations of the oscillators. Solving the Eqs. (14) will be convenient if
we replace the intra-cavity field and the input noise operators by their
quadratures, as well: $\delta x=\frac{1}{\sqrt{2}}\left(\delta
a^{\dagger}+\delta a\right)$, $\delta y=\frac{\iota}{\sqrt{2}}\left(\delta
a^{\dagger}-\delta a\right)$, $\delta x_{in}=\frac{1}{\sqrt{2}}\left(\delta
a^{\dagger}_{in}+\delta a_{in}\right)$, and $\delta
y_{in}=\frac{\iota}{\sqrt{2}}\left(\delta a^{\dagger}_{in}-\delta
a_{in}\right)$. Therefore, the Eqs. (14) takes a simpler form, as given by

$\displaystyle\dot{M}(t)=BM(t)+J(t)\;,$ (16)

where $M(t)^{T}=\left(\delta q_{1},\delta p_{1},\delta q_{2},\delta
p_{2},\delta x,\delta y\right)$, the vector $J(t)$ as given below contains the
noise terms:

$\displaystyle J(t)^{T}$ $\displaystyle=$
$\displaystyle\bigg{(}0,\xi_{1},0,\xi_{2},\sqrt{2\kappa}\delta
x_{in},\sqrt{2\kappa}\delta y_{in}\bigg{)}\;,$ (17)

and

$\displaystyle B=\begin{pmatrix}0&\omega_{m1}&0&0&0&0\\\
-\omega_{m1}-2g_{2}^{(1)}|A|^{2}&-\gamma_{m1}&g_{3}|A|^{2}&0&\sqrt{2}G_{1}\operatorname{Re}\left(A\right)&\sqrt{2}G_{1}\operatorname{Im}\left(A\right)\\\
0&0&0&\omega_{m2}&0&0\\\
g_{3}|A|^{2}&0&-\omega_{m2}-2g_{2}^{(2)}|A|^{2}&-\gamma_{m2}&\sqrt{2}G_{2}\operatorname{Re}\left(A\right)&\sqrt{2}G_{2}\operatorname{Im}\left(A\right)\\\
-\sqrt{2}G_{1}\operatorname{Im}\left(A\right)&0&-\sqrt{2}G_{2}\operatorname{Im}\left(A\right)&0&-\kappa&F\\\
\sqrt{2}G_{1}\operatorname{Re}\left(A\right)&0&\sqrt{2}G_{2}\operatorname{Re}\left(A\right)&0&-F&-\kappa\\\
\end{pmatrix}$ (18)

is a $6\times 6$ time-dependent coefficient matrix.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 2: (a) Limit-cycle trajectories in the $Q_{1}\leftrightharpoons P_{1}$
(red) and $Q_{2}\leftrightharpoons P_{2}$ (blue) spaces, Variation of (b) the
mean values ${Q_{1}}$ (red) and ${Q_{2}}$ (blue), (c) the mean values
${P_{1}}$ (red) and ${P_{2}}$ (blue), (d) synchronization $S_{q}$, (e)
entanglement $E_{D}$, with respect to time (in the units of
$\tau=1/\omega_{m1}$). The parameters chosen are $\omega_{m1}=-\Delta=1$,
$\omega_{m2}=1.005$, $\bar{n}_{m}=0.5$ , $g_{1}=5\times 10^{-5}$,
$g_{2}=g_{1}\times 10^{-2}$, $g_{3}=10^{-6}$, $\gamma_{mj}=0.009$,
$\kappa=0.1$, $E=250$, $\eta_{D}=4$ and $\Omega_{D}=1$. All frequencies are
normalized with respect to $\omega_{m1}$. (f) Evolution of time-averaged
values of synchronization $\bar{S}_{q}$ (blue) and entanglement $\bar{E}_{D}$
(red) with respect to the frequency difference
$\delta_{m}=\omega_{m2}-\omega_{m1}$ of the mechanical oscillators. [all the
other parameters are the same as in (a)-(e) above].

The initial states of the oscillators can be best approximated as Gaussian
states with the peak at the mean position of the respective oscillators, along
with minimum uncertainties in position quadratures. The fluctuation dynamics
of the system is governed by a set of linearized equations, and this ensures
that the evolved states also remain Gaussian. Since the Gaussian state can be
fully characterized by its covariance matrix $C(t)$, we use it to calculate
the correlation between the quantum fluctuations of the system variables and
to compute all the relevant measures, namely, $S_{q}$ and $E_{D}$. The
temporal dynamics of $C(t)$ is governed by the following linear ordinary
differential equation:

$\displaystyle\dot{C}(t)=B(t)C(t)+B(t)C(t)^{T}+\zeta\;,$ (19)

where the elements of $C$ are given by $C_{ij}=\left[\left\langle
M_{i}(\infty)M_{j}(\infty)+M_{j}(\infty)M_{i}(\infty)\right\rangle\right]/2$
and the matrix $\zeta$, as below, describes the diffusion of the system:

$\displaystyle\zeta$ $\displaystyle=$ $\displaystyle{\rm
diag}\big{[}0,(2\bar{n}_{m1}+1)\gamma_{m1},0,(2\bar{n}_{m2}+1)\gamma_{m2},\kappa,\kappa\big{]}\;.$
(20)

In the matrix $C$, the diagonal (off-diagonal) element represents the variance
of the respective mode (covariance of two modes).

The QS quantifier $S_{q}(t)$ can then be expressed as

$\displaystyle
S_{q}(t)=\left\\{\frac{1}{2}\left[C_{11}(t)+C_{33}(t)-C_{13}(t)-C_{31}(t)+C_{22}(t)+C_{44}(t)-C_{24}(t)-C_{42}(t)\right]\right\\}^{-1}\;,$
(21)

and entanglement marker $E_{D}(t)$ as

$\displaystyle
E_{D}(t)=\frac{1}{4}\left[C_{11}(t)+C_{33}(t)-C_{13}(t)-C_{31}(t)\right]\times\left[C_{22}(t)+C_{44}(t)+C_{24}(t)+C_{42}(t)\right]\;.$
(22)

###  II.2 Numerical results

We show in Figs. 2, the time-evolution of the quadratures of the oscillators,
as obtained by simultaneously solving Eqs. (13). In Fig. 2a, we show how the
evolution of $Q_{1}\rightleftharpoons P_{1}$ and $Q_{2}\rightleftharpoons
P_{2}$ of the two oscillators tend to an asymptotic periodic orbit and the two
limit cycles tend to be consistent with each other. This becomes apparent from
the Figs. 2b and 2c, as the mean positions $Q_{1}$ and $Q_{2}$ and also the
mean linear momenta $P_{1}$ and $P_{2}$ are found to oscillate exactly in
phase at the long times. This refers to the complete classical synchronization
between mechanical oscillators.

With the classical limit cycle as a pre-condition, we next analyze how the
$S_{q}$ and $E_{D}$ vary with time. From the Fig. (2d), we can see that
$S_{q}$ reaches a stable value close to unity that refers to complete QS. More
importantly, we find that this is associated with entanglement. From the Fig.
(2e), we can see that $E_{D}$ becomes less than 0.25. This shows that
entanglement criterion (3) is violated by the state of the two oscillators. By
comparing Figs. (2e) and (2d), it can be observed that the generation of
entanglement and synchronization between the oscillators occur at a similar
time scale. This indicates a possible relation between the onset of QS and the
generation of entanglement between the oscillators. In fact, the long-time
behavior, as displayed in Figs. (2d) and (2e), clearly shows that $S_{q}$ and
$E_{D}$ oscillate in similar way at a particular frequency.

We further show in Fig. (2f), how the QS and entanglement behave with increase
in frequency difference $\delta_{m}=\omega_{m2}-\omega_{m1}$. The time-
averaged values of the respective markers $S_{q}$ and $E_{D}$ for QS and the
entanglement deteriorates with rise in frequency difference $\delta_{m}$. One
can observe that the quantum complete synchronization (i.e., $\bar{S}_{q}=1$)
and the maximal violation of entanglement inseparability (i.e.,
$\bar{E}_{D}-1/4$ is minimum) occur around the resonance condition
$\delta_{m}=0$, i.e., for a pair of identical mechanical oscillators (i.e.,
$\omega_{m1}=\omega_{m2}$). Interestingly, they stay robust (i.e., $S_{q}$
remains less than 1 and $E_{D}$ less than 0.25) and therefore maintain both QS
and entanglement, against the asymmetry in the oscillator frequencies, as
large as $\delta_{m}\sim 3.2$. As $\delta_{m}\gtrsim 3.2$, the oscillators no
longer remain entanglement (as $E_{D}$ exceeds 0.25) and quantum-synchronized
(as $S_{q}$ exceeds 1) and both these features vanish together. This further
strengthens our conjecture that QS and entanglement must originate from the
same type of quantum correlation.

##  III Analytic solution of fluctuations

In order to have a deeper insight into entanglement and synchronization
behaviour, we next obtain an analytical solution of the mean square
fluctuations in the relative displacement $q_{-}$, the total momentum $p_{+}$
and relative momentum $p_{-}$. We start with the Eqs. (9), that can be
rewritten for these variables as

$\displaystyle\frac{dq_{-}}{dt}$ $\displaystyle=$
$\displaystyle\omega_{m}p_{-}\;,$ $\displaystyle\frac{dq_{+}}{dt}$
$\displaystyle=$ $\displaystyle\omega_{m}p_{+}\;,$
$\displaystyle\frac{dp_{-}}{dt}$ $\displaystyle=$
$\displaystyle-\omega_{m}q_{-}-2g_{2}a^{\dagger}aq_{-}-g_{3}a^{\dagger}aq_{-}-\gamma_{m}p_{-}+\frac{1}{\sqrt{2}}\left(\xi_{1}-\xi_{2}\right)\;,$
$\displaystyle\frac{dp_{+}}{dt}$ $\displaystyle=$
$\displaystyle-\omega_{m}q_{+}+\sqrt{2}g_{1}a^{\dagger}a-2g_{2}a^{\dagger}aq_{+}+g_{3}a^{\dagger}aq_{+}-\gamma_{m}p_{+}+\frac{1}{\sqrt{2}}\left(\xi_{1}+\xi_{2}\right)\;,$
$\displaystyle\frac{da}{dt}$ $\displaystyle=$
$\displaystyle-(\kappa+\iota\Delta)a+\iota\left[\sqrt{2}g_{1}q_{+}-g_{2}(q_{1}^{2}+q_{2}^{2})+g_{3}q_{1}q_{2}\right]a+E\left[1+\eta_{D}\cos\left(\Omega_{D}t\right)\right]+\sqrt{2\kappa}a_{in}\;,$
(23)

where we have assumed $\omega_{m1}=\omega_{m2}=\omega_{m}$,
$g_{1}^{(j)}=g_{1}$ and $g_{2}^{(j)}=g_{2}$ for all $j$, for simplicity.

###  III.1 Asymptotic solutions of first moments

We next expand the operators $q_{\pm}$ and $p_{\pm}$ as a sum of their mean
values and fluctuations as follows: $q_{\pm}=Q_{\pm}+\delta q_{\pm}$,
$p_{\pm}=P_{\pm}+\delta p_{\pm}$. The Langevin equations (23) can then be
split into two sets: one for the mean values, as listed below, and the other
one for the fluctuation operators, to be discussed later.

$\displaystyle\frac{dQ_{-}}{dt}$ $\displaystyle=$
$\displaystyle\omega_{m}P_{-}\;,$ $\displaystyle\frac{dQ_{+}}{dt}$
$\displaystyle=$ $\displaystyle\omega_{m}P_{+}\;,$
$\displaystyle\frac{dP_{-}}{dt}$ $\displaystyle=$
$\displaystyle-\omega_{m}Q_{-}-2g_{2}|A|^{2}Q_{-}-g_{3}|A|^{2}Q_{-}-\gamma_{m}P_{-}\;,$
$\displaystyle\frac{dP_{+}}{dt}$ $\displaystyle=$
$\displaystyle-\omega_{m}Q_{+}+\sqrt{2}g_{1}|A|^{2}-2g_{2}|A|^{2}Q_{+}+g_{3}|A|^{2}Q_{+}-\gamma_{m}P_{+}\;,$
$\displaystyle\frac{dA}{dt}$ $\displaystyle=$
$\displaystyle-(\kappa+\iota\Delta)A+\iota\left[\sqrt{2}g_{1}Q_{+}-g_{2}(Q_{1}^{2}+Q_{2}^{2})+g_{3}Q_{1}Q_{2}\right]A+E\left[1+\eta_{D}\cos\left(\Omega_{D}t\right)\right]\;.$
(24)

Since the cavity is driven by a modulating field
$E(t)=E[1+\eta_{D}\cos\left(\Omega_{D}t\right)]=E_{0}+E_{1}e^{-i\Omega_{D}t}+E_{-1}e^{i\Omega_{D}t}$,
the amplitudes of the cavity mode and the mechanical modes would also follow
the dynamics of this field at a long-time limit, according to the Floquet
theorem, i.e., $\lim_{t\rightarrow\infty}A(t)=A(t+\tau)$,
$\lim_{t\rightarrow\infty}Q_{1,2}(t)=Q_{1,2}(t+\tau)$ and
$\lim_{t\rightarrow\infty}P_{1,2}(t)=P_{1,2}(t+\tau)$ mari2009gently . We
therefore redefine these classical amplitudes in terms of Fourier components
to the first harmonic of the modulating frequency $\Omega_{D}$, as follows.

$\displaystyle A$ $\displaystyle=$ $\displaystyle
A_{-1}e^{i\Omega_{D}t}+A_{0}+A_{1}e^{-i\Omega_{D}t}\;,$ $\displaystyle Q_{-}$
$\displaystyle=$ $\displaystyle
Q^{-}_{-1}e^{i\Omega_{D}t}+Q^{-}_{0}+Q^{-}_{1}e^{-i\Omega_{D}t}\;,$
$\displaystyle P_{-}$ $\displaystyle=$ $\displaystyle
P^{-}_{-1}e^{i\Omega_{D}t}+P^{-}_{0}+P^{-}_{1}e^{-i\Omega_{D}t}\;,$
$\displaystyle Q_{+}$ $\displaystyle=$ $\displaystyle
Q^{+}_{-1}e^{i\Omega_{D}t}+Q^{+}_{0}+Q^{+}_{1}e^{-i\Omega_{D}t}\;,$
$\displaystyle P_{+}$ $\displaystyle=$ $\displaystyle
P^{+}_{-1}e^{i\Omega_{D}t}+P^{+}_{0}+P^{+}_{1}e^{-i\Omega_{D}t}\;.$ (25)

####  III.1.1 Obtaining the Fourier coefficients of $Q_{-}$ and $P_{-}$

To find the time-independent coefficients $Q_{-1,0,1}^{-}$, we first
substitute $A$, $Q_{-}$, and $P_{-}$ from Eqs. (25) into the Langevin
equations (24) for $Q_{-}$ and $P_{-}$. This leads us to the following set of
three coupled nonlinear equations, when $\Omega_{D}=\omega_{m}$:

$\displaystyle\left(\omega_{m}+W\right)Q_{0}^{-}+W_{0}Q_{1}^{-}+W_{0}^{*}Q_{-1}^{-}$
$\displaystyle=$ $\displaystyle 0\;,$ $\displaystyle
W_{0}Q_{0}^{-}+W_{1}Q_{1}^{-}+\left(\iota\gamma_{m}+W\right)Q_{-1}^{-}$
$\displaystyle=$ $\displaystyle 0\;,$ $\displaystyle
W_{0}^{*}Q_{0}^{-}+\left(-\iota\gamma_{m}+W\right)Q_{1}^{-}+W_{1}^{*}Q_{-1}^{-}$
$\displaystyle=$ $\displaystyle 0\;,$ (26)

where

$\displaystyle W_{0}$ $\displaystyle=$
$\displaystyle(2g_{2}+g_{3})(A_{0}A_{1}^{*}+A_{0}^{*}A_{-1})\;,$
$\displaystyle W_{1}$ $\displaystyle=$
$\displaystyle(2g_{2}+g_{3})(A_{1}^{*}A_{-1})\;,$ $\displaystyle W$
$\displaystyle=$
$\displaystyle(2g_{2}+g_{3})(|A_{0}|^{2}+|A_{1}|^{2}+|A_{-1}|^{2})\;.$ (27)

The trivial solution of the above set of equations (26) is $Q_{-1,0,1}^{-}=0$,
so that $Q_{-}$ vanishes. From the equation of $P_{-}$, we also have
$P_{-}(t)=\exp(-\gamma_{m}t)P_{-}(0)$, which vanishes as $t\rightarrow\infty$.
Both the results $Q_{-}=0,P_{-}=0$ are expected in a synchronized system. This
further implies that $Q_{1}=Q_{2}=Q_{+}/\sqrt{2}$.

####  III.1.2 Obtaining the Fourier coefficients $A_{-1,0,1}$

Putting the above results of $Q_{1,2}$ in the Langevin equation of $A$ in
(24), we find that this reduces to the following equation:

$\displaystyle\frac{dA}{dt}$ $\displaystyle=$
$\displaystyle-(\kappa+\iota\Delta)A+\iota\left[\sqrt{2}g_{1}Q_{+}-\left(\frac{2g_{2}-g_{3}}{2}\right)Q_{+}^{2}\right]A+E\left[1+\eta_{D}\cos\left(\Omega_{D}t\right)\right]\;,$
(28)

which contains only $Q_{+}$. Substituting $A$ and $Q_{+}$ from the Eqs. (25)
into (28) and separating the different Fourier components, we obtain the
following coupled equations of the time-independent coefficients $A_{-1,0,1}$:

$\displaystyle UA_{0}+U_{-1}A_{1}+U_{1}A_{-1}$ $\displaystyle=$ $\displaystyle
E_{0}\;,$ $\displaystyle
U_{1}A_{0}+\left[-\iota\omega_{m}+U\right]A_{1}+\left[\frac{\iota}{2}(2g_{2}-g_{3})(Q_{1}^{+})^{2}\right]A_{-1}$
$\displaystyle=$ $\displaystyle E_{1}\;,$ $\displaystyle
U_{-1}A_{0}+\left[\iota\omega_{m}+U\right]A_{-1}+\left[\frac{\iota}{2}(2g_{2}-g_{3})(Q_{-1}^{+})^{2}\right]A_{1}$
$\displaystyle=$ $\displaystyle E_{-1}\;,$ (29)

where

$\displaystyle U$ $\displaystyle=$
$\displaystyle\left(\kappa+\iota\Delta\right)-\iota\sqrt{2}g_{1}Q_{0}^{+}+\iota(2g_{2}-g_{3})\left(\frac{1}{2}Q_{0}^{+2}+Q_{1}^{+}Q_{-1}^{+}\right)\;,$
$\displaystyle U_{1}$ $\displaystyle=$
$\displaystyle-\iota\sqrt{2}g_{1}Q_{1}^{+}+\iota(2g_{2}-g_{3})\left(Q_{0}^{+}Q_{1}^{+}\right)\;,$
$\displaystyle U_{-1}$ $\displaystyle=$
$\displaystyle-\iota\sqrt{2}g_{1}Q_{-1}^{+}+\iota(2g_{2}-g_{3})\left(Q_{0}^{+}Q_{-1}^{+}\right)\;.$
(30)

Here we consider the weak optomechanical coupling regime, i.e.,
$|g_{i}|\ll\omega_{m},\kappa$ $(i\in 1,2,3)$, such that
$U\approx(\kappa+\iota\Delta)$ and $U_{\pm 1}$ may be neglected. Therefore, in
the long-time limit, we have

$\displaystyle
A_{0}=\frac{E_{0}}{\kappa+\iota\Delta},\;A_{1}=\frac{E_{1}}{\kappa+\iota(\Delta-\omega_{m})},\;A_{-1}=\frac{E_{-1}}{\kappa+\iota(\Delta+\omega_{m})}\;.\;$
(31)

####  III.1.3 Obtaining the Fourier coefficients of $Q_{+}$ and $P_{+}$

So far, we have obtained the expressions of the Fourier components of $Q_{-}$
and $A$. To find out $Q_{+}$, we again substitute $A$, $Q_{+}$, and $P_{+}$
from Eqs. (25) into the Langevin equation (24) for $Q_{+}$ and $P_{+}$, that
leads us to the following set of three coupled linear equations of
$Q^{+}_{-1,0,1}$:

$\displaystyle\left(\omega_{m}+V\right)Q_{0}^{+}+V_{0}Q_{1}^{+}+V_{0}^{*}Q_{-1}^{+}$
$\displaystyle=$
$\displaystyle\sqrt{2}g_{1}(|A_{0}|^{2}+|A_{1}|^{2}+|A_{-1}|^{2})\;,$
$\displaystyle
V_{0}Q_{0}^{+}+V_{1}Q_{1}^{+}+\left(\iota\gamma_{m}+V\right)Q_{-1}^{+}$
$\displaystyle=$
$\displaystyle\sqrt{2}g_{1}(A_{0}A_{1}^{*}+A_{0}^{*}A_{-1})\;,$ $\displaystyle
V_{0}^{*}Q_{0}^{+}+\left(-\iota\gamma_{m}+V\right)Q_{1}^{+}+V_{1}^{*}Q_{-1}^{+}$
$\displaystyle=$
$\displaystyle\sqrt{2}g_{1}(A_{0}A_{-1}^{*}+A_{0}^{*}A_{1})\;,$ (32)

where

$\displaystyle V_{0}$ $\displaystyle=$
$\displaystyle(2g_{2}-g_{3})(A_{0}A_{1}^{*}+A_{0}^{*}A_{-1})\;,$
$\displaystyle V_{1}$ $\displaystyle=$
$\displaystyle(2g_{2}-g_{3})(A_{1}^{*}A_{-1})\;,$ $\displaystyle V$
$\displaystyle=$
$\displaystyle(2g_{2}-g_{3})(|A_{0}|^{2}+|A_{1}|^{2}+|A_{-1}|^{2})\;.$ (33)

Using Eq. (31), we can find the time-independent coefficients $Q^{+}_{-1,0,1}$
by algebraically solving the above equations. The expressions of $P_{+}$ can
then be easily obtained using the equation for $Q_{+}$ in (24).

###  III.2 Frequency spectrum of fluctuations

To investigate the QS and entanglement between mechanical oscillators, we need
to calculate the mean square fluctuations in $q_{-}$ and $p_{\pm}$. The
dynamical behaviour of these fluctuations can be described by the following
Langevin equations, as obtained from Eqs. (23) [see Sec. III.1]:

$\displaystyle\frac{d}{dt}\delta q_{-}$ $\displaystyle=$
$\displaystyle\omega_{m}\delta p_{-}\;,$ $\displaystyle\frac{d}{dt}\delta
p_{-}$ $\displaystyle=$
$\displaystyle-\left[\omega_{m}+\left(2g_{2}+g_{3}\right)\left(|A_{0}|^{2}+|A_{1}|^{2}+|A_{-1}|^{2}\right)\right]\delta
q_{-}-\gamma_{m}\delta p_{-}+\frac{\xi_{1}-\xi_{2}}{\sqrt{2}}\;,$
$\displaystyle\frac{d}{dt}\delta q_{+}$ $\displaystyle=$
$\displaystyle\omega_{m}\delta p_{+}\;,$ $\displaystyle\frac{d}{dt}\delta
p_{+}$ $\displaystyle=$ $\displaystyle-F_{0}\delta q_{+}+F_{1}\delta
a^{\dagger}+F_{2}\delta a-\gamma_{m}\delta
p_{+}+\frac{\xi_{1}+\xi_{2}}{\sqrt{2}}\;,$ $\displaystyle\frac{d}{dt}\delta a$
$\displaystyle=$ $\displaystyle-\left(\kappa+\iota\Delta^{\prime}\right)\delta
a+\iota F_{1}\delta q_{+}+\sqrt{2\kappa}\delta a_{in}\;,$
$\displaystyle\frac{d}{dt}\delta a^{\dagger}$ $\displaystyle=$
$\displaystyle-\left(\kappa-\iota\Delta^{\prime}\right)\delta
a^{\dagger}-\iota F_{2}\delta q_{+}+\sqrt{2\kappa}\delta a_{in}^{\dagger}\;,$
(34)

where

$\displaystyle F_{0}$ $\displaystyle=$
$\displaystyle\omega_{m}+(2g_{2}-g_{3})\left(|A_{0}|^{2}+|A_{1}|^{2}+|A_{-1}|^{2}\right)\;,$
$\displaystyle F_{1}$ $\displaystyle=$
$\displaystyle\sqrt{2}g_{1}A_{0}-(2g_{2}-g_{3})\left(Q_{0}^{+}A_{0}+Q_{1}^{+}A_{-1}+Q_{-1}^{+}A_{1}\right)\;,$
$\displaystyle F_{2}$ $\displaystyle=$
$\displaystyle\sqrt{2}g_{1}A_{0}^{*}-(2g_{2}-g_{3})\left(Q_{0}^{+}A_{0}^{*}+Q_{1}^{+}A_{1}^{*}+Q_{-1}^{+}A_{-1}^{*}\right)\;,$
$\displaystyle\Delta^{\prime}$ $\displaystyle=$
$\displaystyle\Delta-\sqrt{2}g_{1}Q_{0}^{+}+\left(\frac{2g_{2}-g_{3}}{2}\right)\left(Q_{0}^{+2}+2Q_{1}^{+}Q_{-1}^{+}\right)\;.$
(35)

These equation explicitly depend upon the Fourier components of $A$ and
$Q_{+}$. Introducing the cavity field quadratures $\delta x=\frac{\delta
a+\delta a^{\dagger}}{\sqrt{2}}$ and $\delta y=\frac{i(\delta
a^{\dagger}-\delta a)}{\sqrt{2}}$, and the input noise quadratures $\delta
x_{in}=\frac{\delta a_{in}+\delta a_{in}^{\dagger}}{\sqrt{2}}$ and $\delta
y_{in}=\frac{i(\delta a_{in}^{\dagger}-\delta a_{in})}{\sqrt{2}}$, the Eqs.
(34) can be reduced to the following matrix form:

$\displaystyle\dot{f}(t)=Ff(t)+G(t)\;,$ (36)

where $f(t)$ is the fluctuation vector and $G(t)$ is the noise vector, with
the respective transposes given by,

$\displaystyle f(t)^{T}=\left(\delta q_{+},\delta p_{+},\delta x,\delta
y\right),G(t)^{T}=\left(0,\frac{\xi_{1}+\xi_{2}}{\sqrt{2}},\sqrt{2\kappa}\delta
x_{in},\sqrt{2\kappa}\delta y_{in}\right)\;,$ (37)

and the matrix $F$ is given by

$\displaystyle F=\begin{pmatrix}0&\omega_{m}&0&0\\\
-F_{0}&-\gamma_{m}&\frac{F_{1}+F_{2}}{\sqrt{2}}&\frac{i(F_{2}-F_{1})}{\sqrt{2}}\\\
\frac{-i(F_{2}-F_{1})}{\sqrt{2}}&0&-\kappa&\Delta^{\prime}\\\
\frac{F_{1}+F_{2}}{\sqrt{2}}&0&-\Delta^{\prime}&-\kappa\\\ \end{pmatrix}\;.$
(38)

Taking the Fourier transformation of each operator in Eqs. (34) and solving
them in frequency domain, the fluctuations spectrum of the $q_{-}$ and $p_{-}$
can be obtained as

$\displaystyle\delta q_{-}(\omega)$ $\displaystyle=$
$\displaystyle\frac{-\omega_{m}}{d(\omega)}\left(\frac{\xi_{1}(\omega)-\xi_{2}(\omega)}{\sqrt{2}}\right)\;,$
$\displaystyle\delta p_{-}(\omega)$ $\displaystyle=$
$\displaystyle\frac{\iota\omega}{d(\omega)}\left(\frac{\xi_{1}(\omega)-\xi_{2}(\omega)}{\sqrt{2}}\right)\;,$
(39)

where
$d(\omega)=\omega^{2}+\iota\omega\gamma_{m}-\omega_{m}^{2}-\omega_{m}(2g_{2}+g_{3})\left(|A_{0}|^{2}+|A_{1}|^{2}+|A_{-1}|^{2}\right)$.
Similarly, the fluctuation in the total momentum has the following spectrum:

$\displaystyle\delta p_{+}(\omega)$ $\displaystyle=$
$\displaystyle\frac{\iota\omega}{D(\omega)}\left[\sqrt{2\kappa}\left\\{F_{1}[\kappa+\iota(\Delta^{\prime}-\omega)]\delta
a_{in}^{\dagger}(-\omega)+F_{2}[\kappa-\iota(\Delta^{\prime}+\omega)]\delta
a_{in}(\omega)\right\\}+[(\kappa-\iota\omega)^{2}+\Delta^{\prime
2}]\frac{\xi_{1}(\omega)+\xi_{2}(\omega)}{\sqrt{2}}\right]\;,$

where

$D(\omega)=2\Delta^{\prime}\omega_{m}F_{1}F_{2}+[\omega^{2}+\iota\omega\gamma_{m}-\omega_{m}^{2}-\omega_{m}(2g_{2}-g_{3})\left(|A_{0}|^{2}+|A_{1}|^{2}+|A_{-1}|^{2}\right)][(\kappa-\iota\omega)^{2}+\Delta^{\prime
2}]\;.$

We must emphasize here that the stability of such a system is essential to
achieve any synchronization. If all the eigenvalues of the matrix $F$ have
negative real parts, the stability can be ensured at long times. Using the
Routh-Hurwitz criterion dejesus1987routh , the corresponding conditions can be
derived as follows:

$\displaystyle\kappa\gamma_{m}\left[\left(\Delta^{\prime
2}+\kappa^{2}\right)^{2}+\left(F_{0}\omega_{m}+\gamma_{m}\kappa\right)^{2}+2\gamma_{m}\kappa\left(\kappa^{2}+\Delta^{\prime
2}\right)\right.$
$\displaystyle\left.+2F_{0}\omega_{m}\left(\kappa^{2}-\Delta^{\prime
2}\right)+\gamma_{m}^{2}\Delta^{\prime
2}\right]+F_{1}F_{2}\Delta^{\prime}\omega_{m}\left(\gamma_{m}+2\kappa\right)^{2}$
$\displaystyle>$ $\displaystyle 0\;,$ $\displaystyle
F_{0}\omega_{m}\left(\kappa^{2}+\Delta^{\prime
2}\right)-2F_{1}F_{2}\omega_{m}\Delta^{\prime}$ $\displaystyle>$
$\displaystyle 0\;.$ (41)

We assume that the system satisfies these stability conditions.

####  III.2.1 Derivation of mean square fluctuations

The mean square fluctuation of an operator $O(t)$ is determined by
huang2009entangling

$\displaystyle\langle\delta O(t)^{2}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{4\pi^{2}}\iint_{-\infty}^{\infty}d\omega d\Omega
e^{-\iota(\omega+\Omega)t}\langle\delta O(\omega)\delta O(\Omega)\rangle$ (42)

where the $\langle\delta O(\omega)\delta O(\Omega)\rangle$ refers to two-
frequency correlation function of the operator $O$. In the present case, we
use Eq. (39) and (LABEL:f7) into Eq. (42) to obtain the following forms of
mean square fluctuations:

$\displaystyle\langle\delta q_{-}(t)^{2}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{2\pi}\int_{-\infty}^{\infty}[\omega_{m}^{2}\mu]d\omega\;,$
(43) $\displaystyle\langle\delta p_{-}(t)^{2}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{2\pi}\int_{-\infty}^{\infty}[\omega^{2}\mu]d\omega\;,$
(44) $\displaystyle\langle\delta p_{+}(t)^{2}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{2\pi}\int_{-\infty}^{\infty}[\omega^{2}\nu]d\omega\;,$
(45)

where

$\displaystyle\mu$ $\displaystyle=$
$\displaystyle\frac{\gamma_{m}(2\bar{n}_{m}+1)}{d(\omega)d(-\omega)}\;,$
$\displaystyle\nu$ $\displaystyle=$
$\displaystyle\frac{1}{D(\omega)D(-\omega)}\left\\{2\kappa
F_{1}F_{2}[\kappa^{2}+(\Delta^{\prime}+\omega)^{2}]+\gamma_{m}(2\bar{n}_{m}+1)[(\Delta^{\prime
2}+\kappa^{2}-\omega^{2})^{2}+4\kappa^{2}\omega^{2}]\right\\}\;.$ (46)

In the above derivation of $\mu$ and $\nu$, we have used the following non-
vanishing frequency-domain correlation functions of the input noise operators,
as obtained by Fourier transformation of the Eqs. (11) and (12):

$\displaystyle\langle\delta a_{in}(\omega)\delta
a_{in}^{\dagger}(-\Omega)\rangle$ $\displaystyle=$ $\displaystyle
2\pi\delta(\omega+\Omega)\;,$
$\displaystyle\left\langle\xi_{j}(\omega)\xi_{k}\left(\Omega)\right)\right\rangle$
$\displaystyle=$ $\displaystyle
2\pi\delta_{jk}\gamma_{m}\left(2\bar{n}_{m}+1\right)\delta\left(\omega+\Omega\right)\;.$
(47)

Here we have assumed $\bar{n}_{mj}=\bar{n}_{m}$ for all $j$, for simplicity.

According to Wiener-Khintchine theorem mandel , the functions
$\omega_{m}^{2}\mu$, $\omega^{2}\mu$ and $\omega^{2}\nu$ represent the
spectral density of the random fluctuation of $q_{-}$, $p_{-}$, and $p_{+}$,
respectively. We find that they are even and converging functions of $\omega$,
and therefore lead to finite values, when integrated over frequencies.

(a)

(b)

(c)

Figure 3:  Variation of (a) $K$ from Eq. (51), (b) QS $S_{q}$, (c)
entanglement $E_{D}$ between the mechanical oscillators as a function of
modulation amplitude $\eta_{D}$ and modulation frequency $\Omega_{D}$ with
$E=250$ and $\bar{n}_{m}=0.5$. The other parameters are the same as in Fig. 2.

###  III.3 Generic relation between entanglement and QS

As seen in Eq. (I), $S_{q}$ always remains less than unity according to the
uncertainty principle, and therefore, we find a lower limit for $\langle\delta
q_{-}(t)^{2}\rangle$, as given by

$\langle\delta q_{-}(t)^{2}\rangle\geq 1-\langle\delta p_{-}(t)^{2}\rangle\;.$
(48)

that corresponds to a state with QS. In a similar way, we can find out from
(3) that there is an upper limit, as well, for $\langle\delta
q_{-}(t)^{2}\rangle$, to achieve entanglement. This is given by

$\langle\delta q_{-}(t)^{2}\rangle<\frac{1}{4\langle\delta
p_{+}(t)^{2}\rangle}\;.\\\ $ (49)

Therefore, to make the entanglement and QS, we must have

$\frac{1}{4\langle\delta p_{+}(t)^{2}\rangle}>1-\langle\delta
p_{-}(t)^{2}\rangle\;,$ (50)

which can be simplified to the following analytical condition:

$K=\frac{1}{4\langle\delta p_{+}(t)^{2}\rangle}+\langle\delta
p_{-}(t)^{2}\rangle-1>0\;.$ (51)

When $K=0$, the two limits of $\langle\delta q_{-}(t)^{2}\rangle$ become
equal, and it suggests a critical set of parameters, that are required for
simultaneous occurrence of entanglement and QS.

We can finally express the QS $S_{q}$ in terms of $E_{D}$ as

$S_{q}=\frac{\langle\delta p_{+}(t)^{2}\rangle}{E_{D}+\langle\delta
p_{-}(t)^{2}\rangle\langle\delta p_{+}(t)^{2}\rangle}\;.$ (52)

From the above expression, it is clear that when $E_{D}$ is less (a case when
both $\langle\delta q_{-}(t)^{2}\rangle$ and $\langle\delta
p_{+}(t)^{2}\rangle$ are less), the above expression of $S_{q}$ tends to
$\langle\delta p_{-}(t)^{2}\rangle^{-1}$. Thus for smaller $\langle\delta
p_{-}(t)^{2}\rangle$, $S_{q}$ increases. For sufficiently smaller $E_{D}$
(e.g., when it is less than 1/4), $S_{q}$ can increase towards unity [see Eq.
(52) above]. This means when the two oscillators become entangled, they can
tend to near complete synchronization. On the other hand, for larger $E_{D}$
(e.g., for example, when it much larger than 1/4), the value $S_{q}$ becomes
much less than unity. This suggests that in absence of entanglement, the
degree of QS reduces. Our analysis therefore suggests that entanglement
becomes a sufficient condition for QS. In the next Subsection, we will explore
this aspect more critically in our numerical studies.

(a)

(b)

Figure 4:  Variation of (a) QS $S_{q}$ and (b) entanglement $E_{D}$ between
the mechanical oscillators as a function of the driving field intensity $E$
and the average number of thermal phonons $\bar{n}_{m}$ with $\eta_{D}=4$ and
$\Omega_{D}=1$. The other parameters are the same as in Fig. 2.

###  III.4 Numerical results

As noted in Secs. II.2 and III.3, there exists a correlation between QS and
entanglement. To test this conjectiure, we present further numerical results
based on the analytical results in Sec. III. We choose all the parameters such
that they maintain the stability in the system, following Eqs. (41). We first
show in Figs. (3b) and (3c), the variation of the QS $S_{q}$ and entanglement
$E_{D}$, respectively, as a function of modulation amplitude $\eta_{D}$ and
modulation frequency $\Omega_{D}$. It is observed that both the QS and the
entanglement trace the classic ’Arnold tongue’, yielding a range of values of
$\eta_{D}$ and the corresponding values of $\Omega_{D}$ for which
synchronization and entanglement are simultaneously achievable. We further
find that the tongue is symmetric about $\Omega_{D}=\omega_{m}=1$. The
synchronization $S_{q}$ tends to one and the entanglement $E_{D}$ decreases
further below the upper limit of 0.25, as $\eta_{D}$ increases, i.e., as the
modulation strength increases. We note that similar entanglement tongue was
reported for two coupled van der Pol oscillators in lee2014 . However, this
was shown associated only with classical phase synchronization.

We show in Fig. (3a) the variation of $K$ [Eq. (51)] as a function of
$\eta_{D}$ and $\Omega_{D}$. This plot also exhibits the classic ’Arnold
tongue’ around modulation frequency $\Omega_{D}=1$. For larger $\eta_{D}$ at
$\Omega_{D}=1$, the term $K$ increases further. By comparing with the Figs. 3b
and 3c, we can conclude that with increase in $K$ the system approaches to
nearly complete QS (i.e., $S_{q}\lesssim 1$) and larger violation of the
inseparability inequality (3) (i.e., $E_{D}\ll 0.25$). It is also interesting
to see that $K$ remains positive for a large range of $\eta_{D}$ and
$\Omega_{D}$. Thus, both the QS (albeit partial) and entanglement stay robust
over a large variation of modulation parameters.

To get further insight into the correlated behaviour of synchronization and
entanglement, we plot them with respect to the driving field intensity $E$ and
average number of phonons $\bar{n}_{m}$, respectively, in the Figs. (4a) and
(4b), for a fixed modulation amplitude $\eta_{D}=4$ and modulation frequency
$\Omega_{D}=1$. It can be seen, for higher $E$ and at lower temperatures of
the bath ($\bar{n}_{m}=0.5$), the oscillators are nearly complete quantum
synchronized, along with maximal violation of inseparability inequality (3).
With increase in thermal excitation, the synchronization deteriorates and the
entanglement decreases too. However, the synchronization stays more robust to
temperature as compared to entanglement. This indicates that though some
residual synchronization may still exist without entanglement, the
entanglement is always associated with the synchronization and therefore
becomes the sufficient condition for the latter. In fact, this can also be
intuitively understood from the definition (I) of QS, which always remains
less than unity due to uncertainty principle, though the quantum state of the
system may not exhibit entanglement [i.e., may not violate the Eq. (3)]. It
must be borne in mind that the entanglement would indicate the existence of
QS, only in the presence of limit cycle of mean values of the joint
quadratures.

##  IV CONCLUSION

In conclusion, we have explored the interconnection between QS and
entanglement between two mechanical oscillators in an optomechanical system.
In our model, both mechanical oscillators are coupled with the same cavity
mode via linear and quadratic dependence on their displacement from their
respective equilibrium positions. An indirect always-on coupling, proportional
to $g_{3}a^{\dagger}a$ also arises between the oscillators, which is the key
to generation of the synchronization between them. Since the same uncertainty
relation sets the upper limit for QS $S_{q}$ and lower limit for entanglement
marker $E_{D}$, so we expected that there could be a correlation between the
two. In this regard, we have first demonstrated the classical synchronization
via limit cycle trajectories of the mean quadratures at long times, which is a
precondition to achieve QS. Our numerical results show that the two coupled
mechanical oscillators exhibit entanglement and QS, when the cavity is
strongly amplitude-modulated. A nearly complete QS and entanglement between
the coupled oscillators can be achieved for a large range of modulation
parameters. This result leads us to a strong conjecture that they arise from
the same EPR-type correlation. We further show that the entanglement manifests
itself as the sufficient condition for QS. One can characterize the QS as a
function of entanglement, when defined in terms of variances of quadratures.
We have provided all the relevant analysis, supported by numerical results in
this context. Our results open up a newer perspective to interpret the QS and
entanglement on the same footing.

## References

## References

  * (1) S. H. Strogatz, Nonlinear dynamics and chaos, vol. 7.  Levant Books, Kolkata, 2007.
  * (2) C. Huygens, Oeuvres complètes, vol. 7.  M. Nijhoff, 1897.
  * (3) M. J. Leone, B. N. Schurter, B. Letson, V. Booth, M. Zochowski, and C. G. Fink, “Synchronization properties of heterogeneous neuronal networks with mixed excitability type,” Physical Review E, vol. 91, no. 3, p. 032813, 2015\.
  * (4) G. Tang, K. Xu, and L. Jiang, “Synchronization in a chaotic neural network with time delay depending on the spatial distance between neurons,” Physical Review E, vol. 84, no. 4, p. 046207, 2011.
  * (5) J. Cao and R. Li, “Fixed-time synchronization of delayed memristor-based recurrent neural networks,” Science China Information Sciences, vol. 60, no. 3, p. 1, 2017.
  * (6) S. Vaidyanathan, “Dynamics and control of tokamak system with symmetric and magnetically confined plasma,” Int J ChemTech Res, vol. 8, no. 6, p. 795, 2015.
  * (7) T. E. Lee and H. Sadeghpour, “Quantum synchronization of quantum van der pol oscillators with trapped ions,” Physical Review Letters, vol. 111, no. 23, p. 234101, 2013.
  * (8) S. Walter, A. Nunnenkamp, and C. Bruder, “Quantum synchronization of a driven self-sustained oscillator,” Physical Review Letters, vol. 112, no. 9, p. 094102, 2014.
  * (9) S. Shirasaka, N. Watanabe, Y. Kawamura, and H. Nakao, “Optimizing stability of mutual synchronization between a pair of limit-cycle oscillators with weak cross coupling,” Physical Review E, vol. 96, no. 1, p. 012223, 2017.
  * (10) B. Van der Pol and J. Van der Mark, “Frequency demultiplication,” Nature, vol. 120, no. 3019, p. 363, 1927.
  * (11) L. M. Pecora and T. L. Carroll, “Synchronization in chaotic systems,” Physical Review Letters, vol. 64, no. 8, p. 821, 1990.
  * (12) A. Mari, A. Farace, N. Didier, V. Giovannetti, and R. Fazio, “Measures of quantum synchronization in continuous variable systems,” Physical Review Letters, vol. 111, no. 10, p. 103605, 2013.
  * (13) A. Peres, “Separability criterion for density matrices,” Physical Review Letters, vol. 77, no. 8, p. 1413, 1996.
  * (14) L.-M. Duan, G. Giedke, J. I. Cirac, and P. Zoller, “Inseparability criterion for continuous variable systems,” Physical Review Letters, vol. 84, no. 12, p. 2722, 2000.
  * (15) S. Mancini, V. Giovannetti, D. Vitali, and P. Tombesi, “Entangling macroscopic oscillators exploiting radiation pressure,” Physical Review Letters, vol. 88, no. 12, p. 120401, 2002.
  * (16) G. S. Agarwal and A. Biswas, “Inseparability inequalities for higher order moments for bipartite systems,” New Journal of Physics, vol. 7, p. 211, 2005.
  * (17) D. Garg, Manju, S. Dasgupta, and A. Biswas, “Quantum synchronization and entanglement of indirectly coupled mechanical oscillators in cavity optomechanics: A numerical study,” Physics Letters A, vol. 457, p. 128557, 2023.
  * (18) G. Manzano, F. Galve, G. L. Giorgi, E. Hernández-García, and R. Zambrini, “Synchronization, quantum correlations and entanglement in oscillator networks,” Scientific Reports, vol. 3, p. 1439, 2013.
  * (19) D. Witthaut, S. Wimberger, R. Burioni, and M. Timme, “Classical synchronization indicates persistent entanglement in isolated quantum systems,” Nature communications, vol. 8, no. 1, p. 1, 2017.
  * (20) T. E. Lee and M. Cross, “Quantum-classical transition of correlations of two coupled cavities,” Physical Review A, vol. 88, no. 1, p. 013834, 2013.
  * (21) G. L. Giorgi, F. Galve, G. Manzano, P. Colet, and R. Zambrini, “Quantum correlations and mutual synchronization,” Physical Review A, vol. 85, p. 052101, 2012.
  * (22) D. Viennot and L. Aubourg, “Quantum chimera states,” Physics Letters A, vol. 380, p. 678683, 2016.
  * (23) V. Ameri, M. Eghbali-Arani, A. Mari, A. Farace, F. Kheirandish, V. Giovannetti, and R. Fazio, “Mutual information as an order parameter for quantum synchronization,” Physical Review A, vol. 91, no. 1, p. 012301, 2015.
  * (24) V. M. Bastidas, I. Omelchenko, A. Zakharova, E. Schöll, and T. Brandes, “Quantum signatures of chimera states.,” Physical Review E, vol. 92, p. 062924, 2015.
  * (25) Y. Shen, H. Yi Soh, W. Fan, and L.-C. Kwek, “Enhancing quantum synchronization through homodyne measurement and squeezing,” arXiv:2302.13465, 2023.
  * (26) Y. Kato and H. Nakao, “Enhancement of quantum synchronization via continuous measurement and feedback control,” New Journal of Physics, vol. 23, p. 013007, 2021.
  * (27) J. Thompson, B. Zwickl, A. Jayich, F. Marquardt, S. Girvin, and J. Harris, “Strong dispersive coupling of a high-finesse cavity to a micromechanical membrane,” Nature, vol. 452, no. 7183, p. 72, 2008.
  * (28) D. Walls and G. J. Milburn, “Input–output formulation of optical cavities,” in Quantum Optics, p. 127, Springer, 2008.
  * (29) V. Giovannetti and D. Vitali, “Phase-noise measurement in a cavity with a movable mirror undergoing quantum brownian motion,” Physical Review A, vol. 63, no. 2, p. 023812, 2001.
  * (30) A. Mari and J. Eisert, “Gently modulating optomechanical systems,” Physical Review Letters, vol. 103, no. 21, p. 213603, 2009.
  * (31) E. X. DeJesus and C. Kaufman, “Routh-hurwitz criterion in the examination of eigenvalues of a system of nonlinear ordinary differential equations,” Physical Review A, vol. 35, no. 12, p. 5288, 1987.
  * (32) S. Huang and G. Agarwal, “Entangling nanomechanical oscillators in a ring cavity by feeding squeezed light,” New Journal of Physics, vol. 11, no. 10, p. 103044, 2009.
  * (33) L. Mandel and E. Wolf, Optical coherence and quantum optics.  Cambridge University Press, New York, 1995.
  * (34) T. E. Lee, C.-K. Chan, and S. Wang, “Entanglement tongue and quantum synchronization of disordered oscillators,” Physical Review E, vol. 89, p. 022913, 2014.