feos_pets/
parameters.rs

1use feos_core::joback::JobackRecord;
2use feos_core::parameter::{Parameter, PureRecord};
3use ndarray::{Array, Array1, Array2};
4use num_traits::Zero;
5use serde::{Deserialize, Serialize};
6use std::collections::HashMap;
7use std::fmt::Write;
8
9/// PeTS parameter set.
10#[derive(Serialize, Deserialize, Debug, Clone, Default)]
11pub struct PetsRecord {
12    /// Segment diameter in units of Angstrom
13    pub sigma: f64,
14    /// Energetic parameter in units of Kelvin
15    pub epsilon_k: f64,
16    /// Entropy scaling parameters for viscosity
17    #[serde(skip_serializing_if = "Option::is_none")]
18    pub viscosity: Option<[f64; 4]>,
19    /// Entropy scaling parameters for self-diffusion coefficient
20    #[serde(skip_serializing_if = "Option::is_none")]
21    pub diffusion: Option<[f64; 5]>,
22    /// Entropy scaling parameters for thermal conductivity
23    #[serde(skip_serializing_if = "Option::is_none")]
24    pub thermal_conductivity: Option<[f64; 4]>,
25}
26
27impl std::fmt::Display for PetsRecord {
28    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
29        write!(f, "PetsRecord(sigma={}", self.sigma)?;
30        write!(f, ", epsilon_k={}", self.epsilon_k)?;
31        if let Some(n) = &self.viscosity {
32            write!(f, ", viscosity={:?}", n)?;
33        }
34        if let Some(n) = &self.diffusion {
35            write!(f, ", diffusion={:?}", n)?;
36        }
37        if let Some(n) = &self.thermal_conductivity {
38            write!(f, ", thermal_conductivity={:?}", n)?;
39        }
40        write!(f, ")")
41    }
42}
43
44impl PetsRecord {
45    pub fn new(
46        sigma: f64,
47        epsilon_k: f64,
48        viscosity: Option<[f64; 4]>,
49        diffusion: Option<[f64; 5]>,
50        thermal_conductivity: Option<[f64; 4]>,
51    ) -> PetsRecord {
52        PetsRecord {
53            sigma,
54            epsilon_k,
55            viscosity,
56            diffusion,
57            thermal_conductivity,
58        }
59    }
60}
61
62#[derive(Serialize, Deserialize, Clone, Default, Debug)]
63pub struct PetsBinaryRecord {
64    k_ij: f64,
65}
66
67impl From<f64> for PetsBinaryRecord {
68    fn from(k_ij: f64) -> Self {
69        Self { k_ij }
70    }
71}
72
73#[derive(Serialize, Deserialize, Debug, Clone, Default)]
74pub struct PetsParameters {
75    pub molarweight: Array1<f64>,
76    pub sigma: Array1<f64>,
77    pub epsilon_k: Array1<f64>,
78    pub k_ij: Array2<f64>,
79    pub sigma_ij: Array2<f64>,
80    pub epsilon_k_ij: Array2<f64>,
81    pub e_k_ij: Array2<f64>,
82    pub viscosity: Option<Array2<f64>>,
83    pub diffusion: Option<Array2<f64>>,
84    pub thermal_conductivity: Option<Array2<f64>>,
85    pub pure_records: Vec<PureRecord<PetsRecord, JobackRecord>>,
86    pub joback_records: Option<Vec<JobackRecord>>,
87    pub binary_records: Array2<PetsBinaryRecord>,
88}
89
90impl Parameter for PetsParameters {
91    type Pure = PetsRecord;
92    type IdealGas = JobackRecord;
93    type Binary = PetsBinaryRecord;
94
95    fn from_records(
96        pure_records: Vec<PureRecord<Self::Pure, Self::IdealGas>>,
97        binary_records: Array2<PetsBinaryRecord>,
98    ) -> Self {
99        let n = pure_records.len();
100
101        let mut molarweight = Array::zeros(n);
102        let mut sigma = Array::zeros(n);
103        let mut epsilon_k = Array::zeros(n);
104        let mut viscosity = Vec::with_capacity(n);
105        let mut diffusion = Vec::with_capacity(n);
106        let mut thermal_conductivity = Vec::with_capacity(n);
107
108        let mut component_index = HashMap::with_capacity(n);
109
110        for (i, record) in pure_records.iter().enumerate() {
111            component_index.insert(record.identifier.clone(), i);
112            let r = &record.model_record;
113            sigma[i] = r.sigma;
114            epsilon_k[i] = r.epsilon_k;
115            viscosity.push(r.viscosity);
116            diffusion.push(r.diffusion);
117            thermal_conductivity.push(r.thermal_conductivity);
118            molarweight[i] = record.molarweight;
119        }
120
121        let k_ij = binary_records.map(|br| br.k_ij);
122        let mut epsilon_k_ij = Array::zeros((n, n));
123        let mut sigma_ij = Array::zeros((n, n));
124        let mut e_k_ij = Array::zeros((n, n));
125        for i in 0..n {
126            for j in 0..n {
127                e_k_ij[[i, j]] = (epsilon_k[i] * epsilon_k[j]).sqrt();
128                epsilon_k_ij[[i, j]] = (1.0 - k_ij[[i, j]]) * e_k_ij[[i, j]];
129                sigma_ij[[i, j]] = 0.5 * (sigma[i] + sigma[j]);
130            }
131        }
132
133        let viscosity_coefficients = if viscosity.iter().any(|v| v.is_none()) {
134            None
135        } else {
136            let mut v = Array2::zeros((4, viscosity.len()));
137            for (i, vi) in viscosity.iter().enumerate() {
138                v.column_mut(i).assign(&Array1::from(vi.unwrap().to_vec()));
139            }
140            Some(v)
141        };
142
143        let diffusion_coefficients = if diffusion.iter().any(|v| v.is_none()) {
144            None
145        } else {
146            let mut v = Array2::zeros((5, diffusion.len()));
147            for (i, vi) in diffusion.iter().enumerate() {
148                v.column_mut(i).assign(&Array1::from(vi.unwrap().to_vec()));
149            }
150            Some(v)
151        };
152
153        let thermal_conductivity_coefficients = if thermal_conductivity.iter().any(|v| v.is_none())
154        {
155            None
156        } else {
157            let mut v = Array2::zeros((4, thermal_conductivity.len()));
158            for (i, vi) in thermal_conductivity.iter().enumerate() {
159                v.column_mut(i).assign(&Array1::from(vi.unwrap().to_vec()));
160            }
161            Some(v)
162        };
163
164        let joback_records = pure_records
165            .iter()
166            .map(|r| r.ideal_gas_record.clone())
167            .collect();
168
169        Self {
170            molarweight,
171            sigma,
172            epsilon_k,
173            k_ij,
174            sigma_ij,
175            epsilon_k_ij,
176            e_k_ij,
177            viscosity: viscosity_coefficients,
178            diffusion: diffusion_coefficients,
179            thermal_conductivity: thermal_conductivity_coefficients,
180            pure_records,
181            joback_records,
182            binary_records,
183        }
184    }
185
186    fn records(
187        &self,
188    ) -> (
189        &[PureRecord<PetsRecord, JobackRecord>],
190        &Array2<PetsBinaryRecord>,
191    ) {
192        (&self.pure_records, &self.binary_records)
193    }
194}
195
196impl PetsParameters {
197    pub fn to_markdown(&self) -> String {
198        let mut output = String::new();
199        let o = &mut output;
200        write!(
201            o,
202            "|component|molarweight|$\\sigma$|$\\varepsilon$|\n|-|-|-|-|"
203        )
204        .unwrap();
205        for i in 0..self.sigma.len() {
206            let component = self.pure_records[i].identifier.name.clone();
207            let component = component.unwrap_or(format!("Component {}", i + 1));
208            write!(
209                o,
210                "\n|{}|{}|{}|{}|",
211                component, self.molarweight[i], self.sigma[i], self.epsilon_k[i],
212            )
213            .unwrap();
214        }
215
216        output
217    }
218}
219
220impl std::fmt::Display for PetsParameters {
221    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
222        write!(f, "PetsParameters(")?;
223        write!(f, "\n\tmolarweight={}", self.molarweight)?;
224        write!(f, "\n\tsigma={}", self.sigma)?;
225        write!(f, "\n\tepsilon_k={}", self.epsilon_k)?;
226        if !self.k_ij.iter().all(|k| k.is_zero()) {
227            write!(f, "\n\tk_ij=\n{}", self.k_ij)?;
228        }
229        write!(f, "\n)")
230    }
231}
232
233#[cfg(test)]
234pub mod utils {
235    use super::*;
236    use feos_core::joback::JobackRecord;
237    use std::rc::Rc;
238
239    pub fn argon_parameters() -> Rc<PetsParameters> {
240        let argon_json = r#"
241            {
242                "identifier": {
243                    "cas": "7440-37-1",
244                    "name": "argon",
245                    "iupac_name": "argon",
246                    "smiles": "[Ar]",
247                    "inchi": "InChI=1/Ar",
248                    "formula": "Ar"
249                },
250                "model_record": {
251                    "sigma": 3.4050,
252                    "epsilon_k": 119.8,
253                    "viscosity": [0.0, 0.0, 0.0, 0.0],
254                    "thermal_conductivity": [0.0, 0.0, 0.0, 0.0],
255                    "diffusion": [0.0, 0.0, 0.0, 0.0, 0.0]
256                },
257                "molarweight": 39.948
258            }"#;
259        let argon_record: PureRecord<PetsRecord, JobackRecord> =
260            serde_json::from_str(argon_json).expect("Unable to parse json.");
261        Rc::new(PetsParameters::new_pure(argon_record))
262    }
263
264    pub fn krypton_parameters() -> Rc<PetsParameters> {
265        let krypton_json = r#"
266            {
267                "identifier": {
268                    "cas": "7439-90-9",
269                    "name": "krypton",
270                    "iupac_name": "krypton",
271                    "smiles": "[Kr]",
272                    "inchi": "InChI=1S/Kr",
273                    "formula": "Kr"
274                },
275                "model_record": {
276                    "sigma": 3.6300,
277                    "epsilon_k": 163.10
278                },
279                "molarweight": 83.798
280            }"#;
281        let krypton_record: PureRecord<PetsRecord, JobackRecord> =
282            serde_json::from_str(krypton_json).expect("Unable to parse json.");
283        Rc::new(PetsParameters::new_pure(krypton_record))
284    }
285
286    pub fn argon_krypton_parameters() -> Rc<PetsParameters> {
287        let binary_json = r#"[
288            {
289                "identifier": {
290                    "cas": "7440-37-1",
291                    "name": "argon",
292                    "iupac_name": "argon",
293                    "smiles": "[Ar]",
294                    "inchi": "1/Ar",
295                    "formula": "Ar"
296                },
297                "model_record": {
298                    "sigma": 3.4050,
299                    "epsilon_k": 119.8,
300                    "viscosity": [0.0, 0.0, 0.0, 0.0],
301                    "thermal_conductivity": [0.0, 0.0, 0.0, 0.0],
302                    "diffusion": [0.0, 0.0, 0.0, 0.0, 0.0]
303                },
304                "molarweight": 39.948
305            },
306            {
307                "identifier": {
308                    "cas": "7439-90-9",
309                    "name": "krypton",
310                    "iupac_name": "krypton",
311                    "smiles": "[Kr]",
312                    "inchi": "InChI=1S/Kr",
313                    "formula": "Kr"
314                },
315                "model_record": {
316                    "sigma": 3.6300,
317                    "epsilon_k": 163.10,
318                    "viscosity": [0.0, 0.0, 0.0, 0.0],
319                    "thermal_conductivity": [0.0, 0.0, 0.0, 0.0],
320                    "diffusion": [0.0, 0.0, 0.0, 0.0, 0.0]
321                },
322                "molarweight": 83.798
323            }
324        ]"#;
325        let binary_record: Vec<PureRecord<PetsRecord, JobackRecord>> =
326            serde_json::from_str(binary_json).expect("Unable to parse json.");
327        Rc::new(PetsParameters::new_binary(binary_record, None))
328    }
329}