\documentclass[reqno]{amsart}
\usepackage{hyperref}
\usepackage{graphicx}

\AtBeginDocument{{\noindent\small
\emph{Electronic Journal of Differential Equations},
Vol. 2017 (2017), No. 40, pp. 1--10.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu}
\thanks{\copyright 2017 Texas State University.}
\vspace{8mm}}

\begin{document}
\title[\hfilneg EJDE-2017/40\hfil 
Periodic oscillations of the relativistic pendulum]
{Periodic oscillations of the relativistic pendulum with friction}

\author[Q. Liu, L. Huang, G. Jiang \hfil EJDE-2017/40\hfilneg]
{Qihuai Liu, L\"ukai Huang, Guirong Jiang}

\address{Qihuai Liu \newline
School of Mathematics and Computing Sciences,
Guangxi Colleges and Universities Key Laboratory
of Data Analysis and Computation,
Guilin University of Electronic Technology,
Guilin 541002, China}
\email{qhuailiu@gmail.com}

\address{L\"ukai Huang \newline
School of Mathematics and Computing Sciences,
Guilin University of Electronic Technology,
Guilin 541002, China}
\email{756060523@qq.com}

\address{Guirong Jiang \newline
School of Mathematics and Computing Sciences,
Guilin University of Electronic Technology,
Guilin 541002, China}
\email{grjiang9@163.com}

\thanks{Submitted April 6, 2016. Published February 6, 2017.}
\subjclass[2010]{34B15}
\keywords{Relativistic pendulum; Poincar\'e-Miranda theorem; averaging;
\hfill\break\indent periodic solutions}

\begin{abstract}
 We consider the existence and multiplicity of periodic oscillations for
 the forced pendulum model with relativistic effects by using the
 Poincar\'e-Miranda theorem. Some detailed information about the bound
 for the period of forcing term is obtained. To support our analytical work,
 we also consider a forced pendulum oscillator with the special force
 $\gamma_0\sin(\omega t)$ including a sufficiently small parameter.
 The result shows us that for all $\omega\in(0,+\infty)$,
 there exists a $2\pi/\omega$ periodic solution under our settings.
\end{abstract}

\maketitle
\numberwithin{equation}{section}
\newtheorem{theorem}{Theorem}[section]
\newtheorem{lemma}[theorem]{Lemma}
\allowdisplaybreaks


\section{Introduction}\label{s1}

In this article, we consider the existence and multiplicity of periodic
oscillations for the forced pendulum model with relativistic effects
\begin{equation}\label{1.1}
\bigg(\frac{x'}{\sqrt{1-\frac{x'^2}{c^2}}} \bigg)'+kx'+a\sin x=p(t),
\end{equation}
where $c>0$ is the speed of light in the vacuum, $k>0$ is a possible
viscous friction coefficient and $p$ is a continuous and $T$-periodic
forcing term with mean value zero
 \begin{equation*}
 \bar{p}=\frac{1}{T}\int_0^Tp(t)\,\mathrm{d}t=0.
 \end{equation*}

This equation has received much attention as a
prototype of equation with singular $\phi$-Laplacian
(see \cite{mawhin2010periodic} and \cite{Brezis2010,Liu2011,Mawhin2011}).
An essential difference between the relativistic and the newtonian
($c=+\infty$) case has been explained in \cite{torres2008periodic}.
In \cite{torres2008periodic}, Torres proved the following theorem.

\begin{theorem}\label{thm1}
 Let us assume that $2cT\leq1$. For any $a, k$ and any continuous
$T$-periodic function $p(t)$ with mean value zero, \eqref{1.1}
has at least one $T$-periodic solution.
\end{theorem}

The proof of the above theorem is an interesting application of the Schauder fixed
point theorem. The bound in Theorem \ref{thm1} was improved to
$2cT\leq4\sqrt{3}\approx6.9282$ in \cite{torres2011nondegeneracy}
for the general pendulum-type equation was considered,
see also \cite{Cid2013,Liu2014}.

Motivated by \cite{torres2008periodic}, in this paper we give
 detailed information on the bounds for $T$ which depend on the parameters
$a,k$ and the forcing term $p$.
Without loss of generality, we assume $a\geq0$; otherwise we only require
replacing $x$ with $x+\pi$. Let us define
\begin{equation*}
 \|p\|_\infty=\sup_{t\in[0,T]}|p(t)|,
\end{equation*}
and the constant
\begin{equation}\label{1.2}
 c_*=\frac{c\big(kTc_*+3k\pi+2(a+\|p\|_\infty)T\big)}{\sqrt{c^2+\big(kTc_*+3k\pi+2(a+\|p\|_\infty)T\big)^2}} < c.
\end{equation}
Our main result reads as follows.

\begin{theorem}\label{thm2}
 For any values $a,k$ and for any continuous and $T$-periodic function $p(t)$
with mean value zero satisfying $2c_*T\leq2\pi$, \eqref{1.1} has at least
two distinct $T$-periodic solutions.
\end{theorem}

The proof of Theorem \ref{thm2} is an elementary application of a variation of
the Poincar\'e-Miranda theorem (see \cite{liu2015}) which will be given in the
next section. We remark that the two distinct $T$-periodic solutions in
Theorem \ref{thm2} are indeed geometrically different periodic solutions, which
generalizes Theorem \ref{thm1}. Moreover, when $k$ or $\|p\|_\infty$ tend to infinity,
we show that $c_*\to c$ so that $2cT\leq2\pi$. This case does not improve
the previous bound.

To support our analytical work, based on the method of averaging,
we also consider the existence of periodic oscillations for a special
forced pendulum oscillator with a sufficiently small parameter $\varepsilon$,
\begin{equation}\label{1.3}
\bigg(\frac{x'}{\sqrt{1-\frac{x'^2}{c^2}}} \bigg)'+\varepsilon^2 kx'+a\sin x
=\varepsilon^3 \gamma_0\sin(\omega t),
\end{equation}
where $\omega^2=a+\varepsilon^2\beta_0$ with $\beta_0>0$.
We summarize our results as follows.

\begin{theorem}\label{thm3}
 For any $\gamma_0, k$, $\beta_0>0$ and $\omega>0$, \eqref{1.3} has at
least one $2\pi/\omega$-periodic solution when $\varepsilon$ is sufficiently
 small. Moreover, this periodic solution is stable for $k>0$ and is unstable
for $k<0$.
\end{theorem}

Noticed that $T=2\pi/\omega\to+\infty$ as $\omega\to0$.
Thus, in this case, \eqref{1.3} does not meet the hypotheses of Theorem \ref{thm1}.
From \eqref{1.2} we also see that $c_*\to0$ when $\varepsilon\to0$,
 satisfying the hypotheses of Theorem \ref{thm2}.

In Section 2, we introduce a variation of the Poincar\'e-Miranda theorem
in two-dimensional case which is used to prove Theorem \ref{thm2}.
 We prove Theorem \ref{thm2} in Section 3. In the last section, we prove
 Theorem \ref{thm3} using the method of averaging and perform some numerical
simulations.

\section{A variation of the Poincar\'e-Miranda theorem}\label{S2} 

We first introduce a variation of the Poincar\'e-Miranda theorem
(see \cite{rouche1980ordinary,liu2015} for instance) in two-dimensional case,
which goes back to Poincar\'e (1883) and has been used many times in the
study of boundary value problems and the existence of periodic solutions.
For an example see \cite{mawhin2013variations} and the references therein.

Consider the closed rectangle
\begin{equation*}
 D=\{(x,y)\in \mathbb{R}^2:\alpha_1\leq
x\leq\alpha_2, \beta_1\leq y\leq\beta_2\},
\end{equation*}
 where $\alpha_i, \beta_i$ $(i=1,2)$ are constants such that
$\alpha_1<\alpha_2$, $\beta_1<\beta_2$. The
boundary of the rectangle consists of four faces as follows:
\begin{gather*}
V_{-}^1=\{(x,y)\in
\mathbb{R}^2:x=\alpha_1,\beta_1\leq y\leq\beta_2\},\\
V_{+}^1=\{(x,y)\in
\mathbb{R}^2:x=\alpha_2,\beta_1\leq y\leq\beta_2\},\\
V_{-}^2=\{(x,y)\in
\mathbb{R}^2:y=\beta_1,\alpha_1\leq x\leq\alpha_2\},\\
V_{+}^2=\{(x,y)\in
\mathbb{R}^2:y=\beta_2,\alpha_1\leq x\leq\alpha_2\}.
\end{gather*}

We say that a continuous map
$F=(F_1,F_2):
D\to \mathbb{R}^2$ satisfies the \emph{bend-twist
condition} on $D$ provided that
\[
F_1(V_-^1)F_1(V_+^1)\leq0,\quad F_2(V_-^2) F_2(V_+^2)\leq0
\]
or
\[
F_2(V_-^1)F_2(V_+^1)\leq0,\quad F_1(V_-^2) F_1(V_+^2)\leq0,
\]
where $F_j(V_-^i) F_j(V_+^i)\leq0$ means that $F_j(V_-^i)\leq0$ and
$F_j(V_+^i)\geq0$, or $F_j(V_-^i)\geq0$ and $F_j(V_+^i)\leq0$;
$F_j(V_\pm^i)<0$ (resp. $F_j(V_\pm^i)>0$) means that $F_j(x,y)\leq0$
(resp. $F_j(x,y)\geq0$) for all $(x,y)\in V_\pm^i$ and there exists
at least $(x_0,y_0)\in V_\pm^i$
such that $F_j(x_0,y_0)<0$ (resp. $F_j(x_0,y_0)>0$); and $F_j(V_\pm^i)=0$
means that $F_j(x,y)=0$ for all $(x,y)\in V_\pm^i$, $i,j=1,2$.

\begin{theorem}[{See \cite[Theorem 2.1]{liu2015}}]\label{thm2.1}
Assume the continuous map $F:D\to \mathbb{R}^2$ satisfies the bend-twist
condition, then there exists at least one point
$(x_0,y_0)\in D$ such that $F(x_0,y_0)=0$.
\end{theorem}


\section{Proof of Theorem \ref{thm2}}\label{s2}

Equation \eqref{1.1} is equivalent to the plane system
\begin{gather} \label{3.1}
x'=\frac{c(y-kx)}{\sqrt{c^2+(y-kx)^2}}, \\
 \label{3.2}
y'=-a\sin x+p(t).
\end{gather}
Let $\alpha,\beta$ be positive constants and
\begin{equation*}
 \|p\|_\infty=\sup_{t\in[0,T]}|p(t)|,\quad \lambda=kT,\quad
 \mu=\beta+k\alpha+(a+\|p\|_\infty)T.
\end{equation*}
Define $\phi:(-\infty,+\infty)\to(-c,c)$ by
\begin{equation*}
 \phi(u)=\frac{cu}{\sqrt{c^2+u^2}}.
\end{equation*}
It is easy to verify that $\phi$ is an increasing homeomorphism such that
$\phi(-u)=-\phi(u)$.

\begin{lemma}\label{lem1}
 Assume that $p(t)$ is a continuous $T$-periodic function. Then for any values
$a,k$ and any initial value
$(x_0,y_0)\in\big\{(x,y)\big|\,|x|\leq \alpha, |y|\leq \beta \text{ and }
 \alpha,\beta>0\big\}$, the solution $(x(t;x_0,y_0), y(t;x_0,y_0))$ of
\eqref{3.1}-\eqref{3.2} with the initial value $(x_0,y_0)$ satisfies
 \begin{equation*}
 |x'(t)|\leq c_*(\alpha,\beta)<c, \quad \forall t\in[0,T],
 \end{equation*}
 where $c_*$, depending on $\alpha,\beta$, is a solution of $u=\phi(\lambda u+\mu)$.
\end{lemma}

\begin{proof}
First we note that $|x'(t)|\leq c_1:=c$ for all $t\in[0,T]$.
Hence $|x(t)|\leq \alpha+c_1T$ for all $t\in[0,T]$, and by \eqref{3.2} we see
that $|y(t)|\leq \beta+(a+\|p\|_\infty)T$ for all $t\in[0,T]$. Therefore,
\[
 |y(t)-kx(t)|<k\alpha+\beta+(a+c_1k+\|p\|_\infty)T=\lambda c_1+\mu, \quad \forall
 t\in[0,T].
\]
Let $c_2:=\phi(\lambda c_1+\mu)$. By \eqref{3.2}, we have
$|x'(t)|\leq c_2$ for all $t\in[0,T]$. Obviously, $c_2<c_1$.
Repeating this argument we have a sequence $\{c_n\}_{n\in \mathbb{N}}$
defined by $c_n=\phi(\lambda c_{n-1}+\mu)$.

Since $\phi$ is an increasing homeomorphism and $c_2<c_1$, we know that
$c_3=\phi(\lambda c_{2}+\mu)<\phi(\lambda c_{1}+\mu)=c_2$, \dots,
$c_n=\phi(\lambda c_{n-1}+\mu)<\phi(\lambda c_{n-2}+\mu)=c_{n-1}, \dots$.
That is, $\{c_n\}_{n\in \mathbb{N}}$ is a decreasing sequence.
On the other hand, $|c_n|=|\phi(\lambda c_{n-1}+\mu)|<c$. Hence
$\{c_n\}_{n\in \mathbb{N}}$ converges to some value $c_*< \infty$.
Since $\phi$ is continuous, by passing
 we have $c_*=\phi(\lambda c_*+\mu)$, that is
\begin{align*}
 c_*&=\phi(\lambda c_*+\mu)\\
 &=\frac{c\big(kTc_*+\beta+k\alpha+(a+\|p\|_\infty)T\big)}
 {\sqrt{c^2+\big(kTc_*+\beta+k\alpha+(a+\|p\|_\infty)T\big)^2}}.
\end{align*}
\end{proof}


\begin{proof}[Proof of Theorem \ref{thm2}]
Let $\gamma=\frac{3}{2}k\pi+(a+\|p\|_\infty)T$ and $c_*=c_*(3\pi/{2},\gamma)$.
Let us construct a rectangle as follows
\[
D_1=\{(x,y)\in
\mathbb{R}^2:-\frac{\pi}{2}\leq x\leq \frac{\pi}{2}, -\gamma\leq y\leq \gamma\}.
\]
The boundary is of $D_1$ is given by
\begin{gather*}
 V_-^1=\{(x,y)\in D_1: x=-\frac{\pi}{2}\},\\
 V_+^1=\{(x,y)\in D_1: x=\frac{\pi}{2}\}\\
 V_-^2=\{(x,y)\in D_1: y=-\gamma\},\\
 V_+^2=\{(x,y)\in D_1: y=\gamma\}.
\end{gather*}
Let $(x(t;x_0,y_0),y(t;x_0,y_0))$ be the solution of \eqref{3.1} and \eqref{3.2}
with the initial value $(x_0,y_0)\in D_1$. Define the continuous mapping
$F:\mathbb{R}^2\to\mathbb{R}^2$ by
\begin{equation*}
F(x_0,y_0)=\begin{pmatrix}
 F_1(x_0,y_0) \\
 F_2(x_0,y_0)
 \end{pmatrix}
=(P-\mathrm{id})(x_0,y_0),
\end{equation*}
where $P$ denotes the Poincar\'e mapping associated with
 system \eqref{3.1}-\eqref{3.2}.

(i) When $(x_0,y_0)\in V_-^1$, using Lemma \ref{lem1} we know that
\[
 |x'(t)|<c_*,\quad \forall t\in[0,T].
\]
Then it follows that
\[
 -\frac{\pi}{2}-c_*t\leq x(t)\leq-\frac{\pi}{2}+c_*t, \quad \forall t\in[0,T].
\]
When $t\in [0,\pi/c_*]$, we know that
\[
 -\frac{\pi}{2}\leq-\frac{c_*t}{2}\leq\frac{x(t)+\frac{\pi}{2}}{2}
 \leq\frac{c_*t}{2}\leq\frac{\pi}{2}, \quad \forall t\in[0,T].
\]
Then it follows that, for any $t\in [0,\pi/c_*]$,
\[
 \Big(\cos\frac{c_*t}{2}\Big)^2\leq\Big(\cos \frac{x(t)+\frac{\pi}{2}}{2}\Big)^2\leq1.
\]
Therefore,
\begin{align*}
 \int_0^{\pi/c_*}[-\sin x(t)]\,\mathrm{d}t
&=\int_0^{\pi/c_*}\cos\big[x(t)+\frac{\pi}{2}\big]\,\mathrm{d}t\\
&=\int_0^{\pi/c_*}\big[2\big(\cos \frac{x(t)+\frac{\pi}{2}}{2}\big)^2-1\big]
 \,\mathrm{d}t \\
&\geq\int_0^{\pi/c_*}\big[2\big(\cos\frac{c_*t}{2}\big)^2-1\big]\,\mathrm{d}t=0.
\end{align*}
Therefore, we have
\begin{equation} \label{3.3}
\begin{aligned}
y(T)-y(0)&=\int_0^T[-a\sin x(t)]\,\mathrm{d}t\\
 &=\int_0^{\pi/c_*}[-a\sin x(t)]\,\mathrm{d}t
+\int_{\pi/c_*}^T[-a\sin x(t)]\,\mathrm{d}t \\
&\geq-a\big(T-\frac{\pi}{c_*}\big)\geq0.
\end{aligned}
\end{equation}
The above inequality is obtained by the hypothesis $2c_*T\leq2\pi$.

 When $(x_0,y_0)\in V_+^1$, using the same arguments, we have
\[
 \frac{\pi}{2}-c_*t<x(t)<\frac{\pi}{2}+c_*t, \quad \forall t\in[0,T].
\]
When $t\in [0,\pi/c_*]$, we know that
\begin{align*}
 \frac{\pi}{2}\leq\pi-\frac{c_*t}{2}\leq\frac{x(t)+\frac{3\pi}{2}}{2}
 \leq\pi+\frac{c_*t}{2}\leq\frac{3\pi}{2}, \quad \forall t\in[0,T].
\end{align*}
Similarly, for any $t\in [0,\pi/c_*]$, we have
\begin{align*}
 \big(\cos\frac{c_*t}{2}\big)^2\leq\big(\cos \frac{x(t)+\frac{3\pi}{2}}{2}
\big)^2\leq1.
\end{align*}
Therefore,
\begin{align*}
 \int_0^{\pi/c_*}[-\sin x(t)]\,\mathrm{d}t
&=-\int_0^{\pi/c_*}\cos\big[x(t)+\frac{3\pi}{2}\big]\,\mathrm{d}t\\
&=-\int_0^{\pi/c_*}\big[2\big(\cos \frac{x(t)+\frac{3\pi}{2}}{2}\big)^2-1\big]
 \,\mathrm{d}t \\
&\leq-\int_0^{\pi/c_*}\big[2\big(\cos\frac{c_*t}{2}\big)^2-1\big]\,\mathrm{d}t=0.
\end{align*}
Therefore,
\begin{equation} \label{3.4}
\begin{aligned}
y(T)-y(0)
&=\int_0^T[-a\sin x(t)]\,\mathrm{d}t\\
&=\int_0^{\pi/c_*}[-a\sin x(t)]\,\mathrm{d}t
 +\int_{\pi/c_*}^T[-a\sin x(t)]\,\mathrm{d}t \\
&\leq a\big(T-\frac{\pi}{c_*}\big)\leq0.
\end{aligned}
\end{equation}
The last inequality is obtained by the hypothesis $2c_*T\leq2\pi$.
From \eqref{3.3} and \eqref{3.4}, we have that $F_2(V_-^1) F_2(V_+^1)\leq0$.

(ii) When $(x_0,y_0)\in V_-^2$, using the inequality $|x'(t)|\leq c_*$ we know that,
for all $t\in[0,T]$,
\[
 -\frac{3\pi}{2}\leq-c_*T-\frac{\pi}{2}\leq x(t)\leq\frac{\pi}{2}+c_*T\leq\frac{3\pi}{2}.
\]
 Then using \eqref{3.2} we know that
\begin{align*}
 y(t)-kx(t)
&=y_0+\int_0^t(-a\sin x(s)+p(s))\,\mathrm{d}s-kx(t)\\
&\leq -\gamma+(a+\|p\|_\infty)T+\frac{3}{2}k\pi=0, \quad t\in[0,T].
\end{align*}
Since $\phi$ is a continuous
homeomorphism such that $\phi(0)=0$, we have $\phi(u)u\geq0$.
Then it follows that $x'(t)=\phi(y(t)-kx(t))\leq0$ for all $t\in[0,T]$,
 which yields
\begin{equation}\label{3.7}
 x(T)-x(0)=\int_0^Tx'(\tau)\,\mathrm{d}\tau\leq0.
\end{equation}

When $(x_0,y_0)\in V_+^2$, we also know that for all $t\in[0,T]$,
$|x(t)|\leq\frac{3\pi}{2}$, and
\begin{align*}
 y(t)-kx(t)
&=y_0+\int_0^t(-a\sin x(s)+p(s))\,\mathrm{d}s-kx(t)\\
&\geq \gamma-(a+\|p\|_\infty)T-\frac{3}{2}k\pi=0, \quad t\in[0,T].
\end{align*}
With the same arguments we have $x(T)-x(0)\geq0$. Therefore,
\begin{equation*}
 F_1(V_-^2)F_1(V_+^2)\leq0.
\end{equation*}

We have verified that $F$ satisfies the bend-twist condition on $D_1$.
By Theorem \ref{thm1}, there exists at least one point
$(x_1,y_1)\in D_1$ such that $F(x_1,y_1)=0$, which is corresponding
to a fixed point of the Poincar\'e mapping.

Similarly, we can construct the rectangle
\[
D_2=\{(x,y)\in \mathbb{R}^2:\frac{\pi}{2}\leq x\leq \frac{3\pi}{2}, -\gamma\leq
y\leq\gamma\}.
\]
With the same arguments, we can verify that $F$ satisfies the bend-twist condition
on $D_2$ and obtain another fixed point of the Pincar\'e mapping in $D_2$.

Let $V=D_1\cap D_2$. To prove that such two fixed points of $F$ are distinct,
it is sufficient to prove that there is no $T$-periodic solution with the
initial value on $V$. Assume that $(x(t;\frac{\pi}{2},y_0),y(t;\frac{\pi}{2},y_0))$
is a $T$-periodic solution of \eqref{3.1} and \eqref{3.2}. Then we know that
$\{x(t;\frac{\pi}{2},y_0)\big | t\in[0,T]\}$ is contained in $[0, \pi]$,
since the maximum of the derivative of $x(t)$ is $c_*$ and $c_*T\leq \pi$.
Then we have
\[
 y(T)-y(0)=\int_0^T[-a\sin x(t)]\,\mathrm{d}t
\leq \int_0^{\pi/(3c_*)}[-a\sin x(t)]\,\mathrm{d}t <0.
\]
Therefore, we obtain two distinct fixed points, which are corresponding
to two distinct $T$-periodic solutions of equation \eqref{1.1}.
\end{proof}

\section{Numerical examples and proof of Theorem \ref{thm3}}\label{s4}

 First we prove Theorem \ref{thm3} by using the method of averaging. Recall that
\begin{equation}\label{4.1}
\bigg(\frac{x'}{\sqrt{1-\frac{x'^2}{c^2}}} \bigg)'+\varepsilon^2 kx'+a\sin x
=\varepsilon^3 \gamma_0\sin(\omega t),
\end{equation}
where $\omega^2=a+\varepsilon^2\beta_0$ and $\varepsilon$ is a small parameter.

Equation \eqref{4.1} is equivalent to the plane system
\begin{equation} \label{4.2}
\begin{gathered}
 x'=\frac{c(y-\varepsilon^2 kx)}{\sqrt{c^2+(y-\varepsilon^2 kx)^2}}, \\
 y'=-a\sin x+\varepsilon^3\gamma_0\sin(\omega t).
 \end{gathered}
\end{equation}
Let $x=\varepsilon u$, $y=\varepsilon v$ and $\epsilon=\varepsilon^2$.
We expand system \eqref{4.2} into the form of power series by
\begin{equation} \label{4.3}
\begin{gathered}
 u'=v+\epsilon f_1(u,v,t)=v -\Big(k u+\frac{1}{2}c^{-2}v^3\Big) \epsilon
 +O(\epsilon^2 ), \\
 y'=-\omega^2u+\epsilon f_2(u,v,t)=-\omega^2 u+\beta_0u\epsilon
 +\frac{1}{6}\omega^2 u^3\epsilon+\epsilon \gamma_0\sin\omega t+O(\epsilon^2).
 \end{gathered}
\end{equation}
Using the van der Pol transformation
\[
 u=q\sin\omega t+p\cos\omega t,\quad v=\omega(q\cos\omega t-p\sin\omega t),
\]
we obtain
\begin{equation} \label{4.4}
\begin{gathered}
 q'=\epsilon\Big(f_1(u,v,t)\sin\omega t+\frac{\cos\omega t}{\omega}f_2(u,v,t)\Big)
 +O(\epsilon^2), \\[0.2em]
 p'=\epsilon\Big(f_1(u,v,t)\cos\omega t-\frac{\sin\omega t}{\omega}f_2(u,v,t)\Big)
 +O(\epsilon^2).
\end{gathered}
\end{equation}
Then it follows that
\begin{align*}
 q'&=\epsilon F_1(q,p,t,\epsilon) \\
 &=\frac{1}{48c^2\omega}\Big(\omega\Big(9p(p^2+q^2)\omega^3
 +3c^2\big(-8k q+p(p^2+q^2)\omega\big) \\
&\quad + 4\big(-3p^3\omega^3+c^2(6kq+p^3\omega)\big)\cos(2\omega t)
 +p (p^2-3 q^2) \omega (c^2+3 \omega ^2) \cos(4 \omega t) \\
&\quad +2 \big(-3 q (3 p^2+q^2) \omega ^3
 +c^2 (-12 k p+3 p^2 q \omega +q^3 \omega )\big) \sin(2 \omega t) \\
&\quad -q (-3 p^2+q^2) \omega (c^2+3 \omega ^2) \sin(4 \omega t)\Big) \\
&\quad + 24 c^2 \big((p+p \cos(2 \omega t)+q \sin(2 \omega t)) \beta _0
 +\sin(2 \omega t)\gamma _0\big)\Big))+O(\epsilon^2),
\end{align*}
and
\begin{align*}
 p'&=\epsilon F_2(q,p,t,\epsilon) \\
&=\frac{1}{48 c^2 \omega }\Big(\omega \Big(-9 q (p^2+q^2) \omega ^3-3 c^2
\big(8 k p+q (p^2+q^2) \omega \big)\\
&\quad + 4 \big(-3 q^3 \omega ^3+c^2 (-6 k p+q^3 \omega )\big)
 \cos(2 \omega t)-q (-3 p^2+q^2) \omega (c^2+3 \omega ^2)\cos(4 \omega t) \\
&\quad -2\Big(-3 p (p^2+3 q^2) \omega ^3+c^2
\big(p^3 \omega +3 q (4 k+p q \omega )\big)\Big) \sin(2 \omega t) \\
&\quad - p(p^2-3 q^2) \omega (c^2+3 \omega ^2) \sin(4 \omega t )\Big)\\
&\quad - 48 c^2 \sin(\omega t) \big(p \cos( \omega t) \beta _0+\sin(\omega t)
 (q \beta _0+\gamma _0)\big)\Big)+O(\epsilon^2).
\end{align*}

It is not difficult to obtain the averaging system
\begin{equation} \label{4.5}
\begin{aligned}
 \bar{q}'
&=\epsilon G_1(\bar{q},\bar{p}) \\
&=\epsilon\frac{\omega}{2\pi} \int_0^{\frac{2\pi}{\omega}}F_1(\bar{q},\bar{p},t,0)\,\mathrm{d}t\\
&=\epsilon\frac{1}{16} \omega \Big(p (p^2+q^2)-\frac{8 k q}{\omega }
+\frac{3 p (p^2+q^2) \omega ^2}{c^2}+\frac{8 p \beta _0}{\omega ^2}\Big),
\\
\bar{p}'
&=\epsilon G_2(\bar{q},\bar{p}) \\
&=\epsilon\frac{\omega}{2\pi}
 \int_0^{\frac{2\pi}{\omega}}F_2(\bar{q},\bar{p},t,0) \,\mathrm{d}t\\
& =-\epsilon\frac{\omega \big(3 q (p^2+q^2) \omega ^3+c^2
\big(8 k p+q (p^2+q^2) \omega \big)\big)
 +8 c^2 (q \beta _0+\gamma _0)}{16 c^2 \omega }.
 \end{aligned}
\end{equation}
 The other equilibrium points of system \eqref{4.5} correspond to the solutions of
 \begin{gather*}
 G_1(\bar{q},\bar{p})=\frac{1}{16c^2}
\Big(\omega c^2pr^2-8kc^2q+3\omega^3pr^2+\frac{8c^2\beta_0}{\omega}p\Big)=0, \\
 G_2(\bar{q},\bar{p})=-\frac{1}{16c^2}
\Big(\omega c^2qr^2+8kc^2p+3\omega^3qr^2+\frac{8c^2\beta_0}{\omega}q
 +\frac{8c^2\gamma_0}{\omega}\Big)=0,
 \end{gather*}
where $r^2=q^2+p^2$, which is equivalent to
\begin{equation*}
 q=-\Big(\frac{\omega^2(c^2+3\omega^2)}{8c^2\gamma_0}r^2+\frac{\beta_0}{\gamma_0}
\Big)r^2,
 \quad p=-\frac{k\omega}{\gamma_0}r^2.
\end{equation*}
Thus, $r$ satisfies
\begin{equation*}
 \Phi(r)=\Big(\frac{\omega^2(c^2+3\omega^2)}{8c^2\gamma_0}r^2
+\frac{\beta_0}{\gamma_0}\Big)^2r^2
 +\big(\frac{k\omega}{\gamma_0}\big)r^2-1=0.
\end{equation*}
Since $\Phi(0)=-1<0$ and $\Phi(+\infty)=+\infty$, by the intermediate value
theorem we know that there is a $r_*\in(0,+\infty)$ such that $\Phi(r_*)=0$.
The Jacobi determinant of $(G_1,G_2)$ at $r_*$ is
 \begin{align*}
 J&=\frac{\partial(G_1,G_2)}{\partial(q,p)}\Big|_{(q,p)=(q(r_*),p(r_*))}\\
 &=\frac{k^2}{4}+\frac{3 p^4 \omega ^2}{256}+\frac{3}{128} p^2 q^2 \omega ^2
 +\frac{3 q^4 \omega ^2}{256}+\frac{9 p^4 \omega ^4}{128 c^2} \\
&\quad +\frac{9 p^2 q^2 \omega ^4}{64 c^2}+\frac{9 q^4 \omega ^4}{128 c^2}
 +\frac{27 p^4 \omega ^6}{256 c^4}+\frac{27 p^2 q^2 \omega ^6}{128 c^4} \\
&\quad +\frac{27 q^4 \omega ^6}{256 c^4}+\frac{p^2 \beta _0}{8}
 +\frac{q^2 \beta _0}{8}+\frac{3 p^2 \omega ^2 \beta _0}{8 c^2}
 +\frac{3 q^2 \omega ^2 \beta _0}{8 c^2}+\frac{\beta _0^2}{4 \omega ^2}
>0,
 \end{align*}
for $\beta _0>0$. The Jacobi Matrix has a pair of conjugate imaginary eigenvalues
$\lambda_{1,2}$ which satisfy that $\mathfrak{Re}(\lambda_{1,2})=-{k}/{2}$.
With the classical arguments of averaging theory, we know that system \eqref{4.4}
has a $2\pi/\omega$-periodic solution
$(q_\epsilon, p_\epsilon)$ such that $(q_\epsilon, p_\epsilon)\to(q(r_*),p(r_*))$
as $\epsilon\to0$, which yields
a $2\pi/\omega$-periodic solution $x(t)$ of \eqref{4.1}.
The periodic solution $x(t)$ is stable for $k>0$, while it is unstable for $k<0$.
Now we have finished the proof of Theorem \ref{thm3}.


\begin{figure}[htb] 
\begin{center}
 \includegraphics[width=0.4\textwidth]{fig1a}
 \includegraphics[width=0.4\textwidth]{fig1b} \\
(a) \hfil (b)
\end{center}
\caption{Profiles of the $2\pi/\omega$-periodic solution $(u,v)$ of
\eqref{4.3} with $\omega =10^{-3}$, $k=1$, $\beta _0=2$,
$\gamma _0=1$, $c=100$, $\epsilon=10^{-5}$.}
\label{fig1}
\end{figure}

To support our analytical work, we numerically simulate the $2\pi/\omega$-periodic 
solution of \eqref{4.3} for $\omega =0.001,k=1,\beta _0=2,\gamma _0=1,c=100$. 
We obtain that the rest point of \eqref{4.4} is 
$(q(r_*),p(r_*))=(-0.50000,-0.00025)$, the Jacobi determinant of 
$(G_1,G_2)$ at $(q(r_*),p(r_*))$ is $1.0\times10^6$ and the corresponding 
eigenvalues are $\lambda_{1,2}=-0.5\pm1000\mathfrak{i}$. 
We depict the corresponding stable $2\pi/\omega$-periodic solution 
of \eqref{4.3} in figure \ref{fig1}.


\subsection*{Acknowledgements}
This research  is supported by the National Natural Science Foundation (No 11301106), 
the Guangxi Natural Science Foundation (Nos. 2014GXNSFBA118017, 2015GXNSFGA139004) 
and the Guangxi Experiment Center of Information Science (Grant No. YB1410).


\begin{thebibliography}{00}

\bibitem{Brezis2010}    H.~Brezis,  J.~Mawhin;
\emph{Periodic solutions of the forced relativistic pendulum},
  Differential Integral Equations, \textbf{23}(9) (2010), 801-810.

\bibitem{Cid2013}  J.~A.~Cid,  P.~J.~Torres;
\emph{On the existence and stability of periodic solutions for pendulum-like
 equations with friction and ¦Õ-Laplacian},
  Discrete Contin. Dyn. Syst., \textbf{33}(9) (2013), 141-152.

\bibitem{Liu2011}  Q.~Liu,  D.~Qian,  B.~Chu;
\emph{Nonlinear systems with singular vector $\phi$-Laplacian under the 
Hartman-type condition},   Nonlinear Anal. TMA, \textbf{74}(8) (2011), 2880-2886.

\bibitem{Liu2014}  Q. Liu,  C.~Wang,  Z.~Wang;
\emph{On Littlewood's boundedness problem for relativistic oscillators with
 anharmonic potentials},  J. Differential Equations, \textbf{257}(12) (2014), 
4542-4571.

\bibitem{mawhin2010periodic}  J.~Mawhin;
\emph{Periodic solutions of the forced pendulum: classical vs relativistic},
  Matematiche (Catania), \textbf{65}(2) (2010), 97--107.

\bibitem{Mawhin2011}  J.~Mawhin;
\emph{Radial solutions of Neumann problem for periodic perturbations of the
 mean extrinsic curvature operator},
  Milan J. Math., \textbf{79}(1) (2011), 95-112.

\bibitem{mawhin2013variations}  J.~Mawhin;
\emph{Variations on Poincar\'e-Miranda's theorem},
  Adv. Nonlinear Stud., \textbf{13} (2013), 209--217.

\bibitem{rouche1980ordinary}  N.~Rouche.  J. Mawhin;
\emph{Ordinary Differential Equations: Stability and Periodic Solutions},
Pitman  Advanced Publishing Program, Boston, 1980.

\bibitem{torres2008periodic}  P.~Torres;
\emph{Periodic oscillations of the relativistic pendulum with friction},
  Phys. Lett. A, \textbf{372}(42) (2008), 6386--6387.

\bibitem{torres2011nondegeneracy}  P.~Torres;
\emph{Nondegeneracy of the periodically forced li{\'e}nard differential equation with
 $\phi$-laplacian},
  Communications in Contemporary Mathematics, \textbf{13}(02) (2011), 283--292.

\bibitem{liu2015}  X.~Wang,  D.~Qian,  Q.~Liu;
\emph{Existence and multiplicity results for some nonlinear problems withsingular
 phi-laplacian via a geometric approach},
 Bound. Value Probl., \textbf{2016}(47) (2016), 1--27.
 
\end{thebibliography}

\end{document}