Project 2 - Rankine cycle Simulator

RANKINE CYCLE SIMULATOR IN MATLAB

AIM:

                To simulate the Rankine cycle in MATLAB and plot corresponding t-s and h-s diagram.

INTRODUCTION:

                Rankine cycle is the theoretical cycle on which the steam turbine works.

Process 1-2: Reversible adiabatic expansion in the turbine (or steam engine).

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

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

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

(i)For boiler

h4 + Q1 = h1

Q1 = h1 – h4

 

(ii)For turbine

h1 = WT + h2, (WT = turbine work)

WT = h1 – h2

 

(iii)For condenser

h2 = Q2 + h3

Q2 = h2 – h3

 

(iv) For the feed pump,

h3 + WP = h4, WP = Pump work

WP = h4– h3

 

Rankine cycle Efficiency:

                

                Net workdone=Work Extracted-Work input=Wt-Wp. Since we provide the work input to the pump, we must remove the same amount of work from the work output we extract from the system.

             

Back Work Ratio:

                Back work ratio is the amount of work consumed by the pump which is generated by the turbine. Lower the value of the back work ratio higher the performance of the cycle.

               

CODE:

%CODE TO SIMULATE RANKINE CYCLE
clear all
close all
clc

%Rankine cycle processes
fprintf('\t\t\t\tRANKINE CYCLE SIMULATOR\n\n')
fprintf('Process 1-2:Isentropic expansion in turbine\n');
fprintf('Process 2-3:Constant pressure heat rejection in condenser\n');
fprintf('Process 3-4:Isentropic compression in the pump\n');
fprintf('Process 4-1:Constant pressure heat addition in the boiler\n\n');

%Inputs
p1=input('Enter the pressure at turbine Inlet(in bar):');
t1=input('Enter the Temperature at turbine Inlet(in deg C):');
p2=input('Enter the pressure at Condenser Inlet(in bar):');
t3=XSteam('Tsat_p',p2);

fprintf('\n\n\t\t\t\tRESULTS');

%Values at state point 1
h1=XSteam('h_pT',p1,t1);
s1=XSteam('s_pT',p1,t1);
fprintf('\nAt state point 1\n');
fprintf('P1 is :%.2f bar\n',p1);
fprintf('T1 is :%.2f deg C\n',t1);
fprintf('H1 is :%.2f kJ/Kg\n',h1);
fprintf('S1 is :%.2f kJ/KgK\n\n',s1);

%Values at state point 2
s2=s1;
sf2=XSteam('sL_p',p2);
sg2=XSteam('sV_p',p2);
hf2=XSteam('hL_p',p2);
hg2=XSteam('hV_p',p2);
x2=(s2-sf2)/(sg2-sf2);
h2=hf2+(x2*(hg2-hf2));
t2=t3;
fprintf('At state point 2\n');
fprintf('P2 is :%.2f bar\n',p2);
fprintf('T2 is :%.2f deg C\n',t2);
fprintf('H2 is :%.2f kJ/Kg\n',h2);
fprintf('S2 is :%.2f kJ/KgK\n',s2);
fprintf('X2 is :%.2f \n\n',x2);

%Values at state point 3
p3=p2;
h3=XSteam('hL_p',p3);
s3=XSteam('sL_p',p3);
cp3=XSteam('CpL_p',p3);
cv3=XSteam('CvL_p',p3);
fprintf('At state point 3\n');
fprintf('P3 is :%.2f bar\n',p3);
fprintf('T3 is :%.2f deg C\n',t3);
fprintf('H3 is :%.2f kJ/Kg\n',h3);
fprintf('S3 is :%.2f kJ/KgK\n\n',s3);

%Values at state point 4
s4=s3;
p4=p1;
gamma=cp3/cv3;
pow=(1-(1/gamma));
h4=XSteam('h_ps',p4,s4);
t4=((p4/p3)^pow)*t3;
fprintf('At state point 4\n');
fprintf('P4 is :%.2f bar\n',p4);
fprintf('T4 is :%.2f deg C\n',t4);
fprintf('H4 is :%.2f kJ/Kg\n',h4);
fprintf('S4 is :%.2f kJ/KgK\n\n',s4);

%Calcualting the efficiency values
Wt=h1-h2;%Workdone by the turbine
Wp=h4-h3;%Workdone by the Pump
Wnet=Wt-Wp;
Qs=h1-h4;%Heat Supplied
Efficiency=(Wnet/Qs)*100;
SSC=(3600/Wnet);
BWR=Wp/Wt;
fprintf('Wt is:%.2f\n',Wt);
fprintf('Wp is:%.2f\n',Wp);
fprintf('Wnet is:%.2f\n',Wnet);
fprintf('Efficiency is:%.2f\n',Efficiency);
fprintf('S.S.C is:%.2f\n',SSC);
fprintf('Back Work Ratio is:%.2f\n',BWR);

%Points at saturation curve at 35Bar;
sL_saturated=XSteam('sL_p',p4);
tL_saturated=XSteam('T_ps',p4,sL_saturated);
hL_saturated=XSteam('hL_p',p4);

sV_saturated=XSteam('sV_p',p4);
tV_saturated=XSteam('T_ps',p4,sV_saturated);
hV_saturated=XSteam('hV_p',p4);

%Plotting T-s diagram
%Plotting the saturation curve
t_L=linspace(0,374.15,500);
for i=1:500
s_L(i)=XSteam('sL_T',t_L(i));
end
figure()
plot(s_L,t_L,'color','r');
title('t-s Diagram of Rankine cycle');
xlabel('Entropy(s)-[kJ/KgK]');
ylabel('Temperature(t)-[K]');
grid on;
hold on
t_V=linspace(374.15,0,500);
for i=1:500
s_V(i)=XSteam('sV_T',t_V(i));
end
plot(s_V,t_V,'color','r');
hold on;

%Plotting the Process curves
plot(s1,t1,'o','markers',8,'markerfacecolor','r');
text(s1+0.5,t1,'1');
hold on;
plot(s2,t2,'o','markers',8,'markerfacecolor','r');
text(s2+0.5,t2,'2');
hold on;
plot(s3,t3,'o','markers',8,'markerfacecolor','r');
text(s3,t3-10,'3');
hold on;
plot(s4,t4,'o','markers',8,'markerfacecolor','r');
text(s4-0.2,t4+20,'4');
hold on;
plot(sL_saturated,tL_saturated,'o','markers',8,'markerfacecolor','r');
hold on;
plot(sV_saturated,tV_saturated,'o','markers',8,'markerfacecolor','r');
line([s1 s2],[t1 t2],'color','k','linewidth',3);
line([s2 s3],[t2 t3],'color','k','linewidth',3);
line([s3 s4],[t3 t4],'color','k','linewidth',3);
line([s4 sL_saturated],[t4 tL_saturated],'color','k','linewidth',3);
line([sL_saturated sV_saturated],[tL_saturated tV_saturated],'color','k','linewidth',3);
line([sV_saturated s1],[tV_saturated t1],'color','k','linewidth',3);


%Plotting h-s diagram
for i=1:500
hL_saturated(i)=XSteam('hL_T',t_L(i));
end
figure()
plot(s_L,hL_saturated,'color','r');
title('h-s Diagram of Rankine cycle');
xlabel('Entropy(s)-[kJ/KgK]');
ylabel('Enthalpy(h)-[kJ/Kg]');
grid on;
hold on
for i=1:500
hV_saturated(i)=XSteam('hV_T',t_V(i));
end
plot(s_V,hV_saturated,'color','r');
hold on;

%Plotting the process curves
plot(s1,h1,'o','markers',8,'markerfacecolor','r');
text(s1+0.1,h1,'1');
hold on;
plot(s2,h2,'o','markers',8,'markerfacecolor','r');
text(s2+0.1,h2,'2');
hold on;
plot(s3,h3,'o','markers',8,'markerfacecolor','r');
text(s3,h3-150,'3');
hold on;
plot(s4,h4,'o','markers',8,'markerfacecolor','r');
text(s4,h4+200,'4');
hold on;
line([s1 s2],[h1 h2],'color','k','linewidth',3);
line([s2 s3],[h2 h3],'color','k','linewidth',3);
line([s3 s4],[h3 h4],'color','k','linewidth',3);
line([s4 s1],[h4 h1],'color','k','linewidth',3);

CODE EXPLANATION:

Lines 2-4:

clear all

close all

clc

Lines will reset the command window.

Lines: 6-11:

%Rankine cycle processes

fprintf('\t\t\t\tRANKINE CYCLE SIMULATOR\n\n')

fprintf('Process 1-2:Isentropic expansion in turbine\n');

fprintf('Process 2-3:Constant pressure heat rejection in condenser\n');

fprintf('Process 3-4:Isentropic compression in the pump\n');

fprintf('Process 4-1:Constant pressure heat addition in the boiler\n\n');

fprintf command will displays whatever the text we had entered in a single quote.

Lines 13-17:

%Inputs

p1=input('Enter the pressure at turbine Inlet(in bar):');

t1=input('Enter the Temperature at turbine Inlet(in deg C):');

p2=input('Enter the pressure at Condenser Inlet(in bar):');

t3=XSteam('Tsat_p',p2);

Using input command in MATLAB we can able to read the input values from the command window p1,t1,p2 values will be provided as the inputs by the user. T3 value can be found by using the steam table here steam table is provided as a matlab code and we can call the values using the function XSteam.

Syntax is XSTeam(‘Variable’,input value)

Variable refer to the variable for which the value is required.

Input value is the value based on which we need to find the variable.

Line 19:

fprintf('\n\n\t\t\t\tRESULTS');

\\n-Cursor goes to next line in the command window. It is similar to pressing the enter button while typing.

Line 21-28:

%Values at state point 1

h1=XSteam('h_pT',p1,t1);

s1=XSteam('s_pT',p1,t1);

fprintf('\nAt state point 1\n');

fprintf('P1 is :%.2f bar\n',p1);

fprintf('T1 is :%.2f deg C\n',t1);

fprintf('H1 is :%.2f kJ/Kg\n',h1);

fprintf('S1 is :%.2f kJ/KgK\n\n',s1);

Determining the h and s at the state point 1 and print the conditions of the state point 1 using the fprintf command.

 

Lines 30-44:

%Values at state point 2

s2=s1;

sf2=XSteam('sL_p',p2);

sg2=XSteam('sV_p',p2);

hf2=XSteam('hL_p',p2);

hg2=XSteam('hV_p',p2);

x2=(s2-sf2)/(sg2-sf2);

h2=hf2+(x2*(hg2-hf2));

t2=t3;

fprintf('At state point 2\n');

fprintf('P2 is :%.2f bar\n',p2);

fprintf('T2 is :%.2f deg C\n',t2);

fprintf('H2 is :%.2f kJ/Kg\n',h2);

fprintf('S2 is :%.2f kJ/KgK\n',s2);

fprintf('X2 is :%.2f \n\n',x2);

At 2^nd state the entropy value will be same as the state 1. Practically entropy change will occur, but the amount of change will be low that it cannot affect the entire output. Entropy and enthalpy at state point 2 for liquid and the gaseous phase to find the dryness fraction. Dryness fraction is the ratio of the vapor and liquid in the fluid. Total enthalpy is the sum of the Enthalpy at the fluid state and the product of the dryness fraction and the enthalpy at the vaporisation state. Pressure and temperature will be constant during the condensation process.

Printing the values of the entropy, enthalpy, pressure, and temperature values at the state 2 using the fprintf command.

 

Lines 46-56:

%Values at state point 3

p3=p2;

h3=XSteam('hL_p',p3);

s3=XSteam('sL_p',p3);

cp3=XSteam('CpL_p',p3);

cv3=XSteam('CvL_p',p3);

fprintf('At state point 3\n');

fprintf('P3 is :%.2f bar\n',p3);

fprintf('T3 is :%.2f deg C\n',t3);

fprintf('H3 is :%.2f kJ/Kg\n',h3);

fprintf('S3 is :%.2f kJ/KgK\n\n',s3);

Since the temperature and the pressure remains the same in during the condensation process pressure at state point 3 will be same as pressure at state point 2.Saturated Enthalpy and Entropy values at state point 3 will calculated for the liquid phase at the pressure p3. We know that the process 3-4 which is Reversible adiabatic process to form the equation we need the gamma value. Since we know Gamma=Cp/Cv along the process curve. We are finding the values of the Cp and Cv at point 3 using the pressure at the state point 3. Printing the state variables in command window using the fprintf command.

 

Lines 58-69:

%Values at state point 4

s4=s3;

p4=p1;

gamma=cp3/cv3;

pow=(1-(1/gamma));

h4=XSteam('h_ps',p4,s4);

t4=((p4/p3)^pow)*t3;

fprintf('At state point 4\n');

fprintf('P4 is :%.2f bar\n',p4);

fprintf('T4 is :%.2f deg C\n',t4);

fprintf('H4 is :%.2f kJ/Kg\n',h4);

fprintf('S4 is :%.2f kJ/KgK\n\n',s4);

Since the heat addition is the isobaric process pressure at state point 4 will be equal to the pressure at state 1. As we know that process 3-4 is adiabatic process the pressure-temperature relationship between the points can be given by the equation=. Since we know all the variable except the temperature at the state 4 we’ll get the temperature at state 4 using the equation mentioned above. Enthalpy of the state point 4 can be calculated by using the steam table using the pressure and  the enthalpy values at the state 4.

 

Lines 71-84:

%Calcualting the efficiency values

Wt=h1-h2;%Workdone by the turbine

Wp=h4-h3;%Workdone by the Pump

Wnet=Wt-Wp;

Qs=h1-h4;%Heat Supplied

Efficiency=(Wnet/Qs)*100;

SSC=(3600/Wnet);

BWR=Wp/Wt;

fprintf('Wt is:%0.2f\n',Wt);

fprintf('Wp is:%0.2f\n',Wp);

fprintf('Wnet is:%0.2f\n',Wnet);

fprintf('Efficiency is:%0.2f\n',Efficiency);

fprintf('S.S.C is:%0.2f\n',SSC);

fprintf('Back Work Ratio is:%.2f\n',BWR);

 

Here we are calculating the performance parameters of the rankine cycle using the equations in the theory part and printing the values in the command window using the fprintf command.

 

 

 

Lines 86-93:

%Points at saturation curve at 35Bar;

sL_saturated=XSteam('sL_p',p4);

tL_saturated=XSteam('T_ps',p4,sL_saturated);

hL_saturated=XSteam('hL_p',p4);

 

sV_saturated=XSteam('sV_p',p4);

tV_saturated=XSteam('T_ps',p4,sV_saturated);

hV_saturated=XSteam('hV_p',p4);

In the boiler isobaric process will takes place. The temperature and the entropy of the fluid will increase steadily increases up to the saturation line. After reaching the saturation curve the water will completely changes into vapour state with the constant temperature it will reaches the vapour saturation line further heating the temperature and the entropy of the vapour increases and the resulting in the superheated steam. In the above lines we’re calculating the points on the Liquid saturation curve and the vapour saturation curve. Using the pressure values at state point 1 or 4.

 

Line 95-113:

%Plotting T-s diagram

%Plotting the saturation curve

t_L=linspace(0,374.15,500);

for i=1:500

s_L(i)=XSteam('sL_T',t_L(i));

end

figure()

title('t-s Diagram of Rankine cycle');

plot(s_L,t_L,'color','r');

xlabel('Entropy(s)-[kJ/KgK]');

ylabel('Temperature(t)-[K]');

grid on;

hold on

t_V=linspace(374.15,0,500);

for i=1:500

s_V(i)=XSteam('sV_T',t_V(i));

end

plot(s_V,t_V,'color','r');

hold on;

These lines will plot the saturation curve in the t-s diagram. We’re storing a 500 temperature values from 0 to 374.15 using the linspace command. Now we know the temperature Entropy values for the corresponding temperature values is calculated using the for loop. From the steam table Entropy in saturated liquid line for the temperature t will be extracted and the stores in a variable s_L. Now we have the temperature and entropy values along the saturation curve and they were plotted. Similarly vapour line also plotted in the graph.

 

Lines 115-136:

%Plotting the Process curves

plot(s1,t1,'o','markers',8,'markerfacecolor','r');

text(s1+0.5,t1,'1');

hold on;

plot(s2,t2,'o','markers',8,'markerfacecolor','r');

text(s2+0.5,t2,'2');

hold on;

plot(s3,t3,'o','markers',8,'markerfacecolor','r');

text(s3,t3-10,'3');

hold on;

plot(s4,t4,'o','markers',8,'markerfacecolor','r');

text(s4-0.2,t4+20,'4');

hold on;

plot(sL_saturated,tL_saturated,'o','markers',8,'markerfacecolor','r');

hold on;

plot(sV_saturated,tV_saturated,'o','markers',8,'markerfacecolor','r');

line([s1 s2],[t1 t2],'color','k','linewidth',3);

line([s2 s3],[t2 t3],'color','k','linewidth',3);

line([s3 s4],[t3 t4],'color','k','linewidth',3);

line([s4 sL_saturated],[t4 tL_saturated],'color','k','linewidth',3);

line([sL_saturated sV_saturated],[tL_saturated tV_saturated],'color','k','linewidth',3);

line([sV_saturated s1],[tV_saturated t1],'color','k','linewidth',3);

Marking the state points using the temperature and the entropy values of the corresponding points in the t-s diagram. Using the text option state point number is denoted at the corresponding points.Syntax for the text command is text(x,y,’text’) x,y is the coordinate points at where the text to be displayed and whatever we type in place of ‘text’ it will be displayed in the plot in the corresponding coordinates. Line command is used to join the plots.

 

Lines 139-154:

%Plotting h-s diagram

for i=1:500

hL_saturated(i)=XSteam('hL_T',t_L(i));

end

figure()

plot(s_L,hL_saturated,'color','r');

title('h-s Diagram of Rankine cycle');

xlabel('Entropy(s)-[kJ/KgK]');

ylabel('Enthalpy(h)-[kJ/Kg]');

grid on;

hold on

for i=1:500

hV_saturated(i)=XSteam('hV_T',t_V(i));

end

plot(s_V,hV_saturated,'color','r');

hold on;

 

These lines will plot the saturation line in the h-s diagram in a new figure window. Since we have already calculated the entropy values enthalpy values along the saturation line will be calculated using the temperature values along the saturation curve and the curved is plotted using the s and h values.

 

Lines 156-172:

%Plotting the process curves

plot(s1,h1,'o','markers',8,'markerfacecolor','r');

text(s1+0.1,h1,'1');

hold on;

plot(s2,h2,'o','markers',8,'markerfacecolor','r');

text(s2+0.1,h2,'2');

hold on;

plot(s3,h3,'o','markers',8,'markerfacecolor','r');

text(s3,h3-150,'3');

hold on;

plot(s4,h4,'o','markers',8,'markerfacecolor','r');

text(s4,h4+200,'4');

hold on;

line([s1 s2],[h1 h2],'color','k','linewidth',3);

line([s2 s3],[h2 h3],'color','k','linewidth',3);

line([s3 s4],[h3 h4],'color','k','linewidth',3);

line([s4 s1],[h4 h1],'color','k','linewidth',3);

Using these lines, we’re plotting the process curves in the h-s diagram.

 

OUTPUT:

Output plots were shown below:

Command Window outputs:

REFERENCES:

1.Engineering Thermodynamics by R. K. Rajput.

2.https://youtu.be/lhilSmE2Ee0

3.https://youtu.be/HoSsUucDFHI

4.https://www.learnthermo.com/T1-tutorial/ch09/lesson-B/pg07.php