feos_core/parameter/
association.rs

1use super::{BinaryParameters, BinaryRecord, GroupCount, PureParameters};
2use crate::{FeosResult, parameter::PureRecord};
3use nalgebra::DVector;
4use num_traits::Zero;
5use serde::{Deserialize, Serialize};
6use std::collections::HashMap;
7
8/// Pure component association parameters.
9#[derive(Serialize, Deserialize, Clone, Debug)]
10pub struct AssociationRecord<A> {
11    #[serde(skip_serializing_if = "String::is_empty")]
12    #[serde(default)]
13    pub id: String,
14    #[serde(flatten)]
15    pub parameters: Option<A>,
16    /// \# of association sites of type A
17    #[serde(skip_serializing_if = "f64::is_zero")]
18    #[serde(default)]
19    pub na: f64,
20    /// \# of association sites of type B
21    #[serde(skip_serializing_if = "f64::is_zero")]
22    #[serde(default)]
23    pub nb: f64,
24    /// \# of association sites of type C
25    #[serde(skip_serializing_if = "f64::is_zero")]
26    #[serde(default)]
27    pub nc: f64,
28}
29
30impl<A> AssociationRecord<A> {
31    pub fn new(parameters: Option<A>, na: f64, nb: f64, nc: f64) -> Self {
32        Self::with_id(Default::default(), parameters, na, nb, nc)
33    }
34
35    pub fn with_id(id: String, parameters: Option<A>, na: f64, nb: f64, nc: f64) -> Self {
36        Self {
37            id,
38            parameters,
39            na,
40            nb,
41            nc,
42        }
43    }
44}
45
46/// Binary association parameters.
47#[derive(Serialize, Deserialize, Clone, Debug)]
48pub struct BinaryAssociationRecord<A> {
49    /// Identifier of the association site on the first molecule.
50    #[serde(skip_serializing_if = "String::is_empty")]
51    #[serde(default)]
52    pub id1: String,
53    /// Identifier of the association site on the second molecule.
54    #[serde(skip_serializing_if = "String::is_empty")]
55    #[serde(default)]
56    pub id2: String,
57    /// Binary association parameters
58    #[serde(flatten)]
59    pub parameters: A,
60}
61
62impl<A> BinaryAssociationRecord<A> {
63    pub fn new(parameters: A) -> Self {
64        Self::with_id(Default::default(), Default::default(), parameters)
65    }
66
67    pub fn with_id(id1: String, id2: String, parameters: A) -> Self {
68        Self {
69            id1,
70            id2,
71            parameters,
72        }
73    }
74}
75
76/// The definition of an individual association site.
77#[derive(Clone, Debug)]
78pub struct AssociationSite {
79    /// The index of the component (or group) that the association site is on.
80    pub assoc_comp: usize,
81    /// The identifier of the site (if there are multiple sites on one component/group).
82    pub id: String,
83    /// The multplicity of the association site (NA/NB/NC).
84    pub n: f64,
85}
86
87impl AssociationSite {
88    fn new(assoc_comp: usize, id: String, n: f64) -> Self {
89        Self { assoc_comp, id, n }
90    }
91}
92
93/// The combining rule for association parameters that is used as fallback
94/// if no binary association parameters are available.
95pub trait CombiningRule<P> {
96    fn combining_rule(comp_i: &P, comp_j: &P, parameters_i: &Self, parameters_j: &Self) -> Self;
97}
98
99impl<P> CombiningRule<P> for () {
100    fn combining_rule(_: &P, _: &P, _: &Self, _: &Self) {}
101}
102
103/// Parameter set required for the SAFT association Helmoltz energy
104/// contribution and functional.
105#[derive(Clone)]
106pub struct AssociationParameters<A> {
107    pub component_index: DVector<usize>,
108    pub sites_a: Vec<AssociationSite>,
109    pub sites_b: Vec<AssociationSite>,
110    pub sites_c: Vec<AssociationSite>,
111    pub binary_ab: Vec<BinaryParameters<A, ()>>,
112    pub binary_cc: Vec<BinaryParameters<A, ()>>,
113}
114
115impl<A: Clone> AssociationParameters<A> {
116    pub fn new<P, B>(
117        pure_records: &[PureRecord<P, A>],
118        binary_records: &[BinaryRecord<usize, B, A>],
119    ) -> FeosResult<Self>
120    where
121        A: CombiningRule<P>,
122    {
123        let mut sites_a = Vec::new();
124        let mut sites_b = Vec::new();
125        let mut sites_c = Vec::new();
126        let mut pars_a = Vec::new();
127        let mut pars_b = Vec::new();
128        let mut pars_c = Vec::new();
129
130        for (i, record) in pure_records.iter().enumerate() {
131            for site in record.association_sites.iter() {
132                if site.na > 0.0 {
133                    sites_a.push(AssociationSite::new(i, site.id.clone(), site.na));
134                    pars_a.push(&site.parameters);
135                }
136                if site.nb > 0.0 {
137                    sites_b.push(AssociationSite::new(i, site.id.clone(), site.nb));
138                    pars_b.push(&site.parameters);
139                }
140                if site.nc > 0.0 {
141                    sites_c.push(AssociationSite::new(i, site.id.clone(), site.nc));
142                    pars_c.push(&site.parameters);
143                }
144            }
145        }
146
147        let record_map: HashMap<_, _> = binary_records
148            .iter()
149            .flat_map(|br| {
150                br.association_sites.iter().flat_map(|a| {
151                    [
152                        ((br.id1, br.id2, &a.id1, &a.id2), &a.parameters),
153                        ((br.id2, br.id1, &a.id2, &a.id1), &a.parameters),
154                    ]
155                })
156            })
157            .collect();
158
159        let mut binary_ab = Vec::new();
160        for ((a, site_a), pa) in sites_a.iter().enumerate().zip(&pars_a) {
161            for ((b, site_b), pb) in sites_b.iter().enumerate().zip(&pars_b) {
162                if let Some(&record) =
163                    record_map.get(&(site_a.assoc_comp, site_b.assoc_comp, &site_a.id, &site_b.id))
164                {
165                    binary_ab.push(BinaryParameters::new(a, b, record.clone(), ()));
166                } else if let (Some(pa), Some(pb)) = (pa, pb) {
167                    binary_ab.push(BinaryParameters::new(
168                        a,
169                        b,
170                        A::combining_rule(
171                            &pure_records[site_a.assoc_comp].model_record,
172                            &pure_records[site_b.assoc_comp].model_record,
173                            pa,
174                            pb,
175                        ),
176                        (),
177                    ));
178                }
179            }
180        }
181
182        let mut binary_cc = Vec::new();
183        for ((a, site_a), pa) in sites_c.iter().enumerate().zip(&pars_c) {
184            for ((b, site_b), pb) in sites_c.iter().enumerate().zip(&pars_c) {
185                if let Some(&record) =
186                    record_map.get(&(site_a.assoc_comp, site_b.assoc_comp, &site_a.id, &site_b.id))
187                {
188                    binary_cc.push(BinaryParameters::new(a, b, record.clone(), ()));
189                } else if let (Some(pa), Some(pb)) = (pa, pb) {
190                    binary_cc.push(BinaryParameters::new(
191                        a,
192                        b,
193                        A::combining_rule(
194                            &pure_records[site_a.assoc_comp].model_record,
195                            &pure_records[site_b.assoc_comp].model_record,
196                            pa,
197                            pb,
198                        ),
199                        (),
200                    ));
201                }
202            }
203        }
204
205        let component_index = DVector::from_vec((0..pure_records.len()).collect());
206
207        Ok(Self {
208            component_index,
209            sites_a,
210            sites_b,
211            sites_c,
212            binary_ab,
213            binary_cc,
214        })
215    }
216
217    pub fn new_hetero<P, C: GroupCount>(
218        groups: &[PureParameters<P, C>],
219        association_sites: &[Vec<AssociationRecord<A>>],
220        binary_records: &[BinaryParameters<Vec<BinaryAssociationRecord<A>>, ()>],
221    ) -> FeosResult<Self>
222    where
223        A: CombiningRule<P>,
224    {
225        let mut sites_a = Vec::new();
226        let mut sites_b = Vec::new();
227        let mut sites_c = Vec::new();
228        let mut pars_a = Vec::new();
229        let mut pars_b = Vec::new();
230        let mut pars_c = Vec::new();
231
232        for (i, (record, sites)) in groups.iter().zip(association_sites).enumerate() {
233            for site in sites.iter() {
234                if site.na > 0.0 {
235                    let na = site.na * record.count.into_f64();
236                    sites_a.push(AssociationSite::new(i, site.id.clone(), na));
237                    pars_a.push(&site.parameters)
238                }
239                if site.nb > 0.0 {
240                    let nb = site.nb * record.count.into_f64();
241                    sites_b.push(AssociationSite::new(i, site.id.clone(), nb));
242                    pars_b.push(&site.parameters)
243                }
244                if site.nc > 0.0 {
245                    let nc = site.nc * record.count.into_f64();
246                    sites_c.push(AssociationSite::new(i, site.id.clone(), nc));
247                    pars_c.push(&site.parameters)
248                }
249            }
250        }
251
252        let record_map: HashMap<_, _> = binary_records
253            .iter()
254            .flat_map(|br| {
255                br.model_record.iter().flat_map(|a| {
256                    [
257                        ((br.id1, br.id2, &a.id1, &a.id2), &a.parameters),
258                        ((br.id2, br.id1, &a.id2, &a.id1), &a.parameters),
259                    ]
260                })
261            })
262            .collect();
263
264        let mut binary_ab = Vec::new();
265        for ((a, site_a), pa) in sites_a.iter().enumerate().zip(&pars_a) {
266            for ((b, site_b), pb) in sites_b.iter().enumerate().zip(&pars_b) {
267                if let Some(&record) =
268                    record_map.get(&(site_a.assoc_comp, site_b.assoc_comp, &site_a.id, &site_b.id))
269                {
270                    binary_ab.push(BinaryParameters::new(a, b, record.clone(), ()));
271                } else if let (Some(pa), Some(pb)) = (pa, pb) {
272                    binary_ab.push(BinaryParameters::new(
273                        a,
274                        b,
275                        A::combining_rule(
276                            &groups[site_a.assoc_comp].model_record,
277                            &groups[site_b.assoc_comp].model_record,
278                            pa,
279                            pb,
280                        ),
281                        (),
282                    ));
283                }
284            }
285        }
286
287        let mut binary_cc = Vec::new();
288        for ((a, site_a), pa) in sites_c.iter().enumerate().zip(&pars_c) {
289            for ((b, site_b), pb) in sites_c.iter().enumerate().zip(&pars_c) {
290                if let Some(&record) =
291                    record_map.get(&(site_a.assoc_comp, site_b.assoc_comp, &site_a.id, &site_b.id))
292                {
293                    binary_cc.push(BinaryParameters::new(a, b, record.clone(), ()));
294                } else if let (Some(pa), Some(pb)) = (pa, pb) {
295                    binary_cc.push(BinaryParameters::new(
296                        a,
297                        b,
298                        A::combining_rule(
299                            &groups[site_a.assoc_comp].model_record,
300                            &groups[site_b.assoc_comp].model_record,
301                            pa,
302                            pb,
303                        ),
304                        (),
305                    ));
306                }
307            }
308        }
309
310        let component_index =
311            DVector::from_vec(groups.iter().map(|pr| pr.component_index).collect());
312
313        Ok(Self {
314            component_index,
315            sites_a,
316            sites_b,
317            sites_c,
318            binary_ab,
319            binary_cc,
320        })
321    }
322
323    pub fn is_empty(&self) -> bool {
324        (self.sites_a.is_empty() | self.sites_b.is_empty()) & self.sites_c.is_empty()
325    }
326}
feos_core/parameter/association.rs

feos_core/parameter/
association.rs
