CoolProp icon indicating copy to clipboard operation
CoolProp copied to clipboard

Peng-Robinson backend hangs at some P, T conditions

Open infinity77 opened this issue 5 years ago • 1 comments

I am running CoolProp 6.4.2dev with Python 2.7 (I know, it's old but it's a long story) on Windows 10 64 bit.

Running the following simple Python script:


from __future__ import print_function

import CoolProp as MCP
import CoolProp.CoolProp as CP

T0 = 520.0
compound = 'Methane&Nitrogen&CO2&Ethane&Propane&Isobutane&Butane&Isopentane&Pentane&Hexane&Heptane'
composition = [0.9092, 0.0271, 0.0018, 0.0386, 0.011, 0.0037, 0.0037, 0.00135, 0.00135, 0.0008, 0.0014]

print('CoolProp version:', MCP.__version__)
print()
# Change to SRK - "works"
backend = 'PR'
EOS = CP.AbstractState(backend, compound)
EOS.set_mole_fractions(composition)

for P0 in [999999.0, 10132500.0]:
    print('Updating EoS          :', 'T(K) =', T0, ', P(Pa) =', P0)
    EOS.update(CP.PT_INPUTS, P0, T0)
    Mm   = EOS.molar_mass()
    print('Molar mass    (kg/mol):', Mm)
    Dmol = EOS.rhomolar()
    print('Molar density (Kg/m3) :', Dmol)
    print()

It will print the following lines:

CoolProp version: 6.4.2dev

Updating EoS : T(K) = 520.0 , P(Pa) = 999999.0 Molar mass (kg/mol): 0.017904342546 Molar density (Kg/m3) : 231.683091854

Updating EoS : T(K) = 520.0 , P(Pa) = 10132500.0

And it will be stuck forever in trying to update the PR EoS with the second set of temperature and pressure conditions. Using the SRK backend works.

If you let it run for long enough - meaning several minutes on a powerful workstation, you get this message:

ValueError: cubic has three roots, but phase not imposed and guess density not provided

Of course I can try and always force the phase as gas, but I have to deal with thousands of different mixtures with different T,P ranges specified by the user, so it’s not that obvious that at 150 K and 123 bar the mixture XYZ is actually gas...

I am a bit unclear on why it takes about 15 minutes for the solver to decide that the cubic has 3 roots. It seems to me that it might be going through endless loops trying to identify the phase, although the difference in pressure between “it works instantly” and “it needs 15 minutes to error out” are very, very small.

infinity77 avatar Mar 06 '21 11:03 infinity77

I've confirmed this on my side. Problem has to do with calculation of critical points, which is needed when three density roots are found for given P,T, pair. Recommend to specify the phase, although it doesn't solve the problem for mixtures.

ibell avatar Mar 08 '21 00:03 ibell