# Mechanical Systems With Nonideal Constraints: Explicit Equations Without the Use of Generalized Inverses

Firdaus Udwadia, Robert Kalaba, Phailaung Phohomsiri

► To cite this version:

Firdaus Udwadia, Robert Kalaba, Phailaung Phohomsiri. Mechanical Systems With Nonideal Constraints: Explicit Equations Without the Use of Generalized Inverses. Journal of Applied Mechanics, 2004, 71 (5), pp.615-621. 10.1115/1.1767844. hal-01399623

HAL Id: hal-01399623

https://hal.science/hal-01399623

Submitted on 19 Nov 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Mechanical Systems With Nonideal Constraints: Explicit Equations Without the Use of Generalized Inverses

In this paper we obtain the explicit equations of motion for mechanical systems under nonideal constraints without the use of generalized inverses. The new set of equations is shown to be equivalent to that obtained using generalized inverses. Examples demonstrating the use of the general equations are provided.

Firdaus E. Udwadia

Professor of Civil Engineering,
Aerospace and Mechanical Engineering,
Mathematics, and Information and Operations
Management,
University of Southern California,
Los Angeles, CA 90089
e-mail: fudwadia@usc.edu

Robert E. Kalaba

Professor of Biomedical Engineering,
Electrical Engineering, and Economics,
University of Southern California, Los
Angeles, CA 90089

Phailaung Phohomsiri

Graduate Student,
Department of Aerospace and Mechanical
Engineering,
University of Southern California,
Los Angeles, CA 90089

## 1 Introduction

When constraints are applied to mechanical systems, additional forces of constraint are produced that guarantee their satisfaction. The development of the equations of motion for constrained mechanical systems has been pursued by numerous scientists and mathematicians, like Appell [1], Beghin [2], Chetaev [3], Dirac [4], Gauss [5], Gibbs [6], and Hamel [7]. All these investigators have used as their starting point the D'Alembert-Lagrange Principle. This principle, which was enunciated first by Lagrange in his *Mechanique Analytique*, [8], can be presumed as being, at the present time, at the core of classical analytical dynamics.

D'Alembert's principle makes an assumption regarding the nature of constraint forces in mechanical systems, and this assumption seems to work well in many practical situations. It states that the total work done by the forces of constraint under virtual displacements is always zero. In 1992 Udwadia and Kalaba [9] obtained a simple, explicit set of equations of motion, suited for general mechanical systems, with holonomic and/or nonholonomic constraints. Though their equations encompass time dependent constraints that are (1) not necessarily independent, and (2) nonlinear in the generalized velocities, their equations are valid only when D'Alembert's principle is observed by the constraint forces.

However, in many situations in nature, the forces of constraint in mechanical systems do not satisfy D'Alembert's principle. As stated in Pars's *A Treatise on Analytical Dynamics* [10], "There are in fact systems for which the principle enunciated [D'Alembert's principle] . . . does not hold. But such systems will not be considered in this book." Such systems have been considered to lie beyond the scope of Lagrangian mechanics. Recently, Udwadia and Kalaba [11,12] have developed general, explicit equations of motion for constrained mechanical systems that may or may not satisfy D'Alembert's principle. The statement of their result involves the use of generalized inverses of various matrix quantities, and they derive their results by using the special properties of these generalized inverses.

In this paper we give a new, alternative set of explicit equations that describes the motion of constrained mechanical systems that may or may not satisfy D'Alembert's principle. Thus these equations are valid when the forces of constraint may do work under virtual displacements. We show here that there is no need to use any concepts related to generalized inverses in the development of these general equations. The explicit equations developed herein can handle time dependent constraints that are nonlinear in the generalized velocities, as do the equations obtained using generalized inverses. Instead of relying on the properties of generalized inverses, our explicit equations rely on a deeper understanding of virtual displacements as provided in Refs. [13,14].

After obtaining the new equations, we show that they are indeed equivalent to those given earlier by Udwadia and Kalaba [11,12] which make extensive use of generalized inverses. Three illustrative examples are provided showing the use of the new equations. The last example deals with sliding friction.

## 2 Explicit Equations of Motion for Mechanical Systems With Nonideal Constraints

For an unconstrained system of $N$ particles, Lagrange's equation of motion for the system at time $t$ can be written, using generalized coordinates, as

$$M(q,t)\ddot{q} = F(q, \dot{q}, t); \quad q(0) = q_0, \dot{q}(0) = \dot{q}_0, \quad (1)$$

where, $q$ is the generalized coordinate $n$-vector $q = [q_1, q_2, \dots, q_n]^T$; $M$ is an $n$ by $n$ symmetric positive definite matrix; and, $F(q, \dot{q}, t)$ is the $n$-vector of the "given" force which is a known function of $q, \dot{q}$, and time, $t$. The number of degrees-of-freedom of the system is equal to the number of generalized coordinates, $n$, characterizing the configuration of the system at any time, $t$. The acceleration, $a(t)$, of the unconstrained system described by Eq. (1) is then given by $a(t) = M(q,t)^{-1}F(q, \dot{q}, t)$.

Let the system described by Eq. (1) be now further constrained by the $m$ constraint equations

$$\varphi_i(q, \dot{q}, t) = 0, \quad i = 1, 2, \dots, m, \quad (2)$$
in which $k<m$ of these constraint equations are independent. We shall assume that the constraint equations satisfy the initial conditions given in Eq. (1). Equation set (2) includes both holonomic and nonholonomic constraints. Assuming sufficient smoothness, we can differentiate equation set (2) with respect to time $t$ to obtain

$$A(q, \dot{q}, t)\ddot{q} = b(q, \dot{q}, t), \quad (3)$$

where the elements of $A$ and $b$ are known functions of $q, \dot{q}$, and $t$ and the matrix $A$ is an $m$ by $n$ matrix that has rank $k$.

The presence of the constraints causes additional constraint forces to arise at each instant of time to assure that the constraints are satisfied. The equation of motion for the constrained system can be then expressed as

$$M\ddot{q} = F(q, \dot{q}, t) + F^c(q, \dot{q}, t), \quad (4)$$

where, $F$ is the given force, and $F^c$ is the additional force engendered by the presence of the constraints.

Premultiplying Eq. (4) by $M^{-1/2}$, we have

$$M^{1/2}\ddot{q} = M^{-1/2}F(q, \dot{q}, t) + M^{-1/2}F^c(q, \dot{q}, t), \quad (5)$$

which can be written as

$$\ddot{q}_s(t) - F_s^c = a_s(t). \quad (6)$$

Here we have denoted the “scaled” acceleration of the constrained system, $M^{1/2}\ddot{q}$, by $\ddot{q}_s(t)$, the “scaled” force of constraint, $M^{-1/2}F^c$, by $F_s^c$, and the “scaled” acceleration of the unconstrained system, $M^{-1/2}F$, by $a_s(t)$. In the same manner, the constraint Eq. (3) at time $t$ can be expressed as $(AM^{-1/2})(M^{1/2}\ddot{q}) = b$. Denoting $AM^{-1/2}$ by the $m$ by $n$ matrix $B$, we obtain

$$B\ddot{q}_s = b. \quad (7)$$

A virtual displacement (see Refs. [13,14]) is any nonzero vector, $w$, that satisfies the equation

$$A(q, \dot{q}, t)w = 0. \quad (8)$$

When the constraints are nonideal, the work done, $W(t)$, by the constraint force, $F^c$, under virtual displacements, $w$, needs to be specified through knowledge of the $n$-vector $C$, so that, [11]

$$W(t) = w^T F^c(q, \dot{q}, t) = w^T C(q, \dot{q}, t), \quad (9)$$

where $C(q, \dot{q}, t)$ is a known $n$-vector, and characterizes the nature of the nonideal constraint force, $F^c$. This is an extension of D'Alembert's principle.

Equation (8) can be rewritten as

$$(AM^{-1/2})(M^{1/2}w) = 0. \quad (10)$$

Similarly, Eq. (9) can be rewritten as

$$(w^T M^{1/2})(M^{-1/2}F^c) = (w^T M^{1/2})(M^{-1/2}C). \quad (11)$$

Denoting $v = M^{1/2}w$, Eq. (10) becomes

$$Bv = 0. \quad (12)$$

Since $M^{1/2}$ is nonsingular, $v$ is then any nonzero vector such that relation (12) is satisfied. Furthermore, after denoting $C_s = M^{-1/2}C$, Eq. (11) can be written as

$$v^T F_s^c = v^T C_s. \quad (13)$$

Since $B$ has rank $k$, there are $n-k$ linearly independent vectors, $v$, such that $Bv = 0$. Assembling then such vectors $v_1, \dots, v_{n-k}$ in the matrix $V$, we obtain

$$V^T F_s^c = V^T C_s. \quad (14)$$

The matrix $V$ can be constructed by a judicious use of the Gram-Schmidt procedure.

Consider the linear Eqs. (6), (7), and (14). These equations can be expressed as

$$Lr = \begin{bmatrix} [I]_{nxn} & -[I]_{nxn} \\ [B]_{mxn} & [0]_{mxn} \\ [0]_{(n-k)xn} & [V^T]_{(n-k)xn} \end{bmatrix} \begin{bmatrix} \ddot{q}_s \\ F_s^c \end{bmatrix} = \begin{bmatrix} a_s \\ b \\ V^T C_s \end{bmatrix} = s, \quad (15)$$

where $L$ is a $(2n+m-k)$ by $2n$ matrix, $r$ is a $2n$-vector, and $s$ is a $(2n+m-k)$-vector.

The equation set (15) constitutes the fundamental linear set of equations that needs to be solved to obtain both the scaled acceleration, $\ddot{q}_s$, of the constrained system as well as the scaled constraint force, $F_s^c$. In what follows, we shall show that a solution to this linear system of equations exists and is unique.

We premultiply both sides of Eq. (15) by $L^T$ to obtain the equation

$$L^T L r = \begin{bmatrix} I & B^T & 0 \\ -I & 0 & V \\ 0 & 0 & V^T \end{bmatrix} r = \begin{bmatrix} I & B^T & 0 \\ -I & 0 & V \\ 0 & 0 & V^T C_s \end{bmatrix}. \quad (16)$$

Let us denote

$$D = B^T B, \quad (17)$$

and

$$E = VV^T. \quad (18)$$

Equation (16) can be written as

$$Gr = \begin{bmatrix} [I+D]_{nxn} & -[I]_{nxn} \\ -[I]_{nxn} & [I+E]_{nxn} \end{bmatrix} r = \begin{bmatrix} a_s + B^T b \\ -a_s + EC_s \end{bmatrix}, \quad (19)$$

where $G$ is the $2n$ by $2n$ symmetric matrix $L^T L$. We next show that the inverse of the matrix $G$ exists, and we determine it explicitly.

**LEMMA 1.**

**Result 1:** The inverse of the matrix $G$ given in Eq. (19) exists and is

$$G^{-1} = \begin{bmatrix} P & J \\ J & S \end{bmatrix}, \quad (20)$$

where

$$J = (D + E)^{-1} = (B^T B + VV^T)^{-1}, \quad (21)$$

$$P = J(I + E), \quad (22)$$

and

$$S = J(I + D), \quad (23)$$

**Result 2.**

$$SE = I - JD, \quad (24)$$

which is a property that we shall use for the determination of the “scaled” force of constraint, $F_s^c$.

**Proof.**

**Result 1.** For simplicity, let us write $G^{-1}$ as

$$G^{-1} = \begin{bmatrix} P & J \\ J & S \end{bmatrix}. \quad (25)$$

Beginning with the condition $G^{-1}G = I$, we obtain

$$\begin{bmatrix} P & J \\ J & S \end{bmatrix} \begin{bmatrix} I+D & -I \\ -I & I+E \end{bmatrix} = \begin{bmatrix} I & 0 \\ 0 & I \end{bmatrix},$$

which can be written as

$$\begin{bmatrix} P(I+D) - J & -P + J(I+E) \\ J(I+D) - S & -J + S(I+E) \end{bmatrix} = \begin{bmatrix} I & 0 \\ 0 & I \end{bmatrix}. \quad (26)$$

A comparison of the corresponding members on either side of the equality in Eq. (26) shows that
$$P(I+D)-J=I, \qquad (27)$$

$$S=J(I+D), \qquad (28)$$

$$P=J(I+E), \qquad (29)$$

and

$$-J+S(I+E)=I. \qquad (30)$$

Then replacing the matrix P obtained from Eq. (29) in Eq. (27) and simplifying that, we have

$$J(I+E)(I+D)-J=J(D+E+ED)=I. \qquad (31)$$

Since $E=VV^T$, $D=B^TB$, and $BV=0$, we have $ED=VV^TB^TB = V(BV)^TB=0$.

Thus Eq. (31) can be simplified to

$$J(D+E)=I. \qquad (32)$$

From Eqs. (17) and (18), it can be seen that the $n$ by $n$ matrix $D+E=B^TB+VV^T=[B^T V]_{VT}^B$. Since the matrix $[B^T V]$ has full rank, the rank of $[B^T V]_{VT}^B$ is $n$. Hence, $D+E$ has an inverse and from Eq. (32) the matrix $J$ is given by

$$J=(D+E)^{-1}=(B^TB+VV^T)^{-1}. \qquad (33)$$

By Eqs. (25), (28), and (29), the inverse of the matrix $G$ can be then written as

$$G^{-1} = \begin{bmatrix} P & J \\ J & S \end{bmatrix} = \begin{bmatrix} J(I+E) & J \\ J & J(I+D) \end{bmatrix}, \qquad (34)$$

where $J$ is given by Eq. (33). □

**Result 2:** By substituting Eq. (28) in Eq. (30), we obtain $-J+S(I+E)=-J+S+SE=-J+J(I+D)+SE=I$, which can be simplified to

$$SE=I-JD. \qquad (35)$$

From Eqs. (19) and (20), the vector $r=[\dot{q}_s]_{F_c}$, can be uniquely found as

$$\begin{bmatrix} \ddot{q}_s \\ F_s^c \end{bmatrix} = \begin{bmatrix} P & J \\ J & S \end{bmatrix} \begin{bmatrix} a_s + B^T b \\ -a_s + EC_s \end{bmatrix}. \qquad (36)$$

Using Eq. (36), the “scaled” force of constraint can be expanded as

$$F_s^c=Ja_s+JB^Tb-Sa_s+SEC_s. \qquad (37)$$

From Eqs. (23) and (24), Eq. (37) can be expressed as

$$\begin{aligned} F_s^c &= Ja_s + JB^T b - Ja_s - JDa_s + (I - JD)C_s \\ &= JB^T b - JDa_s + (I - JD)C_s. \end{aligned}$$

Noting that $D=B^TB$, the last equation gives a simple form for the constraint force

$$F_s^c=JB^T(b-Ba_s)+(I-JB^T B)C_s. \qquad (38)$$

Since the acceleration of the unconstrained system is defined as $a=M^{-1}F$, we have $Ba_s=(AM^{-1/2})(M^{-1/2}F)=A(M^{-1}F)=Aa$. Using this equality, and substituting $C_s$ by $M^{-1/2}C$ in Eq. (38), we get

$$\begin{aligned} F^c &= M^{1/2} F_s^c = M^{1/2} J B^T (b - Aa) + M^{1/2} (I - J B^T B) M^{-1/2} C \\ &:= F_i^c + F_{ni}^c, \end{aligned} \qquad (39)$$

which gives the force of constraint $F^c$ explicitly for the constrained system. The subscript $i$ is used to describe the force of constraint were all the constraints to be ideal ($C=0$); the subscript $ni$ is used to describe the contribution to the total constraint

force because of the nonideal nature of the constraints. The explicit equation of motion with nonideal constraints can then be written as

$$M\ddot{q}=F+F^c=F+M^{1/2}JB^T(b-Aa)+M^{1/2}(I-JB^T B)M^{-1/2}C. \qquad (40)$$

We emphasize that Eq. (40), which gives explicitly the motion of nonholonomic systems with nonideal constraint forces, does not involve any generalized inverses, or any Lagrange multipliers.

Previous investigators, so far as we know, have not obtained explicit equations of motion for non-ideal constraints. The only other general equation of motion for constrained mechanical systems with nonideal constraints available in the literature to date appears to be the one obtained in Refs. [11,12] and [15,16]. However, the results that have been obtained so far use the concept of the generalized inverse of a matrix, and the derivations are heavily dependent on the properties of generalized inverses. The equation obtained herein is: (1) explicit; (2) applicable to nonideal constraints; and (3) does not use generalized inverses. In the next section we shall compare our result with those obtained in Refs. [11,12].

There are, however, a number of formulations of the equations of motion for constrained mechanical systems under the assumption that the constraints are all ideal, i.e., when $C$ in Eq. (40) in identically zero for all time. It is then perhaps worthwhile comparing Eq. (40) for $C=0$, thereby restricting it to only ideal constraints, with formulations that have been obtained by previous investigators. So, to elucidate our equation further, we compare the form of the equation obtained by us with those obtained previously. Though Eq. (40) is also valid for nonideal constraints, in the next paragraph we restrict ourselves, for purposes of comparison with other formulations of the equations of motion obtained by other researchers only to when all the constraints are ideal.

Unlike the results obtained in Beghin [2], Chataev [3], Hamel [7], and Lagrange [8], Eq. (40) explicitly gives the force of constraint; no Lagrange multipliers are involved. The use of Lagrange multipliers constitutes one approach to solving the problem of constrained motion. We use in this paper a different approach that is innocent of this notion. These multipliers, which were invented by Lagrange, are an intermediary mathematical device for solving the problem of constrained motion. As such, they are not intrinsic (essential) to either the description of the physical problem of constrained motion or to the final equation of motion that is obtained, as witnessed by the fact that we make no mention of Lagrange multipliers in our approach. Another important point of difference is that the constraint equations we use to obtain Eq. (40) are more general than those in Appell [1], Beghin [2], Chataev [3], Gibbs [6], Hamel [7], and Synge [17] because the elements of the matrix $A$ are allowed to be not just functions of $a$ and $t$, but also of $\dot{q}$. This greatly expands the scope of the type of constraints that we use. However, it entails a more delicate interpretation of the concept of virtual displacements (see, Ref. [14]). Furthermore, unlike the formulations of Gibbs [6] and Appell [1] the coordinates we use to describe the constrained motion are the same as those used to describe the unconstrained motion; no quasi-coordinates are used, and no coordinate transformations are needed. Dirac [4] developed a set of equations for the constrained motion of hamiltonian systems in which the constraints are not explicitly dependent on time. Our equation differs from his in that: (1) Eq. (40) (with $C(t)=0$) is also applicable to non-hamiltonian, and dissipative systems, and (2) it allows constraints that contain time explicitly in them. However, Eq. (40) assumes that $M$ is positive definite, while Dirac's method can handle singular Lagrangians; such Lagrangians are more relevant to the field of quantum mechanics (for which Dirac developed his equation) and are seldom found in well-posed problems in classical mechanics.

One consequence of the fact that we use the same set of coordinates to describe the motion of the constrained system as we use to describe the unconstrained system is that our equation provides
the exceptional insight that the total force of constraint is the sum of two forces, as seen from the last two members on the right hand side of the last equality in Eq. (40). The first corresponds to what would result were all the constraints ideal; the second corresponds to the force caused solely by the nonideal nature of the constraints. Our ability to obtain the general equation of motion explicitly gives an additional insight when $C=0$. Nature appears to be acting like a “control engineer,” because the second term on the right-hand side of Eq. (40) may be viewed as a “feedback control force” proportional to the error, $(b-Aa)$, in the satisfaction of the constraint Eq. (3). We observe that the feedback “control gain matrix,” $M^{1/2}JB^T$, which nature uses turns out to be, in general, a highly nonlinear, time-dependent function of $q, \dot{q}$, and $t$. Such insights into the fundamental nature of constrained motion have been unavailable from previous formulations of the equations for constrained mechanical systems, such as those of Appell, Begin, Chataev, Hamel, Gibbs, Jacobi, Lagrange, and Synge.

### 3 Connection of Eq. (40) With Previous Results

In this section we show that the equation of motion obtained above is equivalent to the ones previously obtained in Refs. [11,12].

**LEMMA 2.**

$$JB^T = B^+, \quad (41)$$

where $B^+$ is the Moore-Penrose inverse of the matrix $B$.

**Proof.**

Let us consider a condition $GG^{-1} = G^{-1}G$,

$$\begin{bmatrix} I+D & -I \\ -I & I+E \end{bmatrix} \begin{bmatrix} P & J \\ J & S \end{bmatrix} = \begin{bmatrix} P & J \\ J & S \end{bmatrix} \begin{bmatrix} I+D & -I \\ -I & I+E \end{bmatrix},$$

which can be expanded to

$$\begin{bmatrix} (I+D)P - J & (I+D)J - S \\ -P + (I+E)J & -J + (I+E)S \end{bmatrix} = \begin{bmatrix} P(I+D) - J & -P + J(I+E) \\ J(I+D) - S & -J + S(I+E) \end{bmatrix}. \quad (42)$$

Equating the first element of the second column on either side of Eq. (42), we get

$$(I+D)J - S = -P + J(I+E). \quad (43)$$

After substituting Eqs. (28) and (29) in Eq. (43), we obtain

$$DJ = JD. \quad (44)$$

Similarly, equating the second element of the second column on either side of Eq. (42), we get $-J + (I+E)S = -J + S(I+E)$, which simplifies to

$$ES = SE. \quad (45)$$

As a result of Eqs. (44), (24), and (45), we have

$$DJ = JD = I - SE = I - ES. \quad (46)$$

To show that $JB^T$ is the Moore-Penrose (MP) inverse of the matrix $B$, we need to prove the following conditions:

1. $B(JB^T)B = B$;

2. $(JB^T)B(JB^T) = JB^T$;

3. $(BJB^T)^T = BJB^T$;

and

4. $(JB^T B)^T = JB^T B$.

1. By using the relations obtained from Eqs. (17), (46), and (18), we have $B(JB^T)B = BJD = B(I - ES) = B - BES = B - (BV)V^T S$. Since $BV=0$, $B(JB^T)B = B$. Thus the first MP condition is satisfied.

2. Due to Eqs. (17) and (46), $(JB^T)B(JB^T) = J(B^T B)JB^T = JDJB^T = J(I - ES)B^T = JB^T - JESB^T$. Since $SE = ES$, $E = VV^T$, and $BV=0$, $JESB^T = JSEB^T$ and $EB^T = VV^T B^T = V(BV)^T = 0$; thus $(JB^T)B(JB^T) = JB^T$, and the second MP condition is satisfied.

3. Since the matrices $D$ and $E$ are symmetric, $J = (D+E)^{-1}$, is also symmetric. Hence $(BJB^T)^T = BJ^T B^T = BJB^T$; thus the third MP condition is satisfied

4. Using Eqs. (17) and (44) we get $(JB^T B)^T = (B^T B)J^T = DJ = JD = J(B^T B)$; thus the fourth MP condition is satisfied.

From the result of lemma 2, after substituting $B^+ = JB^T$ in Eq. (39), we obtain

$$F^c = M^{1/2}B^+(b-Aa) + M^{1/2}(I-B^+B)M^{-1/2}C. \quad (47)$$

The first member on the right of Eq. (47) is the force of constraint that would be generated were all the constraints ideal, the second member gives the contribution to total force of constraint because of its non-ideal nature. Since $B = AM^{-1/2}$, Eq. (47) can be rewritten as

$$F^c = M^{1/2}(AM^{-1/2})^+(b-Aa) + M^{1/2}(I-B^+B)M^{-1/2}C. \quad (48)$$

From Eq. (4), we have $\ddot{q} = M^{-1}F + M^{-1}F^c = a + M^{-1}F^c$. Hence, the explicit equation of motion of the constrained system can be expressed as

$$\ddot{q} = a + M^{-1/2}(AM^{-1/2})^+(b - Aa) + M^{-1/2}(I - B^+B)M^{-1/2}C, \quad (49)$$

which is identical to the equation given by Udwadia and Kalaba (Refs. [11,12]). When $C=0$, the constraint forces are ideal and D'Alembert's principle is satisfied. Equation (49) then reduces to the result given in Refs. [9] and [13].

### 4 Examples

In this section, we provide examples that demonstrate the use of the equations of motion (40) for systems with nonideal constraints. The last example deals with a problem of sliding friction.

(a) Consider a particle of unit mass traveling in a three-dimensional configuration space with “given” forces $f_x(x,y,z,t)$, $f_y(x,y,z,t)$ and $f_z(x,y,z,t)$ and satisfying the nonholonomic constraint $\ddot{y} = z^2\dot{x} + \alpha g(x,t)$, where $\alpha$ is a constant and $g(x,t)$ is a given function of $x$ and $t$. The initial conditions are taken to be compatible with the nonholonomic constraint.

Since the mass of particle is unity, the unconstrained acceleration is given by

$$a = \begin{bmatrix} \ddot{x} \\ \ddot{y} \\ \ddot{z} \end{bmatrix} = \begin{bmatrix} f_x(x,y,z,t) \\ f_y(x,y,z,t) \\ f_z(x,y,z,t) \end{bmatrix}. \quad (50)$$

After differentiating the constraint equation with respect to time, we get

$$[-z^2 \quad 1 \quad 0] \begin{bmatrix} \ddot{x} \\ \ddot{y} \\ \ddot{z} \end{bmatrix} = 2z\dot{z}\dot{x} + \alpha g_x \dot{x} + \alpha g_t, \quad (51)$$
where $g_x$ and $g_t$ are partial derivatives of $g(x,t)$ with respect to $x$ and $t$, respectively. A comparison with Eqs. (3) provides us

$$A = [-z^2 \quad 1 \quad 0] \qquad (52)$$

and

$$b = 2z\dot{z}\dot{x} + \alpha g_x \dot{x} + \alpha g_t. \qquad (53)$$

Since $M=I_3$,

$$B = AM^{-1/2} = A. \qquad (54)$$

In addition, the solution vectors $v_1$ and $v_2$ to Eq. (12) are

$$V = [v_1 \ v_2] = \begin{bmatrix} 1 & 1 \\ z^2 & z^2 \\ k_1 & k_2 \end{bmatrix}, \qquad (55)$$

where, $k_1$ and $k_2$ are arbitrarily chosen, with $k_1 \neq k_2$, so that the column vectors $v_1$ and $v_2$ are linearly independent.

As previously shown in lemma 1, $J=(D+E)^{-1}=(B^T B + VV^T)^{-1}$. By Eqs. (54) and (55), we obtain (with $k_1 \neq k_2$)

$$J = \frac{1}{\Delta} \begin{bmatrix} k_1^2 + k_2^2 + z^4(k_1 - k_2)^2 & 2k_1k_2z^2 & -(k_1 + k_2)(z^4 + 1) \\ 2k_1k_2z^2 & z^4(k_1^2 + k_2^2) + (k_1 - k_2)^2 & -z^2(k_1 + k_2)(z^4 + 1) \\ -(k_1 + k_2)(z^4 + 1) & -z^2(k_1 + k_2)(z^4 + 1) & 2(z^4 + 1)^2 \end{bmatrix},$$

where $\Delta=(k_1-k_2)^2(z^4+1)^2$.

This gives

$$JB^T = \frac{1}{(z^4 + 1)} \begin{bmatrix} -z^2 \\ 1 \\ 0 \end{bmatrix}. \qquad (56)$$

We could have, of course, started by choosing, say, $k_1=1$ and $k_2=0$ in Eq. (55); we would then have arrived at relation (56) with much less algebra.

Suppose that the constraint force is nonideal and it does work under virtual displacements. Let us assume that the work done by the constraint force is given, for any virtual displacement, $w$, by

$$w^T F^c = -w^T a_0 (u^T u)^\beta (u/|u|), \qquad (57)$$

where $u=[\dot{x} \ \dot{y} \ \dot{z}]^T$ is the velocity of the particle, $|u|=\sqrt{u^T u}$, and $a_0$ and $\beta$ are constants. In this case, $C$ is a known 3-vector, and can be written as

$$C = -a_0(u^T u)^\beta (u/|u|) = -a_0(\dot{x}^2 + \dot{y}^2 + \dot{z}^2)^{\beta-1/2} \begin{bmatrix} \dot{x} \\ \dot{y} \\ \dot{z} \end{bmatrix}. \qquad (58)$$

After substituting Eqs. (50), and (52) through (58) in Eq. (39), we obtain

$$F^c = \left( \frac{2z\dot{z}\dot{x} + \alpha g_x \dot{x} + \alpha g_t + z^2 f_x - f_y}{z^4 + 1} \right) \begin{bmatrix} -z^2 \\ 1 \\ 0 \end{bmatrix} - a_0 \frac{(\dot{x}^2 + \dot{y}^2 + \dot{z}^2)^{\beta-1/2}}{z^4 + 1} \begin{bmatrix} \dot{x} + z^2 \dot{y} \\ z^2 \dot{x} + z^4 \dot{y} \\ \dot{z}(1 + z^4) \end{bmatrix}. \qquad (59)$$

From Eq. (40), the equation of motion of the constrained system is then

$$\begin{bmatrix} \ddot{x} \\ \ddot{y} \\ \ddot{z} \end{bmatrix} = \begin{bmatrix} f_x \\ f_y \\ f_z \end{bmatrix} + \left( \frac{2z\dot{z}\dot{x} + \alpha g_x \dot{x} + \alpha g_t + z^2 f_x - f_y}{z^4 + 1} \right) \begin{bmatrix} -z^2 \\ 1 \\ 0 \end{bmatrix} - a_0 \frac{(\dot{x}^2 + \dot{y}^2 + \dot{z}^2)^{\beta-1/2}}{z^4 + 1} \begin{bmatrix} \dot{x} + z^2 \dot{y} \\ z^2 \dot{x} + z^4 \dot{y} \\ \dot{z}(1 + z^4) \end{bmatrix}. \qquad (60)$$

The first member on the right-hand side of Eq. (60) is the impressed force. The second member is the constraint force that would be generated had the constraint been ideal, and the third member results from the nonideal nature of the constraint that is

described by Eq. (57). When $\alpha=0$, and $\beta=1$, the equation of motion (60) becomes identical to that given by Udwadia and Kalaba [11]. We note that here the result is obtained without any reference to generalized inverses.

(b) Consider a bead having a mass $m$. Suppose that it moves on a circular ring of radius $R$ as shown in Fig. 1. The motion can be described by the coordinates $(x,y)$. The gravitational acceleration, $g$, is downwards. We assume that the initial conditions on the motion of the bead are compatible with the constraint that it lie on the ring.

Were the bead not constrained to lie on the ring, its unconstrained acceleration would be

$$a = \begin{bmatrix} 0 \\ -g \end{bmatrix}. \qquad (61)$$

In this problem, the constraint equation is $x^2+y^2=R^2$. After differentiating the constraint equation twice, we obtain

$$[x \ y] \begin{bmatrix} \ddot{x} \\ \ddot{y} \end{bmatrix} = -\dot{x}^2 - \dot{y}^2, \qquad (62)$$

so that

$$A = [x \ y], \qquad (63)$$

<figure><img src="image_2.png" /><figcaption>Fig. 1 A bead of mass, $m$, moving on a circular ring of radius, $R$</figcaption></figure>
and

$$b = -\dot{x}^2 - \dot{y}^2. \quad (64)$$

Since the mass matrix $M = mI_2$,

$$B = AM^{-1/2} = m^{-1/2}[x \ y]. \quad (65)$$

For any virtual displacement $w \neq 0$ such that $Aw = 0$, we have $w = [-x]$ so that

$$V = M^{1/2}w = m^{1/2} \begin{bmatrix} y \\ -x \end{bmatrix}. \quad (66)$$

Using Eq. (21), (65) and (66), we obtain

$$ \begin{aligned} \mathcal{J}B^T &= (\mathcal{J}^T B + VV^T)^{-1} \mathcal{J}^T = \frac{1}{m^{3/2}R^4} \begin{bmatrix} x^2m^2 + y^2 & xy(m^2-1) \\ xy(m^2-1) & y^2m^2 + x^2 \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix} \\ &= \frac{m^{1/2}}{R^2} \begin{bmatrix} x \\ y \end{bmatrix}. \end{aligned} \quad (67) $$

Suppose that the nonideal constraint force, due to the rough surface of the ring, is given by

$$w^T F^c = -w^T \frac{h(x, \dot{x}, y, \dot{y}, t)}{\sqrt{\dot{x}^2 + \dot{y}^2}} \begin{bmatrix} \dot{x} \\ \dot{y} \end{bmatrix}, \quad (68)$$

for any virtual displacement $w$, where $h$ is a known function of $x$, $\dot{x}$, $y$, $\dot{y}$, and $t$.

From the calculation in Eq. (39), the force of constraint on the bead can be expressed as

$$F^c = - \frac{m(\dot{x}^2 + \dot{y}^2 - yg)}{R^2} \begin{bmatrix} x \\ y \end{bmatrix} + \frac{h(x, \dot{x}, y, \dot{y}, t)}{\sqrt{\dot{x}^2 + \dot{y}^2}} \cdot \frac{(x\dot{y} - y\dot{x})}{R^2} \begin{bmatrix} y \\ -x \end{bmatrix}. \quad (69)$$

Finally, by Eq. (40), the equation of motion of the constrained system is

$$ \begin{bmatrix} m\ddot{x} \\ m\ddot{y} \end{bmatrix} = \begin{bmatrix} 0 \\ -mg \end{bmatrix} - \frac{m(\dot{x}^2 + \dot{y}^2 - yg)}{R^2} \begin{bmatrix} x \\ y \end{bmatrix} + \frac{h(x, \dot{x}, y, \dot{y}, t)}{\sqrt{\dot{x}^2 + \dot{y}^2}} \cdot \frac{(x\dot{y} - y\dot{x})}{R^2} \begin{bmatrix} y \\ -x \end{bmatrix}. \quad (70) $$

The first member on the right-hand side of Eq. (70) is the given force acting on the unconstrained system; the second is the constraint force that would have been generated had the constraint been ideal; and, the last member accounts for the nonideal nature of the constraint.

(c) Consider a rigid block of mass $m$ sliding on an inclined plane that oscillates in the vertical direction with amplitude $\beta$ and frequency $\omega$, the coefficient of Coulomb friction between the plane and the surface of the block being $\mu$. See Fig. 2. We shall assume that the acceleration of the inclined plane is sufficiently small so that the block does not leave the surface of the plane as it moves under gravity.

In the absence of the inclined plane, the unconstrained equations of motion of the block of mass $m$ and under gravity can be written as

$$\begin{bmatrix} m & 0 \\ 0 & m \end{bmatrix} \begin{bmatrix} \ddot{x} \\ \ddot{y} \end{bmatrix} = \begin{bmatrix} 0 \\ mg \end{bmatrix}, \quad (71)$$

so that the acceleration, $a$, of the unconstrained system is given by $a = [0 \ g]^T$, and $M = mI_2$.

The unconstrained system is then subjected to the constraint, namely that the block must lie on the vibrating inclined plane. Hence, the constraint is given by the kinematic relation $y(t) = x(t)\tan \alpha - \beta \sin \omega t$, which can be expressed after differentiation with respect to time $t$ as

<figure><img src="image_3.png" /><figcaption>Fig. 2 A block sliding under gravity on an inclined plane ($0 < \alpha < \pi/2$) that is vibrating vertically with amplitude $\beta$ and frequency $\omega$. The coefficient of Coulomb friction between the plane and the block is $\mu$.</figcaption></figure>

$$[-\tan \alpha \ 1] \begin{bmatrix} \ddot{x} \\ \ddot{y} \end{bmatrix} = \beta \omega^2 \sin \omega t. \quad (72)$$

Thus, $A = [-\tan \alpha \ 1]$, and $b$ is the scalar $\beta\omega^2 \sin \omega t$. By Eq. (8), we have the virtual displacement

$$w = \delta \begin{bmatrix} 1 \\ \tan \alpha \end{bmatrix}, \quad (73)$$

where $\delta$ is any nonzero constant. Hence, we get

$$V = v = M^{1/2}w = m^{1/2} \delta \begin{bmatrix} 1 \\ \tan \alpha \end{bmatrix}. \quad (74)$$

Since $B = AM^{-1/2} = m^{-1/2}A$, using Eq. (21), we have

$$ \begin{aligned} \mathcal{J}B^T = (\mathcal{J}^T B + VV^T)^{-1} \mathcal{J}^T = & \left( m^{-1} \begin{bmatrix} -\tan \alpha \\ 1 \end{bmatrix} \begin{bmatrix} -\tan \alpha & 1 \end{bmatrix} \right. \\ & \left. + m \delta^2 \begin{bmatrix} 1 \\ \tan \alpha \end{bmatrix} \begin{bmatrix} 1 & \tan \alpha \end{bmatrix} \right)^{-1} \left( m^{-1/2} \begin{bmatrix} -\tan \alpha \\ 1 \end{bmatrix} \right), \end{aligned} $$

which can be simplified to

$$\mathcal{J}B^T = m^{1/2} \cos^2 \alpha \begin{bmatrix} -\tan \alpha \\ 1 \end{bmatrix}. \quad (75)$$

Therefore, the force of constraint, were the constraint to be ideal, would then be given by

$$F_i^c = M^{1/2} \mathcal{J} B^T (b - Aa) = m \cos^2 \alpha \begin{bmatrix} \tan \alpha \\ -1 \end{bmatrix} (g - \beta \omega^2 \sin \omega t). \quad (76)$$

In the presence of Coulomb friction, the magnitude of the frictional force is $\mu|F_i^c|$, where $|z| = +\sqrt{z^T z}$. We note that Coulomb's law of friction is an approximate empirical relation (see Ref. [18]). The relative velocity of the block with respect to the inclined plane is given by $\dot{q} = [\dot{x} \ \dot{x} \tan \alpha]^T$. The frictional force is in a direction opposite that of this relative velocity. The work done by Coulomb friction under a virtual displacement $w$ is then

$$W = -w^T \left( \mu |F_i^c| \frac{\dot{q}}{|\dot{q}|} \right), \quad (77)$$

so that

$$C = -\mu |F_i^c| \frac{\dot{q}}{|\dot{q}|} = \frac{-\mu |F_i^c|}{|\dot{x}| \sec \alpha} \begin{bmatrix} \dot{x} \\ \dot{x} \tan \alpha \end{bmatrix} = -\mu |F_i^c| \cos \alpha \begin{bmatrix} 1 \\ \tan \alpha \end{bmatrix} s, \quad (78)$$

where, $s = \operatorname{sgn}(\dot{x})$.

Relation (76) yields
$$|F_i^c| = m \cos \alpha |(g - \beta \omega^2 \sin \omega t)|. \quad (79)$$

The contribution to the total force of constraint generated by the non-ideal nature of the constraint is then

$$ 
\begin{align} 
F_{ni}^c &= M^{1/2}(I - JB^T B)M^{-1/2}C = m^{1/2} \left( I - m^{1/2} \cos^2 \alpha \begin{bmatrix} -\tan \alpha \\ 1 \end{bmatrix} \right. \nonumber \\
&\qquad \left. \times [-\tan \alpha \ 1] m^{-1/2} \right) m^{-1/2} \left( -\mu |F_i^c| \cos \alpha \begin{bmatrix} 1 \\ \tan \alpha \end{bmatrix} s \right) \nonumber \\
&= -\mu |F_i^c| \cos \alpha \begin{bmatrix} 1 \\ \tan \alpha \end{bmatrix} s \nonumber \\
&= -\mu m \cos^2 \alpha |(g - \beta \omega^2 \sin \omega t)| \begin{bmatrix} 1 \\ \tan \alpha \end{bmatrix} s. \tag{80}
\end{align} 
$$

Note that if the block is to remain in contact with the plane we require $(g - \beta \omega^2 \sin \omega t) \ge 0$. The equation of motion of the block sliding on the plane, by Eqs. (40), is then

$$ 
\begin{bmatrix} m\ddot{x} \\ m\ddot{y} \end{bmatrix} = \begin{bmatrix} 0 \\ mg \end{bmatrix} + m \cos^2 \alpha \begin{bmatrix} \tan \alpha \\ -1 \end{bmatrix} (g - \beta \omega^2 \sin \omega t) \\ - \mu m \cos^2 \alpha \begin{bmatrix} 1 \\ \tan \alpha \end{bmatrix} (g - \beta \omega^2 \sin \omega t)s. \quad (81) 
$$

We note that each of the three members on the right-hand side of Eq. (81) has a simple interpretation. And it is precisely to expose this essential simplicity with which nature seems to operate that we have desisted from simplifying the equation any further. For the first member is the “given” force; the second is the force of constraint were the constraint to be ideal; and, the third is the constraint force engendered by the nonideal nature of the constraint.

## 5 Conclusions

The explicit equations of motion for holonomic and nonholonomic mechanical systems with nonideal constraints have so far been obtained in terms of generalized inverses of matrices. These inverses were first proposed by Moore [19], and their properties were first extensively developed by Penrose [20]. Since the properties of generalized inverses have appeared to be essential in developing these explicit equations, it had been felt that it was because of their relatively recent introduction—in the 1950s—in the scientific literature that the explicit equations of motion for nonholonomic mechanical systems were unavailable until quite recently (see Refs. [11–13], and [15,16]).

In this paper we show that this line of reasoning does not appear to be correct. Rather than reliance on generalized inverses of matrices and their properties, what we may have needed to get the explicit equations of motion is a more refined understanding of the problem, and a further development of concepts that have long since been with us. Among these are: (1) a proper conceptualization of the problem of constrained motion in terms of an unconstrained system, which is then subjected to the imposed constraints, (2) the generalized concept of a virtual displacement vector, described in Ref. [14] and, (3) the use of linear algebra. It is somewhat surprising that though the equations of motion that govern even some of the simplest constrained mechanical systems are nonlinear, it is linear algebra that plays a central role in their development.

We point out that the explicit equations of motion obtained herein, like those obtained earlier (Refs. [11–16]), are completely

innocent of the notion of Lagrange multipliers. Over the last 200 years, Lagrange multipliers have been so widely used in the development of the equations of motion of constrained mechanical systems that it is sometimes tempting to mistakenly believe that they possess an intrinsic presence in the description of constrained motion. This is not true. As shown in this paper, neither in the formulation of the physical problem of the motion of constrained mechanical systems nor in the equations governing their motion are any Lagrange multipliers involved. The use of Lagrange multipliers (a mathematical tool invented by Lagrange [8]) constitutes only one of the several intermediary mathematical devices invented for handling constraints. And, in fact, the direct use of this device appears inapplicable when the constraints are functionally dependent. Lagrange multipliers do not appear in the physical description of constrained motion, and therefore cannot, and do not, ultimately appear in the equations governing such motion.

The explicit equations of motion obtained in this paper apply to general, holonomic, and nonholonomic systems that may have nonideal constraint forces. These constraint forces may, in general, do positive, zero, or negative work under virtual displacements at any time during the motion of the system. The equations given here are the first of their kind that are explicit, and that do not require the use of any generalized inverses, nor use of any of their properties.

## References

[1] Appell, P, 1899, "Sur une Forme Generale des Equations de la Dynamique," C R Acad Sci, Paris, 129, pp 459-460

[2] Beghin, H, 1967, *Mecanique Theorique et Appliquee*, Gauthier-Villars, Paris

[3] Chetaev, N G, 1989, *Theoretical Mechanics*, Mir Publisher, Moscow

[4] Dirac, P A M, 1964, *Lectures in Quantum Mechanics*, Yeshiva University, New York

[5] Gauss, C F, 1829, "Über Ein Neues Allgemeines Grundgsetz der Mechanik," J Reine Angew Math, 4, pp 232-235

[6] Gibbs, J W, 1879, "On the Fundamental Formulas of Dynamics," Am J Math, 2, pp 49-64

[7] Hamel, P G, 1949, *Theoretische Mechanik: Eine Einheitliche Einfuhrung in die Gesamte Mechanik*, Springer-Verlag, Berlin

[8] Lagrange, J L, 1787, *Mechanique Analytique*, Mme Ve Courcier, Paris

[9] Udwadia, F E, and Kalaba, R E, 1992, "A New Perspective on Constrained Motion," Proc R Soc London, Ser A, 439, pp 407-410

[10] Pars, L A, 1979, *A Treatise on Analytical Dynamics*, Oxbow Press, Woodbridge, CT, p 14

[11] Udwadia, F E, and Kalaba, R E, 2001, "Explicit Equations of Motion for Mechanical Systems With Nonideal Constraints," ASME J Appl Mech, 68, pp 462-467

[12] Udwadia, F E, and Kalaba, R E, 2002, "On the Foundations of Analytical Dynamics," Int J Non-Linear Mech, 37, pp 1079-1090

[13] Udwadia, F E, and Kalaba, R E, 1996, *Analytical Dynamics: A New Approach*, Cambridge University Press, Cambridge, UK

[14] Udwadia, F E, Kalaba, R E, and Hee-Chang, E, 1997, "Equation of Motion for Constrained Mechanical Systems and The Extended D'Alembert's Principle," Q Appl Math, LV(2), pp 321-331

[15] Udwadia, F E, 2000, "Fundamental Principles of Lagrangian Dynamics: Mechanical Systems With Non-Ideal, Holonomic, and Non-Holonomic Constraints," J Math Anal Appl, 252, pp 341-355

[16] Udwadia, F E, and Kalaba, R E, 2002, "What is the General Form of the Explicit Equations of Motion for Constrained Mechanical Systems?" ASME J Appl Mech, 69, pp 335-339

[17] Synge, J L, 1927, "On the Geometry of Dynamics," Philos Trans R Soc London, Ser A, 226, pp 31-106

[18] Bowden, F P, and Tabor, D, 1974, *Friction: An Introduction to Tribology*, Heinemann, London

[19] Moore, E H, 1920, "On the Reciprocal of The General Algebraic Matrix," Bull Am Math Soc, 26, pp 394-395

[20] Penrose, R, 1955, "A Generalized Inverse of Matrices," Proc Cambridge Philos Soc, 51, pp 406-413