Simulating a Rankine Cycle using MatLab
Objective:
- To Simulate a Rankine cycle using matlab.
- Extracting Steam and water properties from XSteam data (IAPWS IF - 97 standard)
- Calculating all the state points of the cycle based on User inputs.
- Plotting the T vs s and H vs s plot for the given input.
Rankine Cycle :
- Rankine cycle is an idealized thermodynamic cycle or heat engine in which a fuel is used to produce heat within a boiler , converting water into steam which then expands through a turbine producing useful wok . This process was developed in 1859 by scottish engineer william J.m Rankine . It is a cycle which converts heat into mechanical energy which usually is converted into electrical energy .
( Rankine cycle in a nuclear power plant )
(1) Feed Pump:
Saturated liquid which is outcome of the condensor enters into pump which increases its pressure to the boiler pressure isentropically. This process is work consuming which is given as follows :
(2) Boiler:
In boiler, constant pressure heat addition process takes place which increases the enthalpy of the working fluid to sufficient value. The fluid could be heated upto saturated vapour state or could be superheated. The heat added to the working fluid in the boiler can be given as:
(3) Turbine:
Isentropic expansion of the working fluid takes place in the turbine which results in work output. Thus, this process is work producing. The outlet of this process is mostly a mixture of liquid and vapour but it could be saturated vapour depending upon turbine inlet condition. The work produced in this process can be given as:
(4) Condensor:
In condensor, constant pressure heat rejection process convert all the vapour output of the turbine into saturated liquid. It is then entered into feed pump. The heat rejected in the process can be given as:
Processes and Working :
There are four processes in the rankine cycle (order maybe varied ) :
Process 1-2 : Isentropic Expansion in the turbine.
Process 2-3 : Constant pressure Heat rejection in condenser.
Process 3-4 : Isentropic compression in the pump.
Process 4-1 : Constant pressure heat addition in the boiler.
( Plot obtained from matlab program , T vs S plot )
Working of the rankine cycle simulator
- In process 1-2 , The dry saturated vapour expands through a turbine, generating power. This decreases the temperature and pressure of the vapour, and some condensation may occur.
- Process 2-3 , The vapor-liquid mixture leaving the turbine (2) is condensed at low pressure, usually in a surface condenser using cooling water. In well designed and maintained condensers, the pressure of the vapor is well below atmospheric pressure, approaching the saturation pressure of the operating fluid at the cooling water temperature.
- Process 3-4 , The pressure of the condensate is raised in the feed pump. Because of the low specific volume of liquids, the pump work is relatively small and often neglected in thermodynamic calculations.
- Process 4-1 , High pressure liquid enters the boiler from the feed pump (4) and is heated to the saturation temperature (on the saturation curve). Further addition of energy causes evaporation of the liquid until it is fully converted to saturated steam (1) or superheated (1).
( Rankine cycle in a nuclear power plant )
Elaborating (Rankine cycle paragraph point v ) : As the water condenses in the condenser, waste heat is given off in the form of water vapour which can be seen billowing from a plant's cooling towers. This waste heat is necessary in any thermodynamic cycle. Due to this condensation step, the pressure at the turbine outlet is lowered. This means the pump requires less work to compress the water—resulting in higher overall efficiencies.
Code Algorithm :
(1) Option to provide input values are provided using inputdlg command with default input values.
(2) Since the values received are string, they are converted into numeric values using str2num command.
(3) Parameters at different states are evaluated using XSteam function with the help of input parameters.
(4) At state 3, which is turbine inlet, three possibilities are there. The first is provided temperature at the turbine inlet is less than saturation temperature at the provided boiler pressure. In this case data is impractical and no claculation of the other states can be made. The other two cases are that the turbine inlet condition is equal to vapour saturation temperature or the working fluid is superheated. These condition at state 3 is provided using if-elseif command.
(5) Once the thermodynamic parameters at all states are known, net work done and back work ratio is evaluated.
(6) Output is printed using fprintf command.
(7) To plot the curve, it is necessary to know the extreme limits of the pressure which is possible. Along that pressure limits, saturation temperature, saturation liquid entropy and saturation vapour enthalpy is evaluated.
(8) To known the pressure lower limits while loop is used with a condition that pressure value keeps decreasing till it cease to exist. For that purpose isnan command is used.
(9) Similarly, for pressure upper limits, pressure is increased till it cease to exist.
(10) Once the pressure limit is known saturation temperature, saturation liquid entropy and saturation vapour enthalpy is evaluated and using plot command, the results are plotted.
(11) The Rankine cycle states are also plotted in the same figure using hold on command including the possible variation at state 3.
(12) State 2 to state 2' is not varying linearly similar to the case in superheating from state 3 to state 3'. This variation is captured by evaluating values of entropy at different values of temperature from state 2 to state 2'.
(13) Legend, title, axis labels and text are added to the plot for clear representation.
(15) To plot h-s diagram, enthaply values are evaluated between extreme limits of pressure.
(16) The graph is plotted and labelled in same way as for T-s diagram.
Matlab Program :
clear all
close all
clc
% Input paramters
x = inputdlg({'Turbine Inlet temperature','Turbine Inlet Pressure',...
'Turbine Outlet Pressure'},'Inputs', [1 50; 1 50; 1 50],...
{'400','30','0.05'});
turbine_inlet_temp = str2num(x{1});
turbine_inlet_press = str2num(x{2});
turbine_outlet_press = str2num(x{3});
% Paramters at stage 1 i.e. pump inlet
p1 = turbine_outlet_press;
t1 = XSteam('TSat_p',p1);
h1 = XSteam('hL_p',p1);
s1 = XSteam('sL_p',p1);
% Parameters at stage 2 i.e pump outlet
p2 = turbine_inlet_press;
t2 = XSteam('T_ps',p2,s1);
h2 = XSteam('h_pT',p2,t2);
s2 = s1;
% Parameters at 2' i.e. saturated liquid state reached in boiler
p2_dash = p2;
t2_dash = XSteam('TSat_p',p2);
h2_dash = XSteam('hL_p',p2);
s2_dash = XSteam('sL_p',p2);
% Parameter at 3 i.e turbine inlet
tsat3 = t2_dash;
t3 = turbine_inlet_temp;
if t3tsat3
p3 = turbine_inlet_press;
h3 = XSteam('h_pT',p3,t3);
s3 = XSteam('s_pT',p3,t3);
% parameter at saturated vapour state in boiler
p3_dash = p3;
t3_dash = XSteam('TSat_p',p3);
h3_dash = XSteam('hV_p',p3);
s3_dash = XSteam('sV_p',p3);
end
% parameters at state 4 i.e. condensor inlet
p4 = turbine_outlet_press;
t4 = t1;
s4 = s3;
h4 = XSteam('h_ps',p4,s4);
% Evaluating net work output and back work ratio
Net_work = (h3-h4)-(h2-h1);
back_work_ratio = (h2-h1)/(h3-h4);
% Printing output
fprintf('Rankine Cycle Simulatorn')
fprintf('1-2 is isentropic Compression in Pumpn')
fprintf('2-3 is isobaric heat addition in Boilern')
fprintf('3-4 is isentropic Expansion in Turbinen')
fprintf('4-1 is isobaric heat rejection in Condensornn')
fprintf('Pressure at turbine Inlet(in bar):%dn',p3)
fprintf('Temperature at turbine Inlet(in Degree Celcius):%dn',t3)
fprintf('Pressure at the Condensor(in bar):%dnn',p4)
fprintf('RESULTSn')
fprintf('At State point 1:n')
fprintf('P1 is :%d %sn',p1,'Bar')
fprintf('T1 is :%d %sn',t1,'Deg Celcius')
fprintf('H1 is :%d %sn',h1,'KJ/Kg')
fprintf('S1 is :%d %snn',s1,'KJ/KgK')
fprintf('At State point 2:n')
fprintf('P2 is :%d %sn',p2,'Bar')
fprintf('T2 is :%d %sn',t2,'Deg Celcius')
fprintf('H2 is :%d %sn',h2,'KJ/Kg')
fprintf('S2 is :%d %snn',s2,'KJ/KgK')
fprintf('At State point 3:n')
fprintf('P3 is :%d %sn',p3,'Bar')
fprintf('T3 is :%d %sn',t3,'Deg Celcius')
fprintf('H3 is :%d %sn',h3,'KJ/Kg')
fprintf('S3 is :%d %snn',s3,'KJ/KgK')
fprintf('At State point 4:n')
fprintf('P4 is :%d %sn',p4,'Bar')
fprintf('T4 is :%d %sn',t4,'Deg Celcius')
fprintf('H4 is :%d %sn',h4,'KJ/Kg')
fprintf('S4 is :%d %snn',s4,'KJ/KgK')
fprintf('Net Work done is :%d %sn',Net_work,'KJ/Kg')
fprintf('Back Work ratio is :%dn',back_work_ratio)
% Plotting graphs :
% Finding lower limit of pressure corresponding to temperature
tf=p1;
tt = XSteam('Tsat_p',tf);
while ~isnan(tt)
tf = tf-0.00001;
tt = XSteam('Tsat_p',tf);
end
p_minimum = tf+0.00001;
% Finding Upper limit of pressure corresponding to temperature
tf1=p3;
tt1 = XSteam('Tsat_p',tf1);
while ~isnan(tt1)
tf1 = tf1+1;
tt1 = XSteam('Tsat_p',tf1);
end
p_maximum = tf1-1;
p = linspace(p_minimum,p_maximum,1000);
for i = 1:length(p)
T_saturation(i) = XSteam('Tsat_p',p(i));
S_sat_liquid(i) = XSteam('sL_p',p(i));
S_sat_vapour(i) = XSteam('sV_p',p(i));
end
figure(1)
hold on
plot1 = plot(S_sat_liquid,T_saturation,'r','Linewidth',1.2);
plot2 = plot(S_sat_vapour,T_saturation,'r','Linewidth',1.2);
% state 1 to state2
plot3 = plot([s1,s2],[t1, t2],'Color','k','Linewidth',1.2);
% state 2 to state 2'
t_range = linspace(t2,t2_dash,100);
for i =1:length(t_range)
s_range(i) = XSteam('s_pT',p2,t_range(i));
end
plot(s_range,t_range,'Color','k','Linewidth',1.2);
% state 2' to 3
if t3==tsat3
plot([s2_dash s3],[t2_dash t3],'Color','k','Linewidth',1.2);
elseif t3>tsat3
plot([s2_dash s3_dash],[t2_dash t3_dash],'Color','k','Linewidth',1.2);
%fr0m state 3' to 3
t_range1 = linspace(t3_dash,t3,10000);
for i =1:length(t_range1)
s_range1(i) = XSteam('s_pT',p3,t_range1(i));
end
plot(s_range1,t_range1,'Color','k','Linewidth',1.2);
end
% state 3 to state 4 to state 1
s = [s3 s4 s1];
t = [t3 t4 t1];
plot(s,t,'Color','k','Linewidth',1.2);
legend([plot1 plot3],'Saturation Curve','Rankine Cycle','Location','NorthWest')
text(s1,t1,'1')
text(s2,t2,'2')
text(s2_dash,t2_dash,'2dash','HorizontalAlignment','Right')
if t3==tsat3
text(s3,t3,'3')
elseif t3>tsat3
text(s3_dash,t3_dash,'3dash')
text(s3,t3,'3')
end
text(s4,t4,'4')
xlabel('Entropy(KJ/KgK)')
ylabel('Temperature(K)')
title('T-s plot for Rankine Cycle', 'Color','b','Fontsize',12)
% Plotting h-s curve
p = linspace(p_minimum,p_maximum,1000);
for i = 1:length(p)
h_sat_liquid(i) = XSteam('hL_p',p(i));
h_sat_vapour(i) = XSteam('hV_p',p(i));
end
figure(2)
hold on
axis([0 10 0 3500])
plot1 = plot(S_sat_liquid,h_sat_liquid,'r','Linewidth',1.2);
plot2 = plot(S_sat_vapour,h_sat_vapour,'r','Linewidth',1.2);
h = [h1 h2 h3 h4 h1];
ss = [s1 s2 s3 s4 s1];
plot3 = plot(ss,h,'Color','k','Linewidth',1.2);
legend([plot1 plot3],'Saturation Curve','Rankine Cycle','Location','NorthWest')
text(s1,h1,'1')
text(s2,h2,'2')
text(s3,h3,'3')
text(s4,h4,'4')
xlabel('Entropy(KJ/KgK)')
ylabel('Enthalpy(KJ/Kg)')
title('h-s plot for Rankine Cycle', 'Color','b','Fontsize',12)Results :
References :
https://en.wikipedia.org/wiki/Rankine_cycle
https://energyeducation.ca/encyclopedia/Rankine_cycle
https://in.mathworks.com/matlabcentral/fileexchange/9817-x-steam-thermodynamic-properties-of-water-and-steam
https://www.youtube.com/watch?v=QFZN71MY71o
https://www.sciencedirect.com/topics/earth-and-planetary-sciences/rankine-cycle