From 3b256bf22a2289640b6c51e4ac1eae8c860364b1 Mon Sep 17 00:00:00 2001 From: Gernot Bauer Date: Wed, 7 Dec 2022 16:50:31 +0100 Subject: [PATCH 1/7] add benchmarks --- Cargo.toml | 9 ++ benches/benchmarks.rs | 287 ++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 296 insertions(+) create mode 100644 benches/benchmarks.rs diff --git a/Cargo.toml b/Cargo.toml index aa5ccd923..f7159f935 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -46,6 +46,15 @@ optional = true [dev-dependencies] approx = "0.4" +criterion = "*" + +[[bench]] +name = "benchmarks" +harness = false + +[profile.release-lto] +inherits = "release" +lto = true [features] default = [] diff --git a/benches/benchmarks.rs b/benches/benchmarks.rs new file mode 100644 index 000000000..7234f6cd5 --- /dev/null +++ b/benches/benchmarks.rs @@ -0,0 +1,287 @@ +use criterion::{criterion_group, criterion_main, Criterion}; +use feos::pcsaft::{PcSaft, PcSaftParameters}; +use feos_core::{ + parameter::{IdentifierOption, Parameter}, + Contributions, EquationOfState, HelmholtzEnergy, HelmholtzEnergyDual, State, StateBuilder, + StateHD, +}; +use ndarray::{arr1, Array1}; +use num_dual::{Dual64, DualNum, HyperDual64, StaticMat, StaticVec}; +use quantity::{si::*, QuantityArray1}; +use std::sync::Arc; + +/// State generation using `State::new` +fn state(inp: (&Arc, SINumber, SINumber, &Array1)) -> State { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() +} + +/// Residual Helmholtz energy given a StateHD +fn a_res(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { + let t = inp.1.to_reduced(KELVIN).unwrap(); + let v = inp.2.to_reduced(ANGSTROM.powi(3)).unwrap(); + let m = inp.3.to_reduced(MOL).unwrap(); + let s = StateHD::new(t, v, m); + inp.0.evaluate_residual(&s) +} + +/// Residual Helmholtz energy given a StateHD +fn a_res_shd>(inp: (&Arc, &StateHD)) -> D +where + (dyn HelmholtzEnergy + 'static): HelmholtzEnergyDual, +{ + inp.0.evaluate_residual(inp.1) +} + +fn da_dv(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { + let t = Dual64::from_re(inp.1.to_reduced(KELVIN).unwrap()); + let v = Dual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(); + let m = inp + .3 + .to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual64::from_re(mi)) + .collect(); + let s = StateHD::new(t, v, m); + inp.0.evaluate_residual(&s).eps[0] +} + +fn compressibility(inp: (&Arc, SINumber, SINumber, &Array1)) -> f64 { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() + .compressibility(Contributions::Total) +} + +fn compressibility_(state: &State) -> f64 { + state.compressibility(Contributions::Total) +} + +fn ln_phi(inp: (&Arc, SINumber, SINumber, &Array1)) -> Array1 { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() + .ln_phi() +} + +fn da_dni(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> Array1 { + let t = Dual64::from_re(inp.1.to_reduced(KELVIN).unwrap()); + let v = Dual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()); + let m: Array1<_> = inp + .3 + .to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual64::from_re(mi)) + .collect(); + Array1::from_shape_fn(inp.3.len(), |i| { + let mut n = m.clone(); + n[i] = n[i].derive(); + let s = StateHD::new(t, v, n); + inp.0.evaluate_residual(&s).eps[0] + }) +} + +fn c_v(inp: (&Arc, SINumber, SINumber, &Array1)) -> SINumber { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() + .c_v(Contributions::ResidualNvt) +} + +fn d2a_dt2(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { + let t = HyperDual64::from_re(inp.1.to_reduced(KELVIN).unwrap()).derive2(); + let v = HyperDual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()); + let m: Array1<_> = inp + .3 + .to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| HyperDual64::from_re(mi)) + .collect(); + + let s = StateHD::new(t, v, m); + inp.0.evaluate_residual(&s).eps1eps2[0] +} + +fn molar_volume(inp: (&Arc, SINumber, SINumber, &Array1)) -> SIArray1 { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() + .molar_volume(Contributions::ResidualNvt) +} + +fn benchmark_helmholtz_energy(c: &mut Criterion) { + 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))); + let t = 300.0 * KELVIN; + let density = 71.18 * KILO * MOL / METER.powi(3); + let volume = 100.0 * MOL / density; + let x = arr1(&[1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0]); + let rho_i = &x * density; + let m = &x * 100.0 * MOL; + let s = StateHD::new( + 300.0, + volume.to_reduced(ANGSTROM.powi(3)).unwrap(), + &x * 100.0, + ); + + let mut group = c.benchmark_group("helmholtz_energy"); + group.bench_function("state", |b| b.iter(|| state((&eos, t, density, &x)))); + group.bench_function("a_res", |b| b.iter(|| a_res((&eos, t, volume, &m)))); + group.bench_function("a_res_shd", |b| b.iter(|| a_res_shd((&eos, &s)))); +} + +fn benchmark_first_derivative(c: &mut Criterion) { + 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))); + let t = 300.0 * KELVIN; + let density = 71.18 * KILO * MOL / METER.powi(3); + let volume = 100.0 * MOL / density; + let x = arr1(&[1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0]); + let rho_i = &x * density; + let m = &x * 100.0 * MOL; + + let s = StateBuilder::new(&eos) + .temperature(t) + .partial_density(&rho_i) + .total_moles(m.sum()) + .build() + .unwrap(); + + let mut dual = c.benchmark_group("first_derivative"); + // group.bench_function("state", |b| b.iter(|| state((&eos, t, density, &x)))); + dual.bench_function("compressibility", |b| { + b.iter(|| compressibility((&eos, t, density, &x))) + }); + dual.bench_function("compressibility_c", |b| { + b.iter(|| compressibility((&eos, t, density, &x))) + }); + dual.bench_function("da_dv", |b| b.iter(|| da_dv((&eos, t, volume, &m)))); + dual.bench_function("ln_phi", |b| b.iter(|| ln_phi((&eos, t, density, &x)))); + dual.bench_function("da_dni", |b| b.iter(|| da_dni((&eos, t, volume, &m)))); +} + +fn benchmark_second_derivative(c: &mut Criterion) { + 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))); + let t = 300.0 * KELVIN; + let density = 71.18 * KILO * MOL / METER.powi(3); + let volume = 100.0 * MOL / density; + let x = arr1(&[1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0]); + let rho_i = &x * density; + let m = &x * 100.0 * MOL; + + let s = StateBuilder::new(&eos) + .temperature(t) + .partial_density(&rho_i) + .total_moles(m.sum()) + .build() + .unwrap(); + + let mut hyper_dual = c.benchmark_group("second_derivative"); + hyper_dual.bench_function("c_v", |b| b.iter(|| c_v((&eos, t, density, &x)))); + hyper_dual.bench_function("d2a_dt2", |b| b.iter(|| d2a_dt2((&eos, t, volume, &m)))); + + hyper_dual.bench_function("molar_volume", |b| { + b.iter(|| molar_volume((&eos, t, density, &x))) + }); +} + +criterion_group!( + bench, + benchmark_helmholtz_energy, + benchmark_first_derivative, + benchmark_second_derivative +); +criterion_main!(bench); From 88fd7f622371a089ac8bef59827db196ab9451f2 Mon Sep 17 00:00:00 2001 From: Gernot Bauer Date: Wed, 7 Dec 2022 18:11:45 +0100 Subject: [PATCH 2/7] Add benchmarks for binary mixture of methane/co2 --- Cargo.toml | 4 ++ benches/biogas.rs | 162 ++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 166 insertions(+) create mode 100644 benches/biogas.rs diff --git a/Cargo.toml b/Cargo.toml index f7159f935..8850d32c8 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -52,6 +52,10 @@ criterion = "*" name = "benchmarks" harness = false +[[bench]] +name = "biogas" +harness = false + [profile.release-lto] inherits = "release" lto = true diff --git a/benches/biogas.rs b/benches/biogas.rs new file mode 100644 index 000000000..445dee300 --- /dev/null +++ b/benches/biogas.rs @@ -0,0 +1,162 @@ +use criterion::{criterion_group, criterion_main, Criterion}; +use feos::pcsaft::{PcSaft, PcSaftParameters}; +use feos_core::{ + parameter::{IdentifierOption, Parameter}, + Contributions, EquationOfState, HelmholtzEnergy, HelmholtzEnergyDual, State, StateBuilder, + StateHD, +}; +use ndarray::{arr1, Array1}; +use num_dual::{Dual3, Dual64, DualNum, HyperDual64, StaticMat, StaticVec}; +use quantity::{si::*, QuantityArray1}; +use std::sync::Arc; + +/// State generation using `State::new` +fn state(inp: (&Arc, SINumber, SINumber, &Array1)) -> State { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() +} + +/// Residual Helmholtz energy given a StateHD +fn a_res>(inp: (&Arc, &StateHD)) -> D +where + (dyn HelmholtzEnergy + 'static): HelmholtzEnergyDual, +{ + inp.0.evaluate_residual(inp.1) +} + +fn da_dv(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { + let t = Dual64::from_re(inp.1.to_reduced(KELVIN).unwrap()); + let v = Dual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(); + let m = inp + .3 + .to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual64::from_re(mi)) + .collect(); + let s = StateHD::new(t, v, m); + inp.0.evaluate_residual(&s).eps[0] +} + +fn da_dni(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> Array1 { + let t = Dual64::from_re(inp.1.to_reduced(KELVIN).unwrap()); + let v = Dual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()); + let m: Array1<_> = inp + .3 + .to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual64::from_re(mi)) + .collect(); + Array1::from_shape_fn(inp.3.len(), |i| { + let mut n = m.clone(); + n[i] = n[i].derive(); + let s = StateHD::new(t, v, n); + inp.0.evaluate_residual(&s).eps[0] + }) +} + +fn d2a_dt2(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { + let t = HyperDual64::from_re(inp.1.to_reduced(KELVIN).unwrap()).derive2(); + let v = HyperDual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()); + let m: Array1<_> = inp + .3 + .to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| HyperDual64::from_re(mi)) + .collect(); + + let s = StateHD::new(t, v, m); + inp.0.evaluate_residual(&s).eps1eps2[0] +} + +fn benchmark_helmholtz_energy(c: &mut Criterion) { + let parameters = PcSaftParameters::from_multiple_json( + &[ + (vec!["methane"], "./parameters/pcsaft/gross2001.json"), + ( + vec!["carbon dioxide"], + "./parameters/pcsaft/gross2005_fit.json", + ), + ], + None, + IdentifierOption::Name, + ) + .unwrap(); + let k_ij = -0.0192211646; + let parameters = + PcSaftParameters::new_binary(parameters.pure_records.clone(), Some(k_ij.into())); + let eos = Arc::new(PcSaft::new(Arc::new(parameters))); + + // 230 K, 50 bar, x0 = 0.15 + let t = 230.0 * KELVIN; + let density = 24.16896 * KILO * MOL / METER.powi(3); + let volume = 10.0 * MOL / density; + let x = arr1(&[0.15, 0.85]); + let m = &x * 10.0 * MOL; + + let mut group = c.benchmark_group("helmholtz_energy"); + + // real valued evaluation + let s = StateHD::new( + t.to_reduced(KELVIN).unwrap(), + volume.to_reduced(ANGSTROM.powi(3)).unwrap(), + m.to_reduced(MOL).unwrap(), + ); + group.bench_function("a_res_f64", |b| b.iter(|| a_res((&eos, &s)))); + + // da_dv - dual number + let s = StateHD::new( + Dual64::from_re(t.to_reduced(KELVIN).unwrap()), + Dual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual64::from_re(mi)) + .collect(), + ); + group.bench_function("a_res_dual", |b| b.iter(|| a_res((&eos, &s)))); + + // d2a_dv2 - hyperdual number + let s = StateHD::new( + HyperDual64::from_re(t.to_reduced(KELVIN).unwrap()).derive1(), + HyperDual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive2(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| HyperDual64::from_re(mi)) + .collect(), + ); + group.bench_function("a_res_hyperdual", |b| b.iter(|| a_res((&eos, &s)))); + + // d3a_dv3 - dual3 number + let s = StateHD::new( + Dual3::from_re(t.to_reduced(KELVIN).unwrap()), + Dual3::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual3::from_re(mi)) + .collect(), + ); + group.bench_function("a_res_dual3", |b| b.iter(|| a_res((&eos, &s)))); +} + +criterion_group!(bench, benchmark_helmholtz_energy,); +criterion_main!(bench); From 71c730b867bec0754117f5b90a45eb4783902600 Mon Sep 17 00:00:00 2001 From: Gernot Bauer Date: Tue, 13 Dec 2022 17:37:13 +0100 Subject: [PATCH 3/7] Renamed and restructured benchmarks, added README for benchmarks, added Cargo profile for LTO --- CHANGELOG.md | 1 + Cargo.toml | 6 +- README.md | 29 +++- benches/README.md | 16 ++ benches/benchmarks.rs | 287 ------------------------------------ benches/biogas.rs | 162 -------------------- benches/dual_numbers.rs | 103 +++++++++++++ benches/state_properties.rs | 144 ++++++++++++++++++ 8 files changed, 289 insertions(+), 459 deletions(-) create mode 100644 benches/README.md delete mode 100644 benches/benchmarks.rs delete mode 100644 benches/biogas.rs create mode 100644 benches/dual_numbers.rs create mode 100644 benches/state_properties.rs diff --git a/CHANGELOG.md b/CHANGELOG.md index 24f824759..50d9ce071 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,6 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added - Added SAFT-VRQ Mie equation of state and Helmholtz energy functional for first order Feynman-Hibbs corrected Mie fluids. [#79](https://github.com/feos-org/feos/pull/79) - Added `estimator` module to documentation. [#86](https://github.com/feos-org/feos/pull/86) +- Added benchmarks for the evaluation of the Helmholtz energy and some properties of the `State` object for PC-SAFT. [#89](https://github.com/feos-org/feos/pull/89) ### Changed - Export `EosVariant` and `FunctionalVariant` directly in the crate root instead of their own modules. [#62](https://github.com/feos-org/feos/pull/62) diff --git a/Cargo.toml b/Cargo.toml index 8850d32c8..58f415ba0 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -46,14 +46,14 @@ optional = true [dev-dependencies] approx = "0.4" -criterion = "*" +criterion = "0.4" [[bench]] -name = "benchmarks" +name = "state_properties" harness = false [[bench]] -name = "biogas" +name = "dual_numbers" harness = false [profile.release-lto] diff --git a/README.md b/README.md index 9e561eb04..0bd621457 100644 --- a/README.md +++ b/README.md @@ -96,37 +96,52 @@ If there is no compiled package for your system available from PyPI and you have pip install git+https://github.com/feos-org/feos ``` +This command builds the package without link-time optimization (LTO) which results in considerably worse performance. +See the *Building from source* section for information about building the wheel including LTO. + ### Building from source To compile the code you need the Rust compiler and `maturin` (>=0.13,<0.14) installed. -To install the package directly into the active environment, use +To install the package directly into the active environment (virtualenv or conda), use ``` -maturin develop --release --features python +maturin develop --release ``` -and specify the models that you want to include in the python package as additional features, e.g. +which uses the `python` and `all_models` feature as specified in the `pyproject.toml` file. + +Alternatively, you can specify the models or features that you want to include in the python package explicitly, e.g. ``` maturin develop --release --features "python pcsaft dft" ``` -for the PC-SAFT equation of state and Helmholtz energy functional. If you want to include all available models, use +for the PC-SAFT equation of state and Helmholtz energy functional. + +To build wheels including link-time optimization (LTO), use ``` -maturin develop --release --features "python all_models" +maturin build --profile="release-lto" ``` -To build wheels, use +which will use the `python` and `all_models` features specified in the `pyproject.toml` file. +Use the following command to build a wheel with specific features: ``` -maturin build --release --out dist --features "python ..." +maturin build --profile="release-lto" --features "python ..." ``` +LTO increases compile times drastically but the resulting wheel is way more performant and has a smaller size. +For development however, we recommend using the `--release` flag. + ## Documentation For a documentation of the Python API, Python examples, and a guide to the underlying Rust framework check out the [documentation](https://feos-org.github.io/feos/). +## Benchmarks + +Check out the [benches](https://github.com/feos-org/feos/tree/main/benches) directory for information about provided Rust benchmarks and how to run them. + ## Developers This software is currently maintained by members of the groups of diff --git a/benches/README.md b/benches/README.md new file mode 100644 index 000000000..8f146baae --- /dev/null +++ b/benches/README.md @@ -0,0 +1,16 @@ +# Benchmarks + +This directory contains different benchmarks. +For best performance, use the `release-lto` profile. +Depending on the benchmark, you might have to consider different Cargo `features` as denoted in the table below. + +For example, to run the benchmarks in `dual_numbers`, which uses PC-SAFT, use + +``` +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. Substances: methane + carbon dioxide. Model: PC-SAFT incl. dipole-dipole contributions.|`pcsaft`| +|`state_properties`|Properties of `State`. Including state creation using the natural variables of the Helmholtz energy (no density iteration). Substances: methane + ethane + propane. Model: PC-SAFT (hard-sphere, hard-chain, dispersion).|`pcsaft`| \ No newline at end of file diff --git a/benches/benchmarks.rs b/benches/benchmarks.rs deleted file mode 100644 index 7234f6cd5..000000000 --- a/benches/benchmarks.rs +++ /dev/null @@ -1,287 +0,0 @@ -use criterion::{criterion_group, criterion_main, Criterion}; -use feos::pcsaft::{PcSaft, PcSaftParameters}; -use feos_core::{ - parameter::{IdentifierOption, Parameter}, - Contributions, EquationOfState, HelmholtzEnergy, HelmholtzEnergyDual, State, StateBuilder, - StateHD, -}; -use ndarray::{arr1, Array1}; -use num_dual::{Dual64, DualNum, HyperDual64, StaticMat, StaticVec}; -use quantity::{si::*, QuantityArray1}; -use std::sync::Arc; - -/// State generation using `State::new` -fn state(inp: (&Arc, SINumber, SINumber, &Array1)) -> State { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() -} - -/// Residual Helmholtz energy given a StateHD -fn a_res(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { - let t = inp.1.to_reduced(KELVIN).unwrap(); - let v = inp.2.to_reduced(ANGSTROM.powi(3)).unwrap(); - let m = inp.3.to_reduced(MOL).unwrap(); - let s = StateHD::new(t, v, m); - inp.0.evaluate_residual(&s) -} - -/// Residual Helmholtz energy given a StateHD -fn a_res_shd>(inp: (&Arc, &StateHD)) -> D -where - (dyn HelmholtzEnergy + 'static): HelmholtzEnergyDual, -{ - inp.0.evaluate_residual(inp.1) -} - -fn da_dv(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { - let t = Dual64::from_re(inp.1.to_reduced(KELVIN).unwrap()); - let v = Dual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(); - let m = inp - .3 - .to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual64::from_re(mi)) - .collect(); - let s = StateHD::new(t, v, m); - inp.0.evaluate_residual(&s).eps[0] -} - -fn compressibility(inp: (&Arc, SINumber, SINumber, &Array1)) -> f64 { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() - .compressibility(Contributions::Total) -} - -fn compressibility_(state: &State) -> f64 { - state.compressibility(Contributions::Total) -} - -fn ln_phi(inp: (&Arc, SINumber, SINumber, &Array1)) -> Array1 { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() - .ln_phi() -} - -fn da_dni(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> Array1 { - let t = Dual64::from_re(inp.1.to_reduced(KELVIN).unwrap()); - let v = Dual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()); - let m: Array1<_> = inp - .3 - .to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual64::from_re(mi)) - .collect(); - Array1::from_shape_fn(inp.3.len(), |i| { - let mut n = m.clone(); - n[i] = n[i].derive(); - let s = StateHD::new(t, v, n); - inp.0.evaluate_residual(&s).eps[0] - }) -} - -fn c_v(inp: (&Arc, SINumber, SINumber, &Array1)) -> SINumber { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() - .c_v(Contributions::ResidualNvt) -} - -fn d2a_dt2(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { - let t = HyperDual64::from_re(inp.1.to_reduced(KELVIN).unwrap()).derive2(); - let v = HyperDual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()); - let m: Array1<_> = inp - .3 - .to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| HyperDual64::from_re(mi)) - .collect(); - - let s = StateHD::new(t, v, m); - inp.0.evaluate_residual(&s).eps1eps2[0] -} - -fn molar_volume(inp: (&Arc, SINumber, SINumber, &Array1)) -> SIArray1 { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() - .molar_volume(Contributions::ResidualNvt) -} - -fn benchmark_helmholtz_energy(c: &mut Criterion) { - 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))); - let t = 300.0 * KELVIN; - let density = 71.18 * KILO * MOL / METER.powi(3); - let volume = 100.0 * MOL / density; - let x = arr1(&[1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0]); - let rho_i = &x * density; - let m = &x * 100.0 * MOL; - let s = StateHD::new( - 300.0, - volume.to_reduced(ANGSTROM.powi(3)).unwrap(), - &x * 100.0, - ); - - let mut group = c.benchmark_group("helmholtz_energy"); - group.bench_function("state", |b| b.iter(|| state((&eos, t, density, &x)))); - group.bench_function("a_res", |b| b.iter(|| a_res((&eos, t, volume, &m)))); - group.bench_function("a_res_shd", |b| b.iter(|| a_res_shd((&eos, &s)))); -} - -fn benchmark_first_derivative(c: &mut Criterion) { - 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))); - let t = 300.0 * KELVIN; - let density = 71.18 * KILO * MOL / METER.powi(3); - let volume = 100.0 * MOL / density; - let x = arr1(&[1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0]); - let rho_i = &x * density; - let m = &x * 100.0 * MOL; - - let s = StateBuilder::new(&eos) - .temperature(t) - .partial_density(&rho_i) - .total_moles(m.sum()) - .build() - .unwrap(); - - let mut dual = c.benchmark_group("first_derivative"); - // group.bench_function("state", |b| b.iter(|| state((&eos, t, density, &x)))); - dual.bench_function("compressibility", |b| { - b.iter(|| compressibility((&eos, t, density, &x))) - }); - dual.bench_function("compressibility_c", |b| { - b.iter(|| compressibility((&eos, t, density, &x))) - }); - dual.bench_function("da_dv", |b| b.iter(|| da_dv((&eos, t, volume, &m)))); - dual.bench_function("ln_phi", |b| b.iter(|| ln_phi((&eos, t, density, &x)))); - dual.bench_function("da_dni", |b| b.iter(|| da_dni((&eos, t, volume, &m)))); -} - -fn benchmark_second_derivative(c: &mut Criterion) { - 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))); - let t = 300.0 * KELVIN; - let density = 71.18 * KILO * MOL / METER.powi(3); - let volume = 100.0 * MOL / density; - let x = arr1(&[1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0]); - let rho_i = &x * density; - let m = &x * 100.0 * MOL; - - let s = StateBuilder::new(&eos) - .temperature(t) - .partial_density(&rho_i) - .total_moles(m.sum()) - .build() - .unwrap(); - - let mut hyper_dual = c.benchmark_group("second_derivative"); - hyper_dual.bench_function("c_v", |b| b.iter(|| c_v((&eos, t, density, &x)))); - hyper_dual.bench_function("d2a_dt2", |b| b.iter(|| d2a_dt2((&eos, t, volume, &m)))); - - hyper_dual.bench_function("molar_volume", |b| { - b.iter(|| molar_volume((&eos, t, density, &x))) - }); -} - -criterion_group!( - bench, - benchmark_helmholtz_energy, - benchmark_first_derivative, - benchmark_second_derivative -); -criterion_main!(bench); diff --git a/benches/biogas.rs b/benches/biogas.rs deleted file mode 100644 index 445dee300..000000000 --- a/benches/biogas.rs +++ /dev/null @@ -1,162 +0,0 @@ -use criterion::{criterion_group, criterion_main, Criterion}; -use feos::pcsaft::{PcSaft, PcSaftParameters}; -use feos_core::{ - parameter::{IdentifierOption, Parameter}, - Contributions, EquationOfState, HelmholtzEnergy, HelmholtzEnergyDual, State, StateBuilder, - StateHD, -}; -use ndarray::{arr1, Array1}; -use num_dual::{Dual3, Dual64, DualNum, HyperDual64, StaticMat, StaticVec}; -use quantity::{si::*, QuantityArray1}; -use std::sync::Arc; - -/// State generation using `State::new` -fn state(inp: (&Arc, SINumber, SINumber, &Array1)) -> State { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() -} - -/// Residual Helmholtz energy given a StateHD -fn a_res>(inp: (&Arc, &StateHD)) -> D -where - (dyn HelmholtzEnergy + 'static): HelmholtzEnergyDual, -{ - inp.0.evaluate_residual(inp.1) -} - -fn da_dv(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { - let t = Dual64::from_re(inp.1.to_reduced(KELVIN).unwrap()); - let v = Dual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(); - let m = inp - .3 - .to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual64::from_re(mi)) - .collect(); - let s = StateHD::new(t, v, m); - inp.0.evaluate_residual(&s).eps[0] -} - -fn da_dni(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> Array1 { - let t = Dual64::from_re(inp.1.to_reduced(KELVIN).unwrap()); - let v = Dual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()); - let m: Array1<_> = inp - .3 - .to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual64::from_re(mi)) - .collect(); - Array1::from_shape_fn(inp.3.len(), |i| { - let mut n = m.clone(); - n[i] = n[i].derive(); - let s = StateHD::new(t, v, n); - inp.0.evaluate_residual(&s).eps[0] - }) -} - -fn d2a_dt2(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { - let t = HyperDual64::from_re(inp.1.to_reduced(KELVIN).unwrap()).derive2(); - let v = HyperDual64::from_re(inp.2.to_reduced(ANGSTROM.powi(3)).unwrap()); - let m: Array1<_> = inp - .3 - .to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| HyperDual64::from_re(mi)) - .collect(); - - let s = StateHD::new(t, v, m); - inp.0.evaluate_residual(&s).eps1eps2[0] -} - -fn benchmark_helmholtz_energy(c: &mut Criterion) { - let parameters = PcSaftParameters::from_multiple_json( - &[ - (vec!["methane"], "./parameters/pcsaft/gross2001.json"), - ( - vec!["carbon dioxide"], - "./parameters/pcsaft/gross2005_fit.json", - ), - ], - None, - IdentifierOption::Name, - ) - .unwrap(); - let k_ij = -0.0192211646; - let parameters = - PcSaftParameters::new_binary(parameters.pure_records.clone(), Some(k_ij.into())); - let eos = Arc::new(PcSaft::new(Arc::new(parameters))); - - // 230 K, 50 bar, x0 = 0.15 - let t = 230.0 * KELVIN; - let density = 24.16896 * KILO * MOL / METER.powi(3); - let volume = 10.0 * MOL / density; - let x = arr1(&[0.15, 0.85]); - let m = &x * 10.0 * MOL; - - let mut group = c.benchmark_group("helmholtz_energy"); - - // real valued evaluation - let s = StateHD::new( - t.to_reduced(KELVIN).unwrap(), - volume.to_reduced(ANGSTROM.powi(3)).unwrap(), - m.to_reduced(MOL).unwrap(), - ); - group.bench_function("a_res_f64", |b| b.iter(|| a_res((&eos, &s)))); - - // da_dv - dual number - let s = StateHD::new( - Dual64::from_re(t.to_reduced(KELVIN).unwrap()), - Dual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual64::from_re(mi)) - .collect(), - ); - group.bench_function("a_res_dual", |b| b.iter(|| a_res((&eos, &s)))); - - // d2a_dv2 - hyperdual number - let s = StateHD::new( - HyperDual64::from_re(t.to_reduced(KELVIN).unwrap()).derive1(), - HyperDual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive2(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| HyperDual64::from_re(mi)) - .collect(), - ); - group.bench_function("a_res_hyperdual", |b| b.iter(|| a_res((&eos, &s)))); - - // d3a_dv3 - dual3 number - let s = StateHD::new( - Dual3::from_re(t.to_reduced(KELVIN).unwrap()), - Dual3::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual3::from_re(mi)) - .collect(), - ); - group.bench_function("a_res_dual3", |b| b.iter(|| a_res((&eos, &s)))); -} - -criterion_group!(bench, benchmark_helmholtz_energy,); -criterion_main!(bench); diff --git a/benches/dual_numbers.rs b/benches/dual_numbers.rs new file mode 100644 index 000000000..1373ad364 --- /dev/null +++ b/benches/dual_numbers.rs @@ -0,0 +1,103 @@ +//! Benchmarks for the evaluation of the Helmholtz energy function +//! for a given `StateHD` for different types of dual numbers. +//! These should give an idea about the expected slow-down depending +//! on the dual number type used without the overhead of the `State` +//! creation. +//! +//! The example system is the binary mixture of CH4/CO2 modelled +//! with the PCP-SAFT equation of state. The considered Helmholtz +//! energy contributions are: hard-sphere, hard-chain, dispersion, +//! and polar (dipolar).. +use criterion::{criterion_group, criterion_main, Criterion}; +use feos::pcsaft::{PcSaft, PcSaftParameters}; +use feos_core::{ + parameter::{IdentifierOption, Parameter}, + EquationOfState, HelmholtzEnergy, HelmholtzEnergyDual, StateHD, +}; +use ndarray::arr1; +use num_dual::{Dual3, Dual64, DualNum, HyperDual64}; +use quantity::si::*; +use std::sync::Arc; + +/// Residual Helmholtz energy given a StateHD +fn a_res>(inp: (&Arc, &StateHD)) -> D +where + (dyn HelmholtzEnergy + 'static): HelmholtzEnergyDual, +{ + inp.0.evaluate_residual(inp.1) +} + +fn benchmark_dual_numbers(c: &mut Criterion) { + let parameters = PcSaftParameters::from_multiple_json( + &[ + (vec!["methane"], "./parameters/pcsaft/gross2001.json"), + ( + vec!["carbon dioxide"], + "./parameters/pcsaft/gross2005_fit.json", + ), + ], + None, + IdentifierOption::Name, + ) + .unwrap(); + let k_ij = -0.0192211646; + let parameters = + PcSaftParameters::new_binary(parameters.pure_records.clone(), Some(k_ij.into())); + let eos = Arc::new(PcSaft::new(Arc::new(parameters))); + + // 230 K, 50 bar, x0 = 0.15 + let t = 230.0 * KELVIN; + let density = 24.16896 * KILO * MOL / METER.powi(3); + let volume = 10.0 * MOL / density; + let x = arr1(&[0.15, 0.85]); + let m = &x * 10.0 * MOL; + + let mut group = c.benchmark_group("dual_numbers"); + + // real valued evaluation + let s = StateHD::new( + t.to_reduced(KELVIN).unwrap(), + volume.to_reduced(ANGSTROM.powi(3)).unwrap(), + m.to_reduced(MOL).unwrap(), + ); + group.bench_function("a_f64", |b| b.iter(|| a_res((&eos, &s)))); + + // da_dv - dual number + let s = StateHD::new( + Dual64::from_re(t.to_reduced(KELVIN).unwrap()), + Dual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual64::from_re(mi)) + .collect(), + ); + group.bench_function("a_dual", |b| b.iter(|| a_res((&eos, &s)))); + + // d2a_dv2 - hyperdual number + let s = StateHD::new( + HyperDual64::from_re(t.to_reduced(KELVIN).unwrap()).derive1(), + HyperDual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive2(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| HyperDual64::from_re(mi)) + .collect(), + ); + group.bench_function("a_hyperdual", |b| b.iter(|| a_res((&eos, &s)))); + + // d3a_dv3 - dual3 number + let s = StateHD::new( + Dual3::from_re(t.to_reduced(KELVIN).unwrap()), + Dual3::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual3::from_re(mi)) + .collect(), + ); + group.bench_function("a_dual3", |b| b.iter(|| a_res((&eos, &s)))); +} + +criterion_group!(bench, benchmark_dual_numbers); +criterion_main!(bench); diff --git a/benches/state_properties.rs b/benches/state_properties.rs new file mode 100644 index 000000000..f436ebfb1 --- /dev/null +++ b/benches/state_properties.rs @@ -0,0 +1,144 @@ +use criterion::{criterion_group, criterion_main, Criterion}; +use feos::pcsaft::{PcSaft, PcSaftParameters}; +use feos_core::{ + parameter::{IdentifierOption, Parameter}, + Contributions, State, +}; +use ndarray::{arr1, Array1}; +use quantity::si::*; +use std::sync::Arc; + +fn helmholtz_energy(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { + State::new( + inp.0, + Some(inp.1), + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() + .compressibility(Contributions::Total) +} + +fn compressibility(inp: (&Arc, SINumber, SINumber, &Array1)) -> f64 { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() + .compressibility(Contributions::Total) +} + +fn ln_phi(inp: (&Arc, SINumber, SINumber, &Array1)) -> Array1 { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() + .ln_phi() +} + +fn c_v(inp: (&Arc, SINumber, SINumber, &Array1)) -> SINumber { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() + .c_v(Contributions::ResidualNvt) +} + +fn molar_volume(inp: (&Arc, SINumber, SINumber, &Array1)) -> SIArray1 { + State::new( + inp.0, + Some(inp.1), + None, + Some(inp.2), + None, + None, + None, + Some(inp.3), + None, + None, + None, + None, + feos_core::DensityInitialization::None, + None, + ) + .unwrap() + .molar_volume(Contributions::ResidualNvt) +} + +fn benchmark_properties(c: &mut Criterion) { + 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))); + let t = 300.0 * KELVIN; + let density = 71.18 * KILO * MOL / METER.powi(3); + let volume = 100.0 * MOL / density; + 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("properties"); + group.bench_function("a", |b| b.iter(|| helmholtz_energy((&eos, t, volume, &m)))); + group.bench_function("compressibility", |b| { + b.iter(|| compressibility((&eos, t, density, &x))) + }); + group.bench_function("ln_phi", |b| b.iter(|| ln_phi((&eos, t, density, &x)))); + group.bench_function("c_v", |b| b.iter(|| c_v((&eos, t, density, &x)))); + group.bench_function("molar_volume", |b| { + b.iter(|| molar_volume((&eos, t, density, &x))) + }); +} + +criterion_group!(bench, benchmark_properties,); +criterion_main!(bench); From b79360f51d8d7dfda2b86fbdd4b1d9a6eda5120b Mon Sep 17 00:00:00 2001 From: Gernot Bauer Date: Tue, 13 Dec 2022 17:41:54 +0100 Subject: [PATCH 4/7] Changed profile in workflow for publishing wheels to use LTO --- .github/workflows/release.yml | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/.github/workflows/release.yml b/.github/workflows/release.yml index 835d72831..a00b89e80 100644 --- a/.github/workflows/release.yml +++ b/.github/workflows/release.yml @@ -17,7 +17,7 @@ jobs: with: manylinux: auto command: build - args: --release --out dist + args: --profile release-lto --out dist - name: Upload wheels uses: actions/upload-artifact@v2 with: @@ -42,11 +42,11 @@ jobs: uses: messense/maturin-action@v1 with: target: x86_64 - args: --release --out dist + args: --profile release-lto --out dist - name: Build wheels - universal2 uses: messense/maturin-action@v1 with: - args: --release --universal2 --out dist + args: --profile release-lto --universal2 --out dist - name: Upload wheels uses: actions/upload-artifact@v2 with: @@ -74,7 +74,7 @@ jobs: uses: messense/maturin-action@v1 with: target: ${{ matrix.target }} - args: --release --out dist + args: --profile release-lto --out dist - name: Upload wheels uses: actions/upload-artifact@v2 with: From bb93730fa9cb26bab90297441a7770f021426662 Mon Sep 17 00:00:00 2001 From: Gernot Bauer Date: Tue, 13 Dec 2022 18:34:17 +0100 Subject: [PATCH 5/7] added more systems (methane and water) for dual number evaluations, smaller description for README --- benches/README.md | 4 +- benches/dual_numbers.rs | 148 +++++++++++++++++++++++++++++++++--- benches/state_properties.rs | 6 +- 3 files changed, 144 insertions(+), 14 deletions(-) diff --git a/benches/README.md b/benches/README.md index 8f146baae..bc9ab59e5 100644 --- a/benches/README.md +++ b/benches/README.md @@ -12,5 +12,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. Substances: methane + carbon dioxide. Model: PC-SAFT incl. dipole-dipole contributions.|`pcsaft`| -|`state_properties`|Properties of `State`. Including state creation using the natural variables of the Helmholtz energy (no density iteration). Substances: methane + ethane + propane. Model: PC-SAFT (hard-sphere, hard-chain, dispersion).|`pcsaft`| \ No newline at end of file +|`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`| \ No newline at end of file diff --git a/benches/dual_numbers.rs b/benches/dual_numbers.rs index 1373ad364..d99be75b1 100644 --- a/benches/dual_numbers.rs +++ b/benches/dual_numbers.rs @@ -3,16 +3,11 @@ //! These should give an idea about the expected slow-down depending //! on the dual number type used without the overhead of the `State` //! creation. -//! -//! The example system is the binary mixture of CH4/CO2 modelled -//! with the PCP-SAFT equation of state. The considered Helmholtz -//! energy contributions are: hard-sphere, hard-chain, dispersion, -//! and polar (dipolar).. use criterion::{criterion_group, criterion_main, Criterion}; use feos::pcsaft::{PcSaft, PcSaftParameters}; use feos_core::{ parameter::{IdentifierOption, Parameter}, - EquationOfState, HelmholtzEnergy, HelmholtzEnergyDual, StateHD, + EquationOfState, HelmholtzEnergy, HelmholtzEnergyDual, State, StateHD, }; use ndarray::arr1; use num_dual::{Dual3, Dual64, DualNum, HyperDual64}; @@ -27,7 +22,137 @@ where inp.0.evaluate_residual(inp.1) } -fn benchmark_dual_numbers(c: &mut Criterion) { +fn methane_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))); + + let cp = State::critical_point(&eos, None, None, Default::default()).unwrap(); + + let t = 0.8 * cp.temperature; + let density = cp.density; + let volume = 10.0 * MOL / density; + let x = arr1(&[1.0]); + let m = &x * 10.0 * MOL; + + let mut group = c.benchmark_group("PC-SAFT methane"); + + // real valued evaluation + let s = StateHD::new( + t.to_reduced(KELVIN).unwrap(), + volume.to_reduced(ANGSTROM.powi(3)).unwrap(), + m.to_reduced(MOL).unwrap(), + ); + group.bench_function("a_f64", |b| b.iter(|| a_res((&eos, &s)))); + + // da_dv - dual number + let s = StateHD::new( + Dual64::from_re(t.to_reduced(KELVIN).unwrap()), + Dual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual64::from_re(mi)) + .collect(), + ); + group.bench_function("a_dual", |b| b.iter(|| a_res((&eos, &s)))); + + // d2a_dv2 - hyperdual number + let s = StateHD::new( + HyperDual64::from_re(t.to_reduced(KELVIN).unwrap()).derive1(), + HyperDual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive2(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| HyperDual64::from_re(mi)) + .collect(), + ); + group.bench_function("a_hyperdual", |b| b.iter(|| a_res((&eos, &s)))); + + // d3a_dv3 - dual3 number + let s = StateHD::new( + Dual3::from_re(t.to_reduced(KELVIN).unwrap()), + Dual3::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual3::from_re(mi)) + .collect(), + ); + group.bench_function("a_dual3", |b| b.iter(|| a_res((&eos, &s)))); +} + +fn water_4c_polar_pcsaft(c: &mut Criterion) { + let parameters = PcSaftParameters::from_json( + vec!["water_4C_polar"], + "./parameters/pcsaft/rehner2020.json", + None, + IdentifierOption::Name, + ) + .unwrap(); + let eos = Arc::new(PcSaft::new(Arc::new(parameters))); + + let cp = State::critical_point(&eos, None, None, Default::default()).unwrap(); + + let t = 0.8 * cp.temperature; + let density = cp.density; + let volume = 10.0 * MOL / density; + let x = arr1(&[1.0]); + let m = &x * 10.0 * MOL; + + let mut group = c.benchmark_group("PC-SAFT water (4C, polar)"); + + // real valued evaluation + let s = StateHD::new( + t.to_reduced(KELVIN).unwrap(), + volume.to_reduced(ANGSTROM.powi(3)).unwrap(), + m.to_reduced(MOL).unwrap(), + ); + group.bench_function("a_f64", |b| b.iter(|| a_res((&eos, &s)))); + + // da_dv - dual number + let s = StateHD::new( + Dual64::from_re(t.to_reduced(KELVIN).unwrap()), + Dual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual64::from_re(mi)) + .collect(), + ); + group.bench_function("a_dual", |b| b.iter(|| a_res((&eos, &s)))); + + // d2a_dv2 - hyperdual number + let s = StateHD::new( + HyperDual64::from_re(t.to_reduced(KELVIN).unwrap()).derive1(), + HyperDual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive2(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| HyperDual64::from_re(mi)) + .collect(), + ); + group.bench_function("a_hyperdual", |b| b.iter(|| a_res((&eos, &s)))); + + // d3a_dv3 - dual3 number + let s = StateHD::new( + Dual3::from_re(t.to_reduced(KELVIN).unwrap()), + Dual3::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), + m.to_reduced(MOL) + .unwrap() + .iter() + .map(|&mi| Dual3::from_re(mi)) + .collect(), + ); + group.bench_function("a_dual3", |b| b.iter(|| a_res((&eos, &s)))); +} + +fn methane_co2_pcsaft(c: &mut Criterion) { let parameters = PcSaftParameters::from_multiple_json( &[ (vec!["methane"], "./parameters/pcsaft/gross2001.json"), @@ -52,7 +177,7 @@ fn benchmark_dual_numbers(c: &mut Criterion) { let x = arr1(&[0.15, 0.85]); let m = &x * 10.0 * MOL; - let mut group = c.benchmark_group("dual_numbers"); + let mut group = c.benchmark_group("PC-SAFT methane + co2"); // real valued evaluation let s = StateHD::new( @@ -99,5 +224,10 @@ fn benchmark_dual_numbers(c: &mut Criterion) { group.bench_function("a_dual3", |b| b.iter(|| a_res((&eos, &s)))); } -criterion_group!(bench, benchmark_dual_numbers); +criterion_group!( + bench, + methane_pcsaft, + water_4c_polar_pcsaft, + methane_co2_pcsaft +); criterion_main!(bench); diff --git a/benches/state_properties.rs b/benches/state_properties.rs index f436ebfb1..74467c1df 100644 --- a/benches/state_properties.rs +++ b/benches/state_properties.rs @@ -113,7 +113,7 @@ fn molar_volume(inp: (&Arc, SINumber, SINumber, &Array1)) -> SIArra .molar_volume(Contributions::ResidualNvt) } -fn benchmark_properties(c: &mut Criterion) { +fn properties_pcsaft(c: &mut Criterion) { let parameters = PcSaftParameters::from_json( vec!["methane", "ethane", "propane"], "./parameters/pcsaft/gross2001.json", @@ -128,7 +128,7 @@ fn benchmark_properties(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("properties"); + let mut group = c.benchmark_group("PC-SAFT methane + ethane + propane"); group.bench_function("a", |b| b.iter(|| helmholtz_energy((&eos, t, volume, &m)))); group.bench_function("compressibility", |b| { b.iter(|| compressibility((&eos, t, density, &x))) @@ -140,5 +140,5 @@ fn benchmark_properties(c: &mut Criterion) { }); } -criterion_group!(bench, benchmark_properties,); +criterion_group!(bench, properties_pcsaft); criterion_main!(bench); From 962592b4690048cd32bbc36404acde5adc36ab8b Mon Sep 17 00:00:00 2001 From: Gernot Bauer Date: Wed, 14 Dec 2022 12:25:53 +0100 Subject: [PATCH 6/7] Apply suggestions from code review concerning changes to the text and typos. Co-authored-by: Philipp Rehner <69816385+prehner@users.noreply.github.com> --- README.md | 4 ++-- benches/state_properties.rs | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/README.md b/README.md index 0bd621457..9c06c54bc 100644 --- a/README.md +++ b/README.md @@ -96,7 +96,7 @@ If there is no compiled package for your system available from PyPI and you have pip install git+https://github.com/feos-org/feos ``` -This command builds the package without link-time optimization (LTO) which results in considerably worse performance. +This command builds the package without link-time optimization (LTO) that can be used to increase the performance further. See the *Building from source* section for information about building the wheel including LTO. ### Building from source @@ -131,7 +131,7 @@ Use the following command to build a wheel with specific features: maturin build --profile="release-lto" --features "python ..." ``` -LTO increases compile times drastically but the resulting wheel is way more performant and has a smaller size. +LTO increases compile times measurably but the resulting wheel is more performant and has a smaller size. For development however, we recommend using the `--release` flag. ## Documentation diff --git a/benches/state_properties.rs b/benches/state_properties.rs index 74467c1df..d8db356b7 100644 --- a/benches/state_properties.rs +++ b/benches/state_properties.rs @@ -26,7 +26,7 @@ fn helmholtz_energy(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { None, ) .unwrap() - .compressibility(Contributions::Total) + .helmholtz_energy(Contributions::Total) } fn compressibility(inp: (&Arc, SINumber, SINumber, &Array1)) -> f64 { From a33171544c3da8084d7567051799fa3d0f06d669 Mon Sep 17 00:00:00 2001 From: Gernot Bauer Date: Wed, 14 Dec 2022 14:17:28 +0100 Subject: [PATCH 7/7] added utilities to better extend benchmarks reusing the deriveX-methods from `State` --- benches/dual_numbers.rs | 232 ++++++++++-------------------------- benches/state_properties.rs | 146 ++++++----------------- feos-core/src/lib.rs | 2 +- feos-core/src/state/mod.rs | 18 ++- 4 files changed, 112 insertions(+), 286 deletions(-) diff --git a/benches/dual_numbers.rs b/benches/dual_numbers.rs index d99be75b1..bf9f37f38 100644 --- a/benches/dual_numbers.rs +++ b/benches/dual_numbers.rs @@ -7,22 +7,58 @@ use criterion::{criterion_group, criterion_main, Criterion}; use feos::pcsaft::{PcSaft, PcSaftParameters}; use feos_core::{ parameter::{IdentifierOption, Parameter}, - EquationOfState, HelmholtzEnergy, HelmholtzEnergyDual, State, StateHD, + Derivative, EquationOfState, HelmholtzEnergy, HelmholtzEnergyDual, State, StateHD, }; -use ndarray::arr1; -use num_dual::{Dual3, Dual64, DualNum, HyperDual64}; +use ndarray::{arr1, Array}; +use num_dual::DualNum; use quantity::si::*; use std::sync::Arc; -/// Residual Helmholtz energy given a StateHD -fn a_res>(inp: (&Arc, &StateHD)) -> D +/// Helper function to create a state for given parameters. +/// - temperature is 80% of critical temperature, +/// - volume is critical volume, +/// - molefracs (or moles) for equimolar mixture. +fn state_pcsaft(parameters: PcSaftParameters) -> State { + let n = parameters.pure_records.len(); + let eos = Arc::new(PcSaft::new(Arc::new(parameters))); + let moles = Array::from_elem(n, 1.0 / n as f64) * 10.0 * MOL; + let cp = State::critical_point(&eos, Some(&moles), None, Default::default()).unwrap(); + let temperature = 0.8 * cp.temperature; + State::new_nvt(&eos, temperature, cp.volume, &moles).unwrap() +} + +/// Residual Helmholtz energy given an equation of state and a StateHD. +fn a_res, E: EquationOfState>(inp: (&Arc, &StateHD)) -> D where (dyn HelmholtzEnergy + 'static): HelmholtzEnergyDual, { inp.0.evaluate_residual(inp.1) } -fn methane_pcsaft(c: &mut Criterion) { +/// Benchmark for evaluation of the Helmholtz energy for different dual number types. +fn bench_dual_numbers( + c: &mut Criterion, + group_name: &str, + state: State, +) { + let mut group = c.benchmark_group(group_name); + group.bench_function("a_f64", |b| { + b.iter(|| a_res((&state.eos, &state.derive0()))) + }); + group.bench_function("a_dual", |b| { + b.iter(|| a_res((&state.eos, &state.derive1(Derivative::DV)))) + }); + group.bench_function("a_hyperdual", |b| { + b.iter(|| a_res((&state.eos, &state.derive2(Derivative::DV, Derivative::DV)))) + }); + group.bench_function("a_dual3", |b| { + b.iter(|| a_res((&state.eos, &state.derive3(Derivative::DV)))) + }); +} + +/// Benchmark for the PC-SAFT equation of state +fn pcsaft(c: &mut Criterion) { + // methane let parameters = PcSaftParameters::from_json( vec!["methane"], "./parameters/pcsaft/gross2001.json", @@ -30,64 +66,9 @@ fn methane_pcsaft(c: &mut Criterion) { IdentifierOption::Name, ) .unwrap(); - let eos = Arc::new(PcSaft::new(Arc::new(parameters))); - - let cp = State::critical_point(&eos, None, None, Default::default()).unwrap(); - - let t = 0.8 * cp.temperature; - let density = cp.density; - let volume = 10.0 * MOL / density; - let x = arr1(&[1.0]); - let m = &x * 10.0 * MOL; + bench_dual_numbers(c, "methane", state_pcsaft(parameters)); - let mut group = c.benchmark_group("PC-SAFT methane"); - - // real valued evaluation - let s = StateHD::new( - t.to_reduced(KELVIN).unwrap(), - volume.to_reduced(ANGSTROM.powi(3)).unwrap(), - m.to_reduced(MOL).unwrap(), - ); - group.bench_function("a_f64", |b| b.iter(|| a_res((&eos, &s)))); - - // da_dv - dual number - let s = StateHD::new( - Dual64::from_re(t.to_reduced(KELVIN).unwrap()), - Dual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual64::from_re(mi)) - .collect(), - ); - group.bench_function("a_dual", |b| b.iter(|| a_res((&eos, &s)))); - - // d2a_dv2 - hyperdual number - let s = StateHD::new( - HyperDual64::from_re(t.to_reduced(KELVIN).unwrap()).derive1(), - HyperDual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive2(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| HyperDual64::from_re(mi)) - .collect(), - ); - group.bench_function("a_hyperdual", |b| b.iter(|| a_res((&eos, &s)))); - - // d3a_dv3 - dual3 number - let s = StateHD::new( - Dual3::from_re(t.to_reduced(KELVIN).unwrap()), - Dual3::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual3::from_re(mi)) - .collect(), - ); - group.bench_function("a_dual3", |b| b.iter(|| a_res((&eos, &s)))); -} - -fn water_4c_polar_pcsaft(c: &mut Criterion) { + // water (4C, polar) let parameters = PcSaftParameters::from_json( vec!["water_4C_polar"], "./parameters/pcsaft/rehner2020.json", @@ -95,63 +76,21 @@ fn water_4c_polar_pcsaft(c: &mut Criterion) { IdentifierOption::Name, ) .unwrap(); - let eos = Arc::new(PcSaft::new(Arc::new(parameters))); - - let cp = State::critical_point(&eos, None, None, Default::default()).unwrap(); + bench_dual_numbers(c, "water_4c_polar", state_pcsaft(parameters)); - let t = 0.8 * cp.temperature; - let density = cp.density; - let volume = 10.0 * MOL / density; - let x = arr1(&[1.0]); - let m = &x * 10.0 * MOL; - - let mut group = c.benchmark_group("PC-SAFT water (4C, polar)"); - - // real valued evaluation - let s = StateHD::new( - t.to_reduced(KELVIN).unwrap(), - volume.to_reduced(ANGSTROM.powi(3)).unwrap(), - m.to_reduced(MOL).unwrap(), - ); - group.bench_function("a_f64", |b| b.iter(|| a_res((&eos, &s)))); - - // da_dv - dual number - let s = StateHD::new( - Dual64::from_re(t.to_reduced(KELVIN).unwrap()), - Dual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual64::from_re(mi)) - .collect(), - ); - group.bench_function("a_dual", |b| b.iter(|| a_res((&eos, &s)))); - - // d2a_dv2 - hyperdual number - let s = StateHD::new( - HyperDual64::from_re(t.to_reduced(KELVIN).unwrap()).derive1(), - HyperDual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive2(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| HyperDual64::from_re(mi)) - .collect(), - ); - group.bench_function("a_hyperdual", |b| b.iter(|| a_res((&eos, &s)))); - - // d3a_dv3 - dual3 number - let s = StateHD::new( - Dual3::from_re(t.to_reduced(KELVIN).unwrap()), - Dual3::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual3::from_re(mi)) - .collect(), - ); - group.bench_function("a_dual3", |b| b.iter(|| a_res((&eos, &s)))); + // methane, ethane, propane + let parameters = PcSaftParameters::from_json( + vec!["methane", "ethane", "propane"], + "./parameters/pcsaft/gross2001.json", + None, + IdentifierOption::Name, + ) + .unwrap(); + bench_dual_numbers(c, "methane_ethane_propane", state_pcsaft(parameters)); } +/// Benchmark for the PC-SAFT equation of state. +/// Binary system of methane and co2 used to model biogas. fn methane_co2_pcsaft(c: &mut Criterion) { let parameters = PcSaftParameters::from_multiple_json( &[ @@ -171,63 +110,14 @@ fn methane_co2_pcsaft(c: &mut Criterion) { let eos = Arc::new(PcSaft::new(Arc::new(parameters))); // 230 K, 50 bar, x0 = 0.15 - let t = 230.0 * KELVIN; + let temperature = 230.0 * KELVIN; let density = 24.16896 * KILO * MOL / METER.powi(3); let volume = 10.0 * MOL / density; let x = arr1(&[0.15, 0.85]); - let m = &x * 10.0 * MOL; - - let mut group = c.benchmark_group("PC-SAFT methane + co2"); - - // real valued evaluation - let s = StateHD::new( - t.to_reduced(KELVIN).unwrap(), - volume.to_reduced(ANGSTROM.powi(3)).unwrap(), - m.to_reduced(MOL).unwrap(), - ); - group.bench_function("a_f64", |b| b.iter(|| a_res((&eos, &s)))); - - // da_dv - dual number - let s = StateHD::new( - Dual64::from_re(t.to_reduced(KELVIN).unwrap()), - Dual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual64::from_re(mi)) - .collect(), - ); - group.bench_function("a_dual", |b| b.iter(|| a_res((&eos, &s)))); - - // d2a_dv2 - hyperdual number - let s = StateHD::new( - HyperDual64::from_re(t.to_reduced(KELVIN).unwrap()).derive1(), - HyperDual64::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive2(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| HyperDual64::from_re(mi)) - .collect(), - ); - group.bench_function("a_hyperdual", |b| b.iter(|| a_res((&eos, &s)))); - - // d3a_dv3 - dual3 number - let s = StateHD::new( - Dual3::from_re(t.to_reduced(KELVIN).unwrap()), - Dual3::from_re(volume.to_reduced(ANGSTROM.powi(3)).unwrap()).derive(), - m.to_reduced(MOL) - .unwrap() - .iter() - .map(|&mi| Dual3::from_re(mi)) - .collect(), - ); - group.bench_function("a_dual3", |b| b.iter(|| a_res((&eos, &s)))); + let moles = &x * 10.0 * MOL; + let state = State::new_nvt(&eos, temperature, volume, &moles).unwrap(); + bench_dual_numbers(c, "methane_co2", state); } -criterion_group!( - bench, - methane_pcsaft, - water_4c_polar_pcsaft, - methane_co2_pcsaft -); +criterion_group!(bench, pcsaft, methane_co2_pcsaft); criterion_main!(bench); diff --git a/benches/state_properties.rs b/benches/state_properties.rs index d8db356b7..e9b386c2b 100644 --- a/benches/state_properties.rs +++ b/benches/state_properties.rs @@ -2,115 +2,37 @@ use criterion::{criterion_group, criterion_main, Criterion}; use feos::pcsaft::{PcSaft, PcSaftParameters}; use feos_core::{ parameter::{IdentifierOption, Parameter}, - Contributions, State, + Contributions, EquationOfState, State, }; -use ndarray::{arr1, Array1}; +use ndarray::arr1; use quantity::si::*; use std::sync::Arc; -fn helmholtz_energy(inp: (&Arc, SINumber, SINumber, &SIArray1)) -> f64 { - State::new( - inp.0, - Some(inp.1), - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() - .helmholtz_energy(Contributions::Total) -} - -fn compressibility(inp: (&Arc, SINumber, SINumber, &Array1)) -> f64 { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() - .compressibility(Contributions::Total) -} - -fn ln_phi(inp: (&Arc, SINumber, SINumber, &Array1)) -> Array1 { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() - .ln_phi() -} +type S = State; -fn c_v(inp: (&Arc, SINumber, SINumber, &Array1)) -> SINumber { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() - .c_v(Contributions::ResidualNvt) +/// Evaluate a property of a state given the EoS, the property to compute, +/// temperature, volume, moles, and the contributions to consider. +fn property, Contributions) -> T>( + (eos, property, t, v, n, contributions): ( + &Arc, + F, + SINumber, + SINumber, + &SIArray1, + Contributions, + ), +) -> T { + let state = State::new_nvt(eos, t, v, n).unwrap(); + property(&state, contributions) } -fn molar_volume(inp: (&Arc, SINumber, SINumber, &Array1)) -> SIArray1 { - State::new( - inp.0, - Some(inp.1), - None, - Some(inp.2), - None, - None, - None, - Some(inp.3), - None, - None, - None, - None, - feos_core::DensityInitialization::None, - None, - ) - .unwrap() - .molar_volume(Contributions::ResidualNvt) +/// Evaluate a property with of a state given the EoS, the property to compute, +/// temperature, volume, moles. +fn property_no_contributions) -> T>( + (eos, property, t, v, n): (&Arc, F, SINumber, SINumber, &SIArray1), +) -> T { + let state = State::new_nvt(eos, t, v, n).unwrap(); + property(&state) } fn properties_pcsaft(c: &mut Criterion) { @@ -124,19 +46,25 @@ fn properties_pcsaft(c: &mut Criterion) { let eos = Arc::new(PcSaft::new(Arc::new(parameters))); let t = 300.0 * KELVIN; let density = 71.18 * KILO * MOL / METER.powi(3); - let volume = 100.0 * MOL / density; + let v = 100.0 * MOL / density; 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("PC-SAFT methane + ethane + propane"); - group.bench_function("a", |b| b.iter(|| helmholtz_energy((&eos, t, volume, &m)))); + let mut group = c.benchmark_group("methane_ethane_propane"); + group.bench_function("a", |b| { + b.iter(|| property((&eos, S::helmholtz_energy, t, v, &m, Contributions::Total))) + }); group.bench_function("compressibility", |b| { - b.iter(|| compressibility((&eos, t, density, &x))) + b.iter(|| property((&eos, S::compressibility, t, v, &m, Contributions::Total))) + }); + group.bench_function("ln_phi", |b| { + b.iter(|| property_no_contributions((&eos, S::ln_phi, t, v, &m))) + }); + group.bench_function("c_v", |b| { + b.iter(|| property((&eos, S::c_v, t, v, &m, Contributions::ResidualNvt))) }); - group.bench_function("ln_phi", |b| b.iter(|| ln_phi((&eos, t, density, &x)))); - group.bench_function("c_v", |b| b.iter(|| c_v((&eos, t, density, &x)))); group.bench_function("molar_volume", |b| { - b.iter(|| molar_volume((&eos, t, density, &x))) + b.iter(|| property((&eos, S::molar_volume, t, v, &m, Contributions::ResidualNvt))) }); } diff --git a/feos-core/src/lib.rs b/feos-core/src/lib.rs index c8639edf7..3f2084ba5 100644 --- a/feos-core/src/lib.rs +++ b/feos-core/src/lib.rs @@ -42,7 +42,7 @@ pub use errors::{EosError, EosResult}; pub use phase_equilibria::{ PhaseDiagram, PhaseDiagramHetero, PhaseEquilibrium, SolverOptions, Verbosity, }; -pub use state::{Contributions, DensityInitialization, State, StateBuilder, StateHD, StateVec}; +pub use state::{Contributions, DensityInitialization, State, StateBuilder, StateHD, StateVec, Derivative}; #[cfg(feature = "python")] pub mod python; diff --git a/feos-core/src/state/mod.rs b/feos-core/src/state/mod.rs index 5b061686c..13ff9792a 100644 --- a/feos-core/src/state/mod.rs +++ b/feos-core/src/state/mod.rs @@ -188,7 +188,7 @@ where /// Derivatives of the helmholtz energy. #[derive(Clone, Copy, Eq, Hash, PartialEq, Debug, PartialOrd, Ord)] #[allow(non_camel_case_types)] -pub(crate) enum Derivative { +pub enum Derivative { /// Derivative with respect to system volume. DV, /// Derivative with respect to temperature. @@ -670,7 +670,8 @@ impl State { Err(EosError::NotConverged("State::update_gibbs_energy".into())) } - fn derive0(&self) -> StateHD { + /// Creates a [StateHD] cloning temperature, volume and moles. + pub fn derive0(&self) -> StateHD { StateHD::new( self.reduced_temperature, self.reduced_volume, @@ -678,7 +679,8 @@ impl State { ) } - fn derive1(&self, derivative: Derivative) -> StateHD { + /// Creates a [StateHD] taking the first derivative. + pub fn derive1(&self, derivative: Derivative) -> StateHD { let mut t = Dual64::from(self.reduced_temperature); let mut v = Dual64::from(self.reduced_volume); let mut n = self.reduced_moles.mapv(Dual64::from); @@ -690,7 +692,12 @@ impl State { StateHD::new(t, v, n) } - fn derive2(&self, derivative1: Derivative, derivative2: Derivative) -> StateHD { + /// Creates a [StateHD] taking the first and second (partial) derivatives. + pub fn derive2( + &self, + derivative1: Derivative, + derivative2: Derivative, + ) -> StateHD { let mut t = HyperDual64::from(self.reduced_temperature); let mut v = HyperDual64::from(self.reduced_volume); let mut n = self.reduced_moles.mapv(HyperDual64::from); @@ -707,7 +714,8 @@ impl State { StateHD::new(t, v, n) } - fn derive3(&self, derivative: Derivative) -> StateHD { + /// Creates a [StateHD] taking the first, second, and third derivative with respect to a single property. + pub fn derive3(&self, derivative: Derivative) -> StateHD { let mut t = Dual3_64::from(self.reduced_temperature); let mut v = Dual3_64::from(self.reduced_volume); let mut n = self.reduced_moles.mapv(Dual3_64::from);