Week 2 Air standard Cycle

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.

Once you have done this, create a project in projects.skill-lync.com.

Use the "INSERT CODE SAMPLE" option to paste your code and attach screenshots of the output.

""" 
otto cycle simulator
"""
import math
import matplotlib.pyplot as plt
import numpy as np

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

	#geometric parameters
	#crank pin radius
	a=stroke/2
	R=con_rod/a

	#volume comutation
	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_value=50

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

	V=[]

	for i in range(0,num_value):
		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

bore=0.1
stroke=0.1
con_rod=0.15
cr=13


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

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

#p2v2/t2=p1v1/t1 |  rhs=p1v1/t1 | p2v2/t2=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 | p3v3/t3=rhs |

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=p3v3/t3 | rhs=p3v3/t3  | p4v4/t4=rhs |

rhs=p3*v3/t3 

t4=p4*v4/rhs
plt.figure(1)
plt.plot(V_compression,P_compression,label='Adiabatic compression')
plt.plot([v2,v3],[p2,p3],label='Constant volume heat addition')
plt.plot(V_expansion,P_expansion,label='Adiabatic expansion')
plt.plot([v4,v1],[p4,p1],label='Constant volume heat rejection')
plt.xlabel('volume')
plt.ylabel('pressure')
plt.legend()
plt.savefig('otto cycle pv diagram')
plt.show()


#thermal efficiency

eff=1-(1/pow(cr,(gamma-1)))
eff=str(eff*100)
print('thermal efficiency =' + eff + '%')