feos_core/state/
properties.rs

1use super::{Contributions, State};
2use crate::equation_of_state::{Molarweight, Total};
3use crate::{FeosResult, ReferenceSystem, Residual};
4use nalgebra::allocator::Allocator;
5use nalgebra::{DefaultAllocator, OVector};
6use num_dual::{Dual, DualNum, Gradients, partial, partial2};
7use quantity::*;
8use std::ops::{Div, Neg};
9
10type InvP<T> = Quantity<T, <_Pressure as Neg>::Output>;
11type InvT<T> = Quantity<T, <_Temperature as Neg>::Output>;
12
13impl<E: Total<N, D>, N: Gradients, D: DualNum<f64> + Copy> State<E, N, D>
14where
15    DefaultAllocator: Allocator<N>,
16{
17    /// Chemical potential: $\mu_i=\left(\frac{\partial A}{\partial N_i}\right)_{T,V,N_j}$
18    pub fn chemical_potential(&self, contributions: Contributions) -> MolarEnergy<OVector<D, N>> {
19        let residual = || self.residual_chemical_potential();
20        let ideal_gas = || {
21            quantity::ad::gradient_copy(
22                partial2(
23                    |n: Dimensionless<_>, &t, &v| {
24                        self.eos.lift_total().ideal_gas_helmholtz_energy(t, v, &n)
25                    },
26                    &self.temperature,
27                    &self.molar_volume,
28                ),
29                &Dimensionless::new(self.molefracs.clone()),
30            )
31            .1
32        };
33        Self::contributions(ideal_gas, residual, contributions)
34    }
35
36    /// Partial derivative of chemical potential w.r.t. temperature: $\left(\frac{\partial\mu_i}{\partial T}\right)_{V,N_i}$
37    pub fn dmu_dt(&self, contributions: Contributions) -> MolarEntropy<OVector<D, N>> {
38        let residual = || self.dmu_res_dt();
39        let ideal_gas = || {
40            quantity::ad::partial_hessian_copy(
41                partial(
42                    |(n, t): (Dimensionless<_>, _), &v| {
43                        self.eos.lift_total().ideal_gas_helmholtz_energy(t, v, &n)
44                    },
45                    &self.molar_volume,
46                ),
47                (
48                    &Dimensionless::new(self.molefracs.clone()),
49                    self.temperature,
50                ),
51            )
52            .3
53        };
54        Self::contributions(ideal_gas, residual, contributions)
55    }
56
57    /// Molar isochoric heat capacity: $c_v=\left(\frac{\partial u}{\partial T}\right)_{V,N_i}$
58    pub fn molar_isochoric_heat_capacity(&self, contributions: Contributions) -> MolarEntropy<D> {
59        self.temperature * self.ds_dt(contributions)
60    }
61
62    /// Partial derivative of the molar isochoric heat capacity w.r.t. temperature: $\left(\frac{\partial c_V}{\partial T}\right)_{V,N_i}$
63    pub fn dc_v_dt(
64        &self,
65        contributions: Contributions,
66    ) -> <MolarEntropy<D> as Div<Temperature<D>>>::Output {
67        self.temperature * self.d2s_dt2(contributions) + self.ds_dt(contributions)
68    }
69
70    /// Molar isobaric heat capacity: $c_p=\left(\frac{\partial h}{\partial T}\right)_{p,N_i}$
71    pub fn molar_isobaric_heat_capacity(&self, contributions: Contributions) -> MolarEntropy<D> {
72        match contributions {
73            Contributions::Residual => self.residual_molar_isobaric_heat_capacity(),
74            _ => {
75                self.temperature
76                    * (self.ds_dt(contributions)
77                        - (self.dp_dt(contributions) * self.dp_dt(contributions))
78                            / self.dp_dv(contributions))
79            }
80        }
81    }
82
83    /// Entropy: $S=-\left(\frac{\partial A}{\partial T}\right)_{V,N_i}$
84    pub fn entropy(&self, contributions: Contributions) -> FeosResult<Entropy<D>> {
85        Ok(self.molar_entropy(contributions) * self.total_moles()?)
86    }
87
88    /// Molar entropy: $s=\frac{S}{N}$
89    pub fn molar_entropy(&self, contributions: Contributions) -> MolarEntropy<D> {
90        let residual = || self.residual_molar_entropy();
91        let ideal_gas = || {
92            -quantity::ad::first_derivative(
93                partial2(
94                    |t, &v, n| self.eos.lift_total().ideal_gas_helmholtz_energy(t, v, n),
95                    &self.molar_volume,
96                    &self.molefracs,
97                ),
98                self.temperature,
99            )
100            .1
101        };
102        Self::contributions(ideal_gas, residual, contributions)
103    }
104
105    /// Partial molar entropy: $s_i=\left(\frac{\partial S}{\partial N_i}\right)_{T,p,N_j}$
106    pub fn partial_molar_entropy(&self) -> MolarEntropy<OVector<D, N>> {
107        let c = Contributions::Total;
108        -(self.dmu_dt(c) + self.n_dp_dni(c) * (self.dp_dt(c) / self.dp_dv(c)))
109    }
110
111    /// Partial derivative of the molar entropy w.r.t. temperature: $\left(\frac{\partial s}{\partial T}\right)_{V,N_i}$
112    pub fn ds_dt(
113        &self,
114        contributions: Contributions,
115    ) -> <MolarEntropy<D> as Div<Temperature<D>>>::Output {
116        let residual = || self.ds_res_dt();
117        let ideal_gas = || {
118            -quantity::ad::second_derivative(
119                partial2(
120                    |t, &v, n| self.eos.lift_total().ideal_gas_helmholtz_energy(t, v, n),
121                    &self.molar_volume,
122                    &self.molefracs,
123                ),
124                self.temperature,
125            )
126            .2
127        };
128        Self::contributions(ideal_gas, residual, contributions)
129    }
130
131    /// Second partial derivative of the molar entropy w.r.t. temperature: $\left(\frac{\partial^2 s}{\partial T^2}\right)_{V,N_i}$
132    pub fn d2s_dt2(
133        &self,
134        contributions: Contributions,
135    ) -> <<MolarEntropy<D> as Div<Temperature<D>>>::Output as Div<Temperature<D>>>::Output {
136        let residual = || self.d2s_res_dt2();
137        let ideal_gas = || {
138            -quantity::ad::third_derivative(
139                partial2(
140                    |t, &v, n| self.eos.lift_total().ideal_gas_helmholtz_energy(t, v, n),
141                    &self.molar_volume,
142                    &self.molefracs,
143                ),
144                self.temperature,
145            )
146            .3
147        };
148        Self::contributions(ideal_gas, residual, contributions)
149    }
150
151    /// Enthalpy: $H=A+TS+pV$
152    pub fn enthalpy(&self, contributions: Contributions) -> FeosResult<Energy<D>> {
153        Ok(self.molar_enthalpy(contributions) * self.total_moles()?)
154    }
155
156    /// Molar enthalpy: $h=\frac{H}{N}$
157    pub fn molar_enthalpy(&self, contributions: Contributions) -> MolarEnergy<D> {
158        self.temperature * self.molar_entropy(contributions)
159            + self.molar_helmholtz_energy(contributions)
160            + self.pressure(contributions) * self.molar_volume
161    }
162
163    /// Partial molar enthalpy: $h_i=\left(\frac{\partial H}{\partial N_i}\right)_{T,p,N_j}$
164    pub fn partial_molar_enthalpy(&self) -> MolarEnergy<OVector<D, N>> {
165        let s = self.partial_molar_entropy();
166        let mu = self.chemical_potential(Contributions::Total);
167        s * self.temperature + mu
168    }
169
170    /// Helmholtz energy: $A$
171    pub fn helmholtz_energy(&self, contributions: Contributions) -> FeosResult<Energy<D>> {
172        Ok(self.molar_helmholtz_energy(contributions) * self.total_moles()?)
173    }
174
175    /// Molar Helmholtz energy: $a=\frac{A}{N}$
176    pub fn molar_helmholtz_energy(&self, contributions: Contributions) -> MolarEnergy<D> {
177        let residual = || self.residual_molar_helmholtz_energy();
178        let ideal_gas = || {
179            quantity::ad::zeroth_derivative(
180                partial2(
181                    |t, &v, n| self.eos.lift_total().ideal_gas_helmholtz_energy(t, v, n),
182                    &self.molar_volume,
183                    &self.molefracs,
184                ),
185                self.temperature,
186            )
187        };
188        Self::contributions(ideal_gas, residual, contributions)
189    }
190
191    /// Internal energy: $U=A+TS$
192    pub fn internal_energy(&self, contributions: Contributions) -> FeosResult<Energy<D>> {
193        Ok(self.molar_internal_energy(contributions) * self.total_moles()?)
194    }
195
196    /// Molar internal energy: $u=\frac{U}{N}$
197    pub fn molar_internal_energy(&self, contributions: Contributions) -> MolarEnergy<D> {
198        self.temperature * self.molar_entropy(contributions)
199            + self.molar_helmholtz_energy(contributions)
200    }
201
202    /// Gibbs energy: $G=A+pV$
203    pub fn gibbs_energy(&self, contributions: Contributions) -> FeosResult<Energy<D>> {
204        Ok(self.molar_gibbs_energy(contributions) * self.total_moles()?)
205    }
206
207    /// Molar Gibbs energy: $g=\frac{G}{N}$
208    pub fn molar_gibbs_energy(&self, contributions: Contributions) -> MolarEnergy<D> {
209        self.pressure(contributions) * self.molar_volume
210            + self.molar_helmholtz_energy(contributions)
211    }
212
213    /// Joule Thomson coefficient: $\mu_{JT}=\left(\frac{\partial T}{\partial p}\right)_{H,N_i}$
214    pub fn joule_thomson(&self) -> <Temperature<D> as Div<Pressure<D>>>::Output {
215        let c = Contributions::Total;
216        -(self.molar_volume + self.temperature * self.dp_dt(c) / self.dp_dv(c))
217            / self.molar_isobaric_heat_capacity(c)
218    }
219
220    /// Isentropic compressibility: $\kappa_s=-\frac{1}{V}\left(\frac{\partial V}{\partial p}\right)_{S,N_i}$
221    pub fn isentropic_compressibility(&self) -> InvP<D> {
222        let c = Contributions::Total;
223        -self.molar_isochoric_heat_capacity(c)
224            / (self.molar_isobaric_heat_capacity(c) * self.dp_dv(c) * self.molar_volume)
225    }
226
227    /// Isenthalpic compressibility: $\kappa_H=-\frac{1}{V}\left(\frac{\partial V}{\partial p}\right)_{H,N_i}$
228    pub fn isenthalpic_compressibility(&self) -> InvP<D> {
229        self.isentropic_compressibility() * Dimensionless::new(self.grueneisen_parameter() + 1.0)
230    }
231
232    /// Thermal expansivity: $\alpha_p=-\frac{1}{V}\left(\frac{\partial V}{\partial T}\right)_{p,N_i}$
233    pub fn thermal_expansivity(&self) -> InvT<D> {
234        let c = Contributions::Total;
235        -self.dp_dt(c) / (self.dp_dv(c) * self.molar_volume)
236    }
237
238    /// Grueneisen parameter: $\phi=V\left(\frac{\partial p}{\partial U}\right)_{V,n_i}=\frac{v}{c_v}\left(\frac{\partial p}{\partial T}\right)_{v,n_i}=\frac{\rho}{T}\left(\frac{\partial T}{\partial \rho}\right)_{s, n_i}$
239    pub fn grueneisen_parameter(&self) -> D {
240        let c = Contributions::Total;
241        (self.dp_dt(c) / (self.molar_isochoric_heat_capacity(c) * self.density)).into_value()
242    }
243
244    /// Chemical potential $\mu_i$ evaluated for each contribution of the equation of state.
245    pub fn chemical_potential_contributions(
246        &self,
247        component: usize,
248        contributions: Contributions,
249    ) -> Vec<(&'static str, MolarEnergy<D>)> {
250        let t = Dual::from_re(self.temperature.into_reduced());
251        let v = Dual::from_re(self.density.into_reduced().recip());
252        let mut x = self.molefracs.map(Dual::from_re);
253        x[component].eps = D::one();
254        let mut res = Vec::new();
255        if let Contributions::IdealGas | Contributions::Total = contributions {
256            res.push((
257                self.eos.ideal_gas_model(),
258                self.eos
259                    .lift_total()
260                    .ideal_gas_molar_helmholtz_energy(t, v, &x),
261            ));
262        }
263        if let Contributions::Residual | Contributions::Total = contributions {
264            res.extend(self.eos.lift().helmholtz_energy_contributions(t, v, &x));
265        }
266        res.into_iter()
267            .map(|(s, v)| (s, MolarEnergy::from_reduced(v.eps)))
268            .collect()
269    }
270}
271
272impl<E: Total<N, D> + Molarweight<N, D>, N: Gradients, D: DualNum<f64> + Copy> State<E, N, D>
273where
274    DefaultAllocator: Allocator<N>,
275{
276    /// Specific isochoric heat capacity: $c_v^{(m)}=\frac{C_v}{m}$
277    pub fn specific_isochoric_heat_capacity(
278        &self,
279        contributions: Contributions,
280    ) -> SpecificEntropy<D> {
281        self.molar_isochoric_heat_capacity(contributions) / self.total_molar_weight()
282    }
283
284    /// Specific isobaric heat capacity: $c_p^{(m)}=\frac{C_p}{m}$
285    pub fn specific_isobaric_heat_capacity(
286        &self,
287        contributions: Contributions,
288    ) -> SpecificEntropy<D> {
289        self.molar_isobaric_heat_capacity(contributions) / self.total_molar_weight()
290    }
291
292    /// Specific entropy: $s^{(m)}=\frac{S}{m}$
293    pub fn specific_entropy(&self, contributions: Contributions) -> SpecificEntropy<D> {
294        self.molar_entropy(contributions) / self.total_molar_weight()
295    }
296
297    /// Specific enthalpy: $h^{(m)}=\frac{H}{m}$
298    pub fn specific_enthalpy(&self, contributions: Contributions) -> SpecificEnergy<D> {
299        self.molar_enthalpy(contributions) / self.total_molar_weight()
300    }
301
302    /// Specific Helmholtz energy: $a^{(m)}=\frac{A}{m}$
303    pub fn specific_helmholtz_energy(&self, contributions: Contributions) -> SpecificEnergy<D> {
304        self.molar_helmholtz_energy(contributions) / self.total_molar_weight()
305    }
306
307    /// Specific internal energy: $u^{(m)}=\frac{U}{m}$
308    pub fn specific_internal_energy(&self, contributions: Contributions) -> SpecificEnergy<D> {
309        self.molar_internal_energy(contributions) / self.total_molar_weight()
310    }
311
312    /// Specific Gibbs energy: $g^{(m)}=\frac{G}{m}$
313    pub fn specific_gibbs_energy(&self, contributions: Contributions) -> SpecificEnergy<D> {
314        self.molar_gibbs_energy(contributions) / self.total_molar_weight()
315    }
316}
317
318impl<E: Total<N> + Molarweight<N>, N: Gradients> State<E, N>
319where
320    DefaultAllocator: Allocator<N>,
321{
322    /// Speed of sound: $c=\sqrt{\left(\frac{\partial p}{\partial\rho^{(m)}}\right)_{S,N_i}}$
323    pub fn speed_of_sound(&self) -> Velocity {
324        (self.density * self.total_molar_weight() * self.isentropic_compressibility())
325            .inv()
326            .sqrt()
327    }
328}
feos_core/state/properties.rs

feos_core/state/
properties.rs
