Kinematics of a Holonomic Constraint Rod & Cube System

Haiduke Sarafian

The Pennsylvania State University, University College, York, PA, USA.

DOI: 10.4236/wjm.2018.86018                

Abstract

A loose rod of mass m1 and length l leans against one of the faces of a cube of mass m₂ and side length a. The assembly is placed on a horizontal table with one end of the rod touching the table and its other end leaning against the edge of the cube. We set the rod and the center of mass of the cube on the same vertical plane, and then we release the assembly from the rest. For frictionless contacts, we calculate the separation runtime of the rod from the cube as a function of m₂/m₁ and a/l. This entails forming the equation describing the motion of the system. The equation of motion is analytically unsolvable nonlinear differential equation. Applying a Computer Algebra System, specifically Mathematica [1] [2], we solve the equation numerically. Utilizing the solution, in addition to evaluating the separation runtime, we quantify a list of dynamic quantities, such as the time-dependent interface forces, and, geometric quantities, such as the trajectory of the loose end of the rod. A robust Mathematica code addresses the “what if” scenarios.

Keywords

Lagrangian and Hamiltonian Methods, Computer Algebra System, Mathematica

1. Introduction

Motivation of suggested investigation stems from the fact that the proposed assembly is composed of a sliding cube that acts as a point-like object with only linear kinematics in contrast to the rod that in addition to the former possesses rotational degrees of freedom. The length of the rod adds additional features to the physics of the problem in contrast to the movement of the pair point-like objects in the two-body problem reported in [3] . The leaning rod against the cube constitutes a holonomic constraint reducing the number of degrees of freedom. These two features, namely, 1) rotation of the rod and 2) its constraint movement sparked the investigation. The problem appears to be a 2D two-body problem with numerous degrees of freedom however; the restricted motion of the rod practically reduces to a problem with only one degree of freedom. In so many words, the kinematic information about the entire system comes about from the solution of a single nonlinear differential equation with its dependent variable being the inclination angle of the rod w/floor. The equation of motion is analytically unsolvable, proving the effective usefulness of a Computer Algebra System, specifically Mathematica. Utilizing the numeric solution in hand, we quantify a list of kinematic information primary as noted the separation runtime of the rod-cube system. Additional quantities such as speed, acceleration and contact forces of the bodies are also quantified as well. This article is composed of three sections. In addition to Introduction in Section 2 applying Lagrangian and Hamiltonian methods, we develop the needed formulation conducive to equation of motion. This section includes sub-sections: Rod-Cube contact force and Trajectory of the tip of the rod. Section 3 is the Results. In this section, we utilize Mathematica numeric utilities solving the needed equations conducive to quantifying kinematic and dynamic quantities of interest. We close the report with concluding remarks suggesting ideas furthering the investigation.

2. Physics of the Problem and Its Solution

Figure 1 depicts the schematic of the problem on hand.

As shown in Figure 1, one end of a uniform rod of length ℓ and mass m₁ touches the horizontal surface at x₁, its other end leans against the left vertical face of the cube at C(x₂, y₂). The cube has a mass of m₂ with side length a. The rod and the center of mass (cm) of the cube are on a vertical plane. Relevant forces acting on the rod namely the contact with the floor, N, the weight, W, and the one with the cube, F₂₁, are also shown.

Figure 1. A rod of mass m₁ and length ℓ leans against a cube of mass m₂ and side length a. The assembly is on a frictionless horizontal surface.

Releasing the assembly from rest sets the system in motion. This is because for the frictionless surfaces there is an unbalanced horizontal contact force acting on the cube, not shown in Figure 1; it acts opposite to shown F₂₁. As depicted, the same value of the force F₂₁ acts on the rod pushing it to the left setting it in motion. However, as the rod slides its center of mass falls and the rod rotates about the sliding contact point x₁. While the rod falls its “circulation radius” shortens hence separation is warranted. In a different scenario if the rod was hinged at x₁ so that circulation radius would have stayed constant separation would have occurred. For either scenario one speculates the separation runtime would depend to the masses in motion and respective dimensions of the bodies. It is one of the objectives of this investigation to quantify separation runtime. As we were developing the formulation of the problem also we quantify a list of interesting quantities. The proposed problem is different from the “classic” version [4] [5] where the rod leans against a vertical wall.

To bypass the dynamic quantities, i.e. forces, for time being we formulate the problem following Lagrangian method. Lagrangian of the system is composed of two pieces, one for the cube and the other one for the rod. Kinetic energy of the

cube is $T_{cube}=\frac{1}{2}{m}_2{\stackrel{\dot{}}{x}}_2^2$ and its potential energy with respect to the table top is

$V_{cube}=m_{2}g\frac{a}{2}$ . Kinetic energy of the rod in addition to a 2D motion of the cm includes its rotational energy about the cm. These two terms are,

$T_{rod}=\frac{1}{2}{m}_1{v}_{cm}^2+\frac{1}{2}{I}_{cm}{\omega }^2$

where ω equals $\stackrel{\dot{}}{\theta }$ is the angular velocity of the rod about the cm; $I_{cm}=\frac{1}{12}{m}_2{\mathcal{l}}^2$ . Its potential energy is $V_{rod}=m_{1}g{y}_{1cm}$ ; g is the gravity constant; according to tradition dotted quantities are the derivatives w/time. As shown in Figure 1, it is clear that the angular velocity of the rod about its cm is the same as its angular velocity w/contact point; θ is the inclination angle of the rod w/horizontal.

Conveniently the origin is set as shown. Accordingly,

$\{\begin{array}{l}{x}_2=x_1+\mathcal{l}\cos \left(\theta \right)\\ x_{1cm}=x_1+\frac{1}{2}\mathcal{l}\cos \left(\theta \right)\\ y_{1cm}=\frac{1}{2}\mathcal{l}\sin \left(\theta \right)\end{array}$ (1)

their respective velocities are,

$\{\begin{array}{l}{\stackrel{\dot{}}{x}}_2={\stackrel{\dot{}}{x}}_1-\mathcal{l}\stackrel{\dot{}}{\theta }\sin \left(\theta \right)\\ {\stackrel{\dot{}}{x}}_{1cm}={\stackrel{\dot{}}{x}}_1-\frac{1}{2}\mathcal{l}\stackrel{\dot{}}{\theta }\sin \left(\theta \right)\\ {\stackrel{\dot{}}{y}}_{1cm}=\frac{1}{2}\mathcal{l}\stackrel{\dot{}}{\theta }\cos \left(\theta \right)\end{array}$ (2)

conducive to, $v_{cm}^2={\stackrel{\dot{}}{x}}^2{}_{1cm}+{\stackrel{\dot{}}{y}}^2{}_{1cm}$ .

Utilizing the last equation and Equation (2), Lagrangian of the system, e.g. $L=T_{rod}+T_{cube}-\left(V_{rod}+V_{cube}\right)$ gives,

$\begin{array}{c}L=\frac{1}{2}{m}_1\left[{\stackrel{\dot{}}{x}}_1^2+{\left(\frac{\mathcal{l}}{2}\right)}^2{\stackrel{\dot{}}{\theta }}^2-\mathcal{l}{\stackrel{\dot{}}{x}}_1\stackrel{\dot{}}{\theta }\sin \left(\theta \right)\right]+\frac{1}{2}{I}_{cm}{\stackrel{\dot{}}{\theta }}^2\\ +\frac{1}{2}{m}_2{\left[{\stackrel{\dot{}}{x}}_1-\mathcal{l}\stackrel{\dot{}}{\theta }\sin \left(\theta \right)\right]}^2-m_{1}g\frac{\mathcal{l}}{2}\sin \left(\theta \right)-m_{2}g\frac{a}{2}\end{array}$ (3)

Equation (3) reveals the Lagrangian is independent of x₁, meaning x₁ is a cyclic variable; $\frac{\partial L}{\partial x_1}=0$ . Because $\frac{\partial L}{\partial {\stackrel{\dot{}}{x}}_1}\equiv P_{x_1}$ according to Euler-Lagrange equation $\frac{\text{d}}{\text{d}t}{P}_{x_1}=0$ warrants the conservation of the linear momentum of the system

along the horizontal direction. Explicitly, momentum of the system $P_{x_1}$ is,

$P_{x_1}\equiv \frac{\partial L}{\partial {\stackrel{\dot{}}{x}}_1}=m_1\left[{\stackrel{\dot{}}{x}}_1-\frac{1}{2}\mathcal{l}\stackrel{\dot{}}{\theta }\sin \left(\theta \right)\right]+m_2\left[{\stackrel{\dot{}}{x}}_1-\mathcal{l}\stackrel{\dot{}}{\theta }\sin \left(\theta \right)\right]=0$ (4)

The RHS of Equation (4) is set zero in accordance the initial rest state of the system. Accordingly, Equation (4) sets a relationship between the two variables, namely,

${\stackrel{\dot{}}{x}}_1=\left(\frac{\frac{1}{2}{m}_1+m_2}{m_1+m_2}\right)\mathcal{l}\stackrel{\dot{}}{\theta }\sin \left(\theta \right)$ (5)

Substituting Equation (5) in (3), yields,

$L\left(\theta ,\stackrel{\dot{}}{\theta }\right)=\frac{1}{2}{m}_1{\mathcal{l}}^2\left[-\frac{1}{4}\frac{m_1}{m_1+m_2}{\sin }^2\left(\theta \right)+\frac{1}{3}\right]{\stackrel{\dot{}}{\theta }}^2-m_{1}g\frac{\mathcal{l}}{2}\sin \left(\theta \right)-m_{2}g\frac{a}{2}$ (6)

Therefore substituting Equation (6) to Euler-Lagrange equation, $\frac{\text{d}}{\text{d}t}\frac{\partial L}{\partial \stackrel{\dot{}}{\theta }}-\frac{\partial L}{\partial \theta }=0$

yields the equation of motion. Equation of motion is a challenging second order ODE. Although it can be solved numerically, however, the primary objective is not to identify θ = θ(t), but rather t = t(θ). Therefore, we make a strategy change, namely, because there are no dissipative forces, conservation of energy is warrant we consider the Hamiltonian of the system, namely, H = T + V $\equiv$ Energy, this yields,

$H\left(\theta ,\stackrel{\dot{}}{\theta }\right)\equiv E=\frac{1}{2}{m}_1{\mathcal{l}}^2\left[-\frac{1}{4}\frac{m_1}{m_1+m_2}{\sin }^2\left(\theta \right)+\frac{1}{3}\right]{\stackrel{\dot{}}{\theta }}^2+m_{1}g\frac{\mathcal{l}}{2}\sin \left(\theta \right)+m_{2}g\frac{a}{2}$ (7)

Equating this to the value of the energy of the system at rest, E₀, gives,

$\begin{array}{l}\frac{1}{2}{m}_1{\mathcal{l}}^2\left[-\frac{1}{4}\frac{m_1}{m_1+m_2}{\sin }^2\left(\theta \right)+\frac{1}{3}\right]{\stackrel{\dot{}}{\theta }}^2+m_{1}g\frac{\mathcal{l}}{2}\sin \left(\theta \right)+m_{2}g\frac{a}{2}\\ =E_0\equiv m_{1}g\frac{\mathcal{l}}{2}\left(\frac{a}{\mathcal{l}}\right)+m_{2}g\frac{a}{2}\end{array}$ (8)

where, (a/ℓ), the first term of the RHS equals sin(θ₀) = a/ℓ, is the initial inclination angle of the rod w/horizontal.

Equation (8) simplifies,

${\stackrel{\dot{}}{\theta }}^2=4\frac{g}{\mathcal{l}}\frac{\frac{a}{\mathcal{l}}-\sin \left(\theta \right)}{\frac{4}{3}-\left(\frac{1}{1+\frac{m_2}{m_1}}\right){\sin }^2\left(\theta \right)}$ (9)

Equation (9) explicitly shows the dependence of the angular velocity, $\stackrel{\dot{}}{\theta }$ to the masses and the relevant dimensions. Rearranging Equation (9) its integration gives,

$t=\frac{1}{2}\sqrt{\frac{\mathcal{l}}{g}}\underset{0}{\overset{\arcsin \left(\frac{a}{\mathcal{l}}\right)}{{\displaystyle \int }}}\sqrt{\frac{\frac{4}{3}-\frac{1}{1+n}{\sin }^2\left(\theta \right)}{\frac{a}{\mathcal{l}}-\sin \left(\theta \right)}}\text{d}\theta$ (10)

where we set n = (m₂/m₁). The lower limit of integration indicates the runtime of the rod to the table top. Our primary objective is to evaluate the separation runtime; hence, the lower limit of the integration should be the separation angle. This angle comes about from dynamics aspect of the problem.

The driving force acting on the cube, F₁₂, accelerates the cube according to, $F_{12}=m_2{\ddot{x}}_2$ . At separation instant contact force, F₁₂, vanishes result an equation conducive to the separation angle. First we calculate the needed acceleration. Applying the first equation of set Equation (2) and equation (5) gives,

${\ddot{x}}_2=-\frac{\mathcal{l}}{2\left(1+n\right)}\left[\ddot{\theta }\sin \left(\theta \right)+{\stackrel{\dot{}}{\theta }}^2\cos \left(\theta \right)\right]$ (11)

On the other hand differentiating Equation (9) w/time gives,

$\ddot{\theta }=\frac{g}{2\mathcal{l}}\left\{-\frac{4\cos \left(\theta \right)}{\frac{4}{3}-\frac{1}{1+n}{\sin }^2\left(\theta \right)}+\frac{1}{1+n}\frac{\left(1-4\sin \left(\theta \right)\right)\sin \left(2\theta \right)}{{\left[\frac{4}{3}-\frac{1}{1+n}{\sin }^2\left(\theta \right)\right]}^2}\right\}$ (12)

And finally substituting Equations (9) and (12) in Equation (11) and simplification yields the needed equation. The output is a lengthy fraction, its numerator is,

$-6g\mathcal{l}\left(1+n\right)\cos \left(\theta \right)\left[-8\frac{a}{\mathcal{l}}\left(1+n\right)+12\left(1+n\right)\sin \left(\theta \right)-3{\sin }^3\left(\theta \right)\right]=\text{numerator}$ (13)

Setting Equation (13) zero gives the needed equation whose root is the separation angle, namely,

$3{\sin }^3\left(\theta \right)-12\left(1+n\right)\sin \left(\theta \right)+8\frac{a}{\mathcal{l}}\left(1+n\right)=0$ (14)

As expected terms such as, n $\equiv$ m₂/m₁ and a/ℓ control the value of the separation angle, θ_sep, yielding the needed runtime,

$t=\frac{1}{2}\sqrt{\frac{\mathcal{l}}{g}}\underset{{\theta }_{sep}}{\overset{\arcsin \left(\frac{a}{\mathcal{l}}\right)}{{\displaystyle \int }}}\sqrt{\frac{\frac{4}{3}-\frac{1}{1+n}{\sin }^2\left(\theta \right)}{\frac{a}{\mathcal{l}}-\sin \left(\theta \right)}}\text{d}\theta$ (15)

2.1. Rod-Cube Contact Force

As the rod slides its inclination angle, θ, with respect to the table top decreases, reduces the contact force between the rod and the cube. To quantify the contact force at any given time, we write, $F_{12}=m_2{\ddot{x}}_2$ . Acceleration of the cube, ${\ddot{x}}_2$ , utilizing Equation (11), (12) and (9), yields,

$F_{12}=m_2\left\{-3g\frac{\cos \left(\theta \right)\left[8\frac{a}{\mathcal{l}}\left(1+n\right)-12\left(1+n\right)\sin \left(\theta \right)+3{\sin }^3\left(\theta \right)\right]}{{\left[4\left(1+n\right)-3{\sin }^2\left(\theta \right)\right]}^2}\right\}$ (16)

This is the formal solution, however, because θ = θ(t) there is a need to solve the equation of motion for θ(t). The equation of motion is actually Equation (9). Its differential form is,

$\stackrel{\dot{}}{\theta }+\sqrt{4\left(\frac{g}{\mathcal{l}}\right)\frac{\frac{a}{\mathcal{l}}-\sin \left(\theta \right)}{\frac{4}{3}-\left(\frac{1}{1+n}\right){\sin }^2\left(\theta \right)}}=0$ (17)

As demonstrated, by arranging Equation (17) we calculate the runtime; see Equation (10) and (15). Reversing these equations to obtain θ(t) is not possible; one needs to solve Equation (17) directly, this can’t be done analytically. In the next section for a chosen set of parameters we present its numeric solution.

2.2. Trajectory of the Tip of the Rod

As the contact point of the rod with the table top slides to the left its other end that leans against the mobile cube before breaking away traces a certain curve. It is our interest to plot its trajectory. We offer two solutions. First, a formal analytic solution, and second a straight forward plotting routine. Both solutions hinge upon the know solution of Equation (17).

Analytic solution:

Integrating both sides of Equation (5) w/time yields,

$x_1=-\frac{1+2n}{2\left(1+n\right)}\mathcal{l}\left[\cos \left(\theta \right)-\beta \right]$ (18)

here we set the initial contact point of the rod with the table at the origin yielding the value of β = cos[arcsin(a/ℓ)]. Substituting Equation (18) in the first equation of set of Equation (1) after simplification, gives,

$\cos \left(\theta \right)=\frac{2}{\mathcal{l}}\left(1+n\right)\left[x_2-\frac{1+2n}{2\left(1+n\right)}\mathcal{l}\beta \right]$ (19)

on the other hand Figure 1 shows,

$\sin \left(\theta \right)=\frac{y_2}{\mathcal{l}}$ (20)

Squaring both sides of Equation (19) and (20) and adding them results,

${\left[2\left(1+n\right)\frac{x_2}{\mathcal{l}}-\left(1+2n\right)\beta \right]}^2+{\left(\frac{y_2}{\mathcal{l}}\right)}^2=1$ (21)

Equation (21) deceptively appears as an ellipse; it is not. Coordinates of {x₂, y₂} are time-dependent, they are implicit functions of θ(t). Equation (21) can be utilized tracing the tip of the rod; this is left to the interest of the reader. Alternatively, utilizing Equations (19) and (20) are arranged as,

$\left\{x_2\left(t\right),y_2\left(t\right)\right\}=\left\{\frac{1}{2\left(1+n\right)}\mathcal{l}\cos \left[\theta \left(t\right)\right]+\frac{1+2n}{2\left(1+n\right)}\mathcal{l}\beta ,\mathcal{l}\sin \left[\theta \left(t\right)\right]\right\}$ (22)

where θ(t) is the solution of the equation of motion; Equation (17). For a set of chosen parameters n and ℓ and utilizing the solution of Equation (17) for the same set of parameters along with a parametric plot of Equation (22) yields the seek for trajectory. This is shown in the next section.

3. Results

It is important emphasizing that entire calculation presented in this section is carried out applying one of the most powerful Computer Algebra System, Mathematica [1] . Equation (14) is generic; its root(s) is the separation angle. The input parameters of this equation may be any desired values, for instance for a practical case we set n $\equiv$ m₂/m₁ = 3, and a/ℓ = 1/4, making the combination a/ℓ (1 + n) = 1, simplifies the equation. With this choice of parameters it yields a unique solution, namely, θ = 9.61˚. In other words the rod from its initial position θ(t = 0) = 14.47˚ falls only 4.86˚ and separates. Equation (14) is a cubic trigonometric equation, its cubic term, 3sin³(θ), can be dropped reducing it to a linear trig equation, 6sin(θ) − 1 = 0, with a root of, θ = 9.59˚. This neat approach underestimating the exact solution only by 0.2%!

Utilizing the separation angle and applying Equation (15) we evaluate the associated runtime by numeric integration. For a 1.0 m rod this gives, t_runtime = 0.108 s. This is compared to the runtime when we assume the rod falls all the way to the table, t = 0.187 s; Equation (10). The latter is compared to the runtime

of a point-like mass at the cm of the rod, $t=\sqrt{\frac{a}{g}}=0.159\text{s}$ . As expected the

aforementioned elapsed time overvalues the latter, this is because the rotating rod requires a longer time.

3.1. Solution of the Equation of Motion

Figure 2 depicts the time dependent variation of the inclination angle. The rod begins sliding with the initial inclination angle of arcsin(1/4) = 14.4˚, in 0.108 s at 9.61˚ it separates from the cube.

3.2. Contact Force between the Rod and the Cube

Figure 3 shows the maximum value of the contact force between the rod and the cube occurs at the beginning where the inclination angle is at maximum. As the rod falls and rotates the contact force weakens. By the time the rod separates from the cube its value plunges zero. Comparing Figure 2 and Figure 3 reveals

Figure 2. Time dependent inclination angle θ vs. t.

Figure 3. Display of the contact force per mass vs. t.

despite of the complicated nature of the force, Equation (16), contact force closely traces the time dependent features of the inclination angle vs. time!

3.3. Trajectory (Figure 4)

It appears the tip of the rod as it falls, slides and rotates traces a vertical “straight” line. We have utilized parametric equation of the trajectory; Equation (22).

Figure 4. The solid line is the trajectory of the tip of the rod as it falls to the table.

4. Conclusions and Comments

“Sliding Ladder” is a classic physics problem. Its solution is readily available. For instance, the author referred to reference [4] . This particular reference is cited because its free animated solution is also available as well [5] . Our investigation is a modified version of the “Sliding Ladder” problem. Here we replace the stationary wall with a massive mobile cube.

The author searched for the analysis of the proposed problem in vain ending up crafting its solution in this note. Applying Mathematica, the equation of motion, Equation (17), is solved numerically. Utilizing the solution, we determined time-dependent quantities such as runtime, contact force, and trajectory of the tip of the rod. We have also determined additional information such as separation angle of the rod and the cube. Our investigation embodies long-hand analytic formation conducive to parameter-dependent equations including: side length of the cube and its mass, length of the rod and its mass. Applying Mathematica numeric solver and its powerful graphic features as an example, we presented variety of output for a specific set of parameters. Presented information in this report challenges the interested reader to explore new features applying a different set of parameters. Reference [6] embodies Mathematica codes and information needed to produce the graphs embedded in this article.

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	(2015) Mathematica Is Symbolic Computation Software. V11.0, Wolfram Research Inc.
[2]	Wolfram, S. (1996) Mathematica Book. 3rd Edition, Cambridge University Press, Cambridge.
[3]	Sarafian, H. and Hickey, N. (2017) Characteristics of a Two-Body Holonomic Constraint Mechanical System. World Journal of Mechanics, 7, 161-166. http://www.scirp.org/journal/wjm
[4]	Variations of the Sliding Ladder Problem. https://www.maa.org/sites/default/files/Kapranidis.pdf
[5]	http://demonstrations.wolfram.com/PhysicsOfASlidingLadder/
[6]	Sarafian, H. (2015) Mathematica Graphics Example Book for Beginners. Scientific Research Publishing. http://www.scirp.org