Forward and Inverse Dynamics of a Unicycle-Like Mobile Robot

In this research work, a new method for solving forward and inverse dynamic problems of mechanical systems having an underactuated structure and subjected to holonomic and/or nonholonomic constraints is developed. The method devised in this paper is based on the combination of the Udwadia-Kalaba Equations with the Underactuation Equivalence Principle. First, an analytical method based on the Udwadia-Kalaba Equations is employed in the paper for handling dynamic and control problems of nonlinear nonholonomic mechanical systems in the same computational framework. Subsequently, the Underactuation Equivalence Principle is used for extending the capabilities of the Udwadia-Kalaba Equations from fully actuated mechanical systems to underactuated mechanical systems. The Underactuation Equivalence Principle represents an efficient method recently developed in the field of classical mechanics. The Underactuation Equivalence Principle is used in this paper for mathematically formalizing the underactuation property of a mechanical system considering a particular set of nonholonomic algebraic constraints defined at the acceleration level. On the other hand, in this study, the Udwadia-Kalaba Equations are analytically reformulated in a mathematical form suitable for treating inverse dynamic problems. By doing so, the Udwadia-Kalaba Equations are employed in conjunction with the Underactuation Equivalence Principle for developing a nonlinear control method based on an inverse dynamic approach. As shown in detail in this investigation, the proposed method can be used for analytically solving in an explicit manner the forward and inverse dynamic problems of several nonholonomic mechanical systems. In particular, the tracking control of the unicycle-like mobile robot is considered in this investigation as a benchmark example. Numerical experiments on the dynamic model of the unicycle-like mobile robot confirm the effectiveness of the nonlinear dynamic and control approaches developed in this work.


Introduction
In this section, background information is provided first.Subsequently, the problem of interest for this investigation is formulated and a brief literature survey is reported.Then, the scope and the contributions of this paper are described and the organization of the manuscript is summarized.

Background
In industrial applications, the final design solution of a new product is almost always an engineering approximation that is intrinsically prone to unpredictable uncertainties [1][2][3][4][5][6][7].A viable approach to the analysis of the approximation in a given design solution is based on the use of a reasonably complex mathematical model of the mechanical system at hand [8][9][10][11].A mathematical model covers almost all physical, economical, biological, industrial, and technical phenomena [12][13][14][15][16][17][18].In particular, dynamic models can be used for the development of new control algorithms suitable for solving several engineering problems [19][20][21][22][23][24].This process is particularly effective in the field of robotics for controlling mobile robots that represent the main object of the present research work [25][26][27][28][29][30].

Formulation of the Problem of Interest for This Study
In the last thirty years, the use of mobile robots has gained great attention because of their broad potential in practical applications oriented toward the solution of important engineering problems [31][32][33][34][35][36].For this purpose, several control approaches were developed for obtaining the stabilization and the trajectory tracking control of mobile robots [37][38][39][40][41][42][43].In particular, in the field of robotics, the advantage of using wheeled mobile robots instead of legged mobile robots is widely accepted.When compared with legged mobile robots, wheeled mobile robots have a simple mechanical structure and, consequently, a mathematical model of this class of robots can be easily obtained.Legged mobile robots, on the other hand, are more versatile machines which make use of mechanical limbs for their movement and, therefore, can easily traverse several different terrains.The advantages in terms of dexterity of legged mobile robots over wheeled mobile robots involve a more complex structural design, a different paradigm for the construction of the robot, an increased power consumption, and require the development of a more complex control system.However, in general, both wheeled and legged mobile robots are subject to some degree of complexity in their mechanical behavior and, more importantly, their motion must satisfy some geometric constraints in order to be suitable for widely exploring the external environment [44][45][46][47][48].One of the fundamental challenges associated with dynamic and control problems of both wheeled and legged mobile robots arise from the mathematical representation of their mechanical models which can be classified as underactuated nonholonomic dynamical systems.Underactuated mechanical systems are dynamical systems in which the number of control actuators is lower than the number of degrees of freedom.Moreover, nonholonomic dynamical systems are mechanical systems in which the constraint equations are of holonomic as well as nonholonomic nature, namely, they involve the complete set of generalized coordinates and their first and second time derivatives [49][50][51][52].Thus, underactuated nonholonomic nonlinear mechanical systems represent a broad class of dynamical systems that are particularly challenging to analyze and control.

Literature Review
In the field of robotics, stabilization and tracking problems of mobile robots have grown in importance in recent years, as testified by a large number of research projects devoted to the investigation of these topics [53][54][55].Furthermore, the nonlinear control problem focused on mobile robots have recently attracted significant attention in the control community because of the wide scope of applications of these mechanical systems [56][57][58][59][60][61][62][63].In several engineering applications, mobile robots are often treated as mechanical systems constrained by holonomic and/or nonholonomic algebraic equations, which can be modeled employing nonlinear dynamics techniques such as the multibody approach to the dynamics of mechanical systems.In particular, wheeled mobile robots, which represent the mobile robots of interest for this investigation, must be modeled as nonholonomic mechanical systems to capture the pure rolling conditions of the wheels.Thus, the nonlinear control problem of this family of mechanical systems represents a challenging engineering issue [64][65][66][67].In the literature, the nonlinear control methods employed for this class of mechanical systems are based on non-standard approaches that cannot be easily extended to both holonomic and nonholonomic mechanical systems [68][69][70].For example, neural networks have been effectively used in recent years to approximate the dynamic behavior of mobile robots and for developing advanced nonlinear control strategies [71,72].The control algorithms obtained by integrating kinematic controllers and neural network computed-torque controllers can be used for the solution of the three basic navigation problems associated with nonholonomic mobile robots, namely the tracking of a reference trajectory, the following of a prescribed path, and the stabilization around the desired configuration [73,74].
In this investigation, an effective control approach based on an inverse dynamic paradigm suitable for handling underactuated mechanical systems constrained by holonomic and/or nonholonomic algebraic equations is developed and its use is demonstrated by means of numerical experiments.The mobile robot that represents the object of this work is a wheeled robot having a unicycle-like underactuated structure.As described in detail in the paper, a unicycle-like mobile robot is a wheeled robot having three degrees of freedom that is endowed with a single steering wheel.This is a common design and is of great interest to engineers.The name unicycle for this type of mobile robot is commonly used also to remark the difference between this mobile robot and the bicycle mobile robot [75,76].As mentioned before, the unicycle mobile robot has three degrees of freedom, is controlled by two control actuators, and is subjected to only one nonholonomic constraint equation associated with the pure rolling condition of the single wheel.The bicycle mobile robot, on the other hand, is endowed with four degrees of freedom, is controlled by using two control actuators, and is subjected to two nonholonomic constraint equations associated with the pure rolling conditions of the two wheels.Typically, the dynamic behavior of unicycle mobile robots is studied for testing nonlinear control strategies which can be directly extended to bicycle mobile robots.The performance of these control algorithms applied to bicycle mobile robots is particularly important because these mobile robots are kinematically equivalent to car-like mobile robots which are mechanically balanced.In the literature, the simple kinematic model this mobile robot is often used for investigating the performance of advanced control algorithms aimed at solving the stabilization, tracking, flocking, and rendezvous problems [77][78][79][80][81].

Scope and Contributions of This Research Work
This research work is focused on the dynamics and control of a unicycle-like mobile robot.The forward and inverse dynamic problems associated with the wheeled mobile robot considered in the paper represent archetypical examples pertaining to the general dynamic and control analysis of underactuated nonholonomic mechanical systems.To address these problems, the computational procedure developed in this investigation is based on an effective combination of the Udwadia-Kalaba Equations with the Underactuation Equivalence Principle.The Udwadia-Kalaba Equations allows for obtaining a clear and concise analytical formulation of dynamic problems which involve holonomic and/or nonholonomic algebraic constraints.Furthermore, as shown in this investigation, a general inverse dynamic method suitable for controlling nonlinear mechanical systems based on the Udwadia-Kalaba Equations can be readily devised and implemented.On the other hand, the Underactuation Equivalence Principle is a method recently developed in the field of analytical dynamics which allows for extending the application of the Udwadia-Kalaba Equations from fully actuated mechanical systems to underactuated mechanical systems.Therefore, as illustrated in this paper by means of detailed analytical derivations and extensive numerical computations, the Udwadia-Kalaba Equations can be effectively employed for solving forward and inverse dynamic problems associated with wheeled mobile robots.The Underactuation Equivalence Principle, on the other hand, is suitable for modeling underactuated nonholonomic mechanical systems and can be used for the synthesis of effective nonlinear control strategies.A unicycle-like mobile robot having one steerable drive wheel that is assumed to be always perpendicular to the ground, which represents one of the most common types of mobile robots employed in engineering applications, is used in this work for demonstrating the effectiveness of the inverse dynamic method devised in this investigation [82][83][84][85].

Organization of the Manuscript
This manuscript is organized as follows.In Section 2, the mathematical background related to Udwadia-Kalaba Equations and to the Underactuation Equivalence Principle is described in detail.In Section 3, the wheeled mobile robot having a unicycle-like underactuated structure used in the paper as a benchmark problem is illustrated and the derivation of a nonlinear tracking controller associated with a planar path described in a parametric form specifically designed for this mobile robot is analyzed.In Section 4, a discussion on the numerical results found by using the approach developed in the paper is reported and the directions for future research works are formulated.In Section 5, the summary of the paper and the conclusions obtained in this investigation are provided.

Background Material and Analytical Methods
In this section, the background material and the analytical methods of interest for this investigation are provided.First, the general form of the fundamental problem of constrained motion is recalled and, subsequently, these equations are rewritten in an analytical form suitable for solving inverse dynamics problems.Thereafter, the Udwadia-Kalaba Equations are discussed, and the Underactuation Equivalence Principle is presented.Moreover, the combined use of the Udwadia-Kalaba Equations with the Underactuation Equivalence Principle necessary for the development of nonlinear control laws suitable for controlling underactuated nonholonomic mechanical systems is analyzed.

Fundamental Problem of Constrained Dynamics
In this subsection, the general form of the fundamental problem of constrained motion is briefly recalled.In the central problem of constrained dynamics, the main goal is to predict the dynamical evolution of a mechanical system subjected to a given set of algebraic constraint equations [86][87][88].For this purpose, the system initial conditions, namely the initial generalized coordinates and velocities, are specified data for the problem at hand.The external active forces not related to the kinematic constraints that are applied to the mechanical system are also assumed as known functions of time [89,90].On the other hand, while the algebraic constraint equations represent specified limitations of the motion of the mechanical system and appear explicitly in the fundamental problem of constrained dynamics, the generalized constraint forces corresponding to the algebraic equations are additional unknowns of the dynamic problem to be solved [91,92].
An effective method for computing a closed-form analytical solution of the fundamental problem of constrained motion is based on the approach originally proposed by Udwadia and Kalaba [93].The analytical method named after Udwadia and Kalaba is described in detail in a series of papers and in a book developed by the same authors [94,95].For the purposes of the present study, the formulation and the solution of the general form of the fundamental problem of constrained motion developed by Udwadia and Kalaba can be summarized in the following steps.First, the equations of motion of a nonlinear mechanical system subjected to a general set of holonomic and/or nonholonomic constraints are analytically formulated by using the classical principles of analytical mechanics.Subsequently, the index-three form of the equations of motion is transformed into its index-one counterpart employing a standard index reduction technique.Finally, the general formulation of the Udwadia-Kalaba Equations is applied to the index-one form of the differential algebraic dynamic equations to obtain the system generalized acceleration vector and the generalized force vector associated with the algebraic constraints.In particular, an equivalent alternative form of the Udwadia-Kalaba Equations that is more suitable for describing the dynamic behavior of mechanical systems subjected to nonlinear algebraic constraints, such as multibody mechanical systems and wheeled mobile robots, is used in this paper.
The first step for the formulation and the solution of the fundamental problem of constrained dynamics is the analytical derivation of the equations of motion of a general mechanical system constrained by holonomic and/or nonholonomic constraints.To this end, the classical principles of analytical mechanics, such as the D'Alembert-Lagrange principle of virtual work formulated in conjunction with the Lagrange multiplier technique or the Lagrange equations, can be effectively used.To illustrate this fact, consider a nonlinear mechanical system described by a set of n b generalized coordinates that are used to represent the configuration space of the dynamical system [96].Assume that the nonlinear mechanical system is initially unconstrained, namely the number of generalized coordinates n b of the dynamical system is equal to the number of the system degrees of freedom m f .In this case, it is well known that the equations of motion of a general unconstrained mechanical system can be directly obtained by using the Lagrange equations of the second kind.Thus, one can write: where t is time, q ≡ q(t) ∈ R n b denotes the generalized coordinate vector of the unconstrained mechanical system having dimension n b = m f , L ≡ L(q, q, t) identifies the Lagrangian function associated with the dynamical structure of the mechanical system, and Q e,nc ≡ Q e,nc (q, q, t) ∈ R n b represents the generalized force vector associated with the nonconservative external forces acting on the mechanical system.The generalized force vector associated with the nonconservative external forces Q e,nc can be found by calculating the virtual work of the nonconservative external forces as follows: where δW e,nc ≡ δW e,nc (q, q, δq) denotes the virtual work of the nonconservative external forces applied on the mechanical system.For a given instant of time t, the generalized coordinate vector q and the generalized velocity vector q serve to determine uniquely the positions and the velocities of all the material points and rigid bodies that form the mechanical system [97].Moreover, one can readily express the Lagrangian function of a mechanical systems as follows: where T ≡ T(q, q, t) identifies the kinetic energy associated with the inertial effects of the mechanical system and U ≡ U(q, t) denotes the potential energy relative to the conservative external forces acting on the dynamical system.For a mechanical system composed of rigid bodies, one can write the total kinetic energy of the nonlinear dynamical system in the following general form: where M ≡ M(q, t) ∈ R n b ×n b is a positive-definite symmetric matrix that represents the system mass matrix.Considering the definition of the kinetic energy T given by Equation ( 4) as well as the general structure of the potential energy U, the application of the Lagrange equations of the second kind (1) leads to: or equivalently: where is the inertia quadratic velocity vector which absorbs the generalized inertia effects associated with the squares of the generalized velocities and Q e,c ≡ Q e,c (q, t) ∈ R n b is the generalized force vector relative to the conservative external forces [98].The inertia quadratic velocity vector Q v and the generalized force vector associated with the conservative external forces Q e,c are respectively given by: It is, therefore, apparent that the definition of the Lagrangian function L of a nonlinear mechanical system and the specification of the virtual work associated with the nonconservative external forces δW e,nc allow for calculating the equations of motion of a nonlinear mechanical system by using the Lagrange equations of the second kind given by Equation (1).For simplicity, one can rewrite the general form of the equations of motion of an unconstrained mechanical system given by Equation ( 6) as follows: where Q b ≡ Q b (q, q, t) ∈ R n b represents the total generalized force vector of the body forces relative to the mechanical system given by: where: where Q e ≡ Q e (q, q, t) ∈ R n b identifies the total generalized force vector associated with the external active forces.The equations of motion of a mechanical system written in the form given by Equation ( 8) constitute a nonlinear set of Ordinary Differential Equations (ODEs) in which the effect of the kinematic constraints is intrinsically embedded.This analytical description of the motion of a mechanical system is referred to as Minimal Coordinate Formulation (MCF) since by using this approach the general configuration of the mechanical system is represented employing a minimal set of generalized coordinates q formed by a vector of n b = m f variables [99].In particular, the analytical approach based on the MCF leads to a set of equations of motion given by Equation ( 8) in which the generalized force vector associated with the algebraic constraints does not appear explicitly in accordance also with the basic assumption of ideal kinematic constraints or workless constraint forces.On the other hand, in principle, the same mechanical system can be equally described by using a Redundant Coordinate Formulation (RCF) in which an analytical approach based on a redundant set of generalized coordinates q is employed in order to univocally describe the geometric configuration of the mechanical system under consideration.In this case, the dimension of the generalized coordinate vector n b is larger than the number of the system degrees of freedom m f .In the RCF, the equations of motion can be readily obtained considering the Lagrange equations of the first kind as follows: where q identifies the redundant set of generalized coordinates having dimension n b > m f , L ≡ L(q, q, t) denotes the Lagrangian function associated with the redundant set of generalized coordinates, Q e,nc ≡ Q e,nc (q, q, t) ∈ R n b is the generalized force vector of the nonconservative external forces, and Q c ≡ Q c (q, q, t) ∈ R n b represents the total generalized force vector associated with the algebraic constraints applied to the mechanical system.Following the same mathematical manipulation previously used for the Lagrange equations of the second kind, the equations of motion of a general mechanical system expressed employing the analytical approach based on the RCF can be written as: where M ≡ M(q, t) ∈ R n b ×n b denotes the system mass matrix associated with the redundant generalized coordinate vector and is the total generalized force vector of the body forces described by using the redundant coordinate approach.In the fundamental problem of constrained motion, the generalized force vector associated with the kinematic constraints Q c represents a set of additional unknowns of the problem.However, one can effectively use the Lagrange multiplier technique to express the analytical form of the generalized constraint force vector in terms of the vector of the Lagrange multipliers [100].To this end, the algebraic equations that model the kinematic constraints can be effectively used.

Udwadia-Kalaba Equations in Forward and Inverse Dynamic Problems
In this subsection, the use of the Udwadia-Kalaba Equations for solving dynamic and control problems of nonlinear mechanical systems constrained by a general set of algebraic equations is discussed [101][102][103].To this end, consider the following general form of the differential algebraic equations of motion of a mechanical system represented by using a redundant set of n b generalized coordinates whose motion is restricted by a set of n c algebraic constraints: where t is the time variable, q ≡ q(t) ∈ R n b is the generalized coordinate vector that represents that configuration of the mechanical system, q ≡ q(t) ∈ R n b represents the system generalized velocity vector, q ≡ q(t) ∈ R n b identifies the system generalized acceleration vector, M ≡ M(q, t) ∈ R n b ×n b denotes the system mass matrix, Q b ≡ Q b (q, q, t) ∈ R n b identifies the total generalized force vector that acts on the mechanical system which includes the generalized external forces and the quadratic velocity inertial forces, λ ≡ λ(t) ∈ R n c represents the total vector of Lagrange multipliers associated with the constraint equations, A ≡ A(q, t) ∈ R n c ×n b is the total Jacobian matrix relative to the constraint equations, and b ≡ b(q, q, t) ∈ R n c is the total constraint quadratic velocity vector.The index-three differential algebraic form of the dynamic equations is a general mathematical representation of the equations of motion of a nonlinear mechanical system that includes the effect of holonomic as well as nonholonomic constraints.In fact, considering a set of n f holonomic constraints defined at the position level, a set of n g nonholonomic constraints defined at the velocity level, and a set of n h nonholonomic constraints defined at the acceleration level, one can write the following general form of the algebraic equations: where f ≡ f(q, t) ∈ R n f is a nonlinear vector function associated with the holonomic constraint equations defined at the position level, g ≡ g(q, q, t) ∈ R n g denotes a nonlinear vector function relative to the nonholonomic constraint equations defined at the velocity level, and h ≡ h(q, q, q, t) ∈ R n h identifies a nonlinear vector function corresponding to the nonholonomic constraint equations defined at the acceleration level.This general set of constraint equations is included in the general form of the equations of motion by means of the total constraint Jacobian matrix denoted with A as well as considering the total constraint quadratic velocity vector indicated with b, which can be respectively obtained as follows: where the total number of algebraic constraints is given by n and A h ≡ A h (q, t) ∈ R n h ×n b respectively represent the Jacobian matrices relative to the algebraic constraint equations f, g, and h, while b f ≡ b f (q, q, t) ∈ R n f , b g ≡ b g (q, q, t) ∈ R n g , and b h ≡ b h (q, q, t) ∈ R n h respectively represent the quadratic velocity vectors associated with the algebraic constraint equations f, g, and h.These matrix and vector quantities are respectively defined as follows: where, for simplicity, the following structure of nonholonomic constraint equations defined at the acceleration level is assumed: f q q + f q q q q + 2f qt q + f tt = 0 g q q + g q q + g t = 0 h q q + c = 0 (17) where the partial derivatives and the Jacobian matrices that appear in these equations are respectively defined as: where the term c ≡ c(q, q, t) ∈ R n h is assumed as a known nonlinear vector function.On the other hand, the Udwadia-Kalaba Equations can be written in a general form suitable for handling dynamic problems of nonlinear mechanical systems constrained by a general set of holonomic as well as nonholonomic constraints [104].For this purpose, one can write: where a b ≡ a b (q, q, t) ∈ R n b is the generalized acceleration vector of the mechanical system released from the constraint equations, a c ≡ a c (q, q, t) ∈ R n b is the generalized acceleration vector of the mechanical system induced by the constraint equations, F ≡ F(q, q, t) ∈ R n c ×n c represents the constraint feedback matrix, and e ≡ e(q, q, t) ∈ R n c is the error vector associated with the generalized accelerations.The unconstrained acceleration vector a b can be explicitly computed by ignoring the action of the generalized constraint forces as follows: Furthermore, the acceleration error vector e can be analytically obtained by substituting the generalized acceleration vector corresponding to the unconstrained mechanical system a b into the set of constraint equations written in the standard form defined at the acceleration level: The constraint feedback matrix F, on the other hand, can be obtained by computing the Moore-Penrose pseudoinverse matrix of a special matrix associated with the constrained mechanical system.To this end, one can write: where the + symbol in the superscript indicates the Moore-Penrose pseudoinverse matrix and K ≡ K(q, q, t) ∈ R n c ×n c is called the kinetic matrix of the constrained mechanical system.The kinetic matrix K is defined as follows: The kinetic matrix K of a mechanical system constrained by a general set of holonomic as well as nonholonomic constraint equations completely defines the dynamic structure of a mechanical system as a result of the application of the nonlinear algebraic equations which represent general restrictions of the system motion.By using the Udwadia-Kalaba Equations, the generalized force vector Q c induced by the algebraic constraints can be readily obtained as follows: Finally, one can easily compute the generalized acceleration vector induced by the algebraic constraints a c as follows: The complete set of equations which form the Udwadia-Kalaba Equations can be effectively used for solving forward and inverse dynamic problems of nonholonomic mechanical systems modeled in a seamless analytical framework.

Underactuation Equivalence Principle
In this subsection, the Underactuation Equivalence Principle is presented.For this purpose, the key aspects of this complementary principle of mechanics can be described as follows.Firstly, underactuated mechanical systems are dynamical systems in which the number of the control actions n u is lower than the number of the system degrees of freedom n b .The general structure of the equations of motion of an underactuated mechanical system can be written as follows: where B u ≡ B u (q, t) ∈ R n b ×n u is an influence matrix associated with the input control actions that defines the collocations and the combinations of the control inputs, whereas u ≡ u(t) ∈ R n u is vector of control inputs.If a mechanical system is underactuated, then the rank r B u of the input influence matrix B u is lower than the number of the system degrees of freedom n b [105].According to the Underactuation Equivalence Principle, one can replace the previous set of ordinary differential equations of motion with the following set of differential algebraic equations which involves a zero generalized constraint force vector: where it can be proved that the corresponding generalized constraint force vector Q c is identically equal to the zero vector.If a mechanical system is underactuated, the rank of the input influence matrix B u is lower than the number of system degrees of freedom, namely rank(B u ) < n b .Therefore, it is apparent that the inverse dynamic problem for an underactuated mechanical system is a challenging control problem even if the generalized control force vector Q c can be simply written as a linear combination of the entries of the vector of control actions u.However, for a general class of constrained mechanical systems, the inverse dynamic problem for underactuated mechanical systems can still be addressed using the Udwadia-Kalaba Equations considering the Underactuation Equivalence Principle.This principle states that an unconstrained mechanical system having an underactuated structure of the control actions is dynamically equivalent to a constrained mechanical system subjected to a set of nonholonomic constraints defined at the acceleration level that generates a generalized constraint force vector identically equal to the zero vector.In general, two different mechanical systems are said to be dynamically equivalent when they are characterized by the same generalized acceleration vector.The validity of the Underactuation Equivalence Principle can be readily demonstrated by means of simple mathematical derivations.For this purpose, assume for simplicity the following simplified form of the input influence matrix B u : One can partition the equations of motion of an unconstrained mechanical system having an underactuated structure as follows: where: where u 1 ≡ u 1 (t) represents the effective vector of control actions, q 1 ≡ q 1 (t) is the set of generalized coordinates directly influenced by the control actions, and q 2 ≡ q 2 (t) denotes the set of generalized coordinates not directly influenced by the control actions, whereas M 1,1 ≡ M 1,1 (q, t), M 1,2 ≡ M 1,2 (q, t), and M 2,2 ≡ M 2,2 (q, t) represent the block matrices that form the partition of the mass matrix, while Q b,1 ≡ Q b,1 (q, q, t) and Q b,2 ≡ Q b,2 (q, q, t) identify the partitions of the total generalized force vector.Employing the Underactuation Equivalence Principle, one can consider the underactuated structure of the control actions applied to mechanical system as an equivalent set of nonholonomic constraints defined at the acceleration level.Therefore, one can demonstrate that the following set of equations of motion are associated with a mechanical system that is dynamically equivalent to the system under study: where Q c is the generalized force vector associated with the set of nonholonomic constraints that represent the underactuation property of the original mechanical system, while the constraint matrix A and the constraint quadratic velocity vector b can be readily obtained by using the Underactuation Equivalence Principle applied to the partitioned equations of motion (29).These matrix and vector quantities are respectively given by: It can be easily demonstrated that the nonholonomic constraint equations associated with the Underactuation Equivalence Principle given by Equation (32) lead to an identically zero generalized constraint force vector Q c .To this end, assume that the block matrix M 1,1 (q, t) is not a singular matrix and consider the Woodbury matrix identity.By doing so, the inverse of the mass matrix M can be obtained as follows: where the matrices C 1 ≡ C 1 (q, t), C 2 ≡ C 2 (q, t), and G ≡ G(q, t) are respectively defined as: Considering the analytical formulation of the inverse of the mass matrix M, the central equations of constrained motion can be readily applied to the problem at hand.For this purpose, the generalized acceleration vector of the mechanical system obtained in absence of the algebraic constraints can be computed as follows: Subsequently, using the analytical expression of the generalized acceleration vector associated with the unconstrained mechanical system a b , one can compute the acceleration error vector e as: This equation shows that the acceleration error vector e associated with the mechanical system subjected to the nonholonomic constraints defined at the acceleration level arising from the Underactuation Equivalence Principle is identically equal to zero and, therefore, the corresponding generalized constraint force vector Q c is a zero vector as well.When this condition is met, there is a particular set of generalized accelerations that cannot be altered by the implementation of any vector of control inputs u.Thus, there are some trajectories in the configuration space of the underactuated mechanical system that cannot be imposed to follow by using any external control actions.This analytical result is consistent with the fact that the set of nonholonomic underactuation constraints arises from a subset of the equations of motion of the unconstrained mechanical system that leads to the unconstrained generalized acceleration vector a b .Consequently, the acceleration vector arising from the underactuation constraints is consistent with the space of the unconstrained motion of the mechanical system.Considering the Underactuation Equivalence Principle, the underactuation property of a dynamical system can be converted into a set of acceleration nonholonomic constraints to be able to analytically solve the inverse dynamic problem relative to underactuated mechanical systems employing the Udwadia-Kalaba Equations.Thus, the analytical method based on the Underactuation Equivalence Principle allows for extending the Udwadia-Kalaba nonlinear control approach to underactuated mechanical systems.By using the Underactuation Equivalence Principle, the Udwadia-Kalaba Equations can be applied to underactuated mechanical systems by formulating a proper set of nonholonomic constraints at the acceleration level that consistently reflects the underactuation property of the mechanical system.By doing so, the underactuation property of a mechanical system is recognized as a limitation of the system motion arising from the collocation of the control actuators and this property is correctly taken into account in the formulation of the stabilization or tracking control problem to be solved.To achieve this goal, one can assume that the constraint equations represent the desired behavior of the mechanical system on hand instead of modeling some physical limitations on the motion of the system.Therefore, the control analyst needs to formulate in a proper manner the constraint equations associated with the mechanical joints as well as the constraint equations representing the desired behavior of the mechanical system under study.Consequently, both the forward and inverse dynamic formulations are carried out in the same analytical framework and can be analytically treated by using the Udwadia-Kalaba Equations combined with the Underactuation Equivalence Principle discussed in this section.It is important to note that, by using the approach described in this section, one can formulate some control requirements expressed also in terms of nonholonomic constraint equations instead of using only holonomic algebraic constraints, leading to a more flexible and general formulation of the goals for the control problem to be solved.

Numerical Results
In this section, a set of numerical results and an evaluation of the performance of the method proposed in the paper are provided to demonstrate the application of the Udwadia-Kalaba Equations in conjunction with the Underactuation Equivalence Principle to nonlinear control problems.For this purpose, the problem of the tracking control of a wheeled mobile robot is considered.First, the description of the unicycle-like mobile robot employed as a benchmark problem is illustrated.Subsequently, the derivation of a nonlinear tracking controller is analyzed, and the numerical results found using the approach developed in the paper are presented.To achieve these goals, the equations of motion of the wheeled robot are obtained in a compact analytic form.Finally, the numerical simulations performed in the dynamic analysis are described considering different operative conditions.

Unicycle-Like Mobile Robot Model
In this subsection, the description of the unicycle-like mobile robot used as a demonstrative example is presented.A schematic representation of the wheeled mobile robot considered herein is shown in Figure 1.The dynamical system considered in this subsection is endowed with n b = 3 degrees of freedom.This mechanical system is subjected to a nonlinear force field which mathematically represents the dissipation phenomena.Furthermore, the mobile robot must satisfy a nonholonomic constraint equation which represents the pure rolling condition and, therefore, n c = 1.Thus, the unicycle-like mobile robot is described by the following symbolic set of equations of motion: where q is the configuration vector of the mobile robot, M represents the mass matrix of the mobile robot, Q σ is the generalized force vector associated with the dissipative effects which influence the dynamic behavior of the mobile robot, Q r identifies the generalized force vector arising from the pure rolling nonholonomic constraint equation that can be explicitly obtained by applying the Udwadia-Kalaba Equations, A r is the pure rolling constraint Jacobian matrix, and b r represents the pure rolling constraint quadratic velocity vector.These vector and matrix quantities are respectively given by: and where m is the robot mass, I zz denotes the robot moment of inertia, σ is the viscous damping coefficient that models the dissipation forces acting on the mobile robot, x represents the abscissa of the robot center of mass, y represents the ordinate of the robot center of mass, and θ is the angular orientation of the mobile robot.In particular, the generalized force vector arising from the pure rolling nonholonomic constraint Q r was explicitly obtained by applying the general form of the Udwadia-Kalaba Equations to the pure rolling constraint equations defined in the standard form described by the vector and matrix quantities A r and b r .To this end, Equation ( 19) can be effectively used.To achieve this goal, the pure rolling condition of the mobile robot is modeled as a nonholonomic constraint equation assuming that the velocity v G of the robot center of mass G must be orthogonal to the direction vector j associated with the robot wheels for each instant of time.Therefore, the vectors v G and j must be orthogonal for each instant of time to model the pure rolling of the wheeled mobile robot.Thus, one can write: where: Equation ( 40) represents a nonholonomic constraint equation that involves the generalized configuration vector of the mobile robot q and its time derivative q.By imposing this nonlinear condition, one obtains an algebraic equation which involves both the system generalized positions and velocities.This algebraic equation cannot be rewritten only in terms of the system generalized configuration vector q.Therefore, the pure rolling condition is a nonholonomic constraint equation.The Udwadia-Kalaba Equations, on the other hand, were used for obtaining in a concise form the generalized force vector Q r starting from the analytical definition of the pure rolling nonholonomic constraint equation given by Equation (40).To achieve this goal, one can effectively employ the following step-by-step derivation.First, the kinetic matrix K r associated with the mobile robot subjected to the pure rolling condition can be computed as follows: where A r and b r are the vector and matrix quantities that represent the pure rolling condition in the standard form associated with nonholonomic constraint equations.The constraint feedback matrix F r corresponding to the kinetic matrix K r relative to the pure rolling condition of the mobile robot is given by: The unconstrained generalized acceleration vector a r obtained in absence of external force fields acting on the mobile robot is a zero vector: Therefore, the error vector obtained by introducing the unconstrained generalized acceleration vector a r into the pure rolling nonholonomic constraint equation can be readily calculated to yield: Finally, one can obtain the generalized force vector Q r associated with the pure rolling condition of the mobile robot as follows: The explicit calculation of the generalized force vector Q r is one of the main results of this work which allows for performing dynamical simulation by using standard numerical integration techniques based on a sequential marching of the numerical solution on the time grid.

Nonlinear Trajectory Tracking
In this subsection, the derivation of a nonlinear tracking controller based on the inverse dynamic approach developed in the paper is presented.For this purpose, consider the application of an external control action on the mobile robot having an underactuated structure, which implies the following set of equations of motion: where u is a vector of external control actions and B u denotes an input influence matrix that are respectively given by: where the number of the active control actions is n u = 2, while u 1 is the control torque applied on the right wheel of the mobile robot, u 2 is the control torque applied on the left wheel of the mobile robot, L denotes half the length of the robot axle track, and R corresponds to the radius of the robot wheels.To obtain a nonlinear tracking controller for the mobile robot and, at the same time, to impose the underactuation property of this mechanical system by considering a nonholonomic constraint equation defined at the acceleration level, one can assume the following nonholonomic constraint equations defined in the standard form: where A t is the trajectory tracking constraint Jacobian matrix and b t represents the trajectory tracking constraint quadratic velocity vector.These matrix and vector quantities are respectively defined as: where ξ represents a prescribed function of time for the abscissa x of the robot center of mass and ϕ denotes a prescribed function of time for the angular orientation θ of the robot chassis, whereas the first equation in these matrix expressions is a mathematical representation of the underactuation requirement which can be readily obtained by applying the Underactuation Equivalence Principle given by Equation ( 29).These equations are written at the acceleration level and arise from the following set of nonholonomic constraint equations: By imposing the previous set of constraint equations on the unicycle-like mobile robot, one obtains the generalized force vector associated with the nonlinear controller for the trajectory tracking by means of the Udwadia-Kalaba Equations given by Equation ( 19).This generalized force vector is denoted with Q t and is given by: This generalized force vector allows for imposing on the unicycle-like mobile robot a set of control actions suitable for tracking a prescribed state-space trajectory defined in terms of the arbitrary functions ξ and ϕ.However, it is important to note that the functions ξ and ϕ, which define the desired state-space trajectory, must be also consistent with the pure rolling condition imposed on the unicycle-like mobile robot.

Dynamic Analysis
In this subsection, the numerical results found using the approach developed in the paper are presented.For this purpose, consider the numerical data for the unicycle-like mobile robot reported in Table 1.In the dynamic analysis, a time span equal to T = 10 (s) is considered and a step size equal to ∆t = 10 −3 (s) is assumed.The computational algorithm used for performing the time marching of the numerical solution of the equations of motion is the sixth-order Adams-Bashforth method.To demonstrate the effectiveness of the proposed method for performing the tracking controller of the unicycle-like mobile robot, consider the following parametric form for the system path and time law: where ξ is the assigned horizontal displacement, η is the assigned vertical displacement, ϕ is the assigned angular displacement, and γ is the assigned time law.Thus, the time law chosen for the mobile robot is a uniformly accelerated motion whereas the planar curve considered for the path is an astroid.The numerical values of the constant parameters that appear in these equations are respectively defined in Table 2. To test the performance of the inverse dynamic method developed in this paper in multiple scenarios, three different sets of initial conditions are considered for the unicycle-like mobile robot.For simplicity, the three set of initial conditions considered in this work are respectively labeled with the numbers 1, 2, and 3 as shown in Table 3.
In Table 3, considering a general set of initial conditions labeled with the number i = 1, 2, 3, x i 0 is the initial horizontal displacement, y i 0 represents the initial vertical displacement, θ i 0 denotes the initial angular displacement, u i 0 is the initial horizontal velocity, v i 0 represents the initial vertical velocity, and ω i 0 denotes the initial angular velocity.The resulting horizontal displacement, vertical displacement, and angular displacement of the mechanical system obtained in the dynamical analysis by using the designed tracking control actions are respectively shown in Figure 2a-c.The time law used for performing the dynamic analysis is represented in Figure 3a, while the planar path of the unicycle-like mobile robot resulting from the application of the control inputs obtained by using the proposed inverse dynamic approach is represented in Figure 3b.The resulting horizontal velocity, vertical velocity, and angular velocity of the mechanical system obtained in the dynamical analysis by using the designed tracking control actions are respectively shown in Figure 4a-c.

Descriptions Symbols Data (Units)
Initial horizontal displacement x 1 0 , x 2 0 , x 3  In Figures 2a-c, 3a,b, 4a-c and 5a,b, the solid lines refer to the set of initial conditions labeled with the number 1, the dashed lines refer to the set of initial conditions labeled with the number 2, and the dotted lines refer to the set of initial conditions labeled with the number 3, respectively.

Discussion
In this section, a critical discussion of the numerical results obtained in the paper is reported and some directions for future research are provided.

Performance Analysis
The motion of the mobile robot represented in Figure 2a,b demonstrates that the method proposed in the paper can be effectively used for devising nonlinear control actions for mechanical systems constrained by holonomic as well as nonholonomic constraints.To quantify the performance of the proposed approach, one can compute the deviation of the actual state-space trajectory obtained by using the nonlinear controller from the theoretical state-space trajectory imposed as the reference solution.To this end, the absolute errors between the actual state-space trajectory and the desired state-space trajectory can be used for evaluating the performance of the control method proposed in the paper.The absolute error for the planar path is simply defined as: where e x (t) and e y (t) respectively denote the absolute errors of the horizontal and vertical displacements x and y. Figure 6a shows the time evolution of the absolute error related to the horizontal displacement, while Figure 6b shows the time evolution of the absolute error related to the vertical displacement.In Figure 6a,b, the solid lines refer to the set of initial conditions labeled with the number 1, the dashed lines refer to the set of initial conditions labeled with the number 2, and the dotted lines refer to the set of initial conditions labeled with the number 3, respectively.Furthermore, the root-mean-square deviations between the actual state-space trajectory and the desired state-space trajectory are used in this investigation as a quantitative metric.The root-mean-square deviations are defined as follows: where x RMS and y RMS respectively denote the RMS deviation of the horizontal and vertical displacements x and y.The numerical results found for the horizontal and vertical root-mean-square deviations are reported in Table 4.The numerical results found for the root-mean-square deviations are referred to the three sets of initial conditions considered in this research work that are labeled with the numbers 1, 2, and 3. Table 4 shows that the numerical solution obtained in the dynamic simulation involves small deviations from the designed state-space trajectory.These errors are caused by the numerical approximation used for obtaining the numerical solution of the equations of motion which affect each numerical integration scheme suitable for solving nonlinear sets of ordinary differential equations.Therefore, the numerical results found in this section show a good agreement between the actual state-space trajectory obtained using the proposed inverse dynamic method and the desired state-space trajectory arising from a preliminary motion planning phase, thereby demonstrating the effectiveness of the methodology developed in this paper.

General Considerations
In this paper, the starting point is the analysis of the robot kinematic model, while the complete nonlinear dynamic model of the unicycle robot is subsequently employed for testing the inverse dynamic method proposed in this investigation to obtain an analytical solution for the tracking control problem.Because of the no side-slip condition of the wheel, which is kinematically equivalent to a no-slip condition, the unicycle robot cannot move in the direction of the wheel axis and, therefore, the motion of the robot must take place only in the direction of its heading.This pure rolling condition implies that the velocity vector of the unicycle robot must be parallel to the longitudinal direction of the robot chassis.Therefore, the no side-slip condition must be mathematically modeled as a nonholonomic algebraic equation.In general, mechanical systems subjected to nonholonomic constraints are difficult to control and the main goal of this paper is to provide a general method for addressing this challenge.In this research work, an effective inverse dynamic control method based on the combined application of the Udwadia-Kalaba Equations with the Underactuation Equivalence Principle was developed and tested in the case of a unicycle mobile robot.In principle, the approach developed in the paper can be directly extended to the solution of the control problem of other wheeled mobile robots, such as the bicycle mobile robot, as well as to legged mobile robot.However, these developments will be the objects of future investigations.
The nonlinear control problem of wheeled mobile robots represents a challenging issue.Several effective methods were proposed in the literature for solving this problem [56][57][58][59][60][61][62][63].However, the main difficulty associated with the nonlinear control approaches employed in the literature for wheeled mobile robots is that they rely on non-standard techniques which cannot be readily extended to other control problems [71][72][73][74].Unlike the methods already available in the literature, an inverse dynamic method based on the use of the Udwadia-Kalaba Equations was developed in this work.The proposed method is effective, general, and can be easily applied to different control problems related to mobile robots.Furthermore, in principle, the inverse dynamic method developed in this study can be combined with the conventional optimal control method for obtaining the design of a feedforward plus feedback control architecture.This promising combination will be investigated in future research works.
Future research will be devoted to the application of the nonlinear control algorithm developed in this paper to more complex mechanical systems having an underactuated structure and constrained by holonomic as well as nonholonomic algebraic equations.For instance, the immediate extension of this work, which is focused on the stabilization and tracking problem of a unicycle-like mobile robot, is the development of an inverse dynamic computational procedure applied to a bicycle model of a wheeled mobile robot subjected to a set of nonholonomic constraints that analytically define the pure rolling condition of the robot wheels.Moreover, further directions of future research are related to dynamic and control problems of multibody mechanical systems, namely mechanical systems formed by rigid as well as flexible components interconnected by kinematic joints, which comprise many degrees of freedom and represent dynamical systems particularly challenging to mathematically model and computationally analyze.Another important direction for future research works to be explored is the use of the methods of optimal control theory in conjunction with the inverse dynamic approach based on the Udwadia-Kalaba Equations combined with the Underactuation Equivalence Principle developed in this investigation.

Summary and Conclusions
The scientific research of the authors is focused on three main areas of interest that are closely related to practical applications in the field of mechanical engineering.The general topics addressed by the investigation of the authors are system identification, nonlinear control, and multibody dynamics [106][107][108][109][110]. System identification is aimed at obtaining a good estimation of the parameters that define a dynamical model of a given mechanical system.Nonlinear control deals with the development of nonlinear control strategies for mechanical systems described by nonlinear mathematical models.Multibody dynamics is a branch of analytical mechanics devoted to the study of the dynamic behavior of articulated systems composed of rigid or flexible bodies constrained by mechanical joints.
This paper is concerned with forward and inverse dynamic problems of wheeled mobile robots having an underactuated structure and is part of a wider research plan inspired by the vision of the authors.In particular, this work is focused on the dynamic and control problem of wheeled mobile robots having a unicycle-like structure.This type of mobile robots represents an important example of a large class of underactuated nonholonomic nonlinear mechanical systems which is used in the field of nonlinear control for testing the effectiveness of new computational algorithms.The method developed in the paper is based on the use of the Udwadia-Kalaba Equations in conjunction with the Underactuation Equivalence Principle.It is shown in this paper that the Udwadia-Kalaba Equations can be employed for deriving the generalized force vector necessary for solving both forward and inverse dynamic problems of nonholonomic mechanical systems.Also, this investigation demonstrates that the Underactuation Equivalence Principle allows for extending the effectiveness of the use of the Udwadia-Kalaba Equations from fully actuated mechanical systems to underactuated mechanical systems.To demonstrate this fact, the nonlinear trajectory tracking problem of a unicycle-like mobile robot is considered in this work as an illustrative example.Dynamical simulations and numerical experiments corroborate the analytical developments obtained in this paper.In particular, the numerical results found for the unicycle robot show a small error between the desired trajectory and the actual trajectory followed by the mobile robot.
In future research works, the method developed in this paper will be used as an open-loop (feedforward) control law and will be combined with the optimal control technique that is useful for the design of compensator controllers based on a closed-loop (feedback) architecture.Furthermore, the method developed in this work will be applied to the solution of inverse dynamic problems of more complex mechanical systems like, for example, the bicycle model of a wheeled mobile robot.

Table 1 .
Mobile robot data.

Table 3 .
Sets of initial conditions.

Table 4 .
State-space trajectory root-mean-square deviations for different sets of initial conditions.