Project 2 - Rankine cycle Simulator

Rankine cycle:

The Rankine cycle or Rankine Vapor Cycle is the process widely used by power plants such as coal-fired power plants or nuclear reactors. In this mechanism, a fuel is used to produce heat within a boiler, converting water into steam which then expands through a turbine producing work.

The four processes of an ideal Rankine cycle are as follows;

• Process 1-2: Isentropic compression in a pump. Here the working fluid is pumped from low to high pressure.
• 2-3: Constant pressure heat addition in the boiler. In this process, the high-pressure liquid is added to the boiler and is heated at constant pressure to become a dry saturated vapour.
• 3-4: Isentropic expansion in a turbine. The dry saturated vapour expands through the turbine resulting in the generation of power.
• 4-1: Constant pressure heat rejection in a condenser. Wet vapour is condensed at a constant pressure to become a saturated liquid

  clear all
  close all
  clc
  %Program to simulate rankine cycle
  %Input the desired pressure 
  P1=input('Enter the pressure at turbine inlet (in bar): ')
  T1=input('Enter the temperature at turbine inlet (in degree celsius): ')
  P2=input('Enter the pressure at the condensor inlet (in bar): ')
  %% Printing the properties of first state
  S1=XSteam('S_Pt',P1,T1);
  H1=XSteam('H_pt',P1,T1);
  disp('At state 1:')
  fprintf('P1=%f barn',P1);
  fprintf('T1=%f degree celsiusn',T1);
  fprintf('S1=%f KJ/KgKn',S1);
  fprintf('H1=%f KJ/Kgn',H1);
  %% Finding unknowns at state 2, T2,S2,H2
  %Let Adiabatic exponent be 1.4
  gamma=1.4;
  S2=S1;
  T2=XSteam('T_ps',P2,S2);
  
  disp('At state 2:')
  if S2>XSteam('SV_P',P2)
      %State 2 is in superheated region 
      Tsat2=XSteam('Tsat_p',P2);
      Ssat2=XSteam('Sv_p',P2);
     
      H2=XSteam('H_PT',P2,T2);
  else 
      Tsat2=T2;
      Ssat2=S2;
      H2=XSteam('h_ps',P2,S2);
   
  end
      fprintf('P2=%f barn',P2);
      fprintf('T2=%f deg celn',T2);
      fprintf('S2=%f KJ/KgKn',S2);
      fprintf('H2=%f KJ/Kgn',H2);
      %% FINDING UNKNOWNS AT STATE 3,P3,T3,S3,H3
      T3=XSteam('Tsat_p',P2)
      P3=P2;
      S3=XSteam('Sl_T',T3);
      H3=XSteam('HL_T',T3);
      disp('At state 3:')
      fprintf('P3=%f barn',P3);
      fprintf('T3=%f deg celn',T3);
      fprintf('S3=%f KJ/KgKn',S3);
      fprintf('H3=%f KJ/Kgn',H3);
      %% FINDING UNKNOWNS AT STATE 4 P4,T4,S4,H4
      P4=P1;
      S4=S3;
      T4=XSteam('T_PS',P4,S4);
      H4=XSteam('H_PS',P4,S4);
      disp('At state 4:')
      fprintf('P4=%f barn',P4);
      fprintf('T4=%f deg celn',T4);
      fprintf('S4=%f KJ/KgKn',S4);
      fprintf('H4=%f KJ/Kgn',H4);
  
  
      %% PLOTTING T VS S PLOT
      % LETS PLOT THE DOME FIRST
      T=linspace(0,500,1000);
      %Defining critical temperature
      for i=1:length(T)
          if ~isnan(XSteam('Sl_T',T(i)))
          Sat1(i)=XSteam('Sl_T',T(i));
          end
      end
      for i=1:length(T)
          if ~isnan(XSteam('Sl_T',T(i)))     
          Sat2(i)=XSteam('SV_T',T(i));
          end
     end
      figure(1)
      subplot(2,1,1)
      plot(Sat1,T(1:748),Sat2,T(1:748))
      hold on
  
    %For plotting constant pressure curve from state 4 to state 1
      temp1_4=linspace(T4,T1,250);
      for i=1:250
          s1_4(i)=XSteam('s_pt',P4,temp1_4(i));
      end
       plot([S1 S2 Ssat2 S3 S4 s1_4],[T1 T2 Tsat2 T3 T4 temp1_4])
       xlabel('Entropy')
       ylabel('Temperature')
      
       %% PLOTTING H VS S PLOT
      
     % LETS PLOT THE DOME FIRST
      for i=1:length(T)
          if ~isnan(XSteam('Hl_T',T(i)))
          H_1(i)=XSteam('Hl_T',T(i));
          end
      end
      for i=1:length(T)
          if ~isnan(XSteam('Hl_T',T(i)))     
          H_2(i)=XSteam('HV_T',T(i));
          end
     end
      figure(1)
      subplot(2,1,2)
      plot(Sat1,H_1,Sat2,H_2)
      hold on
      %For plotting constant pressure curve from state 4 to state 1
      temp1_4=linspace(T4,T1,250);
      for i=1:250
          H1_4(i)=XSteam('H_pt',P4,temp1_4(i));
      end
       plot([S1 S2 S3 S4 s1_4],[H1 H2 H3 H4 H1_4])
       xlabel('ENTROPY')
       ylabel('ENTHALPY')
       
       %% EVALUATING NET OUTPUT AND WORK RATIO
  disp(sprintf('nThe net work output = %f KJ/Kg',(H1-H2)-(H4-H3)))
  disp(sprintf('The back work ratio=%f',(H4-H3)/(H1-H2)))
  
  
  
  ​