The Rolling Wheel Equations without Magic: A Combined Slip-Stiction Approach by the Udwadia-Kalaba Formulation with Constraints Relaxation for a Smooth Simulation

Bukoko Ikoki

KDYNAMICS Logiciel Inc., Sherbrooke, Canada.

DOI: 10.4236/jtts.2021.113024                

Abstract

The Udwadia-Kalaba formulation is proposed to model the dynamics of the rolling wheel. A unified approach that addresses both the slip and the stiction in the contact section is considered. Purely rolling constraints are associated with stiction and are suitably lifted as slip occurs. An extended formulation for the Uwadia-Kalaba equations of motion is introduced for that matter. It resorts to the weighted minimum norm and the weighted semi-least-squares solutions of the constraints equations. This not only allows a bias on constraints, by an appropriate description of weight functions based on friction, it also leads to a smooth activation or deactivation of selected constraints without rewriting the equations of motion or upsetting their numerical integration.

Keywords

Rolling Wheel, Friction, Stiction, Slip, Constraints Relaxation

People were aiming cannonballs… for thousands of years before Newton, but without Newton’s laws, we would not have been able to go to the moon [1].

Udwadia

1. Introduction

In 1992, Udwadia and Kalaba [2] proposed a straightforward formulation of the equations of motion for constrained mechanisms, without resorting to Lagrangian multipliers. This formulation has the highly appreciable advantage of handling constraints with a Jacobian that is not necessarily regular. This makes it possible to deal with systems with vanishing constraints. With this formulation, such constraints are indifferently holonomic or non-holonomic. They are ideal or non-ideal, with the possibility of integrating friction, according to Udwadia and Kalaba.

These interesting features of the Udwadia-Kalaba approach perfectly relate to the problem of the rolling wheel when realistic situations such as slipping or braking are encountered. With this approach, realistic modeling of ground vehicles can be envisioned. Among other applications, thanks to an accurately modeled target vehicle, an improved automated HIL test procedure for different controllers can be achieved, for both performance and safety goals. While the otherwise involved simulations are smoothly performed without upsetting the numerical processes with the ability to handle both transient and stationary dynamics. Thanks to the provided explicit form of the equations of motion which do not require to be rewritten when constraints vanish [3] and the physically-inspired model of ground-wheel interaction [4]. So far, to the best of our knowledge, the tools of analytical dynamics [3] have been limited in their ability to bring this type of developments.

First, let’s introduce our notion of friction as we intend to apply it to the current work from a quasi-static point of view. Considering a block of mass m depicted in Figure 1, as it is pulled by a traction force $T=T{u}_1$ on a horizontal plane, the block remains in its condition of rest as long as the applied traction force $T$ is less than or equal in magnitude to the horizontal component $F_c=F_c{u}_1$ of the contact reaction force.

Figure 1. The sliding block.

We designate such a condition of rest under the pulling force as stiction. We extend it to the condition of motion in which the pulling force is exactly balanced by the contact reaction force in a uniform movement. Slip occurs when the traction force overtakes the contact reaction force. Such a contact reaction force in both cases of stiction and slip can suitably be represented by the following expression:

$F_c=-\mu N_c\tanh \left(\frac{T}{\mu N_c}\right)$ (1)

Indeed:

$\begin{array}{l}\text{For smaller traction forces,}T\lesssim \frac{1}{3}\mu N_c:\tanh \left(\frac{T}{\mu N_c}\right)\simeq \left(\frac{T}{\mu N_c}\right)\Rightarrow F_c=-T\\ \text{For greater traction forces,}T\gtrsim 3\mu N_c:\tanh \left(\frac{T}{\mu N_c}\right)\simeq 1\Rightarrow F_c=-\mu N_c\end{array}$

Expression (1) translates into an equation the fact that independently of what actually happens in the contact zone, the contact force is confined between two identifiable values. The applied force and the saturation limit.

2. The Equations of Motion of Constrained Systems

Udwadia and Kalaba [2] provide the explicit equations of motion of constrained systems by the following formulation:

$\{\begin{array}{l}M\left(q\right)\ddot{q}=Q+Q^c(\text{2a})\\ A\ddot{q}=b(\; 2b\; )\end{array}$

M is the n-order positive definite mass matrix, and $q$ is the n dimensional generalized coordinates vector. Matrix A with dimension $m\times n$ and m dimensional vector $b$ are obtained from the m constraints functions of geometric (3a) or kinematical (3b) kinds. The later might be indifferently scleronomic, rheonomic, holonomic or non-holonomic in this approach. Scleronomic constraints yield:

$\{\begin{array}{l}{h}^g\left(q\right)=0(\text{3a})\\ h^{\kappa }\left(q,\stackrel{\dot{}}{q}\right)=0\end{array}$

A and $b$ in (2b) are given by

$A=\{\begin{array}{l}\partial h^g\left(q\right)/\partial q\\ \partial h^{\kappa }\left(q,\stackrel{\dot{}}{q}\right)/\partial \stackrel{\dot{}}{q}\end{array},b=-\frac{\partial A}{\partial t}\stackrel{\dot{}}{q}$

As for Q in (2a), it is the matrix of generalized applied forces, conservative or dissipative forces and other complementary inertial forces which do not depend on $\ddot{q}$ . $Q^c$ represents the reaction forces that are needed to fulfill the constraints. They are of two natures: ideal, $Q_i^c$ , and non-ideal, $Q_{ni}^c$ . Force $Q^c=Q_i^c+Q_{ni}^c$ is added to the unconstrained system. In the presence of non-ideal constraints, we have

$Q_i^c=M{A}_M^{+}\left(b-Aa\right)$ (4)

Vector $a$ is the unconstrained acceleration and is given by $a=M^{-1}Q$ .

The non-ideal contribution is

$Q_{ni}^c=M\left(I-A_M^{+}A\right)M^{-1}C$ (5)

$A_M^{+}$ is the minimum norm weighted Moore-Penrose inverse of the constraints matrix A. It can provide the weighted minimum-norm and the least-squares solution to the constraints Equations (2b). A least-squares solution of a $m\times n$ system

$Ax=b$ (6)

is defined as a vector that minimizes the Euclidean norm of the residual $r=Ax-b$ :

$\underset{x}{\text{min}}{\Vert Ax-b\Vert }^2$ (7)

On the compatibility condition [5], $A{A}^{+}b=b$ , the general solution to (6) is:

$x=A^{+}b+\left(I-A^{+}A\right)z$ (8)

$z$ is an arbitrary vector. Applied to $\left(I-A^{+}A\right)$ , it is projected onto the null space of A. Hence any such vector $z$ would always have (7) fulfilled. Among all least-squares solutions provided by (8), the minimum-norm solution is

$x=A^{+}b$ (9)

For the minimization of the residual $r$ , weights might be considered. The weighted least-squares and weighted minimum-norm problem is then:

$\underset{x}{\text{min}}{\Vert Ax-b\Vert }_N^2$ (10)

with a sought weighted minimum-norm ${\Vert x\Vert }_M={\left(x^{\text{T}}Mx\right)}^{\frac{1}{2}}=\Vert M^{\frac{1}{2}}x\Vert$ while the weighted least-squares is expressed by ${\Vert Ax-b\Vert }_N^2={\left(Ax-b\right)}^{\text{T}}N\left(Ax-b\right)$ . Problem (10) can be reduced to its unweighted equivalent by setting [6]

$\tilde{x}=M^{\frac{1}{2}}x\mathrm{,}\tilde{b}=N^{\frac{1}{2}}b\mathrm{,}\tilde{z}=M^{\frac{1}{2}}z\mathrm{,}\tilde{A}=N^{\frac{1}{2}}A{M}^{-\frac{1}{2}}$

By replacing in (8), the general solution to the weighted least-squares and weighted minimum-norm problem (10) is then

$x=A_{NM}^{+}b+\left(I-A_{NM}^{+}A\right)z$ (11)

The M minimum-norm and N least-squares solution is: $x=Xb=A_{NM}^{+}b$ . Where X is the weighted Moore-Penrose generalized inverse which is defined by the 4 characteristic properties:

$AXA=A\mathrm{,}XAX=X\mathrm{,}{\left(NAX\right)}^{\text{T}}=NAX\text{and}{\left(XAM\right)}^{\text{T}}=XAM$ (12)

For the weighted MP inverse with M and N as positive definite, X is given explicitly by:

$A_{NM}^{+}=M^{-\frac{1}{2}}{\left(N^{\frac{1}{2}}A{M}^{-\frac{1}{2}}\right)}^{+}{N}^{\frac{1}{2}}$ (13)

Following his work on errors in mathematics, Gauss [7] pointed out a perfect analogy between the least-squares problem that he pioneered and the dynamics of constrained systems [8]. Compared with the unconstrained motion with an acceleration $a$ , he suggested in his least-constraint principle [9] that the actual motion of a constrained system is the one that minimizes the function:

$G={\left(\ddot{q}-a\right)}^{\text{T}}M\left(\ddot{q}-a\right)$ (14)

subjected to some constraints $h\left(q\mathrm{,}\stackrel{\dot{}}{q}\mathrm{,}t\right)=0$ . Which after an appropriate number of time differentiations generally result in a linear form identical to (2b):

$A\ddot{q}=b$

The least-constraint principle clearly amounts to a least-square problem (2b) for a weighted minimum-norm solution ${\Vert \ddot{q}\Vert }_M$ from (14).

Building on Gauss’ least-constraint principle, in order to obtain the Udwadia-Kalaba explicit form of the equations of motion [2] subjected to the constraints (2b), we set $\widetilde{\stackrel{¨}{q}}=\tilde{d}+\tilde{a}$ . Constraints (2b) can then be written as

$\tilde{A}\left(\tilde{d}+\tilde{a}\right)=\tilde{b}\iff \tilde{A}\tilde{d}=\tilde{b}-\tilde{A}\tilde{a}$ (15)

The least squares solution to (15) leads to

$d=X\left(b-Aa\right)+\left(I-XA\right)z$ (16)

Since $\ddot{q}=d+a$ , we have $M\ddot{q}=M\left(d+a\right)=Ma+Md$ . Knowing that $a=M^{-1}Q$ , by (16), we obtain the equations of motion:

$M\ddot{q}=Q+MX\left(b-Aa\right)+M\left(I-XA\right)z$ (17)

Unless $z=0$ , the general solution obtained in (16) is not of minimal norm though. Nonetheless, according to Udwadia and Kalaba [9], (17) can be used to address the case of non-ideal constraints in a comprehensive formulation. Noting that vector $z$ is an acceleration, it can conveniently be represented by $M^{-1}C$ . Where C is a n dimensional force vector which describes the non-ideal contribution.

By examining (17), it clearly appears that an additional term to the ideal constraints matrix $Q_i^c$ , would be a major pitfall to the Gauss’ least-constraint principle which unequivocally requires a least-squares and weighted minimum-norm solution. As a matter of course, an alteration of such a minimum norm, with an unfit $\ddot{q}$ as a result, does not render the true dynamics of the system while nevertheless faithfully allowing the fulfilling of the constraints.

The alteration of the motion constraints force by the non-ideal term, in our view, must be related to a relief of such constraints to a degree representative of their violation. As a solution to the problem, instead of a classical least-squares solution to constraints equations, we consider the seminorm weighted least-squares solution to such equations. This allows the mitigation of associated constraints equations whenever the second term in (17) comes into play. A different generalized inverse is chosen:

$X=M^{-\frac{1}{2}}{\left({\left(N^{\frac{1}{2}}\right)}^{+}A{M}^{-\frac{1}{2}}\right)}^{+}{N}^{\frac{1}{2}}$ (18)

N and M are hermitian nonegative definite. The seminorms are defined as ${\Vert x\Vert }_M={\left(x^{\text{T}}Mx\right)}^{\frac{1}{2}}$ and ${\Vert b\Vert }_N={\left(b^{\text{T}}Nb\right)}^{\frac{1}{2}}$ . We can verify that $r\left(A^{\text{T}}XA\right)=r\left(A\right)$ and also

$NAXA=NA,{\left(AX\right)}^{\text{T}}N=NAX,MXAX=MX,{\left(XA\right)}^{\text{T}}M=MXA$

Hence, according to [6], X is a minimum M-seminorm, N-semileast squares inverse of A. If M is positive definite, X is shown to be unique.

3. The Rolling Wheel

A driving torque $T_m$ is applied to a wheel of mass m and moment of inertia $I_2$ . The wheel is rolling on a horizontal plane with an angular velocity $\stackrel{\dot{}}{\theta }$ , as shown in Figure 2. We would like to get the explicit equations of motion with the constraints reactions. We consider the case of a rigid wheel presented by Popp and Schiehlen in ( [10]: p. 116). The generalization of such a wheel [4] allows the extension of the equations to the elastic wheel.

The required matrices to compute the equations of motion are:

$M=\left[\begin{array}{ccc}m& 0& 0\\ 0& m& 0\\ 0& 0& I_3\end{array}\right],Q=\left[\begin{array}{c}{F}_w\\ -mg\\ T_m+T_b\end{array}\right],C=\left[\begin{array}{c}{F}_c\\ 0\\ 0\end{array}\right],N=\left[\begin{array}{cc}1& 0\\ 0& s^2\end{array}\right]$

The generalized coordinates are given by the vector $q={\left[\begin{array}{ccc}x& z& \theta \end{array}\right]}^{\text{T}}$

Constraints equations

The speed of contact point $x_c^K$ relative to the non-rotating reference frame K on the wheel is expressed with respect to the global frame I:

${\stackrel{\dot{}}{x}}_c^I={\stackrel{\dot{}}{x}}_{o_1}^I+{\tilde{\omega }}_{IK}{x}_c^K={\stackrel{\dot{}}{x}}_{o_1}+{\stackrel{\dot{}}{R}}_w{R}_w^{\text{T}}{x}_c^K$ (19)

$R_w=\left[\begin{array}{ccc}\cos \theta & 0& \sin \theta \\ 0& 1& 0\\ -\sin \theta & 0& \cos \theta \end{array}\right],x_{o_1}=\left[\begin{array}{c}x\\ 0\\ z\end{array}\right],x_c^K=\left[\begin{array}{c}0\\ 0\\ -R\end{array}\right]$

The pure rolling condition is obtained by equating (19) to zero.

$\left[\begin{array}{c}\stackrel{\dot{}}{x}\\ 0\\ \stackrel{\dot{}}{z}\end{array}\right]+\stackrel{\dot{}}{\theta }\left[\begin{array}{ccc}0& 0& 1\\ 0& 0& 0\\ -1& 0& 0\end{array}\right]\left[\begin{array}{c}0\\ 0\\ -R\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\end{array}\right]$

A is obtained from the equations of constraints:

$\begin{array}{ccc}\{\begin{array}{l}{h}_1\left(q\mathrm{,}\stackrel{\dot{}}{q}\right)=0\equiv \stackrel{\dot{}}{z}=0\hfill \\ h_2\left(q\mathrm{,}\stackrel{\dot{}}{q}\right)=0\equiv \stackrel{\dot{}}{x}-\stackrel{\dot{}}{\theta }R=0\hfill \end{array}& \Rightarrow & \{\begin{array}{l}\frac{\partial }{\partial t}{h}_1\left(q\mathrm{,}\stackrel{\dot{}}{q}\right)=0\equiv \ddot{z}=0\hfill \\ \frac{\partial }{\partial t}{h}_2\left(q\mathrm{,}\stackrel{\dot{}}{q}\right)=0\equiv \ddot{x}-\ddot{\theta }R=0\hfill \end{array}\end{array}$

Figure 2. The rolling wheel.

In a matrix form, we obtain:

$A\ddot{q}=\left[\begin{array}{ccc}0& 1& 0\\ 1& 0& -R\end{array}\right]\left[\begin{array}{c}\ddot{x}\\ \ddot{z}\\ \ddot{\theta }\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\end{array}\right]=b$

$X=\left[\begin{array}{cc}0& \frac{I_2{s}^2}{m{R}^2+I_2}\\ 1& 0\\ 0& -\frac{mR{s}^2}{m{R}^2+I_2}\end{array}\right]$

For $F_w=0$ and $T_b=0$ , both the ideal and non-ideal constraints reaction forces are then respectively obtained from (4) and (5) as

$Q_i^c=\left[\begin{array}{c}\frac{R{s}^2}{\left(\frac{I_2}{m}+R^2\right)}{T}_m\\ mg\\ -\frac{R^2{s}^2}{\left(\frac{I_2}{m}+R^2\right)}{T}_m\end{array}\right]\mathrm{,}{Q}_{ni}^c=\left[\begin{array}{c}\frac{R^2+m^{-1}{I}_2\left(1-s^2\right)}{\left(\frac{I_2}{m}+R^2\right)}{F}_c\\ 0\\ \frac{m^{-1}{I}_{2}R{s}^2}{\left(\frac{I_2}{m}+R^2\right)}{F}_c\end{array}\right]$ (20)

$Q_i^c\left(1\right)$ is the stiction constraint force on the ground that allows the rolling without slipping when fully enforced. The stiction parameter s alters such a force while allowing the relief of the associated constraint as required. It is given by:

$s=1-{\tanh }^2\left(\frac{1}{3}\frac{k_s{T}_m}{\mu N_{c}R}\right)$ (21)

$Q_i^c\left(2\right)=mg=N_c$ is the reaction to the force of gravity that maintains the center of wheel on the required constant distance R from the ground. $Q_i^c\left(3\right)$ is the required torque for the rolling constraint. The friction force $F_c$ is written as

$F_c=-\mu N_c\tanh \left(\frac{1}{3}\frac{k_F{T}_m}{\mu N_{c}R}\right)\tanh \left(R\stackrel{\dot{}}{\theta }-\stackrel{\dot{}}{x}\right)$ (22)

We set $k_F=k_s=1$ . Different values can be used to fit different contact types and situations. For small values of the stiction parameter s, which correspond to severe slip as $T_m\gg \mu N_{c}R$ by (21) and $F_c\simeq -\mu N_c$ by (22), we have from (20):

$\underset{s\to 0}{lim}{Q}_i^c=\left[\begin{array}{c}0\\ mg\\ 0\end{array}\right]\mathrm{,}{Q}_{ni}^c=\left[\begin{array}{c}-\mu N_c\\ 0\\ 0\end{array}\right]$

When $T_m\ll \mu N_{c}R$ , either because applied torque $T_m$ is relatively small or the reaction of the ground $N_c$ is particularly important, then $F_c\simeq 0$ . And $s\simeq 1$ according to (21). The constraint forces $Q_i^c$ , $Q_{ni}^c$ are then given by:

$Q_i^c\simeq \left[\begin{array}{c}\frac{R}{\left(\frac{I_2}{m}+R^2\right)}{T}_m\\ mg\\ -\frac{R^2}{\left(\frac{I_2}{m}+R^2\right)}{T}_m\end{array}\right]\mathrm{,}{Q}_{ni}^c\simeq \left[\begin{array}{c}0\\ 0\\ 0\end{array}\right]$

Smaller applied torques are associated with stiction. The friction contact reaction force is negligible and the traction force is given by $Q_i^c$ . These results are consistent with those¹ presented in [10] for the rigid wheel.

3.1. Simulation Results

3.1.1. Driving Mode

A wheel of mass $m=25\text{kg}$ and radius $R=0.3\text{m}$ receives a gradually increasing torque on an occasionally slippery ground with $0.3\le \mu \le 0.9$ according to

$T_m\left(t\right)=T_{\max }\cdot \tanh \left(0.5t\right)$ (23)

$\mu \left(t\right)={\mu }_l+\left({\mu }_u-{\mu }_l\right)\cdot \tanh {\left(\frac{1}{s_d}\left(t-t_s\right)\right)}^4$ (24)

A comparison is made between the velocities for the dry phase of the ground and for the slippery one. For a relatively moderate torque with $T_{\max }=20\text{N}\cdot \text{m}$ , for a rolling resistance of $F_w=-1.7{\stackrel{\dot{}}{q}}_1^2$ and a bearing viscous damping effort of $T_b=-0.34{\stackrel{\dot{}}{q}}_3$ , as show in Figure 3, despite the low value of $\mu \left(t\right)$ upon the supposed slippery phase, the wheel shows a stiction behavior since the applied tractive force $T_m/R$ remains lower than the saturation limit $\mu N_c$ .

Figure 3. Rolling with stiction for $T_{\max }=20\text{N}\cdot \text{m}$ .

For a much higher torque, with $T_{\max }=60\text{N}\cdot \text{m}$ , as the wheel enters the slippery phase, Figure 4 shows that it spins noticeably and longitudinally slows down, before accelerating as it regains traction.

With $T_{\max }=220\text{N}\cdot \text{m}$ , Figure 5 shows the effective traction force given by the constraint force for rolling $Q_{ci}$ . The applied traction force $T_m/R$ is also represented against the saturation contact force $\mu N_c$ . As the applied traction overcomes the saturation limit, the contact traction force $Q_{ci}$ is saturated and subsequently lowered to near zero as slip occurs. Thus, fully relaxing the pure rolling constraint.

Experimental tests have to be conducted to determine parameters $k_s$ and $k_F$ that can be associated with a range of different contact behaviors in terms of performance and response. Though the modeling of the rolling wheel presented in this example is about the rigid wheel, it also applies for both the elastic wheel, for

Figure 4. Speed variation and constraint relaxation for $T_{\max }=60\text{N}\cdot \text{m}$ .

Figure 5. Traction saturation and constraint relaxation for $T_{\max }=220\text{N}\cdot \text{m}$ .

longitudinal dynamics and pneumatic tire for both the longitudinal dynamics and the lateral dynamics according to the works of Fufaev and Neĭmark [4].

3.1.2. Braking Mode

A braking torque $T_b$ is concurrently applied with the driving torque $T_m$ on a suddenly slippery ground. For different ground conditions, the longitudinal velocity and the braking distance are observed with the braking force. The braking force and the ground condition are expressed as:

$T_b\left(t\right)=T_{bm}\left(1+\tanh \left(t-t_b\right)\right)\tanh \left(k_v\stackrel{\dot{}}{\theta }\right),k_v\gg 1$ (25)

$\mu \left(t\right)=\frac{1}{2}-\frac{1}{2}\left({\mu }_u-{\mu }_l\right)\tanh \left(k_m\left(t-t_s\right)\right),k_m\gg 1$ (26)

For $T_m=50\text{Nm}$ , $T_{bm}=15\text{s}$ , $t_b=15$ , $k_v=100$ , $k_m=10$ , ${\mu }_l=0.3$ and ${\mu }_u=0.9$ , results are presented in Figure 6.

And then:

For ${\mu }_l=0.1$ , Figure 7 shows the results.

On both Figure 6 and Figure 7, we can observe, as $\mu \left(t\right)$ lowers, the wheel is

Figure 6. Braking for ${\mu }_l=0.3$ .

Figure 7. Braking for ${\mu }_l=0.1$ .

first spinning before braking. For ${\mu }_l=0.3$ , in Figure 6, the braking time is approximately the same as for the reference case with constant friction coefficient $\mu \left(t\right)=0.9$ . At a lower friction coefficient, ${\mu }_l=0.1$ , in Figure 7, the braking time significantly increases as the wheel slides with an almost locked angular velocity.

4. Conclusions

The complex wheel-and-ground interaction is modeled by considering the rolling constraint. Rolling without slipping is associated with stiction. This constraint needs to be lifted as slip kicks in since a classic variational approach cannot handle dry friction. To a certain extent, the Udwadia-Kalaba formulation for non-ideal constraints allows the consideration of friction. However, by considering the general solution of the least-squares problem associated with constraints equations, it fails to preserve the necessary minimum norm. According to Gauss’ least-constraint principle, the accelerations that are found in this way are not the actual ones. As they remain bound to the constraints equations that they persistently fulfill in a least-square sense, it appears that they do not obey the dynamics of the rolling. Nevertheless, the matrix-based Udwadia-Kalaba formulation appeals in its ability to not only allow the interpretation of the physics behind the equations; conversely, it also allows the transposition of an idea into equations. We have suggested a way to relax ideal constraints whenever the additional term of non-ideal constraints comes into play. The weighted semi-least-squares and weighted minimum-norm problem is considered instead of the classic approach by Udwadia and Kalaba. With friction-aware weights on the constraints equations, the residuals are altered accordingly. The bias on constraints is achieved to the extent slip is important by diverting norm minimization efforts to constraints with smaller residuals. By seamlessly setting the computed generalized coordinates free from selected constraints equations as needed, the true dynamics of the system is allowed to take place. A smooth simulation of the rolling is thus achieved. Without needing to rewrite the equations of motion or resort to tricky switching strategies that significantly upset the numerical integration process.

In their presentation of the joints modeling in flexible multibody systems ( [11]: p. 173), Cardona and Géradin address the rolling of the elastic wheel with a slip-stiction approach involving a regularization function proposed by Oden and Martins in ( [12]: p. 587). The authors describe such function as problematic with regard to the numerical integration when a perfect zero slipping situation is encountered. However, the retained approach is considered for its simplicity compared to schemes involving constraints activation and deactivation which are regarded as complicated in terms of the time integration procedure. As, according to the authors, they cause violent oscillations of constraints and velocities in the transition phases from the sliding to the stiction situations.

Through the relaxation of constraints, we have presented a simple method for the simulation of a rolling wheel, which amounts to a smooth activation and deactivation of such constraints, without the inconveniences experienced otherwise. The extended Udwadia-Kalaba equations of motion and the matrix of constraints relaxation we have introduced, both add a new feature to the multibody dynamics formalism. Such a contribution equally applies to the treatment of a variety of multibody systems that involve intermittent constraints.

NOTES

¹Caution must be exercised though with the sign convention of the authors.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Caverley, R. (2001) Constrained or Unconstrained, That Is the Equation. USC News.
[2]	Udwadia, F.E. and Kalaba, R. (1992) A New Perspective on Constrained Motion. Proceedings of the Royal Society, 439, 407-410.
[3]	Soltakhanov, S.K., Yushkov, M.P. and Zeghzda, S. (2009) Mechanics of Nonholonomic Systems. Springer, Berlin Heidelberg. https://doi.org/10.1007/978-3-540-85847-8
[4]	Neimark, J.I. and Fufaev, N.A. (1972) Dynamics of Nonholonomic Systems. American Mathematics Society.
[5]	Rao, C.R. and Mitra, S.K. (1971) Generalized Inverse of Matrices and Its Applications. Wiley, New York, London.
[6]	Mitra, S.K. and Rao, C.R. (1974) Projections under Seminorms and Generalized Moore Penrose Inverses. Linear Algebra and Its Applications, 9, 155-167. https://doi.org/10.1016/0024-3795(74)90034-2
[7]	Gauss, C.F. (1829) über ein neues allgemeines Grundgesetz der Mechanik. Journal für die reine und angewandte Mathematik Journal für die reine und angewandte Mathematik, 1829, 232-235. https://doi.org/10.1515/crll.1829.4.232
[8]	Lanczos, C. (1949) The Variational Principles of Mehanics. University of Toronto Press, Toronto. https://doi.org/10.3138/9781487583057
[9]	Kalaba, R., Natsuyama, H.H. and Udwadia, F.E. (2004) An Extension of Gauss’ Principle of Least Constraint. International Journal of General Systems, 33, 63-69. https://doi.org/10.1080/0308107031000139996
[10]	Popp, K. and Schiehlen, W. (2010) Ground Vehicle Dynamics. Springer, Berlin Heidelberg. https://doi.org/10.1007/978-3-540-68553-1
[11]	Cardona, A. and Géradin, M. (2001) Flexible Multibody Dynamics. Wiley, New York.
[12]	Oden, J. T. and Martins, J. (1985) Models and Computational Methods for Dynamic Friction Phenomena. Computer Methods in Applied Mechanics and Engineering, 52, 527-634. https://doi.org/10.1016/0045-7825(85)90009-X