Visualizing the pV Diagram of otto cycle and finding its thermal efficiency

Aim: To create the PV diagram of an otto cycle using matlab program and output its thermal efficiency.

Theory

Otto cycle is a set of processes the spark ignition engines use to produce power.It involves 4 strokes.At the first stroke the piston moves from Top Dead Cente(TDC) to Bottom Dead cente(BDC) and sucked the air into the cylinder.This stroke is called Suction Stroke.Then the piston will move up from BDC to TDC.This stroke is called compression stroke.At the end of compression stroke,when the piston reaches the TDC, spark plug will ignite a spark into the air fuel mixture and the combustion of fuel occurs .Then the piston expands from TDC to BDC.This is called a expansion stroke.Then the piston again moves up from BDC to TDC to scavenge the exhaust gauses through the exhaust valve.This stroke is called an exhaust stroke.

Figure above represents the PV diagram of petrol engine cycle.This is a constant volume cycle in which the heat is added and heat is removed at constant volume.Compression and expansion stroke are reversible adiabatic /Isentropic processes.

• 5-1 is the suction stroke in which gasoline and air is drawn in.
• 1-2 is the compression stroke.
• 2-3 is the head addition process in which the spark plug ignite the air-fuel mixture.
• 3-4 is the expansion stroke or power stroke.
• 4-1 is the heat rejection process .
• 1-5 is the exhaust stroke.

V2 is the clearance volume.v1 is the cylinder volume and v1=v4(constant volume process).Swept volume is equal to v1-v2.For isentropic process pv^_ϒ = constant.And the gamma, which is the ratio of specific heats is taken as 1.4.At the end of the program we will be finding out the Thermal efficiency of the system.Compression ratio(cr) of the system is equal to v1/v2.

Thermal efficiency of the otto cycle =_{1-1/(cr)^ϒ-1.}

We can find the state variables at each point in ottocycle using ideal gas equation and isentropic equations.But with matlab we can only mark these points on plot intially.After finding out how the volume changes with engine kinematics,we will be able to join those points.So intially we are writing the program to mark the relevant points on otto cycle.

Intial Program

clear all 
close all
clc

%inputs
gamma = 1.4
t3 = 2000
 %state variables 
 p1 = 101325;
 t1 = 500;
 
 %engine parameters 
 
 bore = .1;
 stroke = .1;
 con_rod = .15;
 cr = 12;
 
 %calculating swept volume and clearance volume
 
 v_swept = (pi/4)*bore^2*stroke;
 v_clearance = v_swept/(cr-1);
 v1 = v_clearance + v_swept;
 v2 = v_clearance; 
 
 %state variables at the point 2
 %p1*v1^gamma = p2*v2^gamma;
 
 p2 = p1*cr^gamma;
 
 %p1v1/t1 = p2v2/t2 | t2 = p2v2t2/(p1v1)
 
t2 = p2*v2*t1 /(p1*v1);

%state variables at point 3
v3 = v2;
p3 = p2*t3/t2;

%state variables at point 4
v4 = v1;
p4 = p3*(v3/v4)^gamma;
t4 = p4*v4*t3/(p3*v3);
%THERMAL EFFICIENCY
Thermal_efficiency =  1-(1/cr)^(gamma-1)

figure(1)
hold on
plot(v1,p1,'*','color','r')
plot(v2,p2,'*','color','r')
plot(v3,p3,'*','color','r')
plot(v4,p4,'*','color','r')

• Here the program is starting with clear all, close all, clc to close the plots,clear the workspace and empty the command window.
• When entering into the program we need to type the inputs required which are the engine parameters like bore,stroke,connecting rod length,compression ratio and the state variables at point 1 and the temperature after expansion(t3).
• We have the equations to calculate the swept volume and clearance volume.From Vc and Vs we can find V1.And by using the adiabatic equation, we can find out the P2.With the ideal gas equation we can find out t2.
• Similarily the other state variables can be found out.And the equations needed to find the state variables are to be input in the program.
• Then we need to type the equation of thermal efficiency of ottocycle.And since we need the program to output the thermal efficiency ,we are not gonna put a double colon at end of the code line.
• Then we can use plot command to mark the 4 points.We can use plot command 4 times with the appropriate coordinates to get the points plotted on graph.Since plot command only displays latest plot,Use hold on to displays the 4 points.We can use color to the marker in plots using color inside the plot command.

To join the markers straight line is not right.We need to find out how volume changes with angle of crank.Then only we can find out how p changes with v .Finding that only helps to join the marks of isentriopic processes.So we make use of piston kinematics here.

Basically a equation connecting v and angle of crank is use here.

(V/Vc)=1+(1/2)(cr-1){R+1-cos(theta)-(R^2-sin(theta)^.5}.

Since the volume is the function of theta ,we can make use of function command in matlab.So we can call it when we needed it.

Making the Template for calling the function

Program

 function[V] = engine_kinematics(bore,stroke,con_rod,cr,start_crank,end_crank)
 
 a = stroke/2 ;
  R = con_rod/a ;
 
 v_s = (pi/4)*bore^2*stroke ;
 v_c = v_s/(cr-1);
 theta = linspace(start_crank,end_crank,100);
 
 term1 = .5*(cr-1);
 term2 = R+1-cosd(theta);
 term3 = (R^2-sind(theta).^2).^.5 ;
 
 V = (1+term1*(term2-term3)).*v_c;
 
 end
 

We are creating a function engine_kinematics whose inputs are bore,stroke,connecting rod length,compression ratio,start crank angle,end crank angle and the output is V.Since the bore,stroke etc.. are inputs, we need not type in the program instead we can give the inputs of function when we call it.Then type in the relevant inputs,since we are using equation of v there is parameters like R,a needs to be defined in the program as inputs.So define it at first.Then Vs and Vc are inputted in the program.Then we need to inputs the theta.We can use linspace command to have a range of values of theta from start to end of crank angles.So that we can call the function for any start and end values of crank angles.Then we basically split the long equation of V into three terms and we defined each terms as term1,term2,term3.And then we define the equation in program with these three terms.End should be provided at the end of function.Thus a function template is made and we can call it whenever we needed it.

After entering the inputs in function for 0 to 720 start to end crank angle ,this is the variation of V with theta.

FINAL PROGRAM

clear all 
close all
clc

%inputs
gamma = 1.4
t3 = 2000
 %state variables 
 p1 = 101325;
 t1 = 500;
 
 %engine parameters 
 
 bore = .1;
 stroke = .1;
 con_rod = .15;
 cr = 12;
 
 %calculating swept volume and clearance volume
 
 v_swept = (pi/4)*bore^2*stroke;
 v_clearance = v_swept/(cr-1);
 v1 = v_clearance + v_swept;
 v2 = v_clearance; 
 
 %state variables at the point 2
 %p1*v1^gamma = p2*v2^gamma;
 
 p2 = p1*cr^gamma;
 
 %p1v1/t1 = p2v2/t2 | t2 = p2v2t2/(p1v1)
 
 constant_c = p1*v1^gamma;
t2 = p2*v2*t1 /(p1*v1);
v_compression = engine_kinematics(.1,.1,.15,12,180,360);
p_compression = constant_c./v_compression.^gamma;

%state variables at point 3
v3 = v2;
p3 = p2*t3/t2;

constant_c = p3*v3^gamma;
v_expansion = engine_kinematics(.1,.1,.15,12,0,180);
p_expansion = constant_c./v_expansion.^gamma;


%state variables at point 4
v4 = v1;
p4 = p3*(v3/v4)^gamma;
t4 = p4*v4*t3/(p3*v3);
Thermal_efficiency =  1-(1/cr)^(gamma-1)
figure(1)
hold on
plot(v1,p1,'*','color','r')
plot(v_compression,p_compression,'color','b')
plot(v2,p2,'*','color','r')
plot([v2,v3],[p2,p3],'color','r')
plot(v3,p3,'*','color','r')
plot(v_expansion,p_expansion,'color','b')
plot([v4 v1],[p4 p1],'color','r')
plot(v4,p4,'*','color','r')
xlabel('volume')
ylabel('pressure')

CODE EXPLANATION

• Strating with clear all,close all,clc.
• Entering the inputs and engine parameters.
• Entering the relevant equations such as ideal gas equation and isentropic equation to find out the state variables.
• We know in isentropic process pV^gamma is a constant.So when in compression stroke Pcompression *Vcompression^gamma will be a constant.Vcompression can be finding out by calling the function Engine kinematics by entering the start crank angle and end crank angle as 180 and 360 respectively.Then Pcompression can be found out from isentropic equation.
• Similarly as above, Vexpansion and Pexpansion can be found out by using start crank angle and end crank angle as 0 and 180 respectively.Appropriate Equations are typped in.
• Code for finding state variables at point 4 is entered. And equation for Thermal Efficiency is also typped in .
• Plot command is used to connect points.Plot command with Vcompression and Pcompression as coordinates will connect points 1 and 2.And plot command with Vexpansion and Pexpansion as the cordinates will connect points 3 and 4.Color also given to each plot command.
• To connect between points 1 and 2 and 4 and 1,Since they are constant volume processes,we can use plot command and include the appropriate V and P values at two points (which is to be connected) in pairs.
• Labelling the axis as velocity and pressure using xlabel and ylabel.

ERRORS FACED

I missed the input which is compression ratio and then this error is shown.So clearly entering the all the inputs in correct order will solve the problem.

Here since the v_compression is a matrix,constant_c will be a matrix and so when dividing we should use the dot vector.I forget it here.And when put that dot vector ,the problem is solved.

RESULT

Otto cycle is plotted by matlab programming.

Thermal Efficiency of the emgine is calculated as ~ 63%