Add critical_temperature to VaporPressure::new and Python construct…

…or. Add estimator module to documentation
feos-org · g-bauer · Dec 13, 2022 · Dec 13, 2022 · Dec 13, 2022 · Dec 13, 2022
commit b60f41df58561b5639ca73d4a30c2356fff00a01
diff --git a/CHANGELOG.md b/CHANGELOG.md
@@ -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.

diff --git a/docs/api/eos.md b/docs/api/eos.md
@@ -34,4 +34,18 @@ The `State` and `PhaseEquilibrium` objects are used to define thermodynamic cond
     State
     PhaseEquilibrium
     PhaseDiagram
-```
+```
+
+## Estimator
+
+```{eval-rst}
+.. currentmodule:: feos.eos.estimator
+
+.. autosummary::
+    :toctree: generated/
+
+    DataSet
+    Loss
+    Estimator
+    Phase
+```
diff --git a/src/estimator/python.rs b/src/estimator/python.rs
@@ -11,6 +11,7 @@ impl From<EstimatorError> for PyErr {
 #[macro_export]
 macro_rules! impl_estimator {
     ($eos:ty, $py_eos:ty) => {
+        /// Loss function that is applied to the residuum to calculate costs.
         #[pyclass(name = "Loss")]
         #[derive(Clone)]
         pub struct PyLoss(Loss);
@@ -221,6 +222,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.
@@ -234,12 +239,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>,
@@ -248,6 +260,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()),
                 )?)))
             }
@@ -439,7 +452,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()

diff --git a/src/estimator/vapor_pressure.rs b/src/estimator/vapor_pressure.rs
@@ -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,
@@ -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);