hoomd_geometry/shape/
convex_polytope.rs

1// Copyright (c) 2024-2026 The Regents of the University of Michigan.
2// Part of hoomd-rs, released under the BSD 3-Clause License.
3
4//! N-Dimensional generalization of a convex polyhedron.
5
6use serde::{Deserialize, Serialize};
7
8use crate::{BoundingSphereRadius, Error, SupportMapping};
9use arrayvec::ArrayVec;
10use hoomd_utility::valid::PositiveReal;
11use hoomd_vector::{Cartesian, InnerProduct};
12
13/// A faceted solid defined by the convex hull of a set of points.
14///
15/// [`ConvexPolytope`] stores the given point set without any modification.
16/// Therefore, it can be constructed quickly. The *implicit* convex
17/// hull is formed by [`SupportMapping`] during intersection tests of
18/// `Convex(ConvexPolytope)` with other `Convex(_)` types.
19///
20/// Every vertex in the convex hull is an elements of [`vertices`]. [`vertices`]
21/// may also include duplicate, collinear, coplanar, and/or interior points.
22/// They are exactly the points given at construction.
23///
24/// [`vertices`]: Self::vertices
25///
26/// # Examples
27///
28/// Construction and basic methods:
29/// ```
30/// use approxim::assert_relative_eq;
31/// use hoomd_geometry::{BoundingSphereRadius, shape::ConvexPolyhedron};
32///
33/// # fn main() -> Result<(), hoomd_geometry::Error> {
34/// let tetrahedron = ConvexPolyhedron::with_vertices([
35///     [1.0, 1.0, 1.0].into(),
36///     [1.0, -1.0, -1.0].into(),
37///     [-1.0, 1.0, -1.0].into(),
38///     [-1.0, -1.0, 1.0].into(),
39/// ])?;
40///
41/// let bounding_radius = tetrahedron.bounding_sphere_radius();
42///
43/// assert_relative_eq!(bounding_radius.get(), 3.0_f64.sqrt());
44/// # Ok(())
45/// # }
46/// ```
47///
48/// Intersection tests:
49/// ```
50/// use hoomd_geometry::{Convex, IntersectsAt, shape::ConvexPolygon};
51/// use hoomd_vector::{Angle, Cartesian};
52/// use std::f64::consts::PI;
53///
54/// # fn main() -> Result<(), hoomd_geometry::Error> {
55/// let rectangle = ConvexPolygon::with_vertices([
56///     [-2.0, -1.0].into(),
57///     [2.0, -1.0].into(),
58///     [2.0, 1.0].into(),
59///     [-2.0, 1.0].into(),
60/// ])?;
61/// let rectangle = Convex(rectangle);
62///
63/// assert!(!rectangle.intersects_at(
64///     &rectangle,
65///     &[0.0, 2.1].into(),
66///     &Angle::default()
67/// ));
68/// assert!(rectangle.intersects_at(
69///     &rectangle,
70///     &[0.0, 2.1].into(),
71///     &Angle::from(PI / 2.0)
72/// ));
73/// # Ok(())
74/// # }
75/// ```
76#[derive(Clone, Debug, PartialEq, Serialize, Deserialize)]
77pub struct ConvexPolytope<const N: usize, const MAX_VERTICES: usize = 64> {
78    /// The vertices of the shape.
79    vertices: ArrayVec<Cartesian<N>, MAX_VERTICES>,
80    /// The radius of a bounding sphere of the geometry.
81    bounding_radius: PositiveReal,
82}
83
84/// A faceted convex body in two dimensions.
85///
86/// ```rust
87/// use hoomd_geometry::shape::ConvexPolygon;
88///
89/// # fn main() -> Result<(), hoomd_geometry::Error> {
90/// let hexagon = ConvexPolygon::regular(6);
91/// let square = ConvexPolygon::with_vertices([
92///     [-1.0, -1.0].into(),
93///     [1.0, -1.0].into(),
94///     [1.0, 1.0].into(),
95///     [-1.0, 1.0].into(),
96/// ])?;
97/// # Ok(())
98/// # }
99/// ```
100pub type ConvexPolygon = ConvexPolytope<2, 32>;
101
102/// A faceted convex body in three dimensions.
103///
104/// # Example
105///
106/// ```
107/// use hoomd_geometry::shape::{ConvexPolyhedron, Simplex3};
108/// # fn main() -> Result<(), hoomd_geometry::Error> {
109/// let poly = ConvexPolyhedron::with_vertices([
110///     [1.0, 1.0, 1.0].into(),
111///     [1.0, -1.0, -1.0].into(),
112///     [-1.0, 1.0, -1.0].into(),
113///     [-1.0, -1.0, 1.0].into(),
114/// ])?;
115///
116/// assert_eq!(poly.vertices(), Simplex3::default().vertices());
117/// # Ok(())
118/// # }
119/// ```
120pub type ConvexPolyhedron = ConvexPolytope<3, 32>;
121
122impl<const MAX_VERTICES: usize> ConvexPolytope<2, MAX_VERTICES> {
123    /// Create a regular *n*-gon with *n* vertices and circumradius 0.5.
124    ///
125    /// # Example
126    /// ```
127    /// use hoomd_geometry::shape::ConvexPolygon;
128    ///
129    /// let equilateral_triangle = ConvexPolygon::regular(3);
130    /// ```
131    #[inline]
132    #[must_use]
133    #[expect(
134        clippy::missing_panics_doc,
135        reason = "panic will never occur on a hard-coded constant"
136    )]
137    pub fn regular(n: usize) -> ConvexPolytope<2, MAX_VERTICES> {
138        ConvexPolytope {
139            vertices: (0..n)
140                .map(|x| {
141                    let theta = 2.0 * std::f64::consts::PI * (x as f64) / (n as f64);
142                    Cartesian::from([0.5 * f64::cos(theta), 0.5 * f64::sin(theta)])
143                })
144                .collect(),
145            bounding_radius: 0.5
146                .try_into()
147                .expect("hard-coded constant should be positive"),
148        }
149    }
150}
151
152impl<const N: usize, const MAX_VERTICES: usize> ConvexPolytope<N, MAX_VERTICES> {
153    /// Create an `N`-polytope with the given vertices.
154    ///
155    /// # Example
156    /// ```
157    /// use hoomd_geometry::shape::ConvexPolytope;
158    ///
159    /// # fn main() -> Result<(), hoomd_geometry::Error> {
160    /// let equilateral_triangle = ConvexPolytope::<2>::with_vertices([
161    ///     [1.0, 0.0].into(),
162    ///     [0.5, f64::sqrt(3.0) / 2.0].into(),
163    ///     [-0.5, f64::sqrt(3.0) / 2.0].into(),
164    /// ])?;
165    /// # Ok(())
166    /// # }
167    /// ```
168    /// # Errors
169    ///
170    /// [`Error::DegeneratePolytope`] when no vertices are provided.
171    #[inline]
172    pub fn with_vertices<I>(vertices: I) -> Result<ConvexPolytope<N, MAX_VERTICES>, Error>
173    where
174        I: IntoIterator<Item = Cartesian<N>>,
175    {
176        let mut array_vec: ArrayVec<Cartesian<N>, MAX_VERTICES> = ArrayVec::new();
177        for v in vertices {
178            array_vec.try_push(v).map_err(|_| Error::TooManyVertices)?;
179        }
180
181        if array_vec.is_empty() {
182            return Err(Error::DegeneratePolytope);
183        }
184
185        Ok(ConvexPolytope {
186            bounding_radius: Self::bounding_radius(&array_vec),
187            vertices: array_vec,
188        })
189    }
190
191    /// The vertices of the shape.
192    #[inline]
193    #[must_use]
194    pub fn vertices(&self) -> &[Cartesian<N>] {
195        &self.vertices
196    }
197
198    /// Compute the bounding radius.
199    pub(crate) fn bounding_radius(vertices: &[Cartesian<N>]) -> PositiveReal {
200        vertices
201            .iter()
202            .map(Cartesian::norm_squared)
203            .fold(0.0, f64::max)
204            .sqrt()
205            .try_into()
206            .expect("convex polytope should have a positive bounding radius")
207    }
208}
209
210/// Compute the matrix-vector multiplication of an `ArrayVec` against a `Cartesian<N>`.
211///
212/// This returns an `ExactSizeIterator` of f64 values with `lhs.len()` elements.
213#[inline(always)]
214fn matrix_vector_multiply<const MAX_VERTICES: usize, const N: usize>(
215    lhs: &ArrayVec<Cartesian<N>, MAX_VERTICES>,
216    rhs: Cartesian<N>, // Copy appears to be compiled out, and this lets us elide '_
217) -> impl ExactSizeIterator<Item = f64> + '_ {
218    lhs.iter()
219        .map(move |vertex| (0..N).map(|m| vertex[m] * rhs[m]).sum())
220}
221
222impl<const N: usize, const MAX_VERTICES: usize> SupportMapping<Cartesian<N>>
223    for ConvexPolytope<N, MAX_VERTICES>
224{
225    #[inline]
226    fn support_mapping(&self, n: &Cartesian<N>) -> Cartesian<N> {
227        match N {
228            0 => Cartesian::<N>::default(),
229            1 => self.vertices[0],
230            _ => {
231                let scalars = matrix_vector_multiply(&self.vertices, *n);
232
233                let (mut argmax, mut max_val) = (0, f64::NEG_INFINITY);
234                scalars.enumerate().for_each(|(i, x)| {
235                    if x > max_val {
236                        argmax = i;
237                        max_val = x;
238                    }
239                });
240                self.vertices[argmax]
241            }
242        }
243    }
244}
245
246impl<const N: usize, const MAX_VERTICES: usize> BoundingSphereRadius
247    for ConvexPolytope<N, MAX_VERTICES>
248{
249    #[inline]
250    fn bounding_sphere_radius(&self) -> PositiveReal {
251        self.bounding_radius
252    }
253}
254
255#[cfg(test)]
256mod tests {
257    use super::*;
258    use crate::{Convex, IntersectsAt};
259    use hoomd_vector::{Angle, Cartesian, Rotate, Rotation, Versor};
260
261    use approxim::assert_relative_eq;
262    use rstest::*;
263    use std::f64::consts::{FRAC_1_SQRT_2, PI};
264
265    #[fixture]
266    fn simplex3() -> ConvexPolyhedron {
267        ConvexPolyhedron::with_vertices([
268            [1.0, 1.0, 1.0].into(),
269            [1.0, -1.0, -1.0].into(),
270            [-1.0, 1.0, -1.0].into(),
271            [-1.0, -1.0, 1.0].into(),
272        ])
273        .unwrap()
274    }
275
276    #[fixture]
277    fn equilateral_triangle() -> ConvexPolygon {
278        ConvexPolytope::with_vertices([
279            [1.0, 0.0].into(),
280            [0.5, f64::sqrt(3.0) / 2.0].into(),
281            [-0.5, f64::sqrt(3.0) / 2.0].into(),
282        ])
283        .unwrap()
284    }
285
286    #[rstest]
287    fn test_bounding_radius_computed(
288        simplex3: ConvexPolyhedron,
289        equilateral_triangle: ConvexPolygon,
290    ) {
291        assert_eq!(simplex3.bounding_radius.get(), f64::sqrt(3.0));
292        assert_eq!(equilateral_triangle.bounding_radius.get(), f64::sqrt(1.0));
293    }
294
295    #[rstest]
296    fn test_bounding_radius_regular_polygons(#[values(1, 3, 8, 32)] n: usize) {
297        assert_eq!(ConvexPolygon::regular(n).bounding_radius.get(), 0.5);
298        assert_eq!(
299            ConvexPolytope::<2, 32>::regular(n).bounding_radius.get(),
300            0.5
301        );
302    }
303
304    #[test]
305    fn degenerate_polytope() {
306        let result = ConvexPolytope::<3>::with_vertices([]);
307        assert_eq!(result, Err(Error::DegeneratePolytope));
308    }
309
310    #[test]
311    fn support_mapping_2d() {
312        let cuboid = ConvexPolygon::with_vertices([
313            [-1.0, -2.0].into(),
314            [1.0, -2.0].into(),
315            [1.0, 2.0].into(),
316            [-1.0, 2.0].into(),
317        ])
318        .expect("hard-coded vertices form a polygon");
319
320        assert_relative_eq!(
321            cuboid.support_mapping(&Cartesian::from([1.0, 0.1])),
322            [1.0, 2.0].into()
323        );
324        assert_relative_eq!(
325            cuboid.support_mapping(&Cartesian::from([1.0, -0.1])),
326            [1.0, -2.0].into()
327        );
328        assert_relative_eq!(
329            cuboid.support_mapping(&Cartesian::from([-0.1, 1.0])),
330            [-1.0, 2.0].into()
331        );
332        assert_relative_eq!(
333            cuboid.support_mapping(&Cartesian::from([-0.1, -1.0])),
334            [-1.0, -2.0].into()
335        );
336    }
337
338    #[test]
339    fn support_mapping_3d() {
340        let cuboid = ConvexPolyhedron::with_vertices([
341            [-1.0, -2.0, 3.0].into(),
342            [1.0, -2.0, 3.0].into(),
343            [1.0, 2.0, 3.0].into(),
344            [-1.0, 2.0, 3.0].into(),
345            [-1.0, -2.0, -3.0].into(),
346            [1.0, -2.0, -3.0].into(),
347            [1.0, 2.0, -3.0].into(),
348            [-1.0, 2.0, -3.0].into(),
349        ])
350        .expect("hard-coded vertices form a polygon");
351
352        assert_relative_eq!(
353            cuboid.support_mapping(&Cartesian::from([1.0, 0.1, 0.1])),
354            [1.0, 2.0, 3.0].into()
355        );
356        assert_relative_eq!(
357            cuboid.support_mapping(&Cartesian::from([1.0, 0.1, -0.1])),
358            [1.0, 2.0, -3.0].into()
359        );
360        assert_relative_eq!(
361            cuboid.support_mapping(&Cartesian::from([1.0, -0.1, 0.1])),
362            [1.0, -2.0, 3.0].into()
363        );
364        assert_relative_eq!(
365            cuboid.support_mapping(&Cartesian::from([1.0, -0.1, -0.1])),
366            [1.0, -2.0, -3.0].into()
367        );
368        assert_relative_eq!(
369            cuboid.support_mapping(&Cartesian::from([-1.0, 0.1, 0.1])),
370            [-1.0, 2.0, 3.0].into()
371        );
372        assert_relative_eq!(
373            cuboid.support_mapping(&Cartesian::from([-1.0, 0.1, -0.1])),
374            [-1.0, 2.0, -3.0].into()
375        );
376        assert_relative_eq!(
377            cuboid.support_mapping(&Cartesian::from([-1.0, -0.1, 0.1])),
378            [-1.0, -2.0, 3.0].into()
379        );
380        assert_relative_eq!(
381            cuboid.support_mapping(&Cartesian::from([-1.0, -0.1, -0.1])),
382            [-1.0, -2.0, -3.0].into()
383        );
384    }
385
386    // ConvexPolygon tests from hoomd-blue's test_convex_polygon.cc
387
388    #[fixture]
389    fn square() -> Convex<ConvexPolygon> {
390        Convex(
391            ConvexPolygon::with_vertices([
392                [-0.5, -0.5].into(),
393                [0.5, -0.5].into(),
394                [0.5, 0.5].into(),
395                [-0.5, 0.5].into(),
396            ])
397            .expect("hard-coded vertices form a valid polygon"),
398        )
399    }
400
401    #[fixture]
402    fn triangle() -> Convex<ConvexPolygon> {
403        Convex(
404            ConvexPolygon::with_vertices([
405                [-0.5, -0.5].into(),
406                [0.5, -0.5].into(),
407                [0.5, 0.5].into(),
408            ])
409            .expect("hard-coded vertices form a valid polygon"),
410        )
411    }
412
413    #[rstest]
414    fn square_no_rot(square: Convex<ConvexPolygon>) {
415        let a = Angle::identity();
416        assert!(!square.intersects_at(&square, &[10.0, 0.0].into(), &a));
417        assert!(!square.intersects_at(&square, &[-10.0, 0.0].into(), &a));
418
419        assert!(!square.intersects_at(&square, &[1.1, 0.0].into(), &a));
420        assert!(!square.intersects_at(&square, &[-1.1, 0.0].into(), &a));
421        assert!(!square.intersects_at(&square, &[0.0, 1.1].into(), &a));
422        assert!(!square.intersects_at(&square, &[0.0, -1.1].into(), &a));
423
424        assert!(square.intersects_at(&square, &[0.9, 0.2].into(), &a));
425        assert!(square.intersects_at(&square, &[-0.9, 0.2].into(), &a));
426        assert!(square.intersects_at(&square, &[-0.2, 0.9].into(), &a));
427        assert!(square.intersects_at(&square, &[-0.2, -0.9].into(), &a));
428
429        assert!(square.intersects_at(&square, &[1.0, 0.2].into(), &a));
430    }
431
432    #[rstest]
433    fn square_rot(square: Convex<ConvexPolygon>) {
434        let a = Angle::from(PI / 4.0);
435
436        assert!(!square.intersects_at(&square, &[10.0, 0.0].into(), &a));
437        assert!(!square.intersects_at(&square, &[-10.0, 0.0].into(), &a));
438
439        assert!(!square.intersects_at(&square, &[1.3, 0.0].into(), &a));
440        assert!(!square.intersects_at(&square, &[-1.3, 0.0].into(), &a));
441        assert!(!square.intersects_at(&square, &[0.0, 1.3].into(), &a));
442        assert!(!square.intersects_at(&square, &[0.0, -1.3].into(), &a));
443
444        assert!(!square.intersects_at(&square, &[1.3, 0.2].into(), &a));
445        assert!(!square.intersects_at(&square, &[-1.3, 0.2].into(), &a));
446        assert!(!square.intersects_at(&square, &[-0.2, 1.3].into(), &a));
447        assert!(!square.intersects_at(&square, &[-0.2, -1.3].into(), &a));
448
449        assert!(square.intersects_at(&square, &[1.2, 0.2].into(), &a));
450        assert!(square.intersects_at(&square, &[-1.2, 0.2].into(), &a));
451        assert!(square.intersects_at(&square, &[-0.2, 1.2].into(), &a));
452        assert!(square.intersects_at(&square, &[-0.2, -1.2].into(), &a));
453    }
454
455    fn test_overlap<A, B, R, const N: usize>(
456        r_ab: Cartesian<N>,
457        a: &A,
458        b: &B,
459        o_a: R,
460        o_b: &R,
461    ) -> bool
462    where
463        R: Rotation + Rotate<Cartesian<N>>,
464        A: IntersectsAt<B, Cartesian<N>, R>,
465    {
466        let r_a_inverted = o_a.inverted();
467        let v_ij = r_a_inverted.rotate(&r_ab);
468        let o_ij = o_b.combine(&r_a_inverted);
469        a.intersects_at(b, &v_ij, &o_ij)
470    }
471
472    fn assert_symmetric_overlap<A, B, R, const N: usize>(
473        r_ab: Cartesian<N>,
474        a: &A,
475        b: &B,
476        r_a: R,
477        r_b: R,
478        expected: bool,
479    ) where
480        R: Rotation + Rotate<Cartesian<N>>,
481        A: IntersectsAt<B, Cartesian<N>, R>,
482        B: IntersectsAt<A, Cartesian<N>, R>,
483    {
484        assert_eq!(test_overlap(r_ab, a, b, r_a, &r_b), expected);
485        assert_eq!(test_overlap(-r_ab, b, a, r_b, &r_a), expected);
486    }
487
488    #[rstest]
489    fn square_triangle(square: Convex<ConvexPolygon>, triangle: Convex<ConvexPolygon>) {
490        let r_square = Angle::from(-PI / 4.0);
491        let r_triangle = Angle::from(PI);
492
493        assert_symmetric_overlap(
494            [10.0, 0.0].into(),
495            &square,
496            &triangle,
497            r_square,
498            r_triangle,
499            false,
500        );
501
502        assert_symmetric_overlap(
503            [1.3, 0.0].into(),
504            &square,
505            &triangle,
506            r_square,
507            r_triangle,
508            false,
509        );
510
511        assert_symmetric_overlap(
512            [-1.3, 0.0].into(),
513            &square,
514            &triangle,
515            r_square,
516            r_triangle,
517            false,
518        );
519
520        assert_symmetric_overlap(
521            [0.0, 1.3].into(),
522            &square,
523            &triangle,
524            r_square,
525            r_triangle,
526            false,
527        );
528
529        assert_symmetric_overlap(
530            [0.0, -1.3].into(),
531            &square,
532            &triangle,
533            r_square,
534            r_triangle,
535            false,
536        );
537
538        assert_symmetric_overlap(
539            [1.2, 0.2].into(),
540            &square,
541            &triangle,
542            r_square,
543            r_triangle,
544            true,
545        );
546
547        assert_symmetric_overlap(
548            [-0.7, -0.2].into(),
549            &square,
550            &triangle,
551            r_square,
552            r_triangle,
553            true,
554        );
555
556        assert_symmetric_overlap(
557            [0.4, 1.1].into(),
558            &square,
559            &triangle,
560            r_square,
561            r_triangle,
562            true,
563        );
564
565        assert_symmetric_overlap(
566            [-0.2, -1.2].into(),
567            &square,
568            &triangle,
569            r_square,
570            r_triangle,
571            true,
572        );
573    }
574
575    #[fixture]
576    fn octahedron() -> Convex<ConvexPolyhedron> {
577        Convex(
578            ConvexPolyhedron::with_vertices([
579                [-0.5, -0.5, 0.0].into(),
580                [0.5, -0.5, 0.0].into(),
581                [0.5, 0.5, 0.0].into(),
582                [-0.5, 0.5, 0.0].into(),
583                [0.0, 0.0, FRAC_1_SQRT_2].into(),
584                [0.0, 0.0, -FRAC_1_SQRT_2].into(),
585            ])
586            .expect("hard-coded vertices form a valid polyhedron"),
587        )
588    }
589
590    #[fixture]
591    fn cube() -> Convex<ConvexPolyhedron> {
592        Convex(
593            ConvexPolyhedron::with_vertices([
594                [-0.5, -0.5, -0.5].into(),
595                [0.5, -0.5, -0.5].into(),
596                [0.5, 0.5, -0.5].into(),
597                [-0.5, 0.5, -0.5].into(),
598                [-0.5, -0.5, 0.5].into(),
599                [0.5, -0.5, 0.5].into(),
600                [0.5, 0.5, 0.5].into(),
601                [-0.5, 0.5, 0.5].into(),
602            ])
603            .expect("hard-coded vertices form a valid polyhedron"),
604        )
605    }
606
607    #[rstest]
608    fn overlap_octahedron_no_rot(octahedron: Convex<ConvexPolyhedron>) {
609        let q = Versor::identity();
610
611        assert_symmetric_overlap([0.0, 0.0, 0.0].into(), &octahedron, &octahedron, q, q, true);
612
613        assert_symmetric_overlap(
614            [10.0, 0.0, 0.0].into(),
615            &octahedron,
616            &octahedron,
617            q,
618            q,
619            false,
620        );
621
622        assert_symmetric_overlap(
623            [1.1, 0.0, 0.0].into(),
624            &octahedron,
625            &octahedron,
626            q,
627            q,
628            false,
629        );
630
631        assert_symmetric_overlap(
632            [0.0, 1.1, 0.0].into(),
633            &octahedron,
634            &octahedron,
635            q,
636            q,
637            false,
638        );
639
640        assert_symmetric_overlap(
641            [1.1, 0.2, 0.0].into(),
642            &octahedron,
643            &octahedron,
644            q,
645            q,
646            false,
647        );
648
649        assert_symmetric_overlap(
650            [-1.1, 0.2, 0.0].into(),
651            &octahedron,
652            &octahedron,
653            q,
654            q,
655            false,
656        );
657
658        assert_symmetric_overlap(
659            [-0.2, 1.1, 0.0].into(),
660            &octahedron,
661            &octahedron,
662            q,
663            q,
664            false,
665        );
666
667        assert_symmetric_overlap(
668            [-0.2, -1.1, 0.0].into(),
669            &octahedron,
670            &octahedron,
671            q,
672            q,
673            false,
674        );
675
676        assert_symmetric_overlap([0.9, 0.2, 0.0].into(), &octahedron, &octahedron, q, q, true);
677
678        assert_symmetric_overlap(
679            [-0.9, 0.2, 0.0].into(),
680            &octahedron,
681            &octahedron,
682            q,
683            q,
684            true,
685        );
686
687        assert_symmetric_overlap(
688            [-0.2, 0.9, 0.0].into(),
689            &octahedron,
690            &octahedron,
691            q,
692            q,
693            true,
694        );
695
696        assert_symmetric_overlap(
697            [-0.2, -0.9, 0.0].into(),
698            &octahedron,
699            &octahedron,
700            q,
701            q,
702            true,
703        );
704
705        assert_symmetric_overlap([1.0, 0.2, 0.0].into(), &octahedron, &octahedron, q, q, true);
706    }
707
708    #[rstest]
709    fn overlap_cube_no_rot(cube: Convex<ConvexPolyhedron>) {
710        let q = Versor::identity();
711
712        assert_symmetric_overlap([0.0, 0.0, 0.0].into(), &cube, &cube, q, q, true);
713        assert_symmetric_overlap([10.0, 0.0, 0.0].into(), &cube, &cube, q, q, false);
714
715        assert_symmetric_overlap([1.1, 0.0, 0.0].into(), &cube, &cube, q, q, false);
716        assert_symmetric_overlap([0.0, 1.1, 0.0].into(), &cube, &cube, q, q, false);
717        assert_symmetric_overlap([0.0, 0.0, 1.1].into(), &cube, &cube, q, q, false);
718        assert_symmetric_overlap([1.1, 0.2, 0.0].into(), &cube, &cube, q, q, false);
719
720        assert_symmetric_overlap([-1.1, 0.2, 0.0].into(), &cube, &cube, q, q, false);
721        assert_symmetric_overlap([-0.2, 1.1, 0.0].into(), &cube, &cube, q, q, false);
722        assert_symmetric_overlap([-0.2, -1.1, 0.0].into(), &cube, &cube, q, q, false);
723
724        assert_symmetric_overlap([0.9, 0.2, 0.0].into(), &cube, &cube, q, q, true);
725        assert_symmetric_overlap([-0.9, 0.2, 0.0].into(), &cube, &cube, q, q, true);
726        assert_symmetric_overlap([-0.2, 0.9, 0.0].into(), &cube, &cube, q, q, true);
727        assert_symmetric_overlap([-0.2, -0.9, 0.0].into(), &cube, &cube, q, q, true);
728
729        assert_symmetric_overlap([0.2, 0.0, 0.0].into(), &cube, &cube, q, q, true);
730        assert_symmetric_overlap([0.2, 0.00001, 0.00001].into(), &cube, &cube, q, q, true);
731        assert_symmetric_overlap([0.1, 0.2, 0.1].into(), &cube, &cube, q, q, true);
732        assert_symmetric_overlap([1.0, 0.2, 0.0].into(), &cube, &cube, q, q, true);
733    }
734
735    #[rstest]
736    fn overlap_cube_rot1(cube: Convex<ConvexPolyhedron>) {
737        let q_a = Versor::identity();
738        let q_b = Versor::from_axis_angle(
739            [0.0, 0.0, 1.0]
740                .try_into()
741                .expect("hard-coded vector is non-zero"),
742            PI / 4.0,
743        );
744
745        assert_symmetric_overlap([10.0, 0.0, 0.0].into(), &cube, &cube, q_a, q_b, false);
746
747        assert_symmetric_overlap([1.3, 0.0, 0.0].into(), &cube, &cube, q_a, q_b, false);
748        assert_symmetric_overlap([-1.3, 0.0, 0.0].into(), &cube, &cube, q_a, q_b, false);
749        assert_symmetric_overlap([0.0, 1.3, 0.0].into(), &cube, &cube, q_a, q_b, false);
750        assert_symmetric_overlap([0.0, -1.3, 0.0].into(), &cube, &cube, q_a, q_b, false);
751
752        assert_symmetric_overlap([1.3, 0.2, 0.0].into(), &cube, &cube, q_a, q_b, false);
753        assert_symmetric_overlap([-1.3, 0.2, 0.0].into(), &cube, &cube, q_a, q_b, false);
754        assert_symmetric_overlap([-0.2, 1.3, 0.0].into(), &cube, &cube, q_a, q_b, false);
755        assert_symmetric_overlap([-0.2, -1.3, 0.0].into(), &cube, &cube, q_a, q_b, false);
756
757        assert_symmetric_overlap([1.2, 0.2, 0.0].into(), &cube, &cube, q_a, q_b, true);
758        assert_symmetric_overlap([-1.2, 0.2, 0.0].into(), &cube, &cube, q_a, q_b, true);
759        assert_symmetric_overlap([-0.2, 1.2, 0.0].into(), &cube, &cube, q_a, q_b, true);
760        assert_symmetric_overlap([-0.2, -1.2, 0.0].into(), &cube, &cube, q_a, q_b, true);
761    }
762
763    #[rstest]
764    fn overlap_cube_rot3(cube: Convex<ConvexPolyhedron>) {
765        let q_a = Versor::identity();
766        let q1 = Versor::from_axis_angle(
767            [1.0, 0.0, 0.0]
768                .try_into()
769                .expect("hard-coded vector is non-zero"),
770            PI / 4.0,
771        );
772        let q2 = Versor::from_axis_angle(
773            [0.0, 0.0, 1.0]
774                .try_into()
775                .expect("hard-coded vector is non-zero"),
776            PI / 4.0,
777        );
778        let q_b = q2.combine(&q1);
779
780        assert_symmetric_overlap([10.0, 0.0, 0.0].into(), &cube, &cube, q_a, q_b, false);
781
782        assert_symmetric_overlap([1.4, 0.0, 0.0].into(), &cube, &cube, q_a, q_b, false);
783        assert_symmetric_overlap([-1.4, 0.0, 0.0].into(), &cube, &cube, q_a, q_b, false);
784        assert_symmetric_overlap([0.0, 1.4, 0.0].into(), &cube, &cube, q_a, q_b, false);
785        assert_symmetric_overlap([0.0, -1.4, 0.0].into(), &cube, &cube, q_a, q_b, false);
786
787        assert_symmetric_overlap([1.4, 0.2, 0.0].into(), &cube, &cube, q_a, q_b, false);
788        assert_symmetric_overlap([-1.4, 0.2, 0.0].into(), &cube, &cube, q_a, q_b, false);
789        assert_symmetric_overlap([-0.2, 1.4, 0.0].into(), &cube, &cube, q_a, q_b, false);
790        assert_symmetric_overlap([-0.2, -1.4, 0.0].into(), &cube, &cube, q_a, q_b, false);
791
792        assert_symmetric_overlap([0.0, 1.2, 0.0].into(), &cube, &cube, q_a, q_b, true);
793        assert_symmetric_overlap([0.0, 1.2, 0.1].into(), &cube, &cube, q_a, q_b, true);
794        assert_symmetric_overlap([0.1, 1.2, 0.1].into(), &cube, &cube, q_a, q_b, true);
795        assert_symmetric_overlap([1.2, 0.0, 0.0].into(), &cube, &cube, q_a, q_b, true);
796        assert_symmetric_overlap([1.2, 0.1, 0.0].into(), &cube, &cube, q_a, q_b, true);
797        assert_symmetric_overlap([1.2, 0.1, 0.1].into(), &cube, &cube, q_a, q_b, true);
798
799        assert_symmetric_overlap([-0.9, 0.9, 0.0].into(), &cube, &cube, q_a, q_b, true);
800        assert_symmetric_overlap([-0.9, 0.899, 0.001].into(), &cube, &cube, q_a, q_b, true);
801        assert_symmetric_overlap([0.9, -0.9, 0.0].into(), &cube, &cube, q_a, q_b, true);
802        assert_symmetric_overlap([-0.9, 0.9, 0.1].into(), &cube, &cube, q_a, q_b, true);
803    }
804}
hoomd_geometry/shape/convex_polytope.rs

hoomd_geometry/shape/
convex_polytope.rs
