conspire/geometry/mesh/partition/agglomerate/
mod.rs

1#[cfg(test)]
2mod test;
3
4use super::Partition;
5use crate::{
6    geometry::{
7        Coordinates,
8        mesh::{Connectivities, Connectivity, Mesh},
9    },
10    math::{Set, Tensor, TensorVec},
11};
12use std::collections::HashMap;
13
14impl Partition {
15    pub fn agglomerate(&self, mesh: &Mesh<3>) -> Result<Mesh<3>, &'static str> {
16        let elements_faces = outward_faces(mesh)?;
17        let mut faces_nodes: Vec<Vec<usize>> = Vec::new();
18        let mut indices: HashMap<Vec<usize>, usize> = HashMap::new();
19        let mut cells = Vec::new();
20        for part in 0..self.number_of_parts() {
21            let elements = self.part_elements(part);
22            if elements.is_empty() {
23                continue;
24            }
25            let mut roots = (0..elements.len()).collect::<Vec<_>>();
26            let mut tally: HashMap<Vec<usize>, (usize, usize)> = HashMap::new();
27            let mut order: Vec<(&Vec<usize>, Vec<usize>)> = Vec::new();
28            for (position, &element) in elements.iter().enumerate() {
29                for face in &elements_faces[element] {
30                    let mut key = face.clone();
31                    key.sort_unstable();
32                    match tally.get_mut(&key) {
33                        Some((count, first)) => {
34                            *count += 1;
35                            let (a, b) = (find(&mut roots, *first), find(&mut roots, position));
36                            roots[b] = a;
37                        }
38                        None => {
39                            tally.insert(key.clone(), (1, position));
40                            order.push((face, key))
41                        }
42                    }
43                }
44            }
45            let root = find(&mut roots, 0);
46            if (1..elements.len()).any(|position| find(&mut roots, position) != root) {
47                return Err("a part is not connected through shared faces");
48            }
49            cells.push(
50                order
51                    .into_iter()
52                    .filter(|(_, key)| tally[key].0 == 1)
53                    .map(|(face, key)| {
54                        *indices.entry(key).or_insert_with(|| {
55                            faces_nodes.push(face.clone());
56                            faces_nodes.len() - 1
57                        })
58                    })
59                    .collect::<Vec<usize>>(),
60            );
61        }
62        let mut remap = vec![usize::MAX; mesh.number_of_nodes()];
63        faces_nodes
64            .iter()
65            .flatten()
66            .for_each(|&node| remap[node] = 0);
67        let mut coordinates = Coordinates::new();
68        remap.iter_mut().enumerate().for_each(|(node, new)| {
69            if *new == 0 {
70                *new = coordinates.len();
71                coordinates.push(mesh.coordinates()[node].clone())
72            }
73        });
74        faces_nodes
75            .iter_mut()
76            .flatten()
77            .for_each(|node| *node = remap[*node]);
78        Ok((
79            Connectivities::from(vec![Connectivity::Polyhedral((cells, faces_nodes).into())]),
80            Set::from(coordinates),
81        )
82            .into())
83    }
84}
85
86fn find(roots: &mut [usize], mut node: usize) -> usize {
87    while roots[node] != node {
88        roots[node] = roots[roots[node]];
89        node = roots[node]
90    }
91    node
92}
93
94fn outward_faces(mesh: &Mesh<3>) -> Result<Vec<Vec<Vec<usize>>>, &'static str> {
95    let mut elements_faces = Vec::with_capacity(mesh.number_of_elements());
96    for block in mesh.iter() {
97        match block {
98            Connectivity::Polyhedral(connectivity) => {
99                let mut owner = HashMap::new();
100                connectivity
101                    .elements_faces()
102                    .iter()
103                    .enumerate()
104                    .for_each(|(element, faces)| {
105                        faces.iter().for_each(|&face| {
106                            owner.entry(face).or_insert(element);
107                        })
108                    });
109                connectivity
110                    .elements_faces()
111                    .iter()
112                    .enumerate()
113                    .for_each(|(element, faces)| {
114                        elements_faces.push(
115                            faces
116                                .iter()
117                                .map(|&face| {
118                                    let mut nodes = connectivity.faces_nodes()[face].clone();
119                                    if owner[&face] != element {
120                                        nodes.reverse()
121                                    }
122                                    nodes
123                                })
124                                .collect(),
125                        )
126                    })
127            }
128            Connectivity::Polygonal(_)
129            | Connectivity::Quadrilateral(_)
130            | Connectivity::Triangular(_) => {
131                return Err("agglomeration requires three-dimensional elements");
132            }
133            _ => block
134                .iter()
135                .for_each(|element| elements_faces.push(block.element_faces(element))),
136        }
137    }
138    Ok(elements_faces)
139}
conspire/geometry/mesh/partition/agglomerate/mod.rs

conspire/geometry/mesh/partition/agglomerate/
mod.rs
