Skip to content

Add critical_temperature argument to VaporPressure constructor #86

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
g-bauer merged 4 commits into from
Dec 13, 2022
Merged

Add critical_temperature argument to VaporPressure constructor #86

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
b60f41d
Add critical_temperature to `VaporPressure::new` and Python construct…
g-bauer Dec 13, 2022
fdebf40
improvements to documentation
g-bauer Dec 13, 2022
80410d7
fixed typo in changelog
g-bauer Dec 13, 2022
4cce7d5
better docstring for Loss
g-bauer Dec 13, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions CHANGELOG.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,9 +7,11 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
## [Unreleased]
### 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)

### Changed
- Export `EosVariant` and `FunctionalVariant` directly in the crate root instead of their own modules. [#62](https://github.com/feos-org/feos/pull/62)
- Changed constructors `VaporPressure::new` and `DataSet.vapor_pressure` (Python) to take a new optional argument `critical_temperature`. [#86](https://github.com/feos-org/feos/pull/86)

## [0.3.0] - 2022-09-14
- Major restructuring of the entire `feos` project. All individual models are reunited in the `feos` crate. `feos-core` and `feos-dft` still live as individual crates within the `feos` workspace.
Expand Down
24 changes: 23 additions & 1 deletion docs/api/eos.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,8 @@
The `eos` module contains the `EquationOfState` object that contains all implemented equations of state.
The `State` and `PhaseEquilibrium` objects are used to define thermodynamic conditions and -- once created -- can be used to compute properties.

If you want to adjust parameters of a model to experimental data you can use classes and utilities from the `estimator` module.

## `EquationOfState`

```{eval-rst}
Expand Down Expand Up @@ -34,4 +36,24 @@ The `State` and `PhaseEquilibrium` objects are used to define thermodynamic cond
State
PhaseEquilibrium
PhaseDiagram
```
```

## The `estimator` module

### Import

```python
from feos.eos.estimator import Estimator, DataSet, Loss, Phase
```

```{eval-rst}
.. currentmodule:: feos.eos.estimator

.. autosummary::
:toctree: generated/

Estimator
DataSet
Loss
Phase
```
71 changes: 44 additions & 27 deletions src/estimator/python.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,8 @@ impl From<EstimatorError> for PyErr {
#[macro_export]
macro_rules! impl_estimator {
($eos:ty, $py_eos:ty) => {
/// Collection of loss functions that can be applied to residuals
/// to handle outliers.
#[pyclass(name = "Loss")]
#[derive(Clone)]
pub struct PyLoss(Loss);
Expand Down Expand Up @@ -121,20 +123,26 @@ macro_rules! impl_estimator {
impl PyDataSet {
/// Compute the cost function for each input value.
///
/// The cost function that is used depends on the
/// property. See the class constructors to learn
/// about the cost functions of the properties.
///
/// Parameters
/// ----------
/// eos : PyEos
/// eos : EquationOfState
/// The equation of state that is used.
/// loss : Loss
/// The loss function that is applied to residuals
/// to handle outliers.
///
/// Returns
/// -------
/// numpy.ndarray[Float]
/// The cost function evaluated for each experimental data point.
#[pyo3(text_signature = "($self, eos)")]
///
/// Note
/// ----
/// The cost function that is used depends on the
/// property. For most properties it is the absolute relative difference.
/// See the constructors of the respective properties
/// to learn about the cost functions that are used.
#[pyo3(text_signature = "($self, eos, loss)")]
fn cost<'py>(
&self,
eos: &$py_eos,
Expand All @@ -149,18 +157,12 @@ macro_rules! impl_estimator {
///
/// Parameters
/// ----------
/// eos : PyEos
/// eos : EquationOfState
/// The equation of state that is used.
///
/// Returns
/// -------
/// SIArray1
///
/// See also
/// --------
/// eos_python.saft.estimator.DataSet.vapor_pressure : ``DataSet`` for vapor pressure.
/// eos_python.saft.estimator.DataSet.liquid_density : ``DataSet`` for liquid density.
/// eos_python.saft.estimator.DataSet.equilibrium_liquid_density : ``DataSet`` for liquid density at vapor liquid equilibrium.
#[pyo3(text_signature = "($self, eos)")]
fn predict(&self, eos: &$py_eos) -> PyResult<PySIArray1> {
Ok(self.0.predict(&eos.0)?.into())
Expand All @@ -175,7 +177,7 @@ macro_rules! impl_estimator {
///
/// Parameters
/// ----------
/// eos : PyEos
/// eos : EquationOfState
/// The equation of state that is used.
///
/// Returns
Expand All @@ -198,7 +200,7 @@ macro_rules! impl_estimator {
///
/// Parameters
/// ----------
/// eos : PyEos
/// eos : EquationOfState
/// The equation of state that is used.
///
/// Returns
Expand All @@ -221,6 +223,10 @@ macro_rules! impl_estimator {
/// Use Antoine type equation to extrapolate vapor
/// pressure if experimental data is above critial
/// point of model. Defaults to False.
/// critical_temperature : SINumber, optional
/// Estimate of the critical temperature used as initial
/// value for critical point calculation. Defaults to None.
/// For additional information, see note.
/// max_iter : int, optional
/// The maximum number of iterations for critical point
/// and VLE algorithms.
Expand All @@ -234,12 +240,19 @@ macro_rules! impl_estimator {
/// Returns
/// -------
/// ``DataSet``
///
/// Note
/// ----
/// If no critical temperature is provided, the maximum of the `temperature` input
/// is used. If that fails, the default temperatures of the critical point routine
/// are used.
#[staticmethod]
#[pyo3(text_signature = "(target, temperature, extrapolate)")]
#[pyo3(text_signature = "(target, temperature, extrapolate, critical_temperature=None, max_iter=None, verbosity=None)")]
fn vapor_pressure(
target: &PySIArray1,
temperature: &PySIArray1,
extrapolate: Option<bool>,
critical_temperature: Option<&PySINumber>,
max_iter: Option<usize>,
tol: Option<f64>,
verbosity: Option<Verbosity>,
Expand All @@ -248,6 +261,7 @@ macro_rules! impl_estimator {
target.clone().into(),
temperature.clone().into(),
extrapolate.unwrap_or(false),
critical_temperature.and_then(|tc| Some(tc.clone().into())),
Some((max_iter, tol, verbosity).into()),
)?)))
}
Expand Down Expand Up @@ -439,7 +453,7 @@ macro_rules! impl_estimator {
PySIArray1::from(self.0.target().clone())
}

/// Return `target` as ``SIArray1``.
/// Return number of stored data points.
#[getter]
fn get_datapoints(&self) -> usize {
self.0.datapoints()
Expand Down Expand Up @@ -483,22 +497,25 @@ macro_rules! impl_estimator {

/// Compute the cost function for each ``DataSet``.
///
/// The cost function is:
/// - The relative difference between prediction and target value,
/// - to which a loss function is applied,
/// - and which is weighted according to the number of datapoints,
/// - and the relative weights as defined in the Estimator object.
///
/// Parameters
/// ----------
/// eos : PyEos
/// eos : EquationOfState
/// The equation of state that is used.
///
/// Returns
/// -------
/// numpy.ndarray[Float]
/// The cost function evaluated for each experimental data point
/// of each ``DataSet``.
///
/// Note
/// ----
/// The cost function is:
///
/// - The relative difference between prediction and target value,
/// - to which a loss function is applied,
/// - and which is weighted according to the number of datapoints,
/// - and the relative weights as defined in the Estimator object.
#[pyo3(text_signature = "($self, eos)")]
fn cost<'py>(&self, eos: &$py_eos, py: Python<'py>) -> PyResult<&'py PyArray1<f64>> {
Ok(self.0.cost(&eos.0)?.view().to_pyarray(py))
Expand All @@ -509,7 +526,7 @@ macro_rules! impl_estimator {
///
/// Parameters
/// ----------
/// eos : PyEos
/// eos : EquationOfState
/// The equation of state that is used.
///
/// Returns
Expand All @@ -534,7 +551,7 @@ macro_rules! impl_estimator {
///
/// Parameters
/// ----------
/// eos : PyEos
/// eos : EquationOfState
/// The equation of state that is used.
///
/// Returns
Expand Down Expand Up @@ -562,7 +579,7 @@ macro_rules! impl_estimator {
///
/// Parameters
/// ----------
/// eos : PyEos
/// eos : EquationOfState
/// The equation of state that is used.
///
/// Returns
Expand Down
18 changes: 11 additions & 7 deletions src/estimator/vapor_pressure.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -29,15 +29,18 @@ impl<U: EosUnit> VaporPressure<U> {
target: QuantityArray1<U>,
temperature: QuantityArray1<U>,
extrapolate: bool,
critical_temperature: Option<QuantityScalar<U>>,
solver_options: Option<SolverOptions>,
) -> Result<Self, EstimatorError> {
let datapoints = target.len();
let max_temperature = temperature
.to_reduced(U::reference_temperature())?
.into_iter()
.reduce(|a, b| a.max(b))
.unwrap()
* U::reference_temperature();
let max_temperature = critical_temperature.unwrap_or(
temperature
.to_reduced(U::reference_temperature())?
.into_iter()
.reduce(|a, b| a.max(b))
.unwrap()
* U::reference_temperature(),
);
Ok(Self {
target,
temperature,
Expand Down Expand Up @@ -72,7 +75,8 @@ impl<U: EosUnit, E: EquationOfState> DataSet<U, E> for VaporPressure<U> {
QuantityScalar<U>: std::fmt::Display + std::fmt::LowerExp,
{
let critical_point =
State::critical_point(eos, None, Some(self.max_temperature), self.solver_options)?;
State::critical_point(eos, None, Some(self.max_temperature), self.solver_options)
.or_else(|_| State::critical_point(eos, None, None, self.solver_options))?;
let tc = critical_point.temperature;
let pc = critical_point.pressure(Contributions::Total);

Expand Down