Independent Research Project

The program written is for simulating an Air standard cycle (Otto cycle). The Otto cycle consists of following four processes.

• 1-2 Isentropic compression
• 2-3 Constant volume heat addition
• 3-4 Isentropic expansion
• 4-1 Constant volume heat rejection

This requires inputs like V1, T1, Compression ratio(cr),T3 from which P1,V2,P2,T2,P3 and V3 can be calculated using the formula (P1*V1/T1)=(P2*V2/T2) and (P*V^gamma)=constant. This also requires parameters like Bore, Stroke, Connecting rod length, Compression ratio.

The computer program for plotting PV diagram and thermal efficiency of the Engine is as follows.

"""
  otto cycle simulator
"""

import math
import matplotlib.pyplot as plt

def eng_kinematics(bore,stroke,con_rod,cr,start_crank,end_crank):
    
    #geometric parameters
	a=stroke/2
	R=con_rod/a

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

	num_values=50

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

	v=[]

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

		v.append((1+term1*(term2-term3))*v_c)
	
	return v


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

# eng geometric parameters
bore=0.1
stroke=0.1
con_rod=.15
cr=12

#calculating the swept vol and clearence vol
v_s=(math.pi/4)*pow(bore,2)*stroke
v_c=v_s/(cr-1);

v1=v_s+v_c

#state point 2
v2=v_c

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

#p2*v2/t2=p1*v1/t1 | rhs=p1*v1/t1 | lhs=p2*v2/t2
t2=t1*p2*v2/(p1*v1)
v_compression=eng_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
p3=p2*t3/t2

v_expansion=eng_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
p4=p3*pow((v3/v4),gamma)
t4=t3/pow(cr,gamma-1);


thermal_efficiency=1-(1/pow(cr,gamma-1))
print('Thermal efficiency is=',thermal_efficiency*100)


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('v')
plt.ylabel('p')
plt.show()

plot is as shown below

Output is

Thermodynamic data collected from simulation of a Port Fuel Injected engine is extracted and analysed from a data file using a Python program.

This data file consists of 17 columns of different thermodynamic parameters for which more than 8000 values have been recorded. The objective is to open the file in the directory, read and extract desired data from the file. The program takes the input from the user and plots the graph for the respective parameter and also calculates the Work done per cycle, Power output for specified rpm and Specific Fuel Consumption.

Below is the program for the engine data analysis.

"""
   Engine fata analysis using python
"""

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

#input file from user
filename = input("nEnter the name of CONVERGE output filename with the extension(.out): ")

#compatibility check
try:
 	data=open(filename,'r') #reading the input file
except FileNotFoundError:
	print("File not recognized. Please provide a valid CONVERGE output file")

#skipping the comments
value=[]
for line in data:
	if '#' not in line:
		value.append(line.split())
		
column=np.transpose(value)
column=column.astype(float)

#parameters from the filename
parameter=['Crank Angle','Pressure','Maximum Pressure','Minimum Pressure','Mean Temperature','Maximum Temperature','Minimum Temperature','Volume','Mass','Density','Heat Rejection','Heat Rejection Rate','C_p','C_v','gamma','Kinematic Viscosity','Dynamic Viscosity']
print(parameter)
print('Each parameter from the file is provided with the code as below')
corresponding_num=['crank_angle=0','pressure=1','max_pres=2','min_pres=3','mean_temp=4','max_temp=5','min_temp=6','volume=7','mass=8','density=9','Heat Rejection=10','Heat Rejection Rate=11','C_p=12','C_v=13','gamma=14','Kin Vis=15','Dyn Vis=16']

#taking input for x and y parameter
print(corresponding_num)
print("please enter the code for independent variable then press enter: n")
x=int(input())
print("please enter the code for dependent variable then press enter: n")
y=int(input())

#plotting
plt.plot(column[x,:],column[y,:],linewidth=3,color='blue')
plt.xlabel(parameter[x],fontsize=10)
plt.ylabel(parameter[y],fontsize=10)
plt.title(parameter[y]+'    VS    '+parameter[x],fontsize=10)
plt.savefig(parameter[y]+'VS'+parameter[x]+'.png')
plt.show()

vol=column[7,:]
press=column[1,:]*pow(10,6)  #converting the pressure into Pascal

#area under the curve
work=[]
for i in range(0,8664):
	work.append((vol[i+1]-vol[i])*(press[i+1]+press[i])*0.5)

work_done=0
for j in range(0,len(work)):
	work_done=work_done+work[j]

print('Area under the p-v diagram i.e, Work Done :',work_done,'J')

#Power calculations
speed=1500/60 
power=work_done*speed
print('Power:',power/1000,'KW')

#specific fuel consumption
m_f=20*pow(10,-6)  
bsfc=(m_f*3600)/power
print('Specific Fuel Consumption:',bsfc,'gm/KWhr')

Output plots some parameters entered by me.

Like this we can give various inputs to obtain the required plots.

Area under the p-v diagram i.e, Work Done : 500.58662746800763 J

Power: 12.51466568670019 KW

Specific Fuel Consumption: 5.753249971073308e-06 gm/KWhr

These are the outputs for the input of Pressure and Volume.