MATLAB Program : Air Standard Cycle (Otto Cycle)

Aim: Write a MATLAB program for plotting PV diagram of air standard cycle

Objective:

• Plotting PV diagram(workflow diagram) for Otto cycle
• Determining air standard efficiency (thermal efficiency)

Theory:

Air Standard Cycle: The cycle in which air is used as a working substance is called an air standard cycle.

Air Standard Efficiency: Ratio of work done byt he system to heat input to the system is known as an "Air standard Efficicency". This efficiency is ideal efficiency.

`eta`= Work Done / Heat Supplied

Air standard assumptions:

1. The working fluid is air, which continuously circulates in a closed loop and alwats behaves as an ideal gas.
2. All the processes that make up the cycle are internally reversible.
3. The combustion process is replaced by a heat-addition process from an external source.
4. The exhaust process is replaced by a heat-rejection process that restores the working fluid to its initial state.
5. The physical constants of air are constant.

The Otto Cycle:

The Otto cycle is the Ideal cycle for spark-ignition reciprocating engines. It is named after Nikolaus A. Otto, who built a successful four stroke engine in 1876 in Germany using cycle proposed by Frenchman Beau de Rochas in 1862. Inmost spark-ignition engines, the piston executes four complete strokes (two mechanical cycles) within the cylinder, and the crankshaft completes two revolutions for each thermodynamic cycle. These engines are called Four-stroke Internal Combustion Engines. A schematic of each stroke as well as PV diagram for an actual four-stroke spark-ignition engine is as shown above.

The working of an Ideal Otto Cycle:

• Process 1-2: Isentropic compression of air.
• Process 2-3: Isochoric heat addition.
• Process 3-4: Isentropic expansion of air.
• Process 4-1: Isochoric heat rejection.

Formulae:

1. Governing Equations:

• Isentropic process: `PV^gamma` = Constant
• Isochoric process: `P/T` = Constant

2. Ideal gas equation:

`(PV)/T` = Constant

3. Compression Ratio

CR = `(Vs + Vc) / (Vc)`

4. Change of Volume inside the engine combustion chamber is given by

`(V)/(Vc) = 1+1/2(CR−1)[R+1−cosθ−(R^2−sin^2θ)^(1/2)]`

5. Thermal Efficiency

`eta = 1 - (1)/(RC)^(gamma-1)`

where,

P = Pressure (`N/m^2`)

V = Volume (`m^3`)

T = Temperature (Kelvin)

`gamma` = ratio of specific heats

Vs = Swept Volume

Vc = Clearance Volume

V = Volume when piston moves at an angle

R = Ratio of length of connecting rod to crank pin radius

crank pin radius = stroke length/2

Program:

The above theory when implemented in MATLAB program is as shown

First a function needs to be declared where we will define engine kinematics.

This function was named engine_parameters

function [V] = engine_parameters(bore,stroke,con_rod,cr,start_crank,end_crank)

a = stroke/2;                 %in metres
R = con_rod/a;                %in metres

v_s = (pi/4)*bore^2*stroke;   %swept volume in m^3
v_c = v_s/(cr-1);             %clearance volume in m^3


%a linearly spaced vector of theta
%it would be used to create volume array and then pressure array
theta = linspace(start_crank,end_crank,200);

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;  %Volume at an angle theta
end

This function helps in creating those curves by relating pressure and volume. But still we need a code where we will use this function to determine state variables at the end of each process.

The main program is shown below. We will call our engine_parameters function here.

%Air standard cycle
clear all
close all
clc

%Inputs
gamma = 1.4;
t3 = 2300;                    %in Kelvin

%State variables
p1 = 101325;                  %in Pa
t1 = 500;                     %in Kelvin


%Engine Geometric Parametres 
bore = 0.1;                   %in metres
stroke = 0.1;                 %in metres
con_rod = 0.15;               %in metres
cr = 12;


%Calculating swept volume and clearance volume
v_s = (pi/4)*bore^2*stroke;   %in m^3
v_c = v_s/(cr-1);             %in m^3

v1 = v_s + v_c;
v2 = v_c;


%State variables at state point 2
%Using equation p2v2^gamma = p1v1^gamma
%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);

%creation of voulme and pressure arrays
constant_c = p1*v1^gamma;
V_compression = engine_parameters(bore,stroke,con_rod,cr,180,360);

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;

%creation of voulme and pressure arrays
constant_c = p3*v3^gamma;
V_expansion = engine_parameters(bore,stroke,con_rod,cr,0,180);

P_expansion = constant_c./V_expansion.^gamma;


%State variables at state point 4
v4 = v1;

%p3v3^gamma = p4v4^gamma | p4 = p3(v3/v4)^gamma
p4 = p3*(v3/v4)^gamma;

t4 = p4*v4*t3/(p3*v3);

%Thermal Efficiency
Thermal_Efficiency = 1-((cr)^(1-gamma))

%Plotting
hold on
grid on
colour = [rand rand rand];
plot(v1,p1,'*','color','r')
plot(V_compression,P_compression)
plot(v2,p2,'*','color','b')
plot([v2 v3],[p2 p3],'color',colour)
plot(v3,p3,'*','color','g')
plot(V_expansion,P_expansion)
plot(v4,p4,'*','color','k')
plot([v4 v1],[p4 p1],'color',colour)
title("Air Standard Cycle: Otto Cycle")
legend("Initial State","isentropic compression","before heat addition","isochroric heat addition","after heat addition","isentropic expansion","heat rejection starts","isochoric heat rejection")
xlabel("Volume (in cu. meters)")
ylabel("Pressure (in Pascal)")

Working of Code:

First we define the engine parameters and state variables at initial condition. We apply governing equation and ideal gas equation to compute pressure and volumes at different states. But, in this way we get only four values of pressure and volume which isn't enough. Hence, our function calculates volume as a function of angle, and its corresponding pressure is determined. All this values are plotted sequentially to obtain PV diagram (work flow diagram) of Otto cycle. Various steps have been explained using comments in code (as and when necessary).

Outputs:

Thermal efficiency

PV Diagram

Workspace

Result and Conclusion:

• The PV diagram of Otto Cycle has been plotted successfully.
• The air standard efficiency (thermal efficiency) was found to be 62.989% (for given inputs).
• Hence, air standard efficiency depends only on compression ratio.
• For high compression ratio, air standard efficiency is high.
• But at higher compression ratios knocking occurs, hence an optimum compression ratio must be selected.

References:

• MATLAB for Mechanical Engineers - Course provided by Skill Lync.
• GNU Octave
• Supporting Material: Thermodynamics - an engineering approach (SIE) - Yunus Cengel, Michael Boles, Mehmet Kanoglu
• Supporting Material: Energy Engineering - VP Katekar Sir

Material: Code and Diagrams can be accessed here