From 2104befae9261a20cabb0994dc4060935ac17b62 Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Wed, 8 Nov 2023 19:46:22 +0100 Subject: [PATCH] Some issues with allocation on complex manifolds (#677) * Some issues with allocation on complex manifolds * minor cleanup, export number_of_coordinates --- NEWS.md | 7 +++++++ Project.toml | 2 +- src/Manifolds.jl | 1 + src/groups/general_unitary_groups.jl | 8 +++++++ src/groups/group.jl | 5 +++-- src/groups/unitary.jl | 2 +- src/manifolds/GeneralUnitaryMatrices.jl | 28 ++++++++++++------------- test/groups/general_unitary_groups.jl | 2 ++ 8 files changed, 37 insertions(+), 18 deletions(-) diff --git a/NEWS.md b/NEWS.md index 869178c2f6..a4a7b10ee3 100644 --- a/NEWS.md +++ b/NEWS.md @@ -5,6 +5,13 @@ 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.9.5] - 2023-11-08 + +### Changed + +- `identity_element` now returns a complex matrix for unitary group. +- `number_of_coordinates` is now exported. + ## [0.9.4] - 2023-11-06 ### Added diff --git a/Project.toml b/Project.toml index 97d08c3128..c65d8447f3 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "Manifolds" uuid = "1cead3c2-87b3-11e9-0ccd-23c62b72b94e" authors = ["Seth Axen ", "Mateusz Baran ", "Ronny Bergmann ", "Antoine Levitt "] -version = "0.9.4" +version = "0.9.5" [deps] Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" diff --git a/src/Manifolds.jl b/src/Manifolds.jl index d5f551a0fb..3ca4e4c428 100644 --- a/src/Manifolds.jl +++ b/src/Manifolds.jl @@ -856,6 +856,7 @@ export ×, norm, normal_tvector_distribution, number_eltype, + number_of_coordinates, one, power_dimensions, parallel_transport_along, diff --git a/src/groups/general_unitary_groups.jl b/src/groups/general_unitary_groups.jl index 077497c0fd..a280caba26 100644 --- a/src/groups/general_unitary_groups.jl +++ b/src/groups/general_unitary_groups.jl @@ -45,6 +45,14 @@ function allocate_result(M::Rotations, ::typeof(rand), ::Identity{Multiplication return similar(Matrix{Float64}, representation_size(M)...) end +function allocation_promotion_function( + ::GeneralUnitaryMultiplicationGroup{<:Any,ℂ}, + ::typeof(identity_element), + args::Tuple, +) + return complex +end + decorated_manifold(G::GeneralUnitaryMultiplicationGroup) = G.manifold @doc raw""" diff --git a/src/groups/group.jl b/src/groups/group.jl index ed99cf2ab7..9d138b975e 100644 --- a/src/groups/group.jl +++ b/src/groups/group.jl @@ -244,8 +244,9 @@ end @trait_function identity_element!(G::AbstractDecoratorManifold, p) -function allocate_result(G::AbstractDecoratorManifold, ::typeof(identity_element)) - return zeros(representation_size(G)...) +function allocate_result(G::AbstractDecoratorManifold, f::typeof(identity_element)) + apf = allocation_promotion_function(G, f, ()) + return zeros(apf(Float64), representation_size(G)...) end @doc raw""" diff --git a/src/groups/unitary.jl b/src/groups/unitary.jl index d775f401ef..7b7795d90b 100644 --- a/src/groups/unitary.jl +++ b/src/groups/unitary.jl @@ -31,7 +31,7 @@ function Unitary(n, 𝔽::AbstractNumbers=ℂ; parameter::Symbol=:type) end @doc raw""" - exp_lie(G::Unitary{2,ℂ}, X) + exp_lie(G::Unitary{TypeParameter{Tuple{2}},ℂ}, X) Compute the group exponential map on the [`Unitary(2)`](@ref) group, which is diff --git a/src/manifolds/GeneralUnitaryMatrices.jl b/src/manifolds/GeneralUnitaryMatrices.jl index c31bcf84da..5f81197b8d 100644 --- a/src/manifolds/GeneralUnitaryMatrices.jl +++ b/src/manifolds/GeneralUnitaryMatrices.jl @@ -51,7 +51,7 @@ end check_point(M::GeneralUnitaryMatrices, p; kwargs...) Check whether `p` is a valid point on the [`UnitaryMatrices`](@ref) or [`OrthogonalMatrices`] `M`, -i.e. that ``p`` has an determinante of absolute value one +i.e. that ``p`` has a determinant of absolute value one The tolerance for the last test can be set using the `kwargs...`. """ @@ -79,7 +79,7 @@ end check_point(M::Rotations, p; kwargs...) Check whether `p` is a valid point on the [`UnitaryMatrices`](@ref) `M`, -i.e. that ``p`` has an determinante of absolute value one, i.e. that ``p^{\mathrm{H}}p`` +i.e. that ``p`` has a determinant of absolute value one, i.e. that ``p^{\mathrm{H}}p`` The tolerance for the last test can be set using the `kwargs...`. """ @@ -157,7 +157,7 @@ end 4D rotations can be described by two orthogonal planes that are unchanged by the action of the rotation (vectors within a plane rotate only within the -plane). The cosines of the two angles $α,β$ of rotation about these planes may be +plane). The cosines of the two angles ``α,β`` of rotation about these planes may be obtained from the distinct real parts of the eigenvalues of the rotation matrix. This function computes these more efficiently by solving the system @@ -329,18 +329,18 @@ end get_coordinates(M::OrthogonalMatrices, p, X) get_coordinates(M::UnitaryMatrices, p, X) -Extract the unique tangent vector components $X^i$ at point `p` on [`Rotations`](@ref) -$\mathrm{SO}(n)$ from the matrix representation `X` of the tangent +Extract the unique tangent vector components ``X^i`` at point `p` on [`Rotations`](@ref) +``\mathrm{SO}(n)`` from the matrix representation `X` of the tangent vector. -The basis on the Lie algebra $𝔰𝔬(n)$ is chosen such that -for $\mathrm{SO}(2)$, $X^1 = θ = X_{21}$ is the angle of rotation, and -for $\mathrm{SO}(3)$, $(X^1, X^2, X^3) = (X_{32}, X_{13}, X_{21}) = θ u$ is the -angular velocity and axis-angle representation, where $u$ is the unit vector +The basis on the Lie algebra ``𝔰𝔬(n)`` is chosen such that +for ``\mathrm{SO}(2)``, ``X^1 = θ = X_{21}`` is the angle of rotation, and +for ``\mathrm{SO}(3)``, ``(X^1, X^2, X^3) = (X_{32}, X_{13}, X_{21}) = θ u`` is the +angular velocity and axis-angle representation, where ``u`` is the unit vector along the axis of rotation. -For $\mathrm{SO}(n)$ where $n ≥ 4$, the additional elements of $X^i$ are -$X^{j (j - 3)/2 + k + 1} = X_{jk}$, for $j ∈ [4,n], k ∈ [1,j)$. +For ``\mathrm{SO}(n)`` where ``n ≥ 4``, the additional elements of ``X^i`` are +``X^{j (j - 3)/2 + k + 1} = X_{jk}``, for ``j ∈ [4,n], k ∈ [1,j)``. """ get_coordinates(::GeneralUnitaryMatrices{<:Any,ℝ}, ::Any...) function get_coordinates( @@ -476,7 +476,7 @@ end Convert the unique tangent vector components `Xⁱ` at point `p` on [`Rotations`](@ref) or [`OrthogonalMatrices`](@ref) -to the matrix representation $X$ of the tangent vector. See +to the matrix representation ``X`` of the tangent vector. See [`get_coordinates`](@ref get_coordinates(::GeneralUnitaryMatrices, ::Any...)) for the conventions used. """ get_vector(::GeneralUnitaryMatrices{<:Any,ℝ}, ::Any...) @@ -643,7 +643,7 @@ injectivity_radius(::GeneralUnitaryMatrices) = π @doc raw""" injectivity_radius(G::GeneralUnitaryMatrices{<:Any,ℂ,DeterminantOneMatrices}) -Return the injectivity radius for general complex unitary matrix manifolds, where the determinant is $+1$, +Return the injectivity radius for general complex unitary matrix manifolds, where the determinant is ``+1``, which is[^1] ```math @@ -726,7 +726,7 @@ Compute the logarithmic map on the [`Rotations`](@ref) manifold \log_p q = \operatorname{log}(p^{\mathrm{T}}q) ``` -where $\operatorname{Log}$ denotes the matrix logarithm. For numerical stability, +where ``\operatorname{Log}`` denotes the matrix logarithm. For numerical stability, the result is projected onto the set of skew symmetric matrices. For antipodal rotations the function returns deterministically one of the tangent vectors diff --git a/test/groups/general_unitary_groups.jl b/test/groups/general_unitary_groups.jl index 7168d59651..d681c818bb 100644 --- a/test/groups/general_unitary_groups.jl +++ b/test/groups/general_unitary_groups.jl @@ -89,6 +89,8 @@ include("group_utils.jl") @test manifold_volume(Unitary(3)) ≈ sqrt(3) * 2 * π^6 @test manifold_volume(Unitary(4)) ≈ sqrt(2) * 8 * π^10 / 12 + @test identity_element(U2) isa Matrix{ComplexF64} + for n in [1, 2, 3] Un = Unitary(n) X = zeros(ComplexF64, n, n)