nsys-math-utils 0.2.0

//! 3D volumetric forms

use std::convert::TryFrom;
use derive_more::{From, TryFrom};
use num_traits as num;
#[cfg(feature = "derive_serdes")]
use serde::{Deserialize, Serialize};

use crate::*;

use super::{Aabb3, Capsule3, Cylinder3, Sphere3};

////////////////////////////////////////////////////////////////////////////////////////
//  traits                                                                            //
////////////////////////////////////////////////////////////////////////////////////////

/// A trait for bounded and unbounded shapes.
///
/// A 'Shape' implements the following traits:
///
/// - 'Aabb' -- a shape can compute the axis-aligned bounding volume that encloses it
/// - 'Bsphere' -- TODO: a shape can compute its bounding sphere
///
/// Both of these traits carry the constraint that a shape implements the 'Stereometric'
/// trait which allows a shape to compute its volume (non-zero, possibly infinite).
pub trait Shape <S : Ring> : Aabb <S> /*+ Bsphere <S>*/ { } // TODO: bspheres

/// Trait for computing the axis-aligned bounding volume of a given shape.
///
/// Note for 'Unbounded' shapes some of the components will contain positive or negative
/// infinity.
pub trait Aabb <S : Ring> : Stereometric <S> {
  fn aabb (&self) -> Aabb3 <S>;
}

// TODO: implement bounding spheres
/// A trait for computing bounding spheres
pub trait Bsphere <S : Ring> : Stereometric <S> {
  fn sphere (&self) -> Sphere3 <S>;
}

/// Trait for computing volumes of solid figures.
///
/// Note for 'Unbounded' shapes this will be positive infinitiy.
pub trait Stereometric <S : Ring> {
  fn volume (&self) -> Positive <S>;
}

////////////////////////////////////////////////////////////////////////////////////////
//  enums                                                                             //
////////////////////////////////////////////////////////////////////////////////////////

/// A shape that is either bounded (encloses a finite volume) or unbounded (delimits an
/// infinite volume)
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(From, Clone, Debug, Eq, PartialEq, TryFrom)]
pub enum Variant <S> {
  Bounded   (Bounded   <S>),
  Unbounded (Unbounded <S>)
}

/// A (totally) bounded shape
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(From, Clone, Debug, Eq, PartialEq)]
pub enum Bounded <S> {
  Sphere   (Sphere   <S>),
  Capsule  (Capsule  <S>),
  Cylinder (Cylinder <S>),
  Cone     (Cone     <S>),
  Cube     (Cube     <S>),
  Cuboid   (Cuboid   <S>)
}

/// A shape that is only partly bounded
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(From, Clone, Debug, Eq, PartialEq)]
pub enum Unbounded <S> {
  Orthant   (Orthant),
  Halfspace (Halfspace <S>)
}

////////////////////////////////////////////////////////////////////////////////////////
//  structs                                                                           //
////////////////////////////////////////////////////////////////////////////////////////

/// An axis-aligned halfspace
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Orthant {
  pub normal_axis : SignedAxis3
}

/// A halfspace defined by a normal vector
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Halfspace <S> {
  pub normal_vector : Unit3 <S>
}

/// A sphere defined by a positive radius
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Sphere <S> {
  pub radius : Positive <S>
}

/// A capsule defined by strictly positive radius and half-height
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Capsule <S> {
  pub radius      : Positive <S>,
  pub half_height : Positive <S>
}

/// A cylinder defined by strictly positive radius and half-height
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Cylinder <S> {
  pub radius      : Positive <S>,
  pub half_height : Positive <S>
}

/// A cone defined by strictly positive radius and half-height
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Cone <S> {
  pub radius      : Positive <S>,
  pub half_height : Positive <S>
}

/// A cube defined by a strictly positive half-extent
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Cube <S> {
  pub half_extent : Positive <S>
}

/// A box defined by three strictly positive half extents
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Cuboid <S> {
  pub half_extent_x : Positive <S>,
  pub half_extent_y : Positive <S>,
  pub half_extent_z : Positive <S>
}

////////////////////////////////////////////////////////////////////////////////////////
//  functions                                                                         //
////////////////////////////////////////////////////////////////////////////////////////

pub fn report_sizes() {
  use std::mem::size_of;
  println!("shape report sizes...");

  show!(size_of::<Variant   <f32>>());
  show!(size_of::<Variant   <f64>>());

  show!(size_of::<Bounded   <f32>>());
  show!(size_of::<Bounded   <f64>>());

  show!(size_of::<Unbounded <f32>>());
  show!(size_of::<Unbounded <f64>>());

  show!(size_of::<Sphere    <f32>>());
  show!(size_of::<Capsule   <f32>>());
  show!(size_of::<Cylinder  <f32>>());
  show!(size_of::<Cone      <f32>>());
  show!(size_of::<Cube      <f32>>());
  show!(size_of::<Cuboid    <f32>>());

  show!(size_of::<Orthant>());
  show!(size_of::<Halfspace <f32>>());

  show!(size_of::<Sphere    <f64>>());
  show!(size_of::<Capsule   <f64>>());
  show!(size_of::<Cylinder  <f64>>());
  show!(size_of::<Cone      <f64>>());
  show!(size_of::<Cube      <f64>>());
  show!(size_of::<Cuboid    <f64>>());

  show!(size_of::<Orthant>());
  show!(size_of::<Halfspace <f64>>());

  println!("...shape report sizes");
}

////////////////////////////////////////////////////////////////////////////////////////
//  impls                                                                             //
////////////////////////////////////////////////////////////////////////////////////////

//
//  impl Variant
//
impl <S> Shape <S> for Variant <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{ }
impl <S> Aabb <S> for Variant <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  /// Note that unbounded shapes will return a Cuboid with some infinite half extents
  fn aabb (&self) -> Aabb3 <S> {
    match self {
      Variant::Bounded   (bounded)   => bounded.aabb(),
      Variant::Unbounded (unbounded) => unbounded.aabb()
    }
  }
}
impl <S> Stereometric <S> for Variant <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  /// Note that unbounded shapes will return a Cuboid with some infinite half extents
  fn volume (&self) -> Positive <S> {
    match self {
      Variant::Bounded   (bounded)   => bounded.volume(),
      Variant::Unbounded (unbounded) => unbounded.volume()
    }
  }
}

//
//  impl Bounded
//
impl <S> Shape <S> for Bounded <S> where S : Real + std::fmt::Debug { }
impl <S> Aabb <S> for Bounded <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    match self {
      Bounded::Sphere   (sphere)   => sphere.aabb(),
      Bounded::Capsule  (capsule)  => capsule.aabb(),
      Bounded::Cylinder (cylinder) => cylinder.aabb(),
      Bounded::Cone     (cone)     => cone.aabb(),
      Bounded::Cube     (cube)     => cube.aabb(),
      Bounded::Cuboid   (cuboid)   => cuboid.aabb()
    }
  }
}
impl <S : Real> Stereometric <S> for Bounded <S> {
  /// Volume of a bounded solid is always finite
  fn volume (&self) -> Positive <S> {
    match self {
      Bounded::Sphere   (sphere)   => sphere.volume(),
      Bounded::Capsule  (capsule)  => capsule.volume(),
      Bounded::Cylinder (cylinder) => cylinder.volume(),
      Bounded::Cone     (cone)     => cone.volume(),
      Bounded::Cube     (cube)     => cube.volume(),
      Bounded::Cuboid   (cuboid)   => cuboid.volume()
    }
  }
}

//
//  impl Unbounded
//
impl <S> Shape <S> for Unbounded <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{ }
impl <S> Aabb <S> for Unbounded <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  fn aabb (&self) -> Aabb3 <S> {
    match self {
      Unbounded::Orthant   (orthant)   => orthant.aabb(),
      Unbounded::Halfspace (halfspace) => halfspace.aabb()
    }
  }
}
impl <S> Stereometric <S> for Unbounded <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  /// Volume of an unbounded solid is always infinite
  fn volume (&self) -> Positive <S> {
    if cfg!(debug_assertions) {
      let volume = match self {
        Unbounded::Orthant   (orthant)   => orthant.volume(),
        Unbounded::Halfspace (halfspace) => halfspace.volume()
      };
      debug_assert_eq!(*volume, S::infinity());
    }
    Positive::infinity()
  }
}

//
//  impl Sphere
//
impl <S : OrderedRing> Sphere <S> {
  #[inline]
  /// Sphere with radius 1.0
  pub fn unit() -> Self where S : Field {
    use num::One;
    Sphere { radius: Positive::one() }
  }
  /// Create a new sphere with the absolute value of the given radius.
  ///
  /// Panics if radius is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Sphere;
  /// let s = Sphere::noisy (0.0);
  /// ```
  #[inline]
  pub fn noisy (radius : S) -> Self where S : std::fmt::Debug {
    assert_ne!(radius, S::zero());
    Sphere { radius: Positive::unchecked (radius.abs()) }
  }
  /// Create a new sphere with the absolute value of the given radius.
  ///
  /// Debug panic if radius is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Sphere;
  /// let s = Sphere::noisy (0.0);
  /// ```
  #[inline]
  pub fn unchecked (radius : S) -> Self where S : std::fmt::Debug {
    debug_assert_ne!(radius, S::zero());
    Sphere { radius: Positive::unchecked (radius.abs()) }
  }
  #[inline]
  pub const fn sphere3 (&self, center : Point3 <S>) -> Sphere3 <S> {
    Sphere3 { center, radius: self.radius }
  }
}
impl <S : Real> TryFrom <Bounded <S>> for Sphere <S> {
  type Error = Bounded <S>;
  fn try_from (bounded : Bounded <S>) -> Result <Self, Self::Error> {
    match bounded {
      Bounded::Sphere (sphere) => Ok (sphere),
      _ => Err (bounded)
    }
  }
}
impl <S> Shape <S> for Sphere <S> where S : Real + std::fmt::Debug { }
impl <S> Aabb <S> for Sphere <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    Aabb3::with_minmax ([-*self.radius; 3].into(), [ *self.radius; 3].into())
  }
}
impl <S : Real> Stereometric <S> for Sphere <S> {
  fn volume (&self) -> Positive <S> {
    let four      = Positive::unchecked (S::four());
    let frac_pi_3 = Positive::unchecked (S::frac_pi_3());
    four * frac_pi_3 * self.radius * self.radius * self.radius
  }
}

//
//  impl Capsule
//
impl <S : OrderedRing> Capsule <S> {
  /// Create a new capsule with the absolute values of the given radius and half-height.
  ///
  /// Panics if radius is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Capsule;
  /// let s = Capsule::noisy (0.0, 2.0);
  /// ```
  ///
  /// Panics if half-height is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Capsule;
  /// let s = Capsule::noisy (2.0, 0.0);
  /// ```
  #[inline]
  pub fn noisy (radius : S, half_height : S) -> Self where S : std::fmt::Debug {
    assert_ne!(radius, S::zero());
    assert_ne!(half_height, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Capsule { radius, half_height }
  }
  /// Create a new capsule with the absolute values of the given radius and half-height.
  ///
  /// Debug panic if radius is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Capsule;
  /// let s = Capsule::unchecked (0.0, 2.0);
  /// ```
  ///
  /// Debug panic if half-height is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Capsule;
  /// let s = Capsule::noisy (2.0, 0.0);
  /// ```
  #[inline]
  pub fn unchecked (radius : S, half_height : S) -> Self where
    S : std::fmt::Debug
  {
    debug_assert_ne!(radius, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Capsule { radius, half_height }
  }
  /// Height of the cylinder portion
  #[inline]
  pub fn height (&self) -> NonNegative <S> where S : Field {
    self.half_height * NonNegative::unchecked (S::two())
  }
  #[inline]
  pub const fn capsule3 (&self, center : Point3 <S>) -> Capsule3 <S> {
    Capsule3 { center, radius: self.radius, half_height: self.half_height }
  }
}
impl <S> Shape <S> for Capsule <S> where S : Real + std::fmt::Debug { }
impl <S> Aabb <S> for Capsule <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    let r  = *self.radius;
    let hh = *self.half_height;
    Aabb3::with_minmax ([-r, -r, -r - hh].into(), [ r,  r,  r + hh].into())
  }
}
impl <S : Real> Stereometric <S> for Capsule <S> {
  fn volume (&self) -> Positive <S> {
    let r               = self.radius;
    let h               = self.height();
    let r2              = r * r;
    let r3              = r2 * r;
    let pi              = Positive::unchecked (S::pi());
    let four            = Positive::unchecked (S::four());
    let frac_pi_3       = Positive::unchecked (S::frac_pi_3());
    let cylinder_volume = pi * r2 * h;
    let sphere_volume   = four * frac_pi_3 * r3;
    sphere_volume + cylinder_volume
  }
}

//
//  impl Cylinder
//
impl <S : OrderedRing> Cylinder <S> {
  #[inline]
  /// Cylinder with radius 1.0 and half-height 1.0
  pub fn unit() -> Self where S : Field {
    use num::One;
    Cylinder {
      radius: Positive::one(), half_height: Positive::one()
    }
  }
  /// Create a new cylinder with the absolute values of the given radius and
  /// half-height.
  ///
  /// Panics if radius or half-height are zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cylinder;
  /// let s = Cylinder::noisy (0.0, 0.0);
  /// ```
  #[inline]
  pub fn noisy (radius : S, half_height : S) -> Self where S : std::fmt::Debug {
    assert_ne!(radius,      S::zero());
    assert_ne!(half_height, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Cylinder { radius, half_height }
  }
  /// Create a new cylinder with the absolute values of the given radius and
  /// half-height.
  ///
  /// Debug panic if radius or half-height are zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cylinder;
  /// let s = Cylinder::unchecked (0.0, 0.0);
  /// ```
  #[inline]
  pub fn unchecked (radius : S, half_height : S) -> Self where S : std::fmt::Debug {
    debug_assert_ne!(radius,      S::zero());
    debug_assert_ne!(half_height, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Cylinder { radius, half_height }
  }
  #[inline]
  pub const fn cylinder3 (&self, center : Point3 <S>) -> Cylinder3 <S> {
    Cylinder3 { center, radius: self.radius, half_height: self.half_height }
  }
  #[inline]
  pub fn height (&self) -> Positive <S> where S : Field {
    self.half_height * Positive::unchecked (S::two())
  }
}
impl <S> Shape <S> for Cylinder <S> where S : Real + std::fmt::Debug { }
impl <S> Aabb <S> for Cylinder <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    let r  = *self.radius;
    let hh = *self.half_height;
    Aabb3::with_minmax ([-r, -r, -hh].into(), [ r,  r,  hh].into())
  }
}
impl <S : Real> Stereometric <S> for Cylinder <S> {
  fn volume (&self) -> Positive <S> {
    let pi = Positive::unchecked (S::pi());
    pi * self.radius * self.radius * self.height()
  }
}

//
//  impl Cone
//
impl <S : OrderedRing> Cone <S> {
  /// Create a new cone with the absolute values of the given radius and half-height.
  ///
  /// Panics if radius or half-height are zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cone;
  /// let s = Cone::noisy (0.0, 0.0);
  /// ```
  #[inline]
  pub fn noisy (radius : S, half_height : S) -> Self where S : std::fmt::Debug {
    assert_ne!(radius,      S::zero());
    assert_ne!(half_height, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Cone { radius, half_height }
  }
}
impl <S> Shape <S> for Cone <S> where S : Real + std::fmt::Debug { }
impl <S> Aabb <S> for Cone <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    let r  = *self.radius;
    let hh = *self.half_height;
    Aabb3::with_minmax ([-r, -r, -hh].into(), [ r,  r,  hh].into())
  }
}
impl <S : Real> Stereometric <S> for Cone <S> {
  fn volume (&self) -> Positive <S> {
    let frac_pi_3 = Positive::unchecked (S::frac_pi_3());
    let two       = Positive::unchecked (S::two());
    frac_pi_3 * self.radius * self.radius * two * self.half_height
  }
}

//
//  impl Cube
//
impl <S : OrderedRing> Cube <S> {
  /// Create a new cube with the absolute value of the given half-extent.
  ///
  /// Panics if half-extent is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cube;
  /// let s = Cube::noisy (0.0);
  /// ```
  #[inline]
  pub fn noisy (half_extent : S) -> Self where S : std::fmt::Debug {
    assert_ne!(half_extent, S::zero());
    let half_extent = Positive::unchecked (half_extent.abs());
    Cube { half_extent }
  }
  #[inline]
  pub fn extent (&self) -> Positive <S> where S : Field {
    self.half_extent * Positive::unchecked (S::two())
  }
}
impl <S> Shape <S> for Cube <S> where S : OrderedField + std::fmt::Debug { }
impl <S> Aabb <S> for Cube <S> where S : OrderedField + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    Aabb3::with_minmax (
      [-*self.half_extent; 3].into(),
      [ *self.half_extent; 3].into())
  }
}
impl <S : OrderedField> Stereometric <S> for Cube <S> {
  fn volume (&self) -> Positive <S> {
    let extent = self.extent();
    extent * extent * extent
  }
}

//
//  impl Cuboid
//
impl <S : OrderedRing> Cuboid <S> {
  /// Create a new cuboid with the absolute values of the given half-extents.
  ///
  /// Panics if any half-extent is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cuboid;
  /// let s = Cuboid::noisy ([0.0, 0.0, 0.0]);
  /// ```
  #[inline]
  pub fn noisy (half_extents : [S; 3]) -> Self where S : std::fmt::Debug {
    assert_ne!(half_extents[0], S::zero());
    assert_ne!(half_extents[1], S::zero());
    assert_ne!(half_extents[2], S::zero());
    let half_extent_x = Positive::unchecked (half_extents[0].abs());
    let half_extent_y = Positive::unchecked (half_extents[1].abs());
    let half_extent_z = Positive::unchecked (half_extents[2].abs());
    Cuboid { half_extent_x, half_extent_y, half_extent_z }
  }
  #[inline]
  pub fn extents (&self) -> [Positive <S>; 3] where S : Field {
    let two = Positive::unchecked (S::two());
    [ self.half_extent_x * two,
      self.half_extent_y * two,
      self.half_extent_z * two
    ]
  }
  #[inline]
  pub fn half_extents_vec (&self) -> Vector3 <S> {
    [ *self.half_extent_x,
      *self.half_extent_y,
      *self.half_extent_z
    ].into()
  }
  #[inline]
  pub fn max (&self) -> Point3 <S> {
    Vector3 {
      x: *self.half_extent_x,
      y: *self.half_extent_y,
      z: *self.half_extent_z
    }.into()
  }
  #[inline]
  pub fn min (&self) -> Point3 <S> {
    Vector3 {
      x: -*self.half_extent_x,
      y: -*self.half_extent_y,
      z: -*self.half_extent_z
    }.into()
  }
  #[inline]
  pub fn aabb3 (&self, center : Point3 <S>) -> Aabb3 <S> where S : std::fmt::Debug {
    Aabb3::with_minmax (
      center + self.min().0,
      center + self.max().0
    )
  }
}
impl <S> From <Aabb3 <S>> for Cuboid <S> where S : OrderedField {
  fn from (aabb : Aabb3 <S>) -> Self {
    // TODO: aabbs do not enforce zero width/depth/height so this is unsafe
    Cuboid {
      half_extent_x: Positive::unchecked (*aabb.width()  / S::two()),
      half_extent_y: Positive::unchecked (*aabb.depth()  / S::two()),
      half_extent_z: Positive::unchecked (*aabb.height() / S::two())
    }
  }
}
impl <S> Shape <S> for Cuboid <S> where S : OrderedField + std::fmt::Debug { }
impl <S> Aabb <S>  for Cuboid <S> where S : OrderedField + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    Aabb3::with_minmax (self.min(), self.max())
  }
}
impl <S : OrderedField> Stereometric <S> for Cuboid <S> {
  fn volume (&self) -> Positive <S> {
    let [x, y, z] = self.extents();
    x * y * z
  }
}

//
//  impl Orthant
//
impl Orthant {
  pub fn try_from <S : Real> (halfspace : &Halfspace <S>) -> Option <Self> {
    SignedAxis3::try_from (&halfspace.normal_vector)
      .map (|normal_axis| Orthant { normal_axis })
  }
}
impl From <SignedAxis3> for Orthant {
  fn from (normal_axis : SignedAxis3) -> Self {
    Orthant { normal_axis }
  }
}
impl <S> Shape <S> for Orthant where
  S : OrderedField + num::float::FloatCore + std::fmt::Debug
{ }
impl <S> Aabb <S> for Orthant where
  S : OrderedField + num::float::FloatCore + std::fmt::Debug
{
  /// Computes the containing Aabb3.
  ///
  /// ```
  /// # use math_utils::geometry::shape::Orthant;
  /// # use math_utils::geometry::Aabb3;
  /// # use math_utils::SignedAxis3;
  /// use math_utils::geometry::shape::Aabb;
  /// let s = Orthant::from (SignedAxis3::PosZ);
  /// assert_eq!(s.aabb(), Aabb3::with_minmax (
  ///   [f32::NEG_INFINITY, f32::NEG_INFINITY, f32::NEG_INFINITY].into(),
  ///   [f32::INFINITY,     f32::INFINITY,     0.0].into()
  /// ));
  /// let s = Orthant::from (SignedAxis3::NegZ);
  /// assert_eq!(s.aabb(), Aabb3::with_minmax (
  ///   [f32::NEG_INFINITY, f32::NEG_INFINITY, 0.0].into(),
  ///   [f32::INFINITY,     f32::INFINITY,     f32::INFINITY].into()
  /// ));
  /// ```
  fn aabb (&self) -> Aabb3 <S> {
    use num::float::FloatCore;
    let axis_vec = self.normal_axis.to_vec::<S>();
    let surface  = Vector3::broadcast (S::one()) - axis_vec.map (FloatCore::abs);

    let do_min = |i| if surface[i] == S::one() || axis_vec[i] == S::one() {
      S::neg_infinity()
    } else {
      debug_assert_eq!(axis_vec[i], -S::one());
      S::zero()
    };
    let do_max = |i| if surface[i] == S::one() {
      S::infinity()
    } else if axis_vec[i] == S::one() {
      S::zero()
    } else {
      debug_assert_eq!(axis_vec[i], -S::one());
      S::infinity()
    };
    let min = [ do_min (0), do_min (1), do_min (2) ].into();
    let max = [ do_max (0), do_max (1), do_max (2) ].into();

    Aabb3::with_minmax (min, max)
  }
}
impl <S> Stereometric <S> for Orthant where S : OrderedRing + num::float::FloatCore {
  /// Volume of an unbounded solid is always infinite
  fn volume (&self) -> Positive <S> {
    Positive::infinity()
  }
}

//
//  impl Halfspace
//
impl <S : Real> From <Unit3 <S>> for Halfspace <S> {
  fn from (normal_vector : Unit3 <S>) -> Self {
    Halfspace { normal_vector }
  }
}
impl <S> Shape <S> for Halfspace <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{ }
impl <S> Aabb <S> for Halfspace <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  fn aabb (&self) -> Aabb3 <S> {
    if let Some (orthant) = Orthant::try_from (self) {
      orthant.aabb()
    } else {
      Aabb3::with_minmax ([S::neg_infinity(); 3].into(), [S::infinity(); 3].into())
    }
  }
}
impl <S> Stereometric <S> for Halfspace <S> where S : Real + num::float::FloatCore {
  /// Volume of an unbounded solid is always infinite
  fn volume (&self) -> Positive <S> {
    Positive::infinity()
  }
}

////////////////////////////////////////////////////////////////////////////////////////
//  numcast impls                                                                     //
////////////////////////////////////////////////////////////////////////////////////////

impl <S> Halfspace <S> {
  #[inline]
  pub fn numcast <T> (self) -> Option <Halfspace <T>> where
    S : num::NumCast,
    T : num::NumCast
  {
    Some (Halfspace {
      normal_vector: self.normal_vector.numcast()?
    })
  }
}

impl <S> Sphere <S> {
  #[inline]
  pub fn numcast <T> (self) -> Option <Sphere <T>> where
    S : num::NumCast,
    T : num::NumCast
  {
    Some (Sphere {
      radius: self.radius.numcast()?
    })
  }
}

impl <S> Capsule <S> {
  #[inline]
  pub fn numcast <T> (self) -> Option <Capsule <T>> where
    S : num::NumCast,
    T : num::NumCast
  {
    Some (Capsule {
      radius:      self.radius.numcast()?,
      half_height: self.half_height.numcast()?
    })
  }
}

impl <S> Cylinder <S> {
  #[inline]
  pub fn numcast <T> (self) -> Option <Cylinder <T>> where
    S : num::NumCast,
    T : num::NumCast
  {
    Some (Cylinder {
      radius:      self.radius.numcast()?,
      half_height: self.half_height.numcast()?
    })
  }
}

impl <S> Cone <S> {
  #[inline]
  pub fn numcast <T> (self) -> Option <Cone <T>> where
    S : num::NumCast,
    T : num::NumCast
  {
    Some (Cone {
      radius:      self.radius.numcast()?,
      half_height: self.half_height.numcast()?
    })
  }
}

impl <S> Cube <S> {
  #[inline]
  pub fn numcast <T> (self) -> Option <Cube <T>> where
    S : num::NumCast,
    T : num::NumCast
  {
    Some (Cube {
      half_extent: self.half_extent.numcast()?
    })
  }
}

impl <S> Cuboid <S> {
  #[inline]
  pub fn numcast <T> (self) -> Option <Cuboid <T>> where
    S : num::NumCast,
    T : num::NumCast
  {
    Some (Cuboid {
      half_extent_x: self.half_extent_x.numcast()?,
      half_extent_y: self.half_extent_y.numcast()?,
      half_extent_z: self.half_extent_z.numcast()?
    })
  }
}