% MotionGenesis file: MGWrapParticleAroundHorizontalCylinder2D.txt % Copyright (c) 2016 Motion Genesis LLC. Only for use with MotionGenesis. %-------------------------------------------------------------------- NewtonianFrame N % Cylinder with Nz> along symmetry axis. RigidFrame B % Taut string with By> directed towards Q. Particle Q % Particle at end of string. Point BN( B ) % Point of string B in contact with cylinder N. %-------------------------------------------------------------------- Variable theta'' % Angle from Ny> to By> with sense +Nz>. Variable Tension % Tension in string. Variable L'' % Length of deployed string (between BN and Q). Constant R = 3 cm % Radius of cylinder N. Constant g = 9.8 m/s^2 % Earth's gravitational acceleration. Q.SetMass( m = 100 g ) %-------------------------------------------------------------------- % Rotational and translational kinematics. B.RotateZ( N, theta ) Q.Translate( No, R*Bx> + L*By> ) BN.SetPositionVelocity( Q, -L*By>, B ) %-------------------------------------------------------------------- % Constraint: Total length of string does not change. % Alternate: Dot( v_BN_N>, By> ) = 0. % Let s = Additional length of string wrapped around cylinder after t = 0. % s = R * (theta - theta0) = L0 - L Constant L0 = 70 cm % Initial value of L. Input theta = 0 deg % Initial value of theta. s = R * (theta - Input(theta)) SolveDt( s = L0 - L, L ) % Find L', L'' in terms of theta', theta''. %-------------------------------------------------------------------- % Q's velocity/acceleration simplifies by substituting L', L''. Explicit( v_Q_N>, theta' ) Explicit( a_Q_N>, theta', theta'' ) %-------------------------------------------------------------------- % Add relevant contact/distance forces. Q.AddForce( -Tension*By> - m*g*Ny> ) %-------------------------------------------------------------------- % Form equations of motion via F = m*a and solve for theta'', Tension. Dynamics[1] = Dot( Q.GetDynamics(), Bx> ) Dynamics[2] = Dot( Q.GetDynamics(), By> ) Solve( Dynamics = 0, theta'', Tension ) %-------------------------------------------------------------------- % Use initial velocity to solve for for initial value of theta'. Constant v0 = -4 m/s % Initial Bx> measure of Q's velocity. vx = Dot( Q.GetVelocity(N), Bx> ) SolveSetInput( vx = v0, theta' = 0 rad/sec ) %-------------------------------------------------------------------- % Calculations for plotting Q's trajectory in N. x = Dot( Q.GetPosition(No), Nx> ) y = Dot( Q.GetPosition(No), Ny> ) %-------------------------------------------------------------------- % Calculations for plotting angular momentum and kinetic + potential energy. KineticEnergy = Q.GetKineticEnergy() GravityPotentialEnergy = Q.GetForceGravityPotentialEnergy( -g*Ny>, No ) MechanicalEnergy = KineticEnergy + GravityPotentialEnergy HNo = Dot( Q.GetAngularMomentum( No ), Nz> ) HBN = Dot( Q.GetAngularMomentum( BN ), Nz> ) vSpeed = Q.GetSpeed( N ) %-------------------------------------------------------------------- % List output quantities, provide integration parameters, and solve ODEs. Output t sec, x cm, y cm, L cm, theta deg, Tension N, MechanicalEnergy Joules, vSpeed cm/s, HNo kg*m^2/s^2, HBN kg*m^2/s^2 Input tFinal = 1.536 sec, tStep = 0.002 sec, absError = 1.0E-07 ODE() MGWrapParticleAroundHorizontalCylinder2D %-------------------------------------------------------------------- Save MGWrapParticleAroundHorizontalCylinder2D.html Quit