Project 2 - Rankine cycle Simulator

CODE SIMULATION:

%RANKINE CYCLE SIMULATOR
%CALCULATING THE NETWORK OF THE RANKINE CYCLE
%CALCULATING THE BACK WORK RATIO OF THE RANKINE CYCLE
%TAKING VALUES FROM THE USER
 
clear all
close all
clc
 
%DISPLAYING THE CONTENT 
fprintf("\n                                RANKINE CYCLE SIMULATOR            \n");
fprintf("\n Process 1-2 : Isentropic expansion in the Turbine \n");
fprintf("\n Process 2-3 : Isobaric heat rejection by the Condenser \n");
fprintf("\n Process 3-4 : Isentropic compression in the Pump \n");
fprintf("\n Process 4-1 : Isobaric heat addition by the Boiler \n");
 
%REQUESTING THE USER TO GIVE THE INPUTS
p1 = input('\n Enter the pressure at the Turbine inlet (in bar) = ');
t1 = input('\n Enter the temperature at the Turbine inlet (in Celsius) = ');
p2 = input('\n Enter the pressure at the Condenser (in bar) = ');
 
%TURBINE INLET PARAMETERS
tsat_1 = XSteam('Tsat_p',p1);
if t1 > tsat_1
    h1 = XSteam('h_pt',p1,t1);
    s1 = XSteam('s_pt',p1,t1);
elseif t1 == tsat_1
        h1 = XSteam('hV_p',p1);
        s1 = XSteam('sV_p',p1);
    else 
        fprintf("\n Enter higher value of Turbine inlet temperature (in Celsius)!!! \n\n");
end
 
 
%CONDITIONS AT TURBINE OUTLET
%BASED ON THE LOCATION OF POINT 2 THE FOLLOWING VALUES/PROPERTIES ARE
%EVALUATED
s2 = s1; %isentropic expansion in the turbine
ssat_2 = XSteam('sV_p',p2);
if s2 < ssat_2
    sf_2 = XSteam('sL_p',p2);
    sg_2 = XSteam('sV_p',p2);
    %let x denote the dryness fraction of the steam in the system
    x = (s2 - sf_2)/(sg_2 - sf_2);
    hf_2 = XSteam('hL_p',p2);
    hg_2 = XSteam('hV_p',p2);
    h2  = hf_2 + x*(hg_2 - hf_2);
    t2 = XSteam('Tsat_p',p2);
else
    x = 1;
    h2 = XSteam('hV_p',p2);
    t2 = XSteam('Tsat_p',p2);
end
 
 
%CONDITIONS AT PUMP INLET
p3 = p2;
t3 = XSteam('Tsat_p',p3);
h3 = XSteam('hL_p',p3);
s3 = XSteam('sL_p',p3);
v3 = XSteam('vL_p',p3);
 
 
%CONDITIONS AT BOILER INLET
p4 = p1;
s4 = s3; 
v4 = XSteam('v_ps',p4,s4);
cp4 = XSteam('Cp_ps',p4,s4);
cv4 = XSteam('Cv_ps',p4,s4);
gamma = cp4/cv4;
G = (gamma - 1)/gamma;
t4 = t3*(p2/p1)^G;
 
%WORK REQUIRED BY THE PUMP
W_p = v3*(p4 - p3)*100;
%W_p = h4 - h3
h4 = W_p + h3;
 
 
 
%CONDITION OF INTERMEDIATE STATES: STATE 5 AND STATE 6 IN (4-1) STAGE
t5 = XSteam('Tsat_p',p4);
s5 = XSteam('sL_p',p4);
t6 = t5;
s6 = XSteam('sV_p',p4);
 
 
%WORK DONE BY THE TURBINE
W_t = h1 - h2;
 
%HEAT INPUT
Q_in = h1 - h4;
 
%HEAT OUTPUT
Q_out = h2 - h3;
 
%NET-WORK OBTAINED
W_net = W_t - W_p;
 
%Efficiency of the Rankine Cycle
efficiency = (W_net*100)/Q_in;
 
%BACK WORK RATIO
BWR = W_p/W_t;
 
%SPECIFIC STEAM COSUMPTION
SSC = 3600/W_net;
 
 
%DISPLAYING THE OBTAINED INFORMATION
fprintf("\n                                RESULTS            \n");
fprintf("At State Point 1:                              At State Point 4: \n P1 = %.2f bar                                 P4 = %.2f bar \n T1 = %.2f C                                  T4 = %.2f C  \n H1 = %.2f kJ/kg                             H4 = %.2f kJ/kg  \n S1 = %.2f kJ/kg-C                              S4 = %.2f kJ/kg-C \n",p1,p4,t1,t4,h1,h4,s1,s4);
fprintf("\nAt State Point 2:                              At State Point 3: \n P2 = %.2f bar                                  P3 = %.2f bar \n T2 = %.2f C                                   T3 = %.2f C  \n H2 = %.2f kJ/kg                             H3 = %.2f kJ/kg  \n S2 = %.2f kJ/kg-C                              S3 = %.2f kJ/kg-C \n",p2,p3,t2,t3,h2,h3,s2,s3);
 
fprintf("\n Work required by the pump in the Rankine Cycle = %.2f kJ/kg\n",W_p);
fprintf("\n Work produced by the turbine in the Rankine Cycle = %.2f kJ/kg\n",W_t);
fprintf("\n Net Work obtained in the Rankine Cycle = %.2f kJ/kg\n",W_net);
fprintf("\n Heat input to the Rankine Cycle = %.2f kJ/kg\n",Q_in);
fprintf("\n The efficiency of the Rankine Cycle = %.2f \n",efficiency);
fprintf("\n Back Work Ratio in the Rankine Cycle = %.9f \n",BWR);
fprintf("\n Specific Steam Consumption in the Rankine Cycle = %.2f kg/kWh \n",SSC);
 
 
%PLOTTING THE CURVES
 
%Plotting Rankine Cycle TS DIAGRAM
%Plotting Saturation Curve
T = linspace(0,373.946,100); % critical point of water occurs at 373.946 degreeC
figure(1);
hold on
for i = 1 : length(T)
    s_l(i) = XSteam('sL_T',T(i));
    s_g(i) = XSteam('sV_T',T(i));
    h_l(i) = XSteam('hL_T',T(i));
    h_g(i) = XSteam('hV_T',T(i));
end
 
plot(s_l,T,'linewidth',3,'color','g');
plot(s_g,T,'linewidth',3,'color','g');
plot([s1 s2],[t1 t2],'linewidth',1,'color','k');
text(s1,t1,'1'); 
plot([s2 s3],[t2 t3],'linewidth',1,'color','r');
text(s2,t2,'2');
plot([s3 s4],[t3 t4],'linewidth',1,'color','k');
text(s3,t3,'3');
plot([s4 s5],[t4 t5],'linewidth',1,'color','r');
text(s4,t4,'4');
plot([s5 s6],[t5 t6],'linewidth',1,'color','r');
text(s5,t5,'5');
plot([s6 s1],[t6 t1],'linewidth',1,'color','r');
text(s6,t6,'6');
xlabel('Entropy (kJ/kg-C)');
ylabel('Temperature (C)');
title('Entropy versus Temperature Curve - Rankine Cycle');
legend('Entropy (s)kJ/kg-C','Temperature (t)C');
grid on
hold off
 
 
%Plotting Rankine Cycle HS DIAGRAM
%Plotting Saturation Curve
figure(2);
hold on
plot(s_l,h_l,'linewidth',3,'color','c');
plot(s_g,h_g,'linewidth',3,'color','c');
plot([s1 s2],[h1 h2],'linewidth',1,'color','k');
text(s1,h1,'1');
plot([s2 s3],[h2 h3],'linewidth',1,'color','r');
text(s2,h2,'2');
plot([s3 s4],[h3 h4],'linewidth',1,'color','k');
text(s3,h3,'3');
plot([s4 s1],[h4 h1],'linewidth',1,'color','r');
text(s4,h4,'4');
xlabel('Entropy (kJ/kg-C)');
ylabel('Enthalpy (kJ/kg)');
title('Enthalpy versus Entropy Curve - Rankine Cycle');
legend('Entropy (s)kJ/kg-C','Enthalpy (h)kJ/kg');
grid on
hold off