#Given
di=4*12          #inch, diameter
ri=di/2.0        #Radius
t=0.5           #inch, thickness
sigma=20.0     #Ksi,  stress
#Calculation
#Cylindrical Pressure Vessel
p1=(t*sigma)/ri #sigma = pr/t
#Spherical Vessel
p2=(2*t*sigma)/(ri) #sigma = pr/2t
#Display
print"The maximum internal pressure the cylindrical pressure vessel can sustainis",round(p1*1000,0),"psi"
print"The maximum internal pressure a spherical pressure vessel can sustain is",round(p2*1000,0),"psi"
   
#Determine stress at point B and C
#Given
P = 15000.0     #N, force,
a = 40.0        #mm, length
b = 100.0       #mm, breadth
#CAlculation
#Normal Force
A = a*b         #Area
sigma = P/A
#Bending Moment
I = (a*b**3)/12.0 #I = (1/12)*bh**3
M = P*(b/2.0)  
c = b/2.0
sigma_max =(M*c)/I
#Superposition
x = ((sigma_max-sigma)*b)/((sigma_max+sigma)+(sigma_max-sigma))
sigma_b = (sigma_max-sigma)
sigma_c = (sigma_max + sigma)
#Display
print"The state of stress at B  is(tensile)",sigma_b,"psi"
print"The state of stress at C is (compressive)",sigma_c,"psi"
 
#Calculate the stress
#Given:
ri =24       #inch, radius
t = 0.5      #inch
ro = ri+t
sp_wt_water = 62.4    #lb/ft**3
sp_wt_steel = 490     #lb/ft**3
l_a = 3               #m depth of point A from the top
#Internal Loadings:
import math
v = (math.pi*l_a)*((ro/12.0)**2 - (ri/12.0)**2)
W_st = sp_wt_steel*v
p = sp_wt_water*l_a #lb/ft**2,Pascal's Law
p_=p*0.0069         #psi
#Circumferential Stress:
sigma1 = (p_*ri)/t
#Longitudinal Stress:
A_st = (math.pi)*(ro**2 - ri**2)
sigma2 = W_st/A_st
#Display:
print"The state of stress at A (Circumferential)",round(sigma1,0),"KPa"
print"The state of stress at A (Longitudinal) ",round(sigma2,1),"KPa"
#Determine the state of stress
#Given
y_c = 125/1000.0 #m,  length
x_c = 1.5        #m
y_b = 1.5         #m
x_b = 6.0         #m
udl = 50.0       #kN/m, force per unit length
l_udl = 2.5      #m
l = 250/1000.0   #m
width = 50/1000.0 #m 
#Internal Loadings:
N = 16.45 #kN
V = 21.93 #kN
M = 32.89 #kNm
#Stress Components:
#Normal Force:
A = l*width
sigma1 = N/(A*1000)
#Shear Force:
tou_c = 0
#Bending Moment:
c = y_c
I = (1/12.0)*(width*l**3)
sigma2 = (M*c)/(I*1000)
#Superposition:
sigmaC = sigma1+sigma2
#Display:
print"The stress due to normal force at C   ",round(sigma1,2),"MPa"
print"The stress due to shear force at C    ",tou_c,"MPa"
print"The stress due to bending moment at C ",round(sigma2,1),"MPa"
print"The resultant stress at C             ",round(sigmaC,1),"MPa"
#Given:
r = 0.75*10       #mm, radius
f_x =40000        #N, force along x
f_y =800          #N force along y
l1 = 0.8         #mm
l2 = 0.4        #mm
#Stress Components:
#Normal Force:
A1 =l1*l2
sigma1 = f_x/A1   #stress = P/A
#Bending Moment:
M_y1 = 8000        #N
c1 = l2/2.0
I1 = (1/12.0)*(l1*l2**3)
sigma_A1 = (M_y1*c1)/I1 
M_y2 = 16000        #N
c2 = l2
I2 = (1/12.0)*(l2*l1**3)
sigma_A2 = (M_y2*c2)/I2 
#Resultant:
res_normal= -sigma1-sigma_A1-sigma_A2
#Display:
print"The stress due to normal force at A   ",sigma1/1000,"KPa"
print"The stress due to bending moment 8KN at A  ",sigma_A1/1000,"KPa"
print"The stress due to bending moment 16KN at A  ",sigma_A2/1000,"KPa"
print"The resultant normal stress component at A ",res_normal/1000,"KPa"
#Calculate Stress at A
#Given:
P = 500        #lb, load
r=0.75         #inch, radius
#Stress Components:
#Normal Force:
import math
A = math.pi*r**2
sigma = P/A
#Bendng Moments:
M_x =7000     #lb
cy = r
Ix = (1/4.0)*math.pi*(r**4) 
sigma_max_1 = (M_x*cy)/Ix   
M_y = P*l_bc/2.0
cx = l_bc/2.0
Iy = (1/12.0)*(l_ab*l_bc**3) #I = (1/12)*(bh**3)
sigma_max_2 = (M_y*cx)/Iy #sigma = My/I
#Superposition
sigmaf=round(sigma/1000,3)+round(sigma_max_1/1000,1)
#Display:
print"The normal stress at corner A ",round(sigma/1000,3),"ksi"
print"The normal stress at point A for Bending Moment ",round(sigma_max_1/1000,1),"ksi"
print"The normal stress at point A for Superimposition ",round(sigmaf,1),"ksi"
#Determine the state of stress at point A
#Given
r=0.75           #radius,inch
V=800             #Forca, lb
#Calculation
#shear force
import math
Q=(4*r/(3*math.pi))*(0.5*(math.pi*r**2))
Ix=(1/4.0)*math.pi*(r**4) 
tau=V*Q/(Ix*2*r)
#Since point A is on neutral axis
sigmaA=0
T=11200           #lb inch, force 
Iy=(1/2.0)*math.pi*(r**4) 
sigma_a=T*r/Iy
#Superimposition
sigmayzA=tau+sigma_a
#Result
print"The stress for shear stress distribution is",round(tau/1000,3),"ksi"
print"The stress for Bending moment is",sigmaA,"ksi"
print"The stress for torque",round(sigma_a/1000,2),"ksi"
print"The stress for Superimposition ",round(sigmayzA/1000,1),"ksi"