Skip to content

Commit 950a96e

Browse files
authored
Reimplement Molarweight trait (#258)
1 parent 576b56b commit 950a96e

File tree

28 files changed

+238
-164
lines changed

28 files changed

+238
-164
lines changed

feos-core/src/cubic.rs

+3-1
Original file line numberDiff line numberDiff line change
@@ -4,7 +4,7 @@
44
//! of state - with a single contribution to the Helmholtz energy - can be implemented.
55
//! The implementation closely follows the form of the equations given in
66
//! [this wikipedia article](https://en.wikipedia.org/wiki/Cubic_equations_of_state#Peng%E2%80%93Robinson_equation_of_state).
7-
use crate::equation_of_state::{Components, Residual};
7+
use crate::equation_of_state::{Components, Molarweight, Residual};
88
use crate::parameter::{Identifier, Parameter, ParameterError, PureRecord};
99
use crate::state::StateHD;
1010
use ndarray::{Array1, Array2, ScalarOperand};
@@ -214,7 +214,9 @@ impl Residual for PengRobinson {
214214
self.residual_helmholtz_energy(state),
215215
)]
216216
}
217+
}
217218

219+
impl Molarweight for PengRobinson {
218220
fn molar_weight(&self) -> MolarWeight<Array1<f64>> {
219221
&self.parameters.molarweight * (GRAM / MOL)
220222
}

feos-core/src/equation_of_state/mod.rs

+3-1
Original file line numberDiff line numberDiff line change
@@ -9,7 +9,7 @@ mod ideal_gas;
99
mod residual;
1010

1111
pub use ideal_gas::IdealGas;
12-
pub use residual::{EntropyScaling, NoResidual, Residual};
12+
pub use residual::{EntropyScaling, Molarweight, NoResidual, Residual};
1313

1414
/// The number of components that the model is initialized for.
1515
pub trait Components {
@@ -91,7 +91,9 @@ impl<I: IdealGas, R: Residual> Residual for EquationOfState<I, R> {
9191
) -> Vec<(String, D)> {
9292
self.residual.residual_helmholtz_energy_contributions(state)
9393
}
94+
}
9495

96+
impl<I, R: Molarweight> Molarweight for EquationOfState<I, R> {
9597
fn molar_weight(&self) -> MolarWeight<Array1<f64>> {
9698
self.residual.molar_weight()
9799
}

feos-core/src/equation_of_state/residual.rs

+8-10
Original file line numberDiff line numberDiff line change
@@ -8,7 +8,14 @@ use quantity::*;
88
use std::ops::Div;
99
use typenum::Quot;
1010

11-
/// A reisdual Helmholtz energy model.
11+
/// Molar weight of all components.
12+
///
13+
/// Enables calculation of (mass) specific properties.
14+
pub trait Molarweight {
15+
fn molar_weight(&self) -> MolarWeight<Array1<f64>>;
16+
}
17+
18+
/// A residual Helmholtz energy model.
1219
pub trait Residual: Components + Send + Sync {
1320
/// Return the maximum density in Angstrom^-3.
1421
///
@@ -18,11 +25,6 @@ pub trait Residual: Components + Send + Sync {
1825
/// equation of state anyways).
1926
fn compute_max_density(&self, moles: &Array1<f64>) -> f64;
2027

21-
/// Molar weight of all components.
22-
///
23-
/// Enables calculation of (mass) specific properties.
24-
fn molar_weight(&self) -> MolarWeight<Array1<f64>>;
25-
2628
/// Evaluate the reduced Helmholtz energy of each individual contribution
2729
/// and return them together with a string representation of the contribution.
2830
fn residual_helmholtz_energy_contributions<D: DualNum<f64> + Copy + ScalarOperand>(
@@ -193,8 +195,4 @@ impl Residual for NoResidual {
193195
) -> Vec<(String, D)> {
194196
vec![]
195197
}
196-
197-
fn molar_weight(&self) -> MolarWeight<Array1<f64>> {
198-
panic!("No mass specific properties are available for this model!")
199-
}
200198
}

feos-core/src/lib.rs

+1-1
Original file line numberDiff line numberDiff line change
@@ -33,7 +33,7 @@ pub mod parameter;
3333
mod phase_equilibria;
3434
mod state;
3535
pub use equation_of_state::{
36-
Components, EntropyScaling, EquationOfState, IdealGas, NoResidual, Residual,
36+
Components, EntropyScaling, EquationOfState, IdealGas, Molarweight, NoResidual, Residual,
3737
};
3838
pub use errors::{EosError, EosResult};
3939
pub use phase_equilibria::{

feos-core/src/python/phase_equilibria.rs

+9-7
Original file line numberDiff line numberDiff line change
@@ -943,21 +943,23 @@ macro_rules! impl_phase_equilibrium {
943943
if different_pressures {
944944
dict.insert(String::from("pressure vapor"), p_v);
945945
dict.insert(String::from("pressure liquid"), p_l);
946-
}else {
946+
} else {
947947
dict.insert(String::from("pressure"), p_v);
948948
}
949949
dict.insert(String::from("density liquid"), self.0.liquid().density().convert_to(MOL / METER.powi::<P3>()).into_raw_vec_and_offset().0);
950950
dict.insert(String::from("density vapor"), self.0.vapor().density().convert_to(MOL / METER.powi::<P3>()).into_raw_vec_and_offset().0);
951-
dict.insert(String::from("mass density liquid"), self.0.liquid().mass_density().convert_to(KILOGRAM / METER.powi::<P3>()).into_raw_vec_and_offset().0);
952-
dict.insert(String::from("mass density vapor"), self.0.vapor().mass_density().convert_to(KILOGRAM / METER.powi::<P3>()).into_raw_vec_and_offset().0);
953951
dict.insert(String::from("molar enthalpy liquid"), self.0.liquid().molar_enthalpy(contributions).convert_to(KILO * JOULE / MOL).into_raw_vec_and_offset().0);
954952
dict.insert(String::from("molar enthalpy vapor"), self.0.vapor().molar_enthalpy(contributions).convert_to(KILO * JOULE / MOL).into_raw_vec_and_offset().0);
955953
dict.insert(String::from("molar entropy liquid"), self.0.liquid().molar_entropy(contributions).convert_to(KILO * JOULE / KELVIN / MOL).into_raw_vec_and_offset().0);
956954
dict.insert(String::from("molar entropy vapor"), self.0.vapor().molar_entropy(contributions).convert_to(KILO * JOULE / KELVIN / MOL).into_raw_vec_and_offset().0);
957-
dict.insert(String::from("specific enthalpy liquid"), self.0.liquid().specific_enthalpy(contributions).convert_to(KILO * JOULE / KILOGRAM).into_raw_vec_and_offset().0);
958-
dict.insert(String::from("specific enthalpy vapor"), self.0.vapor().specific_enthalpy(contributions).convert_to(KILO * JOULE / KILOGRAM).into_raw_vec_and_offset().0);
959-
dict.insert(String::from("specific entropy liquid"), self.0.liquid().specific_entropy(contributions).convert_to(KILO * JOULE / KELVIN / KILOGRAM).into_raw_vec_and_offset().0);
960-
dict.insert(String::from("specific entropy vapor"), self.0.vapor().specific_entropy(contributions).convert_to(KILO * JOULE / KELVIN / KILOGRAM).into_raw_vec_and_offset().0);
955+
if self.0.states[0].liquid().eos.residual.has_molar_weight() {
956+
dict.insert(String::from("mass density liquid"), self.0.liquid().mass_density().convert_to(KILOGRAM / METER.powi::<P3>()).into_raw_vec_and_offset().0);
957+
dict.insert(String::from("mass density vapor"), self.0.vapor().mass_density().convert_to(KILOGRAM / METER.powi::<P3>()).into_raw_vec_and_offset().0);
958+
dict.insert(String::from("specific enthalpy liquid"), self.0.liquid().specific_enthalpy(contributions).convert_to(KILO * JOULE / KILOGRAM).into_raw_vec_and_offset().0);
959+
dict.insert(String::from("specific enthalpy vapor"), self.0.vapor().specific_enthalpy(contributions).convert_to(KILO * JOULE / KILOGRAM).into_raw_vec_and_offset().0);
960+
dict.insert(String::from("specific entropy liquid"), self.0.liquid().specific_entropy(contributions).convert_to(KILO * JOULE / KELVIN / KILOGRAM).into_raw_vec_and_offset().0);
961+
dict.insert(String::from("specific entropy vapor"), self.0.vapor().specific_entropy(contributions).convert_to(KILO * JOULE / KELVIN / KILOGRAM).into_raw_vec_and_offset().0);
962+
}
961963
dict
962964
}
963965

feos-core/src/python/state.rs

+17-15
Original file line numberDiff line numberDiff line change
@@ -1254,7 +1254,7 @@ macro_rules! impl_state {
12541254
/// SIArray1
12551255
#[pyo3(signature = (contributions=Contributions::Total), text_signature = "($self, contributions)")]
12561256
fn molar_entropy(&self, contributions: Contributions) -> MolarEntropy<Array1<f64>> {
1257-
StateVec::from(self).molar_entropy(contributions).into()
1257+
StateVec::from(self).molar_entropy(contributions)
12581258
}
12591259

12601260
/// Return mass specific entropy.
@@ -1270,7 +1270,7 @@ macro_rules! impl_state {
12701270
/// SIArray1
12711271
#[pyo3(signature = (contributions=Contributions::Total), text_signature = "($self, contributions)")]
12721272
fn specific_entropy(&self, contributions: Contributions) -> SpecificEntropy<Array1<f64>> {
1273-
StateVec::from(self).specific_entropy(contributions).into()
1273+
StateVec::from(self).specific_entropy(contributions)
12741274
}
12751275

12761276
/// Return molar enthalpy.
@@ -1286,7 +1286,7 @@ macro_rules! impl_state {
12861286
/// SIArray1
12871287
#[pyo3(signature = (contributions=Contributions::Total), text_signature = "($self, contributions)")]
12881288
fn molar_enthalpy(&self, contributions: Contributions) -> MolarEnergy<Array1<f64>> {
1289-
StateVec::from(self).molar_enthalpy(contributions).into()
1289+
StateVec::from(self).molar_enthalpy(contributions)
12901290
}
12911291

12921292
/// Return mass specific enthalpy.
@@ -1302,18 +1302,18 @@ macro_rules! impl_state {
13021302
/// SIArray1
13031303
#[pyo3(signature = (contributions=Contributions::Total), text_signature = "($self, contributions)")]
13041304
fn specific_enthalpy(&self, contributions: Contributions) -> SpecificEnergy<Array1<f64>> {
1305-
StateVec::from(self).specific_enthalpy(contributions).into()
1305+
StateVec::from(self).specific_enthalpy(contributions)
13061306
}
13071307

13081308

13091309
#[getter]
13101310
fn get_temperature(&self) -> Temperature<Array1<f64>> {
1311-
StateVec::from(self).temperature().into()
1311+
StateVec::from(self).temperature()
13121312
}
13131313

13141314
#[getter]
13151315
fn get_pressure(&self) -> Pressure<Array1<f64>> {
1316-
StateVec::from(self).pressure().into()
1316+
StateVec::from(self).pressure()
13171317
}
13181318

13191319
#[getter]
@@ -1323,12 +1323,12 @@ macro_rules! impl_state {
13231323

13241324
#[getter]
13251325
fn get_density(&self) -> Density<Array1<f64>> {
1326-
StateVec::from(self).density().into()
1326+
StateVec::from(self).density()
13271327
}
13281328

13291329
#[getter]
13301330
fn get_moles<'py>(&self, py: Python<'py>) -> Moles<Array2<f64>> {
1331-
StateVec::from(self).moles().into()
1331+
StateVec::from(self).moles()
13321332
}
13331333

13341334
#[getter]
@@ -1337,13 +1337,13 @@ macro_rules! impl_state {
13371337
}
13381338

13391339
#[getter]
1340-
fn get_mass_density(&self) -> MassDensity<Array1<f64>> {
1341-
StateVec::from(self).mass_density().into()
1340+
fn get_mass_density(&self) -> Option<MassDensity<Array1<f64>>> {
1341+
self.0[0].eos.residual.has_molar_weight().then(|| StateVec::from(self).mass_density())
13421342
}
13431343

13441344
#[getter]
1345-
fn get_massfracs<'py>(&self, py: Python<'py>) -> Bound<'py, PyArray2<f64>> {
1346-
StateVec::from(self).massfracs().into_pyarray_bound(py)
1345+
fn get_massfracs<'py>(&self, py: Python<'py>) -> Option<Bound<'py, PyArray2<f64>>> {
1346+
self.0[0].eos.residual.has_molar_weight().then(|| StateVec::from(self).massfracs().into_pyarray_bound(py))
13471347
}
13481348

13491349
/// Returns selected properties of a StateVec as dictionary.
@@ -1385,11 +1385,13 @@ macro_rules! impl_state {
13851385
dict.insert(String::from("temperature"), states.temperature().convert_to(KELVIN).into_raw_vec_and_offset().0);
13861386
dict.insert(String::from("pressure"), states.pressure().convert_to(PASCAL).into_raw_vec_and_offset().0);
13871387
dict.insert(String::from("density"), states.density().convert_to(MOL / METER.powi::<P3>()).into_raw_vec_and_offset().0);
1388-
dict.insert(String::from("mass density"), states.mass_density().convert_to(KILOGRAM / METER.powi::<P3>()).into_raw_vec_and_offset().0);
13891388
dict.insert(String::from("molar enthalpy"), states.molar_enthalpy(contributions).convert_to(KILO * JOULE / MOL).into_raw_vec_and_offset().0);
13901389
dict.insert(String::from("molar entropy"), states.molar_entropy(contributions).convert_to(KILO * JOULE / KELVIN / MOL).into_raw_vec_and_offset().0);
1391-
dict.insert(String::from("specific enthalpy"), states.specific_enthalpy(contributions).convert_to(KILO * JOULE / KILOGRAM).into_raw_vec_and_offset().0);
1392-
dict.insert(String::from("specific entropy"), states.specific_entropy(contributions).convert_to(KILO * JOULE / KELVIN / KILOGRAM).into_raw_vec_and_offset().0);
1390+
if states.0[0].eos.residual.has_molar_weight() {
1391+
dict.insert(String::from("mass density"), states.mass_density().convert_to(KILOGRAM / METER.powi::<P3>()).into_raw_vec_and_offset().0);
1392+
dict.insert(String::from("specific enthalpy"), states.specific_enthalpy(contributions).convert_to(KILO * JOULE / KILOGRAM).into_raw_vec_and_offset().0);
1393+
dict.insert(String::from("specific entropy"), states.specific_entropy(contributions).convert_to(KILO * JOULE / KELVIN / KILOGRAM).into_raw_vec_and_offset().0);
1394+
}
13931395
dict
13941396
}
13951397
}

feos-core/src/python/user_defined.rs

+3-1
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
11
#![allow(non_snake_case)]
2-
use crate::{Components, IdealGas, Residual, StateHD};
2+
use crate::{Components, IdealGas, Molarweight, Residual, StateHD};
33
use ndarray::{Array1, ScalarOperand};
44
use num_dual::*;
55
use numpy::convert::IntoPyArray;
@@ -203,7 +203,9 @@ macro_rules! impl_residual {
203203
) -> Vec<(String, D)> {
204204
vec![("Python".to_string(), self.residual_helmholtz_energy(state))]
205205
}
206+
}
206207

208+
impl Molarweight for PyResidual {
207209
fn molar_weight(&self) -> MolarWeight<Array1<f64>> {
208210
Python::with_gil(|py| {
209211
let py_result = self.0.bind(py).call_method0("molar_weight").unwrap();

feos-core/src/state/properties.rs

+41-39
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
11
use super::{Contributions, Derivative::*, PartialDerivative, State};
2-
use crate::equation_of_state::{IdealGas, Residual};
2+
use crate::equation_of_state::{IdealGas, Molarweight, Residual};
33
use crate::ReferenceSystem;
44
use ndarray::Array1;
55
use quantity::*;
@@ -74,14 +74,6 @@ impl<E: Residual + IdealGas> State<E> {
7474
self.temperature * self.ds_dt(contributions) / self.total_moles
7575
}
7676

77-
/// Specific isochoric heat capacity: $c_v^{(m)}=\frac{C_v}{m}$
78-
pub fn specific_isochoric_heat_capacity(
79-
&self,
80-
contributions: Contributions,
81-
) -> SpecificEntropy {
82-
self.molar_isochoric_heat_capacity(contributions) / self.total_molar_weight()
83-
}
84-
8577
/// Partial derivative of the molar isochoric heat capacity w.r.t. temperature: $\left(\frac{\partial c_V}{\partial T}\right)_{V,N_i}$
8678
pub fn dc_v_dt(
8779
&self,
@@ -103,11 +95,6 @@ impl<E: Residual + IdealGas> State<E> {
10395
}
10496
}
10597

106-
/// Specific isobaric heat capacity: $c_p^{(m)}=\frac{C_p}{m}$
107-
pub fn specific_isobaric_heat_capacity(&self, contributions: Contributions) -> SpecificEntropy {
108-
self.molar_isobaric_heat_capacity(contributions) / self.total_molar_weight()
109-
}
110-
11198
/// Entropy: $S=-\left(\frac{\partial A}{\partial T}\right)_{V,N_i}$
11299
pub fn entropy(&self, contributions: Contributions) -> Entropy {
113100
Entropy::from_reduced(
@@ -120,11 +107,6 @@ impl<E: Residual + IdealGas> State<E> {
120107
self.entropy(contributions) / self.total_moles
121108
}
122109

123-
/// Specific entropy: $s^{(m)}=\frac{S}{m}$
124-
pub fn specific_entropy(&self, contributions: Contributions) -> SpecificEntropy {
125-
self.molar_entropy(contributions) / self.total_molar_weight()
126-
}
127-
128110
/// Partial molar entropy: $s_i=\left(\frac{\partial S}{\partial N_i}\right)_{T,p,N_j}$
129111
pub fn partial_molar_entropy(&self) -> MolarEntropy<Array1<f64>> {
130112
let c = Contributions::Total;
@@ -160,11 +142,6 @@ impl<E: Residual + IdealGas> State<E> {
160142
self.enthalpy(contributions) / self.total_moles
161143
}
162144

163-
/// Specific enthalpy: $h^{(m)}=\frac{H}{m}$
164-
pub fn specific_enthalpy(&self, contributions: Contributions) -> SpecificEnergy {
165-
self.molar_enthalpy(contributions) / self.total_molar_weight()
166-
}
167-
168145
/// Partial molar enthalpy: $h_i=\left(\frac{\partial H}{\partial N_i}\right)_{T,p,N_j}$
169146
pub fn partial_molar_enthalpy(&self) -> MolarEnergy<Array1<f64>> {
170147
let s = self.partial_molar_entropy();
@@ -184,11 +161,6 @@ impl<E: Residual + IdealGas> State<E> {
184161
self.helmholtz_energy(contributions) / self.total_moles
185162
}
186163

187-
/// Specific Helmholtz energy: $a^{(m)}=\frac{A}{m}$
188-
pub fn specific_helmholtz_energy(&self, contributions: Contributions) -> SpecificEnergy {
189-
self.molar_helmholtz_energy(contributions) / self.total_molar_weight()
190-
}
191-
192164
/// Internal energy: $U=A+TS$
193165
pub fn internal_energy(&self, contributions: Contributions) -> Energy {
194166
self.temperature * self.entropy(contributions) + self.helmholtz_energy(contributions)
@@ -199,11 +171,6 @@ impl<E: Residual + IdealGas> State<E> {
199171
self.internal_energy(contributions) / self.total_moles
200172
}
201173

202-
/// Specific internal energy: $u^{(m)}=\frac{U}{m}$
203-
pub fn specific_internal_energy(&self, contributions: Contributions) -> SpecificEnergy {
204-
self.molar_internal_energy(contributions) / self.total_molar_weight()
205-
}
206-
207174
/// Gibbs energy: $G=A+pV$
208175
pub fn gibbs_energy(&self, contributions: Contributions) -> Energy {
209176
self.pressure(contributions) * self.volume + self.helmholtz_energy(contributions)
@@ -214,11 +181,6 @@ impl<E: Residual + IdealGas> State<E> {
214181
self.gibbs_energy(contributions) / self.total_moles
215182
}
216183

217-
/// Specific Gibbs energy: $g^{(m)}=\frac{G}{m}$
218-
pub fn specific_gibbs_energy(&self, contributions: Contributions) -> SpecificEnergy {
219-
self.molar_gibbs_energy(contributions) / self.total_molar_weight()
220-
}
221-
222184
/// Joule Thomson coefficient: $\mu_{JT}=\left(\frac{\partial T}{\partial p}\right)_{H,N_i}$
223185
pub fn joule_thomson(&self) -> <Temperature as Div<Pressure>>::Output {
224186
let c = Contributions::Total;
@@ -278,6 +240,46 @@ impl<E: Residual + IdealGas> State<E> {
278240
}
279241
res
280242
}
243+
}
244+
245+
impl<E: Residual + Molarweight + IdealGas> State<E> {
246+
/// Specific isochoric heat capacity: $c_v^{(m)}=\frac{C_v}{m}$
247+
pub fn specific_isochoric_heat_capacity(
248+
&self,
249+
contributions: Contributions,
250+
) -> SpecificEntropy {
251+
self.molar_isochoric_heat_capacity(contributions) / self.total_molar_weight()
252+
}
253+
254+
/// Specific isobaric heat capacity: $c_p^{(m)}=\frac{C_p}{m}$
255+
pub fn specific_isobaric_heat_capacity(&self, contributions: Contributions) -> SpecificEntropy {
256+
self.molar_isobaric_heat_capacity(contributions) / self.total_molar_weight()
257+
}
258+
259+
/// Specific entropy: $s^{(m)}=\frac{S}{m}$
260+
pub fn specific_entropy(&self, contributions: Contributions) -> SpecificEntropy {
261+
self.molar_entropy(contributions) / self.total_molar_weight()
262+
}
263+
264+
/// Specific enthalpy: $h^{(m)}=\frac{H}{m}$
265+
pub fn specific_enthalpy(&self, contributions: Contributions) -> SpecificEnergy {
266+
self.molar_enthalpy(contributions) / self.total_molar_weight()
267+
}
268+
269+
/// Specific Helmholtz energy: $a^{(m)}=\frac{A}{m}$
270+
pub fn specific_helmholtz_energy(&self, contributions: Contributions) -> SpecificEnergy {
271+
self.molar_helmholtz_energy(contributions) / self.total_molar_weight()
272+
}
273+
274+
/// Specific internal energy: $u^{(m)}=\frac{U}{m}$
275+
pub fn specific_internal_energy(&self, contributions: Contributions) -> SpecificEnergy {
276+
self.molar_internal_energy(contributions) / self.total_molar_weight()
277+
}
278+
279+
/// Specific Gibbs energy: $g^{(m)}=\frac{G}{m}$
280+
pub fn specific_gibbs_energy(&self, contributions: Contributions) -> SpecificEnergy {
281+
self.molar_gibbs_energy(contributions) / self.total_molar_weight()
282+
}
281283

282284
/// Speed of sound: $c=\sqrt{\left(\frac{\partial p}{\partial\rho^{(m)}}\right)_{S,N_i}}$
283285
pub fn speed_of_sound(&self) -> Velocity {

feos-core/src/state/residual_properties.rs

+3-1
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
11
use super::{Contributions, Derivative::*, PartialDerivative, State};
2-
use crate::equation_of_state::{EntropyScaling, Residual};
2+
use crate::equation_of_state::{EntropyScaling, Molarweight, Residual};
33
use crate::errors::EosResult;
44
use crate::phase_equilibria::PhaseEquilibrium;
55
use crate::ReferenceSystem;
@@ -484,7 +484,9 @@ impl<E: Residual> State<E> {
484484
pub fn residual_molar_gibbs_energy(&self) -> MolarEnergy {
485485
self.residual_gibbs_energy() / self.total_moles
486486
}
487+
}
487488

489+
impl<E: Residual + Molarweight> State<E> {
488490
/// Total molar weight: $MW=\sum_ix_iMW_i$
489491
pub fn total_molar_weight(&self) -> MolarWeight {
490492
(self.eos.molar_weight() * Dimensionless::new(&self.molefracs)).sum()

0 commit comments

Comments
 (0)