Air standard Cycle using Python

OBJECTIVE:

To write a python code for an Otto cycle using the defined parameters and make plots.

1. To create a P-V diagram (Pressure Vs. Volume diagram)

2. To print the thermal efficiency of the cycle.

PROBLEM STATEMENT:

 An Otto cycle is an idealized thermodynamic cycle that describes the functioning of a typical spark igniton piston engine. It is the thermodynamic cycle most commonly found in automobile engines. The Otto cycle consists of isentropic compression, heat addition at constant volume, isentropic expansion, and rejection of heat at constant volume. This is shown in the figure below

The processes are described below:

Process 0–1: a mass of air is drawn into piston/cylinder arrangement at constant pressure.

Process 1–2: is an adiabatic (isentropic) compression of the charge as the piston moves from bottom dead center (BDC) to top dead center (TDC).

Process 2–3: is a constant-volume heat transfer to the working gas from an external source while the piston is at top dead center. This process is intended to represent the ignition of the fuel-air mixture and the subsequent rapid burning.

Process 3–4: is an adiabatic (isentropic) expansion (power stroke).

Process 4–1: completes the cycle by a constant-volume process in which heat is rejected from the air while the piston is at bottom dead center.

Process 1–0: the mass of air is released to the atmosphere in a constant pressure process.

Engine Kinematics Parameters:

1. During the compression and expansions processes the volume change with respect to the crank angle is given by the formula below:

where, v = volume traced at particular crank angle

V_c = Clearance Volume

r_c = Compression ratio

R = Length of connecting rod

θ = Angle traced by crank during rotation

This equation is used to get a smooth curve during compression and expansion process in the cycle.

2. Thermal Efficiency of the cycle is given by the formula below:

r = compression ratio

 = Gamma which is the ratio of specific heats

The problem statement is defined by plotting the P-V diagram for different values of pressure and volume values and to find the thermal efficiency of the cycle. The volume array is defined for different crank angles and corresponding pressure is calculated and plotted using the above formulas.

SOLUTION PROCEDURE:

#Otto Cycle Simulator
#By Ashwen Venkatesh

import math
from matplotlib import pyplot as plt
import numpy as np


def engine_kinmeatics(bore, stroke, con_rod, cr, start_crank, end_crank):
	#Geometric Parameters
	a= stroke/2
	R= con_rod/a

	#Volume parameters
	v_s=(math.pi/4)*pow(bore,2)*stroke
	v_c=v_s/(cr-1)

	sc=math.radians(start_crank)
	ec=math.radians(end_crank)

	#Defining the number of iterations for different crank angles
	num_values=50
	dtheta=(ec-sc)/(num_values-1)
	V=[]

	for i in range(0,num_values):
		theta=sc+i*dtheta
		term1=0.5*(cr-1)
		term2=R+1-math.cos(theta)
		term3=pow(pow(R,2)-pow(math.sin(theta),2),0.5)
		V.append((1 + term1*(term2-term3))*v_c)

	return(V)

#Function to compare the colour of the line
def get_color():
    r = g = b = 0
    return r, g, b

red, green, blue = get_color()


def add_arrow(line, position=None, direction='right', size=15, color=None):
    """
    add an arrow to a line.

    line:       Line2D object
    position:   x-position of the arrow. If None, mean of xdata is taken
    direction:  'left' or 'right'

    size:       size of the arrow in fontsize points
    color:      if None, line color is taken.
    """
    if color is None:
        color = line.get_color()

    xdata = line.get_xdata()
    ydata = line.get_ydata()

    if position is None:
        position = xdata.mean()
    # find closest index to create an arrow
    start_ind = np.argmin(np.absolute(xdata - position))
    if direction == 'right':
        end_ind = start_ind + 1
    else:
        end_ind = start_ind - 1

    #To create an arrow over the line 
    line.axes.annotate('',
        xytext=(xdata[start_ind], ydata[start_ind]),
        xy=(xdata[end_ind], ydata[end_ind]),
        arrowprops=dict(arrowstyle="->", color=color),
        size=size
    )
#Inputs
p1=101325 
t1=500
gamma = 1.4
t3=2300

#Geometric Parameters
bore=0.1
stroke=0.1
con_rod=0.15
cr=12

#Volume calculation
v_s=(math.pi/4)*pow(bore,2)*stroke
v_c=v_s/(cr-1)
v1=v_c+v_s


#State point 2
v2=v_c

#p1v1^gamma=p2v2^gamma
p2=p1*pow(v1/v2,gamma)

#p2v2/t2=p1v1/t1
t2=p2*v2*t1/(p1*v1)

#Computing volume for different crank angles during compression
V_Compression = engine_kinmeatics(bore, stroke, con_rod, cr, 180, 0)

constant = p1*pow(v1,gamma)

#Computing pressure for different volumes during compression
P_Compression=[]

for v in V_Compression:
	P_Compression.append(constant/pow(v,gamma))

#State point 3
v3=v2
p3=(t3/t2)*p2

#Computing volume for different crank angles during expansion
V_Expansion = engine_kinmeatics(bore, stroke, con_rod, cr, 0, 180)

constant = p3*pow(v3,gamma)

#Computing volume for different crank angles during expansions
P_Expansion=[]

for v in V_Expansion:
	P_Expansion.append(constant/pow(v,gamma))

#State point 4
v4=v1
p4=pow(v3/v4,gamma)*p3
t4=pow(v3/v4,gamma-1)*t3

#Thermal Efficiency of the cycle
Efficiency= (1-(1/pow(cr,gamma-1)))*100
print(Efficiency)

#Assigning labels for the axes
plt.xlabel("Volume (in m^3)")
plt.ylabel("Pressure (in Pa)")

#Plotting various values of pressure and volumes
HeatAddition_Plot = plt.plot([v2,v3],[p2,p3])[0]
Compression_Plot = plt.plot(V_Compression,P_Compression)[0]
Expansion_Plot = plt.plot(V_Expansion,P_Expansion)[0]
HeatRejection_Plot = plt.plot([v4,v1],[p4,p1])[0]

#Adding arrows for the curves
add_arrow(Compression_Plot)
add_arrow(HeatAddition_Plot)
add_arrow(Expansion_Plot)
add_arrow(HeatRejection_Plot)

#Assigning title for the plot
plt.title('Otto Cycle', size=18)
plt.show()

Note: The arrows are created using the function add_arrow in the code. Also, plt.plot is used for assigning title for the plot

RESULTS AND CONCLUSION:

Therefore, a python code is created for the Otto cycle by defining the parameters.

Efficiency of the cycle is found to be 62.98928%

The P-V diagram for the cycle after running the code is shown below: