Mixture Properties using Peng-Robinson

122 views
Skip to first unread message

Friederike Boehm

unread,
Feb 4, 2022, 11:03:23 AM2/4/22
to coolprop-users
Hi,
I am trying to switch a thermodynamic property generation routine from Aspen Plus to CoolProp. This uses the Peng-Robinson EoS and it does work from a technical standpoint. The problem is that when I compare the results to those from Aspen only the density is the same, everything else (enthalpy-differences, heat capacity, isentropic expansion coefficient) is totally off. For pure fluids it works totally fine. I looked through a lot of the CoolProp code, but I still can't figure out why I get weird values like an expansion coefficient of 300 when it should be something like 1.3.
The testing fluid is a 50/50 mixture of methane and ethane but this will become more complex.
Any ideas what's going on?

Also I noticed that in the solver_rho_Tp function in CubicBackend.cpp always sets the quality to -1 and does't allow the two-phase region.

Best regards,
Friederike

Ian Bell

unread,
Feb 4, 2022, 5:45:20 PM2/4/22
to coolpro...@googlegroups.com
Can you please provide a more comprehensive example? Hard to say otherwise. 

--
You received this message because you are subscribed to the Google Groups "coolprop-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to coolprop-user...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/coolprop-users/8767bf9d-f162-4042-be1f-1dc1a4f87470n%40googlegroups.com.

Friederike Boehm

unread,
Feb 7, 2022, 12:49:01 PM2/7/22
to coolprop-users
Sure thing. Apologies for the insufficient description ealier. I was trying to keep it short but after spending another day with this problem I realised this is way more complex :D
There seem to be two (possibly separate) problems:
  1. mixture properties using Peng-Robinson EoS
  2. mixture properties in the two-phase region
Regarding the first problem (1): I set all property methods in Aspen to the standard Peng-Robinson method which should be the same as the one implemented in CoolProp. For pure fluids (methane, ethane, propane and nitrogen) I get property values which are close enough to assume the methods are similar enough. As far as I know I turned off all non-standard alpha-functions or volume corrections. When I do the same calculations for a 50/50 mixture of methane and ethane the densities are practically the same, but everything else deviates way more. I calculated entropy, enthalpy, heat capacity, isentropic expansion coefficient and density for temperatures from 100 to 300 K and pressures from 1 to 21 bar in Aspen and CoolProp. The tables show the mean and maximum absolute deviation from the Aspen values in percent. In regard to the second problem, the data for the mixture is cleaned from all data points that are in the two-phase region according to the Aspen data, as those data points were even further off or could not be computed.

Methane:
property                                  mean deviation / %        max  deviation / %
____________________________________________________________________
totalEntropy                             0.1054                              0.328845
vaporEntropy                           0.0788673                       0.0865843
liquidEntropy                           0.184282                          0.328845
totalEnthalpy                           0.128219                         0.342046
vaporEnthalpy                         0.106786                         0.149505
liquidEnthalpy                         0.19194                           0.342046
totalHeatCapacity                  0.115254                         0.380907
vaporHeatCapacity                0.098544                         0.380907
liquidHeatCapacity                0.164933                         0.21431
totalExpansionCoefficient    0.155593                         0.301243
vaporExpansionCoefficient  0.148849                         0.161466
liquidExpansionCoefficient  0.240105                         0.653065
totalDensity                            0.100435                         0.356001
vaporDensity                          0.0203411                       0.171851
liquidDensity                          0.338551                          0.356001


Methan[0.5]&Ethane[0.5]
property                                 mean deviation / %        max  deviation / %
___________________________________________________________________
totalEntropy                             54.1622                         67.7097
vaporEntropy                           60.6791                         67.7097
liquidEntropy                           45.9937                         52.8483
totalEnthalpy                           32.8982                         49.2161
vaporEnthalpy                         26.5181                         31.7697
liquidEnthalpy                         40.8952                         49.2161
totalHeatCapacity                  62.8988                         80.908
vaporHeatCapacity                74.8701                         80.908
liquidHeatCapacity                47.8936                         49.3701
totalExpansionCoefficient    3192.22                         33342.3
vaporExpansionCoefficient  5673.56                         33342.3
liquidExpansionCoefficient  81.9964                         158.258
totalDensity                            0.127352                       0.297458
vaporDensity                          0.0105603                     0.0445371
liquidDensity                          0.273741                       0.297458

The deviations for the other pure fluids (ethane, propane and nitrogen) are bigger but nowhere near the mixture deviations.
The AbstractState is initialized like this:

methaneAndEthane = AbstractStateFactory["PR", "Methane&Ethane"];

AbstractStateSetFractions[methaneAndEthane, {0.5, 0.5}]

AbstractStateSetBinaryInteraction[methaneAndEthane, 0, 1, "kij", -0.0026]


The second problem (2) concerns the properties in the two-phase region. Using temperature and pressure as inputs for states with a low quality (<0.2-0.3) CoolProp will fail to calculate the "total" value, which means that no phase is imposed and it should be in the two phase region. For bigger quality values, not imposing a phase will return the value for the gas phase. From what I can gather in the source code for CoolProp::AbstractCubicBackend::solver_rho_Tp "two-phase" cannot be imposed and Q is always set to -1 in the function.
Trying to use the inputs that contain the quality (QT_INPUTS or PQ_INPUTS) also does not always work. QT will fail for all quality inputs that are neither 0 or 1 (Error: Inputs in Brent [50.000000,10000.000000] do not bracket the root.), PQ has worked so far and at least the resulting densities are the same as in Aspen.
One other weird thing I noticed: After building the phase envelope, updating the State with QT_INPUTS will result in an Assertion Failure, other inputs will still work. What exactly does building the phase envelope do? I haven't really figured that out yet...

I hope this makes the problem a bit clearer. If not I'll try again :D

Friederike Boehm

unread,
Feb 15, 2022, 12:15:15 PM2/15/22
to coolprop-users
Hi,

so I found a work around for the second problem. It seems the PT Inputs don't work correctly in the two phase region. Using PQ Inputs and root-finding works to update the state. The densities come out perfectly matching those in Aspen. There are a few problems that remain though.
  • parameter outputs for the saturated states (SatL and SatV) only work for very few parameters (density , temperature, pressure, residual enthalpy, residual entropy, residual Gibbs energy). Most importantly, they don't work for molar specific enthalpy and entropy and I cannot figure out why
  • this leads to neither enthalpy nor entropy being calculated with keyed_output in the two-phase region (they are derived from hmolar and smolar values for the SatV and SatL states)
  • The results are still very different from Aspen-PR values apart from the density
  • After building the phase envelope for mixtures QT_INPUTS and PQ_INPUTS will result in an Assertion Failure
I can access the values for SatV and SatL by constructing new states with the compositions from get_mole_fractions_liquid and get_mole_fractions_vapor and setting the temperature and molar density, but this seems a bit roundabout.

Any ideas?

And does this still belong here or should I open an Issue on GitHub?

Lynn McGuire

unread,
Feb 15, 2022, 5:38:26 PM2/15/22
to coolpro...@googlegroups.com
Hi Friederike,

If you want to try PR in another thermodynamic software, our Win32 programming interface for Design II for Windows handles VB, VBA, C++, F77, etc.
    https://www.winsim.com

We are not freeware but I do give out two week passwords for free.

Sincerely,
Lynn McGuire
President
WinSim Inc.
To view this discussion on the web visit https://groups.google.com/d/msgid/coolprop-users/69754ad9-1b2e-4e4a-b862-e7bffe917208n%40googlegroups.com.

-- 
Sincerely,
Michael Lynn McGuire, P.E.
President
WinSim Inc.
USA 281-545-9200 x102 office
USA 832-202-3918 cell
USA 281-545-8820 fax
l...@winsim.com
www.winsim.com

P.S. Please sign up for our email notification list at 
     https://secure.campaigner.com/CSB/Public/Form.aspx?fid=1502895

Friederike Boehm

unread,
Feb 16, 2022, 5:20:20 AM2/16/22
to coolprop-users
Hi Lynn,

thanks for the offer! We're not using windows machines exclusively in our group and the open-source nature of CoolProp makes things easier for us as well. So as long as I can figure out this problem with PR I am not really willing to switch (we are already using CoolProp in other projects in the group and are trying to keep/move everything to CoolProp).

Best regards,
Friederike

Friederike Boehm

unread,
Feb 16, 2022, 12:35:24 PM2/16/22
to coolprop-users

Finally got around to put my problem into python code (I am usually using Mathematica, but that is not the problem)
PengRobinson.ipynb
Reply all
Reply to author
Forward
0 new messages