Index
Clapeyron.AlyLeeIdeal
Clapeyron.BasicIdeal
Clapeyron.JobackIdeal
Clapeyron.LJRefIdeal
Clapeyron.MonomerIdeal
Clapeyron.PPDSIdeal
Clapeyron.ReidIdeal
Clapeyron.WalkerIdeal
Ideal Models
All Clapeyron.jl
models can be separated between an ideal and a residual contribution. The ideal contribution can be obtained via integration of the ideal isobaric heat capacity:
$\frac{A_{\mathrm{ideal}}}{Nk_\mathrm{B}T} = \sum_{i=1}^{N_{\mathrm{Component}}} x_i\left[\ln{\frac{\rho_i}{\rho_0}} + \frac{1}{Nk_\mathrm{B}T} \int_{T_0}^T \!\!C_{p,i}^0 dT + \frac{H_{0,i}}{Nk_\mathrm{B}T}- \frac{1}{Nk_{B}}\!\!\int_{T_0}^T \frac{C_{p,i}^0}{T} dT -\ln{\frac{T}{T_0}}-\frac{S_{0,i}}{Nk_\mathrm{B}} - 1\right]$
Normally, EoS models contain an ideal model. this model can be accessed by using Clapeyron.idealmodel
.
Clapeyron.BasicIdeal
— TypeBasicIdeal <: IdealModel
BasicIdeal(components;
userlocations = String[],
reference_state = nothing,
verbose = false)
Input parameters
None
Description
Default Ideal Model. Constant specific heat capacity equal to 5R/2
. it's Helmholtz energy is equal to:
a₀ = A₀/nRT = ∑(xᵢlog(nxᵢ/V)) - 1 - 1.5log(T)
Model Construction Examples
# Because this model does not have parameters, all those constructors are equivalent:
idealmodel = BasicIdeal()
idealmodel = BasicIdeal("water")
idealmodel = BasicIdeal(["water","carbon dioxide"])
Clapeyron.ReidIdeal
— TypeReidIdeal <: IdealModel
ReidIdeal(components;
userlocations = String[],
reference_state = nothing,
verbose = false)
Input parameters
a
: Single Parameter (Float64
) - polynomial coefficientb
: Single Parameter (Float64
) - polynomial coefficientc
: Single Parameter (Float64
) - polynomial coefficientd
: Single Parameter (Float64
) - polynomial coefficiente
: Single Parameter (optional) (Float64
) - polynomial coefficient
Model parameters
coeffs
: Single Parameter (NTuple{5,Float64}
)
Description
Reid Ideal Model. Helmholtz energy obtained via integration of specific heat capacity:
Cpᵢ(T) = aᵢ + bᵢT + cᵢT^2 + dᵢT^3 + eᵢT^4
Cp(T) = ∑Cpᵢxᵢ
Model Construction Examples
# Using the default database
idealmodel = ReidIdeal("water") #single input
idealmodel = ReidIdeal(["water","ethanol"]) #multiple components
# Using user-provided parameters
# Passing files or folders
idealmodel = ReidIdeal(["neon","hydrogen"]; userlocations = ["path/to/my/db","reid.csv"])
# Passing parameters directly
idealmodel = ReidIdeal(["water","butane"];
userlocations = (a = [32.24, 9.487],
b = [0.00192, 0.3313],
c = [1.06e-5, -0.0001108],
d = [-3.6e-9, -2.822e-9])
) #e is not used
Clapeyron.JobackIdeal
— TypeJobackIdeal <: JobackIdealModel
JobackIdeal(components;
userlocations = String[],
reference_state = nothing,
verbose = false)
Input parameters
Mw
: Single Parameter (Float64
) - Molecular Weight[g/mol]
N_a
: Single Parameter (Float64
)T_c
: Single Parameter (Float64
)P_c
: Single Parameter (Float64
)V_c
: Single Parameter (Float64
)T_b
: Single Parameter (Float64
)T_m
: Single Parameter (Float64
)H_form
: Single Parameter (Float64
)G_form
: Single Parameter (Float64
)a
: Single Parameter (Float64
)b
: Single Parameter (Float64
)c
: Single Parameter (Float64
)d
: Single Parameter (Float64
)H_fusion
: Single Parameter (Float64
)H_vap
: Single Parameter (Float64
)eta_a
: Single Parameter (Float64
)eta_b
: Single Parameter (Float64
)
Description
Joback Group Contribution Ideal Model. GC version of ReidIdeal
. Helmholtz energy obtained via integration of specific heat capacity:
aᵢ = ∑(νᵢₖbₖ) - 37.93
bᵢ = ∑(νᵢₖbₖ) + 0.210
cᵢ = ∑(νᵢₖcₖ) - 3.91e-4
dᵢ = ∑(νᵢₖbₖ) + 2.06e-7
Cpᵢ(T) = aᵢ + bᵢT + cᵢT^2 + dᵢT^3
The GC-averaged Reid Model is available by using ReidIdeal(model::JobackIdeal)
.
The estimated critical point of a single component can be obtained via crit_pure(model::JobackIdeal)
References
- Joback, K. G., & Reid, R. C. (1987). Estimation of pure-component properties from group-contributions. Chemical Engineering Communications, 57(1–6), 233–243. doi:10.1080/00986448708960487
List of available groups
Name | Description |
---|---|
-CH3 | Methyl |
-CH2- | Methylene |
>CH- | |
>C< | |
CH2=CH- | |
-CH=CH- | |
=C< | |
=C= | |
CH | |
C | |
ring-CH2- | Cyclic alkane |
ring>CH- | |
ring>C< | |
ring=CH- | Aromatic group |
ring=C< | |
-F | Fluoride |
-Cl | Chloride |
-Br | Bromide |
-I | Iodide |
-OH (alcohol) | Hydroxyl group |
-OH (phenol) | |
-O- (non-ring) | |
-O- (ring) | |
>C=O (non-ring) | Ketone |
>C=O (ring) | |
O=CH- (aldehyde) | Aldehyde |
-COOH (acid) | Carboxylic acid |
-COO- (ester) | Ester |
O (other than above) | Ketone |
-NH2 | Amine |
>NH (non-ring) | |
>NH (ring) | |
>N- (non-ring) | |
-N= (non-ring) | |
-N= (ring) | |
=NH | |
-CN | Nitrile |
-NO3 | Nitroxide |
-SH | |
-S- (non-ring) | |
-S- (ring) |
Clapeyron.MonomerIdeal
— TypeMonomerIdeal <: MonomerIdealModel
MonomerIdeal(components;
userlocations = String[],
reference_state = nothing,
verbose = false)
Input parameters
Mw
: Single Parameter (Float64
) - Molecular Weight[g/mol]
Model Parameters
None
Description
Monomer Ideal Model, result obtained from statistical mechanics Λ
Λᵢ = h/√(kᵦTMwᵢ/Nₐ)
a₀ = A₀/nRT = ∑xᵢlog(ρᵢΛᵢ^3)
Model Construction Examples
# Using the default database
idealmodel = MonomerIdeal("water") #single input
idealmodel = MonomerIdeal(["water","ethanol"]) #multiple components
# Using user-provided parameters
# Passing files or folders
idealmodel = MonomerIdeal(["neon","hydrogen"]; userlocations = ["path/to/my/db","mw.csv"])
# Passing parameters directly
idealmodel = MonomerIdeal(["neon","hydrogen"];userlocations = (;Mw = [20.17, 2.]))
Clapeyron.WalkerIdeal
— TypeWalkerIdeal <: WalkerIdealModel
WalkerIdeal(components;
userlocations = String[],
group_userlocations = String[]
verbose = false)
Input parameters
Mw
: Single Parameter (Float64
) - Molecular Weight[g/mol]
Nrot
: Single Parameter (Int
)theta1
: Single Parameter (Float64
)theta2
: Single Parameter (Float64
)theta3
: Single Parameter (Float64
)theta4
: Single Parameter (Float64
)deg1
: Single Parameter (Int
)deg2
: Single Parameter (Int
)deg3
: Single Parameter (Int
)deg4
: Single Parameter (Int
)
Description
Walker [1] Group Contribution Ideal Model.
Cpᵢ(T)/R = (5+NRot)/2 ∑νᵢₖ∑gₖᵥ(θₖᵥ/T)^2*exp(θₖᵥ/T)/(1-exp(θₖᵥ/T)) , v ∈ 1:4
References
- Walker, P. J., & Haslam, A. J. (2020). A new predictive group-contribution ideal-heat-capacity model and its influence on second-derivative properties calculated using a free-energy equation of state. Journal of Chemical and Engineering Data, 65(12), 5809–5829. doi:10.1021/acs.jced.0c00723
Clapeyron.LJRefIdeal
— TypeLJRefIdeal <: IdealModel
LJRef(components;
userlocations = String[],
verbose = false)
Input parameters
sigma
: Single Parameter (Float64
) - particle size [Å]epsilon
: Single Parameter (Float64
) - dispersion energy [K
]Mw
: Single Parameter (Float64
) - Molecular Weight[g/mol]
Description
Lennard-Jones Reference equation of state. Ideal Part. valid from 0.5 < T/Tc < 7 and pressures up to p/pc = 500.
τᵢ = 1.32ϵᵢ/T
δᵢ = n(Nₐσᵢ^3)/0.31V
a⁰ᵢ(δ,τ) = log(δᵢ) + 1.5log(τᵢ) - 1.515151515τᵢ + 6.262265814
a⁰(δ,τ,z) = ∑xᵢ(a⁰ᵢ + log(xᵢ))
LJRefIdeal
acts as a wrapper of LJRef
model, you can access it with LJRef(model::LJRefIdeal)
.
The original model was done with only one component in mind. to support multiple components, a VDW 1-fluid mixing rule (shown above) is implemented, but it is not tested.
References
- Thol, M., Rutkai, G., Köster, A., Lustig, R., Span, R., & Vrabec, J. (2016). Equation of state for the Lennard-Jones fluid. Journal of physical and chemical reference data, 45(2), 023101. doi:10.1063/1.4945000
Clapeyron.AlyLeeIdeal
— TypeAlyLeeIdeal <: AlyLeeIdealModel
AlyLeeIdeal(components;
userlocations = String[],
reference_state = nothing,
verbose = false)
Input parameters
A
: Single Parameter (Float64
) - Model CoefficientB
: Single Parameter (Float64
) - Model CoefficientC
: Single Parameter (Float64
) - Model CoefficientD
: Single Parameter (Float64
) - Model CoefficientE
: Single Parameter (Float64
) - Model CoefficientF
: Single Parameter (Float64
) - Model CoefficientG
: Single Parameter (Float64
) - Model CoefficientH
: Single Parameter (Float64
) - Model CoefficientI
: Single Parameter (Float64
) - Model Coefficient
Description
Aly-Lee Ideal Model (extended):
Cpᵢ(T)/R = A + B(CT⁻¹/sinh(CT⁻¹))² + D(ET⁻¹/cosh(ET⁻¹))² + F(GT⁻¹/sinh(GT⁻¹))² + H(IT⁻¹/cosh(IT⁻¹))²
Model Construction Examples
# Using the default database
idealmodel = AlyLeeIdeal("water") #single input
idealmodel = AlyLeeIdeal(["water","ethanol"]) #multiple components
# Using user-provided parameters
# Passing files or folders
idealmodel = AlyLeeIdeal(["neon","hydrogen"]; userlocations = ["path/to/my/db","alylee.csv"])
# Passing parameters directly
idealmodel = AlyLeeIdeal(["water","carbon dioxide"];
userlocations = (A = [4.004, 3.5],
B = [0.01, 2.044],
C = [268.8, 919.3],
D = [0.99, -1.06],
E = [1141.4, -865.1],
F = [3.07, 2.034],
G = [2507.37, 483.55],
H = [0.0, 0.0139],
I = [0.0, 341.11])
)
References
- Aly, F. A., & Lee, L. L. (1981). Self-consistent equations for calculating the ideal gas heat capacity, enthalpy, and entropy. Fluid Phase Equilibria, 6(3–4), 169–179. doi:10.1016/0378-3812(81)85002-9
Clapeyron.PPDSIdeal
— TypePPDSIdeal <: PPDSIdealModel
PPDSIdeal(components;
userlocations = String[],
reference_state = nothing,
verbose = false)
Input parameters
A
: Single Parameter (Float64
) - Model CoefficientB
: Single Parameter (Float64
) - Model CoefficientC
: Single Parameter (Float64
) - Model CoefficientD
: Single Parameter (Float64
) - Model CoefficientE
: Single Parameter (Float64
) - Model CoefficientF
: Single Parameter (Float64
) - Model CoefficientG
: Single Parameter (Float64
) - Model Coefficient
Description
PPDS Ideal Model:
Cpᵢ(T)/R = B + (C - B)y²[1 + (y − 1)(D + Ey + Fy² + Gy³)]
y = T/(A + T)
Model Construction Examples
# Using the default database
idealmodel = PPDSIdeal("water") #single input
idealmodel = PPDSIdeal(["water","ethanol"]) #multiple components
# Using user-provided parameters
# Passing files or folders
idealmodel = PPDSIdeal(["neon","hydrogen"]; userlocations = ["path/to/my/db","alylee.csv"])
# Passing parameters directly
idealmodel = PPDSIdeal(["water","carbon dioxide"];
userlocations = (A = [4.004, 3.5],
B = [0.01, 2.044],
C = [268.8, 919.3],
D = [0.99, -1.06],
E = [1141.4, -865.1],
F = [3.07, 2.034],
G = [2507.37, 483.55],
H = [0.0, 0.0139],
I = [0.0, 341.11])
)
References
- Gmehling, J., Kleiber, M., Kolbe, B., & Rarey, J. (2019). Chemical thermodynamics for process simulation (2nd ed.). Berlin, Germany: Blackwell Verlag.