Air standard cycle of an otto cycle simulating in python

Objective:

Based on the concepts you have learnt during the forward kinematics routine, go ahead and write code that can solve an otto cycle and make plots for it.

Here are the requirements:

1. Your code should create a PV diagram.
2. You should output the thermal efficiency of the engine.

Theoretical Explanation: Now, with the help of P-V and T-S diagram, we can easily understand all the process of Otto cycle.

cycle. cycle.

Note: For better explanation, while reading watch the P-V and T-S diagram every time.

The Otto cycle was given by Dr. Nikolaus August Otto. It is a gas power cycle that is used in spark ignition system (i.e. petrol engine) for its working. The entire modern petrol engine works on Otto cycle. It consists of four processes, Two isentropic (reversible adiabatic) processes and two isochoric (constant volume) processes. It has low compression ratio ranges from 7:1 to 10:1. Here we will try to understand this cycle with help of its P-V and T-S diagram.

The four processes of this cycle are as follows:

1. Isentropic (reversible adiabatic) compression
   2. Constant volume (Isochoric) heat addition
   3. Isentropic (reversible adiabatic) Expansion
   4. Constant volume heat rejection.

Before we begin to understand these four processes, lets us first understand about isentropic and isochoric process.

1. Isentropic Process: It is a thermodynamic process in which the entropy of the system remains unchanged (i.e. entropy remains constant). There is no dissipation of heat takes place during isentropic process, so sometime isentropic process called as reversible adiabatic process.
2. Isochoric Process: The process which takes place at constant volume is called isochoric process.

Based on the above two process we can now understand the four processes to begin with.

1. Process 1-2: Isentropic Compression

This process involves the motion of piston from TDC to BDC. The air that is sucked into cylinder during suction stroke undergoes reversible adiabatic (isentropic) compression. Since the air is compressed, the pressure increases from P1 to P2, the volume decreases from V1 to V2, temperature rises from T1 to T2, and entropy remains constant.

1. Process 2-3: Constant Volume Heat Addition

This process is an isochoric process i.e. the heat is added to the air at constant volume. The piston in this process rest for a moment at TDC and during this time heat is added to the air through external source. Due to the heat addition, the pressure increases from P2 to P3, pressure, volume remains constant (i.e. V2=V3), temperature increases from T2 to T3 and entropy increases from S2 to S3.

In this process,

Heat Supplied = mC_v(T₃ – T₂)

where,

m → Mass

C_v → Specific heat at constant volume

The amount of heat added is given by

Q = mC_v (T₃ – T₂)

1. Process 3-4: Isentropic Expansion

In this process, the isentropic (reversible adiabatic) expansion of air takes place. The piston moves from TDC to BDC. Power is obtained in this process which is used to do some work. Since this process involves expansion of air, so the pressure decreases from P3 to P4, volume increases from V3 to V4, temperature falls from T3 to T4 and entropy remains unchanged (i.e. S3=S4).

1. Constant Volume Heat Rejection

In this process, the piston rest for a moment at BDC and rejection of heat takes place at constant volume. The pressure decreases from P4 to P1, Volume remains constant (i.e. V4=V1), temperature falls from T4 to T1.

The amount of heat rejected in this process is given by

Q = mC_v (T₃ – T₄)

When this cycle is used in four stroke petrol engines than two process increases. one is suction of air fuel mixture inside the cylinder which takes place at constant atmospheric pressure and other one is exhaust of gases out of the engine cylinder at constant atmospheric pressure. These two processes are not shown in ideal otto cycle that we have discussed above.

The PV plot is designed when the variables at different stages is known and in order to calculate all state variables of points 1-4, we have different formulas to start with. Below is the calculation of all sate variables putting including engine geomatric paramters.

After attaining formula for swept and clearance volume we can calculate v1 and v2 from the equation.

After the state variables at point 1 are defined we are moving on to the point 2 to determine the pressure and temperature value as at this stage we called it as reversible adiabatic process because volume here is constant.

Calculation of state variables at point 3 can be called as a isochoric process where volume is constant.

Same can be said for calculating of state variables at point 4 which is again can be called as a isochoric process where volume is constant.

Now that we have all state variables we an write a code for the same to compute all the variables to plot PV diagram of Otto cycle in Python.

Heat Rejected = mC_v(T₄ – T₁)

Thermal efficiency (air-standard efficiency) of Otto Cycle:

In order to see how the volume changes inside the IC engines geometric relation we are using a an equation to determine that. It is as below:

Where denominator Vc is clearance volume.

Step by step explanation:

Since there are two program for Otto cycle, we are going to proceed with the main program and function program.

# Air standard cycle
# code for otto cycle simulator

import math
import matplotlib.pyplot as plt

We have started using math fucntion and matplotlib function as we are going to use math function and then we are goin to plot the PV diagram.

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

this here is the engine_kinematics function which we wrote it seperately and then copied here in the main program to get the required path curve for PV plot.

#geomatric paramters

a = stroke/2
R = con_rod/a

a and R are the crank pin radius and constant value

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

Here we are defining the clearance volume and swept volume

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

sc and ec are the start crank and ending crank angle output in radians.

num_values = 50

this num_value represent the how many values and point we needin the curve path.

dtheta = (ec-sc)/(num_values-1)

V=[] is an empty array where we can get the values from V engine kinematics.

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)

here we are putting the terms to evalute the curve function we mentioned earlier and give below the equation.

return V

we are returning the value from V and ending the output.

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

Here again are mentioning the given and assume values both.

# engine geomatric parameters
bore = 0.1
stroke = 0.1
con_rod = 0.15
cr = 12

here we have mentioned the predefined values of bore, stroke, connecting rod and compression ratio as these are the engine geometric parameters.

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

Here we are finding the values of swept and clearance volume to determine v1 and v2.

# state point 2
v2 = v_c

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

Since we know p1 v1 and v2 from above formula, we can determine p2 as v1/v2 = compression ratio.

# p2v2/t2 = p1v1/t1 | rhs = p1v1/t1 | t2 = p2v2/rhs
rhs = p1*v1/t1
t2 = p2*v2/rhs

Now we have p1,v1,v2 and t1 we are determining the value of t2.

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))

This particular expression is completely focused on using for loop where there is an empty array of V_compression and P_compresion, using the ideal gas eqn, to create a curve path along joining the points where pressure and volume are compressed in isentropic process.

# state point 3
v3 = v2

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

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

At stage 3 we have p2 t2 and t3 so we can find p3

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))

Here again this particular expression is completely focused on using for loop where there is an empty array of V_expansion and P_expansion, using the ideal gas eqn, to create a curve path along joining the points where pressure and volume are compressed in isentropic process.

# state point 4

v4 = v1

#p4v4^gamma = p3v3^gamma

p4 = p3*pow(v3,gamma)/pow(v4,gamma)

Same as 2-3 stage in isochoric process, we have p3 v3 v4 to find p4.

#rhs = p4v4/t4

t4 = p4*v4/rhs

# thermal efficiency

n_otto = (1 - (1/pow(cr,gamma-1)))*100
print (n_otto)

Here we need to find what is the thermal efficiency of Otto cycle by applying formula as mentioned in theoretical explanation.

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()

Here we are plotting the points and making the required path for the otto cycle.

Original Coding:

Now here we are going to writing the coding where def engine_kinematics(bore, stroke, con_rod, cr, start_crank, end_crank): is a function program that we are going to calling in the main program.

Main Program :

# Air standard cycle
# code for otto cycle simulator

import math
import matplotlib.pyplot as plt

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

	#geomatric paramters

	a = stroke/2
	R = con_rod/a

	#volumetric 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

# engine geomatric parameters
bore = 0.1
stroke = 0.1
con_rod = 0.15
cr = 12

# volume computation
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 = p1v1/t1 | rhs = p1v1/t1 | 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 3
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)

#rhs = p4v4/t4

t4 = p4*v4/rhs

# thermal efficiency

n_otto = (1 - (1/pow(cr,gamma-1)))*100
print (n_otto)

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()

The engine kinematics has been writtrn separately only to copy later here in the main program named as def engine_kinematics(bore, stroke, con_rod, cr, start_crank, end_crank):

Output : Now the output of the above coding is a follows:

As we can see that our PV diagram has been plotted in the above graph.

And consequently, we are able to obtain our Thermal efficiency of the Otto cycle with the compression ratio 12. As we can see below.

As we can see here in command window our thermal efficiency is 62.98%.