Skip to content

Commit b084cd9

Browse files
prehnerprehner
authored
Move CombiningRule to core and provide info on eos parameters (#304)
1 parent 6b2a106 commit b084cd9

File tree

23 files changed

+678
-515
lines changed

23 files changed

+678
-515
lines changed

‎crates/feos-core/src/cubic.rs‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -73,7 +73,7 @@ impl PengRobinsonParameters {
7373
/// A simple version of the Peng-Robinson equation of state.
7474
pub struct PengRobinson {
7575
/// Parameters
76-
parameters: PengRobinsonParameters,
76+
pub parameters: PengRobinsonParameters,
7777
/// Critical temperature in Kelvin
7878
tc: DVector<f64>,
7979
a: DVector<f64>,

‎crates/feos-core/src/parameter/association.rs‎

Lines changed: 145 additions & 137 deletions
Original file line numberDiff line numberDiff line change
@@ -1,9 +1,9 @@
11
use super::{BinaryParameters, BinaryRecord, GroupCount, PureParameters};
2-
use crate::{FeosError, FeosResult, parameter::PureRecord};
2+
use crate::{FeosResult, parameter::PureRecord};
33
use nalgebra::DVector;
44
use num_traits::Zero;
55
use serde::{Deserialize, Serialize};
6-
use std::collections::{HashMap, HashSet};
6+
use std::collections::HashMap;
77

88
/// Pure component association parameters.
99
#[derive(Serialize, Deserialize, Clone, Debug)]
@@ -74,32 +74,34 @@ impl<A> BinaryAssociationRecord<A> {
7474
}
7575

7676
#[derive(Clone, Debug)]
77-
pub struct AssociationSite<A> {
77+
pub struct AssociationSite {
7878
pub assoc_comp: usize,
7979
pub id: String,
8080
pub n: f64,
81-
pub parameters: A,
8281
}
8382

84-
impl<A> AssociationSite<A> {
85-
fn new(assoc_comp: usize, id: String, n: f64, parameters: A) -> Self {
86-
Self {
87-
assoc_comp,
88-
id,
89-
n,
90-
parameters,
91-
}
83+
impl AssociationSite {
84+
fn new(assoc_comp: usize, id: String, n: f64) -> Self {
85+
Self { assoc_comp, id, n }
9286
}
9387
}
9488

89+
pub trait CombiningRule<P> {
90+
fn combining_rule(comp_i: &P, comp_j: &P, parameters_i: &Self, parameters_j: &Self) -> Self;
91+
}
92+
93+
impl<P> CombiningRule<P> for () {
94+
fn combining_rule(_: &P, _: &P, _: &Self, _: &Self) {}
95+
}
96+
9597
/// Parameter set required for the SAFT association Helmoltz energy
9698
/// contribution and functional.
9799
#[derive(Clone)]
98100
pub struct AssociationParameters<A> {
99101
pub component_index: DVector<usize>,
100-
pub sites_a: Vec<AssociationSite<Option<A>>>,
101-
pub sites_b: Vec<AssociationSite<Option<A>>>,
102-
pub sites_c: Vec<AssociationSite<Option<A>>>,
102+
pub sites_a: Vec<AssociationSite>,
103+
pub sites_b: Vec<AssociationSite>,
104+
pub sites_c: Vec<AssociationSite>,
103105
pub binary_ab: Vec<BinaryParameters<A, ()>>,
104106
pub binary_cc: Vec<BinaryParameters<A, ()>>,
105107
}
@@ -108,96 +110,92 @@ impl<A: Clone> AssociationParameters<A> {
108110
pub fn new<P, B>(
109111
pure_records: &[PureRecord<P, A>],
110112
binary_records: &[BinaryRecord<usize, B, A>],
111-
) -> FeosResult<Self> {
113+
) -> FeosResult<Self>
114+
where
115+
A: CombiningRule<P>,
116+
{
112117
let mut sites_a = Vec::new();
113118
let mut sites_b = Vec::new();
114119
let mut sites_c = Vec::new();
120+
let mut pars_a = Vec::new();
121+
let mut pars_b = Vec::new();
122+
let mut pars_c = Vec::new();
115123

116124
for (i, record) in pure_records.iter().enumerate() {
117125
for site in record.association_sites.iter() {
118-
let par = &site.parameters;
119126
if site.na > 0.0 {
120-
sites_a.push(AssociationSite::new(
121-
i,
122-
site.id.clone(),
123-
site.na,
124-
par.clone(),
125-
));
127+
sites_a.push(AssociationSite::new(i, site.id.clone(), site.na));
128+
pars_a.push(&site.parameters);
126129
}
127130
if site.nb > 0.0 {
128-
sites_b.push(AssociationSite::new(
129-
i,
130-
site.id.clone(),
131-
site.nb,
132-
par.clone(),
133-
));
131+
sites_b.push(AssociationSite::new(i, site.id.clone(), site.nb));
132+
pars_b.push(&site.parameters);
134133
}
135134
if site.nc > 0.0 {
136-
sites_c.push(AssociationSite::new(
137-
i,
138-
site.id.clone(),
139-
site.nc,
140-
par.clone(),
141-
));
135+
sites_c.push(AssociationSite::new(i, site.id.clone(), site.nc));
136+
pars_c.push(&site.parameters);
142137
}
143138
}
144139
}
145140

146-
let indices_a: HashMap<_, _> = sites_a
141+
let record_map: HashMap<_, _> = binary_records
147142
.iter()
148-
.enumerate()
149-
.map(|(i, site)| ((site.assoc_comp, &site.id), i))
150-
.collect();
151-
152-
let indices_b: HashMap<_, _> = sites_b
153-
.iter()
154-
.enumerate()
155-
.map(|(i, site)| ((site.assoc_comp, &site.id), i))
156-
.collect();
157-
158-
let indices_c: HashMap<_, _> = sites_c
159-
.iter()
160-
.enumerate()
161-
.map(|(i, site)| ((site.assoc_comp, &site.id), i))
162-
.collect();
163-
164-
let index_set: HashSet<_> = indices_a
165-
.keys()
166-
.chain(indices_b.keys())
167-
.chain(indices_c.keys())
168-
.copied()
143+
.flat_map(|br| {
144+
br.association_sites.iter().flat_map(|a| {
145+
[
146+
((br.id1, br.id2, &a.id1, &a.id2), &a.parameters),
147+
((br.id2, br.id1, &a.id2, &a.id1), &a.parameters),
148+
]
149+
})
150+
})
169151
.collect();
170152

171153
let mut binary_ab = Vec::new();
172-
let mut binary_cc = Vec::new();
173-
for br in binary_records {
174-
let i = br.id1;
175-
let j = br.id2;
176-
for record in &br.association_sites {
177-
let a = &record.id1;
178-
let b = &record.id2;
179-
if !index_set.contains(&(i, a)) {
180-
return Err(FeosError::IncompatibleParameters(format!(
181-
"No association site {a} on component {i}"
182-
)));
183-
}
184-
if !index_set.contains(&(j, b)) {
185-
return Err(FeosError::IncompatibleParameters(format!(
186-
"No association site {b} on component {j}"
187-
)));
188-
}
189-
if let (Some(x), Some(y)) = (indices_a.get(&(i, a)), indices_b.get(&(j, b))) {
190-
binary_ab.push(BinaryParameters::new(*x, *y, record.parameters.clone(), ()));
191-
}
192-
if let (Some(y), Some(x)) = (indices_b.get(&(i, a)), indices_a.get(&(j, b))) {
193-
binary_ab.push(BinaryParameters::new(*x, *y, record.parameters.clone(), ()));
154+
for ((a, site_a), pa) in sites_a.iter().enumerate().zip(&pars_a) {
155+
for ((b, site_b), pb) in sites_b.iter().enumerate().zip(&pars_b) {
156+
if let Some(&record) =
157+
record_map.get(&(site_a.assoc_comp, site_b.assoc_comp, &site_a.id, &site_b.id))
158+
{
159+
binary_ab.push(BinaryParameters::new(a, b, record.clone(), ()));
160+
} else if let (Some(pa), Some(pb)) = (pa, pb) {
161+
binary_ab.push(BinaryParameters::new(
162+
a,
163+
b,
164+
A::combining_rule(
165+
&pure_records[site_a.assoc_comp].model_record,
166+
&pure_records[site_b.assoc_comp].model_record,
167+
pa,
168+
pb,
169+
),
170+
(),
171+
));
194172
}
195-
if let (Some(x), Some(y)) = (indices_c.get(&(i, a)), indices_c.get(&(j, b))) {
196-
binary_cc.push(BinaryParameters::new(*x, *y, record.parameters.clone(), ()));
197-
binary_cc.push(BinaryParameters::new(*y, *x, record.parameters.clone(), ()));
173+
}
174+
}
175+
176+
let mut binary_cc = Vec::new();
177+
for ((a, site_a), pa) in sites_c.iter().enumerate().zip(&pars_c) {
178+
for ((b, site_b), pb) in sites_c.iter().enumerate().zip(&pars_c) {
179+
if let Some(&record) =
180+
record_map.get(&(site_a.assoc_comp, site_b.assoc_comp, &site_a.id, &site_b.id))
181+
{
182+
binary_cc.push(BinaryParameters::new(a, b, record.clone(), ()));
183+
} else if let (Some(pa), Some(pb)) = (pa, pb) {
184+
binary_cc.push(BinaryParameters::new(
185+
a,
186+
b,
187+
A::combining_rule(
188+
&pure_records[site_a.assoc_comp].model_record,
189+
&pure_records[site_b.assoc_comp].model_record,
190+
pa,
191+
pb,
192+
),
193+
(),
194+
));
198195
}
199196
}
200197
}
198+
201199
let component_index = DVector::from_vec((0..pure_records.len()).collect());
202200

203201
Ok(Self {
@@ -211,90 +209,100 @@ impl<A: Clone> AssociationParameters<A> {
211209
}
212210

213211
pub fn new_hetero<P, C: GroupCount>(
214-
pure_records: &[PureParameters<P, C>],
212+
groups: &[PureParameters<P, C>],
215213
association_sites: &[Vec<AssociationRecord<A>>],
216214
binary_records: &[BinaryParameters<Vec<BinaryAssociationRecord<A>>, ()>],
217-
) -> FeosResult<Self> {
215+
) -> FeosResult<Self>
216+
where
217+
A: CombiningRule<P>,
218+
{
218219
let mut sites_a = Vec::new();
219220
let mut sites_b = Vec::new();
220221
let mut sites_c = Vec::new();
222+
let mut pars_a = Vec::new();
223+
let mut pars_b = Vec::new();
224+
let mut pars_c = Vec::new();
221225

222-
for (i, (record, sites)) in pure_records.iter().zip(association_sites).enumerate() {
226+
for (i, (record, sites)) in groups.iter().zip(association_sites).enumerate() {
223227
for site in sites.iter() {
224-
let par = &site.parameters;
225228
if site.na > 0.0 {
226229
let na = site.na * record.count.into_f64();
227-
sites_a.push(AssociationSite::new(i, site.id.clone(), na, par.clone()));
230+
sites_a.push(AssociationSite::new(i, site.id.clone(), na));
231+
pars_a.push(&site.parameters)
228232
}
229233
if site.nb > 0.0 {
230234
let nb = site.nb * record.count.into_f64();
231-
sites_b.push(AssociationSite::new(i, site.id.clone(), nb, par.clone()));
235+
sites_b.push(AssociationSite::new(i, site.id.clone(), nb));
236+
pars_b.push(&site.parameters)
232237
}
233238
if site.nc > 0.0 {
234239
let nc = site.nc * record.count.into_f64();
235-
sites_c.push(AssociationSite::new(i, site.id.clone(), nc, par.clone()));
240+
sites_c.push(AssociationSite::new(i, site.id.clone(), nc));
241+
pars_c.push(&site.parameters)
236242
}
237243
}
238244
}
239245

240-
let indices_a: HashMap<_, _> = sites_a
241-
.iter()
242-
.enumerate()
243-
.map(|(i, site)| ((site.assoc_comp, &site.id), i))
244-
.collect();
245-
246-
let indices_b: HashMap<_, _> = sites_b
247-
.iter()
248-
.enumerate()
249-
.map(|(i, site)| ((site.assoc_comp, &site.id), i))
250-
.collect();
251-
252-
let indices_c: HashMap<_, _> = sites_c
246+
let record_map: HashMap<_, _> = binary_records
253247
.iter()
254-
.enumerate()
255-
.map(|(i, site)| ((site.assoc_comp, &site.id), i))
256-
.collect();
257-
258-
let index_set: HashSet<_> = indices_a
259-
.keys()
260-
.chain(indices_b.keys())
261-
.chain(indices_c.keys())
262-
.copied()
248+
.flat_map(|br| {
249+
br.model_record.iter().flat_map(|a| {
250+
[
251+
((br.id1, br.id2, &a.id1, &a.id2), &a.parameters),
252+
((br.id2, br.id1, &a.id2, &a.id1), &a.parameters),
253+
]
254+
})
255+
})
263256
.collect();
264257

265258
let mut binary_ab = Vec::new();
266-
let mut binary_cc = Vec::new();
267-
for br in binary_records {
268-
let i = br.id1;
269-
let j = br.id2;
270-
for record in &br.model_record {
271-
let a = &record.id1;
272-
let b = &record.id2;
273-
if !index_set.contains(&(i, a)) {
274-
return Err(FeosError::IncompatibleParameters(format!(
275-
"No association site {a} on component {i}"
276-
)));
277-
}
278-
if !index_set.contains(&(j, b)) {
279-
return Err(FeosError::IncompatibleParameters(format!(
280-
"No association site {b} on component {j}"
281-
)));
282-
}
283-
if let (Some(x), Some(y)) = (indices_a.get(&(i, a)), indices_b.get(&(j, b))) {
284-
binary_ab.push(BinaryParameters::new(*x, *y, record.parameters.clone(), ()));
285-
}
286-
if let (Some(y), Some(x)) = (indices_b.get(&(i, a)), indices_a.get(&(j, b))) {
287-
binary_ab.push(BinaryParameters::new(*x, *y, record.parameters.clone(), ()));
259+
for ((a, site_a), pa) in sites_a.iter().enumerate().zip(&pars_a) {
260+
for ((b, site_b), pb) in sites_b.iter().enumerate().zip(&pars_b) {
261+
if let Some(&record) =
262+
record_map.get(&(site_a.assoc_comp, site_b.assoc_comp, &site_a.id, &site_b.id))
263+
{
264+
binary_ab.push(BinaryParameters::new(a, b, record.clone(), ()));
265+
} else if let (Some(pa), Some(pb)) = (pa, pb) {
266+
binary_ab.push(BinaryParameters::new(
267+
a,
268+
b,
269+
A::combining_rule(
270+
&groups[site_a.assoc_comp].model_record,
271+
&groups[site_b.assoc_comp].model_record,
272+
pa,
273+
pb,
274+
),
275+
(),
276+
));
288277
}
289-
if let (Some(x), Some(y)) = (indices_c.get(&(i, a)), indices_c.get(&(j, b))) {
290-
binary_cc.push(BinaryParameters::new(*x, *y, record.parameters.clone(), ()));
291-
binary_cc.push(BinaryParameters::new(*y, *x, record.parameters.clone(), ()));
278+
}
279+
}
280+
281+
let mut binary_cc = Vec::new();
282+
for ((a, site_a), pa) in sites_c.iter().enumerate().zip(&pars_c) {
283+
for ((b, site_b), pb) in sites_c.iter().enumerate().zip(&pars_c) {
284+
if let Some(&record) =
285+
record_map.get(&(site_a.assoc_comp, site_b.assoc_comp, &site_a.id, &site_b.id))
286+
{
287+
binary_cc.push(BinaryParameters::new(a, b, record.clone(), ()));
288+
} else if let (Some(pa), Some(pb)) = (pa, pb) {
289+
binary_cc.push(BinaryParameters::new(
290+
a,
291+
b,
292+
A::combining_rule(
293+
&groups[site_a.assoc_comp].model_record,
294+
&groups[site_b.assoc_comp].model_record,
295+
pa,
296+
pb,
297+
),
298+
(),
299+
));
292300
}
293301
}
294302
}
295303

296304
let component_index =
297-
DVector::from_vec(pure_records.iter().map(|pr| pr.component_index).collect());
305+
DVector::from_vec(groups.iter().map(|pr| pr.component_index).collect());
298306

299307
Ok(Self {
300308
component_index,

0 commit comments

Comments
 (0)