1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103
use core::ops::{Neg, Add, Sub, Mul, Div}; use num_traits::{Zero, One, NumCast}; use crate::{*, transform::*}; #[derive(Copy, Clone, Debug, PartialEq)] pub struct Moebius<U> { data: [U; 4], } impl<U> From<[U; 4]> for Moebius<U> { fn from(array: [U; 4]) -> Self { Self { data: array } } } impl<U> Into<[U; 4]> for Moebius<U> { fn into(self) -> [U; 4] { self.data } } impl<U> Moebius<U> { pub fn new(a: U, b: U, c: U, d: U) -> Self { Self::from([a, b, c, d]) } pub fn a_ref(&self) -> &U { &self.data[0] } pub fn b_ref(&self) -> &U { &self.data[1] } pub fn c_ref(&self) -> &U { &self.data[2] } pub fn d_ref(&self) -> &U { &self.data[3] } pub fn a_mut(&mut self) -> &mut U { &mut self.data[0] } pub fn b_mut(&mut self) -> &mut U { &mut self.data[1] } pub fn c_mut(&mut self) -> &mut U { &mut self.data[2] } pub fn d_mut(&mut self) -> &mut U { &mut self.data[3] } } impl<U: Clone> Moebius<U> { pub fn a(&self) -> U { self.data[0].clone() } pub fn b(&self) -> U { self.data[1].clone() } pub fn c(&self) -> U { self.data[2].clone() } pub fn d(&self) -> U { self.data[3].clone() } } impl<U: Zero + One> Identity for Moebius<U> { fn identity() -> Self { Self::new(U::one(), U::zero(), U::zero(), U::one()) } } impl<U> Chain<U> for Moebius<U> where U: Add<Output=U> + Mul<Output=U> + Div<Output=U> + Clone { fn chain(self, other: Self) -> Self { Self::new( self.a()*other.a() + self.b()*other.c(), self.a()*other.b() + self.b()*other.d(), self.c()*other.a() + self.d()*other.c(), self.c()*other.b() + self.d()*other.d(), ) } } impl<U> Transform<U> for Moebius<U> where U: Add<Output=U> + Mul<Output=U> + Div<Output=U> + Clone { fn apply(&self, x: U) -> U { (self.a()*x.clone() + self.b())/(self.c()*x + self.d()) } } impl<T: Algebra + Clone, U: Algebra<T> + Clone> Transform<Construct<T, Construct<T, U>>> for Moebius<Construct<T, U>> { fn apply(&self, x: Construct<T, Construct<T, U>>) -> Construct<T, Construct<T, U>> { (self.a()*x.clone() + self.b())/(self.c()*x + self.d()) } } impl<U: Neg<Output=U> + Mul<Output=U> + Div<Output=U> + Sub<Output=U> + Clone> Moebius<U> { pub fn det(&self) -> U { self.a()*self.d() - self.b()*self.c() } pub fn normalize(mut self) -> Self { let det = self.det(); self.data.iter_mut().for_each(|x| *x = x.clone() / det.clone()); self } } impl<T: Algebra + Clone> Deriv<Complex<T>> for Moebius<Complex<T>> { fn deriv(&self, p: Complex<T>) -> Complex<T> { let u: Complex<T> = self.a() * p.clone() + self.b(); let d: Complex<T> = self.c() * p + self.d(); return (self.a() * d.clone() - u * self.c()) / (d.clone() * d); } } impl<T: NumCast + Algebra + Dot<Output=T> + Clone> DerivDir<Quaternion<T>> for Moebius<Complex<T>> { fn deriv_dir(&self, p: Quaternion<T>, v: Quaternion<T>) -> Quaternion<T> { let u = self.a() * p.clone() + self.b(); let d = self.c() * p + self.d(); let d2 = d.clone().abs_sqr(); let g1 = (self.a() * v.clone()) / d.clone(); let g21 = (self.c() * v.clone()).conj(); let g22 = d.clone().conj() * (d.dot(self.c() * v) * T::from(2).unwrap() / d2.clone()); let g2 = u * ((g21 - g22) / d2); return g1 + g2; } }