Better HeatCapacityLiquid Extrapolation for Supercritical Temperatures

[BenPortner]

Hello @CalebBell

I feel like the existing HeatCapacityLiquid.extrapolation methods are somewhat insufficient for most common gases when reaching supercritical temperatures. For example, trying to calculate the enthalpy of "liquid" nitrogen at T>126.2 K, the numbers quickly go through the ceiling:

from thermo import Chemical
import matplotlib.pyplot as plt
from numpy import arange

c = Chemical("N2")
c.HeatCapacityGas.method = "COOLPROP"
c.HeatCapacityLiquid.method = "COOLPROP"

def H_l(T):
    return c.H_ref - c.H_int_Tb_to_T_ref_g - c.Hvap_Tb + c.HeatCapacityLiquid.T_dependent_property_integral(c.Tb, T)

def H_g(T):
    return c.H_ref + c.HeatCapacityGas.T_dependent_property_integral(c.T_ref, T)

# calculate liquid and gas enthalpies in temperature range Tc +- 5K
Ts = arange(c.Tc-5,c.Tc+5,0.2)
Hgs = [H_g(T) for T in Ts]
Hls = [H_l(T) for T in Ts]

Blue: gas heat capacity, orange: liquid heat capacity, black: critical temperature

The "explosion" is somewhat dampened when choosing c.HeatCapacityLiquid.extrapolation = 'constant' but the enthalpy still grows steadily and surpasses the gas enthapy eventually:

Expected Behavior

No liquid phase can exist at temperatures above Tc. There are still situations where it may be useful to calculate the liquid enthalpy at these temperatures. However, the steadily growing enthalpy values exhibited by current extrapolation methods can easily lead to numerical instabilities (I can provide an example upon request). The stability could be greatly enhanced if there was an extrapolation method that could guarantee that H_l(T) = H_g(T) for T > Tc. The rationale for this is that the enthalpy difference between both phases disappears beyond the critical point. It contains the implicit assumption that HeatCapacityGas(T) is more accurate for T>Tc (which may not be the case if P>Pc?). The next best thing would be a method that at least makes sure that H_l(T) < H_g(T) for T>Tc. Is there any way to implement either one in thermo?

[CalebBell]

Hi Ben, The whole HeatCapacityLiquid isn't a good choice for a pure component above the critical point. I don't recommend you use it. I am happy with the current extrapolation algorithms, although it is not well suited to CoolProp. Sincerely, Caleb

[BenPortner]

Dear @CalebBell,

The whole HeatCapacityLiquid isn't a good choice for a pure component above the critical point. I don't recommend you use it.

Unfortunately, not using it is not an option sometimes. E.g. when you only know the final temperature after an iterative calculation (which btw. is not uncommon in engineering applications). Such iterations are easily tripped by the current HeatCapacityLiquid extrapolation implementations.

I am happy with the current extrapolation algorithms

H(T=126.4 K) = -24.2e3 J/mol
H(T=126.6 K) = 407.4e3 J/mol
H(T=126.8 K) = 1099.4e3 J/mol

You are happy with exponentially growing enthalpy values? I find that hard to believe.

although it is not well suited to CoolProp

The problem is not limited to CoolProp. It occurs for all methods, except DADGOSTAR_SHAW (which is rather inaccurate).

[CalebBell]

As noted above, gas heat capacity should be used above the supercritical point. Reading the documentation will show how to change the extrapolation methods used. The integration with CoolProp offers some features but others offer awkward behavior. Thermo can't offer answers that make everyone happy out of the box.