Breaking changes for Release v0.10.0

.github/workflows/test.yml

@@ -4,7 +4,7 @@ on:
   push:
     branches: [main]
   pull_request:
-    branches: [main]
+    branches: [main, development]
 
 env:
   CARGO_TERM_COLOR: always

.github/workflows/wheels.yml

@@ -3,7 +3,7 @@ on:
   push:
     branches: [main]
   pull_request:
-    branches: [main]
+    branches: [main, development]
 jobs:
   linux:
     runs-on: ubuntu-latest

CHANGELOG.md

@@ -4,10 +4,37 @@ All notable changes to this project will be documented in this file.
 The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
+## [Breaking]
+### Added
+- Extended tp-flash algorithm to static numbers of components and enabled automatic differentiation for binary systems. [#336](https://github.com/feos-org/feos/pull/336)
+- Rewrote `PhaseEquilibrium::pure_p` to mirror `pure_t` and enabled automatic differentiation. [#337](https://github.com/feos-org/feos/pull/337)
+- Added `boiling_temperature` to the list of properties for parallel evaluations of gradients. [#337](https://github.com/feos-org/feos/pull/337)
+- Added the `Composition` trait to allow more flexibility in the creation of states and phase equilibria. [#330](https://github.com/feos-org/feos/pull/330)
+- Added `PhaseEquilibrium::ph_flash` and `PhaseEquilibrium::ps_flash`. [#338](https://github.com/feos-org/feos/pull/338)
+- Added getters for `vapor_phase_fraction`, `molar_enthalpy`, `molar_entropy`, `total_moles`, `enthalpy`, and `entropy` to `PhaseEquilibrium`. [#338](https://github.com/feos-org/feos/pull/338)
+- Added `PropertyAD` trait in `feos_core::ad` with one struct per property for uniform evaluation with or without parameter derivatives, including parallel variants. [#358](https://github.com/feos-org/feos/pull/358)
+- Added `feos_core::ad::dataset` module with `PureDataset` and `BinaryDataset` types, constructible from records, CSV files, or readers, for use in parameter fits. [#358](https://github.com/feos-org/feos/pull/358)
+- Exposed `Property`, `PureDataset`, and `BinaryDataset` in `py-feos`. [#358](https://github.com/feos-org/feos/pull/358)
+
+### Changed
+- Removed any assumptions about the total number of moles in a `State` or `PhaseEquilibrium`. Evaluating extensive properties now returns a `Result`. [#330](https://github.com/feos-org/feos/pull/330)
+- Redesigned the `IdealGas` trait and added `IdealGasAD` in analogy to `ResidualDyn` and `Residual`. [#330](https://github.com/feos-org/feos/pull/330)
+- Replaced the `PropertiesAD` blanket-impl trait with per-property `PropertyAD` types. [#358](https://github.com/feos-org/feos/pull/358)
+- Replaced `ParametersAD::named_derivatives` with a `build` constructor and `seed_derivatives(&values, names)`. [#358](https://github.com/feos-org/feos/pull/358)
+- Removed the per-property `*_derivatives` free functions in Python in favour of static methods on the new `Property` class. [#358](https://github.com/feos-org/feos/pull/358)
+
+### Removed
+- Removed the `StateBuilder` struct, because it is mostly obsolete with the addition of the `Composition` trait. [#330](https://github.com/feos-org/feos/pull/330)
+
+### Packaging
+- Updated `quantity` dependency to 0.13 and removed the `typenum` dependency. [#328](https://github.com/feos-org/feos/pull/328)
+- Added `csv` as a `feos-core` dependency for the new dataset module. [#358](https://github.com/feos-org/feos/pull/358)
+
 ## [Unreleased]
 ### Added
 - Add Rayon global thread pool control via `FEOS_MAX_THREADS` and `set_num_threads()`/ `get_num_threads()` to Python. [#346](https://github.com/feos-org/feos/pull/346)
 - Added DIPPR107 parameterization for ideal gas heat capacities of Burkhardt et al. [#344](https://github.com/feos-org/feos/pull/344)
+- Implemented `IdealGasAD<D>` for `Dippr`. [#357](https://github.com/feos-org/feos/pull/357)
 
 ### Fixed
 - Fixed the calculation of temperature and pressure derivatives of dew and bubble points. [#347](https://github.com/feos-org/feos/pull/347)

Cargo.toml

@@ -22,7 +22,7 @@ keywords = [
 categories = ["science"]
 
 [workspace.dependencies]
-quantity = "0.12"
+quantity = "0.13"
 num-dual = "0.13"
 ndarray = "0.17"
 nalgebra = "0.34"
@@ -33,7 +33,6 @@ serde = "1.0"
 serde_json = "1.0"
 indexmap = "2.0"
 itertools = "0.14"
-typenum = "1.16"
 rayon = "1.11"
 petgraph = "0.8"
 rustdct = "0.7"
@@ -43,6 +42,7 @@ gauss-quad = "0.2"
 approx = "0.5"
 criterion = "0.8"
 paste = "1.0"
+csv = "1.0"
 
 feos-core = { version = "0.9", path = "crates/feos-core" }
 feos-dft = { version = "0.9", path = "crates/feos-dft" }

crates/feos-core/Cargo.toml

@@ -25,7 +25,7 @@ serde = { workspace = true, features = ["derive"] }
 serde_json = { workspace = true, features = ["preserve_order"] }
 indexmap = { workspace = true, features = ["serde"] }
 rayon = { workspace = true, optional = true }
-typenum = { workspace = true }
+csv = { workspace = true }
 itertools = { workspace = true }
 
 [dev-dependencies]

crates/feos-core/src/ad/dataset/binary.rs

@@ -0,0 +1,304 @@
+use std::{io, path::Path};
+
+use nalgebra::U2;
+use ndarray::{Array1, Array2, ArrayView1, ArrayView2};
+use serde::{Deserialize, Serialize};
+
+use crate::Residual;
+use crate::ad::Gradient;
+use crate::ad::properties::{BubblePointPressure, DewPointPressure, PropertyAD};
+
+use super::{Dataset, DatasetAD, DatasetRecord, DatasetStorage, ParametersAD};
+
+/// The pressure column doubles as the initial guess passed to the VLE solver.
+#[derive(Deserialize, Serialize)]
+pub struct BubblePointRecord {
+    pub temperature_k: f64,
+    pub liquid_molefrac_1: f64,
+    pub bubble_pressure_pa: f64,
+}
+
+impl DatasetRecord for BubblePointRecord {
+    const N_INPUTS: usize = 3;
+
+    fn input(&self, column: usize) -> f64 {
+        match column {
+            0 => self.temperature_k,
+            1 => self.liquid_molefrac_1,
+            2 => self.bubble_pressure_pa,
+            _ => unreachable!("invalid bubble point input column"),
+        }
+    }
+
+    fn target(&self) -> f64 {
+        self.bubble_pressure_pa
+    }
+}
+
+/// The pressure column doubles as the initial guess passed to the VLE solver.
+#[derive(Deserialize, Serialize)]
+pub struct DewPointRecord {
+    pub temperature_k: f64,
+    pub vapor_molefrac_1: f64,
+    pub dew_pressure_pa: f64,
+}
+
+impl DatasetRecord for DewPointRecord {
+    const N_INPUTS: usize = 3;
+
+    fn input(&self, column: usize) -> f64 {
+        match column {
+            0 => self.temperature_k,
+            1 => self.vapor_molefrac_1,
+            2 => self.dew_pressure_pa,
+            _ => unreachable!("invalid dew point input column"),
+        }
+    }
+
+    fn target(&self) -> f64 {
+        self.dew_pressure_pa
+    }
+}
+
+/// Expand a list of binary-mixture property entries into:
+/// - the [`BinaryProperty`] enum, metadata and dispatch methods,
+/// - constructors,
+/// - [`BinaryDataset::from_csv`] and [`BinaryDataset::from_reader`] match arms.
+macro_rules! binary_properties {
+    ($(
+        $variant:ident {
+            record:       $record:ty,
+            property:     $prop:ty,
+            default_name: $default:expr,
+            input_names:  $inputs:expr,
+            target_name:  $target:expr,
+            constructor:  $ctor:ident,
+        }
+    ),* $(,)?) => {
+        /// Binary-mixture properties.
+        #[derive(Debug, Clone, Copy, PartialEq, Eq)]
+        pub enum BinaryProperty {
+            $($variant,)*
+        }
+
+        impl BinaryProperty {
+            pub fn default_name(self) -> &'static str {
+                match self { $(Self::$variant => $default,)* }
+            }
+
+            pub fn input_names(self) -> &'static [&'static str] {
+                match self { $(Self::$variant => $inputs,)* }
+            }
+
+            pub fn target_name(self) -> &'static str {
+                match self { $(Self::$variant => $target,)* }
+            }
+
+            fn evaluate_ad<T: ParametersAD<U2>, const P: usize>(
+                self,
+                names: [String; P],
+                parameters: &[f64],
+                inputs: ArrayView2<f64>,
+            ) -> (Array1<f64>, Array2<f64>, Array1<bool>)
+            where
+                T::Lifted<Gradient<P>>: Sync,
+            {
+                match self {
+                    $(Self::$variant => <$prop>::evaluate_parallel_derivatives::<T, P>(names, parameters, inputs),)*
+                }
+            }
+
+            fn evaluate<E>(self, eos: &E, inputs: ArrayView2<f64>) -> (Array1<f64>, Array1<bool>)
+            where
+                E: Residual + Sync,
+            {
+                match self {
+                    $(Self::$variant => <$prop>::evaluate_parallel(eos, inputs),)*
+                }
+            }
+        }
+
+        impl BinaryDataset {
+            $(
+                pub fn $ctor(records: Vec<$record>) -> Self {
+                    Self {
+                        property: BinaryProperty::$variant,
+                        storage: DatasetStorage::from_records(records),
+                    }
+                }
+            )*
+
+            pub fn from_csv(property: BinaryProperty, path: &Path) -> Result<Self, csv::Error> {
+                let storage = match property {
+                    $(BinaryProperty::$variant => DatasetStorage::from_csv::<$record>(path)?,)*
+                };
+                Ok(Self { property, storage })
+            }
+
+            pub fn from_reader(
+                property: BinaryProperty,
+                reader: impl io::Read,
+            ) -> Result<Self, csv::Error> {
+                let storage = match property {
+                    $(BinaryProperty::$variant => DatasetStorage::from_reader::<$record>(reader)?,)*
+                };
+                Ok(Self { property, storage })
+            }
+        }
+    };
+}
+
+binary_properties! {
+    BubblePointPressure {
+        record:       BubblePointRecord,
+        property:     BubblePointPressure,
+        default_name: "bubble point pressure",
+        input_names:  &["temperature_k", "liquid_molefrac_1"],
+        target_name:  "bubble_pressure_pa",
+        constructor:  bubble_point_pressure,
+    },
+    DewPointPressure {
+        record:       DewPointRecord,
+        property:     DewPointPressure,
+        default_name: "dew point pressure",
+        input_names:  &["temperature_k", "vapor_molefrac_1"],
+        target_name:  "dew_pressure_pa",
+        constructor:  dew_point_pressure,
+    },
+}
+
+/// Binary-mixture dataset: shared data storage plus a property tag.
+#[derive(Clone)]
+pub struct BinaryDataset {
+    property: BinaryProperty,
+    storage: DatasetStorage,
+}
+
+impl BinaryDataset {
+    pub fn with_name(mut self, name: impl Into<String>) -> Self {
+        self.storage.set_name(name.into());
+        self
+    }
+
+    pub fn property(&self) -> BinaryProperty {
+        self.property
+    }
+
+    pub fn inputs(&self) -> ArrayView2<'_, f64> {
+        self.storage.inputs()
+    }
+
+    pub fn target(&self) -> ArrayView1<'_, f64> {
+        self.storage.target()
+    }
+
+    pub fn name(&self) -> &str {
+        self.storage.name().unwrap_or(self.property.default_name())
+    }
+
+    pub fn input_names(&self) -> &'static [&'static str] {
+        self.property.input_names()
+    }
+
+    pub fn target_name(&self) -> &'static str {
+        self.property.target_name()
+    }
+}
+
+impl Dataset for BinaryDataset {
+    fn inputs(&self) -> ArrayView2<'_, f64> {
+        self.inputs()
+    }
+
+    fn target(&self) -> ArrayView1<'_, f64> {
+        self.target()
+    }
+
+    fn name(&self) -> &str {
+        self.name()
+    }
+
+    fn input_names(&self) -> &'static [&'static str] {
+        self.input_names()
+    }
+
+    fn target_name(&self) -> &'static str {
+        self.target_name()
+    }
+
+    fn evaluate<E: Residual + Sync>(&self, eos: &E) -> (Array1<f64>, Array1<bool>) {
+        self.property.evaluate(eos, self.inputs())
+    }
+}
+
+impl DatasetAD<2> for BinaryDataset {
+    fn evaluate_ad_const<T: ParametersAD<U2>, const P: usize>(
+        &self,
+        names: [String; P],
+        parameters: &[f64],
+        inputs: ArrayView2<f64>,
+    ) -> (Array1<f64>, Array2<f64>, Array1<bool>)
+    where
+        T::Lifted<Gradient<P>>: Sync,
+    {
+        self.property.evaluate_ad::<T, P>(names, parameters, inputs)
+    }
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use std::io::Cursor;
+
+    fn csv(s: &str) -> Cursor<&[u8]> {
+        Cursor::new(s.as_bytes())
+    }
+
+    #[test]
+    fn bubble_point_from_reader() {
+        let data = "\
+temperature_k,liquid_molefrac_1,bubble_pressure_pa
+300.0,0.3,500000.0
+320.0,0.5,800000.0
+";
+        let ds =
+            BinaryDataset::from_reader(BinaryProperty::BubblePointPressure, csv(data)).unwrap();
+
+        assert_eq!(ds.inputs().ncols(), 3);
+        assert_eq!(ds.inputs().nrows(), 2);
+        assert_eq!(ds.inputs()[[0, 0]], 300.0);
+        assert_eq!(ds.inputs()[[0, 1]], 0.3);
+        assert_eq!(ds.inputs()[[0, 2]], 500000.0);
+        assert_eq!(ds.target()[0], 500000.0);
+        assert_eq!(ds.target()[1], 800000.0);
+        assert_eq!(ds.name(), "bubble point pressure");
+    }
+
+    #[test]
+    fn dew_point_from_reader() {
+        let data = "\
+temperature_k,vapor_molefrac_1,dew_pressure_pa
+310.0,0.7,400000.0
+330.0,0.9,700000.0
+";
+        let ds = BinaryDataset::from_reader(BinaryProperty::DewPointPressure, csv(data)).unwrap();
+
+        assert_eq!(ds.inputs().ncols(), 3);
+        assert_eq!(ds.inputs()[[0, 0]], 310.0);
+        assert_eq!(ds.inputs()[[0, 1]], 0.7);
+        assert_eq!(ds.inputs()[[0, 2]], 400000.0);
+        assert_eq!(ds.target()[1], 700000.0);
+        assert_eq!(ds.name(), "dew point pressure");
+    }
+
+    #[test]
+    fn dew_point_pressure_doubles_as_initial_guess() {
+        let records = vec![DewPointRecord {
+            temperature_k: 310.0,
+            vapor_molefrac_1: 0.7,
+            dew_pressure_pa: 400000.0,
+        }];
+        let ds = BinaryDataset::dew_point_pressure(records);
+        assert_eq!(ds.inputs()[[0, 2]], ds.target()[0]);
+    }
+}