Air standard Cycle

#Otto cycle simulator

import math
import matplotlib.pyplot as plt

def piston_kinematics(bore,stroke,con_rod,cr,start_crank,end_crank):

	#Geometry
	a = stroke/2
	R = con_rod/a

	#End volumes
	V_s = math.pi*(1/4)*pow(bore,2)*stroke
	V_c = V_s/(cr-1)


	sc = math.radians(start_crank) #Starting crank angle
	ec = math.radians(end_crank) #Ending crank angle

	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


#State inputs
p1 = 101325
t1 = 500
t3 = 2300

#Constants
gamma = 1.4

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


#Swept and clearance volume
v_s = (math.pi/4)*pow(bore,2)*stroke
v_c = v_s/(cr-1)

#State point 1
v1 = v_c+v_s

#state point 2
v2 = v_c

#pv^gamma = constant
p2 = p1*pow(v1,gamma)/pow(v2,gamma)

#pv/t = constant
const1 = p1*v1/t1
t2 = p2*v2/const1
V_compression = piston_kinematics(bore,stroke,con_rod,cr,180,0)

const2 = p1*pow(v1,gamma)

P_compression = []

for v in V_compression:
	P_compression.append(const2/pow(v,gamma))


#state point 3
v3 = v2

#pv/t = constant
const1 = p2*v2/t2
p3 = const1*t3/v3

V_expansion = piston_kinematics(bore,stroke,con_rod,cr,0,180)

const2 = p3*pow(v3,gamma)

P_expansion = []

for v in V_expansion:
	P_expansion.append(const2/pow(v,gamma))

#state point 4
v4 = v1

#pv^gamma = constant
p4 = p3*pow(v3,gamma)/pow(v4,gamma)

#pv/t = constant
t4 = p4*v4/const1

#Thermal efficiency of otto cycle
n_therm = 1-(1/pow(cr,gamma-1))
print(n_therm)

plt.plot([v2,v3],[p2,p3])
plt.plot(V_compression,P_compression)
plt.plot(V_expansion,P_expansion)
plt.plot([v4,v1],[p4,p1])
plt.show()

First we create a function (‘piston_kinematics’) which creates an array of chamber volume as a function of crank angle.

The inputs of the function are bore, stroke, connecting rod length (‘con_rod’), compression ratio (‘cr’), crank start angle, crank end angle.

Calculate engine parameters crank pin radius and ratio of connecting rod and crank pin radius.

Find the swept volume and clearance volume.

Calculate the volume at the specific crank angle as per the equation:

`V=(1+1/2 (r_c-1)[R+1-cosθ-(R^2-sin^2⁡θ )^(1/2) ])*V_\"clearance\"`

End of function.

The above created function will be used in the following program:

This program will plot the PV diagram and calculate the thermal efficiency of the engine.

Input the state variables p1, t1 and t3.

Input the engine geometric parameters such as bore, stroke, connecting rod length, compression ratio.

Input the value of gamma which is a constant.

Calculate the swept volume and clearance volume.

Then state variables v1 and v2 are calculated.

State 1 to 2 is isentropic compression. Using the ideal gas equation t2 is calculated. And using the isentropic process equation (pv^gamma = constant), p2 is calculated.

Using the ‘piston_kinematics’ function, array of chamber volume during compression is calculated.

Corresponding pressure during the compression is calculated using the isentropic process equation.

State 2 to 3 is constant volume heat addition.

p3 is calculated using the ideal gas equation.

State 3 to 4 is isentropic expansion

v4 is equal to v1, because process 4 to 1 is constant volume process.

Using the isentropic process equation (pv^gamma = constant), p4 is calculated.

Using the ‘piston_kinematics’ function, array of chamber volume during expansion is calculated.

Corresponding pressure during the expansion is calculated using the isentropic process equation.

Thermal efficiency is calculated using the equation:

`η=1-1/(cr^(γ-1) )`

Individual plots of the 4 processes are plotted in one figure.