Skip to main content

tess2_rust/mesh/
delaunay.rs

1// Copyright 2025 Lars Brubaker
2// Delaunay refinement methods for Mesh.
3
4use super::{EdgeIdx, Mesh, F_HEAD};
5use crate::geom::Real;
6
7impl Mesh {
8    /// Compute the in-circle predicate for Delaunay refinement.
9    pub fn in_circle(
10        v_s: Real,
11        v_t: Real,
12        v0_s: Real,
13        v0_t: Real,
14        v1_s: Real,
15        v1_t: Real,
16        v2_s: Real,
17        v2_t: Real,
18    ) -> Real {
19        let adx = v0_s - v_s;
20        let ady = v0_t - v_t;
21        let bdx = v1_s - v_s;
22        let bdy = v1_t - v_t;
23        let cdx = v2_s - v_s;
24        let cdy = v2_t - v_t;
25
26        let ab_det = adx * bdy - bdx * ady;
27        let bc_det = bdx * cdy - cdx * bdy;
28        let ca_det = cdx * ady - adx * cdy;
29
30        let a_lift = adx * adx + ady * ady;
31        let b_lift = bdx * bdx + bdy * bdy;
32        let c_lift = cdx * cdx + cdy * cdy;
33
34        a_lift * bc_det + b_lift * ca_det + c_lift * ab_det
35    }
36
37    /// Check if an edge is locally Delaunay.
38    pub fn edge_is_locally_delaunay(&self, e: EdgeIdx) -> bool {
39        let e_sym = e ^ 1;
40        let e_sym_lnext = self.edges[e_sym as usize].lnext;
41        let e_sym_lnext_lnext = self.edges[e_sym_lnext as usize].lnext;
42        let e_lnext = self.edges[e as usize].lnext;
43        let e_lnext_lnext = self.edges[e_lnext as usize].lnext;
44
45        let v = self.edges[e_sym_lnext_lnext as usize].org;
46        let v0 = self.edges[e_lnext as usize].org;
47        let v1 = self.edges[e_lnext_lnext as usize].org;
48        let v2 = self.edges[e as usize].org;
49
50        Self::in_circle(
51            self.verts[v as usize].s,
52            self.verts[v as usize].t,
53            self.verts[v0 as usize].s,
54            self.verts[v0 as usize].t,
55            self.verts[v1 as usize].s,
56            self.verts[v1 as usize].t,
57            self.verts[v2 as usize].s,
58            self.verts[v2 as usize].t,
59        ) < 0.0
60    }
61
62    /// Refine a valid triangulation into a Constrained Delaunay Triangulation.
63    pub fn refine_delaunay(&mut self) {
64        let mut stack: Vec<EdgeIdx> = Vec::new();
65
66        let mut f = self.faces[F_HEAD as usize].next;
67        while f != F_HEAD {
68            if self.faces[f as usize].inside {
69                let e_start = self.faces[f as usize].an_edge;
70                let mut e = e_start;
71                loop {
72                    let is_internal = self.edge_is_internal(e);
73                    self.edges[e as usize].mark = is_internal;
74                    if is_internal && !self.edges[(e ^ 1) as usize].mark {
75                        stack.push(e);
76                    }
77                    e = self.edges[e as usize].lnext;
78                    if e == e_start {
79                        break;
80                    }
81                }
82            }
83            f = self.faces[f as usize].next;
84        }
85
86        let max_iter = stack.len() * stack.len() + 1;
87        let mut iter = 0;
88
89        while let Some(e) = stack.pop() {
90            if iter >= max_iter {
91                break;
92            }
93            iter += 1;
94            self.edges[e as usize].mark = false;
95            self.edges[(e ^ 1) as usize].mark = false;
96
97            if !self.edge_is_locally_delaunay(e) {
98                let neighbors = [
99                    self.edges[e as usize].lnext,
100                    self.lprev(e),
101                    self.edges[(e ^ 1) as usize].lnext,
102                    self.lprev(e ^ 1),
103                ];
104                self.flip_edge(e);
105                for &nb in &neighbors {
106                    if !self.edges[nb as usize].mark && self.edge_is_internal(nb) {
107                        self.edges[nb as usize].mark = true;
108                        self.edges[(nb ^ 1) as usize].mark = true;
109                        stack.push(nb);
110                    }
111                }
112            }
113        }
114    }
115}