# Kinematics of Spherical Robots Rolling over 3D Terrains

Authors: Saeed Moazami; Hassan Zargarzadeh; Srinivas Palanki
Journal: Complexity
 (2019)

---

## Abstract

Although the kinematics and dynamics of spherical robots (SRs) on flat horizontal and inclined 2D surfaces are thoroughly investigated, their rolling behavior on generic 3D terrains has remained unexplored. This paper derives the kinematics equations of the most common SR configurations rolling over 3D surfaces. First, the kinematics equations for a geometrical sphere rolling over a 3D surface are derived along with the characterization of the modeling method. Next, a brief review of current mechanical configurations of SRs is presented as well as a novel classification for SRs based on their kinematics. Then, considering the mechanical constraints of each category, the kinematics equations for each group of SRs are derived. Afterward, a path-tracking method is utilized for a desired 3D trajectory. Finally, simulations are carried out to validate the developed models and the effectiveness of the proposed control scheme.


---

## Body

## 1. Introduction

Spherical robots are a class of mobile robots that are generally recognized by their ball-shaped shell and internal driving components that provide torques required for their rolling motion. From the ball-shaped exterior, SRs inherit multiple advantages over other types of mobile robots, in which skidding, tipping over, falling, or friction with the surface makes them vulnerable or inefficient [1, 2]. Despite all the unique features, the complicated nonlinear behavior of SRs has remained as a hurdle to fully comprehend their dynamics, motion kinematics, and unveil their maneuverability capabilities.Although the earliest efforts to analytically capture the kinematics and dynamics of rolling geometrical balls on mathematical surfaces back to more than a century ago, e.g., in E. Routh and S. A. Chaplygin’s works, this topic is still an open discussion that is being investigated in more recent papers such as [3–5]. Numerous research studies have been conducted on the mathematical modeling of kinematics and dynamics of SRs with a variety of mechanical configurations [6, 7]. Among them, a widespread assumption is that SRs roll over an ideally flat horizontal plane. In a number of available works, where SRs are studied on nonhorizontal surfaces, simplifying assumptions are made. For instance, in [2, 8], two SR designs are investigated that climb obstacles, assuming the condition to be static. From a different view, rolling of SRs is studied where the desired path is assumed to be a straight line with constant slope or a single step obstacle [9–11] and a 2D curved path with variable slope [12] respectively. In [13], the authors have investigated the dynamics of Martian tumbleweed rovers while this special type of SR rolls in its heading direction and the turning action is not considered for them. In fact, while several pieces of research have been done on 3D kinematics of other types of mobile robots such as legged [14] and wheeled robots [15, 16], to the best of the authors’ knowledge, the general problem of kinematics of SRs rolling on 3D terrains has not been investigated in the literature. The motivation to address this problem is that, while many applications of the SRs are on flat surfaces such as indoor [17] and paved roads [18], for a variety of applications such as agriculture [19], surveillance [20], environmental monitoring [21], and even planetary explorations [22], they would get exposed to uneven terrains.In this work, prior to deriving the kinematics of SRs on 3D terrains, a general method for modeling a geometrical sphere rolling over a mathematically known 3D surface is developed. Then, the derived equations are expanded in order to be applied to SRs considering their specifications. Concretely, a variety of mechanisms are utilized in SRs to provide the required propelling torques and forces for their rolling action. Each configuration imposes its own kinematical constraints on the SRs’ rolling motion. Therefore, to study the kinematics of SRs, it is essential to classify current and feasible designs of SRs accordingly.There are a few SR classifications available in the literature. In survey [6], SRs are classified based on their mechanical driving principles as follows: (1) barycenter offset (BCO), (2) conservation of angular momentum (COAM), and (3) shell transformation (OST). In another review, SRs are classified based on their mechanical configurations [23], e.g., single wheel, hamster wheel, pendulum driven, gimbal mechanism, single ball, mass movement, and a set of designs that use flywheels. Again, each class of the above robots can have different kinematics. Therefore, to investigate the rolling behavior of SRs based on the constraints imposed by their kinematical configurations, presenting a new, although brief, classification seems inevitable.The classification proposed in this work divides SRs into two major categories: (1) Continuous Rolling Spherical Robots (CR-SR), with triple axes rolling (3R-SR), dual axes rolling (2R-SR), and rolling and turning (RT-SR) spherical robots as subcategories and (2) Rolling and Steering Spherical Robots (RS-SR). Based on the constraints of each different category, the kinematics equations for SRs rolling over 3D terrains are derived.Afterward, a modified pure pursuit method [24] is utilized for the path-tracking control problem of the proposed models. The path-tracking method indicates the required kinematic actuation states that should be provided in each SR category to track a desired 3D trajectory. Finally, simulations are carried out to verify the kinematics model’s accuracy and its controller’s tracking efficiency in different models along with the comparison between their behavior.The remainder of this paper is organized as follows: Section2 is devoted to the model description. Section 3 elaborates on the problem of kinematics modeling of a sphere rolling over a 3D surface. In Section 4, the objectives of the paper are addressed, including presenting a new classification for SRs by their kinematics properties and then deriving kinematics equations of SRs based on specifications of each selected SR category. In Section 5, a 3D path tracking kinematic control scheme is presented. The simulation results in MATLAB are represented in Section 6, followed by the conclusion and future work.


## 2. Model Description

Consider a geometrical sphere rolling over a 3D surfaceS as shown in Figure 1. To derive the kinematics equations of the rolling motion, the first step is to define the required reference frames, as ℱW=OW,XW−,YW−,ZW−, ℱTr=OTr,xTr−,yTr−,zTr−, ℱT0=OT0,xT0−,yT0−,zT0−, ℱT=OT,xT−,yT−,zT−, and ℱL=OL,xL−,yL−,zL−, defined as world, translated, surface-side tangent, robot-side tangent, and local reference frames, respectively. In each reference frame set, the first element shows its origin and the rest are the reference frame’s axes. Accordingly, as it is illustrated in Figure 1, I^W,J^W,K^W, i^Tr,j^Tr,k^Tr, i^T0,j^T0,k^T0, i^T,j^T,k^T, and h^,l^,n^ represent unit vectors of ℱW, ℱTr, ℱT0, ℱT, and ℱL respectively.Figure 1
Schematic diagram of reference frames and kinematic parameters that are used to derive kinematics equation of the rolling sphere on a 3D surfaceS.In the sequel, vector quantities are defined inℱW; if otherwise, their reference frames are specified in left superscript, e.g., V˗ℱA. Vectors are represented by an overhead right arrow while unit vectors are accented by a caret-shaped hat symbol. Homogeneous transformations, e.g., from arbitrary reference frames ℱA to ℱB, consisting of a rotation matrix of ℛAB3×3 and a translation vector ℳAB3×1, are defined in the form of [25](1)TAB=ℛABℳAB01×31.Through utilizing transformation matrices, a fourth element of0 is added to velocity vectors as translation does not apply to them. A fourth element of 1 in position vectors guarantees both rotation and translation.In Figure1, the position of the contact point of the sphere with the 3D surface S is represented as P0:(2)P0=ℱWP0=x0t,y0t,z0tT,which instantaneously indicates the position coordinates of OTr, OT0, OT, and OL. The imaginary plane tangent to S passing through P0 is called the tangent plane denoted as TS. As it is shown in Figure 1, normal axis of the sphere is defined to be perpendicular to TS coinciding with the center of the sphere OSp, while lateral and longitudinal axes are defined to be parallel to TS; i.e., they do not rotate through rotation of the sphere about axes parallel to TS, but they rotate accordingly when the sphere rotates about the normal axis. ℱTr is resulted from the translation of ℱW to P0 such that xTr−, yTr−, and zTr− remain parallel to XW−, YW−, and ZW−. Obviously, axes of ℱTr do not necessarily lie on TS. In ℱT0, zT0− is along the normal vector of S, nS˗. xT0− is locally tangent to the surface and therefore lies on TS. yT0− is mutually perpendicular to xT0− and zT0− following the right-hand rule. Furthermore, as shown in Figure 1, ℱT is resulted from rotation of ℱT0 about xT0− by 180 degrees. Consequently, xT− coincides with xT0−. yT− and zT− are in the opposite direction to yT0− and zT0− with respect to OT0. In ℱL, xL− and yL− are also tangent to TS. xL− is along the sphere’s longitudinal direction. zL− is normal to TS coinciding with zT− and normal axis of the sphere, and yL− is mutually perpendicular to xL− and zL− following the right-hand rule. One can note that yL− can be considered as the projection of the sphere’s lateral axis on TS.The next step is to define rotation angles of the sphere. According to Figure1, θ is the sphere’s rolling angle, i.e., the rotation angle of the sphere about its lateral axis. ϕ is the sphere’s tilting angle, i.e., the rotation angle of the sphere about its longitudinal axis, and ψ is the sphere’s turning angle about the robot-side normal axis, n^. It can be noted that ℱL can be considered as ℱT rotated through an angle of ψ about zT−. Finally, S is a 3D surface defined in ℱW as S:Fx,y,z=0, in the form of(3)fx,y−z=0,implying that for any x and y, the elevation of the surface S along the vertical axis of ZW− is z=fx,y.Assumptions: the mathematical model is derived based on the following assumptions:(1)
the sphere rolls over a 3D surface without slipping, skidding, or falling; i.e., it does not lose traction with the surface.(2)
The surfaceS is smooth; i.e., z=fx,y is a continuous function of class C1 or higher [26].(3)
The sphere is rigid, and the radius of curvature of the surface is never smaller than the radius of the sphere. Therefore, the sphere always remains in a single point contact with the surface.


## 3. Kinematics Modeling of Spheres on 3D Surfaces

The next phase to derive the kinematics equations of the sphere rolling overS is to calculate the utilized reference frames’ unit vectors analytically along with deriving transformation matrices between reference frames. In this process, first, we will determine the vector quantities in ℱL, and then, by sequentially utilizing transformation matrices, we can calculate the desired quantities in ℱW.As mentioned in the previous section,ℱTr is parallel to ℱW and OTr coincides with P0; therefore, TTrW performs no rotation, i.e., ℛTrW=I3 and a translation of ℳTrW=P0. Thus, TTrW can be written as follows:(4)TTrW=I3P001,where I3×3 is the identity matrix. From (3), the normal vector of surface S at position P0 is given by(5)nS˗=∇F=fx,fy,−1T,where fx=∂fx,y/∂x and fy=∂fx,y/∂y are the partial derivative functions. Now, as nS˗ coincides with zT0−, its unit vector can be written as(6)k^T0=n^S=nS˗nS˗.Definingsn as the following:(7)sn=nS˗−1=1+fx2+fy2−1/2,we can write(8)n^S=snnS˗Seemingly, the direction ofn^S is perpendicular to the surface and always downward, i.e., Zn^S=−1. As robot-site normal, n^, is in the opposite direction to the surface-side normal and based on assumption (1), there is no relative velocity between ℱT0 and ℱT; therefore, we can write(9)TTT0=10000−10000−100001.In fact,n^ and k^T can be calculated by rotating nS˗ by 180 degrees about i^T0; therefore,(10)n^=k^T=−n^S=sn−fx,−fy,1TNext, we calculate the unit vectors ofℱT axes. Let us consider a sphere that rolls over S in a way that ψ=0, which means only θ and ϕ contribute in the rolling action. Due to this rolling, the direction of sphere’s longitudinal and lateral axes changes accordingly to match the geometry of S. Based on Euler’s rotation theorem, any Cartesian coordinate system in 3D space with a common origin are related by a rotation about a unique axis called Euler axis denoted as e^. Intuitively, this rotation maps axes of ℱT to ℱTr. We know that the rotation maps n^ to k^Tr and it is assumed that ψ=0; therefore, as it is depicted in Figure 2, e^ should lie in the tangent plane and be perpendicular to both n^ and k^Tr, leading to(11)e˗=k^Tr×n^.Figure 2
Schematic diagram of reference frames, Euler axise^, and transformations between reference frames.Consideringk^Tr=0,0,1T and n^ from (10), e˗ is calculated as(12)e˗=snfy,−fx,0T.Consequently,(13)e^=e˗e˗=sfy,−fx,0Twiths=e˗−1=1fx2+fy2.Moreover, the angle of rotation can be calculated as(14)γ=atan2k^Tr×n^,k^Tr⋅n^,where atan2 is the four-quadrant inverse tangent:(15)k^Tr×n^=sinγ=sns−1,k^Tr⋅n^=cosγ=sn.Therefore, to avoid singularity, from (14) and (15), γ is written in the following conditional statement form:(16)γ=atan2s−1,1,ifk^Tr≠n^,0,ifk^Tr=n^.Now, to calculatei^T, we should apply the same rotation to i^Tr=1,0,0T. There are several algorithms to transform a vector in space, given a rotation axis and an angle of rotation such as the Rodrigues’ rotation method that is explained in the following.According to Rodrigues’ rotation formula [27], as i^Tr is a vector in ℝ3 and e^=ex,ey,0T is a unit vector representing the axis of rotation about which i^T rotates by an angle of γ resulting in i^Tr, it can be shown that(17)i^T=i^Trcosγ+e^×i^Trsinγ+e^e^⋅i^Tr1−cosγ.Equation (17) can be written in the matrix form to calculate TTTr using Rodrigues’ rotation formula. Let E× denote cross-product matrix for e^; i.e., a matrix when multiplied from the left with a vector gives the cross product. It can be shown that(18)E×=0−ezeyez0−ex−eyex0=00−sfx00−sfysfxsfy0.Now, it can be shown that the rotation matrixℛTTr can be written in the following form:(19)ℛTTr=I3+sinγE×+1−cosγE×2,which, using (18), can be written as the following:(20)ℛTTr=1−s2fx21−Cγ−s2fxfy1−Cγ−sfxSγ−s2fxfy1−Cγ1−s2fy21−Cγ−sfySγsfxSγsfySγ1−s2fx2+fy21−Cγ,where, for the sake of brevity, sin⋅ and cos⋅ are abbreviated to S⋅ and C⋅, respectively. Considering Sγ=sns−1 from (15), we can write:(21)sfxSγ=snfx,sfySγ=snfy.Also, from (13) and (15), it can be concluded that(22)1−s2fx2+fy21−Cγ=sn.AsℳTTr=03×1, and by substituting (21) and (22) in (20), we can write:(23)TT0Tr=ℛTTr001=1−s2fx21−sn−s2fxfy1−sn−snfx0−s2fxfy1−sn1−s2fy21−sn−snfy0snfxsnfysn00001.Alternatively, to calculateTTTr, rotation quaternions can be utilized which is explained in detail in Appendix. Depending on the application, ℛTTr can be calculated from (19) and (23), or through quaternions. To calculate j^T, we can write:(24)j^T=n^×i^TNext, asℱL is resulted from the rotation of ℱT through ψ about zT−, then ℳLT=03×1 and TLT can be written as follows:(25)TLT=CψSψ00−SψCψ0000100001.Finally, from (4), (23), and (25), it can be concluded that(26)TLW=TTrWTTTrTLT=a11a12a13a14a21a22a23a24a31a32a33a34a41a42a43a44,where the elements of TLW are given as(27)a11=1−s2fx21−snCψ+s2fxfy1−snSψ,a12=1−s2fx21−snSψ−s2fxfy1−snCψ,a13=−snfx,a14=x0,a21=−s2fxfy1−snCψ−1−s2fy21−snSψ,a22=−s2fxfy1−snSψ+1−s2fy21−snCψ,a23=−snfy,a24=y0,a31=snfxCψ−snfySψ,a32=snfxSψ+snfyCψ,a33=sn,a34=z0,a41=a42=a43=0,a44=1.Now, consider the sphere that instantaneously rolls over the tangent plane and letV˗ℱL=u,v,w,0T denote the velocity of the sphere in ℱL; then, we can write:(28)V˗ℱW=TLWV˗ℱL,where V˗ℱW=x˙,y˙,z˙,0T. Consequently, the position of the sphere can be calculated by integrating (28) over time. To calculate V˗ℱL, it can be written as a set of functions of kinematics parameters such that(29)u=g1θ,ϕ,ψ,θ˙,ϕ˙,ψ˙,v=g2θ,ϕ,ψ,θ˙,ϕ˙,ψ˙,w=0.Obviously,g1 and g2 functions vary depending on the configuration of the SR. Transforming V˗ℱL to ℱW by plugging (29) into (28) and using terms presented in (27), we have(30)x˙y˙z˙0T=TLWg1g200T.resulting in(31)x˙=a11g1+a12g2,y˙=a21g1+a22g2,z˙=a31g1+a32g2.The next section, first, presents a classification based on kinematics specifications of SRs and then provides calculation steps ofV˗ℱL for each class of SR.


## 4. Kinematics Modeling of Spherical Robots Rolling over 3D Surfaces

This section is dedicated to deriving the kinematics of SRs rolling over 3D surfaces. First, we present a brief classification in which the most popular configurations of SRs are divided into two main categories and a number of subcategories. Secondly, by applying constraints of each category (or type), the kinematics of SRs on 3D surfaces is derived.


### 4.1. Continuous Rolling Spherical Robots (CR-SRs)

CR-SRs are capable of performing rolling maneuvers continuously; i.e., when the robot rolls in any direction, there is no limitation for its rotation angle. However, based on their mechanical configurations, these robots can roll in different directions based on which this category can be divided into three subcategories explained in the following.


#### 4.1.1. Triple-Axis Rolling Spherical Robot (3R-SR)

This subcategory is the most general type of CR-SR. As illustrated in Figure3, in 3R-SRs the internal components are mechanically capable of providing torques required for the robot to rotate about all of its three main axes, namely, angular velocities of θ˙, ϕ˙, and ψ˙. Therefore, 3R-SRs are omnidirectional. Several mechanical configurations are available in the literature that relies on this actuation method [28–30].Figure 3
Schematic diagram of rotation directions of continuous rolling spherical robots.Figure4 illustrates differential rolling lengths for a 3R-SR in an infinitesimally small time dt along h^ and l^ directions. Based on the assumption (1) and using the equal-arc-length rule for rolling without slipping [31], we have(32)dxL=Rdθ,dyL=−Rdϕ,where R denotes the sphere’s radius. Integrating (32) over time, V˗ℱL can be calculated as the following:(33)V˗ℱL=Rθ˙,−Rϕ˙,0T.Figure 4
Differential rolling length of a spherical robot using triple-axis rolling during differential of time,dt, along h^ and l^ directions.Now, from (28), (31), and (33), we can calculate elements of V˗ℱW as follows:(34)x˙=a11Rθ˙−a12Rϕ˙,y˙=a21Rθ˙−a22Rϕ˙,z˙=a31Rθ˙−a32Rϕ˙.Substituting terms presented in (27) in (34), and considering ψ˙ as the angular velocity of ℱL in accordance with ℱT, elements of V˗ℱW for a 3R-SR are calculated as the following:(35)x˙=R1−s2fx21−snCψ+s2fxfy1−snSψθ˙−R1−s2fx21−snSψ−s2fxfy1−snCψϕ˙,y˙=−Rs2fxfy1−snCψ+1−s2fy21−snSψθ˙+Rs2fxfy1−snSψ−1−s2fy21−snCψϕ˙,z˙=snRfxCψ−fySψθ˙−snRfxSψ+fyCψϕ˙.


#### 4.1.2. Dual-Axis Rolling Spherical Robot (2R-SR)

In this subcategory of SRs, the robot provides torques and, consequently, the ability to roll about its two main axes which are parallel toTS. There are a number of SR designs that relate to this type [23, 32–34] such as hamster wheel SRs. Comparing to 3R-SRs, similar steps should be taken in order to derive kinematics equation of 2R-SRs. The only distinction is that, in 2R-SRs, the robot is not capable of rotating about its vertical axis; i.e., it is assumed that ψ˙=0; therefore, we can write(36)x˙=a11Rθ˙−a12Rϕ˙,y˙=a21Rθ˙−a22Rϕ˙,z˙=a31Rθ˙−a32Rϕ˙,ψ˙=0.The termψ˙=0 renders 2R-SRs as robots with completely different behavior. Assuming ψ˙=ψ˙0=0, for 2R-SRs, we can write(37)x˙=R1−s2fx21−snθ˙+Rs2fxfy1−snϕ˙,y˙=−Rs2fxfy1−snθ˙−R1−s2fy21−snϕ˙,z˙=snRfxθ˙−snRfyϕ˙.


#### 4.1.3. Rolling and Turning Spherical Robot (RT-SR)

In addition to rolling action, RT-SRs are able to turn about their vertical axis in order to change their moving direction. In other words, internal parts can provide angular velocities ofθ˙ and ψ˙, while it is assumed that ϕ˙=0. RT-SRs are relatively less common but still some designs can be found in the literature [35–38]. To calculate kinematics of RT-SRs, one can write(38)dxL=Rdθ,dyL=0,x˙=a11Rθ˙,y˙=a21Rθ˙,z˙=a31Rθ˙.Angular velocity ofℱL in accordance with ℱT equals to ψ˙. Substituting terms from (27) in (38), we have(39)x˙=R1−s2fx21−snCψ+s2fxfy1−snSψθ˙,y˙=−Rs2fxfy1−snCψ+1−s2fy21−snSψθ˙,z˙=snRfxCψ−fySψθ˙.


### 4.2. Rolling and Steering Spherical Robots (RS-SRs)

As it is shown in Figure5, consider α˙ as angular velocity of a RS-SR around its transverse axis that provides forward and backward rolling motion of the robot. In addition, the robot can tilt about its longitudinal axis as its steering mechanism provides turning ability for the robot through an arc shape rolling motion. Unlike CR-SRs, in RS-SRs, tilting angle ϕ is limited. Another difference is that in RS-SRs, the transverse axis is defined in a way that tilts with the robot, so it is not always parallel to the tangent plane. This type of SRs uses an underactuated mechanical propulsion mechanism that allows the robot to roll in a nonholonomic manner [39–46].Figure 5
Schematic diagram of a spherical robot using rolling and steering (RS) action kinematics.Figure5 illustrates a RS-SR with tilting angle of ϕ with respect to the tangent plane. The turning action of the sphere can be modeled as a cone with an apex angle of 2ϕ, purely rolling over the surface S, for which the instantaneous axis of rotation is the imaginary line passing through the cone’s vertex V and P0 that lies in the tangent plane TS. Considering the instantaneous differential rolling arc of the cone’s base circle, we have(40)rcdα=ρdψ,where rc is the cone’s base circle radius and dα and dψ are the differential rotation of the sphere during differential of time, dt, about the transverse and vertical axes respectively. ρ is the instantaneous radius of curvature of the turning motion that is the distance between the vertex of the cone and the contact point. These two quantities can be calculated as(41)rc=Rcosϕand(42)ρ=−Rcotϕ.Substituting (41) and (42) into (40) results in(43)dψ=−dαsinϕ.Dividing both sides of (43) by dt, we have(44)ψ˙=−α˙sinϕ.Moving on, the position vector of the center point of the sphere inℱL in accordance with P0 can be written as(45)r˗P0OSpℱL=Rn^,and angular velocity of the sphere in ℱL can be written as(46)Ω˗SpℱL=ϕ˙h^+α˙cosϕl^−α˙sinϕn^.Then, the linear velocity of the sphere’s center point can be calculated as follows:(47)V˗SpℱL=Ω˗SpℱL×r˗P0OSpℱL=Rα˙cosϕh^−Rϕ˙l^.Forthwith, we are able to formV˗ℱL as the following:(48)V˗ℱL=Rα˙cosϕ,−Rϕ˙,0T.Eventually, from (31) and (48), we can calculate components of V˗ℱW as follows:(49)x˙=a11Rα˙cosϕ−a12Rϕ˙,y˙=a21Rα˙cosϕ−a22Rϕ˙,z˙=a31Rα˙cosϕ−a32Rϕ˙.Equations in (49), using terms presented in (27), can be utilized to calculate the position of the SR, P0, based on input parameters of robot, α˙ and ϕ, as the following:(50)x˙=R1−s2fx21−snCψ+s2fxfy1−snSψα˙Cϕ−R1−s2fx21−snSψ−s2fxfy1−snCψϕ˙,y˙=−Rs2fxfy1−snCψ+1−s2fy21−snSψα˙Cϕ+Rs2fxfy1−snSψ−1−s2fy21−snCψϕ˙,z˙=snRfxCψ−fySψα˙Cϕ−snRfxSψ+fyCψϕ˙.In all types of SR, to calculate the position of the center point of the sphere denoted asOSp=XSp,YSp,ZSpT, we have(51)OSp=P0+Rn^.Alternatively, utilizing (26) and using OSpℱL=0,0,R,1T, it can be calculated as(52)OSp=TLWOSpℱL=a13R+a14a23R+a24a33R+a341.Observing terms in (27), one can notice that P0=a14,a24,a34T and n^=a13,a23,a33T; therefore, (52) is analogous to (51).


## 4.1. Continuous Rolling Spherical Robots (CR-SRs)

CR-SRs are capable of performing rolling maneuvers continuously; i.e., when the robot rolls in any direction, there is no limitation for its rotation angle. However, based on their mechanical configurations, these robots can roll in different directions based on which this category can be divided into three subcategories explained in the following.


### 4.1.1. Triple-Axis Rolling Spherical Robot (3R-SR)

This subcategory is the most general type of CR-SR. As illustrated in Figure3, in 3R-SRs the internal components are mechanically capable of providing torques required for the robot to rotate about all of its three main axes, namely, angular velocities of θ˙, ϕ˙, and ψ˙. Therefore, 3R-SRs are omnidirectional. Several mechanical configurations are available in the literature that relies on this actuation method [28–30].Figure 3
Schematic diagram of rotation directions of continuous rolling spherical robots.Figure4 illustrates differential rolling lengths for a 3R-SR in an infinitesimally small time dt along h^ and l^ directions. Based on the assumption (1) and using the equal-arc-length rule for rolling without slipping [31], we have(32)dxL=Rdθ,dyL=−Rdϕ,where R denotes the sphere’s radius. Integrating (32) over time, V˗ℱL can be calculated as the following:(33)V˗ℱL=Rθ˙,−Rϕ˙,0T.Figure 4
Differential rolling length of a spherical robot using triple-axis rolling during differential of time,dt, along h^ and l^ directions.Now, from (28), (31), and (33), we can calculate elements of V˗ℱW as follows:(34)x˙=a11Rθ˙−a12Rϕ˙,y˙=a21Rθ˙−a22Rϕ˙,z˙=a31Rθ˙−a32Rϕ˙.Substituting terms presented in (27) in (34), and considering ψ˙ as the angular velocity of ℱL in accordance with ℱT, elements of V˗ℱW for a 3R-SR are calculated as the following:(35)x˙=R1−s2fx21−snCψ+s2fxfy1−snSψθ˙−R1−s2fx21−snSψ−s2fxfy1−snCψϕ˙,y˙=−Rs2fxfy1−snCψ+1−s2fy21−snSψθ˙+Rs2fxfy1−snSψ−1−s2fy21−snCψϕ˙,z˙=snRfxCψ−fySψθ˙−snRfxSψ+fyCψϕ˙.


### 4.1.2. Dual-Axis Rolling Spherical Robot (2R-SR)

In this subcategory of SRs, the robot provides torques and, consequently, the ability to roll about its two main axes which are parallel toTS. There are a number of SR designs that relate to this type [23, 32–34] such as hamster wheel SRs. Comparing to 3R-SRs, similar steps should be taken in order to derive kinematics equation of 2R-SRs. The only distinction is that, in 2R-SRs, the robot is not capable of rotating about its vertical axis; i.e., it is assumed that ψ˙=0; therefore, we can write(36)x˙=a11Rθ˙−a12Rϕ˙,y˙=a21Rθ˙−a22Rϕ˙,z˙=a31Rθ˙−a32Rϕ˙,ψ˙=0.The termψ˙=0 renders 2R-SRs as robots with completely different behavior. Assuming ψ˙=ψ˙0=0, for 2R-SRs, we can write(37)x˙=R1−s2fx21−snθ˙+Rs2fxfy1−snϕ˙,y˙=−Rs2fxfy1−snθ˙−R1−s2fy21−snϕ˙,z˙=snRfxθ˙−snRfyϕ˙.


### 4.1.3. Rolling and Turning Spherical Robot (RT-SR)

In addition to rolling action, RT-SRs are able to turn about their vertical axis in order to change their moving direction. In other words, internal parts can provide angular velocities ofθ˙ and ψ˙, while it is assumed that ϕ˙=0. RT-SRs are relatively less common but still some designs can be found in the literature [35–38]. To calculate kinematics of RT-SRs, one can write(38)dxL=Rdθ,dyL=0,x˙=a11Rθ˙,y˙=a21Rθ˙,z˙=a31Rθ˙.Angular velocity ofℱL in accordance with ℱT equals to ψ˙. Substituting terms from (27) in (38), we have(39)x˙=R1−s2fx21−snCψ+s2fxfy1−snSψθ˙,y˙=−Rs2fxfy1−snCψ+1−s2fy21−snSψθ˙,z˙=snRfxCψ−fySψθ˙.


## 4.1.1. Triple-Axis Rolling Spherical Robot (3R-SR)

This subcategory is the most general type of CR-SR. As illustrated in Figure3, in 3R-SRs the internal components are mechanically capable of providing torques required for the robot to rotate about all of its three main axes, namely, angular velocities of θ˙, ϕ˙, and ψ˙. Therefore, 3R-SRs are omnidirectional. Several mechanical configurations are available in the literature that relies on this actuation method [28–30].Figure 3
Schematic diagram of rotation directions of continuous rolling spherical robots.Figure4 illustrates differential rolling lengths for a 3R-SR in an infinitesimally small time dt along h^ and l^ directions. Based on the assumption (1) and using the equal-arc-length rule for rolling without slipping [31], we have(32)dxL=Rdθ,dyL=−Rdϕ,where R denotes the sphere’s radius. Integrating (32) over time, V˗ℱL can be calculated as the following:(33)V˗ℱL=Rθ˙,−Rϕ˙,0T.Figure 4
Differential rolling length of a spherical robot using triple-axis rolling during differential of time,dt, along h^ and l^ directions.Now, from (28), (31), and (33), we can calculate elements of V˗ℱW as follows:(34)x˙=a11Rθ˙−a12Rϕ˙,y˙=a21Rθ˙−a22Rϕ˙,z˙=a31Rθ˙−a32Rϕ˙.Substituting terms presented in (27) in (34), and considering ψ˙ as the angular velocity of ℱL in accordance with ℱT, elements of V˗ℱW for a 3R-SR are calculated as the following:(35)x˙=R1−s2fx21−snCψ+s2fxfy1−snSψθ˙−R1−s2fx21−snSψ−s2fxfy1−snCψϕ˙,y˙=−Rs2fxfy1−snCψ+1−s2fy21−snSψθ˙+Rs2fxfy1−snSψ−1−s2fy21−snCψϕ˙,z˙=snRfxCψ−fySψθ˙−snRfxSψ+fyCψϕ˙.


## 4.1.2. Dual-Axis Rolling Spherical Robot (2R-SR)

In this subcategory of SRs, the robot provides torques and, consequently, the ability to roll about its two main axes which are parallel toTS. There are a number of SR designs that relate to this type [23, 32–34] such as hamster wheel SRs. Comparing to 3R-SRs, similar steps should be taken in order to derive kinematics equation of 2R-SRs. The only distinction is that, in 2R-SRs, the robot is not capable of rotating about its vertical axis; i.e., it is assumed that ψ˙=0; therefore, we can write(36)x˙=a11Rθ˙−a12Rϕ˙,y˙=a21Rθ˙−a22Rϕ˙,z˙=a31Rθ˙−a32Rϕ˙,ψ˙=0.The termψ˙=0 renders 2R-SRs as robots with completely different behavior. Assuming ψ˙=ψ˙0=0, for 2R-SRs, we can write(37)x˙=R1−s2fx21−snθ˙+Rs2fxfy1−snϕ˙,y˙=−Rs2fxfy1−snθ˙−R1−s2fy21−snϕ˙,z˙=snRfxθ˙−snRfyϕ˙.


## 4.1.3. Rolling and Turning Spherical Robot (RT-SR)

In addition to rolling action, RT-SRs are able to turn about their vertical axis in order to change their moving direction. In other words, internal parts can provide angular velocities ofθ˙ and ψ˙, while it is assumed that ϕ˙=0. RT-SRs are relatively less common but still some designs can be found in the literature [35–38]. To calculate kinematics of RT-SRs, one can write(38)dxL=Rdθ,dyL=0,x˙=a11Rθ˙,y˙=a21Rθ˙,z˙=a31Rθ˙.Angular velocity ofℱL in accordance with ℱT equals to ψ˙. Substituting terms from (27) in (38), we have(39)x˙=R1−s2fx21−snCψ+s2fxfy1−snSψθ˙,y˙=−Rs2fxfy1−snCψ+1−s2fy21−snSψθ˙,z˙=snRfxCψ−fySψθ˙.


## 4.2. Rolling and Steering Spherical Robots (RS-SRs)

As it is shown in Figure5, consider α˙ as angular velocity of a RS-SR around its transverse axis that provides forward and backward rolling motion of the robot. In addition, the robot can tilt about its longitudinal axis as its steering mechanism provides turning ability for the robot through an arc shape rolling motion. Unlike CR-SRs, in RS-SRs, tilting angle ϕ is limited. Another difference is that in RS-SRs, the transverse axis is defined in a way that tilts with the robot, so it is not always parallel to the tangent plane. This type of SRs uses an underactuated mechanical propulsion mechanism that allows the robot to roll in a nonholonomic manner [39–46].Figure 5
Schematic diagram of a spherical robot using rolling and steering (RS) action kinematics.Figure5 illustrates a RS-SR with tilting angle of ϕ with respect to the tangent plane. The turning action of the sphere can be modeled as a cone with an apex angle of 2ϕ, purely rolling over the surface S, for which the instantaneous axis of rotation is the imaginary line passing through the cone’s vertex V and P0 that lies in the tangent plane TS. Considering the instantaneous differential rolling arc of the cone’s base circle, we have(40)rcdα=ρdψ,where rc is the cone’s base circle radius and dα and dψ are the differential rotation of the sphere during differential of time, dt, about the transverse and vertical axes respectively. ρ is the instantaneous radius of curvature of the turning motion that is the distance between the vertex of the cone and the contact point. These two quantities can be calculated as(41)rc=Rcosϕand(42)ρ=−Rcotϕ.Substituting (41) and (42) into (40) results in(43)dψ=−dαsinϕ.Dividing both sides of (43) by dt, we have(44)ψ˙=−α˙sinϕ.Moving on, the position vector of the center point of the sphere inℱL in accordance with P0 can be written as(45)r˗P0OSpℱL=Rn^,and angular velocity of the sphere in ℱL can be written as(46)Ω˗SpℱL=ϕ˙h^+α˙cosϕl^−α˙sinϕn^.Then, the linear velocity of the sphere’s center point can be calculated as follows:(47)V˗SpℱL=Ω˗SpℱL×r˗P0OSpℱL=Rα˙cosϕh^−Rϕ˙l^.Forthwith, we are able to formV˗ℱL as the following:(48)V˗ℱL=Rα˙cosϕ,−Rϕ˙,0T.Eventually, from (31) and (48), we can calculate components of V˗ℱW as follows:(49)x˙=a11Rα˙cosϕ−a12Rϕ˙,y˙=a21Rα˙cosϕ−a22Rϕ˙,z˙=a31Rα˙cosϕ−a32Rϕ˙.Equations in (49), using terms presented in (27), can be utilized to calculate the position of the SR, P0, based on input parameters of robot, α˙ and ϕ, as the following:(50)x˙=R1−s2fx21−snCψ+s2fxfy1−snSψα˙Cϕ−R1−s2fx21−snSψ−s2fxfy1−snCψϕ˙,y˙=−Rs2fxfy1−snCψ+1−s2fy21−snSψα˙Cϕ+Rs2fxfy1−snSψ−1−s2fy21−snCψϕ˙,z˙=snRfxCψ−fySψα˙Cϕ−snRfxSψ+fyCψϕ˙.In all types of SR, to calculate the position of the center point of the sphere denoted asOSp=XSp,YSp,ZSpT, we have(51)OSp=P0+Rn^.Alternatively, utilizing (26) and using OSpℱL=0,0,R,1T, it can be calculated as(52)OSp=TLWOSpℱL=a13R+a14a23R+a24a33R+a341.Observing terms in (27), one can notice that P0=a14,a24,a34T and n^=a13,a23,a33T; therefore, (52) is analogous to (51).


## 5. Control of Spherical Robots over 3D Surfaces

In this section, the pure pursuit method is utilized to address the path tracking problem of the four different types of spherical robots on a 3D surface. ConsiderP:xdt,ydt,zdtT as a 3D desired trajectory lying on S, which is assumed to be smooth; i.e., x˙d,y˙d,z˙dT is bounded for t>0. It is desired that the location of P0=x0t,y0t,z0tT converges to P with a stable error. et is the tracking error defined as(53)et=xdt,ydt,zdtT−P0.Additionally, defining deviation angleζ as the angle between the heading direction, h^, and et, we have(54)ζ=∠et,h^.The objective of this algorithm is to converge the tracking erroret along with ζ to zero for all four proposed categories of SRs. To that end, in each category, the path-planning algorithm computes values for actuated kinematics states, namely, θ˙, ϕ˙, α˙, and ψ˙.


### 5.1. 3R-SR Kinematics Control

As 3R-SRs are capable of providingθ˙, ϕ˙, and ψ˙, the controller is designed to provide the required values accordingly as follows:(55)θ˙=kθetke+etcosζ+1Rx˙d,y˙d,z˙dT,ϕ˙=−kϕ1etke+etsinζ+kϕ2sinζ,ψ˙=kψζ.In the design ofθ˙, et/ke+et is used to normalize the error gain to provide a smooth transition for control efforts from longer distances to the vicinity of the target point; thus, to compensate the magnitude of et, the first term of θ˙ uses cosζ multiplied with a normalized error gain. cosζ is positive or negative when the SR is behind or ahead of the target, respectively. Additionally, one can easily recognize that x˙d,y˙d,z˙dT/R is the angular velocity that is required for the sphere to pursue the target based on its speed. Therefore, this term is used as the bias value for θ˙. For ϕ˙, if ζ≫0, i.e., the robot is not heading towards the target point; the robot tries to compensate its lateral distance to the target by rolling laterally using ϕ˙ in parallel to other actuators. If the robot is heading towards the target ζ≈0, the robot relies mostly on θ˙ to approach the target. Finally, ψ˙ is designed to be directly proportional to ζ as a measure of required turning to compensate the deviation angle.


### 5.2. 2R-SR Kinematics Control

Considering the fact that 2R-SRs are not able to turn, i.e.,ψ˙=0, it can be implied that a 2R-SR is not capable of following the target by changing its heading direction and consequently the robot can only use θ˙ and ϕ˙ to move towards the target. Therefore, the controller is designed such that(56)θ˙=kθ1etke+etcosζ+kθ2cosζ,ϕ˙=−kϕ1etke+etsinζ+kϕ2sinζ,ψ˙=0.


### 5.3. RT-SR Kinematics Control

RT-SRs can move towards the target point by compensating their deviation error,ζ, through turning action and then using θ˙ to approach and follow the target. This method is used to design the controller to provide the required angular velocities as the following:(57)θ˙=kθetke+etcosζ+1Rx˙d,y˙d,z˙dT,ϕ˙=0,ψ˙=kψζ.


### 5.4. RS-SR Kinematics Control

Calculated values forα˙ and ϕ˙ are as follows for RS-SRs:(58)α˙=kαetke+etcosζ+1Rx˙d,y˙d,z˙dT,ϕ˙=kϕζ.It should be noted thatψ˙ is driven and cannot be used as the actuated value in RS-SRs.


## 5.1. 3R-SR Kinematics Control

As 3R-SRs are capable of providingθ˙, ϕ˙, and ψ˙, the controller is designed to provide the required values accordingly as follows:(55)θ˙=kθetke+etcosζ+1Rx˙d,y˙d,z˙dT,ϕ˙=−kϕ1etke+etsinζ+kϕ2sinζ,ψ˙=kψζ.In the design ofθ˙, et/ke+et is used to normalize the error gain to provide a smooth transition for control efforts from longer distances to the vicinity of the target point; thus, to compensate the magnitude of et, the first term of θ˙ uses cosζ multiplied with a normalized error gain. cosζ is positive or negative when the SR is behind or ahead of the target, respectively. Additionally, one can easily recognize that x˙d,y˙d,z˙dT/R is the angular velocity that is required for the sphere to pursue the target based on its speed. Therefore, this term is used as the bias value for θ˙. For ϕ˙, if ζ≫0, i.e., the robot is not heading towards the target point; the robot tries to compensate its lateral distance to the target by rolling laterally using ϕ˙ in parallel to other actuators. If the robot is heading towards the target ζ≈0, the robot relies mostly on θ˙ to approach the target. Finally, ψ˙ is designed to be directly proportional to ζ as a measure of required turning to compensate the deviation angle.


## 5.2. 2R-SR Kinematics Control

Considering the fact that 2R-SRs are not able to turn, i.e.,ψ˙=0, it can be implied that a 2R-SR is not capable of following the target by changing its heading direction and consequently the robot can only use θ˙ and ϕ˙ to move towards the target. Therefore, the controller is designed such that(56)θ˙=kθ1etke+etcosζ+kθ2cosζ,ϕ˙=−kϕ1etke+etsinζ+kϕ2sinζ,ψ˙=0.


## 5.3. RT-SR Kinematics Control

RT-SRs can move towards the target point by compensating their deviation error,ζ, through turning action and then using θ˙ to approach and follow the target. This method is used to design the controller to provide the required angular velocities as the following:(57)θ˙=kθetke+etcosζ+1Rx˙d,y˙d,z˙dT,ϕ˙=0,ψ˙=kψζ.


## 5.4. RS-SR Kinematics Control

Calculated values forα˙ and ϕ˙ are as follows for RS-SRs:(58)α˙=kαetke+etcosζ+1Rx˙d,y˙d,z˙dT,ϕ˙=kϕζ.It should be noted thatψ˙ is driven and cannot be used as the actuated value in RS-SRs.


## 6. Simulation Result

This section provides the simulation results through MATLAB/Simulink to verify the presented method and the controllers that are proposed in Section5. The simulation is carried out, where for the terrain surface, S is defined as the following:(59)z=fx,y=acosωx+cosωy−2.Two desired path types have been selected for the analysis by defining the path’s projection on theXY plane: a smooth path of circular shape with the following mathematical expression:(60)xt=2cos5t100,yt=2cos5t100,and a continuous but nonsmooth diamond-shaped path with sharp vertices and linear edges. The length of the diagonals of the diamond is selected to be four meters which is equal to the diameter of the circular path.Additionally, the following numerical parameters are used:a=0.2m, ω=2, and R=0.2m; and controller design parameters are tuned by trial-and-error.It can be observed from Figures6 through 9 that the SRs can track both desired trajectories successfully, and the tracking error converges to the vicinity of the origin in approximately 30 seconds. Figures 6 and 7 show 3D representations of all types of SRs (CR (3R, 2R, RT) and RS), rolling over the surface defined in (59) for smooth and nonsmooth trajectories. It can be implied from Figure 7 that 3R and 2R-SRs outperform two other types of SR in tracking the nonsmooth trajectory. The reason stems from the fact that RT and RS-SRs rely on adjusting their heading direction for tracking the target that results in lesser agility at sharp vertices.Figure 6
Robots’ actual and desired trajectories on a 3D surface for a smooth circular path. Scan the QR-code to watch an animation of the same simulation result.Figure 7
Robots’ actual and desired trajectories on a 3D surface for a path with sharp vertices.Figure 8
Robots’ actual and desired smooth trajectories alongX, Y, and Z directions versus time.
(a)
(b)
(c)Figure 9
Robots’ actual and desired nonsmooth trajectories alongX, Y, and Z directions versus time.
(a)
(b)
(c)Figures8 and 9 represent the actual and desired values for X, Y, and Z against time for the robots in smooth and nonsmooth trajectories, respectively. It can be seen that the tracking error is bounded due to the fact that x˙d,y˙d,z˙dT is bounded. Since the desired trajectory is unknown to the pure pursuit, the tracking error remains bounded, while the convergence bound can be reduced to the desired level by tuning the controller’s parameters.Figures10 through 12 show the control signals for 2R, RT, and RS-SRs. 3R-SRs rely on all three rotations; therefore, at any situation, all three control signals try to compensate for the errors. Consequently, the control signals for 3R-SRs do not follow a meaningful trend. It can be seen from Figure 12(b) that a limitation of π/3 has been set for the tilting angle of RS-SRs for the simulation. According to the simulation results presented in this section, all types of SRs can effectively track the desired 3D trajectory by utilizing the designed controllers.Figure 10
2R-SR control signals for smooth and nonsmooth trajectories.
(a)
(b)Figure 11
RT-SR control signals for smooth and nonsmooth trajectories.
(a)
(b)Figure 12
(a) RS-SR control signals for smooth and nonsmooth trajectories. (b) RS-SR tilting angles.
(a)
(b)


## 7. Conclusion

This paper concerns the problem of modeling the kinematics and path-tracking control of popular SR configurations rolling over 3D terrains. It is noteworthy that the aim of this paper is not to compare the performance of different types of spherical robots. In fact, it provides a classification of SRs based on their kinematics behavior to derive kinematics equations on 3D surfaces. Utilizing the proposed kinematics, one can use the Euler–Lagrange method to derive the dynamics of robots on uneven terrains by keeping in mind that the gravity direction is not always perpendicular to the tangent plane.One should also note that any SR has its mechanical configurations and consequently different internal kinematics and dynamics. In this regard, reference frames are selected to be parallel to the tangent plane to provide a general form of kinematics equations. Internal kinematic states can be transformed to the utilized reference frames based on the SR’s composition of internal reference frames. The dynamics of an SR defines the constraints and limitations of the robot. Not all kinematic states presented in this paper are feasible for any SR. Here are some examples, though widely used in the literature; Coriolis terms can violate decoupling kinematic assumptions such as simultaneous rolling and turning about two axes in conditions other than quasi-static. Different RS-SRs can have different maximum tilting angles that should be considered based on their mechanical design. The gravity can make rolling in some directions impossible, and the robot cannot always roll with the angular velocity that the controller suggests, and so on. These are the topics that should be addressed in the integration of kinematics and dynamics.



---

*Source: 7543969-2019-12-16.xml*