Skip to main content

xtalforge_core/
groups.rs

1//! Functional-group catalogue and matcher (Task 1.13).
2
3use crate::{BondGraph, Error, Structure};
4
5/// A functional group the matcher knows.
6#[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord)]
7pub enum GroupKind {
8    /// O-H.
9    Hydroxyl,
10    /// C with two oxygen neighbours.
11    Carboxylate,
12    /// N with one to three neighbours.
13    Amine,
14    /// C-F.
15    Fluoride,
16}
17
18impl GroupKind {
19    /// Stable name.
20    #[must_use]
21    pub const fn name(self) -> &'static str {
22        match self {
23            Self::Hydroxyl => "hydroxyl",
24            Self::Carboxylate => "carboxylate",
25            Self::Amine => "amine",
26            Self::Fluoride => "fluoride",
27        }
28    }
29}
30
31/// One matched group.
32#[derive(Debug, Clone, PartialEq, Eq)]
33pub struct GroupMatch {
34    /// Group type.
35    pub kind: GroupKind,
36    /// Atom indices in the P1 cell.
37    pub atoms: Vec<usize>,
38}
39
40/// Finds every occurrence of the catalogue on the perceived bond graph.
41///
42/// # Errors
43///
44/// Returns the errors of the P1 expansion or bond perception.
45pub fn find_groups(structure: &Structure) -> Result<Vec<GroupMatch>, Error> {
46    let graph = structure.bond_graph()?;
47    let p1 = structure.p1()?;
48    let element_of = |index: usize| {
49        p1.atoms
50            .get(index)
51            .map(|atom| atom.element.clone())
52            .unwrap_or_default()
53    };
54    let mut matches = Vec::new();
55    for index in 0..graph.n_atoms() {
56        let element = element_of(index);
57        let neighbours = neighbours(&graph, index);
58        match element.as_str() {
59            "O" if neighbours
60                .iter()
61                .any(|neighbour| element_of(*neighbour) == "H") =>
62            {
63                if let Some(hydrogen) = neighbours
64                    .iter()
65                    .find(|neighbour| element_of(**neighbour) == "H")
66                {
67                    matches.push(GroupMatch {
68                        kind: GroupKind::Hydroxyl,
69                        atoms: sorted(vec![index, *hydrogen]),
70                    });
71                }
72            }
73            "C" => {
74                let oxygens: Vec<usize> = neighbours
75                    .iter()
76                    .copied()
77                    .filter(|neighbour| element_of(*neighbour) == "O")
78                    .collect();
79                if oxygens.len() >= 2
80                    && let (Some(first), Some(second)) = (oxygens.first(), oxygens.get(1))
81                {
82                    matches.push(GroupMatch {
83                        kind: GroupKind::Carboxylate,
84                        atoms: sorted(vec![index, *first, *second]),
85                    });
86                }
87            }
88            "N" if (1..=3).contains(&neighbours.len()) => {
89                matches.push(GroupMatch {
90                    kind: GroupKind::Amine,
91                    atoms: vec![index],
92                });
93            }
94            "F" if neighbours.len() == 1 => {
95                if let Some(carbon) = neighbours.first() {
96                    matches.push(GroupMatch {
97                        kind: GroupKind::Fluoride,
98                        atoms: sorted(vec![index, *carbon]),
99                    });
100                }
101            }
102            _ => {}
103        }
104    }
105    matches.sort_by(|first, second| (first.kind, &first.atoms).cmp(&(second.kind, &second.atoms)));
106    matches.dedup();
107    Ok(matches)
108}
109
110/// Counts per group kind.
111///
112/// # Errors
113///
114/// Returns the errors of [`find_groups`].
115pub fn group_counts(structure: &Structure) -> Result<Vec<(GroupKind, usize)>, Error> {
116    let mut counts: Vec<(GroupKind, usize)> = Vec::new();
117    for group in find_groups(structure)? {
118        if let Some(entry) = counts.iter_mut().find(|(kind, _)| *kind == group.kind) {
119            entry.1 += 1;
120        } else {
121            counts.push((group.kind, 1));
122        }
123    }
124    Ok(counts)
125}
126
127fn neighbours(graph: &BondGraph, index: usize) -> Vec<usize> {
128    let mut out = Vec::new();
129    for bond_index in graph.bonds_of(index) {
130        if let Some(bond) = graph.bonds.get(*bond_index) {
131            let other = if bond.i == index { bond.j } else { bond.i };
132            if !out.contains(&other) {
133                out.push(other);
134            }
135        }
136    }
137    out
138}
139
140fn sorted(mut atoms: Vec<usize>) -> Vec<usize> {
141    atoms.sort_unstable();
142    atoms.dedup();
143    atoms
144}