Index
Clapeyron.COSMOSAC02
Clapeyron.COSMOSAC10
Clapeyron.COSMOSACdsp
Clapeyron.NRTL
Clapeyron.PSRKUNIFAC
Clapeyron.UNIFAC
Clapeyron.UNIFACFV
Clapeyron.UNIFACFVPoly
Clapeyron.UNIQUAC
Clapeyron.VTPRUNIFAC
Clapeyron.Wilson
Clapeyron.aspenNRTL
Clapeyron.ogUNIFAC
Clapeyron.tcPRWilsonRes
Activity Models
There are two alternatives on the definition of an activity model:
- Defining an excess gibbs energy function
- Defining an activity coefficient function
those two can be converted between one form to another via:
$\gamma_{i} = \frac{\partial{G^E}}{\partial{n_i}}$
$\ln{\gamma_{i}} = \frac{1}{RT}\frac{\partial{G^E}}{\partial{n_i}}$
When defining one form, the other is derived automatically.
Those functions can also be derived from any arbitrary equation of state:
$\frac{\partial{G^E}}{\partial{n_i}}= \mu_i - \mu^{0}_i$
Where $\mu_i$ and $\mu^{0}_i$ are the mixture and pure chemical potentials of component i
. in this case, those potentials are dependent of the pressure. whereas activity models are usually only temperature and composition dependent.
Clapeyron.Wilson
— TypeWilson <: ActivityModel
Wilson(components;
puremodel = PR,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
Tc
: Single Parameter (Float64
) - Critical Temperature[K]
Pc
: Single Parameter (Float64
) - Critical Pressure[Pa]
acentricfactor
: Single Parameter (Float64
) - Acentric FactorMw
: Single Parameter (Float64
) - Molecular Weight[g/mol]
g
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parameter
model parameters
Tc
: Single Parameter (Float64
) - Critical Temperature[K]
Pc
: Single Parameter (Float64
) - Critical Pressure[Pa]
ZRA
: Single Parameter (Float64
) - Rackett compresibility factorMw
: Single Parameter (Float64
) - Molecular Weight[g/mol]
g
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parameter
Input models
puremodel
: model to calculate pure pressure-dependent properties
Description
Wilson activity model, with Rackett correlation for liquid volume:
Gᴱ = nRT∑xᵢlog(∑xⱼjΛᵢⱼ)
Λᵢⱼ = exp(-gᵢⱼ/T)*Vⱼ/Vᵢ
ZRAᵢ = 0.29056 - 0.08775ωᵢ
Vᵢ = (RTcᵢ/Pcᵢ)ZRAᵢ^(1 + (1-T/Tcᵢ)^2/7)
Model Construction Examples
# Using the default database
model = Wilson(["water","ethanol"]) #Default pure model: PR
model = Wilson(["water","ethanol"],puremodel = BasicIdeal) #Using Ideal Gas for pure model properties
model = Wilson(["water","ethanol"],puremodel = PCSAFT) #Using Real Gas model for pure model properties
# Passing a prebuilt model
my_puremodel = AbbottVirial(["water","ethanol"]; userlocations = ["path/to/my/db","critical.csv"])
mixing = Wilson(["water","ethanol"],puremodel = my_puremodel)
# Using user-provided parameters
# Passing files or folders
model = Wilson(["water","ethanol"];userlocations = ["path/to/my/db","wilson.csv"])
# Passing parameters directly
model = Wilson(["water","ethanol"],
userlocations = (g = [0.0 3988.52; 1360.117 0.0],
Tc = [647.13, 513.92],
Pc = [2.19e7, 6.12e6],
acentricfactor = [0.343, 0.643],
Mw = [18.015, 46.069])
)
References
- Wilson, G. M. (1964). Vapor-liquid equilibrium. XI. A new expression for the excess free energy of mixing. Journal of the American Chemical Society, 86(2), 127–130. doi:10.1021/ja01056a002
Clapeyron.NRTL
— TypeNRTL <: ActivityModel
function NRTL(components;
puremodel=PR,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
a
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parameterb
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parameterc
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction ParameterMw
: Single Parameter (Float64
) (Optional) - Molecular Weight[g/mol]
Input models
puremodel
: model to calculate pure pressure-dependent properties
Description
NRTL (Non Random Two Fluid) activity model:
Gᴱ = nRT∑[xᵢ(∑τⱼᵢGⱼᵢxⱼ)/(∑Gⱼᵢxⱼ)]
Gᵢⱼ exp(-cᵢⱼτᵢⱼ)
τᵢⱼ = aᵢⱼ + bᵢⱼ/T
Model Construction Examples
# Using the default database
model = NRTL(["water","ethanol"]) #Default pure model: PR
model = NRTL(["water","ethanol"],puremodel = BasicIdeal) #Using Ideal Gas for pure model properties
model = NRTL(["water","ethanol"],puremodel = PCSAFT) #Using Real Gas model for pure model properties
# Passing a prebuilt model
my_puremodel = AbbottVirial(["water","ethanol"]; userlocations = ["path/to/my/db","critical.csv"])
mixing = NRTL(["water","ethanol"],puremodel = my_puremodel)
# Using user-provided parameters
# Passing files or folders
model = NRTL(["water","ethanol"];userlocations = ["path/to/my/db","nrtl_ge.csv"])
# Passing parameters directly
model = NRTL(["water","ethanol"],
userlocations = (a = [0.0 3.458; -0.801 0.0],
b = [0.0 -586.1; 246.2 0.0],
c = [0.0 0.3; 0.3 0.0])
)
References
- Renon, H., & Prausnitz, J. M. (1968). Local compositions in thermodynamic excess functions for liquid mixtures. AIChE journal. American Institute of Chemical Engineers, 14(1), 135–144. doi:10.1002/aic.690140124
Clapeyron.aspenNRTL
— TypeaspenNRTL <: ActivityModel
function aspenNRTL(components;
puremodel=PR,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
a0
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parametera1
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parametert0
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parametert1
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parametert2
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parametert3
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parameter
Input models
puremodel
: model to calculate pure pressure-dependent properties
Description
NRTL (Non Random Two Fluid) activity model:
Gᴱ = nRT∑[xᵢ(∑τⱼᵢGⱼᵢxⱼ)/(∑Gⱼᵢxⱼ)]
Gᵢⱼ exp(-αᵢⱼτᵢⱼ)
αᵢⱼ = αᵢⱼ₀ + αᵢⱼ₁T
τᵢⱼ = tᵢⱼ₀ + tᵢⱼ₁/T + tᵢⱼ₂*ln(T) + tᵢⱼ₃*T
Model Construction Examples
# Using the default database
model = aspenNRTL(["water","ethanol"]) #Default pure model: PR
model = aspenNRTL(["water","ethanol"],puremodel = BasicIdeal) #Using Ideal Gas for pure model properties
model = aspenNRTL(["water","ethanol"],puremodel = PCSAFT) #Using Real Gas model for pure model properties
# Passing a prebuilt model
my_puremodel = AbbottVirial(["water","ethanol"]; userlocations = ["path/to/my/db","critical.csv"])
mixing = aspenNRTL(["water","ethanol"],puremodel = my_puremodel)
# Using user-provided parameters
# Passing files or folders
model = aspenNRTL(["water","ethanol"];userlocations = ["path/to/my/db","nrtl.csv"])
# Passing parameters directly
model = aspenNRTL(["water","acetone"],
userlocations = (a0 = [0.0 0.228632; 0.228632 0.0],
a1 = [0.0 0.000136; 0.000136 0.0],
t0 = [0.0 13.374756; 16.081848 0.0],
t1 = [0.0 -415.527344; -88.8125 0.0],
t2 = [0.0 -1.913689; -2.930901 0.0],
t3 = [0.0 0.00153; 0.005305 0.0])
)
References
- Renon, H., & Prausnitz, J. M. (1968). Local compositions in thermodynamic excess functions for liquid mixtures. AIChE journal. American Institute of Chemical Engineers, 14(1), 135–144. doi:10.1002/aic.690140124
Clapeyron.UNIQUAC
— TypeUNIQUACModel <: ActivityModel
UNIQUAC(components;
puremodel = PR,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
a
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary Interaction Energy Parameterr
: Single Parameter (Float64
) - Normalized Van der Vals volumeq
: Single Parameter (Float64
) - Normalized Surface Areaq_p
: Single Parameter (Float64
) - Modified Normalized Surface AreaMw
: Single Parameter (Float64
) - Molecular Weight[g/mol]
Input models
puremodel
: model to calculate pure pressure-dependent properties
UNIQUAC (Universal QuasiChemical Activity Coefficients) activity model:
Gᴱ = nRT(gᴱ(comb) + gᴱ(res))
gᴱ(comb) = ∑[xᵢlog(Φᵢ/xᵢ) + 5qᵢxᵢlog(θᵢ/Φᵢ)]
gᴱ(res) = -∑xᵢqᵖᵢlog(∑θᵖⱼτⱼᵢ)
θᵢ = qᵢxᵢ/∑qᵢxᵢ
θᵖ = qᵖᵢxᵢ/∑qᵖᵢxᵢ
Φᵢ = rᵢxᵢ/∑rᵢxᵢ
τᵢⱼ = exp(-aᵢⱼ/T)
Model Construction Examples
# Using the default database
model = UNIQUAC(["water","ethanol"]) #Default pure model: PR
model = UNIQUAC(["water","ethanol"],puremodel = BasicIdeal) #Using Ideal Gas for pure model properties
model = UNIQUAC(["water","ethanol"],puremodel = PCSAFT) #Using Real Gas model for pure model properties
# Passing a prebuilt model
my_puremodel = AbbottVirial(["water","ethanol"]; userlocations = ["path/to/my/db","critical.csv"])
mixing = UNIQUAC(["water","ethanol"],puremodel = my_puremodel)
# Using user-provided parameters
# Passing files or folders
model = UNIQUAC(["water","ethanol"];userlocations = ["path/to/my/db","uniquac_ge.csv"])
# Passing parameters directly
model = UNIQUAC(["water","ethanol"],
userlocations = (a = [0.0 378.1; 258.4 0.0],
r = [0.92, 2.11],
q = [1.4, 1.97],
q_p = [1.0, 0.92],
Mw = [18.015, 46.069])
)
References
- Abrams, D. S., & Prausnitz, J. M. (1975). Statistical thermodynamics of liquid mixtures: A new expression for the excess Gibbs energy of partly or completely miscible systems. AIChE journal. American Institute of Chemical Engineers, 21(1), 116–128. doi:10.1002/aic.690210115
Clapeyron.ogUNIFAC
— TypeogUNIFACModel <: UNIFACModel
ogUNIFAC(components;
puremodel = PR,
userlocations = String[],
group_userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
R
: Single Parameter (Float64
) - Normalized group Van der Vals volumeQ
: Single Parameter (Float64
) - Normalized group Surface AreaA
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy Parameter
Input models
puremodel
: model to calculate pure pressure-dependent properties
UNIFAC (UNIQUAC Functional-group Activity Coefficients) activity model.
Original formulation.
The Combinatorial part corresponds to an GC-averaged modified UNIQUAC model. The residual part iterates over groups instead of components.
Gᴱ = nRT(gᴱ(comb) + gᴱ(res))
Combinatorial part:
gᴱ(comb) = ∑[xᵢlog(Φᵢ/xᵢ) + 5qᵢxᵢlog(θᵢ/Φᵢ)]
θᵢ = qᵢxᵢ/∑qᵢxᵢ
Φᵢ = rᵢxᵢ/∑rᵢxᵢ
rᵢ = ∑Rₖνᵢₖ for k ∈ groups
qᵢ = ∑Qₖνᵢₖ for k ∈ groups
Residual Part:
gᴱ(residual) = -v̄∑XₖQₖlog(∑ΘₘΨₘₖ)
v̄ = ∑∑xᵢνᵢₖ for k ∈ groups, for i ∈ components
Xₖ = (∑xᵢνᵢₖ)/v̄ for i ∈ components
Θₖ = QₖXₖ/∑QₖXₖ
Ψₖₘ = exp(-(Aₖₘ/T)
References
- Fredenslund, A., Gmehling, J., Michelsen, M. L., Rasmussen, P., & Prausnitz, J. M. (1977). Computerized design of multicomponent distillation columns using the UNIFAC group contribution method for calculation of activity coefficients. Industrial & Engineering Chemistry Process Design and Development, 16(4), 450–462. doi:10.1021/i260064a004
Clapeyron.UNIFAC
— TypeUNIFACModel <: ActivityModel
UNIFAC(components;
puremodel = PR,
userlocations = String[],
group_userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
R
: Single Parameter (Float64
) - Normalized group Van der Vals volumeQ
: Single Parameter (Float64
) - Normalized group Surface AreaA
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy ParameterB
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy ParameterC
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy Parameter
Input models
puremodel
: model to calculate pure pressure-dependent properties
Description
UNIFAC (UNIQUAC Functional-group Activity Coefficients) activity model. Modified UNIFAC (Dortmund) implementation. The Combinatorial part corresponds to an GC-averaged modified UNIQUAC
model. The residual part iterates over groups instead of components.
Gᴱ = nRT(gᴱ(comb) + gᴱ(res))
Combinatorial part:
gᴱ(comb) = ∑[xᵢlog(Φ'ᵢ) + 5qᵢxᵢlog(θᵢ/Φᵢ)]
θᵢ = qᵢxᵢ/∑qᵢxᵢ
Φᵢ = rᵢxᵢ/∑rᵢxᵢ
Φ'ᵢ = rᵢ^(0.75)/∑xᵢrᵢ^(0.75)
rᵢ = ∑Rₖνᵢₖ for k ∈ groups
qᵢ = ∑Qₖνᵢₖ for k ∈ groups
Residual Part:
gᴱ(residual) = -v̄∑XₖQₖlog(∑ΘₘΨₘₖ)
v̄ = ∑∑xᵢνᵢₖ for k ∈ groups, for i ∈ components
Xₖ = (∑xᵢνᵢₖ)/v̄ for i ∈ components
Θₖ = QₖXₖ/∑QₖXₖ
Ψₖₘ = exp(-(Aₖₘ + BₖₘT + CₖₘT²)/T)
References
- Fredenslund, A., Gmehling, J., Michelsen, M. L., Rasmussen, P., & Prausnitz, J. M. (1977). Computerized design of multicomponent distillation columns using the UNIFAC group contribution method for calculation of activity coefficients. Industrial & Engineering Chemistry Process Design and Development, 16(4), 450–462. doi:10.1021/i260064a004
- Weidlich, U.; Gmehling, J. A modified UNIFAC model. 1. Prediction of VLE, hE, and.gamma..infin. Ind. Eng. Chem. Res. 1987, 26, 1372–1381.
List of groups available
Name | Description |
---|---|
CH3 | Methyl group |
CH2 | Methylene group |
CH | |
C | |
CH2=CH | |
CH=CH | |
CH2=C | |
CH=C | |
ACH | Aromatic CH |
AC | |
ACCH3 | |
ACCH2 | |
ACCH | |
OH (P) | Primary alcohol |
CH3OH | Methanol |
H2O | Water |
ACOH | |
CH3CO | Methyl ketone |
CH2CO | Methylene ketone |
CHO | |
CH3COO | Acetate group |
CH2COO | |
HCOO | |
CH3O | Methoxy group |
CH2O | |
CHO | |
THF | Tetrahydrofuran |
CH3NH2 | Methylamine |
CH2NH2 | |
CHNH2 | |
CH3NH | |
CH2NH | |
CHNH | |
CH3N | |
CH2N | |
ACNH2 | |
AC2H2N | |
AC2HN | |
AC2N | |
CH3CN | Acetonitrile |
CH2CN | |
COO | Ester group |
COOH | Carboxylate group |
HCOOH | Formic acid |
CH2CL | |
CHCL | |
CCL | |
CH2CL2 | Dichloromethane |
CHCL2 | |
CCL2 | |
CHCL3 | Chloroform |
CCL3 | |
CCL4 | Carbon tetrachloride |
ACCL | |
CH3NO2 | Nitromethane |
CH2NO2 | |
CHNO2 | |
ACNO2 | |
CS2 | Carbon disulfide |
CH3SH | Methanethiol |
CH2SH | |
FURFURAL | Furfural |
DOH | |
I | Iodine group |
BR | Bromine group |
CH=-C | |
C=-C | |
DMSO | Dimethyl sulfoxide |
ACRY | Acrylate |
CL-(C=C) | |
C=C | |
ACF | |
DMF | Dimethylformamide |
HCON(CH2)2 | |
CF3 | |
CF2 | |
CF | |
CY-CH2 | Cycloalkane group |
CY-CH | |
CY-C | |
OH (S) | Second hydroxyl group |
OH (T) | Tertiary hydroxyl group |
CY-CH2O | |
TRIOXAN | Trioxane |
CNH2 | |
NMP | N-Methylpyrrolidone |
NEP | |
NIPP | |
NTBP | |
CONH2 | |
CONHCH3 | |
CONHCH2 |
Clapeyron.PSRKUNIFAC
— TypePSRKUNIFAC(components;
puremodel = BasicIdeal,
userlocations = String[],
group_userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
R
: Single Parameter (Float64
) - Normalized group Van der Vals volumeQ
: Single Parameter (Float64
) - Normalized group Surface AreaA
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy ParameterB
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy ParameterC
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy Parameter
Input models
puremodel
: model to calculate pure pressure-dependent properties
Description
UNIFAC (UNIQUAC Functional-group Activity Coefficients) activity model.
Modified UNIFAC (Dortmund) implementation, with parameters tuned to the Predictive Soave-Redlich-Kwong (PSRK) EoS.
References
- Fredenslund, A., Gmehling, J., Michelsen, M. L., Rasmussen, P., & Prausnitz, J. M. (1977). Computerized design of multicomponent distillation columns using the UNIFAC group contribution method for calculation of activity coefficients. Industrial & Engineering Chemistry Process Design and Development, 16(4), 450–462. doi:10.1021/i260064a004
- Weidlich, U.; Gmehling, J. A modified UNIFAC model. 1. Prediction of VLE, hE, and.gamma..infin. Ind. Eng. Chem. Res. 1987, 26, 1372–1381.
- Horstmann, S., Jabłoniec, A., Krafczyk, J., Fischer, K., & Gmehling, J. (2005). PSRK group contribution equation of state: comprehensive revision and extension IV, including critical constants and α-function parameters for 1000 components. Fluid Phase Equilibria, 227(2), 157–164. doi:10.1016/j.fluid.2004.11.002"
Clapeyron.VTPRUNIFAC
— TypeVTPRUNIFACModel <: UNIFACModel
VTPRUNIFAC(components;
puremodel = BasicIdeal,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
Q
: Single Parameter (Float64
) - Normalized group Surface AreaA
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy ParameterB
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy ParameterC
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy Parameter
Input models
puremodel
: model to calculate pure pressure-dependent properties
Description
UNIFAC (UNIQUAC Functional-group Activity Coefficients) activity model.
Modified UNIFAC (Dortmund) implementation, only residual part, activity model used for the Volume-Translated Peng-Robinson (VTPR) EoS.
The residual part iterates over groups instead of components.
Gᴱ = nRT(gᴱ(res))
gᴱ(res) = -v̄∑XₖQₖlog(∑ΘₘΨₘₖ)
v̄ = ∑∑xᵢνᵢₖ for k ∈ groups, for i ∈ components
Xₖ = (∑xᵢνᵢₖ)/v̄ for i ∈ components
Θₖ = QₖXₖ/∑QₖXₖ
Ψₖₘ = exp(-(Aₖₘ + BₖₘT + CₖₘT²)/T)
References
- Fredenslund, A., Gmehling, J., Michelsen, M. L., Rasmussen, P., & Prausnitz, J. M. (1977). Computerized design of multicomponent distillation columns using the UNIFAC group contribution method for calculation of activity coefficients. Industrial & Engineering Chemistry Process Design and Development, 16(4), 450–462. "doi:10.1021/i260064a004"
- Weidlich, U.; Gmehling, J. A modified UNIFAC model. 1. Prediction of VLE, hE, and.gamma..infin. Ind. Eng. Chem. Res. 1987, 26, 1372–1381.
- Ahlers, J., & Gmehling, J. (2001). Development of an universal group contribution equation of state. Fluid Phase Equilibria, 191(1–2), 177–188. doi:10.1016/s0378-3812(01)00626-4
Clapeyron.UNIFACFV
— TypeUNIFACFVModel <: ActivityModel
UNIFACFV(components;
puremodel = PR,
userlocations = String[],
group_userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
volume
: Single Parameter (Float64
) - specific volume of species[g/cm^3]
R
: Single Parameter (Float64
) - Normalized group Van der Vals volumeQ
: Single Parameter (Float64
) - Normalized group Surface AreaA
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy ParameterMw
: Single Parameter (Float64
) - Molecular weight of groups
Input models
puremodel
: model to calculate pure pressure-dependent properties
Description
UNIFACFV (UNIFAC Free Volume) activity model. specialized for solvent-polymer mixtures
The Combinatorial part corresponds to an GC-averaged modified UNIQUAC
model.
Gᴱ = nRT(gᴱ(comb) + gᴱ(res) + gᴱ(FV))
References
Clapeyron.UNIFACFVPoly
— TypeUNIFACFVPolyModel <: ActivityModel
UNIFACFVPoly(components;
puremodel = PR,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
volume
: Single Parameter (Float64
) - specific volume of species[g/cm^3]
c
Single Parameter (Float64
) - number of external degrees of freedom per solvent moleculeR
: Single Parameter (Float64
) - Normalized group Van der Vals volumeQ
: Single Parameter (Float64
) - Normalized group Surface AreaA
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Binary group Interaction Energy ParameterMw
: Single Parameter (Float64
) - Molecular weight of groups
Input models
puremodel
: model to calculate pure pressure-dependent properties
Description
UNIFAC-FV (polymer) (UNIFAC Free Volume) activity model. specialized for polymer blends
The Combinatorial part corresponds to an GC-averaged modified UNIQUAC
model.
Gᴱ = nRT(gᴱ(comb) + gᴱ(res) + gᴱ(FV))
Clapeyron.tcPRWilsonRes
— TypetcPRWilsonRes <: WilsonModel
tcPRWilsonRes(components;
puremodel = BasicIdeal,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters
g
: Pair Parameter (Float64
, asymetrical, defaults to0
) - Interaction Parameterv
: Single Parameter (optional) (Float64
) - individual volumes.
Input models
puremodel
: model to calculate pure pressure-dependent properties
Description
tc-PR-Wilson residual activity model, meant to used in combination with the tc-PR cubic EoS:
Gᴱᵣ = nRT( ∑xᵢlog(∑xⱼjΛᵢⱼ) - ∑xᵢ*(log(vᵢ/v) )
Λᵢⱼ = exp(-gᵢⱼ/T)*Vⱼ/Vᵢ
Model Construction Examples
# Using the default database
model = tcPRWilsonRes(["water","ethanol"]) #Default pure model: PR
model = tcPRWilsonRes(["water","ethanol"],puremodel = BasicIdeal) #Using Ideal Gas for pure model properties
model = tcPRWilsonRes(["water","ethanol"],puremodel = PCSAFT) #Using Real Gas model for pure model properties
# Passing a prebuilt model
my_puremodel = AbbottVirial(["water","ethanol"]; userlocations = ["path/to/my/db","critical.csv"])
mixing = tcPRWilsonRes(["water","ethanol"],puremodel = my_puremodel)
# Using user-provided parameters
# Passing files or folders
model = tcPRWilsonRes(["water","ethanol"];userlocations = ["path/to/my/db","tcpr_wilson.csv"])
# Passing parameters directly
model = tcPRWilsonRes(["water","ethanol"],userlocations = (;g = [0.0 3988.52; 1360.117 0.0]))
References
- Wilson, G. M. (1964). Vapor-liquid equilibrium. XI. A new expression for the excess free energy of mixing. Journal of the American Chemical Society, 86(2), 127–130. doi:10.1021/ja01056a002
- Piña-Martinez, A., Privat, R., Nikolaidis, I. K., Economou, I. G., & Jaubert, J.-N. (2021). What is the optimal activity coefficient model to be combined with the translated–consistent Peng–Robinson equation of state through advanced mixing rules? Cross-comparison and grading of the Wilson, UNIQUAC, and NRTL aE models against a benchmark database involving 200 binary systems. Industrial & Engineering Chemistry Research, 60(47), 17228–17247. doi:10.1021/acs.iecr.1c03003
Clapeyron.COSMOSAC02
— TypeCOSMOSAC02(components;
puremodel = PR,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters:
Pi
:Single Parameter{String}V
: Single Parameter{Float64}A
: Single Parameter{Float64}
Description
An activity coefficient model using molecular solvation based on the COSMO-RS method.
References
- Klamt, A. (1995). Conductor-like screening model for real solvents: A new approach to the quantitative calculation of solvation phenomena. Journal of Physical Chemistry, 99(7), 2224–2235. doi:10.1021/j100007a062
- Lin, S-T. & Sandler, S.I. (2002). A priori phase equilibrium prediction from a segment contribution solvation model. Industrial & Engineering Chemistry Research, 41(5), 899–913. doi:10.1021/ie001047w
Clapeyron.COSMOSAC10
— TypeCOSMOSAC10(components;
puremodel = PR,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters:
Pnhb
:Single Parameter{String}POH
:Single Parameter{String}POT
:Single Parameter{String}V
: Single Parameter{Float64}A
: Single Parameter{Float64}
Description
An activity coefficient model using molecular solvation based on the COSMO-RS method. Sigma profiles are now split by non-hydrogen bonding, hydrogen acceptor and hydrogen donor.
References
- Klamt, A. (1995). Conductor-like screening model for real solvents: A new approach to the quantitative calculation of solvation phenomena. Journal of Physical Chemistry, 99(7), 2224–2235. doi:10.1021/j100007a062
- Lin, S-T. & Sandler, S.I. (2002). A priori phase equilibrium prediction from a segment contribution solvation model. Industrial & Engineering Chemistry Research, 41(5), 899–913. doi:10.1021/ie001047w
- Hsieh, C-H., Sandler, S.I., & Lin, S-T. (2010). Improvements of COSMO-SAC for vapor–liquid and liquid–liquid equilibrium predictions. Fluid Phase Equilibria, 297(1), 90-97. doi:10.1016/j.fluid.2010.06.011
Clapeyron.COSMOSACdsp
— TypeCOSMOSACdsp(components;
puremodel = PR,
userlocations = String[],
pure_userlocations = String[],
verbose = false)
Input parameters:
Pnhb
:Single Parameter{String}POH
:Single Parameter{String}POT
:Single Parameter{String}V
: Single Parameter{Float64}A
: Single Parameter{Float64}epsilon
: Single Parameter{Float64}COOH
: Single Parameter{Float64}water
: Single Parameter{Float64}hb_acc
: Single Parameter{Float64}hb_don
: Single Parameter{Float64}
Description
An activity coefficient model using molecular solvation based on the COSMO-RS method. Sigma profiles are now split by non-hydrogen bonding, hydrogen acceptor and hydrogen donor. A dispersive correction is included.
References
- Klamt, A. (1995). Conductor-like screening model for real solvents: A new approach to the quantitative calculation of solvation phenomena. Journal of Physical Chemistry, 99(7), 2224–2235. doi:10.1021/j100007a062
- Lin, S-T. & Sandler, S.I. (2002). A priori phase equilibrium prediction from a segment contribution solvation model. Industrial & Engineering Chemistry Research, 41(5), 899–913. doi:10.1021/ie001047w
- Hsieh, C-H., Sandler, S.I., & Lin, S-T. (2010). Improvements of COSMO-SAC for vapor–liquid and liquid–liquid equilibrium predictions. Fluid Phase Equilibria, 297(1), 90-97. doi:10.1016/j.fluid.2010.06.011
- Hsieh, C-H., Lin, S-T. & Vrabec, J. (2014). Considering the dispersive interactions in the COSMO-SAC model for more accurate predictions of fluid phase behavior. Fluid Phase Equilibria, 367, 109-116. doi:10.1016/j.fluid.2014.01.032