MGBeamOnTwoCablesDynamicsViaKaneLagrangeAugmentEmbedConstraints.html  (MotionGenesis input/output).
   (1) %    File: MGBeamOnTwoCablesDynamicsViaKaneLagrangeAugmentEmbedConstraints.txt
   (2) % Problem: Dynamics of a beam supported by two cables.
   (3) %--------------------------------------------------------------------
   (4) NewtonianFrame N         % Ceiling with Nx> horizontally right, Ny> vertically down
   (5) RigidBody      B         % Beam with Bx> pointing from Bo to Bc and Bz> = Nz>.
   (6) Point          Nc( N )   % Point of N attached to cable C.
   (7) Point          Bc( B )   % Point of B attached to cable C.
   (8) %--------------------------------------------------------------------
   (9) Constant  LN = 6 m       % Distance between No and NC.
   (10) Constant  LB = 4 m       % Distance between Bo and BC.
   (11) Constant  LA = 2.7 m     % Length of cable A (connects Bo to No).
   (12) Constant  LC = 3.7 m     % Length of cable C (connects BC to NC).
   (13) Constant  g = 9.8 m/s^2  % Earth's gravitational acceleration.
   (14) Constant  b = 230 N*s/m  % Aerodynamic damping constant.
   (15) B.SetMass( m = 120 kg )
   (16) B.SetInertia( Bcm, Izz = m*LB^2/12, 0, Izz )
-> (17) Izz = 0.08333333*m*LB^2

   (18) %--------------------------------------------------------------------
   (19) Variable  x'', y''       % Nx> and Ny> measures of Bo's position vector from No.
   (20) Variable  q''            % Angle from Nx> to Bx> with positive sense +Nz>.
   (21) Variable  TC             % Tension in cable C (tension in cable A is non-contributing).
   (22) %--------------------------------------------------------------------
   (23) %   Rotational kinematics relating  Bx>, By>, Bz>  to  Nx>, Ny>, Nz>.
   (24) B.RotateZ( N,  q )
-> (25) B_N = [cos(q), sin(q), 0;  -sin(q), cos(q), 0;  0, 0, 1]
-> (26) w_B_N> = q'*Bz>
-> (27) alf_B_N> = q''*Bz>

   (28) %--------------------------------------------------------------------
   (29) %   Translational kinematics.
   (30) Bo.Translate( No,  x*Nx> + y*Ny> )
-> (31) p_No_Bo> = x*Nx> + y*Ny>
-> (32) v_Bo_N> = x'*Nx> + y'*Ny>
-> (33) a_Bo_N> = x''*Nx> + y''*Ny>

   (34) Bcm.Translate( Bo,  0.5*LB*Bx> )
-> (35) p_Bo_Bcm> = 0.5*LB*Bx>
-> (36) v_Bcm_N> = 0.5*LB*q'*By> + x'*Nx> + y'*Ny>
-> (37) a_Bcm_N> = -0.5*LB*q'^2*Bx> + 0.5*LB*q''*By> + x''*Nx> + y''*Ny>

   (38) BC.SetPosition( Bo,  LB*Bx> )
-> (39) p_Bo_Bc> = LB*Bx>

   (40) NC.SetPosition( No,  LN*Nx> )
-> (41) p_No_Nc> = LN*Nx>

   (42) %--------------------------------------------------------------------
   (43) %   Constraint relating cable lengths to magnitudes of position vectors.
   (44) LengthConstraint[1] = LA^2 - Bo.GetDistanceSquared( No )
-> (45) LengthConstraint[1] = LA^2 - x^2 - y^2

   (46) LengthConstraint[2] = LC^2 - Bc.GetDistanceSquared( Nc )
-> (47) LengthConstraint[2] = LC^2 + 2*LB*cos(q)*(LN-x) - LB^2 - y^2 - 2*LB*y*sin(q)
        - (LN-x)^2

   (48) %--------------------------------------------------------------------
   (49) %   Relevant forces on beam (to calculate generalized forces).
   (50) %   Note: Cable A's tension TA does not contribute to generalized forces
   (51) %   because LengthConstraint[1] is embedded, whereas cable C's tension TC does
   (52) %   contribute to generalized forces because LengthConstraint[2] is augmented.
   (53) unitVectorBcToNc> = NC.GetPositionVector(Bc) / LC
-> (54) unitVectorBcToNc> = -LB/LC*Bx> + (LN-x)/LC*Nx> - y/LC*Ny>

   (55) BC.AddForce( TC*unitVectorBcToNc> )
-> (56) Force_Bc> = -LB*TC/LC*Bx> + TC*(LN-x)/LC*Nx> - TC*y/LC*Ny>

   (57) Bcm.AddForce( m*g*Ny> )
-> (58) Force_Bcm> = m*g*Ny>

   (59) %--------------------------------------------------------------------
   (60) %   Use LengthConstraint[1] to solve for y in terms of x, q.
   (61) %   Differentiate LengthConstraint[1] to solve for y' in terms of x', q'.
   (62) SolveQuadraticPositiveRootDt( LengthConstraint[1],  y )
-> (63) y = sqrt(LA^2-x^2)
-> (64) y' = -x*x'/y
-> (65) y'' = -(x'^2+y'^2+x*x'')/y

   (66) %--------------------------------------------------------------------
   (67) %   Dynamics via Kane's augmented/embedded constraint method.
   (68) %   Embed LengthConstraint[1] and augment LengthConstraint[2].
   (69) SetGeneralizedSpeed( x', q' )
   (70) Dynamics = System.GetDynamicsKane()
-> (71) Dynamics[1] = m*g*x/y + TC*(LB*cos(q)-LN-LB*x*sin(q)/y)/LC + m*(1+x^2/y^2)*x''
        - 0.5*m*(LB*cos(q)*q'^2-x*(LB*sin(q)*q'^2+2*(1+x^2/y^2)*x'^2/y)/y)
        - 0.5*m*LB*(sin(q)+x*cos(q)/y)*q''

-> (72) Dynamics[2] = 0.5*LB*(2*TC*(y*cos(q)+sin(q)*(LN-x))/LC-m*g*cos(q)-m*cos
        (q)*(1+x^2/y^2)*x'^2/y) + 0.25*(m*LB^2+4*Izz)*q'' - 0.5*m*LB*(sin(q)+x*
        cos(q)/y)*x''

   (73) %--------------------------------------------------------------------
   (74) %   Integration parameters and initial values.
   (75) Input  tFinal = 9 sec, tStep = 0.02 sec,  absError = 1.0E-07
   (76) Input  x = 2 m,  x' = 0 m/s
   (77) %--------------------------------------------------------------------
   (78) %   Use LengthConstraint[2] to set initial value for q (with guess).
   (79) SolveSetInput( LengthConstraint[2] = 0,   q = 15 deg )

->   %  INPUT has been assigned as follows:
->   %   q                         27.67399398437716       deg

   (80) %--------------------------------------------------------------------
   (81) %   Differentiate length constraints to set initial value for q' (with guess).
   (82) SolveSetInput( Dt(LengthConstraint[2]) = 0,  q' = 0 m/s )

->   %  INPUT has been assigned as follows:
->   %   q'                        0                       m/s

   (83) %--------------------------------------------------------------------
   (84) %   Quantities to be output by ODE command.
   (85) Output  t sec,  x m,  y m,  q degrees,  TC Newtons
   (86) %--------------------------------------------------------------------
   (87) %   Numerically integrate to solve for motion.
   (88) ODE( [Dynamics; DtDt(LengthConstraint[2])] = 0,  x'', q'', TC )  MGBeamOnTwoCablesDynamics

   (89) %--------------------------------------------------------------------
   (90) %   Save input and output responses.
Saved by Motion Genesis LLC.   Copyright (c) 2009-2018 Motion Genesis LLC on command names and syntax. All rights reserved.