Skip to main content

conspire/geometry/mesh/tessellation/base/
mod.rs

1#[cfg(test)]
2mod test;
3
4use crate::{
5    geometry::{
6        Coordinates,
7        bvh::BoundingVolumeHierarchy,
8        mesh::{
9            Connectivity, Mesh,
10            smooth::Smoothing,
11            tessellation::{D, Normals, Tessellation, features::Features},
12        },
13    },
14    math::{Quantity, Scalar, Tensor, TensorVec},
15    units::Length,
16};
17use std::{array::from_fn, cell::OnceCell, collections::HashMap};
18
19const WELD_TOLERANCE: Scalar = 1e-6;
20
21impl Tessellation {
22    pub fn mesh(&self) -> &Mesh<D> {
23        &self.mesh
24    }
25    pub fn normals(&self) -> &Normals {
26        &self.normals
27    }
28    pub fn bvh(&self) -> &BoundingVolumeHierarchy<D> {
29        self.bvh
30            .get_or_init(|| BoundingVolumeHierarchy::from(&self.mesh))
31    }
32    pub fn features(&self) -> &Features {
33        self.features.get_or_init(|| Features::of(self))
34    }
35    pub fn smooth(&mut self, smoothing: Smoothing) -> Result<(), &'static str> {
36        self.mesh.smooth(smoothing)?;
37        self.refresh();
38        Ok(())
39    }
40    pub fn smooth_welded(&mut self, smoothing: Smoothing) -> Result<(), &'static str> {
41        let mut min = [f64::INFINITY; D];
42        let mut max = [f64::NEG_INFINITY; D];
43        for point in self.mesh.coordinates() {
44            (0..D).for_each(|axis| {
45                min[axis] = min[axis].min(point[axis].value());
46                max[axis] = max[axis].max(point[axis].value());
47            });
48        }
49        let diagonal = Quantity::new(
50            (0..D)
51                .map(|axis| (max[axis] - min[axis]).powi(2))
52                .sum::<Scalar>()
53                .sqrt(),
54        );
55        self.smooth_welded_with_tolerance(smoothing, diagonal * WELD_TOLERANCE)
56    }
57    pub(crate) fn smooth_welded_with_tolerance(
58        &mut self,
59        smoothing: Smoothing,
60        tolerance: Quantity<Length>,
61    ) -> Result<(), &'static str> {
62        let mut representatives = Vec::with_capacity(self.mesh.number_of_nodes());
63        let mut anchors: HashMap<[i64; D], Vec<usize>> = HashMap::new();
64        let mut welded = Coordinates::new();
65        for point in self.mesh.coordinates() {
66            let cell = from_fn(|axis| (point[axis] / tolerance).floor().value() as i64);
67            let mut representative = None;
68            'search: for dz in -1i64..=1 {
69                for dy in -1i64..=1 {
70                    for dx in -1i64..=1 {
71                        if let Some(indices) =
72                            anchors.get(&[cell[0] + dx, cell[1] + dy, cell[2] + dz])
73                        {
74                            for &index in indices {
75                                if (point - &welded[index]).norm() <= tolerance {
76                                    representative = Some(index);
77                                    break 'search;
78                                }
79                            }
80                        }
81                    }
82                }
83            }
84            representatives.push(representative.unwrap_or_else(|| {
85                let index = welded.len();
86                welded.push(point.clone());
87                anchors.entry(cell).or_default().push(index);
88                index
89            }));
90        }
91        let triangles = match &self.mesh.connectivities()[0] {
92            Connectivity::Triangular(triangles) => triangles
93                .iter()
94                .map(|&[a, b, c]| [representatives[a], representatives[b], representatives[c]])
95                .collect::<Vec<_>>(),
96            _ => panic!(),
97        };
98        let mut mesh = Mesh::from((vec![Connectivity::Triangular(triangles.into())], welded));
99        mesh.smooth(smoothing)?;
100        let smoothed = mesh.coordinates();
101        self.mesh
102            .coordinates
103            .iter_mut()
104            .zip(&representatives)
105            .for_each(|(point, &representative)| *point = smoothed[representative].clone());
106        self.refresh();
107        Ok(())
108    }
109    fn refresh(&mut self) {
110        self.normals = self.mesh.normals();
111        self.bvh = OnceCell::new();
112        self.features = OnceCell::new();
113    }
114}