From 0ec6e29005b51fc1cd86bf140801c417a63f85f5 Mon Sep 17 00:00:00 2001 From: mattatz Date: Fri, 20 Dec 2024 14:05:46 +0900 Subject: [PATCH 01/14] Create SurfaceBoundingBoxTree. --- src/bounding_box/mod.rs | 21 +++- src/bounding_box/surface_bounding_box_tree.rs | 97 +++++++++++++++++++ src/surface/mod.rs | 10 ++ src/surface/nurbs_surface.rs | 13 +++ 4 files changed, 140 insertions(+), 1 deletion(-) create mode 100644 src/bounding_box/surface_bounding_box_tree.rs diff --git a/src/bounding_box/mod.rs b/src/bounding_box/mod.rs index 6c5ea89..8253cd7 100644 --- a/src/bounding_box/mod.rs +++ b/src/bounding_box/mod.rs @@ -1,17 +1,22 @@ pub mod bounding_box_traversal; pub mod bounding_box_tree; pub mod curve_bounding_box_tree; +pub mod surface_bounding_box_tree; pub use bounding_box_traversal::*; pub use bounding_box_tree::*; pub use curve_bounding_box_tree::*; +pub use surface_bounding_box_tree::*; use nalgebra::{ allocator::Allocator, DefaultAllocator, DimName, DimNameDiff, DimNameSub, OPoint, OVector, U1, }; use simba::scalar::SupersetOf; -use crate::{curve::nurbs_curve::NurbsCurve, misc::FloatingPoint, region::CompoundCurve}; +use crate::{ + curve::nurbs_curve::NurbsCurve, misc::FloatingPoint, region::CompoundCurve, + surface::NurbsSurface, +}; /// A struct representing a bounding box in D space. #[derive(Clone, Debug)] @@ -214,3 +219,17 @@ where Self::from_iter(pts) } } + +impl<'a, T: FloatingPoint, D: DimName> From<&'a NurbsSurface> + for BoundingBox> +where + DefaultAllocator: Allocator, + D: DimNameSub, + DefaultAllocator: Allocator>, +{ + fn from(value: &'a NurbsSurface) -> Self { + let pts = value.dehomogenized_control_points(); + let flatten = pts.into_iter().flatten(); + Self::new_with_points(flatten) + } +} diff --git a/src/bounding_box/surface_bounding_box_tree.rs b/src/bounding_box/surface_bounding_box_tree.rs new file mode 100644 index 0000000..2d27d92 --- /dev/null +++ b/src/bounding_box/surface_bounding_box_tree.rs @@ -0,0 +1,97 @@ +use std::borrow::Cow; + +use nalgebra::{allocator::Allocator, DefaultAllocator, DimName, DimNameDiff, DimNameSub, U1}; + +use crate::{ + misc::FloatingPoint, + split::{Split, SplitSurfaceOption}, + surface::{NurbsSurface, UVDirection}, +}; + +use super::{BoundingBox, BoundingBoxTree}; + +/// A struct representing a bounding box tree with surface in D space. +#[derive(Clone)] +pub struct SurfaceBoundingBoxTree<'a, T: FloatingPoint, D: DimName> +where + DefaultAllocator: Allocator, +{ + surface: Cow<'a, NurbsSurface>, + tolerance: (T, T), + direction: UVDirection, +} + +impl<'a, T: FloatingPoint, D: DimName> SurfaceBoundingBoxTree<'a, T, D> +where + D: DimNameSub, + DefaultAllocator: Allocator, + DefaultAllocator: Allocator>, +{ + /// Create a new bounding box tree from a surface. + pub fn new( + surface: &'a NurbsSurface, + direction: UVDirection, + tolerance: Option<(T, T)>, + ) -> Self { + let tol = tolerance.unwrap_or_else(|| { + let (u, v) = surface.knots_domain_interval(); + let div = T::from_usize(64).unwrap(); + (u / div, v / div) + }); + Self { + surface: Cow::Borrowed(surface), + tolerance: tol, + direction, + } + } + + pub fn surface(&self) -> &NurbsSurface { + self.surface.as_ref() + } + + pub fn surface_owned(self) -> NurbsSurface { + self.surface.into_owned() + } +} + +impl<'a, T: FloatingPoint, D: DimName> BoundingBoxTree for SurfaceBoundingBoxTree<'a, T, D> +where + D: DimNameSub, + DefaultAllocator: Allocator, + DefaultAllocator: Allocator>, +{ + /// Check if the surface is dividable or not. + fn is_dividable(&self) -> bool { + let interval = self.surface.knots_domain_interval(); + let (utol, vtol) = self.tolerance; + interval.0 > utol || interval.1 > vtol + } + + /// Try to divide the surface into two parts. + fn try_divide(&self) -> anyhow::Result<(Self, Self)> { + let div = T::from_f64(0.5).unwrap(); + let t = self.surface().knots_domain_at(self.direction); + let mid = (t.0 + t.1) * div; + let (l, r) = self + .surface() + .try_split(SplitSurfaceOption::new(mid, self.direction))?; + let next = self.direction.opposite(); + Ok(( + Self { + surface: Cow::Owned(l), + tolerance: self.tolerance, + direction: next, + }, + Self { + surface: Cow::Owned(r), + tolerance: self.tolerance, + direction: next, + }, + )) + } + + /// Get the bounding box of the surface. + fn bounding_box(&self) -> BoundingBox> { + self.surface.as_ref().into() + } +} diff --git a/src/surface/mod.rs b/src/surface/mod.rs index bf9e885..9c124e7 100644 --- a/src/surface/mod.rs +++ b/src/surface/mod.rs @@ -10,6 +10,16 @@ pub enum UVDirection { V, } +impl UVDirection { + /// Get the opposite direction. + pub fn opposite(&self) -> Self { + match self { + Self::U => Self::V, + Self::V => Self::U, + } + } +} + /// The direction to flip the surface. #[derive(Debug, Clone, Copy, PartialEq, Eq)] pub enum FlipDirection { diff --git a/src/surface/nurbs_surface.rs b/src/surface/nurbs_surface.rs index 88fdea5..53bce9b 100644 --- a/src/surface/nurbs_surface.rs +++ b/src/surface/nurbs_surface.rs @@ -93,11 +93,24 @@ where self.v_knots.domain(self.v_degree) } + pub fn knots_domain_at(&self, direction: UVDirection) -> (T, T) { + match direction { + UVDirection::U => self.u_knots_domain(), + UVDirection::V => self.v_knots_domain(), + } + } + /// Get the u and v domain of the knot vector by degree pub fn knots_domain(&self) -> ((T, T), (T, T)) { (self.u_knots_domain(), self.v_knots_domain()) } + /// Get the u and v domain interval of the knot vector by degree + pub fn knots_domain_interval(&self) -> (T, T) { + let (u, v) = self.knots_domain(); + (u.1 - u.0, v.1 - v.0) + } + pub fn control_points(&self) -> &Vec>> { &self.control_points } From 1f886e0a4f31d9c1adc32f89c5ec9ddcd89d19a0 Mon Sep 17 00:00:00 2001 From: mattatz Date: Fri, 20 Dec 2024 14:20:42 +0900 Subject: [PATCH 02/14] Refactor Intersection trait. --- src/boolean/boolean_compound_curve.rs | 2 +- src/boolean/boolean_curve.rs | 2 +- src/boolean/degeneracies.rs | 2 +- src/curve/tests.rs | 2 +- .../compound_curve_intersection.rs | 5 +- src/intersection/curve_intersection.rs | 35 ++-- src/intersection/has_intersection.rs | 12 +- .../intersection_compound_curve.rs | 6 +- src/intersection/intersection_curve.rs | 4 +- .../intersection_surface_curve.rs | 165 ++++++++++++++++++ src/intersection/mod.rs | 10 +- .../surface_curve_intersection.rs | 1 - 12 files changed, 208 insertions(+), 38 deletions(-) create mode 100644 src/intersection/intersection_surface_curve.rs diff --git a/src/boolean/boolean_compound_curve.rs b/src/boolean/boolean_compound_curve.rs index 9f8991a..e1d9385 100644 --- a/src/boolean/boolean_compound_curve.rs +++ b/src/boolean/boolean_compound_curve.rs @@ -1,7 +1,7 @@ use argmin::core::ArgminFloat; use nalgebra::U3; -use crate::prelude::Intersection; +use crate::prelude::Intersects; use crate::{ curve::NurbsCurve, misc::FloatingPoint, prelude::CurveIntersectionSolverOptions, region::CompoundCurve, diff --git a/src/boolean/boolean_curve.rs b/src/boolean/boolean_curve.rs index 7d29911..c393e71 100644 --- a/src/boolean/boolean_curve.rs +++ b/src/boolean/boolean_curve.rs @@ -7,7 +7,7 @@ use argmin::core::ArgminFloat; use itertools::Itertools; use nalgebra::U3; -use crate::prelude::{CompoundCurveIntersection, HasIntersectionParameter, Intersection}; +use crate::prelude::{CompoundCurveIntersection, HasIntersectionParameter, Intersects}; use crate::region::CompoundCurve; use crate::{curve::NurbsCurve, misc::FloatingPoint, prelude::CurveIntersectionSolverOptions}; diff --git a/src/boolean/degeneracies.rs b/src/boolean/degeneracies.rs index 349e695..76a0436 100644 --- a/src/boolean/degeneracies.rs +++ b/src/boolean/degeneracies.rs @@ -14,7 +14,7 @@ pub enum Degeneracy { } impl Degeneracy { - pub fn new>( + pub fn new>( it: &I, a: &NurbsCurve>, b: &NurbsCurve>, diff --git a/src/curve/tests.rs b/src/curve/tests.rs index 9869991..6a18cc2 100644 --- a/src/curve/tests.rs +++ b/src/curve/tests.rs @@ -3,7 +3,7 @@ use nalgebra::{Point2, Rotation2, Translation2}; use crate::{ curve::NurbsCurve2D, misc::Transformable, - prelude::{CurveIntersectionSolverOptions, Intersection}, + prelude::{CurveIntersectionSolverOptions, Intersects}, }; use super::KnotStyle; diff --git a/src/intersection/compound_curve_intersection.rs b/src/intersection/compound_curve_intersection.rs index 9ac18f2..4ab6eb4 100644 --- a/src/intersection/compound_curve_intersection.rs +++ b/src/intersection/compound_curve_intersection.rs @@ -53,7 +53,7 @@ where } } -impl<'a, T: FloatingPoint, D: DimName> HasIntersectionParameter +impl<'a, T: FloatingPoint, D: DimName> HasIntersectionParameter for CompoundCurveIntersection<'a, T, D> where D: DimNameSub, @@ -69,7 +69,8 @@ where } } -impl<'a, T: FloatingPoint, D: DimName> HasIntersection, T> +impl<'a, T: FloatingPoint, D: DimName> + HasIntersection, Intersection<'a, T, D>, T, T> for CompoundCurveIntersection<'a, T, D> where D: DimNameSub, diff --git a/src/intersection/curve_intersection.rs b/src/intersection/curve_intersection.rs index 65f9c73..c69ffb6 100644 --- a/src/intersection/curve_intersection.rs +++ b/src/intersection/curve_intersection.rs @@ -2,39 +2,46 @@ use super::{has_intersection::HasIntersection, HasIntersectionParameter}; /// A struct representing the intersection of two curves. #[derive(Debug, Clone)] -pub struct CurveIntersection { - /// The point & parameter of the first curve at the intersection. - a: (P, T), - /// The point & parameter of the second curve at the intersection. - b: (P, T), +pub struct Intersection { + /// The point & parameter of the first object at the intersection. + a: (P, T0), + /// The point & parameter of the second object at the intersection. + b: (P, T1), } -impl CurveIntersection { - pub fn new(a: (P, T), b: (P, T)) -> Self { +impl Intersection { + pub fn new(a: (P, T0), b: (P, T1)) -> Self { Self { a, b } } - pub fn as_tuple(self) -> ((P, T), (P, T)) { + pub fn as_tuple(self) -> ((P, T0), (P, T1)) { (self.a, self.b) } } -impl HasIntersectionParameter for CurveIntersection { - fn a_parameter(&self) -> T { +impl HasIntersectionParameter + for Intersection +{ + fn a_parameter(&self) -> T0 { self.a.1 } - fn b_parameter(&self) -> T { + fn b_parameter(&self) -> T1 { self.b.1 } } -impl HasIntersection<(P, T), T> for CurveIntersection { - fn a(&self) -> &(P, T) { +impl HasIntersection<(P, T0), (P, T1), T0, T1> + for Intersection +{ + fn a(&self) -> &(P, T0) { &self.a } - fn b(&self) -> &(P, T) { + fn b(&self) -> &(P, T1) { &self.b } } + +/// A struct representing the intersection of two curves. +pub type CurveIntersection = Intersection; diff --git a/src/intersection/has_intersection.rs b/src/intersection/has_intersection.rs index d355ab8..31fef1f 100644 --- a/src/intersection/has_intersection.rs +++ b/src/intersection/has_intersection.rs @@ -1,11 +1,11 @@ /// A trait for types that have an intersection. -pub trait HasIntersection: HasIntersectionParameter { - fn a(&self) -> &V; - fn b(&self) -> &V; +pub trait HasIntersection: HasIntersectionParameter { + fn a(&self) -> &VA; + fn b(&self) -> &VB; } /// A trait for types that have intersection parameters. -pub trait HasIntersectionParameter { - fn a_parameter(&self) -> T; - fn b_parameter(&self) -> T; +pub trait HasIntersectionParameter { + fn a_parameter(&self) -> A; + fn b_parameter(&self) -> B; } diff --git a/src/intersection/intersection_compound_curve.rs b/src/intersection/intersection_compound_curve.rs index 1fa7e83..7b67a2d 100644 --- a/src/intersection/intersection_compound_curve.rs +++ b/src/intersection/intersection_compound_curve.rs @@ -6,10 +6,10 @@ use num_traits::Float; use crate::{curve::NurbsCurve, misc::FloatingPoint, region::CompoundCurve}; use super::{ - CompoundCurveIntersection, CurveIntersectionSolverOptions, HasIntersection, Intersection, + CompoundCurveIntersection, CurveIntersectionSolverOptions, HasIntersection, Intersects, }; -impl<'a, T, D> Intersection<'a, &'a NurbsCurve> for CompoundCurve +impl<'a, T, D> Intersects<'a, &'a NurbsCurve> for CompoundCurve where T: FloatingPoint + ArgminFloat, D: DimName + DimNameSub, @@ -66,7 +66,7 @@ where } } -impl<'a, T, D> Intersection<'a, &'a CompoundCurve> for CompoundCurve +impl<'a, T, D> Intersects<'a, &'a CompoundCurve> for CompoundCurve where T: FloatingPoint + ArgminFloat, D: DimName + DimNameSub, diff --git a/src/intersection/intersection_curve.rs b/src/intersection/intersection_curve.rs index 1f4a139..bc136cd 100644 --- a/src/intersection/intersection_curve.rs +++ b/src/intersection/intersection_curve.rs @@ -16,10 +16,10 @@ use crate::{ use super::{ CurveIntersection, CurveIntersectionBFGS, CurveIntersectionProblem, - CurveIntersectionSolverOptions, HasIntersection, Intersection, + CurveIntersectionSolverOptions, HasIntersection, Intersects, }; -impl<'a, T, D> Intersection<'a, &'a NurbsCurve> for NurbsCurve +impl<'a, T, D> Intersects<'a, &'a NurbsCurve> for NurbsCurve where T: FloatingPoint + ArgminFloat, D: DimName + DimNameSub, diff --git a/src/intersection/intersection_surface_curve.rs b/src/intersection/intersection_surface_curve.rs new file mode 100644 index 0000000..e4447f0 --- /dev/null +++ b/src/intersection/intersection_surface_curve.rs @@ -0,0 +1,165 @@ +use argmin::core::{ArgminFloat, State}; +use nalgebra::{ + allocator::Allocator, DefaultAllocator, DimName, DimNameDiff, DimNameSub, OPoint, U1, +}; + +use crate::{ + curve::NurbsCurve, + misc::FloatingPoint, + prelude::{BoundingBoxTraversal, CurveBoundingBoxTree, SurfaceBoundingBoxTree}, + surface::{NurbsSurface, UVDirection}, +}; + +use super::{CurveIntersection, CurveIntersectionSolverOptions, Intersects}; + +impl<'a, T, D> Intersects<'a, &'a NurbsCurve> for NurbsSurface +where + T: FloatingPoint + ArgminFloat, + D: DimName + DimNameSub, + DefaultAllocator: Allocator, + DefaultAllocator: Allocator>, +{ + type Output = anyhow::Result>, T>>>; + type Option = Option>; + + /// + #[allow(clippy::type_complexity)] + fn find_intersections( + &'a self, + other: &'a NurbsCurve, + option: Self::Option, + ) -> Self::Output { + let options = option.unwrap_or_default(); + + let div = T::one() / T::from_usize(options.knot_domain_division).unwrap(); + let interval = self.knots_domain_interval(); + let ta = SurfaceBoundingBoxTree::new( + self, + UVDirection::U, + Some((interval.0 * div, interval.1 * div)), + ); + let tb = CurveBoundingBoxTree::new(other, Some(other.knots_domain_interval() * div)); + + let traversed = BoundingBoxTraversal::try_traverse(ta, tb)?; + let a_domain = self.knots_domain(); + let b_domain = other.knots_domain(); + todo!(); + + /* + let intersections = traversed + .into_pairs_iter() + .filter_map(|(a, b)| { + let ca = a.surface_owned(); + let cb = b.curve_owned(); + + let problem = CurveIntersectionProblem::new(&ca, &cb); + + // let inv = T::from_f64(0.5).unwrap(); + // let d0 = ca.knots_domain(); + // let d1 = cb.knots_domain(); + + // Define initial parameter vector + let init_param = Vector2::::new( + ca.knots_domain().0, + cb.knots_domain().0, + // (d0.0 + d0.1) * inv, + // (d1.0 + d1.1) * inv, + ); + + // Set up solver + let solver = CurveIntersectionBFGS::::new() + .with_step_size_tolerance(options.step_size_tolerance) + .with_cost_tolerance(options.cost_tolerance); + + // Run solver + let res = Executor::new(problem, solver) + .configure(|state| { + state + .param(init_param) + .inv_hessian(Matrix2::identity()) + .max_iters(options.max_iters) + }) + .run(); + + match res { + Ok(r) => { + // println!("{}", r.state().get_termination_status()); + r.state().get_best_param().and_then(|param| { + if (a_domain.0..=a_domain.1).contains(¶m[0]) + && (b_domain.0..=b_domain.1).contains(¶m[1]) + { + let p0 = self.point_at(param[0]); + let p1 = other.point_at(param[1]); + Some(CurveIntersection::new((p0, param[0]), (p1, param[1]))) + } else { + None + } + }) + } + Err(_e) => { + // println!("{}", e); + None + } + } + }) + .filter(|it| { + // filter out intersections that are too close + let p0 = &it.a().0; + let p1 = &it.b().0; + let d = (p0 - p1).norm(); + d < options.minimum_distance + }) + .collect_vec(); + + let sorted = intersections + .into_iter() + .sorted_by(|x, y| x.a().1.partial_cmp(&y.a().1).unwrap_or(Ordering::Equal)) + .collect_vec(); + + // println!("sorted: {:?}", sorted.iter().map(|it| it.a()).collect_vec()); + + // group near parameter results & extract the closest one in each group + let parameter_minimum_distance = T::from_f64(1e-3).unwrap(); + let groups = sorted + .into_iter() + .map(|pt| vec![pt]) + .coalesce(|x, y| { + let x0 = &x[x.len() - 1]; + let y0 = &y[y.len() - 1]; + let da = Float::abs(x0.a().1 - y0.a().1); + let db = Float::abs(x0.b().1 - y0.b().1); + if da < parameter_minimum_distance || db < parameter_minimum_distance { + // merge near parameter results + let group = [x, y].concat(); + Ok(group) + } else { + Err((x, y)) + } + }) + .collect::>, T>>>>() + .into_iter() + .collect_vec(); + + let pts = groups + .into_iter() + .filter_map(|group| match group.len() { + 1 => Some(group[0].clone()), + _ => { + // find the closest intersection in the group + group + .iter() + .map(|it| { + let delta = &it.a().0 - &it.b().0; + let norm = delta.norm_squared(); + (it, norm) + }) + .min_by(|a, b| a.1.partial_cmp(&b.1).unwrap_or(Ordering::Equal)) + .map(|closest| closest.0.clone()) + } + }) + .collect_vec(); + + Ok(pts) + */ + } +} diff --git a/src/intersection/mod.rs b/src/intersection/mod.rs index a2bc660..aee60de 100644 --- a/src/intersection/mod.rs +++ b/src/intersection/mod.rs @@ -6,8 +6,8 @@ pub mod curve_intersection_solver_options; pub mod has_intersection; pub mod intersection_compound_curve; pub mod intersection_curve; -// pub mod intersection_surface_curve; -// pub mod surface_curve_intersection; +pub mod intersection_surface_curve; +pub mod surface_curve_intersection; pub use compound_curve_intersection::*; pub use curve_intersection::*; @@ -15,11 +15,9 @@ pub use curve_intersection_bfgs::*; pub use curve_intersection_problem::*; pub use curve_intersection_solver_options::*; pub use has_intersection::*; -// pub use intersection_surface_curve::*; -// pub use surface_curve_intersection::*; -/// Intersection between two curves trait -pub trait Intersection<'a, T> { +/// Intersection between two objects trait +pub trait Intersects<'a, T> { type Output; type Option; diff --git a/src/intersection/surface_curve_intersection.rs b/src/intersection/surface_curve_intersection.rs index 0c6958a..8b13789 100644 --- a/src/intersection/surface_curve_intersection.rs +++ b/src/intersection/surface_curve_intersection.rs @@ -1,2 +1 @@ -use super::{has_intersection::HasIntersection, HasIntersectionParameter}; From a62e60a0110843f312236bfd8fdc677b37c90e19 Mon Sep 17 00:00:00 2001 From: mattatz Date: Fri, 20 Dec 2024 14:21:54 +0900 Subject: [PATCH 03/14] Rename. --- src/intersection/{curve_intersection.rs => intersection.rs} | 0 src/intersection/mod.rs | 4 ++-- 2 files changed, 2 insertions(+), 2 deletions(-) rename src/intersection/{curve_intersection.rs => intersection.rs} (100%) diff --git a/src/intersection/curve_intersection.rs b/src/intersection/intersection.rs similarity index 100% rename from src/intersection/curve_intersection.rs rename to src/intersection/intersection.rs diff --git a/src/intersection/mod.rs b/src/intersection/mod.rs index aee60de..d77b112 100644 --- a/src/intersection/mod.rs +++ b/src/intersection/mod.rs @@ -1,5 +1,5 @@ pub mod compound_curve_intersection; -pub mod curve_intersection; +pub mod intersection; pub mod curve_intersection_bfgs; pub mod curve_intersection_problem; pub mod curve_intersection_solver_options; @@ -10,7 +10,7 @@ pub mod intersection_surface_curve; pub mod surface_curve_intersection; pub use compound_curve_intersection::*; -pub use curve_intersection::*; +pub use intersection::*; pub use curve_intersection_bfgs::*; pub use curve_intersection_problem::*; pub use curve_intersection_solver_options::*; From 90543338a9464111f920d5fb5be874ecc7f18a69 Mon Sep 17 00:00:00 2001 From: mattatz Date: Fri, 20 Dec 2024 15:17:03 +0900 Subject: [PATCH 04/14] WIP. --- Cargo.toml | 9 +- examples/intersect_surface_curve.rs | 310 ++++++++++++++++++ src/boolean/degeneracies.rs | 2 +- src/intersection/mod.rs | 25 -- .../surface_curve_intersection.rs | 1 - .../compound_curve_intersection.rs | 10 +- .../curve_curve}/curve_intersection_bfgs.rs | 30 +- .../curve_intersection_problem.rs | 0 .../curve_intersection_solver_options.rs | 0 .../intersection_compound_curve.rs | 11 +- .../curve_curve/intersection_curve_curve.rs} | 12 +- src/intersects/curve_curve/mod.rs | 16 + .../has_intersection.rs | 0 .../intersection.rs | 6 +- src/intersects/mod.rs | 16 + .../intersection_surface_curve.rs | 85 +++-- src/intersects/surface_curve/mod.rs | 11 + .../surface_curve_intersection_bfgs.rs | 208 ++++++++++++ .../surface_curve_intersection_problem.rs | 81 +++++ src/lib.rs | 4 +- 20 files changed, 732 insertions(+), 105 deletions(-) create mode 100644 examples/intersect_surface_curve.rs delete mode 100644 src/intersection/mod.rs delete mode 100644 src/intersection/surface_curve_intersection.rs rename src/{intersection => intersects/curve_curve}/compound_curve_intersection.rs (89%) rename src/{intersection => intersects/curve_curve}/curve_intersection_bfgs.rs (89%) rename src/{intersection => intersects/curve_curve}/curve_intersection_problem.rs (100%) rename src/{intersection => intersects/curve_curve}/curve_intersection_solver_options.rs (100%) rename src/{intersection => intersects/curve_curve}/intersection_compound_curve.rs (96%) rename src/{intersection/intersection_curve.rs => intersects/curve_curve/intersection_curve_curve.rs} (94%) create mode 100644 src/intersects/curve_curve/mod.rs rename src/{intersection => intersects}/has_intersection.rs (100%) rename src/{intersection => intersects}/intersection.rs (83%) create mode 100644 src/intersects/mod.rs rename src/{intersection => intersects/surface_curve}/intersection_surface_curve.rs (62%) create mode 100644 src/intersects/surface_curve/mod.rs create mode 100644 src/intersects/surface_curve/surface_curve_intersection_bfgs.rs create mode 100644 src/intersects/surface_curve/surface_curve_intersection_problem.rs diff --git a/Cargo.toml b/Cargo.toml index 47f173e..cc66db6 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -44,8 +44,8 @@ approx = { version = "0.5", default-features = false } serde_json = "1.0.133" [features] -default = [] -# default = ["bevy"] # for debugging example +# default = [] +default = ["bevy"] # for debugging example bevy = [ "dep:easer", "dep:rand_distr", @@ -136,6 +136,11 @@ name = "intersect_curves" path = "examples/intersect_curves.rs" required-features = ["bevy"] +[[example]] +name = "intersect_surface_curve" +path = "examples/intersect_surface_curve.rs" +required-features = ["bevy"] + [[example]] name = "iso_curves" path = "examples/iso_curves.rs" diff --git a/examples/intersect_surface_curve.rs b/examples/intersect_surface_curve.rs new file mode 100644 index 0000000..36714ff --- /dev/null +++ b/examples/intersect_surface_curve.rs @@ -0,0 +1,310 @@ +use std::f32::consts::FRAC_PI_2; + +use bevy::{ + color::palettes::{css::TOMATO, tailwind::LIME_500}, + prelude::*, + render::mesh::{PrimitiveTopology, VertexAttributeValues}, +}; +use bevy_infinite_grid::{InfiniteGridBundle, InfiniteGridPlugin, InfiniteGridSettings}; + +use bevy_panorbit_camera::{PanOrbitCamera, PanOrbitCameraPlugin}; +use bevy_points::{plugin::PointsPlugin, prelude::PointsMaterial}; +use nalgebra::{Point2, Rotation2, Translation2, Vector2}; + +use curvo::prelude::*; + +mod materials; + +use materials::*; + +fn main() { + App::new() + .add_plugins(DefaultPlugins.set(WindowPlugin { + primary_window: Some(Window { + resolution: (640., 480.).into(), + ..Default::default() + }), + ..Default::default() + })) + .add_plugins(LineMaterialPlugin) + .add_plugins(InfiniteGridPlugin) + .add_plugins(PanOrbitCameraPlugin) + .add_plugins(PointsPlugin) + .add_plugins(AppPlugin) + .run(); +} +struct AppPlugin; + +impl Plugin for AppPlugin { + fn build(&self, app: &mut bevy::prelude::App) { + app.add_systems(Startup, setup).add_systems(Update, update); + } +} + +#[derive(Component)] +struct FirstCurve(pub NurbsCurve2D); + +#[derive(Component)] +struct SecondCurve(pub NurbsCurve2D); + +fn setup( + mut commands: Commands, + mut meshes: ResMut>, + mut line_materials: ResMut>, + _points_materials: ResMut>, +) { + let points = vec![ + Point2::new(-1.0, -1.0), + Point2::new(1.0, -1.0), + Point2::new(1.0, 0.0), + Point2::new(-1.0, 0.0), + Point2::new(-1.0, 1.0), + Point2::new(1.0, 1.0), + ]; + let curve = NurbsCurve2D::try_interpolate(&points, 3).unwrap(); + + let line_vertices = curve + .tessellate(Some(1e-8)) + .iter() + .map(|p| p.cast::()) + .map(|p| [p.x, p.y, 0.]) + .collect(); + let line = Mesh::new(PrimitiveTopology::LineStrip, default()).with_inserted_attribute( + Mesh::ATTRIBUTE_POSITION, + VertexAttributeValues::Float32x3(line_vertices), + ); + commands + .spawn(( + Mesh3d(meshes.add(line)), + MeshMaterial3d(line_materials.add(LineMaterial { + color: TOMATO.into(), + ..Default::default() + })), + )) + .insert(FirstCurve(curve.clone())) + .insert(Name::new("curve")); + + let circle = + NurbsCurve2D::try_circle(&Point2::origin(), &Vector2::x(), &Vector2::y(), 1.).unwrap(); + + let line_vertices = circle + .tessellate(Some(1e-8)) + .iter() + .map(|p| p.cast::()) + .map(|p| [p.x, p.y, 0.]) + .collect(); + let line = Mesh::new(PrimitiveTopology::LineStrip, default()).with_inserted_attribute( + Mesh::ATTRIBUTE_POSITION, + VertexAttributeValues::Float32x3(line_vertices), + ); + commands + .spawn(( + Mesh3d(meshes.add(line)), + MeshMaterial3d(line_materials.add(LineMaterial { + color: LIME_500.into(), + ..Default::default() + })), + )) + .insert(SecondCurve(circle.clone())) + .insert(Name::new("circle")); + + /* + let traversed = BoundingBoxTraversal::try_traverse(&curve, &circle, None, None); + if let Ok(traversed) = traversed { + let n = traversed.pairs_iter().count(); + traversed + .pairs_iter() + .enumerate() + .for_each(|(idx, (a, b))| { + let t = (idx as f32) / (n as f32); + let hue = t * 360. * 1e2 % 360.; + let color = Color::hsl(hue, 0.5, 0.5); + let b0 = a.bounding_box(); + let b1 = b.bounding_box(); + let vertices = b0 + .lines() + .iter() + .chain(b1.lines().iter()) + .flat_map(|(a, b)| { + [a, b] + .iter() + .map(|p| p.cast::()) + .map(|p| [p.x, p.y, 0.]) + .collect::>() + }) + .collect(); + let line = Mesh::new(PrimitiveTopology::LineList, default()) + .with_inserted_attribute( + Mesh::ATTRIBUTE_POSITION, + VertexAttributeValues::Float32x3(vertices), + ); + commands + .spawn(MaterialMeshBundle { + mesh: meshes.add(line), + material: line_materials.add(LineMaterial { + // color: Color::WHITE, + color, + opacity: 0.5, + alpha_mode: AlphaMode::Blend, + }), + // visibility: Visibility::Hidden, + ..Default::default() + }) + .insert(Name::new("bounding box")); + + [a.curve(), b.curve()].iter().for_each(|curve| { + let line_vertices = curve + .tessellate(Some(1e-8)) + .iter() + .map(|p| p.cast::()) + .map(|p| [p.x, p.y, 0.]) + .collect(); + let line = Mesh::new(PrimitiveTopology::LineStrip, default()) + .with_inserted_attribute( + Mesh::ATTRIBUTE_POSITION, + VertexAttributeValues::Float32x3(line_vertices), + ); + commands + .spawn(MaterialMeshBundle { + mesh: meshes.add(line), + material: line_materials.add(LineMaterial { + color: Color::WHITE, + ..Default::default() + }), + // visibility: Visibility::Hidden, + ..Default::default() + }) + .insert(Name::new("segment")); + }); + }); + } + */ + + /* + let intersections = curve.find_intersections(&circle, None); + if let Ok(intersections) = intersections { + commands + .spawn(MaterialMeshBundle { + mesh: meshes.add(PointsMesh::from_iter(intersections.iter().flat_map(|it| { + [ + Vec3::from(it.a().0.cast::().coords.to_homogeneous()), + Vec3::from(it.b().0.cast::().coords.to_homogeneous()), + ] + }))), + material: points_materials.add(PointsMaterial { + settings: PointsShaderSettings { + point_size: 0.025, + color: Color::WHITE, + ..Default::default() + }, + circle: true, + ..Default::default() + }), + // visibility: Visibility::Hidden, + ..Default::default() + }) + .insert(Name::new("intersection points")); + } + */ + + let center = Vec3::ZERO; + commands.spawn(( + Transform::from_translation(center + Vec3::new(0., 0., 8.)).looking_at(center, Vec3::Y), + PanOrbitCamera::default(), + )); + commands.spawn(InfiniteGridBundle { + settings: InfiniteGridSettings { + x_axis_color: Color::BLACK, + z_axis_color: Color::BLACK, + ..Default::default() + }, + transform: Transform::from_rotation(Quat::from_rotation_x(FRAC_PI_2)), + ..Default::default() + }); +} + +#[allow(clippy::type_complexity)] +fn update( + time: Res