(* Problem 4.8 *)(* VELOCITY ANALYSIS - input angle phi *)Apply[Clear,Names["Global`*"]];Off[General::spell];Off[General::spell1];n = 100 ; (* rpm *)omega = n*N[Pi]/30 ; (* rad/s *)initdata = {AB->0.12, AC->0.06, BD->0.24, DE->0.33,EF->0.19,La->.3, Lb->.07, phi[t]->N[Pi]150/180, phi'[t]->omega, phi''[t]->0};(* Position of joint A *)xA = yA = 0;(* Position of joint C *)xC = -AC ;yC = 0 ;(* Position of joint F *)xF = La ;yF = Lb ;(* Position of joint B *)xB = AB Cos[phi[t]] ;yB = AB Sin[phi[t]] ;Print["xB = ", xB ," = ", xB/.initdata, " m" ];Print["yB = ", yB ," = ", yB/.initdata, " m" ];(* Linear velocity of joint B *)vBx = D[xB,t];vBy = D[yB,t];Print["vBx = ", vBx ," = ", vBx/.initdata, " m/s" ];Print["vBy = ", vBy ," = ", vBy/.initdata, " m/s" ];(* Linear acceleration of joint B *)aBx = D[vBx,t];aBy = D[vBy,t];Print["aBx = ", aBx ," = ", aBx/.initdata, " m/s^2" ];Print["aBy = ", aBy ," = ", aBy/.initdata, " m/s^2" ];(* Position of joint D *)eqnD1 = ( xDsol - xB )^2 + ( yDsol - yB )^2 - BD^2 == 0 ;eqnD2 = ( yDsol - yB ) / ( xDsol - xB ) == ( yB - yC ) / ( xB - xC );solutionD = Solve [ { eqnD1 , eqnD2 } , { xDsol , yDsol } ];(* Two solutions for D *)xD1 = xDsol /. solutionD[[1]];yD1 = yDsol /. solutionD[[1]];xD2 = xDsol /. solutionD[[2]];yD2 = yDsol /. solutionD[[2]];If [ (yD1/.initdata) <= 0 , xD = xD1 ; yD = yD1 , xD = xD2 ; yD=yD2 ] ;Print["xD = ",xD/.initdata," m"];Print["yD = ",yD/.initdata," m"];(* Linear velocity of joint D *)vDx = D[xD,t];vDy = D[yD,t];Print["vDx = ", vDx/.initdata, " m/s" ];Print["vDy = ", vDy/.initdata, " m/s" ];(* Linear acceleration of joint D *)aDx = D[vDx,t];aDy = D[vDy,t];Print["aDx = ", aDx/.initdata, " m/s^2" ];Print["aDy = ", aDy/.initdata, " m/s^2" ];(* Angular velocity and acceleration of the link BD (link 2)*)phi2 = ArcTan[ (yB-yC)/(xB-xC) ] ;omega2 = D[ phi2 , t ] ;alpha2 = D[ omega2, t ] ;Print["phi2 = phi3 = ", phi2/.initdata , " rad = ",(phi2/.initdata)*180/N[Pi]," deg" ];Print["omega2 = omega3 = ", omega2/.initdata , " rad/s" ];Print["alpha2 = alpha3 = ", alpha2/.initdata , " rad/s^2" ];(* Position of joint E *)eqnE1 = ( xEsol - xD )^2 + ( yEsol - yD )^2 - DE^2 == 0 ;eqnE2 = ( xEsol - xF )^2 + ( yEsol - yF )^2 - EF^2 == 0 ;solutionE = Solve [ { eqnE1 , eqnE2 } , { xEsol , yEsol } ];(* Two solutions for E *)xE1 = xEsol /. solutionE[[1]];yE1 = yEsol /. solutionE[[1]];xE2 = xEsol /. solutionE[[2]];yE2 = yEsol /. solutionE[[2]];If [ (yE1/.initdata) <= 0 , xE = xE1 ; yE = yE1 , xE = xE2 ; yE=yE2 ] ;Print["xE = ",xE/.initdata," m"];Print["yE = ",yE/.initdata," m"];(* Linear velocity of joint E *)vEx = D[xE,t];vEy = D[yE,t];Print["vEx = ", vEx/.initdata, " m/s" ];Print["vEy = ", vEy/.initdata, " m/s" ];(* Linear acceleration of joint D *)aEx = D[vEx,t];aEy = D[vEy,t];Print["aEx = ", aEx/.initdata, " m/s^2" ];Print["aEy = ", aEy/.initdata, " m/s^2" ];(* Angular velocity and acceleration of the link DE (link 4)*)phi4 = ArcTan[ (yD-yE)/(xD-xE) ] ;omega4 = D[ phi4 , t ] ;alpha4 = D[ omega4, t ] ;Print["phi4 = ", phi4/.initdata , " rad = ",(phi4/.initdata)*180/N[Pi]," deg" ];Print["omega4 = ", omega4/.initdata , " rad/s" ];Print["alpha4 = ", alpha4/.initdata , " rad/s^2" ];(* Angular velocity and acceleration of the link EF (link 5)*)phi5 = ArcTan[ (yE-yF)/(xE-xF) ] ;omega5 = D[ phi5 , t ] ;alpha5 = D[ omega5, t ] ;Print["phi5 = ", phi5/.initdata , " rad = ",(phi5/.initdata)*180/N[Pi]," deg" ];Print["omega5 = ", omega5/.initdata , " rad/s" ];Print["alpha5 = ", alpha5/.initdata , " rad/s^2" ];rB={xB/.initdata,yB/.initdata,0};rC={xC/.initdata,yC/.initdata,0};rBC=rC-rB;Print["rB = ", rB, " m" ];Print["rC = ", rC, " m" ];Print["rBC = rC - rB = ", rBC, " m" ];vB={vBx/.initdata,vBy/.initdata,0};aB={aBx/.initdata,aBy/.initdata,0};Print["vB = ", vB , " m/s" ];Print["aB = ", aB , " m/s" ];w2={0,0,omega2/.initdata};a2={0,0,alpha2/.initdata};Print["w2 = ", w2, " rad/s" ];Print["a2 = ", a2, " rad/s^2" ];vC2=vB+Cross[w2,rBC];Print["vC2 = vB + w2 x rBC = ",vC2, " m/s" ];aC2=aB+Cross[a2,rBC]-w2.w2 rBC;Print["aC2 = aB + a2 x rBC - w2^2 rBC = ", aC2, " m/s^2" ];vC32=-vC2;Print["vC32 = vC3 - vC2 = ",vC32, " m/s" ];aC32cor=2 Cross[w2,vC32];Print["aC32cor = 2 w2 x vC32 = ",aC32cor, " m/s^2" ];aC32=-aC2-aC32cor;Print["aC32 = aC3 - aC2 - aC32cor = ",aC32, " m/s^2" ];rD={xD/.initdata,yD/.initdata,0};rBD=rD-rB;Print["rD = ", rD, " m" ];Print["rBD = rD - rB = ", rBD, " m" ];vD=vB+Cross[w2,rBD];Print["vD = vB + w2 x rBD = ",vD, " m/s" ];aD=aB+Cross[a2,rBD]-w2.w2 rBD;Print["aD = aB + a2 x rBD - w2^2 rBD = ", aD, " m/s^2" ];Problem4_8_I.nb 1(* Problem 4.8 *)(* VELOCITY ANALYSIS - input angle phi *)Apply[Clear,Names["Global`*"]];Off[General::spell];Off[General::spell1];n = 100 ; (* rpm *)omega = n*N[Pi]/30 ; (* rad/s *)initdata = {AB->0.12, AC->0.06, BD->0.24, DE->0.33,EF->0.19,La->.3, Lb->.07, phi[t]->N[Pi]150/180, phi'[t]->omega, phi''[t]->0};(* Position of joint A *)xA = yA = 0;(* Position of joint C *)xC = -AC ;yC = 0 ;(* Position of joint F *)xF = La ;yF = Lb ;(* Position of joint B *)xB = AB Cos[phi[t]] ;yB = AB Sin[phi[t]] ;Print["xB = ", xB ," = ", xB/.initdata, " m" ];Print["yB = ", yB ," = ", yB/.initdata, " m" ];(* Linear velocity of joint B *)vBx = D[xB,t];vBy = D[yB,t];Print["vBx = ", vBx ," = ", vBx/.initdata, " m/s" ];Print["vBy = ", vBy ," = ", vBy/.initdata, " m/s" ];(* Linear acceleration of joint B *)aBx = D[vBx,t];aBy = D[vBy,t];Print["aBx = ", aBx ," = ", aBx/.initdata, " m/s^2" ];Print["aBy = ", aBy ," = ", aBy/.initdata, " m/s^2" ];(* Position of joint D *)eqnD1 = ( xDsol - xB )^2 + ( yDsol - yB )^2 - BD^2 == 0 ;eqnD2 = ( yDsol - yB ) / ( xDsol - xB ) == ( yB - yC ) / ( xB - xC );solutionD = Solve [ { eqnD1 , eqnD2 } , { xDsol , yDsol } ];(* Two solutions for D *)xD1 = xDsol /. solutionD[[1]];yD1 = yDsol /. solutionD[[1]];xD2 = xDsol /. solutionD[[2]];yD2 = yDsol /. solutionD[[2]];If [ (yD1/.initdata) <= 0 , xD = xD1 ; yD = yD1 , xD = xD2 ; yD=yD2 ] ;Print["xD = ",xD/.initdata," m"];Print["yD = ",yD/.initdata," m"];(* Linear velocity of joint D *)vDx = D[xD,t];vDy = D[yD,t];Print["vDx = ", vDx/.initdata, " m/s" ];Print["vDy = ", vDy/.initdata, " m/s" ];(* Linear acceleration of joint D *)aDx = D[vDx,t];aDy = D[vDy,t];Print["aDx = ", aDx/.initdata, " m/s^2" ];Print["aDy = ", aDy/.initdata, " m/s^2" ];(* Angular velocity and acceleration of the link BD (link 2)*)phi2 = ArcTan[ (yB-yC)/(xB-xC) ] ;omega2 = D[ phi2 , t ] ;alpha2 = D[ omega2, t ] ;Print["phi2 = phi3 = ", phi2/.initdata , " rad = ",(phi2/.initdata)*180/N[Pi]," deg" ];Print["omega2 = omega3 = ",