Analysis of otto cycle and simultaneously generating P-V plot and calculating thermal efficiency of an otto cycle using python

The code for generating p-v plot for otto cycle is as follows:

# otto cycle simulator
# by aashish rohila

import math
import matplotlib.pyplot as plt

# inputs for temp.(in kelvin) and pressure(in N/m^2)

p1 = 101325*2
t1 = 300
gamma=1.4
t3 = 2400

# geometric parameters of cylinder (in metre), compression ratio = cr

bore = 0.12
stroke = 0.2
con_rod =  0.3
cr = 12
a = stroke/2
R = con_rod/a

# volume calculations, v_s = swept volume, v_c = clearance volume
v_s = (math.pi/4)*pow(bore,2)*stroke
v_c = v_s/(cr-1)
v1 = v_c + v_s

v_comp = []

def engine(sc,ec):
    sc_rad = math.radians(sc)
    ec_rad = math.radians(ec)
    
    n = 50
    dtheta = (ec_rad-sc_rad)/(n-1)
    for i in range (0,n):
        theta = sc_rad+ (i*dtheta)
        term1 =  0.5*(cr-1)
        term2 = R+1-math.cos(theta)
        term3_temp = pow(R,2)-pow(math.sin(theta),2)
        term3 = pow(term3_temp,0.5)
        v_temp_comp = v_c*(1+term1*(term2-term3))
        v_comp.append(v_temp_comp)
    return v_comp
    
v_comp = engine(180,0)   



# p*pow(v,gamma)= constant will  be used to generate pressure array for isentropic compression

p_comp=[]
for vol in v_comp:
    constant = p1*pow(v1,gamma)
    p_temp_comp = constant/pow(vol,gamma)
    p_comp.append(p_temp_comp)


# state point 2
# p*pow(v,gamma)=constant for isentropic compression
# p*v/t = constant for ideal fluid for closed system     
    
         
v2 = v_c
p2 = p1*pow(v1,gamma)/pow(v2,gamma)
rhs = p1*v1/t1
t2 = p2*v2/rhs


# state point 3
# rhs=p1*v1/t1 is valid for all 4 points of otto cycle as it is a closed cycle,so no mass is entering or leaving a system
v3 = v2 = v_c
p3 = rhs*t3/v3

v_exp = []

def engine_exp(sc,ec):
    sc_rad = math.radians(sc)
    ec_rad = math.radians(ec)
    
    n = 60
    dtheta = (ec_rad-sc_rad)/(n-1)
    for i in range (0,n):
        theta = sc_rad+ (i*dtheta)
        term1 =  0.5*(cr-1)
        term2 = R+1-math.cos(theta)
        term3_temp = pow(R,2)-pow(math.sin(theta),2)
        term3 = pow(term3_temp,0.5)
        v_temp_exp = v_c*(1+term1*(term2-term3))
        v_exp.append(v_temp_exp)
    return v_exp
    
v_exp = engine_exp(0,180)   



# p*pow(v,gamma)= constant will  be used to generate pressure array for isentropic compression

p_exp=[]
for vol in v_exp:
    constant = p3*pow(v3,gamma)
    p_temp_exp = constant/pow(vol,gamma)
    p_exp.append(p_temp_exp)



    

# state oint 4
v4 = v1
p4 = p3*pow(v3,gamma)/pow(v4,gamma)
t4 = p4*v4/rhs

plt.plot(v_comp,p_comp)
plt.plot(v_exp,p_exp)
plt.plot((v2,v3),(p2,p3))
plt.plot((v1,v4),(p1,p4))
plt.xlabel(\"volume\")
plt.ylabel(\"pressure\")
plt.show()

# thermal efficiency of an otto cycle = 1-(1/pow(cr,(gamma-1)))
thermal_efficiency = 1-(1/pow(cr,(gamma-1)))
print(thermal_efficiency)

And the output generated for thermal eficiency is as follows:

And the plot for P-V diagram is as follows:

In the above code we have simply used conditions for ideal fluid for isentropic compression as expansion as p(v^gamma)= constant and pv/t=constant for closed system.

we have used above conditions to find various pressure and volume points.

Then we taken an empty array for volume and used for loop and append function to append the array.

we have also used def function to append the volume array.

At the end we have used plt.plot(), plt.xlabel(), plt.ylabel and plt.show() command