Implement binary association parameters

CHANGELOG.md

@@ -8,6 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 ### Added
 - Added `IdealGasModel` enum that collects all implementors of the `IdealGas` trait. [#158](https://github.com/feos-org/feos/pull/158)
 - Added `feos.ideal_gas` module in Python from which (currently) `Joback` and `JobackParameters` are available. [#158](https://github.com/feos-org/feos/pull/158)
+- Added binary association parameters to PC-SAFT. [#167](https://github.com/feos-org/feos/pull/167)
 
 ### Changed
 - Changed the internal implementation of the association contribution to accomodate more general association schemes. [#150](https://github.com/feos-org/feos/pull/150)

feos-core/CHANGELOG.md

@@ -8,6 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 ### Added
 - Added `Components`, `Residual`, `IdealGas` and `DeBroglieWavelength` traits to decouple ideal gas models from residual models. [#158](https://github.com/feos-org/feos/pull/158)
 - Added `JobackParameters` struct that implements `Parameters` including Python bindings. [#158](https://github.com/feos-org/feos/pull/158)
+- Added `Parameter::from_model_records` as a simpler interface to generate parameters.  [#169](https://github.com/feos-org/feos/pull/169)
 
 ### Changed
 - Changed `EquationOfState` from a trait to a `struct` that is generic over `Residual` and `IdealGas` and implements all necessary traits to be used as equation of state including the ideal gas contribution. [#158](https://github.com/feos-org/feos/pull/158)
@@ -19,6 +20,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Moved `StateVec` into own file and module. [#158](https://github.com/feos-org/feos/pull/158)
 - Ideal gas and residual Helmholtz energy models can now be separately implemented in Python via the `PyIdealGas` and `PyResidual` structs. [#158](https://github.com/feos-org/feos/pull/158)
 - Bubble and dew point iterations will not attempt a second iteration if no solution is found for the given initial pressure. [#166](https://github.com/feos-org/feos/pull/166)
+- Made the binary records in the constructions and getters of the `Parameter` trait optional. [#169](https://github.com/feos-org/feos/pull/169)
+- Changed the second argument of `new_binary` in Python from a `BinaryRecord` to the corresponding binary model record (analogous to the Rust implementation). [#169](https://github.com/feos-org/feos/pull/169)
 
 ### Removed
 - Removed `EquationOfState` trait. [#158](https://github.com/feos-org/feos/pull/158)

feos-core/src/cubic.rs

@@ -100,10 +100,7 @@ impl PengRobinsonParameters {
                 PureRecord::new(id, molarweight[i], record)
             })
             .collect();
-        Ok(PengRobinsonParameters::from_records(
-            records,
-            Array2::zeros([pc.len(); 2]),
-        ))
+        Ok(PengRobinsonParameters::from_records(records, None))
     }
 }
 
@@ -114,7 +111,7 @@ impl Parameter for PengRobinsonParameters {
     /// Creates parameters from pure component records.
     fn from_records(
         pure_records: Vec<PureRecord<Self::Pure>>,
-        binary_records: Array2<Self::Binary>,
+        binary_records: Option<Array2<Self::Binary>>,
     ) -> Self {
         let n = pure_records.len();
 
@@ -133,19 +130,21 @@ impl Parameter for PengRobinsonParameters {
             kappa[i] = 0.37464 + (1.54226 - 0.26992 * r.acentric_factor) * r.acentric_factor;
         }
 
+        let k_ij = binary_records.unwrap_or_else(|| Array2::zeros([n; 2]));
+
         Self {
             tc,
             a,
             b,
-            k_ij: binary_records,
+            k_ij,
             kappa,
             molarweight,
             pure_records,
         }
     }
 
-    fn records(&self) -> (&[PureRecord<PengRobinsonRecord>], &Array2<f64>) {
-        (&self.pure_records, &self.k_ij)
+    fn records(&self) -> (&[PureRecord<PengRobinsonRecord>], Option<&Array2<f64>>) {
+        (&self.pure_records, Some(&self.k_ij))
     }
 }
 
@@ -294,7 +293,7 @@ mod tests {
         let propane = mixture[0].clone();
         let tc = propane.model_record.tc;
         let pc = propane.model_record.pc;
-        let parameters = PengRobinsonParameters::from_records(vec![propane], Array2::zeros((1, 1)));
+        let parameters = PengRobinsonParameters::new_pure(propane);
         let pr = Arc::new(PengRobinson::new(Arc::new(parameters)));
         let options = SolverOptions::new().verbosity(Verbosity::Iter);
         let cp = State::critical_point(&pr, None, None, options)?;

feos-core/src/joback.rs

@@ -72,7 +72,6 @@ pub struct JobackParameters {
     d: Array1<f64>,
     e: Array1<f64>,
     pure_records: Vec<PureRecord<JobackRecord>>,
-    binary_records: Array2<JobackBinaryRecord>,
 }
 
 impl Parameter for JobackParameters {
@@ -81,11 +80,10 @@ impl Parameter for JobackParameters {
 
     fn from_records(
         pure_records: Vec<PureRecord<Self::Pure>>,
-        _binary_records: Array2<Self::Binary>,
+        _binary_records: Option<Array2<Self::Binary>>,
     ) -> Self {
         let n = pure_records.len();
 
-        let binary_records = Array::from_elem((n, n), JobackBinaryRecord);
         let mut a = Array::zeros(n);
         let mut b = Array::zeros(n);
         let mut c = Array::zeros(n);
@@ -108,12 +106,11 @@ impl Parameter for JobackParameters {
             d,
             e,
             pure_records,
-            binary_records,
         }
     }
 
-    fn records(&self) -> (&[PureRecord<Self::Pure>], &Array2<Self::Binary>) {
-        (&self.pure_records, &self.binary_records)
+    fn records(&self) -> (&[PureRecord<Self::Pure>], Option<&Array2<Self::Binary>>) {
+        (&self.pure_records, None)
     }
 }
 
@@ -188,11 +185,7 @@ impl Components for Joback {
         component_list
             .iter()
             .for_each(|&i| records.push(self.parameters.pure_records[i].clone()));
-        let n = component_list.len();
-        Self::new(Arc::new(JobackParameters::from_records(
-            records,
-            Array::from_elem((n, n), JobackBinaryRecord),
-        )))
+        Self::new(Arc::new(JobackParameters::from_records(records, None)))
     }
 }
 

feos-core/src/lib.rs

@@ -202,7 +202,6 @@ mod tests {
     use crate::EosResult;
     use crate::StateBuilder;
     use approx::*;
-    use ndarray::Array2;
     use quantity::si::*;
     use std::sync::Arc;
 
@@ -248,7 +247,7 @@ mod tests {
     fn validate_residual_properties() -> EosResult<()> {
         let mixture = pure_record_vec();
         let propane = mixture[0].clone();
-        let parameters = PengRobinsonParameters::from_records(vec![propane], Array2::zeros((1, 1)));
+        let parameters = PengRobinsonParameters::new_pure(propane);
         let residual = Arc::new(PengRobinson::new(Arc::new(parameters)));
         let joback_parameters = Arc::new(JobackParameters::new_pure(PureRecord::new(
             Identifier::default(),

feos-core/src/parameter/mod.rs

@@ -35,13 +35,12 @@ where
     /// Creates parameters from records for pure substances and possibly binary parameters.
     fn from_records(
         pure_records: Vec<PureRecord<Self::Pure>>,
-        binary_records: Array2<Self::Binary>,
+        binary_records: Option<Array2<Self::Binary>>,
     ) -> Self;
 
     /// Creates parameters for a pure component from a pure record.
     fn new_pure(pure_record: PureRecord<Self::Pure>) -> Self {
-        let binary_record = Array2::from_elem([1, 1], Self::Binary::default());
-        Self::from_records(vec![pure_record], binary_record)
+        Self::from_records(vec![pure_record], None)
     }
 
     /// Creates parameters for a binary system from pure records and an optional
@@ -50,19 +49,31 @@ where
         pure_records: Vec<PureRecord<Self::Pure>>,
         binary_record: Option<Self::Binary>,
     ) -> Self {
-        let binary_record = Array2::from_shape_fn([2, 2], |(i, j)| {
-            if i == j {
-                Self::Binary::default()
-            } else {
-                binary_record.clone().unwrap_or_default()
-            }
+        let binary_record = binary_record.map(|br| {
+            Array2::from_shape_fn([2, 2], |(i, j)| {
+                if i == j {
+                    Self::Binary::default()
+                } else {
+                    br.clone()
+                }
+            })
         });
         Self::from_records(pure_records, binary_record)
     }
 
+    /// Creates parameters from model records with default values for the molar weight,
+    /// identifiers, and binary interaction parameters.
+    fn from_model_records(model_records: Vec<Self::Pure>) -> Self {
+        let pure_records = model_records
+            .into_iter()
+            .map(|r| PureRecord::new(Default::default(), Default::default(), r))
+            .collect();
+        Self::from_records(pure_records, None)
+    }
+
     /// Return the original pure and binary records that were used to construct the parameters.
     #[allow(clippy::type_complexity)]
-    fn records(&self) -> (&[PureRecord<Self::Pure>], &Array2<Self::Binary>);
+    fn records(&self) -> (&[PureRecord<Self::Pure>], Option<&Array2<Self::Binary>>);
 
     /// Helper function to build matrix from list of records in correct order.
     ///
@@ -73,7 +84,11 @@ where
         pure_records: &Vec<PureRecord<Self::Pure>>,
         binary_records: &[BinaryRecord<Identifier, Self::Binary>],
         search_option: IdentifierOption,
-    ) -> Array2<Self::Binary> {
+    ) -> Option<Array2<Self::Binary>> {
+        if binary_records.is_empty() {
+            return None;
+        }
+
         // Build Hashmap (id, id) -> BinaryRecord
         let binary_map: HashMap<(String, String), Self::Binary> = {
             binary_records
@@ -86,7 +101,7 @@ where
                 .collect()
         };
         let n = pure_records.len();
-        Array2::from_shape_fn([n, n], |(i, j)| {
+        Some(Array2::from_shape_fn([n, n], |(i, j)| {
             let id1 = pure_records[i]
                 .identifier
                 .as_string(search_option)
@@ -106,7 +121,7 @@ where
                 .or_else(|| binary_map.get(&(id2, id1)))
                 .cloned()
                 .unwrap_or_default()
-        })
+        }))
     }
 
     /// Creates parameters from substance information stored in json files.
@@ -250,7 +265,7 @@ where
             }
         }
 
-        Ok(Self::from_records(pure_records, binary_records))
+        Ok(Self::from_records(pure_records, Some(binary_records)))
     }
 
     /// Creates parameters from segment information stored in json files.
@@ -331,8 +346,10 @@ where
             .map(|&i| pure_records[i].clone())
             .collect();
         let n = component_list.len();
-        let binary_records = Array2::from_shape_fn([n, n], |(i, j)| {
-            binary_records[(component_list[i], component_list[j])].clone()
+        let binary_records = binary_records.map(|br| {
+            Array2::from_shape_fn([n, n], |(i, j)| {
+                br[(component_list[i], component_list[j])].clone()
+            })
         });
 
         Self::from_records(pure_records, binary_records)
@@ -484,24 +501,24 @@ mod test {
 
     struct MyParameter {
         pure_records: Vec<PureRecord<MyPureModel>>,
-        binary_records: Array2<MyBinaryModel>,
+        binary_records: Option<Array2<MyBinaryModel>>,
     }
 
     impl Parameter for MyParameter {
         type Pure = MyPureModel;
         type Binary = MyBinaryModel;
         fn from_records(
             pure_records: Vec<PureRecord<MyPureModel>>,
-            binary_records: Array2<MyBinaryModel>,
+            binary_records: Option<Array2<MyBinaryModel>>,
         ) -> Self {
             Self {
                 pure_records,
                 binary_records,
             }
         }
 
-        fn records(&self) -> (&[PureRecord<MyPureModel>], &Array2<MyBinaryModel>) {
-            (&self.pure_records, &self.binary_records)
+        fn records(&self) -> (&[PureRecord<MyPureModel>], Option<&Array2<MyBinaryModel>>) {
+            (&self.pure_records, self.binary_records.as_ref())
         }
     }
 
@@ -555,7 +572,7 @@ mod test {
 
         assert_eq!(p.pure_records[0].identifier.cas, Some("123-4-5".into()));
         assert_eq!(p.pure_records[1].identifier.cas, Some("678-9-1".into()));
-        assert_eq!(p.binary_records[[0, 1]].b, 12.0)
+        assert_eq!(p.binary_records.unwrap()[[0, 1]].b, 12.0)
     }
 
     #[test]
@@ -608,8 +625,9 @@ mod test {
 
         assert_eq!(p.pure_records[0].identifier.cas, Some("123-4-5".into()));
         assert_eq!(p.pure_records[1].identifier.cas, Some("678-9-1".into()));
-        assert_eq!(p.binary_records[[0, 1]], MyBinaryModel::default());
-        assert_eq!(p.binary_records[[0, 1]].b, 0.0)
+        let br = p.binary_records.as_ref().unwrap();
+        assert_eq!(br[[0, 1]], MyBinaryModel::default());
+        assert_eq!(br[[0, 1]].b, 0.0)
     }
 
     #[test]
@@ -672,11 +690,12 @@ mod test {
         assert_eq!(p.pure_records[0].identifier.cas, Some("000-0-0".into()));
         assert_eq!(p.pure_records[1].identifier.cas, Some("123-4-5".into()));
         assert_eq!(p.pure_records[2].identifier.cas, Some("678-9-1".into()));
-        assert_eq!(p.binary_records[[0, 1]], MyBinaryModel::default());
-        assert_eq!(p.binary_records[[1, 0]], MyBinaryModel::default());
-        assert_eq!(p.binary_records[[0, 2]], MyBinaryModel::default());
-        assert_eq!(p.binary_records[[2, 0]], MyBinaryModel::default());
-        assert_eq!(p.binary_records[[2, 1]].b, 12.0);
-        assert_eq!(p.binary_records[[1, 2]].b, 12.0);
+        let br = p.binary_records.as_ref().unwrap();
+        assert_eq!(br[[0, 1]], MyBinaryModel::default());
+        assert_eq!(br[[1, 0]], MyBinaryModel::default());
+        assert_eq!(br[[0, 2]], MyBinaryModel::default());
+        assert_eq!(br[[2, 0]], MyBinaryModel::default());
+        assert_eq!(br[[2, 1]].b, 12.0);
+        assert_eq!(br[[1, 2]].b, 12.0);
     }
 }

feos-core/src/python/cubic.rs

@@ -4,7 +4,6 @@ use crate::parameter::{
 };
 use crate::python::parameter::PyIdentifier;
 use crate::*;
-use ndarray::Array2;
 use numpy::{PyArray2, PyReadonlyArray2, ToPyArray};
 use pyo3::exceptions::PyTypeError;
 use pyo3::prelude::*;
@@ -56,7 +55,12 @@ impl_binary_record!();
 #[derive(Clone)]
 pub struct PyPengRobinsonParameters(pub Arc<PengRobinsonParameters>);
 
-impl_parameter!(PengRobinsonParameters, PyPengRobinsonParameters);
+impl_parameter!(
+    PengRobinsonParameters,
+    PyPengRobinsonParameters,
+    PyPengRobinsonRecord,
+    f64
+);
 
 #[pymethods]
 impl PyPengRobinsonParameters {

feos-core/src/python/joback.rs

@@ -7,7 +7,6 @@ use crate::{
     impl_binary_record, impl_json_handling, impl_parameter, impl_parameter_from_segments,
     impl_pure_record, impl_segment_record,
 };
-use ndarray::Array2;
 use numpy::{PyArray2, PyReadonlyArray2, ToPyArray};
 use pyo3::exceptions::PyTypeError;
 use pyo3::prelude::*;
@@ -62,6 +61,7 @@ impl_segment_record!(JobackRecord, PyJobackRecord);
 pub struct PyJobackBinaryRecord(pub JobackBinaryRecord);
 
 impl_binary_record!(JobackBinaryRecord, PyJobackBinaryRecord);
+
 /// Create a set of Joback parameters from records.
 ///
 /// Parameters
@@ -83,7 +83,12 @@ impl_binary_record!(JobackBinaryRecord, PyJobackBinaryRecord);
 #[derive(Clone)]
 pub struct PyJobackParameters(pub Arc<JobackParameters>);
 
-impl_parameter!(JobackParameters, PyJobackParameters);
+impl_parameter!(
+    JobackParameters,
+    PyJobackParameters,
+    PyJobackRecord,
+    PyJobackBinaryRecord
+);
 impl_parameter_from_segments!(JobackParameters, PyJobackParameters);
 
 #[pymethods]