conspire/geometry/mesh/quality/improve/untangle/
mod.rs

1#[cfg(test)]
2mod test;
3
4use super::Incidence;
5use crate::{
6    geometry::{
7        Coordinate,
8        mesh::{Mesh, Tessellation},
9    },
10    math::{Quantity, Scalar, Tensor},
11};
12use std::mem::transmute_copy;
13
14const PROBES: usize = 32;
15
16impl<const D: usize> Mesh<D> {
17    pub fn untangle(&mut self, iterations: usize, margin: Scalar, surface: Option<&Tessellation>) {
18        let number_of_nodes = self.number_of_nodes();
19        let neighbors = self.node_node_connectivity().to_vec();
20        let incidence = Incidence::of(self);
21        let constrained = surface.filter(|_| D == 3);
22        let elements: Vec<&[usize]> = constrained
23            .map(|surface| surface.mesh().connectivities().iter().flatten().collect())
24            .unwrap_or_default();
25        let mut boundary = vec![false; number_of_nodes];
26        if constrained.is_some() {
27            self.exterior_faces()
28                .iter()
29                .flatten()
30                .for_each(|&node| boundary[node] = true);
31        }
32        let coordinates = self.coordinates.members_mut();
33        for _ in 0..iterations {
34            for node in 0..number_of_nodes {
35                if neighbors[node].is_empty() {
36                    continue;
37                }
38                let mut current = incidence.inversion(node, coordinates, margin);
39                if current <= 0.0 {
40                    continue;
41                }
42                let mut step = 0.5
43                    * neighbors[node]
44                        .iter()
45                        .map(|&neighbor| {
46                            (&coordinates[node] - &coordinates[neighbor]).norm().value()
47                        })
48                        .sum::<Scalar>()
49                    / (neighbors[node].len() as Scalar);
50                for _ in 0..PROBES {
51                    let mut improved = false;
52                    for axis in 0..D {
53                        for sign in [-1.0, 1.0] {
54                            let original = coordinates[node].clone();
55                            coordinates[node][axis] += Quantity::new(sign * step);
56                            if boundary[node] {
57                                coordinates[node] = project_to_surface(
58                                    constrained.unwrap(),
59                                    &elements,
60                                    &coordinates[node],
61                                );
62                            }
63                            let trial = incidence.inversion(node, coordinates, margin);
64                            if trial < current {
65                                current = trial;
66                                improved = true;
67                            } else {
68                                coordinates[node] = original;
69                            }
70                        }
71                    }
72                    if !improved {
73                        step *= 0.5;
74                    }
75                }
76            }
77        }
78    }
79}
80
81fn project_to_surface<const D: usize>(
82    surface: &Tessellation,
83    elements: &[&[usize]],
84    point: &Coordinate<D>,
85) -> Coordinate<D> {
86    let query: &Coordinate<3> = unsafe { &*(point as *const Coordinate<D>).cast() };
87    let projected = surface
88        .bvh()
89        .closest_point(query, surface.mesh().coordinates(), elements)
90        .map_or_else(|| query.clone(), |(projected, _)| projected);
91    unsafe { transmute_copy(&projected) }
92}
conspire/geometry/mesh/quality/improve/untangle/mod.rs

conspire/geometry/mesh/quality/improve/untangle/
mod.rs
