OTTO CYCLE SIMULATOR USING PYTHON

INTRODUCTION

An Otto cycle is an idealized thermodynamic cycle that describes the functioning of a typical spark ignition piston engine. The Otto cycle is a description of what happens to a mass of gas as it is subjected to changes of pressure, temperature, volume, the addition of heat, and removal of heat. An Otto cycle is constructed by a pair of quasi-parallel isentropic processes and by a pair of parallel isochoric processes. An Otto cycle is represented by a P-V diagram.  

EXPLANATION OF PROGRAM FOR PLOTTING P-V DIAGRAM:

Math library and matplot library are imported to the program using import command. matplotlib.pyplt is imported as plt to plot the graph.

`(p_1*v_1)/T_1 = (p_2*v_2)/T_2`

The program is used to calculate all the state variables(P, V, T) of four state points and make a plot of these values on the P-V diagram. So, to calculate all the state variables at each state points we already know the state variables at state point 1 as we take these values as input values to calculate the remaining state variables of three state points. 

For state variables at state point 1, the values of P1 and T1 are taken arbitrarily and V1 is calculated from Engine geometric parameters(Bore, Stroke, Connecting rod length, Compression ratio) by calculating volume swept by piston and clearance volume by formulae

volume swept = (pi/4)*bore*stroke;

clearance volume = volume swept/(compression ratio - 1);

V1 = volume swept + clearance volume

For state variables at state point 2, the value of V2 is equal to clearance volume and P2 is calculated from `P_1*V_1^gamma = P_2*V_2^gamma`| `p2 = p1*r_c^gamma`, T2 is calculated by ideal gas equation `(p_1*v_1)/T_1 = (p_2*v_2)/T_2` | `T_2 = (p_2*v_2*T_1)/(p_1*v_1)` . I took values of gamma and T3 as input values.

For state variables at state point 3, the value of V3 is equal to V2, P3 is calculated using the ideal-gas equation.

For state variables  at state point 4, the value of V4 is equal to V1, P4 is calculated.

We Know that Volume inside of engine changes with rotation of the crank. So, we need to calculate the volume and pressure for each angle of the crank.

So, I have created a function to calculate Volume inside of engine with rotation of crank using formula `V/V_c=1+(1/2*(r_c-1))*{(R+1-cos(theta))-[(R^2-sin(theta)^2)]^(1/2)}`

So, basically, we just need to give the right inputs to calculate volume using the above formula.

def is used to define the function and engine_kinematics is the name of the function defined which takes the inputs values of bore, stroke, con_rod, cr, start_crank, end_crank. In this function we need to calculate the volume for different angles of the crank depending on the starting crank angle and ending crank angle. This function returns the values of volume in an array. 

This function is used or called to calculate volume during compression and volume during expansion, using thermodynamic formulas pressure during compression is calculated and pressure during expansion is calculated.

\"\"\"
otto cycle simulator
\"\"\"
import math
import matplotlib.pyplot as plt


def engine_kinematics(bore, stroke, con_rod, cr, start_crank, end_crank):
	\"\"\"
	   Engine kinematics
	\"\"\"

	# Geometric parameters

	a = stroke/2
	R= con_rod/a

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

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

	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(R,2) - pow(math.sin(theta),2)
		term3 = pow(term3,0.5)

		V.append((1 + term1*(term2 - term3))*V_c)
	return V

#inputs
p1 = 101325
t1 = 500
gamma = 1.4
t3 = 2300

#geometric parameters
bore = 0.1
stroke = 0.1
con_rod = 0.15
cr = 15

#volume computations
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

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

#p2v2/t2 = p11v1/t1 | rhs = p1v1/t1 |p2v2/t2 = rhs | t2 = p2v2/rhs
rhs = p1*v1/t1
t2 = p2*v2/rhs
v_compression = engine_kinematics(bore, stroke, con_rod, cr, 180, 0)

constant = p1*pow(v1,gamma)

p_compression = []

for v in v_compression:
	p_compression.append(constant/pow(v, gamma))

#state point
v3 = v2

#p3v3/t3 = p2v2/t2 | rhs = p2v2/t2 | p3 = rhs*t3/v3

rhs = p2*v2/t2
p3 = rhs*t3/v3

v_expansion = engine_kinematics(bore, stroke, con_rod, cr, 0, 180)

constant = p3*pow(v3,gamma)

p_expansion = []

for v in v_expansion:
	p_expansion.append(constant/pow(v, gamma))


#state point 4 
v4 = v1

#p4v4^gamma = p3v3^gamma
p4 = p3*pow(v3,gamma)/pow(v4,gamma)

#p4v4/t4 = rhs
t4 = p4*v4/rhs

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.xlabel(\'Volume\')
plt.ylabel(\'Pressure\')
plt.show()

PLOTTING THE VALUES

Using plt.plot command [v2,p2],[v3,p3] and [v4,p4],[v1,p1] as the joining of these points implies change in pressure at constant volume.

Finally to make plot we need to plot V_compressed and P_compressed, And also we need to plot V_expansion and P_expansion.

plt.show is used to display the plot.

RESULT