Julia gives "MethodError: rand!" when trying to optimize over a certain manifold using Manopt

Question

My goal is to find the minimum of a certain function over the manifold of doubly stochastic matrices. I use the package Manopt (https://manoptjl.org/stable/) for this. The relevant part of my code looks like this:

using Random
using Graphs
using Manopt
using Manifolds
using LinearAlgebra

n = 8
G = DiGraph(n)
for i = 1:n
    add_edge!(G, i%n +1, (i+1)%n +1)
    add_edge!(G, (i+1)%n +1, i%n +1)
end

A = adjacency_matrix(G)
P_0 = A / 2
M = MultinomialMatrices(n, n)
f(M, p) = 0 #sum((P_0 - p).^ 2)

particle_swarm(M, f)

Of course I need to change f to a more complex function, however, the code is not even working in this case of f=0. I also tried using the Nelder Mead Method instead of the partical swarm, yielding the same error. The error is the following: "MethodError: rand!(::TaskLocalRNG, ::AbstractManifold, ::SubArray{Float64, 1, Matrix{Float64}, Tuple{Base.Slice{Base.OneTo{Int64}}, Int64}, true}) is ambiguous. Candidates: [...]" and it keeps on listing 50-70 candidates.

Does anyone know what I am doing wrong?

I tried using different manifolds to optimize over and the problem didn't occur in the case of "simpler" one (eg. the sphere). I tried using other methods for optimizing, none of them worked.

Answer 1

The reason behind the initial error (and the follow ups in the comments) is relatively simple.

The algorithms from Manopt.jl follow the ManifoldsBase.jl “interface” of a typical manifold and expect a certain set of functions for a specific manifold to be implemented. For often includes for example a retraction. Details are usually documented on the Technical Details section for every solver, like the one in the question here particle swarm. For now the technical details do not mention rand/rand!, since that is not necessarily needed for the algorithm itself, but just for the initialisation of default values (in any algorithm actually). One could maybe improve that.

Sometimes these functions are not all implemented in Manifolds.jl, either because no one came along yet needing them, or (more likely) because the people implementing a certain manifold were not aware of a formula / paper / reference for there. For example the exponential map is not available in closed form on certain manifolds (though then several retractions usually still are).

For the specific manifolds at hand, I took a peek into literature and noticed that the rand functions could actually be provided from what is given in literature, e.g.

• for MultinomialSymmetric
• for MultinomialDoubleStochastic
• and even for a new manifold I could easily add, MultinomialSymmetricPositiveDefinite, where we could maybe move the doc string also and attach it to rand (which by default uses rand! internally)

All these are now thoroughly tested, documented, and available in the new Manifolds.jl version 0.9.13.