diff --git a/RELEASES.md b/RELEASES.md index 33a459d..02701c7 100644 --- a/RELEASES.md +++ b/RELEASES.md @@ -1,3 +1,8 @@ +# Version 0.4.0 (2024-02-27) + + * Make Wigner rotation n-dimensional. + * Negate `FrameN` axis instead of beta. + # Version 0.3.0 (2024-02-21) * Update dependencies. diff --git a/src/lib.rs b/src/lib.rs index d914de2..864ba11 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -36,8 +36,8 @@ use nalgebra::{ ShapeConstraint, }, storage::{Owned, RawStorage, RawStorageMut, Storage, StorageMut}, - DefaultAllocator, Dim, DimName, DimNameDiff, DimNameSub, Matrix, Matrix4, MatrixView, - MatrixViewMut, OMatrix, OVector, Rotation3, Scalar, SimdRealField, Unit, Vector3, VectorView, + DefaultAllocator, Dim, DimName, DimNameDiff, DimNameSub, Matrix, MatrixView, MatrixViewMut, + OMatrix, OVector, Scalar, SimdRealField, Unit, VectorView, }; use num_traits::{real::Real, sign::Signed}; use std::ops::{Add, Neg, Sub}; @@ -996,25 +996,39 @@ where { frame.velocity().boost(&-self.clone()).frame() } -} -impl FrameN -where - ShapeConstraint: DimEq, - N: SimdRealField + Signed + Real, - D: DimNameSub, - DefaultAllocator: Allocator>, -{ + /// Wigner rotation angle $\epsilon$ of the boost composition `self`$\oplus$`frame`. + /// + /// The angle between the forward and backward composition: + /// + /// $$ + /// \epsilon = \angle (\vec \beta_u \oplus \vec \beta_v, \vec \beta_v \oplus \vec \beta_u) + /// $$ + /// + /// See [`Self::rotation`] for further details. + #[must_use] + pub fn angle(&self, frame: &Self) -> N + where + DefaultAllocator: Allocator, + { + let (u, v) = (self, frame); + let ucv = u.compose(v); + let vcu = v.compose(u); + ucv.axis().angle(&vcu.axis()) + } + /// $ /// \gdef \Bu {B^{\mu'}\_{\phantom {\mu'} \mu} (\vec \beta_u)} /// \gdef \Bv {B^{\mu''}\_{\phantom {\mu''} \mu'} (\vec \beta_v)} /// \gdef \Puv {u \oplus v} /// \gdef \Buv {B^{\mu'}\_{\phantom {\mu'} \mu} (\vec \beta_{\Puv})} /// \gdef \Ruv {R^{\mu''}\_{\phantom {\mu''} \mu'} (\epsilon)} + /// \gdef \R {R (\epsilon)} + /// \gdef \Kuv {K(\epsilon)} /// \gdef \Luv {\Lambda^{\mu''}\_{\phantom {\mu''} \mu} (\vec \beta_{\Puv})} /// $ - /// Wigner rotation axis $\widehat {\vec \beta_u \times \vec \beta_v}$ and angle $\epsilon$ of - /// the boost composition `self`$\oplus$`frame`. + /// Wigner rotation matrix $R(\widehat {\vec \beta_u \wedge \vec \beta_v}, \epsilon)$ of the + /// boost composition `self`$\oplus$`frame`. /// /// The composition of two pure boosts, $\Bu$ to `self` followed by $\Bv$ to `frame`, results in /// a composition of a pure boost $\Buv$ and a *non-vanishing* spatial rotation $\Ruv$ for @@ -1024,6 +1038,77 @@ where /// \Luv = \Ruv \Buv = \Bv \Bu /// $$ /// + /// The returned homogeneous rotation matrix + /// + /// $$ + /// \R = \begin{pmatrix} + /// 1 & \vec{0}^T \\\\ + /// \vec{0} & \Kuv \end{pmatrix} + /// $$ + /// + /// embeds the spatial rotation matrix + /// + /// $$ + /// \Kuv = I + \sin \epsilon (\hat v \hat u^T - \hat u \hat v^T) + /// + (\cos \epsilon - 1) (\hat u \hat u^T + \hat v \hat v^T) + /// $$ + /// + /// rotating in the plane spanned by the orthonormalized pair $\hat u = \hat \beta_u$ and + /// $\hat v = \widehat{\hat \beta_v - (\hat \beta_u \cdot \hat \beta_v) \hat \beta_u}$. The + /// rotation angle $\epsilon$ is found according to [`Self::angle`]. + /// + /// See [`Self::axis_angle`] for $4$-dimensional specialization. + /// + /// ``` + /// use approx::{assert_ulps_eq, assert_ulps_ne}; + /// use nalgebra::{Matrix4, Vector3}; + /// use nalgebra_spacetime::{Frame4, LorentzianN}; + /// + /// let u = Frame4::from_beta(Vector3::new(0.18, 0.73, 0.07)); + /// let v = Frame4::from_beta(Vector3::new(0.41, 0.14, 0.25)); + /// let ucv = u.compose(&v); + /// let vcu = v.compose(&u); + /// + /// let boost_u = Matrix4::new_boost(&u); + /// let boost_v = Matrix4::new_boost(&v); + /// let boost_ucv = Matrix4::new_boost(&ucv); + /// let boost_vcu = Matrix4::new_boost(&vcu); + /// + /// let (matrix_ucv, angle_ucv) = u.rotation(&v); + /// + /// assert_ulps_eq!(angle_ucv, u.angle(&v)); + /// assert_ulps_ne!(boost_ucv, boost_v * boost_u); + /// assert_ulps_ne!(boost_vcu, boost_u * boost_v); + /// assert_ulps_eq!(matrix_ucv * boost_ucv, boost_v * boost_u); + /// assert_ulps_eq!(boost_vcu * matrix_ucv, boost_v * boost_u); + /// ``` + #[must_use] + pub fn rotation(&self, frame: &Self) -> (OMatrix, N) + where + DefaultAllocator: Allocator + + Allocator + + Allocator> + + Allocator, DimNameDiff>, + { + let [u, v] = [self, frame]; + let ang = u.angle(v); + let (sin, cos) = ang.sin_cos(); + let u = u.axis().into_inner(); + let v = v.axis().into_inner(); + let v = (&v - &u * u.dot(&v)).normalize(); + let ut = u.transpose(); + let vt = v.transpose(); + let rot = (&v * &ut - &u * &vt) * sin + (&u * &ut + &v * &vt) * (cos - N::one()); + let mut mat = OMatrix::::identity(); + let (r, c) = mat.shape_generic(); + let mut sub = mat.generic_view_mut((1, 1), (r.sub(U1), c.sub(U1))); + sub += rot; + (mat, ang) + } + + /// Wigner rotation axis $\widehat {\vec \beta_u \times \vec \beta_v}$ and angle $\epsilon$ of + /// the boost composition `self`$\oplus$`frame`. + /// /// $$ /// \epsilon = \arcsin \Bigg({ /// 1 + \gamma + \gamma_u + \gamma_v @@ -1036,6 +1121,8 @@ where /// \gamma = \gamma_u \gamma_v (1 + \vec \beta_u \cdot \vec \beta_v) /// $$ /// + /// See [`Self::rotation`] for $n$-dimensional generalization. + /// /// ``` /// use approx::{assert_abs_diff_ne, assert_ulps_eq}; /// use nalgebra::Vector3; @@ -1055,7 +1142,13 @@ where /// assert_ulps_eq!(axis, ucv.cross(&vcu).normalize(), epsilon = 1e-15); /// ``` #[must_use] - pub fn axis_angle(&self, frame: &Self) -> (Unit>>, N) { + pub fn axis_angle(&self, frame: &Self) -> (Unit>>, N) + where + N: SimdRealField + Signed + Real, + D: DimNameSub, + ShapeConstraint: DimEq, + DefaultAllocator: Allocator>, + { let (u, v) = (self, frame); let (axis, sin) = Unit::new_and_get(u.axis().cross(&v.axis())); let uv = u.axis().dot(&v.axis()); @@ -1067,50 +1160,6 @@ where let prod = (N::one() + cg) * (N::one() + ug) * (N::one() + vg); (axis, (sum / prod * bg * sin).asin()) } - - /// Wigner rotation matrix $R(\widehat {\vec \beta_u \times \vec \beta_v}, \epsilon)$ of the - /// boost composition `self`$\oplus$`frame`. - /// - /// See [`Self::axis_angle`] for further details. - /// - /// ``` - /// use approx::{assert_ulps_eq, assert_ulps_ne}; - /// use nalgebra::{Matrix4, Vector3}; - /// use nalgebra_spacetime::{Frame4, LorentzianN}; - /// - /// let u = Frame4::from_beta(Vector3::new(0.18, 0.73, 0.07)); - /// let v = Frame4::from_beta(Vector3::new(0.41, 0.14, 0.25)); - /// let ucv = u.compose(&v); - /// let vcu = v.compose(&u); - /// - /// let boost_u = Matrix4::new_boost(&u); - /// let boost_v = Matrix4::new_boost(&v); - /// let boost_ucv = Matrix4::new_boost(&ucv); - /// let boost_vcu = Matrix4::new_boost(&vcu); - /// - /// let rotation_ucv = u.rotation(&v); - /// - /// assert_ulps_ne!(boost_ucv, boost_v * boost_u); - /// assert_ulps_ne!(boost_vcu, boost_u * boost_v); - /// assert_ulps_eq!(rotation_ucv * boost_ucv, boost_v * boost_u); - /// assert_ulps_eq!(boost_vcu * rotation_ucv, boost_v * boost_u); - /// ``` - #[must_use] - pub fn rotation(&self, frame: &Self) -> Matrix4 - where - N::Element: SimdRealField, - DefaultAllocator: Allocator, - { - let mut m = Matrix4::::identity(); - let (axis, angle) = self.axis_angle(frame); - let axis = axis.into_inner(); - let axis = Unit::new_unchecked(Vector3::new(axis[0], axis[1], axis[2])); - let r = Rotation3::from_axis_angle(&axis, angle); - let (rows, cols) = m.shape_generic(); - m.generic_view_mut((1, 1), (rows.sub(U1), cols.sub(U1))) - .copy_from(r.matrix()); - m - } } impl From> for FrameN