feos_core/ad/properties/
dew_point_pressure.rs1use super::PropertyAD;
2use crate::ad::Gradient;
3use crate::{Composition, Contributions, FeosResult, PhaseEquilibrium, ReferenceSystem, Residual};
4use nalgebra::allocator::Allocator;
5use nalgebra::{DefaultAllocator, U1};
6use num_dual::{DualNum, DualStruct, Gradients};
7use quantity::{_Pressure, KELVIN, PASCAL, Pressure, Temperature};
8
9pub struct DewPointPressure(pub Temperature, pub f64, pub Option<Pressure>);
15
16impl<'a> From<&'a [f64]> for DewPointPressure {
17 fn from(value: &'a [f64]) -> Self {
18 Self(value[0] * KELVIN, value[1], Some(value[2] * PASCAL))
19 }
20}
21
22impl<N: Gradients> PropertyAD<N> for DewPointPressure
23where
24 DefaultAllocator: Allocator<N> + Allocator<U1, N> + Allocator<N, N>,
25 f64: Composition<f64, N>,
26{
27 type Unit = _Pressure;
28 const REFERENCE: Pressure = PASCAL;
29
30 fn evaluate<E: Residual<N, D>, D: DualNum<f64, Inner = f64> + Copy>(
31 &self,
32 eos: &E,
33 ) -> FeosResult<Pressure<D>>
34 where
35 DefaultAllocator: Allocator<N> + Allocator<U1, N> + Allocator<N, N>,
36 {
37 let t = Temperature::from_inner(&self.0);
38 let p = Option::from_inner(&self.2);
39 let (y, _) = self.1.into_molefracs(&eos.re())?;
40 let y = y.map(D::from);
41 let vle = PhaseEquilibrium::dew_point(eos, t, y, p, None, Default::default())?;
42 Ok(vle.vapor().pressure(Contributions::Total))
43 }
44
45 fn evaluate_gradient<E: Residual<N, Gradient<P>>, const P: usize>(
46 &self,
47 eos: &E,
48 ) -> FeosResult<quantity::Quantity<Gradient<P>, Self::Unit>>
49 where
50 DefaultAllocator: Allocator<N> + Allocator<U1, N> + Allocator<N, N>,
51 {
52 let eos_f64 = eos.re();
53 let (vapor_molefracs, _) = self.1.into_molefracs(&eos_f64)?;
54 let vle = PhaseEquilibrium::dew_point(
55 &eos_f64,
56 self.0,
57 &vapor_molefracs,
58 self.2,
59 None,
60 Default::default(),
61 )?;
62
63 let v_l = 1.0 / vle.liquid().density.to_reduced();
64 let v_v = 1.0 / vle.vapor().density.to_reduced();
65 let x = &vle.liquid().molefracs;
66 let t = self.0.into_reduced();
67 let (a_l, a_v, v_l, v_v) = {
68 let t = Gradient::from(t);
69 let v_l = Gradient::from(v_l);
70 let v_v = Gradient::from(v_v);
71 let x = x.map(Gradient::from);
72 let y = vapor_molefracs.map(Gradient::from);
73
74 let a_l = eos.residual_helmholtz_energy(t, v_l, &x);
75 let (p_v, mu_res_v, dp_v, dmu_v) = eos.dmu_dv(t, v_v, &y);
76 let vi_v = dmu_v / dp_v;
77 let v_v = vi_v.dot(&x);
78 let a_v = (mu_res_v - vi_v * p_v).dot(&x);
79 (a_l, a_v, v_l, v_v)
80 };
81 let rho_l = vle.liquid().partial_density().to_reduced();
82 let rho_l = [rho_l[0], rho_l[1]];
83 let rho_v = vle.vapor().partial_density().to_reduced();
84 let rho_v = [rho_v[0], rho_v[1]];
85 let p = -(a_l - a_v
86 + t * (x[0] * (rho_l[0] / rho_v[0]).ln() + x[1] * (rho_l[1] / rho_v[1]).ln() - 1.0))
87 / (v_l - v_v);
88 Ok(Pressure::from_reduced(p))
89 }
90}