Week 2 Air standard Cycle
PV DIAGRAM OF AIR STANDARD OTTO CYCLE
for plotting the PV diagram and obtaining values of state variables at all four points, we will take a case study of a standard gasoline engine.
ENGINE PARAMETERS (SI units):
- bore = 0.1 m
- stroke length = 0.1 m
- connecting rod length = .15 m
- compression ratio = 12
along with engine parameters, state variables at initial state are defined. and values of state variables at successive stages are defined with various calculations performed in Matlab programming.
- initial pressure (p1): atmospheric pressure
- initial volume (v1): total volume of cylinder (swept volume + clearance volume)
- initial temperature (t1): 500 K
INPUT:
# PV diagram for air standard otto cycle
import math as M
import numpy as np
import matplotlib.pyplot as plt
# engine parameters (all dimensions are in SI system)
bore = 0.1
stroke = 0.1
con_rod = 0.15
comp_ratio = 12
# calculating swept and clearance volume
v_swept = (M.pi/4)*(pow(bore,2))*stroke
v_clearence = v_swept/(comp_ratio - 1) #comp_ratio = (v_swept + v_clearance)/(v_clearence)
# initial state values
p1 = 101325 # pressure at state 1 in N/m^2
t1 = 500 # temperature at state 1 in kelvin
G = 1.4 # adiabatic constatnt for air
v1 = v_swept + v_clearence # volume at state 1 in m^3
# calculating state variables at point2 (after isentropic compression)
v2 = v_clearence
p2 = p1*pow((v1/v2),G) # p1(v1^G) = p2(v2^G); law of adiabatic process
t2 = (t1*p2*v2)/(p1*v1) # treating air as an ideal gas; for same amount of air: (p*v)/t is constant
# calculating state variables at point3(after constant volume heat addition
v3 = v2 # volume is constant during entire process
t3 = 2300 # temperature is measured practically
p3 = (p2*t3)/t2 # following the gay-lussaic's law; p3/t3 = p2/t2
# calculating state variables at point4 (after isentropic expansion)
v4 = v1 # total(initial) volume of cylinder
p4 = p3*pow((v3/v4),G) # p4(v4^G) = p3(v3^G); law of reversible adiabatic(isentropic process)
# plotting all four state points on PV diagram
plt.figure()
plt.plot([v1,v2,v3,v4],[p1,p2,p3,p4],'o',color='black')
plt.axvline(x=0.0000714,color= 'red',label='TDC',linestyle='dashed')
plt.axvline(x=0.0008568,color= 'green',label='BDC',linestyle='dashed')
plt.xlabel('volume in m^3')
plt.ylabel('pressure in N/m^2')
plt.title('PV diagram of air standard otto cycle')
plt.legend()
plt.grid()
plt.show()
OUTPUT:
above plotting shows only representation of four process transition points. to plot the actual processes, we have to join those four points according to the process taking place between those points.
to plot the accurate graph, we have to determine relations between both variables for those particular processes.
The following MatLab program further calculates the values of each variable for states between transition points and plots the actual values between them.
PLOTTING THERMODYNAMIC PROCESSES BETWEEN TRANSITION STATE POINTS
for plotting each thermodynamic process, we have to calculate values of pressure and volume at each successive point between all four transition state points. for doing so, we will derive relation between pressure and volume for reversible adiabatic and isochoric processes.
PLOTTING ISENTROPIC PROCESS ( REVERSIBLE ADIABATIC )
for plotting isentropic process, a function is defined by the Matlab program which calculates the value of piston volume for each crank angle the engine rotates.
1) VOLUME - CRANK ANGLE RELATION
INPUT:
def cylinder_volume(bore,stroke,con_rod,comp_ratio,Cstart_angle,Cend_angle):
C_angle = np.linspace(Cstart_angle,Cend_angle,120) # defining crank angle
C_radius = stroke/2 # relation between crank angle and stroke
R = con_rod/C_radius # defining ratio used in the formula
v_swept = (M.pi/4)*pow(bore,2)*stroke
v_clearance = v_swept/(comp_ratio-1)
# putting volume-crank angle relation in python format
v_cylinder = [ ]
for i in range(0,120):
part1 = 0.5*(comp_ratio - 1)
part2 = R+1-M.cos(np.radians(C_angle[i]))
part3 =pow(pow(R,2)-pow(M.sin(np.radians(C_angle[i])),2),0.5)
# final equation calculating cylinder volume at crank angle
V = (1 + part1*(part2 - part3)) * v_clearance
v_cylinder.append(V)
return(v_cylinder)
2) JOINING STATE POINTS
INPUT:
# process plotting for PV diagram
import math as M
import numpy as np
import matplotlib.pyplot as plt
# engine parameters (all dimensions are in SI system)
bore = 0.1
stroke = 0.1
con_rod = 0.15
comp_ratio = 12
# calculating swept and clearance volume
v_swept = (M.pi/4)*(pow(bore,2))*stroke
v_clearence = v_swept/(comp_ratio - 1) #comp_ratio = (v_swept + v_clearance)/(v_clearence)
# initial state values
p1 = 101325 # pressure at state 1 in N/m^2
t1 = 500 # temperature at state 1 in kelvin
G = 1.4 # adiabatic constatnt for air
v1 = v_swept + v_clearence # volume at state 1 in m^3
# calculating state variables at point2 (after isentropic compression)
v2 = v_clearence
p2 = p1*pow((v1/v2),G) # p1(v1^G) = p2(v2^G); law of adiabatic process
t2 = (t1*p2*v2)/(p1*v1) # treating air as an ideal gas; for same amount of air: (p*v)/t is constant
# calculating state variables at point3(after constant volume heat addition
v3 = v2 # volume is constant during entire process
t3 = 2300 # temperature is measured practically
p3 = (p2*t3)/t2 # following the gay-lussaic's law; p3/t3 = p2/t2
# calculating state variables at point4 (after isentropic expansion)
v4 = v1 # total(initial) volume of cylinder
p4 = p3*pow((v3/v4),G) # p4(v4^G) = p3(v3^G); law of reversible adiabatic(isentropic process)
# volume-crank angle function
def cylinder_volume(bore,stroke,con_rod,comp_ratio,Cstart_angle,Cend_angle):
C_angle = np.linspace(Cstart_angle,Cend_angle,120) # defining crank angle
C_radius = stroke/2 # relation between crank angle and stroke
R = con_rod/C_radius # defining ratio used in the formula
v_swept = (M.pi/4)*pow(bore,2)*stroke
v_clearance = v_swept/(comp_ratio-1)
# putting volume-crank angle relation in python format
v_cylinder = [ ]
for i in range(0,120):
part1 = 0.5*(comp_ratio - 1)
part2 = R+1-M.cos(np.radians(C_angle[i]))
part3 =pow(pow(R,2)-pow(M.sin(np.radians(C_angle[i])),2),0.5)
# final equation calculating cylinder volume at crank angle
V = (1 + part1*(part2 - part3)) * v_clearance
v_cylinder.append(V)
return(v_cylinder)
# plotting isentropic compression processes
v_comp = cylinder_volume(bore,stroke,con_rod,comp_ratio,180,0) # calling functon for cylinder volume
v_exp = cylinder_volume(bore,stroke,con_rod,comp_ratio,0,180) # calling functon for cylinder volume
p_comp = [ ]
p_exp = [ ]
for i in range(0,120):
p_c = p1*pow((v1/v_comp[i]),G)
p_comp.append(p_c)
# plotting isentropic expansion processes
p_e = p3*pow((v3/v_exp[i]),G)
p_exp.append(p_e)
# plotting all four state points on PV diagram
plt.figure()
plt.plot (v_comp,p_comp); # plotting isentropic compresion process
plt.plot ([v2,v3],[p2,p3]); # plotting isochoric heat addition process
plt.plot (v_exp,p_exp); # plotting isentropic expansion process
plt.plot ([v1,v4],[p1,p4]); # plotting isochoric heat rejection process
plt.xlabel('volume in m^3')
plt.ylabel('pressure in N/m^2')
plt.title('PV diagram of air standard otto cycle')
plt.legend({'isentropic compression','isochoric heat addition','isentropic expansion','isochoric heat rejection'})
plt.show()
OUTPUT:
THERMAL EFFICIENCY OF AIR STANDARD OTTO CYCLE:
as we all know, thermal efficiency of an air standard otto cycle is purely a function of compression ratio and adiabatic constant.
Notto = 1 - 1 / (comp_ratio)^(Gamma-1)
with the help of following python program, we will plot thermal efficiency against compression ratio for various values of adiabatic constant gamma.
INPUT:
import math as M
import numpy as np
import matplotlib.pyplot as plt
comp_ratio = np.arange(1,50,0.5)
G = [1.4 ,1.41 ,1.66]
def efficiency_otto(comp_ratio,G):
comp_ratio = np.arange(1,50,0.5) # range of compression ratio
G = [1.4 ,1.41 ,1.66] # array of gamma for three gases
N1 = 1 - 1 / pow((comp_ratio),(G[0]-1)) # efficiency for air standard otto cycle
N2 = 1 - 1 / pow((comp_ratio),(G[1]-1)) # efficiency for hydrogen otto cycle
N3 = 1 - 1 / pow((comp_ratio),(G[2]-1)) #efficiency for helium otto cycle
# plotting section
plt.plot(comp_ratio,N1) # plotting air standard efficiency curve
plt.plot(comp_ratio,N2) # plotting hydrogen efficiency curve
plt.plot(comp_ratio,N3) # plotting helium efficiency curve
plt.xlabel('compression ratio')
plt.ylabel('Thermal Efficiency')
plt.title('EFFICIENCY CURVE FOR OTTO CYCLE')
plt.legend({'air standard curve','helium curve','hydrogen curve'})
plt.show()
a = efficiency_otto(comp_ratio,G)
OUTPUT: