Estimation of kinematic parameters of Otto cycle using MATLAB

AIM:- To create a program for solving the Otto cycle.

PROCEDURE:-

1.Load the input values provided (π(pie),γ(gamma), temperature, pressure, and other values provided )

2.Load the engine parameters(bore, stroke, connecting rod, and compression ratio)

3.Load the formulas that are required to find state point values, Like volume, pressure, and temperature in the first two conditions with the help of engine parameters(area of the bore, swept volume, clearance volume).

4.After loading data for finding out volume at two-state points, we automatically get all the volume for other state points due to constant heat addition and heat rejection.

5. With the help of PV=nRT and PV^gamma(only for the adiabatic process) relation we can get the remaining unknown pressure values.

6. As we need to Trace the PV graph for points that don't come under state point, we will create a volume trace function.

7. The trace will be a function of stroke, bore length, compression ratio, connecting rod, and angle of connecting rod.

8. Plot the volume to theta graph to check if ur function handle is working. Then add the 'end' command at the bottom of the handle to create the function

9. Use this handle to plot all the Adiabatic(Compression and expansion) and Isochoric(Heat addition and heat rejection) processes.

10. For the Adiabatic process the relation of PV^γ is used to equate values of pressure values. During the compression process, the piston moves from BDC to TDC so the starting and ending angles are 180 and 0, whereas in the expansion process it's the reverse(TDC to BDC), and angles would be 0 and 180 respectively.

11.Whereas for the isochoric processes, the piston is assumed to be in a static state (At TDC during heat addition and at BDC during heat rejection). The starting and ending angles for both would be 0,0 and 180,180.

12.Since the Volume remains constant, the pressure growth will be linear in nature and so we can determine the values using the linspace command by dividing it between the two pressure points. To match it with the array of volume trace we divide it into 100 parts.

13. Now That we have all the values, Determine the efficiency of the cycle by loading the formula " 1-1/compression ratio^gamma ".

14. Plot the Pressure vs Volume plot by giving suitable labels and legends.

15. Result is stored in Workspace and we get the efficiency as 62.99%.

%Clearing the previous data
clear all
close all
clc

%Otto_cycle_calculations
%Inputs
G=1.4;%Gamma
pi = 3.14;%π
t3 = 2300;%temperature at state point 3

%State_point_1
p1 = 101325;%Pressure at state point 1
t1 = 500;%Temperature at state point 1

%Engine parameters
bore = 0.1;%diameter of bore
stroke = 0.1;%stroke length
c_rod = 0.15;%connecting rod
cR = 12;%compression ratio

%Calculations
area_of_bore = (pi/4)*bore^2;%(π/4)*D^2
swept_vol = area_of_bore*stroke;%((π/4)*D^2)*L
clearance_vol = swept_vol/(cR-1);

%State 1
v1 = swept_vol + clearance_vol;%Volume at state point 1
v2 = clearance_vol;%Volume at state point 2

%State 2
%p1/p2 = (v1/v2)^G, cR=(v1/v2)
p2 = (p1*cR^G);%Pressure at state point 2
t2 = p2*v2*t1/(p1*v1);%Temperature at state point 2

%Adiabatic compression
constant_c = p1*v1^G; %Term for simplifying calculation
%Compression process piston moves from BDC to TDC
%Starting angle = 180, Ending angle = 0 
v_compression = engine_kinematics2(bore,stroke,c_rod,cR,180,0);
p_compression = constant_c./(v_compression.^G);

%State 3 
v3 = v2;%constant volume addition
p3 = p2*t3/t2;%PV=nRT, V,n,R constant

%Constant volume heat addition 
%No movement of Connecting rod Starting angle= 0,Ending angle = 0
v_isochoric1 = engine_kinematics2(bore,stroke,c_rod,cR,0,0);
p_isochoric1 = linspace(p2,p3,100);%Dividing pressure into equal numbers of hundred between 2 points

%State 4 
v4 = v1;%Constant volume heat rejection
%P1V1^G = P2V2^G
p4 = (p3*v3^G)/(v4^G);
%PV=nRT (n,R)= constant
l = (p4*v4)/(p3*v3);%use of constant for simplification 
t4 = l*t3;

%Adiabatic expansion
constant_c2 = p3*v3^G;%Term for simplifying calculation
%Expansion process piston moves from TDC to BDC
%Starting angle = 0, Ending angle = 180 
v_expansion = engine_kinematics2(bore,stroke,c_rod,cR,0,180);
p_expansion = constant_c2./(v_expansion.^G);%PaVa^G = PbVb^G

%Constant volume heat rejection
%No movement of Connecting rod Starting angle= 180 , Ending angle = 180
v_isochoric2 = engine_kinematics2(bore,stroke,c_rod,12,180,180);
p_isochoric2 = linspace(p4,p1,100);%Dividing pressure into equal numbers of hundred between 2 points


%Efficiency 
%n = 1-1/cR^(G-1)
term1 = cR^(G-1);%denomenator term
n = 1-(1/term1);%substitution

%Plotting
figure(1)
hold on
plot(v1,p1,'*','color','r')
plot(v_compression,p_compression)
plot(v2,p2,'*','color','r')
plot(v_isochoric1,p_isochoric1);
plot(v3,p3,'*','color','r')
plot(v_expansion,p_expansion)
plot(v4,p4,'*','color','r')
plot(v_isochoric2,p_isochoric2);
xlabel('Volume')
ylabel('Pressure')
title('OTTO CYCLE')
grid on
legend('Point 1','Adiabatic compression','Point 2','Isochoric heat addition','Point 3','Adiabatic expansion','Point 4','Isochoric heat rejection')
---------------------------------------------------------------------------------------
%Function handle for Volume trace
function [V] = engine_kinematics2(bore,stroke,c_rod,cR,start_crank,end_crank)

a = stroke/2;%crankpin radius
R = c_rod/a;%Ratio of connecting rod to crankpin radius

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

theta = linspace(start_crank,end_crank,100);%Angle of connnecting rod
%Substitution
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 any given point
end

OUTPUT:-