Stabilization of Steady Rotations of a Spherical Robot on a Vibrating Base Using Feedback

Abstract

This paper treats the problem of a spherical robot with an axisymmetric pendulum drive rolling without slipping on a vibrating plane. The main purpose of the paper is to investigate the stabilization of the upper vertical rotations of the pendulum using feedback (additional control action). For the chosen type of feedback, regions of asymptotic stability of the upper vertical rotations of the pendulum are constructed and possible bifurcations are analyzed. Special attention is also given to the question of the stability of periodic solutions arising as the vertical rotations lose stability.

1 INTRODUCTION

In this paper, we consider the dynamics of a spherical pendulum-type robot rolling on a plane performing vertical periodic oscillations. We assume that there is no slipping at the point of contact of the shell with the plane. The main goal of the paper is to investigate the stabilization of vertical rotations of the pendulum using feedback.

Spherical robots are interesting for developing new approaches to the control of complex nonlinear systems. Detailed reviews of various models of spherical robots moving both by displacing the center of mass and by using other principles of motion and reviews of their technical implementations can be found in the recent papers [1, 2, 3, 4].

This paper investigates a spherical robot moving by displacing the center of mass, more precisely, a robot of pendulum type which is a spherical shell with an axisymmetric pendulum fastened at its center. The dynamics of a spherical robot of this type rolling on fixed surfaces is fairly well understood [5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15].

One of the research directions that are important from a practical point of view is analysis of the motion of a spherical robot on moving surfaces. To date, the dynamics of spherical robots on a vibrating plane remains largely unexplored. There exist several studies examining the free dynamics of spherical bodies on a vibrating plane, for example, the dynamics of a homogeneous sphere [16, 17, 18, 19], a Chaplygin sphere [20], and a spherical top with axisymmetric mass distribution on a smooth [21, 22] and rough [23] horizontal plane. In [24] an analysis is made of the controlled motion of a Chaplygin sphere with three internal gyrostats on a plane that performs horizontal oscillations.

The problem of stabilizing the equilibrium points of a dynamical system by means of vibrations (for example, a suspension point) is a classical problem and was addressed by many authors [25, 26, 27, 28, 29, 30, 31, 32, 33]. However, when it comes to spherical robots (and any robotic devices in general), of more current interest is the problem of stabilization of equilibrium points (or steady motions) by means of internal mechanisms in the presence of external vibrations.

For a spherical pendulum-type robot this problem is solved in [34]. There, the dynamics of a spherical pendulum-type robot rolling without slipping on a plane performing vertical periodic oscillations is considered. In particular, reference [34] examines the influence of the rotation of an internal pendulum on the stabilization of partial solutions where the axis of the pendulum is vertical. It is shown that the dynamics of a linearized system is governed by the classical Mathieu equation, in which one of the coefficients depends on the angular velocity of rotation of the pendulum. The regions of stability can be represented using the Ince – Strutt diagram [31, 35, 36, 37]. In addition, it is shown that, in the presence of rotations of the pendulum, a vibrostabilization without loss of contact with the plane is possible.

Steady motions of a spherical robot can be stabilized not only by means of constant rotations, but also by using feedback. By feedback we mean an additional control action that stabilizes the required motion of the system. The feedback allows an adjustment of the current state of the system, increases the efficiency of control and the stability of motion. For example, it can ensure the asymptotic behavior of stability of the required solution. This approach to the implementation of feedback for the stabilization of steady motions of a spherical pendulum-type robot on a fixed plane is discussed in [38, 39].

In this paper, which continues the study [34], we solve a similar problem, that of developing the feedback of the spherical robot for the stabilization of its steady motions on a vibrating plane. In particular, we consider the stabilization of the upper vertical rotation, a partial solution where the pendulum rotates about its symmetry axis, which is vertical, and the center of mass lies above the suspension point.

In the approach suggested, the feedback depends on the phase variables (the current position and velocities) and does not depend on the specific type of trajectory. It is shown that, for the chosen type of the control function, the linearized equations of pendulum motion near vertical rotations are the damped Mathieu equation [35, 37, 40, 41, 42, 43]. For the above-mentioned partial solution, regions of asymptotic stability are constructed and possible bifurcations are analyzed. Special attention is given to the problem of the stability of periodic solutions emerging as the vertical rotations lose stability.

2 EQUATIONS OF MOTION

Consider a spherical robot rolling without slipping on a horizontal plane. The spherical robot is modeled by a spherical shell of radius $R$ and mass $M$ and an axisymmetric spherical pendulum (Lagrange top) of mass $m$ fastened at its center. The center of mass of the pendulum is displaced along the symmetry axis from the center of the shell by distance $\rho$ (see Fig. 1).

Also, we assume that

• the underlying plane performs vertical periodic oscillations according to the law $\xi(t)$;

• at the point of attachment of the pendulum to the shell $C$, the control torque $\boldsymbol{U}$ is applied. In this work we will assume that the control torque $\boldsymbol{U}$ models a feedback and depends only on the angular velocity of the pendulum and its orientation.

We will describe the position of the system by the coordinates of the geometric center of the shell $\boldsymbol{r}_{C}=(x_{C},y_{C},z_{C})$ in the fixed reference frame. Orientation of the pendulum and the shell will be defined by two orthogonal matrices of rotation ${\bf Q}$ and ${\bf S}$, respectively.

Fig. 1

A schematic model of a spherical robot.

Following [34], one can show that the dynamical equations of motion for the system with additional control torque can be represented as ¹²

$$\begin{aligned} &\displaystyle{\bf\hat{I}}\left(\dot{\boldsymbol{\Omega}}+\boldsymbol{\omega}\times\boldsymbol{\Omega}\right)=m\rho R\left(\dot{\boldsymbol{\omega}}\times\boldsymbol{e}_{3}-(\boldsymbol{\omega}\times\boldsymbol{e}_{3})\times\boldsymbol{\omega}\right)\times\boldsymbol{\gamma}+\boldsymbol{U},\\ &\displaystyle{\bf{I}}\dot{\boldsymbol{\omega}}+\boldsymbol{\omega}\times{\bf{I}}\boldsymbol{\omega}=m\rho R\left((\dot{\boldsymbol{\Omega}}+\boldsymbol{\omega}\times\boldsymbol{\Omega})\times\boldsymbol{\gamma}\right)\times\boldsymbol{e}_{3}+m\rho({\rm g}+\ddot{\xi}(t))\boldsymbol{\gamma}\times\boldsymbol{e}_{3}-\boldsymbol{U},\\ &\displaystyle{\bf\hat{I}}=\tilde{I}{\bf E}-(m+M)R^{2}\boldsymbol{\gamma}\otimes\boldsymbol{\gamma},\quad\tilde{I}=I_{0}+(m+M)R^{2},\quad{\bf{I}}=\mathop{\rm diag}\nolimits(i_{1}+m\rho^{2},i_{1}+m\rho^{2},i_{3}), \end{aligned}$$

(2.1)

where $\boldsymbol{\omega}$, $\boldsymbol{\Omega}$ are the angular velocities of the pendulum and the shell, respectively, $\boldsymbol{\gamma}$ is the normal vector to the supporting plane, ${\bf{I}}$ is the tensor of inertia of the pendulum relative to its point of attachment $C$, $i_{1}$, $i_{3}$ are the principal moments of inertia of the pendulum, ${\bf E}$ is a $3\times 3$ identity matrix, $I_{0}$ is the principal moment of inertia of the shell, $\boldsymbol{e}_{3}=(0,0,1)$, and $\rm g$ is the free-fall acceleration. We will assume that the control torque $\boldsymbol{U}$ depends only on the angular velocity of the pendulum and its orientation, i. e. $\boldsymbol{U}=\boldsymbol{U}(\boldsymbol{\omega},\boldsymbol{\gamma})$.

To form a closed system of equations of motion, Eqs. (2.1) should be supplemented with the equation for evolution of the vector $\boldsymbol{\gamma}$ [34]:

$$\dot{\boldsymbol{\gamma}}+\boldsymbol{\omega}\times\boldsymbol{\gamma}=0.$$

(2.2)

Equations (2.1), (2.2) describe the dynamics of the system in the moving reference frame. For a complete description of the system’s dynamics in absolute space, Eqs. (2.1), (2.2) should be supplemented with kinematic relations for the matrices ${\bf Q}$ and ${\bf S}$ and the position of the center of the shell $\boldsymbol{r}_{C}$ (see [34]).

For a more convenient description of the system we introduce new variables $\boldsymbol{K}$ and $\boldsymbol{M}$ which are related to the angular velocities $\boldsymbol{\omega}$ and $\boldsymbol{\Omega}$ by

$$\boldsymbol{K}=\dfrac{j_{3}}{i_{3}}{\bf\hat{I}}\boldsymbol{\Omega}-\dfrac{\boldsymbol{\gamma}\times(\boldsymbol{e}_{3}\times\boldsymbol{\omega})}{\chi},\quad\boldsymbol{M}={\bf D}\boldsymbol{\omega},$$

where

$$\begin{gathered}\displaystyle{\bf D}={\bf{J}}+\boldsymbol{e}_{3}\times\boldsymbol{\gamma}\otimes\boldsymbol{e}_{3}\times\boldsymbol{\gamma},\qquad{\bf{J}}=\mathop{\rm diag}\nolimits(j_{1},j_{1},j_{3}),\\ \displaystyle j_{1}=\dfrac{i_{1}+m\rho^{2}}{\chi m\rho R}-1,\quad j_{3}=\dfrac{i_{3}}{\chi m\rho R},\quad\chi=\dfrac{m\rho R}{\tilde{I}}.\end{gathered}$$

The vectors $\boldsymbol{M}$ and $\boldsymbol{K}$ are dimensionless ³ analogs of the angular momentum and the projections of Chaplygin’s vector integral onto the axes of the moving coordinate system [34, 44, 45].

The equations for $\boldsymbol{M}$ and $\boldsymbol{\gamma}$ decouple from the complete system and form the reduced system

$$\begin{gathered}\displaystyle\dot{\boldsymbol{M}}+\boldsymbol{\omega}\times\boldsymbol{M}=(\boldsymbol{e}_{3},\boldsymbol{\gamma}\times\boldsymbol{\omega})\boldsymbol{\gamma}\times(\boldsymbol{e}_{3}\times\boldsymbol{\omega})+(\tilde{{\rm g}}+\ddot{\tilde{\xi}}(t))\boldsymbol{\gamma}\times\boldsymbol{e}_{3}+\chi(\boldsymbol{u}\times\boldsymbol{\gamma})\times\boldsymbol{e}_{3}-\boldsymbol{u},\\ \displaystyle\dot{\boldsymbol{\gamma}}+\boldsymbol{\omega}\times\boldsymbol{\gamma}=0,\end{gathered}$$

(2.3)

where $\boldsymbol{\omega}={\bf D}^{-1}\boldsymbol{M}$, and the following notation is introduced:

$$\boldsymbol{u}=\dfrac{j_{3}}{i_{3}}\boldsymbol{U},\quad\tilde{{\rm g}}=\dfrac{{\rm g}}{\ell},\quad\tilde{\xi}(t)=\dfrac{\xi(t)}{\ell},\quad\ell=\chi R.$$

The evolution of the vector $\boldsymbol{K}$ is determined by the equation

$$\dot{\boldsymbol{K}}+\boldsymbol{\omega}\times\boldsymbol{K}=\boldsymbol{u}.$$

(2.4)

The reduced equations (2.3) contain four independent parameters:

$$\tilde{{\rm g}}>0,\quad 0<\chi<1,\quad j_{1}>0,\quad 0<j_{3}<2j_{1},$$

and depend on the function $\tilde{\xi}(t).$

3 DYNAMICS OF A FREE SYSTEM

Before investigating the system with feedback satisfying the above-mentioned conditions, we present some results of [34] for a free system (with $\boldsymbol{u}=0$) which will be used in what follows.

The reduced system (2.3) with $\boldsymbol{u}=0$ admits the following integrals of motion:

• a geometric integral $\boldsymbol{\gamma}^{2}=1$;

• a Lagrange integral (in dimensionless form)

  $$\mathcal{K}_{1}=(\boldsymbol{M},\boldsymbol{e}_{3})=M_{3},$$

  (3.1)

  which coincides, up to a constant multiplier, with the projection of the angular velocity of the pendulum onto the axis of its symmetry

  $$\mathcal{K}_{1}=j_{3}\omega_{3};$$

• a dimensionless analog of the area integral

  $$\mathcal{K}_{2}=(\boldsymbol{M},\boldsymbol{\gamma}).$$

  (3.2)

Remark 1

Note that in the case of free motion, the equations of motion also admit the Chaplygin vector integral [44, 45] referred to the fixed coordinate system

$${\boldsymbol{K}}^{abs}={\bf Q}^{\top}{\boldsymbol{K}}={\rm const}.$$

(3.3)

In particular, the integrals (3.3) allow one to calculate the angular velocity of the spherical shell $\boldsymbol{\Omega}$ using explicit algebraic expressions of $\boldsymbol{\omega}$ and ${\bf Q}$ instead of the Eq. (2.4).

Furthermore, this system possesses a standard invariant measure. Thus, Eqs. (2.3) can be represented in Hamiltonian form. For details on the Hamiltonian nature of these equations and the analogy of the dynamics of the top fastened at the center of the spherical shell rolling without slipping on a plane with the dynamics of the spherical top on a smooth plane, see, for example, [22, 34, 46].

In addition, when $\boldsymbol{u}=0$ the system (2.3) has the symmetry field

$$\begin{gathered}\displaystyle\hat{\boldsymbol{\varsigma}}=M_{1}\dfrac{\partial}{\partial M_{2}}-M_{2}\dfrac{\partial}{\partial M_{1}}+\gamma_{1}\dfrac{\partial}{\partial\gamma_{2}}-\gamma_{2}\dfrac{\partial}{\partial\gamma_{1}},\end{gathered}$$

(3.4)

which corresponds to rotation of the pendulum about its symmetry axis $Ox_{3}$.

In [34], a reduction of the free system (2.3) by the first integrals (3.1), (3.2) and by the symmetry field (3.4) was performed. It was shown that, despite the existence of nonholonomic constraints, the free system can be reduced to a Hamiltonian system with one and a half degrees of freedom.

Also, the following integrals of the symmetry field were chosen as variables of the reduced system:

$$k_{3}=\dfrac{M_{2}\gamma_{1}-M_{1}\gamma_{2}}{k(j_{1}+1-\gamma_{3}^{2})}\quad\text{and}\quad\gamma_{3},$$

where

$$k=\sqrt{\dfrac{1-\gamma_{3}^{2}}{j_{1}+1-\gamma_{3}^{2}}}.$$

On the fixed level set of the integrals (3.1) and (3.2)

$$\mathcal{K}_{1}=k_{1},\quad\mathcal{K}_{2}=k_{2}$$

the equations of motion take the form

$$\begin{array}[]{l}\dot{k}_{3}=-\dfrac{k(k_{2}-k_{1}\gamma_{3})(k_{2}\gamma_{3}-k_{1})}{j_{1}(1-\gamma_{3}^{2})^{2}}-k(\tilde{{\rm g}}+\ddot{\tilde{\xi}}(t)),\\ \dot{\gamma}_{3}=kk_{3}.\end{array}$$

If it is necessary to reconstruct the dynamics of the complete system, the expressions for the first integrals (3.1), (3.2), (3.3) and the kinematic relations [34] need to be used.

4 DYNAMICS OF THE SYSTEM WITH FEEDBACK

4.1 Choice of Control for the Stabilization of Vertical Rotations of the Pendulum

In the absence of controls, Eqs. (2.3) admit two one-parameter families of partial solutions (fixed points)

$$\sigma^{\pm}:\boldsymbol{M}=(0,0,M_{3}),\boldsymbol{\gamma}=(0,0,\pm 1).$$

(4.1)

These solutions correspond to rotations of the pendulum with the constant angular velocity $\omega_{3}=M_{3}/j_{3}$ about its symmetry axis, which is directed vertically. In what follows, we will call such rotations upper ($\sigma^{+}$) or lower ($\sigma^{-}$) vertical rotations. The spherical shell on these solutions moves uniformly in a straight line in absolute space and rotates with constant angular velocity. The stability of these solutions in the absence of feedback is analyzed in detail in [34].

The purpose of this paper is to address the problem of stabilizing vertical rotations $\sigma^{\pm}$ using feedback. In what follows, we will be concerned primarily with the upper vertical rotation $\sigma^{+}$. For the lower vertical rotation $\sigma^{-}$, all results can be obtained in a similar way.

As stated above, we implement the feedback by means of the external control torque $\boldsymbol{u}$, which depends on the variables $\boldsymbol{M}$ and $\boldsymbol{\gamma}$ of the reduced system (2.3) or, which is the same, on the variables $\boldsymbol{\omega}$ and $\boldsymbol{\gamma}$. Thus, the feedback depends on the orientation and the angular velocity of the internal pendulum. The position and the velocities of the external shell are not taken into account in this case.

We also require that the control $\boldsymbol{u}$ satisfy the following conditions.

• the vertical rotations $\sigma^{\pm}$ must persist after adding the control, i.e., the control must vanish on the solutions (4.1).

• After adding the control $\boldsymbol{u}$, the type of fixed points of the system must change from the center (without feedback) to a stable focus or node (with feedback).

To stabilize the vertical rotations of the pendulum, we choose a feedback (control torque $\boldsymbol{u}$) in Eqs. (4.5) in such a way that it damps the horizontal components of the angular velocity of the pendulum. This can be achieved in different ways. In this paper, we choose a feedback in the following form:

$$\boldsymbol{u}=\mu(\boldsymbol{\omega},\boldsymbol{\tau})\boldsymbol{\tau}+\eta\dfrac{(\boldsymbol{M},\boldsymbol{e}_{3}-\boldsymbol{\gamma})}{\sqrt{1-\gamma_{3}^{2}}}\boldsymbol{\nu}.$$

(4.2)

Here $\mu$ and $\eta$ are positive constant coefficients, and $\boldsymbol{\nu}$ and $\boldsymbol{\tau}$ are unit vectors lying in the horizontal plane and directed along the projection of the symmetry axis of the top onto the plane and perpendicularly to it (see Fig. 2). They are given by the equalities

$$\boldsymbol{\tau}=\dfrac{\boldsymbol{e}_{3}\times\boldsymbol{\gamma}}{\sqrt{1-\gamma_{3}^{2}}},\quad\boldsymbol{\nu}=\dfrac{\boldsymbol{\gamma}\times(\boldsymbol{e}_{3}\times\boldsymbol{\gamma})}{\sqrt{1-\gamma_{3}^{2}}}.$$

(4.3)

Fig. 2.

A spherical robot with a top on a plane. In this figure, vectors $\boldsymbol{e}_{1}$, $\boldsymbol{e}_{2}$ and $\boldsymbol{e}_{3}$ are the unit vectors of the moving reference frame $Cx_{1}x_{2}x_{3}$ rigidly attached to the pendulum.

In this system, the feedback (4.2) artificially generates the viscous rolling resistance torque, which decreases the horizontal component of the angular velocity of the pendulum. It is easy to verify that the control (4.2) vanishes on the solutions $\sigma^{\pm}$. Thus, the vertical rotations will persist after adding such a feedback. In addition, its “dissipative” nature leads, as we will see below, to the required asymptotic behavior near the vertical rotations.

4.2 Reduction to a System with One and a Half Degrees of Freedom

1. Reduction by the symmetry field. By analogy with the free system [34], we perform for the system (2.3) a reduction by the symmetry field (3.4). To do so, we pass from the vector $\boldsymbol{M}$ to variables $k_{1}$, $k_{2}$ and $k_{3}$ which are expressed in terms of the vectors $\boldsymbol{M}$ and $\boldsymbol{\gamma}$ as follows:

$$k_{1}=(\boldsymbol{M},\boldsymbol{e}_{3}),\quad k_{2}=(\boldsymbol{M},\boldsymbol{\gamma}),\quad k_{3}=\dfrac{(\boldsymbol{\gamma}\times{\bf D}^{-1}\boldsymbol{M},\boldsymbol{e}_{3})}{k}.$$

The inverse transformation has the form

$$\begin{array}[]{c}M_{1}=\dfrac{(k_{2}-k_{1}\gamma_{3})\gamma_{1}-(j_{1}+1-\gamma_{3}^{2})kk_{3}\gamma_{2}}{1-\gamma_{3}^{2}},\\ M_{2}=\dfrac{(k_{2}-k_{1}\gamma_{3})\gamma_{2}+(j_{1}+1-\gamma_{3}^{2})kk_{3}\gamma_{1}}{1-\gamma_{3}^{2}},\\ M_{3}=k_{1}.\end{array}$$

(4.4)

In the new variables the reduced system of equations can be written as

$$\begin{array}[]{l}\dot{k}_{1}=-(\boldsymbol{u},\boldsymbol{e}_{3}),\\ \dot{k}_{2}=(\chi\left(\boldsymbol{u}\times\boldsymbol{\gamma})\times\boldsymbol{e}_{3}-\boldsymbol{u},\boldsymbol{\gamma}\right),\\ \dot{k}_{3}=-\dfrac{k(k_{2}-k_{1}\gamma_{3})(k_{2}\gamma_{3}-k_{1})}{j_{1}(1-\gamma_{3}^{2})^{2}}-k(\tilde{{\rm g}}+\ddot{\tilde{\xi}}(t))+\dfrac{k(\boldsymbol{u}\times\boldsymbol{\gamma},\boldsymbol{e}_{3})(1+\chi\gamma_{3})}{1-\gamma_{3}^{2}},\\ \dot{\gamma}_{3}=kk_{3},\end{array}$$

(4.5)

where we can represent the control (4.2), using (4.3) and (4.4), in the form

$$\boldsymbol{u}=\dfrac{\mu kk_{3}}{\sqrt{1-\gamma_{3}^{2}}}\boldsymbol{\tau}+\eta\dfrac{k_{1}-k_{2}}{\sqrt{1-\gamma_{3}^{2}}}\boldsymbol{\nu}.$$

(4.6)

The time dependence of the remaining two components $\gamma_{1}$ and $\gamma_{2}$ is defined by the algebraic relations

$$\gamma_{1}=\sqrt{1-\gamma_{3}^{2}}\cos\varphi,\quad\gamma_{2}=\sqrt{1-\gamma_{3}^{2}}\sin\varphi,$$

where $\varphi$ is the angle of proper rotation which is defined by the quadrature

$$\dot{\varphi}=\dfrac{\gamma_{3}(k_{2}-k_{1}\gamma_{3})}{j_{1}(1-\gamma_{3}^{2})}-\dfrac{k_{1}}{j_{3}}.$$

2. Reduction to a system of three equations. Equations (4.5) with a control of the form (4.6) admits an additional integral of motion:

$$\mathcal{C}=\dfrac{\chi k_{1}+k_{2}}{1+\chi}.$$

We perform a reduction of the equations of motion to the level set of this integral. To do so, we introduce a new variable

$$q=\dfrac{(\boldsymbol{M},\boldsymbol{e}_{3}-\boldsymbol{\gamma})}{1+\chi}=\dfrac{k_{1}-k_{2}}{1+\chi}.$$

On the fixed level set of the integral $\mathcal{C}=c$ the dynamics of the system is described by the system of three differential equations

$$\displaystyle\begin{aligned} \displaystyle\dot{q}&\displaystyle=-\eta(1+\chi)q,\\ \displaystyle\dot{k}_{3}&\displaystyle=\dfrac{k(c(1-\gamma_{3})-q(\chi+\gamma_{3}))(c(1-\gamma_{3})+q(1+\chi\gamma_{3}))}{j_{1}(1-\gamma_{3}^{2})^{2}}-k(\tilde{{\rm g}}+\ddot{\tilde{\xi}}(t))-\dfrac{\mu k^{2}k_{3}(1+\chi\gamma_{3})}{1-\gamma_{3}^{2}},\\ \displaystyle\dot{\gamma}_{3}&\displaystyle=kk_{3}.\end{aligned}$$

(4.7)

3. Invariant manifold. It follows from the first equation of (4.7) that the following proposition holds.

Proposition 1.

The system (4.7) admits an asymptotically stable invariant manifold

$$\mathcal{M}_{q}=\{(q,k_{3},\gamma_{3})\mid q=0\}.$$

The attraction to the invariant manifold is exponential in time and is defined by the value of the feedback parameter $\eta$. This implies, in particular, that for sufficiently large $\eta$ all trajectories are attracted sufficiently fast to the manifold $\mathcal{M}_{q}$. Further dynamics of the system is defined by the equations of motion on this manifold.

Remark 2

In a similar way, using the second term in the control (4.6) in the form

$$\eta^{*}\dfrac{k_{1}+k_{2}}{\sqrt{1-\gamma_{3}^{2}}}\boldsymbol{\nu},$$

one can stabilize the invariant manifold $k_{1}=-k_{2}$, on which the solution corresponding to the lower vertical rotation of the pendulum lies.

Setting $q=0$ in (4.7), we obtain a time-dependent system of equations with one degree of freedom which governs the dynamics on the manifold $\mathcal{M}_{q}$:

$$\begin{array}[]{l}\dot{k}_{3}=\dfrac{kc^{2}}{j_{1}(1+\gamma_{3})^{2}}-k(\tilde{{\rm g}}+\ddot{\tilde{\xi}}(t))-\dfrac{\mu k^{2}k_{3}(1+\chi\gamma_{3})}{1-\gamma_{3}^{2}},\\ \dot{\gamma}_{3}=kk_{3}.\end{array}$$

(4.8)

In contrast to the free system [34], the system we consider here is not Hamiltonian and differs from that discussed in [34] by the last term in the first equation of (4.8).

To investigate the system (4.8), we pass from the variables $k_{3}$ and $\gamma_{3}$ to the nutation angle $\vartheta$ and the corresponding momentum $p$ by making the change of variables

$$\gamma_{3}=\cos\vartheta,\quad k_{3}=-\dfrac{p}{\sqrt{j_{1}+\sin^{2}\vartheta}}.$$

Furthermore, in this paper, as in [34], we will assume that the the plane oscillates according to the law

$$\xi(t)=A_{0}\cos\Omega t,$$

(4.9)

where $A_{0}$ and $\Omega$ are the amplitude and the frequency of oscillations of the plane, respectively.

Introducing a dimensionless time variable and defining new dimensionless parameters:

$$t=\dfrac{\tau}{\Omega},\quad p=j_{1}\Omega p_{\vartheta},\quad c=j_{1}\Omega c_{\vartheta},\quad\tilde{{\rm g}}=j_{1}\Omega^{2}\gamma,\quad A_{0}=j_{1}\ell\delta,\quad\mu=\dfrac{j_{1}\Omega\mu_{0}}{1+\chi},$$

we can write the system (4.8) as

$$\begin{array}[]{l}\dot{\vartheta}=\dfrac{j_{1}p_{\vartheta}}{j_{1}+\sin^{2}\vartheta},\\ \dot{p}_{\vartheta}=\dfrac{j_{1}p_{\vartheta}^{2}\sin\vartheta\cos\vartheta}{(j_{1}+\sin^{2}\vartheta)^{2}}-\dfrac{c_{\vartheta}^{2}\sin\vartheta}{(1+\cos\vartheta)^{2}}+(\gamma-\delta\cos\tau)\sin\vartheta-\dfrac{j_{1}\mu_{0}p_{\vartheta}(1+\chi\cos\vartheta)}{(1+\chi)(j_{1}+\sin^{2}\vartheta)},\end{array}$$

(4.10)

where the dot denotes differentiation with respect to the new time variable $\tau$.

For $\mu_{0}=0$ the system (4.10) becomes Hamiltonian. The corresponding Hamiltonian is presented in [34].

Remark 3.

The system (4.10) is defined in the strip $\vartheta\in[0,\pi]$, $p_{\vartheta}\in\mathbb{R}$. However, to facilitate the analysis of the zero solution, we will consider it on the cylinder $\vartheta\in[-\pi,\pi]$, $p_{\vartheta}\in\mathbb{R}$. Also, it should be remembered that the points with opposite values of $\vartheta,p_{\vartheta}$ in the complete system correspond to the same solutions rotated by the angle $\pi$ about the axis $Cx_{3}$.

Next, we address the question of stability of the upper equilibrium point depending on the system parameters, feed back parameters and oscillations of the plane.

4.3 Stability of Vertical Rotations

The system (4.10) admits the obvious partial solution (equilibrium point)

$$\vartheta=0,\quad p_{\vartheta}=0.$$

(4.11)

In the complete system, this solution corresponds to the upper vertical rotations $\sigma^{+}$, under which the axis of the pendulum is directed along the vertical and the center of mass of the pendulum lies above the suspension point ($\boldsymbol{\gamma}=\boldsymbol{e}_{3}$).

To investigate the linear stability of the equilibrium point (4.11), we expand Eqs. (4.10) as a Taylor series near it up to linear terms in the variables $\vartheta$ and $p_{\vartheta}$.

We represent the resulting equations as a linear homogeneous equation with periodic coefficients:

$$\dfrac{d^{2}\vartheta}{d\tau^{2}}+\mu_{0}\dfrac{d\vartheta}{d\tau}+\left(\kappa+\delta\cos\tau\right)\vartheta=0,$$

(4.12)

where

$$\kappa=\dfrac{c_{\vartheta}^{2}}{4}-\gamma.$$

(4.13)

Equation (4.12) is the damped Mathieu equation [35, 37, 42, 43] and contains three independent parameters

$$\kappa,\quad\delta\geqslant 0,\quad\mu_{0}\geqslant 0.$$

The value $\kappa$ depends on the parameter $\gamma>0$ and the level set of the integral $c_{\vartheta}$, which can take arbitrary values (the change of sign corresponds to a change in the direction of rotation of the pendulum about its symmetry axis). Thus, for the upper vertical rotations $\sigma^{+}$ the parameter $\kappa$ can take arbitrary values. However, when the value of $\gamma$ is fixed, the parameter $\kappa$ changes in the interval $[-\gamma,+\infty)$.

Next, we define the regions of linear stability of the solution (4.11) to Eq. (4.12) on the plane $(\kappa,\delta)$ for fixed values of the parameter $\mu_{0}$.

1. Free system. In the case $\mu_{0}=0$ (without control) the stability of the solution (4.11) to the Mathieu equation (4.12) has been well studied [31, 35, 36, 37, 47]. We present some results we will need in what follows.

Let ${\bf{X}}(\tau)$ denote a fundamental solution matrix of the linearized system (4.10) which satisfies the condition ${\bf{X}}(0)={\bf E}$ (in this case, ${\bf E}$ is a $2\times 2$ identity matrix). As is well known [48], the problem of investigating the stability of the zero solution to Eq. (4.12) reduces to estimating the multipliers, the eigenvalues of the monodromy matrix ${\bf{X}}(2\pi)$, which we denote as $\varrho_{i},i=1,2$.

The corresponding characteristic equation for calculation of the multipliers has the form

$$\varrho^{2}-2A(\kappa,\delta)\varrho+1=0,$$

(4.14)

where $A(\kappa,\delta)=x_{11}(2\pi)=x_{22}(2\pi)$ are the diagonal elements of the matrix ${\bf{X}}(2\pi)$.

From (4.14) we find

$$\varrho_{1,2}=A(\kappa,\delta)\pm\sqrt{A(\kappa,\delta)^{2}-1},$$

(4.15)

whence it follows that the boundaries of the stability regions are determined by the values of the parameters $\kappa$ and $\delta$ for which $|A(\kappa,\delta)|=1$, and are represented in this case by the Ince – Strutt diagram (Fig. 3).

• For $|A(\kappa,\delta)|<1$ the multipliers (4.15) are complex conjugate and $|\varrho_{1}|=|\varrho_{2}|=1$. The fixed point (4.11) is a stable point of center type. The corresponding regions (of stability) on the parameter plane $(\kappa,\delta)$ are shown in white in Fig. 3.

• For $|A(\kappa,\delta)|>1$ the absolute value of one of the multipliers is greater than unity and the fixed point (4.11) is an unstable saddle point. The corresponding regions (of instability) on the parameter plane $(\kappa,\delta)$ are shown in dark gray in Fig. 3.

For $|A(\kappa,\delta)|=1$, i.e., on the bifurcation curves separating the regions of stability and instability, the multipliers are equal to $\varrho_{1}=\varrho_{2}=\pm 1$. For values of the parameters $\kappa$ and $\delta$ lying on these curves, the elements of the fundamental matrix ${\bf{X}}(\tau)$ are even ${\rm{ce}}_{n}(\tau/2,-2\delta)$ $(n=0,1,2,\ldots)$ or odd ${\rm{se}}_{n}(\tau/2,-2\delta)$ $(n=1,2,\ldots)$ Mathieu functions of the first kind ($2\pi$-periodic for even $n$, $\varrho_{1}=\varrho_{2}=1$, or $4\pi$-periodic for odd $n$, $\varrho_{1}=\varrho_{2}=-1$) and the corresponding Mathieu functions of the second kind (aperiodic and unbounded). Thus, in the linear approximation, the upper vertical rotations $\sigma^{+}$ on these boundaries are unstable (for more details, see [31, 47]).

The boundaries of these regions emanate from the points $\kappa_{l}$ and $\kappa_{r}$ (with $\delta=0$):

$$\kappa_{l}=\dfrac{(n-1)^{2}}{4},\quad\kappa_{r}=\dfrac{n^{2}}{4},\quad n=1,2,\ldots$$

Fig. 3.

Regions of stability (white) and instability (dark gray) of the pendulum’s vertical rotations $\sigma^{+}$ on the parameter plane $(\kappa,\delta)$.

In the absence of vibrations $(\delta=0)$ and control torque $(\mu_{0}=0)$ the stability of vertical rotations is defined by the Maievskii condition $c_{\vartheta}^{2}\geqslant 4\gamma$. In the chosen notation it takes the form

$$\kappa\geqslant 0.$$

2. System with feedback. To analyze Eq. (4.12) for $\mu_{0}>0$, we make a change of function

$$\vartheta(\tau)=e^{-\mu_{0}\tau/2}\xi(\tau).$$

For the new variable, Eq. (4.12) is reduced to the standard Mathieu equation:

$$\dfrac{d^{2}\xi}{d\tau^{2}}+\left(\tilde{\kappa}+\delta\cos\tau\right)\xi=0,\quad\tilde{\kappa}=\kappa-\dfrac{\mu_{0}^{2}}{4}.$$

(4.16)

Thus, the characteristic equation corresponding to the fundamental solution matrix of the system (4.12) has the form

$$\varrho^{2}-2e^{-\pi\mu_{0}}A_{\mu}\varrho+e^{-2\pi\mu_{0}}=0,$$

(4.17)

where $A_{\mu}=A(\tilde{\kappa},\delta)$. From (4.17) we find

$$\varrho_{1,2}=e^{-\pi\mu_{0}}(A_{\mu}\pm\sqrt{A_{\mu}^{2}-1}),\quad\varrho_{1}\varrho_{2}=e^{-2\pi\mu_{0}}.$$

(4.18)

This implies that, when $\mu_{0}>0$, on the plane $(\kappa,\delta)$ there exist three regions corresponding to different types of the fixed point (4.11):

• When $|A_{\mu}|<1$ the multipliers (4.18) are complex conjugate and $|\varrho_{1}|=|\varrho_{2}|=e^{-\pi\mu_{0}}<1$. The fixed point (4.11) is a stable focus. The corresponding regions (of stability) on the parameter plane $(\kappa,\delta)$ are shown in white in Fig. 4.

• When $1<|A_{\mu}|<\rm{cosh}\pi\mu_{0}$ the multipliers (4.18) are real and $|\varrho_{1}|<1,|\varrho_{2}|<1$. The fixed point (4.11) is a stable node. The corresponding regions (of stability) on the parameter plane $(\kappa,\delta)$ are shown in light gray in Fig. 4.

• When $|A_{\mu}|>\rm{cosh}\pi\mu_{0}$ the multipliers (4.18) are real, but the absolute value of one of them is greater than unity. The fixed point (4.11) is an unstable saddle point. The corresponding regions (of instability) on the parameter plane $(\kappa,\delta)$ are shown in dark gray in Fig. 4.

When $|A_{\mu}|=1$ (see the dashed curves in Fig. 4) the multipliers are real and equal to each other ($\varrho_{1}=\varrho_{2}=\pm e^{-\pi\mu_{0}}<1$). For values of the parameters $\kappa$ and $\delta$ lying on these curves, the fixed point (4.11) is a stable node. According to (4.16), on the parameter plane $(\kappa,\delta)$ these curves form the Ince – Strutt diagram with displacement along the horizontal axis by $\mu_{0}^{2}/4$ to the right relative to the case of the diagram with $\mu_{0}=0$ (see Fig. 3).

When $|A_{\mu}|=\rm{cosh}\pi\mu_{0}$ (the solid curves in Fig. 4), i.e., on the bifurcation curves separating the regions of stability and instability, the multipliers are real and one of them lies on the unit circle ($\varrho_{1}=\pm 1$, $\varrho_{2}=\pm e^{-2\pi\mu_{0}}$). As in the case $\mu_{0}=0$, for values of the parameters $\kappa$ and $\delta$ lying on these curves, the upper vertical rotations $\sigma^{+}$ are unstable in the linear approximation (one of the linearly independent solutions to equations (4.12) is aperiodic and unbounded).

Fig. 4.

Regions of stability (white and light gray) and instability (dark gray) of vertical rotations of the pendulum on a vibrating plane and with the control torque (4.2) on the parameter plane $(\kappa,\delta)$ for $\mu_{0}=1$. The sign “$+$” denotes bifurcation curves corresponding to $\varrho_{1}=+1$, and the sign “$-$” denotes bifurcation curves corresponding to $\varrho_{2}=-1$.

As an example, Fig. 5 shows the real and imaginary parts of the multipliers versus the parameter $\kappa$ for $\mu_{0}=1$ and $\delta={\rm const}=1.3$ (the horizontal dashed straight line in Fig. 4). Figure 5b is an enlarged fragment of Fig. 5a. Let us consider the resulting dependences in more detail, starting with some negative $\kappa$ from the instability region and moving towards its increase.

Fig. 5.

Characteristic dependences of the real $\rm{Re}\varrho$ (a), (b) and imaginary $\rm{Im}\varrho$ (c) parts of the multipliers on $\kappa$ for constant $\delta=1.3$.

When $\kappa<\kappa_{1}\approx-0.396$ both multipliers are positive, but one of them is larger than unity ($\varrho_{1}>1$). In this case, the fixed point (4.11) is a saddle. The corresponding curve in Fig. 5a is shown as a heavy line. When $\kappa=\kappa_{1}$ a pitch-fork bifurcation occurs, and when $\kappa_{1}<\kappa<\kappa_{2}\approx-0.322$ the fixed point is a stable node ($0<\varrho_{1,2}<1$). For $\kappa=\kappa_{2}$ the multipliers become complex conjugate, and on the interval $\kappa_{2}<\kappa<\kappa_{3}\approx-0.304$ the fixed point (4.11) is a stable focus. As $\kappa$ changes from $\kappa_{2}$ to $\kappa_{3}$, the argument of the multipliers changes by $\pi$. When $\kappa=\kappa_{3}$ a reverse bifurcation occurs, and on the interval $\kappa_{3}<\kappa<\kappa_{4}\approx-0.167$ the fixed point is a stable node, but now with negative multipliers ($\varrho_{1,2}<0,|\varrho_{1,2}|<1$). When $\kappa=\kappa_{4}$ the fixed point (4.11) loses stability via a period-doubling bifurcation, and on the interval $\kappa_{4}<\kappa<\kappa_{5}\approx 0.334$ it is a saddle. When $\kappa=\kappa_{5}$ a reverse period-doubling bifurcation occurs. As $\kappa>\kappa_{5}$ increases further, both multipliers lie all the time inside a unit circle, and the regions corresponding to a stable node ($\rm{Im}\varrho=0$) and to a stable focus ($\rm{Im}\varrho\not=0$) alternate.

The analysis presented above has shown that, after the control (4.2) is added, the regions of stability of the fixed point (4.11) on the parameter plane $(\kappa,\delta)$ expand compared to those for the free system (with $\mu_{0}=0$). In addition, as the plane oscillates in a prescribed manner, the vertical rotations of the pendulum of the spherical robot on the vibrating surface can be stabilized by changing the angular velocity $\omega_{3}$ of rotation of the pendulum (or the parameter $\kappa$ with fixed $\Omega$). Moreover, such rotations are asymptotically stable.

4.4 Analysis of Reaction Forces and the Possibility of Loss of Contact with the Plane

At high frequencies or large oscillation amplitudes of the plane, the supporting reaction can vanish, resulting in the shell losing contact with the plane. In this section we consider vertical rotations for which the shell remains contact with the supporting plane at all times.

For the system considered, the vertical component of the reaction force has the form [34, 49]

$$N_{z}=(M+m)({\rm g}+\ddot{\xi}(t))+m\rho\ddot{\gamma}_{3}.$$

After nondimensionalization and transformation to variables $\vartheta$ and $p_{\vartheta}$, the reaction force $N_{z}$ becomes

$$\begin{aligned} \displaystyle N_{z}=(M+m)(\gamma-\delta\cos\tau)&+\dfrac{m\rho\sin^{2}\vartheta}{\ell(j_{1}+\sin^{2}\vartheta)}\left(\dfrac{4(\kappa+\gamma)}{(1+\cos\vartheta)^{2}}-\gamma+\delta\cos\tau\right)\\ \displaystyle&-\dfrac{m\rho}{\ell(j_{1}+\sin^{2}\vartheta)}\left(\dfrac{j_{1}^{2}\cos\vartheta p_{\vartheta}^{2}}{(j_{1}+\sin^{2}\vartheta)^{2}}-\dfrac{\mu_{0}j_{1}(1+\chi\cos\vartheta)\sin\vartheta p_{\vartheta}}{(1+\chi)(j_{1}+\sin^{2}\vartheta)}\right). \end{aligned} $$

(4.19)

In the general case, the reaction $N_{z}$ (4.19) depends on the phase variables $\vartheta$ and $p_{\vartheta}$ and the system parameters. For steady motions corresponding to vertical rotations $\sigma^{+}$ ($\vartheta=p_{\vartheta}=0$), expression for the reaction force simplifies to

$$N_{st}=N_{z}|_{\vartheta=p_{\vartheta}=0}=(M+m)(\gamma-\delta\cos\tau).$$

The reaction $N_{st}$ is a periodic function of time, and its minimum is defined by the simple expression

$$N_{min}=\min\limits_{\tau}N_{st}=(M+m)(\gamma-\delta).$$

The passage of the minimal reaction through zero corresponds to the loss of contact with the plane. Thus, the condition for steady rotations that there be no loss of contact with the plane can be written as

$$\delta<\delta_{*}=\gamma.$$

(4.20)

Condition (4.20) is identical to the one for vertical rotations of the free system (with $\mu_{0}$=0) and does not depend on the velocity of rotation of the pendulum $\omega_{3}$ (or $\kappa$).

For a fixed value of $\Omega={\rm const}$ condition (4.20) defines a straight line $\delta=\delta_{*}$ on the plane $(\kappa,\delta)$ that corresponds to $N_{min}=0$. This straight line divides the parameter plane $(\kappa,\delta)$ into a region of vertical rotation of the pendulum of the spherical robot [34] with loss of contact with the plane ($\delta>\delta_{*}$) and the one without loss of contact with the plane ($\delta<\delta_{*}$). For the value $\delta=\delta_{*}$ (on the boundary of the region of motion without loss of contact) it is necessary to carry out an additional study. Figure 7b shows the straight line $\delta=\delta_{*}$ as a dashed line for $\kappa>\kappa_{1}$.

For other (non-steady-state) solutions the critical value of $\delta_{*}$, which separates motions with and without loss of contact with the plane, now depends on the velocity of rotation of the pendulum $\omega_{3}$ (or $\kappa$). In order to find this dependence $\delta_{*}(\kappa)$ for some solution $\vartheta(\tau),p_{\vartheta}(\tau)$, it is necessary to solve the equation

$$\min\limits_{\tau}N_{z}\left|{}_{\vartheta=\vartheta(\tau),p_{\vartheta}=p_{\vartheta}(\tau)}\right.=0$$

for $\delta$. As a rule, this equation is solved numerically.

An example of the dependence $\delta_{*}(\kappa)$ for periodic solutions emerging as the zero solution (4.11) loses stability is shown in Fig. 7b as a dashed line in the region $\kappa<\kappa_{1}$ (i.e., in the region of existence of these solutions). We consider these solutions in more detail in the next section.

4.5 Stability of Other Partial Solutions

As discussed above, when the zero solution (4.11) passes across the boundaries of the stability region (the solid curves in Fig. 4) it loses stability and becomes a saddle. The loss of stability occurs via pitch-fork bifurcations or period-doubling bifurcations [50, 51]. In the former case, one of the multipliers becomes equal to unity and stable nodes emerge near the saddle. The corresponding bifurcation curves are denoted in Fig. 4 by the sign “$+$”. In the latter case, one of the multipliers becomes equal to minus unity and stable foci emerge near the saddle. The corresponding bifurcation curves in Fig. 4 are denoted by the sign “$-$”. Numerical analysis has shown that no other scenarios of the loss of stability of the zero fixed point are observed in the system considered [51].

We consider the emerging periodic solutions. Specifically, we examine the motion of a spherical robot on a plane performing periodic oscillations in a prescribed manner (4.9) with amplitude $A_{0}$ and frequency $\Omega$

$$\Omega=30\text{s}^{-1},\quad A_{0}=5\cdot 10^{-3}\text{m}.$$

(4.21)

For calculations we will use the mass-geometric characteristics of a full-scale specimen of the spherical robot with a pendulum drive (see Fig. 6):

$$\begin{array}[]{c}m=0.306\text{kg},\quad M=0.078\text{kg},\quad R=0.06\text{m},\quad\rho=0.028\text{m},\\ i_{1}=0.3214\cdot 10^{-3}\text{kg}\cdot\text{m}^{2},\quad i_{3}=0.5169\cdot 10^{-3}\text{kg}\cdot\text{m}^{2},\quad I_{0}=0.183\cdot 10^{-3}\text{kg}\cdot\text{m}^{2}.\\ \end{array}$$

Fig. 6.

A full-scale specimen of the spherical robot with an axisymmetric pendulum drive.

The dimensionless parameters of the problem for this full-scale specimen of a spherical robot and the oscillation parameters (4.21) have the following values:

$$\ell\approx 0.0197,\quad j_{1}\approx 2.325,\quad j_{3}\approx 3.062,\quad\chi\approx 0.33,\quad\delta\approx 0.11.$$

(4.22)

The chosen value of the parameter $\delta$ lies in the region of motion without loss of contact with the plane for the upper vertical rotations (with $\kappa>\kappa_{1}$) and the emerging periodic solutions (with $\kappa<\kappa_{1}$). The corresponding line on the plane $(\kappa,\delta)$ in Fig. 7b is shown as a dash-dotted line.

Furthermore, as pointed out above, when the values of $\Omega$ and $\gamma$ are fixed, it follows from (4.13) that there exists the smallest possible value of the parameter $\kappa$:

$$\kappa_{*}=-\gamma\approx-0.273.$$

This value corresponds to the absence of rotation of the internal pendulum.

Consider a bifurcation of the zero solution (4.11) that arises as the parameter $\kappa$ is varied. As stated above, when $\kappa=\kappa_{1}$ a pitch-fork bifurcation occurs and two new stable periodic solutions ($\sigma_{1}^{+}$ and $\sigma_{1}^{-}$) arise near the zero solution. In the complete system, these solutions correspond to identical periodic nutation axes turned relative to each other by angle $\pi$ about the vertical (see Remark 3 in Section 4.2).

To analyze the stability of these solutions, we construct a bifurcation diagram on the plane $(\kappa,\vartheta)$. We will lay off the minimal and maximal values of $\vartheta$ for the motion along the solutions $\sigma^{\pm}$ on the ordinate. It is obvious that for the steady-state solution (4.11) these values coincide and are equal to zero. The corresponding bifurcation diagram is presented in Fig. 7a.

Fig. 7.

(a) Bifurcation diagram on the plane $(\kappa,\vartheta)$ for the parameters (4.22) ($\kappa_{1}\approx-0.03$, $\kappa_{2}\approx 0.244$); (b) regions of stability of the zero solution $\sigma^{+}$ (to the right of the heavy solid curve) and the periodic solutions $\sigma_{1}^{\pm}$ (to the left of it) on the parameter plane $(\kappa,\delta)$.

It can be seen from the figure that, as $\kappa$ decreases (i.e., as the angular velocity of rotation of the pendulum $\omega_{3}$ decreases), the periodic solutions $\sigma_{1}^{\pm}$ move farther and farther away from the upper vertical rotation. And when $\kappa=\kappa_{*}$, these solutions tend to the lower vertical rotation $\sigma^{-}$, which at this value of $\kappa$ corresponds to a pendulum hanging vertically without rotation.

Figure 7b shows a stability diagram of the solution $\sigma^{+}$, supplemented by regions of stability of the periodic solutions $\sigma_{1}^{\pm}$. As in Fig. 4, light gray denotes regions in which the solution $\sigma^{+}$ is a stable node, white corresponds to regions where it is a stable focus, and dark gray indicates regions in which it is a saddle. The heavy curve in this figure bounds the region of stability of the zero solution. To the right of the heavy curve, the solution $\sigma^{+}$ is stable (a node in the light gray region and a focus in the white region). To the left of the heavy curve, the solution $\sigma^{+}$ is unstable, and the periodic solutions emerging from it are stable. In the shaded regions, the solutions $\sigma_{1}^{\pm}$ are stable nodes, and in the unshaded ones they are stable foci. It can be seen from the figure that the periodic solutions $\sigma_{1}^{\pm}$ remain asymptotically stable for all $\kappa_{*}<\kappa<\kappa_{1}$, changing their type sequentially from node to focus and vice versa.

Thus, at values of the pendulum’s angular velocity of rotation, $\omega_{3}$, corresponding to values of $\kappa$ smaller than those of $\kappa_{1}$, the feedback (4.2) will stabilize periodic nutations and not vertical rotations. In the example considered, the corresponding condition for stabilization of nutations has the form

$$\omega_{3}<\omega_{3}^{*}\approx 22.35\text{rad/s}.$$

As the angular velocity of the pendulum decreases, its average deviation from the vertical increases. As $\omega_{3}\rightarrow 0$ ($\kappa\rightarrow\kappa_{*}$), a stabilization of periodic solutions (nutations) close to the lower vertical rotation $\sigma^{-}$ occurs. When $\omega_{3}=0$, a stabilization of the lower equilibrium point of the pendulum occurs.

In the example considered above, the loss of stability of the zero solution occurs only via a pitch-fork bifurcation ($\varrho_{1}=+1$). For other parameters of the spherical robot and oscillations of the supporting plane, the loss of stability can occur also via a period-doubling bifurcation. Emerging stable solutions can also lose stability via a similar bifurcation. An interesting open question is that of existence of a cascade of period-doubling bifurcations in this system and of emergence of chaotic attracting invariant sets.

5 CONCLUSION

In this paper, we have explored the possibility of stabilizing the upper vertical rotations of the pendulum of the spherical robot on a vibrating supporting plane by using feedback.

We implement feedback by means of a control torque which depends only on the angular velocity of the pendulum and the orientation of its symmetry axis. In this case, similarly to the free system, the equations governing the dynamics of a pendulum decouple from the complete system of equations. For this reason, the stability of the equilibrium points of the pendulum does not depend on the rotation of the shell.

The feedback chosen in this paper artificially generates in the system the viscous resistance torque, which decreases the horizontal component of the angular velocity of the pendulum. Its “dissipative” nature leads to asymptotic behavior near vertical rotations of the pendulum.

We have shown that, in the linear approximation, the dynamics of the system near vertical rotations is described by the damped Mathieu equation. For this equation, we have numerically constructed and analyzed the regions of stability and instability on the parameter plane. We have shown that, as the plane oscillates in a prescribed manner with an amplitude and a frequency corresponding to the region of rotations without loss of contact with the plane, the vertical rotations of the pendulum of the spherical robot on the vibrating surface can be stabilized by changing the angular velocity of rotation of the pendulum. The regions of stability expand (compared to a system without feedback), and the stability becomes asymptotic.

Also, we have numerically analyzed the stability of periodic solutions arising when the vertical rotations of the pendulum lose stability.

A question that remains open is that of finding and analyzing the domains of attraction in the phase space of vertical rotations of the pendulum of the spherical robot and other solutions which arise when the zero solution becomes unstable. In the future, it would also be interesting to experimentally verify the stabilization (destabilization) of the equilibrium points of the spherical robot as the angular velocity of rotation of the pendulum is varied with and without feedback.