Skip to main content

cliffy_core/
transforms.rs

1//! Geometric transformations for state manipulation
2//!
3//! This module provides explicit geometric transformation types:
4//! - `Rotor`: Represents rotations in geometric algebra
5//! - `Versor`: General geometric transformations (rotations + reflections)
6//! - `Motor`: Translations combined with rotations (for CGA)
7//!
8//! These allow users to apply explicit geometric operations to state,
9//! rather than hiding all geometry behind scalar updates.
10
11use crate::geometric::GA3;
12use amari_core::{Bivector, Vector};
13
14/// A rotor represents a rotation in geometric algebra.
15///
16/// Rotors are even-grade multivectors with unit magnitude that
17/// implement rotations via the sandwich product: R v R†
18///
19/// In 3D, a rotor encodes rotation by angle θ around an axis n as:
20/// R = cos(θ/2) + sin(θ/2) * B
21/// where B is the unit bivector representing the rotation plane.
22#[derive(Clone, Debug)]
23pub struct Rotor {
24    /// The internal multivector representation
25    inner: GA3,
26}
27
28impl Rotor {
29    /// Create a rotor from a multivector (assumes it's already a valid rotor)
30    pub fn from_multivector(mv: GA3) -> Self {
31        Self { inner: mv }
32    }
33
34    /// Create the identity rotor (no rotation)
35    pub fn identity() -> Self {
36        Self {
37            inner: GA3::scalar(1.0),
38        }
39    }
40
41    /// Create a rotor for rotation by angle (in radians) in the XY plane
42    ///
43    /// This is equivalent to rotation around the Z axis.
44    pub fn xy(angle: f64) -> Self {
45        Self::from_bivector_angle(angle, 1.0, 0.0, 0.0)
46    }
47
48    /// Create a rotor for rotation by angle (in radians) in the XZ plane
49    ///
50    /// This is equivalent to rotation around the Y axis.
51    pub fn xz(angle: f64) -> Self {
52        Self::from_bivector_angle(angle, 0.0, 1.0, 0.0)
53    }
54
55    /// Create a rotor for rotation by angle (in radians) in the YZ plane
56    ///
57    /// This is equivalent to rotation around the X axis.
58    pub fn yz(angle: f64) -> Self {
59        Self::from_bivector_angle(angle, 0.0, 0.0, 1.0)
60    }
61
62    /// Create a rotor from an angle and bivector components
63    ///
64    /// The bivector (xy, xz, yz) defines the rotation plane.
65    /// Components are normalized internally.
66    pub fn from_bivector_angle(angle: f64, xy: f64, xz: f64, yz: f64) -> Self {
67        let half_angle = angle / 2.0;
68
69        // Create unit bivector
70        // Note: Negate to match standard right-hand rotation convention
71        let biv = Bivector::<3, 0, 0>::from_components(-xy, -xz, -yz);
72        let biv_mv = GA3::from_bivector(&biv);
73        let mag = biv_mv.magnitude();
74
75        if mag < 1e-10 {
76            // Degenerate case - return identity
77            return Self::identity();
78        }
79
80        let biv_unit = &biv_mv * (1.0 / mag);
81
82        // R = cos(θ/2) + sin(θ/2) * B
83        let cos_part = GA3::scalar(half_angle.cos());
84        let sin_part = &biv_unit * half_angle.sin();
85
86        Self {
87            inner: &cos_part + &sin_part,
88        }
89    }
90
91    /// Create a rotor for rotation around an axis vector by an angle
92    ///
93    /// The axis does not need to be normalized.
94    pub fn from_axis_angle(axis_x: f64, axis_y: f64, axis_z: f64, angle: f64) -> Self {
95        // The bivector for rotation around axis (x,y,z) is proportional to
96        // x*(e2^e3) + y*(e3^e1) + z*(e1^e2)
97        // which in our basis order (e12, e13, e23) is:
98        // xy = z, xz = -y, yz = x
99        Self::from_bivector_angle(angle, axis_z, -axis_y, axis_x)
100    }
101
102    /// Get the internal multivector
103    pub fn as_multivector(&self) -> &GA3 {
104        &self.inner
105    }
106
107    /// Apply this rotor to transform a multivector (sandwich product)
108    ///
109    /// Returns R * v * R†
110    pub fn transform(&self, v: &GA3) -> GA3 {
111        let rev = self.inner.reverse();
112        self.inner.geometric_product(v).geometric_product(&rev)
113    }
114
115    /// Compose two rotors: self followed by other
116    ///
117    /// The result applies self first, then other.
118    pub fn then(&self, other: &Rotor) -> Rotor {
119        Rotor {
120            inner: other.inner.geometric_product(&self.inner),
121        }
122    }
123
124    /// Get the inverse rotor (reverse rotation)
125    pub fn inverse(&self) -> Rotor {
126        // For unit rotors, inverse = reverse
127        Rotor {
128            inner: self.inner.reverse(),
129        }
130    }
131
132    /// Normalize the rotor to unit magnitude
133    pub fn normalize(&self) -> Rotor {
134        match self.inner.normalize() {
135            Some(normalized) => Rotor { inner: normalized },
136            None => Self::identity(),
137        }
138    }
139
140    /// Get the rotation angle (in radians)
141    pub fn angle(&self) -> f64 {
142        // For a rotor R = cos(θ/2) + sin(θ/2)*B
143        // The scalar part is cos(θ/2)
144        let scalar = self.inner.get(0);
145        2.0 * scalar.clamp(-1.0, 1.0).acos()
146    }
147
148    /// Spherical linear interpolation between identity and this rotor
149    ///
150    /// t=0 gives identity, t=1 gives self
151    pub fn slerp(&self, t: f64) -> Rotor {
152        let angle = self.angle();
153        let new_angle = angle * t;
154
155        // Extract the bivector part and renormalize
156        // Note: Internal representation has negated bivector, so negate when extracting
157        let biv_xy = -self.inner.get(3); // e12
158        let biv_xz = -self.inner.get(5); // e13
159        let biv_yz = -self.inner.get(6); // e23
160        let biv_mag = (biv_xy * biv_xy + biv_xz * biv_xz + biv_yz * biv_yz).sqrt();
161
162        if biv_mag < 1e-10 {
163            // No rotation - return identity
164            Self::identity()
165        } else {
166            Self::from_bivector_angle(
167                new_angle,
168                biv_xy / biv_mag,
169                biv_xz / biv_mag,
170                biv_yz / biv_mag,
171            )
172        }
173    }
174
175    /// Spherical linear interpolation between two rotors
176    pub fn slerp_to(&self, other: &Rotor, t: f64) -> Rotor {
177        // Compute relative rotation: other = relative * self
178        // So relative = other * self.inverse()
179        let relative = other.then(&self.inverse());
180        let interpolated = relative.slerp(t);
181        interpolated.then(self)
182    }
183}
184
185/// A versor is a general geometric transformation (rotation, reflection, or their composition).
186///
187/// Every versor can be written as a product of vectors.
188/// Versors with an even number of vectors are rotors.
189/// Versors with an odd number of vectors include reflections.
190#[derive(Clone, Debug)]
191pub struct Versor {
192    /// The internal multivector representation
193    inner: GA3,
194    /// Whether this is an even versor (rotor) or odd versor (includes reflection)
195    is_even: bool,
196}
197
198impl Versor {
199    /// Create a versor from a multivector
200    pub fn from_multivector(mv: GA3, is_even: bool) -> Self {
201        Self { inner: mv, is_even }
202    }
203
204    /// Create the identity versor
205    pub fn identity() -> Self {
206        Self {
207            inner: GA3::scalar(1.0),
208            is_even: true,
209        }
210    }
211
212    /// Create a reflection through a plane with normal vector (x, y, z)
213    ///
214    /// Reflects points through the plane passing through origin with the given normal.
215    pub fn reflection(normal_x: f64, normal_y: f64, normal_z: f64) -> Self {
216        // A reflection is represented by the unit vector normal to the plane
217        let vec = Vector::<3, 0, 0>::from_components(normal_x, normal_y, normal_z);
218        let mv = GA3::from_vector(&vec);
219        let normalized = mv
220            .normalize()
221            .unwrap_or_else(|| GA3::from_vector(&Vector::from_components(1.0, 0.0, 0.0)));
222
223        Self {
224            inner: normalized,
225            is_even: false,
226        }
227    }
228
229    /// Convert from a Rotor
230    pub fn from_rotor(rotor: Rotor) -> Self {
231        Self {
232            inner: rotor.inner,
233            is_even: true,
234        }
235    }
236
237    /// Get the internal multivector
238    pub fn as_multivector(&self) -> &GA3 {
239        &self.inner
240    }
241
242    /// Apply this versor to transform a multivector
243    ///
244    /// For even versors: V * v * V†
245    /// For odd versors: V * v * V† (with grade involution handling)
246    pub fn transform(&self, v: &GA3) -> GA3 {
247        let rev = self.inner.reverse();
248        if self.is_even {
249            self.inner.geometric_product(v).geometric_product(&rev)
250        } else {
251            // For odd versors, we need to handle the sign change
252            // This is simplified - full implementation would use grade involution
253            let result = self.inner.geometric_product(v).geometric_product(&rev);
254            &result * -1.0
255        }
256    }
257
258    /// Compose two versors
259    pub fn then(&self, other: &Versor) -> Versor {
260        Versor {
261            inner: other.inner.geometric_product(&self.inner),
262            is_even: self.is_even == other.is_even, // Even * Even = Even, Odd * Odd = Even
263        }
264    }
265
266    /// Check if this is an even versor (rotor)
267    pub fn is_rotor(&self) -> bool {
268        self.is_even
269    }
270
271    /// Try to convert to a Rotor (only succeeds for even versors)
272    pub fn to_rotor(&self) -> Option<Rotor> {
273        if self.is_even {
274            Some(Rotor {
275                inner: self.inner.clone(),
276            })
277        } else {
278            None
279        }
280    }
281}
282
283/// A translation represented geometrically
284///
285/// In standard GA3 (Euclidean), translations are not directly representable
286/// as versors. This type provides a convenient API that internally uses
287/// vector addition.
288#[derive(Clone, Debug)]
289pub struct Translation {
290    /// Translation vector components
291    x: f64,
292    y: f64,
293    z: f64,
294}
295
296impl Translation {
297    /// Create a new translation
298    pub fn new(x: f64, y: f64, z: f64) -> Self {
299        Self { x, y, z }
300    }
301
302    /// Create a translation along the X axis
303    pub fn x(amount: f64) -> Self {
304        Self::new(amount, 0.0, 0.0)
305    }
306
307    /// Create a translation along the Y axis
308    pub fn y(amount: f64) -> Self {
309        Self::new(0.0, amount, 0.0)
310    }
311
312    /// Create a translation along the Z axis
313    pub fn z(amount: f64) -> Self {
314        Self::new(0.0, 0.0, amount)
315    }
316
317    /// Apply this translation to a multivector
318    ///
319    /// For vectors, this adds the translation. For other grades,
320    /// behavior depends on the geometric interpretation.
321    pub fn transform(&self, v: &GA3) -> GA3 {
322        let trans_vec = Vector::<3, 0, 0>::from_components(self.x, self.y, self.z);
323        let trans_mv = GA3::from_vector(&trans_vec);
324        v + &trans_mv
325    }
326
327    /// Compose two translations
328    pub fn then(&self, other: &Translation) -> Translation {
329        Translation {
330            x: self.x + other.x,
331            y: self.y + other.y,
332            z: self.z + other.z,
333        }
334    }
335
336    /// Get the inverse translation
337    pub fn inverse(&self) -> Translation {
338        Translation {
339            x: -self.x,
340            y: -self.y,
341            z: -self.z,
342        }
343    }
344
345    /// Linear interpolation of translation
346    pub fn lerp(&self, t: f64) -> Translation {
347        Translation {
348            x: self.x * t,
349            y: self.y * t,
350            z: self.z * t,
351        }
352    }
353
354    /// Linear interpolation to another translation
355    pub fn lerp_to(&self, other: &Translation, t: f64) -> Translation {
356        Translation {
357            x: self.x + (other.x - self.x) * t,
358            y: self.y + (other.y - self.y) * t,
359            z: self.z + (other.z - self.z) * t,
360        }
361    }
362}
363
364/// A general geometric transformation combining rotation and translation
365#[derive(Clone, Debug)]
366pub struct Transform {
367    /// Rotation component
368    rotor: Rotor,
369    /// Translation component (applied after rotation)
370    translation: Translation,
371}
372
373impl Transform {
374    /// Create a new transform with rotation and translation
375    pub fn new(rotor: Rotor, translation: Translation) -> Self {
376        Self { rotor, translation }
377    }
378
379    /// Create identity transform
380    pub fn identity() -> Self {
381        Self {
382            rotor: Rotor::identity(),
383            translation: Translation::new(0.0, 0.0, 0.0),
384        }
385    }
386
387    /// Create a pure rotation transform
388    pub fn rotation(rotor: Rotor) -> Self {
389        Self {
390            rotor,
391            translation: Translation::new(0.0, 0.0, 0.0),
392        }
393    }
394
395    /// Create a pure translation transform
396    pub fn translation(translation: Translation) -> Self {
397        Self {
398            rotor: Rotor::identity(),
399            translation,
400        }
401    }
402
403    /// Apply this transform to a multivector
404    pub fn transform(&self, v: &GA3) -> GA3 {
405        let rotated = self.rotor.transform(v);
406        self.translation.transform(&rotated)
407    }
408
409    /// Compose two transforms: self followed by other
410    pub fn then(&self, other: &Transform) -> Transform {
411        // First apply self's rotation, then self's translation
412        // Then apply other's rotation, then other's translation
413        //
414        // Combined rotation: other.rotor * self.rotor
415        // Translation is more complex due to rotation interaction
416        let combined_rotor = self.rotor.then(&other.rotor);
417
418        // Transform self's translation by other's rotation, then add other's translation
419        let self_trans_vec = Vector::<3, 0, 0>::from_components(
420            self.translation.x,
421            self.translation.y,
422            self.translation.z,
423        );
424        let self_trans_mv = GA3::from_vector(&self_trans_vec);
425        let rotated_trans = other.rotor.transform(&self_trans_mv);
426
427        let combined_translation = Translation::new(
428            rotated_trans.get(1) + other.translation.x,
429            rotated_trans.get(2) + other.translation.y,
430            rotated_trans.get(4) + other.translation.z,
431        );
432
433        Transform {
434            rotor: combined_rotor,
435            translation: combined_translation,
436        }
437    }
438
439    /// Get the inverse transform
440    pub fn inverse(&self) -> Transform {
441        let inv_rotor = self.rotor.inverse();
442        let inv_trans_base = self.translation.inverse();
443
444        // Rotate the inverse translation by the inverse rotation
445        let trans_vec = Vector::<3, 0, 0>::from_components(
446            inv_trans_base.x,
447            inv_trans_base.y,
448            inv_trans_base.z,
449        );
450        let trans_mv = GA3::from_vector(&trans_vec);
451        let rotated = inv_rotor.transform(&trans_mv);
452
453        Transform {
454            rotor: inv_rotor,
455            translation: Translation::new(rotated.get(1), rotated.get(2), rotated.get(4)),
456        }
457    }
458
459    /// Interpolate between identity and this transform
460    pub fn interpolate(&self, t: f64) -> Transform {
461        Transform {
462            rotor: self.rotor.slerp(t),
463            translation: self.translation.lerp(t),
464        }
465    }
466
467    /// Interpolate to another transform
468    pub fn interpolate_to(&self, other: &Transform, t: f64) -> Transform {
469        Transform {
470            rotor: self.rotor.slerp_to(&other.rotor, t),
471            translation: self.translation.lerp_to(&other.translation, t),
472        }
473    }
474}
475
476#[cfg(test)]
477mod tests {
478    use super::*;
479    use std::f64::consts::PI;
480
481    #[test]
482    fn test_rotor_identity() {
483        let r = Rotor::identity();
484        let v = Vector::<3, 0, 0>::from_components(1.0, 2.0, 3.0);
485        let mv = GA3::from_vector(&v);
486
487        let rotated = r.transform(&mv);
488
489        // Identity should not change the vector
490        assert!((rotated.get(1) - 1.0).abs() < 1e-10);
491        assert!((rotated.get(2) - 2.0).abs() < 1e-10);
492        assert!((rotated.get(4) - 3.0).abs() < 1e-10);
493    }
494
495    #[test]
496    fn test_rotor_xy_90() {
497        // 90 degree rotation in XY plane
498        let r = Rotor::xy(PI / 2.0);
499
500        // Rotate (1, 0, 0) should give (0, 1, 0)
501        let v = GA3::from_vector(&Vector::from_components(1.0, 0.0, 0.0));
502        let rotated = r.transform(&v);
503
504        assert!((rotated.get(1) - 0.0).abs() < 1e-10); // x -> 0
505        assert!((rotated.get(2) - 1.0).abs() < 1e-10); // y -> 1
506        assert!((rotated.get(4) - 0.0).abs() < 1e-10); // z -> 0
507    }
508
509    #[test]
510    fn test_rotor_preserves_magnitude() {
511        let r = Rotor::from_bivector_angle(1.23, 1.0, 2.0, 3.0);
512        let v = GA3::from_vector(&Vector::from_components(3.0, 4.0, 5.0));
513
514        let original_mag = v.magnitude();
515        let rotated = r.transform(&v);
516        let rotated_mag = rotated.magnitude();
517
518        assert!(
519            (original_mag - rotated_mag).abs() < 1e-10,
520            "Magnitude changed: {} -> {}",
521            original_mag,
522            rotated_mag
523        );
524    }
525
526    #[test]
527    fn test_rotor_composition() {
528        // Two 90-degree rotations should equal one 180-degree rotation
529        let r1 = Rotor::xy(PI / 2.0);
530        let r2 = Rotor::xy(PI / 2.0);
531        let r_combined = r1.then(&r2);
532
533        let v = GA3::from_vector(&Vector::from_components(1.0, 0.0, 0.0));
534        let rotated = r_combined.transform(&v);
535
536        // (1, 0, 0) rotated 180 degrees should give (-1, 0, 0)
537        assert!((rotated.get(1) + 1.0).abs() < 1e-10);
538        assert!((rotated.get(2) - 0.0).abs() < 1e-10);
539    }
540
541    #[test]
542    fn test_rotor_inverse() {
543        let r = Rotor::from_bivector_angle(0.7, 1.0, 1.0, 1.0);
544        let v = GA3::from_vector(&Vector::from_components(1.0, 2.0, 3.0));
545
546        let rotated = r.transform(&v);
547        let back = r.inverse().transform(&rotated);
548
549        assert!((back.get(1) - 1.0).abs() < 1e-10);
550        assert!((back.get(2) - 2.0).abs() < 1e-10);
551        assert!((back.get(4) - 3.0).abs() < 1e-10);
552    }
553
554    #[test]
555    fn test_rotor_slerp() {
556        let r = Rotor::xy(PI);
557
558        // Halfway interpolation should give PI/2 rotation
559        let half = r.slerp(0.5);
560
561        let v = GA3::from_vector(&Vector::from_components(1.0, 0.0, 0.0));
562        let rotated = half.transform(&v);
563
564        // (1, 0, 0) rotated 90 degrees should give (0, 1, 0)
565        assert!((rotated.get(1) - 0.0).abs() < 1e-10);
566        assert!((rotated.get(2) - 1.0).abs() < 1e-10);
567    }
568
569    #[test]
570    fn test_translation() {
571        let t = Translation::new(1.0, 2.0, 3.0);
572        let v = GA3::from_vector(&Vector::from_components(0.0, 0.0, 0.0));
573
574        let translated = t.transform(&v);
575
576        assert!((translated.get(1) - 1.0).abs() < 1e-10);
577        assert!((translated.get(2) - 2.0).abs() < 1e-10);
578        assert!((translated.get(4) - 3.0).abs() < 1e-10);
579    }
580
581    #[test]
582    fn test_transform_composition() {
583        let rot = Rotor::xy(PI / 2.0);
584        let trans = Translation::new(1.0, 0.0, 0.0);
585        let transform = Transform::new(rot, trans);
586
587        let v = GA3::from_vector(&Vector::from_components(1.0, 0.0, 0.0));
588        let result = transform.transform(&v);
589
590        // First rotate (1,0,0) by 90deg -> (0,1,0)
591        // Then translate by (1,0,0) -> (1,1,0)
592        assert!((result.get(1) - 1.0).abs() < 1e-10);
593        assert!((result.get(2) - 1.0).abs() < 1e-10);
594    }
595}