CHAPTER 7:INTRODUCTION TO THREE-DIMENSIONAL DYNAMICS OF RIGID BODIES

SAMPLE PROBLEM 7/1,PAGE NUMBER:519

In [1]:
import math
from numpy.linalg import det
# Variable declaration
L=0.8;# m
N=60;# rev/min
betadot=4;# rad/s
beta=30;# degree

# Solution
# (a)
omega_x=betadot;# (i) rad/s
omega_z=(2*math.pi*N/60);# (k) rad/s
omega=[omega_x,0,omega_z];# (i,j,k) rad/s
print"(a)The angular velocity of OA,omega=%1.0fi+%1.2fk rad/s"%(omega[0],omega[2]);
# (b)
omegadot_z=0;# (k) rad/s
omegadot_x=omega_z*omega_x;# (i) rad/s
alpha=omegadot_x+omegadot_z;# (j) rad/s**2
alpha=[0,alpha,0];# (i,j,k) rad/s**2
print"(b)The angular acceleration of OA,alpha=%2.1fj rad/s**2"%(alpha[1]);
# (c)
r=[0,0.693,0.4];# m
# v=omega*r;
v_1=det([(omega[0],omega[2]),(r[1],r[2])]);# m/s
v_2=-det([(omega[0],omega[2]),(r[0],r[2])]);# m/s
v_3=det([(omega[0],omega[1]),(r[0],r[1])]);# m/s
v=[v_1,v_2,v_3];# m/s
print"(c)The velocity of point A,v=%1.2fi+(%1.2f)j+%1.2fk m/s"%(v[0],v[1],v[2]);
# (d)
a_1=det([(alpha[1],alpha[2]),(r[1],r[2])])+det([(omega[1],omega[2]),(v[1],v[2])]);# m/s**2
a_2=-det([(alpha[0],alpha[2]),(r[0],r[2])])+(-det([(omega[0],omega[2]),(v[0],v[2])]));# m/s**2
a_3=det([(alpha[0],alpha[1]),(r[0],r[1])])+det([(omega[0],omega[1]),(v[0],v[1])]);# m/s**2
a=[a_1,a_2,a_3];# m/s**2
print"(d)The acceleration of point A,v=%2.1fi+(%2.1f)j+(%1.2f)k m/s**2"%(a[0],a[1],a[2]);
(a)The angular velocity of OA,omega=4i+6.28k rad/s
(b)The angular acceleration of OA,alpha=25.1j rad/s**2
(c)The velocity of point A,v=-2.75i+(-1.60)j+2.77k m/s
(d)The acceleration of point A,v=20.1i+(-28.4)j+(-6.40)k m/s**2

SAMPLE PROBLEM 7/2,PAGE NUMBER:520

In [2]:
import math
from numpy.linalg import det
# Variable declaration
N_0=120;# rev/min
N=60;# rev/min
gamma=30;# degree
OCbar=10;# inch
CAbar=5;# inch
theta=30;# degree

# Calculation
# (a)
omega_0=(2*math.pi*N_0)/60;# rad/sec
omega_1=(2*math.pi*N)/60;# rad/sec
omega=[0,(omega_1*math.cos(gamma*math.pi/180)),(omega_0+(omega_1*math.sin(theta*math.pi/180)))];# rad/sec
print"(a)The angular velocity,omega=%1.2fj+%2.2fk rad/s"%(omega[1],omega[2]);
alpha=[(omega_1*omega_0*math.cos(theta*math.pi/180)),0,0];# rad/sec**2
print"(b)The angular acceleration,alpha=%2.1fi rad/s**2"%alpha[1];   
r=[0,5,10];# inch
# (c)
# v=omega*r;
v_1=det([(omega[1],omega[2]),(r[1],r[2])]);# in/sec
v_2=-det([(omega[0],omega[2]),(r[0],r[2])]);# in/sec
v_3=det([(omega[0],omega[1]),(r[0],r[1])]);# in/sec
v=[v_1,v_2,v_3];# in/sec
print"(c)The velocity of point A,v=%2.1fi+(%1.0f)j+%1.fk in/sec"%(v[0],v[1],v[2]);
# a=(alpha*r)+(omega*v)
a_1=det([(alpha[1],alpha[2]),(r[1],r[2])])+det([(omega[1],omega[2]),(v[1],v[2])]);# in/sec**2
a_2=-det([(alpha[0],alpha[2]),(r[0],r[2])])+(-det([(omega[0],omega[2]),(v[0],v[2])]));# in/sec**2
a_3=det([(alpha[0],alpha[1]),(r[0],r[1])])+det([(omega[0],omega[1]),(v[0],v[1])]);# in/sec**2
a=[a_1,a_2,a_3];# in/sec**2
print"   The acceleration of point A,a=%1.0fi+(%1.0f)j+%3.0fk in/sec**2"%(a[0],a[1],a[2]);
(a)The angular velocity,omega=5.44j+15.71k rad/s
(b)The angular acceleration,alpha=0.0i rad/s**2
(c)The velocity of point A,v=-24.1i+(-0)j+0k in/sec
   The acceleration of point A,a=0i+(-1063)j+473k in/sec**2

SAMPLE PROBLEM 7/3,PAGE NUMBER:530

In [3]:
from numpy.linalg import norm
from scipy.optimize import fsolve
# Variable declaration
omega_1=6;# rad/s
r_x=50;# mm
r_y=100;# mm
r_z=100;# mm

# Calculation
# v_A=r_x*omega_2;
v_B=r_y*omega_1;# (i) mm/s
# v_A=v_B+(omega_n*r_A/B);
# Expanding the determinant and equating the coefficients of the i, j, k terms give
def equations(p):
    omega_nx,omega_ny,omega_nz,omega_2=p
    return((-6-(omega_ny-omega_nz)),(omega_2-((-2*omega_nx)+omega_nz)),(0-((2*omega_nx)-omega_ny)),(((r_x*omega_nx)+(r_y*omega_ny)+(r_z*omega_nz))))
omega_nx,omega_ny,omega_nz,omega_2=fsolve(equations,(1,1,1,1))
omega_n=[omega_nx,omega_ny,omega_nz];# rad/s
omega_n=norm(omega_n);# rad/s
print"\nThe angular velocity of crank DA,omega_2=%1.0f rad/s\nThe angular velocity of link AB,omega_n=%1.3f rad/s"%(omega_2,omega_n);
The angular velocity of crank DA,omega_2=6 rad/s
The angular velocity of link AB,omega_n=4.472 rad/s

SAMPLE PROBLEM 7/4,PAGE NUMBER:531

In [4]:
import math
from numpy.linalg import norm
from scipy.optimize import fsolve
# Variable declaration
# From sample problem 7/3
omega_1=6;# rad/s
omega_2=6;# rad/s
r_x=50;# mm
r_y=100;# mm
r_z=100;# mm
omega_n=2*math.sqrt(5);# rad/s

# Calculation
r_AB=[r_x,r_y,r_z];# mm
# a_A=[r_x*omega_2**2]i+[r_x*omegadot]j;
# a_B=[r_y*omega_1**2]k+[0]i;
omegadot=[(omega_n**2*r_AB[0]),(omega_n**2*r_AB[1]),(omega_n**2*r_AB[2])];# rad/s**2
# omegadot*r_A/B=(100*omegadot_ny-100*omegadot_nz)i+(50*omegadot_nz-100*omegadot_nx)j+(100*omegadot_nx-50omegadot_ny)k
def equations(p):
    omegadot_nx,omegadot_ny,omegadot_nz,omegadot_2=p
    return((28-(omegadot_ny-omegadot_nz)),((omegadot_2+40)-((-2*omegadot_nx)+omegadot_nz)),(-32-((2*omegadot_nx)-omegadot_ny)),(((2*omegadot_nx)+(4*omegadot_ny)+(4*omegadot_nz))))
omegadot_nx,omegadot_ny,omegadot_nz,omegadot_2=fsolve(equations,(1,1,1,1))
omegadot_nv=[omegadot_nx,omegadot_ny,omegadot_nz];# rad/s**2
omegadot_n=norm(omegadot_nv);# rad/s**2
print"The angular acceleration of crank AD,omegadot_2=%2.0f rad/s\nomegadot_n=%1.0fi+%1.0fj+(%1.0f)z rad/s**2\nThe angular acceleration of link AB,omegadot_n=%2.2f rad/s"%(omegadot_2,omegadot_nv[0],omegadot_nv[1],omegadot_nv[2],omegadot_n);
The angular acceleration of crank AD,omegadot_2=-36 rad/s
omegadot_n=-8i+16j+(-12)z rad/s**2
The angular acceleration of link AB,omegadot_n=21.54 rad/s

SAMPLE PROBLEM 7/5,PAGE NUMBER:532

In [5]:
import math
from numpy.linalg import det,norm
# Variable declaration
omega=3;# rad/s
p=8;# rad/s
gamma=30;# degree
y=0.300;# m
z=0.120;# m

# Calculation
# Velocity
omega=[0,0,3];# rad/s
r_B=[0,0.350,0];# m
v_B1=det([(omega[1],omega[2]),(r_B[1],r_B[2])]);# m/s
v_B2=-det([(omega[0],omega[2]),(r_B[0],r_B[2])]);# m/s
v_B3=det([(omega[0],omega[1]),(r_B[0],r_B[1])]);# m/s
v_B=[v_B1,v_B2,v_B3];# m/s
# Note that k*i=J=jcos(gamma*math.pi/180)-ksin(gamma),K*j=-i*cos(gamma*math.pi/180) and K*k=i*sin(gamma)
r_AB=[0,y,z];# m
# omega*r_AB=3K*(yj+zk);
omegaintor_AB=(-(omega[2]*(y*math.cos(gamma*math.pi/180))))+(omega[2]*(z*math.sin(gamma*math.pi/180)));# m/s
p=[0,8,0];# rad/s
v_rel1=det([(p[1],p[2]),(r_AB[1],r_AB[2])]);# m/s
v_rel2=-det([(p[0],p[2]),(r_AB[0],r_AB[2])]);# m/s
v_rel3=det([(p[0],p[1]),(r_AB[0],r_AB[1])]);# m/s
v_rel=[v_rel1,v_rel2,v_rel3];# m/s
v_A=v_B[0]+omegaintor_AB+v_rel[0];# m/s
print"The velocity of point A,v_A=%0.4fi m/s"%v_A;
# Acceleration
a_B1=det([(omega[1],omega[2]),(v_B[1],v_B[2])]);# m/s**2
a_B2=-det([(omega[0],omega[2]),(v_B[0],v_B[2])]);# m/s**2
a_B3=det([(omega[0],omega[1]),(v_B[0],v_B[1])]);# m/s**2
a_B=[a_B1,a_B2,a_B3];# m/s**2
a_B=[0,((a_B[1]*(math.cos(gamma*math.pi/180)))),-(a_B[1]*(math.sin(gamma*math.pi/180)))];# m/s**2
omegadot=0;# m/s**2
# Assume O=omega*(omega*r_A/B)
O=[0,((omega[2]*omegaintor_AB*(math.cos(gamma*math.pi/180)))),-omega[2]*(omegaintor_AB*(math.sin(gamma*math.pi/180)))];# m/s**2
# Assume O_1=2*omega*v_rel
O_1=[0,((2*omega[2]*v_rel[0]*(math.cos(gamma*math.pi/180)))),-2*omega[2]*(v_rel[0]*(math.sin(gamma*math.pi/180)))];# m/s**2
a_rel1=det([(p[1],p[2]),(v_rel[1],v_rel[2])]);# m/s**2
a_rel2=-det([(p[0],p[2]),(v_rel[0],v_rel[2])]);# m/s**2
a_rel3=det([(p[0],p[1]),(v_rel[0],v_rel[1])]);# m/s**2
a_rel=[a_rel1,a_rel2,a_rel3];# m/s**2
a_A=[(a_B[0]+(omegadot*r_AB[0])+O[0]+O_1[0]+a_rel1),(a_B[1]+(omegadot*r_AB[1])+O[1]+O_1[1]+a_rel2),(a_B[2]+(omegadot*r_AB[2])+O[2]+O_1[2]+a_rel3)];# m/s**2
a_A=norm(a_A);# m/s**2
print"The acceleration of point A,a_A=%1.2f m/s"%a_A;
# Angular Acceleration
# Note that k*i=J=jcos(gamma*math.pi/180)-ksin(gamma),K*j=-i*cos(gamma*math.pi/180) and K*k=i*sin(gamma)
omega=[3,8];# rad/s (K,j)(k*j=-i*cos(gamma*math.pi/180))
alpha=(0+(-omega[0]*omega[1]*math.cos(gamma*math.pi/180)));# (i) rad/s**2
print"The angular acceleration of the disk,alpha=%2.1fi rad/s**2"%alpha;
The velocity of point A,v_A=-0.6894i m/s
The acceleration of point A,a_A=8.12 m/s
The angular acceleration of the disk,alpha=-20.8i rad/s**2

SAMPLE PROBLEM 7/6,PAGE NUMBER:544

In [6]:
import math
# Variable declaration
m=70;# The mass of bent plate in kg
omega=30;# rad/s
x_A=0.125;# m
y_A=0.100;# m
x_B=0.075;# m
y_B=.150;# m
d_x=0.0375;# m
d_y=0.125;# m
d_z=0.075;# m

# Calculation
# Part A
m_A=x_A*y_A*m;# kg
m_B=x_B*y_B*m;# kg
I_xxA=((m_A/12)*(y_A**2+x_A**2))+(m_A*((x_A/2)**2+(y_A/2)**2));# kg.m**2
I_yyA=(m_A/3)*(y_A)**2;# kg.m**2
I_zzA=(m_A/3)*(x_A)**2;# kg.m**2
I_xyA=0;# kg.m**2
I_xzA=0;# kg.m**2
I_yzA=0+(m_A*(x_A/2)*(y_A/2));# kg.m**2
# Part B
I_xxB=((m_B/12)*(y_B**2))+((m_B)*(d_y**2+d_z**2));# kg.m**2
I_yyB=((m_B/12)*(x_B**2+y_B**2))+(m_B*(d_x**2+d_z**2));# kg.m**2
I_zzB=((m_B/12)*(x_B**2))+(m_B*((x_A)**2+(d_x**2)));# kg.m**2
I_xyB=0+(m_B*d_x*d_y);# kg.m**2
I_xzB=0+(m_B*d_x*d_z);# kg.m**2
I_yzB=0+(m_B*d_y*d_z);# kg.m**2
I_xx=I_xxA+I_xxB;# kg.m**2
I_yy=I_yyA+I_yyB;# kg.m**2
I_zz=I_zzA+I_zzB;# kg.m**2
I_xy=I_xyA+I_xyB;# kg.m**2
I_xz=I_xzA+I_xzB;# kg.m**2
I_yz=I_yzA+I_yzB;# kg.m**2
# (a)
H_o=[-(omega*I_xz),-(omega*I_yz),(omega*I_zz)];# The angular momentum of the body in N.m.s
# (b)
T=((omega)*(H_o[2]))/2;#(k.i=0,k.j=0,k.k=1) The kinetic energy in J
print"(a)The angular momentum H of the plate about point O,H_O=%0.4fi+(%0.4f)j+%0.4fk\n(b)The kinetic energy of the plate,T=%1.2f J"%(H_o[0],H_o[1],H_o[2],T);
(a)The angular momentum H of the plate about point O,H_O=-0.0664i+(-0.3035)j+0.5502k
(b)The kinetic energy of the plate,T=8.25 J

SAMPLE PROBLEM 7/8,PAGE NUMBER:566

In [7]:
import math
# Variable declaration
m=1000;# The mass of turbine rotor in kg
k=0.200;# m
N=500;# rev/min
rho=400;# The radius of gyration in m
v=25*0.514;# m/s
d_AG=0.6;# m
d_GB=0.9;# m
d_AB=d_AG+d_GB;# m
g=9.81;# The acceleration due to gravity in m/s**2

# Calculation
# The moment principle from statics easily gives
W=m*g;# N
R_1=(m*g)*d_AG;# N
R_2=W-R_1;# N
omega=(v/rho);# rad/s
I=m*k**2;# kg-m**2
deltaR=(I*omega*((2*math.pi*N)/60))/d_AB;
R_A=R_1-deltaR;# N
R_B=R_2+deltaR;# N
print"The vertical components of the bearing reactions at A and B,R_A=%4.0f N and R_B=%4.0f N"%(R_A,R_B);
# The answer provided in the textbook is wrong
The vertical components of the bearing reactions at A and B,R_A=5841 N and R_B=3969 N

SAMPLE PROBLEM 7/9,PAGE NUMBER:567

In [8]:
# SAMPLE PROBLEM 7/9
import math
# Variable declaration
t=4;# s
theta=20;# degree
p=(2*math.pi)/4;# rad/s

# Calculation
# (a)
# I_zz=(56/3)*mr**2;
# I_xx=(32/3)*mr**2;
# By using coefficient of I_xx and I_zz
I=56/3;# The moment of inertia
I_0=32/3;# The moment of inertia
costheta=1;# radian
n=I/((I_0-I)*costheta);# The ratio of angular rates
# (b)
tau=((2*math.pi)/p)*abs(((I_0-I)/I)*math.cos(theta*math.pi/180));# s
beta=math.atan((I/I_0)*math.tan(theta));# degree
print"(a)The number of complete cycles,n=%1.2f \n   The minus sign indicates retrograde precession where, in the present case,and p are essentially of opposite sense. Thus, the station will make seven wobbles for every three revolutions.\n(b)The period of precession,tau=%1.3f s"%(n,tau);
(a)The number of complete cycles,n=-3.00 
   The minus sign indicates retrograde precession where, in the present case,and p are essentially of opposite sense. Thus, the station will make seven wobbles for every three revolutions.
(b)The period of precession,tau=3.759 s