Cubic Equations of State¶
Introduction¶
CoolProp (as of version 6) comes with two standard cubic equations of state: Soave-Redlich-Kwong (SRK) and Peng-Robinson (PR). These two equations of state can be expressed in a common, generalized form:
where for pure fluids, \(a\) and \(b\) are not composition dependent, whereas for mixtures, they have composition dependence. These cubic EOS can be converted to a form that is compatible with the multi-fluid model used in CoolProp according to the analysis in Bell and Jager [168].
The motivations for the use of cubic EOS are twofold:
The only required information for the EOS are \(T_c\), \(p_c\), and the acentric factor of the pure fluids
They are much more computationally efficient (see below)
Caveats¶
Warning
NOT ALL PROPERTIES ARE AVAILABLE AS INPUTS/OUTPUTS
Only a limited subset of properties are available currently. You can do:
Flash calculations with TP, PQ, DT, QT inputs
Calculation of mixture critical point(s)
Calculation of some mixture flashes
Pure Fluids¶
Usage¶
As a user, in the high-level interface, all that you have to do to evaluate properties using the SRK or PR backends is to change the backend name.
In [1]: import CoolProp.CoolProp as CP
# The multi-parameter Helmholtz backend
In [2]: CP.PropsSI("T","P",101325,"Q",0,"HEOS::Propane")
Out[2]: 231.03621464431768
# SRK
In [3]: CP.PropsSI("T","P",101325,"Q",0,"SRK::Propane")
Out[3]: 231.2704159248594
# Peng-Robinson
In [4]: CP.PropsSI("T","P",101325,"Q",0,"PR::Propane")
Out[4]: 230.95724200364114
The same holds for the low-level interface:
In [5]: import CoolProp.CoolProp as CP
In [6]: AS = CP.AbstractState("SRK", "Propane"); AS.update(CP.QT_INPUTS, 0, 300); print(AS.p())
1008752.0286813603
In [7]: AS = CP.AbstractState("PR", "Propane"); AS.update(CP.QT_INPUTS, 0, 300); print(AS.p())
997556.0816415057
All the fluids available in CoolProp are also available through the cubic equations of state. The fluids can be extended according to the analysis shown below
Speed¶
The increase in speed for evaluating properties from a cubic EOS is one of the primary motivations. Here we show an example of the speedup to VLE calculations:
In [8]: import CoolProp.CoolProp as CP
# The multi-parameter Helmholtz backend
In [9]: AS = CP.AbstractState("HEOS", "Propane")
In [10]: %timeit AS.update(CP.QT_INPUTS, 0, 300)
15.3 us +- 166 ns per loop (mean +- std. dev. of 7 runs, 100,000 loops each)
# The cubic SRK backend
In [11]: AS = CP.AbstractState("SRK", "Propane")
In [12]: %timeit AS.update(CP.QT_INPUTS, 0, 300)
2.21 us +- 12.5 ns per loop (mean +- std. dev. of 7 runs, 100,000 loops each)
And here, we run the PT flash
In [13]: import CoolProp.CoolProp as CP
# The multi-parameter Helmholtz backend
In [14]: AS = CP.AbstractState("HEOS", "Propane")
In [15]: AS.specify_phase(CP.iphase_gas)
In [16]: %timeit AS.update(CP.PT_INPUTS, 101325, 300)
4.53 us +- 11.3 ns per loop (mean +- std. dev. of 7 runs, 100,000 loops each)
# The cubic SRK backend
In [17]: AS = CP.AbstractState("SRK", "Propane")
In [18]: AS.specify_phase(CP.iphase_gas)
In [19]: %timeit AS.update(CP.PT_INPUTS, 101325, 300)
368 ns +- 0.986 ns per loop (mean +- std. dev. of 7 runs, 1,000,000 loops each)
As you can see, the speed difference is quite significant
Mixtures¶
Interaction Parameters¶
For mixtures, cubic EOS (and their modifications) are heavily used. Cubic EOS allow for reasonable prediction of VLE with only one adjustable parameter, at least for reasonable mixtures. CoolProp does not come with any values for the \(k_{ij}\) parameter, but it is straightforward to add it yourself.
Warning
The ability to adjust interaction parameters is ONLY available in the low-level interface
Warning
When you call set_binary_interaction_double, it only applies to the given instance of the AbstractState
In [20]: import CoolProp.CoolProp as CP
In [21]: AS = CP.AbstractState("SRK", "Methane&Ethane")
In [22]: AS.set_mole_fractions([0.5, 0.5])
In [23]: AS.update(CP.QT_INPUTS, 0, 120); print(AS.p())
104467.46237356304
In [24]: AS.set_binary_interaction_double(0,1,"kij",-0.05)
In [25]: AS.update(CP.QT_INPUTS, 0, 120); print(AS.p())
82034.10345545325
Critical Points¶
According to a forthcoming paper from Bell et al., it is possible to calculate critical points of mixtures from cubic EOS. For instance, here is how to calculate all the critical points (there can be more than one) that are found for an equimolar mixture of methane and ethane:
In [26]: import CoolProp.CoolProp as CP
In [27]: AS = CP.AbstractState("SRK", "Methane&Ethane")
In [28]: AS.set_mole_fractions([0.5, 0.5])
In [29]: pts = AS.all_critical_points()
In [30]: [(pt.T, pt.p, pt.rhomolar, pt.stable) for pt in pts]
Out[30]: [(266.3194710014874, 6859304.778393156, 7895.129640980038, True)]
Cubics in multi-fluid model¶
The cubic equations of state can also be used to replace a single fluid in a multi-fluid model (the GERG-like model) (see Bell and Jager [168]), by appending either -SRK or -PengRobinson to the fluid name. For instance, you could calculate the NBP with both the conventional multi-fluid model (as in GERG), or translating both to cubic equations of state
In [31]: import CoolProp.CoolProp as CP
# With the normal multi-fluid model
In [32]: CP.PropsSI('T', 'P', 101325, 'Q', 0, 'HEOS::Methane-SRK[0.4]&Ethane-SRK[0.6]')
Out[32]: 122.01168651100696
# With both fluids using cubic translations in the multi-fluid model
In [33]: CP.PropsSI('T', 'P', 101325, 'Q', 0, 'HEOS::Methane[0.4]&Ethane[0.6]')
Out[33]: 121.51910871460049
Detailed Example¶
Here we plot phase envelopes and critical points for an equimolar methane/ethane mixture
(Source code, png, .pdf)
Adding your own fluids¶
The cubic fluids in CoolProp are defined based on a JSON format, which could yield something like this for a FAKE (illustrative) fluid. For instance if we had the file fake_fluid.json (download it here: fake_fluid.json):
[
{
"CAS": "000-0-00",
"Tc": 400.0,
"Tc_units": "K",
"acentric": 0.1,
"aliases": [
],
"molemass": 0.04,
"molemass_units": "kg/mol",
"name": "FAKEFLUID",
"pc": 4000000.0,
"pc_units": "Pa"
}
]
The JSON-formatted fluid information is validated against the JSON schema, which can be obtained by calling the get_global_param_string function:
In [34]: import CoolProp.CoolProp as CP
# Just the first bit of the schema (it's a bit large; see below for the complete schema)
In [35]: CP.get_global_param_string("cubic_fluids_schema")[0:60]
Out[35]: '{\n "$schema": "http://json-schema.org/draft-04/schema#",\n '
A JSON schema defines the structure of the data, and places limits on the values, defines what values must be provided, etc.. There are a number of online tools that can be used to experiment with and validate JSON schemas: http://www.jsonschemavalidator.net , http://jsonschemalint.com/draft4/.
In [36]: import CoolProp.CoolProp as CP, json
# Normally you would read it from a file!
In [37]: fake_fluids = [
....: {
....: "CAS": "000-0-00",
....: "Tc": 400.0,
....: "Tc_units": "K",
....: "acentric": 0.1,
....: "aliases": [
....: ],
....: "molemass": 0.04,
....: "molemass_units": "kg/mol",
....: "name": "FAKEFLUID",
....: "pc": 4000000.0,
....: "pc_units": "Pa"
....: }
....: ]
....:
# Adds fake fluid for both SRK and PR backends since they need the same parameters from the same internal library
In [38]: CP.add_fluids_as_JSON("PR", json.dumps(fake_fluids))
In [39]: CP.PropsSI("T","P",101325,"Q",0,"SRK::FAKEFLUID")
Out[39]: 247.0950826605273
# Put in a placeholder interaction parameter (all beta and gamma values are 1.0)
In [40]: CP.apply_simple_mixing_rule("000-0-00", CP.get_fluid_param_string("Ethane","CAS"), "Lorentz-Berthelot")
# Once a fluid has been added to the cubic library, it can be used in the multi-fluid model
In [41]: CP.PropsSI("T","P",101325,"Q",0,"HEOS::FAKEFLUID-SRK[0.3]&Ethane-SRK[0.7]")
Out[41]: 190.93809486697978
Complete fluid schema¶
For completeness, here is the whole JSON schema used for the cubic backends:
In [42]: import CoolProp.CoolProp as CP
In [43]: print(CP.get_global_param_string("cubic_fluids_schema"))
{
"$schema": "http://json-schema.org/draft-04/schema#",
"title": "Fluid for cubic EOS in CoolProp",
"items": {
"properties": {
"name": {
"description": "Name of the fluid",
"type": "string"
},
"CAS": {
"description": "CAS registry number of the fluid",
"type": "string"
},
"Tc": {
"description": "Critical temperature (K)",
"type": "number",
"minimum": 0.1,
"maximum": 20000
},
"pc": {
"description": "Critical pressure (Pa)",
"type": "number",
"minimum": 0,
"maximum": 500000000
},
"rhomolarc": {
"description": "Critical density (mol/m^3)",
"type": "number",
"minimum": 0.1,
"maximum": 2000000
},
"rhomolarc_units": {
"description": "Units of the critical density provided",
"enum": [
"mol/m^3",
"kg/m^3"
]
},
"acentric": {
"description": "Acentric factor (-)",
"type": "number",
"minimum": -10,
"maximum": 10
},
"molemass": {
"description": "Molar mass (kg/mol)",
"type": "number",
"minimum": 0,
"maximum": 1
},
"molemass_units": {
"description": "Units of the molar mass provided",
"enum": [
"kg/mol"
]
},
"pc_units": {
"description": "Units of the critical pressure",
"enum": [
"Pa"
]
},
"Tc_units": {
"description": "Units of the critical temperature",
"enum": [
"K"
]
},
"aliases": {
"type": "array",
"items": {
"type": "string"
},
"minItems": 0
},
"alpha": {
"description": "The alpha function being used",
"properties": {
"type":{
"enum": [
"Twu",
"Mathias-Copeman",
"default"
]
},
"c":{
"type": "array",
"items": {
"type": "number"
},
"minItems": 3,
"maxItems": 3
}
}
}
},
"required": [
"name",
"CAS",
"Tc",
"Tc_units",
"pc",
"pc_units",
"acentric",
"molemass",
"molemass_units"
]
},
"type": "array"
}
