otto cycle pv diagram

%Inputs
gamma=1.4
t3=2300
% state variables
p1=101325
t1=500
% Engine geometric parameters
bore=0.1;
stroke=0.1;
con_rod=0.15;
cr=12;
% calculating the 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*(v1/v2)^gamma=p1*cr^gamma
p2=p1*cr^gamma;
constant_c=p1*v1^gamma;

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

V_compression= engine_kinematics(bore,stroke,con_rod,cr,180,0)
P_compression=constant_c./V_compression.^gamma;
%p1v1/t1=p2v2/t2 | t2=p2*v2*t1/(p1*v1)
t2= p2*v2*t1/(p1*v1);
% 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
%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)
plot(v2,p2,\'*\',\'color\',\'r\')
plot(v3,p3,\'*\',\'color\',\'r\')
plot(V_expansion,P_expansion)
plot(v4,p4,\'*\',\'color\',\'r\')
figure(2)
plot([v1 v2 v3 v4 v1],[p1 p2 p3 p4 p1],\'color\',\'b\')