diff --git a/NEWS.md b/NEWS.md index 8cdb6e7a..41eb7873 100644 --- a/NEWS.md +++ b/NEWS.md @@ -5,6 +5,12 @@ 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.15.17] 04/10/2024 + +### Changed + +* **Mildly breaking**: the number system parameter now corresponds to the coefficients standing in front of basis vectors in a linear combination instead of components of a vector. For example, `DefaultOrthonormalBasis() == DefaultOrthonormalBasis(ℝ)` of `DefaultManifold(3, field=ℂ)` now has 6 vectors, and `DefaultOrthonormalBasis(ℂ)` of the same manifold has 3 basis vectors. + ## [0.15.16] 13/09/2024 ### Changed diff --git a/Project.toml b/Project.toml index bf7b878d..ad144c1a 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "ManifoldsBase" uuid = "3362f125-f0bb-47a3-aa74-596ffd7ef2fb" authors = ["Seth Axen ", "Mateusz Baran ", "Ronny Bergmann ", "Antoine Levitt "] -version = "0.15.16" +version = "0.15.17" [deps] LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" diff --git a/ext/ManifoldsBaseQuaternionsExt.jl b/ext/ManifoldsBaseQuaternionsExt.jl index 742ec168..be45760d 100644 --- a/ext/ManifoldsBaseQuaternionsExt.jl +++ b/ext/ManifoldsBaseQuaternionsExt.jl @@ -2,13 +2,13 @@ module ManifoldsBaseQuaternionsExt if isdefined(Base, :get_extension) using ManifoldsBase - using ManifoldsBase: ℍ + using ManifoldsBase: ℍ, QuaternionNumbers using Quaternions else # imports need to be relative for Requires.jl-based workflows: # https://github.com/JuliaArrays/ArrayInterface.jl/pull/387 using ..ManifoldsBase - using ..ManifoldsBase: ℍ + using ..ManifoldsBase: ℍ, QuaternionNumbers using ..Quaternions end @@ -20,4 +20,12 @@ end return QuaternionF64 end +@inline function ManifoldsBase.coordinate_eltype( + ::AbstractManifold, + p, + 𝔽::QuaternionNumbers, +) + return quat(number_eltype(p)) +end + end diff --git a/src/DefaultManifold.jl b/src/DefaultManifold.jl index cbaa0029..6bf2f41d 100644 --- a/src/DefaultManifold.jl +++ b/src/DefaultManifold.jl @@ -65,23 +65,22 @@ embed!(::DefaultManifold, Y, p, X) = copyto!(Y, X) exp!(::DefaultManifold, q, p, X) = (q .= p .+ X) -function get_basis_orthonormal(::DefaultManifold, p, N::RealNumbers) - return CachedBasis( - DefaultOrthonormalBasis(N), - [_euclidean_basis_vector(p, i) for i in eachindex(p)], - ) -end -function get_basis_orthogonal(::DefaultManifold, p, N::RealNumbers) - return CachedBasis( - DefaultOrthogonalBasis(N), - [_euclidean_basis_vector(p, i) for i in eachindex(p)], - ) -end -function get_basis_default(::DefaultManifold, p, N::RealNumbers) - return CachedBasis( - DefaultBasis(N), - [_euclidean_basis_vector(p, i) for i in eachindex(p)], +for fname in [:get_basis_orthonormal, :get_basis_orthogonal, :get_basis_default] + BD = Dict( + :get_basis_orthonormal => :DefaultOrthonormalBasis, + :get_basis_orthogonal => :DefaultOrthogonalBasis, + :get_basis_default => :DefaultBasis, ) + BT = BD[fname] + @eval function $fname(::DefaultManifold{ℝ}, p, N::RealNumbers) + return CachedBasis($BT(N), [_euclidean_basis_vector(p, i) for i in eachindex(p)]) + end + @eval function $fname(::DefaultManifold{ℂ}, p, N::ComplexNumbers) + return CachedBasis( + $BT(N), + [_euclidean_basis_vector(p, i, real) for i in eachindex(p)], + ) + end end function get_basis_diagonalizing(M::DefaultManifold, p, B::DiagonalizingOrthonormalBasis) vecs = get_vectors(M, p, get_basis(M, p, DefaultOrthonormalBasis())) @@ -89,10 +88,10 @@ function get_basis_diagonalizing(M::DefaultManifold, p, B::DiagonalizingOrthonor return CachedBasis(B, DiagonalizingBasisData(B.frame_direction, eigenvalues, vecs)) end -# Complex manifold, real basis -> coefficients c are complesx -> reshape +# Complex manifold, real basis -> coefficients c are complex -> reshape # Real manifold, real basis -> reshape -function get_coordinates_orthonormal!(M::DefaultManifold, c, p, X, ::RealNumbers) - return copyto!(c, reshape(X, number_of_coordinates(M, ℝ))) +function get_coordinates_orthonormal!(M::DefaultManifold, c, p, X, N::AbstractNumbers) + return copyto!(c, reshape(X, number_of_coordinates(M, N))) end function get_coordinates_diagonalizing!( M::DefaultManifold, @@ -103,11 +102,11 @@ function get_coordinates_diagonalizing!( ) return copyto!(c, reshape(X, number_of_coordinates(M, ℝ))) end -function get_coordinates_orthonormal!(::DefaultManifold, c, p, X, ::ComplexNumbers) +function get_coordinates_orthonormal!(::DefaultManifold{ℂ}, c, p, X, ::RealNumbers) m = length(X) return copyto!(c, [reshape(real(X), m); reshape(imag(X), m)]) end -function get_vector_orthonormal!(M::DefaultManifold, Y, p, c, ::RealNumbers) +function get_vector_orthonormal!(M::DefaultManifold, Y, p, c, ::AbstractNumbers) return copyto!(Y, reshape(c, representation_size(M))) end function get_vector_diagonalizing!( @@ -119,7 +118,7 @@ function get_vector_diagonalizing!( ) return copyto!(Y, reshape(c, representation_size(M))) end -function get_vector_orthonormal!(M::DefaultManifold{ℂ}, Y, p, c, ::ComplexNumbers) +function get_vector_orthonormal!(M::DefaultManifold{ℂ}, Y, p, c, ::RealNumbers) n = div(length(c), 2) return copyto!(Y, reshape(c[1:n] + c[(n + 1):(2n)] * 1im, representation_size(M))) end diff --git a/src/bases.jl b/src/bases.jl index 53302a8b..aa37118e 100644 --- a/src/bases.jl +++ b/src/bases.jl @@ -48,7 +48,7 @@ Abstract type that represents a basis of vector space of type `VST` on a manifol a subset of it. The type parameter `𝔽` denotes the [`AbstractNumbers`](@ref) that will be used -for the vectors elements. +as coefficients in linear combinations of the basis vectors. # See also @@ -63,7 +63,7 @@ An arbitrary basis of vector space of type `VST` on a manifold. This will usuall be the fastest basis available for a manifold. The type parameter `𝔽` denotes the [`AbstractNumbers`](@ref) that will be used -for the vectors elements. +as coefficients in linear combinations of the basis vectors. # See also @@ -89,7 +89,7 @@ Abstract type that represents an orthonormal basis of vector space of type `VST` manifold or a subset of it. The type parameter `𝔽` denotes the [`AbstractNumbers`](@ref) that will be used -for the vectors elements. +as coefficients in linear combinations of the basis vectors. # See also @@ -104,7 +104,7 @@ An arbitrary orthogonal basis of vector space of type `VST` on a manifold. This be the fastest orthogonal basis available for a manifold. The type parameter `𝔽` denotes the [`AbstractNumbers`](@ref) that will be used -for the vectors elements. +as coefficients in linear combinations of the basis vectors. # See also @@ -137,7 +137,7 @@ Abstract type that represents an orthonormal basis of vector space of type `VST` manifold or a subset of it. The type parameter `𝔽` denotes the [`AbstractNumbers`](@ref) that will be used -for the vectors elements. +as coefficients in linear combinations of the basis vectors. # See also @@ -153,7 +153,7 @@ An arbitrary orthonormal basis of vector space of type `VST` on a manifold. This be the fastest orthonormal basis available for a manifold. The type parameter `𝔽` denotes the [`AbstractNumbers`](@ref) that will be used -for the vectors elements. +as coefficients in linear combinations of the basis vectors. # See also @@ -184,7 +184,7 @@ of the ambient space projected onto the subspace representing the tangent space at a given point. The type parameter `𝔽` denotes the [`AbstractNumbers`](@ref) that will be used -for the vectors elements. +as coefficients in linear combinations of the basis vectors. Available methods: - `:gram_schmidt` uses a modified Gram-Schmidt orthonormalization. @@ -205,6 +205,7 @@ end An orthonormal basis obtained from a basis. # Constructor + GramSchmidtOrthonormalBasis(𝔽::AbstractNumbers = ℝ) """ struct GramSchmidtOrthonormalBasis{𝔽} <: AbstractOrthonormalBasis{𝔽,TangentSpaceType} end @@ -218,7 +219,7 @@ An orthonormal basis `Ξ` as a vector of tangent vectors (of length determined b tensor ``R(u,v)w`` and where the direction `frame_direction` ``v`` has curvature `0`. The type parameter `𝔽` denotes the [`AbstractNumbers`](@ref) that will be used -for the vectors elements. +as coefficients in linear combinations of the basis vectors. # Constructor @@ -302,15 +303,15 @@ function allocate_result( f::typeof(get_coordinates), p, X, - basis::AbstractBasis, -) - T = allocate_result_type(M, f, (p, X)) + basis::AbstractBasis{𝔽}, +) where {𝔽} + T = coordinate_eltype(M, p, 𝔽) return allocate_coordinates(M, p, T, number_of_coordinates(M, basis)) end @inline function allocate_result_type( M::AbstractManifold, - f::Union{typeof(get_coordinates),typeof(get_vector)}, + f::typeof(get_vector), args::Tuple{Any,Vararg{Any}}, ) apf = allocation_promotion_function(M, f, args) @@ -361,18 +362,17 @@ function combine_allocation_promotion_functions(::typeof(identity), ::typeof(com end """ - coordinate_eltype(M::AbstractManifold{M𝔽}, p, 𝔽::AbstractNumbers) where {M𝔽} + coordinate_eltype(M::AbstractManifold, p, 𝔽::AbstractNumbers) Get the element type for 𝔽-field coordinates of the tangent space at a point `p` from manifold `M`. This default assumes that usually complex bases of complex manifolds have real coordinates but it can be overridden by a more specific method. """ -@inline function coordinate_eltype(::AbstractManifold{M𝔽}, p, 𝔽::AbstractNumbers) where {M𝔽} - if M𝔽 === 𝔽 - return real(number_eltype(p)) - else - return number_eltype(p) - end +@inline function coordinate_eltype(::AbstractManifold, p, 𝔽::ComplexNumbers) + return complex(float(number_eltype(p))) +end +@inline function coordinate_eltype(::AbstractManifold, p, ::RealNumbers) + return real(float(number_eltype(p))) end @doc raw""" @@ -406,15 +406,17 @@ function _dual_basis( return DefaultOrthonormalBasis{𝔽}(TangentSpaceType()) end -function _euclidean_basis_vector(p::StridedArray, i) - X = zero(p) +# if `p` has complex eltype but you'd like to have real basis vectors, +# you can pass `real` as a third argument to get that +function _euclidean_basis_vector(p::StridedArray, i, eltype_transform = identity) + X = zeros(eltype_transform(eltype(p)), size(p)...) X[i] = 1 return X end -function _euclidean_basis_vector(p, i) +function _euclidean_basis_vector(p, i, eltype_transform = identity) # when p is for example a SArray - X = similar(p) - copyto!(X, zero(p)) + X = similar(p, eltype_transform(eltype(p))) + fill!(X, zero(eltype(X))) X[i] = 1 return X end @@ -559,8 +561,8 @@ function _get_coordinates(M::AbstractManifold, p, X, B::DefaultOrthonormalBasis) return get_coordinates_orthonormal(M, p, X, number_system(B)) end function get_coordinates_orthonormal(M::AbstractManifold, p, X, N) - Y = allocate_result(M, get_coordinates, p, X, DefaultOrthonormalBasis(N)) - return get_coordinates_orthonormal!(M, Y, p, X, N) + c = allocate_result(M, get_coordinates, p, X, DefaultOrthonormalBasis(N)) + return get_coordinates_orthonormal!(M, c, p, X, N) end function _get_coordinates(M::AbstractManifold, p, X, B::DiagonalizingOrthonormalBasis) @@ -572,8 +574,8 @@ function get_coordinates_diagonalizing( X, B::DiagonalizingOrthonormalBasis, ) - Y = allocate_result(M, get_coordinates, p, X, B) - return get_coordinates_diagonalizing!(M, Y, p, X, B) + c = allocate_result(M, get_coordinates, p, X, B) + return get_coordinates_diagonalizing!(M, c, p, X, B) end function _get_coordinates(M::AbstractManifold, p, X, B::CachedBasis) @@ -585,7 +587,7 @@ function get_coordinates_cached( p, X, B::CachedBasis, - ::RealNumbers, + ::ComplexNumbers, ) return map(vb -> conj(inner(M, p, X, vb)), get_vectors(M, p, B)) end @@ -595,7 +597,7 @@ function get_coordinates_cached( p, X, C::CachedBasis, - ::𝔽, + ::RealNumbers, ) where {𝔽} return map(vb -> real(inner(M, p, X, vb)), get_vectors(M, p, C)) end @@ -644,7 +646,7 @@ function get_coordinates_cached!( p, X, B::CachedBasis, - ::RealNumbers, + ::ComplexNumbers, ) map!(vb -> conj(inner(M, p, X, vb)), Y, get_vectors(M, p, B)) return Y @@ -656,7 +658,7 @@ function get_coordinates_cached!( p, X, C::CachedBasis, - ::𝔽, + ::RealNumbers, ) where {𝔽} map!(vb -> real(inner(M, p, X, vb)), Y, get_vectors(M, p, C)) return Y @@ -939,26 +941,22 @@ For array manifolds, this converts a vector representation of the tangent vector to an array representation. The [`vee`](@ref) map is the `hat` map's inverse. """ -@inline hat(M::AbstractManifold, p, X) = - get_vector(M, p, X, VeeOrthogonalBasis(number_system(M))) -@inline hat!(M::AbstractManifold, Y, p, X) = - get_vector!(M, Y, p, X, VeeOrthogonalBasis(number_system(M))) +@inline hat(M::AbstractManifold, p, X) = get_vector(M, p, X, VeeOrthogonalBasis(ℝ)) +@inline hat!(M::AbstractManifold, Y, p, X) = get_vector!(M, Y, p, X, VeeOrthogonalBasis(ℝ)) """ - number_of_coordinates(M::AbstractManifold{𝔽}, B::AbstractBasis) - number_of_coordinates(M::AbstractManifold{𝔽}, ::𝔾) + number_of_coordinates(M::AbstractManifold, B::AbstractBasis) + number_of_coordinates(M::AbstractManifold, ::𝔾) Compute the number of coordinates in basis of field type `𝔾` on a manifold `M`. This also corresponds to the number of vectors represented by `B`, or stored within `B` in case of a [`CachedBasis`](@ref). """ -function number_of_coordinates(M::AbstractManifold{𝔽}, ::AbstractBasis{𝔾}) where {𝔽,𝔾} +function number_of_coordinates(M::AbstractManifold, ::AbstractBasis{𝔾}) where {𝔾} return number_of_coordinates(M, 𝔾) end -function number_of_coordinates(M::AbstractManifold{𝔽}, f::𝔾) where {𝔽,𝔾} - # for odd manifolds this first case has to match. - (real_dimension(𝔽) == real_dimension(f)) && return manifold_dimension(M) - return div(manifold_dimension(M), real_dimension(𝔽)) * real_dimension(f) +function number_of_coordinates(M::AbstractManifold, f::𝔾) where {𝔾} + return div(manifold_dimension(M), real_dimension(f)) end """ @@ -1084,8 +1082,7 @@ For array manifolds, this converts an array representation of the tangent vector to a vector representation. The [`hat`](@ref) map is the `vee` map's inverse. """ -vee(M::AbstractManifold, p, X) = - get_coordinates(M, p, X, VeeOrthogonalBasis(number_system(M))) +vee(M::AbstractManifold, p, X) = get_coordinates(M, p, X, VeeOrthogonalBasis(ℝ)) function vee!(M::AbstractManifold, Y, p, X) - return get_coordinates!(M, Y, p, X, VeeOrthogonalBasis(number_system(M))) + return get_coordinates!(M, Y, p, X, VeeOrthogonalBasis(ℝ)) end diff --git a/test/ManifoldsBaseTestUtils.jl b/test/ManifoldsBaseTestUtils.jl index fcdc9130..ce9fe2aa 100644 --- a/test/ManifoldsBaseTestUtils.jl +++ b/test/ManifoldsBaseTestUtils.jl @@ -306,7 +306,7 @@ function ManifoldsBase.exp!( end function ManifoldsBase.get_basis_orthonormal( - ::DefaultManifold, + ::DefaultManifold{ℝ}, p::NonBroadcastBasisThing, 𝔽::RealNumbers, ) @@ -319,7 +319,7 @@ function ManifoldsBase.get_basis_orthonormal( ) end function ManifoldsBase.get_basis_orthogonal( - ::DefaultManifold, + ::DefaultManifold{ℝ}, p::NonBroadcastBasisThing, 𝔽::RealNumbers, ) @@ -332,7 +332,7 @@ function ManifoldsBase.get_basis_orthogonal( ) end function ManifoldsBase.get_basis_default( - ::DefaultManifold, + ::DefaultManifold{ℝ}, p::NonBroadcastBasisThing, N::ManifoldsBase.RealNumbers, ) diff --git a/test/allocation.jl b/test/allocation.jl index 3838a20d..fdbfe172 100644 --- a/test/allocation.jl +++ b/test/allocation.jl @@ -95,4 +95,7 @@ ManifoldsBase.representation_size(::AllocManifold4) = (2, 3) @test default_type(M2) === Matrix{Float64} @test default_type(M2, TangentSpaceType()) === Matrix{Float64} + + @test ManifoldsBase.coordinate_eltype(M, fill(true, 2, 3), ℝ) === Float64 + @test ManifoldsBase.coordinate_eltype(AllocManifold4(), a4, ℍ) === QuaternionF64 end diff --git a/test/bases.jl b/test/bases.jl index 49e2de2b..c40cf55e 100644 --- a/test/bases.jl +++ b/test/bases.jl @@ -219,14 +219,15 @@ using ManifoldsBaseTestUtils p = [1.0, 2.0im, 3.0] X = [1.2, 2.2im, 2.3im] b = [Matrix{Float64}(I, 3, 3)[:, i] for i in 1:3] - Bℝ = CachedBasis(DefaultOrthonormalBasis{ℝ}(), b) + bℂ = [b..., (b .* 1im)...] + + Bℝ = CachedBasis(DefaultOrthonormalBasis{ℂ}(), b) aℝ = get_coordinates(M, p, X, Bℝ) Yℝ = get_vector(M, p, aℝ, Bℝ) @test Yℝ ≈ X @test ManifoldsBase.number_of_coordinates(M, Bℝ) == 3 - bℂ = [b..., (b .* 1im)...] - Bℂ = CachedBasis(DefaultOrthonormalBasis{ℂ}(), bℂ) + Bℂ = CachedBasis(DefaultOrthonormalBasis{ℝ}(), bℂ) aℂ = get_coordinates(M, p, X, Bℂ) Yℂ = get_vector(M, p, aℂ, Bℂ) @test Yℂ ≈ X @@ -380,23 +381,23 @@ using ManifoldsBaseTestUtils X = [2.0, 1.0im] Bc = DefaultOrthonormalBasis(ManifoldsBase.ℂ) CBc = get_basis(Mc, p, Bc) - @test CBc.data == [[1.0, 0.0], [0.0, 1.0], [1.0im, 0.0], [0.0, 1.0im]] + @test CBc.data == [[1.0, 0.0], [0.0, 1.0]] B = DefaultOrthonormalBasis(ManifoldsBase.ℝ) CB = get_basis(Mc, p, B) - @test CB.data == [[1.0, 0.0], [0.0, 1.0]] - @test get_coordinates(Mc, p, X, CBc) == [2.0, 0.0, 0.0, 1.0] - @test get_coordinates(Mc, p, X, CB) == [2.0, 1.0im] + @test CB.data == [[1.0, 0.0], [0.0, 1.0], [1.0im, 0.0], [0.0, 1.0im]] + @test get_coordinates(Mc, p, X, CB) == [2.0, 0.0, 0.0, 1.0] + @test get_coordinates(Mc, p, X, CBc) == [2.0, 1.0im] # ONB cc = zeros(4) - @test get_coordinates!(Mc, cc, p, X, Bc) == [2.0, 0.0, 0.0, 1.0] + @test get_coordinates!(Mc, cc, p, X, B) == [2.0, 0.0, 0.0, 1.0] @test cc == [2.0, 0.0, 0.0, 1.0] c = zeros(ComplexF64, 2) - @test get_coordinates!(Mc, c, p, X, B) == [2.0, 1.0im] + @test get_coordinates!(Mc, c, p, X, Bc) == [2.0, 1.0im] @test c == [2.0, 1.0im] # Cached - @test get_coordinates!(Mc, cc, p, X, CBc) == [2.0, 0.0, 0.0, 1.0] + @test get_coordinates!(Mc, cc, p, X, CB) == [2.0, 0.0, 0.0, 1.0] @test cc == [2.0, 0.0, 0.0, 1.0] - @test get_coordinates!(Mc, c, p, X, CB) == [2.0, 1.0im] + @test get_coordinates!(Mc, c, p, X, CBc) == [2.0, 1.0im] @test c == [2.0, 1.0im] end diff --git a/test/default_manifold.jl b/test/default_manifold.jl index ef1c0c86..44bc97d7 100644 --- a/test/default_manifold.jl +++ b/test/default_manifold.jl @@ -636,11 +636,19 @@ using ManifoldsBaseTestUtils MC = ManifoldsBase.DefaultManifold(3; field = ManifoldsBase.ℂ) p = [1.0im, 2.0im, -1.0im] CB = get_basis(MC, p, DefaultOrthonormalBasis(ManifoldsBase.ℂ)) - @test CB.data isa Vector{Vector{ComplexF64}} - @test ManifoldsBase.coordinate_eltype(MC, p, ManifoldsBase.ℂ) === Float64 - @test ManifoldsBase.coordinate_eltype(MC, p, ManifoldsBase.ℝ) === ComplexF64 + @test CB.data == [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]] + @test CB.data isa Vector{Vector{Float64}} + @test ManifoldsBase.coordinate_eltype(MC, p, ManifoldsBase.ℂ) === ComplexF64 + @test ManifoldsBase.coordinate_eltype(MC, p, ManifoldsBase.ℝ) === Float64 CBR = get_basis(MC, p, DefaultOrthonormalBasis()) - @test CBR.data == [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]] + @test CBR.data == [ + [1.0 + 0.0im, 0.0 + 0.0im, 0.0 + 0.0im], + [0.0 + 0.0im, 1.0 + 0.0im, 0.0 + 0.0im], + [0.0 + 0.0im, 0.0 + 0.0im, 1.0 + 0.0im], + [0.0 + 1.0im, 0.0 + 0.0im, 0.0 + 0.0im], + [0.0 + 0.0im, 0.0 + 1.0im, 0.0 + 0.0im], + [0.0 + 0.0im, 0.0 + 0.0im, 0.0 + 1.0im], + ] end @testset "Show methods" begin @@ -671,7 +679,7 @@ using ManifoldsBaseTestUtils end @testset "scalars" begin - M = DefaultManifold() + M = ManifoldsBase.DefaultManifold() p = 1.0 X = 2.0 @test copy(M, p) === p