Skip to content

Commit c33df59

Browse files
committed
include phase fraction in PhaseEquilibrium
1 parent d97cb98 commit c33df59

35 files changed

+573
-620
lines changed

.github/workflows/test.yml

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -2,7 +2,7 @@ name: Test
22

33
on:
44
push:
5-
branches: [main, development]
5+
branches: [main]
66
pull_request:
77
branches: [main, development]
88

.github/workflows/wheels.yml

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,7 +1,7 @@
11
name: Build Wheels
22
on:
33
push:
4-
branches: [main, development]
4+
branches: [main]
55
pull_request:
66
branches: [main, development]
77
jobs:

CHANGELOG.md

Lines changed: 8 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -7,8 +7,15 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
77
## [Breaking]
88
### Added
99
- 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)
10-
- Rewrote `PhaseEquilibrium::pure_p` to mirror `pure_t` and enable automatic differentiation. [#337](https://github.com/feos-org/feos/pull/337)
10+
- Rewrote `PhaseEquilibrium::pure_p` to mirror `pure_t` and enabled automatic differentiation. [#337](https://github.com/feos-org/feos/pull/337)
1111
- Added `boiling_temperature` to the list of properties for parallel evaluations of gradients. [#337](https://github.com/feos-org/feos/pull/337)
12+
- 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)
13+
14+
### Changed
15+
- 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)
16+
17+
### Removed
18+
- 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)
1219

1320
### Packaging
1421
- Updated `quantity` dependency to 0.13 and removed the `typenum` dependency. [#323](https://github.com/feos-org/feos/pull/323)

crates/feos-core/src/ad/mod.rs

Lines changed: 13 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
use crate::DensityInitialization::Liquid;
22
use crate::density_iteration::density_iteration;
3-
use crate::{FeosResult, PhaseEquilibrium, ReferenceSystem, Residual};
3+
use crate::{Composition, FeosResult, PhaseEquilibrium, ReferenceSystem, Residual};
44
use nalgebra::{Const, SVector, U1, U2};
55
#[cfg(feature = "rayon")]
66
use ndarray::{Array1, Array2, ArrayView2, Zip};
@@ -215,20 +215,21 @@ pub trait PropertiesAD {
215215
)
216216
}
217217

218-
fn bubble_point_pressure<const P: usize>(
218+
fn bubble_point_pressure<const P: usize, X: Composition<f64, U2>>(
219219
&self,
220220
temperature: Temperature,
221221
pressure: Option<Pressure>,
222-
liquid_molefracs: SVector<f64, 2>,
222+
liquid_molefracs: X,
223223
) -> FeosResult<Pressure<Gradient<P>>>
224224
where
225225
Self: Residual<U2, Gradient<P>>,
226226
{
227227
let eos_f64 = self.re();
228+
let (liquid_molefracs, _) = liquid_molefracs.into_molefracs(&eos_f64)?;
228229
let vle = PhaseEquilibrium::bubble_point(
229230
&eos_f64,
230231
temperature,
231-
&liquid_molefracs,
232+
liquid_molefracs,
232233
pressure,
233234
None,
234235
Default::default(),
@@ -265,20 +266,21 @@ pub trait PropertiesAD {
265266
Ok(Pressure::from_reduced(p))
266267
}
267268

268-
fn dew_point_pressure<const P: usize>(
269+
fn dew_point_pressure<const P: usize, X: Composition<f64, U2>>(
269270
&self,
270271
temperature: Temperature,
271272
pressure: Option<Pressure>,
272-
vapor_molefracs: SVector<f64, 2>,
273+
vapor_molefracs: X,
273274
) -> FeosResult<Pressure<Gradient<P>>>
274275
where
275276
Self: Residual<U2, Gradient<P>>,
276277
{
277278
let eos_f64 = self.re();
279+
let (vapor_molefracs, _) = vapor_molefracs.into_molefracs(&eos_f64)?;
278280
let vle = PhaseEquilibrium::dew_point(
279281
&eos_f64,
280282
temperature,
281-
&vapor_molefracs,
283+
vapor_molefracs,
282284
pressure,
283285
None,
284286
Default::default(),
@@ -329,12 +331,8 @@ pub trait PropertiesAD {
329331
parameters,
330332
input,
331333
|eos: &Self::Lifted<Gradient<P>>, inp| {
332-
eos.bubble_point_pressure(
333-
inp[0] * KELVIN,
334-
Some(inp[2] * PASCAL),
335-
SVector::from([inp[1], 1.0 - inp[1]]),
336-
)
337-
.map(|p| p.convert_into(PASCAL))
334+
eos.bubble_point_pressure(inp[0] * KELVIN, Some(inp[2] * PASCAL), inp[1])
335+
.map(|p| p.convert_into(PASCAL))
338336
},
339337
)
340338
}
@@ -353,12 +351,8 @@ pub trait PropertiesAD {
353351
parameters,
354352
input,
355353
|eos: &Self::Lifted<Gradient<P>>, inp| {
356-
eos.dew_point_pressure(
357-
inp[0] * KELVIN,
358-
Some(inp[2] * PASCAL),
359-
SVector::from([inp[1], 1.0 - inp[1]]),
360-
)
361-
.map(|p| p.convert_into(PASCAL))
354+
eos.dew_point_pressure(inp[0] * KELVIN, Some(inp[2] * PASCAL), inp[1])
355+
.map(|p| p.convert_into(PASCAL))
362356
},
363357
)
364358
}

crates/feos-core/src/errors.rs

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -24,6 +24,10 @@ pub enum FeosError {
2424
InvalidState(String, String, f64),
2525
#[error("Undetermined state: {0}")]
2626
UndeterminedState(String),
27+
#[error(
28+
"Extensive properties can only be evaluated for states that are initialized with extensive properties."
29+
)]
30+
IntensiveState,
2731
#[error("System is supercritical.")]
2832
SuperCritical,
2933
#[error("No phase split according to stability analysis.")]

crates/feos-core/src/lib.rs

Lines changed: 11 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -299,8 +299,8 @@ mod tests {
299299

300300
// residual properties
301301
assert_relative_eq!(
302-
s.helmholtz_energy(Contributions::Residual),
303-
sr.residual_helmholtz_energy(),
302+
s.helmholtz_energy(Contributions::Residual)?,
303+
sr.residual_helmholtz_energy()?,
304304
max_relative = 1e-15
305305
);
306306
assert_relative_eq!(
@@ -309,8 +309,8 @@ mod tests {
309309
max_relative = 1e-15
310310
);
311311
assert_relative_eq!(
312-
s.entropy(Contributions::Residual),
313-
sr.residual_entropy(),
312+
s.entropy(Contributions::Residual)?,
313+
sr.residual_entropy()?,
314314
max_relative = 1e-15
315315
);
316316
assert_relative_eq!(
@@ -319,8 +319,8 @@ mod tests {
319319
max_relative = 1e-15
320320
);
321321
assert_relative_eq!(
322-
s.enthalpy(Contributions::Residual),
323-
sr.residual_enthalpy(),
322+
s.enthalpy(Contributions::Residual)?,
323+
sr.residual_enthalpy()?,
324324
max_relative = 1e-15
325325
);
326326
assert_relative_eq!(
@@ -329,8 +329,8 @@ mod tests {
329329
max_relative = 1e-15
330330
);
331331
assert_relative_eq!(
332-
s.internal_energy(Contributions::Residual),
333-
sr.residual_internal_energy(),
332+
s.internal_energy(Contributions::Residual)?,
333+
sr.residual_internal_energy()?,
334334
max_relative = 1e-15
335335
);
336336
assert_relative_eq!(
@@ -339,12 +339,12 @@ mod tests {
339339
max_relative = 1e-15
340340
);
341341
assert_relative_eq!(
342-
s.gibbs_energy(Contributions::Residual)
343-
- s.total_moles()
342+
s.gibbs_energy(Contributions::Residual)?
343+
- s.total_moles()?
344344
* RGAS
345345
* s.temperature
346346
* s.compressibility(Contributions::Total).ln(),
347-
sr.residual_gibbs_energy(),
347+
sr.residual_gibbs_energy()?,
348348
max_relative = 1e-15
349349
);
350350
assert_relative_eq!(

crates/feos-core/src/phase_equilibria/bubble_dew.rs

Lines changed: 36 additions & 27 deletions
Original file line numberDiff line numberDiff line change
@@ -4,7 +4,7 @@ use crate::state::{
44
Contributions,
55
DensityInitialization::{InitialDensity, Liquid, Vapor},
66
};
7-
use crate::{ReferenceSystem, Residual, SolverOptions, State, Verbosity};
7+
use crate::{Composition, ReferenceSystem, Residual, SolverOptions, State, Verbosity};
88
use nalgebra::allocator::Allocator;
99
use nalgebra::{DMatrix, DVector, DefaultAllocator, Dim, Dyn, OVector, U1};
1010
#[cfg(feature = "ndarray")]
@@ -141,10 +141,10 @@ where
141141
{
142142
/// Calculate a phase equilibrium for a given temperature
143143
/// or pressure and composition of the liquid phase.
144-
pub fn bubble_point<TP: TemperatureOrPressure<D>>(
144+
pub fn bubble_point<TP: TemperatureOrPressure<D>, X: Composition<D, N>>(
145145
eos: &E,
146146
temperature_or_pressure: TP,
147-
liquid_molefracs: &OVector<D, N>,
147+
liquid_molefracs: X,
148148
tp_init: Option<TP::Other>,
149149
vapor_molefracs: Option<&OVector<f64, N>>,
150150
options: (SolverOptions, SolverOptions),
@@ -162,10 +162,10 @@ where
162162

163163
/// Calculate a phase equilibrium for a given temperature
164164
/// or pressure and composition of the vapor phase.
165-
pub fn dew_point<TP: TemperatureOrPressure<D>>(
165+
pub fn dew_point<TP: TemperatureOrPressure<D>, X: Composition<D, N>>(
166166
eos: &E,
167167
temperature_or_pressure: TP,
168-
vapor_molefracs: &OVector<D, N>,
168+
vapor_molefracs: X,
169169
tp_init: Option<TP::Other>,
170170
liquid_molefracs: Option<&OVector<f64, N>>,
171171
options: (SolverOptions, SolverOptions),
@@ -181,35 +181,43 @@ where
181181
)
182182
}
183183

184-
pub(super) fn bubble_dew_point<TP: TemperatureOrPressure<D>>(
184+
pub(super) fn bubble_dew_point<TP: TemperatureOrPressure<D>, X: Composition<D, N>>(
185185
eos: &E,
186186
temperature_or_pressure: TP,
187-
vapor_molefracs: &OVector<D, N>,
187+
vapor_molefracs: X,
188188
tp_init: Option<TP::Other>,
189189
liquid_molefracs: Option<&OVector<f64, N>>,
190190
bubble: bool,
191191
options: (SolverOptions, SolverOptions),
192192
) -> FeosResult<Self> {
193-
let (temperature, pressure, iterate_p) =
194-
temperature_or_pressure.temperature_pressure(tp_init);
195-
Self::bubble_dew_point_tp(
196-
eos,
197-
temperature,
198-
pressure,
199-
vapor_molefracs,
200-
liquid_molefracs,
201-
bubble,
202-
iterate_p,
203-
options,
204-
)
193+
if eos.components() == 1 {
194+
let mut vle = Self::pure(eos, temperature_or_pressure, None, options.1)?;
195+
if bubble {
196+
vle.phase_fractions = [D::from(0.0), D::from(1.0)];
197+
}
198+
Ok(vle)
199+
} else {
200+
let (temperature, pressure, iterate_p) =
201+
temperature_or_pressure.temperature_pressure(tp_init);
202+
Self::bubble_dew_point_tp(
203+
eos,
204+
temperature,
205+
pressure,
206+
vapor_molefracs,
207+
liquid_molefracs,
208+
bubble,
209+
iterate_p,
210+
options,
211+
)
212+
}
205213
}
206214

207215
#[expect(clippy::too_many_arguments)]
208-
fn bubble_dew_point_tp(
216+
fn bubble_dew_point_tp<X: Composition<D, N>>(
209217
eos: &E,
210218
temperature: Option<Temperature<D>>,
211219
pressure: Option<Pressure<D>>,
212-
molefracs_spec: &OVector<D, N>,
220+
composition: X,
213221
molefracs_init: Option<&OVector<f64, N>>,
214222
bubble: bool,
215223
iterate_p: bool,
@@ -218,6 +226,7 @@ where
218226
let eos_re = eos.re();
219227
let mut temperature_re = temperature.map(|t| t.re());
220228
let mut pressure_re = pressure.map(|p| p.re());
229+
let (molefracs_spec, total_moles) = composition.into_molefracs(eos)?;
221230
let molefracs_spec_re = molefracs_spec.map(|x| x.re());
222231
let (v1, rho2) = if iterate_p {
223232
// temperature is specified
@@ -315,7 +324,7 @@ where
315324
Self::newton_step_t(
316325
eos,
317326
t,
318-
molefracs_spec,
327+
&molefracs_spec,
319328
&mut p,
320329
&mut molar_volume,
321330
&mut rho2,
@@ -325,7 +334,7 @@ where
325334
Self::newton_step_p(
326335
eos,
327336
&mut t,
328-
molefracs_spec,
337+
&molefracs_spec,
329338
p,
330339
&mut molar_volume,
331340
&mut rho2,
@@ -348,11 +357,11 @@ where
348357
x2,
349358
)?;
350359

351-
Ok(PhaseEquilibrium(if bubble {
352-
[state2, state1]
360+
Ok(if bubble {
361+
PhaseEquilibrium::with_vapor_phase_fraction(state2, state1, D::from(0.0), total_moles)
353362
} else {
354-
[state1, state2]
355-
}))
363+
PhaseEquilibrium::with_vapor_phase_fraction(state1, state2, D::from(1.0), total_moles)
364+
})
356365
}
357366

358367
fn newton_step_t(

0 commit comments

Comments
 (0)