P-V Diagram Visualizer for Air Standard Cycle Using MATLAB

AIM: To create a PV diagram and output the thermal efficiency of the engine using MATLAB.

INTRODUCTION:

The entire modern petrol engine works on Otto cycle. It consist of four processes, Two isentropic (reversible adiabatic) processes and two isochoric (constant volume) processes. The four processes of this cycle is as follows:

Process 1-2:Isentropic compression ( reversible adiabatic)
Process 2-3:Constant volume (Isochoric) heat addition
Process 3-4:Isentropic Expansion (reversible adiabatic)
Process 4-1:Constant volume heat rejection

PROCEDURE:

• The basic known parameters such as p1 and t1 ae given as input and t3 is assumed.
• Then the engine geometric parameters are entered.
• The Swept volume and Clearance Volume of the engine is then calculated using the formula given below,

`V_s = (pi/4)*(bore^2)*stroke` and `V_c = V_s/(cr-1)`

• Considering the state variable at each point and applying ideal gas law, we can determine p2, p3, p4, and t2 using the following formula.

For state variables, `p_2*v_2^gamma= p_1* v_1^gamma`

Applying Ideal Gas Law, `p_1*v_1/t_1=p_2*v_2/t_2`

• Create a function for engine kinematics using the engine kinematics formula seperated and designated to various terms,

`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`

• The function is called into the already existing function using the calling fucntion V_compression= engine_kinematics (bore,stroke,con_rod,cr,180,0)
  P_compression= constant_c./V_compression.^gamma for compession stroke and V_expansion= engine_kinematics (bore,stroke,con_rod,cr,0,180).
  P_expansion= constant_c./V_expansion.^gamma for expansion stroke.
• The constant_c is calcualted using the formula, `constant_c= p1*v1^gamma`and`constant_c= p3*v3^gamma`.
• The plotting is done using the plot command which is represented as figure 1.
• Finally, calculate the efficiency using the formula,

`E=1-(t_1/t_2)`

CODE 1 : Otto Cycle

clear all 
close all
clc

%Input
gamma= 1.4;
%Assume
t3= 2300;

% State Variables
p1= 101325;
t1= 500;

%Engine Geometric Parameters
bore= 0.1;
con_rod= 0.15;
stroke= 0.1;
cr=12;

%Calculate Swept Volume and Clearance Volume
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

%p2v2^gamma=p1v1^gamma
%p2= p1*cr^gamma
p2= p1*cr^gamma;

%Applying Ideal Gas Law
%p1v1/t1 =p2v2/t2 | t2= p2*v2*t1/(p1*v1)
 t2= p2*v2*t1/(p1*v1);
 
 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;
 
 %Applying Ideal Gas Law
 %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;
 
 %p3v3^gamma=p4v4^gamma  | p4= p3*(v3/v4)^gamma
 p4= p3*(v3/v4)^gamma;
 
 figure (1)
 hold on
 plot(v1,p1,'*','color','r');
 plot(V_compression,P_compression,'color','g');
 plot(v2,p2,'*','color','r');
 plot(v3,p3,'*','color','r');
 plot(V_expansion,P_expansion,'color','g');
 plot(v4,p4,'*','color','r');
 plot([v1 v4],[p1 p4],'color','g');
 plot([v2 v3],[p2 p3],'color','g');
 grid on
 xlabel('Volume')
 ylabel('Pressure')
 title('Otto Cycle P-V Diagram') 
 
 % To find efficieny
E=1-(t1/t2)
 

CODE 2: Engine Kinematics Function

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= 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

ERRORS:

In the above image, the function was incooect due to a misatke in the formula used in the engine kinematics function. In the formula pi/2 was used instead of pi/4.

OUTPUT:

Efficiency:

nth= 62.99%