From 536cc72663a2908090eb69a4f4257472713804c9 Mon Sep 17 00:00:00 2001
From: Johannes Terblanche <Affie@users.noreply.github.com>
Date: Wed, 11 Aug 2021 17:28:59 +0200
Subject: [PATCH 1/2] Mutating function in _solveLambdaNumeric cost

---
 Project.toml                 | 4 ++--
 src/NumericalCalculations.jl | 9 +++++++--
 2 files changed, 9 insertions(+), 4 deletions(-)

diff --git a/Project.toml b/Project.toml
index 024a7f775..0579d88be 100644
--- a/Project.toml
+++ b/Project.toml
@@ -55,8 +55,8 @@ FunctionalStateMachine = "0.2.9"
 JLD2 = "0.3, 0.4"
 JSON2 = "0.3"
 KernelDensityEstimate = "0.5.6"
-Manifolds = "0.6"
-ManifoldsBase = "0.11, 0.12"
+Manifolds = "0.6.3"
+ManifoldsBase = "0.12.6"
 MetaGraphs = "0.6.4"
 NLSolversBase = "7.6"
 NLsolve = "3, 4"
diff --git a/src/NumericalCalculations.jl b/src/NumericalCalculations.jl
index 956e2869c..1767d6933 100644
--- a/src/NumericalCalculations.jl
+++ b/src/NumericalCalculations.jl
@@ -94,8 +94,9 @@ function _solveLambdaNumeric( fcttype::Union{F,<:Mixture{N_,F,S,T}},
   # norm(M, p, X) == distance(M, p, X)
   #TODO fix closure for performance
   fM = getManifold(fcttype)
-  function cost(Xc)
-    p = exp(M, ϵ, hat(M, ϵ, Xc))  
+  function cost(p, X, Xc)
+    hat!(M, X, ϵ, Xc)
+    exp!(M, p, ϵ, X)  
     # X = objResX(p)
     # return norm(fM, p, X)^2 #TODO the manifold of p and X are not always the same
     #options getPointIdentity or leave it to factor 
@@ -110,6 +111,10 @@ function _solveLambdaNumeric( fcttype::Union{F,<:Mixture{N_,F,S,T}},
     Optim.optimize(cost, X0c, Optim.NelderMead())
   end
 
+  X0 = hat(M, ϵ, X0c)
+  p0 = exp(M, ϵ, X0)
+  r = Optim.optimize(Xc->cost(p0, X0, Xc), X0c, alg)
+  
   return exp(M, ϵ, hat(M, ϵ, r.minimizer)) 
 
 end

From 3955af877fcadac96a780d3368525ac2c5ea8487 Mon Sep 17 00:00:00 2001
From: Johannes Terblanche <Affie@users.noreply.github.com>
Date: Wed, 11 Aug 2021 20:42:36 +0200
Subject: [PATCH 2/2] Fix code messed up by merge

---
 src/NumericalCalculations.jl | 14 +++++++-------
 1 file changed, 7 insertions(+), 7 deletions(-)

diff --git a/src/NumericalCalculations.jl b/src/NumericalCalculations.jl
index 1767d6933..8544d6a26 100644
--- a/src/NumericalCalculations.jl
+++ b/src/NumericalCalculations.jl
@@ -104,13 +104,13 @@ function _solveLambdaNumeric( fcttype::Union{F,<:Mixture{N_,F,S,T}},
     return sum(residual.^2)
   end
 
-  # separate statements to try improve type-stability 
-  r = if islen1
-    Optim.optimize(cost, X0c, Optim.BFGS())
-  else
-    Optim.optimize(cost, X0c, Optim.NelderMead())
-  end
-
+  # # separate statements to try improve type-stability 
+  # r = if islen1
+  #   Optim.optimize(cost, X0c, Optim.BFGS())
+  # else
+  #   Optim.optimize(cost, X0c, Optim.NelderMead())
+  # end
+  alg = islen1 ? Optim.BFGS() : Optim.NelderMead()
   X0 = hat(M, ϵ, X0c)
   p0 = exp(M, ϵ, X0)
   r = Optim.optimize(Xc->cost(p0, X0, Xc), X0c, alg)