S=1.

091; %(dm)
B=0.915; %(dm)
R=S/2; %(dm)
TSN=13;
Vh=pi*S*(B^2)/4;
Vc=Vh/(TSN-1);
ro = 1;
Vz = ro*Vc;
Va=Vh+Vc;
Vb = Va;
lamda = 1/3;
n1 = 1.37;
n2 = 1.3376;
pa = 0.096235;
pc = 3.23;
pb = 0.3833;
pr = 0.11;
pz = 11;847;
Fp = (pi*(B)^2)/4;
n = 6600;
w = 2*pi*n/60;
mnp = 120.75 ;
mtt = 161.;
mk = 161 ;
mA = mtt/3;
mB = 2*mtt/3;
mt = mnp + mA;
mr = mk + mB;
% do thi PV
% Doan tu 0 den 10
ahc=[Link];
phc=[0.12 (pr+pa)/2 0.135];
a0= linspace (0,10,1000);
x0=R*(1-cosd(a0)+lamda/4*(1-cosd(2*a0)));
v0= x0*(pi*B^2)/4 +Vc;
p0=interp1(ahc, phc, a0,'spline');
J0 = R*w^2.*(cosd(a0)+lamda.*cosd(2.*a0));

plot(v0,p0,'r','linewidth',1.2 )
hold on

%Qua trinh nap

a1 = linspace (10,180,1000);
x1=R*(1-cosd(a1)+lamda/4*(1-cosd(2*a1)));
v1= x1*(pi*B^2)/4 +Vc;
p1=linspace(pa,pa,1000);
J1 = R*w^2.*(cosd(a1)+lamda.*cosd(2.*a1));
plot(v1,p1,'r','linewidth',1.2);
hold on
%Qua trinh nen
a2= linspace (180,350,1000);
x2=R*(1-cosd(a2)+lamda/4*(1-cosd(2*a2)));
v2= x2*(pi*B^2)/4 +Vc;
p2=p1.*((Va./v2).^n1);
J2 = R*w^2.*(cosd(a2)+lamda.*cosd(2.*a2));
plot(v2,p2,'r','linewidth',1.2 );
hold on

%Doan 350 den 362.5

pcc= (pz*0.85-pc)/3+pc;
ahc=[350 356.25 362.5];
phc=[max(p2) (pcc+max(p2))/2 (pcc)];
a22= linspace (350,362.5,1000);
x22=R*(1-cosd(a22)+lamda/4*(1-cosd(2*a22)));
v22= x22*(pi*B^2)/4 +Vc;
p22=interp1(ahc,phc,a22,'spline');
J22 = R*w^2.*(cosd(a22)+lamda.*cosd(2.*a22));
plot(v22,p22,'r','linewidth',1.2 );
hold on

%Qua trinh gian no

a3=linspace (375,500,1000);
x3=R*(1-cosd(a3)+lamda/4*(1-cosd(2*a3)));
v3= x3*(pi*B^2)/4 +Vc;
p3=pz.*((Vc./v3).^n2);
J3 = R*w^2.*(cosd(a3)+lamda.*cosd(2.*a3));
plot (v3,p3,'r','linewidth',1.2 );
hold on
%Doan 362.5 den 375
ahc=[362.5 370 375];
phc=[(pcc) 0.85*pz max(p3)];
a23=linspace (362.5,375,1000);
x23=R*(1-cosd(a23)+lamda/4*(1-cosd(2*a23)));
v23= x23*(pi*B^2)/4 +Vc;
p23=interp1(ahc,phc,a23,'spline');
J23 = R*w^2.*(cosd(a23)+lamda.*cosd(2.*a23));
plot(v23,p23,'r ','linewidth',1.2 );
hold on
%Doan 500 den 580
ahc=[500 520 540];
phc=[min(p3) (min(p3)+(min(p3)+pa)/2)/2 (min(p3)+pa)/2];
a33= linspace (500,540,1000);
x33=R*(1-cosd(a33)+lamda/4*(1-cosd(2*a33)));
v33= x33*(pi*B^2)/4 +Vc;
p33=interp1(ahc,phc,a33,'spline');
J33 = R*w^2.*(cosd(a33)+lamda.*cosd(2.*a33));
plot(v33,p33,'r ','linewidth',1.2 )
hold on
ahc=[540 560 580];
phc=[(min(p3)+pa)/2 ((min(p3)+pa)/2+pr)/2 pr];
a34= linspace (540,580,1000);
x34=R*(1-cosd(a34)+lamda/4*(1-cosd(2*a34)));
v34= x34*(pi*B^2)/4 +Vc;
p34=interp1(ahc,phc,a34,'spline');
J34 = R*w^2.*(cosd(a34)+lamda.*cosd(2.*a34));
plot(v34,p34,'r ','linewidth',1.2 )
hold on

%Qua trinh xa
a4=linspace (580,720,1000);
x4=R*(1-cosd(a4)+lamda/4*(1-cosd(2*a4)));
v4= x4*(pi*B^2)/4 +Vc;
p4= linspace(pr,pr,1000);
J4 = R*w^2.*(cosd(a4)+lamda.*cosd(2.*a4));
plot(v4,p4,'r','linewidth',1.2 );
hold on

%Chu thich
title('DO THI CONG P-V');
xlabel('The tich (dm3)')
ylabel('Ap suat (MN/m2)');
grid on

a=[a0 a1 a2 a22 a23 a3 a33 a34 a4];

x=[ x0 x1 x2 x22 x23 x3 x33 x34 x4 ];
v=[ v0 v1 v2 v22 v23 v3 v33 v34 v4 ];pkt=[ p0 p1 p2 p22 p23 p3 p33 p34 p4
];
J=[ J0 J1 J2 J22 J23 J3 J33 J34 J4 ];
% do thi van toc
figure
V = R/10*w*(sind(a)+lamda*sind(2*a)/2);
plot(a,V,'g','linewidth',1.2);
hold on
xlabel('Alpha (Do)');
ylabel('V (m/s)');
title('DO THI VAN TOC CUA PISTON');
legend('V');
% do thi pkt pj pt
figure
pj = -mt*R*w^2*(cosd(a)+lamda*cosd(2*a))*10^-7;
plot(a,pj,'g','linewidth',1.2)
hold on
plot(a,pkt,'b','linewidth',1.2)
grid on
hold on
pt = pkt +pj ;
plot(a,pt,'r','linewidth',1.2 )
title('DO THI pkt pj p1');
xlabel('Goc quay truc khuyu (Do)');
ylabel('Pkt(MN/m2) pj(MN/m2) p1(MN/m2)');
legend('pj','pkt','p1');
% do thi gia toc
figure
j = R/10*w^2*(cosd(a)+lamda*cosd(2*a));
plot(a,j,'g','linewidth',1.2);
hold on;
xlabel('GQTK (Do)');
ylabel('jp (m/s^2)');
title('DO THI GIA TOC CUA PISTON');
legend('jp');
% do thi s
figure
a0=linspace(0,360,1000);
s=R/10*(1-cosd(a0)+lamda/4*(1-cosd(2*a0)));
plot(a0,s,'g','LineWidth',1.2);
hold on;
xlabel('GQTK (Do)');
ylabel('s (m)');
title('DO CHUYEN VI CUA PISTON');
legend('s');
%% xuat bang excel
%lap bang ap suat MCCT
Vxn = linspace(Va,Vc,50);
Pxn = pa.*(Va./Vxn).^n1;
Vxg = linspace(Vz,Vb,50);
Pxg = pz.*(Vz./Vxg).^n2;
xlswrite('[Link]',[Vxn(:),Pxn(:)]);
xlswrite('[Link]',[Vxg(:),Pxg(:)]);
%lap bang gia tri ket qua tinh toan dong luc hoc
xlswrite('[Link]',[a(:),pkt(:),pj(:),pt(:)]);
%lap bang gia tri ket qua tinh toan dong hoc
xlswrite('[Link]',[a(:),x(:),v(:),J(:)]);