search

usnistgov / REFPROP-issues

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

Inconsistencies with temperature calculation and phase calculation in mixtures #687

Closed daveknippers closed 1 month ago

daveknippers commented 1 month ago

Description

I'm attempting to use REFPROP's API to predict the phase of mixtures. There's a discrepancy between what I'm seeing in the REFPROP GUI and the return values from the API.

The mixture is a simple hydrocarbon, 0.85 mole fraction Methane and 0.15 mole fraction Ethane. This mixture should be in a supercritical state. The REFPROP GUI states the phase as supercritical as expected, as shown in this attachment:

image For each row, the highlighted portions were the inputs.

This mixture should be a single phase at that temperature and pressure. Instead, the REFPROP API predicts it's a two phase mixture.

We use the REFPROP api to calculate the density based on the temperature and pressure (TP -> D), then we feed that density back into refprop's api and call it again to get the temperature (DP -> T)

(T  246.133831 , P  7414750.528 , no phase forced) D:  107.4275631434379
(D  107.4275631434379 , P  7414750.528 , no phase forced) T:  212.1822146332417

That doesn't jive with the GUI returns.

Using either density or the temperature as inputs along side pressure, the phase is predicted to be two-phase:

(T  246.133831 , P  7414750.528 , no phase forced) phase:  Two-phase
(D  107.4275631434379 , P  7414750.528 , no phase forced) phase:  Two-phase
(D  107.4275631434379 , P  7414750.528 , no phase forced) qmole:  0.7122786145413702 0.2877213854586298

When vapor is forced (or liquid, for that matter), the calculated phase and temperature agrees with the GUI (and Coolprop):

(D  107.4275631434379 , P  7414750.528 , forcing "V") T:  246.12279422157877
(T  246.133831 , P  7414750.528 , forcing "V") D:  107.40755988984448
(T  246.133831 , P  7414750.528 , forcing "V") phase:  Supercritical

This bug is similar to https://github.com/usnistgov/REFPROP-issues/issues/626 but it also shows an issue with calculation of temperature from density.

Steps to Reproduce

import os
os.environ['RPPREFIX'] = "fill in here"
from ctREFPROP.ctREFPROP import REFPROPFunctionLibrary
RP = REFPROPFunctionLibrary(os.environ['RPPREFIX'])
RP.SETPATHdll(os.environ['RPPREFIX'])
MASS_BASE_SI = RP.GETENUMdll(0, "MASS BASE SI").iEnum

species = 'Methane;Ethane'
mole_frac = [0.85,0.15]
mixture_temperature_K = 246.133831
mixture_pressure_Pa = 7414750.528
mixture_density_kg_m3 = 107.40852

print(species,mole_frac)

print('Mixture temperature (K): ',mixture_temperature_K)
print('Mixture pressure (Pa): ',mixture_pressure_Pa)

crit_temp_K = RP.REFPROPdll(species, '', 'TC', MASS_BASE_SI,0,0,0,0,mole_frac).Output[0]
print('Critical temperature (K): ',crit_temp_K)
crit_pres_Pa = RP.REFPROPdll(species, '', 'PC', MASS_BASE_SI,0,0,0,0,mole_frac).Output[0]
print('Critical pressure (Pa): ',crit_pres_Pa)

print('mixture_temperature_K > crit_temp_K ? ',mixture_temperature_K > crit_temp_K)
print('mixture_pressure_Pa > crit_pres_Pa ? ',mixture_pressure_Pa > crit_pres_Pa)

D = RP.REFPROPdll(species, 'TP', 'D', MASS_BASE_SI,0,0,mixture_temperature_K,mixture_pressure_Pa,mole_frac).Output[0]
print('(T ',mixture_temperature_K,', P ',mixture_pressure_Pa, ', no phase forced) D: ',D)

phase = RP.REFPROPdll(species, 'TP', 'PHASE', MASS_BASE_SI,0,0,mixture_temperature_K,mixture_pressure_Pa,mole_frac).hUnits
print('(T ',mixture_temperature_K,', P ',mixture_pressure_Pa, ', no phase forced) phase: ',phase)

T = RP.REFPROPdll(species, 'DP', 'T', MASS_BASE_SI,0,0,D,mixture_pressure_Pa,mole_frac).Output[0]
print('(D ',D,', P ',mixture_pressure_Pa, ', no phase forced) T: ',T)

phase = RP.REFPROPdll(species, 'DP', 'PHASE', MASS_BASE_SI,0,0,D,mixture_pressure_Pa,mole_frac).hUnits
print('(D ',D,', P ',mixture_pressure_Pa, ', no phase forced) phase: ',phase)

qmole_0,qmole_1 = RP.REFPROPdll(species, 'DP', 'QMOLE', MASS_BASE_SI,0,0,D,mixture_pressure_Pa,mole_frac).x[:2]
print('(D ',D,', P ',mixture_pressure_Pa, ', no phase forced) qmole: ',qmole_0,qmole_1)

print('-----------------------------------')

T = RP.REFPROPdll(species, 'DPV', 'T', MASS_BASE_SI,0,0,D,mixture_pressure_Pa,mole_frac).Output[0]
print('(D ',D,', P ',mixture_pressure_Pa, ', forcing "V") T: ',T)

D = RP.REFPROPdll(species, 'TPV', 'D', MASS_BASE_SI,0,0,mixture_temperature_K,mixture_pressure_Pa,mole_frac).Output[0]
print('(T ',mixture_temperature_K,', P ',mixture_pressure_Pa, ', forcing "V") D: ',D)

phase = RP.REFPROPdll(species, 'TPV', 'PHASE', MASS_BASE_SI,0,0,mixture_temperature_K,mixture_pressure_Pa,mole_frac).hUnits
print('(T ',mixture_temperature_K,', P ',mixture_pressure_Pa, ', forcing "V") phase: ',phase)

Expected behavior:

Methane;Ethane [0.85, 0.15]
Mixture temperature (K):  246.133831
Mixture pressure (Pa):  7414750.528
Critical temperature (K):  218.01563065614843
Critical pressure (Pa):  6223003.54090857
mixture_temperature_K > crit_temp_K ?  True
mixture_pressure_Pa > crit_pres_Pa ?  True
(T  246.133831 , P  7414750.528 , no phase forced) D:  107.4275631434379
(T  246.133831 , P  7414750.528 , no phase forced) phase:  Two-phase
(D  107.4275631434379 , P  7414750.528 , no phase forced) T:  212.1822146332417
(D  107.4275631434379 , P  7414750.528 , no phase forced) phase:  Two-phase
(D  107.4275631434379 , P  7414750.528 , no phase forced) qmole:  0.7122786145413702 0.2877213854586298
-----------------------------------
(D  107.4275631434379 , P  7414750.528 , forcing "V") T:  246.12279422157877
(T  246.133831 , P  7414750.528 , forcing "V") D:  107.40755988984448
(T  246.133831 , P  7414750.528 , forcing "V") phase:  Supercritical

Actual behavior:

Methane;Ethane [0.85, 0.15]
Mixture temperature (K):  246.133831
Mixture pressure (Pa):  7414750.528
Critical temperature (K):  218.01563065614843
Critical pressure (Pa):  6223003.54090857
mixture_temperature_K > crit_temp_K ?  True
mixture_pressure_Pa > crit_pres_Pa ?  True
(T  246.133831 , P  7414750.528 , no phase forced) D:  107.4275631434379
(T  246.133831 , P  7414750.528 , no phase forced) phase:  Supercritical
(D  107.4275631434379 , P  7414750.528 , no phase forced) T:  246.12279422157877
(D  107.4275631434379 , P  7414750.528 , no phase forced) phase:  Supercritical
(D  107.4275631434379 , P  7414750.528 , no phase forced) qmole: 0.85 0.15 (herr would be "State is single phase, quality not calculated")
-----------------------------------
(D  107.4275631434379 , P  7414750.528 , forcing "V") T:  246.12279422157877
(T  246.133831 , P  7414750.528 , forcing "V") D:  107.40755988984448
(T  246.133831 , P  7414750.528 , forcing "V") phase:  Supercritical
(D  107.40755988984448 , P  7414750.528 , forcing "V") herr:  State is single phase, quality not calculated

Versions

REFPROP Version: 10.0
Operating System and Version: Linux, Windows Access Method: python 3.11

ianhbell commented 1 month ago

Probably the saturation splines are on in the GUI and not enabled in your code.

daveknippers commented 1 month ago

Well, credit's due, adding a single flag to my call fixed my issue. Thank you!

renjie-github commented 1 month ago

Hi @daveknippers Could you share how you set the flag? I'm not familiar with this API usage and any suggestion will be much appreciated!

daveknippers commented 1 month ago

from what i remember setting the 6th argument to 1 turns splines on.

calc = rp.REFPROPdll(rp_species_string, k, tout, MASS_BASE_SI,0,1,v1,v2,rp_moles)

after splines are on, you can manually deactivate them with a FLAGSdll call for future calculations.

rp.FLAGSdll('Splines off',1)
calc = REFPROPdll(rp_species_string, k, tout, MASS_BASE_SI,0,0,v1,v2,rp_moles)

in certain cases my solutions failed to converge with splines on, so it was necessary to manage.

renjie-github commented 1 month ago

@daveknippers thanks a lot!