import math
#calculate pressure and for longitudinal stress
di=0.3 ##m
de=0.4 ##m
v=0.3
sigmathetamax=250*10**6 ##Pa
p0=0.
pi=0.
##solution a:
a=0.15
b=0.2
r=a
##sigmathetamax=pi*((b**2+a**2)/(b**2-a**2))
pi=sigmathetamax*((b**2-a**2.)/(b**2+a**2.))
print'%s %.2f %s'%("in Pa is= ",pi,"")
##solution b:
r=a
##sigmathetamax=-2*p0*(b**2/(b**2-a**2))
p0=-(-sigmathetamax)*((b**2.-a**2.)/(2.*b**2.))
print'%s %.2f %s'%("in Pa is= ",p0,"")
##solution c:
u=((a**3*pi)/(b**2-a**2))*(0.7+1.3*(b**2./a**2.))
print'%s %.2f %s'%("in per E meter is= ",u,"")
sigmaz=(pi*a**2-p0*b**2)/(b**2-a**2)
print'%s %.2f %s'%("for longitudinal stress is",sigmaz,"")
import math
#calculate pressure
sigmayp=340. ##MPa
tauyp=sigmayp/2. ##MPa
print'%s %.2f %s'%("in MPa is=",tauyp,"")
a=0.1 ##m
b=0.15 ##m
v=0.3
##pi=4*p0
##sigmatheta=(pi*(a**2+b**2)-2*p0*b**2)/(b**2-a**2)
##sigmatheta=1.7*pi
##sloution a: maxi principal stress theory
sigmatheta=1.7
pi=sigmayp/sigmatheta
print'%s %.2f %s'%("in MPa is= ",pi,"")
##sloution b: maxi shearing stress theory
##(sigmatheta-sigmar)/2=1.35*pi
pi=tauyp/1.35
print'%s %.2f %s'%("in MPa is= ",pi,"")
##solution c: energy of distortion theory
sigmar=-1
sigmayp1=math.sqrt(sigmatheta**2+sigmar**2-sigmatheta*sigmar)##*pi
print(sigmayp1)
pi=sigmayp/sigmayp1
print'%s %.2f %s'%("in MPa is=",pi,"")
##solution d: maxi principal strain theory
##(sigmatheta-v*sigmar)/E=sigmayp/E
pi=sigmayp/(sigmatheta-v*sigmar)
print'%s %.2f %s'%("in MPa is= ",pi,"")
##solution e: octahedral shearing stress theory:
pi=(math.sqrt(2.)*sigmayp)/math.sqrt((sigmatheta-sigmar)**2.+sigmar**2.+(-sigmatheta)**2)
print'%s %.2f %s'%("in MPa is= ",pi,"")
import math
#calculate contact pressure and tangential stresses in the outer cylinder
a=0.15 ##m
b=0.2 ##m
c=0.25 ##m
E=200*10**9. ##Pa
delta=0.0001 ##m
140 ##MPa
p=((E*delta)/8.)*(((b**2.-a**2.)*(c**2.-b**2.))/(2.*(b**2.)*(c**2.-a**2.)))
print'%s %.2f %s'%("the contact pressure in Pa is= ",p,"") ## textbook ans is wrong
p=12.3*10**6
sigmatheta=p*((b**2+c**2.)/(c**2.-b**2.)) ## where r=0.2
print'%s %.2f %s'%("tangential stresses in the outer cylinder in Pa is= ",sigmatheta,"")
sigmatheta1=(2*p*b**2)/(c**2-b**2) ## where r=0.25
print'%s %.2f %s'%("tangential stresses in the outer cylinder in Pa is= ",sigmatheta1,"")
sigmatheta3=-(2*p*b**2)/(b**2-a**2) ## where r=0.15
print'%s %.2f %s'%("tangential stresses in the inner cylinder in Pa is= ",sigmatheta3,"")
sigmatheta4=-p*((b**2.+a**2.)/(b**2.-a**2.)) ## where r=0.2
print'%s %.2f %s'%("tangential stresses in the inner cylinder in Pa is= ",sigmatheta4,"")
import math
#calculate radial displacement of disk and shaft
dn=0.1 ##m
do=0.5 ##m
t=0.08 ##m
w=6900*(2*math.pi/60.) ##rpm
row=7.8*10**3##Ns**2/m**4
E=200*10**9 ##Pa
v=0.3
b=0.05
c=0.25
##solution a:
##ud=((0.05*3.3*0.7)*(0.0025+0.0625-(1.3/3.3)*0.0025+(1.3/0.7)*0.0625)*row*w**2)/(8*E)
ud=((0.05*3.3*0.7)*(b**2+c**2-(1.3/3.3)*b**2+(1.3/0.7)*c**2))/(8)
print'%s %.4f %s'%("radial displacement of the disk in meter is= ",ud,"")
##us=((0.05*0.7)*(3.3*0.0025-1.3*0.0025)*row*w**2)/(8*E)
us=((0.05*0.7)*(3.3*b**2-1.3*b**2))/(8)
print'%s %.6f %s'%("radial displacement of the shaft in meter is= ",us,"")
delta=(ud-us)*row*w**2./E
print(delta)
##solution b:
##p=E*delta*(c**2-b**2)/(2*b*c**2)
p=E*delta*(c**2-b**2)/(2*b*c**2)
print'%s %.2f %s'%("in Pa is= ",p,"")
sigmathetamax=p*(c**2+b**2)/(c**2-b**2)
print'%s %.2f %s'%("in Pa is= ",sigmathetamax,"")
##solution c:
sigmathetamax=3.3*(b**2.+c**2.-(1.9/3.3)*b**2.+c**2.)*row*w**2./8.
print'%s %.2f %s'%("in Pa is= ",sigmathetamax,"")
import math
#calculate sigma theta and sigma
ti=0.075 ##m
to=0.015##m
a=0.05##m
b=0.25##m
delta=0.05 ##mm
w=6900*(2*math.pi/60.) ##rpm
s=1.
row=7.8*10**3##Ns**2/m**4
E=200. ##GPa
##solution a:
t1=ti*a**s
print'%s %.4f %s'%("t1 is=",t1,"")
t1=to*b**2.
print'%s %.4f %s'%("t1 is=",t1,"")
##(ti/to)=(t1*a**-s)/(t1*b**-s)=(b/a)**s
c=(b/a)**s
(ti/to)==c
print'%s %.2f %s'%("ti/t0 is=",c,"")
m1=-0.5+math.sqrt((0.5)**2+(1+0.3*1))
print'%s %.2f %s'%("m1 is=",m1,"")
m2=-0.5-math.sqrt((0.5)**2.+(1.+0.3*1.))
print'%s %.2f %s'%("m2 is=",m2,"")
##sigmar=0=(c1/t1)*(0.05)**m1+(c2/t1)*(0.05)**(m2)-0.00176*row*w**2 ## r=0.05
##sigmar=0=(c1/t1)*(0.25)**m1+(c2/t1)*(0.25)**(m2)-0.0439*row*w**2 ## r=0.25
c1=t1*0.12529*row*w**2.
print'%s %.2f %s'%("c1 is=",c1,"")
c2=t1*-6.272*10**-5*row*w**2
print'%s %.2f %s'%("c2 is=",c2,"")
r=0.05
sigmar=(0.12529*r**0.745-6.272*10**-5*r**(-1.745)-0.70*r**2)##*row*w**2
print'%s %.5f %s'%("sigmar is= ",sigmar,"")
sigmatheta=(0.09334*r**0.745+1.095*10**-4*r**(-1.745)-0.40*r**2)##*row*w**2
print'%s %.2f %s'%("sigmatheta is= ",sigmatheta,"")
##solution b:
r=0.05
##ur=(r*sigmatheta)/E
ur=(r*sigmatheta)
print'%s %.7f %s'%("ur is= ",ur,"")
import math
#calculate pressure
b=0.25 ##m
w=6900.*(2*math.pi/60.) ##rpm
t1=0.075 ##m
t2=0.015 ##m
row=7.8*10**3##Ns**2/m**4
c1=t1
x=t2/t1
print(x)
##(t2/t1)==(c1*exp(-(row*w**2/2*sigma)*b**2))/c1
##exp(-(row*w**2/2*sigma)*b**2)=x
##log(x)=-(row*w**2*b**2/2*sigma)
y=2.*math.log(x)
print(y)
sigma=-(row*w**2.*b**2.)/y
print'%s %.2f %s'%("in Pa is= ",sigma,"")
##t=c1*exp(-row*(w**2/2*sigma)*r**2)
z=row*(w**2./(2.*sigma))
print(z)