Project 2 - Rankine cycle Simulator

Aim: Write a code to simulate Rankine cycle.

The Rankine cycle is a model used to predict the performance of steam turbine systems. It was also used to study the performance of reciprocating steam engines. The Rankine cycle is an idealized thermodynamic cycle of a heat engine that converts heat into mechanical work while undergoing phase change. It is an idealized cycle in which friction losses in each of the four components are neglected.

There are four processes in the Rankine cycle. The states are identified by numbers (in brown) in the T–s diagram.

• Process 1–2: The working fluid is pumped from low to high pressure. As the fluid is a liquid at this stage, the pump requires little input energy.

In other words Process 1-2 is [Isentropic compression]

• Process 2–3: The high-pressure liquid enters a boiler, where it is heated at constant pressure by an external heat source to become a dry saturated vapour. The input energy required can be easily calculated graphically, using an enthalpy–entropy chart (h–s chart, or Mollier diagram), or numerically, using steam tables.

In other words Process 2-3 is [Constant pressure heat addition in boiler]

• Process 3–4: The dry saturated vapour expands through a turbine, generating power. This decreases the temperature and pressure of the vapour, and some condensation may occur. The output in this process can be easily calculated using the chart or tables noted above.

In other words Process 3-4 is [Isentropic expansion]

• Process 4–1: The wet vapour then enters a condenser, where it is condensed at a constant pressure to become a saturated liquid.

In other words Process 4-1 is [Constant pressure heat rejection in condenser]

Matlab code:

clear all
close all
clc

disp('Rankine cycle simulator')
disp('1-2 is isentropic expansion in turbine')
disp('2-3 is constant pressure heat rejection in condenser')
disp('3-4 is isentropic compression in pump')
disp('4-1 is constant pressure heat addition in boiler')

p1 = input('enter the pressure at the inlet of turbine (in bar)=   ');
t1 = input('enter the temperature at the inlet of the turbine (in degree C)=   ');
p2 = input('enter the pressure at the inlet of the condenser (in bar)=   ');


%steam condition at the inlet of turbine
h1 = XSteam('h_pt',p1,t1);
s1 = XSteam('s_pt',p1,t1);

%steam condition at the outlet of turbine
sg2 = XSteam('sV_p',p2);
sf2 = XSteam('sL_p',p2);
sfg2 = sg2-sf2;
x2 = (s1-sf2)/sfg2;
s2 = sf2+x2*sfg2;

hf2 = XSteam('hL_p',p2);
hg2 = XSteam('hV_p',p2);
hfg2 = hg2-hf2;
h2 = hf2+x2*hfg2;

t2 = XSteam('Tsat_p',p2);

%Workdone by the Turbine
W_t = h1-h2;

%Steam condition at the outlet of condenser
p3=p2;
t3 = XSteam('Tsat_p',p2);
h3 = XSteam('hL_p',p2);
s3 = XSteam('sL_p',p2);
vf = XSteam('vL_p',p2);
t22 = t3;

%Calculation of pump work
W_p = vf*(p1-p2)*100;

%steam condition at state 4
p4 = p1;
h4 = W_p+h3;
s4 = s3;
%h4 = XSteam('h_ps',p4,s4);
t4 = XSteam('T_ph',p4,h4);

%work_net
W_net = W_t-W_p;

%Back work ratio
back_work_ratio = (h2-h1)/(h3-h4);

%heat supplied
Q_s = h1-h4;

%heat rejected
Q_r = h2-h3;

%efficiency
efficiency = (W_net/Q_s);

%specific steam consumption
SSC = 3600/W_net;
tf_sat = XSteam('Tsat_p',p1);
sf_sat = XSteam('sL_p',p1);
hf_sat = XSteam('hL_p',p1);
sg_sat = XSteam('sV_p',p1);
hg_sat = XSteam('hV_p',p1);

%Plotting T-s diagram
hold on

%for plotting saturation curve
tsat1 = linspace(0,374,1000);
for i=1:length(tsat1)
    ss1(i) = XSteam('sL_T',tsat1(i));
end

for j=1:length(tsat1);
    ss2(j) = XSteam('sV_T',tsat1(j));
end

plot (ss1,tsat1,'b');
plot (ss2,tsat1,'r');

%Saturation temp at p1
temp = XSteam('Tsat_p',p1);
t = linspace(t4,temp,30);
for i=1:length(t)
    S(i)=XSteam('s_pt',p1,t(i));
    hold on
end

t0 = linspace(temp,t1,300);
for i=1:length(t0);
    S0(i)=XSteam('s_pt',p1,t0(i));
end

%Plotting turbine work
plot([s1 s2], [t1 t2],'g');

%Plotting heat rejection in condenser
plot([s2 s3], [t2,t3],'b');

%Plotting work done by pump
plot([s3,s4],[t3,t4],'r');

%plotting the curved line from outlet of the pump to the point where the
%phase change starts
plot(S,t,'y');

%plotting the phase change region
plot([sf_sat sg_sat],[tf_sat tf_sat],'y');

%plotting the curved constant pressure line fro where the phase change ends
%to the superheated point
plot(S0,t0,'y');
legend('Saturation curve L', 'Saturation curve V', 'Turbine', 'condenser', 'pump', 'boiler', 'location','northwest');
title('T-s Diagram');
xlabel('Entropy');
ylabel('Temperature');

%Plotting H-s diagram
figure(2)
hold on

%Plotting the saturation curve
tsat1 = linspace(0,374,1000);

for i=1:length(tsat1)
    hh1(i) = XSteam('hL_T',tsat1(i));
end

for j = 1:length(tsat1);
    hh2(j)=XSteam('hV_T',tsat1(j));
end

plot(ss1,hh1,'b')
plot(ss2,hh2,'r')

%plotting turbine
plot([s1 s2], [h1 h2],'g')

%plotting condenser
plot([s2 s3], [h2 h3],'b')

%plotting pump
plot([s3 s4], [h3 h4],'r')

s = linspace(s3,s1,1000);
for i=1:length(s)
    h(i) = XSteam('h_ps',p1,s(i));
end

%plotting boiler
plot([s4 s1],[h4 h1],'y');
%plot(s,h,'y');
legend('Saturation curve L', 'Saturation curve V', 'Turbine', 'condenser', 'pump', 'boiler', 'location','northwest');
title('H-s diagram')
xlabel('Entropy')
ylabel('Enthalpy')

%output at state 1
fprintf('At state point 1:\r\n')
fprintf('P1 is : %f Bar\r\n',p1)
fprintf('T1 is : %f Deg Celsius\r\n',t1)
fprintf('H1 is : %f KJ/Kg\r\n',h1)
fprintf('S1 is : %f KJ/KgK\r\n',s1)

%output at state 2
fprintf('At state point 2:\r\n')
fprintf('P2 is : %f Bar\r\n',p2)
fprintf('T2 is : %f Deg Celsius\r\n',t2)
fprintf('H2 is : %f KJ/Kg\r\n',h2)
fprintf('S2 is : %f KJ/KgK\r\n',s2)
fprintf('X is : %f \r\n',x2)

%output at state 3
fprintf('At state point 3:\r\n')
fprintf('P3 is : %f Bar\r\n',p3)
fprintf('T3 is : %f Deg Celsius\r\n',t3)
fprintf('H3 is : %f KJ/Kg\r\n',h3)
fprintf('S3 is : %f KJ/KgK\r\n',s3)

%output at state 4
fprintf('At state point 4:\r\n')
fprintf('P4 is : %f Bar\r\n',p4)
fprintf('T4 is : %f Deg Celsius\r\n',t4)
fprintf('H4 is : %f KJ/Kg\r\n',h4)
fprintf('S4 is : %f KJ/KgK\r\n',s4)

%Calculations
fprintf('Work done by turbine: %f KJ/Kg\r\n',W_t)
fprintf('Work done by pump: %f KJ/Kg\r\n',W_p)
fprintf('Net work done: %f KJ/Kg\r\n',W_net)
fprintf('Back work ratio: %f \r\n',back_work_ratio)
fprintf('Heat added in the boiler: %f KJ/kg\r\n',Q_s)
fprintf('Heat rejection in the condenser: %f KJ/Kg\r\n',Q_r)
fprintf('Thermal Efficiency is: %f %%\r\n',efficiency*100)
fprintf('Specific steam consumption: %f KJ/KWh\r\n',SSC)