More rigorous treatment of extensive states

feos-org · prehner · Feb 14, 2026 · Jan 27, 2026 · Jan 27, 2026 · Feb 12, 2026
commit 33ececbfcb8546b86cfdfad3e3074f47ad9c29bf
diff --git a/crates/feos-core/src/ad/mod.rs b/crates/feos-core/src/ad/mod.rs
@@ -57,8 +57,8 @@ pub trait PropertiesAD {
             let v2 = Gradient::from(v2);
             let x = Self::pure_molefracs();
 
-            let a1 = self.residual_molar_helmholtz_energy(t, v1, &x);
-            let a2 = self.residual_molar_helmholtz_energy(t, v2, &x);
+            let a1 = self.residual_helmholtz_energy(t, v1, &x);
+            let a2 = self.residual_helmholtz_energy(t, v2, &x);
             (a1, a2)
         };
 
@@ -93,7 +93,7 @@ pub trait PropertiesAD {
             let residual_entropy = |v| {
                 let (a, s) = first_derivative(
                     partial2(
-                        |t, &v, x| self.lift().residual_molar_helmholtz_energy(t, v, x),
+                        |t, &v, x| self.lift().residual_helmholtz_energy(t, v, x),
                         &v,
                         &x,
                     ),
@@ -248,16 +248,16 @@ pub trait PropertiesAD {
             let y = y.map(Gradient::from);
             let x = liquid_molefracs.map(Gradient::from);
 
-            let a_v = self.residual_molar_helmholtz_energy(t, v_v, &y);
+            let a_v = self.residual_helmholtz_energy(t, v_v, &y);
             let (p_l, mu_res_l, dp_l, dmu_l) = self.dmu_dv(t, v_l, &x);
             let vi_l = dmu_l / dp_l;
             let v_l = vi_l.dot(&y);
             let a_l = (mu_res_l - vi_l * p_l).dot(&y);
             (a_l, a_v, v_l, v_v)
         };
-        let rho_l = vle.liquid().partial_density.to_reduced();
+        let rho_l = vle.liquid().partial_density().to_reduced();
         let rho_l = [rho_l[0], rho_l[1]];
-        let rho_v = vle.vapor().partial_density.to_reduced();
+        let rho_v = vle.vapor().partial_density().to_reduced();
         let rho_v = [rho_v[0], rho_v[1]];
         let p = -(a_v - a_l
             + t * (y[0] * (rho_v[0] / rho_l[0]).ln() + y[1] * (rho_v[1] / rho_l[1]).ln() - 1.0))
@@ -298,16 +298,16 @@ pub trait PropertiesAD {
             let x = x.map(Gradient::from);
             let y = vapor_molefracs.map(Gradient::from);
 
-            let a_l = self.residual_molar_helmholtz_energy(t, v_l, &x);
+            let a_l = self.residual_helmholtz_energy(t, v_l, &x);
             let (p_v, mu_res_v, dp_v, dmu_v) = self.dmu_dv(t, v_v, &y);
             let vi_v = dmu_v / dp_v;
             let v_v = vi_v.dot(&x);
             let a_v = (mu_res_v - vi_v * p_v).dot(&x);
             (a_l, a_v, v_l, v_v)
         };
-        let rho_l = vle.liquid().partial_density.to_reduced();
+        let rho_l = vle.liquid().partial_density().to_reduced();
         let rho_l = [rho_l[0], rho_l[1]];
-        let rho_v = vle.vapor().partial_density.to_reduced();
+        let rho_v = vle.vapor().partial_density().to_reduced();
         let rho_v = [rho_v[0], rho_v[1]];
         let p = -(a_l - a_v
             + t * (x[0] * (rho_l[0] / rho_v[0]).ln() + x[1] * (rho_l[1] / rho_v[1]).ln() - 1.0))

diff --git a/crates/feos-core/src/cubic.rs b/crates/feos-core/src/cubic.rs
@@ -221,7 +221,7 @@ mod tests {
         let parameters = PengRobinsonParameters::new_pure(propane)?;
         let pr = PengRobinson::new(parameters);
         let options = SolverOptions::new().verbosity(Verbosity::Iter);
-        let cp = State::critical_point(&&pr, None, None, None, options)?;
+        let cp = State::critical_point(&&pr, (), None, None, options)?;
         println!("{} {}", cp.temperature, cp.pressure(Contributions::Total));
         assert_relative_eq!(cp.temperature, tc * KELVIN, max_relative = 1e-4);
         assert_relative_eq!(

diff --git a/crates/feos-core/src/density_iteration.rs b/crates/feos-core/src/density_iteration.rs
@@ -87,7 +87,7 @@ where
     let (a_res, da_res) = first_derivative(
         |molar_volume| {
             eos.lift()
-                .residual_molar_helmholtz_energy(t, molar_volume, &x)
+                .residual_helmholtz_energy(t, molar_volume, &x)
         },
         molar_volume,
     );

diff --git a/crates/feos-core/src/equation_of_state/mod.rs b/crates/feos-core/src/equation_of_state/mod.rs
@@ -1,9 +1,10 @@
-use crate::{ReferenceSystem, StateHD};
+use crate::ReferenceSystem;
+use crate::state::StateHD;
 use nalgebra::{
     Const, DVector, DefaultAllocator, Dim, Dyn, OVector, SVector, U1, allocator::Allocator,
 };
 use num_dual::DualNum;
-use quantity::{Energy, MolarEnergy, Moles, Temperature, Volume};
+use quantity::{Dimensionless, MolarEnergy, MolarVolume, Temperature};
 use std::ops::Deref;
 
 mod residual;
@@ -164,17 +165,17 @@ where
     fn ideal_gas_helmholtz_energy<D2: DualNum<f64, Inner = D> + Copy>(
         &self,
         temperature: Temperature<D2>,
-        volume: Volume<D2>,
-        moles: &Moles<OVector<D2, N>>,
-    ) -> Energy<D2> {
+        volume: MolarVolume<D2>,
+        moles: &OVector<D2, N>,
+    ) -> MolarEnergy<D2> {
         let total_moles = moles.sum();
         let molefracs = moles / total_moles;
-        let molar_volume = volume / total_moles;
+        let molar_volume = volume.into_reduced() / total_moles;
         MolarEnergy::from_reduced(self.ideal_gas_molar_helmholtz_energy(
             temperature.into_reduced(),
-            molar_volume.into_reduced(),
+            molar_volume,
             &molefracs,
-        )) * total_moles
+        )) * Dimensionless::new(total_moles)
     }
 }