Skip to content
Merged
4 changes: 4 additions & 0 deletions Cargo.toml
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,10 @@ criterion = "0.4"
name = "state_properties"
harness = false

[[bench]]
name = "state_creation"
harness = false

[[bench]]
name = "dual_numbers"
harness = false
Expand Down
3 changes: 2 additions & 1 deletion benches/README.md
Original file line number Diff line number Diff line change
Expand Up @@ -13,4 +13,5 @@ cargo bench --profile=release-lto --features=pcsaft --bench=dual_numbers
|Name|Description|Cargo features|
|--|--|--|
|`dual_numbers`|Helmholtz energy function evaluated using `StateHD` with different dual number types.|`pcsaft`|
|`state_properties`|Properties of `State`. Including state creation using the natural variables of the Helmholtz energy (no density iteration).|`pcsaft`|
|`state_properties`|Properties of `State`. Including state creation using the natural variables of the Helmholtz energy (no density iteration).|`pcsaft`|
|`state_creation`|Different constructors of `State` and `PhaseEquilibrium` including critical point calculations. For pure substances and mixtures.|`pcsaft`|
13 changes: 8 additions & 5 deletions benches/dual_numbers.rs
Original file line number Diff line number Diff line change
Expand Up @@ -48,8 +48,11 @@ fn bench_dual_numbers<E: EquationOfState>(
group.bench_function("a_dual", |b| {
b.iter(|| a_res((&state.eos, &state.derive1(Derivative::DV))))
});
group.bench_function("a_dual2", |b| {
b.iter(|| a_res((&state.eos, &state.derive2(Derivative::DV))))
});
group.bench_function("a_hyperdual", |b| {
b.iter(|| a_res((&state.eos, &state.derive2(Derivative::DV, Derivative::DV))))
b.iter(|| a_res((&state.eos, &state.derive2_mixed(Derivative::DV, Derivative::DV))))
Comment thread
g-bauer marked this conversation as resolved.
});
group.bench_function("a_dual3", |b| {
b.iter(|| a_res((&state.eos, &state.derive3(Derivative::DV))))
Expand All @@ -66,7 +69,7 @@ fn pcsaft(c: &mut Criterion) {
IdentifierOption::Name,
)
.unwrap();
bench_dual_numbers(c, "methane", state_pcsaft(parameters));
bench_dual_numbers(c, "dual_numbers_pcsaft_methane", state_pcsaft(parameters));

// water (4C, polar)
let parameters = PcSaftParameters::from_json(
Expand All @@ -76,7 +79,7 @@ fn pcsaft(c: &mut Criterion) {
IdentifierOption::Name,
)
.unwrap();
bench_dual_numbers(c, "water_4c_polar", state_pcsaft(parameters));
bench_dual_numbers(c, "dual_numbers_pcsaft_water_4c_polar", state_pcsaft(parameters));

// methane, ethane, propane
let parameters = PcSaftParameters::from_json(
Expand All @@ -86,7 +89,7 @@ fn pcsaft(c: &mut Criterion) {
IdentifierOption::Name,
)
.unwrap();
bench_dual_numbers(c, "methane_ethane_propane", state_pcsaft(parameters));
bench_dual_numbers(c, "dual_numbers_pcsaft_methane_ethane_propane", state_pcsaft(parameters));
}

/// Benchmark for the PC-SAFT equation of state.
Expand Down Expand Up @@ -116,7 +119,7 @@ fn methane_co2_pcsaft(c: &mut Criterion) {
let x = arr1(&[0.15, 0.85]);
let moles = &x * 10.0 * MOL;
let state = State::new_nvt(&eos, temperature, volume, &moles).unwrap();
bench_dual_numbers(c, "methane_co2", state);
bench_dual_numbers(c, "dual_numbers_pcsaft_methane_co2", state);
}

criterion_group!(bench, pcsaft, methane_co2_pcsaft);
Expand Down
171 changes: 171 additions & 0 deletions benches/state_creation.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,171 @@
use criterion::{criterion_group, criterion_main, Criterion};
use feos::pcsaft::{PcSaft, PcSaftParameters};
use feos_core::{
parameter::{IdentifierOption, Parameter},
Contributions, DensityInitialization, EquationOfState, PhaseEquilibrium, State,
};
use ndarray::{Array, Array1};
use quantity::si::*;
use std::sync::Arc;

/// Evaluate NPT constructor
fn npt<E: EquationOfState>(
(eos, t, p, n, rho0): (
&Arc<E>,
SINumber,
SINumber,
&SIArray1,
DensityInitialization<SIUnit>,
),
) {
State::new_npt(eos, t, p, n, rho0).unwrap();
}

/// Evaluate critical point constructor
fn critical_point<E: EquationOfState>((eos, n): (&Arc<E>, Option<&SIArray1>)) {
State::critical_point(eos, n, None, Default::default()).unwrap();
}

/// VLE for pure substance for given temperature or pressure
fn pure<E: EquationOfState>((eos, t_or_p): (&Arc<E>, SINumber)) {
PhaseEquilibrium::pure(eos, t_or_p, None, Default::default()).unwrap();
}

/// Evaluate temperature, pressure flash.
fn tp_flash<E: EquationOfState>((eos, t, p, feed): (&Arc<E>, SINumber, SINumber, &SIArray1)) {
PhaseEquilibrium::tp_flash(eos, t, p, feed, None, Default::default(), None).unwrap();
}

fn bubble_point<E: EquationOfState>((eos, t, x): (&Arc<E>, SINumber, &Array1<f64>)) {
PhaseEquilibrium::bubble_point(
eos,
t,
x,
None,
None,
(Default::default(), Default::default()),
)
.unwrap();
}

fn dew_point<E: EquationOfState>((eos, t, y): (&Arc<E>, SINumber, &Array1<f64>)) {
PhaseEquilibrium::dew_point(
eos,
t,
y,
None,
None,
(Default::default(), Default::default()),
)
.unwrap();
}

fn bench_states<E: EquationOfState>(c: &mut Criterion, group_name: &str, eos: &Arc<E>) {
let ncomponents = eos.components();
let x = Array::from_elem(ncomponents, 1.0 / ncomponents as f64);
let n = &x * 100.0 * MOL;
let crit = State::critical_point(&eos, Some(&n), None, Default::default()).unwrap();
let vle = if ncomponents == 1 {
PhaseEquilibrium::pure(&eos, crit.temperature * 0.95, None, Default::default()).unwrap()
} else {
PhaseEquilibrium::tp_flash(
&eos,
crit.temperature,
crit.pressure(Contributions::Total) * 0.95,
&crit.moles,
None,
Default::default(),
None,
)
.unwrap()
};

let mut group = c.benchmark_group(group_name);
group.bench_function("new_npt_liquid", |b| {
b.iter(|| {
npt((
&eos,
vle.liquid().temperature,
vle.liquid().pressure(Contributions::Total) * 1.01,
&n,
DensityInitialization::Liquid,
))
})
});
group.bench_function("new_npt_vapor", |b| {
b.iter(|| {
npt((
&eos,
vle.vapor().temperature,
vle.vapor().pressure(Contributions::Total) * 0.99,
&n,
DensityInitialization::Vapor,
))
})
});
group.bench_function("critical_point", |b| {
b.iter(|| critical_point((&eos, Some(&n))))
});
if ncomponents != 1 {
group.bench_function("tp_flash", |b| {
b.iter(|| {
tp_flash((
&eos,
crit.temperature,
crit.pressure(Contributions::Total) * 0.99,
&n,
))
})
});

group.bench_function("bubble_point", |b| {
b.iter(|| bubble_point((&eos, vle.liquid().temperature, &vle.liquid().molefracs)))
});

group.bench_function("dew_point", |b| {
b.iter(|| dew_point((&eos, vle.vapor().temperature, &vle.vapor().molefracs)))
});
} else {
group.bench_function("pure_t", |b| {
b.iter(|| pure((&eos, vle.vapor().temperature)))
});
group.bench_function("pure_p", |b| {
b.iter(|| pure((&eos, vle.vapor().pressure(Contributions::Total))))
});
}
}

fn pcsaft(c: &mut Criterion) {
let parameters = PcSaftParameters::from_json(
vec!["methane"],
"./parameters/pcsaft/gross2001.json",
None,
IdentifierOption::Name,
)
.unwrap();
let eos = Arc::new(PcSaft::new(Arc::new(parameters)));
bench_states(c, "state_creation_pcsaft_methane", &eos);

let parameters = PcSaftParameters::from_json(
vec!["methane", "ethane", "propane"],
"./parameters/pcsaft/gross2001.json",
None,
IdentifierOption::Name,
)
.unwrap();
let eos = Arc::new(PcSaft::new(Arc::new(parameters)));
bench_states(c, "state_creation_pcsaft_methane_ethane", &eos);

let parameters = PcSaftParameters::from_json(
vec!["methane", "ethane", "propane"],
"./parameters/pcsaft/gross2001.json",
None,
IdentifierOption::Name,
)
.unwrap();
let eos = Arc::new(PcSaft::new(Arc::new(parameters)));
bench_states(c, "state_creation_pcsaft_methane_ethane_propane", &eos);
}

criterion_group!(bench, pcsaft);
criterion_main!(bench);
2 changes: 1 addition & 1 deletion benches/state_properties.rs
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,7 @@ fn properties_pcsaft(c: &mut Criterion) {
let x = arr1(&[1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0]);
let m = &x * 100.0 * MOL;

let mut group = c.benchmark_group("methane_ethane_propane");
let mut group = c.benchmark_group("state_properties_pcsaft_methane_ethane_propane");
group.bench_function("a", |b| {
b.iter(|| property((&eos, S::helmholtz_energy, t, v, &m, Contributions::Total)))
});
Expand Down
5 changes: 5 additions & 0 deletions feos-core/CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,11 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0

### Changed
- Added `Sync` and `Send` as supertraits to `EquationOfState`. [#57](https://github.com/feos-org/feos/pull/57)
- Added `Dual2_64` dual number to `HelmholtzEnergy` trait to facilitate efficient non-mixed second order derivatives. [#94](https://github.com/feos-org/feos/pull/94)
- Renamed `PartialDerivative::Second` to `PartialDerivative::SecondMixed`. [#94](https://github.com/feos-org/feos/pull/94)
- Added `PartialDerivative::Second` to enable non-mixed second order partial derivatives. [#94](https://github.com/feos-org/feos/pull/94)
- Changed `dp_dv_` and `ds_dt_` to use `Dual2_64` instead of `HyperDual64`. [#94](https://github.com/feos-org/feos/pull/94)
- Added `get_or_insert_with_d2_64` to `Cache`. [#94](https://github.com/feos-org/feos/pull/94)

## [0.3.1] - 2022-08-25
### Added
Expand Down
6 changes: 5 additions & 1 deletion feos-core/src/equation_of_state.rs
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ use crate::errors::{EosError, EosResult};
use crate::state::StateHD;
use crate::EosUnit;
use ndarray::prelude::*;
use num_dual::{Dual, Dual3, Dual3_64, Dual64, DualNum, DualVec64, HyperDual, HyperDual64};
use num_dual::{Dual, Dual3, Dual3_64, Dual64, DualNum, DualVec64, HyperDual, HyperDual64, Dual2_64};
use num_traits::{One, Zero};
use quantity::{QuantityArray1, QuantityScalar};
use std::fmt;
Expand All @@ -28,6 +28,7 @@ pub trait HelmholtzEnergy:
+ HelmholtzEnergyDual<Dual64>
+ HelmholtzEnergyDual<Dual<DualVec64<3>, f64>>
+ HelmholtzEnergyDual<HyperDual64>
+ HelmholtzEnergyDual<Dual2_64>
+ HelmholtzEnergyDual<Dual3_64>
+ HelmholtzEnergyDual<HyperDual<Dual64, f64>>
+ HelmholtzEnergyDual<HyperDual<DualVec64<2>, f64>>
Expand All @@ -46,6 +47,7 @@ impl<T> HelmholtzEnergy for T where
+ HelmholtzEnergyDual<Dual64>
+ HelmholtzEnergyDual<Dual<DualVec64<3>, f64>>
+ HelmholtzEnergyDual<HyperDual64>
+ HelmholtzEnergyDual<Dual2_64>
+ HelmholtzEnergyDual<Dual3_64>
+ HelmholtzEnergyDual<HyperDual<Dual64, f64>>
+ HelmholtzEnergyDual<HyperDual<DualVec64<2>, f64>>
Expand Down Expand Up @@ -99,6 +101,7 @@ pub trait IdealGasContribution:
+ IdealGasContributionDual<Dual64>
+ IdealGasContributionDual<Dual<DualVec64<3>, f64>>
+ IdealGasContributionDual<HyperDual64>
+ IdealGasContributionDual<Dual2_64>
+ IdealGasContributionDual<Dual3_64>
+ IdealGasContributionDual<HyperDual<Dual64, f64>>
+ IdealGasContributionDual<HyperDual<DualVec64<2>, f64>>
Expand All @@ -115,6 +118,7 @@ impl<T> IdealGasContribution for T where
+ IdealGasContributionDual<Dual64>
+ IdealGasContributionDual<Dual<DualVec64<3>, f64>>
+ IdealGasContributionDual<HyperDual64>
+ IdealGasContributionDual<Dual2_64>
+ IdealGasContributionDual<Dual3_64>
+ IdealGasContributionDual<HyperDual<Dual64, f64>>
+ IdealGasContributionDual<HyperDual<DualVec64<2>, f64>>
Expand Down
17 changes: 15 additions & 2 deletions feos-core/src/python/statehd.rs
Original file line number Diff line number Diff line change
@@ -1,7 +1,9 @@
use crate::StateHD;
use ndarray::Array1;
use num_dual::python::{PyDual3Dual64, PyDual3_64, PyDual64, PyHyperDual64, PyHyperDualDual64};
use num_dual::{Dual, Dual3, Dual3_64, Dual64, DualVec64, HyperDual, HyperDual64};
use num_dual::python::{
PyDual2_64, PyDual3Dual64, PyDual3_64, PyDual64, PyHyperDual64, PyHyperDualDual64,
};
use num_dual::{Dual, Dual2_64, Dual3, Dual3_64, Dual64, DualVec64, HyperDual, HyperDual64};
use pyo3::prelude::*;
use std::convert::From;

Expand Down Expand Up @@ -75,6 +77,16 @@ impl From<StateHD<Dual<DualVec64<3>, f64>>> for PyStateDDV3 {
}
}

#[pyclass(name = "StateD2")]
#[derive(Clone)]
pub struct PyStateD2(StateHD<Dual2_64>);

impl From<StateHD<Dual2_64>> for PyStateD2 {
fn from(s: StateHD<Dual2_64>) -> Self {
Self(s)
}
}

#[pyclass(name = "StateD3")]
#[derive(Clone)]
pub struct PyStateD3(StateHD<Dual3_64>);
Expand Down Expand Up @@ -167,6 +179,7 @@ impl_state_hd!(
HyperDual<Dual64, f64>
);
impl_state_hd!(PyStateD, PyDual64, StateHD<Dual64>, Dual64);
impl_state_hd!(PyStateD2, PyDual2_64, StateHD<Dual2_64>, Dual2_64);
impl_state_hd!(PyStateD3, PyDual3_64, StateHD<Dual3_64>, Dual3_64);
impl_state_hd!(
PyStateD3D,
Expand Down
7 changes: 5 additions & 2 deletions feos-core/src/python/user_defined.rs
Original file line number Diff line number Diff line change
@@ -1,8 +1,10 @@
use crate::python::statehd::*;
use crate::*;
use ndarray::prelude::*;
use num_dual::python::{PyDual3Dual64, PyDual3_64, PyDual64, PyHyperDual64, PyHyperDualDual64};
use num_dual::{Dual, Dual3, Dual3_64, Dual64, DualVec64, HyperDual, HyperDual64};
use num_dual::python::{
PyDual2_64, PyDual3Dual64, PyDual3_64, PyDual64, PyHyperDual64, PyHyperDualDual64,
};
use num_dual::{Dual, Dual2_64, Dual3, Dual3_64, Dual64, DualVec64, HyperDual, HyperDual64};
use numpy::convert::IntoPyArray;
use pyo3::prelude::*;
use quantity::python::PySIArray1;
Expand Down Expand Up @@ -215,6 +217,7 @@ impl_helmholtz_energy!(
PyHyperDualDualVec64_3,
HyperDual<DualVec64<3>, f64>
);
impl_helmholtz_energy!(PyStateD2, PyDual2_64, Dual2_64);
impl_helmholtz_energy!(PyStateD3, PyDual3_64, Dual3_64);
impl_helmholtz_energy!(PyStateD3D, PyDual3Dual64, Dual3<Dual64, f64>);
impl_helmholtz_energy!(PyStateD3DV2, PyDual3DualVec64_2, Dual3<DualVec64<2>, f64>);
Expand Down
Loading