Air standard Otto Cycle analysis using python
AIM
To analyse Otto cycle and plot its P-V diagram
GOVERNING EQUATIONS
- Ideal gas equation: `(P*V)/T` = constant. This equation is true for ideal gases and since Otto cycle is an air standard cycle, this equation holds true.
Here, P = pressure
T = temperature
V = volume
- For adiabatic process, `P*V^gamma` = constant
`gamma` = `C_P/C_v`
Value of `gamma` for air = 1.4
- Swept volume, `V_(swept)` = `pi/4*d^2*l`
Here, d = bore or diameter of piston cylinder
l = stroke or distance traversed by piston
- Clearance volume = `V_(cl earance)` = `V_(swept)/(r-1)`
- `V/V_c = 1 + 1/2(r-1)[R+1-costheta-(R^2-sin^2theta)^0.5]`
- Efficiency of Otto cycle: `eta` = work done/heat added
= `(Q_(i n) - Q_(out))/Q_(i n) = 1 - Q_(out)/Q_(i n)`
= `1-(m*C_v*(T_4-T_1))/(m*C_v*(T_3-T_2)`
= 1 – `(T_4-T_1)/(T_3-T_2)`
The above equation can be written in the terms of compression ratio, r as
`eta = 1-1/(r^(gamma-1))`
OBJECTIVE
The objective of this program is to analyse the Otto cycle, plot its P-V curve and determine its thermal efficiency. An Otto cyclee is an air standard thermodynamic cycle that describes working of a typical spark ignition piston engine in which a mass of gas (usually air) is subjected to changes in pressure, temperature, volume, heat addition and rejection. There are usually 4 processes involved in a 2 stroke Otto cycle.
The 4 processes involved in an Otto cycle
1) Process 1-2: It is the adiabatic or isentropic (constant entropy) compression of charge in which the piston moves from the Bottom Dead centre (BDC) to Top Dead Centre (TDC). For this process, `P*V^gamma`= constant. Therefore, `P_1*V_1^gamma = P_2*V_2^gamma`
2) Process 2-3: It is isochoric (constant volume) heat addition while piston remains at TDC. It causes fuel ignition and rapid burning take place. Therefore, `V_2 = V_3`.
3) Process 3-4: It is the adiabatic (isentropic) expansion of charge in which piston moves from TDC to BDC. Again, `P*V^gamma`=constant. This implies `P_3*V_3^gamma = P_4*V_4^gamma`.
4) Process 4-1: It is the isochoric heat rejection while the piston remains at BDC. Thus, `V_4 = V_1`.
Given the state variables at point 1 and the temperature after the heat addition process, the state variables at all the states are determined using the suitable governing equations.
The nature of volume change during the adiabatic processes is determined using the governing equation and then the plot between pressure and volume is generated.
At last, thermal efficiency of the cycle is calculated using the formula.
PROGRAM
Initially, the required modules are imported. In order to determine the nature of curve for the adiabatic (isentropic) processes, a function is defined which is then called in the main program to get the accurate plot. A function named engine_kinematics is created with volume (v) as output for various inputs like bore, stroke, connecting rod length, compression ratio, start crank angle and end crank angle. Additional parameters like crank pin radius (a) and R (ratio of connecting rod lngth and crank pin radius) are defined.
An equation is formed to calculate actual volume change using a geometric relation. Whenever this function would be called in the main program, it would calculate the volume (`V_(compression) and V_(expansion)`) in this case.
The required inputs are provided which include the geometric parameters (stroke, bore, connecting rod length), compression ratio and required state variables.
Then state variables for different state points are calculated using suitable equations. After `V_(compresion), P_(compression), V_(expansion) and P_(expansion)` are calculated, P-V diagram is plotted to get the actual P-V curve for Otto cycle using maplotlib module.
Thermal efficiency of the cycle is calculated and displayed.
# Program for analysis of otto cycle
import math
import matplotlib.pyplot as plt
import numpy as num
def engine_kinematics(d,l,cr,r,start_crank,end_crank):
a = l/2 #crank pin radius
R = cr/a
V_s = (math.pi/4)*pow(d,2)*l
V_c = V_s/(r-1)
theta = num.linspace(start_crank,end_crank,180)
v = []
# v/V_c = 1 + 0.5(r-1)(R + 1 - cos(theta) - (R^2 - sin^2(theta))^0.5)
#splitting the terms
for i in theta:
term1 = 0.5*(r-1)
term2 = R + 1 - math.cos(i)
term3 = pow(pow(R,2) - pow(math.sin(i),2),0.5)
v.append((1 + term1*(term2 - term3))*V_c)
return v
# inputs
gamma = 1.4
T3 = 2500
d = 0.2 #bore
l = 0.4 #stroke
cr = 0.25 #connecting rod length (m)
r = 6 #compression ratio
V_s = (math.pi/4)*pow(d,2)*l # Swept volume
V_c = V_s/(r-1) # clearance volume
# State 1
P1 = 101325
T1 = 500
V1 = V_s + V_c
# State 2
P2 = P1*pow(r,gamma)
V2 = V_c
T2 = P2*V2*T1/(P1*V1)
c1 = P1*pow(V1,gamma)
V_compression = engine_kinematics(d,l,cr,r,math.pi,0)
P_compression = []
for vol in V_compression:
P_compression.append(c1/pow(vol,gamma))
# State 3
V3 = V2
P3 = (P2*T3)/T2
c2 = P3*pow(V3,gamma)
V_expansion = engine_kinematics(d,l,cr,r,0,math.pi)
P_expansion = []
for volume in V_expansion:
P_expansion.append(c2/pow(volume,gamma))
# State 4
V4 = V1
P4 = P3/pow(r,gamma)
T4 = P4*V4*T3/(P3*V3)
# Plotting
plt.plot(V1,P1,'*')
plt.text(V1,P1,'1')
plt.plot(V_compression,P_compression,color='red')
plt.plot(V2,P2,'*')
plt.text(V2,P2,'2')
plt.plot([V2,V3],[P2,P3],color='green')
plt.plot(V3,P3,'*')
plt.text(V3,P3,'3')
plt.plot(V_expansion,P_expansion,color='blue')
plt.plot(V4,P4,'*')
plt.text(V4,P4,'4')
plt.plot([V4,V1],[P4,P1],color='yellow')
plt.xlabel('Volume')
plt.ylabel('Pressure')
plt.title('P-V diagram')
plt.grid('on')
plt.show()
"""
Efficiency
eta = work done/ heat added
eta = (heat added - heat rejected)/ Heat added
"""
eta = 1 - (T4 - T1)/(T3 - T2)
print('Efficiency of cycle, Eta =',eta)
Output
CONCLUSION
The P-V diagram for the otto cycle is generated and the thermal efficiency for the given state parameters comes out to be 0.5116