From e49cf2a1a7ea22cad800bab509cf03a2c738f78e Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Mon, 15 Jul 2024 12:05:07 +0200 Subject: [PATCH 01/13] Typed instance-free allocation --- NEWS.md | 6 ++++ Project.toml | 2 +- src/ManifoldsBase.jl | 65 +++++++++++++++++++++++++++++++++++++++++--- src/bases.jl | 7 +++++ test/allocation.jl | 19 +++++++++++-- 5 files changed, 92 insertions(+), 7 deletions(-) diff --git a/NEWS.md b/NEWS.md index b457a004..caf32c89 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.11] 16/07/2024 + +### Added + +* Function `allocate_as` to generically allocate point and tangent vectors on a manifold without a pre-existing instance but of a particular type. + ## [0.15.10] 19/05/2024 ### Added diff --git a/Project.toml b/Project.toml index a2f5afa3..2aa09697 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.10" +version = "0.15.11" [deps] LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" diff --git a/src/ManifoldsBase.jl b/src/ManifoldsBase.jl index 43b30f4a..c2ceff1e 100644 --- a/src/ManifoldsBase.jl +++ b/src/ManifoldsBase.jl @@ -59,7 +59,7 @@ Type `T` is the new number element type [`number_eltype`](@ref), if it is not gi the element type of `a` is retained. The `dims` argument can be given for non-nested allocation and is forwarded to the function `similar`. -It's behavior can be overriden by a specific manifold, for example power manifold with +It's behavior can be overridden by a specific manifold, for example power manifold with nested replacing representation can decide that `allocate` for `Array{<:SArray}` returns another `Array{<:SArray}` instead of `Array{<:MArray}`, as would be done by default. """ @@ -86,6 +86,62 @@ function allocate(::AbstractManifold, a, T::Type, dim1::Integer, dims::Integer.. end allocate(::AbstractManifold, a, T::Type, dims::Tuple) = allocate(a, T, dims) +""" + allocate_as(M::AbstractManifold, [T:::Type]) + allocate_as(M::AbstractManifold, F::FiberType, [T:::Type]) + +Allocate a new point on manifold `M` with optional type given by `T`. Note that `T` is not +number element type as in [`allocate`](@ref) but rather the type of the entire point to be +returned. + +If `F` is provided, then an element of the corresponding fiber is allocated, assuming it is +independent of the base point. + +To allocate a tangent vector, use `` + +# Example + +```julia-repl +julia> using ManifoldsBase + +julia> M = ManifoldsBase.DefaultManifold(4) +DefaultManifold(4; field = ℝ) + +julia> allocate_as(M) +4-element Vector{Float64}: + 0.0 + 0.0 + 0.0 + 0.0 + +julia> allocate_as(M, Array{Float64}) +4-element Vector{Float64}: + 0.0 + 0.0 + 0.0 + 0.0 + +julia> allocate_as(M, TangentSpaceType()) +4-element Vector{Float64}: + 0.0 + 0.0 + 0.0 + 0.0 + +julia> allocate_as(M, TangentSpaceType(), Array{Float64}) +4-element Vector{Float64}: + 0.0 + 0.0 + 0.0 + 0.0 + +``` +""" +allocate_as(M::AbstractManifold) = similar(Array{Float64}, representation_size(M)) +function allocate_as(M::AbstractManifold, T::Type{<:AbstractArray}) + return similar(T, representation_size(M)) +end + """ _pick_basic_allocation_argument(::AbstractManifold, f, x...) @@ -477,7 +533,7 @@ embed!(M::AbstractManifold, Y, p, X) = copyto!(M, Y, p, X) Embed `p` from manifold `M` an project it back to `M`. For points from `M` this is identity but in case embedding is defined for points outside of `M`, this can serve as a way -to for example remove numerical innacuracies caused by some algorithms. +to for example remove numerical inaccuracies caused by some algorithms. """ function embed_project(M::AbstractManifold, p) return project(M, embed(M, p)) @@ -487,8 +543,8 @@ end Embed vector `X` tangent at `p` from manifold `M` an project it back to tangent space at `p`. For points from that tangent space this is identity but in case embedding is -defined for tagent vectors from outside of it, this can serve as a way to for example remove -numerical innacuracies caused by some algorithms. +defined for tangent vectors from outside of it, this can serve as a way to for example remove +numerical inaccuracies caused by some algorithms. """ function embed_project(M::AbstractManifold, p, X) return project(M, p, embed(M, p, X)) @@ -1219,6 +1275,7 @@ export ×, ℝ, ℂ, allocate, + allocate_as, angle, base_manifold, base_point, diff --git a/src/bases.jl b/src/bases.jl index 75800df3..eb237129 100644 --- a/src/bases.jl +++ b/src/bases.jl @@ -280,6 +280,13 @@ const all_uncached_bases{T} = Union{ DefaultOrthonormalBasis{<:Any,T}, } +function allocate_as(M::AbstractManifold, ::TangentSpaceType) + return similar(Array{Float64}, representation_size(M)) +end +function allocate_as(M::AbstractManifold, ::TangentSpaceType, T::Type{<:AbstractArray}) + return similar(T, representation_size(M)) +end + """ allocate_coordinates(M::AbstractManifold, p, T, n::Int) diff --git a/test/allocation.jl b/test/allocation.jl index d9b8e820..eccee8f8 100644 --- a/test/allocation.jl +++ b/test/allocation.jl @@ -67,9 +67,24 @@ ManifoldsBase.representation_size(::AllocManifold3) = (2, 3) @test number_eltype(Any[[2.0], [3.0]]) === Float64 @test number_eltype(typeof([[1.0, 2.0]])) === Float64 - alloc2 = ManifoldsBase.allocate_result(AllocManifold2(), rand) + M2 = AllocManifold2() + alloc2 = ManifoldsBase.allocate_result(M2, rand) @test alloc2 isa Matrix{Float64} - @test size(alloc2) == representation_size(AllocManifold2()) + @test size(alloc2) == representation_size(M2) @test ManifoldsBase.allocate_result(AllocManifold3(), rand) isa Matrix{ComplexF64} + + an = allocate_as(M2) + @test an isa Matrix{Float64} + @test size(an) == representation_size(M2) + an = allocate_as(M2, Array{Float32}) + @test an isa Matrix{Float32} + @test size(an) == representation_size(M2) + + an = allocate_as(M2, TangentSpaceType()) + @test an isa Matrix{Float64} + @test size(an) == representation_size(M2) + an = allocate_as(M2, TangentSpaceType(), Array{Float32}) + @test an isa Matrix{Float32} + @test size(an) == representation_size(M2) end From ce14cf4bb59a8a312caed19d8b8388c0a9365c54 Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Mon, 15 Jul 2024 12:32:09 +0200 Subject: [PATCH 02/13] specialize for ProductManifold --- .../ManifoldsBaseRecursiveArrayToolsExt.jl | 2 ++ .../ProductManifoldRecursiveArrayToolsExt.jl | 17 +++++++++++++++++ test/product_manifold.jl | 19 +++++++++++++++++++ 3 files changed, 38 insertions(+) diff --git a/ext/ManifoldsBaseRecursiveArrayToolsExt/ManifoldsBaseRecursiveArrayToolsExt.jl b/ext/ManifoldsBaseRecursiveArrayToolsExt/ManifoldsBaseRecursiveArrayToolsExt.jl index 6e4dc7dc..89e327f5 100644 --- a/ext/ManifoldsBaseRecursiveArrayToolsExt/ManifoldsBaseRecursiveArrayToolsExt.jl +++ b/ext/ManifoldsBaseRecursiveArrayToolsExt/ManifoldsBaseRecursiveArrayToolsExt.jl @@ -9,6 +9,7 @@ if isdefined(Base, :get_extension) import ManifoldsBase: allocate, + allocate_as, allocate_result, default_inverse_retraction_method, default_retraction_method, @@ -40,6 +41,7 @@ else import ..ManifoldsBase: allocate, + allocate_as, allocate_result, default_inverse_retraction_method, default_retraction_method, diff --git a/ext/ManifoldsBaseRecursiveArrayToolsExt/ProductManifoldRecursiveArrayToolsExt.jl b/ext/ManifoldsBaseRecursiveArrayToolsExt/ProductManifoldRecursiveArrayToolsExt.jl index e519fd4c..914aee6e 100644 --- a/ext/ManifoldsBaseRecursiveArrayToolsExt/ProductManifoldRecursiveArrayToolsExt.jl +++ b/ext/ManifoldsBaseRecursiveArrayToolsExt/ProductManifoldRecursiveArrayToolsExt.jl @@ -5,6 +5,23 @@ function allocate(x::ArrayPartition, T::Type) return ArrayPartition(map(t -> allocate(t, T), submanifold_components(x))...) end +allocate_as(M::ProductManifold) = ArrayPartition(map(N -> allocate_as(N), M.manifolds)...) +function allocate_as(M::ProductManifold, ::Type{ArrayPartition{T,U}}) where {T,U} + return ArrayPartition(map((N, V) -> allocate_as(N, V), M.manifolds, U.parameters)...) +end +function allocate_as(M::ProductManifold, ft::TangentSpaceType) + return ArrayPartition(map(N -> allocate_as(N, ft), M.manifolds)...) +end +function allocate_as( + M::ProductManifold, + ft::TangentSpaceType, + ::Type{ArrayPartition{T,U}}, +) where {T,U} + return ArrayPartition( + map((N, V) -> allocate_as(N, ft, V), M.manifolds, U.parameters)..., + ) +end + @inline function allocate_result(M::ProductManifold, f) return ArrayPartition(map(N -> allocate_result(N, f), M.manifolds)) end diff --git a/test/product_manifold.jl b/test/product_manifold.jl index 67c84266..8030d7b4 100644 --- a/test/product_manifold.jl +++ b/test/product_manifold.jl @@ -1,3 +1,4 @@ +using Revise using Test using ManifoldsBase using ManifoldsBase: DefaultManifold, submanifold_component, submanifold_components @@ -132,6 +133,24 @@ using ManifoldsBaseTestUtils @test q[1].x[1] isa Vector end + @testset "allocate_as" begin + p1 = allocate_as(M) + @test p1 isa ArrayPartition{Float64,Tuple{Vector{Float64},Matrix{Float64}}} + + p1 = allocate_as(M, ArrayPartition{Float32,Tuple{Vector{Float32},Matrix{Float32}}}) + @test p1 isa ArrayPartition{Float32,Tuple{Vector{Float32},Matrix{Float32}}} + + X1 = allocate_as(M, TangentSpaceType()) + @test X1 isa ArrayPartition{Float64,Tuple{Vector{Float64},Matrix{Float64}}} + + X1 = allocate_as( + M, + TangentSpaceType(), + ArrayPartition{Float32,Tuple{Vector{Float32},Matrix{Float32}}}, + ) + @test X1 isa ArrayPartition{Float32,Tuple{Vector{Float32},Matrix{Float32}}} + end + p1 = ArrayPartition([1.0, 0.0, 0.0], [4.0 5.0; 6.0 7.0]) p2 = ArrayPartition([0.0, 1.0, 0.0], [4.0 8.0; 3.0 7.5]) X1 = ArrayPartition([0.0, 1.0, 0.2], [4.0 0.0; 2.0 7.0]) From e62fb72b25a1ebdb6974c9b2c9aa6787cf445e6f Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Mon, 15 Jul 2024 12:46:08 +0200 Subject: [PATCH 03/13] specialize for nested PowerManifold --- src/PowerManifold.jl | 48 +++++++++++++++++++++++++++++++++++++++- test/power.jl | 20 +++++++++++++++++ test/product_manifold.jl | 1 - 3 files changed, 67 insertions(+), 2 deletions(-) diff --git a/src/PowerManifold.jl b/src/PowerManifold.jl index c188c713..95ca9e01 100644 --- a/src/PowerManifold.jl +++ b/src/PowerManifold.jl @@ -172,6 +172,52 @@ function Base.:^( return PowerManifold(M, size...) end +function allocate_as( + M::PowerManifold{ + 𝔽, + TM, + TSize, + <:Union{NestedPowerRepresentation,NestedReplacingPowerRepresentation}, + }, +) where {𝔽,TM<:AbstractManifold{𝔽},TSize} + return [allocate_as(M.manifold) for _ in get_iterator(M)] +end +function allocate_as( + M::PowerManifold{ + 𝔽, + TM, + TSize, + <:Union{NestedPowerRepresentation,NestedReplacingPowerRepresentation}, + }, + ::Type{<:Array{U}}, +) where {𝔽,TM<:AbstractManifold{𝔽},TSize,U} + return [allocate_as(M.manifold, U) for _ in get_iterator(M)] +end + +function allocate_as( + M::PowerManifold{ + 𝔽, + TM, + TSize, + <:Union{NestedPowerRepresentation,NestedReplacingPowerRepresentation}, + }, + ft::TangentSpaceType, +) where {𝔽,TM<:AbstractManifold{𝔽},TSize} + return [allocate_as(M.manifold, ft) for _ in get_iterator(M)] +end +function allocate_as( + M::PowerManifold{ + 𝔽, + TM, + TSize, + <:Union{NestedPowerRepresentation,NestedReplacingPowerRepresentation}, + }, + ft::TangentSpaceType, + ::Type{<:Array{U}}, +) where {𝔽,TM<:AbstractManifold{𝔽},TSize,U} + return [allocate_as(M.manifold, ft, U) for _ in get_iterator(M)] +end + """ _allocate_access_nested(M::PowerManifoldNested, y, i) @@ -194,7 +240,7 @@ function allocate_result(M::PowerManifoldNested, f, x...) ] end end -# avoid ambituities - though usually not used +# avoid ambiguities - though usually not used function allocate_result( M::PowerManifoldNested, f::typeof(get_coordinates), diff --git a/test/power.jl b/test/power.jl index 390efa72..9d8359db 100644 --- a/test/power.jl +++ b/test/power.jl @@ -90,6 +90,26 @@ end allocate(M, p) isa Vector{SMatrix{2,2,Float64,4}} end + @testset "allocate_as" begin + M = ManifoldsBase.DefaultManifold(2, 2) + N = PowerManifold(M, NestedReplacingPowerRepresentation(), 2) + p = allocate_as(N) + @test p isa Vector{Matrix{Float64}} + @test size(p) == (2,) + + p = allocate_as(N, Vector{Matrix{Float32}}) + @test p isa Vector{Matrix{Float32}} + @test size(p) == (2,) + + X = allocate_as(N, TangentSpaceType()) + @test X isa Vector{Matrix{Float64}} + @test size(X) == (2,) + + X = allocate_as(N, TangentSpaceType(), Vector{Matrix{Float32}}) + @test X isa Vector{Matrix{Float32}} + @test size(X) == (2,) + end + for PowerRepr in [NestedPowerRepresentation, NestedReplacingPowerRepresentation] @testset "PowerManifold with $(PowerRepr)" begin M = ManifoldsBase.DefaultManifold(3) diff --git a/test/product_manifold.jl b/test/product_manifold.jl index 8030d7b4..d1dbe785 100644 --- a/test/product_manifold.jl +++ b/test/product_manifold.jl @@ -1,4 +1,3 @@ -using Revise using Test using ManifoldsBase using ManifoldsBase: DefaultManifold, submanifold_component, submanifold_components From 67d5511d0fff01cae40f127ab11aa7bd31885f5d Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Thu, 18 Jul 2024 09:36:24 +0200 Subject: [PATCH 04/13] add default_type --- NEWS.md | 1 + src/ManifoldsBase.jl | 16 ++++++++++++++++ test/allocation.jl | 2 ++ 3 files changed, 19 insertions(+) diff --git a/NEWS.md b/NEWS.md index caf32c89..c95b783d 100644 --- a/NEWS.md +++ b/NEWS.md @@ -10,6 +10,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added * Function `allocate_as` to generically allocate point and tangent vectors on a manifold without a pre-existing instance but of a particular type. +* Function `default_type` to get the default type of points and tangent vectors for a manifold. ## [0.15.10] 19/05/2024 diff --git a/src/ManifoldsBase.jl b/src/ManifoldsBase.jl index c2ceff1e..a14b5fff 100644 --- a/src/ManifoldsBase.jl +++ b/src/ManifoldsBase.jl @@ -407,6 +407,21 @@ function at the level, where also information from `p` and `M` can be accessed. """ copyto!(::AbstractManifold, Y, p, X) = copyto!(Y, X) +""" + default_type(M::AbstractManifold) + +Get the default type of points on manifold `M`. +""" +default_type(M::AbstractManifold) = typeof(allocate_as(M)) +""" + default_type(M::AbstractManifold, ft::FiberType) + +Get the default type of points from the fiber `ft` of the fiber bundle based on manifold `M`. +For example, call `default_type(MyManifold(), TangentSpaceType())` to get the default type +of a tangent vector. +""" +default_type(M::AbstractManifold, ft::FiberType) = typeof(allocate_as(M, ft)) + @doc raw""" distance(M::AbstractManifold, p, q) @@ -1293,6 +1308,7 @@ export ×, default_approximation_method, default_inverse_retraction_method, default_retraction_method, + default_type, default_vector_transport_method, distance, exp, diff --git a/test/allocation.jl b/test/allocation.jl index eccee8f8..e7aac766 100644 --- a/test/allocation.jl +++ b/test/allocation.jl @@ -87,4 +87,6 @@ ManifoldsBase.representation_size(::AllocManifold3) = (2, 3) an = allocate_as(M2, TangentSpaceType(), Array{Float32}) @test an isa Matrix{Float32} @test size(an) == representation_size(M2) + + @test default_type(M2) === Matrix{Float64} end From 9037377968ab98a32cd6a212ef693603de07d50c Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Thu, 18 Jul 2024 16:58:17 +0200 Subject: [PATCH 05/13] add a test --- test/allocation.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/test/allocation.jl b/test/allocation.jl index e7aac766..9b6501b6 100644 --- a/test/allocation.jl +++ b/test/allocation.jl @@ -89,4 +89,5 @@ ManifoldsBase.representation_size(::AllocManifold3) = (2, 3) @test size(an) == representation_size(M2) @test default_type(M2) === Matrix{Float64} + @test default_type(M2, TangentSpaceType()) === Matrix{Float64} end From c967487fe0a3b01b98588bfc1c71a02a5fac5adc Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Tue, 23 Jul 2024 14:57:19 +0200 Subject: [PATCH 06/13] rename allocate_as to allocate_on --- NEWS.md | 4 ++-- .../ManifoldsBaseRecursiveArrayToolsExt.jl | 4 ++-- .../ProductManifoldRecursiveArrayToolsExt.jl | 14 ++++++------ src/ManifoldsBase.jl | 22 +++++++++---------- src/PowerManifold.jl | 16 +++++++------- src/bases.jl | 4 ++-- test/allocation.jl | 8 +++---- test/power.jl | 10 ++++----- test/product_manifold.jl | 10 ++++----- 9 files changed, 46 insertions(+), 46 deletions(-) diff --git a/NEWS.md b/NEWS.md index c95b783d..fa2b5aae 100644 --- a/NEWS.md +++ b/NEWS.md @@ -5,11 +5,11 @@ 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.11] 16/07/2024 +## [0.15.11] 23/07/2024 ### Added -* Function `allocate_as` to generically allocate point and tangent vectors on a manifold without a pre-existing instance but of a particular type. +* Function `allocate_on` to generically allocate point and tangent vectors on a manifold without a pre-existing instance but of a particular type. * Function `default_type` to get the default type of points and tangent vectors for a manifold. ## [0.15.10] 19/05/2024 diff --git a/ext/ManifoldsBaseRecursiveArrayToolsExt/ManifoldsBaseRecursiveArrayToolsExt.jl b/ext/ManifoldsBaseRecursiveArrayToolsExt/ManifoldsBaseRecursiveArrayToolsExt.jl index 89e327f5..12cd73f8 100644 --- a/ext/ManifoldsBaseRecursiveArrayToolsExt/ManifoldsBaseRecursiveArrayToolsExt.jl +++ b/ext/ManifoldsBaseRecursiveArrayToolsExt/ManifoldsBaseRecursiveArrayToolsExt.jl @@ -9,7 +9,7 @@ if isdefined(Base, :get_extension) import ManifoldsBase: allocate, - allocate_as, + allocate_on, allocate_result, default_inverse_retraction_method, default_retraction_method, @@ -41,7 +41,7 @@ else import ..ManifoldsBase: allocate, - allocate_as, + allocate_on, allocate_result, default_inverse_retraction_method, default_retraction_method, diff --git a/ext/ManifoldsBaseRecursiveArrayToolsExt/ProductManifoldRecursiveArrayToolsExt.jl b/ext/ManifoldsBaseRecursiveArrayToolsExt/ProductManifoldRecursiveArrayToolsExt.jl index 914aee6e..86ba2956 100644 --- a/ext/ManifoldsBaseRecursiveArrayToolsExt/ProductManifoldRecursiveArrayToolsExt.jl +++ b/ext/ManifoldsBaseRecursiveArrayToolsExt/ProductManifoldRecursiveArrayToolsExt.jl @@ -5,20 +5,20 @@ function allocate(x::ArrayPartition, T::Type) return ArrayPartition(map(t -> allocate(t, T), submanifold_components(x))...) end -allocate_as(M::ProductManifold) = ArrayPartition(map(N -> allocate_as(N), M.manifolds)...) -function allocate_as(M::ProductManifold, ::Type{ArrayPartition{T,U}}) where {T,U} - return ArrayPartition(map((N, V) -> allocate_as(N, V), M.manifolds, U.parameters)...) +allocate_on(M::ProductManifold) = ArrayPartition(map(N -> allocate_on(N), M.manifolds)...) +function allocate_on(M::ProductManifold, ::Type{ArrayPartition{T,U}}) where {T,U} + return ArrayPartition(map((N, V) -> allocate_on(N, V), M.manifolds, U.parameters)...) end -function allocate_as(M::ProductManifold, ft::TangentSpaceType) - return ArrayPartition(map(N -> allocate_as(N, ft), M.manifolds)...) +function allocate_on(M::ProductManifold, ft::TangentSpaceType) + return ArrayPartition(map(N -> allocate_on(N, ft), M.manifolds)...) end -function allocate_as( +function allocate_on( M::ProductManifold, ft::TangentSpaceType, ::Type{ArrayPartition{T,U}}, ) where {T,U} return ArrayPartition( - map((N, V) -> allocate_as(N, ft, V), M.manifolds, U.parameters)..., + map((N, V) -> allocate_on(N, ft, V), M.manifolds, U.parameters)..., ) end diff --git a/src/ManifoldsBase.jl b/src/ManifoldsBase.jl index a14b5fff..59e4c093 100644 --- a/src/ManifoldsBase.jl +++ b/src/ManifoldsBase.jl @@ -87,8 +87,8 @@ end allocate(::AbstractManifold, a, T::Type, dims::Tuple) = allocate(a, T, dims) """ - allocate_as(M::AbstractManifold, [T:::Type]) - allocate_as(M::AbstractManifold, F::FiberType, [T:::Type]) + allocate_on(M::AbstractManifold, [T:::Type]) + allocate_on(M::AbstractManifold, F::FiberType, [T:::Type]) Allocate a new point on manifold `M` with optional type given by `T`. Note that `T` is not number element type as in [`allocate`](@ref) but rather the type of the entire point to be @@ -107,28 +107,28 @@ julia> using ManifoldsBase julia> M = ManifoldsBase.DefaultManifold(4) DefaultManifold(4; field = ℝ) -julia> allocate_as(M) +julia> allocate_on(M) 4-element Vector{Float64}: 0.0 0.0 0.0 0.0 -julia> allocate_as(M, Array{Float64}) +julia> allocate_on(M, Array{Float64}) 4-element Vector{Float64}: 0.0 0.0 0.0 0.0 -julia> allocate_as(M, TangentSpaceType()) +julia> allocate_on(M, TangentSpaceType()) 4-element Vector{Float64}: 0.0 0.0 0.0 0.0 -julia> allocate_as(M, TangentSpaceType(), Array{Float64}) +julia> allocate_on(M, TangentSpaceType(), Array{Float64}) 4-element Vector{Float64}: 0.0 0.0 @@ -137,8 +137,8 @@ julia> allocate_as(M, TangentSpaceType(), Array{Float64}) ``` """ -allocate_as(M::AbstractManifold) = similar(Array{Float64}, representation_size(M)) -function allocate_as(M::AbstractManifold, T::Type{<:AbstractArray}) +allocate_on(M::AbstractManifold) = similar(Array{Float64}, representation_size(M)) +function allocate_on(M::AbstractManifold, T::Type{<:AbstractArray}) return similar(T, representation_size(M)) end @@ -412,7 +412,7 @@ copyto!(::AbstractManifold, Y, p, X) = copyto!(Y, X) Get the default type of points on manifold `M`. """ -default_type(M::AbstractManifold) = typeof(allocate_as(M)) +default_type(M::AbstractManifold) = typeof(allocate_on(M)) """ default_type(M::AbstractManifold, ft::FiberType) @@ -420,7 +420,7 @@ Get the default type of points from the fiber `ft` of the fiber bundle based on For example, call `default_type(MyManifold(), TangentSpaceType())` to get the default type of a tangent vector. """ -default_type(M::AbstractManifold, ft::FiberType) = typeof(allocate_as(M, ft)) +default_type(M::AbstractManifold, ft::FiberType) = typeof(allocate_on(M, ft)) @doc raw""" distance(M::AbstractManifold, p, q) @@ -1290,7 +1290,7 @@ export ×, ℝ, ℂ, allocate, - allocate_as, + allocate_on, angle, base_manifold, base_point, diff --git a/src/PowerManifold.jl b/src/PowerManifold.jl index 95ca9e01..037b9032 100644 --- a/src/PowerManifold.jl +++ b/src/PowerManifold.jl @@ -172,7 +172,7 @@ function Base.:^( return PowerManifold(M, size...) end -function allocate_as( +function allocate_on( M::PowerManifold{ 𝔽, TM, @@ -180,9 +180,9 @@ function allocate_as( <:Union{NestedPowerRepresentation,NestedReplacingPowerRepresentation}, }, ) where {𝔽,TM<:AbstractManifold{𝔽},TSize} - return [allocate_as(M.manifold) for _ in get_iterator(M)] + return [allocate_on(M.manifold) for _ in get_iterator(M)] end -function allocate_as( +function allocate_on( M::PowerManifold{ 𝔽, TM, @@ -191,10 +191,10 @@ function allocate_as( }, ::Type{<:Array{U}}, ) where {𝔽,TM<:AbstractManifold{𝔽},TSize,U} - return [allocate_as(M.manifold, U) for _ in get_iterator(M)] + return [allocate_on(M.manifold, U) for _ in get_iterator(M)] end -function allocate_as( +function allocate_on( M::PowerManifold{ 𝔽, TM, @@ -203,9 +203,9 @@ function allocate_as( }, ft::TangentSpaceType, ) where {𝔽,TM<:AbstractManifold{𝔽},TSize} - return [allocate_as(M.manifold, ft) for _ in get_iterator(M)] + return [allocate_on(M.manifold, ft) for _ in get_iterator(M)] end -function allocate_as( +function allocate_on( M::PowerManifold{ 𝔽, TM, @@ -215,7 +215,7 @@ function allocate_as( ft::TangentSpaceType, ::Type{<:Array{U}}, ) where {𝔽,TM<:AbstractManifold{𝔽},TSize,U} - return [allocate_as(M.manifold, ft, U) for _ in get_iterator(M)] + return [allocate_on(M.manifold, ft, U) for _ in get_iterator(M)] end """ diff --git a/src/bases.jl b/src/bases.jl index eb237129..53302a8b 100644 --- a/src/bases.jl +++ b/src/bases.jl @@ -280,10 +280,10 @@ const all_uncached_bases{T} = Union{ DefaultOrthonormalBasis{<:Any,T}, } -function allocate_as(M::AbstractManifold, ::TangentSpaceType) +function allocate_on(M::AbstractManifold, ::TangentSpaceType) return similar(Array{Float64}, representation_size(M)) end -function allocate_as(M::AbstractManifold, ::TangentSpaceType, T::Type{<:AbstractArray}) +function allocate_on(M::AbstractManifold, ::TangentSpaceType, T::Type{<:AbstractArray}) return similar(T, representation_size(M)) end diff --git a/test/allocation.jl b/test/allocation.jl index 9b6501b6..b99983ef 100644 --- a/test/allocation.jl +++ b/test/allocation.jl @@ -74,17 +74,17 @@ ManifoldsBase.representation_size(::AllocManifold3) = (2, 3) @test ManifoldsBase.allocate_result(AllocManifold3(), rand) isa Matrix{ComplexF64} - an = allocate_as(M2) + an = allocate_on(M2) @test an isa Matrix{Float64} @test size(an) == representation_size(M2) - an = allocate_as(M2, Array{Float32}) + an = allocate_on(M2, Array{Float32}) @test an isa Matrix{Float32} @test size(an) == representation_size(M2) - an = allocate_as(M2, TangentSpaceType()) + an = allocate_on(M2, TangentSpaceType()) @test an isa Matrix{Float64} @test size(an) == representation_size(M2) - an = allocate_as(M2, TangentSpaceType(), Array{Float32}) + an = allocate_on(M2, TangentSpaceType(), Array{Float32}) @test an isa Matrix{Float32} @test size(an) == representation_size(M2) diff --git a/test/power.jl b/test/power.jl index 9d8359db..c21f644f 100644 --- a/test/power.jl +++ b/test/power.jl @@ -90,22 +90,22 @@ end allocate(M, p) isa Vector{SMatrix{2,2,Float64,4}} end - @testset "allocate_as" begin + @testset "allocate_on" begin M = ManifoldsBase.DefaultManifold(2, 2) N = PowerManifold(M, NestedReplacingPowerRepresentation(), 2) - p = allocate_as(N) + p = allocate_on(N) @test p isa Vector{Matrix{Float64}} @test size(p) == (2,) - p = allocate_as(N, Vector{Matrix{Float32}}) + p = allocate_on(N, Vector{Matrix{Float32}}) @test p isa Vector{Matrix{Float32}} @test size(p) == (2,) - X = allocate_as(N, TangentSpaceType()) + X = allocate_on(N, TangentSpaceType()) @test X isa Vector{Matrix{Float64}} @test size(X) == (2,) - X = allocate_as(N, TangentSpaceType(), Vector{Matrix{Float32}}) + X = allocate_on(N, TangentSpaceType(), Vector{Matrix{Float32}}) @test X isa Vector{Matrix{Float32}} @test size(X) == (2,) end diff --git a/test/product_manifold.jl b/test/product_manifold.jl index d1dbe785..8635a110 100644 --- a/test/product_manifold.jl +++ b/test/product_manifold.jl @@ -132,17 +132,17 @@ using ManifoldsBaseTestUtils @test q[1].x[1] isa Vector end - @testset "allocate_as" begin - p1 = allocate_as(M) + @testset "allocate_on" begin + p1 = allocate_on(M) @test p1 isa ArrayPartition{Float64,Tuple{Vector{Float64},Matrix{Float64}}} - p1 = allocate_as(M, ArrayPartition{Float32,Tuple{Vector{Float32},Matrix{Float32}}}) + p1 = allocate_on(M, ArrayPartition{Float32,Tuple{Vector{Float32},Matrix{Float32}}}) @test p1 isa ArrayPartition{Float32,Tuple{Vector{Float32},Matrix{Float32}}} - X1 = allocate_as(M, TangentSpaceType()) + X1 = allocate_on(M, TangentSpaceType()) @test X1 isa ArrayPartition{Float64,Tuple{Vector{Float64},Matrix{Float64}}} - X1 = allocate_as( + X1 = allocate_on( M, TangentSpaceType(), ArrayPartition{Float32,Tuple{Vector{Float32},Matrix{Float32}}}, From b248f495530cd57dd370547bf87edc4737432c29 Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Sun, 28 Jul 2024 12:50:23 +0200 Subject: [PATCH 07/13] improved quaternion allocation --- NEWS.md | 5 +++++ Project.toml | 3 +++ ext/ManifoldsBaseQuaternionsExt.jl | 23 +++++++++++++++++++++++ src/ManifoldsBase.jl | 4 ++-- test/allocation.jl | 7 ++++++- 5 files changed, 39 insertions(+), 3 deletions(-) create mode 100644 ext/ManifoldsBaseQuaternionsExt.jl diff --git a/NEWS.md b/NEWS.md index fa2b5aae..abd7a980 100644 --- a/NEWS.md +++ b/NEWS.md @@ -11,6 +11,11 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 * Function `allocate_on` to generically allocate point and tangent vectors on a manifold without a pre-existing instance but of a particular type. * Function `default_type` to get the default type of points and tangent vectors for a manifold. +* Package extension for the `Quaternions.jl` package that handles allocation. + +### Changed + +* Default allocation method was made more robust to custom promotion functions. ## [0.15.10] 19/05/2024 diff --git a/Project.toml b/Project.toml index 2aa09697..adb92a5a 100644 --- a/Project.toml +++ b/Project.toml @@ -12,11 +12,13 @@ Requires = "ae029012-a4dd-5104-9daa-d747884805df" [weakdeps] Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" +Quaternions = "94ee1d12-ae83-5a48-8b1c-48b8ff168ae0" RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" [extensions] ManifoldsBasePlotsExt = "Plots" +ManifoldsBaseQuaternionsExt = "Quaternions" ManifoldsBaseRecursiveArrayToolsExt = "RecursiveArrayTools" ManifoldsBaseStatisticsExt = "Statistics" @@ -35,6 +37,7 @@ Requires = "1" ReverseDiff = "1" StaticArrays = "1" Statistics = "1.6" +Quaternions = "0.7" Test = "1.6" julia = "1.6" diff --git a/ext/ManifoldsBaseQuaternionsExt.jl b/ext/ManifoldsBaseQuaternionsExt.jl new file mode 100644 index 00000000..742ec168 --- /dev/null +++ b/ext/ManifoldsBaseQuaternionsExt.jl @@ -0,0 +1,23 @@ +module ManifoldsBaseQuaternionsExt + +if isdefined(Base, :get_extension) + using ManifoldsBase + using ManifoldsBase: ℍ + 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 ..Quaternions +end + +@inline function ManifoldsBase.allocate_result_type( + ::AbstractManifold{ℍ}, + f::TF, + args::Tuple{}, +) where {TF} + return QuaternionF64 +end + +end diff --git a/src/ManifoldsBase.jl b/src/ManifoldsBase.jl index 59e4c093..6e6a4b59 100644 --- a/src/ManifoldsBase.jl +++ b/src/ManifoldsBase.jl @@ -186,8 +186,8 @@ Return type of element of the array that will represent the result of function ` @inline eti_to_one(eti) = one(number_eltype(eti)) return typeof(sum(map(eti_to_one, args))) end -@inline function allocate_result_type(::AbstractManifold, f::TF, args::Tuple{}) where {TF} - return Float64 +@inline function allocate_result_type(M::AbstractManifold, f::TF, args::Tuple{}) where {TF} + return allocation_promotion_function(M, f, ())(Float64) end """ diff --git a/test/allocation.jl b/test/allocation.jl index b99983ef..3838a20d 100644 --- a/test/allocation.jl +++ b/test/allocation.jl @@ -1,7 +1,8 @@ using ManifoldsBase using Test using ManifoldsBase: - combine_allocation_promotion_functions, allocation_promotion_function, ℝ, ℂ + combine_allocation_promotion_functions, allocation_promotion_function, ℝ, ℂ, ℍ +using Quaternions struct AllocManifold{𝔽} <: AbstractManifold{𝔽} end AllocManifold() = AllocManifold{ℝ}() @@ -19,6 +20,9 @@ ManifoldsBase.representation_size(::AllocManifold2) = (2, 3) struct AllocManifold3 <: AbstractManifold{ℂ} end ManifoldsBase.representation_size(::AllocManifold3) = (2, 3) +struct AllocManifold4 <: AbstractManifold{ℍ} end +ManifoldsBase.representation_size(::AllocManifold4) = (2, 3) + @testset "Allocation" begin a = [[1.0], [2.0], [3.0]] b = [[2.0], [3.0], [-3.0]] @@ -73,6 +77,7 @@ ManifoldsBase.representation_size(::AllocManifold3) = (2, 3) @test size(alloc2) == representation_size(M2) @test ManifoldsBase.allocate_result(AllocManifold3(), rand) isa Matrix{ComplexF64} + @test ManifoldsBase.allocate_result(AllocManifold4(), rand) isa Matrix{QuaternionF64} an = allocate_on(M2) @test an isa Matrix{Float64} From 83452118da6200d951a1957753ff32492edea695 Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Sun, 28 Jul 2024 12:53:25 +0200 Subject: [PATCH 08/13] fix CI? --- Project.toml | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index adb92a5a..1799d424 100644 --- a/Project.toml +++ b/Project.toml @@ -47,6 +47,7 @@ DoubleFloats = "497a8b3b-efae-58df-a0af-a86822472b78" ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" +Quaternions = "94ee1d12-ae83-5a48-8b1c-48b8ff168ae0" RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd" ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" @@ -54,4 +55,4 @@ Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] -test = ["Test", "Aqua", "DoubleFloats", "ForwardDiff", "OrdinaryDiffEq", "Plots", "ReverseDiff", "StaticArrays", "Statistics", "RecursiveArrayTools"] +test = ["Test", "Aqua", "DoubleFloats", "ForwardDiff", "OrdinaryDiffEq", "Plots", "Quaternions", "ReverseDiff", "StaticArrays", "Statistics", "RecursiveArrayTools"] From a84ad3cf0798393ae353702df8c738c069780c92 Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Sun, 28 Jul 2024 13:02:50 +0200 Subject: [PATCH 09/13] fix CI on 1.6 too --- src/ManifoldsBase.jl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/ManifoldsBase.jl b/src/ManifoldsBase.jl index 6e6a4b59..b5a43ef7 100644 --- a/src/ManifoldsBase.jl +++ b/src/ManifoldsBase.jl @@ -1186,6 +1186,9 @@ function __init__() @require Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" begin include("../ext/ManifoldsBasePlotsExt.jl") end + @require Quaternions = "94ee1d12-ae83-5a48-8b1c-48b8ff168ae0" begin + include("../ext/ManifoldsBaseQuaternionsExt.jl") + end @require Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" begin include("../ext/ManifoldsBaseStatisticsExt.jl") end From 9fea4931865dd35570bb31b7d399e691d28dfc34 Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Sun, 28 Jul 2024 13:11:42 +0200 Subject: [PATCH 10/13] update Julia versions on CI --- .github/workflows/ci.yml | 2 +- .github/workflows/documenter.yml | 2 +- .github/workflows/format.yml | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index dbb73679..1f5aace2 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -11,7 +11,7 @@ jobs: runs-on: ${{ matrix.os }} strategy: matrix: - julia-version: ["1.6", "1.9", "~1.10.0-0"] + julia-version: ["1.6", "1.10", "~1.11.0-0"] os: [ubuntu-latest, macOS-latest, windows-latest] steps: - uses: actions/checkout@v4 diff --git a/.github/workflows/documenter.yml b/.github/workflows/documenter.yml index 0090d0fb..62d0c5ef 100644 --- a/.github/workflows/documenter.yml +++ b/.github/workflows/documenter.yml @@ -16,7 +16,7 @@ jobs: version: 1.3.353 - uses: julia-actions/setup-julia@latest with: - version: 1.9 + version: 1.10 - name: Julia Cache uses: julia-actions/cache@v2 - name: Cache Quarto diff --git a/.github/workflows/format.yml b/.github/workflows/format.yml index f8eb9a2c..f487f6c2 100644 --- a/.github/workflows/format.yml +++ b/.github/workflows/format.yml @@ -13,7 +13,7 @@ jobs: - uses: actions/checkout@v4 - uses: julia-actions/setup-julia@v2 with: - version: 1.9 + version: 1.10 - name: Install JuliaFormatter and format run: | julia -e 'using Pkg; Pkg.add(PackageSpec(name="JuliaFormatter"))' From 771aa4b812da3ae6d002114428d63f0f1aae2735 Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Sun, 28 Jul 2024 14:34:46 +0200 Subject: [PATCH 11/13] 1.10 == 1.1 --- .github/workflows/documenter.yml | 4 ++-- .github/workflows/format.yml | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/.github/workflows/documenter.yml b/.github/workflows/documenter.yml index 62d0c5ef..0c5e136f 100644 --- a/.github/workflows/documenter.yml +++ b/.github/workflows/documenter.yml @@ -13,10 +13,10 @@ jobs: - uses: actions/checkout@v4 - uses: quarto-dev/quarto-actions/setup@v2 with: - version: 1.3.353 + version: "1.4.550" - uses: julia-actions/setup-julia@latest with: - version: 1.10 + version: "1.10" - name: Julia Cache uses: julia-actions/cache@v2 - name: Cache Quarto diff --git a/.github/workflows/format.yml b/.github/workflows/format.yml index f487f6c2..67a9ff6c 100644 --- a/.github/workflows/format.yml +++ b/.github/workflows/format.yml @@ -13,7 +13,7 @@ jobs: - uses: actions/checkout@v4 - uses: julia-actions/setup-julia@v2 with: - version: 1.10 + version: "1.10" - name: Install JuliaFormatter and format run: | julia -e 'using Pkg; Pkg.add(PackageSpec(name="JuliaFormatter"))' From bcd7cc92b396ca9b39dbb1773ee817ccd020d2fa Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Sun, 28 Jul 2024 14:41:42 +0200 Subject: [PATCH 12/13] quarto --- .github/workflows/documenter.yml | 2 +- tutorials/_quarto.yml | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/.github/workflows/documenter.yml b/.github/workflows/documenter.yml index 0c5e136f..3079c7d1 100644 --- a/.github/workflows/documenter.yml +++ b/.github/workflows/documenter.yml @@ -13,7 +13,7 @@ jobs: - uses: actions/checkout@v4 - uses: quarto-dev/quarto-actions/setup@v2 with: - version: "1.4.550" + version: "1.5.55" - uses: julia-actions/setup-julia@latest with: version: "1.10" diff --git a/tutorials/_quarto.yml b/tutorials/_quarto.yml index 50bbc704..e36da3f6 100644 --- a/tutorials/_quarto.yml +++ b/tutorials/_quarto.yml @@ -24,4 +24,4 @@ format: variant: -raw_html+tex_math_dollars wrap: preserve -jupyter: julia-1.9 \ No newline at end of file +jupyter: julia-1.10 \ No newline at end of file From 5dfb16c0fdd39d39e595c47881f1bc3c2017d102 Mon Sep 17 00:00:00 2001 From: Mateusz Baran Date: Sun, 28 Jul 2024 16:02:34 +0200 Subject: [PATCH 13/13] update date --- NEWS.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/NEWS.md b/NEWS.md index abd7a980..f303ecfc 100644 --- a/NEWS.md +++ b/NEWS.md @@ -5,7 +5,7 @@ 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.11] 23/07/2024 +## [0.15.11] 28/07/2024 ### Added