Fix a few bugs.

src/plans/vectorial_plan.jl

@@ -708,7 +708,7 @@ function get_hessian!(
     for (j, f) in zip(1:n, vhf.hessians!![i])
         copyto!(M, _write(pM, rep_size, Y, (j,)), f(M, p, X))
     end
-    return X
+    return Y
 end
 function get_hessian!(
     M::AbstractManifold,
@@ -721,10 +721,11 @@ function get_hessian!(
 ) where {FT,JT}
     n = _vgf_index_to_length(i, vgf.range_dimension)
     pM = PowerManifold(M, range, n)
+    rep_size = representation_size(M)
     for (j, f) in enumerate(vgf.hessians!!)
         copyto!(M, _write(pM, rep_size, Y, (j,)), p, f(M, p, X))
     end
-    return X
+    return Y
 end
 # Part I(c) A single gradient function
 function get_hessian!(
@@ -803,7 +804,7 @@ function get_hessian!(
     y = zero_vector(pM, P)
     vhf.hessians!!(M, y, p, X)
     copyto!(M, Y, p, y[pM, i])
-    return X
+    return Y
 end
 function get_hessian!(
     M::AbstractManifold,
@@ -823,7 +824,7 @@ function get_hessian!(
     vhf.hessians!!(M, y, p, X)
     # Luckily all documented access functions work directly on `x[pM_temp,...]`
     copyto!(pM_out, Y, P[pM_temp, i], y[pM_temp, i])
-    return X
+    return Y
 end
 
 get_hessian_function(vgf::VectorGradientFunction, recursive=false) = vgf.hessians!!

test/plans/test_vectorial_plan.jl

@@ -5,8 +5,8 @@ using Manopt: get_value, get_value_function, get_gradient_function
     g(M, p) = [p[1] - 1, -p[2] - 1]
     # # Function
     grad_g(M, p) = [[1.0, 0.0, 0.0], [0.0, -1.0, 0.0]]
-    hess_g(M, p, X) = [copy(X), copy(X)]
-    hess_g!(M, Y, p, X) = (Y .= [copy(X), copy(X)])
+    hess_g(M, p, X) = [copy(X), -copy(X)]
+    hess_g!(M, Y, p, X) = (Y .= [copy(X), -copy(X)])
     # since the ONB of M is just the identity in coefficients, JF is gradients'
     jac_g(M, p) = [1.0 0.0; 0.0 -1.0; 0.0 0.0]'
     jac_g!(M, J, p) = (J .= [1.0 0.0; 0.0 -1.0; 0.0 0.0]')
@@ -24,7 +24,7 @@ using Manopt: get_value, get_value_function, get_gradient_function
     g2(M, p) = -p[2] - 1
     grad_g2(M, p) = [0.0, -1.0, 0.0]
     grad_g2!(M, X, p) = (X .= [0.0, -1.0, 0.0])
-    hess_g2(M, p, X) = -copy(X)
+    hess_g2(M, p, X) = copy(-X)
     hess_g2!(M, Y, p, X) = copyto!(Y, -X)
     # verify a few case
     vgf_fa = VectorGradientFunction(g, grad_g, 2)
@@ -69,6 +69,7 @@ using Manopt: get_value, get_value_function, get_gradient_function
         2;
         function_type=ComponentVectorialType(),
         jacobian_type=ComponentVectorialType(),
+        hessian_type=ComponentVectorialType(),
     )
     vhf_fi = VectorHessianFunction(g, grad_g!, hess_g!, 2; evaluation=InplaceEvaluation())
     vhf_vi = VectorHessianFunction(
@@ -78,6 +79,7 @@ using Manopt: get_value, get_value_function, get_gradient_function
         2;
         function_type=ComponentVectorialType(),
         jacobian_type=ComponentVectorialType(),
+        hessian_type=ComponentVectorialType(),
         evaluation=InplaceEvaluation(),
     )
 
@@ -103,22 +105,28 @@ using Manopt: get_value, get_value_function, get_gradient_function
         @test Z == gg[2]
     end
 
+    X = [1.0, 0.5, 0.25]
+    gh = [X, -X]
     # Hessian
-    for vgf in
-        [vgf_fa, vgf_va, vgf_fi, vgf_vi, vgf_ja, vgf_ji, vhf_fa, vhf_fi, vhf_va, vhf_vi]
-        @test get_gradient(M, vgf, p) == gg
-        @test get_gradient(M, vgf, p, :) == gg
-        @test get_gradient(M, vgf, p, 1:2) == gg
-        @test get_gradient(M, vgf, p, [1, 2]) == gg
-        @test get_gradient(M, vgf, p, 1) == gg[1]
-        @test get_gradient(M, vgf, p, 2) == gg[2]
+    for vhf in [
+        vhf_fa,
+        vhf_va,
+        vhf_fi,
+        #vhf_vi,
+    ]
+        @test get_hessian(M, vhf, p, X) == gh
+        @test get_hessian(M, vhf, p, X, :) == gh
+        @test get_hessian(M, vhf, p, X, 1:2) == gh
+        @test get_hessian(M, vhf, p, X, [1, 2]) == gh
+        @test get_hessian(M, vhf, p, X, 1) == gh[1]
+        @test get_hessian(M, vhf, p, X, 2) == gh[2]
         Y = [zero_vector(M, p), zero_vector(M, p)]
-        get_gradient!(M, Y, vgf, p, :)
-        @test Y == gg
+        get_hessian!(M, Y, vhf, p, X, :)
+        @test Y == gh
         Z = zero_vector(M, p)
-        get_gradient!(M, Z, vgf, p, 1)
-        @test Z == gg[1]
-        get_gradient!(M, Z, vgf, p, 2)
-        @test Z == gg[2]
+        get_hessian!(M, Z, vhf, p, X, 1)
+        @test Z == gh[1]
+        get_hessian!(M, Z, vhf, p, X, 2)
+        @test Z == gh[2]
     end
 end