Skip to main content

geop_core_sketch/space/
wire.rs

1//! A 3-D sketch's drawn curves as the edges of a wire. A point two drawn
2//! curves pass through — an end of one, a point a spline is drawn through,
3//! an arc's middle point, or points made one by a
4//! [`Constraint3d::Coincident`](super::Constraint3d::Coincident) — is one
5//! vertex of the wire, and every curve passing through it is cut there: the
6//! wire is a graph whose edges meet only at its vertices, which is what the
7//! kernel requires of a body's edges. Curves that meet anywhere else are
8//! refused by name, since a crossing that is no point of the sketch is one
9//! the user cannot see.
10
11use std::collections::{BTreeMap, BTreeSet};
12
13use geop_core_geometry::{
14    intersection::curve_curve::curve_curve_intersect,
15    nurb_curve::{NurbCurve, NurbCurve3D},
16};
17use geop_core_math::{
18    geop_error::{GeopError, GeopResult, WithContext},
19    scalars::Scalar,
20    vector::Vector3,
21    with_context,
22};
23
24use super::{
25    CurveKind3d, Enclosure3d, Sketch3d,
26    curves::{arc_pieces, halvings},
27    geometry::Arc3,
28};
29use crate::{CurveId, PointId};
30
31/// Budgets for the crossing search of every pair of drawn curves: how many
32/// crossings one pair may have, how hard a single search tries, and the
33/// subdivision size it stops isolating at — one part in `MIN_SUBDIVISION`
34/// of a curve's domain.
35const MAX_CROSSINGS: usize = 16;
36const MAX_NODES: usize = 4000;
37const MIN_SUBDIVISION: i64 = 1_000_000;
38
39/// An edge of the wire a 3-D sketch is built as: the drawn curve `curve`
40/// from its point `start` to its point `end` — the whole curve, or, `cut`,
41/// one of the pieces it is cut into at the points it shares with other
42/// drawn curves. Points are the lowest id of their class (see
43/// [`Sketch3d::point_classes`]).
44#[derive(Clone, Copy, Debug, PartialEq, Eq)]
45pub struct WireEdge {
46    pub curve: CurveId,
47    pub start: PointId,
48    pub end: PointId,
49    pub cut: bool,
50}
51
52impl Sketch3d {
53    /// The indices of the points of `curve` between its ends that another
54    /// drawn curve passes through too: where the wire cuts it.
55    fn junctions(
56        &self,
57        curve: CurveId,
58        class: &BTreeMap<PointId, PointId>,
59    ) -> GeopResult<Vec<usize>> {
60        let points = self.curve(curve)?.points();
61        // What defines a helix between its ends is its axis, off the curve.
62        if matches!(self.curve(curve)?.kind, CurveKind3d::Helix { .. }) {
63            return Ok(Vec::new());
64        }
65        let others: BTreeSet<PointId> = self
66            .drawn_curves()
67            .filter(|(id, _)| **id != curve)
68            .flat_map(|(_, c)| c.points())
69            .map(|p| class[&p])
70            .collect();
71        Ok((1..points.len().saturating_sub(1))
72            .filter(|&k| others.contains(&class[&points[k]]))
73            .collect())
74    }
75
76    /// The edges of the wire (see the module docs): oldest curve first, each
77    /// curve's from its start to its end.
78    pub fn wire_edges(&self) -> GeopResult<Vec<WireEdge>> {
79        let class = self.point_classes();
80        let mut edges = Vec::new();
81        for (&curve, c) in self.drawn_curves() {
82            let points = c.points();
83            if points.len() < 2 {
84                return Err(GeopError::new(format!(
85                    "the drawn curve {curve} has no ends"
86                )));
87            }
88            let junctions = self.junctions(curve, &class)?;
89            let cut = !junctions.is_empty();
90            let stations: Vec<PointId> = (std::iter::once(0).chain(junctions))
91                .chain(std::iter::once(points.len() - 1))
92                .map(|k| class[&points[k]])
93                .collect();
94            edges.extend(stations.windows(2).map(|w| WireEdge {
95                curve,
96                start: w[0],
97                end: w[1],
98                cut,
99            }));
100        }
101        Ok(edges)
102    }
103
104    /// The edges of `curve` (see [`Sketch3d::wire_edges`]) as NURBS on
105    /// `[0, 1]`, in order, built from `geometry`: the whole curve, or its
106    /// pieces between the points it is cut at — a spline split at the knot
107    /// of each such point, an arc rebuilt from its middle point to each end.
108    /// An edge of several pieces (an arc's quarters) is joined into one
109    /// curve; an edge of one piece is that piece as it is, its weights
110    /// included: a spline's weights of one come out of each knot inserted
111    /// an ulp wide, which is honest and bars nothing.
112    pub fn wire_edge_nurbs(
113        &self,
114        curve: CurveId,
115        geometry: &Enclosure3d,
116    ) -> GeopResult<Vec<NurbCurve3D>> {
117        let ctx = with_context!("3-D sketch curve {curve} as the wire's edges");
118        let edge = |pieces: &[NurbCurve3D]| match pieces {
119            [one] => one.with_unit_domain(),
120            many => NurbCurve::join(many),
121        };
122        let junctions = self.junctions(curve, &self.point_classes())?;
123        let whole = self.curve_nurbs(curve, geometry).with_context(ctx)?;
124        if junctions.is_empty() {
125            return Ok(vec![edge(&whole).with_context(ctx)?]);
126        }
127        let at = |p: &PointId| geometry.points[p];
128        let pieces: Vec<Vec<NurbCurve3D>> = match &self.curve(curve)?.kind {
129            CurveKind3d::Spline { .. } => {
130                let [spline] = whole.as_slice() else {
131                    return Err(GeopError::new("a spline is one NURBS")).with_context(ctx);
132                };
133                // Knots `0 0 0 0 t1 .. t(m-2) 1 1 1 1`: the point `k` is
134                // passed at the knot `3 + k` (see `NurbCurve::cubic_spline`).
135                let knots = spline.knot_vector.clone();
136                let mut rest = spline.clone();
137                let mut pieces = Vec::new();
138                for k in junctions {
139                    let (left, right) = rest.split(knots[3 + k]).with_context(ctx)?;
140                    pieces.push(vec![left]);
141                    rest = right;
142                }
143                pieces.push(vec![rest]);
144                pieces
145            }
146            CurveKind3d::Arc {
147                start,
148                through,
149                end,
150            } => {
151                let drawn = Arc3 {
152                    s: self.points[start].at,
153                    m: self.points[through].at,
154                    e: self.points[end].at,
155                };
156                let arc = Arc3 {
157                    s: at(start),
158                    m: at(through),
159                    e: at(end),
160                };
161                let normal = arc.normal().with_context(ctx)?;
162                let halves = [
163                    ((&arc.s, &arc.m), (&drawn.s, &drawn.m)),
164                    ((&arc.m, &arc.e), (&drawn.m, &drawn.e)),
165                ];
166                halves
167                    .into_iter()
168                    .map(|((p, q), (drawn_p, drawn_q))| {
169                        arc_pieces(
170                            p,
171                            q,
172                            &normal,
173                            arc.half_sweep_between(p, q)?,
174                            halvings(drawn.half_sweep_between(drawn_p, drawn_q)?.0)?,
175                        )
176                    })
177                    .collect::<GeopResult<_>>()
178                    .with_context(ctx)?
179            }
180            CurveKind3d::Line { .. } | CurveKind3d::Helix { .. } | CurveKind3d::Reference => {
181                return Err(GeopError::new(format!(
182                    "curve {curve} has no points between its ends to be cut at"
183                )))
184                .with_context(ctx);
185            }
186        };
187        (pieces.iter().map(|p| edge(p)))
188            .collect::<GeopResult<_>>()
189            .with_context(ctx)
190    }
191
192    /// Every point of the sketch that is one — the lowest id of its class
193    /// (see [`Sketch3d::point_classes`]) — and where it is: the enclosure of
194    /// all the points made one with it, at `geometry`.
195    pub fn point_enclosures(&self, geometry: &Enclosure3d) -> BTreeMap<PointId, Vector3> {
196        let mut out: BTreeMap<PointId, Vector3> = BTreeMap::new();
197        for (p, &rep) in &self.point_classes() {
198            let at = geometry.points[p];
199            out.entry(rep)
200                .and_modify(|e| *e = e.union(&at))
201                .or_insert(at);
202        }
203        out
204    }
205
206    /// Checks the wire's edges meet only at its vertices (see the module
207    /// docs): two points of the sketch are nowhere at one place unless they
208    /// are one point, and two drawn curves meet only at points both pass
209    /// through. Says which points or curves do otherwise, and where.
210    pub fn check_wire(&self, geometry: &Enclosure3d) -> GeopResult<()> {
211        let enclosures = self.point_enclosures(geometry);
212        let points: Vec<(&PointId, &Vector3)> = enclosures.iter().collect();
213        for (i, (p, at)) in points.iter().enumerate() {
214            for (q, other) in &points[i + 1..] {
215                if at.could_be_equal(other) {
216                    return Err(GeopError::new(format!(
217                        "points {p} and {q} are both at {at:?}, but are two points of the sketch: \
218                         draw through one of them, or make them coincident"
219                    )));
220                }
221            }
222        }
223        let class = self.point_classes();
224        let min_subdivision = Scalar::from_ratio(1, MIN_SUBDIVISION)?;
225        let curves: Vec<(CurveId, BTreeSet<PointId>, NurbCurve3D)> = (self.drawn_curves())
226            .map(|(&id, c)| {
227                let nurbs = NurbCurve::join(&self.curve_nurbs(id, geometry)?)?;
228                Ok((id, c.points().iter().map(|p| class[p]).collect(), nurbs))
229            })
230            .collect::<GeopResult<_>>()?;
231        for (i, (a, a_points, a_curve)) in curves.iter().enumerate() {
232            for (b, b_points, b_curve) in &curves[i + 1..] {
233                let ctx = with_context!("curves {a} and {b} of the 3-D sketch");
234                let hits = curve_curve_intersect(
235                    a_curve,
236                    b_curve,
237                    MAX_CROSSINGS,
238                    MAX_NODES,
239                    min_subdivision,
240                )
241                .with_context(ctx)?;
242                if hits.is_coincident() {
243                    return Err(GeopError::new(format!(
244                        "curves {a} and {b} run along each other"
245                    )));
246                }
247                for (ta, tb) in hits.into_vec() {
248                    let at = (a_curve.evaluate(ta).with_context(ctx)?)
249                        .union(&b_curve.evaluate(tb).with_context(ctx)?);
250                    let passes = |p: &PointId| enclosures[p].could_be_equal(&at);
251                    if a_points.intersection(b_points).any(passes) {
252                        continue;
253                    }
254                    let of =
255                        |points: &BTreeSet<PointId>| points.iter().copied().find(|p| passes(p));
256                    return Err(match (of(a_points), of(b_points)) {
257                        (Some(p), _) => GeopError::new(format!(
258                            "curves {a} and {b} cross at point {p} of curve {a}, which curve {b} does \
259                             not pass through: draw curve {b} through it"
260                        )),
261                        (None, Some(p)) => GeopError::new(format!(
262                            "curves {a} and {b} cross at point {p} of curve {b}, which curve {a} does \
263                             not pass through: draw curve {a} through it"
264                        )),
265                        (None, None) => GeopError::new(format!(
266                            "curves {a} and {b} cross at {at:?}, where neither has a point: draw both \
267                             through a point there"
268                        )),
269                    });
270                }
271            }
272        }
273        Ok(())
274    }
275}