OTTO CYCLE & PLOTS USING MATLAB

Otto Cycle- An otto cycle is an idealized thermodynamic cycle that describes the functioning of a typical S.I. engine. It incorporates mainly 4 processes. 2 isentropic process & 2 constant vol. process.

The shaded region in the above fig. gives the work-output.

Thermodynamic relations for various process-

Constant vol. (Isochoric) process - V = constant, dV = 0

constant press. (Isobaric)process - P = constant, dP = 0

constant temp. (Isothermal) process - T = constant, dT = 0, PV = constant

constant entropy (Isentropic)/ reversible adiabatic process - S = 0, dS = 0, PV^gamma = constant

Code explaination -

1. Define the objective of the program by comment.

2. Define the inputs-

• gamma (a constant) = 1.4
• Generally, the temp. or press. after combustion is specified, let we know temp. T3 = 2300 K

3. Engine geometric parameters -

• bore (diametere of engine cylinder).
• stroke/ stroke length (the dictance travelled by piston from TDC to BDC or within the engine cylinder)
• cr = connecting rod (Part that connects piston to crankshaft, It converts reciprocation motion of piston to rotary motion of crankshaft).
• v_s = swept volume (Vol. within the the engine cylinder from TDC to BDC).
• v_c = clearance volume (Vol. of engine cylinder when piston is at TDC).
• v_s = (pi/4)*(bore)^2*stroke
• v_c = v_s/(cr-1)

From the above fig. it is clear that for one complete cycle, there are 4 state points (1,2,3&4),each having a set of variables. And to plot a P-V curve, we must have the values of all these variables for all the 4 state points.

i. State variables at state point 1 -

• Let press. at state point 1 be atmospheric press.(1atm. or 101325 Pa.), p1 = 101325
• And temp. be 500 K, T1 = 500
• volume, V1 = v_s +v_c

ii. State variables at state point 2 -

• process 1-2 is an isentropic process, therefore

p1V1^gamma = p2V2^gamma

• p2 = p1 (V1/V2)^gamma, sinnce V1/V2 = cr, therefore,

p2 = p1*(cr)^gamma

• From ideal gas equation - p1V1/T1 = p2V2/T2, therefore

T2 = p2V2T1/(p1V1)

• V2 = v_c
• Here, we call our function "engine_kinematics "to get volume trace for compression stroke. Here 'theta' varies from 180-0 as piston moves from BDC to TDC

wkt pV^gamma = constant, let p1V1^gamma = constant_c

Values for vol. trace get stored in V_compression.

Therefore, p_compression = constant_c/V_compression

iii. State variables at state point 3 -

• process 2-3 is a constant vol. process, therefore

V3 = V2 = v_c

• From ideal gas equation - p2V2/T2 = p3V3/T3, therefore

p3 = p3T3/T2

• Here again, we call our function "engine_kinematics" to get volume trace for expansion stroke. Here 'theta' varies from 0-180 as piston moves from TDC to BDC

wkt pV^gamma = constant, let p3V3^gamma = constant_c

Values for vol. trace get stored in V_expansion.

Therefore, p_expansion = constant_c/V_expansion

iv. State variables at state point 4 -

• process 3-4 is an isentropic process, therefore

p3V3^gamma = p4V4^gamma

• p4 = p3(V3/V4)^gamma
• process 4-1 is a constant vol. process, therefore

V4 = V1

Now, our process generally depends on how our vol. changes from state point1 to state point2.

To know it(volume trace) , we have to apply engine kinematics.

4. Engine Kinematics :-

The rectilinear motion of piston is converted into rotary motion of crank shaft by crank- slider mechanism. Engine kinematics deals with the motion of the piston which prodces the velocity& acceleration.

Piston kinematics equation -

v/v_c = 1 + (1/2)(cr-1)[R+1-cos(theta) - (R^2-sin(theta)^2)^(1/2)

where,

v = cylinder vol. at any crank angle.

v_c = clearance vol.

cr = compression ratio

R = ratio of connecting rod length to crank pin radius(a)

a = stroke/2

The above equation shows how piston moves from TDC to BDC & vice-versa and how vol. of combustion chamber changes w.r.t. crank angle (volume trace).

Assumption :- theta =0 when piston is at TDC & 180 at BDC.

theta = crank angle

We have to use the above equation multiple tmes in our program. Therefore instead of writing the whole stuff again & again, we can create a functon of engine_kinematics and then call it into our main program when required.

engine_kinematics function code explaination-

1. initiate the program by typing function followed by entity in square bracket whose value will be obtained by this function.
2. In RHS , first write the function name followed by parameters on which it depends in paranthesis.
3. Define the variables, a, R, v_s & v_c.
4. theta =changes in between 0-180, therefore we use linspace comand, & get 100 values in between them.
5. starting value = starting crank angle (start_crank), final value = end crank angle (ebdcrank).
6. Now we use our piston kinematics equation.Though the equation is too long, therefore we break it into 3 terms as follows-

term1 = 0.5*(cr-1)
term2 = R+1- cosd(theta)
term3 = (R^2-(sind(theta)).^2).^0.5, d(theta) as our values are in degrees.

g. And finally merge them acc. to eqn.-

V = (1+ term1*(term2- term3)).*v_c

h. To multiply arrays we have to use element-wise operators.

i. To finish the function , we write end.

engine_kinematics function code :-

 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);
 %  v/v_c = 1+ 0.5(cr-1)(R+1-cosd(theta)-((R^2 - sind(theta)^2))^0.5)
 
 term1 = 0.5*(cr-1);
 term2 = R+1- cosd(theta);
 term3 = (R^2-(sind(theta)).^2).^0.5;
 
 V = (1+ term1*(term2- term3)).*v_c;
 
end

As we defined 100 values for crank angle,& for each value we are going to have volume. crank angle mainly changes during compression & expansin stroke.Therefore, we used this function twice. one for compression stroke & another for expansion stroke.

MainProgram code :-

 % P-V diagram for ideal otto cycle
clear all
close all
clc

% Inputs
gamma = 1.4;
T3 = 2300;
%state variables at state point 1 
p1 = 101325;
T1 = 500
 %Engine Geometric Parameters
 bore = 0.1;
 stroke = 0.1;
 con_rod = 0.15;
 cr = 10;
 % swept vol. & clearance vol. calculation
 v_swept = (pi/4)*bore^2*stroke
 v_clearance = (v_swept)/(cr-1)
 
 V1 = v_swept + v_clearance
 V2 = v_clearance
 
 % state variables at state point 2
 % p2 = p1*(v1/v2)^gamma = p1*cr^gamma
 p2 = p1*cr^gamma;
 %p1v1/t1 = p2v2/t2 | T2 = (p2*v2*T1)/(p1*v1)
 T2 = (p2*V2*T1)/(p1*V1);
 %for isentropic compression pv^gamma = constant
 constant_c = p1*V1^gamma;
 
 V_compression = engine_kinematics (bore,stroke,con_rod,cr,180,0);

 P_compression = constant_c./V_compression.^gamma;
 
 %state variables at state point 3
V3 = V2

%p3v3/t3 = p2v2/t2 | p3 = (p2*T3)/T2
 p3 = (p2*T3)/T2

 constant_c = p3*V3^gamma;
 
 V_expansion = engine_kinematics (bore,stroke,con_rod,cr,0,180);

 P_expansion = constant_c./V_expansion.^gamma;
 
 
 
 %state variables at state point 4
 V4 = V1
 %p3*v3^gamma = p4*v4^gamma | p4 =p3*(v3/v4)^gamma
 p4 =p3*(V3/V4)^gamma
 %p3V3/T3 = p4V4/T4 | T4 = p4V4T3/(p3V3)
 T4 = p4*V4*T3/(p3*V3)

 %Thermal efficiency
 nth1 = (1-(1/(cr^(gamma-1))))*100
nth2 = (1-((T4-T1)/(T3-T2)))*100
 % plotting
 figure (1)
 hold on 

 axis ([-1e-04,10e-04,-1e+06 7e+06])
 xlabel ('Volume')
 ylabel ('pressure')

 plot (V1, p1,'*','color', 'r')
 plot (V_compression, P_compression)
 plot (V2, p2,'*','color', 'r')
 plot (V3, p3,'*','color', 'r')
 plot (V_expansion, P_expansion)
 plot (V4, p4,'*','color', 'r')
 plot ([V2 V3], [p2 p3],'color','b')
 plot ([V4 V1], [p4 p1],'color','b')
 

Thermal efficiency-

nth1 = (1-(1/cr^(gamma-1)))*100

nth2 = (1-((T4-T1)/(T3-T2)))*100

Both formule are for thermal efficiency of otto cycle. Just to cross-check , I used both of them.

nth1 = nth2

Plotting -

1. Define figure
2. hold on , as multiple curves are going to be plotted in a single figure.
3. Adjust axis, so that our diagram will be clearly visible.
4. Enter x-axis & y-axis parameter names.
5. First we plot the state points.To plot he points simply define point coordinates with symbol & color under plot command.
6. Then we connect them.
7. For connecting 1-2, we use (V_compression & p_compression) and for 3-4 we use (V_expansion & p_expansion), as it gives the actual curve & not just the straight line.
8. Point 2-3 and 4-1 get connected as a straight line as processes 2-3 & 4-1 are costant volume process.

error :-

error1 - volume at state point 4 was defined by V4, and I wrote v4 while plotting.It was corrected by chainging v4 to V4.

error2 - figure (Ideal otto cycle), but MATLAB only consider numeric values in parenthesis after figure command, so, I changed it to figure (1).

Output :-

P-V curve - We got req. P-V diagram with 2 isentropic & 2 constant vol. curves.

Thermal efficiency -

1. nth1 =nth2 = 60.1893%
2. In general otto cycle engines have cr in between 8 to 11, and efficiency in between 56% to 61%.
3. For our otto cycle with cr=10, efficiency is 60.1893% that is completely in range.
4. Hence our assumptions & calculations are correct.