Week 6 - Data analysis
"""
Task to write a script basic performance of Engine with given engine_data
Basic performance calculation
Calculate the area under the P-V diagram. P is in column 2 and V is in column 8.
Calculate the power output of this engine. Assume that RPM is 1500
This engine consumed 20 micro grams of fuel (per 1 stroke cycle). Calculate its specific fuel consumption.
"""
import math
import matplotlib.pyplot as plt
import numpy as np
from scipy.integrate import simps
#Compatibility check
fname = input("Enter the filename:..")
if fname == 'engine_data.out':
print('file is valid')
else:
print('file is invalid')
exit()
def compatibility_check(a):
if 'CONVERGE' not in a:
print('Please provide converge output file to process')
exit()
else:
print('Compatibility check is successful')
pass
#try:
open('engine_data.out')
line_count = 1
# Creating an empty array to store the values of Volume and pressure respectively
Volume =[]
pressure =[]
v = 8
p = 2
# Creating null array and storing all the values of data file in to it
data =[]
for line in open('engine_data.out'):
if '#' not in line:
data.append(line.split())
if line_count >= 6:
Volume.append(float(line.split()[v-1]))
pressure.append(float(line.split()[p-1]))
line_count = line_count + 1
#plotting p-v diagram
plt.plot(Volume, pressure, linewidth = 2, color = 'b')
plt.xlabel('Volume')
plt.ylabel('pressure')
plt.title("p-v diagram of engine")
plt.grid()
plt.show()
area = simps(pressure,Volume)
print('Area under the curve using Simpsons Rule :' +str(area) + ' units')
print('Work done by the engine : ' + str(area) + ' MJ')
area = area * 1e6
print('Work done by engine : ' + str(area) + ' J')
# Inputs
# Speed of an Engine
rpm = 1500
rps = 1500/60
stroke_cycle = rps/2
# Fuel consumption in micrograms per 1 stroke cycle
fuel_consumption = 20
# Engine specifications
# Number of cylinder
n = 1
# 4 stroke engine
# Calculation of power output and Specific fuel consumption
power_output = n*area*stroke_cycle/(1000)
print('power output from engine : ' + str(power_output) + 'KW')
SFC = (fuel_consumption*pow(10,-6)*stroke_cycle*3600)/(power_output*1000)
print('Specific Fuel consumption, SFC = ' + str(SFC) + ' (kg/hr)/KW')
print('t')
# Calculation of Area under PV diagram
Area = np.trapz(pressure,Volume)
print('Area under curve using trapezoidal rule : ' + str(Area) + 'units')
print('Work done by the engine : ' + str(Area) + ' MJ' )
Area = Area * 1e6
print('Work done by engine : ' + str(Area) + ' J')
power_output = n*Area*stroke_cycle/(1000)
print('power output from the engine : ' + str(power_output) + 'KW')
SFC = (fuel_consumption*pow(10,-6)*stroke_cycle*3600)/(power_output*1000)
print('Specific Fuel consumption, SFC = ' + str(SFC) + ' (kg/hr)/KW')
print('t')
# Transposing columns in to row matrix and converting them in to float type from string type
transpose = np.transpose(data)
float_f = transpose.astype(np.float)
# Array for header and their units
parameters = ['Crank','Pressure','Maximum Pressure','Minimum Pressure','Mean Temperature','Maximum Temperature','Minimum Temperature','Volume','Mass','Density','Integrated HR','HR rate','Cp','Cv','gamma','Kinematic Viscosity','Dynamic Viscosity']
units = ['(DEG)','(Mpa)','(Mpa)','(Mpa)','(K)','(K)','(K)','pow(m,3)','(kg)','(kg/(pow(m,3)))','(J)','(J)','(J/time)','(J/kgK)','(J/kgK)','','pow(m,2)/s','(Ns/pow(m,2)']
# Providing user choice of picking 2 parameters for plotting
print('Choose from following options and press enter key')
print('1 Crank \n 2 Pressure \n 3 Maximum Pressure \n 4 Minimum Pressure \n 5 Mean Temperature \n 6 Maximum Temperature \n 7 Minimum Temperature \n 8 Volume \n 9 Mass \n 10 Density \n 11 Integrated HR \n 12 HR rate \n 13 Cp \n 14 Cv \n 15 gamma \n 16 Kinematic Viscosity \n 17 Dynamic Viscosity')
X = int(input('Enter your choice on X-axis \n'))
Y = int(input('Enter your choice on Y-axis \n'))
# Plotting based on user choice
plt.figure(1)
plt.plot(float_f[X-1],float_f[Y-1],linewidth=3)
plt.title(' plot for ' + parameters[Y-1]+ ' versus ' +parameters[X-1],color='k')
plt.xlabel(parameters[X-1] + ' ' +units[X-1])
plt.ylabel(parameters[Y-1] + ' ' +units[Y-1])
plt.grid()
plt.show()