Check number of polar and associating segments in HomoGC method

Cargo.toml

@@ -14,15 +14,17 @@ exclude = ["/.github/*", "*.ipynb", "/examples/*", "/parameters/*"]
 
 [dependencies]
 quantity = "0.5"
-feos-core = "0.2"
-feos-dft = "0.2"
+feos-core = { path = "../feos-core" }
+feos-dft = { path = "../feos-dft", optional = true }
+feos-campd = { path = "../feos-campd", optional = true }
 num-dual = "0.5"
 num-traits = "0.2"
 ndarray = { version = "0.15", features=["approx"] }
 serde = { version = "1.0", features = ["derive"] }
 serde_json = "1.0"
 indexmap = "1.8"
 lazy_static = "1.4"
+conv = "0.3"
 numpy = { version = "0.16", optional = true }
 pyo3 = { version = "0.16", optional = true }
 
@@ -31,4 +33,6 @@ approx = "0.4"
 
 [features]
 default = []
-python = ["pyo3", "quantity/python", "feos-core/python", "feos-dft/python", "numpy"]
+dft = ["feos-dft"]
+campd = ["feos-campd"]
+python = ["pyo3", "numpy", "feos-core/python", "feos-dft?/python"]

build_wheel/Cargo.toml

@@ -10,8 +10,8 @@ crate-type = ["cdylib"]
 
 [dependencies]
 quantity = "0.5"
-feos-core = "0.2"
-feos-dft = "0.2"
-feos-pcsaft = { path = "..", features = ["python"] }
+feos-core = { path = "../../feos-core" }
+feos-dft = { path = "../../feos-dft" }
+feos-pcsaft = { path = "..", features = ["dft", "python"] }
 pyo3 = { version = "0.16", features = ["extension-module", "abi3", "abi3-py37"] }
 numpy = "0.16"

build_wheel/src/dft.rs

@@ -8,7 +8,7 @@ use feos_dft::*;
 use feos_pcsaft::python::*;
 use feos_pcsaft::{PcSaftFunctional, PcSaftOptions};
 use numpy::*;
-use pyo3::exceptions::PyValueError;
+use pyo3::exceptions::{PyIndexError, PyValueError};
 use pyo3::prelude::*;
 use quantity::python::*;
 use quantity::si::*;

build_wheel/src/eos.rs

@@ -3,7 +3,7 @@ use feos_pcsaft::python::PyPcSaftParameters;
 use feos_pcsaft::{PcSaft, PcSaftOptions};
 use numpy::convert::ToPyArray;
 use numpy::{PyArray1, PyArray2};
-use pyo3::exceptions::PyValueError;
+use pyo3::exceptions::{PyIndexError, PyValueError};
 use pyo3::prelude::*;
 use quantity::python::*;
 use quantity::si::*;

src/campd.rs

@@ -0,0 +1,48 @@
+use crate::eos::PcSaft;
+use crate::parameters::{PcSaftParameters, PcSaftRecord};
+use feos_campd::{PropertyModel, SegmentAndBondCount};
+use feos_core::joback::JobackRecord;
+use feos_core::parameter::{Parameter, ParameterError, PureRecord, SegmentRecord};
+use serde::{Deserialize, Serialize};
+use std::path::Path;
+use std::rc::Rc;
+
+#[derive(Clone, Serialize, Deserialize)]
+pub struct PcSaftPropertyModel(Vec<SegmentRecord<PcSaftRecord, JobackRecord>>);
+
+impl PcSaftPropertyModel {
+    pub fn new<P: AsRef<Path>>(file: P) -> Result<Self, ParameterError> {
+        Ok(Self(SegmentRecord::from_json(file)?))
+    }
+}
+
+impl PropertyModel<SegmentAndBondCount> for PcSaftPropertyModel {
+    type Eos = PcSaft;
+
+    fn build_eos(&self, chemical_record: SegmentAndBondCount) -> Result<PcSaft, ParameterError> {
+        Ok(PcSaft::new(Rc::new(PcSaftParameters::from_segments(
+            vec![chemical_record],
+            self.0.clone(),
+            None,
+        )?)))
+    }
+}
+
+#[derive(Clone, Serialize, Deserialize)]
+pub struct PcSaftFixedPropertyModel(PureRecord<PcSaftRecord, JobackRecord>);
+
+impl PcSaftFixedPropertyModel {
+    pub fn new(pure_record: PureRecord<PcSaftRecord, JobackRecord>) -> Self {
+        Self(pure_record)
+    }
+}
+
+impl PropertyModel<()> for PcSaftFixedPropertyModel {
+    type Eos = PcSaft;
+
+    fn build_eos(&self, _: ()) -> Result<PcSaft, ParameterError> {
+        Ok(PcSaft::new(Rc::new(PcSaftParameters::new_pure(
+            self.0.clone(),
+        ))))
+    }
+}

src/lib.rs

@@ -1,9 +1,16 @@
 #![warn(clippy::all)]
 #![allow(clippy::too_many_arguments)]
+
+#[cfg(feature = "campd")]
+mod campd;
+#[cfg(feature = "dft")]
 mod dft;
 mod eos;
 mod parameters;
 
+#[cfg(feature = "campd")]
+pub use campd::{PcSaftFixedPropertyModel, PcSaftPropertyModel};
+#[cfg(feature = "dft")]
 pub use dft::PcSaftFunctional;
 pub use eos::{PcSaft, PcSaftOptions};
 pub use parameters::{PcSaftParameters, PcSaftRecord};

src/parameters.rs

@@ -1,5 +1,8 @@
+use conv::ValueInto;
 use feos_core::joback::JobackRecord;
-use feos_core::parameter::{FromSegments, FromSegmentsBinary, Parameter, PureRecord};
+use feos_core::parameter::{
+    FromSegments, FromSegmentsBinary, Parameter, ParameterError, PureRecord,
+};
 use ndarray::{Array, Array1, Array2};
 use num_traits::Zero;
 use quantity::si::{JOULE, KB, KELVIN};
@@ -45,41 +48,41 @@ pub struct PcSaftRecord {
     pub thermal_conductivity: Option<[f64; 4]>,
 }
 
-impl FromSegments for PcSaftRecord {
-    fn from_segments(segments: &[(Self, f64)]) -> Self {
+impl FromSegments<f64> for PcSaftRecord {
+    fn from_segments(segments: &[(Self, f64)]) -> Result<Self, ParameterError> {
         let mut m = 0.0;
         let mut sigma3 = 0.0;
         let mut epsilon_k = 0.0;
 
         segments.iter().for_each(|(s, n)| {
-            m += s.m * *n;
-            sigma3 += s.m * s.sigma.powi(3) * *n;
-            epsilon_k += s.m * s.epsilon_k * *n;
+            m += s.m * n;
+            sigma3 += s.m * s.sigma.powi(3) * n;
+            epsilon_k += s.m * s.epsilon_k * n;
         });
 
         let q = segments
             .iter()
-            .filter_map(|(s, n)| s.q.map(|q| q * *n))
+            .filter_map(|(s, n)| s.q.map(|q| q * n))
             .reduce(|a, b| a + b);
         let mu = segments
             .iter()
-            .filter_map(|(s, n)| s.mu.map(|mu| mu * *n))
+            .filter_map(|(s, n)| s.mu.map(|mu| mu * n))
             .reduce(|a, b| a + b);
         let kappa_ab = segments
             .iter()
-            .filter_map(|(s, n)| s.kappa_ab.map(|k| k * *n))
+            .filter_map(|(s, n)| s.kappa_ab.map(|k| k * n))
             .reduce(|a, b| a + b);
         let epsilon_k_ab = segments
             .iter()
-            .filter_map(|(s, n)| s.epsilon_k_ab.map(|e| e * *n))
+            .filter_map(|(s, n)| s.epsilon_k_ab.map(|e| e * n))
             .reduce(|a, b| a + b);
         let na = segments
             .iter()
-            .filter_map(|(s, n)| s.na.map(|na| na * *n))
+            .filter_map(|(s, n)| s.na.map(|na| na * n))
             .reduce(|a, b| a + b);
         let nb = segments
             .iter()
-            .filter_map(|(s, n)| s.nb.map(|nb| nb * *n))
+            .filter_map(|(s, n)| s.nb.map(|nb| nb * n))
             .reduce(|a, b| a + b);
 
         // entropy scaling
@@ -110,19 +113,19 @@ impl FromSegments for PcSaftRecord {
 
         let n_t = segments.iter().fold(0.0, |acc, (_, n)| acc + n);
         segments.iter().for_each(|(s, n)| {
-            let s3 = s.m * s.sigma.powi(3) * *n;
+            let s3 = s.m * s.sigma.powi(3) * n;
             if let Some(p) = viscosity.as_mut() {
                 let [a, b, c, d] = s.viscosity.unwrap();
                 p[0] += s3 * a;
                 p[1] += s3 * b / sigma3.powf(0.45);
-                p[2] += *n * c;
-                p[3] += *n * d;
+                p[2] += n * c;
+                p[3] += n * d;
             }
             if let Some(p) = thermal_conductivity.as_mut() {
                 let [a, b, c, d] = s.thermal_conductivity.unwrap();
-                p[0] += *n * a;
-                p[1] += *n * b;
-                p[2] += *n * c;
+                p[0] += n * a;
+                p[1] += n * b;
+                p[2] += n * c;
                 p[3] += n_t * d;
             }
             // if let Some(p) = diffusion.as_mut() {
@@ -136,7 +139,7 @@ impl FromSegments for PcSaftRecord {
         // correction due to difference in Chapman-Enskog reference between GC and regular formulation.
         viscosity = viscosity.map(|v| [v[0] - 0.5 * m.ln(), v[1], v[2], v[3]]);
 
-        Self {
+        Ok(Self {
             m,
             sigma: (sigma3 / m).cbrt(),
             epsilon_k: epsilon_k / m,
@@ -149,7 +152,49 @@ impl FromSegments for PcSaftRecord {
             viscosity,
             diffusion,
             thermal_conductivity,
+        })
+    }
+}
+
+impl FromSegments<usize> for PcSaftRecord {
+    fn from_segments(segments: &[(Self, usize)]) -> Result<Self, ParameterError> {
+        // We do not allow more than a single segment for q, mu, kappa_ab, epsilon_k_ab
+        let quadpole_comps = segments.iter().filter_map(|(s, _)| s.q).count();
+        if quadpole_comps > 1 {
+            return Err(ParameterError::IncompatibleParameters(format!(
+                "{quadpole_comps} segments with quadrupole moment."
+            )));
+        };
+        let dipole_comps = segments.iter().filter_map(|(s, _)| s.mu).count();
+        if dipole_comps > 1 {
+            return Err(ParameterError::IncompatibleParameters(format!(
+                "{dipole_comps} segment with dipole moment."
+            )));
+        };
+        let faulty_assoc_comps = segments
+            .iter()
+            .filter_map(|(s, _)| s.kappa_ab.xor(s.epsilon_k_ab))
+            .count();
+        if faulty_assoc_comps > 0 {
+            return Err(ParameterError::IncompatibleParameters(format!(
+                "Incorrectly specified association sites on {faulty_assoc_comps} segment(s)"
+            )));
         }
+        let assoc_comps = segments
+            .iter()
+            .filter_map(|(s, _)| s.kappa_ab.and(s.epsilon_k_ab))
+            .count();
+        if assoc_comps > 1 {
+            return Err(ParameterError::IncompatibleParameters(format!(
+                "{assoc_comps} segments with association sites."
+            )));
+        }
+        let segments: Vec<_> = segments
+            .iter()
+            .cloned()
+            .map(|(s, c)| (s, c as f64))
+            .collect();
+        Self::from_segments(&segments)
     }
 }
 
@@ -232,10 +277,25 @@ impl From<f64> for PcSaftBinaryRecord {
     }
 }
 
-impl FromSegmentsBinary for PcSaftBinaryRecord {
-    fn from_segments_binary(segments: &[(Self, f64, f64)]) -> Self {
-        let k_ij = segments.iter().map(|(br, n1, n2)| br.k_ij * n1 * n2).sum();
-        Self { k_ij }
+impl From<PcSaftBinaryRecord> for f64 {
+    fn from(binary_record: PcSaftBinaryRecord) -> Self {
+        binary_record.k_ij
+    }
+}
+
+impl<T: Copy + ValueInto<f64>> FromSegmentsBinary<T> for PcSaftBinaryRecord {
+    fn from_segments_binary(segments: &[(Self, T, T)]) -> Result<Self, ParameterError> {
+        let k_ij = segments
+            .iter()
+            .map(|(br, n1, n2)| br.k_ij * (*n1).value_into().unwrap() * (*n2).value_into().unwrap())
+            .sum();
+        Ok(Self { k_ij })
+    }
+}
+
+impl std::fmt::Display for PcSaftBinaryRecord {
+    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
+        write!(f, "PcSaftBinaryRecord(k_ij={})", self.k_ij)
     }
 }
 

src/python.rs

@@ -1,15 +1,17 @@
 use crate::eos::polar::DQVariants;
-use crate::parameters::{PcSaftParameters, PcSaftRecord};
+use crate::parameters::{PcSaftBinaryRecord, PcSaftParameters, PcSaftRecord};
 use feos_core::joback::JobackRecord;
 use feos_core::parameter::{
-    BinaryRecord, IdentifierOption, Parameter, ParameterError, PureRecord, SegmentRecord,
+    BinaryRecord, Identifier, IdentifierOption, Parameter, ParameterError, PureRecord,
+    SegmentRecord,
 };
 use feos_core::python::joback::PyJobackRecord;
 use feos_core::python::parameter::{PyBinarySegmentRecord, PyChemicalRecord, PyIdentifier};
 use feos_core::*;
-use numpy::{PyArray2, ToPyArray};
+use numpy::{PyArray2, PyReadonlyArray2, ToPyArray};
+use pyo3::exceptions::PyTypeError;
 use pyo3::prelude::*;
-use std::convert::TryFrom;
+use std::convert::{TryFrom, TryInto};
 use std::rc::Rc;
 
 impl From<&str> for DQVariants {
@@ -133,6 +135,14 @@ impl_json_handling!(PyPcSaftRecord);
 impl_pure_record!(PcSaftRecord, PyPcSaftRecord, JobackRecord, PyJobackRecord);
 impl_segment_record!(PcSaftRecord, PyPcSaftRecord, JobackRecord, PyJobackRecord);
 
+#[pyclass(name = "PcSaftBinaryRecord", unsendable)]
+#[pyo3(
+    text_signature = "(pure_records, binary_records=None, substances=None, search_option='Name')"
+)]
+#[derive(Clone)]
+pub struct PyPcSaftBinaryRecord(PcSaftBinaryRecord);
+impl_binary_record!(PcSaftBinaryRecord, PyPcSaftBinaryRecord);
+
 /// Create a set of PC-SAFT parameters from records.
 ///
 /// Parameters
@@ -163,15 +173,6 @@ impl_parameter_from_segments!(PcSaftParameters, PyPcSaftParameters);
 
 #[pymethods]
 impl PyPcSaftParameters {
-    #[getter]
-    fn get_pure_records(&self) -> Vec<PyPureRecord> {
-        self.0
-            .pure_records
-            .iter()
-            .map(|r| PyPureRecord(r.clone()))
-            .collect()
-    }
-
     #[getter]
     fn get_k_ij<'py>(&self, py: Python<'py>) -> &'py PyArray2<f64> {
         self.0.k_ij.view().to_pyarray(py)

tests/dft.rs

@@ -1,4 +1,5 @@
 #![allow(clippy::excessive_precision)]
+#![cfg(feature = "dft")]
 use approx::assert_relative_eq;
 use feos_core::parameter::{IdentifierOption, Parameter};
 use feos_core::{Contributions, PhaseEquilibrium, State};