Project 2 - Rankine cycle Simulator

%Printing Introduction
fprintf ('RANKINE CYCLE SIMULATOR \n')
fprintf('\n1-2 is isentropic expansion in the turbine \n')
fprintf('2-3 is constant pressure heat rejection by the condenser \n')
fprintf('3-4 is isentropic compression in pump \n')
fprintf('4-1 is constant pressure heat addition by the boiler \n')

% User Inputs
P1 = input('\n enter the pressure at the turbine inlet (in bar):');
T1 = input('enter the temprature at the turbine inlet (in degree celcius):');
P2 = input('enter the pressure at the turbine outlet(in bar):');

% State point2
S2 = S1;
Sf2 = XSteam('sL_p',P2);
Sg2 = XSteam('sV_p',P2);
x = (S2-Sf2)/(sg2-sf2);

% Dryness fracion
Hf2 = XSteam('hL_p',P2);
Hg2 = XSteam('hV_p',P2);
H2 = Hf2 + x*(Hg2-Hf2);
T2 = XSteam('T-ph',P2,H2);
Ts2 = Xsteam('Tsat_p',P2);

% State Point 3
H3 = Hf2;
S3 = Sf2;
T3 = T2;
P3 = P2;

% State point 4
P4 = P1;
S4 = S3;
H4 = Xsteam('h_ps',P4,S4);
T4 = Xsteam('T_ps',P4,S4);

t=linspace(0,374,500); %temp space for plotting saturation cutves

%Plot T-S Diagram
figure(1)
hold on
% Plot Saturation Curve
fori=1:length(t)
Plot (XSteam('sL_T',t(i)),t(i);:,'color','b')
Plot (XSteam('sV_T',t(i)),t(i);:,'color','b')
end
%Plot Cycle
plot([S1 S2],[T1 T2],'LineWidth',2,'Color','r')
if x>1
T3=Ts2;
plot([S2 Sg2 S3 S4],[T2 Ts2 T3 T4],'linewide',2,'Color','r')
T2=Ts2;
else
plot([S2 S3 S4],[T2 T3 T4],'linewidth',2,'Color','r')
end
l=linspace(T1,T2,1000);
for i=1:length(length)
plot(XSteam('s_pT',P1,l(i)),l(i);:,'color','r')
end
plot([Sf1 Sg1],[Ts1 Ts1],'linewidth',2,'color','r')
xlabel('Entropy[Kj/Kg-K]')
ylabel('Temparature[K]')
title('T-S Diagram')

% Plot Saturation Curve
for i=1:length(t)
plot(XSteam('sL_T',t(i)),XSteam('hL_T',t(i));:,'color','b')
plot(XSteam('sV_T',t(i)),XSteam('hV_T',t(i));:,'color','b')
end
%Plot Cycle
plot([S1 S2 S3 S4 S1],[H1 H2 H3 H4 H1],'linewidth',2,'Color','r')
xlable('Entropy[Kj/Kg-k]')
ylable('Enthalpy[Kj/Kg]')
title('H-S Diagram')

%calcualtions

%Turbine Work
Wt=H1-H2;

%Pump Work
Wp=H4-H3;

%Heat Addition
Qin=H1-H4;

%Heat Rejection
Qout=H2-H3;

%Net Work
Wnet=Wt-Wp;

%Efficiency
Ntherm=(Wnet/Qin)*100;

%Specific steam consumption
SSC=3600/Wnet;

%Back Work Ratio
BWR=Wp/Wt;

%Print Results
fprintf('\n\nResults\n')
fprintf('At State Point 1:\n')
fprintf('P1 is:%.2f bar\n',P1)
fprintf('T1 is:%.2f Deg Celcius\n',T1)
fprintf('H1 is:%.2f KJ/kg\n',H1)
fprintf('S1 is:%.2f KJ/kg-K\n',S1)

fprintf('\nAt State Point 2:\n')
fprintf('P2 is:%.2f bar\n',P2)
fprintf('T2 is:%.2f Deg Celcius\n',T2)
fprintf('H2 is:%.2f KJ/kg\n',H2)
fprintf('S2 is:%.2f KJ/kg-K\n',S2)

fprintf('\nAt State Point 3:\n')
fprintf('P3 is:%.2f bar\n',P3)
fprintf('T3 is:%.2f Deg Celcius\n',T3)
fprintf('H3 is:%.2f KJ/kg\n',H3)
fprintf('S3 is:%.2f KJ/kg-K\n',S3)

fprintf('\nAt State Point 4:\n')
fprintf('P4 is:%.2f bar\n',P4)
fprintf('T4 is:%.2f Deg Celcius\n',T4)
fprintf('H4 is:%.2f KJ/kg\n',H4)
fprintf('S4 is:%.2f KJ/kg-K\n',S4)

fprintf('\nWt is:%.2f Kj/Kg\n',Wt)
fprintf('Wp is:%.2f Kj/Kg\n',Wp)
fprintf('Wnet is:%.2f Kj/Kg\n',Wnet)
fprintf('Thermal Efficiency is:%.2f percent\n',Ntherm)
fprintf('S.S.C. is:%.2f Kg/KWh\n',SSC)
fprintf('Back Work Ratio is:%.3f\n',BWR)