Analysis of Rigid Body Motion in 2D using Matlab - Prof. Dan Marghitu, Study notes of Mechanical Engineering

Download Analysis of Rigid Body Motion in 2D using Matlab - Prof. Dan Marghitu and more Study notes Mechanical Engineering in PDF only on Docsity! A θ L O x y A θ x y C A θ x y C 0 1 OO Figure 1  ı α ω an ta vC C C F01x F01y G (a) (b) (c) aC n t n t Figure 2 G1 G2 F F F F F F (b) (a) moment equation: IO alpha = sum M wrt O = rC x G 1 ÅÅÅÅ 6 L m H-3 g Cos@theta@tDD + 2 L theta££@tDL=0 Solution: theta''@tD= 3 g Cos@theta@tDDÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅ 2 L Method II force equation: m aC = sum F = FO + G projection on x: 1 ÅÅÅÅ 2 H-2 FOx - L m Cos@theta@tDD theta£@tD2 - L m Sin@theta@tDD theta££@tDL=0 H1L projection on y: 1 ÅÅÅÅ 2 H-2 HFOy + g mL - L m Sin@theta@tDD theta£@tD2 + L m Cos@theta@tDD theta££@tDL=0 H2L moment eqation: IC alpha = sum M wrt C = -rC x FO projection on z: 1 ÅÅÅÅÅÅÅ 12 L H6 FOy Cos@theta@tDD - 6 FOx Sin@theta@tDD + L m theta££@tDL=0 H3L from Eq.H1L => FOx = 1ÅÅÅÅ 2 H-L m Cos@theta@tDD theta£@tD2 - L m Sin@theta@tDD theta££@tDL from Eq.H2L => FOy = 1ÅÅÅÅ 2 H-2 g m - L m Sin@theta@tDD theta£@tD2 + L m Cos@theta@tDD theta££@tDL from Eqs.H1LH2LH3L => theta''@tD = 3 g Cos@theta@tDDÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅ 2 L Initial Conditions at t=0: theta@0D=0, theta'@0D=omega@0D=0 numerical data: m=12ê32.2,L=3,g=32.2 theta''@0D=alpha@0D= 3 gÅÅÅÅÅÅÅÅ 2 L =16.1 radês^2 V=F0y@0D = -g m + 1ÅÅÅÅ 2 L m theta££@0D=- g mÅÅÅÅÅÅÅÅ 4 =-3. lb H=F0x@0D = 0 lb 0.5 1 1.5 2 t@sD 25 50 75 100 125 150 175 theta@degD Appendix1.nb 3 0.5 1 1.5 2 t@sD -4 -2 2 4 omega@radêsD 0.5 1 1.5 2 t@sD -10 10 20 alpha@radêsD Solve::ifun : Inverse functions are being used by Solve, so some solutions may not be found; use Reduce for complete solution information . More… 99theta@tD Ø 2 JacobiAmplitudeA- 1ÅÅÅÅ 2 "##############################################################################H-2 C - C@1DL H-t2 - 2 t C@2D - C@2D2L , 4 CÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅ 2 C + C@1D E=, 9theta@tD Ø 2 JacobiAmplitudeA 1ÅÅÅÅ 2 "##############################################################################H-2 C - C@1DL H-t2 - 2 t C@2D - C@2D2L , 4 CÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅÅ 2 C + C@1D E== Appendix1.nb 4 Apply@Clear, Names@"Global`*"DD; Off@General::spellD; Off@General::spell1D; H*Input data*L data = 8m1 Ø 1., m2 Ø 1., L1 Ø 1., L2 Ø 1., g Ø 10<; H*data=8m1Øm, m2Øm,L1ØL,L2ØL<;*L IC1 = m1ê 12* HL1^2L; IA = IC1 + m1 HL1ê 2L^2; IC2 = m2ê 12* HL2^2L; H*Position,velocity and acceleration vectors*L xB = L1* Cos@q1@tDD; yB = L1* Sin@q1@tDD; rB = 8xB, yB, 0<; rC1 = rB ê2.; vC1 = D@rC1, tD; aC1 = D@vC1, tD; xC = xB + L2* Cos@q2@tDD; yC = yB + L2* Sin@q2@tDD; rC = 8xC, yC, 0<; rC2 = HrB + rCLê 2.; vC2 = D@rC2, tD; aC2 = D@vC2, tD; Print@"rC1=", rC1D; Print@"rB=", rBD; Print@"rC2=", rC2D; Print@"rC=", rCD; Print@"aC1=", aC1D; Print@"aC2=", aC2D; H*Angular velocities and accelerations*L omega1 = 80, 0, q1'@tD<; omega2 = 80, 0, q2'@tD<; alpha1 = 80, 0, q1''@tD<; alpha2 = 80, 0, q2''@tD<; Print@"alpha1=", alpha1D; Print@"alpha2=", alpha2D; F01 = 8F01x, F01y, 0<; F21 = 8F21x, F21y, 0<; Print@"Joint reaction at B: F21=8F21x,F21y,0<"D; G1 = 80, -m1 g, 0<; G2 = 80, -m2 g, 0<; H*Newton equations*L Appendix2.nb 1 1 2 3 4 5 t@sD -10 -8 -6 -4 -2 q2@radD Appendix2.nb 4 Figure 3 Figure 4 t5= (Evaluate[theta[t]]/.sol/.t->.5)[[1]]; Print["r[.5] = ",r5, " m"]; Print["theta[.5] = ",t5, " rad"]; Method I: cartesian coordinates rC = 9 1€€€€ 2 L Cos@theta@tDD, 1€€€€ 2 L Sin@theta@tDD, 0= rA = 8Cos@theta@tDD r@tD, r@tD Sin@theta@tDD, 0< rP = 8p Cos@theta@tDD, p Sin@theta@tDD, 0< F21 = 8−f21 Sin@theta@tDD, f21 Cos@theta@tDD, 0< F21 joint forces at Q H f21 and p unknowns L IO alpha − rC x G1 − rP x F21 = 0 , H1L 1 €€€€ 6 H−6 f21 p + 3 g L m Cos@theta@tDD + 2 L2 m theta′′@tDL == 0 aA = d^2HrALêdt^2 =8−2 Sin@theta@tDD r′@tD theta′@tD + Cos@theta@tDD r′′@tD + r@tD H−Cos@theta@tDD theta′@tD2 − Sin@theta@tDD theta′′@tDL, 2 Cos@theta@tDD r′@tD theta′@tD + Sin@theta@tDD r′′@tD + r@tD H−Sin@theta@tDD theta′@tD2 + Cos@theta@tDD theta′′@tDL, 0< m2 aA = − F21 + G2 =>HxL: m2 aAx + F21x = 0 H2L −f21 Sin@theta@tDD + m H−2 Sin@theta@tDD r′@tD theta′@tD + Cos@theta@tDD r′′@tD + r@tD H−Cos@theta@tDD theta′@tD2 − Sin@theta@tDD theta′′@tDLL == 0HyL: m2 aAy + F21y − G2 = 0 H3L g m + f21 Cos@theta@tDD + m H2 Cos@theta@tDD r′@tD theta′@tD + Sin@theta@tDD r′′@tD + r@tD H−Sin@theta@tDD theta′@tD2 + Cos@theta@tDD theta′′@tDLL == 0 IA alpha − HrP−rAL x H−F21L = 0 , H4L f21 p − f21 r@tD + IA theta′′@tD == 0 From Eq.H4L => p = f21 r@tD − IA theta′′@tD €€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€€ f21 From Eq.H2L => f21 = m Csc@theta@tDD H−2 Sin@theta@tDD r′@tD theta′@tD + Cos@theta@tDD r′′@tD + r@tD H−Cos@theta@tDD theta′@tD2 − Sin@theta@tDD theta′′@tDLL From Eqs. H1L and H3L => two ODE rt_cartesian.nb 2 1 €€€€ 6 H3 g L m Cos@theta@tDD + 2 L2 m theta′′@tD − 6 H−IA theta′′@tD + m Csc@theta@tDD r@tD H−2 Sin@theta@tDD r′@tD theta′@tD + Cos@theta@tDD r′′@tD + r@tD H−Cos@theta@tDD theta′@tD2 − Sin@theta@tDD theta′′@tDLLLL == 0 g m + m H2 Cos@theta@tDD r′@tD theta′@tD + Sin@theta@tDD r′′@tD + r@tD H−Sin@theta@tDD theta′@tD2 + Cos@theta@tDD theta′′@tDLL + m Cot@theta@tDD H−2 Sin@theta@tDD r′@tD theta′@tD + Cos@theta@tDD r′′@tD + r@tD H−Cos@theta@tDD theta′@tD2 − Sin@theta@tDD theta′′@tDLL == 0 numerical application 5. Cos@theta@tDD + r@tD H2. r′@tD theta′@tD − 1. Cot@theta@tDD r′′@tDL + 1.33333 theta′′@tD + r@tD2 H1. Cot@theta@tDD theta′@tD2 + 1. theta′′@tDL == 0 10. + 0. Cos@theta@tDD r′@tD theta′@tD + 1. Csc@theta@tDD r′′@tD == 1. Csc@theta@tDD r@tD theta′@tD2 0.2 0.4 0.6 0.8 1 t@sD 0.2 0.4 0.6 0.8 1.2 1.4 r@mD 0.2 0.4 0.6 0.8 1 t@sD -1.5 -1.25 -1 -0.75 -0.5 -0.25 theta@radD … Graphics … r@.5D = 0.0971902 m theta@.5D = −0.444956 rad rt_cartesian.nb 3 (* RT *) Apply [Clear,Names["Global`*"]]; Off[General::spell]; Off[General::spell1]; "Method II: polar coordinates " rC= {L/2, 0, 0} ; rA= {r[t], 0, 0} ; rP= {p, 0, 0} ; Print["rC = ",rC]; Print["rA = ",rA]; Print["rP = ",rP]; F21= {0, f21, 0} ; Print["F21 = ", F21]; " F21 joint forces at Q ( f21 and p unknowns ) " G1= {-m g Sin[theta[t]], -m g Cos[theta[t]], 0}; G2= {-m g Sin[theta[t]], -m g Cos[theta[t]], 0}; Print["G1 = G2 = ", G1]; alpha= {0, 0, theta''[t]}; IO=m L^2/3; "IO alpha - rC x G1 - rP x F21 = 0 , (1)" e1=IO alpha-Cross[rC,G1]-Cross[rP,F21]; e1z=Simplify[e1[[3]]]==0 aAr=r''[t] - r[t] (theta'[t])^2; aAt=r[t]theta''[t] + 2 r'[t] theta'[t]; aA={aAr,aAt,0}; Print["aA = ", aA]; "m2 aA = - F21 + G2 => " e2= m aA + F21 - G2; "(r): Equation (2) on r" e2x=e2[[1]]==0 "(t): Equation (3) on t" e2y=e2[[2]]==0 "IA alpha - (rP-rA) x (-F21) = 0 , (4)" e2=IA alpha-Cross[(rP-rA),-F21]; e2z=Simplify[e2[[3]]]==0 "From Eq.(4) => " solP=Solve[e2z,p]; ps=p/.solP[[1]]; Print["p = ",ps]; "From Eq.(3) => " solf=Solve[e2y,f21]; f21s=f21/.solf[[1]]; Print["f21 = ",f21s]; "From Eqs. (1) and (3) => two ODE " eq1=e1z/.solP[[1]]/.solf[[1]] eq2=e2x/.solP[[1]]/.solf[[1]] "numerical application" rule={m→1.,L→1,IA→1.,g→10.}; equation1=Simplify[eq1/.rule] equation2=Simplify[eq2/.rule] sol=NDSolve[{equation1,equation2,r[0]Š.1,theta[0]Š.1, r'[0]Š0.,theta'[0]Š0.},{r,theta},{t,0.,1.}]; Plot[Evaluate[r[t]]/.sol,{t,0.,1.},PlotRange→All, AxesLabel→{"t[s]","r[m]"}]; Plot[Evaluate[theta[t]]/.sol,{t,0.,1.},PlotRange→All, rt_polar.nb 1 Figure 5 Figure 6 O0 C1 q A B C2  ı 1 2 Figure 7 y x 0 Figure 10 O 1 A B Figure 11 1 A B α β 0 1 r B2 0 α A Figure 13 Figure 16 0 O A 1 x β 0 O A 1 β G 01F β Figure 17 0 O A 1 B 2 C 0 A 1 B C A B 2 O β β β 21F 01F 12F 1G 2G x r Figure 18 O0 A  ı 2 y x 0 B C φ O0 A 1 2 0 B C 3 φ C2 F32 01y F12 G2 Figure 21 A α 0 1 2 B Figure 22 1 20 Figure 23 0 O A 1 C 0 B a b D k