usnistgov / REFPROP-issues

A repository solely used for reporting issues with NIST REFPROP
26 stars 13 forks source link

REFPROP 10.0 - Erratic P/T Curve for LNG Mixture #146

Closed CallieAYYC closed 4 years ago

CallieAYYC commented 5 years ago

I have been trying to plot a P/T phase diagram for LNG and are getting some erratic results at the lower pressures.

My mol% mixture is as follows Nitrogen 2.32% Carbon Dioxide 0.01% Methane 97.66% Ethane 0.01%

I was expecting a smoother curve than the one provided and the dew point not crossing the bubble point. I tried various other EOS including default, GERG and PR. AGA8 did not come out at all.

REFPROP Version: [10.0]
Operating System and Version: [Win10]
Access Method: [Standard GUI]

Default: image

GERG: image

Peng-Robinson: image

ianhbell commented 5 years ago

Indeed, this seems to be a problem with the algorithm for calculation of the isopleth of the phase envelope for this mixture in REFPROP. I tried the same calculation with CoolProp and it worked better (though still had something funny at the end):

import CoolProp.CoolProp as CP
import matplotlib.pyplot as plt
AS = CP.AbstractState('HEOS','Nitrogen&CO2&Methane&Ethane')
z = [_/100 for _ in [2.32, 0.01, 97.66, 0.01]]
AS.set_mole_fractions(z)
AS.build_phase_envelope('')
PE = AS.get_phase_envelope_data()
plt.plot(PE.T, PE.p)
plt.xlabel('$T$ / K')
plt.ylabel('$p$ / Pa')
plt.yscale('log')
plt.show()

image

EricLemmon commented 4 years ago

I've tried several times over the past year to fix this problem (and several like it), but I could not come up with a solution that had the potential to damage calculations for other mixtures. This is a very tricky problem in that one algorithm is used for the calculation of the saturation states for any mixture, and thus changes can (and will) seriously impact the calculations elsewhere. However, this week I think I've found something that works, and while plotting the phase diagram for all 150 predefined mixtures, I did not notice any problems.

The fix was applied to the subroutine SATSPLN. That routine works by calculating the entire phase diagram starting from a low density in the vapor phase. Once that point is obtained, the density is slowly incremented until it fails in the liquid phase at high densities, generally at states that would be in the solid phase. The calculations at each state point are saved in arrays. These are then used later in spline fits to obtain very good initial estimates for use in the VLE solver in Refprop. For the situation here, the initial state point first calculated was actually wrong (although the criteria for VLE was met, but the root it found was in the two-phase and caused by the erratic behavior of the pure fluid equations between the saturation boundaries). After about 10 steps, the SATSPLN routine jumped to the right solution and then all proceeded as it should. The solution to this problem was to step backwards through the calculations for a composition in the equilibrium phase that radically changed. If found, the preceding points were used to estimate values at the state where the solver failed. The VLE routine was then called again but now with these initial estimates to obtain the correct point. This continued until all incorrect points were replaced.

If you would like an updated DLL with the new changes, please let me know (eric.lemmon@nist.gov). However, keep in mind that you will need to watch carefully when working with other mixtures to make sure they have not been affected. This is something that will require many beta testers to ensure that what I implemented is working properly.