From cdfbe7f0a21129e55fffc070ba249faf96256cdc Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Fri, 25 Oct 2024 14:47:44 +0200 Subject: [PATCH] New rotation actions (#765) * New rotation actions * documentation improvements * export two symbols --- NEWS.md | 9 ++- src/Manifolds.jl | 2 + src/groups/rotation_action.jl | 105 +++++++++++++++++++++++++++++++++ test/groups/rotation_action.jl | 39 ++++++++++++ 4 files changed, 154 insertions(+), 1 deletion(-) diff --git a/NEWS.md b/NEWS.md index 741b7be6b..447a74546 100644 --- a/NEWS.md +++ b/NEWS.md @@ -5,7 +5,14 @@ All notable changes to this project will be documented in this file. The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/), and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html). -## [0.10.5] - 2024-10-24 +## [0.10.6] – 2024-10-25 + +### Added + +* Two new actions: `ComplexPlanarRotation`, `QuaternionRotation`. +* New function `quaternion_rotation_matrix` for converting quaternions to rotation matrices. + +## [0.10.5] – 2024-10-24 ### Added diff --git a/src/Manifolds.jl b/src/Manifolds.jl index 5c81913a5..2b36221c2 100644 --- a/src/Manifolds.jl +++ b/src/Manifolds.jl @@ -987,6 +987,7 @@ export AbstractGroupAction, ActionDirection, AdditionOperation, CircleGroup, + ComplexPlanarRotation, GeneralLinear, GroupManifold, GroupOperationAction, @@ -1001,6 +1002,7 @@ export AbstractGroupAction, PowerGroup, ProductGroup, ProductOperation, + QuaternionRotation, RealCircleGroup, RightAction, RightInvariantMetric, diff --git a/src/groups/rotation_action.jl b/src/groups/rotation_action.jl index 35e1bc86a..5cfb8e16c 100644 --- a/src/groups/rotation_action.jl +++ b/src/groups/rotation_action.jl @@ -304,3 +304,108 @@ function optimal_alignment( Ostar = det(UVt) ≥ 0 ? UVt : F.U * Diagonal([i < L ? 1 : -1 for i in 1:L]) * F.Vt return convert(typeof(Xmul), Ostar) end + +@doc raw""" + ComplexPlanarRotation() + +Action of the circle group [`CircleGroup`](@ref) on ``ℝ^2`` by left multiplication. +""" +struct ComplexPlanarRotation <: AbstractGroupAction{LeftAction} end + +base_group(::ComplexPlanarRotation) = CircleGroup() + +group_manifold(::ComplexPlanarRotation) = Euclidean(2) + +@doc raw""" + apply(A::ComplexPlanarRotation, g::Complex, p) + +Rotate point `p` from [`Euclidean(2)`](@ref) manifold by the group element `g`. +The formula reads +````math +p_{rot} = \begin{bmatrix} +\cos(θ) & \sin(θ)\\ +-\sin(θ) & \cos(θ) +\end{bmatrix} p, +```` +where `θ` is the argument of complex number `g`. +""" +function apply(::ComplexPlanarRotation, g::Complex, p) + sinθ, cosθ = g.im, g.re + return (@SMatrix [cosθ -sinθ; sinθ cosθ]) * p +end +apply(::ComplexPlanarRotation, ::Identity{MultiplicationOperation}, p) = p + +function apply!(A::ComplexPlanarRotation, q, g::Complex, p) + return copyto!(q, apply(A, g, p)) +end +function apply!(::ComplexPlanarRotation, q, ::Identity{MultiplicationOperation}, p) + return copyto!(q, p) +end + +@doc raw""" + QuaternionRotation + +Action of the unit quaternion group [`Unitary`](@ref)`(1, ℍ)` on ``ℝ^3``. +""" +struct QuaternionRotation <: AbstractGroupAction{LeftAction} end + +base_group(::QuaternionRotation) = Unitary(1, ℍ) + +group_manifold(::QuaternionRotation) = Euclidean(3) + +@doc raw""" + apply(A::QuaternionRotation, g::Quaternion, p) + +Rotate point `p` from [`Euclidean`](@ref)`(3)` manifold through conjugation by the group +element `g`. +The formula reads +````math +(0, p_{rot,x}, p_{rot,y}, p_{rot,z}) = g ⋅ (0, p_x, p_y, p_z) ⋅ g^{\mathrm{H}} +```` +where ``(0, p_x, p_y, p_z)`` is quaternion with non-real coefficients from encoding +the point `p` and ``g^{\mathrm{H}}`` is quaternion conjugate of ``g``. +""" +function apply(::QuaternionRotation, g::Quaternions.Quaternion, p::SVector) + p_quat = Quaternions.Quaternion(0, p[1], p[2], p[3]) + p_conj = g * p_quat * conj(g) + return @SVector [p_conj.v1, p_conj.v2, p_conj.v3] +end +apply(::QuaternionRotation, ::Identity{MultiplicationOperation}, p) = p + +function apply!(::QuaternionRotation, q, g::Quaternions.Quaternion, p) + p_quat = Quaternions.Quaternion(0, p[1], p[2], p[3]) + p_conj = g * p_quat * conj(g) + q[1] = p_conj.v1 + q[2] = p_conj.v2 + q[3] = p_conj.v3 + return q +end +function apply!(::QuaternionRotation, q, ::Identity{MultiplicationOperation}, p) + return copyto!(q, p) +end + +""" + quaternion_rotation_matrix(g::Quaternions.Quaternion) + +Compute rotation matrix for [`RotationAction`](@ref) corresponding to +[`QuaternionRotation`](@ref) by `g`. + +See https://www.songho.ca/opengl/gl_quaternion.html for details. +""" +function quaternion_rotation_matrix(g::Quaternions.Quaternion) + r11 = 1 - 2 * (g.v2^2 + g.v3^2) + r12 = 2 * (g.v1 * g.v2 - g.v3 * g.s) + r13 = 2 * (g.v1 * g.v3 + g.v2 * g.s) + r21 = 2 * (g.v1 * g.v2 + g.v3 * g.s) + r22 = 1 - 2 * (g.v1^2 + g.v3^2) + r23 = 2 * (g.v2 * g.v3 - g.v1 * g.s) + r31 = 2 * (g.v1 * g.v3 - g.v2 * g.s) + r32 = 2 * (g.v2 * g.v3 + g.v1 * g.s) + r33 = 1 - 2 * (g.v1^2 + g.v2^2) + + return @SMatrix [ + r11 r12 r13 + r21 r22 r23 + r31 r32 r33 + ] +end diff --git a/test/groups/rotation_action.jl b/test/groups/rotation_action.jl index 81db444ce..bcd046e22 100644 --- a/test/groups/rotation_action.jl +++ b/test/groups/rotation_action.jl @@ -182,3 +182,42 @@ end apply!(A, q, a, p1) @test isapprox(q, apply(A, a, p1)) end + +@testset "ComplexPlanarRotation" begin + M = Euclidean(2) + G = CircleGroup() + A = ComplexPlanarRotation() + + @test group_manifold(A) === M + @test base_group(A) === G + + p = [1.0, 2.0] + @test apply(A, 1.0im, p) ≈ [-2.0, 1.0] + @test apply(A, Identity(G), p) == p + q = similar(p) + apply!(A, q, 1.0im, p) + @test isapprox(q, [-2.0, 1.0]) + apply!(A, q, Identity(G), p) + @test q == p +end + +@testset "QuaternionRotation" begin + M = Euclidean(3) + G = Unitary(1, ℍ) + A = QuaternionRotation() + + @test group_manifold(A) === M + @test base_group(A) === G + + p = @SVector [1.0, 2.0, 5.0] + g = sign(Quaternion(1.0, 3.0, 5.0, 7.0)) + @test apply(A, g, p) ≈ [2.7142857142857135, 3.571428571428571, 3.142857142857142] + @test apply(A, Identity(G), p) == p + q = similar(p) + apply!(A, q, g, p) + @test isapprox(q, [2.7142857142857135, 3.571428571428571, 3.142857142857142]) + apply!(A, q, Identity(G), p) + @test q == p + + @test Manifolds.quaternion_rotation_matrix(g) * p ≈ apply(A, g, p) +end