Project 2 - Rankine cycle Simulator

clear all
close all
clc

% step 1
% print description of all process

fprintf('RANKINE CYCLE SIMULATOR');
fprintf('1-2 Isentropic expansion in the turbine');
fprintf('Constant pressure heat rejection by the condenser');
fprintf('Isentropic compression in the pump');
fprintf('Constant pressure heat addition by the boiler');

% Step 2
% Taking all input arguments

P1 = input('Enter the pressure at the turbine inlet(in bar):');
T1 = input('Enter the temperature at the turbime inlet(in degree celcius):');
P2 = input('Enter the pressure at the condenser(in bar):');

fprintf('RESULT');

% Calculation for general properties
% Step 3
% Calculation of cp value at p1 & p2 for saturated liquid & vapour
CP1_L=XSteam('CpL_p',P1);
CP1_V=XSteam('CpV_p',P1);
CP2_L=XSteam('CpL_p',P2);
CP2_V=XSteam('CpV_p',P2);

% Step 4
% Calculation of sauration temperature at P1 & P2
Tsat_P1=XSteam('Tsat_p',P1);
Tsat_P2=XSteam('Tsat_p',P2);

% Calculation of enthalpy & entropy at saturation points
% Step 5
% Calculation of saturation properties at state 1
Hg_1=XSteam('hv_p',P1);
Sg_1=XSteam('sv_p',P1);

% Step 6
% Calculation of saturation properties at state 2
Hg_2=XSteam('hv_p',P2);
Sg_2=XSteam('sv_p',P2);

% Step 7
% Calculation of saturation properties at state 3
Hf_3=XSteam('hl_p',P2);
Sf_3=XSteam('sl_p',P2);

% Step 8
% Calculation of saturation properties at state 4
Hf_4=XSteam('hl_p',P1);
Sf_4=XSteam('sl_p',P1);

% Calculation of process properties at state 1
% Step 9
fprintf('at state point 1');
fprintf('P1 is:%f Bar', P1);
fprintf('T1 is:%f degC', T1);
H1=Xsteam('h_pT', P1, T1);
fprintf('H1 is :%f kj/kg',H1);
S1=Xsteam('s_pT', P1, T1);
fprintf('S1 is :%f kj/kg',S1);

% Calculation of process properties at state 2
% Step 10
fprintf('at state point 2');
fprintf('P2 is:%f Bar', P2);
S2=S1;
fprintf('S2 is:%f kj/kgK', S2);   

% Calculating value of X
if(S2<=Sg_2)

    X2=(S2-Sf_3)/(Sg_2-Sf_3);
else
    X2=1;
end
fprintf('X2 is:%f n',X2)

% Calculation of value of T2
if(X2<1);
    T2=Tsat_P2;
else
%S2=Sg_2+CP2_V*(log(T/Tsat_P2))
T2=(Tsat_P2+273.15)*exp((S2-Sg_2)/CP2_V);
T2=T2-273.15;
end
fprintf('T2 is:%f degC',T2);
if(X2<1);
    H2=Hf_3 + X2*(Hf_3-Hg_2);
else
    H2=Xsteam('h_pT',P2,T2);
end
fprintf('H2 is:%f kj/kg',H2)


% Calculation of process properties at state 3
% Step 11

fprintf('at state point 3');
P3=P2;
fprintf('P3 is:%f Bar',P3)
V3 = XSteam('vL_p', P3);
T3=Tsat_P2;
fprintf('T3 is:%f degC',T3)
H3=Hf_3;
fprintf('H3 is : %f kj/kg',H3)
S3=Sf_3;
fprintf('S3 is : %f kj/kgK',S3) 

% Calculation of process properties at state 4
% Step 12

fprintf('at state point 4');
P4=P1;
fprintf('P4 is:%f Bar',P4)
S4=S3;
fprintf('S4 is:%f kj/kgK',S4);
Wp=V3*(P4-P3)*100;
H4=Wp + H3;
fprintf('H4 is: %f kj/kg',H4);
T4=(Tsat_P1+273.15)/(exp((Sf_4-S4)/CP1_L));
T4=T4-273.15;
fprintf('T4 is: %f degC',T4)

% Calculation for plotting saturation liquid & vapour lines
% Step 13

temp=linspace(1,T1,1000);
for i=1:length(temp)
    H_L(i)=Xsteam('hL_T',temp(i));
    H_V(i)=Xsteam('hV_T',temp(i));
    S_L(i)=Xsteam('sL_T',temp(i));
    S_V(i)=Xsteam('sV_T',temp(i));
end

% Plotting T-S & H-S Diagram
% Step 14
% Plotting T-S Plot

figure(1)
hold on
plot([S1 S2],[T1 T2],'r','linewidth',2)
text(S1,T1,'1')
text(S2,T2,'2')
plot([S2 Sg_2],[T2 Tsat_P2],'r','linewidth',2)
text(Sg_2,Test_P2,'2s')
plot([Sg_2 S3],[Tsat_P2 T3],'r','linewidth',2)
text(S3,T3,'3')
plot([S3 S4],[T3 T4],'r','linewidth',2)
text(S4,T4,'4')
plot([S4 Sf_4],[T4 Tsat_P1],'r','linewidth',2)
text(Sf_4,Tsat_P1,'4s')
plot([Sf_4 Sg_1],[Tsat_P1 Tsat_P1],'r','linewidth',2)
text(Sg_1,Tsat_P1,'1s')
plot([Sg_1 S1],[Tsat_P1 T1],'r','linewidth',2)
plot(S_L,temp,'b','linestyle','--','linewidth', 1)
plot(S_V,temp,'b','linestyle','--','linewidth', 1)

title('T-S Diagram of ideal rankine cycle')
ylabel('Temperature (degC)')
xlabel('specific entropy (Kj/kgK)')

% Step 15
% Plotting H-S Plot

figure(2)
hold on
plot([S1 S2],[H1 H2],'r','linewidth',2)
text(S1,H1,'1')
text(S2,H2,'2')
plot([S2 Sg_2],[H2 Hg_2],'r','linewidth',2)
text(Sg_2,Hg_2,'2s')
plot([Sg_2 S3],[Hg_2 H3],'r','linewidth',2)
text(S3,H3,'3')
plot([S3 S4],[H3 H4],'r','linewidth',2)
text(S4,H4,'4')
plot([S4 Sf_4],[H4 Hf_4],'r','linewidth',2)
text(Sf_4,Hf_4,'4s')
plot([Sf_4 Sg_1],[Hf_4 Hg_1],'r','linewidth',2)
text(Sg_1,Hg_1,'1s')
plot([Sg_1 S1],[Hg_1 H1],'r','linewidth',2)
plot(S_L,H_L,'b','linestyle','--','linewidth', 1.5)
plot(S_V,H_V,'b','linestyle','--','linewidth', 1.5)

title('H-S Diagram of ideal rankine cycle')
ylabel('Specific enthalpy (kj/kg)')
xlabel('specific entropy (Kj/kgK)')

% workdone, heat, SSC & efficiency calculations
%step 16

Wt = H1-H2;
fprintf('Work by turbine is: %f kj/kg',Wt)
fprintf('Work by compressor is: %f kj/kg',Wp)
Wnet= Wt-Wp;
fprintf('Net workdone is: %f kj/kg',Wnet)
Qin=H1-H4;
eff_thermal=Wnet*100/Qin;
fprintf('Thermal efficiency is: %f',eff_thermal);
SSC=3600/Wt
fprintf('SSC is: %f kg/kWh',SSC)