Project 2 - Rankine cycle Simulator
OBJECTIVE:
To create a Rankine Cycle Simulator using OCTAVE/MATLAB, which calculates the state points of the cycle based on user inputs and plots the corresponding T-s and h-s plots.
SOLUTION:
The Rankine cycle is an idealized thermodynamic cycle, which is used in steam turbines. The principle of the cycle is the conversion of heat energy into mechanical work, with the working fluid undergoing phase changes. The cycle consists of four different processes; Isentropic expansion in Turbine, Constant pressure heat rejection in condenser, Isentropic compression in Compressor and Constant pressure heat addition in heater.
- Isentropic compression: The working fluid is compressed or pumped from low to high pressure. This requires an input work to be given to the compressor or pump.
- Constant pressure heat addition: The high pressure fluid is heated, to a high temperature at constant pressure, in a heater. The working fluid undergoes phase change from liquid to vapour.
- Isentropic expansion: The high temperature fluid is expanded in a turbine, which extracts output work and in the process the fluid releases heat and pressure.
- Constant pressure heat rejection: The fluid releases heat and condense to liquid form in the condenser.
The working fluid used in the cycle is usually water and common heat sources include nuclear fission or the combustion of fossil fuels such as coal, natural gas and oil or concentrated solar power.
The following program takes input from the user and calculates the different parameters at various state points. It evaluates the performance of the cycle and displays various results, T-s and h-s plots. A sample of the calculation performed by the program is shown:
Sample Calculation
The following calculations has been done by taking the processes as below:
1-2 : Isentropic expansion in Turbine; Work output takes place
2-3 : Constant pressure heat rejection in condenser
3-4 : Isentropic compression in Compressor; Work input takes place
4-1 : Constant pressure heat addition in heater
Assuming, p1 = 30 bar, p2 = 0.05 bar, t1 = 400 deg C
p1 = p4 = 30 bar
p2 = p3 = 0.05 bar
For given p1 and t1, s1 = 6.9234 kJ/kg K, h1 = 3231.7 kJ/kg
s2 = s1
s2 = sf + x (sg – sf) = 0.47792 + x (7.9134)
x2 = 0.814
h2 = hf + x (hg – hf) = 138.27 + 0.814 x 2422.7 = 2110.34 kJ/kg
h3 = hf @ 0.05 bar = 138.27 kJ/kg
s3 = sf @ 0.05 bar = 0.47792 kJ/kg K
s4 = s3
At state point 3, nf = 0.00100537 m3/kg
h4 = nf (p4 – p3) + h3 = 0.00100537 x (3000 – 5) + 138.27 = 141.28 kJ/kg
Total work, Wnet = (h1 – h2) – (h4 – h3) = 1121.36 – 3.01 = 1118.35 kJ/kg
For a rate of steam flow, m = 1 kg/s
Wnet = 1118.35 kJ/s = 1118.35 kW
WT = h1 – h2 = 1121.36 kJ/kg
Back work ratio, WB = Wnet / WT = 1118.35/1121.36 = 0.997
Wc = h4 – h3 = 3.011 kJ/kg
Specific steam consumption, SSC = 3600/WT = 3600/1121.36 = 3.212 kg/kWh
Thermal efficiency = Wnet/Qin = 1118.35/(h1 – h4) = 0.3617 = 36.17 %
The program is given below:
clear all
close all
clc
%Intro Statements
fprintf("\n RANKINE CYCLE SIMULATOR \n");
fprintf("\n 1-2 : Isentropic expansion in Turbine; Work output takes place");
fprintf("\n 2-3 : Constant pressure heat rejection in condenser");
fprintf("\n 3-4 : Isentropic compression in Compressor; Work input takes place");
fprintf("\n 4-1 : Constant pressure heat addition in heater");
%Request for user input
t1 = input("\n\n INPUT TURBINE INLET TEMPERATURE IN DEG CELSIUS = ");
p1 = input("\n INPUT TURBINE INLET PRESSURE IN BARS = ");
p2 = input("\n INPUT TURBINE OUTLET PRESSURE IN BARS = ");
%known conditions
p4 = p1;
p3 = p2;
%calculation of state points
s1 = XSteam('s_pT',p1,t1);
h1 = XSteam('h_pT',p1,t1);
s2 =s1;
x2 = (s2 - XSteam('sL_p',p2))/(XSteam('sV_p',p2)- XSteam('sL_p',p2));
h2 = XSteam('hL_p',p2)+ x2*(XSteam('hV_p',p2) - XSteam('hL_p',p2));
h3 = XSteam('hL_p',p2);
s3 = XSteam('sL_p',p2);
s4 = s3;
h4 = (XSteam('vL_p',p2)*(p4 - p3)*100) + h3;
t3 = XSteam('Tsat_p',p3);
%calculation of Gamma
cp4 = XSteam('Cp_ps',p4,s4);
cv4 = XSteam('Cv_ps',p4,s4);
gamma4 = cp4/cv4;
const4 = (gamma4-1)/gamma4;
t4 = t3/((p3/p4)^const4);
t2 = t3;
%Calculation of ouputs
Wnet = (h1 - h2) - (h4 - h3);
Wt = h1 - h2;
Wc = h4 - h3;
eeta = (Wnet/(h1-h4))*100;
Wb = Wnet/Wt;
SSC = 3600/Wt;
%Display the results
fprintf("\n RESULTS \n");
fprintf("\n At state point 1 \n\n p1 = %d bar \n t1 = %d deg C \n h1 = %d kJ/kg \n s1 = %d kJ/kg K \n",p1,t1,h1,s1);
fprintf("\n At state point 2 \n\n p2 = %d bar \n t2 = %d deg C \n h2 = %d kJ/kg \n s2 = %d kJ/kg K \n",p2,t2,h2,s2);
fprintf("\n At state point 3 \n\n p3 = %d bar \n t3 = %d deg C \n h3 = %d kJ/kg \n s3 = %d kJ/kg K \n",p3,t3,h3,s3);
fprintf("\n At state point 4 \n\n p4 = %d bar \n t4 = %d deg C \n h4 = %d kJ/kg \n s4 = %d kJ/kg K \n",p4,t4,h4,s4);
fprintf("\n Turbine Work = %d kJ/kg",Wt);
fprintf("\n Compressor Work = %d kJ/kg",Wc);
fprintf("\n Net Work Output = %d kJ/kg",Wnet);
fprintf("\n Back Work Ratio = %d",Wb);
fprintf("\n Specific Steam Consumption = %d kg/kWh",SSC);
fprintf("\n Thermal Efficiency = %d Percent\n",eeta);
%Plotting T-s and h-s diagrams
temperature = linspace(0,1000,1000);
temp = linspace(XSteam('Tsat_p',p4), t1, 500);
temp1 = linspace(s3,s2,500);
for i = 1: length(temperature)
a(i) = XSteam('sL_T',temperature(i));
b(i) = XSteam('sV_T',temperature(i));
c(i) = XSteam('hL_T',temperature(i));
d(i) = XSteam('hV_T',temperature(i));
end
for j = 1: length(temp)
e(j) = XSteam('s_pT',p4,temp(j));
end
for k = 1: length(temp1)
f(k) = XSteam('h_ps',p4,temp1(k));
g(k) = XSteam('h_ps',p3,temp1(k));
end
t_SL = XSteam('Tsat_p',p4);
s_SL = XSteam('sL_p',p4);
t_SV = XSteam('Tsat_p',p4);
s_SV = XSteam('sV_p',p4);
%T-s diagram
figure(1)
hold on
plot(a,temperature,'color','b');
plot(b,temperature,'color','b');
plot([s3 s4],[t3 t4],'color','r');
plot([s2 s3],[t2 t3],'color','r');
plot([s1 s2],[t1 t2],'color','r');
plot(e,temp,'color','r');
plot([s_SL s_SV],[t_SL t_SV], 'color', 'r');
plot([s4 s_SL],[t4 t_SL],'color','r');
text(s1,t1,'1');
text(s2,t2,'2');
text(s3,t3,'3');
text(s4,t4,'4');
axis([0 10 0 500]);
xlabel('Entropy (kJ/kg K)');
ylabel('Temperature (deg C)');
title('T-s Diagram');
%h-s diagram
figure (2)
hold on
plot(a,c,'color','b');
plot(b,d,'color','b');
plot([s3 s4],[h3 h4],'color','r');
plot([s1 s2],[h1 h2],'color','r');
plot(temp1,f,'color','r');
plot(temp1,g,'color','r');
text(s1,h1,'1');
text(s2,h2,'2');
text(s3,h3,'3');
text(s4,h4,'4');
xlabel('Entropy (kJ/kg K)');
ylabel('Enthalpy (kJ/kg)');
title('h-s Diagram');
Output window:
T-s diagram
h-s diagram: