Rankine cycle (with units)

Calculate the efficiency of a Rankine vapor power cycle using a pure fluid model for water. Includes the units of quantities in the calculations.

Requires: Cantera >= 3.0.0, pint

Tags: Python thermodynamics thermodynamic cycle non-ideal fluid units

***************** State 1 ******************

  water:

       temperature   300 K
          pressure   3528.2 Pa
           density   996.59 kg/m^3
  mean mol. weight   18.016 kg/kmol
    vapor fraction   0
   phase of matter   liquid-gas-mix

                          1 kg             1 kmol
                     ---------------   ---------------
          enthalpy       -1.5858e+07        -2.857e+08  J
   internal energy       -1.5858e+07        -2.857e+08  J
           entropy            3913.2             70500  J/K
    Gibbs function       -1.7032e+07       -3.0685e+08  J
 heat capacity c_p            4181.3             75330  J/K
 heat capacity c_v              4131             74425  J/K


***************** State 2 ******************

  water:

       temperature   300.14 K
          pressure   8e+05 Pa
           density   996.91 kg/m^3
  mean mol. weight   18.016 kg/kmol
    vapor fraction   0
   phase of matter   liquid

                          1 kg             1 kmol
                     ---------------   ---------------
          enthalpy       -1.5857e+07       -2.8568e+08  J
   internal energy       -1.5858e+07       -2.8569e+08  J
           entropy              3915             70532  J/K
    Gibbs function       -1.7032e+07       -3.0685e+08  J
 heat capacity c_p            4178.6             75282  J/K
 heat capacity c_v            4127.9             74368  J/K


***************** State 3 ******************

  water:

       temperature   443.62 K
          pressure   8e+05 Pa
           density   4.1587 kg/m^3
  mean mol. weight   18.016 kg/kmol
    vapor fraction   1
   phase of matter   liquid-gas-mix

                          1 kg             1 kmol
                     ---------------   ---------------
          enthalpy       -1.3202e+07       -2.3784e+08  J
   internal energy       -1.3394e+07       -2.4131e+08  J
           entropy             10183        1.8346e+05  J/K
    Gibbs function       -1.7719e+07       -3.1922e+08  J
 heat capacity c_p            2464.4             44399  J/K
 heat capacity c_v            1764.5             31790  J/K


***************** State 4 ******************

  water:

       temperature   300 K
          pressure   3528.2 Pa
           density   0.030558 kg/m^3
  mean mol. weight   18.016 kg/kmol
    vapor fraction   0.83516
   phase of matter   liquid-gas-mix

                          1 kg             1 kmol
                     ---------------   ---------------
          enthalpy       -1.3822e+07       -2.4902e+08  J
   internal energy       -1.3938e+07        -2.511e+08  J
           entropy             10700        1.9277e+05  J/K
    Gibbs function       -1.7032e+07       -3.0685e+08  J
 heat capacity c_p   inf               inf              J/K
 heat capacity c_v   nan               nan              J/K

efficiency =  0.23320852221836424 dimensionless

import cantera.with_units as ctu

# parameters
eta_pump = 0.6 * ctu.units.dimensionless  # pump isentropic efficiency
eta_turbine = 0.8 * ctu.units.dimensionless # turbine isentropic efficiency
p_max = 116.03 * ctu.units.psi  # maximum pressure


def pump(fluid, p_final, eta):
    """Adiabatically pump a fluid to pressure p_final, using
    a pump with isentropic efficiency eta."""
    h0 = fluid.h
    s0 = fluid.s
    fluid.SP = s0, p_final
    h1s = fluid.h
    isentropic_work = h1s - h0
    actual_work = isentropic_work / eta
    h1 = h0 + actual_work
    fluid.HP = h1, p_final
    return actual_work


def expand(fluid, p_final, eta):
    """Adiabatically expand a fluid to pressure p_final, using
    a turbine with isentropic efficiency eta."""
    h0 = fluid.h
    s0 = fluid.s
    fluid.SP =s0, p_final
    h1s = fluid.h
    isentropic_work = h0 - h1s
    actual_work = isentropic_work * eta
    h1 = h0 - actual_work
    fluid.HP = h1, p_final
    return actual_work


def print_state(n, fluid):
    print('\n***************** State {0} ******************'.format(n))
    print(fluid.report())


if __name__ == '__main__':
    # create an object representing water
    w = ctu.Water()

    # start with saturated liquid water at 80.33 degrees Fahrenheit
    w.TQ = ctu.Q_(80.33, "degF"), 0.0 * ctu.units.dimensionless
    h1 = w.h
    p1 = w.P
    print_state(1, w)

    # pump it adiabatically to p_max
    pump_work = pump(w, p_max, eta_pump)
    h2 = w.h
    print_state(2, w)

    # heat it at constant pressure until it reaches the saturated vapor state
    # at this pressure
    w.PQ = p_max, 1.0 * ctu.units.dimensionless
    h3 = w.h
    heat_added = h3 - h2
    print_state(3, w)

    # expand back to p1
    turbine_work = expand(w, p1, eta_turbine)
    print_state(4, w)

    # efficiency
    eff = (turbine_work - pump_work)/heat_added

    print('efficiency = ', eff)