MGWrapParticleAroundHorizontalCylinder2D.html  (MotionGenesis input/output).
```   (1) % MotionGenesis file:  MGWrapParticleAroundHorizontalCylinder2D.txt
(2) % Copyright (c) 2016 Motion Genesis LLC.  Only for use with MotionGenesis.
(3) %--------------------------------------------------------------------
(4) NewtonianFrame  N                % Cylinder with Nz> along symmetry axis.
(5) RigidFrame      B                % Taut string with By> directed towards Q.
(6) Particle        Q                % Particle at end of string.
(7) Point           BN( B )          % Point of string B in contact with cylinder N.
(8) %--------------------------------------------------------------------
(9) Variable    theta''              % Angle from Ny> to By> with sense +Nz>.
(10) Variable    Tension              % Tension in string.
(11) Variable    L''                  % Length of deployed string (between BN and Q).
(12) Constant    R = 3 cm             % Radius of cylinder N.
(13) Constant    g = 9.8 m/s^2        % Earth's gravitational acceleration.
(14) Q.SetMass( m = 100 g )
(15) %--------------------------------------------------------------------
(16) %   Rotational and translational kinematics.
(17) B.RotateZ( N, theta )
-> (18) B_N = [cos(theta), sin(theta), 0;  -sin(theta), cos(theta), 0;  0, 0, 1]
-> (19) w_B_N> = theta'*Bz>
-> (20) alf_B_N> = theta''*Bz>

(21) Q.Translate( No, R*Bx> + L*By> )
-> (22) p_No_Q> = R*Bx> + L*By>
-> (23) v_Q_N> = -L*theta'*Bx> + (L'+R*theta')*By>
-> (24) a_Q_N> = (-2*L'*theta'-R*theta'^2-L*theta'')*Bx> + (L''+R*theta''-L*theta'^2)*By>

(25) BN.SetPositionVelocity( Q,  -L*By>, B )
-> (26) p_Q_BN> = -L*By>
-> (27) v_BN_N> = (L'+R*theta')*By>

(28) %--------------------------------------------------------------------
(29) %   Constraint: Total length of string does not change.
(30) %   Alternate:  Dot( v_BN_N>, By> ) = 0.
(31) %   Let  s = Additional length of string wrapped around cylinder after t = 0.
(32) %        s = R * (theta - theta0) =  L0 - L
(33) Constant    L0 = 70 cm           % Initial value of L.
(34) Input       theta = 0 deg        % Initial value of theta.
(35) s = R * (theta - Input(theta))
-> (36) s = R*theta

(37) SolveDt( s = L0 - L,   L )       % Find L', L'' in terms of theta', theta''.
-> (38) L = L0 - s
-> (39) L' = -R*theta'
-> (40) L'' = -R*theta''

(41) %--------------------------------------------------------------------
(42) %   Q's velocity/acceleration simplifies by substituting L', L''.
(43) Explicit( v_Q_N>,  theta' )
-> (44) v_Q_N> = -L*theta'*Bx>

(45) Explicit( a_Q_N>,  theta', theta'' )
-> (46) a_Q_N> = (R*theta'^2-L*theta'')*Bx> - L*theta'^2*By>

(47) %--------------------------------------------------------------------
(48) %   Add relevant contact/distance forces.
(49) Q.AddForce(  -Tension*By> - m*g*Ny> )
-> (50) Force_Q> = -Tension*By> - m*g*Ny>

(51) %--------------------------------------------------------------------
(52) %   Form equations of motion via F = m*a and solve for theta'', Tension.
(53) Dynamics[1] = Dot(  Q.GetDynamics(),  Bx>  )
-> (54) Dynamics[1] = m*(g*sin(theta)+R*theta'^2-L*theta'')

(55) Dynamics[2] = Dot(  Q.GetDynamics(),  By>  )
-> (56) Dynamics[2] = Tension + m*g*cos(theta) - m*L*theta'^2

(57) Solve( Dynamics = 0,   theta'',  Tension )
-> (58) theta'' = (g*sin(theta)+R*theta'^2)/L
-> (59) Tension = -m*(g*cos(theta)-L*theta'^2)

(60) %--------------------------------------------------------------------
(61) %   Use initial velocity to solve for for initial value of theta'.
(62) Constant  v0 =  -4 m/s          % Initial Bx> measure of Q's velocity.
(63) vx = Dot( Q.GetVelocity(N), Bx> )
-> (64) vx = -L*theta'

(65) SolveSetInput( vx = v0,   theta' = 0 rad/sec )

->   %  INPUT has been assigned as follows:

(66) %--------------------------------------------------------------------
(67) %   Calculations for plotting Q's trajectory in N.
(68) x = Dot(  Q.GetPosition(No),  Nx>  )
-> (69) x = R*cos(theta) - L*sin(theta)

(70) y = Dot(  Q.GetPosition(No),  Ny>  )
-> (71) y = R*sin(theta) + L*cos(theta)

(72) %--------------------------------------------------------------------
(73) %   Calculations for plotting angular momentum and kinetic + potential energy.
(74) KineticEnergy = Q.GetKineticEnergy()
-> (75) KineticEnergy = 0.5*m*L^2*theta'^2

(76) GravityPotentialEnergy = Q.GetForceGravityPotentialEnergy( -g*Ny>, No )
-> (77) GravityPotentialEnergy = m*g*(R*sin(theta)+L*cos(theta))

(78) MechanicalEnergy = KineticEnergy + GravityPotentialEnergy
-> (79) MechanicalEnergy = GravityPotentialEnergy + KineticEnergy

(80) HNo = Dot(  Q.GetAngularMomentum( No ),  Nz>  )
-> (81) HNo = m*L^2*theta'

(82) HBN = Dot(  Q.GetAngularMomentum( BN ),  Nz>  )
-> (83) HBN = m*L^2*theta'

(84) vSpeed = Q.GetSpeed( N )
-> (85) vSpeed = abs(L*theta')

(86) %--------------------------------------------------------------------
(87) %   List output quantities, provide integration parameters, and solve ODEs.
(88) 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
(89) Input  tFinal = 1.536 sec,  tStep = 0.002 sec,  absError = 1.0E-07
(90) ODE()  MGWrapParticleAroundHorizontalCylinder2D

(91) %--------------------------------------------------------------------
```