CalebBell / thermo

Thermodynamics and Phase Equilibrium component of Chemical Engineering Design Library (ChEDL)
MIT License
594 stars 114 forks source link

Flash calculations: Differences in entropy #119

Closed AndreasMatthias closed 1 year ago

AndreasMatthias commented 1 year ago

In the following example I'm using two different flash calculation to get the entropy. The first one is a T,P-flash:

S1 = flasher.flash(T=T, P=P).S()

The second one is a T,V-flash:

S2 = flash.flash(T=T, V=V).S()

I do not understand why S1 and S2 are not always equal. There are always differences in the liquid phase right below the saturation temperature.

This is the output of the example code below showing the differences in entropy dS = S1 - S1 at certain temperatures T. Just as the differences in enthaply dH.

Tc:    305.3 K
Tsat:  220.4 K
Tm:     90.3 K
--------------------------------------
     T       dS       dH
--------------------------------------
   100      0.0      0.0   L  L
   110      0.0      0.0   L  L
   120      0.0      0.0   L  L
   130      0.0      0.0   L  L
   140      0.0     -0.0   L  L
   150    -29.3  -6760.7   L  L
   160    -27.3  -6444.0   L  L
   170    -25.4  -6127.9   L  L
   180    -23.5  -5809.9   L  L
   190    -21.8  -5487.5   L  L
   200    -20.1  -5158.0   L  L
   210    -18.4  -4819.5   L  L
   220    -16.8  -4471.0   L  L
   230      0.0     -0.0   V  V
   240      0.0      0.0   V  V
   250      0.0      0.0   V  V
   260      0.0      0.0   V  V
   270      0.0      0.0   V  V
   280      0.0     -0.0   V  V
   290      0.0      0.0   V  V

This is the example code:

import thermo

constants, properties = thermo.ChemicalConstantsPackage.from_IDs(['ethane'])
CAS = constants.CASs[0]

eos_kwargs = dict(
    Tcs = constants.Tcs,
    Pcs = constants.Pcs,
    omegas = constants.omegas)

liquid = thermo.CEOSLiquid(
    thermo.PRMIX,
    eos_kwargs=eos_kwargs,
    HeatCapacityGases = [thermo.HeatCapacityLiquid(CAS)])

gas = thermo.CEOSGas(
    thermo.PRMIX,
    eos_kwargs=eos_kwargs,
    HeatCapacityGases = [thermo.HeatCapacityGas(CAS)])

flasher = thermo.FlashPureVLS(
    constants, properties,
    gas=gas, liquids=[liquid], solids=[])

P = 5e5
Tc = constants.Tcs[0]
Tsat = flasher.flash(P=P, VF=0).T
Tm = constants.Tms[0]

print(f'Tc:   {Tc:6.1f} K')
print(f'Tsat: {Tsat:6.1f} K')
print(f'Tm:   {Tm:6.1f} K')

print('--------------------------------------')
print('     T       dS       dH')
print('--------------------------------------')
for T in range(100, 300, 10):
    res1 = flasher.flash(T=T, P=P)
    S1 = res1.S()
    H1 = res1.H()
    V = res1.V()

    res2 = flasher.flash(T=T, V=V)
    S2 = res2.S()
    H2 = res2.H()

    dS = S1 - S2
    dH = H1 - H2

    print(f'{T:6.0f} {dS:8.1f} {dH:8.1f}   {res1.phase}  {res2.phase}')
CalebBell commented 1 year ago

Hi Andreas, The error in your code is to provide the CEOSLiquid a HeatCapacityLiquid object instead of the HeatCapacityGas object. The EOS already predicts the difference between the idea-gas and liquid state. Providing the HeatCapacityLiquid object results in a thermodynamic inconsistency. Sincerely, Caleb

AndreasMatthias commented 1 year ago

Thanks Caleb! I see it's working now. But I'm afraid I do not understand why. Could you please provide some more insight into the difference between HeatCapacityLiquid and HeatCapacityGas? I chose HeatCapacityLiquid on purpose, because I thought it should match to the liquid state.

CalebBell commented 1 year ago

No, it shouldn't take the HeatCapacityLiquid; that's why the argument name is HeatCapacityGases. Cubic EOSs are real fluid models. They are implemented starting with an idea-gas reference state; and then to obtain the properties of the liquid phase, those are calculated using the ideal-gas heat capacity and the departure functions of the cubic EOS. The Peng-robinson equation of state predicts things like the heat of vaporization and vapor pressure directly, mathematically, as the result of the solution of a set of equations. A real fluid's heat capacity comes in two parts - the ideal-gas part, and the departure from ideal gas part. The fact the Peng-Robinson equation of state can have a liquid root doesn't change the fact the heat capacity comes from those two parts. The equation for departure from ideal gas part is the same whether a liquid or gas is present.

AndreasMatthias commented 1 year ago

Thank you for the clarification. Much appreciated.