\documentclass[reqno]{amsart}
\usepackage{hyperref}
\AtBeginDocument{{\noindent\small
{\em Electronic Journal of Differential Equations},
Vol. 2007(2007), No. 133, pp. 1--20.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu
\newline ftp ejde.math.txstate.edu (login: ftp)}
\thanks{\copyright 2007 Texas State University - San Marcos.}
\vspace{9mm}}
\begin{document}
\title[\hfilneg EJDE-2007/133\hfil General solution]
{General solution of overdamped Josephson junction equation in the case of
phase-lock}
\author[S. I. Tertychniy\hfil EJDE-2007/133\hfilneg]
{Sergey I. Tertychniy}
\address{Sergey I. Tertychniy \newline
VNIIFTRI, Mendeleevo, Moscow Region, 141570, Russia}
\email{bpt97@mail.ru}
\thanks{Submitted May 25, 2007. Published October 12, 2007.}
\subjclass[2000]{33E30, 34A05, 34A25, 34B30, 34B60, 34M05, \hfill\break\indent
34M35, 70K40}
\keywords{Overdamped Josephson junction;
nonlinear first order ODE; \hfill\break\indent
linear second order ODE;
double confluent Heun equation; phase lock; general solution}
\begin{abstract}
The first order nonlinear ordinary differential equation
$\dot \varphi(t) + \sin\varphi(t)=B+A\cos\omega t$,
which is commonly used as a simple model of an overdamped
Josephson junction in superconductors is investigated.
Its general solution is obtained in the case known as
phase-lock where all but one solution converge to a
common `essentially periodic' attractor.
The general solution is represented in explicit form
in terms of the Floquet solution of a double confluent Heun equation.
In turn, the latter
solution is represented through the Laurent series
which defines an analytic function on the Riemann sphere with
punctured poles.
The series coefficients are given in terms of infinite products
of $2\times2$ matrices with a single zero element.
The closed form of the phase-lock condition is obtained and
represented as the condition for existence of a real root of
a transcendental function.
The efficient phase-lock criterion is conjectured, and
its plausibility is confirmed in numerical tests.
\end{abstract}
\maketitle
\numberwithin{equation}{section}
\newtheorem{theorem}{Theorem}[section]
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{corollary}[theorem]{Corollary}
\newtheorem{conjecture}[theorem]{Conjecture}
\section{Introduction}
The nonlinear first order ordinary differential equation
\begin{equation}
\dot \varphi(t) + \sin\varphi(t)=q(t), \label{eq1}
\end{equation}
is commonly used in applied physics
as the simple mathematical model describing the
basic electric properties of
Josephson junction (JJ) in
superconductors \cite{R1,R2}.
Here the right-hand side function $q(t)$,
assumed to be known, specifies the external impact
representing the appropriately normalized
\textit{bias current} (or simply \textit{bias})
supplied by an external current source.
The unknown real valued function
$\varphi(t)$ called {\it the phase\/}
describes the macroscopic quantum state of JJ.
In particular, it is connected with the instantaneous
voltage $V$ applied across JJ in accord with
the equation $V=(\hbar/2e)\mathrm{d}\varphi/\mathrm{d}\tau$,
where $\hbar$ is the Plank constant, $e$ is the electron charge,
$\tau$ is the (dimensional) current time.
The dimensionless variable $t$ entering \eqref{eq1}
is defined as $t=\omega_c\tau$, where
$\omega_c$ is a constant parameter
depending on the junction properties and
named JJ characteristic frequency.
See \cite{R5} for more details of JJ physics.
Equation \eqref{eq1} arises as the limiting case
of the second order ODE utilized in
more general Resistively Shunted Junction (RSJ)
model \cite{R3,R4}. The reduction is legitimate
if the role of the
junction capacitance proves negligible.
In practice, if JJ can be
described by \eqref{eq1} it is named {\it overdamped}.
Summarizing the aforementioned relationships, we shall name \eqref{eq1},
for brevity,
{\em overdamped Josephson junction equation} (\eqref{eq1}--\eqref{eq2}).
Under concordant conditions,
the theoretical modelling applying \eqref{eq1}--\eqref{eq2}
is in excellent agreement with experiments.
It is also worth noting that nowadays
electronic devices based on the Josephson effect
play the important role
in various branches of measurement technology.
In particular, JJ arrays serve
the heart element of the modern DC voltage standards \cite{R6}.
The development of JJ-based synthesizers of AC voltage waveforms
is currently in progress \cite{R7,R8,R9}.
These and other successful applications stimulate
the growing interest to
the theoretical study and the modelling of JJ properties including
investigation of capabilities of RSJ model
and its limiting cases
and the predictions they lead to.
To be more specific, the case most important
from viewpoint of applications and simultaneously distinguished by
the wealthiness of the underlaid mathematics
is definitely the one of the bias function representing harmonic oscillations.
Without loss of generality, it can be represented
in the form
\begin{equation}\label{eq2}
q(t)=B+A\cos\omega t,
\end{equation}
where $A,B,\omega>0$ are some real constant parameters.
Hereinafter, the abbreviation \eqref{eq1}--\eqref{eq2} introduced above
will refer to the couple (\ref{eq1}), (\ref{eq2}).
In spite of apparent simplicity of \eqref{eq1}--\eqref{eq2},
few facts of its specific analytic theory
had been available until recently.
Some preliminary results concerning
the problem of derivation of analytic solutions
of \eqref{eq1}--\eqref{eq2} in the general setting had been obtained in \cite{R11a}.
The approach put forward therein
is elaborated in the present work.
The focus is made on the case of manifestation
by \eqref{eq1} of the {\em phase-lock\/} property
which is one of its most important features from viewpoint of applications.
The phase-lock property is formalized as follows:
in the case of phase-lock
any solution $\varphi=\varphi(t)$ of
\eqref{eq1} either yields a
periodic exponent $e^{\mathrm{i}\varphi}$, exactly
two such distinguished solutions existing,
or, as the time parameter grows, the exponent
$e^{\mathrm{i}\varphi}$ exponentially converges
to the similar exponent for the one (common for all
$\varphi$'s) of the periodic functions just noted
(another one plays the role of repeller).
(Here and in what follows we shall not distinguish
phase functions which differ by a constant equal to $2\pi$ times an integer.)
The corresponding period coincides
with one of $q(t)$; i.e., in the case (\ref{eq2}),
$2\pi\omega^{-1}$.
This behavior is stable with respect to weak parameter perturbations;
i.e., the subset of parameter values leading to phase-lock is open.
It is worth noting for completeness that in the opposite
(no phase-lock) case
no {\it stable\/} periodicity
in the behavior of $e^{\mathrm{i}\varphi}$
is observed. There is also
a third, intermediate, type of the phase
behavior, where the asymptotic
attracting and repelling solutions are, in a sense, merged \cite{R11}.
It is realized on the lower-dimensional subset of the space of
parameter values.
In the present work, the complete analytic solution
of \eqref{eq1}--\eqref{eq2} is obtained
under assumption of the parameter choice ensuring phase-lock.
The closed form of the phase-lock criterion
in the form of the constraint imposed on the
problem parameters
is conjectured.
\section{Overdamped Josephson junction equation against
reduced double confluent Heun equation}\label{s1}
The analytic theory of \eqref{eq1}--\eqref{eq2}
can be based on its reduction to the following
system of two {\em linear\/} first order
ODEs, \cite{R11a},
\begin{equation}\label{eq3}
\begin{gathered}
4\mathrm{i}\omega z^2 \;x'(z)
= 2z x(z) +\left[2B z+A\left(z^2+1\right)\right] y(z),\\
-4\mathrm{i}\omega z^2\; y'(z)
=\left[2B z+A\left(z^2+1\right)\right] x(z)+ 2z y(z),
\end{gathered}
\end{equation}
where $z$ is the free complex variable.
Indeed, on the universal
covering $\Omega_1\simeq\mathbb{R}\ni t $
of the unit circle $S^1$ in $\mathbb{C}$;
i.e., for $z=\exp(\mathrm{i}\omega t)$,
any non-trivial solution of \eqref{eq3} determines
a solution of \eqref{eq1}--\eqref{eq2} in accordance with the
equation
\begin{equation}\label{eq4}
\exp(\mathrm{i}\varphi)=\frac{\Re x-\mathrm{i} \Re y}{\Re x+\mathrm{i} \Re y}
\end{equation}
($\Re$ means the real part)
supplemented with the continuity requirement.
[The fulfillment of \eqref{eq1}--\eqref{eq2}
follows from a straightforward computation
taking into account the equality $\overline z=z^{-1}$
holding true on $\Omega_1$
which is utilized for the demonstration that, for real $A,B,\omega$,
the functions
$\Re x,\Re y$ also verify \eqref{eq3} on $\Omega_1$.]
Conversely,
any real valued solution $\varphi(t)$ of \eqref{eq1}--\eqref{eq2}
induces through \eqref{eq4}
some `initial data' $x(0),y(0)$ (with arbitrary norm
$(x^2(0)+y^2(0))^{1/2}>0$)
for \eqref{eq3}. Having solved the latter
on $\Omega_1$ for $x,y$,
one obtains, applying \eqref{eq4},
another phase function obeying \eqref{eq1}--\eqref{eq2}.
It however must coincide,
due to the identical initial value assumed at $t=0$,
with the original $\varphi(t)$. Finally, since
these $x,y$ defined on $\Omega_1$ simultaneously
obey the linear ODEs (\ref{eq3})
with meromorphic coefficients which have
the only singular points $z=0$ and $z^{-1}=0$,
they can be extended, integrating (\ref{eq3}), say,
along radial directions, to analytic functions
defined on the whole $\Omega=\Omega_1\times \mathbb{R}_+$,
the universal covering of the Riemann sphere with the punctured
poles $z=0$ and $z^{-1}=0$.
It is worth noting that on $\Omega_1$
the real valued functions
$\tilde{x} =\tilde{x}(t)=\Re x(e^{\mathrm{i}\omega t}),
\tilde{y}=\tilde{y}(t)=\Re y(e^{\mathrm{i}\omega t})$
satisfy the equations
\begin{equation}\label{eq5}
2\mathrm{d} \tilde{x}/\mathrm{d} t = \tilde{x}+q \tilde{y},\quad
-2\mathrm{d} \tilde{y}/\mathrm{d} t = q \tilde{x}+ \tilde{y},
\end{equation}
leading, together with \eqref{eq4}, to the equation
$({\mathrm{d}/\mathrm{d} t})\log(\tilde{x}^2+\tilde{y}^2)=\cos\varphi$
which implies
\begin{equation}\label{eq5a}
\mathrm{const}_1\,e^{-t}\le
|\tilde{x}+\mathrm{i}\tilde{y}|^2\le \mathrm{const}_2\,e^t
\quad \mbox{for }t>0.
\end{equation}
Note that if $x,y\not\equiv 0$ then
$\mathrm{const}_1,\mathrm{const}_2$ may be assumed to be strictly positive.
We shall refer to these bounds later on.
The key observation enabling one
to radically simplify the problem of
description of the space of
solutions of \eqref{eq3} is as follows \cite{R11a}.
Let us introduce the analytic function $v(z)$
which satisfies the equation
\begin{equation} \label{eq6}
\Big[
z^2\frac{\mathrm{d}^2}{ \mathrm{d}^2 z}
+\left(
\mu (z^2+1)
-n z\right)
\frac{\mathrm{d}}{\mathrm{d} z}
+(2\omega)^{-2}
\Big]
v =0,
\end{equation}
where
\begin{equation} \label{eq7}
n= -\big(\frac{B}{\omega}+1\big),\quad
\mu=\frac{A}{2\omega},
\end{equation}
are the constants replacing original $A,B$
which will be used below whenever it proves convenient.
Then a straightforward calculation shows that
the functions $x,y$ determined by the equations
\begin{gather} \label{eq8}
v = \hphantom{(2\omega z)^{-1}}
\mathrm{i} z^{{n+1\over2}}
\exp \Big( \frac{1}{2}\mu \left(-z+z^{-1}\right) \Big) (x-\mathrm{i} y), \\
\label{eq9}
v' = (2\omega z)^{-1} z^{{n+1\over2}}
\exp\Big(\frac{1}{2}\mu \left(-z+z^{-1}\right)\Big)(x+\mathrm{i} y)
\end{gather}
verify \eqref{eq3}. Conversely, defining the
function $v(z)$ through the solution $x,y$ of \eqref{eq3}
in accordance with \eqref{eq8},
a straightforward computation proves satisfaction of
\eqref{eq9} and, then, \eqref{eq6} follows.
Thus, \eqref{eq3} are equivalent to (\ref{eq6}) and
\eqref{eq8}, \eqref{eq9} represent the corresponding one-to-one
transformation.
Equation \eqref{eq6} coincides,
after appropriate identification of the constant parameters,
with \cite[Eq.\ 1.4.40]{R14}.
It represents therefore a
particular instance of the double confluent Heun equation which can
be shown to be non-degenerated
in all cases where \eqref{eq1}--\eqref{eq2} is non-trivial.
It is also worth reproducing here the {\it canonical form\/}
of the ``generic'' double confluent Heun equation as it is given in
\cite[Eq. 4.5.1]{R14}. It reads
\begin{equation}\label{eq10}
z^2\frac{\mathrm{d}^2y}{\mathrm{d} z^2}
+(-z^2+c z+t)\frac{\mathrm{d} y}{\mathrm{d} z}+(-a z+\lambda)y=0
\end{equation}
where $a,c,t,\lambda$ are some constants. To adjust it to our case,
a single term has to be eliminated
setting $a=0$. For brevity, we shall name this subclass of double confluent
Heun equation's
{\it reduced.}
Besides, some obvious rescaling of the free variable $z$ is
to be carried out. After these, the three residuary constant parameters
exactly correspond to our constant parameters
$n,\mu,\omega$ (we shall not need and so omit the reproducing of
the concrete form of this transformation).
The general analytic theory of double confluent Heun equation is given
in the chapter 8 of the treatise \cite{R13}.
Its solutions are there represented,
up to nonzero factors given in explicit form,
through the Laurent series whose coefficients are assumed to be
computable through the `endless' chain of
3-term linear homogeneous equations (`recurrence relations').
In the present work, we derive the solution of reduced double
confluent Heun equation in a cognate but more explicit form.
As a technical limitation, we also stipulate
in the present work for the additional
condition to be imposed on the free constant problem parameters
claiming of them the ensuring of the phase-lock property.
On the base of practice of numerical computations,
it can be conjectured that
such parameter values fill up a non-empty open subset
({\em phase-lock area}) in the whole
parameter space (see also the Conjecture \ref{conjA} below).
The case where the parameters
belong to its complement
is left beyond the scope of the present work.
\section{Formal solution of reduced double confluent Heun equation}
Let us introduce yet another unknown function $E(z)$ replacing $v(z)$
by means of the transformation
\begin{equation} \label{eq11}
v(z)= z^{{n+\epsilon\over2}-\mathrm{i}{\kappa}} e^{-{\mu}z} E(z),
\end{equation}
where the discrete {\em parity parameter $\epsilon$\/}
may assume
one of the two values, either $\epsilon=0$ or $\epsilon=1$ (i.e.\ obeys
the equation $\epsilon^2=\epsilon$),
and $\kappa$
is some real positive constant
which will be determined latter.
For $v$ obeying (\ref{eq6}),
$E(z)$ verifies the equation
\begin{equation} \label{eq12}
\begin{aligned}
0&= z^3 E'' +z\left[ \left(\epsilon- 2\mathrm{i}{\kappa}\right) z
- \mu \left(z^2-1\right) \right] E' \\
&\quad +\Big[ \mu \big( {n-\epsilon\over2}+\mathrm{i} {\kappa}\big) z^2
+\Big((1-\epsilon)\big({1\over 4} + \mathrm{i} {\kappa}\big)
-\kappa^2 -\big({n+1\over2}\big)^2 +\lambda\Big) z \\
&\quad +\mu\big({n+\epsilon\over 2}-\mathrm{i} {\kappa}\big) \Big]E,
\end{aligned}
\end{equation}
where $\lambda=(2\omega)^{-2}-\mu^2$.
Conversely, (\ref{eq11}),(\ref{eq12}) imply
the fulfillment of \eqref{eq6}.
At first glance, \eqref{eq12} seems `much worse'
than the original
double confluent Heun equation representation.
Nevertheless,
it is this equation which
we shall
attempt to solve
searching for its solution in the form of Laurent series
\begin{equation} \label{eq13}
E=\sum^\infty_{k=-\infty} a_k z^{k}
\end{equation}
`centered' in the points $z=0$ and $z^{-1}=0$
(which are the only singular points for \eqref{eq12})
with unknown $z$-independent coefficients $a_k$.
Then, carrying out
straightforward substitution, one gets a sequence of 3-term
recurrence relations
binding triplets of neighboring $a$'s %series coefficients
which can be written down %represented
either as
\begin{equation}
\begin{aligned}
0&=- \mu \left(k-1-{n-\epsilon\over2}- \mathrm{i} {\kappa}\right) a_{k-1}\\
&\quad +(Z_{k}+\lambda) a_k
+ \mu \left(k+1+{n+\epsilon\over2} - \mathrm{i} {\kappa}\right) a_{k+1},
\end{aligned}
\label{eq14}
\end{equation}
where
\begin{equation}
Z_k= \left(k+{\epsilon-1\over2}
-\mathrm{i} {\kappa}\right)^2-\left({n+1\over2}\right)^2,
\label{eq15}
\end{equation}
or as
\begin{equation}
\begin{aligned}
0&= - \mu \left(k-1-{n+\epsilon\over2}+\mathrm{i} {\kappa}\right)a_{-(k-1)} \\
&\quad +(\tilde {Z_{k}}+\lambda) a_{-k}
+ \mu \left(k+1+{n-\epsilon\over2} +\mathrm{i} {\kappa}\right) a_{-(k+1)},
\end{aligned} \label{eq16}
\end{equation}
where
\begin{equation}
\tilde Z_k=\left(k+{1-\epsilon\over2}+
\mathrm{i} {\kappa}\right)^2-\left({n+1\over2}\right)^2
\label{eq17}
\end{equation}
and $k=0,\pm1,\pm2\dots$.
The sets of
\eqref{eq14}) and \eqref{eq16} are exactly equivalent
and each of them separately
covers the whole set of equations
the coefficients $a_k$ have to obey.
However, in the approach utilized in the present work,
we shall consider them both in conjunction, employing (\ref{eq14})
for coefficients $a_k$ with indices $k\ge -1$ and (\ref{eq16})
for $a_k$ with indices $k\le 1$.
Thus, \eqref{eq14} and (\ref{eq16})
will be considered separately but on the common index
variation `half-interval' $k\ge 0$
(remaining legitimate, in principle, for arbitrary integer $k$).
Obviously, these two equation sets
cover the complete set of conditions imposed to the coefficients $a_k$
and are `almost disjoined' coupling % intersecting
in their `boundary' $k=\pm0$-members alone.
Let us further consider for $k\ge 0$ the following {\it formal
infinite products\/}
\begin{equation} \label{eq18}
R_{k}=\prod_{j=k}^\infty M_j,\quad
\tilde R_{k}=\prod_{j=k}^\infty \tilde M_j
\end{equation}
of the $2\times2$ matrices
\begin{gather} \label{eq19}
M_j= \begin{pmatrix}
1+\lambda/ Z_{j} & \mu^2/ Z_{j} \\
1 & 0 \end{pmatrix}, \\ \label{eq20}
\tilde M_j= \begin{pmatrix}
1+\lambda/\tilde Z_{j} & \mu^2/\tilde Z_{j}
\\
{\tilde Z_{j-1}}/{\tilde Z_{j}}
& 0
\end{pmatrix}.
\end{gather}
It is assumed throughout that
the matrices $M_j,\tilde M_j$ with {\it larger\/} indices $j$
are situated in products {\it to the right\/} with respect to
ones labelled with {\it lesser\/} index values.
Notice that in the case $\kappa=0$ and integer $n$, zero may
appear in denominators of $M$-factors,
making the above definitions meaningless.
This apparent fault admits a simple
resolution (see \eqref{eq34} and the discussion following it).
For a while, we temporary leave out
consideration of such specific parameter choices.
It is also worth
noting that the above definitions of $R_{k},\,\tilde R_{k}$
may be understood as a concise form of representation of
the `descending' recurrence relations
\begin{gather}\label{eq21}
R_k = M_{k} R_{k+1},\\ \label{eq21a}
\tilde R_k = \tilde M_{k} \tilde R_{k+1},\quad k=0,1,\dots
\end{gather}
among the neighboring $R$'s.
These are the only dependencies which
will be actually
used below in derivations involving $R_{k},\tilde R_k$.
The formulas (\ref{eq18}) are `formal' since
neither the issue of the convergence of such sequences nor
how one should understand the `initial values'
$R_\infty, \tilde R_\infty$ necessary for the actual
determination of the `finite index value' $R$-matrices
by means of \eqref{eq21}--\eqref{eq21a}
are here addressed.
Now a straightforward calculation
applying \eqref{eq21}--\eqref{eq21a}
allows one to show that
the following formulas
\begin{gather} \label{eq22}
\tilde a_k
= \mu^k
{\Gamma\left(1+{n+\epsilon\over2}- \mathrm{i} {\kappa}\right)
\over
\Gamma\left(k+1+{n+\epsilon\over2}- \mathrm{i} {\kappa}\right)
}
(0,1)\cdot R_k \cdot
\begin{pmatrix}1\\ 0\end{pmatrix} \\
\label{eq23}
\tilde a_{-k}=
{\mu^k \over {\tilde Z_{k-1}}}
{\Gamma\left(1+{n-\epsilon\over2}+ \mathrm{i} {\kappa}\right)
\over \Gamma\left(k+1+{n-\epsilon\over2}+ \mathrm{i} {\kappa}\right) }
(0,1)\cdot \tilde R_k \cdot
\begin{pmatrix}1\\0 \end{pmatrix}
\end{gather}
(where $\Gamma$ is the Euler gamma-function)
yield {\it the formal solutions\/}
to \eqref{eq14} and \eqref{eq16},
respectively.
\section{Validation of the formal solution}
In this section we show that
the formal solution of \eqref{eq12}
presented in the form of expansion (\ref{eq13})
with coefficients
given by \eqref{eq22},\eqref{eq23}
represents a well defined analytic function of $z$.
This means, first of all,
that the infinite matrix
products $R_k, \tilde R_k$ it involves converge.
Moreover, the convergence takes place
for any constant parameter values.
The key auxiliary result
which may be utilized for the proof of
this assertion is as follows.
\begin{theorem} \label{thm1}
Let us consider the sequences of complex numbers
$\alpha_j,\beta_j,\gamma_j,\delta_j$ satisfying
the following `ascending'
matrix recurrence relation
\begin{equation} \label{eq24}
\begin{pmatrix}
\alpha_j &\beta_j\\
\gamma_j & \delta_j
\end{pmatrix}
= \begin{pmatrix}
\alpha_{j-1} &\beta_{j-1}\\
\gamma_{j-1} & \delta_{j-1}
\end{pmatrix}
\begin{pmatrix}
1+\lambda/ Z_{j} & \mu^2/ Z_{j} \\
\sigma_j & 0
\end{pmatrix},
\end{equation}
where
\begin{equation} \label{eq25}
\mbox{either $\sigma_j= 1$ or $\sigma_j= Z_{j-1}/ Z_{j}$.}
\end{equation}
Then all they converge as $j\to \infty$. Moreover,
$\lim\beta_j=0=\lim\delta_j$ whereas for $\alpha,\gamma$-sequences there exist
positive quantities $N_\alpha$, $N_\gamma$
and
positive integer $j_0$, all depending at most on
$n,\mu,\lambda,\kappa,\epsilon$,
such that for all $j> j_0$
\begin{equation} \label{eq26}
\begin{gathered}
|\alpha_{j}-\lim_{j'\to \infty} \alpha_{j'}|<
{N_\alpha
\max_{j'>j_0}(|\alpha_{j'}|) j^{-1}},\\
|\gamma_{j}-lim_{j'\to \infty} \gamma_{j'}|
< {N_\gamma \max_{j'>j_0}(|\gamma_{j'}|) j^{-1}},
\end{gathered}
\end{equation}
where the maxima are finite.
\end{theorem}
The outline of the proof can be found in Appendix.
\smallskip
\noindent{\it Remark}:
Formally, we need not include in the theorem stipulation the requirement that
either $\kappa\not =0$ or $n$ is non-integer (which would a priori evade
possibility arising of
contribution with vanishing $Z_{*}$ in
the denominator) because with fixed constant parameters
and sufficiently large $j_0$
no zero $Z_{*}$ may appear.
Let us return to \eqref{eq18} and consider the
four double-indexed sets of complex numbers
$\{\alpha,\beta,\gamma,\delta\}^{(j_0)}_j$
defined as follows:
\begin{equation} \label{eq26a}
\begin{pmatrix}
\alpha_{j}^{(j_0)} &\beta_{j}^{(j_0)}\\
\gamma_{j}^{(j_0)} & \delta_{j}^{(j_0)}
\end{pmatrix}
=R_{j_0}^{(j)}
=\prod_{k=j_0}^{j} M_k.
\end{equation}
It is straightforward to verify that the
sequences obtained by the picking out the elements with common
value of the upper index $j_0$ obey
the recurrence relations (\ref{eq24})
for the first choice in (\ref{eq25}). Hence it follows from the theorem that
all they converge.
We denote the corresponding limits as $\alpha^{(j_0)}$ etc.
We have therefore the consistent definition for the infinite matrix
products (\ref{eq18}),
\begin{equation} \label{eq26b}
R_{j_0}=
\begin{pmatrix}
\alpha^{(j_0)} &\beta^{(j_0)}\\
\gamma^{(j_0)} & \delta^{(j_0)}
\end{pmatrix}.
\end{equation}
Let us further note that, increasing $j_0$, the `starting' sequence
elements
($\alpha_{j_0}^{(j_0)},\alpha_{j_0+1}^{(j_0)}$ for $\alpha$ etc) tend
to the corresponding elements of the
idempotent matrix
\begin{equation} \label{eq26c}
M_{\infty}=
\begin{pmatrix}1&0\\1&0
\end{pmatrix},\quad M_{\infty}^2=M_{\infty},
\end{equation}
the discrepancy decreasing like
$O(j_{0}^{-2})$. On the other hand,
in accordance with (\ref{eq26}) the elements of
$R_{j_0}^{(j)}$ differ from their $j$-limits constituting $R_{j_0}$
by $O(j_{0}^{-1})$-order quantities. This means that
$R_{j_0}$ tends to $M_{\infty}$ as $j_0$ goes to
infinity with the difference going to zero as $O(j_{0}^{-1})$.
In other words, we have the following result.
\begin{corollary} \label{cor0}
1. $R_{j_0}-M_{\infty}=O(j_{0}^{-1})$.
Furthermore, in view of the convergence,
2. The modules of the elements of all the matrices $R_{j_0}$
are bounded from above in total
--- provided `no-zeroes-in~denominators' condition
for the parameter choice
is respected, of course.
\end{corollary}
In according with the above relations,
decomposing $R_{j}$ into the product of the two factors,
$R_{j}=R_j^{j_0}\cdot R_{j_0}$,
and approximating $R_{j_0} $ by $M_\infty$,
one obtains a simple but important algorithm of computation
of the products (\ref{eq18})
\begin{equation} \label{eq18a}
R_{j}\approx\prod_{k=j}^{j_0} M_k\cdot M_\infty.\;
\end{equation}
The approximation
is the better, the larger $j_0>j$ is selected.
In the limit, one gets
\begin{equation} \label{eq28a}
R_j= \begin{pmatrix}
\alpha^{(j)} &0\\
\gamma^{(j)} &0
\end{pmatrix}.
\end{equation}
This interpretation resolves
(in a quite obvious way, though)
the aforementioned uncertainty in %problem of
specification of the `initial value'
for the sequence $R_j$
treated through the `descending' recurrence
relation
`$[R_\infty]\dots\to R_j\to R_{j-1}\to \dots$'
implied by (\ref{eq21}).
After introducing a consistent representation of
the matrices (\ref{eq18}), it is straightforward
to do the same for the matrices $\tilde R_j$ (\ref{eq19}).
The above speculation applies to them with minor modifications as well.
The only distinction is the making use of the
second choice in (\ref{eq25})
and the operating with complex conjugated quantities
(equivalent in our case
to the replacing
$\kappa$ by $-\kappa$) throughout.
We shall mark
the elements of $\tilde R_j$ obtained in this way
with tildes over
the corresponding $\alpha$'s and $\gamma$'s.
Now, with well defined $R_j$ and $\tilde R_j$
(\ref{eq18}) in hand, one is able to
calculate the coefficients $a_k,a_{-k}$ for $k=0,1,2\dots$
in accordance with \eqref{eq22} \eqref{eq23}.
The triple matrix products reduce to
separate elements of $R$'s (or $\tilde R$'s)
denoted above as $\gamma^{(k)}$ (for (\ref{eq22})),
and $\tilde\gamma^{(k)}$(for \eqref{eq23}, respectively)
which are the functions of the parameters $n,\mu,\lambda,\kappa,\epsilon$.
Therefore, the sequences
\begin{gather} \label{eq31}
\tilde a_k = \mu^k \gamma^{(k)}
{
\Gamma\left(1+{n+\epsilon\over2}- \mathrm{i} {\kappa}\right)
\over
\Gamma\left(k+1+{n+\epsilon\over2}- \mathrm{i} {\kappa}\right)
} \\
\label{eq32}
\tilde a_{-k}
={\mu^k \tilde \gamma^{(k)}
\over {\tilde Z_{k-1}}}
{
\Gamma\left(1+{n-\epsilon\over2}+ \mathrm{i} {\kappa}\right)
\over
\Gamma\left(k+1+{n-\epsilon\over2}+ \mathrm{i} {\kappa}\right)
},\quad k=0,1,2,\dots
\end{gather}
are well defined and
solve \eqref{eq14}, \eqref{eq16}, respectively,
for $k=1,2,\dots$
The important feature of the expressions (\ref{eq31}), (\ref{eq32})
which is used below
is their
asymptotic behavior for
large values of the index
$k$ which is easy to derive in explicit form.
Specifically, in accordance with inequalities
(\ref{eq26}), the set of modules of
$\gamma$- and $\tilde\gamma$-factors
involved in \eqref{eq31}--\eqref{eq32}
is bounded in total from above and
each of these sequences
converge to a finite limit.
These imply in particular the validity of the estimates
\begin{equation} \label{eq31a}
|\tilde a_k| \propto {\zeta_1^k \over k!},\quad
|\tilde a_{-k}| \propto {\zeta_2^k \over k^2 k!},
\end{equation}
asymptotically, in the leading order,
for some $k$-independent
$\zeta_1,\zeta_2$.
\section{Matching condition and phase-lock criterion}
By construction, the
$\tilde a$-coefficients
defined by
\eqref{eq22} and \eqref{eq23}
obey the linear homogeneous equations
(\ref{eq14}) and (\ref{eq16}), respectively,
which are `the same', essentially, differing only
in the associated intervals of the variation of the index,
consisting of the
positive integers for the former and negative ones for the latter.
However, these two sequences
cannot be joined, automatically, since they
are `differently normalized'.
This means, in particular, that their `edge elements'
indexed with zeroes,
generally speaking, differ.
We will denote
them as $\tilde a_0$ and $\tilde a_{-0}$, respectively,
distinguishing here, in notations, the index `$0$' from the index `$-0$'.
Now, referring to \eqref{eq31},(\ref{eq32}),
one notes that
in view of the factors of $\Gamma$-functions present in denominators
and
leading to asymptotic behaviors (\ref{eq31a}),
the
following power series in $z$
\begin{equation} \label{eq33}
E_+(z)=1+\sum_{k=1}^\infty {\tilde a_k\over\tilde a_0}z^k,\quad
E_-(z)=1+\sum_{k=1}^\infty {\tilde a_{-k}\over\tilde a_{-0}}z^{-k},
\end{equation}
admit absolutely converging majorants.
(We assumed above $\tilde a_{\pm0}\not=0$.
Otherwise, i.e.\ if $\tilde a_{\pm0}=0$,
$\tilde a_{\pm1}$ may not vanish
and the series with the coefficients $\tilde a_{\pm k}/\tilde a_{\pm1}$
can be utilized instead.)
Indeed,
the Maclaurin series for the exponent function can play this role.
Therefore, the series $E_+(z)$ and $E_-(z^{-1})$ define some {\it entire
functions} of $z$.
As a consequence, the expression
\begin{equation} \label{eq34}
E(z)=
{{4\over\pi^2}\sin\left({\pi\over2}(n+\epsilon- 2\mathrm{i} {\kappa})\right)
\sin\left({\pi\over2}(n+\epsilon+ 2\mathrm{i} {\kappa})\right)
\over
(n+\epsilon- 2\mathrm{i} {\kappa})
(n+2-\epsilon+2\mathrm{i} {\kappa})
}
(E_+(z)+E_-(z)-1)
\end{equation}
represents a single-valued function of $z$
analytic everywhere on the Riemann sphere except the poles $z=0,z^{-1}=0$.
They are the essential singular points for $E$.
It has to be noted that the additional
$z$-independent fractional factor in (\ref{eq34})
may be regarded as a specific
common
`normalization' of the
(\ref{eq33})-type series which may be, in principle, arbitrary.
However, its given form is,
essentially, unique being fixed (up to
an insignificant nonzero numerical factor) in view of the following reasons.
The two sine-factors in the numerator
regarded as holomorphic functions of $n+\epsilon\pm 2\mathrm{i} \kappa$
are introduced for the cancelling out
zeroes in denominators arising
due to the poles in the factors $Z_{j}^{-1}$ and $\tilde Z_{j}^{-1}$
involved in the products
(\ref{eq21}) and regarded as the functions of the same parameters.
The set of these (vicious, essentially) singularities
constitute a homogeneous grid which is just covered by the grid of
roots of the sine-factors in (\ref{eq34}) --- with the two exceptions.
These two `superfluous' sine-factor roots are, in turn, `neutralized'
by the two linear factors in the denominator in (\ref{eq34})
which are therefore also uniquely determined.
As the result, in vicinity of
any zero in denominators in coefficients
of the power series defining
$E(z)$ considered as the function of $n+\epsilon\pm 2\mathrm{i} \kappa$
(a root of some $Z_{*}$ or $\tilde Z_{*}$),
the resulting expression takes the form of
the ratio $\sin x/x\;(x\simeq 0)$ and is
not now associated with any irregular behavior.
Thus, as a matter of fact,
the fractional factor involved in (\ref{eq34})
is distinguished (up to a numerical factor)
by the claims (i) to cancel out the poles in the
original expressions of the $\tilde a$-coefficients
(\ref{eq31}), (\ref{eq32}) considered as the functions of
$n+\epsilon\pm 2\mathrm{i} \kappa$ and
(ii) to introduce neither more zeroes nor more poles
as a result of such a `renormalization'.
Now, when plugging the function (\ref{eq34})
in \eqref{eq12}
in order to verify its fulfillment,
we may provisionally drop out, sparing the space,
$z$-independent fractional factor (restoring it afterwards).
It is important to emphasize that
the expressions (\ref{eq31}), (\ref{eq32})
verify, by construction, all the 3-term
recurrence relations (\ref{eq14}),(\ref{eq16})
which bind the $a$-coefficients with indices
{\em of a common sign}, either
non-negative or non-positive.
The only equation which does not fall into these categories,
and, accordingly,
has not been automatically fulfilled,
is the `central'
one binding the coefficients $a_{-1},a_0=1=a_{-0},a_1$; i.e.,
the equation
\begin{equation} \label{eq35}
\mu \left(1+{n-\epsilon\over2}+\mathrm{i} {\kappa}\right) a_{-1}
+(Z_{0}+\lambda) a_0
+ \mu \left(1+{n+\epsilon\over2}-\mathrm{i} {\kappa}\right)
a_{1}=0.
\end{equation}
With normalization adopted in
(\ref{eq33}), one has $a_0=1, a_1=\tilde a_1/\tilde a_0,
a_{-1}=\tilde a_{-1}/\tilde a_{-0}$.
Further,
in accordance with (\ref{eq31}),(\ref{eq32}):
\begin{gather*} % \label{eq36}
\tilde a_0 = \gamma^{(0)},\quad
\tilde a_1= {\mu\over 1+{n+\epsilon\over2}-\mathrm{i} {\kappa}}\gamma^{(1)} \\
\tilde a_{-0} = {\tilde\gamma^{(0)}\over\tilde Z_{-1}},\quad
\tilde a_{-1}= {\mu\over 1+{n-\epsilon\over2}+\mathrm{i}
{\kappa}}{\tilde\gamma^{(1)}\over\tilde Z_0}
\end{gather*}
Besides, one has $\gamma^{(0)}=\alpha^{(1)} $,
$\tilde\gamma^{(0)}=\tilde\alpha^{(1)}\tilde Z_{-1}/\tilde Z_{0} $.
Combining these dependencies,
the following representation of (\ref{eq35}) arises
\begin{equation}\label{eq37}
0 = \mu^2{\gamma^{(1)} \over\alpha^{(1)} }
+(Z_{0}+\lambda) + \mu^2
{\tilde\gamma^{(1)}\over \tilde\alpha^{(1)}}.
\end{equation}
Accordingly, it is convenient to
introduce the following function of
the parameters $\kappa$ and $n,\lambda,\mu,\epsilon$
\begin{equation} \label{eq38}
\begin{aligned}
\Xi&= {{4\over\pi^2}
\sin\left({\pi\over2}(n+\epsilon- 2\mathrm{i} {\kappa})\right)
\sin\left({\pi\over2}(n+\epsilon+ 2\mathrm{i} {\kappa})\right)
\over
(n+\epsilon- 2\mathrm{i} {\kappa})
(n+2-\epsilon+2\mathrm{i} {\kappa})
}\\
&\quad \times
\left( \mu^2{\gamma^{(1)}\tilde\alpha^{(1)}}
+(Z_{0}+\lambda)\tilde\alpha^{(1)}\alpha^{(1)}
+ \mu^2
{\tilde\gamma^{(1)}\alpha^{(1)}}
\right)
\end{aligned}
\end{equation}
where $Z_0=\left({(\epsilon-1)/2}
-\mathrm{i} {\kappa}\right)^2-\left({(n+1)/2}\right)^2$ (see (\ref{eq15}))
and $\alpha$'s, $\gamma$'s
are defined as the elements of the convergent matrix products
as follows
\begin{equation} \label{eq39}
\begin{pmatrix}
\alpha^{(1)} &0\\
\gamma^{(1)} &0
\end{pmatrix}
= \prod_{j=1}^\infty M_j,
\quad \begin{pmatrix}
\tilde \alpha^{(1)} &0\\
\tilde\gamma^{(1)} &0
\end{pmatrix}
= \prod_{j=1}^\infty \tilde M_j
\end{equation}
(see \eqref{eq18}-(\ref{eq20})).
The fractional multiplier in the first line of \eqref{eq38}
coincides with the one entering \eqref{eq34} and plays the
identical role: It eliminates the vicious singularities
arising for specific values of the parameters $n,\kappa$.
We shall name $\Xi=\Xi(\kappa,n,\mu,\lambda,\epsilon)\equiv\Xi(\kappa;\dots)$
the {\em discriminant function\/} for brevity. The following
statement holds true.
\begin{proposition} \label{prop1}
Restricting $\kappa$ to real values, the
vanishing of the discriminant function is necessary and sufficient
for the single valued analytic function \eqref{eq34}
to verify \eqref{eq12}
everywhere on the Riemann sphere except its poles $z=0,z^{-1}=0$.
\end{proposition}
Indeed, the vanishing of $\Xi$ implies the fulfillment
of \eqref{eq35} (where $a$'s are expressed through $\tilde a$'s),
the last equation among those binding coefficients of the expansion
(\ref{eq13}) which has not been fulfilled
as the result of the very coefficients definition.
Now all the 3-term recurrence relations for $a$-coefficients,
which \eqref{eq12} is equivalent to, are satisfied and the analytic
function (\ref{eq34}) verifies \eqref{eq12} everywhere
except of its own singular points $z=0$ and $z^{-1}=0$.
The equation
\begin{equation} \label{eq40}
\Xi(\kappa;\dots) =0
\end{equation}
referred to in the above proposition
can be named {\it the matching condition\/} since it enforces
the sequences of the coefficients $a_k$,
$a_{-k}$, separately obeying the corresponding `halves'
of the equation chain (\ref{eq14}) (equivalently, (\ref{eq16}))
to be `matched'
in their `edge' elements $a_{\pm0}=1,a_{\pm1}$.
It is worth emphasizing that,
up to this point, the parameter $\kappa$ (absent
in \eqref{eq6})
has not been restricted in
any way (apart of the claim to be real).
Now and in what follows
we regard the condition (\ref{eq40})
as $\kappa$ definition eliminating
this odd `degree of freedom'.
Hereinafter, it is a well defined function of the other parameters.
It seems interesting enough that
the addition of unspecified constant
$\kappa$ to the transformation (\ref{eq11}) and its subsequent
`fine tuning' by means of the claim of fulfillment of
\eqref{eq40}) is necessary
for the representation of
solution of double confluent Heun equation (\ref{eq6}) in terms of
convergent Laurent series.
More precisely, it is clear that \eqref{eq15} can be solved
{\em for any\/} (including trivial zero) choice
of $\kappa$,
choosing loosely $a_{j_0},a_{j_0+1}$ for
arbitrary $j_0$ and then
calculating, term by term, all the other coefficients
$a_j$, advancing in parallel in both
directions of $j$-index variation `off $j_0$ towards $\pm\infty$'.
Then (\ref{eq13}) immediately yields a ($\kappa$-dependent!)
formal solution to \eqref{eq12} and hence,
through transformation (\ref{eq11}),
to \eqref{eq6}.
However, {\em it can be only formal and
will necessarily diverge\/} for any $z$ unless
the matching condition (\ref{eq40}) is fulfilled -- just in view of the
uniqueness of solution with the analytic properties presupposed.
On the other hand,
considering separately the `halves' of the set of \eqref{eq15}
and resolving them `in index variation directions'
opposite to the ones assumed above
(in a sense, `off $\pm\infty$ towards $\pm0$'),
we obtain the {\em always converging\/}
series (\ref{eq33}). However, as we have seen,
we again have no solution, in this case even formal,
unless the matching condition
fixing $\kappa$
is fulfilled. Obviously, the uniqueness property
implies that the introduction of the
`branching'
power function factor,
as in (\ref{eq11}), is the only way to obtain a solution to \eqref{eq6}
admitting representation in terms of convergent power series.
Now,
tracking back the relationship connecting \eqref{eq12}
with the primary \eqref{eq1}
and invoking the general theory of the latter
applicable in the case of arbitrary continuous periodic $q(t)$ \cite{R11},
one can infer the following statements
which however, in the present context, are only of
the status of {\em conjectures\/} in view of the lack of
their proof `from the first principles';
i.e., on the base of the properties of the discriminant function $\Xi$
following from its explicit definition.
\begin{conjecture} \label{conjA}
1. There exists an open non-empty subset of the space of the
problem parameters $n,\mu,\omega,\epsilon$ where
\eqref{eq40}
admits a real valued positive solution.
2. If real $\kappa$ solves \eqref{eq40}, $-\kappa$
is the solution as well. No more real roots exist.
3. Real roots of \eqref{eq40} obey the condition
$|\kappa|\le (2\omega)^{-1}$.
\end{conjecture}
{\it Remark}:
The last statement
which is very important from viewpoint of applications,
is nothing else
but the form of the limitation on the rate of
growth or decreasing of the functions $\tilde{x},\tilde
y=x,y|_{z=e^{\mathrm{i}\omega t}}$ of the real variable $t$
implied by the inequalities (\ref{eq5a}) and the equation
$$
x-\mathrm{i} y=-iz^{{\epsilon-1\over2}-\mathrm{i} \kappa}
e^{-{\mu\over2}\left(z+{1\over z}\right)} E,
$$
following from definitions.
Notice that the latter clarifies the role of the {\em parity\/}
parameter $\epsilon\; (=0\;\mbox{or}\; 1)$
which
determines the multiplicity of the inverse to the map
$S^1\ni z\to (x+\mathrm{i} y)/|x+\mathrm{i} y|\in S^1$ induced
by the solution (\ref{eq34}).
If $\epsilon=0$, the revolution along the circle $|z|=1$
leads to the reversing of the direction of the vector
with components $(x,y)$
whereas in the case $\epsilon=1$ its direction is preserved.
The above inverse map is double-valued in the former case and one-to-one
in the latter one.
More properties of the discriminant function
can be inferred from the
numerical experiments although, as opposed to the assertions
of the Conjecture \ref{conjA}, they
have, to date, no analytic arguments in their favor yet, even indirect.
Nevertheless, the first item below is
important enough from viewpoint of applications
(seeming also plausible enough)
to be
{\em explicitly formulated\/} here.
\begin{conjecture} \label{conjB}
1. {\bf Phase-lock criterion.}
Equation \eqref{eq40} admits
a real non-zero solution if and only
if
\begin{equation}
\Xi(0;\dots)>0.
\end{equation}
This means in particular that $\Xi(0;\dots)$ is real. Moreover,
the numerical study makes evidence that
2. $\Xi(\kappa;\dots)$ is real for real $\kappa$
(assuming the other
parameters to be also real, of course).
\end{conjecture}
\section{Floquet solutions of double confluent Heun equation and involutive
solution maps}
Let us assume now that there exists a real positive solution $\kappa$ of
\eqref{eq40}. With this $\kappa$, the function $E(z)$
defined by \eqref{eq34} verifies \eqref{eq12}.
Let us consider the function $E_\#(z)$ defined through $E(z)$
as follows
\begin{equation}
\label{eq41}
E_\#(z)= z^{2\mathrm{i}{\kappa}-\epsilon}
\Big[ E'\big({1\over z}\big)
+\Big(\big({n+\epsilon\over2}-\mathrm{i}{\kappa}\big)z-\mu\Big)
E\big({1\over z}\big)\Big].
\end{equation}
Then a straightforward computation shows that it verifies
\eqref{eq12}, provided $E(z)$ does.
It is worthwhile to note that the repetition of the
transformation (\ref{eq41}) yields no more solutions to \eqref{eq12}.
As a matter of fact, one has $_{\#}\circ_{\#}=(2\omega)^{-2} Id$.
Thus $(2\omega)^{-1}{}_\#$ is the involutive map on the space
of its solutions.
Next, the functions $E(z)$
and $E_\#(z)$ are linearly independent for nonzero
$\kappa$. Indeed, utilizing (\ref{eq12}), one finds
%
\begin{equation} \label{eq42}
\begin{aligned}
E_\# '(z) &= z^{ 2\mathrm{i}{\kappa}-\epsilon-1}
\Big[ \Big(-{n+\epsilon\over2}+\mathrm{i}{\kappa} +\mu z\Big)
E'({1\over z}) \\
&\quad +\Big(
\mu\big({n+\epsilon\over2}-\mathrm{i}{\kappa}\big)(z^2+1)
+
\big( -({n+\epsilon\over2}-\mathrm{i}{\kappa})^2
+\lambda \big)z
\Big) E({1\over z}) \Big].
\end{aligned}
\end{equation}
This reduction allows one to calculate
the determinant of the linear transformation
binding the pairs of functions $E_\#,E_\# '(z)$ and $E,E'(1/z)$
which proves
equal to $(2\omega)^{-2}z^{2(-\epsilon + 2\mathrm{i}\kappa)}$ and is therefore
nonzero. Hence $E_\#(z)$ is not zero, identically,
(and may vanish at isolated points, at most, as well as
$E(z)$). Finally,
$E(z)$ is periodic on the unit circle centered at zero whereas $E_\#(z)$,
for real $\kappa\not=0$, is not. Hence they are linear independent.
The functions
$E(z)$ and $E_\#(z)$ constitute therefore the fundamental
system for \eqref{eq12} and any its solution %of the latter
can be expanded in this basis with constant expansion coefficients.
Consequently, for $\kappa\not=0$,
the domain of general solution to \eqref{eq12}
is $\Omega$, the universal covering of Riemann sphere with punctured poles.
The analytic properties of $E,E_\#$ identify these solutions
as (the unique pair of) the
{\it Floquet solutions\/} of the reduced double confluent Heun equation
under consideration. See \cite[section 2.4]{R13}.
The function $E(z)$ obeys the important functional equation
which can be derived as follows.
A straightforward calculation
shows that the right-hand side expression of \eqref{eq41}
with the `branched' factor
$z^{2\mathrm{i}{\kappa}}$ removed (i.e.,
$z^{-2\mathrm{i}{\kappa}} E_\#(z)$)
satisfies the ODE which coincides with (\ref{eq12})
{\em up to the opposite sign of the parameter $\kappa$}. This means
that, for real $\kappa$, $n,\mu,\omega$,
the analytic function
\begin{equation} \label{eq43}
\hat E(z)=z^{-\epsilon}
\Big[ \bar E'\big({1\over z}\big)
+\Big(\big({n+\epsilon\over2}+\mathrm{i}{\kappa}
\big)z-\mu\Big)\bar E\big({1\over z}\big) \Big],
\end{equation}
where $\bar E(z)=\overline{E(\overline z)}$,
is the solution of \eqref{eq12} itself.
As opposed to $E_\#$, this `yet another' solution has
{\it the same\/} analytic
properties as $E(z)$ and
hence must coincide with it
up to some numerical factor $C_C$; i.e.,
\begin{equation} \label{eq44}
C_C E(z)=
z^{-\epsilon}
\Big[ \bar E'\big({1\over z}\big)
+\Big(\Big({n+\epsilon\over2}+\mathrm{i}{\kappa}
\big)z-\mu\Big)\bar E\big({1\over z}\big) \Big].
\end{equation}
($C_C$ may not vanish since otherwise $E_\#$
would also be identical zero.) This is the generalization of the
similar property of the so
called `Heun polynomials' established in Ref.~\cite{R11a}.
The complex valued constant $C_C$
actually reduces to a single real constant.
To show that, let us
notice at first that if $E(z)$ verifies \eqref{eq12},
it follows from the latter and (\ref{eq44}),
\begin{equation} \label{eq45}
\begin{aligned}
C_C E'(z)&= z^{-\epsilon-1}
\Big[ \big( -{n+\epsilon\over2}-\mathrm{i}{\kappa}
+\mu z \big)\bar E'(\frac{1}{z}) \\
&\quad + \Big( \mu {\epsilon+n\over2}(z^2+1)
-\mathrm{i}\mu{\kappa}(z^2-1)\\
&\quad +\big( -{(n+\epsilon)^2\over4}
-{\kappa^2}+\lambda \big)z \Big) \bar E(\frac{1}{z})
\Big].
\end{aligned}
\end{equation}
Evaluating \eqref{eq44}--\eqref{eq45}
together with their complex conjugated
versions at the point $z=z^{-1}=1$, one
obtains four linear homogeneous equations binding
the quantities $E(1),E'(1),\bar E(1),\bar E'(1)$ which may not
vanish simultaneously. The corresponding consistency condition
reads $ |C_C|^2=(2\omega)^{-2}$ implying
\begin{equation} \label{eq46}
C_C=(2\omega)^{-1}e^{\mathrm{i} C_c},
\end{equation}
where $C_c$ is some {\em real\/} constant (actually, the function of the
parameters $n,\mu,\omega,\epsilon$). It
encodes all the monodromy data for \eqref{eq12}, essentially.
It is straightforward to
show that the transformation
(\ref{eq44})
is also involutive.
It manifests the specific symmetry
in the behaviors of the function $E(z)$ in vicinities
of the essentially singular points $z=0$ and $z^{-1}=0$.
Remarkably, the possessing of
this symmetry suffices itself for the fulfillment of \eqref{eq12}.
Indeed,
differentiating (\ref{eq44}) and taking into account (\ref{eq46}),
one arrives at \eqref{eq12}. In a sense, (\ref{eq44})
together with stipulation for the analyticity of $E(z)$
can be considered as the equivalent to \eqref{eq12}.
Additionally, (\ref{eq44})
implies anti-linear (involving complex conjugation)
dependencies among the `distant'
Laurent series
coefficients $a_{-k}$ and $a_{k},a_{k+1}$ (whereas
(\ref{eq14}) (or (\ref{eq16}))
binds `nearby' $a_k,a_{k\pm1}$).
In particular, it suffices to find all $a_k$ for $k>0$
and then $a_{-k}$ can be computed
from the latter by means of a simple (anti-)linear transformation.
\section{Essentially periodic and general
solutions of overdamped Josephson junction equations}
The connection between the functions $E(z),E_\#(z),\hat E(z)$
pointed out above is important for the consistent projecting
the results concerning solutions of \eqref{eq12}
to the level of original \eqref{eq1}--\eqref{eq2}.
This procedure applying equations
\eqref{eq4}, (\ref{eq8}), (\ref{eq9}), (\ref{eq11})
is straightforward and leads to the following conclusions.
At first, the representation of the two special (and the most important)
solutions to \eqref{eq1}--\eqref{eq2} for which the
exponents $\exp(\mathrm{i}\varphi)$ are {\it periodic\/}
(for brevity, we shall call such phase functions {\em essentially periodic})
follows. It reads
\begin{gather} \label{eq47}
\exp(-\mathrm{i}\varphi)
= 2\mathrm{i}\omega \Big(z {E'(z)\over E(z)}
+{n+\epsilon\over2}-\mathrm{i}{\kappa}-\mu z\Big),\\
\label{eq48}
\exp(\mathrm{i}\varphi)=-2\mathrm{i}\omega
\Big(z^{-1}{E'(z^{-1})\over E(z^{-1})}+{n+\epsilon\over2}-\mathrm{i}{\kappa}
-\mu z^{-1}\Big),
\end{gather}
where
\begin{equation}
z=\exp(\mathrm{i}\omega t). \label{eq48a}
\end{equation}
For $\kappa>0$, the first of these formulas
determines the asymptotic
limit (the attractor) of a generic solution to \eqref{eq1}--\eqref{eq2}
whereas the second solution
is unstable (the repeller).
It is important to emphasize that the functions $\varphi(t)$
defined by \eqref{eq47}, \eqref{eq48},
and \eqref{eq48a} are {\em real\/} and
\eqref{eq44}, \eqref{eq46} are the crucial ones utilized in
the calculation establishing this fact.
At second, it is straightforward to carry out
the `nonlinear superposition' of solutions (\ref{eq47}), (\ref{eq48})
operating with linear problem
-related counterparts.
The result is represented by the formula
\begin{equation} \label{eq49}
\begin{aligned}
&\exp(\mathrm{i}\varphi) \\
&= -{\mathrm{i}\over2} \Big\{\cos\psi\cdot E(z)
+\sin\psi\cdot z^{-\epsilon+2\mathrm{i}{\kappa}}
\Big[ E'(z^{-1})
+\Big( \big({n+\epsilon\over2}-\mathrm{i}{\kappa} \big)z -\mu \Big)E(z^{-1}) \Big]
\Big\} \\
&\times
\Big\{\omega\cos\psi\cdot
\Big[ z E'(z)+\big({n+\epsilon\over2}-\mathrm{i}{\kappa}-\mu z\big)E(z)\Big]
+
{1\over4\omega} \sin\psi\cdot
z^{-\epsilon+1+2\mathrm{i}{\kappa}}
E(z^{-1})\Big\}^{-1},
\end{aligned}
\end{equation}
where $\psi$ is an arbitrary real constant.
More exactly, the set of all functions $\varphi$ %solutions
described by (\ref{eq49}) is parameterized by points of
$S^1$. %the unit circle.
As opposed to (\ref{eq47}) and (\ref{eq48})
defined on the Riemann sphere with punctured poles,
the function (\ref{eq45}) is defined on its universal covering, $\Omega$.
Continuous (and then real analytic) function $\varphi(t)$ determined
by this equation on $\Omega_1\in\Omega$, where $z$ is understood
as $e^{\mathrm{i}\omega t}$, is just
the general solution of \eqref{eq1}--\eqref{eq2}
in the case of phase-lock.
In particular,
\eqref{eq47}, \eqref{eq48} represent the special
instances of \eqref{eq49} arising
for $\psi=0$ and $\psi=\pi/2$, respectively.
As a consequence,
asymptotic properties of the general solution mentioned above
immediately follow.
Indeed,
as $t$ increases, the exponent (\ref{eq49})
is converging to (\ref{eq47}) and is taking off
(\ref{eq48}) (unless it coincides with the latter).
The two solutions described by
(\ref{eq47}),(\ref{eq48}) are the only
ones which are not affected by the translations %shifts
$t\to {}t+2\pi\omega^{-1}$
(in the sense %that
the exponents (\ref{eq47}),(\ref{eq48}) are kept unchanged)
and preserve their functional form in asymptotics.
Finally, at third, considering $\varphi$ defined by (\ref{eq47})
as analytic function of
$z$ and taking in account \eqref{eq12},
one obtains
\begin{equation} \label{eq50}
\begin{aligned}
{\mathrm{d}\varphi\over\mathrm{d} z}
&= -\mathrm{i} z^{-2} \Big\{
z^3 \Big({E'(z)\over E(z)}\Big)^2
+z\big((1-z^2)\mu +z\big((\epsilon-1)
-2\mathrm{i}{\kappa}\big)\big) {E'(z)\over E(z)} \\
&\quad -\Big(z+z\big({n-\epsilon\over2}+\mathrm{i}{\kappa}\big)-\mu\Big)
\Big(\big({n+\epsilon\over2}-\mathrm{i}{\kappa}\big)-\mu z\Big)
+{z\over4\omega^2}\Big\} \\
&\quad \times \Big\{
z {E'(z)\over E(z)}
+\big({n+\epsilon\over2}-\mathrm{i}{\kappa}\big)
-\mu z \Big\}^{-1}
\end{aligned}
\end{equation}
On the unit circle, this $\varphi$ verifies \eqref{eq1}--\eqref{eq2}.
It is therefore smooth (even real analytic) therein.
Hence (\ref{eq50}) is continuous on $S^1$.
Finally, since $\exp\mathrm{i}\varphi|_{z=\exp\mathrm{i}\omega t }$ is
periodic, the following proposition holds.
\begin{proposition} \label{prop2}
The quantity
\begin{equation} \label{eq51}
k=(2\pi)^{-1} \oint {\mathrm{d}\varphi\over\mathrm{d} z} \mathrm{d} z,
\end{equation}
where ${\mathrm{d}\varphi/\mathrm d\, z}$ denotes the right-hand side
in \eqref{eq50} and the
integration is carried out over the circle $|z|=1$,
is well defined and amounts to an
integer.
\end{proposition}
This integer is the degree
of the map ${S}^1\Rightarrow {S}^1$ induced by the function (\ref{eq47}).
In physical applications, it is called {\it the phase-lock order\/}
and is considered as the integer-valued function of the parameters
which is obviously continuous and, thus, locally constant.
Phase-lock order is involved in the formula representing the
property of being `essentially periodic' for
the phase function defined by (\ref{eq47})
(and asymptotically for a generic phase function)
which reads
$$
\varphi(t+2\pi\omega^{-1})=\varphi(t)+2\pi k \quad \forall t.
$$
In a phase-lock state of JJ,
the uniformly distributed discrete
levels of averaged voltage equal to
$k\cdot (\hbar\omega/2e)$ for some $k=0,\pm 1,\pm 2\dots$
are observed.
\begin{conjecture} \label{conjC}
Any integer map degree \eqref{eq51} is realized on some
non-empty open subset of the space of the
problem parameters $n,\mu,\omega,\epsilon$ (or, equivalently,
$A,B,\omega,\epsilon$).
\end{conjecture}
This assertion is closely cognate to the item 1 of the above
Conjecture \ref{conjA}.
\section{Summary}
It the present work, the general solution of the
overdamped Josephson junction equation
(\ref{eq1}) was obtained
for the (co)sinusoidal right-hand side function (\ref{eq2})
in the case of one of three possible asymptotic
behaviors known as the phase-lock mode.
The solution is represented in explicit form
in terms of the Floquet solution
of the particular instance (arising in case of the vanishing
of one of the four free constant parameters)
of the
double confluent Heun equation.
The Floquet solution of double confluent Heun equation is represented in terms of
the Laurent series whose coefficients are determined by
explicit formulas
involving
convergent
infinite products of $2\times 2$ matrices with a single
zero element tending to the idempotent matrix (\ref{eq26c}).
The derivation % method
presupposes the existence
of a real solution of the transcendental
equation
(\ref{eq40})
which is equivalent
to the claim of realization of the phase-lock mode
for the given parameter values.
The plausible criterion
of its existence (i.e., {\em the phase-lock criterion})
is conjectured.
The derivation of general solution of \eqref{eq1}--\eqref{eq2}
consists of a number of steps.
In a concise form they can be described as follows.
\begin{itemize}
\item The investigation of the basic properties
of \eqref{eq1}
for arbitrary periodic (sufficiently regular)
$q(t)$ allows one to establish the division of the
space of the problem parameters
into the two open subsets of which one corresponds
to the phase-lock property
of the \eqref{eq1}--\eqref{eq2} solutions
whereas another corresponds to their pseudo-chaotic
behavior revealing no %apparent
stable periodicity.
In the (lower-dimensional) complement to these areas the
intermediate behavior is observed.
The corresponding results are discussed in
sufficient details
in Ref.~\cite{R11}.
\item The next important point
is the intimate connection
(first mentioned by Buchstaber; see, e.g.,
Ref.~\cite{R12})
between \eqref{eq1}
and a simple {\em linear} system
of the two first order ODEs
(\ref{eq5}). For the (co)sinusoidal right-hand side function (\ref{eq2}),
the latter
takes the form (\ref{eq3}).
\item At the next step, the transformation
(\ref{eq8}) was found which converts the linear system (\ref{eq3})
to a
particular instance of the double confluent Heun equation{} (\ref{eq10}).
\end{itemize}
Generally speaking,
it could be solved by means of the expansion
of the (modified) unknown function in Laurent series
\cite{R13},\cite{R11a}
centered at the singular points $z=0,z^{-1}=0$
but preliminarily
the additional simple but important transformation
has to be carried out:
\begin{itemize}
\item the `branched' power factor involving unspecified constant
($\kappa$-dependent contribution in (\ref{eq11}))
is introduced.
\item Additionally, the
discrete `parity' parameter $\epsilon$,
assuming either the value 0 or the value 1 (also absent in the
equation to be solved),
is introduced in the aforementioned power factor.
This modification
proves necessary for the subsequent exhaustive
`indexing' of the solution space.
\item After that, the standard technique of the
power expansion leads to the `endless'
sequence of the 3-term constraints (\ref{eq14})
(or, which is the same, (\ref{eq16}))
imposed on the
unknown series coefficients.
\item The next step is the devision
of the set of power series coefficients into the two subsets.
The
non-negative-index-value coefficients and
non-positive-index-value ones
are treated separately, solving the
separate subsets of the equations (\ref{eq14}) and (\ref{eq16}),
respectively, for $k\ge 1$.
Analyzing them, the application of
the continued fraction technique leads, after some transformations,
to the `explicit' formulas (\ref{eq22}),\eqref{eq23}
for the series coefficients
involving infinite products of
$2\times 2$ matrices converging for large index values to the
idempotent matrix (\ref{eq26c}).
This convergence is sufficiently fast to imply the
convergence of the matrix products and, accordingly, the finiteness of
the series coefficients. Moreover,
the associated estimates
make evident the existence of
the absolutely converging
majorants for the resulting Laurent series.
Therefore, they actually determine
an analytic function serving
the Floquet solution of double confluent Heun equation.
It proves
representable as the sum of the two entire functions of the
arguments $z$ and $z^{-1}$, respectively.
\item The procedure producing Laurent series coefficients
noted above proves however suffering
from the improper introduction of a kind of
vicious singularities arising as zeroes in denominator
which appear for some special parameter values.
They are eliminated my means of multiplication
of the `raw' coefficient expressions
by some $z$-independent (but parameter dependent)
factors given in explicit form.
\end{itemize}
One more {\em remark\/} is here in order.
The case where the above `regularizing' factor actually enters
the play exactly corresponds to the situation
where the problem solution can be built in terms of polynomials
(it was discovered
in Ref.~\cite{R11a}).
In brief, the `vanishing' factors `neutralizing' diverging
terms with vanishing
denominators, eliminate, in pass, the terms which do not
diverge. Ultimately, only finite number of `originally diverging'
terms survive which lead to
the reduction of the transcendental function
describing generic solution to a
polynomial. (The details of the realization of
this reduction have not been properly elaborated yet.)
\begin{itemize}
\item Now the `solution candidate' for \eqref{eq12}
can be represented as the analytic function (\ref{eq34})
which is well-defined for any parameter values.
However, at the price of automatic convergence
of the power series it has been built upon,
generally speaking,
it does not verify \eqref{eq12}.
The equation is fulfilled if and only if the
fulfillment of \eqref{eq40},
which is the transcendental equation
for the still unspecified parameter $\kappa$,
is stipulated.
\end{itemize}
At this stage, having solved \eqref{eq40},
a single solution (the Floquet solution) of double confluent Heun
equation can be regarded as been
explicitly constructed.
\begin{itemize}
\item The invariance of the space of solutions of double confluent Heun equation
under consideration
with respect to the transformation (\ref{eq41})
allows one to immediately obtain the fundamental system of
its solutions in terms of the single Floquet solution
noted above.
\item The existence of automorphism represented by
(\ref{eq44}) expresses the important intrinsic property
of the Floquet solution of double confluent Heun equation.
It is used
for the derivation of
the explicit representation of
the exponent $\exp(\mathrm{i}\varphi)$
specifying the {\em real valued\/}
phase function $\varphi$.
It is obtained by means of
the restriction
of the analytic function (see (\ref{eq49}))
from the universal covering
of the Riemann sphere with punctured poles
to the embedded universal covering of $S^1$.
This trick yields
the general solution of
\eqref{eq1}, \eqref{eq2} {\em in the case of phase-lock}.
\item Finally, employing the analytic properties and the
periodicity (on $S^1$) of $\exp(\mathrm{i}\varphi)(z)$,
the formula (\ref{eq51}) involving Floquet solution of
double{} confluent Heun equation follows which gives the
degree of the map ${S}^1\Rightarrow {S}^1$
it induces (the winding number)
also known in application fields as the phase-lock order.
\end{itemize}
It is worth noting in conclusion that
all the major constructions derived above
admit a straightforward
algorithmic implementation
which have been used for the numeric verification
of the relevant relationships.
\section{Appendix: Outline of the proof of Theorem \ref{thm1}}
Equation (\ref{eq24}) implies
$ \beta_j=\mu^2\alpha_{j-1} Z_{j}^{-1}$,
$\delta_j=\mu^2\gamma_{j-1} Z_{j}^{-1}$
and hence the asserted properties of the sequences $\beta_j,\delta_j$
follow from
the existence of finite limits for the sequences $\alpha_j,\gamma_j$.
As to the sequences
$\alpha_j$ and $\gamma_j$, their elements have
to obey the identical decoupled 3-term recurrence relations
which, for $\alpha$'s, read
\begin{equation} \label{eq52}
\alpha_j=
(1+\lambda Z_{j+(\epsilon-1)/2}^{-1})\alpha_{j-1}
+\mu^2 Z_{j-\tilde\epsilon+(\epsilon-1)/2}^{-1}\alpha_{j-2},
\end{equation}
where $\tilde\epsilon=1$ for the upper choice of $\sigma$
in (\ref{eq25}) and
$\tilde\epsilon=0$ for the lower $\sigma$ choice therein.
It suffices to consider the $\alpha$-sequence case.
Evidently, for every integer $j_0>0$ and $l>0$
any solution of equations (\ref{eq52})
can be represented in terms of the
decomposition
\begin{equation}
\label{eq53}
\alpha_{j_0+l}=(1+p_{j_0,l})\alpha_{j_0-1}+q_{j_0,l}\alpha_{j_0-2},
\end{equation}
for some coefficients
$p_{j,l},q_{j,l}$
independent on the `starting' terms $\alpha_{j_0-1},\alpha_{j_0-2}$.
Applying (\ref{eq52}),
it is straightforward to derive the reduction
\begin{equation} \label{eq:53a}
q_{j_0,l+1}= (1+p_{j_0+1,l} )\mu^2 Z_{j_0-\tilde\epsilon}^{-1},
\end{equation}
whereas for $p_{*,*}$ one gets the following recurrence relation:
\begin{equation} \label{eq54}
p_{j_0,l+1}= p_{j_0+1,l}+
\lambda(1+p_{j_0+1,l}) Z_{j}^{-1}
+\mu^2(1+p_{j_0+2,l-1})Z_{j+1-\tilde\epsilon}^{-1}.
\end{equation}
Equation (\ref{eq54}) is equivalent to
(\ref{eq52}), essentially, but it possesses
the advantage of being endowed
with the standard `initial conditions'
\begin{equation} \label{eq55}
p_{j_0,-1}= 0, p_{j_0,-2}=-1
\end{equation}
independent of the specific `starting' values
$\alpha_{j_0-1},\alpha_{j_0-2}$.
Furthermore, one gets
\begin{equation} \label{eq56}
q_{j_0,-1}= 0, q_{j_0,-2}=1.
\end{equation}
It proves convenient to carry out one more rearrangement
of unknowns introducing the differences
\begin{equation} \label{eq57}
\Delta p_{j_0,l}=p_{j_0,l+1}-p_{j_0,l}.
\end{equation}
which obey the own `initial conditions'
\begin{gather} \label{eq58}
\Delta p_{j_0,-2}= 1, \\ \label{eq58a}
\Delta p_{j_0,-1}= \lambda Z^{-1}_{j_0},
\end{gather}
and analogous recurrence relations
\begin{equation} \label{eq59}
\Delta p_{j_0,l+1}= \Delta p_{j_0+1,l}+
\lambda \Delta p_{j_0+1,l} Z_{j_0}^{-1}
+\mu^2 \Delta p_{j_0+2,l-1}Z_{j_0+1-\tilde\epsilon}^{-1}.
\end{equation}
Now, summing up the subset of all the equations (\ref{eq59})
distinguished by
the common sum of indices
at the left and taking into account (\ref{eq58a}),
all but two `free' $\Delta p$-terms cancel out and
one obtains
\begin{equation} \label{eq60}
\begin{aligned}
\Delta p_{j_0,l+1}
&= \lambda Z^{-1}_{j_0+2+l} +\lambda \sum_{m=0}^{l+1}
\Delta p_{j_0+1+m,l-m} Z_{j_0+m}^{-1} \\
&\quad +\mu^2 \sum_{m=0}^{l+1}
\Delta p_{j_0+2+m,l-1-m}Z_{j_0+m+1-\tilde\epsilon}^{-1}.
\end{aligned}
\end{equation}
In the sums, the second index of $\Delta p_{*,*}$ is
everywhere less than the same index at the left
that allows to apply the method of mathematical induction.
For the `starting' values $-1,0$ of the second index one has
\begin{gather*}
Z^{}_{j_0}\Delta p_{j_0,-1} = \lambda, \\
Z_{j_0+1} \Delta p_{j_0,0} =
\lambda(1+\lambda Z^{-1}_{j_0})
+\mu^2 Z_{j_0+1} Z^{-1}_{j_0+2-\tilde\epsilon}.
\end{gather*}
Therefore, for $l=-1,0$ there exist the finite limits
$\lim_{j_0\to \infty}|Z^{}_{j_0+l+1}\Delta p_{j_0,l}|$.
As a consequence, for these $l$'s,
\begin{equation} \label{eq61}
|\Delta p_{j_0,l}|< \tilde N |Z_{j_0+l+1}|^{-1}
\end{equation}
for appropriate constant $\tilde N$ which is convenient to choose $>1$.
Let us consider this fact as the starting point
of mathematical induction and assume that for some integer $l_0\ge 0$
and any integer $l$ from the interval $[-1,l_0]$ (\ref{eq61})
holds true. We may apply it for the estimating from above of the quantity
$|\Delta p_{j_0,l_0+1}|$. This can be realized making use of
the `decomposition' (\ref{eq60}) and the following
elementary inequalities
\begin{gather}\label{eq62}
\sum_{m=j_0}^{L+j_0+1}
|Z_{m}|^{-1}<{1+|n+1|^{-1} \over j_0-|n+1|/2+(\epsilon-1)/2 },
\\
\sum_{m=j_0+1-\tilde\epsilon}^{L+j_0+2-\tilde\epsilon}
|Z_{m}|^{-1}< {1+|n+1|^{-1}\over j_0-|n+1|/2+(\epsilon-1)/2 },
\end{gather}
where $L>0$ (and $n\not=-1$) .
These imply the inequalities
\begin{align*}
| \Delta p_{j_0,l_0+1}|
&\le |\lambda| |Z_{j_0+l_0+2}|^{-1} \\
&\quad +\tilde N |Z_{j_0+l_0+1}|^{-1}
\Big( |\lambda|
\sum_{m=0}^{l_0+1} |Z_{j_0+m}|^{-1}
+ |\mu^2| \sum_{m=0}^{l_0+1}
|Z_{j_0+m+1-\tilde\epsilon}|^{-1} \Big) \\
&< |Z_{j_0+l_0+2}|^{-1}
\Big(1+\tilde N{|Z_{j_0+l_0+2}|\over |Z_{j_0+l_0+1}|}
{(|\lambda|+|\mu^2|)(1+|n+1|^{-1})\over (j_0-|n+1|/2+(\epsilon-1)/2)}
\Big).
\end{align*}
Since we assumed $\tilde N>1$, there exists the bound of
lower index values such that for any $j_0$ exceeding
it the factor in brackets is less than $\tilde N $, and then
the above inequalities imply
$|\Delta p_{j_0,l_0+1}|<\tilde N |Z_{j_0+l_0+2}|^{-1}$.
The induction step is accomplished.
(\ref{eq61}) is therefore established
for sufficiently large $j_0$ and arbitrary
$l\ge0$. Increasing $\tilde N$ if necessary, (\ref{eq61})
proves valid for arbitrary $j_0$.
In view of this property, one sees that the
sum $\sum_{l=0}^\infty\Delta p_{j_0,l}$
has the majorant const$\times\sum_l|Z_{j_0+l+1}|^{-1}$
and thus converges itself.
The sequence of its partial sums
$\sum_{m=0}^l\Delta p_{j_0,m}=p_{j_0,l+1}-p_{j_0,-1}=p_{j_0,l+1}$
also converges as $l\to \infty$.
Moreover, in view of (\ref{eq55}), (\ref{eq57}), (\ref{eq62}),
(\ref{eq61}) one has the $l$-uniform bound
\begin{equation}\label{eq62a}
|p_{j_0,l+1}|
<\tilde N\sum_{m=0}^l |Z_{j_0+m+2}|^{-1}
<{\tilde N(1+|n+1|^{-1}) \over j_0+1+(\epsilon-1)/2}.
\end{equation}
It follows from \eqref{eq53}, \eqref{eq:53a}
\begin{equation} \label{eq53b}
\alpha_j= \alpha_{j_0+l}=(1+p_{j_0,l})\alpha_{j_0-1}
+ (1+p_{j_0+1,l-1} )\mu^2 Z_{j_0-\tilde\epsilon}^{-1}\alpha_{j_0-2}
\end{equation}
and the convergence of $\alpha$-sequence follows from the
convergence of $p_{*,l}$. % as $l\to \infty$.
Moreover, one has
\begin{equation} \label{eq53c}
\lim\alpha_j-\alpha_{j_0-1}=
\lim_l p_{j_0,l}\cdot \alpha_{j_0-1}
+(1+\lim_l p_{j_0+1,l} )\mu^2 Z_{j_0-\tilde\epsilon}^{-1}\cdot\alpha_{j_0-2}.
\end{equation}
The factors in front of
the first and the second terms to the right scales as $j_0^{-1}$ and
$j_0^{-2}$, respectively.
We may therefore write down the inequality
\begin{equation*}
|\lim\alpha_j-\alpha_{j_0-1}|=N\max(|\alpha_{j_0-1}|,
|\alpha_{j_0-2}|)j_0^{-1},
\end{equation*}
where $N$ may depend on the parameters $n,\lambda,\mu,\kappa,\epsilon$
but not on the specific specimen of $\alpha$-sequence.
This obviously implies the inequality (\ref{eq26}).
It has also to be noted in conclusion that the case $n=-1$
formally falling off the above speculation
does not actually correspond to an exceptional situation.
Although inequalities (\ref{eq62})
formally fail, similar ones differing from \eqref{eq62}
in the values of `constant' ($j_0$-independent) terms alone
can be derived. The further reasoning holds true
and leads to the same conclusions.
\begin{thebibliography}{00}
\bibitem{R1} B. D. Josephson,
\emph{Possible new effects in superconducting tunelling},
{Phys.\ Lett.} {1}(1962) 251-253.
\bibitem{R2} B. D. Josephson,
Coupled superconductors, { Rev.\ Mod.\ Phys.} {36} (1964) 216-220.
\bibitem{R3} D. E. McCumber,
\emph{Effect of ac Impedance on dc Voltage-Current Characteristic
of Superconductor Weak-Link Junctions},
{J.\ Appl.\ Phys.} {39}(1968) 3113-3118.
\bibitem {R4} W. C. Stewart,
\emph{Current-voltage characteristics of Josephson junctions},
{Appl.\ Phys\ Lett.} {12},(1968) 277-280.
\bibitem{R5}
A.\ Barone and G.\ Paterno,
\emph{Physics and Applications of the Josephson Effect}, Willey, N.Y., 1982.
\bibitem{R6} A. H. Hamilton,
\emph{Josephson voltage standards},
{Rev.\ of Sci.\ Instr.}, {71} (2000) 3611-3623.
\bibitem{R7}
S. P. Benz and C. A. Hamilton,
\emph{Application of the Josephson effect to Voltage Metrology},
{Proc IEEE} {92} (2004) 1617-1629.
\bibitem{R8} J. Niemeyer,
\emph{Josephson Voltage Metrology},
in: {Abstract Booklet EUCAS 2005}, Vienna (2005) 83-84.
\bibitem {R9}
A. Kemooinen, J. Nissil\"a, K. Ojasalo, J. Hassel, A. Manninen,
P. Helist\"o and H. Sepp\"a,
\emph{AC Voltage Standard based on an Externally-Shunted SIS Josephson
Junction Array},
in: Abstract Booklet EUCAS 2005, Vienna (2005) 336.
\bibitem{R10}
R. L. Kautz;
\emph{Noise, chaos, and Josephson voltage standard},
Rep. Prog. Phys., {59} (1996) 935-992.
\bibitem{R11}
S. I. Tertychniy,
\emph{Long-term behavior of solutions to the equation
$\dot \phi + \sin\phi=f$
with periodic $f$ and the modeling of dynamics
of overdamped Josephson junctions}. Preprint
{Arxiv:math-ph/0512058} (2005).
\bibitem{R11a}
S. I. Tertychniy,
\emph{The modeling of a Josephson junction and Heun polynomials}.
Preprint {Arxiv:math-ph/0601064} (2006).
\bibitem{R12}
V. M. Buchstaber, O. V. Karpov and S. I. Tertychniy,
\emph{Quantum Josephson D/A converter driven by trains of short $2\pi$-pulses},
in: {Conference Digest CPEM 2002}, Ottawa, (2002) 502-503
(2002).
\bibitem {R13}
D. Schmidt, G. Wolf,
\emph{Double confluent Heun equation, in: Heun's differential equations},
Ronveaux (Ed.) Oxford Univ. Press, Oxford, N.Y., (1995), Part C.
\bibitem{R14}
S. Yu. Slavyanov and W. Lay,
\emph{Special functions: A Unified Theory Based on Singularities},
Nevskiy dialect, SPb, (2002), in Russian;
English edition: Oxford Univ. Press, Oxford, N.Y., (2000).
\end{thebibliography}
\end{document}