Refrigeration and liquefaction

Yeast works best at about 32 °C, and fermentation releases a great deal of heat. In winter, cooling water from Kingfisher's towers holds the fermenters at temperature easily, but on hot summer days the tower water arrives at 28 °C or more, too warm to remove heat fast enough. The plant needs about 3 MW of mechanical refrigeration on those days. Next door, Lakeshore Gas needs refrigeration on a larger scale and at lower temperature to condense carbon dioxide at −20 °C. Both use ammonia vapor-compression refrigeration, the subject of most of this lesson, which then finishes with the related problem of liquefying gases such as nitrogen. These are the last applications of the course, and they draw on almost everything before them: equations of state, residual properties, diagrams, flow processes and the second law.

The vapor-compression cycle

A refrigerator is a heat engine run backwards (Lesson 7). In the vapor-compression cycle, a working fluid, the refrigerant, circulates through four steps. Low-pressure liquid evaporates in the evaporator, absorbing heat from the space or stream being cooled. A compressor raises the vapor to a pressure at which it can condense at a temperature above the surroundings. The condenser rejects heat to cooling water or air, turning the vapor back to liquid. An expansion valve drops the liquid's pressure back to the evaporator pressure, flashing part of it to vapor and chilling the rest.

On a pressure-enthalpy diagram (Lesson 9), the cycle is almost a rectangle: two horizontal lines for evaporation and condensation at constant pressure, a vertical line for the expansion valve at constant enthalpy, and a sloped line for compression. The refrigeration effect per kilogram is the enthalpy rise in the evaporator, the compressor work is the enthalpy rise in the compressor, and their ratio is the coefficient of performance.

The expansion valve is the one deliberately irreversible step. Replacing it with an expander would recover some work and raise the COP, but for most refrigerants the work available from expanding a liquid is small and the valve is far cheaper. Large carbon dioxide and natural gas liquefaction plants, where the recoverable work is larger, do use liquid expanders.

Choosing a refrigerant

A good refrigerant has a vapor pressure above atmospheric at the evaporator temperature, so that air cannot leak in, and a moderate pressure at the condenser temperature. It has a large heat of vaporization, so that little fluid needs to circulate. It should be chemically stable, compatible with the equipment, and as safe and environmentally benign as possible. Ammonia scores very well on thermodynamic grounds and has no ozone depletion or direct global warming effect, so it dominates industrial refrigeration in food, beverage and chemical plants. It is toxic and somewhat flammable, which means ammonia plants are built to strict codes with leak detection, ventilation and emergency procedures. Commercial and domestic systems have moved through several generations of synthetic refrigerants as ozone and climate rules have tightened; carbon dioxide itself is increasingly used as a refrigerant in supermarkets and heat pumps.

How operating temperatures affect performance

Every kelvin matters in refrigeration. Lowering the evaporator temperature or raising the condenser temperature widens the temperature lift, reduces the COP and raises the compressor power for the same duty. As a rule of thumb, each degree of extra lift costs roughly 2 to 4 percent in power. Designers therefore keep the evaporator as warm as the process allows and the condenser as cool as the cooling medium allows, and they size heat exchangers generously to keep approach temperatures small. In winter, when the condenser can run cooler, the same plant uses noticeably less power.

# Vapor-compression refrigeration with ammonia and Linde liquefaction, using Peng-Robinson properties
import math
R = 8.314
def pr_state(T, P, Tc, Pc, w, phase="vap"):
    """Peng-Robinson: returns Z, residual enthalpy H^R (J/mol), residual entropy S^R (J/mol K), ln(fugacity coefficient)."""
    k = 0.37464 + 1.54226 * w - 0.26992 * w ** 2
    sa = 1 + k * (1 - math.sqrt(T / Tc))
    ac = 0.45724 * R ** 2 * Tc ** 2 / Pc
    a, b = ac * sa ** 2, 0.07780 * R * Tc / Pc
    dadT = -ac * k * sa / math.sqrt(T * Tc)
    A, B = a * P / (R * T) ** 2, b * P / (R * T)
    c2, c1, c0 = -(1 - B), A - 3 * B ** 2 - 2 * B, -(A * B - B ** 2 - B ** 3)
    p = c1 - c2 ** 2 / 3; q = 2 * c2 ** 3 / 27 - c2 * c1 / 3 + c0
    disc = (q / 2) ** 2 + (p / 3) ** 3
    if disc > 0:
        u, v = -q / 2 + math.sqrt(disc), -q / 2 - math.sqrt(disc)
        roots = [math.copysign(abs(u) ** (1 / 3), u) + math.copysign(abs(v) ** (1 / 3), v) - c2 / 3]
    else:
        r = math.sqrt(-(p / 3) ** 3); phi = math.acos(max(-1.0, min(1.0, -q / (2 * r)))); m = 2 * math.sqrt(-p / 3)
        roots = [m * math.cos((phi + 2 * math.pi * j) / 3) - c2 / 3 for j in range(3)]
    roots = sorted(z for z in roots if z > B)
    Z = roots[-1] if phase == "vap" else roots[0]
    L = math.log((Z + (1 + math.sqrt(2)) * B) / (Z + (1 - math.sqrt(2)) * B))
    HR = R * T * (Z - 1) + (T * dadT - a) / (2 * math.sqrt(2) * b) * L
    SR = R * math.log(Z - B) + dadT / (2 * math.sqrt(2) * b) * L
    lnphi = Z - 1 - math.log(Z - B) - A / (2 * math.sqrt(2) * B) * L
    return Z, HR, SR, lnphi

# fluid: Tc (K), Pc (Pa), omega, M (g/mol), ideal gas Cp/R constants (A, B, C, D) - approximate values
FL = {"NH3": (405.7, 112.8e5, 0.253, 17.03, (3.578, 3.020e-3, 0.0, -0.186e5)),
      "CO2": (304.2, 73.8e5, 0.224, 44.01, (5.457, 1.045e-3, 0.0, -1.157e5)),
      "N2":  (126.2, 34.0e5, 0.038, 28.01, (3.280, 0.593e-3, 0.0, 0.040e5))}
def H(f, T, P, ph="vap"):
    Tc, Pc, w, M, (A, B, C, D) = FL[f]; T0 = 298.15
    hig = R * (A * (T - T0) + B / 2 * (T ** 2 - T0 ** 2) + C / 3 * (T ** 3 - T0 ** 3) - D * (1 / T - 1 / T0))
    return (hig + pr_state(T, P, Tc, Pc, w, ph)[1]) / M       # kJ/kg
def S(f, T, P, ph="vap"):
    Tc, Pc, w, M, (A, B, C, D) = FL[f]; T0 = 298.15
    sig = R * (A * math.log(T / T0) + B * (T - T0) + C / 2 * (T ** 2 - T0 ** 2) - D / 2 * (1 / T ** 2 - 1 / T0 ** 2) - math.log(P / 1e5))
    return (sig + pr_state(T, P, Tc, Pc, w, ph)[2]) / M       # kJ/(kg K)
def psat(f, T):
    Tc, Pc, w, M, _ = FL[f]; Tr = T / Tc
    P = Pc * math.exp(5.92714 - 6.09648 / Tr - 1.28862 * math.log(Tr) + 0.169347 * Tr ** 6
                      + w * (15.2518 - 15.6875 / Tr - 13.4721 * math.log(Tr) + 0.43577 * Tr ** 6))
    for _ in range(200):
        d = pr_state(T, P, Tc, Pc, w, "liq")[3] - pr_state(T, P, Tc, Pc, w, "vap")[3]
        P *= math.exp(d)
        if abs(d) < 1e-10: break
    return P
def bisect(fun, target, lo, hi):
    for _ in range(80):
        mid = (lo + hi) / 2
        if fun(mid) < target: lo = mid
        else: hi = mid
    return mid

def ammonia_cycle(t_evap, t_cond, eta=0.80):
    Te, Tk = t_evap + 273.15, t_cond + 273.15
    Pe, Pk = psat("NH3", Te), psat("NH3", Tk)
    h1, s1 = H("NH3", Te, Pe), S("NH3", Te, Pe)                  # saturated vapor leaving the evaporator
    T2s = bisect(lambda T: S("NH3", T, Pk), s1, Tk, Tk + 200)
    h2 = h1 + (H("NH3", T2s, Pk) - h1) / eta                     # actual compressor outlet
    T2 = bisect(lambda T: H("NH3", T, Pk), h2, Tk, Tk + 250)
    h3 = H("NH3", Tk, Pk, "liq")                                 # saturated liquid leaving the condenser
    q_evap, w = h1 - h3, h2 - h1                                 # h4 = h3 across the expansion valve
    return Pe, Pk, h1, h2, h3, T2, q_evap, w

Pe, Pk, h1, h2, h3, T2, q, w = ammonia_cycle(-10, 35)
print(f"Ammonia, -10 C / 35 C: Pe = {Pe/1e5:.2f} bar, Pk = {Pk/1e5:.2f} bar, discharge {T2-273.15:.0f} C")
print(f"   refrigeration effect {q:.0f} kJ/kg, compressor work {w:.0f} kJ/kg, COP = {q/w:.2f} (Carnot 5.85)")
load = 3.0e3                                                     # kW of fermenter cooling in summer
m = load / q
print(f"   for {load/1000:.1f} MW of cooling: ammonia {m:.2f} kg/s, compressor power {m*w:.0f} kW, condenser duty {m*(h2-h3):.0f} kW")

# Condensing the carbon dioxide stream: 226.5 mol/s from 20 bar, 35 C vapor to 20 bar, -20 C liquid
m_co2 = 226.5 * 44.01 / 1000                                     # kg/s
q_co2 = m_co2 * (H("CO2", 308.15, 20e5) - H("CO2", 253.15, 20e5, "liq"))
Pe2, Pk2, a1, a2, a3, T22, q2, w2 = ammonia_cycle(-28, 35)
print(f"CO2 condenser duty {q_co2:.0f} kW; ammonia at -28 C / 35 C: COP {q2/w2:.2f}, compressor power {q_co2 * w2 / q2:.0f} kW")

# Linde liquefaction of nitrogen: feed 200 bar and 300 K, vapor leaves the exchanger at 1 bar and 295 K
Tsat = bisect(lambda T: psat("N2", T), 1.013e5, 60, 120)
h_hp, h_vap_out = H("N2", 300, 200e5), H("N2", 295, 1.013e5)
h_liq = H("N2", Tsat, 1.013e5, "liq")
y = (h_vap_out - h_hp) / (h_vap_out - h_liq)
print(f"Nitrogen normal boiling point by Peng-Robinson: {Tsat:.1f} K")
print(f"Linde process from 200 bar and 300 K: fraction liquefied y = {y:.3f}")

Ammonia, -10 C / 35 C: Pe = 2.88 bar, Pk = 13.41 bar, discharge 131 C
   refrigeration effect 1107 kJ/kg, compressor work 286 kJ/kg, COP = 3.87 (Carnot 5.85)
   for 3.0 MW of cooling: ammonia 2.71 kg/s, compressor power 775 kW, condenser duty 3775 kW
CO2 condenser duty 3380 kW; ammonia at -28 C / 35 C: COP 2.37, compressor power 1425 kW
Nitrogen normal boiling point by Peng-Robinson: 77.3 K
Linde process from 200 bar and 300 K: fraction liquefied y = 0.072

The program builds every state from the Peng-Robinson equation: vapor pressures by equal fugacities, enthalpies and entropies from the ideal gas state plus residuals, the isentropic compressor outlet by bisection on entropy, and the actual outlet from the efficiency. It is a compact version of what a property package does inside a process simulator.

Liquefying gases: the Linde process

Gases such as nitrogen, oxygen and natural gas have critical temperatures far below ambient, so no ordinary refrigerant can condense them. Liquefaction instead uses the Joule-Thomson effect of Lesson 8. In the Linde process, gas compressed to high pressure and cooled to ambient temperature passes through a counter-current heat exchanger and then through a throttle valve. Throttling cools it, and the cold, low-pressure gas that does not liquefy flows back through the heat exchanger, precooling the incoming high-pressure gas. Each pass gets colder, until part of the gas condenses after the valve. At steady state, an energy balance around the exchanger, valve and separator gives the fraction liquefied.

For liquefaction to work, the high-pressure feed must have a lower enthalpy than the low-pressure vapor at the same temperature, which is exactly the condition of a positive Joule-Thomson effect, a negative residual enthalpy at high pressure. Hydrogen and helium have negative Joule-Thomson coefficients at room temperature and must be precooled, with liquid nitrogen for hydrogen, before a Linde stage can liquefy them. The Claude process improves on Linde by sending part of the high-pressure gas through an expander, which cools it much more than a throttle valve and recovers work.

Pulling the course together

The three Kingfisher projects show how the parts of this course fit together. Lakeshore Gas's liquefaction unit needed phase diagrams (Lesson 2), equations of state for gas and liquid (Lessons 3 to 5), residual properties (Lesson 8), diagrams and tables (Lesson 9), compressor calculations (Lesson 10) and refrigeration (this lesson). The steam turbine study needed the second law (Lesson 7), steam tables and the Rankine cycle (Lessons 9 to 11). The fermenter chiller needed heat effects (Lesson 6) and refrigeration cycles. In every case, the first law gave the energy, the second law gave the limit, and property methods supplied the numbers. Thermodynamics II takes the same tools into mixtures, phase equilibrium and chemical reaction equilibrium.