import math
# Variables
P_1 = 30; #[bar]
P_2 = 0.04; #[bar]
#(1).Carnot cycle
#It has been reported in the book that at 30 bar pressure (saturated) :
H_liq_1 = 1008.42; #[kJ/kg]
H_vap_1 = 2804.2; #[kJ/kg]
S_liq_1 = 2.6457; #[kJ/kg-K]
S_vap_1 = 6.1869; #[kJ/kh-K]
#Therefore, H_1 = H_liq_1, H_2 = H_vap_1, S_1 = S_liq_1 and S_2 = S_vap_1
H_1 = H_liq_1;
H_2 = H_vap_1;
S_1 = S_liq_1;
S_2 = S_vap_1;
#At 0.04 bar pressure (saturated) :
H_liq_2 = 121.46; #[kJ/kg]
H_vap_2 = 2554.4; #[kJ/kg]
S_liq_2 = 0.4226; #[kJ/kg-K]
S_vap_2 = 8.4746; #[kJ/kh-K]
# Calculations and Results
#Dryness fraction at state 3 can be found the fact that S_3 = S_2
x_3 = (S_2 - S_liq_2)/(S_vap_2 - S_liq_2);
H_3 = H_liq_2*(1 - x_3) + x_3*H_vap_2; #[kJ/kg]
#Dryness fraction at state 4 can be found the fact that S_4 = S_1
x_4 = (S_1 - S_liq_2)/(S_vap_2 - S_liq_2);
H_4 = H_liq_2*(1 - x_4) + x_4*H_vap_2; #[kJ/kg]
#Work done by turbine W_tur = -delta_H = -(H_3 - H_2)
W_tur = H_2 - H_3; #[kJ/kg]
#Work supplied by boiler,
q_H = H_2 - H_1; #[kJ/kg]
#Work transfer in compressor is given by
W_com = -(H_1 - H_4); #[kJ/kg]
#Efficiency can now be calculated as
#n = (Net work done/Work supplied by boiler)
n_carnot = (W_tur + W_com)/q_H;
#Efficiency of the Carnot cycle can also be determined from the formula
# n = 1 - (T_L/T_H), Where T_L is saturated temperature at 0.04 bar and T_H is saturated temperature at 30 bar
print "1.Carnot cycle";
print "The work done by the turbine is %f kJ/kg"%(W_tur);
print "The heat transfer in the boiler is %f kJ/kg"%(q_H);
print "The cycle efficiency is %f"%(n_carnot);
#(2).Rankine cycle
#The enthalpies at state 2 and 3 remain as in the Carnot cycle
#Saturated liquid enthalpy at 0.04 bar is
H_4_prime = H_liq_2;
#Saturated liquid volume at 0.04 bar as reported in the book is
V_liq = 0.001004; #[m**(3)/kg]
#Work transfer in pump can be calculated as
W_pump = -V_liq*(P_1 - P_2)*100; #[kJ/kg]
#Work transfer around pump gives, W_pump = -delta_H = -(H_1_prime - H_4_prime);
H_1_prime = H_4_prime - W_pump; #[kJ/kg]
#Heat supplied to boiler is
q_H_prime = H_2 - H_1_prime; #[kJ/kg]
#Work done by turbine is
W_tur_prime = H_2 - H_3; #[kJ/kg]
#Efficiency can now be calculated as
#n = (Net work done/Heat input)
n_rankine = (W_tur_prime + W_pump)/q_H_prime; #
print "2.Rankine cycle";
print "The work done by the turbine is %f kJ/kg"%(W_tur_prime);
print "The heat transfer in the boiler is %f kJ/kg"%(q_H_prime);
print "The cycle efficiency is %f"%(n_rankine);
T_max = 700+273.15; #[K] - Maximum temperature.
P_boiler = 10*10**(6); #[Pa] - Constant pressure in the boiler
P_condenser = 10*10**(3); #[Pa] - Constant pressure in the condenser
#At state 2 i.e, at 700 C and 10 MPa,it has been reported in the book that from steam table
S_2 = 7.1687; #[kJ/kg-K] - Entropy
H_2 = 3870.5; #[kJ/kg] - Enthalpy
#At state 3 i.e, at 700 C and 10 KPa,
S_3 = S_2; #[kJ/kg-K]- Entropy
#For sturated steam at 10 kPa, it has been reported in the book that from steam table
S_liq = 0.6493; #[kJ/kg-K]- Entropy of saturated liquid
S_vap = 8.1502; #[kJ/kg-K] - Enthalpy of saturated liquid
#Therefore steam is saturated and its dryness factor can be calculated as
x = (S_2 - S_liq)/(S_vap - S_liq);
# Calculations and Results
#The enthalpy at state 3 is now calculated. For steam at 10 kPa,it has been reported in the book that from steam table
H_liq = 191.83; #[kJ/kg]
H_vap = 2584.7; #[kJ/kg]
#Therefore enthalpy at state 3 is
H_3 = H_liq*(1-x) + H_vap*x; #[kJ/kg]
#Work done by the turbine
W_tur = -(H_3 - H_2); #[kJ/kg]
#Now we have to calculate work input to the pump
#State 4:Saturated liquid at 10 kPa
#State 4:Compressed liquid at 10 MPa
#Since volume of liquid does not get affected by pressure we take volume of saturated liquid at 10 kPa,
V_liq = 0.001010; #[m**(3)/kg]
#Work transfer in the pump is
W_pump = -V_liq*(P_boiler - P_condenser)*10**(-3); #[kJ/kg]
#Energy balance around pump gives, W_pump = -delta_H = -(H_1 - H_4)
H_4 = H_liq; # Enthalpy at state 4 (saturated liquid at 10 kPa)
H_1 = H_4 - W_pump; #[kJ/kg]
#Heat supplied to boiler is
q_H = H_2 - H_1; #[kJ/kg]
#Efficiency can now be calculated as
#n = (Net work done/Heat input)
n_rankine = (W_tur + W_pump)/q_H;
print "The efficiency of the Rankine cycle is found to be %f"%(n_rankine);
#Now let us determine the efficiency of Carnot cycle. The maximun temperature is 700 C and minimum temperature is that of saturated steam at 10 kPa,
T_min = 45.81 + 273.15; #[K] - From steam table as reported in the book
n_carnot = 1-(T_min/T_max);
#Note that the efficiency of Rankine cycle is less than that of carnot cycle
W = 1.1; #[kW] - Work done per ton of refrigeration
#1 ton refrigeration = 3.517 kW, therefore
H = 3.517; #[kW] - Heat absorbed
T_low = -30 + 273.15; #[K] - Low temperature maintained
# Calculations
#COP can be calculated as
#COP = (Heat absorbed/Work done)
COP = H/W;
#For reversed carnot cycle, COP = T_low/(T_high - T_low). Solving this we get
T_high = (T_low/COP) + T_low; #[K] - Higher temperature
#Heat rejected is
H_rej = W + H; #[kW];
# Results
print "The COP is %f"%(COP);
print "The higher temperature of the cycle is %f K"%(T_high);
print "The heat rejected per ton of refrigeration is %f kW"%(H_rej);
T_high = 20 + 273.15; #[K] - High temperature
T_low = 0 + 273.15; #[K] - Low temperature
Q_H = 10; #[kW] - Heat supplied
# Calculations
#If 'Q_H' is the rate at which heat is taken from surrounding and 'W' is the rate at which work is done,then
# Q_H = W + Q_L
#(Q_H/Q_L) = (T_high/T_low)
#Also for a reversible cycle, (Q_H/Q_L) = 1 + (W/Q_L). Solving we get,
Q_L = (T_low/T_high)*Q_H; #[kW]
W = (Q_H - Q_L) ; #[kW]
# Results
print "The minimum power required is %f kW"%(W);
T_high = 40 + 273.15; #[K] - High temperature
T_low = -20 + 273.15; #[K] - Low temperature
C = 10; #[tons of refrigeration] - Capacity
#1 ton refrigeration = 3.517 kW, therefore
H = C*3.517; #[kW] - Heat absorbed
# Calculations
#For reversed carnot cycle, COP = T_low/(T_high - T_low)
COP = T_low/(T_high - T_low);
# COP = (Refrigerating effect)/(Work input), therefore power required is given by
P = (H/COP); #[kW]
# Results
print "The COP is %f"%(COP);
print "The power required is %f kW"%(P);
COP = 4; #Coefficient of performance
P = 10; #[kW] - Work done on the cycle
# Calculations
#For reversed carnot cycle, COP = T_low/(T_high - T_low)
#ratio = (T_high/T_low),therefore
ratio = -1/(COP + 1);
# Refrigerating effect = (COP)*Work input, therefore refrigeration is given by
H = COP*P; #[kW]
#Maximum refrigearation in tons is given by
H_max = (H/3.517);
# Results
print "The maximum refrigeration value is %f ton"%(H_max);
m = 0.6; #[kg/s] - mass flow rate
T_low = -20+273.15; #[K] - Temperature at which vapour enters the compressor
T_high = 30+273.15; #[K] - Temperature at which vapour leaves the condenser
#From saturated refrigeration-12 tables we get,at -20 C
H_1 = 178.74; #[kJ/kg] - (H_1 = H_vap)
P_1 = 0.15093; #[MPa] - (P_1 = P_sat)
P_4 = P_1;
S_1 = 0.7087; #[kJ/kg-K] - (S_1 = S_vap)
S_2 = S_1;
#At 30 C
P_2 = 0.7449; #[MPa] - (P_2 = P_sat)
P_3 = P_2;
H_3 = 64.59; #[kJ/kg] - (H_3 = H_liq)
H_4 = H_3;
S_3 = 0.24; #[kJ/kg-K] - (S_3 = S_liq)
# Calculations and Results
#It is assumed that presssure drop in the evaporator and condenser are negligible. The heat transfer rate in the evaporator is
Q_L = m*(H_1 - H_4);
print "The heat transfer rate in the evaporator is %f kW"%(Q_L);
#At state 2 (P = 0.7449 MPa and S = 0.7087 kJ/kg-K) and looking in the superheated tables we have to calculate the enthalpy at state 2
#At P = 0.7 MPa and S = 0.6917 kJ/kg-K,
H_11 = 200.46; #[kJ/kg]
#At P = 0.7 MPa and S = 0.7153 kJ/kg-K,
H_12 = 207.73; #[kJ/kg]
#Thus at P = 0.7 MPa and S = 0.7087 kJ/kg-K, enthalpy is given by
H_13 = ((S_2 -0.6917)/(0.7153 - 0.6917))*(H_12 - H_11) + H_11; #[kJ/kg]
#At P = 0.8 MPa and S = 0.7021 kJ/kg-K,
H_21 = 206.07; #[kJ/kg]
#At P = 0.8 MPa and S = 0.7253 kJ/kg-K,
H_22 = 213.45; #[kJ/kg]
#Thus at P = 0.8 MPa and S = 0.7087 kJ/kg-K, enthalpy is given by
H_23 = ((S_2 -0.7021)/(0.7253 - 0.7021))*(H_22 - H_21) + H_21; #[kJ/kg]
#At P = 0.7449 MPa, S = 0.7087 kJ/kg-K, the enthalpy is
H_2 = ((0.7449 - 0.7)/(0.8 - 0.7))*(H_23 - H_13) + H_13; #[kJ/kg]
#Power consumed by the compressor is
W_comp = m*(H_2 - H_1); #[kW]
print "The power consumed by the compressor is %f kW"%(W_comp);
#Heat removed in evaporator/work done on compressor
COP_R = Q_L/W_comp;
print "The COP the refrigerator is %f kW"%(COP_R);
#At -20 C,saturated conditions
H_liq = 17.82; #[kJ/kg]
H_vap = 178.74; #[kJ/kg]
x_4 = (H_4 - H_liq)/(H_vap - H_liq);
print "The dryness factor of refrigerant after the expansion valve is %f"%(x_4);
#The heat transfer rate in the condenser is
Q_H = m*(H_3 - H_2); #[kW]
print "The heat transfer rate in the condenser is %f kW"%(Q_H);
#If the cycle would have worked as a pump then,
#COP_HP = (Heat supplied from condenser/Work done on compressor)
COP_HP = (-Q_H)/W_comp;
print "The COP if cycle would work as a heat pump is %f kW"%(COP_HP);
#If the cycle would have been a reversed Carnot cycle then
COP_C = T_low/(T_high - T_low);
print "The COP if cycle would run as reversed Carnot cycle is %f kW"%(COP_C);
from scipy.optimize import fsolve
import math
# Variables
H_1 = 310.38; #[kJ/kg]
H_2 = 277.7; #[kJ/kg]
H_5 = -122.6; #[kJ/kg]
H_6 = 77.8; #[kJ/kg]
#The enthalpy at point 3 is same at point 4 as the expansion is isenthalpic
# Calculations and Results
#The mass condensed is 1 kg and therefore m_1 = m+6 + 1
#Enthalpy balance around heat exchanger
#m_2*H_2 + m_2*H_6 = m_3*H_3 + m_7*H_7
#Enthalpy balance around separator
#m_4*H_4 = m_5*H_5 + m_6*H_6
#It can be seen that m_1 = m_2 = m_3 = m_4
#and m_6 = m_7 = m_1 - 1
#Substituting the values for enthalpy balance around heat exchanger we get,
#m_1*H_2 + (m_1 - 1)*(H_6) = m_1*H_3 + (m_1 - 1)*H_1
#and substituting the values for enthalpy balance around seperator we get
#m_1*H_3 = (1)*(-122.6) + (m_1 - 1)*77.8
#H_3 = ((1)*(-122.6) + (m_1 - 1)*77.8)/m_1
#Substituting the expression for 'H_3' in the above equation and then solving for m_1, we get
def f(m_1):
return m_1*H_2+(m_1-1)*(H_6)-m_1*(((1)*(-122.6) + (m_1 - 1)*77.8)/m_1)-(m_1-1)*H_1
m_1 = fsolve(f,4)
#Thus to liquify 1 kg of air compression of m_1 kg of air is carried out.
#Now substituting this value of m_1 to get the value of H_3,
H_3 = ((1)*(-122.6) + (m_1 - 1)*77.8)/m_1; #[kJ/kg]
#From given compressed air table we see at 200 bar and 160 K,
H_3_1 = 40.2; #[kJ/kg]
#At 200 bar and 180 K,
H_3_2 = 79.8; #[kJ/kg]
#By interpolation we get,
T_3 = ((H_3 - H_3_1)*(180 - 160))/(H_3_2 - H_3_1) + 160; #[K]
print "Temperature before throttling is %f K"%(T_3);
# Variables
#At 1 bar, 310 K
H_1 = 310.38; #[kJ/kg]
#At 200 bar, 310 K
H_2 = 277.7; #[kJ/kg]
#At 1 bar, Saturated liquid
H_7 = -122.6; #[kJ/kg]
#At 1 bar, Saturated vapour
H_8 = 77.8; #[kJ/kg]
#At 200 bar, 200 K
H_3 = 117.6; #[kJ/kg]
#At 1 bar, 100 K
H_11 = 98.3; #[kJ/kg]
# Calculations and Results
#For 1 kg of liquid air obtained,the overall enthalpy balance is
#m_2*H_2 = W - 122.6 + (m_2 - 1)*H_1
#W = - 0.8*m_2*(H_11 - H_3)
#Overall enthalpy balance equation becomes
#H_2*m_2 = 15.44*m_2 - H_7 + (m_2 - 1)*H_1, solving
m_2_prime = (H_7 - H_1)/(H_2 - 15.44 - H_1);
print "The number of kimath.lograms of air compressed per kg of liquid air produced is %f kg"%(m_2_prime);
#(2)
#Enthalpy balance around separator is
#0.2*m_2*H_5 = -H_7 + (0.2*m_2 - 1)*H_8, solving
m_2 = m_2_prime;
H_5_prime = ((0.2*m_2-1)*H_8 - H_7)/(0.2*m_2);
#At point 5, P = 200 bar and enthalpy is
H_5_1 = -33.53; #[kJ/kg]
#From compressed air tables at 200 bar and 140 K,
H_5_2 = 0.2; #[kJ/kg]
#At 200 bar and 120 K,
H_5_3 = -38.0; #[kJ/kg]
#Solving by interpolation we get
T_5 = ((H_5_1 - H_5_3)*(140 - 120))/(H_5_2 - H_5_3) + 120; #[K]
print "The temperature of air before throttling is %f K"%(T_5);
#(3)
#During mixing of streams 8 and 11 to produce stream 9, the enthalpy balance is
# (0.2*m_2 - 1)*H_8 + 0.8*m_2*H_11 = (m_2 - 1)*H_9,Solving for H_9
H_9_prime = ((0.2*m_2-1)*H_8+0.8*m_2*H_11)/(m_2 - 1);
#From given compressed air tables at 1 bar and 100 K,
H_9_1 = H_11;
#At 1 bar and 90 K
H_9_2 = 87.9; #[kJ/kg]
#Solving by interpolation we get
T_9 = ((H_9_prime - H_9_2)*(100 - 90))/(H_9_1 - H_9_2) + 90; #[K]
print "The temperature of stream entering second heat exchanger is %f K"%(T_9);
#(4)
#Enthalpy balance around first heat exchanger is
#H_2*m_2 + (m_2 - 1)*H_10 = H_3*m-2 + (m-2 - 1)*H_1, solving for H_10
H_10_prime = ((m_2 - 1)*H_1 + H_3*m_2 - H_2*m_2)/(m_2 - 1);
#From given compressed air tables at 1 bar and 140 K,
H_10_1 = 139.1; #[kJ/kg]
#At 1 bar and 120 K
H_10_2 = 118.8; #[kJ/kg]
#Solving by interpolation we get
T_10 = ((H_10_prime - H_10_2)*(140 - 120))/(H_10_1 - H_10_2) + 120; #[K]
print "The temperature of stream exiting second heat exchanger is %f K"%(T_10);
P_high = 40; #[bar]
P_low = 5; #[bar]
m_1 = 0.5; #[kg/s] - Rate of mass moving through the expander
m_2 = 0.1; #[kg/s] - Rate of mass of vapour mixing with air
e = 0.7; #Efficiency
#At state 3,(40 bar and 200 K),enthalpy and entropy is given by
H_3 = 179.7; #[kJ/kg]
S_3 = 5.330; #[kJ/kg-K]
#If isentropic conditions exits in the turbine then state 11 is at 5 bar
S_11 = 5.330; #[kJ/kg-K]
#From given compressed air tables at 5 bar and 120 K,
H_11_1 = 113.6; #[kJ/kg]
S_11_1 = 5.455; #[kJ/kg-K]
#At 5 bar and 100 K
H_11_2 = 90.6; #[kJ/kg]
S_11_2 = 5.246; #[kJ/kg-K]
#The enthalpy has to be determined when S = S_3
# Calculations and Results
#Solving by interpolation we get
H_11_s = ((H_11_1 - H_11_2)*(S_3 - S_11_2))/(S_11_1 - S_11_2) + H_11_2; #[kJ/kg]
#The adiabatic efficiency of tyrbine is given by
#(H_3 - H_11_a)/(H_3 - H_11_s) = e
H_11_a = H_3 - e*(H_3 - H_11_s); #[kJ/kg] - Actual enthalpy
#At 5 bar,the saturated enthalpy is given to be
H_8 = 88.7; #[kJ/kg]
#From enthalpy balance during mixing we get,
#0.1*H_8 + 0.5*H_11_a = 0.6*H_9
H_9 = (m_2*H_8 + m_1*H_11_a)/(m_1 + m_2); #[kJ/kg]
#From given compressed air tables at 5 bar and 140 K,
H_9_1 = 135.3; #[kJ/kg]
#At 5 bar and 120 K
H_9_2 = 113.6; #[kJ/kg]
#By interpolation we get
T_9 = ((H_9 - H_11_1)*(140 - 120))/(H_9_1 - H_11_1) + 120; #[K]
print " The temperature of air entering the second heat exchanger is %f K"%(T_9);