math_utils/
traits.rs

1//! Abstract traits
2
3use either::Either;
4#[cfg(feature = "derive_serdes")]
5use serde;
6
7use crate::{approx, num};
8use crate::types::{
9  Affinity, LinearIso, NonZero, NonNegative, Normalized, Projectivity, Sign
10};
11
12/// Projective completion (homogeneous coordinates)
13pub trait ProjectiveSpace <S : Field> : AffineSpace <S> {
14  /// Affine subspace
15  type Patch : AffineSpace <S>;
16  /// Construct an augmented matrix for the affinity
17  fn homography (affinity :
18    Affinity <S, Self::Patch, Self::Patch,
19      <<Self::Patch as AffineSpace <S>>::Translation as Module <S>>::LinearEndo>
20  ) -> Projectivity <S, Self, Self, <Self::Translation as Module <S>>::LinearEndo> where
21    <Self::Translation as Module <S>>::LinearEndo :
22      From <<<Self::Patch as AffineSpace <S>>::Translation as Module <S>>::LinearEndo>
23  {
24    Projectivity::new (LinearIso::new ((*affinity.linear_iso).into()).unwrap())
25  }
26  /// Return the projective completion (homogeneous coordinates) of the affine point or
27  /// vector
28  fn homogeneous (point_or_vector :
29    Either <Self::Patch, <Self::Patch as AffineSpace <S>>::Translation>
30  ) -> Self where Self : From <(Self::Patch, S)> {
31    match point_or_vector {
32      Either::Left  (point)  => (point, S::one()).into(),
33      Either::Right (vector) => (Self::Patch::from_vector (vector), S::zero()).into()
34    }
35  }
36}
37
38/// `AffineSpace` with translations in a (Euclidean) real inner product space
39pub trait EuclideanSpace <S : Real> : AffineSpace <S> + MetricSpace <S> { }
40
41/// Space of `Point`s (positions) and `Vector`s (displacements)
42pub trait AffineSpace <S : Field> : Point <Self::Translation> {
43  type Translation : VectorSpace <S> + GroupAction <Self>;
44}
45
46/// Point types convertible to and from a vector type, with difference function that
47/// follows from the free and transitive group action
48pub trait Point <V> : Sized + std::ops::Sub <Self, Output=V> where
49  V : AdditiveGroup + GroupAction <Self>
50{
51  fn to_vector (self) -> V;
52  fn from_vector (vector : V) -> Self;
53  fn origin() -> Self {
54    Self::from_vector (V::zero())
55  }
56}
57
58/// (Right) action of a group on a set
59pub trait GroupAction <X> : Group {
60  fn action (self, x : X) -> X;
61}
62
63impl <G> GroupAction <G> for G where G : Group {
64  #[expect(clippy::renamed_function_params)]
65  fn action (self, g : Self) -> Self {
66    Self::operation (g, self)
67  }
68}
69
70/// Set of points with distance function
71pub trait MetricSpace <S : Field> : NormedVectorSpace <S> {
72  fn distance_squared (self, other : Self) -> NonNegative <S> {
73    (self - other).norm_squared()
74  }
75  fn distance (self, other : Self) -> NonNegative <S> where S : Sqrt {
76    self.distance_squared (other).sqrt()
77  }
78}
79
80/// `VectorSpace` with vector length/magnitude function
81pub trait NormedVectorSpace <S : Field> : InnerProductSpace <S> {
82  type Unit : Into <Self>;
83  // required
84  fn norm_squared (self) -> NonNegative <S>;
85  /// Infinity norm or uniform norm
86  fn norm_max (self) -> NonNegative <S> where S : SignedExt;
87  #[must_use]
88  fn normalize (self) -> Self::Unit where S : Sqrt;
89  // provided
90  fn norm (self) -> NonNegative <S> where S : Sqrt {
91    self.norm_squared().sqrt()
92  }
93  /// Returns zero vector if input is the zero vector, otherwise a normalized vector of
94  /// each axis sign.
95  fn unit_sigvec (self) -> Self where S : SignedExt + Sqrt {
96    let v = self.sigvec();
97    if v.is_zero() {
98      v
99    } else {
100      v.normalize().into()
101    }
102  }
103}
104
105/// Bilinear form on a `VectorSpace`
106pub trait InnerProductSpace <S : Field> : VectorSpace <S> + Dot <S> {
107  fn inner_product (self, other : Self) -> S {
108    self.dot (other)
109  }
110  fn outer_product (self, other : Self) -> Self::LinearEndo;
111  fn orthogonal (self) -> Self where S : approx::AbsDiffEq <Epsilon = S>;
112}
113
114/// Module with scalars taken from a `Field`
115pub trait VectorSpace <S : Field> : Module <S> + std::ops::Div <S, Output=Self> + Copy {
116  type NonZero;
117  // required
118  #[must_use]
119  fn map <F> (self, f : F) -> Self where F : FnMut (S) -> S;
120  // provided
121  /// Map `signum_or_zero` over each element of the given vector
122  fn sigvec (self) -> Self where S : SignedExt {
123    self.map (SignedExt::signum_or_zero)
124  }
125  /// Linear interpolation
126  fn interpolate (a : Self, b : Self, t : Normalized <S>) -> Self {
127    a + (b - a) * *t
128  }
129  /// Linear extrapolation
130  fn extrapolate (a : Self, b : Self, t : S) -> Self {
131    a + (b - a) * t
132  }
133}
134
135/// Additional `Vector2` methods
136pub trait Vector2Ext <S> {
137  /// Exterior or wedge product: $(a,b) \wedge (c,d) = ad - bc$
138  fn exterior_product (self, rhs : Self) -> S;
139}
140
141/// Scalar product (bilinear form) on a `Module`
142pub trait Dot <S : Ring> : Module <S> {
143  fn dot (self, other : Self) -> S;
144  /// Self dot product
145  fn self_dot (self) -> NonNegative <S> where S : OrderedRing {
146    NonNegative::unchecked (self.dot (self))
147  }
148}
149
150/// Additive group combined with scalar multiplication
151pub trait Module <S : Ring> : AdditiveGroup + std::ops::Mul <S, Output=Self> + Copy {
152  /// Linear endomorphism represented by a square matrix type
153  type LinearEndo : LinearMap <S, Self, Self> + Ring;
154}
155
156/// Module homomorphism
157pub trait LinearMap <S, V, W> : std::ops::Mul <V, Output=W> + Copy where
158  S : Ring,
159  V : Module <S>,
160  W : Module <S>
161{
162  fn determinant (self) -> S;
163  fn transpose   (self) -> Self;
164}
165
166/// Additional matrix methods
167pub trait Matrix <S> : Copy {
168  type Rows;
169  type Submatrix;
170  /// Returns the row-major matrix
171  fn rows (self) -> Self::Rows;
172  /// Returns the submatrix formed by removing the given (col,row)
173  fn submatrix (self, i : usize, j : usize) -> Self::Submatrix;
174  /// Returns the matrix formed by filling additional row and column with zeros
175  fn fill_zeros (submatrix : Self::Submatrix) -> Self where S : num::Zero;
176  /// Computes the determinant of the (i,j)-submatrix
177  #[inline]
178  fn minor <V, W> (self, i : usize, j : usize) -> S where
179    S : Ring,
180    V : Module <S>,
181    W : Module <S>,
182    Self::Submatrix : LinearMap <S, V, W>
183  {
184    self.submatrix (i, j).determinant()
185  }
186}
187
188/// Householder transformation
189pub trait ElementaryReflector {
190  type Vector;
191  /// Construct a Householder transformation
192  fn elementary_reflector (v : Self::Vector, index : usize) -> Self;
193}
194
195/// `OrderedField` with special functions
196pub trait Real : OrderedField + Exp + Powf + Sqrt + Trig {
197  // TODO: more constants
198  fn pi() -> Self;
199  fn frac_pi_3() -> Self;
200  fn sqrt_3() -> Self;
201  fn frac_1_sqrt_3() -> Self;
202}
203
204/// A `Field` with an ordering
205pub trait OrderedField : Field + OrderedRing { }
206
207/// A (commutative) `Ring` where $1 \neq 0$ and all non-zero elements are invertible
208pub trait Field : Ring + MultiplicativeGroup + Powi + Rational {
209  fn half() -> Self {
210    Self::one() / Self::two()
211  }
212}
213
214/// Some additional methods for fractional parts of numbers
215pub trait Rational {
216  fn floor (self) -> Self;
217  fn ceil (self) -> Self;
218  fn trunc (self) -> Self;
219  fn fract (self) -> Self;
220}
221
222/// Interface for a group with identity represented by `one`, operation defined by `*`
223/// and `/`
224pub trait MultiplicativeGroup : MultiplicativeMonoid +
225  std::ops::Div <Self, Output=Self> + std::ops::DivAssign + num::Inv <Output=Self>
226{ }
227
228/// Ring of integers
229pub trait Integer : OrderedRing + SignedExt + num::PrimInt { }
230
231/// A Ring with an ordering
232pub trait OrderedRing : Ring + SignedExt + MinMax + PartialOrd + std::ops::Rem { }
233
234/// Interface for a ring as a combination of an additive group and a distributive
235/// multiplication operation.
236///
237/// This is basically the base "scalar" trait (a [`Module`] is defined over a `Ring`).
238/// Here we also add the `Debug` trait as a constraint because basically all scalars
239/// should implement `Debug`, and this avoids us needing to explicitly give a `Debug`
240/// constraint when using functions or macros like `assert_eq` that require it.
241pub trait Ring : AdditiveGroup + MultiplicativeMonoid
242  + num::MulAdd <Self, Self, Output=Self> + num::MulAddAssign <Self, Self>
243  + std::fmt::Debug
244{
245  fn two() -> Self {
246    Self::one() + Self::one()
247  }
248  fn three() -> Self {
249    Self::one() + Self::one() + Self::one()
250  }
251  fn four() -> Self {
252    Self::two() * Self::two()
253  }
254  fn five() -> Self {
255    Self::two() + Self::three()
256  }
257  fn six() -> Self {
258    Self::two() * Self::three()
259  }
260  fn seven() -> Self {
261    Self::three() + Self::four()
262  }
263  fn eight() -> Self {
264    Self::two() * Self::two() * Self::two()
265  }
266  fn nine() -> Self {
267    Self::three() * Self::three()
268  }
269  fn ten() -> Self {
270    Self::two() * Self::five()
271  }
272}
273
274/// Interface for a group with identity represented by `zero`, and operation defined by
275/// `+` and `-`
276pub trait AdditiveGroup : AdditiveMonoid +
277  std::ops::Sub <Self, Output=Self> + std::ops::SubAssign + std::ops::Neg <Output=Self>
278{ }
279
280/// Set with identity, inverses, and (associative) binary operation
281pub trait Group : PartialEq + std::ops::Neg <Output=Self> {
282  fn identity() -> Self;
283  fn operation (a : Self, b : Self) -> Self;
284}
285
286/// Set with identity represented by `one` and (associative) binary operation defined by
287/// `*`
288pub trait MultiplicativeMonoid : Copy + PartialEq +
289  std::ops::Mul <Self, Output=Self> + std::ops::MulAssign + num::One
290{
291  fn squared (self) -> Self {
292    self * self
293  }
294  fn cubed (self) -> Self {
295    self * self * self
296  }
297}
298
299/// Set with identity represented by `zero` and (associative) binary operation defined
300/// by `+`
301pub trait AdditiveMonoid : Sized + PartialEq + std::iter::Sum
302  + std::ops::Add <Self, Output=Self> + std::ops::AddAssign + num::Zero
303{ }
304
305/// Interface for angle units
306pub trait Angle <S : OrderedField> : Clone + Copy + PartialEq + PartialOrd + Sized +
307  AdditiveGroup + std::ops::Div <Self, Output=S> + std::ops::Mul <S, Output=Self> +
308  std::ops::Div <S, Output=Self> + std::ops::Rem <Self, Output=Self>
309{
310  // required
311  /// Full rotation
312  fn full_turn() -> Self;
313  // provided
314  /// Half rotation
315  fn half_turn() -> Self {
316    Self::full_turn() / S::two()
317  }
318  /// Restrict to `(-half_turn, half_turn]`
319  fn wrap_signed (self) -> Self {
320    if self > Self::half_turn() || self <= -Self::half_turn() {
321      let out = (self + Self::half_turn()).wrap_unsigned() - Self::half_turn();
322      if out == -Self::half_turn() {
323        Self::half_turn()
324      } else {
325        out
326      }
327    } else {
328      self
329    }
330  }
331  /// Restrict to `[0, full_turn)`
332  fn wrap_unsigned (self) -> Self {
333    if self >= Self::full_turn() {
334      self % Self::full_turn()
335    } else if self < Self::zero() {
336      self + Self::full_turn() * ((self / Self::full_turn()).trunc().abs() + S::one())
337    } else {
338      self
339    }
340  }
341}
342
343/// Unsigned integer power function
344pub trait Pow {
345  fn pow (self, exp : u32) -> Self;
346}
347/// Signed integer power function
348pub trait Powi {
349  fn powi (self, n : i32) -> Self;
350}
351/// Fractional power function
352pub trait Powf {
353  fn powf (self, n : Self) -> Self;
354}
355/// Exponential function
356pub trait Exp {
357  fn exp (self) -> Self;
358}
359/// Square root function
360pub trait Sqrt {
361  fn sqrt (self) -> Self;
362}
363/// Cube root function
364pub trait Cbrt {
365  fn cbrt (self) -> Self;
366}
367/// Trigonometric functions
368pub trait Trig : Sized {
369  fn sin      (self) -> Self;
370  fn sin_cos  (self) -> (Self, Self);
371  fn cos      (self) -> Self;
372  fn tan      (self) -> Self;
373  fn asin     (self) -> Self;
374  fn acos     (self) -> Self;
375  fn atan     (self) -> Self;
376  fn atan2    (self, other : Self) -> Self;
377}
378
379/// Provides `min`, `max`, and `clamp` that are not necessarily consistent with those
380/// from `Ord`. This is provided because `f32` and `f64` do not implement `Ord`, so this
381/// trait is defined to give a uniform interface with `Ord` types.
382pub trait MinMax {
383  fn min   (self, other : Self) -> Self;
384  fn max   (self, other : Self) -> Self;
385  fn clamp (self, min : Self, max : Self) -> Self;
386}
387
388/// Function returning number representing sign of self
389pub trait SignedExt : num::Signed {
390  #[inline]
391  fn sign (self) -> Sign {
392    if self.is_zero() {
393      Sign::Zero
394    } else if self.is_positive() {
395      Sign::Positive
396    } else {
397      debug_assert!(self.is_negative());
398      Sign::Negative
399    }
400  }
401  /// Maps `0.0` to `0.0`, otherwise equal to `S::signum` (which would otherwise map
402  /// `+0.0 -> 1.0` and `-0.0 -> -1.0`)
403  #[inline]
404  fn signum_or_zero (self) -> Self where Self : num::Zero {
405    if self.is_zero() {
406      Self::zero()
407    } else {
408      self.signum()
409    }
410  }
411
412  #[inline]
413  fn signum_or_zero_approx (self) -> Self where
414    Self : OrderedRing + approx::AbsDiffEq <Epsilon = Self>
415  {
416    let one = Self::one();
417    if self.abs() < Self::default_epsilon() * (one + one + one + one) {
418      Self::zero()
419    } else {
420      self.signum()
421    }
422  }
423}
424
425/// Adds serde `Serialize` and `DeserializeOwned` constraints.
426///
427/// This makes it easier to conditionally add these constraints to type definitions when
428/// `derive_serdes` feature is enabled.
429#[cfg(not(feature = "derive_serdes"))]
430pub trait MaybeSerDes { }
431#[cfg(feature = "derive_serdes")]
432pub trait MaybeSerDes : serde::Serialize + serde::de::DeserializeOwned { }
433
434impl <S, T> MetricSpace <S> for T where S : Field, T : NormedVectorSpace <S> { }
435impl <T> OrderedField for T where T : Field + OrderedRing { }
436impl <T> Field        for T where
437  T : Ring + MultiplicativeGroup + Powi + Rational { }
438impl <T> Integer      for T where T : OrderedRing + num::PrimInt { }
439impl <T> OrderedRing  for T where
440  T : Ring + SignedExt + MinMax + PartialOrd + std::ops::Rem
441{ }
442impl <T> Ring         for T where
443  T : AdditiveGroup + MultiplicativeMonoid
444    + num::MulAdd<Self, Self, Output=Self> + num::MulAddAssign <Self, Self>
445    + std::fmt::Debug
446{ }
447impl <T> AdditiveGroup for T where
448  T : AdditiveMonoid + std::ops::Sub <Self, Output=Self> + std::ops::SubAssign
449    + std::ops::Neg <Output=Self>
450{ }
451impl <T> MultiplicativeGroup for T where
452  T : MultiplicativeMonoid + std::ops::Div <Self, Output=Self> + std::ops::DivAssign
453    + num::Inv <Output=Self>
454{ }
455impl <T> AdditiveMonoid for T where
456  T : Sized + PartialEq + std::iter::Sum + std::ops::Add <Self, Output=Self>
457    + std::ops::AddAssign + num::Zero
458{ }
459impl <T> MultiplicativeMonoid for T where
460  T : Copy + PartialEq + std::ops::Mul <Self, Output=Self> + std::ops::MulAssign
461    + num::One
462{ }
463impl <T : num::Signed> SignedExt for T { }
464
465impl <
466  #[cfg(not(feature = "derive_serdes"))]
467  T,
468  #[cfg(feature = "derive_serdes")]
469  T : serde::Serialize + serde::de::DeserializeOwned
470> MaybeSerDes for T { }
471
472macro impl_integer ($type:ty) {
473  impl MinMax         for $type {
474    fn min (self, other : Self) -> Self {
475      Ord::min (self, other)
476    }
477    fn max (self, other : Self) -> Self {
478      Ord::max (self, other)
479    }
480    fn clamp (self, min : Self, max : Self) -> Self {
481      Ord::clamp (self, min, max)
482    }
483  }
484  impl Pow for $type {
485    fn pow (self, exp : u32) -> Self {
486      self.pow (exp)
487    }
488  }
489}
490impl_integer!(i8);
491impl_integer!(i16);
492impl_integer!(i32);
493impl_integer!(i64);
494
495macro impl_real_float ($type:ident) {
496  impl AffineSpace <Self> for $type {
497    type Translation = $type;
498  }
499  impl Point <Self> for $type {
500    fn to_vector (self) -> $type {
501      self
502    }
503    fn from_vector (vector : $type) -> Self {
504      vector
505    }
506  }
507  impl VectorSpace <Self> for $type {
508    type NonZero = NonZero <$type>;
509    fn map <F> (self, mut f : F) -> Self where F : FnMut (Self) -> Self {
510      f (self)
511    }
512  }
513  impl Module <Self> for $type {
514    type LinearEndo = Self;
515  }
516  impl Group for $type {
517    fn identity() -> Self {
518      1.0
519    }
520    fn operation (a : Self, b : Self) -> Self {
521      a * b
522    }
523  }
524  impl LinearMap <Self, Self, Self> for $type {
525    fn determinant (self) -> Self {
526      self
527    }
528    fn transpose (self) -> Self {
529      self
530    }
531  }
532  impl Real for $type {
533    fn pi() -> Self {
534      std::$type::consts::PI
535    }
536    fn frac_pi_3() -> Self {
537      std::$type::consts::FRAC_PI_3
538    }
539    #[expect(clippy::excessive_precision)]
540    #[allow(clippy::allow_attributes)]
541    #[allow(clippy::cast_possible_truncation)]
542    #[allow(trivial_numeric_casts)]
543    fn sqrt_3() -> Self {
544      1.732050807568877293527446341505872366942805253810380628055f64 as $type
545    }
546    #[expect(clippy::excessive_precision)]
547    #[allow(clippy::allow_attributes)]
548    #[allow(clippy::cast_possible_truncation)]
549    #[allow(trivial_numeric_casts)]
550    fn frac_1_sqrt_3() -> Self {
551      (1.0f64 / 1.732050807568877293527446341505872366942805253810380628055f64) as $type
552    }
553  }
554  impl Rational for $type {
555    fn floor (self) -> Self {
556      self.floor()
557    }
558    fn ceil (self) -> Self {
559      self.ceil()
560    }
561    fn trunc (self) -> Self {
562      self.trunc()
563    }
564    fn fract (self) -> Self {
565      self.fract()
566    }
567  }
568  impl MinMax for $type {
569    fn min (self, other : Self) -> Self {
570      self.min (other)
571    }
572    fn max (self, other : Self) -> Self {
573      self.max (other)
574    }
575    fn clamp (self, min : Self, max : Self) -> Self {
576      self.clamp (min, max)
577    }
578  }
579  impl Powi for $type {
580    fn powi (self, n : i32) -> Self {
581      self.powi (n)
582    }
583  }
584  impl Powf for $type {
585    fn powf (self, n : Self) -> Self {
586      self.powf (n)
587    }
588  }
589  impl Exp for $type {
590    fn exp (self) -> Self {
591      self.exp()
592    }
593  }
594  impl Sqrt for $type {
595    fn sqrt (self) -> Self {
596      self.sqrt()
597    }
598  }
599  impl Cbrt for $type {
600    fn cbrt (self) -> Self {
601      self.cbrt()
602    }
603  }
604  impl Trig for $type {
605    fn sin (self) -> Self {
606      self.sin()
607    }
608    fn sin_cos (self) -> (Self, Self) {
609      self.sin_cos()
610    }
611    fn cos (self) -> Self {
612      self.cos()
613    }
614    fn tan (self) -> Self {
615      self.tan()
616    }
617    fn asin (self) -> Self {
618      self.asin()
619    }
620    fn acos (self) -> Self {
621      self.acos()
622    }
623    fn atan (self) -> Self {
624      self.atan()
625    }
626    fn atan2 (self, other : Self) -> Self {
627      self.atan2 (other)
628    }
629  }
630}
631
632impl_real_float!(f32);
633impl_real_float!(f64);
math_utils/traits.rs

math_utils/
traits.rs
