Introduce fill and fill! on Power manifolds. (#190)

* Implement and test fill(M,p)
* adapt tests.
* bump version.
* Remove the ambiguity.
* `fill!` for array power manifold

---------

Co-authored-by: Mateusz Baran <[email protected]>

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.10] 15/05/2024
+
+### Added
+
+* Functions `fill(p, N)` and `fill!(P, p, N)` to fill values into a point on a power manifold `N`.
+
 ## [0.15.9] 02/05/2024
 
 ### Added

Project.toml

@@ -1,7 +1,7 @@
 name = "ManifoldsBase"
 uuid = "3362f125-f0bb-47a3-aa74-596ffd7ef2fb"
 authors = ["Seth Axen <[email protected]>", "Mateusz Baran <[email protected]>", "Ronny Bergmann <[email protected]>", "Antoine Levitt <[email protected]>"]
-version = "0.15.9"
+version = "0.15.10"
 
 [deps]
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"

src/ManifoldsBase.jl

@@ -9,6 +9,8 @@ import Base:
     copyto!,
     angle,
     eltype,
+    fill,
+    fill!,
     isempty,
     length,
     similar,
@@ -1241,6 +1243,8 @@ export ×,
     embed!,
     embed_project,
     embed_project!,
+    fill,
+    fill!,
     geodesic,
     geodesic!,
     get_basis,

src/PowerManifold.jl

@@ -72,6 +72,9 @@ can be used to nest manifold, i.e. `PowerManifold(M, NestedPowerRepresentation()
 represents vectors of length 2 whose elements are vectors of length 3 of points on N
 in a nested array representation.
 
+The third signature `M^(...)` is equivalent to the first one, and hence either yields
+a combination of power manifolds to _one_ larger power manifold, or a power manifold with the default representation.
+
 Since there is no default [`AbstractPowerRepresentation`](@ref) within this interface, the
 `^` operator is only available for `PowerManifold`s and concatenates dimensions.
 
@@ -525,6 +528,45 @@ function exp!(M::PowerManifoldNestedReplacing, q, p, X)
     return q
 end
 
+@doc raw"""
+    fill(p, M::AbstractPowerManifold)
+
+Create a point on the [`AbstractPowerManifold`](@ref) `M`, where every entry is set to the
+point `p`.
+
+!!! note
+    while usually the manifold is a first argument in all functions in `ManifoldsBase.jl`,
+    we follow the signature of `fill`, where the power manifold serves are the size information.
+"""
+function fill(p, M::AbstractPowerManifold)
+    P = allocate_result(M, rand) # rand finds the right way to allocate our point usually
+    return fill!(P, p, M)
+end
+
+@doc raw"""
+    fill!(P, p, M::AbstractPowerManifold)
+
+Fill a point `P` on the [`AbstractPowerManifold`](@ref) `M`, setting every entry to `p`.
+
+!!! note
+    while usually the manifold is the first argument in all functions in `ManifoldsBase.jl`,
+    we follow the signature of `fill!`, where the power manifold serves are the size information.
+"""
+function fill!(P, p, M::PowerManifoldNestedReplacing)
+    for i in get_iterator(M)
+        P[M, i] = p
+    end
+    return P
+end
+function fill!(P, p, M::AbstractPowerManifold)
+    rep_size = representation_size(M.manifold)
+    for i in get_iterator(M)
+        copyto!(M.manifold, _write(M, rep_size, P, i), p)
+    end
+    return P
+end
+
+
 function get_basis(M::AbstractPowerManifold, p, B::AbstractBasis)
     rep_size = representation_size(M.manifold)
     vs = [get_basis(M.manifold, _read(M, rep_size, p, i), B) for i in get_iterator(M)]

test/power.jl

@@ -1,7 +1,12 @@
 using Test
 using ManifoldsBase
 using ManifoldsBase:
-    AbstractNumbers, ℝ, ℂ, NestedReplacingPowerRepresentation, VectorSpaceType
+    AbstractNumbers,
+    DefaultManifold,
+    ℝ,
+    ℂ,
+    NestedReplacingPowerRepresentation,
+    VectorSpaceType
 using StaticArrays
 using LinearAlgebra
 using Random
@@ -24,6 +29,23 @@ end
 
 struct TestArrayRepresentation <: AbstractPowerRepresentation end
 
+const TestPowerManifoldMultidimensional =
+    AbstractPowerManifold{𝔽,<:AbstractManifold{𝔽},TestArrayRepresentation} where {𝔽}
+
+function ManifoldsBase.representation_size(M::TestPowerManifoldMultidimensional)
+    return (representation_size(M.manifold)..., ManifoldsBase.get_parameter(M.size)...)
+end
+
+@inline function ManifoldsBase._write(
+    ::TestPowerManifoldMultidimensional,
+    rep_size::Tuple,
+    x::AbstractArray,
+    i::Tuple,
+)
+    return view(x, ManifoldsBase.rep_size_to_colons(rep_size)..., i...)
+end
+
+
 @testset "Power Manifold" begin
 
     @testset "Power Manifold with a test representation" begin
@@ -334,7 +356,7 @@ struct TestArrayRepresentation <: AbstractPowerRepresentation end
         @test norm(N, P, Z .- Zc) ≈ 0
     end
 
-    @testset "Other stuff" begin
+    @testset "Curvature" begin
         M1 = TestSphere(2)
         @testset "Weingarten" begin
             Mpr = PowerManifold(M1, NestedPowerRepresentation(), 2)
@@ -388,4 +410,36 @@ struct TestArrayRepresentation <: AbstractPowerRepresentation end
         p = rand(N)
         @test zero_vector(N, p) == 0 .* p
     end
+
+    @testset "fill" begin
+        M = ManifoldsBase.DefaultManifold(3)
+        N = PowerManifold(M, NestedPowerRepresentation(), 2)
+        p = [1.0, 2.0, 3.0]
+        P1 = fill(p, N)
+        @test P1[N, 1] == p
+        @test P1[N, 2] == p
+        P2 = [zeros(3), zeros(3)]
+        fill!(P2, p, N)
+        @test P2[N, 1] == p
+        @test P2[N, 2] == p
+
+        M = ManifoldsBase.DefaultManifold(3)
+        NR = PowerManifold(M, NestedReplacingPowerRepresentation(), 2)
+        P1 = fill(p, NR)
+        @test P1[NR, 1] === p
+        @test P1[NR, 2] === p
+        P2 = [zeros(3), zeros(3)]
+        fill!(P2, p, NR)
+        @test P2[NR, 1] === p
+        @test P2[NR, 2] === p
+
+        NAR = PowerManifold(M, TestArrayRepresentation(), 2)
+        P1 = fill(p, NAR)
+        @test P1 isa Matrix{Float64}
+        @test P1 == [1.0 1.0; 2.0 2.0; 3.0 3.0]
+
+        P2 = similar(P1)
+        fill!(P2, p, NAR)
+        @test P2 == [1.0 1.0; 2.0 2.0; 3.0 3.0]
+    end
 end