Rankine cycle Simulator

AIM:

Create a Rankine Cycle Simulator using OCTAVE/MATLAB.

OBJECTIVE:

• Create a Rankine Cycle Simulator using OCTAVE/MATLAB.
• calculate the state points of the Rankine Cycle (any type of your choice) based on user inputs.
• plot the corresponding T-s and h-s plots for the given set of inputs.

LIBRARAY USED:

The link for the standard Xsteam library to get the steam table data for water, is given below. It must be used to calculate the points in between the state variables.

Click here to download the Xsteam Library.

Provided Inputs:

• Turbine Inlet temperature = 400° C
• Turbine Inlet Pressure = 30 bars
• Turbine Outlet Pressure = 0.05 bars = Condenser Inlet Pressure

THEORY:

Process 1-2: Reversible adiabatic expansion in the turbine.

Process 2-3: Constant-pressure transfer of heat in the condenser.

Process 3-4: Reversible adiabatic pumping process in the feed pump.

Process 4-1: Constant-pressure transfer of heat in the boiler.

STEPS & ALGORITHM EXPLANATION:

• Taking all inputs arguments.
• Calculation of Cp value at P1 & P2 for saturated liquid & vapor
• Calculation of Saturation Temperature at P1 & P2
• Calculation of saturation properties at state 1, 2, 3 and 4
• Calculation of Process properties at state 1, 2, 3, and 4
• Calculation & Plotting saturation Liquid and Vapour lines.
• Plotting T-s Plot
• Plotting h-s Plot
• Calculation of Work done, heat, SSC, and efficiency

MATLAB Code:

clear all
close all
clc

% Print discription of all process

P1 = input('Pressure at the Turbine Inlet (in bar): ');
T1 = input('Temperature at the Turbine inlet (in Degree celsius): ');
P2 = input('Pressure at the Condenser (in bar): ');

fprintf('RESULT\n');

% 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);

Tsat_P1=XSteam('Tsat_p',P1);
Tsat_P2=XSteam('Tsat_p',P2);

% Calculation of saturation properties at state 1
Hg_1=XSteam('hV_p',P1);
Sg_1=XSteam('sV_p',P1);

% Calculation of saturation properties at state 2
Hg_2=XSteam('hV_p',P2);
Sg_2=XSteam('sV_p',P2);

% Calculation of saturation properties at state 3
Hf_3=XSteam('hL_p',P2);
Sf_3=XSteam('sL_p',P2);

% 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
fprintf('At State point 1\n');
fprintf('P1 is :%f Bar\n',P1);
fprintf('T1 is :%f degC\n',T1);
H1=XSteam('h_pT',P1,T1);
fprintf('H1 is :%f kJ/kg\n',H1);
S1=XSteam('s_pT',P1,T1);
fprintf('S1 is :%f kJ/kgn\n',S1);

% Calculation of Process properties at state 2

fprintf('At state point 2\n');
fprintf('P2 is : %f Bar\n', P2);
S2 = S1;
fprintf('S2 is : %f kJ/kgK\n', S2);

%Calculating value of X
if(S2<=Sg_2)
    X2 = (S2 - Sf_3) / (Sg_2-Sf_3);   % Dryness Fraction at state point 2
else
    X2=1;    
end
fprintf('X2 is :%f \n', X2);

%calculating value of T2
if(X2<1)
    T2=Tsat_P2;
else 
    %S2=Sg_2+CP2_V*(log(T2/Tsat_P2))
    T2=(Tsat_P2+273.15)*exp((S2-Sg_2)/CP2_V) ;
    T2=T2-273.15;
end    
fprintf('T2 is : %f degC \n', 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/kgn\n', H2)

% Calculation of Process properties at state 3
fprintf('At state point 3\n')
P3 = P2;
fprintf('P3 is : %f Bar\n', P3)
V3 = XSteam('vL_p', P3);
T3 = Tsat_P2;    % No subcooling assumed
fprintf('T3 is : %f degC\n', T3)
H3 =Hf_3;
fprintf('H3 is : %f kJ/kg\n', H3)
S3 =Sf_3;
fprintf('S3 is : %f kJ/kgK\n\n', S3)

% Calculation of Process properties at state 4
fprintf('At state point 4\n')
P4 = P1;
fprintf('P4 is : %f Bar\n', P4)
S4 = S3;
fprintf('S4 is : %f kJ/kgK\n', S4)
Wp = V3 * (P4 - P3) * 100;
H4 = Wp + H3;
fprintf('H4 is : %f kJ/kg\n', H4)
T4= (Tsat_P1+273.15)/(exp((Sf_4-S4)/CP1_L));
T4=T4-273.15;
fprintf('T4 is : %f Deg celsius \n\n', T4)

% Calculation & plotting saturation Liquid and Vapour lines

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 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,Tsat_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('Temprature (degC)')
xlabel('Specific Entropy (kJ/kgK)')

% 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)')

mkdir('C:\Users\ASUS\Desktop\matlab\rankIne cYcle\');
cd('C:\Users\ASUS\Desktop\matlab\rankIne cYcle\')
saveas(1,'1');
saveas(2,'2');


% Calculation of Work done, heat, SSC and efficiency
Wt = H1 - H2;
fprintf('Work by turbine is : %f kJ/kg\n', Wt)
fprintf('Work by compressor is : %f kJ/kg\n', Wp)
Wnet = Wt - Wp;
fprintf('Net Wok done is : %f kJ/kg\n', Wnet)
Qin = H1 - H4;
eff_thermal = Wnet*100 / Qin ;    % Effieciency
fprintf('Thermal efficiency is : %f', eff_thermal); disp('%')
SSC = 3600 / Wt;    % Specific steam consuption
fprintf('S.S.C is : %f kg/kWh ', SSC)

RESULTS:

ERRORS :

• Wrong file directory:

• Typing errors
• Difficulty in understanding the 'XSteam.m' library function at first glance